Dear Henrik 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.
http://cran.r-project.org/doc/Rnews/Rnews_2005-1.pdf Hope this helps Christoph -------------------------------------------------------------- Christoph Buser <[EMAIL PROTECTED]> Seminar fuer Statistik, LEO C13 ETH Zurich 8092 Zurich SWITZERLAND phone: x-41-44-632-4673 fax: 632-1228 http://stat.ethz.ch/~buser/ -------------------------------------------------------------- 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 > lmer. > > 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): > http://www.bio.ic.ac.uk/research/mjcraw/statcomp/data/rats.txt > // > 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: > > attach(rats) > rats$Glycogen <- as.numeric(Glycogen) > rats$Treatment <- as.factor(Treatment) > rats$Rat <- as.factor(Rat) > rats$Liver <- as.factor(Liver) > str(rats) > > model1 <- lmer(Glycogen ~ Treatment + (1|Liver) + (1|Rat), data=rats) > summary(model1) > > 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) > summary(model2) > > 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)) > str(rats) > > ...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! > > P.S. > 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. > > Regards, > > Henrik > > -- > ************************ > Henrik Pärn > Department of Biology > NTNU > 7491 Trondheim > Norway > > +47 735 96282 (office) > +47 909 89 255 (mobile) > +47 735 96100 (fax) > > ______________________________________________ > [email protected] 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. ______________________________________________ [email protected] 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.
