[R] p value in MANOVA

[Marc Schwartz]

On Thu, 2003-11-20 at 11:12, Tu Yu-Kang wrote:
> Dear R users,
> 
> Can anyone tell me how to get the p value out of the output from 
> summary.manova?
> 
> I tried all the methods I can think of, but failed.
> 
> Many thanks
> 
> Yu-Kang


Using the example in ?summary.manova:

## Example on producing plastic film from Krzanowski (1998, p. 381)
tear <- c(6.5, 6.2, 5.8, 6.5, 6.5, 6.9, 7.2, 6.9, 6.1, 6.3,
          6.7, 6.6, 7.2, 7.1, 6.8, 7.1, 7.0, 7.2, 7.5, 7.6)
gloss <- c(9.5, 9.9, 9.6, 9.6, 9.2, 9.1, 10.0, 9.9, 9.5, 9.4,
           9.1, 9.3, 8.3, 8.4, 8.5, 9.2, 8.8, 9.7, 10.1, 9.2)
opacity <- c(4.4, 6.4, 3.0, 4.1, 0.8, 5.7, 2.0, 3.9, 1.9, 5.7,
             2.8, 4.1, 3.8, 1.6, 3.4, 8.4, 5.2, 6.9, 2.7, 1.9)
Y <- cbind(tear, gloss, opacity)
rate <- factor(gl(2,10), labels=c("Low", "High"))
additive <- factor(gl(2, 5, len=20), labels=c("Low", "High"))
fit <- manova(Y ~ rate * additive)

summary(fit) results in:

> summary(fit)
              Df Pillai approx F num Df den Df   Pr(>F)   
rate           1 0.6181   7.5543      3     14 0.003034 **
additive       1 0.4770   4.2556      3     14 0.024745 * 
rate:additive  1 0.2229   1.3385      3     14 0.301782   
Residuals     16                                          
---
Signif. codes:  0 `***' 0.001 `**' 0.01 `*' 0.05 `.' 0.1 ` ' 1 


To get just the p values, you would use:

> summary(fit)$stats[1:3, "Pr(>F)"]
         rate      additive rate:additive 
  0.003034045   0.024745281   0.301781645 

and if you just want just the numbers as a vector:

> as.numeric(summary(fit)$stats[1:3, "Pr(>F)"])
[1] 0.003034045 0.024745281 0.301781645



To gain some insight into how to get to the above, ?summary.manova tells
you that the values returned by the method are in a list with one
component being 'stats'. Alternatively, you can use:

str(summary(fit))

which will give you the details of the internal structure of the list
object returned by summary.manova(). See ?str for more information.

The output of the above includes the following:

...
 $ stats      : num [1:4, 1:6]  1.000  1.000  1.000 16.000  0.618 ...
  ..- attr(*, "dimnames")=List of 2
  .. ..$ : chr [1:4] "rate" "additive" "rate:additive" "Residuals"
  .. ..$ : chr [1:6] "Df" "Pillai" "approx F" "num Df" ...
...

This tells you that summary(fit)$stats is a 4 x 6 matrix in this
example.

Then to get the complete list of dimnames for 'stats' use:

> dimnames(summary(fit)$stats) 
[[1]]
[1] "rate"          "additive"      "rate:additive" "Residuals"    

[[2]]
[1] "Df"       "Pillai"   "approx F" "num Df"   "den Df"   "Pr(>F)"  


So to get the p values, you want the first three rows (in this case) and
the "Pr(>F)" column, hence:

summary(fit)$stats[1:3, "Pr(>F)"]

HTH,

Marc Schwartz