wilcox_R: What the package does (short line)

ddep <-
function(x,est=onestep,alpha=.05,grp=NA,nboot=500,plotit=TRUE,SEED=TRUE,pr=TRUE,WT=TRUE,...){
#
#   Do ANOVA on dependent groups
#   using the partially centered method plus
#   depth of zero among  bootstrap values.
#
#   An improved version of ddep that better handles heteroscedasticity
#   (A weighted grand mean is used in this version.)
#
#   The data are assumed to be stored in x in list mode
#   or in a matrix. In the first case
#   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.
#   If stored in a matrix, columns correspond to groups.
#
#   grp is used to specify some subset of the groups, if desired.
#   By default, all J groups are used.
#
#   The default number of bootstrap samples is nboot=2000
#
if(pr)print("Warning: Might not be level robust if the number of groups is relatively large and n is small")
if(pr)print("Currently seems that rmmismcp is preferable")
if(is.list(x)){
nv<-NA
for(j in 1:length(x))nv[j]<-length(x[[j]])
if(var(nv) !=0){
stop("The groups are stored in list mode and appear to have different sample sizes")
}
temp<-matrix(NA,ncol=length(x),nrow=nv[1])
for(j in 1:length(x))temp[,j]<-x[[j]]
x<-temp
}
J<-ncol(x)
if(!is.na(grp[1])){ #Select the groups of interest
J<-length(grp)
for(j in 1:J)temp[,j]<-x[,grp[j]]
x<-temp
}
x<-elimna(x) # Remove any rows with missing values.
bvec<-matrix(0,ncol=J,nrow=nboot)
hval<-vector("numeric",J)
if(SEED)set.seed(2) # set seed of random number generator so that
#             results can be duplicated.
print("Taking bootstrap samples. Please wait.")
n<-nrow(x)
totv<-apply(x,2,est,na.rm=TRUE,...)
data<-matrix(sample(n,size=n*nboot,replace=TRUE),nrow=nboot)
for(ib in 1:nboot)bvec[ib,]<-apply(x[data[ib,],],2,est,na.rm=TRUE,...) #nboot by J matrix
if(!WT){
gv<-rep(mean(totv),J) #Grand mean
#m1<-rbind(bvec,gv)
}
bplus<-nboot+1
center<-totv
cmat<-var(bvec)
if(WT){
wt=1/diag(cmat)
ut=sum(wt)
gv<-rep(sum(wt*totv)/ut,J) #Grand mean
}
m1<-rbind(bvec,gv)
discen<-mahalanobis(m1,totv,cmat)
#print("Bootstrap complete; computing significance level")
if(plotit && ncol(x)==2){
plot(bvec,xlab="Group 1",ylab="Group 2")
temp.dis<-order(discen[1:nboot])
ic<-round((1-alpha)*nboot)
xx<-bvec[temp.dis[1:ic],]
xord<-order(xx[,1])
xx<-xx[xord,]
temp<-chull(xx)
lines(xx[temp,])
lines(xx[c(temp[1],temp[length(temp)]),])
abline(0,1)
}
sig.level<-sum(discen[bplus]<=discen)/bplus
list(p.value=sig.level,center=totv,grand.mean=gv)
}