R语言年度截面回归结果写入数据框的问题及报错排查

问题1：提取回归系数与p值并整合到结果数据框 #

可以用broom包的tidy()函数快速提取回归模型的系数、标准误、p值，再结合glance()提取的R²、F统计量等指标，按年份合并成结构化数据框。如果不想用第三方包，也可以直接从summary(model)的输出中手动提取系数和p值。

问题2：解决BP/BG检验的报错 #

报错Error in y - z$residuals : non-numeric argument to binary operator通常是两个原因：

• 因变量或自变量不是数值型，导致模型残差类型异常
• 调用检验函数时参数传递错误，手动计算残差时类型不匹配

解决核心：确保拟合模型前所有变量都是数值型，调用检验函数时直接传入模型对象，不要手动处理残差。

完整可运行代码示例 #

# 加载必要包
library(readxl)
library(broom)
library(dplyr)
library(lmtest)

# 初始化结果存储列表
result_list <- list()

# 循环处理1995-2018年数据
for(year in 1995:2018){
  # 读取对应年份的工作表（假设工作表名为年份字符串，如"1995"）
  df <- read_excel("你的数据文件路径.xlsx", sheet = as.character(year))
  
  # 强制转换核心变量为数值型（解决检验报错的关键步骤）
  df <- df %>% mutate(across(c(y, x1, x2), as.numeric)) # 替换成你的因变量和自变量名
  
  # 拟合年度截面回归
  model <- lm(y ~ x1 + x2, data = df) # 替换为你的实际回归公式
  
  # 提取系数与p值
  coeff_data <- tidy(model) %>% 
    mutate(year = year) %>%
    select(year, term, estimate, p.value) %>%
    rename(coef = estimate, coef_p = p.value)
  
  # 提取模型整体指标（R²、F统计量等）
  model_metrics <- glance(model) %>%
    mutate(year = year) %>%
    select(year, r.squared, adj.r.squared, statistic, p.value) %>%
    rename(f_stat = statistic, f_p = p.value)
  
  # 执行BP和BG检验
  bp_result <- tidy(bptest(model)) %>%
    mutate(year = year, test_type = "BP") %>%
    select(year, test_type, statistic, p.value) %>%
    rename(test_stat = statistic, test_p = p.value)
  
  bg_result <- tidy(bgtest(model)) %>%
    mutate(year = year, test_type = "BG") %>%
    select(year, test_type, statistic, p.value) %>%
    rename(test_stat = statistic, test_p = p.value)
  
  # 合并当前年份的所有结果
  year_all_data <- coeff_data %>%
    left_join(model_metrics, by = "year") %>%
    bind_rows(bp_result, bg_result)
  
  # 存入结果列表
  result_list[[as.character(year)]] <- year_all_data
}

# 合并所有年份数据为最终结构化数据框
final_output <- bind_rows(result_list)

# 查看结果
print(final_output)

输出说明 #

最终的final_output是结构化数据框，包含字段：

• year：年份
• term/test_type：变量名或检验类型（BP/BG）
• coef/test_stat：回归系数或检验统计量
• coef_p/test_p：系数p值或检验p值
• r.squared、adj.r.squared：决定系数
• f_stat、f_p：F统计量及对应p值

注意事项 #

• 替换代码中的y ~ x1 + x2为你的实际回归公式，变量名要与Excel工作表中的列名完全一致
• 如果Excel工作表名不是纯数字年份（如带有前缀），需调整sheet参数的取值逻辑
• 若不想使用broom包，手动提取系数和p值可使用summary(model)$coefficients[, c(1,4)]

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