[R] 2-way Factorial with random factors

Eric Fuchs e.j.fuchs at gmail.com
Fri Oct 26 18:02:23 CEST 2007


  I am using R mainly on windows XP, version 2.5.  I´m a biologist,
with a medium level statistics background.  I have a problem stating a
two-way factorial design where both factors are random.  I´m using the
lmer() function implemented in the Matrix package version 0.99.

  My design is as follows:  Two species were randomly selected from a
total of 4 species.  This species are present in five different
populations, again selected at random.  Five trees per species were
selected randomly in each population and fruit number was estimated.
The main question is to estimate variance components for species,
population and their interaction.  Since both species are present in
each population, I figured they are not nested, but are crossed.
There is no fixed factor.

I conducted a normal aov() to estimate variance components manually
via EMS approximation:

fm1 <- aov(fruto~pop*spp, data=p1)

to test the significance of the variance components I did:

fm2 <- aov(fruto~pop+spp+Error(pop:spp), data=p1)

Manually, the variance component for the interaction is negative,
therefore I turned to lmer() to try to get a better estimate as well
as confidence intervals.  After invoking the Matrix library, I typed:

fm3 <- lmer(fruto~(1|pop) + (1|spp) + (1|pop:spp), data=p1)

  During calculations I get an error which indicates that the variance
components for pop and pop:spp are almost zero.  The summary output
presents variance components that do not resemble the manually
estimated values (eg. sigma_pop = 3650 vs sigma_pop(reml) = 7840).  If
I conduct the analysis in SAS, with GLM and the random option, I get
correct estimates of variance components.

  After reading a couple of help files, it looks like I need a fixed
component, which I don´t know how to code (dummy coding?).
 My specific questions are:

1.  Is my lmer formula correct?  If not, how should it be stated?
2.  Should I use another formula or command to estimate the
interaction component?
3.  Do I have to create a fixed factor in order to correctly use lmer()?

I know I should have consulted Pinhero´s book, but I´m currently in
Costa Rica, and the book is taking a long time to get here.   Thanks
in advance for all your help.  I´m new in the mailing list, and I hope
I´ve made myself clear enough.


Eric Fuchs
Universidad de Costa Rica

Eric J. Fuchs
lex parsimoniae:
entia non sunt multiplicanda praeter necessitatem

More information about the R-help mailing list