MW<-function(x,y){ n=table(y) k=length(n) N=sum(n) mu<-function(l){ N1=sum(n[1:l]) N2=sum(n[l:k]) m=(N1^2+N2^2-sum(n^2)-n[l]^2)/4 sigma=(2*(N1^3+N2^3)+3*(N1^2+N2^2)-sum((n^2*(2*n+3)))-n[l]^2*(2*n[l]+3)+12*n[l]*N1*N2-12*n[l]^2*N)/72 return(c(m,sigma)) } U<-function(x1,x2){ n1<-length(x1) n2<-length(x2) u<-0 for(j in 1:n2){ u<-u+sum(x1<=x2[j]) } return(u) } A=array(0,dim=c(k,1)) m=array(0,dim=c(k,1)) s=array(0,dim=c(k,1)) for(l in 1:k){ t=mu(l) m[l]=t[1] s[l]=sqrt(t[2]) if(l>1){ for(i in 1:(l-1)){ for(j in (i+1):l){ A[l]=A[l]+U(x[y==i],x[y==j]) } } } if(l