R语言那些事 统计系列（一） rstatix and ggpubr 做 t 检验与可视化之添加P值

缘由

统计是生信的基石之一

碎碎念

要爱一个具体的人，不要爱一个抽象的人； 要爱生活，不要爱生活的意义。 ——《卡拉马佐夫兄弟》

Ref

How to Add P-Values onto Basic GGPLOTS - Datanovia

解读

数据加载

R 自带的 ToothGrowth 数据集：

library(ggpubr)
library(rstatix)
 
df <- ToothGrowth
df$dose <- as.factor(df$dose)
head(df, 3)

数据结构：

  len supp dose
1  4.2   VC  0.5
2 11.5   VC  0.5
3  7.3   VC  0.5

• len：牙齿长度
• supp：补充剂类型（VC / OJ）
• dose：剂量

---

2. 两组独立样本 t 检验

stat.test <- df %>%
  t_test(len ~ supp) %>%
  add_significance()
stat.test

绘制箱线图并标注 p 值：

bxp <- ggboxplot(df, x = "supp", y = "len", fill = "#00AFBB")
stat.test <- stat.test %>% add_xy_position(x = "supp")
 
bxp +
  stat_pvalue_manual(stat.test, label = "p") +
  scale_y_continuous(expand = expansion(mult = c(0.05, 0.1)))

可以灵活调整 p 值标签和括号位置：

bxp + stat_pvalue_manual(
  stat.test, label = "T-test, p = {p}",
  vjust = -1,
  bracket.nudge.y = 1
) +
  scale_y_continuous(expand = expansion(mult = c(0.05, 0.15)))

---

3. 分组比较 + 多重检验校正

按照不同 dose 分组，比较 supp：

stat.test <- df %>%
  group_by(dose) %>%
  t_test(len ~ supp) %>%
  adjust_pvalue() %>%
  add_significance("p.adj")

分面箱线图：

stat.test <- stat.test %>% add_xy_position(x = "supp")
 
bxp <- ggboxplot(df, x = "supp", y = "len", fill = "#00AFBB",
                 facet.by = "dose")
 
bxp + stat_pvalue_manual(stat.test, label = "p.adj")

还可以自定义 p 值格式：

stat.test$p.format <- p_format(stat.test$p.adj, accuracy = 0.01, leading.zero = FALSE)
 
bxp + stat_pvalue_manual(stat.test, label = "p.format")

或者科学计数法：

stat.test$p.scient <- format(stat.test$p.adj, scientific = TRUE)
 
bxp + stat_pvalue_manual(stat.test, label = "p.scient")

---

4. 配对样本 t 检验

如果数据是配对设计（如同一组受试者在两种处理条件下的结果），可以设置 paired = TRUE：

stat.test <- df %>%
  t_test(len ~ supp, paired = TRUE) %>%
  add_significance()

绘制配对箱线图：

bxp <- ggpaired(df, x = "supp", y = "len", fill = "#E7B800",
                line.color = "gray", line.size = 0.4)
stat.test <- stat.test %>% add_xy_position(x = "supp")
 
bxp + stat_pvalue_manual(stat.test, label = "p.signif")

显示 p 值 + 显著性星号：

bxp + stat_pvalue_manual(stat.test, label = "{p}{p.signif}")

---

5. 多组比较

对不同剂量下的 len 做多组比较：

bxp <- ggboxplot(df, x = "dose", y = "len", fill = "dose",
                 palette = c("#00AFBB", "#E7B800", "#FC4E07"))
 
stat.test <- df %>% t_test(len ~ dose)
stat.test <- stat.test %>% add_xy_position(x = "dose")
 
bxp + stat_pvalue_manual(stat.test, label = "p.adj.signif")

也可以配合柱状图（均值 ± SD）：

bp <- ggbarplot(df, x = "dose", y = "len", add = "mean_sd", fill = "dose",
                palette = c("#00AFBB", "#E7B800", "#FC4E07"))
 
stat.test <- stat.test %>% add_xy_position(fun = "mean_sd", x = "dose")
 
bp + stat_pvalue_manual(stat.test, label = "p.adj.signif")