Programmers,
I'm trying specify in R that the error term in aov is the interaction of block and my treatment. Can you help me? I constructed an example below.
-Rik
d <- data.frame("a"=as.factor(rep(1:3,12)),
"b"=as.factor(rep(1:3,each=3)),
"block"=as.factor(rep(1:4,each=9)));
d$y <- as.numeric(d$a)*rnorm(36, mean=.75, sd=1) +
as.numeric(d$b)*rnorm(36, mean=1.2, sd=1) +
as.numeric(d$block)*rnorm(36, mean=1.2, sd=1)+
rnorm(36);
anova1=aov(y~a+b+Error(a*b*block),data=d)
summary(anova1)
I did the above but my output contains no p-value:
Error: a
Df Sum Sq Mean Sq
a 2 50.38 25.19
Error: b
Df Sum Sq Mean Sq
b 2 16.61 8.307
Error: block
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 3 138.2 46.06
Error: a:b
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 4 46.39 11.6
Error: a:block
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 6 69.71 11.62
Error: b:block
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 6 110 18.33
Error: a:b:block
Df Sum Sq Mean Sq F value Pr(>F)
Residuals 12 289.5 24.12
As help for summary.aov states, the p-values are only returned "if there are non-zero residual degrees of freedom." This is sensible, because with out the residual df, there is no randomness left to estimate in the model.
Your model is 'fully saturated' meaning there is no replication within the combinations of factor levels it attempts to estimate. See this cross-validated answer for further discussion.
To diagnose this problem in your data, try the code below. It demonstrates there is only one observation within each block/a*b combination:
> xtabs( ~ a + b + block, d)
, , block = 1
b
a 1 2 3
1 1 1 1
2 1 1 1
3 1 1 1
, , block = 2
b
a 1 2 3
1 1 1 1
2 1 1 1
3 1 1 1
, , block = 3
b
a 1 2 3
1 1 1 1
2 1 1 1
3 1 1 1
, , block = 4
b
a 1 2 3
1 1 1 1
2 1 1 1
3 1 1 1