如何批量执行多元回归、提取交互项结果并生成可读表格

需求背景 #

需要对数据集内多个变量批量执行数百个多元回归模型，手动操作效率极低。仅需提取p:wt交互项的统计结果，同时为每行添加对应变量（血液蛋白）标识，最终整理为易读的汇总表格。以下以mtcars数据集为示例（真实场景中p对应约300种血液蛋白，示例中对应mtcars第2至5列）。

原始代码及输出如下：

library(tidyverse)

models <- lapply(mtcars[2:5], function(p) 
  (summary(lm(mpg ~ p * wt + qsec, data = mtcars))))
models_stats <- lapply(models, tidy) 
models_table <- do.call(rbind, models_stats)
models_table

原始输出：

# A tibble: 20 × 5
   term        estimate std.error statistic     p.value
 * <chr>          <dbl>     <dbl>     <dbl>       <dbl>
 1 (Intercept)  41.6      8.20        5.08  0.0000247  
 2 p            -3.39     0.963      -3.52  0.00154    
 3 wt          -10.8      2.40       -4.51  0.000113   
 4 qsec          0.757    0.349       2.17  0.0390     
 5 p:wt          0.977    0.317       3.08  0.00473    
 6 (Intercept)  30.4      5.82        5.22  0.0000167  
 7 p            -0.0378   0.0138     -2.73  0.0109     
 8 wt           -7.61     1.26       -6.04  0.00000188 
 9 qsec          0.780    0.291       2.68  0.0123     
10 p:wt          0.0106   0.00298     3.55  0.00143    
11 (Intercept)  40.3      7.68        5.25  0.0000156  
12 p            -0.106    0.0263     -4.04  0.000395   
13 wt           -8.68     1.29       -6.72  0.000000328
14 qsec          0.503    0.361       1.39  0.174      
15 p:wt          0.0278   0.00730     3.81  0.000733   
16 (Intercept)  -7.79    11.3        -0.692 0.495      
17 p             7.52     2.81        2.68  0.0124     
18 wt            2.80     3.21        0.873 0.390      
19 qsec          0.874    0.246       3.55  0.00142    
20 p:wt         -2.12     0.926      -2.29  0.0301   

需要解决三个问题：

• a) 移除主效应行，仅保留p:wt交互项行
• b) 为每行添加对应变量（蛋白名）标识
• c) 整理为易读的汇总表格

解决方案 #

使用tidyverse工具链优化代码，一次性完成模型拟合、交互项提取、变量标识添加及表格整理：

library(tidyverse)

# 批量处理：拟合模型、提取交互项、添加变量名并整理
interactions_summary <- mtcars[2:5] %>%
  imap_dfr(function(p_data, var_name) {
    # 拟合多元回归模型
    lm(mpg ~ p_data * wt + qsec, data = mtcars) %>%
      # 提取模型统计量
      tidy() %>%
      # 仅保留p:wt交互项
      filter(term == "p:wt") %>%
      # 添加当前变量（蛋白）名称
      mutate(protein_name = var_name) %>%
      # 调整列顺序，将变量名列放在最前面
      relocate(protein_name, .before = everything())
  })

# 查看最终汇总表格
interactions_summary

最终输出（易读表格） #

关键代码说明 #

• imap_dfr：同时遍历变量的值和名称，自动将结果行绑定为tibble，替代繁琐的lapply+do.call(rbind)组合
• filter(term == "p:wt")：精准筛选出仅有的交互项行，剔除主效应和截距项
• mutate(protein_name = var_name)：将当前遍历的变量名（即蛋白名）添加为新列，明确对应关系
• relocate：调整列的顺序，让变量名列处于最前端，提升表格可读性

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