[R] extracting p-values from Anova objects (from the car library)

David Winsemius dwinsemius at comcast.net
Tue Aug 24 00:19:42 CEST 2010


On Aug 23, 2010, at 5:35 PM, Johan Steen wrote:

> Thanks for your replies,
>
> but unfortunately none of them seem to help.
> I do get p-values in the output, but can't seem to locate them  
> anywhere in these objects via the str() function.

That is because in the case of Anova.mlm objects ... they aren't there.

> I also get very different output using str() than you obtained from  
> the lm help page
>
> Here's my output:
>
> > A <- factor( rep(1:2,each=3) )
> > B <- factor( rep(1:3,times=2) )
> > idata <- data.frame(A,B)
> > idata
>  A B
> 1 1 1
> 2 1 2
> 3 1 3
> 4 2 1
> 5 2 2
> 6 2 3
> >
> > fit <- lm( cbind(a1_b1,a1_b2,a1_b3,a2_b1,a2_b2,a2_b3) ~ sex,  
> data=Data.wide)
> > result <- Anova(fit, type="III", test="Wilks", idata=idata,  
> idesign=~A*B)
> > result
>
> Type III Repeated Measures MANOVA Tests: Wilks test statistic
>            Df test stat approx F num Df den Df    Pr(>F)
> (Intercept)  1   0.02863   610.81      1     18 2.425e-15
> sex          1   0.76040     5.67      1     18   0.02849
> A            1   0.91390     1.70      1     18   0.20925
> sex:A        1   0.99998     0.00      1     18   0.98536
> B            1   0.26946    23.05      2     17 1.443e-05
> sex:B        1   0.98394     0.14      2     17   0.87140
> A:B          1   0.27478    22.43      2     17 1.704e-05
> sex:A:B      1   0.98428     0.14      2     17   0.87397

Remember I suggested looking for print methods?   methods(print)  
produces a print.Anova.mlm citation. The above output is produced by  
that function which is not visible unless you do:

 > getAnywhere(print.Anova.mlm)
A single object matching ‘print.Anova.mlm’ was found
It was found in the following places
   registered S3 method for print from namespace car
   namespace:car
with value

function (x, ...)
{
     test <- x$test
     repeated <- x$repeated
     ntests <- length(x$terms)
     tests <- matrix(NA, ntests, 4)
     if (!repeated)
         SSPE.qr <- qr(x$SSPE)
     for (term in 1:ntests) {
         eigs <- Re(eigen(qr.coef(if (repeated) qr(x$SSPE[[term]])  
else SSPE.qr,
             x$SSP[[term]]), symmetric = FALSE)$values)
         tests[term, 1:4] <- switch(test, Pillai = stats:::Pillai(eigs,
             x$df[term], x$error.df), Wilks = stats:::Wilks(eigs,
             x$df[term], x$error.df), `Hotelling-Lawley` =  
stats:::HL(eigs,
             x$df[term], x$error.df), Roy = stats:::Roy(eigs,
             x$df[term], x$error.df))
     }
     ok <- tests[, 2] >= 0 & tests[, 3] > 0 & tests[, 4] > 0
     ok <- !is.na(ok) & ok
     tests <- cbind(x$df, tests, pf(tests[ok, 2], tests[ok, 3],
         tests[ok, 4], lower.tail = FALSE))
     rownames(tests) <- x$terms
     colnames(tests) <- c("Df", "test stat", "approx F", "num Df",
         "den Df", "Pr(>F)")
     tests <- structure(as.data.frame(tests), heading = paste("\nType ",
         x$type, if (repeated)
             " Repeated Measures", " MANOVA Tests: ", test, " test  
statistic",
         sep = ""), class = c("anova", "data.frame"))
     print(tests)
     invisible(x)
}
<environment: namespace:car>

Notice that after printing "tests",  those results are basically  
thrown away, and the function invisibly returns its argument. Hacking  
that function would be simple matter if you wanted to build a list of   
results.


> > summary(result)
>
> snipped long and uninformative output.

David Winsemius, MD
West Hartford, CT



More information about the R-help mailing list