1
$\begingroup$

I'm just starting to experiment with the mgcv package in r. My problem is this - I'm modelling the count of a bird survey in space, with a number of different habitat predictor variables. I have a GAM that has a spatial smooth s(x,y) and then ten parametric linear terms. The GAM takes the form:

~gam(count ~ s(X,Y) A + B + C + D + E + F + G + H + I + J,
             offset = log(d),
             family = poisson,
             link = log)

I have 60 different data sets, each one representing a specific species counts in space over a specific period of time, together with the predictor variables. My objective is to fit a GAM for each dataset, automatically. I want to identify the "best" fit. Having read around a bit, stepwise model selection seems to not be the done thing. Also I'm aware of how smooth terms can be penalised, effectively selecting them out of the fit. But it's not clear how to do this with linear parametric terms. So I'd welcome any methodological advice on how to automate model selection for my particular problem.

$\endgroup$
1
$\begingroup$

I think the easiest way to do this with mgcv would be to not force linearity for the parametric effects, but to admit potential for a small amount of non-linearity and then you can use the double penalty on all terms with select = TRUE.

gam(count ~ s(X,Y) +
      s(A, k = 3) + s(B, k = 3) + s(C, k = 3) + s(D, k = 3) + s(E, k = 3) + 
      s(F, k = 3) + s(G, k = 3) + s(H, k = 3) + s(I, k = 3) + s(J, k = 3) +
      offset(log(d)),
    method = "REML",
    family = poisson(link = "log"),
    select = TRUE)
| cite | improve this answer | |
$\endgroup$
  • $\begingroup$ Thank you. I had tried this after my post, and it does seem to work with some toy data. But when I run it on my actual data I get Error in if (abs(old.score - score) > score.scale * conv.tol) { : missing value where TRUE/FALSE needed - is this a problem with my data? $\endgroup$ – Anthony W Jun 24 at 8:09
  • $\begingroup$ the problem seems to be that a poisson distribution only works with integers. If I change family to family = quasipoisson(link="log") the problem goes a way, presumably because this distribution works with continuous numbers. $\endgroup$ – Anthony W Jun 24 at 8:25
  • $\begingroup$ @AnthonyW Why do you have none integers? You start with counts, right? $\endgroup$ – Gavin Simpson Jun 24 at 15:03
  • $\begingroup$ Yes. But I am taking a mean across several years - so I end up with non-integers. $\endgroup$ – Anthony W Jun 25 at 15:40
  • $\begingroup$ Why not build that into the offset (so you can model the sum rather than the mean), or use all the data and include Year & Year related terms in the model formula? $\endgroup$ – Gavin Simpson Jun 25 at 18:58

Your Answer

By clicking “Post Your Answer”, you agree to our terms of service, privacy policy and cookie policy

Not the answer you're looking for? Browse other questions tagged or ask your own question.