如何通过for循环提取lm模型的系数、P值及显著性？

问题描述 #

假设已有如下线性回归模型：

set.seed(416)
N = 100
X  = data.frame(
  x1 = rnorm(N,0,1),
  x2 = rnorm(N,0,1)
)
X$y <- 2*X$x1 + 3*X$x2 + rnorm(N)

model <- lm(y ~ x1 + x2, data=X)

需求是为每个自变量创建包含Estimate（估计值）、**P-Value（P值）和Sig（判断p<0.05）**的命名列表，原手动方式效率较低：

# 低效的手动方式
model_data <- c(
         est_x1= coef(summary(model))['x1', 'Estimate'], 
         p_x1=coef(summary(model))['x1', 'Pr(>|t|)'], 
         sig_x1=coef(summary(model))['x1', 'Pr(>|t|)'] < 0.05, 
         est_x2= coef(summary(model))['x2', 'Estimate'], 
         p_x2= coef(summary(model))['x2', 'Pr(>|t|)'], 
         sig_x2=coef(summary(model))['x2', 'Pr(>|t|)'] < 0.05
         )

输出结果：

est_x1         p_x1       sig_x1       est_x2         p_x2       sig_x2 
1.947646e+00 1.204075e-34 1.000000e+00 3.008251e+00 1.158695e-51 1.000000e+00

提问：是否可以通过for循环实现该功能，以便后续轻松添加更多变量？

解决方案 #

方法1：使用for循环实现 #

先提取模型摘要的系数表，再遍历每个自变量（排除截距项），自动生成所需指标及命名：

# 提取模型系数摘要表
coef_table <- coef(summary(model))
# 获取所有自变量名称（剔除截距项）
vars <- rownames(coef_table)[-1]

# 初始化空向量用于存储结果
model_data_loop <- c()

# 遍历每个自变量
for(var in vars){
  # 提取当前变量的三个指标
  est_val <- coef_table[var, "Estimate"]
  p_val <- coef_table[var, "Pr(>|t|)"]
  sig_flag <- as.integer(p_val < 0.05)
  
  # 为指标命名并添加到结果向量
  model_data_loop <- c(model_data_loop,
                       setNames(est_val, paste0("est_", var)),
                       setNames(p_val, paste0("p_", var)),
                       setNames(sig_flag, paste0("sig_", var))
                       )
}

# 查看最终结果
model_data_loop

运行后输出与手动方式完全一致，后续新增变量时只需修改模型公式，循环代码无需调整。

方法2：高效向量化实现（无需循环） #

如果追求更简洁的代码，可使用向量化操作直接生成结果，效率更高：

# 提取自变量对应的系数行（排除截距项）
coef_table <- coef(summary(model))[-1, ]
# 生成所有指标的命名
name_list <- c(paste0("est_", rownames(coef_table)),
               paste0("p_", rownames(coef_table)),
               paste0("sig_", rownames(coef_table)))
# 提取所有指标的数值
value_list <- c(coef_table[, "Estimate"],
                coef_table[, "Pr(>|t|)"],
                as.integer(coef_table[, "Pr(>|t|)"] < 0.05))
# 组合成命名向量
model_data_vec <- setNames(value_list, name_list)

model_data_vec

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