[prev in list] [next in list] [prev in thread] [next in thread]
List: r-sig-mixed-models
Subject: [R-sig-ME] Res: Random effects, autocorrelation and nesting in lme
From: walmes zeviani <walmeszeviani () yahoo ! com ! br>
Date: 2010-02-16 19:04:42
Message-ID: 293066.22516.qm () web32105 ! mail ! mud ! yahoo ! com
[Download RAW message or body]
Anders,
I think that the correct sintax or model specification is different. I think you \
should use form=~1|Person instead of form=~1|Person_replicate at correlation \
argument, this because repeated measures on the same person along time induces \
correlation. More about it in Pinheiro & Bates (2000). The code is the following:
> Data <- data.frame(Person=factor(rep(1:3,each=12)),
+ Person_replicate=factor(rep(1:12,each=3)),
+ Time=rep(1:4,9),
+ response=round(5+rnorm(36),2))
>
> lme(response~Time, random=list(~Time|Person),
+ correlation=corGaus(form=~1|Person), data=Data)
Linear mixed-effects model fit by REML
Data: Data
Log-restricted-likelihood: -41.12866
Fixed: response ~ Time
(Intercept) Time
5.3595834 -0.2110061
Random effects:
Formula: ~Time | Person
Structure: General positive-definite, Log-Cholesky parametrization
StdDev Corr
(Intercept) 2.270560e-05 (Intr)
Time 1.236894e-08 0
Residual 7.507016e-01
Correlation Structure: Gaussian spatial correlation
Formula: ~1 | Person
Parameter estimate(s):
range
0.8230017
Number of Observations: 36
Number of Groups: 3
>
Walmes Zeviani. Brasil.
________________________________
De: Anders Christian Jensen <m03acj@math.ku.dk>
Para: r-sig-mixed-models@r-project.org
Enviadas: Terça-feira, 16 de Fevereiro de 2010 13:57:15
Assunto: [R-sig-ME] Random effects, autocorrelation and nesting in lme
Dear R - users
I have a problem with lme when trying to specify both "random" and
"correlation" arguments to my model.
The problem seems to be related to the nesting structure of my data:
For each person (p) I have multiple replicates (r) of a time series of
observations. Observations within a time series is indexed by t. That is,
my outcome (y) is uniquely determined by indices y_{prt}.
I would like to specify a model that includes
1) a random slope and intercept for each person, and
2) some sort of correlation on the time series, within each
(person,replicate); possibly an AR(1) process.
In other words, I would like to specify the following model (disregarding
the fixed effects part of the model):
y_{prt} = a_p + b_p * t + e_{prt},
where a_p and b_p are the (Gaussian) random slope and intercept on person
level, and e_{prt} should be such that
cor(e_{prt},e_{p'r't'}) = f(|t-t'|), for some function f
when p=p' and r=r', and 0 otherwise.
When I try to fit the model with lme I get an error message:
"Incompatible formulas for groups in "random" and "correlation""
A very small toy data set (without any autocorrelation) illustrates my
problem:
> Data <- data.frame(Person = factor(rep(1:3,each=12)),
+ Person_replicate = factor(rep(1:12,each=3)),
+ Time = rep(1:4,9),
+ response = round(5+rnorm(36),2))
> str(Data)
'data.frame': 36 obs. of 4 variables:
$ Person : Factor w/ 3 levels "1","2","3": 1 1 1 1 1 1 1 1 1 1 ...
$ Person_replicate: Factor w/ 12 levels "1","2","3","4",..: 1 1 1 2 2 2 3
3 3 4 ...
$ Time : int 1 2 3 4 1 2 3 4 1 2 ...
$ response : num 5.21 4.95 3.32 4.86 6.18 5.68 5.14 3.81 6.17
5.08 ...
> head(Data)
Person Person_replicate Time response
1 1 1 1 5.21
2 1 1 2 4.95
3 1 1 3 3.32
4 1 2 4 4.86
5 1 2 1 6.18
6 1 2 2 5.68
I then try to fit the model, and get an error:
> library(splines)
> library(nlme)
> lme(response ~ Time,
+ random = list(~Time|Person),
+ correlation=corGaus(form=~1|Person_replicate),
+ data = Data)
Error in lme.formula(response ~ Time, random = list(~Time | Person),
correlation = corGaus(form = ~1 | :
Incompatible formulas for groups in "random" and "correlation"
Basically I would like the "random"-argument for lme to work on person
level, but the "correlation"-argument should work on (person,replicate)
level. I found old posts on the mailing list with similar problems but
none of them include an answer.
Any help would be much appreciated! Thanks
_______________________________________________
R-sig-mixed-models@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
____________________________________________________________________________________
[[elided Yahoo spam]]
[[alternative HTML version deleted]]
_______________________________________________
R-sig-mixed-models@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
[prev in list] [next in list] [prev in thread] [next in thread]
Configure |
About |
News |
Add a list |
Sponsored by KoreLogic