[R] autocorrelation in gams

[Kingsford Jones]

Keeping Gavin's advice in mind, you may also want to look at ?acf (and
see section 14.1 of MASS) and help(ACF, package=nlme) (see section 5.3
of MEMSS).  These are useful functions for exploring the 1d empirical
autocorrelation structure of model residuals.

hth,
Kingsford Jones

On Fri, Aug 15, 2008 at 1:18 AM, Gavin Simpson <gavin.simpson at ucl.ac.uk> wrote:
> On Thu, 2008-08-14 at 16:12 +0100, Abigail McQuatters-Gollop wrote:
>> Hi,
>>
>> I am looking at the effects of two explanatory variables on chlorophyll.
>> The data are an annual time-series (so are autocorrelated) and the
>> relationships are non-linear. I want to account for autocorrelation in
>> my model.
>>
>>
>>
>> The model I am trying to use is this:
>>
>>
>>
>> Library(mgcv)
>>
>>
>>
>> gam1 <-gam(Chl~s(wintersecchi)+s(SST),family=gaussian,
>> na.action=na.omit, correlation=corAR1(form =~ Year))
>>
>
> There is no correlation argument in mgcv::gam you need gamm(). gam() has
> a ... argument which I suspect is silently mopping up your correlation
> argument so that no error/warning is raised.
>
> Note that gamm() uses lme from the nlme package (in the Gaussian case)
> and works that function very hard (see Wood 2006 GAM book). In my
> experience with this function separating trend from the correlation is
> quite difficult when also estimating the degree of smoothness and one
> has to work hard with the options.
>
> As an alternative you might also take a look at the paper by Ferguson et
> al (2007):
>
> http://www3.interscience.wiley.com/journal/119392062/abstract?CRETRY=1&SRETRY=0
>
> Which has R code using the sm package to do a very similar thing.
>
> HTH
>
> G
>
>>
>> the result I get is this:
>>
>>
>>
>> Family: gaussian
>>
>> Link function: identity
>>
>>
>>
>> Formula:
>>
>> CPRChl ~ s(wintersecchi) + s(SST)
>>
>>
>>
>> Parametric coefficients:
>>
>>             Estimate Std. Error t value Pr(>|t|)
>>
>> (Intercept)  3.57000    0.05061   70.54   <2e-16 ***
>>
>> ---
>>
>> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>>
>>
>>
>> Approximate significance of smooth terms:
>>
>>                   edf Est.rank     F p-value
>>
>> s(wintersecchi) 2.445        5 4.672 0.00887 **
>>
>> s(SST)          2.408        5 4.301 0.01237 *
>>
>> ---
>>
>> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
>>
>>
>>
>> R-sq.(adj) =  0.676   Deviance explained = 75.4%
>>
>> GCV score = 0.074563   Scale est. = 0.053781  n = 21
>>
>>
>>
>> The result look good - significant, with a lot of deviance explained,
>> but I am not convinced the model is actually accounting for the
>> autocorrelation (based on the formula in the results). How can I tell?
>>
>>
>>
>> Many thanks,
>>
>>
>>
>>
>>
>>
>>
>> Dr Abigail McQuatters-Gollop
>>
>> Sir Alister Hardy Foundation for Ocean Science (SAHFOS)
>>
>> The Laboratory
>>
>> Citadel Hill
>>
>> Plymouth UK PL1 2PB
>>
>> tel: +44 1752 633233
>>
>>
>>
>>
>>       [[alternative HTML version deleted]]
>>
>> ______________________________________________
>> R-help at r-project.org mailing list
>> https://stat.ethz.ch/mailman/listinfo/r-help
>> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
>> and provide commented, minimal, self-contained, reproducible code.
>
> ______________________________________________
> R-help at r-project.org mailing list
> https://stat.ethz.ch/mailman/listinfo/r-help
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
>