简单柱形图

Choosing the Right Statistical Test

Statistical Hypothesis Testing Overview

Nonparametric Tests vs. Parametric Tests

Multiple comparisons

Parametric tests of means	Nonparametric tests of medians
1-sample t-test	1-sample Sign, 1-sample Wilcoxon
2-sample t-test	Mann-Whitney test
One-Way ANOVA	Kruskal-Wallis, Mood’s median test
Factorial DOE with a factor and a blocking variable	Friedman test

Parametric analyses	Sample size requirements for nonnormal data
1-sample t-test	Greater than 20
2-sample t-test	Each group should have more than 15 observations
One-Way ANOVA	For 2-9 groups, each group should have more than 15 observations For 10-12 groups, each group should have more than 20 observations

“Once you have determined that differences exist among the means, post hoc range tests and pairwise multiple comparisons can determine which means differ. Range tests identify homogeneous subsets of means that are not different from each other. Pairwise multiple comparisons test the difference between each pair of means, and yield a matrix where asterisks indicate significantly different group means at an alpha level of 0.05” (SPSS, Inc.).

Kruskal–Wallis Test
Kruskal-Wallis Test in R
Primer on multiple testing
下述数据样品共有5组，每组有4个平行，所以选用Nonparametric的Kruskal-Wallis test计算差异，这里选择了pairwise multiple comparisons去看哪两组间存在差异，具体使用的是Mann-Whitney test，并使用Benjamini–Hochberg procedure去校正FDR（一般使用Bonferroni校正FWER）

library(dplyr)
library(tidyr)
library(ggplot2)
library(ggsci)
library(multcompView)

calc_data <- data.frame(x=c("WT","CK","Treatment A","Treatment B","Treatment C"),val1=c(21,20,11,10,13),val2=c(21.2,20.3,11.2,10.2,13.1),val3=c(20.8,20.1,10.9,10.3,14),val4=c(22,19,12,9,12))

long_data <- gather(calc_data,replicates, measurement, val1:val4, factor_key=TRUE)

tmp <- group_by(long_data, x) %>%
  summarise(
    count = n(),
    mean = mean(measurement, na.rm = TRUE),
    sd = sd(measurement, na.rm = TRUE),
    median = median(measurement, na.rm = TRUE),
    IQR = IQR(measurement,na.rm = TRUE),
    max_value = max(measurement)
  )
  
kruskal_test_res <- kruskal.test(measurement ~ x, data = long_data)
# pairwise.wilcox.test(measurement, x, data = long_data, p.adjust.method = 'BH')
# dunn.test::dunn.test(signal, group, method = 'bh')  

calc_data <- as.data.frame(t(calc_data))
colnames(calc_data) <- calc_data[1,]
calc_data <- calc_data[-1,]
wilcox_test_res <- combn(colnames(calc_data), 2, function(x) list(data = c(paste(x[1], x[2])), p = wilcox.test(x = as.numeric(calc_data[[x[1]]]), y =  as.numeric(calc_data[[x[2]]]))$p.value), simplify = FALSE)
wilcox_test_res <- as.data.frame(do.call(rbind,wilcox_test_res))
wilcox_test_res$p <- p.adjust(unlist(wilcox_test_res$p),method="BH")
wilcox_test_res$data <- gsub("Treatment ","Treatment_",wilcox_test_res$data)
wilcox_test_res$data <- gsub(" ","-",wilcox_test_res$data)
signif_letter_input <- vector()
signif_letter_input <- unlist(wilcox_test_res$p)
names(signif_letter_input) <- unlist(wilcox_test_res$data)
signif_letters <- multcompLetters(signif_letter_input)
anno <- data.frame(x=names(signif_letters$Letters),letter=signif_letters$Letters)
anno$x <- gsub("_"," ",anno$x)
plot_data <- merge(anno,tmp)
ggplot(plot_data,aes(x=x,y=mean,fill=x))+geom_bar(stat="identity",color="black",width=.75)+geom_errorbar(aes(ymin=mean, ymax=mean+sd), width=.2)+geom_text(aes(x=x,label = letter,y=max_value+1) , hjust=0.5,inherit.aes = FALSE,size=4) +xlab(NULL)+ylab("The content of drug (mg/g)")+scale_y_continuous(limits = c(0,25), expand = c(0, 0))+theme_classic()+theme(axis.text.x=element_text(angle=30,hjust=1,size=10),axis.text.y=element_text(size=10),axis.title.y=element_text(size=12,margin = margin(t = 0, r = 5, b = 0, l = 0)),text=element_text(family="sans"),legend.position="none",plot.title = element_text(size=16,colour = "black",hjust = 0.5))+ggtitle("The Content of Drug")+scale_fill_npg() + annotate("text",  x=Inf, y = Inf, label = paste0("Kruskal-Wallis, p = ",formatC(kruskal_test_res$p.value, format = "e", digits = 2)), vjust=1, hjust=1,size = 4)