Использование объекта 'mlm` в выходных данных `mtable`
Есть ли способ работать с mlm объекты в mtable от memisc пакет?
Без использования матрицы множественных ответов я хочу получить что-то вроде:
library(car)
library(memisc)
lm1 = lm(Sepal.Length ~ Petal.Length + Petal.Width + Species, data=iris)
lm2 = lm(Sepal.Width ~ Petal.Length + Petal.Width + Species, data=iris)
mtable(lm1, lm2)
который производит
Calls:
lm1: lm(formula = Sepal.Length ~ Petal.Length + Petal.Width + Species,
data = iris)
lm2: lm(formula = Sepal.Width ~ Petal.Length + Petal.Width + Species,
data = iris)
===============================================
lm1 lm2
-----------------------------------------------
(Intercept) 3.683*** 3.048***
(0.107) (0.094)
Petal.Length 0.906*** 0.155*
(0.074) (0.065)
Petal.Width -0.006 0.623***
(0.156) (0.136)
Species: versicolor/setosa -1.598*** -1.764***
(0.206) (0.180)
Species: virginica/setosa -2.113*** -2.196***
(0.304) (0.265)
-----------------------------------------------
R-squared 0.837 0.551
adj. R-squared 0.832 0.539
sigma 0.339 0.296
F 185.769 44.496
p 0.000 0.000
Log-likelihood -48.116 -27.711
Deviance 16.681 12.708
AIC 108.231 67.423
BIC 126.295 85.486
N 150 150
===============================================
Но:
mlmIris = lm(cbind(Sepal.Length, Sepal.Width) ~ Petal.Length + Petal.Width + Species, data=iris)
mtable(mlmIris)
производит
Error in qt(p = alpha/2, df = dendf) :
Non-numeric argument to mathematical function
Я не собираюсь воспроизводить способы, которые я пытался извлечь lm объект, который я могу использовать в mtable, Достаточно сказать, что ни один из них не работал.
1 ответ
Как предлагается в комментариях, вам нужно написать getSummary метод для mlm объекты, потому что там в настоящее время нет никаких. lm метод, который вы получаете по наследству, не работает.
Я сделал быстрый удар в этом и предоставить getSummary.mlm ниже вместе с подходящим setSummaryTemplate, Коэффициенты и стандартные ошибки теперь правильно извлекаются. Что требует дополнительной работы - это извлечение сводной статистики (R-квадрат, остаточная сумма квадратов, F-статистика и т. Д.), Над которой я не работал. Это должно дать вам хорошее начало, хотя. Если вы улучшите метод дальше, пожалуйста, рассмотрите возможность предоставления его Мартину Эльфу (memisc сопровождающий), чтобы он был непосредственно доступен в memisc,
После поиска функций, представленных ниже, это работает:
R> mtable(mlmIris)
Calls:
mlmIris: lm(formula = cbind(Sepal.Length, Sepal.Width) ~ Petal.Length +
Petal.Width + Species, data = iris)
=====================================================
Sepal.Length Sepal.Width
-----------------------------------------------------
(Intercept) 3.683*** 3.048***
(0.107) (0.094)
Petal.Length 0.906*** 0.155*
(0.074) (0.065)
Petal.Width -0.006 0.623***
(0.156) (0.136)
Species: versicolor/setosa -1.598*** -1.764***
(0.206) (0.180)
Species: virginica/setosa -2.113*** -2.196***
(0.304) (0.265)
-----------------------------------------------------
N 150
=====================================================
Исходный код:
getSummary.mlm <- function(obj, alpha = 0.05, ...)
{
## extract coefficient summary
cf <- lapply(summary(mlmIris), "[[", "coefficients")
## augment with confidence intervals
cval <- qnorm(1 - alpha/2)
for(i in seq_along(cf)) cf[[i]] <- cbind(cf[[i]],
cf[[i]][, 1] - cval * cf[[i]][, 2],
cf[[i]][, 1] + cval * cf[[i]][, 2])
## collect in array
nam <- unique(unlist(lapply(cf, rownames)))
acf <- array(dim = c(length(nam), 6, length(cf)),
dimnames = list(nam, c("est", "se", "stat", "p", "lwr", "upr"), colnames(coef(obj))))
for(i in seq_along(cf)) acf[rownames(cf[[i]]), , i] <- cf[[i]]
## return everything
return(list(
coef = acf,
sumstat = c(
"N" = nobs(obj)
),
contrasts = obj$contrasts,
xlevels = obj$xlevels,
call = obj$call
))
}
setSummaryTemplate("mlm" = c(
"N" = "($N:d)"
))