[R-sig-ME] exemples, tutorial on lmer
Roberts, J. Kyle
jkrobert at bcm.tmc.edu
Tue May 29 18:45:56 CEST 2007
This is very short, but it is an excerpt from a chapter that I just wrote about software options for multilevel analysis. Hope it helps,
Those who are familiar with the R programming language probably know it as the free open-sourced version of S-PLUS. In fact, both packages are based on the S programming language and use commands that are remarkably similar in both form and execution. R was originally written by Gentleman
and Ihaka at the University of Auckland with the first version released in 2000. The original goal of R was to make a "leaner" S language-based system which did not have the same memory hampering problems that S-PLUS had in its early stages of development. As R has evolved, its strength has come
in that it is a powerful free software package with a network of users who write code which is openly shared among other users. Although S-PLUS supports more of a GUI (Graphical User Interface) design, recent programs in R such as R-Commander developed by John Fox (Fox, 2005) give it the same feel
as the S-PLUS interface.
As the code for running a linear multilevel model in R is exactly the same as running the model in S-PLUS, it does not bear repeating here. The nlme package run in R uses the same set of commands used by the nlme library in S-PLUS (available at http://cran.r-project.org/). There are some
subtle differences, however. First, R does not allow the user to run the plot(augPred(lme.object)) command to produce the graph seen in Figure 6. This graph can still be produced, but it must be done with a groupedData object. Second, the commands for running more complicated models like binary
response models are different in R. In S-PLUS, the correlatedData library is used to fit binary response models, whereas the lme4 package in R which is used to fit these models.
Recently, new routines have been introduced to R that also allow for the fitting of multilevel models. In addition to the nlme library, the Matrix package contains the routine "lmer" which is also capable of running these models. This library is very similar to nlme in terms of its
functionality, but the main benefit that it has is that the computations are based on analysis of sparse matrices. This allows lmer to analyze data at much faster rates than previous versions of nlme (Bates, 2006a).
The lmer command syntax is slightly different than the nlme syntax. In lmer the command line takes the form:
lmer(formula, data, family, method, control , start, subset, weights, na.action, offset, model, x, y, ...)
A thorough description of all of the arguments in this command line can be obtained by typing ?lmer at any command prompt in R. The format of the command line is mostly similar to lme, with the exception of the fact that lme splits the commands for the fixed and random parts of the command, where
lmer simply has a single "formula" command line. Part of the reasoning behind this change is because multilevel models are sometimes better thought of as single linear (in this case) equation in stead of a series of two equations (or levels) estimated separately (D. Bates, personal communication,
April 6, 2006). Therefore, the syntax necessary to run the Roberts (2004) data would be:
ex1<-lmer(SCIENCE ~ URBAN + (URBAN|GROUP), example.data)
The above command line may be read as the variable "SCIENCE" being estimated by "urban" plus the random effects for both "URBAN" and for the intercept defined by the level-2 grouping variable "GROUP". Just as with S-PLUS, if a three level model were to be specified with level-2 and level-3 grouping
structures of classrooms and schools, respectively (these are not real variables in the dataset, but meant as heuristic), the command line for this model would read:
ex1<-lmer(SCIENCE ~ URBAN + (URBAN|SCHOOL/CLASSROOM), example.data)
Typing the command summary(ex1) will produce the following:
Linear mixed-effects model fit by REML
Formula: SCIENCE ~ URBAN + (URBAN | GROUP)
AIC BIC logLik MLdeviance REMLdeviance
424.1713 442.6223 -206.0857 413.2216 412.1713
Groups Name Variance Std.Dev. Corr
GROUP (Intercept) 113.60403 10.65852
URBAN 0.25200 0.50200 -0.625
Residual 0.27066 0.52025
# of obs: 160, groups: GROUP, 16
Estimate Std. Error DF t value Pr(>|t|)
(Intercept) 22.39124 2.71703 158 8.2411 6.152e-14 ***
URBAN -0.86700 0.12981 158 -6.6790 3.884e-10 ***
Signif. codes: 0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
Correlation of Fixed Effects:
This output is also remarkably similar to the lme output. One notable difference is that the df has changed from previous versions of the package (see the S-PLUS output from above). The reasoning behind this df change is a discussion that could easily take up another volume, and a problem
that probably won't be resolved in the near future. For an interesting discussion of the problems related to df and p-values in multilevel models, see Bates (2006b).
The only disadvantage to using R, as opposed to S-PLUS, is that it doesn't support a graphical user interface (GUI) design. However, even this drawback has been partially overcome with the development of the R-Commander program by Fox (2005). Regardless, R is an extremely powerful package
using object oriented language. Couple this with the fact that the software is free and R makes a very inviting package. The one caveat is that this package, like S-PLUS, has a steep learning curve. Learning new code is like learning a new language; but the payoff may well be worth the effort.
J. Kyle Roberts, Ph.D.
Baylor College of Medicine
Center for Educational Outreach
One Baylor Plaza, MS: BCM411
Houston, TX 77030-3411
713-798-6672 - 713-798-8201 Fax
jkrobert at bcm.edu
From: r-sig-mixed-models-bounces at r-project.org [mailto:r-sig-mixed-models-bounces at r-project.org] On Behalf Of Olivier MARTIN
Sent: Tuesday, May 29, 2007 9:36 AM
To: r-help at stat.math.ethz.ch; r-sig-mixed-models at r-project.org
Subject: [R-sig-ME] exemples, tutorial on lmer
I have some difficulties to work with the function lmer from lme4 Does somebody have a tutorial or different examples to use this function?
R-sig-mixed-models at r-project.org mailing list https://stat.ethz.ch/mailman/listinfo/r-sig-mixed-models
More information about the R-sig-mixed-models