如何在gtsummary中为事后检验所有两两比较生成对应p值列

结论 #

gtsummary完全支持输出4组分组的整体检验+所有两两比较独立p值列的需求，可通过内置的add_p()搭配自定义add_stat()函数实现，也可对接emmeans等事后检验工具输出校正后的p值。

实现步骤 #

• 第一步：用tbl_summary()生成基线汇总表，通过add_p()添加整体检验p值（连续变量可选anova，分类变量可选卡方检验）
• 第二步：自定义两两比较p值计算函数，通过add_stat()循环追加所有两两比较的p值列
• 第三步：自定义列头和p值格式化规则，输出最终表格

代码示例 #

首先加载依赖包：

library(gtsummary)
library(dplyr)

自定义两两比较p值计算函数：

# 入参pair指定要比较的两组，返回对应p值列
calc_pair_p <- function(data, variable, by, pair, ...) {
  # 筛选目标两组数据
  comp_data <- filter(data, .data[[by]] %in% pair)
  comp_data[[by]] <- droplevels(comp_data[[by]])
  # 自动匹配检验方法：连续变量用t检验，分类变量用卡方检验
  test_method <- ifelse(is.numeric(comp_data[[variable]]), "t.test", "chisq.test")
  # 计算p值
  p_val <- calculate_p(
    data = comp_data, variable = variable, by = by, test = test_method
  )$p.value
  # 返回格式化的单列数据
  return(tibble(!!paste0("p_", pair[1], "_vs_", pair[2]) := p_val))
}

生成最终表格：

# 示例用trial数据集，生成分组为a/b/c/d的4组样本
trial_4group <- trial %>%
  mutate(group = sample(c("a", "b", "c", "d"), n(), replace = TRUE)) %>%
  select(group, age, grade, response)

# 生成汇总表+所有p值列
tbl_result <- trial_4group %>%
  tbl_summary(by = group) %>%
  # 添加整体检验p值
  add_p(test = list(
    all_continuous() ~ "anova",
    all_categorical() ~ "chisq.test"
  )) %>%
  # 逐一添加所有两两比较p值列
  add_stat(fns = everything() ~ ~calc_pair_p(.x, pair = c("a","b"))) %>%
  add_stat(fns = everything() ~ ~calc_pair_p(.x, pair = c("a","c"))) %>%
  add_stat(fns = everything() ~ ~calc_pair_p(.x, pair = c("a","d"))) %>%
  add_stat(fns = everything() ~ ~calc_pair_p(.x, pair = c("b","c"))) %>%
  add_stat(fns = everything() ~ ~calc_pair_p(.x, pair = c("b","d"))) %>%
  add_stat(fns = everything() ~ ~calc_pair_p(.x, pair = c("c","d"))) %>%
  # 自定义列头
  modify_header(
    p.value = "**整体检验p值**",
    p_a_vs_b = "**a vs b p值**",
    p_a_vs_c = "**a vs c p值**",
    p_a_vs_d = "**a vs d p值**",
    p_b_vs_c = "**b vs c p值**",
    p_b_vs_d = "**b vs d p值**",
    p_c_vs_d = "**c vs d p值**"
  ) %>%
  # 统一格式化p值
  modify_fmt_fun(c(p.value, starts_with("p_")), ~ style_pvalue(.x))

补充说明 #

如果需要输出经多重检验校正的p值，只需修改calc_pair_p中的p值计算逻辑，对接emmeans、multcomp等工具的输出即可，列追加的逻辑无需改动。

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