Как происходит регрессия экстракта из lme, lmer, glmer в Latex? - PullRequest
7 голосов
/ 23 февраля 2012

Подгоняю модели с 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))

Любая помощь?

спасибо

Ответы [ 4 ]

7 голосов
/ 13 мая 2013

Я только что обнаружил, что существует coef метод для summary.mer объектов, который предоставляет все необходимые данные (для фиксированных эффектов).Возвращенный объект (после приведения к data.frame) может быть легко передан выбранному пакету форматирования (например, xtable или ascii).
См. Следующий пример (который дает только пригодное для использования 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
6 голосов
/ 26 февраля 2012

Редактировать:

Во время редактирования пакет lme4 обновился, и memisc больше не работает с этими объектами.Пакет техрег является альтернативой.Я оставил этот ответ на тот случай, если memisc обновится и он снова начнет работать.

Пакет memisc выполняет таблицы lme4:

Вот фрагмент кода, который яwrote:

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 - используются по умолчаниюих в преамбуле латекса или отключите их с помощью команд из файла справки (useBooktabs = FALSE, useDcolumn = FALSE).

3 голосов
/ 23 февраля 2012

Для lme моя личная версия ниже; Вы можете скачать его с другими подобными дополнениями, например, извлечь \Sexpr{} строку для p-значений таблиц lme / lm / glm как Dmisc из

http://www.menne -biomed.de / загрузить

Это очень персонализировано, но если мне нравится округление до действительно значащих цифр. Извините, пакет 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,...)
}
0 голосов
/ 09 марта 2012

Вот мое решение: предположим, fit - это результат вашей модели, например. fit <- lme(...). Если вы хотите, чтобы все переменные отображались summary(fit), вы можете просто набрать

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

и вы увидите результат, подобный структуре. Затем вы можете сохранить определенные компоненты сводного отчета в текстовый файл или файл Rdata.

Добро пожаловать на сайт PullRequest, где вы можете задавать вопросы и получать ответы от других членов сообщества.
...