如何在dplyr的map函数中结合relevel提取回归系数

核心思路 #

针对e的每个水平（E、M、H）分别重新设置因子参考水平，拟合m ~ g*e的线性回归，再批量提取系数估计值和标准误。关键是要在每个模型拟合前动态重置e的参考水平，避免因因子水平固定导致的错误。

具体实现代码 #

先加载所需工具包，构造示例数据（替换成你的真实数据集即可）：

# 加载包
library(dplyr)
library(purrr)
library(broom)
library(tibble)

# 构造示例数据集（替换为你的真实数据）
set.seed(123)
df <- tibble(
  g = factor(sample(c("F", "M"), 100, replace = TRUE)),  # 性别因子
  e = factor(sample(c("E", "M", "H"), 100, replace = TRUE)),  # 分类变量
  m = rnorm(100, 
            mean = 5 + (g == "M")*0.5 + (e == "M")*1 + (e == "H")*2 + 
              (g == "M" & e == "M")*0.3 + (g == "M" & e == "H")*0.6, 
            sd = 1)  # 模拟因变量
)

然后批量处理每个参考水平：

# 定义e的所有参考水平
ref_levels <- c("E", "M", "H")

# 批量拟合模型并提取系数信息
result <- tibble(ref_e = ref_levels) %>%
  # 针对每个参考水平，重置e的因子水平后拟合模型
  mutate(model = map(ref_e, ~ lm(m ~ g*e, data = df %>% mutate(e = relevel(e, ref = .x)))),
         # 提取模型系数的估计值和标准误
         coef_details = map(model, tidy)) %>%
  # 展开系数详情为表格
  unnest(coef_details) %>%
  # 保留关键列
  select(ref_e, term, estimate, std.error)

代码说明 #

1. relevel(e, ref = .x)：在每个模型拟合前，将e的参考水平动态设置为当前循环的.x值（即E/M/H中的一个），确保每个模型的参考组正确。
2. map()：批量处理每个参考水平，避免重复写三次拟合代码。
3. tidy()：来自broom包，自动将模型系数整理成整洁的表格，直接提取estimate（估计值）和std.error（标准误）。

结果示例 #

输出的result表格会按每个参考水平分组，展示对应模型的所有系数（包括截距、g的主效应、e的其他水平效应、g与e的交互项）及其估计值和标准误，和手动逐个拟合的结果完全一致。

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