R语言非参数检验多重比较
医学和生信笔记,专注R语言在临床医学中的使用,R语言数据分析和可视化。
之前介绍了多个样本均数的多重比较,今天说说kruskal-Wallis H检验后的多重比较,Friedman M检验后的多重比较。
也是和课本对照着来,孙振球,徐勇勇《医学统计学》第四版。
非参数检验后的多重比较,我们也是用这个宝藏R包:PMCMRplus。
kruskal-Wallis H检验及多重比较
使用课本例8-5的数据。
rm(list = ls()) death_rate <- c(32.5,35.5,40.5,46,49,16,20.5,22.5,29,36,6.5,9.0,12.5,18,24) drug <- rep(c("Drug_A","drug_B","drug_C"),each=5) mydata <- data.frame(death_rate,drug) # 分类变量因子化 mydata$drug <- factor(mydata$drug) str(mydata) ## data.frame: 15 obs. of 2 variables: ## $ death_rate: num 32.5 35.5 40.5 46 49 16 20.5 22.5 29 36 ... ## $ drug : Factor w/ 3 levels "Drug_A","drug_B",..: 1 1 1 1 1 2 2 2 2 2 ...
进行kruskal-Wallis H 检验:
fit <- kruskal.test(death_rate ~ drug, data = mydata) fit ## ## Kruskal-Wallis rank sum test ## ## data: death_rate by drug ## Kruskal-Wallis chi-squared = 9.74, df = 2, p-value = 0.007673
多重检验,课本上用的是Nemenyi检验,我们通过多重比较的全能R包PMCMRplus实现。
library(PMCMRplus)
也是提供两种输入方式,直接提供kruskal-Wallis H检验的结果,或者formula形式,都可以。
res <- kwAllPairsNemenyiTest(fit) res <- kwAllPairsNemenyiTest(death_rate ~ drug, data = mydata) summary(res) ## ## Pairwise comparisons using Tukey-Kramer-Nemenyi all-pairs test with Tukey-Dist approximation ## data: death_rate by drug ## alternative hypothesis: two.sided ## P value adjustment method: single-step ## H0 ## q value Pr(>|q|) ## drug_B - Drug_A == 0 2.5 0.1805089 ## drug_C - Drug_A == 0 4.4 0.0052932 ** ## drug_C - drug_B == 0 1.9 0.3710425 ## --- ## Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
除此之外,还提供了kwAllPairsConoverTest()和kwAllPairsDunnTest()方法。
Friedman M检验及多重比较
使用课本本例8-9的数据,这个方式适用于随机区组设计资料的多样本比较。
进行Friedman M检验需要矩阵形式的数据(这个是R语言里为数不多的不支持formula形式的统计检验函数之一),可以自己输入,也可以直接读取spss格式数据,然后变成矩阵即可。
df <- matrix( c(8.4, 11.6, 9.4, 9.8, 8.3, 8.6, 8.9, 7.8, 9.6, 12.7, 9.1, 8.7, 8, 9.8, 9, 8.2, 9.8, 11.8, 10.4, 9.9, 8.6, 9.6, 10.6, 8.5, 11.7, 12, 9.8, 12, 8.6, 10.6, 11.4, 10.8 ), byrow = F, nrow = 8, dimnames = list(1:8,LETTERS[1:4]) ) print(df) ## A B C D ## 1 8.4 9.6 9.8 11.7 ## 2 11.6 12.7 11.8 12.0 ## 3 9.4 9.1 10.4 9.8 ## 4 9.8 8.7 9.9 12.0 ## 5 8.3 8.0 8.6 8.6 ## 6 8.6 9.8 9.6 10.6 ## 7 8.9 9.0 10.6 11.4 ## 8 7.8 8.2 8.5 10.8
进行Friedman M检验:
fit <- friedman.test(df) fit ## ## Friedman rank sum test ## ## data: df ## Friedman chi-squared = 15.152, df = 3, p-value = 0.001691
使用q检验(quade test)进行多重比较:
res <- quadeAllPairsTest(df) summary(res) ## ## Pairwise comparisons using Quades test with TDist approximation ## data: df ## P value adjustment method: holm ## H0 ## t value Pr(>|t|) ## B - A == 0 1.126 0.27307148 ## C - A == 0 3.236 0.01583866 * ## D - A == 0 5.093 0.00028885 *** ## C - B == 0 2.110 0.14099045 ## D - B == 0 3.967 0.00351141 ** ## D - C == 0 1.857 0.15475978 ## --- ## Signif. codes: 0 *** 0.001 ** 0.01 * 0.05 . 0.1 1
可以看到非常简单,函数名称清晰易懂,结果也是非常直观,直接给出了两两比较的P值和统计量。
除此之外,还提供了非常多其他方法:
-
frdAllPairsConoverTest() frdAllPairsExactTest() frdAllPairsConoverTest() frdAllPairsMillerTest() frdAllPairsSiegelTest()
大家选择自己需要的即可。
医学和生信笔记,专注R语言在临床医学中的使用,R语言数据分析和可视化。