【问题标题】:How extract regression results from lme, lmer, glmer to Latex?如何从lme、lmer、glmer中提取回归结果到Latex?
【发布时间】:2012-02-23 04:36:51
【问题描述】:

我正在使用 lme、lmer 和 glmer 拟合模型。我需要使用 summary() 对象构建表格并导出到 Latex 以显示我的结果。 xtable、mtable 和 apsrtable 不起作用。我看到了一篇关于 lme4 对象的解决方案的前一篇文章(下面的链接),但不适用于这些对象。

http://leftcensored.skepsi.net/2011/03/13/code-latex-tables-for-lme4-models/

这是我正在拟合的模型的两个示例:

lme(y ~  time, data, na.action=na.omit, method="REML", random = ~ 1 | subject, control=lmeControl(msMaxIter = 200, msVerbose = TRUE))

glmer(y ~ time + (time | subject), data, family=binomial(link = "logit"), REML=T, control=list(maxIter = 800, maxFN=1000, msVerbose = TRUE))

有什么帮助吗?

谢谢

【问题讨论】:

    标签: r latex


    【解决方案1】:

    我刚刚发现summary.mer 对象存在一个coef 方法,它提供了所有必要的数据(用于固定效果)。返回的对象(在强制转换为 data.frame 之后)可以轻松地移交给选择的格式化包(例如,xtableascii)。
    请参阅以下示例(仅生成可用的data.frame):

    require(lme4)
    
    gm1 <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd),
                  family = binomial, data = cbpp)
    
    (res.table <- as.data.frame(coef(summary(gm1))))
    ##             Estimate Std. Error z value        Pr(>|z|)
    ## (Intercept)  -1.3985     0.2279  -6.137 0.0000000008416
    ## period2      -0.9923     0.3054  -3.249 0.0011562741408
    ## period3      -1.1287     0.3260  -3.462 0.0005368285553
    ## period4      -1.5804     0.4288  -3.686 0.0002282168737
    

    【讨论】:

      【解决方案2】:

      编辑:

      在编辑时 lme4 包已更新,memisc 不再适用于这些对象。包 texreg 是一个替代方案。我已经留下了这个答案,以防 memisc 得到更新并再次开始工作。

      memisc 包执行 lme4 表:

      这是我写的一些代码的sn-p:

      GPusenonMH=lmer(GPEtc_c~Age.y+Measure+Gender+Marital2+Work2+(1|NHS), family="poisson", data=subset(lemurdata, Measure %in% c(1,3)))
      
      model1=mtable(GPusetotal, GPuseMH, GPusenonMH, summary.stats=FALSE)
      
      toLatex(model1)
      

      如果你想要这些东西,显然你可以将 summary.stats=TRUE 设置为。

      请注意,dcolumn 和 booktabs Latex 包都是默认使用的,因此要么将它们放在你的 Latex 序言中,要么使用帮助文件中的命令(useBooktabs=FALSE,useDcolumn=FALSE)关闭它们。

      【讨论】:

      • 谢谢! mtable 与我的 GLMER 完美配合。帮助页面 (?mtable) 也非常有用,它展示了如何以比您可能选择在 R 中使用的更正式的术语重新标记变量和模型
      【解决方案3】:

      lme,下面是我的个人版本;您可以使用其他类似的插件下载它,例如从

      中提取 lme/lm/glm 表的 p 值的 \Sexpr{} 字符串作为 Dmisc

      http://www.menne-biomed.de/download

      这是非常个性化的,但如果我非常喜欢四舍五入到真正有效的数字。抱歉,nlme 包可以满足我的所有需求(而且不仅仅是 lme/gaussian),所以还没有 lme4 版本。

      "latex.summary.lme" <-
      function(object, title="",parameter=NULL, file="",
        shadep=0.05,caption=NULL,label=NULL,ctable=FALSE,form=NULL,
        interceptp = FALSE, moredec=0, where="!htbp", ...) {
        # This function can be mis-used for gls models when an explicit
        # form is given
        options(Hverbose=FALSE)
        require('Hmisc')
        require('nlme')
        dd <- object$dims
        method <- object$method
        fixF <- object$call$fixed
        xtTab <- as.data.frame(object$tTable)
        sigp <- xtTab[,"p-value"]< shadep # cells that will be shaded
        if (!interceptp){
          sigp[1] <- FALSE # intercept will never be shaded
          # Replace small significances, discarding p-value for (Intercept)
          xtTab[1,"p-value"] = 1 # we do not show it anyway, easier formatting
        }
        pval <- format(zapsmall(xtTab[, "p-value"],4))
        pval[as.double(pval) < 0.0001] <- "$< .0001$"
        xtTab[, "p-value"] <- pval
        xtTab[,"t-value"] <- round(xtTab[,"t-value"],1)
        if (ncol(xtTab) == 5) # not for gls
          xtTab[,"DF"] <- as.integer(xtTab[,"DF"])
        # extract formula
        if (is.null(form)) {
          if (!is.null(object$terms)) {
            form=object$terms
          } else {
            form = formula(object)
          }
        }
        if (is.null(parameter)) {
          parameter=as.character(form[[2]])
        }
        if (any(wchLv <- (as.double(levels(xtTab[, "p-value"])) == 0))) {
            levels(xtTab[, "p-value"])[wchLv] <- "<.0001"
        }
        if (is.null(label))
          label <- lmeLabel("contr",form)
        form <- deparse(removeFormFunc(as.formula(form)),width.cutoff=500)
      
        form <- paste(sub('~','$\\\\sim$ ',form),sep="")
        # All I( in factors are replaced with (This could be improved)
        row.names(xtTab) <- 
          gsub("I\\(","(",dimnames(object$tTable)[[1]])
        row.names(xtTab) <-  gsub("\\^2","\\texttwosuperior",row.names(xtTab))
      
        # Determine base level  
        levs <- lapply(object$contrasts,function(object) {dimnames(object)[[1]][1]})
        levnames <- paste(names(levs),levs,sep=" = ",collapse=", ")
        # Try to locate numeric covariables
      #  v1 <- all.vars(formula(object))[-1]
      ## Changed 8.10.2008, not regression-tested
        v1 <- all.vars(form)[-1]
        numnames <- v1[is.na(match(v1,names(levs)))]
        if (length(numnames > 0)) {
          numnames <- paste(numnames," = 0",collapse=", ")
          levnames <- paste(levnames,numnames,sep=", ")
        }
        if (is.null(caption)){ # TODO: Allow %s substitution
          if (inherits(object,"lme"))
            md = "Mixed model (lme)" else
          if (inherits(object,"gls"))
            md = "Extended linear model (gls)" else
            md = "Linear model"
          caption <- paste(md," contrast table for \\emph{",
             parameter, "} (model ",form,
          "). The value in row (Intercept) gives the reference value for ",
            levnames,".",sep='')
        }
        caption.lot <- paste("Contrast table for ",parameter, " by ",
            levnames)
        ndec <- pmax(round(1-log10(xtTab[,2]+0.000001)+moredec),0)
        xtTab[,1] <- formatC(round(xtTab[,1],ndec))
        xtTab[,2] <- formatC(round(xtTab[,2],ndec))
        if (ncol(xtTab) == 5) {
          names(xtTab) <- c("Value","StdErr","DF","t","p")
          pcol = 5
        } else {# gls misuse
          names(xtTab) <- c("Value","StdErr","t","p")
          pcol = 4
        }
        # Only show intercept p/t when explicitely required
        if (!interceptp){
          xtTab[1,pcol-1] <- NA
          xtTab[1,pcol] <- ''
        }
        cellTex <- matrix(rep("", NROW(xtTab) * NCOL(xtTab)), nrow=NROW(xtTab))
        cellTex[sigp,pcol] <- "cellcolor[gray]{0.9}"
        rowlabel <- ifelse(nchar(parameter) >9,"",parameter)
        latex(xtTab, title=title, file=file, caption=caption,caption.lot=caption.lot,
          caption.loc="bottom", label=label, cellTexCmds = cellTex,
          rowlabel=rowlabel, ctable=ctable, where=where,
          booktabs = !ctable, numeric.dollar=FALSE,col.just=rep("r",5),...)
      }
      
      "latex.lme" <-
      function(object, title="",parameter=NULL,file="",shadep=0.05,
        caption=NULL,label=NULL,ctable=FALSE,form=NULL,
        interceptp=FALSE,  moredec= 0, where="!htbp",...) {
        options(Hverbose=FALSE)
        require('Hmisc')
        require('nlme')
        latex.summary.lme(summary(object),title=title,parameter=parameter, 
          file=file, shadep=shadep, caption=caption,
          label=label, ctable=ctable, form=form, moredec=moredec, where=where,...)
      }
      

      【讨论】:

      • 谢谢,迪特!这真的很有帮助。
      • 我在使用您的功能时遇到错误。你可以帮帮我吗?错误是'[.data.frame(xtTab, , "p-value") 中的错误:选择了未定义的列'
      • 会不会是语言问题?检查在您的语言环境中是否有一个列“p-value”当你这样一个摘要(lme(....))
      【解决方案4】:

      这是我的解决方案:假设 fit 是您的 lme 模型的结果,例如fit &lt;- lme(...)。如果您想让summary(fit) 显示所有变量,您只需键入

      > fit_text <- unclass(fit)
      > attributes(fit_text)
      

      你会看到类似结构的结果。然后,您可以将汇总报告的某些组件保存到 txt 文件或 Rdata 文件中。

      【讨论】:

        猜你喜欢
        • 1970-01-01
        • 2015-09-27
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        • 1970-01-01
        • 2018-01-01
        • 1970-01-01
        • 1970-01-01
        相关资源
        最近更新 更多