[R] How to Code Random Nested Variables within Two-way Fixed Model in lmer or lme

[Helios de Rosario]

Dear Dave, there are some inconsistencies in your explanation of the
problem. You said your variables are:

> CO is a continuous response variable,
> 
> Week is a fixed categorical factor,
> 
> Habitat is a fixed categorical factor, and
> 
> Location is a random categorical factor nested within Habitat.

What does this last statement mean? Are the Locations identified with
the same names in both Habitats (e.g. Location=1,2,3,... for
Habitat=Control, and then Location=1,2,3... for Habitat=Treatment,
although the Locations of both Habitats have nothing to do with each
other? Or do all 13 Locations receive different names?

If Location is really nested within Habitat, you would be in the former
case, and then your random terms should include the interaction
"Location:Habitat". In the latter case, the random term would just be
"Location".

But then, your model with aov is:

> mCO = aov(CO ~ Week * Habitat + Error(Location/Week))

Since you don't include Habitat in the Error term, I would say that
Location is not really nested within Habitat. But then, why is Week
nested within Location? Do you mean that the effect of Week may be
affected by the random Location?

Anyway, your interpretation of the ANOVA table is misleading:

> Error: Location
> 
>            Df Sum Sq Mean Sq F value   Pr(>F)   
> 
> Habitat     1 182566   182566   8.6519 0.01341 *
> 
> Residuals 11 232115    21101                   
> 
>   
> 
> Error: Location:Week
> 
>               Df Sum Sq Mean Sq F value     Pr(>F)     
> 
> Week          10 596431    59643 11.0534    7.5e-13 ***
> 
> Week:Habitat 10 196349    19635   3.6389 0.0003251 ***
> 
> Residuals    110 593551     5396               

This actually means:

For the F test of Habitat, the denominator MS is that for Location
 
For the F test of Week, the denominator MS is that for the Location x
Week
 
For the F test of Habitat x Week, the denominator MS is that for
Location x Week


And then, you wrote your attempt with lmer:

> m. = lmer(CO ~ Week * Habitat + (1|Habitat/Location))

The random term here (1|Habitat/Location) has nothing to do with the
Error term you used in aov.

If location is really nested within Habitat, perhaps you meant

m. = lmer(CO ~ Week * Habitat + (1|Habitat:Location)) 

(Habitat/Location means that Habitat has a random effect per se as
well, and I guess you don't mean that!)
Or if Location is not really nested,

m. = lmer(CO ~ Week * Habitat + (1|Location))

or if you really wanted the same model as with aov:

m. = lmer(CO ~ Week * Habitat + (1|Location/Week)) 

Please clarify your model. Otherwise it would be impossible to make any
comparison.

Helios





>>> El día 29/09/2011 a las 6:30, Dave Robichaud <drobichaud at lgl.com>
escribió:
> Hi All,
> 
> I am frustrated by mixed-effects model! I have searched the web for 
> hours, and found lots on the nested anova, but nothing useful on my 
> specific case, which is: a random factor (C) is nested within one of
the 
> fixed-factors (A), and a second fixed factor (B) is crossed with the

> first fixed factor:
> 
> C/A
> 
> A
> 
> B
> 
> A x B
> 
> My question: I have a functioning model using the aov command (see 
> below), and I would now would like to recode it, using a more
flexible 
> command such as lme or lmer. Once I have the equivalent syntax down,
I 
> would ideally like to re-run my analysis using "family = poisson", as
CO 
> is actually count data.
> 
> I have a dataset including a response variable CO, measured once per

> Week (for 11 weeks) at 13 Locations. The 13 Locations are divided
into 2 
> habitat types (Control and Treatment).
> 
> Thus:
> 
> CO is a continuous response variable,
> 
> Week is a fixed categorical factor,
> 
> Habitat is a fixed categorical factor, and
> 
> Location is a random categorical factor nested within Habitat.
> 
> Here is my model in R:
> 
> mCO = aov(CO ~ Week * Habitat + Error(Location/Week))
> 
> summary(mCO)
> 
> And the output:
> 
> Error: Location
> 
>            Df Sum Sq Mean Sq F value   Pr(>F)   
> 
> Habitat     1 182566   182566   8.6519 0.01341 *
> 
> Residuals 11 232115    21101                   
> 
>   
> 
> Error: Location:Week
> 
>               Df Sum Sq Mean Sq F value     Pr(>F)     
> 
> Week          10 596431    59643 11.0534    7.5e-13 ***
> 
> Week:Habitat 10 196349    19635   3.6389 0.0003251 ***
> 
> Residuals    110 593551     5396                       
> 
> Given that this is a mixed model, I believe the appropriate error
terms 
> are as follows:
> 
> For the F test of Habitat, the denominator MS is that for
location/habitat;
> 
> For the F test of Week, the denominator MS is the residual; and
> 
> For the F test of Habitat x Week, the denominator MS is the
residual.
> 
> My tinkering with lmer and lme have not produced results similar to
the 
> above
> 
> For example,
> 
> m. = lmer(CO ~ Week * Habitat + (1|Habitat/Location))
> 
> anova(m.)
> 
> produces:
> 
> Analysis of Variance Table
> 
>                Df Sum Sq Mean Sq F value
> 
> Week           10 596431    59643 11.0534
> 
> Habitat        1   28652    28652   5.3100
> 
> Week:Habitat 10 196349    19635   3.6389
> 
> Any coding advice would be greatly appreciated!
> 
> Thanks for your consideration,
> 
> Dave Robichaud
> 
> 
> 	[[alternative HTML version deleted]]

INSTITUTO DE BIOMECÁNICA DE VALENCIA
Universidad Politécnica de Valencia • Edificio 9C
Camino de Vera s/n • 46022 VALENCIA (ESPAÑA)
Tel. +34 96 387 91 60 • Fax +34 96 387 91 69
www.ibv.org

  Antes de imprimir este e-mail piense bien si es necesario hacerlo.
En cumplimiento de la Ley Orgánica 15/1999 reguladora de la Protección
de Datos de Carácter Personal, le informamos de que el presente mensaje
contiene información confidencial, siendo para uso exclusivo del
destinatario arriba indicado. En caso de no ser usted el destinatario
del mismo le informamos que su recepción no le autoriza a su divulgación
o reproducción por cualquier medio, debiendo destruirlo de inmediato,
rogándole lo notifique al remitente.