[R-sig-ME] Retrieve autocorrelation-corrected errors from gls (nlme) or, gamm (mgcv) (Daniel Malter)

[Highland Statistics Ltd]

> ----------------------------------------------------------------------
>
> Message: 1
> Date: Tue, 21 Aug 2012 16:23:43 -0400
> From: Daniel Malter <daniel at umd.edu>
> To: "r-sig-mixed-models at r-project.org"
> 	<r-sig-mixed-models at r-project.org>
> Subject: [R-sig-ME] Retrieve autocorrelation-corrected errors from gls
> 	(nlme) or gamm (mgcv)
> Message-ID:
> 	<DF5B15061A61C2409088F43AE2C87EC6018C3F473B00 at OITMXCMS02VI.AD.UMD.EDU>
> Content-Type: text/plain; charset="us-ascii"
>
> Hi,
>
> In the example below, I am modeling the dependent variable Y as a function of X when the errors are autoregressive at the first lag.  There is a number of ways/functions with which to model this. arima (tseries), gls (nlme), and gamm should produce similar results in the simulated example below, and they do.
>
> However, I need the residuals from this analysis and both gls and gamm seem to return the errors before correction for autocorrelation, whereas arima returns the corrected errors (i.e., the estimate of the innovation). My question is whether there is an easy way to retrieve the corrected errors from gamm or gls. My real data are panel data with AR3 errors for which I model nonlinear effects of the independent variables on the dependent variable. Hand-computing the appropriate errors would be painful.
>
> #simulate data with AR1 errors
>
> set.seed(394857395)
> e<-rnorm(101)
> e<-e[2:101]+0.5*e[1:100]
> x<-rnorm(100)
> y<-x+e
>
> #OLS for comparison
> reg<-lm(y~x)
> summary(reg)
> pacf(residuals(reg))


It is better to use rstandard instead of resid

> #arima
> require(tseries)
> reg1<-arima(y,order=c(1,0,0),xreg=x)
> reg1
>
> #gls
> require(nlme)
> reg2<-gls(y~x,correlation=corAR1())
> summary(reg2)
>
> #gamm
> require(mgcv)
> reg3<-gamm(y~s(x),correlation=corAR1())
> summary(reg3$lme)
>
> par(mfcol=c(1,3))
> pacf(residuals(reg1),main="ARIMA",lag.max=10)
> pacf(residuals(reg2),main="GLS",lag.max=10)
> pacf(residuals(reg3$lme),main="GAMM",lag.max=10)
>


Use:  resid(reg1, type = "n")


Alain

> Thanks,
> Daniel
>


-- 

Dr. Alain F. Zuur
First author of:

1. Analysing Ecological Data (2007).
Zuur, AF, Ieno, EN and Smith, GM. Springer. 680 p.
URL: www.springer.com/0-387-45967-7


2. Mixed effects models and extensions in ecology with R. (2009).
Zuur, AF, Ieno, EN, Walker, N, Saveliev, AA, and Smith, GM. Springer.
http://www.springer.com/life+sci/ecology/book/978-0-387-87457-9


3. A Beginner's Guide to R (2009).
Zuur, AF, Ieno, EN, Meesters, EHWG. Springer
http://www.springer.com/statistics/computational/book/978-0-387-93836-3


4. Zero Inflated Models and Generalized Linear Mixed Models with R. (2012) Zuur, Saveliev, Ieno.
http://www.highstat.com/book4.htm

Other books: http://www.highstat.com/books.htm


Statistical consultancy, courses, data analysis and software
Highland Statistics Ltd.
6 Laverock road
UK - AB41 6FN Newburgh
Tel: 0044 1358 788177
Email: highstat at highstat.com
URL: www.highstat.com
URL: www.brodgar.com