# Modeling academic seasonality with hierarchical GAMs

I am trying to model the seasonality of daily pageviews to calculus-related Wikipedia articles using a hierarchical GAM, assuming that there is a shared 'academic calendar' seasonality and that each page may deviate from this in a smooth manner depending on where the topic tends to fall in the calculus syllabus.

I've tried using a variety of conventional time series models (things in the seasonal ARIMA family, mostly) as well as Facebook's Prophet, and this is the only model I've tried that's been able to give me an estimate useful for deseasonalizing the daily pageview time series -- Wikipedia traffic data is full of outliers, and most models don't do well with that.

Unfortunately, it doesn't appear to scale well. Modeling with 10 pages takes three or four minutes. With 20 pages, about half an hour. I'd like to model this with at least 200 time series.

The data looks more or less like this, with features for page title, week-of-year and day-of-week:

    Weekly Daily           wikipage value
1:     26     4  Improper_integral   182
2:     27     5  Improper_integral   161
3:     27     6  Improper_integral   126
4:     27     7  Improper_integral   125
5:     27     1  Improper_integral   108


My model is set up like this:

bam(value ~ te(Weekly,Daily,
bs=c("cc", "cc"),
k=c(52, 7), m=c(2, 2)) +
te(Weekly, wikipage,
bs=c("cc", "re"),
k=c(52, 50), m=c(2, 2)),
data=matrix_wiki, method="fREML", family="poisson",
knots = list(Weekly = c(0, 52), Daily=c(0,7)))


I have 1217 observations per time series.

I am following the examples from this paper. https://peerj.com/preprints/27320.pdf

Does anybody know how to speed up this model? Alternatively, is there a simpler way I can obtain the clean seasonality estimates that I'm looking for, given that ARIMA-type models don't seem to work for me?

1. Using discrete = TRUE and running some of the code in parallel threads (see arguments to bam()
2. Change to t2() smooths (see the example at the end of ?t2 for a way to specify something exactly equivalent to the fs basis type, and then run the model through gamm4::gamm4() or mgcv::gamm(), which should be quicker.
• Thanks for the response! When I try running it with discrete=TRUE, it tells me that bam can not discretize with this nesting structure. I'm not able to find any explanation of this online. I switched to t2 smooths and tried running the model through gamm4 but fitting the model doesn't seem to be appreciably faster. Perhaps I'm missing something? I'm not exactly sure which part of the ?t2 docs you were referring to. – dakotalk Dec 21 '18 at 16:07
• Oh, sorry; the latest version has borked t2 smooths (see the ChangeLog) I don't think you'll get massive speed-ups as these are complex models. I should add to this answer that most of the issue is likely in the use of 52 basis function for each marginal smooth of Weekly. Do you really need to use that many basis functions for these marginals? – Gavin Simpson Dec 21 '18 at 16:39
• Using (20,4) instead of (52,7) seems to give me somewhat acceptable results with substantial speedup. I don't have a good intuition for how many to use. My previous experience trying to use Fourier components for this problem led me to think that even a lot of components wouldn't be enough. What would you suggest regarding t2? I'd like to run the code in parallel if possible. – dakotalk Dec 21 '18 at 21:13