1.npmc包的最后更新版本（1.0.7）无法正确运行在R 3.0以上的版本。HK Zhang 在Rstudio中做了调试和编译，对npmc.R的原代码做了一些细微的更改以支持最新的mvtnorm包。npmc包最后编译运行成功。（注：R>4.0 目前也不行了，3.4可行）
最新原程序共享在GitHub:
https://github.com/HK-Zhang/Rice/tree/master/npmc
也可以通过下面链接直接下载这个包
http://pan.baidu.com/s/1nuqdXcX
2.单因素方差分析（单因素ANOVA）后接多重比较（Tukey HSD检验），是比较多组数据两两差异的常规做法。不过由于单因素ANOVA的前提假设条件比较严格，要求数据必须满足正态性、方差齐性等，而很多情况下我们的数据并不符合ANOVA的条件，因此多重比较也无法直接使用。此时除了考虑转换数据外，另一种思路就是使用非参数的方法。
对于非参数的方法，有两种

1. 此前提到过可以首先执行Kruskal-Wallis检验代替单因素ANOVA比较多组数据的整体差异，然后再通过Wilcoxon检验代替Tukey HSD继续确定两组间差异水平，并通过p值校正控制I类错误率。

2. 另一种常见做法则是使用非参数的多重比较。《R语言实战 第一版》（154页）提到了一个R包，npmc，提供了一种非参数多组比较程序npmc()，可在控制I类错误的前提下，进行所有组之间的成对比较。

1.示例数据
npmc的内置数据集brain，来自某项对大脑疾病的研究，左半球（LH）或右半球（RH）病变对反应时间的影响。

library(npmc)
#数据集，详情 ?brain
data(brain)
head(brain)
tail(brain)

7SH5LX{7SZFY78Z38S{%~YQ.png

npmc的非参数多重比较

在进行非参数多重比较之前，可首先使用Kruskal-Wallis检验查看整体差异。

#Kruskal-Wallis 检验，整体差异
fit <- kruskal.test(var~class, data = brain)
fit$p.value

image.png

结果显示3组数据在整体上存在显著差异。

接下来，继续使用npmc包中的方法进行两组间的比较。

npmc()函数中提供了两种非参数多重比较方法，Behrens-Fisher和Steel，均可用于成对和多对一的比较情况，并同时计算相对效应的置信区间（1-alpha）。
npmc()函数接受的输入为一个两列的数据框，其中一列名为class（分组变量），另一列名为var（因变量）。因此在输入数据前，##更改变量名称保证能够被npmc()识别 (一定要更改，否则报错)。##
npmc包的功能实现需要mvtnorm包计算检验统计量的p值，如果mvtnorm包不可用，则结果将包含NA的p值。mvtnorm包中的函数使用随机值进行积分计算，因此npmc()关于p值和置信区间的结果因调用而异，并且只能被视为近似解。

#npmc 的非参数多重比较
npmc_test <- npmc(brain, df = 2, alpha = 0.05)
#概要
summary(npmc_test, type = 'both', short = FALSE)
#或者
npmc_test$test

image.png

对于结果输出概要的解释：

BF/Steel，分别为Behrens-Fisher和Steel的检验结果；
cmp，比较的组名；
gn，两组样本量之和；
effect，相对效应估计值；
variance，方差估计值；
std，标准偏差；
statistic，检验统计量；
p-value 1s，单侧检验p值，与双侧检验相比它只考虑一端情况，若a-b比较组中a和b存在显著差异则代表a显著小于b；
p-value 2s，双侧检验p值；
zero，TRUE代表0方差或近似0方差，否则为FALSE。

一般情况下，关注双侧检验的结果就可以了。
对于本示例，尽管Kruskal-Wallis检验显示数据整体存在显著差异，但p值比较接近0.05显著性阈值了，暗示两组间的差异可能不明显。
非参数多重比较结果中，Behrens-Fisher显示l组（LH脑损伤的患者）和r组（RH脑损伤的患者）的反应时间存在区别，但不明显；Steel则显示无差异。
最后画个箱线图描述数据分布特征。

#ggplot2 箱线图
library(ggplot2)
 
p <- ggplot(data = brain, aes(x = class, y = var, fill = class)) +
geom_boxplot(outlier.size = 1) +
theme(panel.grid = element_blank(), panel.background = element_rect(fill = 'transparent', color = 'black'),
    legend.title = element_blank(), legend.key = element_blank(), plot.title = element_text(hjust = 0.5)) +
labs(x = '', y = '', title = 'Behrens-Fisher\n')
 
p +
annotate('segment', x = 1, xend = 2, y = 50, yend = 50) +
annotate('segment', x = 1, xend = 3, y = 55, yend = 55) +
annotate('segment', x = 2, xend = 3, y = 60, yend = 60) +
annotate('text', x = 1.5, y = 52, label = 'p = 0.147') +
annotate('text', x = 2, y = 57, label = 'p = 0.754') +
annotate('text', x = 2.5, y = 62, label = 'p = 0.044')

image.png

参考： 小白鱼的生统笔记（微信公众号）