0

之前已经询问过这个问题的变体(例如,有没有办法从 `glmer` 对象中获得“边际效应”),其中大多数建议使用 ggeffects(或 sjPlot)。然而,我们学院的一位统计学家在同意 ggeffects 的输出时遇到了一些麻烦。

我对统计和 R 都是新手(所以我可能不应该使用 glmer ......),而且我在理解我做错了什么时遇到了一些麻烦。

我的模型最简单的形式是来自 lme4 的 glmer:结果(二进制)~ FE(二进制)+(1|RE)。固定效应是对我的随机效应中的一些但不是所有个体进行的测试。

模型输出

Intercept: 1.2654

FE_2: -0.2305

RE Std. Dev.: 2.896

gg预测:

ggpredict(model, type = "fe", terms = "FE")

1: 0.780

2: 0.738

现在,据我所知,我可以得到测试 2 的概率边际效应,如下所示:

y <- 1.2654 - 0.2305
prob <- exp(y) / (1 + exp(y))

这与 ggpredict 的输出完全相同:0.738。

然而,他说这是每个人的条件概率,我需要“插入一些我不理解的东西”来获得我可以推广到人群的概率。他的解决方案快速示例:

y <- 1.2654 - 0.2305 + rnorm(100000000) * 2.896 
prob_trus <- exp(y) / (1 + exp(y))
mean(prob_trus)
0.62

使用与前面的 ggpredict 具有相同参数的 ggaverage 几乎给了我相同的概率,0.639。所以,最后我的问题:

  • ggaverage 是否与他的解决方案相同,只是有一些模拟差异和/或舍入误差?
  • 如果他的方法是要走的路,我如何从 ggeffects 获得相同的结果?
4

2 回答 2

0

实际上,您从具有非线性链接函数的混合模型(例如,混合效应逻辑回归)中获得的系数通常具有以随机效应为条件的解释。

赫德克尔等人。(2018 年)最近提出了一种通过边际/总体解释获得回归系数的新想法。这是在marginal_coefs()R 包GLMMadaptive的函数中实现的,该包使用自适应高斯正交拟合混合模型。混合效应逻辑回归的一般示例是:

library("GLMMadaptive")
fm <- mixed_model(fixed = y ~ x1 + x2 + x3, random = ~ x1 | id, data = DF, family = binomial())

marginal_coefs(fm)
marginal_coefs(fm, std_errors = TRUE)

有关其他示例,例如,如何获得边际预测或产生效果图,请查看小插图:https ://drizopoulos.github.io/GLMMadaptive/articles/Methods_MixMod.html

于 2018-08-27T19:00:18.663 回答
0

ggaverage()与“我不明白的东西”-方法不同,类似的结果是偶然的。您可以在小插图ggpredict()中阅读有关和之间的差异的信息。ggaverage()

我不会说应该将随机效应方差(您的同事也希望以预测为条件的影响)添加到预测值中以计算预测概率。相反,这种不确定性应该在预测的标准误差(以及由此产生的置信区间)中考虑。

ggpredict(model, type = "fe", terms = "FE")与 比较时,您会看到差异ggpredict(model, type = "re", terms = "FE")。注意 changed type = "re",它现在考虑了随机效应,导致不同的预测值和更大的置信区间。

这是包中示例数据的示例:

data(efc_test)

fit <- glmer(
  negc7d ~ c12hour + e42dep + c161sex + c172code + (1 | grp),
    data = efc_test,
    family = binomial(link = "logit")
  )

ggpredict(fit, "c12hour", type = "fe")

#> # Predicted probabilities for Negative impact with 7 items 
#> # x = average number of hours of care per week 
#> 
#>   x predicted conf.low conf.high
#>   0     0.342    0.249     0.449
#>   5     0.344    0.251     0.450
#>  10     0.346    0.254     0.452
#>  15     0.348    0.256     0.453
#>  20     0.350    0.258     0.455
#>  25     0.352    0.260     0.456
#>  30     0.354    0.261     0.458
#>  35     0.356    0.263     0.460
#>  40     0.357    0.265     0.463
#>  45     0.359    0.266     0.465
#>  ... and 25 more rows.
#> 
#> Adjusted for:
#> *   e42dep = 2.92
#> *  c161sex = 1.76
#> * c172code = 1.97

ggpredict(fit, "c12hour", type = "re")

#> # Predicted probabilities for Negative impact with 7 items 
#> # x = average number of hours of care per week 
#> 
#>   x predicted conf.low conf.high
#>   0     0.472    0.107     0.870
#>   5     0.475    0.108     0.871
#>  10     0.477    0.109     0.872
#>  15     0.479    0.110     0.873
#>  20     0.481    0.111     0.873
#>  25     0.483    0.111     0.874
#>  30     0.485    0.112     0.875
#>  35     0.487    0.113     0.876
#>  40     0.489    0.114     0.877
#>  45     0.491    0.115     0.878
#>  ... and 25 more rows.
#> 
#> Adjusted for:
#> *   e42dep = 2.92
#> *  c161sex = 1.76
#> * c172code = 1.97
于 2018-08-27T09:45:23.123 回答