Introduccion

Este es un test alternativo a la U Wilcoxon-Mann-Whitney y como estadístico, veo muchos procedimientos y los pruebo a ver su efectividad, en este caso probé el procedimiento brunnermunzel de la libreria LAWSTAT y encontré algunas diferencias y por eso produje el procedimiento BMTest de la librería ANALITICA (por ahora disponible en github) pero pronto en cran.

A continuación les muestro el ejemplo de la libreria lawstat

desarrollamos los pasos de forma manual en R

los datos

Y <- c(1, 2, 1, 1, 1, 1, 1, 1, 1, 1, 2, 4, 1, 1)
N <- c(3, 3, 4, 3, 1, 2, 3, 1, 1, 5, 4)

La fórmula para el cálculo del valor crítico de Brunner Munzel es

\(T_{BM} = \frac{\bar{R}_X - \bar{R}_Y}{\sqrt{\frac{S_X^2}{m} + \frac{S_Y^2}{n}}}\)

primero calcular los rangos entonces unir los vectores

datos<-c(Y,N)


asignar rangos, debiendo quedar de esta forma

##  [1]  7.5 16.0  7.5  7.5  7.5  7.5  7.5  7.5  7.5  7.5 16.0 23.0  7.5  7.5 19.5
## [16] 19.5 23.0 19.5  7.5 16.0 19.5  7.5  7.5 25.0 23.0
## Ry 7.5 16 7.5 7.5 7.5 7.5 7.5 7.5 7.5 7.5 16 23 7.5 7.5
## Rn 19.5 19.5 23 19.5 7.5 16 19.5 7.5 7.5 25 23

Calcular promedio y varianza de rangos por grupo

  Ry_m <- mean(Ry)
  Rn_m <- mean(Rn)

  cat ("mean(Ry)=",Ry_m)
## mean(Ry)= 9.821429
  cat ("mean(Rn)=",Rn_m)
## mean(Rn)= 17.04545
  Sy_sq <- sum((Ry - Ry_m)^2) / (14 - 1)
  Sn_sq <- sum((Rn - Rn_m)^2) / (11 - 1)
  
  cat ("var(y)=",Sy_sq)
## var(y)= 23.79258
  cat ("var(n)=",Sn_sq)
## var(n)= 43.27273
  S_sq <- (Sy_sq / 14) + (Sn_sq / 11)
  
  cat ("var total =",S_sq)
## var total = 5.633354
  T_stat <- (Ry_m - Rn_m) / sqrt(S_sq)
  
  cat ("Tbm =",T_stat)
## Tbm = -3.043657

Aquí la primera diferencia, en el Tbm

3.13 vs -3.04

los grados de libertad

gl= \(\frac{\left( \frac{S_1^2}{m} + \frac{S_2^2}{n} \right)^2}{\frac{S_1^4}{m^2 (m - 1)} + \frac{S_2^4}{n^2 (n - 1)}}\)

df<-(S_sq^2)/((Sy_sq^2 / (14^2 * (14 - 1))) + (Sn_sq^2 / (11^2 * (11 - 1))))

cat("df=",df)
## df= 17.9321

segunda diferencia

df 17.68 vs 17.93

pv2t<-2 * (1 - pt(abs(T_stat), df))

pvtr<-1 - pt(T_stat, df)

pvtl<- pt(T_stat, df)

cat("p-value 2t=",pv2t)
## p-value 2t= 0.00701071
cat("p-value greater =",pvtr)
## p-value greater = 0.9964946
cat("p-value less =",pvtl)
## p-value less = 0.003505355

tercera diferencia

ninguno de los p valor coincide

rb<-BMTest(Y,N)
summary(rb)
## ========================================
##   Summary of Pairwise Comparison Test
## ========================================
## Method: Brunner-Munzel (two.sided) 
## 
## Comparison          : Grupo1 - Grupo2
## Mean difference     : -1.3701
## Degrees of freedom  : 17.93
## Standard error      : 2.3735
## t critical value    : 2.1015
## p-value             : 0.007
## P(X < Y) + 0.5 P(X = Y): 0.789
## Significance        : **
## 
## Group means (ordered from highest to lowest):
##   Grupo2: 2.727
##   Grupo1: 1.357