ウィルコクソンの符号順位和検定(概要とpython実装)

ウィルコクソンの符号順位和検定について

  • ノンパラメトリック検定の一つ
  • サインランク検定ともいう
  • 順位尺度を用いており、データ間における代表値(中央値)の差を検定する方法
  • 対応のあるデータにおいて、データに正規性を仮定できない(データが正規分布してない)時に用いられる。(下図参照)
  • ウィルコクソンの順位和検定とは異なる検定法なので注意 (代表値の差を検定する部分は一緒だが、ウィルコクソンの順位和検定は対応のないデータに対する検定法)

コード

違いに対応のあるデータであり、それぞれ正規分布が仮定できないとする。
帰無仮説を「2群の代表値に差がないこと」と設定し、有意水準0.05で検定を行う。

import numpy as np
from scipy import stats
import seaborn as sns

#用いるデータ(互いに対応のあるデータであり、それぞれ正規性がないと仮定)
A = np.array([1.83, 1.50, 1.62, 2.48, 1.68, 1.88, 1.55, 3.06, 1.30])
B = np.array([0.88, 0.65, 0.60, 1.05, 1.06, 1.29, 1.06, 2.14, 1.29])

#試しにバイオリンプロットを作成
sns.violinplot(x=['A' for i in range(len(A))]+['B' for i in range(len(B))],
y=np.concatenate([A, B]))
result = stats.wilcoxon(A, B)
print('p-value: ', round(result.pvalue, 4))
# >>> p-value:  0.0039

p < 0.05なので帰無仮説は棄却され、2群間の代表値には差があると言える。

 


バイオリンプロット


検定の方法詳細

手順

  1. 仮説の設定
    (H0):「2群間の代表値に差がない」と仮定する。
    対立仮説(H1):「2群間の代表値に差がある」と仮定する。
  2. 有意水準 $\alpha$ を設定する。
  3. 統計量 T を算出する。
  4. 自由度 $\nu$ を算出する。
  5. サンプル数 $n$ に応じて判定を行う。
    $n \leqq 25$ : 統計量Tがウィルコクソンの符号順位検定表における棄却限界値以下の時に帰無仮説を棄却する。
    $n > 25$:統計量Zを計算する。統計量Zは近似的に標準正規分布に従うため、これが棄却限界値以下の時に帰無仮説を棄却する。

$$ Z = \frac{|T-\frac{n(n+1)}{4} |}{\sqrt{\frac{n(n+1)(2n+1)}{24} }} $$

統計量Tの求め方

  1. 2グループ間の差を$ d_i = x_{1i} - x_{2i} $ とおく。
  2. 同順位を持つ要素が存在する場合は、平均順位を割り当てる。
    (例:2つの値が違いに10位タイの時は,10.5(11と10の平均)を割り当てる)
  3. $ x_{1i} = x_{2i} $ (差が0)には順位を付けない。
  4. 絶対値 $ |d_i| $ の小さい順に順位をつける。
  5. $d_i$ が負のもの($W^-$)、正のもの($W^+$)の順位を加算し、小さい法の値をTとおく。

$$ T = min(W^+, W^-) $$

詳細コード

n ≦ 25のとき

import matplotlib.pyplot as plt

#用いるデータ(互いに対応のあるデータであり、それぞれ正規性がないと仮定)
A = np.array([9, 7, 8, 5, 7, 6, 4, 5])
B = np.array([6, 5, 8, 7, 6, 7, 7, 8])

#差を求めて絶対値の小さい順位並び替え
delta = [round(a-b,2) for a, b in zip(A,B)]
delta = sorted(delta, key=abs)

#差が0のものは順位からはずす
delta.remove(0)
delta = np.array(delta)
print(delta)
# >>> [ 1 -1  2 -2  3 -3 -3]

#順位を調べる絶対値が同じものは平均順位とする。
rank_list = np.zeros(len(delta))
i = 0
while i < len(delta):
number = abs(delta[i])

#絶対値が同じ値がないか調べる
same_value = len(np.where(abs(delta)==number)[0])

#絶対値が同じ値がある場合
if same_value >= 2:
    rank_list[i:i+same_value] = np.arange(i+1,i+1+same_value).sum()/same_value
    i += same_value
else:
    rank_list[i] = i+1
    i += 1
print(rank_list)
# >>> [1.5 1.5 3.5 3.5 6.  6.  6. ]

plus_list = rank_list[np.where(delta>0)[0]]
minus_list = rank_list[np.where(delta<0)[0]]
print('差がプラスの順位:',plus_list)
print('差がマイナスの順位:',minus_list)
"""output
差がプラスの順位: [1.5 3.5 6. ]
差がマイナスの順位: [1.5 3.5 6.  6. ]
"""

W_plus = plus_list.sum()
W_minus = minus_list.sum()
print('W+:',W_plus)
print('W-:',W_minus)
"""output
W+: 11.0
W-: 17.0
"""

#よって統計量Tは
T = min([W_plus, W_minus])
print('T:',T)
# >>> T: 11.0

#サンプル数(差が0のデータは外す)
n = len(delta)
print('n:',n)
# >>> n: 7
"""

n≦25なのでウィルコクソンの符号順位和検定の数表より、n=2の時の両側検定の棄却限界値2が求まる。
Tが2より大きいので、帰無仮説は棄却されず二つの代表値に差があるとは言えない

n > 25のとき

import matplotlib.pyplot as plt

#用いるデータ(違いに対応のないデータ、かつ不等分散であり、それぞれ正規分布に従うと仮定)
A = np.array([9, 7, 8, 5, 7, 6, 4, 5, 6, 7,
9, 9, 7, 9, 7, 6, 7, 8, 4, 8,
9, 6, 8, 7, 5, 6, 7, 9])
B = np.array([6, 5, 8, 7, 6, 7, 7, 8, 5, 8,
7, 8, 9, 8, 9, 9, 9, 7, 7, 6,
9, 8, 7, 8, 8, 9, 6, 8])

#差を求めて絶対値の小さい順位並び替え
delta = [round(a-b,2) for a, b in zip(A,B)]
delta = sorted(delta, key=abs)

#差が0のものは順位からはずす
delta.remove(0)
delta = np.array(delta)

#順位を調べる絶対値が同じものは平均順位とする。
rank_list = np.zeros(len(delta))
i = 0
while i < len(delta):
    number = abs(delta[i])

#絶対値が同じ値がないか調べる
same_value = len(np.where(abs(delta)==number)[0])

#絶対値が同じ値がある場合
if same_value >= 2:
    rank_list[i:i+same_value] = np.arange(i+1,i+1+same_value).sum()/same_value
    i += same_value
else:
    rank_list[i] = i+1
    i += 1

plus_list = rank_list[np.where(delta>0)[0]]
minus_list = rank_list[np.where(delta<0)[0]]
W_plus = plus_list.sum()
W_minus = minus_list.sum()
print('W+:',W_plus)
print('W-:',W_minus)
"""output
W+: 129.5
W-: 247.5
"""

#よって統計量Tとサンプル数nは
T = min([W_plus, W_minus])
n = len(delta)
print('T:',T)
print('n:',n)
"""output
T: 129.5
n: 27
"""

#n>25なので統計量を計算する
Z = abs(T - (n*(n+1)/4)) / (np.sqrt(n*(n+1)*(2*n+1)/24))
print('Z:',round(Z,3))
#>>> Z: 1.429

# 結果の判定と図示 =====
# x軸の等差数列を生成
X = np.linspace(start=-6, stop=6, num=100)

# 標準正規分布作成
norm = stats.norm.pdf(x=X, loc=0)

#パーセント点関数(ppf)でα=0.05の棄却域を計算する
upper_005 = stats.norm.ppf(q=0.95, loc=0)
lower_005 = stats.norm.ppf(q=0.05, loc=0)

#可視化
plt.plot(X, norm)
X_lower = np.linspace(start=-6, stop=lower_005, num=100)
X_upper = np.linspace(start=upper_005, stop=6, num=100)
lower = stats.norm.pdf(x=X_lower, loc=0)
upper = stats.norm.pdf(x=X_upper, loc=0)
base = np.zeros(100)
plt.fill_between(X_lower, base, lower, facecolor='lime', alpha=0.5)
plt.fill_between(X_upper, base, upper, facecolor='lime', alpha=0.5)

#求めたZをプロット
plt.plot(Z, 0.005, 's',label='Z value')
plt.title('norm distribution', fontsize=15)
plt.ylim(0,0.5)
plt.legend(fontsize=15)
plt.savefig('ウィルコクソンの符号順位和検定.png')
plt.show()
print('下限点:',round(lower_005,3))
print('上限点:',round(upper_005,3))
"""output
下限点: -1.645
上限点: 1.645
"""


Z値が有意水準0.05の時の棄却域(黄緑)に入っていないので、帰無仮説は棄却されず2群間の代表値には差があるとは言えない

参考