[R] How to construct confidence bands from a gls fit?

Maik Renner maikrenner at googlemail.com
Thu Apr 23 19:42:54 CEST 2009


Dear R-list,

I would like to show the implications of estimating a linear trend to
time series,
which contain significant serial correlation.

I want to demonstrate this, comparing lm() and an gls() fits, using
the LakeHuron
data set, available in R.
Now in my particular case I would like to draw confidence bands on the plot and
show that there are differences. Unfortunately, I do not know how to
construct those bands
for the gls fit object.
Does anyone know how to do that?

Here is the code:

library(nlme)
library(car)
lm.hu  <- lm(LakeHuron~time(LakeHuron))
## use durbin-watson test to estimate the AR(p) model
durbin.watson(lm.hu,max.lag=7)
# I decided to use a AR(2) model
gls.hu <- gls(LakeHuron ~ time(LakeHuron), correlation=corARMA(p=2),
method="ML")

plot(LakeHuron, main="Lake level time series (Lake Huron)")
lines(as.numeric(time(LakeHuron)),predict(lm.hu,interval=c("conf"))[,1],lty=1,col=2)
lines(as.numeric(time(LakeHuron)),predict(lm.hu,interval=c("conf"))[,2],lty=2,col=2)
lines(as.numeric(time(LakeHuron)),predict(lm.hu,interval=c("conf"))[,3],lty=2,col=2)
## predict of the gls class does not contain confidence intervals
lines(as.numeric(time(LakeHuron)),predict(gls.hu,interval=c("conf")),lty=1,col=4)

Thank you,

Maik

--
Maik Renner
Institut für Hydrologie und Meteorologie
Technische Universität Dresden




More information about the R-help mailing list