R/twopcor.R In WRS2: A Collection of Robust Statistical Methods

```twopcor <- function(x1, y1, x2, y2, nboot = 599){
#
#   Compute a .95 confidence interval for
#   the difference between two Pearson
#   correlations corresponding to two independent
#   goups.
#
#   This function uses an adjusted percentile bootstrap method that
#   gives good results when the error term is heteroscedastic.
#
#   WARNING: If the number of boostrap samples is altered, it is
#   unknown how to adjust the confidence interval when n1+n2 < 250.
#
cl <- match.call()
if (nboot < 599) warning("It is unknown how to adjust the confidence interval when n1+n2 < 250.")
X<-elimna(cbind(x1,y1))
x1<-X[,1]
y1<-X[,2]
X<-elimna(cbind(x2,y2))
x2<-X[,1]
y2<-X[,2]

data1<-matrix(sample(length(y1),size=length(y1)*nboot,replace=TRUE),nrow=nboot)
bvec1 <- apply(data1, 1, function(xx) cor(x1[xx], y1[xx]))

data2<-matrix(sample(length(y2),size=length(y2)*nboot,replace=TRUE),nrow=nboot)
bvec2<-apply(data2,1,function(xx) cor(x2[xx], y2[xx]))

bvec<-bvec1-bvec2
ilow<-15
ihi<-584
if(length(y1)+length(y2) < 250){
ilow<-14
ihi<-585
}
if(length(y1)+length(y2) < 180){
ilow<-11
ihi<-588
}
if(length(y1)+length(y2) < 80){
ilow<-8
ihi<-592
}
if(length(y1)+length(y2) < 40){
ilow<-7
ihi<-593
}
bsort<-sort(bvec)
r1<-cor(x1,y1)
r2<-cor(x2,y2)
ci<-c(bsort[ilow],bsort[ihi])
result <- list(r1 = r1, r2 = r2, ci = ci, call = cl)
class(result) <- "twocor"
result
}
```

Try the WRS2 package in your browser

Any scripts or data that you put into this service are public.

WRS2 documentation built on May 2, 2019, 4:46 p.m.