在數(shù)據(jù)分析過(guò)程中劈猪，常常需要把組間的顯著性添加到圖形中丙唧，但是在ggplot2中實(shí)現(xiàn)起來(lái)略顯麻煩办桨，幸運(yùn)的是搭综，有很多R包可以幫助我們實(shí)現(xiàn)這一操作垢箕，比如ggsignif和ggpubr。

1. 統(tǒng)計(jì)方法

一般根據(jù)數(shù)據(jù)是否符合正態(tài)分布兑巾，選擇合適的統(tǒng)計(jì)方法：

2. ggsignif

2.1 ggsignif介紹

ggsignif包主要函數(shù)為：geom_signif()和stat_signif()条获，常用geom_signif()。

# geom_signif參數(shù)
geom_signif(mapping = NULL, data = NULL, stat = "signif",
position = "identity", na.rm = FALSE, show.legend = NA,
inherit.aes = TRUE, comparisons = NULL, test = "wilcox.test",
test.args = NULL, annotations = NULL, map_signif_level = FALSE,
y_position = NULL, xmin = NULL, xmax = NULL, margin_top = 0.05,
step_increase = 0, tip_length = 0.03, size = 0.5, textsize = 3.88,
family = "", vjust = 0, ...)

2.2 ggsignif使用：

安裝

# 安裝穩(wěn)定版本：
install.packages("ggsignif")
#安裝最新開發(fā)版本：
devtools::install_github("const-ae/ggsignif")

使用iris數(shù)據(jù)集做演示

library(ggplot2)
library(ggthemr)     #載入主題配置包
library(ggsignif)    #載入ggsignif
library(patchwork) #載入拼圖包
head(iris)
#   Sepal.Length Sepal.Width Petal.Length Petal.Width Species
# 1          5.1         3.5          1.4         0.2  setosa
# 2          4.9         3.0          1.4         0.2  setosa
# 3          4.7         3.2          1.3         0.2  setosa
# 4          4.6         3.1          1.5         0.2  setosa
# 5          5.0         3.6          1.4         0.2  setosa
# 6          5.4         3.9          1.7         0.4  setosa

Species的三組兩兩分別作差異性檢驗(yàn)府框，提前設(shè)定好配對(duì)分析的list

compaired <- list(c("versicolor", "virginica"), 
                  c("versicolor","setosa"), 
                  c("virginica","setosa"))

繪制geom_boxplot()和小提琴圖geom_violin()

ggthemr("flat")
p1 <- ggplot(iris, aes(Species, Sepal.Width, fill = Species)) +
    geom_boxplot() +
    ylim(1.5, 6.5) +
    geom_signif(comparisons = compaired,
                step_increase = 0.3,
                map_signif_level = F,
                test = wilcox.test)
p2 <- ggplot(iris, aes(Species, Sepal.Width, fill = Species)) +
    geom_violin() +
    ylim(1.5, 6.5) +
    geom_signif(comparisons = compaired,
                step_increase = 0.3,
                map_signif_level = T, #修改參數(shù)map_signif_level=TRUE
                test = wilcox.test)
p1|p2

3. ggpubr

3.1 介紹

ggpubr添加p-value主要使用ggpubr包中的兩個(gè)函數(shù)：compare_means()和stat_compare_mean()吱窝。

• compare_means():可以進(jìn)行一組或多組間的比較。

compare_means(formula, data, method = "wilcox.test", paired = FALSE,
group.by = NULL, ref.group = NULL, ...)

• stat_compare_mean():自動(dòng)添加p-value盹牧、顯著性標(biāo)記到ggplot圖中

stat_compare_means(mapping = NULL, comparisons = NULL hide.ns = FALSE,
label = NULL, label.x = NULL, label.y = NULL, ...)

3.2 比較獨(dú)立的兩組

#先加載包
library(ggpubr)
library(patchwork)
data("ToothGrowth") #加載演示數(shù)據(jù)集ToothGrowth
head(ToothGrowth)
#    len supp dose
# 1  4.2   VC  0.5
# 2 11.5   VC  0.5
# 3  7.3   VC  0.5
# 4  5.8   VC  0.5
# 5  6.4   VC  0.5
# 6 10.0   VC  0.5

compare_means(len~supp, data=ToothGrowth)
## 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 len   OJ     VC     0.0645 0.064 0.064    ns       Wilcox…

y:測(cè)試中使用的y變量
p:p-value
p.adj:調(diào)整后的p-value。默認(rèn)為p.adjust.method=“holm”
p.format:四舍五入后的p-value
p.signif:顯著性水平
method:用于統(tǒng)計(jì)檢驗(yàn)的方法

繪制箱線圖

p1 <- ggboxplot(ToothGrowth, x="supp", y="len", color = "supp",
               palette = "lancet", add = "jitter")#添加p-valuep+stat_compare_means()
#使用其他統(tǒng)計(jì)檢驗(yàn)方法
p2 <- p1+stat_compare_means(method = "t.test")
p1|p2

上述顯著性標(biāo)記可以通過(guò)label.x踩寇、label.y啄清、hjust及vjust來(lái)調(diào)整
顯著性標(biāo)記可以通過(guò)aes()映射來(lái)更改：
aes(label=…p.format…)或aes(lebel=paste0(“p=”,…p.format…))：只顯示p-value，不顯示統(tǒng)計(jì)檢驗(yàn)方法
aes(label=…p.signif…)：僅顯示顯著性水平
aes(label=paste0(…method…,"\n", “p=”,…p.format…))：p-value與顯著性水平分行顯示
也可以將標(biāo)簽指定為字符向量俺孙，不要映射辣卒，只需將p.signif兩端的…去掉即可

3.3 配對(duì)樣本

compare_means(len~supp, data=ToothGrowth, paired = TRUE)
## A tibble: 1 x 8
#  .y.   group1 group2       p  p.adj p.format p.signif
#  <chr> <chr>  <chr>    <dbl>  <dbl> <chr>    <chr>   
#1 len   OJ     VC     0.00431 0.0043 0.0043   **      
## … with 1 more variable: method <chr>

利用ggpaired()進(jìn)行可視化

ggpaired(ToothGrowth, x="supp", y="len", color = "supp", line.color = "gray",
line.size = 0.4, palette = "jco")+ stat_compare_means(paired = TRUE)

3.4 多組比較 Global test

compare_means(len~dose, data=ToothGrowth, method = "anova")
## A tibble: 1 x 6
#  .y.          p   p.adj p.format p.signif method
#  <chr>    <dbl>   <dbl> <chr>    <chr>    <chr> 
#1 len   9.53e-16 9.5e-16 9.5e-16  ****     Anova 

繪圖

p1=ggboxplot(ToothGrowth, x="dose", y="len", color = "dose", palette = "jco")+
stat_compare_means()
#使用其他的方法
p2=ggboxplot(ToothGrowth, x="dose", y="len", color = "dose", palette = "jco")+
stat_compare_means(method = "anova")
p1|p2

Pairwise comparisons:如果分組變量中包含兩個(gè)以上的水平，那么會(huì)自動(dòng)進(jìn)行pairwise test,默認(rèn)方法為”wilcox.test”

compare_means(len~dose, data=ToothGrowth)
## A tibble: 3 x 8
#  .y.   group1 group2           p    p.adj p.format p.signif
#  <chr> <chr>  <chr>        <dbl>    <dbl> <chr>    <chr>   
#1 len   0.5    1          7.02e-6   1.4e-5 7.0e-06  ****    
#2 len   0.5    2          8.41e-8   2.5e-7 8.4e-08  ****    
#3 len   1      2          1.77e-4   1.8e-4 0.00018  ***     
## … with 1 more variable: method <chr>

#可以指定比較哪些組
my_comparisons <- list(c("0.5", "1"), c("1", "2"), c("0.5", "2"))
p1=ggboxplot(ToothGrowth, x="dose", y="len", color = "dose",palette = "jco")+
  stat_compare_means(comparisons=my_comparisons)+ # Add pairwisecomparisons p-value 
  stat_compare_means(label.y = 50) # Add global p-value
#可以通過(guò)修改參數(shù)label.y來(lái)更改標(biāo)簽的位置
p2=ggboxplot(ToothGrowth, x="dose", y="len", color = "dose",palette = "jco")+
  stat_compare_means(comparisons=my_comparisons, label.y = c(29, 35, 40))+
  # Add pairwise comparisons p-value
  stat_compare_means(label.y = 45) # Add global p-value
p1|p2

參考：https://blog.csdn.net/xj4math/article/details/115448669