5

过去,我曾使用sjPlotsjp.glmer包中的sjPlot来可视化广义混合效应模型中的不同斜率。但是,使用新软件包,我无法弄清楚如何绘制各个斜率,如图所示(随机)组级别的固定效应概率,位于此处

这是我认为应该允许生成图形的代码。我似乎无法在新版本的sjPlot.

library(lme4)
library(sjPlot)
data(efc)
# create binary response
efc$hi_qol = 0
efc$hi_qol[efc$quol_5 > mean(efc$quol_5,na.rm=T)] = 1
# prepare group variable
efc$grp = as.factor(efc$e15relat)
# data frame for 2nd fitted model
mydf <- na.omit(data.frame(hi_qol = as.factor(efc$hi_qol),
                           sex = as.factor(efc$c161sex),
                           c12hour = as.numeric(efc$c12hour),
                           neg_c_7 = as.numeric(efc$neg_c_7),
                           grp = efc$grp))
# fit 2nd model
fit2 <- glmer(hi_qol ~ sex + c12hour + neg_c_7 + (1|grp),
              data = mydf,
              family = binomial("logit"))

我尝试使用以下代码绘制模型。

plot_model(fit2,type="re")
plot_model(fit2,type="prob")
plot_model(fit2,type="eff") 

我认为我可能缺少一个标志,但是在阅读了文档之后,我无法找出那个标志可能是什么。

4

2 回答 2

6

看起来这可能会做你想要的:

(pp <- plot_model(fit2,type="pred",
       terms=c("c12hour","grp"),pred.type="re"))
  • type="pred":绘制预测值
  • terms=c("c12hour", "grp"):包括c12hour(作为x轴变量)和grp预测
  • pred.type="re": 随机效应

我还没有得到置信区间丝带(试过了ci.lvl=0.9,但没有运气......)

pp+facet_wrap(~group)更接近链接的博客文章中显示的情节(每个随机效应级别都有自己的方面......)

在此处输入图像描述

于 2018-12-19T17:05:45.237 回答
3

本已经发布了正确答案。sjPlot 使用ggeffects-package来绘制边际效应图,因此另一种方法是直接使用 ggeffects:

ggpredict(fit2, terms = c("c12hour", "grp"), type="re") %>% plot()

一个新的小插图描述了如何获得混合模型/随机效应的边际效应。但是,此绘图类型目前不提供置信区间。

链接的博客文章中的type = "ri.prob"选项没有针对协变量进行调整,这就是为什么我首先删除了该选项,然后在 ggeffects / sjPlot 中(正确地)重新实现了它。链接的博客文章中显示的置信区间也不正确。一旦我想出了如何获得 CI 或预测区间的方法,我也会添加这个选项。

于 2018-12-22T14:44:38.893 回答