按TIME_SLOT分组执行多组Wilcoxon检验遇报错，求解决方案

错误原因分析 #

• 原代码的group_by(LANE_TYPE, BARRIER_POS, TIME_SLOT)会将数据拆分为单一BARRIER_POS取值的小组，无法满足R_CO_GC ~ BARRIER_POS的两组对比需求，这是触发'x' must be numeric报错的核心原因。
• 在dplyr管道中直接调用wilcoxon.test需要结合列表存储结果，否则会传递错误的参数类型。

修正后的代码实现 #

1. 数据预处理 #

先过滤出符合条件的数据集：

library(dplyr)
library(purrr) # 用于处理列表结果

filtered_data <- lr2_analysis_CO_FIL %>%
  filter(BARRIER_POS != 'no_street' & HOLIDAY == 'N' & HOUR >= 6)

2. 按时间段执行三组Wilcoxon检验 #

针对每个TIME_SLOT，分别执行指定的三组对比，并提取关键结果：

stat.test <- filtered_data %>%
  group_by(TIME_SLOT) %>%
  summarize(
    # 组1：footway receptor vs footway no_barrier
    test1 = list(wilcox.test(
      x = R_CO_GC[LANE_TYPE == "footway" & BARRIER_POS == "receptor"],
      y = R_CO_GC[LANE_TYPE == "footway" & BARRIER_POS == "no_barrier"],
      paired = FALSE # 根据实验设计选择是否配对，默认非配对
    )),
    # 组2：footway receptor vs roadway emisor
    test2 = list(wilcox.test(
      x = R_CO_GC[LANE_TYPE == "footway" & BARRIER_POS == "receptor"],
      y = R_CO_GC[LANE_TYPE == "roadway" & BARRIER_POS == "emisor"],
      paired = FALSE
    )),
    # 组3：roadway emisor vs roadway no_barrier
    test3 = list(wilcox.test(
      x = R_CO_GC[LANE_TYPE == "roadway" & BARRIER_POS == "emisor"],
      y = R_CO_GC[LANE_TYPE == "roadway" & BARRIER_POS == "no_barrier"],
      paired = FALSE
    ))
  ) %>%
  # 提取检验结果中的p值、统计量，并生成可视化标注文本
  mutate(
    test1_p = map_dbl(test1, ~ .x$p.value),
    test1_stat = map_dbl(test1, ~ .x$statistic),
    test2_p = map_dbl(test2, ~ .x$p.value),
    test2_stat = map_dbl(test2, ~ .x$statistic),
    test3_p = map_dbl(test3, ~ .x$p.value),
    test3_stat = map_dbl(test3, ~ .x$statistic),
    # 格式化标注文本，保留3位小数
    test1_label = paste0("footway receptor vs no_barrier\np = ", round(test1_p, 3)),
    test2_label = paste0("footway receptor vs roadway emisor\np = ", round(test2_p, 3)),
    test3_label = paste0("roadway emisor vs no_barrier\np = ", round(test3_p, 3))
  )

3. 箱线图添加检验结果标注 #

使用ggplot2和ggsignif包将检验结果标注到分面箱线图中：

library(ggplot2)
library(ggsignif)

# 绘制基础分面箱线图
base_plot <- ggplot(filtered_data, aes(x = interaction(LANE_TYPE, BARRIER_POS), y = R_CO_GC)) +
  geom_boxplot(fill = "#69b3a2") +
  facet_wrap(~TIME_SLOT, ncol = 2) + # 按时间段分面，调整列数
  labs(x = "组别", y = "R_CO_GC数值", title = "不同时间段下的R_CO_GC分布与Wilcoxon检验结果") +
  theme_bw()

# 添加显著性标注
final_plot <- base_plot +
  geom_signif(
    data = stat.test,
    aes(
      xmin = c(1, 1, 2), # 每组对比的左侧x位置
      xmax = c(2, 3, 3), # 每组对比的右侧x位置
      annotations = c(test1_label, test2_label, test3_label),
      facet_group = TIME_SLOT
    ),
    y_position = max(filtered_data$R_CO_GC) * 1.1, # 标注位置在箱线图上方10%处
    tip_length = 0.01 # 调整标注线的尖端长度
  )

print(final_plot)

关键说明 #

• 若你的实验是配对设计，需将wilcox.test中的paired参数改为TRUE。
• 可根据数据分布调整y_position参数，确保标注不遮挡箱线图。
• 若不需要保留原始检验对象，可在mutate后用select(-test1, -test2, -test3)删除冗余列。

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