[R] AOV in MASS not the same??

[Peter Ho]

I would appreciate it if someone could explain the results of the 
example from the aov() help file. The output given below is different 
from book
Venables and Ripley -  MASS

The results In R1.5.1 (Under windows) is as follows:
 > N <- c(0,1,0,1,1,1,0,0,0,1,1,0,1,1,0,0,1,0,1,0,1,1,0,0)
 > P <- c(1,1,0,0,0,1,0,1,1,1,0,0,0,1,0,1,1,0,0,1,0,1,1,0)
 > K <- c(1,0,0,1,0,1,1,0,0,1,0,1,0,1,1,0,0,0,1,1,1,0,1,0)
 > yield <- c(49.5,62.8,46.8,57.0,59.8,58.5,55.5,56.0,62.8,55.8,69.5,55.0,
+            62.0,48.8,45.5,44.2,52.0,51.5,49.8,48.8,57.2,59.0,53.2,56.0)
 > npk <- data.frame(block=gl(6,4), N=factor(N), P=factor(P),
+                   K=factor(K), yield=yield)
 >
 > ( npk.aov <- aov(yield ~ block + N*P*K, npk) )
Call:
   aov(formula = yield ~ block + N * P * K, data = npk)

Terms:
                           block        N                P        K      
        N:P              N:K      P:K
Sum of Squares  343.2950 189.2817   8.4017  95.2017  21.2817  33.1350   
0.4817
Deg. of Freedom        5        1                    1        1          
      1                1                1
                Residuals
Sum of Squares   185.2867
Deg. of Freedom        12

Residual standard error: 3.929447
1 out of 13 effects not estimable
Estimated effects may be unbalanced
 > summary(npk.aov)
                     Df     Sum Sq   Mean Sq   F value        Pr(>F)  
block             5      343.30    68.66        4.4467        0.015939 *
N                  1      189.28    189.28      12.2587      0.004372 **
P                   1       8.40       8.40          0.5441       
 0.474904  
K                  1       95.20     95.20        6.1657        0.028795 *
N:P               1      21.28      21.28        1.3783        0.263165  
N:K               1      33.13     33.13         2.1460        0.168648  
P:K                1       0.48      0.48           0.0312       
 0.862752  
Residuals   12 185.29   15.44                   
---
Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1
 > coefficients(npk.aov)
(Intercept)      block2      block3      block4      block5      block6
 51.8250000   3.4250000   6.7500000  -3.9000000  -3.5000000   2.3250000
         N1          P1          K1       N1:P1       N1:K1       P1:K1
  9.8500000   0.4166667  -1.9166667  -3.7666667  -4.7000000   0.5666667


The output for the ANOVA table is exactly the same as in Venables and 
Ripley MASS page 177, but the values of the coefficients are different. 
Also if
alias(npk.aov) is used,as in the book, a different result is also obtained :

Model :
yield ~ block + N + P + K + N:P + N:K + P:K + N:P:K

Complete :
                     (Intercept) block2 block3 block4 block5 block6     
N1        P1        K1        N1:P1     N1:K1     P1:K1
N1:P1:K1                          0.25   0.25   0.25                  
          -0.25     -0.25     -0.25      0.50          0.50      0.50



Can someone explain why this is different from the book MASS.

I also have a more general question in relation to the coefficient 
terms. Why are there more than one coefficint for the blocking factor. I 
want to construct a normal probability plot of effects and since the 
ANOVA table gives me one block term, why are are more than one 
coefficient term for blocks. Is there another way to extract the 
regression coefficients? I don't understand why there are more than one 
blocking coefficient .

Also in this analysis the block coefficients are from block2 to block 6, 
whereas we have block1 to block5 in MASS.

Another point is that at the end of the ANOVA (summary) table  the 
warning "Estimated effects may be unbalanced" (aslo different from the 
book) . In this case, should aov() be used or should I refit the model 
with , say lme() ?


Thanks

Peter






-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-.-
r-help mailing list -- Read http://www.ci.tuwien.ac.at/~hornik/R/R-FAQ.html
Send "info", "help", or "[un]subscribe"
(in the "body", not the subject !)  To: r-help-request at stat.math.ethz.ch
_._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._._