ウィルコクソンの順位和検定、マン・ホイットニーのU検定

正規近似による方法
マン・ホイットニーのU検定

正規近似による方法

データの準備と要約

data <- c(
  12, 4, 9, 6, 2, 10, 12, 11, 7, 10,
  10, 10, 25, 13, 18, 7, 15, 9, 18, 19)
group <- c(rep("X",10), rep("Y",10))
dat_xy <- data.frame(data, group)
#View(dat_xy)

順位を追加

#dataの順位列 (rank) を追加
dat_xy2 <- dat_xy %>%
  mutate(rank = rank(data, ties.method = "average"))
#順位でソート
dat_xy3 <- dat_xy2 %>%
  arrange(rank)
#View(dat_xy3)

片側検定（帰無仮説Y群>X群）

ここからは、coinパッケージのwilcox_testを使用します

#groupをファクター変数に変更してYを基準に変更
dat_xy$group <- factor(dat_xy$group)

yを基準にして片側検定

dat_xy$group <- relevel(dat_xy$group, ref="Y") 
coin::wilcox_test(
  data ~ group,
  data = dat_xy,
  alternative = "greater"
)

有意水準5％のもとで、YはXより優位に高いことが分かりました

wilcox_testの Z値の求め方

Z=2.2798 の算出方法について解説します

X群とY群の順位和 $(W)$

rank_sum <- dat_xy3 %>%<br>  group_by(group) %>%<br>  summarise(sum_rank = sum(rank))<br>print(rank_sum)

X群の順位和: $W_x=75$
Y群の順位和: $W_y=135$
X群のサイズ: $n_x$
Y群のサイズ: $n_y$
サンプルサイズ: $N=n_x+n_y$

Z値の求め方

Y群を基準として記載

順位和の期待値 $E[W]=\dfrac{n_y(n+1)}{2}$

順位和の分散 $V[W]=\dfrac{n_xn_y(ｎ+1)}{12}$

Z値

$z=\dfrac{W-E[W]}{\sqrt{V[W]}}$

# グループ別のランク和の計算
Wx = 75
Wy = 135

# サンプルサイズ
nx <- sum(group == "X")
ny <- sum(group == "Y")
n <- nx + ny

# U統計量の計算
# Y群を基準にしているので、ここではWyを使用します


#期待値
EW <- ny*(n+1)/2
#分散
VW <-  nx*ny*(n+1)/12
# Z値の計算
Zw <- (Wy - EW) / sqrt(VW)

# 結果の出力
Zw

#p値
pnorm(Zw, lower.tail = FALSE)

coinパッケージのwilcox_testの答えと微妙に違うようです・・・？
タイデータを考慮したマン・ホイットニーU検定を実行してみましょう

マン・ホイットニーのU検定

タイデータのがあるので補正したZ値を求めます

マン・ホイットニーのU統計量　$U=n_xn_y+\dfrac{n_y(n_y+1)}{2}-W_y$

順位和の期待値 $E[U]=\dfrac{n_xn_y}{2}$

補正項 $T=\sum{t_i^3-t_i}$ （t=タイデータの個数）

分散　$V[U]=\dfrac{n_xn_y(N+1)}{12}-\dfrac{n_xn_y}{12N(N-1)}T$

検定統計量

$z=\dfrac{|U-E[U]|}{\sqrt{V[U]}}$

# グループ別のランク和の計算
Wx = 75
Wy = 135

# サンプルサイズ
nx <- sum(group == "X")
ny <- sum(group == "Y")
n <- nx + ny

# U統計量の計算
# Y群を基準にしているので、ここではWyを使用します
U <- nx*ny+ny*(ny+1)/2-Wy

#期待値
EU <- nx*ny/2

# タイデータの個数
t　<- table(dat_xy3$rank) 
# タイの修正項
T <- sum(t^3 - t) 
#分散
VU <- (nx*ny/12)*((n+1)-T/(n*(n-1)))

# Z値の計算
Zu <- abs(U - EU) / sqrt(VU)

# 結果の出力
Zu

#p値
pnorm(Zu, lower.tail = FALSE)

これでwilcox_testのZ値、p値と同値になりました

前半で使用した　wilcox.test　を使用する場合には、correct =FALSE とすることで連続性の補正が省かれてマン・ホイットニーU検定を実行することができます

wilcox.test(
  data ~ group,
  data = dat_xy,
  alternative = "greater",
  correct =FALSE
)

ｐ値は同値となりました。タイがあるので正確検定ができず”警告メッセージ”が記載されますが問題ありません（correct =FALSEにしているので、このメッセージは不要かと思うのですが、、、）。