glmmTMB零膨胀负二项混合模型的R平方计算方法咨询

问题背景 #

你基于glmmTMB构建了零膨胀负二项混合模型：

M2<- glmmTMB(psychological100~ (1|ID) + time*MNM01, data=mnmlong, ziformula=~ (1|ID) + time*MNM01, family=nbinom2())
summary(M2)

尝试用以下两种方法计算R²均失败：

MuMIn::r.squaredGLMM(M2)

报错：

Error in r.squaredGLMM.glmmTMB(M2) : r.squaredGLMM cannot (yet) handle 'glmmTMB' object with zero-inflation

performance::r2_zeroinflated(M2)

报错：

Error in residuals.glmmTMB(model, type = "pearson") : pearson residuals are not implemented for models with zero-inflation or variable dispersion

可行解决方案 #

1. 手动近似计算边际/条件R²

对于零膨胀混合模型，我们可以分别针对计数部分和零膨胀部分计算近似的R²，或者计算整体的拟合相关系数平方：

# 提取观测值和计数部分的拟合值（包含随机效应）
y <- M2$frame$psychological100
mu_full <- predict(M2, type = "response")

# 计算条件R²（固定+随机效应解释的变异）
r2_conditional <- cor(y, mu_full)^2

# 计算边际R²（仅固定效应解释的变异）
mu_fixed <- predict(M2, type = "response", re.form = NA)
r2_marginal <- cor(y, mu_fixed)^2

# 针对零膨胀部分，计算观测零值与预测零概率的相关平方
y_zi <- as.numeric(y == 0)
zi_prob <- predict(M2, type = "zprob")
r2_zi <- cor(y_zi, zi_prob)^2

# 输出结果
cat("条件R²（计数部分）：", round(r2_conditional, 3), "\n")
cat("边际R²（计数部分）：", round(r2_marginal, 3), "\n")
cat("零膨胀部分R²：", round(r2_zi, 3), "\n")

这种方法基于皮尔逊相关系数的平方，是一种直观的近似，能快速反映模型对数据的解释能力。

2. 使用DHARMa包通过模拟残差评估拟合优度

虽然无法直接得到传统R²，但DHARMa可以通过模拟残差全面评估零膨胀混合模型的拟合情况，间接反映模型的表现：

library(DHARMa)
# 生成模拟残差
sim_res <- simulateResiduals(M2)

# 可视化残差分布，检查拟合情况
plot(sim_res)

# 执行分散性检验和零膨胀检验
testDispersion(sim_res)
testZeroInflation(sim_res)

通过残差的QQ图、分布直方图，以及检验结果，可以判断模型是否充分捕捉了数据的变异和零膨胀特征，这是零膨胀模型拟合评估的常用手段。

3. 等待工具包版本更新

从报错信息可以看到，MuMIn和performance目前尚未支持带零膨胀的glmmTMB模型的R²计算。你可以关注这两个包的GitHub仓库或CRAN更新日志，后续版本可能会添加相关功能支持。

内容的提问来源于stack exchange，提问作者Sisay