[R] Specifying an appropriate error term in a hierarchical regression

[Chris Bergstresser]

Hi all --

   So I obtained a copy of _Statistical Models in S_, which I hoped
would answer my question, but it doesn't quite.  On page 151, it
discusses nested terms, which is what I'm working with.  So I tried
using "value ~ a/b" which apparently should fit factor b within factor
a.  But the summary is the exact same summary that "value ~ a + b" and
"value ~ a * b" produces, which doesn't seem right.  And the F value
reported for A is wrong; it should compare the mean square of A with
the mean square of B(A), not the residual mean square as it does.
   What am I doing wrong?

On 4/11/06, Chris Bergstresser <chris at subtlety.com> wrote:
> Hi all --
>
>    So I'm working through my statistics homework again, and trying to
> reproduce the examples in the book (Kirk's _Experimental Design_,
> third edition) in R.  This is a completely randomized hierarchical
> design (CRH-28(A)).  The B factor is completely nested within the A
> factor.  Pages 480-482, for those playing along at home.
>
>    I can use:
>
> summary(aov(value ~ a + Error(b), data = ex));
>
>    to get the correct F value for the main effect of A.  I can use
>
> summary(aov(value ~ b, data = ex));
>
>    to get the correct values for B(A) and the within cell SS.  But I
> can't find any documentation about constructing the Error term to get
> this output in a single analysis (except for
> http://www.psych.upenn.edu/~baron/rpsych/rpsych.html, but Kirk doesn't
> talk about these tests in terms of Error strata, so it's a little hard
> to figure out the correspondence).
>    Also, the documentation on the Error term in ?aov is rather
> perfunctory.  There's no mention of the "/" operator, for example.
>
> ex = scan()
>  3  6  3  3
>  1  2  2  2
>  5  6  5  6
>  2  3  4  3
>  7  8  7  6
>  4  5  4  3
>  7  8  9  8
> 10 10  9 11
>
> a = factor(rep(paste("a", 1:2, sep = ""), each = 16));
> b = factor(rep(paste("b", 1:8, sep = ""), each = 4));
>
> ex = data.frame(value = ex, a, b);
>
> summary(aov(value ~ a + Error(b), data = ex));
> summary(aov(value ~ b, data = ex));
>