等分散じゃないからKruskal-Wallis検定ってダメだから(続き2)

あ前回、前々回で、どうやらKruskal-Wallis検定はt検定の代替手法としてはダメっぽいことを実験で示しました。KWは平均の差の検定でも中央値の差の検定でもありません

 じゃあ、もうこのテストが全然使えないかというと、一応使い道はあります……がかなり用途が限られます。結局、もともとの話に戻るのですが、「複数のグループで順位(ランク)の平均値が等しい」ことを示したいときは有力です。でも研究を行っていて、そんな場面はあまりなかったのですが。

 例によって、この仮説に従うサンプルをRで生成し、P値の分布がどうなるかを試してみましょう。
 擬似的に、順位平均が等しい2つのグループを作成します。この方法は、長いベクトルを作成して、最大と最小の値を順に組み合わせていけば可能です(Geminiに聞きました)。
 こうしてできた2グループについて、並べ替え検定(permutation test)を使います。それぞれのシミュレーション反復でKW検定のP値をストックしておき、それの累積分布関数(ecdf)を調べます。うまくいけば直線になります。また、第一種の過誤(差がないのにあると結論する誤り)は5%におさまるはずです。

# Kruskal-Wallis test 
# 並べ替えpermutationした場合のP値の分布
# 2026-01-29

# x, y サンプルサイズ
nx <- 100
ny <- 50
N <- nx + ny


# 母集団 x_ranks, y_ranks

all_ranks <- 1:N # 1からNまでのランクを用意

x_ranks <- c(1:(nx/2), (N - nx/2 + 1):N)  # x: 大きい・小さいほうから(nx/2)個とる
                       # xのランク平均はall_ranksの平均
mean_x <- mean(x_ranks)

y_ranks <- setdiff(all_ranks, x_ranks)  # y: x以外の要素を抽出
mean_y <- mean(y_ranks)					 # ランク平均はやはりall_ranksの平均

mean_x; mean_y   # 確認: x, y とも (1+N)/2


# Simulation

nrep <- 1000
res <- rep (NA, nrep)		# P.valuesの格納用
res_mean_x <- numeric (nrep) # Xの平均ランク

for (i in 1:nrep){
	
	shuffled <- sample(all_ranks, replace = FALSE)	# ランクをランダムに並べ替え

	x <- shuffled[1:nx]
	y <- shuffled[(nx + 1):N]
	
	res_mean_x[i] <- mean(x)		
	res[i] <- kruskal.test(list(x, y))$p.value
}

mean(res < 0.05)		# Type-I error


# Plots

par(mfrow = c(2,1))

hist(res, freq = FALSE, breaks = 20, col = "orange",
		 main = "P-value Frequency Distribution", xlab = "P-value")
abline (h = 1)

# CDFのプロット

plot(ecdf(res), 
		 main="CDF of P-value", 
		 xlab="P-value", 
		 ylab="Cumulative Probability",
		 col="blue",
		 cex= 0.2
)
abline(0, 1, col="red", lty=2) # 理論上の一様分布(赤線)

par(mfrow = c(1,1))

P値の累積分布は直線によくあてはまるようです。

累積分布関数(青:実測値、赤:理論値)

P値のヒストグラムは、予想どおり一様分布になっているようです。(前々回参照)

Pのヒストグラム
mean(res < 0.05)		# Type-I error

第一種の過誤は、ちょうど 0.05 になりました(結果は示していません)。
調べた範囲では、2つのグループの標本数が変わっても問題はありませんでした。そういうわけで、第一種の過誤はうまく調整されており、優秀な方法といえそうです。

なお、このコードはランクの分散について制御していません。
最初に作成したx, y のランクの分散は、
  V(x) = 2735.6
  V(y) =  212.5
と全く違う値になっています。にもかかわらず、KW検定はうまく機能しているようです。このコードでは permutation test で要素をグチャグチャに並べ替えているのだから、当たり前なのかもしれません。


いやはや、難しいですね。KW検定は使い道はあるが、検定の背景を正しく理解しなければいけないようです。「分散が違う、だめだ、KW検定だ」というのは戒めなければいけません。




いいなと思ったら応援しよう!