wilcox_R: What the package does (short line)

cid <-
function(x,y,alpha=.05,plotit=FALSE,pop=0,fr=.8,rval=15,xlab="",ylab=""){
#
#  Compute a confidence interval for delta using the method in
#  Cliff, 1996, p. 140, eq 5.12
#
#  The null hypothesis is that for two independent group, P(X<Y)=P(X>Y).
#  This function reports a 1-alpha confidence interval for
#  P(X>Y)-P(X<Y)
#
#  plotit=TRUE creates a plot of the difference scores.
#  pop=0 adaptive kernel density estimate
#  pop=1 results in the expected frequency curve.
#  pop=2 kernel density estimate (Rosenblatt's shifted histogram)
#  pop=3 boxplot
#  pop=4 stem-and-leaf
#  pop=5 histogram
#  pop=6  kernel density estimate
#
x<-x[!is.na(x)]
y<-y[!is.na(y)]
m<-outer(x,y,FUN="-")
msave<-m
m<-sign(m)
d<-mean(m)
phat<-(1-d)/2
flag=T
if(phat==0 || phat==1)flag=F
q0<-sum(msave==0)/length(msave)
qxly<-sum(msave<0)/length(msave)
qxgy<-sum(msave>0)/length(msave)
c.sum<-matrix(c(qxly,q0,qxgy),nrow=1,ncol=3)
dimnames(c.sum)<-list(NULL,c("P(X<Y)","P(X=Y)","P(X>Y)"))
if(flag){
sigdih<-sum((m-d)^2)/(length(x)*length(y)-1)
di<-NA
for (i in 1:length(x))di[i]<-sum(x[i]>y)/length(y)-sum(x[i]<y)/length(y)
dh<-NA
for (i in 1:length(y))dh[i]<-sum(y[i]>x)/length(x)-sum(y[i]<x)/length(x)
sdi<-var(di)
sdh<-var(dh)
sh<-((length(y)-1)*sdi+(length(x)-1)*sdh+sigdih)/(length(x)*length(y))
zv<-qnorm(alpha/2)
cu<-(d-d^3-zv*sqrt(sh)*sqrt((1-d^2)^2+zv^2*sh))/(1-d^2+zv^2*sh)
cl<-(d-d^3+zv*sqrt(sh)*sqrt((1-d^2)^2+zv^2*sh))/(1-d^2+zv^2*sh)
}
if(!flag){
sh=NULL
nm=max(c(length(x),length(y)))
if(phat==1)bci=binomci(nm,nm,alpha=alpha)
if(phat==0)bci=binomci(0,nm,alpha=alpha)
}
if(plotit){
if(pop==1 || pop==0){
if(length(x)*length(y)>2500){
print("Product of sample sizes exceeds 2500.")
print("Execution time might be high when using pop=0 or 1")
print("If this is case, might consider changing the argument pop")
#print("pop=2 might be better")
}}
if(pop==0)akerd(as.vector(msave),xlab=xlab,ylab=ylab)
if(pop==1)rdplot(as.vector(msave),fr=fr,xlab=xlab,ylab=ylab)
if(pop==2)kdplot(as.vector(msave),rval=rval,xlab=xlab,ylab=ylab)
if(pop==3)boxplot(as.vector(msave))
if(pop==4)stem(as.vector(msave))
if(pop==5)hist(as.vector(msave),xlab=xlab)
if(pop==6)skerd(as.vector(msave))
}
if(flag)pci=c((1-cu)/2,(1-cl)/2)
if(!flag){
pci=bci$ci
cl=1-2*pci[2]
cu=1-2*pci[1]
}
list(n1=length(x),n2=length(y),cl=cl,cu=cu,d=d,sqse.d=sh,phat=phat,summary.dvals=c.sum,ci.p=pci)
}