wilcox_R: What the package does (short line)

depthg2 <-
function(x,y,alpha=.05,nboot=500,MD=FALSE,plotit=TRUE,op=FALSE,fast=FALSE,SEED=TRUE,
xlab="VAR 1",ylab="VAR 2"){
#
#   Compare two independent groups based on p measures
#   for each group.
#
#   The method is based on Tukey's depth if MD=F;
#   otherwise the Mahalanobis depth is used.
#   If p>2, then Mahalanobis depth is used automatically
#
#   The method is designed to be sensitive to differences in scale
#
if(SEED)set.seed(2) # set seed of random number generator so that
#             results can be duplicated.
x=elimna(x)
y=elimna(y)
x=as.matrix(x)
y=as.matrix(y)
if(is.matrix(x) && is.matrix(y)){  # YES, code is odd.
nv1<-nrow(x)
nv2<-nrow(y)
if(ncol(x)!=ncol(y))stop("Number of columns of x is not equal to number for y")
if(ncol(x) >2)MD<-T
if(ncol(x)==2 && plotit){
plot(rbind(x,y),type="n",xlab=xlab,ylab=ylab)
points(x,pch="*")
points(y,pch="o")
temp<-NA
for(i in 1:nrow(x)){
temp[i]<-depth(x[i,1],x[i,2],x)
}
flag<-(temp>=median(temp))
xx<-x[flag,]
xord<-order(xx[,1])
xx<-xx[xord,]
temp<-chull(xx)
xord<-order(xx[,1])
xx<-xx[xord,]
temp<-chull(xx)
lines(xx[temp,])
lines(xx[c(temp[1],temp[length(temp)]),])
temp<-NA
for(i in 1:nrow(y)){
temp[i]<-depth(y[i,1],y[i,2],y)
}
flag<-(temp>=median(temp))
xx<-y[flag,]
xord<-order(xx[,1])
xx<-xx[xord,]
temp<-chull(xx)
flag<-(temp>=median(temp))
xord<-order(xx[,1])
xx<-xx[xord,]
temp<-chull(xx)
lines(xx[temp,],lty=2)
lines(xx[c(temp[1],temp[length(temp)]),],lty=2)
}
print("Taking bootstrap samples. Please wait.")
data1<-matrix(sample(nv1,size=nv1*nboot,replace=TRUE),nrow=nboot)
data2<-matrix(sample(nv2,size=nv2*nboot,replace=TRUE),nrow=nboot)
qhatd<-NA
dhatb<-NA
for(ib in 1:nboot){
if(op)print(paste("Bootstrap sample ",ib," of ",nboot, "is complete."))
if(!fast)temp<-lsqs2(x[data1[ib,],],y[data2[ib,],],plotit=FALSE,MD=MD)
if(fast)temp<-lsqs2.for(x[data1[ib,],],y[data2[ib,],],plotit=FALSE,MD=MD)
qhatd[ib]<-temp[[1]]-temp[[2]]
}
temp<-sort(qhatd)
lv<-round(alpha*nboot/2)
uv<-nboot-lv
difci<-c(temp[lv+1],temp[uv])
}
#
if(!is.matrix(x) && !is.matrix(y)){
nv1<-length(x)
nv2<-length(y)
print("Taking bootstrap samples. Please wait.")
data1<-matrix(sample(nv1,size=nv1*nboot,replace=TRUE),nrow=nboot)
data2<-matrix(sample(nv2,size=nv2*nboot,replace=TRUE),nrow=nboot)
qhatd<-NA
dhatb<-NA
for(ib in 1:nboot){
if(!fast)temp<-lsqs2(x[data1[ib,]],y[data2[ib,]],plotit=FALSE,MD=MD)
if(fast)temp<-lsqs2.for(x[data1[ib,]],y[data2[ib,]],plotit=FALSE,MD=MD)
qhatd[ib]<-temp[[1]]-temp[[2]]
dhatb[ib]<-(temp[[1]]+temp[[2]])/2
#print(paste("Bootstrap sample ",ib," of ",nboot, "is complete."))
}}
temp<-sort(qhatd)
temp2<-sort(dhatb)
lv<-round(alpha*nboot/2)
uv<-nboot-lv
difci<-c(temp[lv+1],temp[uv])
list(difci=difci)
}