R/twoDcorR.R

twoDcorR <-
function(x,y,corfun=wincor,alpha=.05,nboot=500,SEED=TRUE,MC=FALSE){
#
# Comparing two robust dependent correlations: Overlapping case
# Winsorized correlation is used by default.
#
# x is assumed to be a matrix with 2 columns
#
#  Compare correlation of x[,1] with y to x[,2] with y
#
m1=cbind(x,y)
m1<-elimna(m1)  # Eliminate rows with missing values
nval=nrow(m1)
x<-m1[,1:2]
y=m1[,3]
est<-cor2xy(x,y,corfun=corfun)$cor 
r12=est[1]
r13=est[2]
if(SEED)set.seed(2) # set seed of random number generator so that
#             results can be duplicated.
#
#  If you use corfun=scor, set plotit=F
#
data<-matrix(sample(length(y),size=length(y)*nboot,replace=TRUE),nrow=nboot)
data=listm(t(data))
if(MC){
library(parallel)
bvec<-mclapply(data,twoDcorR_sub,x,y,corfun) 
}
if(!MC)bvec<-lapply(data,twoDcorR_sub,x,y,corfun) 
mat=matrix(NA,nrow=nboot,ncol=2)
for(i in 1:nboot)mat[i,]=bvec[[i]]
ihi<-floor((1-alpha/2)*nboot+.5)
ilow<-floor((alpha/2)*nboot+.5)
corhat=cor(mat^2)[1,2]
bsort<-sort(mat[,1]-mat[,2])
ci12<-1
ci12[1]<-bsort[ilow]
ci12[2]<-bsort[ihi]
pv=mean(bsort<0)
pv=2*min(c(pv,1-pv))
list(est.rho1=r12,est.rho2=r13,est.dif=r12-r13,ci=ci12,p.value=pv)
}
musto101/wilcox_R documentation built on May 23, 2019, 10:52 a.m.