多重插补后拟合多水平广义线性模型报错求助

核心问题根源 #

1. 模型语法误用：survey包的svyglm不支持(fld5|locations)这种混合效应（随机效应）语法，R会将|识别为逻辑或运算，导致‘|’ not meaningful for factors警告。
2. 插补方法错误：你对所有变量统一使用pmm（预测均值匹配）插补，该方法仅适用于连续变量。分类变量用此方法会导致类型异常、类别水平丢失，触发contrasts can be applied only to factors with 2 or more levels错误。
3. 插补后变量类型丢失：pmm插补会将分类因子转为数值型，或在部分插补数据集里让分类变量只剩单一水平，无法生成模型所需的对比矩阵。

具体解决步骤 #

1. 修正多重插补方法

针对不同变量类型指定对应插补策略：

• 连续变量（dm2、age）：保留pmm
• 二分类变量（gender）：用logreg
• 多分类变量（fld5、race、locations、types）：用polyreg

修改后的插补代码：

library(survey)
library(mitools)
library(mice)

# 为每个变量指定插补方法
method_vec <- rep("pmm", ncol(mydata))
names(method_vec) <- colnames(mydata)

# 二分类变量用logreg
method_vec["gender"] <- "logreg"
# 多分类变量用polyreg
multi_cat_vars <- c("fld5", "race", "locations", "types")
method_vec[multi_cat_vars] <- "polyreg"

# 执行插补
imp <- mice(mydata, seed = 123, m = 5, method = method_vec)
implog <- mice::complete(imp, action = "long", include = TRUE)

# 强制保留分类变量的原始因子水平
for(var in c(multi_cat_vars, "gender")) {
  implog[[var]] <- factor(implog[[var]], levels = levels(mydata[[var]]))
}

2. 修正模型语法（二选一）

根据你的分析需求选择对应方案：

方案A：同时处理多水平随机效应与调查权重

如果需要保留(fld5|locations)这类随机效应结构，改用lme4的glmer结合调查权重，再用MIcombine合并插补结果：

# 转换为插补列表
imp_list <- imputationList(split(implog, implog$.imp)[-1])

# 对每个插补数据集运行混合效应模型
models <- with(imp_list, 
               lme4::glmer(dm2 ~ fld5 + age + gender + race + (fld5|locations) + (fld5|types), 
                           family = binomial(link = "logit"),
                           weights = weight))

# 合并插补结果
model <- MIcombine(models)

方案B：仅处理调查设计（无随机效应）

如果(fld5|locations)是交互项的笔误，改用:表示交互项，继续使用svyglm：

# 设置调查设计
designs <- svydesign(id =~ locations, weights =~ weight, data = imp_list)

# 运行调查GLM（交互项用:表示）
model <- MIcombine(with(designs, 
                        svyglm(dm2 ~ fld5 + age + gender + race + fld5:locations + fld5:types, 
                               family = binomial(link = "logit"))))

3. 前置检查（可选）

运行模型前，检查每个插补数据集的分类变量水平是否完整：

# 遍历所有插补数据集，输出分类变量的水平数
for(i in 1:length(imp_list$imputations)) {
  cat("插补数据集", i, ":\n")
  for(var in c("fld5", "gender", "race", "locations", "types")) {
    cat(var, "水平数:", length(levels(imp_list$imputations[[i]][[var]])), "\n")
  }
}

若某变量在某数据集里水平数不足2，需调整插补方法或检查原始数据的类别分布。

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