# te( ) interactions and AIC model selection with GAM in R

I'm working with a time-series of several years and to analyze it, I’m using GAM smoothers from the package mgcv. I’m constructing models where zooplankton biomass (bm) is the dependent variable and the continuous explanatory variables are:

-time in Julian days (t), to creat a long-term linear trend

-Julian days of the year (t_year) to create an annual cycle

-Mean temperature of Winter (temp_W), Temperature of September (temp_sept) or Chla.

Questions:

1) To introduce a tensor product modifying the annual cycle in my model, I tried 2 different approaches:

a) gam( bm ~ t + te (t_year, temp_W, temp_sept, k = c( 5,30 ), d = ( 1,2), bs = c( “cc”,”cr” ) ), data = data )

b) gam( bm ~ t + te ( t_year, temp_W, temp_sept, k = 5, bs = c( “cc”,”cr”,”cr” ) ), data = data )

Here is my problem: when I’m using just 2 variables (e.g., t_year and temp_W) for the tensor product, I can understand pretty well how the interpolation works and visualize it with vis.gam() as a 3d plot or a contour one. But with 3 variables it is difficult for me to understand how it works. Besides, I don’t know which one is the proper way to construct it, a) or b). Finally, when I plot a) or b) as vis.gam (model_name , view= c(“t_year”, “temp_W”)), How should I interpret the plot? The effect of temp_W on the annual cycle after considering already the effect of temp_sept or just the individual effect of temp_W on the annual cycle?

2) I’m trying to do a model selection using AIC criteria. I have several questions about it:

Should I use always the same type of smoothing basis (bs), the same type of smoother ( e.g te) and the same dimension of the basis (k)? Example:

Option 1:

a) mod1 <- gam(bm ~ t, data = data )

b)mod2 <- gam( bm ~ te ( t, k = 5, bs = “cr” ), data = data )

c) mod3 <- gam( bm ~ te ( t_year, k = 5, bs = “cc”), data = data )

d) mod4 <- gam( bm ~ te ( t_year, temp_W, k = 5, bs = c( “cc”,”cr” ) ), data = data )

e) mod5 <- gam( bm ~ te ( t_year, temp_W, temp_sept, k = 5, bs = c( “cc”,”cr”,”cr” ) ), data = data ).

Here the limitation for k = 5, is due to mod5, I don’t use s () because in mod4 and mod5 te () is used and finally, I always use “cr” and “cc”.

Option 2:

a) mod1 <- gam( bm ~ t, data = data )

b) mod2 <- gam( bm ~ s ( t, k = 13, bs = “cr” ), data = data )

c) mod3 <- gam( bm ~ s( t_year, k = 13, bs = “cc” ), data = data )

d) mod4 <- gam( bm ~ te( t_year, temp_W, k = 11, bs = c( “cc”,”cr” ) ), data = data)

e) mod5 <- gam( bm ~ te( t_year, temp_W, temp_sept, k = 5, bs = c( “cc”,”cr”,”cr” ) ), data = data )

I can get lower AIC for each of the models with Option 2, but are they comparable when I use AIC criteria? Is it therefore the proper way to do it as in Option 1?

AIC (mod1, mod2, mod3, mod4, mod5).

As an example of how the data frame looks like:

> time_series_data

bm      Chla year month       t t_year temp_W temp_sept
1  2.1680335 54.718891 1993     1 2449009     20   12.1      19.3
2  4.6180770 29.372938 1993     2 2449043     54   12.1      19.3
3  4.6871990 99.198623 1993     3 2449064     75   12.1      19.3
4  4.9862020 59.835987 1993     4 2449094    105   12.1      19.3
5  3.4977156 79.990143 1993     5 2449120    131   12.1      19.3
6  3.1030763 68.018739 1993     6 2449148    159   12.1      19.3
7  2.0312841 70.850406 1993     7 2449181    192   12.1      19.3
8  1.2477797 62.381780 1993     8 2449211    222   12.1      19.3
9  2.1445538 99.094776 1993     9 2449254    265   12.1      19.3
10 6.7026438 82.397907 1993    10 2449282    293   12.1      19.3
11 1.6524655 44.977256 1993    11 2449303    314   12.1      19.3
12 2.1627389 52.624779 1993    12 2449342    353   12.1      19.3
13 3.0981200 58.274128 1994     1 2449374     20   11.3      18.6
14 2.4342291 14.733698 1994     2 2449402     48   11.3      18.6
15 4.8691345 51.508774 1994     3 2449431     77   11.3      18.6
16 3.7366294 38.206928 1994     4 2449458    104   11.3      18.6
17 3.3565706 72.763028 1994     5 2449500    146   11.3      18.6
18 2.7869220 81.265662 1994     6 2449520    166   11.3      18.6
19 2.6971469 50.692921 1994     7 2449540    186   11.3      18.6
20 1.3758862 94.396013 1994     8 2449569    215   11.3      18.6
21 5.7578197 59.357898 1994     9 2449620    266   11.3      18.6
22 2.8941763 21.974925 1994    10 2449645    291   11.3      18.6
23 0.9530070  7.781981 1994    11 2449673    319   11.3      18.6
24 0.3713342 84.950835 1994    12 2449697    343   11.3      18.6

-
Can you explain what temp_W and temp_sept really are. You say mean but does this mean they are scalar (same for each observation in same year??) –  Gavin Simpson Jul 31 '12 at 13:27
@Gavin Simpson Yes, that's right. There is one value of temp_W for each zooplankton biomass value in a particular year, and the same happens with temp_sept. –  Ricardo González-Gil Jul 31 '12 at 14:33