Langbehn, Douglas
2012-Mar-09 00:36 UTC
[R] pdMat class in LME to mimic SAS proc mixed group option? Group-specific random slopes
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-langbehn at 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.
