PH525x series - Rank tests

Wilcoxon Rank Sum Test

由于样本均数、标准差对离群值很敏感，而t检验是基于这些统计量进行的，所以t检验也对离群值比较敏感。此时，可以选择使用Wilcoxon rank test进行统计检验。Wilcoxon rank test检验的基本思想是：

1）合并所有数据；
2）将数据转换为秩次；
3）将秩次放回原组中；
4）计算秩次之和或平均值并进行检验。

具体过程：

set.seed(779) ##779 picked for illustration purposes
#生成满足正态分布的2组数据，每组25个值
N=25
x<- rnorm(N,0,1)
y<- rnorm(N,0,1)
#更改x的某个值，模拟离群值
x[1] <- 5
x[2] <- 7

library(rafalib)
mypar(1,2)
#绘制点带图，就是把点绘制在一条线段
stripchart(list(x,y),vertical=TRUE,ylim=c(-7,7),ylab="Observations",pch=21,bg=1)
abline(h=0)
#将x和y混在一起进行排序，返回每个值得秩次后放回原来的组
xrank <- rank(c(x,y))[seq(along=x)]
yrank <- rank(c(x,y))[-seq(along=y)]
#绘制点带图
stripchart(list(xrank,yrank),vertical=TRUE,ylab="Ranks",pch=21,bg=1,cex=1.25)
#将x中的元素一个一个拿出来与y放一起进行排序，求x的秩次，但是为啥减一？？？？？
ws <- sapply(x,function(z) rank(c(z,y))[1]-1)
#在点带图上添加文字，文字内容是ws
text( rep(1.05,length(ws)), xrank, ws, cex=0.8)

stripchart.png

W <- sum(ws)

W是第一组中每个元素相对第二组的秩次之和，我们可以基于组合数学去计算W的p值，还可以利用CLT理论去计算，因为根据CLT理论，W大致满足正态分布。我们可以构建一个z-score，如下：

n1<-length(x);n2<-length(y)
Z <- (mean(ws)-n2/2)/ sqrt(n2*(n1+n2+1)/12/n1)
print(Z)
## [1] 1.523124

由于Z不够大所以得出的p值大于0.05，上述过程便是R中wilcox.test函数的计算过程的一部分。

参考文章