2

我在 R 中有一个glmer模型,我想为其绘制预测。我从 sjPlot 库中找到了 plot_model 函数,它工作正常。

这是一个MWE:

library(lattice)

cbpp$response <- sample(c(0,1), replace=TRUE, size=nrow(cbpp))
gm1 <- glmer(response ~ size + incidence + (1 | herd),
              data = cbpp, family = binomial)

例如,调用plot_model(gm1, type = "pred", show.data = TRUE)会产生下图:

在此处输入图像描述

但是,我对 R 并不熟悉,而且我很难控制情节美学并将多个模型绘制到同一个图中(已经在此处提出了关于该问题的问题)。我熟悉 Python 和 matplotlib,让这些数字在 Python 环境中工作对我来说要简单得多。

我猜想实现这一点的一种方法是从 R 中获取 y 值(预测的火灾概率)并导出它们,以便我可以在 Python 中读取它们,以便在本例中针对每个协变量(evi prev)绘制它们。但是,我不知道该怎么做。此外,我试图阅读sjPlot源代码以弄清楚它是如何绘制预测的,但也无法弄清楚。

4

2 回答 2

2

ggpredict()实际上返回更多的值(并且沿着x轴,即对于有问题的术语 -size在你的例子中 - 这些是偶数间隔的),但只打印更少的值。

library(lme4)
#> Loading required package: Matrix
library(ggeffects)

cbpp$response <- sample(c(0,1), replace=TRUE, size=nrow(cbpp))
gm1 <- glmer(response ~ size + incidence + (1 | herd), data = cbpp, family = binomial)

pr1 <- ggpredict(gm1, term = "size")

pr1
#> 
#> # Predicted probabilities of response
#> # x = size
#> 
#>   x predicted std.error conf.low conf.high
#>   2     0.632     0.717    0.297     0.875
#>   6     0.610     0.550    0.347     0.821
#>  10     0.587     0.407    0.390     0.759
#>  14     0.563     0.321    0.407     0.708
#>  18     0.539     0.339    0.376     0.695
#>  22     0.515     0.448    0.306     0.719
#>  26     0.491     0.601    0.229     0.758
#>  34     0.444     0.951    0.110     0.837
#> 
#> Adjusted for:
#> * incidence = 1.77
#> *      herd = 0 (population-level)
#> Standard errors are on link-scale (untransformed).

as.data.frame(pr1)
#>     x predicted std.error  conf.low conf.high group
#> 1   2 0.6323758 0.7168742 0.2967912 0.8751705     1
#> 2   4 0.6211339 0.6316777 0.3221952 0.8497229     1
#> 3   6 0.6097603 0.5501862 0.3470481 0.8212222     1
#> 4   8 0.5982662 0.4743133 0.3701925 0.7904902     1
#> 5  10 0.5866630 0.4072118 0.3898523 0.7592017     1
#> 6  12 0.5749627 0.3539066 0.4033525 0.7302266     1
#> 7  14 0.5631779 0.3213384 0.4071542 0.7076259     1
#> 8  16 0.5513213 0.3159857 0.3981187 0.6953669     1
#> 9  18 0.5394060 0.3391396 0.3759558 0.6947993     1
#> 10 20 0.5274456 0.3857000 0.3438768 0.7038817     1
#> 11 22 0.5154536 0.4484344 0.3063836 0.7192510     1
#> 12 24 0.5034437 0.5215385 0.2672889 0.7380720     1
#> 13 26 0.4914299 0.6012416 0.2292244 0.7584368     1
#> 14 28 0.4794260 0.6852450 0.1938167 0.7791488     1
#> 15 30 0.4674458 0.7721464 0.1619513 0.7994688     1
#> 16 32 0.4555030 0.8610687 0.1339908 0.8189431     1
#> 17 34 0.4436111 0.9514457 0.1099435 0.8373008     1

reprex 包(v0.2.1)于 2019 年 5 月 6 日创建

有一些小插曲显示了包的不同特性,这里演示了如何在焦点项的特定值/级别计算边际效应。

Ben 发布的配方展示了如何计算置信区间(以随机效应为条件或不以随机效应为条件)在 中实现,这里ggpredict()有一个简短的小插曲来解释差异。

于 2019-05-06T13:45:06.200 回答
2

最简单的方法可能是使用ggeffects::ggpredict().

就像是

library(ggeffects)
pred_frame <- ggpredict(myModel, term="evi_prev")

应该产生一个具有预测、上下置信水平的数据框。我不确定它是否会预测沿 x 轴均匀分布的值(这会很好),或者如何欺骗它这样做。(如果你提供一个可重现的例子,我可能会试一试。)

玩弄您发布的 MWE 确实表明很难获得均匀间隔值的预测(或更一般地说,对于不在原始数据中的值);我尝试了类似的东西terms="size [1:35]",但这限制了预测值的范围,而不是填充它们。

更简单地说,可以使用对象的内置predict()方法(可能用于指定例如均匀间隔的值)来获得预测[用于获得关于概率而不是对数赔率尺度的预测];置信区间更难,但可以使用此处显示的配方生成merModnewdatatype="response"

于 2019-05-04T00:42:10.983 回答