I would like to be able to use lme to fit random effect models In which some
but not all of the random effects are constrained to be independent. It seems
as thought the pdMat options in lme are a promising avenue. However, none of
the existing pdMat classes seem to allow what I want.
As a specific example, I would like to fit a random intercept/slope mixed model
to longitudinal observations in which the subjects are divided into 5
subgroups. There is good reason to suspect that the intercept/slope variance
will be heteroscedastic among the 5 subgroups (as well as the residual R-side
variances).
If I'm willing to assume that each subject's random slope and intercepts is
independent (a dangerous assumption), then I can get close to the model using
the following code:
#R CODE MODEL1
G.Data<-groupedData(y~1|Subject,data=OriginalData)
G.Data$subgroup<-as.factor(G.Data$subgroup)
pd<-pdDiag(~subgroup+subgroup:time-1)
# subgroup is the subgroup-specific intercept and subgroup:time is the
subgroup-specific intercept.
model1<-lme(fixed=y~x*time,
random=pd,
weights=varIdent(form=~1|subgroup),
na.action=na.omit,
data=G.Data)
# END R CODE
Unfortunately, the pdDiag command constrains the random intercepts and slopes
among all groups to be independent.
The default class for pdMat, pdSymm, assumes that ALL random intercept and
slope parameters are correlated. The reality is that the intercepts and slopes
should be correlated within each subgroup, but these random parameters are
independent between the subgroups. Use of the pdSymmm class leads to parameter
non-identifiability.
What I need is something like the 'value' argument for pdSymm, by which one can
specify initial covariance parameter estimates for the random effects. If there
were a way to specify that initial covariance values of 0 are to be fixed at 0,
then the problem would be solved (I think). For example, with 2 subgroups, the
"value" matrix would look like
1 0 r 0
0 1 0 r
r 0 1 0
0 r 0 1
Where the corresponding row/column order is
Intercept,1 intercept2, slope 1, slope 2.
I've been digging through Pinheiro and Bates, as well as these help forums and
files, but couldn't find info on writing or modifying pdMat classes.
In SAS proc mixed, such a model is specified fairly easily by using the group
option:
/* SAS CODE FRAGMENT WITHIN PROC MIXED */
random int visit/subject=subject group=subgroup type=un;
repeated /subject=subject group=subgroup;
/*Random statement, with group and type options, accomplishes what I want to do
in R
/* The repeat statement, with its own group option, allows the r-side variance
to also differ by subgroup.
/* This is accomplished in lme by using the
/* weights=varIdent(form=~1|subgroup)
/* argument.
Doug Langbehn MD, PhD
Professor
Dept. of Psychiatry
University of Iowa
Iowa City, Iowa
douglas-langb...@uiowa.edu
________________________________
Notice: This UI Health Care e-mail (including attachments) is covered by the
Electronic Communications Privacy Act, 18 U.S.C. 2510-2521, is confidential and
may be legally privileged. If you are not the intended recipient, you are
hereby notified that any retention, dissemination, distribution, or copying of
this communication is strictly prohibited. Please reply to the sender that you
have received the message in error, then delete it. Thank you.
______________________________________________
R-help@r-project.org mailing list
https://stat.ethz.ch/mailman/listinfo/r-help
PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
and provide commented, minimal, self-contained, reproducible code.
