GLMM与零膨胀混合模型中随机效应显著性检验方法咨询

1. 纵向研究纳入随机效应的研究依据 #

• 方法学层面：你使用的分道路、分年度事故数据属于典型的纵向面板计数数据，同一道路不同年份的观测不满足独立同分布假设，Cameron & Trivedi（1990，即你做过离散检验引用的经典方法提出者）在计数数据回归的权威专著中明确提出，这类存在组内相关的计数数据必须通过随机效应捕捉组水平未观测异质性，否则会导致标准误估计偏误、推断结果失效。广义线性混合模型框架下纳入组水平随机效应的方法，自Breslow & Clayton（1993）的经典研究发表后，就已经成为非独立计数数据建模的标准方案。
• 领域应用层面：道路交通事故分析中，不同道路的几何条件、区位属性、常年通行特征等不随时间变化的未观测因素会持续影响事故风险，这类异质性无法通过可观测协变量完全覆盖，纳入道路水平随机效应是交通领域事故计数建模的通用规范，相关应用在事故分析领域的主流研究中被广泛采用，不存在方法合理性争议。

2. 随机效应显著性的有效检验方法（AIC/BIC之外） #

注意：你之前尝试的Hausman检验仅用于对比固定效应与随机效应模型的系数一致性，无法用于检验随机效应本身是否需要纳入，确实不适用于当前场景
可用的规范检验方法包括：

• 修正边界似然比检验（LRT）：检验“随机效应方差=0”属于边界检验（原假设下方差参数落在参数空间边界），普通LRT的p值会偏保守，使用修正后的LRT对比严格嵌套的“无随机效应基准模型”和“带随机效应的混合模型”，是目前该场景下最常用的检验方法，若检验p<0.05，说明纳入随机效应后模型拟合度有显著提升。
• 参数自助法似然比检验：通过自助抽样得到经验p值，完全规避边界检验的p值偏误问题，结果比修正LRT更稳健，适合对结果可靠性要求高的场景。
• 随机效应方差置信区间检验：计算随机效应标准差/方差的剖面似然置信区间，若区间不包含0，即可判定随机效应统计显著，该方法不依赖大样本渐近假设，结果稳定。
  注意：不要使用普通Wald检验判断随机效应显著性，Wald检验在参数空间边界上的统计性质极差，结果不可靠

3. R中的实现工具包与调用方式 #

以下方法均适配你已经在使用的lme4、glmmTMB拟合流程，不需要额外调整模型设定：

3.1 基础修正似然比检验 #

直接使用你已安装的lme4、glmmTMB包内置的anova()方法即可，前提是待对比的两个模型严格嵌套（除待检验的随机效应外，固定效应、分布族、零膨胀设定、偏移项完全一致）。

• 针对lme4拟合的GLMM模型示例代码：

# 拟合无随机效应的基准NB模型（固定效应、偏移项与GLMM完全一致）
m_glm_nb <- MASS::glm.nb(accident_cnt ~ speed_limit + lane_num + offset(log(road_length)), data = crash_dat)
# 拟合带道路随机效应的GLMM NB模型
m_glmm_nb <- lme4::glmer.nb(accident_cnt ~ speed_limit + lane_num + offset(log(road_length)) + (1|road_id), data = crash_dat)
# 执行LRT检验
anova(m_glm_nb, m_glmm_nb, test = "LRT")

• 针对glmmTMB拟合的零膨胀混合模型示例代码：

# 拟合无随机效应的基准ZINB模型
m_zinb <- glmmTMB(accident_cnt ~ speed_limit + lane_num + offset(log(road_length)),
                  ziformula = ~1, family = nbinom2, data = crash_dat)
# 拟合带道路随机效应的ZINB混合模型
m_zinb_mm <- glmmTMB(accident_cnt ~ speed_limit + lane_num + offset(log(road_length)) + (1|road_id),
                     ziformula = ~1, family = nbinom2, data = crash_dat)
# 执行LRT检验
anova(m_zinb, m_zinb_mm, test = "Chisq")

3.2 参数自助法LRT（稳健检验） #

使用pbkrtest包实现，该包可直接适配lme4、glmmTMB拟合的模型对象，通过参数自助抽样得到经验p值，示例代码：

# 加载包，首次使用需先安装
library(pbkrtest)
# 执行自助法LRT，nsim为自助抽样次数，建议设置为1000及以上
boot_test_res <- PBmodcomp(m_zinb_mm, m_zinb, nsim = 1000)
# 输出检验结果，查看bootstrap经验p值即可
summary(boot_test_res)

3.3 随机效应置信区间检验 #

直接使用confint()方法计算随机效应的剖面似然置信区间即可，不需要额外安装工具包，示例代码：

# 检验glmer拟合模型的随机效应，parm指定提取随机效应参数，method指定用剖面似然法计算区间
confint(m_glmm_nb, parm = "theta_", method = "profile")
# 检验glmmTMB拟合模型的随机效应
confint(m_zinb_mm, parm = "ranef", method = "profile")

若输出的道路随机效应标准差的置信区间下限大于0，即可判定随机效应显著。

内容的提问来源于stack exchange，提问作者나수빈