[R] nested mixed-effect model: variance components

[Christoph Buser]

Dear Eric

Do you really have habitats nested within lagoons or are they
partially crossed (meaning that you have the same habitats in
different lagoons)?

If you have them perfectly nested, I think that you cannot
calculate both a fixed effect for habitats and a random effect
for lagoon (see the example below, lme and aov).

You can compare e.g. two lagoons by defining a contrast of the
habitats of one lagoon against the habitats of the other (if you
think that this is a meaningful test to interpret), but you
cannot estimate a random effect lagoon in presence of a nested
FIXED effect habitat. 

aov() will not return you the test and warn you about the
singular model. 

lme() will estimate a variance component for lagoon, but does
not provide you a test for the fixed factor.

Regards,

Christoph Buser



set.seed(1)
dat <- data.frame(y = rnorm(100), lagoon = factor(rep(1:4,each = 25)),
                  habitat = factor(rep(1:20, each = 5)))

summary(aov(y~habitat + Error(lagoon/habitat), data = dat))

library(nlme)
summary(lme(y~habitat, random = ~1|lagoon/habitat, data = dat))


--------------------------------------------------------------
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
http://stat.ethz.ch/~buser/
--------------------------------------------------------------



Eric Pante writes:
 > Dear listers,
 > 
 > I am trying to assess variance components for a nested, mixed-effects 
 > model. I think I got an answer that make sense from R, but I have a 
 > warning message and I wanted to check that what I am looking at is 
 > actually what I need:
 > 
 > my data are organized as transects within stations, stations within 
 > habitats, habitats within lagoons.
 > lagoons: random, habitats: fixed
 > the question is: how much variation is due to lagoons? habitats? 
 > lagoons*habitat? transects?
 > 
 > Here is my code:
 > 
 > res <- aov(COVER ~ HABITAT + Error(HABITAT+LAGOON+LAGOON/HABITAT), 
 > data=cov)
 > summary(res)
 > 
 > and I get Sum Sq for each to calculate variance components:
 > 
 > Error: STRATE
 >         Df Sum Sq Mean Sq
 > STRATE  5 4493.1   898.6
 > 
 > Error: ATOLL
 >            Df Sum Sq Mean Sq F value Pr(>F)
 > Residuals  5 3340.5   668.1
 > 
 > Error: STRATE:ATOLL
 >            Df  Sum Sq Mean Sq F value Pr(>F)
 > Residuals 18 2442.71  135.71
 > 
 > Error: Within
 >             Df Sum Sq Mean Sq F value Pr(>F)
 > Residuals 145 6422.0    44.3
 > 
 > My error message seems to come from the LAGOON/HABITAT, the Error is 
 > computed.
 > Warning message: Error() model is singular in: aov(COVER ~ HABITAT + 
 > Error(HABITAT+LAGOON+LAGOON/HABITAT), data=cov),
 > 
 > THANKS !!!
 > eric
 > 
 > 
 > 
 > Eric Pante
 > ----------------------------------------------------------------
 > College of Charleston, Grice Marine Laboratory
 > 205 Fort Johnson Road, Charleston SC 29412
 > Phone: 843-953-9190 (lab)  -9200 (main office)
 > ----------------------------------------------------------------
 > 
 > 	"On ne force pas la curiosite, on l'eveille ..."
 > 	Daniel Pennac
 > 
 > ______________________________________________
 > R-help at 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