[R] FW: Re: Doubt about nested aov output

[Douglas Bates]

On 9/7/05, John Wilkinson (pipex) <wilks at dial.pipex.com> wrote:
> Ronaldo,
> 
> Further to my previous posting  on your Glycogen nested aov model.
> 
> Having read Douglas Bates' response and  Reflected on his lmer analysis
> output of your aov nested model example as given.The Glycogen treatment has
> to be a Fixed Effect.If a 'treatment' isn't a Fixed Effect what is ? If
> Douglas Bates' lmer model is modified to treat Glycogen Treatment as a
> purely Fixed Effect,with Rat and the interaction Rat:Liver as random effects
> then--
> 
> > model.lmer<-lmer(Glycogen~Treatment+(1|Rat)+(1|Rat:Liver))
> > summary(model.lmer)
> Linear mixed-effects model fit by REML
> Formula: Glycogen ~ Treatment + (1 | Rat) + (1 | Rat:Liver)
>      AIC      BIC    logLik MLdeviance REMLdeviance
>  239.095 248.5961 -113.5475   238.5439      227.095
> Random effects:
>  Groups    Name        Variance   Std.Dev.
>  Rat:Liver (Intercept) 2.1238e-08 0.00014573
>  Rat       (Intercept) 2.0609e+01 4.53976242
>  Residual              4.2476e+01 6.51733769
> # of obs: 36, groups: Rat:Liver, 6; Rat, 2
> 
> Fixed effects:
>             Estimate Std. Error DF t value  Pr(>|t|)
> (Intercept) 140.5000     3.7208 33 37.7607 < 2.2e-16 ***
> Treatment2   10.5000     2.6607 33  3.9463 0.0003917 ***
> Treatment3   -5.3333     2.6607 33 -2.0045 0.0532798 .
> ---
> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> 
> Correlation of Fixed Effects:
>            (Intr) Trtmn2
> Treatment2 -0.358
> Treatment3 -0.358  0.500
> 
> > anova(model.lmer)
> Analysis of Variance Table
>           Df  Sum Sq Mean Sq   Denom F value    Pr(>F)
> Treatment  2 1557.56  778.78   33.00  18.335 4.419e-06 ***
>  --------------------------------------------------------
> which agrees with the aov model below.
> 
> >  model <- aov(Glycogen~Treatment+Error(Rat/Liver))
> > summary(model)
> 
> 
> John
> 
> -----Original Message-----
> From: John Wilkinson (pipex) [mailto:wilks at dial.pipex.com]
> Sent: 07 September 2005 12:04 PM
> To: "Ronaldo Reis-Jr."
> Cc: r-help
> Subject: Re: [R] Doubt about nested aov output
> 
> 
> Ronaldo ,
> 
> It looks as though you have specified your model incorrectly.
> 
> In the Rats example ,the Treatment is the only fixed effect,Rat and Liver
> are random effects
> 
> In aov testing for sig of 'Means' of Random Effects is pointless and that is
> why 'p' values are not given.Further more the interaction between a Random
> Effect and a Fixed Effect is also a Random Effect. The 'aov' with error
> structure terms output reflects this by only giving 'p' values to
> Fixed Effects and their interactions
> 
> >  model <- aov(Glycogen~Treatment+Error(Rat/Liver))
> > summary(model)
> 
> Error: Rat
>           Df Sum Sq Mean Sq F value Pr(>F) #Rat is random effect
> Residuals  1 413.44  413.44
> 
> Error: Rat:Liver                                      #Rat:Liver is Random effect
>           Df  Sum Sq Mean Sq F value Pr(>F)
> Residuals  4 164.444  41.111
> 
> Error: Within
>           Df  Sum Sq Mean Sq F value    Pr(>F)
> Treatment  2 1557.56  778.78  18.251 8.437e-06 *** #Fixed effect
> Residuals 28 1194.78   42.67
> ---
> Signif. codes:  0 '***' 0.001 '**' 0.01 '*' 0.05 '.' 0.1 ' ' 1
> 
> 
> 
> I hope that this is of help.
> 
> John

The difference between models like
  lmer(Glycogen~Treatment+(1|Rat)+(1|Rat:Liver))
and
  lmer(Glycogen~Treatment+(1|Treatment:Rat)+(1|Treatment:Rat:Liver))

is more about the meaning of the levels of "Rat" than about the
meaning of "Treatment".  As I understood it there are three different
rats labelled 1.  There is a rat 1 on treatment 1 and a rat 1 on
treatment 2 and a rat 1 on treatment 3.  Thus the levels of Rat do not
designate the "experimental unit", it is the levels of Treatment:Rat
that do this.