account for temporal correlation

Thanks a lot,

just to make sure i got it right,

if (using the real dataset) from the LogLikelihood ratio test model1 isn't "better" than model,
means that temporal auto correlation isn't seriously affecting the model?

and, I shouldn't be nesting time within subject because is implicit
that observation from the same subject are the repeated measures?

The need of nesting would for example be:
an experiment ?also having spatial correlation, to say, if subjects are are also grouped by their geographical position,
should i in this case be nesting location within subject?

cheers
m.

Il giorno 18 Nov 2011, alle ore 18:26, Ben Bolker ha scritto:

matteo dossena <matteo.dossena at ...> writes:

Dear list,

I have a data frame like this:

set.seed(5)
mydata <- data.frame(var = rnorm(100,20,1),
? ? ? ? ? ? ? ? ? ? ?temp = sin(sort(rep(c(1:10),10))),
? ? ? ? ? ? ? ? ? ? ?subj = as.factor(rep(c(1:10),5)),
? ? ? ? ? ? ? ? ? ? ?time = sort(rep(c(1:10),10)),
? ? ? ? ? ? ? ? ? ? ?trt = rep(c("A","B"), 50))

I need to model the response of var as a function of temp*trt
and to do so I'm using the following model:

library(nlme)
model <- lme(var~temp*trt,random=~1|subj,mydata)

however, since i have repeated measurement of the same subject,
i.e. I measured var in subj1 at time1,2,3..
I must consider the non independence of the residuals.
moreover, temp is also a function of time, but i'm not sure how
to include this in my model.

I'm following the approach in Zuur et al 2009, so I checked for
temporal auto-correlation using the
function afc()
In fact the residuals follow the temporal patter of temperature with time.
However, here I'm not interested in the temporal dependence of temperature
and consequently the effect of
this on var.
Instead what i need to do is to account for the
correlation between consecutive measures (made on the same
subject) in the error term of the model.

here is my attempt to do so:

model1 <- lme(var~temp*trt,random=~1|subj,
? ? ? ? ? ? ? ?correlation=corAR1(form=~1|subj),mydata)

model1$modelStruct$corStruct

Correlation structure of class corAR1 representing
? ? ? ?Phi
-0.05565362

?You shouldn't be nesting time within subject. ?'subject' is all the grouping
you need here.

intervals(model1)

gives (approximate!!) CI for the correlation structure parameter
(-0.27,0.77) in this case

?Of course, in this case we don't expect to find anything interesting
because these are simulated data without any correlation built in.

?These plots are useful.

plot(ACF(model),alpha=0.05)
plot(ACF(model1),alpha=0.05) ?## should be ALMOST identical to the one above
## taking correlation into account:
plot(ACF(model1,resType="normalized"),alpha=0.05)