[R] ANOVA in Randomized-complete blocks design

[Simon Blomberg]

January Weiner wrote:
> Ops, sorry.
>
> OK, as you see, I am going through S&R and doing the examples using R.
> I have a further question, on the nested ANOVA example from Box 10.1
> (mosquito female wing lengths).  I made sure that the data is correct
> :-) (data below).
>
> I am not sure how to create in R a Model II nested ANOVA.
> aov( wing ~ cage / female ) which is, I believe, a fixed effect nested
> ANOVA where females are nested within cages produced the correct sum
> of squares, but computates a different F value for the cage effect:
>
>             Df  Sum Sq Mean Sq F value    Pr(>F)
> cage         2  665.68  332.84  255.70 1.452e-10 ***
> cage:female  9 1720.68  191.19  146.88 6.981e-11 ***
> Residuals   12   15.62    1.30
>
> (whereas in S&R, among cages F = 1.741; the rest is same).  The F
> value for the cage:female effect is the same as in S&R.  Why do I get
> a higly significant cage effect?
>
> In S&R, the significance of the cage effect is done by F = MSamong /
> MSsubgr = 332.84 / 191.19 = 1.74. In the R model, F = 255.70, and I do
> not understand where this value comes from (at first I thought that it
> is MSamong / MSerror = 332.84 / 1.3 = 256.0308).
>
> I am confused...
>   
aov(wing ~ cage/female) gives you a test of cage against the residual 
variance and cage:female against the residual variance. As you noticed, 
the residual variance is the wrong error stratum for the cage effect. To 
get the correct test of the cage effect, try:

summary(aov(wing ~ cage + Error(female)))

You will get the correct F test from that.

To get the variance components, it is easier to use the lme function in 
the nlme package:

fit <- lme(wing ~ cage, random=~1|female)
summary(fit)
where cage is fixed and female is random.

Compare:

fit <- lme(wing ~ 1, random=~1|cage/female)
summary(fit)

This treats both cage and female as random effects.

VarCorr(fit)

will give you the variance components, or you can read off the standard 
deviations (sqrt of variance components) from summary(fit). Yet another 
way to analyse the problem is with the lmer function (package lme4).

Cheers,

Simon.

> How should I modify the model?
>
> Another question: is there a way to automatically estimate the
> variance components, or do I have to take the respective MS' and
> calculate it myself?
>
> Thanks,
> January
>
> OK, here is the data:
>
> cage female wing
> cage1 f1 58.5
> cage1 f1 59.5
> cage1 f2 77.8
> cage1 f2 80.9
> cage1 f3 84.0
> cage1 f3 83.6
> cage1 f4 70.1
> cage1 f4 68.3
> cage2 f1 69.8
> cage2 f1 69.8
> cage2 f2 56.0
> cage2 f2 54.5
> cage2 f3 50.7
> cage2 f3 49.3
> cage2 f4 63.8
> cage2 f4 65.8
> cage3 f1 56.6
> cage3 f1 57.5
> cage3 f2 77.8
> cage3 f2 79.2
> cage3 f3 69.9
> cage3 f3 69.2
> cage3 f4 62.1
> cage3 f4 64.5
>
> Cheers,
> January
>
> --
>
> ______________________________________________
> 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
> and provide commented, minimal, self-contained, reproducible code.
>
>   


-- 
Simon Blomberg, B.Sc.(Hons.), Ph.D, M.App.Stat. 
Centre for Resource and Environmental Studies
The Australian National University              
Canberra ACT 0200                               
Australia                                       
T: +61 2 6125 7800 email: Simon.Blomberg_at_anu.edu.au
F: +61 2 6125 0757
CRICOS Provider # 00120C

The combination of some data and an aching desire for 
an answer does not ensure that a reasonable answer 
can be extracted from a given body of data.
- John Tukey.