Экспортировать сводные результаты mlogit в Latex с помощью toLatex() или xtable()
Я пытаюсь экспортировать результаты mlogit () в латексную таблицу, но ни одна из моих попыток не удалась!
1) Сначала я попробовал с пакетом xtable ():
> library(xtable)
> s<-summary(mx1)
> tab<-xtable(s, caption= "RPL results")
Errore in UseMethod("xtable") :
no applicable method for 'xtable' applied to an object of class "c('summary.mlogit', 'mlogit')"
2) Затем я попытался с toLatex() из пакета memsic():
> library("memisc")
> s<-summary(mx1)
> toLatex(mtable(s))
Errore in UseMethod("getSummary") :
no applicable method for 'getSummary' applied to an object of class "c('summary.mlogit', 'mlogit')"
Любая идея? Похоже, что в mlogit () отсутствует метод getSummary()
4 ответа
Как сказал @JakobR xtable не умеет обращаться с объектом класса mlogit или же summary.mlogit, Но с тех пор xtable полагаться на S3 Система ООП просто добавить такой метод (используя, например, xtable.summary.lm как шаблон)
require(mlogit)
require(xtable)
### from help page
data(Fishing)
Fish <- mlogit.data(Fishing, varying = c(2:9), shape = "wide", choice = "mode")
modelsum <- summary(mlogit(mode ~ price + catch, data = Fish))
modelsum$CoefTable
## Estimate Std. Error t-value Pr(>|t|)
## boat:(intercept) 0.87137 0.1140428 7.6408 2.1538e-14
## charter:(intercept) 1.49889 0.1329328 11.2755 0.0000e+00
## pier:(intercept) 0.30706 0.1145738 2.6800 7.3627e-03
## price -0.02479 0.0017044 -14.5444 0.0000e+00
## catch 0.37717 0.1099707 3.4297 6.0420e-04
Теперь мы можем написать наш собственный метод:
## check the class first
class(modelsum)
[1] "summary.mlogit" "mlogit"
### write a method from summary.mlogit
xtable.summary.mlogit <- function (x, caption = NULL, label = NULL, align = NULL, digits = NULL,
display = NULL, ...)
{
x <- data.frame(x$CoefTable, check.names = FALSE)
class(x) <- c("xtable", "data.frame")
caption(x) <- caption
label(x) <- label
align(x) <- switch(1 + is.null(align), align, c("r", "r",
"r", "r", "r"))
digits(x) <- switch(1 + is.null(digits), digits, c(0, 4,
4, 2, 4))
display(x) <- switch(1 + is.null(display), display, c("s",
"f", "f", "f", "f"))
return(x)
}
Давайте сделаем простой тест
xtable(modelsum, digits = 2)
## % latex table generated in R 2.15.1 by xtable 1.7-0 package
## % Thu Aug 9 09:09:26 2012
## \begin{table}[ht]
## \begin{center}
## \begin{tabular}{rrrrr}
## \hline
## & Estimate & Std. Error & t-value & Pr($>$$|$t$|$) \\
## \hline
## boat:(intercept) & 0.87 & 0.11 & 7.64 & 0.00 \\
## charter:(intercept) & 1.50 & 0.13 & 11.28 & 0.00 \\
## pier:(intercept) & 0.31 & 0.11 & 2.68 & 0.01 \\
## price & -0.02 & 0.00 & -14.54 & 0.00 \\
## catch & 0.38 & 0.11 & 3.43 & 0.00 \\
## \hline
## \end{tabular}
## \end{center}
## \end{table}
Небольшое редактирование, так как ОП просит поддержки значимости звезд (asterisk функция не выглядит элегантно, я знаю)
## function to add star...
asterisk <- function(y) ifelse(y < 0.001, "***",
ifelse(y < 0.01, "**" ,
ifelse(y < 0.05, "*",
ifelse(y < 0.1, ".", ""))))
DF <- read.table(text = capture.output(data.frame(modelsum$CoefTable)))
DF$V6 <- asterisk(DF[,4])
names(DF) <- c(colnames(modelsum$CoefTable), " ")
xtable(DF)
## % latex table generated in R 2.15.1 by xtable 1.7-0 package
## % Thu Aug 9 11:46:31 2012
## \begin{table}[ht]
## \begin{center}
## \begin{tabular}{rrrrrl}
## \hline
## & Estimate & Std. Error & t-value & Pr($>$$|$t$|$) & \\
## \hline
## boat:(intercept) & 0.87 & 0.11 & 7.64 & 0.00 & *** \\
## charter:(intercept) & 1.50 & 0.13 & 11.28 & 0.00 & *** \\
## pier:(intercept) & 0.31 & 0.11 & 2.68 & 0.01 & ** \\
## price & -0.02 & 0.00 & -14.54 & 0.00 & *** \\
## catch & 0.38 & 0.11 & 3.43 & 0.00 & *** \\
## \hline
## \end{tabular}
## \end{center}
## \end{table}
Решение, вдохновленное этой темой
Проблема в том, что xtable не сейчас, как справиться с чем-то вроде summary.mlogitОднако вы можете, например, извлечь таблицу коэффициентов с помощью s$CoefTable и поэтому xtable(s$CoefTable) буду работать.
Вы также можете получить хорошую сводную таблицу без написания функции, если вы просто используете функцию latex от Hmisc пакет. Пытаться
library(Hmisc)
latex(modelsum$CoefTable, digits=3) # using @dickoa's example
Как вы можете видеть, это дает вам нечто похожее на то, что было получено с помощью решения @dickoa.
# With caption
latex(modelsum$CoefTable, digits=3,
caption='A mlogit summary table')
Вы можете прочитать файл справки, где вы можете получить много вариантов для игры (?latex).
Что касается функции mtable() пакета memisc, одним из решений является написание специального метода getSummary, как предлагается здесь для функции lme4 (): https://stat.ethz.ch/pipermail/r-sig-mixed-models/2009q1/002058.html
library(lme4)
library(memisc)
### create three models
fm1 <- lmer(Reaction ~ 1 + (Days|Subject), sleepstudy)
fm1.1 <- lmer(Reaction ~ Days + (Days|Subject), sleepstudy)
fm1.2 <- lmer(Reaction ~ as.factor(Days) + (Days|Subject), sleepstudy)
### note: need to run the code below fro setCoefTemplate and
### getSummary.lmer first
mtable("Model 1"=fm1, "Model 2"=fm1.1, "Model 3"=fm1.2,
coef.style = "est.ci", # using "homegrown" est.ci, specified above
summary.stats=c("AIC","BIC"),
getSummary = "getSummary.lmer")#,
setCoefTemplate(
est.ci=c(
est = "($est:#)($p:*)",
ci = "[($lwr:#),($upr:#)]"))
getSummary.lmer <- function (obj, alpha = 0.05, ...)
{
require(lme4)
smry <- summary(obj)
#N <- if (length(weights(obj))) ### NOTE: how to deal with groups/samp size?
# sum(weights(obj))
#else sum(smry$df[1:2])
coef <- smry at coefs
lower <- qnorm(p = alpha/2, mean = coef[, 1], sd = coef[,2])
upper <- qnorm(p = 1 - alpha/2, mean = coef[, 1], sd = coef[,2])
if (ncol(smry at coefs) == 3) {
p <- (1 - pnorm(smry at coefs[,3]))*2 # NOTE: no p-values for lmer() due to
# unclear dfs; calculate p-values based on z
coef <- cbind(coef, p, lower, upper)
} else {
coef <- cbind(coef, lower, upper) # glmer will have 4 columns with p-values
}
colnames(coef) <- c("est", "se", "stat", "p", "lwr", "upr")
#phi <- smry$dispersion
#LR <- smry$null.deviance - smry$deviance
#df <- smry$df.null - smry$df.residual
ll <- smry at AICtab[3][,1]
deviance <- smry at AICtab[4][,1]
#if (df > 0) {
# p <- pchisq(LR, df, lower.tail = FALSE)
# L0.pwr <- exp(-smry$null.deviance/N)
# McFadden <- 1 - smry$deviance/smry$null.deviance
# Cox.Snell <- 1 - exp(-LR/N)
# Nagelkerke <- Cox.Snell/(1 - L0.pwr)
#}
#else {
# LR <- NA
# df <- NA
# p <- NA
# McFadden <- NA
# Cox.Snell <- NA
# Nagelkerke <- NA
#}
AIC <- smry at AICtab[1][,1] # NOTE: these are both data.frames? not sure why...
BIC <- smry at AICtab[2][,1]
### NOTE: don't see a similar slot for "xlevels" to get levels of
### factor variables used as predictors; for time being, force
### user to specify explicitly; nope that didn't work...
#if (fac != NULL) {
# n <- length(fac)
# xlevels <- vector(n, mode = "list")
# for (i in 1:n) {
# xlevels[i] <- levels(obj at frame[,fac[i]])
# }
# }
#sumstat <- c(phi = phi, LR = LR, df = df, p = p, logLik = ll,
# deviance = deviance, McFadden = McFadden, Cox.Snell = Cox.Snell,
# Nagelkerke = Nagelkerke, AIC = AIC, BIC = BIC, N = N)
sumstat <- c(logLik = ll, deviance = deviance, AIC = AIC, BIC = BIC)
list(coef = coef, sumstat = sumstat,
contrasts = attr(model.matrix(obj), "contrasts"),
xlevels = NULL, call = obj at call)
}