在R中针对多重插补数据集按分组重复执行混合效应回归

针对你需要在50个多重插补数据集上按cyl分组拟合lmer模型的需求，这里提供两种可行的实现方式，解决tidyverse与with.mice的兼容问题：

方法1：用tidyverse遍历插补数据集并分组建模 #

直接用purrr遍历每个插补数据集，结合dplyr分组嵌套后拟合模型，最后按Rubin规则合并多重插补结果：

# 加载必要包
library(mice)
library(lme4)
library(tidyverse)
library(broom.mixed)

# 模拟带缺失值的mtcars并生成多重插补数据集（对应你的implist）
set.seed(123)
mtcars_miss <- mtcars
mtcars_miss$mpg[sample(1:nrow(mtcars), 10)] <- NA
implist <- mice(mtcars_miss, m = 50, printFlag = FALSE)

# 定义单插补数据集的分组建模函数
fit_grouped <- function(data) {
  data %>%
    group_by(cyl) %>%
    nest() %>%
    mutate(
      # 对每个cyl分组拟合lmer模型
      model = map(data, ~ lmer(mpg ~ disp + hp + (1|gear), data = .x)),
      # 提取固定效应结果
      tidied = map(model, tidy, effects = "fixed")
    ) %>%
    unnest(tidied) %>%
    select(cyl, term, estimate, std.error)
}

# 遍历所有插补数据集，得到每个插补的分组模型结果
imputed_res <- map(implist$imp, fit_grouped)

# 合并所有插补结果并按Rubin规则计算合并统计量
pooled_res <- bind_rows(imputed_res, .id = "imp") %>%
  group_by(cyl, term) %>%
  summarise(
    pooled_est = mean(estimate),
    # Rubin规则：合并插补内方差与插补间方差
    pooled_se = sqrt(mean(std.error^2) + var(estimate) * (1 + 1/n())),
    statistic = pooled_est / pooled_se,
    p_value = 2 * pt(abs(statistic), df = n() - 1, lower.tail = FALSE)
  )

# 查看最终合并结果
print(pooled_res)

方法2：结合mice的with函数与by分组 #

如果习惯用with.mice，可以通过by函数在每个插补数据集中直接按cyl分组建模：

# 用with+by在每个插补数据集上按cyl分组拟合模型
grouped_models <- with(implist, 
                       by(data = .data, INDICES = cyl, 
                          FUN = function(x) lmer(mpg ~ disp + hp + (1|gear), data = x)))

# 整理所有插补的模型结果
tidied_res <- map(grouped_models, function(imp_mods) {
  map_dfr(imp_mods, tidy, effects = "fixed", .id = "cyl")
})

# 合并并pool结果（同方法1的Rubin规则计算）
pooled_res2 <- bind_rows(tidied_res, .id = "imp") %>%
  group_by(cyl, term) %>%
  summarise(
    pooled_est = mean(estimate),
    pooled_se = sqrt(mean(std.error^2) + var(estimate) * (1 + 1/n())),
    statistic = pooled_est / pooled_se,
    p_value = 2 * pt(abs(statistic), df = n() - 1, lower.tail = FALSE)
  )

关键注意事项 #

• 确保每个cyl分组内有足够样本（尤其是gear的水平数），避免混合效应模型无法收敛
• 多重插补结果必须用Rubin规则合并，不能直接取系数均值，否则会低估标准误
• broom.mixed包专门用于整理混合效应模型的结果，比基础broom更适配lmer

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