[R] lmer applied to a wellknown (?) example
buser at stat.math.ethz.ch
Thu Aug 31 12:28:05 CEST 2006
There is an article in the R-News "Fitting linear mixed models
in R" in which you can find some examples for the syntax of
nested and non-nested design.
Hope this helps
Christoph Buser <buser at stat.math.ethz.ch>
Seminar fuer Statistik, LEO C13
ETH Zurich 8092 Zurich SWITZERLAND
phone: x-41-44-632-4673 fax: 632-1228
Henrik Parn writes:
> Dear all,
> During my pre-R era I tried (yes, tried) to understand mixed models by
> working through the 'rat example' in Sokal and Rohlfs Biometry (2000)
> 3ed p 288-292. The same example was later used by Crawley (2002) in his
> Statistical Computing p 363-373 and I have seen the same data being used
> elsewhere in the litterature.
> Because this example is so thoroughly described, I thought it would be
> very interesting to analyse it also using lmer and to see how the
> different approaches and outputs differs - from the more or less manual
> old-school (?) approach in Sokal, aov in Crawley, and to mixed models by
> In the example, three treatments (Treatment) with two rats (Rat) each
> (i.e six unique rats in total). Three liver preparations (Liver) are
> taken from each rat (i.e 18 unique liver preparations), and two glycogen
> readings (Glycogen) are taken from each liver preparation (36 readings).
> We want to test if treatments has affected the glycogen levels. The
> readings are nested in preparation and the preparations nested in rats.
> The data can be found here (or on p. 289 in Sokal):
> I was hoping to use the rat example as some kind of reference on my way
> to understand mixed models and using lmer. However, first I wish someone
> could check my suggested models!
> My suggestions:
> rats$Glycogen <- as.numeric(Glycogen)
> rats$Treatment <- as.factor(Treatment)
> rats$Rat <- as.factor(Rat)
> rats$Liver <- as.factor(Liver)
> model1 <- lmer(Glycogen ~ Treatment + (1|Liver) + (1|Rat), data=rats)
> Was that it?
> I also tried to make the 'liver-in-rat' nesting explicit (as suggested
> in 'Examples from...')
> model2 <- lmer(Glycogen ~ Treatment + (1|Rat:Liver) + (1|Rat), data=rats)
> but then the random effects differs from model1.
> Does the non-unique coding of rats and preparations in the original data
> set mess things up? Do I need to recode the ids to unique levels...
> rats$rat2 <- as.factor(rep(1:6, each=6))
> rats$liver2 <- as.factor(rep(1:18, each=2))
> ...and then:
> model3 <- lmer(Glycogen ~ Treatment + (1|liver2) + (1|rat2), data=rats)
> # or maybe
> model3 <- lmer(Glycogen ~ Treatment + (1|rat2:liver2) + (1|rat2), data=rats)
> Can anyone help me to get this right! Thanks in advance!
> Thanks to all contributors to lme/lmer topic on the list (yes, I have
> searched for 'lmer nested'...) and also the examples provided by John
> Fox' 'Linear mixed models, Appendix to An R and S-PLUS companion...' and
> Douglas Bates' 'Examples from Multilevel Software...' and R-news 5/1.
> Very helpful, but as usually I bet I missed something...Sorry.
> Henrik Pärn
> Department of Biology
> 7491 Trondheim
> +47 735 96282 (office)
> +47 909 89 255 (mobile)
> +47 735 96100 (fax)
> R-help at stat.math.ethz.ch mailing list
> PLEASE do read the posting guide http://www.R-project.org/posting-guide.html
> and provide commented, minimal, self-contained, reproducible code.
More information about the R-help