glmmTMB拟合零膨胀混合GLM执行Bootstrap报错如何解决

• glmmTMB 模型的coef()、fixef()返回值不是普通lm/glm模型输出的一维数值向量，而是分为条件分布子模型、零膨胀子模型、离散参数子模型三个部分的嵌套列表。初始化结果矩阵时直接用length(coef(m_F))、names(coef(m_F))，拿到的是列表顶层的元素属性，和实际需要提取的系数长度、名称完全不匹配，赋值时会触发维度不匹配错误。
• 有放回的Bootstrap重抽样很容易抽到随机效应水平缺失、零值占比异常的子样本，导致重拟合模型失败、返回系数长度和预期不符，也会触发同类下标报错。

1. 自定义统一的系数提取函数，将不同子模型的固定效应拼接为标准一维向量，同时给系数名加前缀区分所属子模型，避免名称混淆。
2. 用自定义提取函数的输出初始化Bootstrap结果矩阵，保证列数、列名和实际提取的系数完全对应。
3. 循环中增加容错判断，遇到拟合失败、系数长度不匹配的重抽样样本直接跳过，避免整个循环中断。
4. 标准非参数Bootstrap需要将重抽样样本量设置为和原数据集一致，若做子抽样可保留原设置的bootsize=100。

修正后的可运行代码如下：

library(glmmTMB)
library(boot)
my.ds <- read.csv("https://raw.githubusercontent.com/Leprechault/trash/main/ds.desenvol.csv")

# 拟合原模型
m_F <- glmmTMB(development ~ poly(temp,2) + (1 | storage), data = my.ds,
               family = ziGamma(link = "log"),
               ziformula = ~ 1)

# 自定义固定效应提取函数：拼接三个子模型的固定效应为一维向量
get_fixef <- function(model) {
  fe_list <- fixef(model)
  c(
    setNames(fe_list$cond, paste0("cond_", names(fe_list$cond))),
    setNames(fe_list$zi, paste0("zi_", names(fe_list$zi))),
    setNames(fe_list$disp, paste0("disp_", names(fe_list$disp)))
  )
}

# 初始化bootstrap结果矩阵
nboot <- 1000
true_fe <- get_fixef(m_F)
bres <- matrix(NA, nrow = nboot, ncol = length(true_fe),
               dimnames = list(rep = seq(nboot), coef = names(true_fe)))

set.seed(1000)
# 标准bootstrap抽样量和原数据一致，若需子抽样可修改为100
bootsize <- nrow(my.ds) 

for (i in seq(nboot)) {
  # 重抽样本
  bdat <- my.ds[sample(nrow(my.ds), size = bootsize, replace = TRUE),]
  # 容错拟合，遇到报错不中断循环
  bfit <- try(update(m_F, data = bdat), silent = TRUE)
  # 判断是否拟合成功，再提取系数赋值
  if (!inherits(bfit, "try-error")) {
    current_fe <- try(get_fixef(bfit), silent = TRUE)
    # 只有系数长度匹配时才写入结果
    if (!inherits(current_fe, "try-error") && length(current_fe) == ncol(bres)) {
      bres[i,] <- current_fe
    }
  }
}

# 去除拟合失败的空行，得到最终bootstrap结果
bres_final <- na.omit(bres)

如果需要提取包含随机效应的分组水平系数（即coef()默认返回的各随机效应水平对应系数），只需要修改自定义提取函数的逻辑，把对应输出拼接为长度固定的一维向量，保证和结果矩阵列数一一对应即可。

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