R/Dancovapb.R

Dancovapb <-
function(x1,y1,x2,y2,fr1=1,fr2=1,est=hd,alpha=.05,plotit=TRUE,pts=NA,sm=FALSE,xout=FALSE,outfun=out,DIF=FALSE,...){
#
# Compare two dependent  groups using the ancova method 
# (a method similar to the one used by the  R function ancova).
# No parametric assumption is made about the form of
# the regression lines--a running interval smoother is used.
#
# percentile bootstrap method is used. 
#
# est indicates estimator to be used; Harrell-Davis median estimator is default.
#
#  Assume data are in x1 y1 x2 and y2
#
#  sm=T will create smooths using bootstrap bagging.
#  pts can be used to specify the design points where the regression lines
#  are to be compared.
#  pts=NA means five points will be picked empirically. 
#
#
if(ncol(as.matrix(x1))>1)stop('One covariate only is allowed with this function')
if(length(x1)!=length(y1))stop('x1 and y1 have different lengths')
if(length(x1)!=length(x2))stop('x1 and y2 have different lengths')
if(length(x2)!=length(y2))stop('x2 and y2 have different lengths')
if(length(y1)!=length(y2))stop('y1 and y2 have different lengths')
xy=elimna(cbind(x1,y1,x2,y2))
x1=xy[,1]
y1=xy[,2]
x2=xy[,3]
y2=xy[,4]
if(is.na(pts[1])){
npt<-5
isub<-c(1:5)  # Initialize isub
test<-c(1:5)
xorder<-order(x1)
y1<-y1[xorder]
x1<-x1[xorder]
xorder<-order(x2)
y2<-y2[xorder]
x2<-x2[xorder]
n1<-1
n2<-1
vecn<-1
for(i in 1:length(x1))n1[i]<-length(y1[near(x1,x1[i],fr1)])
for(i in 1:length(x1))n2[i]<-length(y2[near(x2,x1[i],fr2)])
for(i in 1:length(x1))vecn[i]<-min(n1[i],n2[i])
sub<-c(1:length(x1))
isub[1]<-min(sub[vecn>=12])
isub[5]<-max(sub[vecn>=12])
isub[3]<-floor((isub[1]+isub[5])/2)
isub[2]<-floor((isub[1]+isub[3])/2)
isub[4]<-floor((isub[3]+isub[5])/2)
mat<-matrix(NA,5,7)
dimnames(mat)<-list(NULL,c('X','n','DIF','ci.low','ci.hi','p.value','p.crit'))
for (i in 1:5){
t1=near(x1,x1[isub[i]],fr1)
t2=near(x2,x1[isub[i]],fr2)
pick=as.logical(t1*t2)
test=rmmcppb(y1[pick],y2[pick],est=est,dif=DIF,plotit=FALSE,alpha=alpha,pr=FALSE,SEED=FALSE,...)
mat[i,1]<-x1[isub[i]]
mat[i,2]<-length(y1[pick])
mat[i,3]<-test$output[,2]
mat[i,3]<-test$output[,2]
mat[i,4]<-test$output[,5]
mat[i,5]<-test$output[,6]
mat[i,6]<-test$output[,3]
}
temp2<-order(0-mat[,6])
bot=c(1:nrow(mat))
dvec=sort(alpha/bot,decreasing=TRUE)
mat[temp2,7]=dvec
}
if(!is.na(pts[1])){
n1<-1
n2<-1
vecn<-1
for(i in 1:length(pts)){
n1[i]<-length(y1[near(x1,pts[i],fr1)])
n2[i]<-length(y2[near(x2,pts[i],fr2)])
}
# First check sample size
#
flage=rep(TRUE,length(pts))
for (i in 1:length(pts)){
t1<-near(x1,pts[i],fr1)
t2<-near(x2,pts[i],fr2)
pick=as.logical(t1*t2)
if(sum(pick)<=5){print(paste('Warning: there are',sum(pick),' points corresponding to the design point X=',pts[i]))
flage[i]=FALSE
}}
pts=pts[flage]
mat<-matrix(NA,length(pts),7)
dimnames(mat)<-list(NULL,c('X','n','DIF','ci.low','ci.hi',
'p.value','p.crit'))
for (i in 1:length(pts)){
t1<-near(x1,pts[i],fr1)
t2<-near(x2,pts[i],fr2)
pick=as.logical(t1*t2)
#print(y1[pick])
test=rmmcppb(y1[pick],y2[pick],est=est,dif=DIF,plotit=FALSE,alpha=alpha,pr=FALSE,SEED=FALSE,...)
mat[i,3]<-test$output[,2]
mat[i,1]<-pts[i]
mat[i,2]<-length(y1[pick])
mat[i,4]<-test$output[,5]
mat[i,5]<-test$output[,6]
mat[i,6]<-test$output[,3]
}
temp2<-order(0-mat[,6])
bot=c(1:nrow(mat))
dvec=sort(alpha/bot,decreasing=TRUE)
mat[temp2,7]=dvec
}
if(plotit){
runmean2g(x1,y1,x2,y2,fr=fr1,est=est,sm=sm,xout=xout,outfun=outfun,,...)
}
list(output=mat)
}
musto101/wilcox_R documentation built on May 23, 2019, 10:52 a.m.