如何使用rms::ols函数生成多模型对象的估计值与置信区间

由于broom::tidy暂不支持rms::ols对象，我们可以通过直接提取ols模型的系数、置信区间信息，或者自定义适配函数来实现类似tidy(lm)的批量输出效果。

方法一：直接在分组流程中提取信息 #

library(tidyverse)
library(rms)

mtcars %>%
  # 将预测变量转换为长格式，和lm示例逻辑一致
  gather(predictor, measure, -mpg) %>%
  group_by(predictor) %>%
  do({
    # 对每组数据拟合ols模型
    ols_model <- ols(mpg ~ measure, data = .)
    # 提取置信区间并整理列名
    ci_table <- as_tibble(confint(ols_model), rownames = "term") %>%
      rename(conf.low = `2.5%`, conf.high = `97.5%`)
    # 提取系数、标准误、p值信息
    coef_table <- as_tibble(coef(summary(ols_model)), rownames = "term") %>%
      rename(estimate = Estimate, std.error = S.E., p.value = Pr(>|t|))
    # 合并信息并返回
    left_join(coef_table, ci_table, by = "term")
  })

方法二：自定义适配ols的tidy函数 #

如果需要更简洁的调用方式，可以模仿broom::tidy的格式自定义一个处理ols对象的函数：

library(tidyverse)
library(rms)

# 自定义tidy.ols函数，支持提取置信区间
tidy.ols <- function(model, conf.int = TRUE, conf.level = 0.95) {
  # 提取系数、标准误、p值
  base_table <- as_tibble(coef(summary(model)), rownames = "term") %>%
    rename(estimate = Estimate, std.error = S.E., p.value = Pr(>|t|))
  
  # 如果需要置信区间则合并
  if (conf.int) {
    ci_data <- as_tibble(confint(model, level = conf.level), rownames = "term") %>%
      rename(conf.low = names(.)[2], conf.high = names(.)[3])
    base_table <- left_join(base_table, ci_data, by = "term")
  }
  
  return(base_table)
}

# 批量拟合并输出结果
mtcars %>%
  gather(predictor, measure, -mpg) %>%
  group_by(predictor) %>%
  do(tidy.ols(ols(mpg ~ measure, data = .), conf.int = TRUE))

两种方法最终都会输出和tidy(lm)结构一致的结果，包含分组的预测变量、系数估计值、标准误、置信区间和p值。

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