[R-sig-ME] what kind of analysis should be used when there are several correlated dependent variables?
datkins at u.washington.edu
Fri Dec 18 19:51:30 CET 2009
One route would be to fit a multivariate mixed-effects model, in which
the multiple outcomes are jointly treated as a DV. Neither lme() nor
lmer() have a way of explicitly specifying a multivariate model (ie,
MCMCglmm can "cbind" multiple outcomes together), but it is possible to
fit at least some basic multivariate models with some data manipulation.
The main "trick" here is to reshape your data by stacking all DVs into a
single DV column and creating a factor that identifies which values in
the DV correspond to your 4 separate DVs.
I recently used lme() to do this with a study of daily positive mood,
negative mood, and fatigue in women with breast cancer (each woman had
30 days of data).
Here was the code to reshape the data (where "id" identifies the
individual, and "day" identifies -- shockingly -- day):
tmp.long.df <- reshape(my.df,
varying = c("posemo","negemo","fatigue"),
v.names = "dv",
times = 1:3,
names(tmp.long.df) <- "dv.f" # rename times
tmp.long.df$dv.f <- factor(tmp.long.df$dv.f, 1:3,
So, your separate columns of DVs are specified on "varying", and the new
data.frame has a column called "times" that identifies which DV is which
in your new "dv" column. I change times (which just happened to be the
8th column in my data, you'll likely need to change that) to "dv.f" and
then change that to factor with levels identifying separate DVs.
Then, here is an example of a call to lme() that I used:
mult.lme3.2 <- lme(dv ~ -1 + dv.f/(supp.sat.c*(emot.iss.c + phys.iss.c +
data = mult.df, na.action = na.omit, control = list(msVerbose = TRUE),
random = list(id = ~ -1 + dv.f),
weights = varIdent(form = ~ 1 | dv.f),
correlation = corAR1(form = ~ 1 | id/dv.f))
So, a couple comments:
-- With this set-up, getting rid of the intercept makes the models a bit
easier to interpret (for me) as you get separate intercepts for each DV.
-- Similarly, the "dv.f/..." fits the covariates nested within each DV,
which again makes things a bit easier to interpret.
-- The random statement fits separate random-intercepts by DV that are
-- The weights statement allows different residual error terms by DV.
-- For our data, we fit an AR1 correlation model because of 30 days of data.
Hope that helps; I think David Afshartous has some postings about
multivariate models in lme which you can probably find in the archive
(if memory serves).
[And, I am again impressed at how flexible the mixed-effects tools are
in R -- kudos to Jose and Doug for the great software.]
Dave Atkins, PhD
Research Associate Professor
Center for the Study of Health and Risk Behaviors
Department of Psychiatry and Behavioral Science
University of Washington
1100 NE 45th Street, Suite 300
Seattle, WA 98105
datkins at u.washington.edu
More information about the R-sig-mixed-models