2012-08-09 309 views
2

我試圖將mlogit()結果 導出到乳膠表中,但是我的任何嘗試都沒有成功!使用toLatex()或xtable()輸出mlogit摘要結果到乳膠

1)首先,我試圖與包xtable():從包MEMSIC()

> 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():

> 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()方法

+0

請指明可以找到'mlogit()'的位置。另外,xtable非常好地轉換數據幀,所以一個簡單的方法就是在彙總結果中使用'str',提取所需的組件,然後調用'xtable'。此外,可複製的代碼(包含數據)將使人們更容易地爲您提供幫助 – richiemorrisroe 2012-08-09 08:39:03

+0

@richiemorrisroe:我在這裏的超鏈接mlogit:http://cran.r-project.org/web/packages/mlogit/index.html – danfreak 2012-08-09 09:18:11

回答

4

的問題是,xtable不是現在如何處理類似summary.mlogit 但是例如,您可以用s$CoefTable提取係數中的表,從而xtable(s$CoefTable)將工作做。

+0

這真的很快,而@dickoa解決方案提供更個性化的結果(即數字等) – danfreak 2012-08-09 10:47:36

4

由於@JakobR表示xtable不知道如何處理類mlogitsummary.mlogit的對象。 但由於xtable依靠S3 OOP系統是簡單的添加這樣的方法(例如使用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} 

小編輯自OP詢問意義明星支持(The asterisk function doesn't look elegant I know

## 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} 

溶液通過此thread

+0

感謝偉大的解決方案!我怎麼能將p值轉換爲星號以獲得更快的重要讀數? – danfreak 2012-08-09 10:48:52

+0

@danfreak:我做了一個小小的編輯,加入了意義明星 – dickoa 2012-08-09 12:22:15

+0

很棒:非常感謝! – danfreak 2012-08-09 15:08:56

0

相對於memisc包一種解決方案是編寫自定義getSummary方法,如下建議功能lme4()的mtable()函數的啓發: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) 
} 
1

您還可以得到一個很好的總結表,而無需編寫,如果你只是使用功能latexHmisc包的功能。嘗試

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)玩。