在R中按分组对数据框各列执行配对样本Wilcoxon检验求助

需求说明 #

将给定的数据框按Group（Group1和Group2）拆分后，针对每个微生物丰度列，以ID为配对依据，对Time（T0和T2）的配对数据执行配对样本Wilcoxon检验。

输入数据 #

df <- structure(list(ID = c(1L, 2L, 3L, 4L, 5L, 10L, 11L, 12L, 13L, 
14L, 1L, 2L, 3L, 4L, 5L, 10L, 11L, 12L, 13L, 14L), Group = c("Group1", 
"Group1", "Group1", "Group1", "Group1", "Group2", "Group2", "Group2", 
"Group2", "Group2", "Group1", "Group1", "Group1", "Group1", "Group1", 
"Group2", "Group2", "Group2", "Group2", "Group2"), Time = c("T0", 
"T0", "T0", "T0", "T0", "T0", "T0", "T0", "T0", "T0", "T2", "T2", 
"T2", "T2", "T2", "T2", "T2", "T2", "T2", "T2"), Alistipes_putredinis = c(20.54997, 
0, 0, 0, 0.00226, 0, 12.83106, 0.38555, 16.45834, 0, 0.6405, 
0, 0, 0, 0, 16.32105, 0, 0, 0, 0), Bidobacterium_l = c(14.43141, 
0.25318, 0.83121, 0.54282, 8.50687, 0.10432, 0.0051, 0.00139, 
0.11766, 0.57905, 0.74302, 0.79018, 0.20329, 5.34884, 1.75612, 
0.25502, 0, 0.60675, 0.59414, 0.85824), Dialister_invisus = c(12.45032, 
7.88459, 0.99032, 1.84241, 0.93828, 1.96771, 5.16734, 3.82884, 
1.55069, 0.97391, 1.43845, 9.77174, 2.59287, 4.70876, 2.30655, 
3.08246, 10.56866, 0.49946, 0.95196, 0.14022), Clostridia_bacterium = c(7.5127, 
0, 0, 0, 0, 6.65269, 0, 0, 4.11219, 0, 5.34908, 0.00794, 0, 0, 
0, 0, 0, 4.33676, 6.2422, 0), Ruminococcus_sp_NSJ_71 = c(6.57903, 
1.45815, 0.28668, 1.66816, 1.66008, 1.85348, 1.80051, 0.91537, 
3.12064, 3.00647, 2.91748, 1.97839, 0.0726, 3.7829, 1.59076, 
1.05453, 0.03881, 3.84824, 0.01241, 2.88977)), class = "data.frame", row.names = c(NA, 
-20L))

解决方案 #

方法1：使用tidyverse包（简洁易读）

# 加载tidyverse包
library(tidyverse)

# 将宽格式数据转为长格式，方便分组处理
long_df <- df %>%
  pivot_longer(
    cols = -c(ID, Group, Time),  # 排除非丰度列
    names_to = "Microbe",        # 微生物名称列
    values_to = "Abundance"      # 丰度值列
  )

# 按Group和Microbe分组，执行配对Wilcoxon检验
wilcox_results <- long_df %>%
  group_by(Group, Microbe) %>%
  summarise(
    wilcox_stat = wilcox.test(Abundance ~ Time, paired = TRUE)$statistic,
    p_value = wilcox.test(Abundance ~ Time, paired = TRUE)$p.value,
    .groups = "drop"
  )

# 输出结果
print(wilcox_results)

方法2：使用Base R（无需额外加载包）

# 按Group拆分数据框
group_data <- split(df, df$Group)

# 定义检验函数：对单个分组的所有丰度列执行配对Wilcoxon检验
run_paired_wilcox <- function(group_df) {
  # 获取需要检验的丰度列名
  target_cols <- setdiff(names(group_df), c("ID", "Group", "Time"))
  
  # 遍历每列执行检验
  result_list <- lapply(target_cols, function(col) {
    # 提取T0和T2的配对数据
    t0_data <- group_df[group_df$Time == "T0", col]
    t2_data <- group_df[group_df$Time == "T2", col]
    
    # 执行配对Wilcoxon检验
    test_result <- wilcox.test(t0_data, t2_data, paired = TRUE)
    
    # 整理结果为数据框
    data.frame(
      Microbe = col,
      wilcox_stat = test_result$statistic,
      p_value = test_result$p.value
    )
  })
  
  # 合并所有列的结果
  do.call(rbind, result_list)
}

# 对两个分组分别执行检验
group1_results <- run_paired_wilcox(group_data$Group1)
group2_results <- run_paired_wilcox(group_data$Group2)

# 查看结果
cat("Group1检验结果：\n")
print(group1_results)
cat("\nGroup2检验结果：\n")
print(group2_results)

结果说明 #

两种方法都会输出每个分组下，各微生物丰度列的Wilcoxon检验统计量和p值，可用于判断T0和T2时间点的丰度是否存在显著差异。

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