wilcox_R: What the package does (short line)

medr <-
function(x,est=median,alpha=.05,nboot=500,grp=NA,op=1,MM=FALSE,cop=3,pr=TRUE,
SEED=TRUE,...){
#
#   Test the hypothesis that the distribution for each pairwise
#   difference has a measure of location = 0
#   By default, the  median estimator is used
#
#   Independent groups are assumed.
#
#   The data are assumed to be stored in x in list mode or in a matrix.
#   If stored in list mode,
#   x[[1]] contains the data for the first group, x[[2]] the data
#   for the second group, etc. Length(x)=the number of groups = J, say.
#   If stored in a matrix, columns correspond to groups.
#
#   By default, all pairwise differences are used, but contrasts
#   can be specified with the argument con.
#   The columns of con indicate the contrast coefficients.
#   Con should have J rows, J=number of groups.
#   For example, con[,1]=c(1,1,-1,-1,0,0) and con[,2]=c(,1,-1,0,0,1,-1)
#   will test two contrasts: (1) the sum of the first
#   two measures of location is
#   equal to the sum of the second two, and (2) the difference between
#   the first two is equal to the difference between the
#   measures of location for groups 5 and 6.
#
#   The default number of bootstrap samples is nboot=500
#
#   op controls how depth is measured
#   op=1, Mahalanobis
#   op=2, Mahalanobis based on MCD covariance matrix
#   op=3, Projection distance
#   op=4, Projection distance using FORTRAN version
#
#   for arguments MM and cop, see pdis.
#
if(is.matrix(x)){
xx<-list()
for(i in 1:ncol(x)){
xx[[i]]<-x[,i]
}
x<-xx
}
if(!is.list(x))stop("Data must be stored in list mode or in matrix mode.")
if(!is.na(grp)){  # Only analyze specified groups.
xx<-list()
for(i in 1:length(grp))xx[[i]]<-x[[grp[1]]]
x<-xx
}
J<-length(x)
mvec<-NA
for(j in 1:J){
temp<-x[[j]]
temp<-temp[!is.na(temp)] # Remove missing values.
x[[j]]<-temp
mvec[j]<-est(temp,...)
}
Jm<-J-1
d<-(J^2-J)/2
data<-list()
bvec<-matrix(NA,ncol=d,nrow=nboot)
if(SEED)set.seed(2) # set seed of random number generator so that
#             results can be duplicated.
if(pr)print("Taking bootstrap samples. Please wait.")
for(it in 1:nboot){
for(j in 1:J)data[[j]]<-sample(x[[j]],size=length(x[[j]]),replace=TRUE)
dval<-0
for(j in 1:J){
for(k in 1:J){
if(j<k){
dval<-dval+1
bvec[it,dval]<-loc2dif(data[[j]],data[[k]],est=est,...)
}}}}
output<-matrix(NA,nrow=d,ncol=3)
dimnames(output)<-list(NULL,c("Group","Group","psihat"))
tvec<-NA
dval<-0
for(j in 1:J){
for(k in 1:J){
if(j<k){
dval<-dval+1
output[dval,1]<-j
output[dval,2]<-k
tvec[dval]<-loc2dif(x[[j]],x[[k]],est=est,...)
output[dval,3]<-tvec[dval]
}}}
tempcen<-apply(bvec,1,mean)
vecz<-rep(0,d)
smat<-var(bvec-tempcen+tvec)
temp<-bvec-tempcen+tvec
bcon<-rbind(bvec,vecz)
if(op==1)dv<-mahalanobis(bcon,tvec,smat)
if(op==2){
smat<-cov.mcd(temp)$cov
dv<-mahalanobis(bcon,tvec,smat)
}
if(op==3){
print("Computing p-value. Might take a while with op=3")
dv<-pdis(bcon,MM=MM,cop=cop,center=tvec)
}
if(op==4)dv<-pdis.for(bcon,MM=MM,cop=cop,pr=FALSE,center=tvec)
bplus<-nboot+1
sig.level<-1-sum(dv[bplus]>=dv[1:nboot])/nboot
if(op==4)print(sig.level)
list(sig.level=sig.level,output=output)
}