额问问 发表于 2016-7-1 18:02 
运行错误,说没有BM。test这个函数。是不是需要装软件包啊?
BM.test 这个函数是前面那些代码,这个函数是自己编写的。
# Brown-Mood检验##################################################
BMq.test = function(x,y,q,alt)
{
xy = c(x,y)
quantile.xy = quantile(xy,q)
t = sum(xy>quantile.xy)
lx = length(x[x!=quantile.xy])
ly = length(y[y!=quantile.xy])
lxy = lx + ly
A = sum(x>quantile.xy) #检验统计量A
z = (A-lx*t/(lx+ly))/(lx*ly*t*(lx+ly-t)/(lx+ly)^3)^0.5
if(A>(lx*t/(lx+ly))){
z1 = (A+0.5-lx*t/(lx+ly))/(lx*ly*t*(lx+ly-t)/(lx+ly)^3)^0.5
#正态近似时的标准化统计量
}
else
{z1 = (A-0.5-lx*t/(lx+ly))/(lx*ly*t*(lx+ly-t)/(lx+ly)^3)^0.5}
if(alt == 'greater')
{
pv1 = 1-phyper(A,lx,ly,t)
pv2 = 1-pnorm(z)
pv3 = 1-pnorm(z1)
}
if(alt == 'less')
{
pv1 = phyper(A,lx,ly,t)
pv2 = pnorm(z)
pv3 = pnorm(z1)
}
if(alt == 'two.sided')
{
pv1 = 2*min(1-phyper(A,lx,ly,t),phyper(A,lx,ly,t))
pv2 = 2*min(1-pnorm(z),pnorm(z))
pv3 = 2*min(1-pnorm(z1),pnorm(z1))
}
conting.table = matrix(c(A,lx-A,lx,t-A,ly-(t-A),ly,t,lxy-t,lxy),3,3)
col.name = c('X','Y','X+Y')
row.name = c('>MQXY','<MQXY','TOTAL')
dimnames(conting.table) = list(row.name,col.name)
list(contingency.table = conting.table,p.value = pv1,pvnorm = pv2,pvnr = pv3)
}
x = c(10,8,12,16,5,9,7,11,6)
y = c(12,15,20,18,13,14,9,16)
xy = c(x,y)
median(xy)
BMq.test(x,y,0.25,'two.sided')