# Fitting multivariate linear mixed model in R

I'm wondering how to fit multivariate linear mixed model and finding multivariate BLUP in R. I'd appreciate if someone come up with example and R code.

Edit

I wonder how to fit multivariate linear mixed model with lme4. I fitted univariate linear mixed models with the following code:

library(lme4)
lmer.m1 <- lmer(Y1 ~ A*B + (1|Block) + (1|Block:A), data=Data)
summary(lmer.m1)
anova(lmer.m1)

lmer.m2 <- lmer(Y2 ~ A*B + (1|Block) + (1|Block:A), data=Data)
summary(lmer.m2)
anova(lmer.m2)


I'd like to know how to fit multivariate linear mixed model with lme4. The data is below:

Block A B    Y1    Y2
1 1 1 135.8 121.6
1 1 2 149.4 142.5
1 1 3 155.4 145.0
1 2 1 105.9 106.6
1 2 2 112.9 119.2
1 2 3 121.6 126.7
2 1 1 121.9 133.5
2 1 2 136.5 146.1
2 1 3 145.8 154.0
2 2 1 102.1 116.0
2 2 2 112.0 121.3
2 2 3 114.6 137.3
3 1 1 133.4 132.4
3 1 2 139.1 141.8
3 1 3 157.3 156.1
3 2 1 101.2  89.0
3 2 2 109.8 104.6
3 2 3 111.0 107.7
4 1 1 124.9 133.4
4 1 2 140.3 147.7
4 1 3 147.1 157.7
4 2 1 110.5  99.1
4 2 2 117.7 100.9
4 2 3 129.5 116.2

• Unfortunately, your question cannot motivate any decent answer as it stands. You might consider adding information about your design (crossed effects--random or not--vs. nested ones), what you qualify as multivariate (is this on the response variable(s) or the number of predictors), etc.
– chl
May 6, 2011 at 19:37
• I've updated my question. Oct 21, 2011 at 0:19

Fitting multivariate models with lme4 or nlmeis a bit fiddly, but solutions can be found in this document by Ben Bolker.

Else if you want to stay in a frequentist framework, the mcglm package can handle multivariate models, even with non-normal distributions. A detailed tutorial should be published soon. If you are not familiar with design matrices, designing the matrix of random effects can be a bit tricky though.

In a Bayesian framework, the MCMCglmm package is also very good at modelling multivariate traits incl. non-normal distributions, and handles random effects in a simpler way than mcglm does. Its use for multivariate models is also well described by Ben Bolker. But you first have to make yourself familiar with Markov chain Monte Carlo, as well as with the principles of Bayesian statistics.

So either way, there might be a slow learning curve at the beginning, depending on your familiarity with either methods!

Try the R package nlme

You can find some examples, theory and further documentation in: http://cran.r-project.org/doc/contrib/Fox-Companion/appendix-mixed-models.pdf

The nlme package is able to calculate pooled estimates [or the so called BLUP= best linear unbiased predictor].

Once you've downloaded the package, type in R console: help(predict.lme)

For more information, look at page 17 in Fox's paper. There you can find an example on how to pool information across subjects.

Hope this helps :)

• These are univariate BLUP. But I'm looking for multivariate BLUP. May 7, 2011 at 16:32