R/TWOpNOVPV.R

TWOpNOVPV <-
function(x,y,HC4=TRUE,alpha=.05){
#
# Comparing two dependent correlations: Non-overlapping case
#
#   Compute a .95 confidence interval
#   for the difference between two dependent Pearson correlations,
#   non-overlapping case.
#
#    Both x and y are assumed to be matrices with two columns.
#   The function compares the correlation between x[,1] and x[,2]
#   to the correlation between y[,1] and y[,2].
#
#  For simulation results, see Wilcox (2009).
#  COMPARING PEARSON CORRELATIONS: DEALING WITH
#  HETEROSCEDASTICITY AND NON-NORMALITY, Communications in Statistics--Simulations
#   and Computations, 38, 2220-2234.
#
# This function is exactly like TWOpNOV, only it returns a p-value as well.
#
#  Note: To get a p-value, HC4=TRUE must be used. 
#
alph<-c(1:99)/100
for(i in 1:99){
irem<-i
chkit<-TWOpNOV(x,y,alpha=alph[i],HC4=TRUE)
chkit=c(chkit$ci.lower,chkit$ci.upper)
if(sign(chkit[1]*chkit[2])==1)break
}
p.value<-irem/100
if(p.value<=.1){
iup<-(irem+1)/100
alph<-seq(.001,iup,.001)
for(i in 1:length(alph)){
p.value<-alph[i]
alph<-c(1:99)/100
for(i in 1:99){
irem<-i
chkit<-TWOpNOV(x,y,alpha=alph[i],HC4=TRUE)
chkit=c(chkit$ci.lower,chkit$ci.upper)                 
if(sign(chkit[1]*chkit[2])==1)break   
}}}
p.value<-irem/100
if(p.value<=.1){
iup<-(irem+1)/100
alph<-seq(.001,iup,.001)
for(i in 1:length(alph)){
p.value<-alph[i]
chkit<-TWOpNOV(x,y,alpha=alph[i],HC4=TRUE)
chkit=c(chkit$ci.lower,chkit$ci.upper)        
if(sign(chkit[1]*chkit[2])==1)break   
}}
if(p.value<=.001){
alph<-seq(.0001,.001,.0001)
for(i in 1:length(alph)){
p.value<-alph[i]
chkit<-TWOpNOV(x,y,alpha=alph[i],HC4=TRUE)
chkit=c(chkit$ci.lower,chkit$ci.upper)        
if(sign(chkit[1]*chkit[2])==1)break  
}}
if(p.value<=.001){
alph<-seq(.0001,.001,.0001)
for(i in 1:length(alph)){
p.value<-alph[i]
chkit<-TWOpNOV(x,y,alpha=alph[i],HC4=TRUE)
chkit=c(chkit$ci.lower,chkit$ci.upper)        
if(sign(chkit[1]*chkit[2])==1)break  
}}
res=TWOpNOV(x,y,alpha=alpha,HC4=TRUE)
ci=c(res$ci.lower,res$ci.upper)
list(p.value=p.value,est.1=res$est.1,est.2=res$est.2,ci=ci) #ci.lower=res$ci.lower,ci.upper=res$ci.upper)
}
musto101/wilcox_R documentation built on May 23, 2019, 10:52 a.m.