2014-09-17 7 views
3

Интересно, как извлечь часть Multivariate Tests: Site с выхода fm1 в следующем MWE. Любая помощь будет высоко оценена. БлагодаряИзвлечение многомерных тестов с выхода функции Anova или Manova из пакета автомобилей

library(car) 
fm1 <- summary(Anova(lm(cbind(Al, Fe, Mg, Ca, Na) ~ Site, data=Pottery))) 
fm1 

Type II MANOVA Tests: 

Sum of squares and products for error: 
      Al   Fe   Mg   Ca   Na 
Al 48.2881429 7.08007143 0.60801429 0.10647143 0.58895714 
Fe 7.0800714 10.95084571 0.52705714 -0.15519429 0.06675857 
Mg 0.6080143 0.52705714 15.42961143 0.43537714 0.02761571 
Ca 0.1064714 -0.15519429 0.43537714 0.05148571 0.01007857 
Na 0.5889571 0.06675857 0.02761571 0.01007857 0.19929286 

------------------------------------------ 

Term: Site 

Sum of squares and products for the hypothesis: 
      Al   Fe   Mg   Ca   Na 
Al 175.610319 -149.295533 -130.809707 -5.8891637 -5.3722648 
Fe -149.295533 134.221616 117.745035 4.8217866 5.3259491 
Mg -130.809707 117.745035 103.350527 4.2091613 4.7105458 
Ca -5.889164 4.821787 4.209161 0.2047027 0.1547830 
Na -5.372265 5.325949 4.710546 0.1547830 0.2582456 

Multivariate Tests: Site 
       Df test stat approx F num Df den Df  Pr(>F)  
Pillai   3 1.55394 4.29839  15 60.00000 2.4129e-05 *** 
Wilks    3 0..08854  15 50.09147 1.8404e-12 *** 
Hotelling-Lawley 3 35.43875 39.37639  15 50.00000 < 2.22e-16 *** 
Roy    3 34.16111 136.64446  5 20.00000 9.4435e-15 *** 
--- 
Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

ответ

2

Я также не мог найти, как извлечь таблицу тестов, но в качестве обходного пути вы можете вычислить результаты, выполнив команду Anova по всем типам тестов.

Однако метод печати print.Anova.mlm не возвращает результаты, поэтому это нужно немного подкорректировать.

library(car) 

# create new print function 
outtests <- car:::print.Anova.mlm 

# allow the function to return the results and disable print 
body(outtests)[[16]] <- quote(invisible(tests)) 
body(outtests)[[15]] <- NULL 

# Now run the regression 
mod <- lm(cbind(Al, Fe, Mg, Ca, Na) ~ Site, data=Pottery) 

# Run the Anova over all tests 
tab <- lapply(c("Pillai", "Wilks", "Hotelling-Lawley", "Roy"), 
        function(i) outtests(Anova(mod, test.statistic=i))) 

tab <- do.call(rbind, tab) 
row.names(tab) <- c("Pillai", "Wilks", "Hotelling-Lawley", "Roy") 
tab 

    # Type II MANOVA Tests: Pillai test statistic 
    #    Df test stat approx F num Df den Df Pr(>F)  
#Pillai   3  1.554 4.298  15 60.000 2.413e-05 *** 
#Wilks    3  0.012 13.089  15 50.091 1.840e-12 *** 
#Hotelling-Lawley 3 35.439 39.376  15 50.000 < 2.2e-16 *** 
#Roy    3 34.161 136.644  5 20.000 9.444e-15 *** 
#--- 
#Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

В выходной таблице класса anova и data.frame вы можете использовать xtable на нем.

xtable:::xtable(tab) 
+0

(+1): Хорошее решение. – MYaseen208

+0

+1. Этот метод даже никогда не приходил мне в голову. Я делаю так много манипуляций с данными, что забыл о статистике. –

3

fm1$multivariate.tests получает вас к Site части fm1 продукции.

Тогда вы можете использовать комбинацию cat и capture.output для приятной печати или просто capture.output для символа.

> cat(capture.output(fm1$multivariate.tests)[18:26], sep = "\n") 
# 
# Multivariate Tests: Site 
#     Df test stat approx F num Df den Df  Pr(>F)  
# Pillai   3 1.55394 4.29839  15 60.00000 2.4129e-05 *** 
# Wilks    3 0..08854  15 50.09147 1.8404e-12 *** 
# Hotelling-Lawley 3 35.43875 39.37639  15 50.00000 < 2.22e-16 *** 
# Roy    3 34.16111 136.64446  5 20.00000 9.4435e-15 *** 
# --- 
# Signif. codes: 0 ‘***’ 0.001 ‘**’ 0.01 ‘*’ 0.05 ‘.’ 0.1 ‘ ’ 1 

Update: Из результата

unlist(fm1$multivariate.tests, recursive = FALSE) 

это не выглядит как результаты легко доступны как числовые значения. Итак, как вы просили, вот что нужно, чтобы манипулировать результатами в матрицу. Сделав это, а затем увидев ответ пользователя20650, я рекомендую вам следовать этому предложению и получить значения через таблицу ANOVA.

co <- capture.output(fm1$multivariate.tests)[20:24] 
s <- strsplit(gsub("([*]+$)|[<]", "", co[-1]), "\\s+") 
dc <- do.call(rbind, lapply(s, function(x) as.numeric(x[-1]))) 
row.names(dc) <- sapply(s, "[", 1) 
s2 <- strsplit(co[1], " ")[[1]] 
s2 <- s2[nzchar(s2)] 
s3 <- s2[-c(1, length(s2))] 
colnames(dc) <- c(s2[1], paste(s3[c(TRUE, FALSE)], s3[c(FALSE, TRUE)]), s2[10]) 
dc 
#     Df test stat approx F num Df den Df  Pr(>F) 
# Pillai   3 1.55394 4.29839  15 60.00000 2.4129e-05 
# Wilks    3 0..08854  15 50.09147 1.8404e-12 
# Hotelling-Lawley 3 35.43875 39.37639  15 50.00000 2.2200e-16 
# Roy    3 34.16111 136.64446  5 20.00000 9.4435e-15 

Если кому-то хочется улучшить мой второй кусок кода, не стесняйтесь.

+0

(+1): Ницца обходной путь. Может ли это быть преобразовано в data.frame для использования для xtable? Спасибо – MYaseen208

+0

Хм, хорошо, я мог бы попытаться получить его в матрицу с некоторыми манипуляциями с строкой –

+0

Это было бы хорошо для меня. – MYaseen208