(takes about half an hour to run without cache, most of which is evaluating chunks dosomepvals and dosomepvalsNULL).

To cite the hyper2 package in publications, please use Hankin (2017). Here I discuss the hyper2 generalization of the Mann-Whitney-Wilcoxon test, applying it to the chorioamnion dataset at wilcox.test.Rd:

x <- c(0.80, 0.83, 1.89, 1.04, 1.45, 1.38, 1.91, 1.64, 0.73, 1.46)
y <- c(1.15, 0.88, 0.90, 0.74, 1.21)

First do the Mann-Whitney-Wilcox test:

wilcox.test(x,y)
## 
##  Wilcoxon rank sum exact test
## 
## data:  x and y
## W = 35, p-value = 0.25
## alternative hypothesis: true location shift is not equal to 0

Now use hyper3:

names(x) <- rep("x",length(x))
names(y) <- rep("y",length(y))
equalp.test(ordervec2supp3(names(sort(c(x,y)))))
## 
##  Constrained support maximization
## 
## data:  ordervec2supp3(names(sort(c(x, y))))
## null hypothesis: x = y
## null estimate:
##   x   y 
## 0.5 0.5 
## (argmax, constrained optimization)
## Support for null:  -27.899 + K
## 
## alternative hypothesis:  sum p_i=1 
## alternative estimate:
##       x       y 
## 0.24016 0.75984 
## (argmax, free optimization)
## Support for alternative:  -26.484 + K
## 
## degrees of freedom: 1
## support difference = 1.4148
## p-value: 0.092537

1 Some permutation tests under the null:

n <- 100
a <- matrix(0,n,2)
colnames(a) <- c("mww","hyper3")
for(i in seq_len(n)){
  d <- rrace3(pn=c(a=10,b=10),ps=c(a=0.5,b=0.5))  #null 
  jj <- seq_along(d)
  a[i,1] <- wilcox.test(jj[d=="a"],jj[d=="b"])$p.value
  a[i,2] <- equalp.test(ordervec2supp3(d))$p.value
}
head(a)
##          mww  hyper3
## [1,] 0.39305 0.32161
## [2,] 0.85343 0.42308
## [3,] 0.52885 0.64111
## [4,] 0.52885 0.25089
## [5,] 0.85343 0.53734
## [6,] 0.79594 0.83824
par(pty="s")
plot(log10(a),asp=1,xlim=c(-2,0),ylim=c(-2,0),main="p-values under the NULL (log scale)",
xlab="log(p),Wilcox",ylab="log(p), hyper3")
abline(0,1)

plot(a,asp=1,xlim=c(0,1),ylim=c(0,1),main="p-values under the NULL (linear scale)",xlab="p,Wilcox",ylab="p,hyper3")
abline(0,1)

(size_wilcox <- sum(a[,1]<0.05)/nrow(a))
## [1] 0.04
(size_hyper3 <- sum(a[,2]<0.05)/nrow(a))
## [1] 0.03
table(type_I_wilcox=a[,1]<0.05,type_I_hyper3=a[,2]<0.05)
##              type_I_hyper3
## type_I_wilcox FALSE TRUE
##         FALSE    94    2
##         TRUE      3    1

2 Some permutation tests under the alternative

n <- 100
a <- matrix(0,n,2)
colnames(a) <- c("mww","hyper3")
for(i in seq_len(n)){
  d <- rrace3(pn=c(a=10,b=10),ps=c(a=0.2,b=0.8))
  jj <- seq_along(d)
  a[i,1] <- wilcox.test(jj[d=="a"],jj[d=="b"])$p.value
  a[i,2] <- equalp.test(ordervec2supp3(d))$p.value
}
head(a)
##           mww    hyper3
## [1,] 0.279861 0.0600666
## [2,] 0.011496 0.0026632
## [3,] 0.105122 0.0352050
## [4,] 0.005196 0.0027129
## [5,] 0.005196 0.0053327
## [6,] 0.005196 0.0017941
par(pty="s")
plot(log10(a),asp=1,xlim=c(-4,0),ylim=c(-4,0))
abline(0,1)

(power_wilcox <- sum(a[,1]<0.05)/nrow(a))
## [1] 0.64
(power_hyper3 <- sum(a[,2]<0.05)/nrow(a))
## [1] 0.78
table(wilcox=a[,1]<0.05,hyper3=a[,2]<0.05)
##        hyper3
## wilcox  FALSE TRUE
##   FALSE    21   15
##   TRUE      1   63

References

Hankin, R. K. S. 2017. “Partial Rank Data with the hyper2 Package: Likelihood Functions for Generalized Bradley-Terry Models.” The R Journal 9 (2): 429–39.