分别绘制 glmmTMB 对零截断计数分量和二进制分量的影响

Question

我拟合了一个 Hurdle 混合模型（glmmTMB 包中的 glmmTMB 函数）来同时探索感染流行率（模型的二进制部分，零和非零数据）和感染强度（零截断负二项式模型，计数数据）如何响应环境变化。 我拟合了模型，使用 Anova 获得显着的固定效应，然后使用 emmeans 和对比度进行成对比较。对于每个步骤，我都能够获得二进制组件和计数组件的单独输出。 如何根据模型获得这两个组件的单独绘图？

作为示例，我使用 Owls 数据。

library("glmmTMB") 
library("emmeans")

data("Owls")

# fit model
mod <- glmmTMB(SiblingNegotiation ~ FoodTreatment + Nest + offset(log(BroodSize)), 
                       zi= ~ FoodTreatment + Nest, 
                       family=truncated_nbinom2, data=Owls)

# Anova for significant fixed effects
#count part
car::Anova(mod, type="II", component="cond") 
#binary part
car::Anova(mod, type="II", component="zi")

# pairwise comparisons
#binary part
emm.nest.zi <- emmeans::emmeans(mod, c("Nest"), component="zi")
contrast(emm.nest.zi, alpha=0.05, method="pairwise", adjust="bh")
#count part
emm.nest.nb <- emmeans::emmeans(mod, c("Nest"), component="cond")
contrast(emm.nest.nb, alpha=0.05, method="pairwise", adjust="bh")

现在，我想获得两个基于模型的图。一个用于模型的二进制部分，另一个用于模型的计数部分。这是我的尝试。

# plot the count part of the glmmTMB model
predict_nest.nb <- ggeffect(mod,"Nest", component="cond")
ggplot(predict_nest.nb,aes(x=x, y=predicted), fill=group, group=group, color=group)+
  geom_errorbar(data=predict_nest.nb, 
                mapping=aes(x=x, ymin=conf.low, ymax=conf.high), 
                width=0.1, size=1) +
  geom_point(size=4, aes()) +
  theme(axis.text.x=element_text(size=10, angle=90, vjust=1, hjust=1, face="plain"))

# plot the binary part of the glmmTMB model
predict_nest.zi <- ggeffect(mod,"Nest", component="zi")
ggplot(predict_nest.zi,aes(x=x, y=predicted), fill=group, group=group, color=group)+
  geom_errorbar(data=predict_nest.nb, 
                mapping=aes(x=x, ymin=conf.low, ymax=conf.high), 
                width=0.1, size=1) +
  geom_point(size=4, aes()) +
  theme(axis.text.x=element_text(size=10, angle=90, vjust=1, hjust=1, face="plain"))

但是我得到了两个相同的图，尽管指定了组件=“zi”或“cond”。 在我原来的 df 中，成对对比显示零部分存在一些显着差异，这与模型的计数部分不同。所以，这两个图应该看起来不同。

为什么不起作用？如何分别绘制 glmmTMB 模型的两个组成部分？

感谢您的任何意见=）

Answer 1

emmeans 包本身的绘图函数可以很好地完成这种事情......

require("emmeans")

plot(pairs(emm.nest.zi))
plot(pairs(emm.nest.nb))