glmer 输出的优势比和置信区间答案

【问题标题】：Odds ratio and confidence intervals from glmer outputglmer 输出的优势比和置信区间
【发布时间】：2014-12-12 13:32:13
【问题描述】：

我制作了一个模型，该模型查看了许多变量以及对妊娠结果的影响。结果是一个分组的二进制文件。一群动物将有 34 个怀孕和 3 个空，接下来将有 20 个怀孕和 4 个空，依此类推。

我使用 glmer 函数对这些数据进行了建模，其中 y 是妊娠结果（怀孕或空腹）。

mclus5 <- glmer(y~adg + breed + bw_start + year + (1|farm),
                data=dat, family=binomial)

我得到了所有带有系数等的常用输出，但为了解释，我想将其转换为每个系数的优势比和置信区间。

在过去的逻辑回归模型中，我使用了以下代码

round(exp(cbind(OR=coef(mclus5),confint(mclus5))),3)

这会很好地提供我想要的，但它似乎不适用于我运行的模型。

有谁知道我可以通过 R 为我的模型获取此输出的方法吗？

【问题讨论】：

键入：?fixef 和 class(mclus5) 和 showMethods("fixef", includeDefs=TRUE)。 glmer 是一个 S4 函数。

标签： r confidence-interval

【解决方案1】：

我认为其他选择是只使用包emmeans：

library(emmeans)
data.frame(confint(pairs(emmeans(fit, ~ factor_name,type="response"))))

【讨论】：

不过，无需将呼叫包装在 data.frame 中。事实上，这样做会抑制可能有用的消息。
哈，是的，这很好——不过，我通常会用它制作一个数据框以将其保存到 csv 文件中。但可能最好分两步完成...
请注意，这可能要切换到现在使用emmeans。
是的，当然现在也可以使用更新的 packae emmeans，好点！

【解决方案2】：

唯一真正的区别是您必须使用fixef() 而不是coef() 来提取固定效应系数（coef() 为您提供估计的系数每个组）。

我将使用 lme4 包中的内置示例进行说明。

library("lme4")
gm1 <- glmer(cbind(incidence, size - incidence) ~ period + (1 | herd),
               data = cbpp, family = binomial)

固定效应系数和置信区间，对数优势尺度：

cc <- confint(gm1,parm="beta_")  ## slow (~ 11 seconds)
ctab <- cbind(est=fixef(gm1),cc)

（如果您想要更快但不太准确的 Wald 置信区间，您可以改用 confint(gm1,parm="beta_",method="Wald")；这将等同于 @Gorka 的答案，但更方便一些。）

指数获得优势比：

rtab <- exp(ctab)
print(rtab,digits=3)
##               est 2.5 % 97.5 %
## (Intercept) 0.247 0.149  0.388
## period2     0.371 0.199  0.665
## period3     0.324 0.165  0.600
## period4     0.206 0.082  0.449

稍微简单/更通用的解决方案：

library(broom.mixed)
tidy(gm1,conf.int=TRUE,exponentiate=TRUE,effects="fixed")

对于 Wald 区间，或添加 conf.method="profile" 作为配置文件置信区间。

【讨论】：

【解决方案3】：

我相信还有另一种更快的方法（如果您可以接受不太准确的结果）。

发件人：http://www.ats.ucla.edu/stat/r/dae/melogit.htm

首先我们得到估计的置信区间

se <- sqrt(diag(vcov(mclus5)))
# table of estimates with 95% CI
tab <- cbind(Est = fixef(mclus5), LL = fixef(mclus5) - 1.96 * se, UL = fixef(mclus5) + 1.96 * se)

然后是 95% CI 的优势比

print(exp(tab), digits=3)

【讨论】：