export mlogit summary results to Latex with toLatex() or xtable()

Question

I'm trying to export mlogit() results into a latex table but none of my attempts succeeded!

1) First I tried with the package 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) Then I tried with toLatex() from the package 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')"

Any idea? It seems that mlogit() is missing a getSummary() method

Answer 1

You can also obtain a nice summary table without writing a function if you just use function latex from Hmisc package. Try

library(Hmisc)
latex(modelsum$CoefTable, digits=3) # using @dickoa's example

As you can see this gives you something similar to that obtained using @dickoa's solution.

# With caption
latex(modelsum$CoefTable, digits=3, 
      caption='A mlogit summary table')

You can read the help file where you can get a lot of options to play with (?latex).

Answer 2

As @JakobR said xtable doesn't know how to deal with object of class mlogit or summary.mlogit. But since xtable rely on S3 OOP system is simple to add such method (using for example xtable.summary.lm as template)

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

Now we can write our own method :

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

Let's make a simple test

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}

Small edit since the OP ask for significance stars support (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}

Solution inspired by this thread

Answer 3

with regard to the mtable() function of the memisc package one solution is to write a custom getSummary method as suggested here for the function 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)
}

Answer 4

The problem is, that xtable does not now how to handle something like summary.mlogit However you can for example extract the coefficent table with s$CoefTable and thus xtable(s$CoefTable) will work.