[R-sig-ME] Question about autocorrelation in nlme

[Thierry Onkelinx]

Dear Paul,

Please keep the mailing list in cc.

I'd recommend the book of Pinheiro and Bates on nlme.

Best regards,
ir. Thierry Onkelinx
Instituut voor natuur- en bosonderzoek / Research Institute for Nature
and Forest
team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
Kliniekstraat 25
1070 Anderlecht
Belgium

To call in the statistician after the experiment is done may be no
more than asking him to perform a post-mortem examination: he may be
able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher
The plural of anecdote is not data. ~ Roger Brinner
The combination of some data and an aching desire for an answer does
not ensure that a reasonable answer can be extracted from a given body
of data. ~ John Tukey


2015-09-16 12:23 GMT+02:00 Paul Buerkner <paul.buerkner op gmail.com>:
> Thanks for the quick answer! That makes sense.
>
> Do you know where I can find the complete model formulation used by nlme so
> that I can simulate data from it?
>
> Thanks again!
>  Paul
>
>
>
> 2015-09-16 11:08 GMT+02:00 Thierry Onkelinx <thierry.onkelinx op inbo.be>:
>>
>> Dear Paul,
>>
>> The correlation functions in nlme work on the residuals, not on the
>> response.
>>
>> Best regards,
>> ir. Thierry Onkelinx
>> Instituut voor natuur- en bosonderzoek / Research Institute for Nature
>> and Forest
>> team Biometrie & Kwaliteitszorg / team Biometrics & Quality Assurance
>> Kliniekstraat 25
>> 1070 Anderlecht
>> Belgium
>>
>> To call in the statistician after the experiment is done may be no
>> more than asking him to perform a post-mortem examination: he may be
>> able to say what the experiment died of. ~ Sir Ronald Aylmer Fisher
>> The plural of anecdote is not data. ~ Roger Brinner
>> The combination of some data and an aching desire for an answer does
>> not ensure that a reasonable answer can be extracted from a given body
>> of data. ~ John Tukey
>>
>>
>> 2015-09-16 10:19 GMT+02:00 Paul Buerkner <paul.buerkner op gmail.com>:
>> > Dear list,
>> >
>> > I have a question concerning autocorrelation models as estimated by
>> > nlme.
>> >
>> > When I try to recover an autocorrelative (AR1) effect from simulated
>> > data,
>> > I do not get the expected results, that is nlme gives different
>> > estimations
>> > than I would expect based on the parameters used to simulate the data
>> > (reproducible example below).
>> >
>> > ### autocor with random effects
>> > set.seed(12345)
>> > phi <- 0.5  # autocorrelation
>> > b <- 0.3  # effect of predictor x
>> > x <- sample(1:20, 1000, TRUE)  # predictor x
>> > ID <- rep(1:100, each = 10)  # person IDs
>> > r <- rnorm(100, sd = 1)  # random Intercepts with sd = 1
>> > e <- rnorm(1000, sd = 2)  # residuals with sd = 2
>> > y <- rep(0,1000)  # initialize y
>> > y[1] <- b*x[1] + r[ID[1]] + e[1]
>> > for (i in 2:1000) {
>> >   if (ID[i] == ID[i-1])
>> >     y[i] <- phi * y[i-1] + b*x[i] + r[ID[i]] + e[i]  # not the first
>> > observation for this ID
>> >   else y[i] <- b*x[i] + r[ID[i]] + e[i]  # first observation for this ID
>> > }
>> >
>> > data <- data.frame(y = y, x = x, ID = ID, time = rep(1:10, 100))
>> > head(data)
>> > p <- ggplot(data, aes(x=time, y=y, group=ID))
>> > p + geom_line()
>> >
>> > ### fit with AR1 effect
>> > library(nlme)
>> > fit_lmeA <- lme(y ~ x, data = data, random = ~ 1 | ID, cor = corAR1(form
>> > =
>> > ~1|ID))
>> > summary(fit_lmeA)
>> >
>> >
>> > Accordingly, the model I think of as AR1 does not seem to be the one
>> > nlme
>> > fits to the data. Does someone know what type of model nlme actually
>> > assumes in this situation?
>> >
>> > Thanks for your help!
>> >
>> > Best
>> >  Paul
>> >
>> >         [[alternative HTML version deleted]]
>> >
>> > _______________________________________________
>> > R-sig-mixed-models op r-project.org mailing list
>> > https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
>
>