[R] Translating lme code into lmer was: Mixed effect model in R

From: Doran, Harold <HDoran_at_air.org>
Date: Fri 20 Oct 2006 - 12:22:55 GMT


This question comes up periodically, probably enough to give it a proper thread and maybe point to this thread for reference (similar to the 'conservative anova' thread not too long ago).

Moving from lme syntax, which is the function found in the nlme package, to lmer syntax (found in lme4) is not too difficult. It is probably useful to first explain what the differences are between the two functions for some background.

The lme function in the nlme package was originally designed to handle models with nested random effects. For example, students nested in a single school or employees nested in a single organization. It was possible to fit models with crossed random effects using lme, but it was inefficient and the algorithms were not optimized for such situations.

The lmer function, on the other hand, can handle models with crossed (or partially crossed) random effects as well as models with nested random effects. Doug Bates has also significantly improved the computational algorithms such that the speed by which very large data problems can be resolved is impressively fast.

A model fit using lmer is typically of the form

lmer(response ~ covariate(s) + (covariates|group), data)

Notice that the random terms in lmer are always of the form (a|b) denoting that a is a random effect conditional on the grouping factor b. Also notice that random effects terms must be enclosed within parentheses. In lme the random term is of the form

random = ~ covariate | group

To illustrate the difference, here is an example using the sleepstudy data where the first is the lmer model and the second is the lme model. Parameter estimates from both models are the same.

fm1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy, method='ML') fm2 <- lme(Reaction ~ Days, random=~Days|Subject, sleepstudy, method='ML')

Now, assume you want to extend this to include additional levels of random variation. The egsingle data set in the mlmRev package has repeated measures on students and those students are nested in schools. We can fit a model with random intercepts and slopes for students and students within schools. The first model is fit using lme and the second is fit using lmer.

fm3 <- lme(math ~ year, random=~year|schoolid/childid, egsingle) fm4 <- lmer(math ~ year +(year|schoolid:childid) +(year|schoolid), egsingle, control=list(gradient = FALSE, niterEM = 0))

Both result in parameter estimates that are exactly the same. The newest release of lme4 allows for an even easier transition between lme and lmer. Model syntax in lmer can be specified as:

fm5 <- lmer(math ~ year +(year|schoolid/childid), egsingle, control=list(gradient = FALSE, niterEM = 0))

Notice how the grouping structure is used in this case, which resembles the lme syntax for random terms.

The reason that the control statement is used is because Doug Bates told me to. I hope he pipes in and explains it, but from what I recall, the way he wrote the program has some fancy math. As it turns out, it is better to skip some of that math by turning off the EM iterations.

Just for documention, I mentioned the CPU time is much faster using lmer than lme.

> system.time(fm3 <- lme(math~year, random=~year|schoolid/childid,
egsingle))

[1] 63.08 0.03 63.69 NA NA

> system.time(fm4 <- lmer(math ~ year +(year|schoolid:childid)
+(year|schoolid), egsingle, control=list(gradient = FALSE, niterEM = 0)))
[1] 2.22 0.00 2.22 NA NA

On my Windows machine, the first model fit using lme took just over a minute whereas the same model using lmer was estimated in 2.22 seconds.

There are other examples of how to use the lmer function in the mlmRev package found using

vignette("MlmSoftRev")

This is anything but comprehenssive, but I hope others can chime in and share their experiences, or correct any errors I may have made above.

Best,

Harold

> -----Original Message-----
> From: r-help-bounces@stat.math.ethz.ch
> [mailto:r-help-bounces@stat.math.ethz.ch] On Behalf Of Lina Jansen
> Sent: Friday, October 20, 2006 4:58 AM
> To: R-help@stat.math.ethz.ch
> Subject: Re: [R] Mixed effect model in R
>
> Thanks for the helping links. Now, I worked out that I have
> to use the lme4 package (with the lmer function) for my
> analysis. But now I do not understand the input to the lmer function.
>
> In the lme function (of the nlme package) the correct input
> would in my case
> be:
>
> lme(fixed=Ac_LC~cond_ind,random=~img_cond|sub_ind/cond_ind)
>
> but also after reading the help and the R news I do not
> understand the formula I have to use for the lmer function.
> Could someone help me translating the lme input to a lmer
> input? And does someone know of a good explanation of the
> kinds of formulas you can input? In the books they only
> explain the lme input.
>
> Thanks,
> Lina
>
> 2006/10/17, Stefan Grosse <singularitaet@gmx.net>:
> >
> > Please always reply to the list as well as there always might be
> > someone faster/better answering (or it could be that I am wrong, so
> > someone might correct me)
> >
> > Indeed Pinheiro/Bates assume gaussian error terms... but I am not
> > really sure whether you meant that with "non normally distributed
> > respond variable" resp. "with non-normal data"
> >
> > however:
> > "/ Mixed-effects models: / The recommended nlme
> > <http://cran.r-project.org/src/contrib/Descriptions/nlme.html>
> > package, associated with Pinheiro and Bates, /
> Mixed-Effects Models in
> > S and S-PLUS / (Springer, 2000), fits linear and nonlinear
> > mixed-effects models, commonly used in the social sciences for
> > hierarchical and longitudinal data. Generalized linear
> mixed-effects
> > models may be fit by the glmmPQL function in the MASS
> package, and by
> > the lmer function in the Matrix
> > <http://cran.r-project.org/src/contrib/Descriptions/Matrix.html>
> > package (related to the lme4
> > <http://cran.r-project.org/src/contrib/Descriptions/lme4.html>
> > package, which largely supersedes nlme
> >
> <http://cran.r-project.org/src/contrib/Descriptions/nlme.html> for /
> > linear / mixed models). Also see the lmeSplines
> >
> <http://cran.r-project.org/src/contrib/Descriptions/lmeSplines.html>
> > and lmm
> <http://cran.r-project.org/src/contrib/Descriptions/lmm.html>
> > packages." [
> > http://cran.r-project.org/src/contrib/Views/SocialSciences.html ]
> >
> > Lina Jansen schrieb:
> > >
> > >
> > > 2006/10/17, Stefan Grosse <singularitaet@gmx.net
> > > <mailto:singularitaet@gmx.net>>:
> > >
> > > Interesting packages for you might be the nlme and
> lme4 packages
> > > and as
> > > a book Pinheiro/Bates, "Mixed-Effects Models in S and S-Plus"
> > >
> > >
> > > Thank you for the answer. I am always unsure concerning the
> > > non-normality. Can I use the nlme and lme4 with non-normal data?
> > > First, I thought they would work like an ANOVA but with
> random and
> > > fixed effects.
> >
> >
>
> [[alternative HTML version deleted]]
>
> ______________________________________________
> R-help@stat.math.ethz.ch 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.
>



R-help@stat.math.ethz.ch 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. Received on Fri Oct 20 22:44:06 2006

Archive maintained by Robert King, hosted by the discipline of statistics at the University of Newcastle, Australia.
Archive generated by hypermail 2.1.8, at Fri 20 Oct 2006 - 19:30:11 GMT.

Mailing list information is available at https://stat.ethz.ch/mailman/listinfo/r-help. Please read the posting guide before posting to the list.