R语言使用lm模型拟合多组数据、预测新值及提取模型参数问询

先修正原代码的分组逻辑 #

你原代码中group_by(group, numb)会导致每个嵌套子集仅1行数据，无法正常拟合线性模型，首先需要调整为仅按group分组。 #

完整实现代码 #

1. 嵌套建模、提取模型参数、生成历史+未来预测值 #

你需要的两个需求都可以在嵌套数据结构内完成，不需要拆分表格：

library(tidyverse)
library(broom)

# 导入你的原始数据
test_data <- structure(list(group = c("Group_1", "Group_1", "Group_1", "Group_1", 
"Group_2", "Group_2", "Group_2", "Group_2", "Group_3", "Group_3", 
"Group_3", "Group_3"), numb = c(1, 2, 3, 4, 1, 2, 3, 4, 1, 2, 
3, 4), total = c(616597744.48, 516080403.54476, 990894258.72, 
923317167.70895, 3850620416.96513, 3823237639.55, 4150030206.48, 
4317861944.93, 6403590027.27012, 6078175252.18719, 6951291610.00877, 
6432993298.93888)), class = c("grouped_df", "tbl_df", "tbl", 
"data.frame"), row.names = c(NA, -12L), groups = structure(list(
    group = c("Group_1", "Group_2", "Group_3"), .rows = structure(list(
        1:4, 5:8, 9:12), ptype = integer(0), class = c("vctrs_list_of", 
    "vctrs_vctr", "list"))), class = c("tbl_df", "tbl", "data.frame"
), row.names = c(NA, -3L), .drop = TRUE)

# 嵌套处理：同一张表存储原始数据、模型、系数、模型指标、预测值
result <- test_data %>%
  # 仅按分组变量分组
  group_by(group) %>%
  nest() %>%
  # 拟合每个组的lm模型
  mutate(model = map(data, ~ lm(total ~ numb, data = .x))) %>%
  # 提取模型系数（截距、斜率、系数p值等）
  mutate(coef = map(model, ~ tidy(.x))) %>%
  # 提取模型整体指标（R平方、调整R平方、模型p值、F统计量等）
  mutate(model_metrics = map(model, ~ glance(.x))) %>%
  # 生成历史数据的拟合值
  mutate(fit_hist = map(model, ~ augment(.x))) %>%
  # 生成后续3个新值（numb=5、6、7）的预测值
  mutate(fit_new = map(model, ~ predict(.x, newdata = tibble(numb = 5:7)) %>%
                         as_tibble() %>%
                         rename(total = value) %>%
                         mutate(numb = 5:7, .before = 1))) %>%
  ungroup()

2. 提取你需要的模型参数 #

• 查看所有组的回归系数：

result %>% unnest(coef)

• 查看所有组的R平方等模型整体指标：

result %>% unnest(model_metrics)

3. 绘制包含原始值、拟合趋势、未来预测的分面图 #

# 合并绘图所需的所有数据
plot_data <- bind_rows(
  # 原始观测值
  test_data %>% mutate(type = "原始观测值"),
  # 历史拟合值
  result %>% unnest(fit_hist) %>% select(group, numb, total = .fitted) %>% mutate(type = "历史拟合值"),
  # 未来预测值
  result %>% unnest(fit_new) %>% mutate(type = "未来预测值")
)

# 绘图
ggplot(plot_data, aes(x = numb, y = total, color = type)) +
  geom_point(size = 2) +
  # 为原始观测值添加拟合趋势线
  stat_smooth(data = plot_data %>% filter(type == "原始观测值"),
              method = "lm", se = FALSE, color = "grey50", linetype = "dashed") +
  facet_wrap(~group, scales = "free_y") +
  labs(x = "numb", y = "total", color = "数据类型") +
  theme_bw()

疑问明确答复 #

• 可以直接在嵌套的data列体系下完成所有预测操作，不需要拆分嵌套结构
• 可以通过broom包的tidy()和glance()函数直接在同一张嵌套表中提取所有组的系数、R平方等参数，不需要单独逐个提取模型结果

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