compare_means()与stat_compare_means()同数据结果差异问题咨询

问题描述 #

使用compare_means()计算Treatment组间Protein表达的差异显著性时，得到Wilcoxon检验p值为0.00794，对应显著性标记为**，执行代码与返回结果如下：

compare_means(Protein ~ Treatment, data = patient_data[])

# A tibble: 1 x 8
  .y.       group1 group2       p  p.adj p.format p.signif method  
  <chr>     <chr>  <chr>    <dbl>  <dbl> <chr>    <chr>    <chr>   
1 Protein Treated   Ctrl   0.00794 0.0079 0.0079   **       Wilcoxon

但使用如下ggplot2代码调用stat_compare_means()添加统计标注时，得到的检验结果与上述计算结果差异极大：

ggplot(patient_data[], aes(x = Treatment, y = Protein, color = Experiment)) +
  geom_jitter(width = 0.1) +
  stat_compare_means(method = "wilcox.test")

偏差原因 #

核心问题是ggplot2全局美学映射的继承机制导致检验分组错误：
你在ggplot()全局参数中设置了color = Experiment，这个分组映射会被后续所有几何对象、统计计算函数默认继承。stat_compare_means()会自动按照映射的分组变量（此处为3个水平的Experiment）拆分数据集，在每个Experiment子集内分别做Treatment组间的Wilcoxon检验，而非预期的、跨Experiment合并所有Ctrl和Treated样本做整体组间比较。
单独运行compare_means()时没有引入Experiment分组，是对全量样本按Treatment分组计算，因此和绘图时拆分3个子集算出的检验结果存在极大差异。从提供的数据集可以看到，每个Experiment子集内仅包含2-4个样本，部分子集单组样本量仅1个，无法输出稳定可靠的检验结果。

修复方案 #

选择以下任意一种方式调整代码即可得到和单独计算一致的检验结果：

• 方案1：将color = Experiment的美学映射放到geom_jitter()层内部，不放在全局ggplot映射中，避免统计函数继承该分组规则：

ggplot(patient_data[], aes(x = Treatment, y = Protein)) +
  geom_jitter(aes(color = Experiment), width = 0.1) +
  stat_compare_means(method = "wilcox.test")

• 方案2：在stat_compare_means()层显式指定分组规则，覆盖全局的分组映射，强制按x轴的Treatment类别做全样本比较：

ggplot(patient_data[], aes(x = Treatment, y = Protein, color = Experiment)) +
  geom_jitter(width = 0.1) +
  stat_compare_means(method = "wilcox.test", aes(group = 1))

说明：aes(group = 1)的作用是告知ggplot在统计计算环节忽略其他分组变量，将x轴每个类别下的所有观测值视为同一组进行检验。

附：使用的patient_data数据集子集 #

Sample.Name     Protein Experiment  Treatment
Sample_name_01  1070    01          Ctrl
Sample_name_02  222     02          Ctrl
Sample_name_03  243     02          Ctrl
Sample_name_04  264     03          Ctrl
Sample_name_05  285     03          Ctrl
Sample_name_06  187     01          Treated
Sample_name_07  55      02          Treated
Sample_name_08  72      02          Treated
Sample_name_09  52      03          Treated
Sample_name_10  55      03          Treated

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