I'm fitting a GAMM with correlation structure, using a non-Gaussian family. Here's an example of my global model: M0 <- gamm(response ~ var1*var2 + var3 + s(var4) + s(var5) + s(var6,var7), random=list(placeID= ~1), correlation= corAR1(form= ~ year | placeID), data=data, family=quasipoisson)

I'd like to do term selection on the global model to drop any nonrelevant variables, but I'm not sure how it should be done. Based on the help text of gam.selection - and hoping it applies to GAMM as well! - I've used backward stepwise selection so far. The text says "It is perfectly possible to perform backwards selection using p-values in the usual way". Thus, I've started from the smoothers (removing smoothers from any linear and/or nonsignificant variable) and then moved on to linear variables. I would be interested to hear, whether you think this method is useful in GAMM or not. Please, also tell if you have something else to suggest (e.g. GCV mentioned below).

I would've also liked to compare models that differ in their structure regarding coordinates (smoothed or not, interaction included or not). I suspect that s(latitude,longitude) might be overfitting, and I would like to check this somehow. Since the model uses PQL, I gather that AIC is not recommended (although included in lme part).

The gam.selection text mentiones different score: "In general the most logically consistent method to use for deciding which terms to include in the model is to compare GCV/UBRE/ML scores for models with and without the term." Is this method only for comparing models that differ by a single variable? If so, do you have any suggestions on how to compare models suchs as:

• M1 <- gamm(response ~ var1*var2 + s(lat,long)+ s(var5), random=list(placeID= ~1), correlation= corAR1(form= ~ year | placeID), data=data, family=quasipoisson)
• M2 <- gamm(response ~ var1*var2 + lat + s(var5), random=list(placeID= ~1), correlation= corAR1(form= ~ year | placeID), data=data, family=quasipoisson)

I had been under the impression my global GAMM-model M0 uses GCV (the default), but when I call gam.check(M0$gam), I only get the four plots and no GCV score or other text output. Is there a way to acquire a GCV score of a GAMM object?

I'm fitting a GAMM with correlation structure, using a non-Gaussian family. Here's an example of my global model: M0 <- gamm(response ~ var1*var2 + var3 + s(var4) + s(var5) + s(var6,var7), random=list(placeID= ~1), correlation= corAR1(form= ~ year | placeID), data=data, family=quasipoisson)

I'd like to do term selection on the global model to drop any nonrelevant variables, but I'm not sure how it should be done. Based on the help text of gam.selection - and hoping it applies to GAMM as well! - I've used backward stepwise selection so far. The text says "It is perfectly possible to perform backwards selection using p-values in the usual way". Thus, I've started from the smoothers (removing smoothers from any linear and/or nonsignificant variable) and then moved on to linear variables. I would be interested to hear, whether you think this method is useful in GAMM or not. Please, also tell if you have something else to suggest (e.g. GCV mentioned below).

I would've also liked to compare models that differ in their structure regarding coordinates (smoothed or not, interaction included or not). I suspect that s(latitude,longitude) might be overfitting, and I would like to check this somehow. Since the model uses PQL, I gather that AIC is not recommended (although included in lme part).

The gam.selection text mentions different scores: "In general the most logically consistent method to use for deciding which terms to include in the model is to compare GCV/UBRE/ML scores for models with and without the term." Is this method only for comparing models that differ by a single variable? If so, do you have any suggestions on how to compare models such as:

• M1 <- gamm(response ~ var1*var2 + s(lat,long)+ s(var5), random=list(placeID= ~1), correlation= corAR1(form= ~ year | placeID), data=data, family=quasipoisson)
• M2 <- gamm(response ~ var1*var2 + lat + s(var5), random=list(placeID= ~1), correlation= corAR1(form= ~ year | placeID), data=data, family=quasipoisson)

I had been under the impression that my global GAMM-model M0 uses GCV (the default), but when I call gam.check(M0$gam), I only get the four plots and no GCV score or other text output. Is there a way to acquire a GCV score of a GAMM object?