R help - Sep 2011 - How to Code Random Nested Variables within Two-way Fixed Model in lmer or lme

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]]

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.