Se construye una simulación de la potencia de los tres tests mencionados anteriormente. Lo haremos considerando que los coeficientes de correlación toman valores en el intervalo \((-1,0)\), esto es, una cola a izquierda. Antes de enunciar el código que nos permite evidenciar los resultados gráficos, observemos unos ejemplos de los valores de la correlación bajo la hipótesis alterna.

Algunos ejemplos de como varía la correlación bajo la hipótesis alterna

set.seed(1)
x=rnorm(100,0,1)
y1= 2- 10*(x/1)+rnorm(50, 0, 4)
y2= 2- 10*(x/5)+rnorm(50, 0, 4)
y3= 2- 10*(x/100)+rnorm(50, 0, 4)

par(mfrow=c(1,3))

plot(x,y1)
plot(x,y2)
plot(x,y3)

Con lo cual vemos que los datos inician con una asociación negativa bastante pronunciada y luego esta va disminuyendo hasta aproximarse a una asociación nula. Lo confirmamos de la siguiente forma:

cor.test(x,y1,method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  x and y1
## t = -24.696, df = 98, p-value < 2.2e-16
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.9512013 -0.8949470
## sample estimates:
##        cor 
## -0.9282022
cor.test(x,y2,method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  x and y2
## t = -5.8523, df = 98, p-value = 6.437e-08
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.6412214 -0.3471828
## sample estimates:
##        cor 
## -0.5088977
cor.test(x,y3,method = "pearson")
## 
##  Pearson's product-moment correlation
## 
## data:  x and y3
## t = 0.069881, df = 98, p-value = 0.9444
## alternative hypothesis: true correlation is not equal to 0
## 95 percent confidence interval:
##  -0.1896222  0.2031952
## sample estimates:
##         cor 
## 0.007058857

Con esto en mente, en el siguiente código hacemos una simulación de tamaño 10.000 que nos permitirá determinar cuál test resulta más potente bajo el supuesto de normalidad.

prob_rechazo_resul1=NULL
prob_rechazo_resul2=NULL
prob_rechazo_resul3=NULL
z=seq(10,1,length=100)
tsim=10000
r_xy=seq(-1,-0.1,length=100)   ####coeficiente corr bajo la alternativa
for (j in 1:100) 
{
  valorp_pearson=NULL
  valorp_spearman=NULL
  valorp_kendall=NULL
  for (i in 1:tsim) 
  {
    x=rnorm(100,0,1)
    y= 2- 10*(x/z[j])+rnorm(50, 0, 4)
    resul1=cor.test(x,y, method="pearson",alternative = "less", conf.level = 0.95)
    resul2=cor.test(x,y,method="spearman",alternative = "less",conf.level=0.95)
    resul3=cor.test(x,y,method = "kendall",alternative = "less",conf.level = 0.95)
    valorp_pearson[i]=resul1$p.value
    valorp_spearman[i]=resul2$p.value
    valorp_kendall[i]=resul3$p.value
    
  }
  prob_rechazo_resul1[j]=sum(ifelse(valorp_pearson<0.05,1,0))/tsim
  prob_rechazo_resul2[j]=sum(ifelse(valorp_spearman<0.05,1,0))/tsim
  prob_rechazo_resul3[j]=sum(ifelse(valorp_kendall<0.05,1,0)/tsim)
}

plot(r_xy, prob_rechazo_resul1, type="l", col=2,
     ylab="Probailidad de rechazo", xlab="Hipótesis Alterna",
     main=" Pearson vs Spearman vs Kendall")
lines (r_xy, prob_rechazo_resul2, type="l", col=4)
lines(r_xy,prob_rechazo_resul3,type = "l",col=1)
legend("topleft",bty = "n",legend = c("Pearson","Spearman","Kendall"),lty = 1,col = c(2,4,1),cex=0.8)

Referencias