カイ二乗検定(χ²検定)

統計的仮説検定を行う際に、カイ二乗(χ²)値とχ²分布を用いて2変数間の独立・適合度を検定する方法です。

クロス分割表の関連性の検定、比率の差の検定、分布の適合性の検定などを行う際に使用します。クロス分割表に入る期待値が5以下のものがある場合にはχ²検定を用いるのは不適当になるため、他の方法を適用します。

【χ²値算出法】
O:観測値
E:期待値
χ² = Σ{(O-E)²/E}

〔例〕
血液型と疾患Aの間の関連性検定

	疾患A	健常者	合計
A型	16	36	52
B型	12	16	28
O型	25	26	51
AB型	11	4	15
合計	64	82	146

帰無仮説 H₀:「疾患と血液型は独立である(関連がない)」
対立仮説 H₁:「疾患と血液型は独立ではない(関連がある)」

統計的仮説検定の基本は「ある仮定のもとで起こりにくいことが起きたときは、その仮定を棄てる」なので、否定したい仮説をH₀としてこの仮説下で「起こりにくいことが起きる」ことを示します。

分割表の各セルの値をNij(i=列、j=行)、期待値を Eij (i=列、j=行)で表すと周辺度数より

E_ij = (N_{i total} × N_{total j}) / N_{total total}

の式にて各セルの期待値が計算され、

χ²値 = Σ{(O_ij – E_ij)²/_Eij}

の式にてχ²値が計算されます。
分割表の自由度は F = (行のセル数-1)×(列のセル数-1)
にて計算されます。

上の例題では、χ²値≒9.47489となります。
このχ²値は自由度が3のχ²分布に従います。

自由度3のχ²分布において有意確率0.05を取る値 = 7.8147です。 これより χ²値(=9.4749) > 7.8147 なので、「起こりにくいことが起きた」 として、帰無仮説を棄却します。 上記の場合の正確な有意確率:p-value = 0.0236004になります。

χ2値 < 7.8147 のとき、帰無仮説を採択することになるので
「2変数は独立でないとはいえない(関連があるとはいえない)」

χ2値 > 7.8147 のとき、帰無仮説を棄却することになるので
「2変数は独立ではない(関連がある)」

よってこの場合、「血液型と疾患の間に有意な関連がある」という結論に至ります。 有意確率は慣例的に0.05(95%の確率)もしくは0.01(99%の確率)を使用します。

有意確率 (p-value)

統計的仮説検定において、帰無仮説を棄却し対立仮説を採択する時の基準です。有意水準とも言います。一般的に5%、1%が使用されます。

例えば、検定した結果が有意確率5%以下であるならば、最初に立てた帰無仮説を「滅多に起こらないこと」として棄却し対立仮説を採択します。

また、「有意確率αで検定すると有意な差が認められた」ということと、「危険率αのもとで有意な差があるといえる」は同じような意味で使用されます。

有意確率5% で検定を行うということは、第1種の過誤をおかす危険率が5%であることを意味しており、同様の調査・検定を行うと、20 回に 1 回は得られた結論が誤っていることを表します。

連鎖不平衡係数

Dについて

【算出法】
遺伝子座Aのアレル頻度 = P_A P_a
遺伝子座Bのアレル頻度 = P_B P_b
各ハプロタイプ頻度 = P_AB P_Ab P_aB P_ab

D = P_AB-P_AP_BP_AB = P_AP_B+D
D = P_aP_B-P_aBP_aB = P_aP_B-D
D = P_AP_b-P_AbP_Ab = P_AP_b-D
D = P_ab-P_aP_bP_ab = P_aP_b+D
D = P_ABP_ab-P_AbP_aB
(値の範囲 = -0.25 ～ 0.25)

ある遺伝子Aと遺伝子Bの間に関連が無い(連鎖平衡である)場合、アレルAとアレルBは独立であるため、ハプロタイプ頻度はAとBのアレル頻度を掛け合わせたものとなります(積の法則)。

しかしAとBの間に関連がある(連鎖不平衡が存在する)場合、ハプロタイプ頻度は各アレル頻度の積の値からずれます。このズレを連鎖不平衡の尺度として定義しているのが連鎖不平衡係数:Dです。

Dは、D’、r²を導く際に必要となる重要な値です。

D’について

【算出法】
Dmin = max (-P_AP_B, -P_aP_b) = (-1)×min(P_AP_B, P_aP_b)
Dmax = min (P_AP_b, P_aP_B)
D’ = D / D_max(D>0の場合)
D’ = D / D_min(D<0の場合)
(値の範囲 = -1 ～ 1 (絶対値を取る場合は0～1))

すべてのハプロタイプ頻度とアレル頻度は確率なので、0から1の間の値を取ります。故に連鎖不平衡尺度Dは限られた範囲の値を取る値です。

Dがとりうる範囲を限定されていることより、連鎖不平衡尺度を規格化した値がD’です。D>0の場合にはDの最大値で、D<0の場合にはDの最小値で尺度Dを割ることで算出します。なお、D’の定義は2種類あり、どちらのアレルを規準にするかを考えて正負の符号付きで表示する場合と、絶対値をとって0-1の範囲に収める定義があります。SNPAlyzeでは、前者で計算しておりますが、絶対値をとる|D’|も表示してあります。

D’は主にマーカーと遺伝子座との遺伝的距離に依存することから、尺度として適切な性質を持つとも言われます。D’= 1の場合には2つの遺伝子が組換えを起こしておらず、完全に連鎖不平衡にある状態を指します。
しかしD’= 1の場合の問題点として、以下のような場合の時、

ハプロタイプ数

		遺伝子座1
		c	d
遺伝子座2	a	50	0
	b	49	1

アレルdの遺伝子多型は偶然サンプルに入ってきたとみなせるので、この結果を連鎖不平衡と考えるのは不自然です。

しかしながら、この場合はD’=1となってしまいます。片方あるいは両方の遺伝子多型頻度が極めて低い場合には、高い値が算出される場合があるので、D’を使うときは注意が必要となります。このようにD’はサンプルサイズに依存する性質を持ちます。

また、ある1遺伝子について考えるとすると、D’はその遺伝子の中で組換えが起こったかどうかの判断に用いることが出来ます。しかし原因変異を追及しようとすると、組換えが起こっていないとほとんどのSNPでD’の値は高くなってしまいます。そのため、遺伝子内のどのSNPに原因があるかを推定する場合にはr²などの、他の尺度も考慮する必要があります。

D’は組み換えの歴史を現す尺度です。

r²について

【算出法】
r² = D² / (P_AP_BP_aP_b) ( r²=Δ²と表される場合もあり)
(値の範囲 = 0 ～ 1)

r²はD’よりもシビアな値をとるため、1遺伝子内でも連鎖不平衡にあるSNPの組合せを検出しやすい尺度です。そのためD’を補完する性質を持ちます。

また、r²はD’に比べてサンプル数が少ない場合の値上昇率は低いです。r²が1/3以上の値を取るならば、十分に強い連鎖不平衡があるとみなすことが出来るとのことです。

r²は各遺伝子座の2つの対立遺伝子を成分とする2×2表を対象としたχ²検定を行い、χ²/n とした値です。よってこの値を用いて連鎖不平衡を評価すると、χ²検定を用いた場合と似た結果になります。

r²は変異の成り立ちを表す尺度です。

EMアルゴリズム

不完全データからの最尤推定アルゴリズムです。最尤法の直接応用が困難な場合、完全データが得られたという仮想的状況での尤度を用いて反復的手順にて不完全データに対する最尤推定量を計算します。

EMアルゴリズムは2段階のステップから成り、E-step(Expectation step)にて、理論量を計算し、M-step(Maximization step)にてその理論量を最大化します。これらを交互に繰り返してパラメータを更新します。これにより、最尤推定あるいは尤度関数の極大点を得ることが出来ます。

〔例〕
血液型データからアレル頻度を推定

A型	B型	AB型	O型
160	80	40	120

このような血液型データが得られている場合、
表現型 = A、B、AB、O の4つであるのに対して
遺伝子型 = a/a、a/o、b/b、b/o、a/b、o/o の6パターンが存在します。

表現型の実測データから遺伝子型を知ることは不可能ですが、EMアルゴリズムを用いて遺伝子型を推定することが出来ます。

アレル頻度をそれぞれ P(a)、P(b)、P(o)
集団がハーディー・ワインベルグ平衡にあると仮定します。
EMアルゴリズムでは最初に推定する値を代入する必要があるので、
P(a)=0.2、P(b)=0.2、P(o)=0.6と仮定します(不完全データ)。

—1回目計算開始——————————————————–
E-Step
  アレルaの期待値:
    E(a) = 160×{P(a)(P(a)+P(O)) / P(a)(P(a)+2P(o))}+80×0 + 40×{P(a)P(b) / 2P(a)P(b)}+120×0 = 111
  アレルbの期待値:
    E(a) = 160×0+80×{P(b)(P(b)+P(O)) / P(b)(P(b)+2P(o))} + 40×{P(a)P(b) / 2P(a)P(b)}+120×0 = 66
  アレルoの期待値:
    E(a) = 160×{P(a)+P(O) / P(a)(P(a)+2P(o))}+80×{P(b)+P(O) / P(b)(P(b)+2P(o))} + 40×0 +120 = 223
  〔ここでは全アレル数に占める各アレル数を算出しています〕

M-Step
  アレルaの頻度を再計算:
    P(a) = 111/ (111+66+223) = 0.2775
  アレルbの頻度を再計算:
    P(b) = 66/ (111+66+223) = 0.165
  アレルoの頻度を再計算:
    P(o) = 223/ (111+66+223) = 0.5575
  〔ここではE-Stepの結果を元にアレル頻度の再計算を行います〕
—1回目計算終了——————————————————–

—2回目計算開始——————————————————–
E-Step
    E(a) = 116 E(b) = 65 E(o) = 219

M-Step
    P(a) = 0.29 P(b) = 0.1625 P(o) = 0.5475
—2回目計算終了——————————————————–

—3回目計算開始——————————————————–
E-Step
    E(a) = 118 E(b) = 65 E(o) = 217

M-Step
    P(a) = 0.295 P(b) = 0.1625 P(o) = 0.5425
—3回目計算終了——————————————————–
・・・
このようにE-StepとM-Stepを繰り返し、アレル頻度を更新し続けると、頻度はある値まで収束します。
この例ではP(a) = 0.2923 P(b) = 0.1630 P(o) = 0.5447が収束値となり、表現型から遺伝子型頻度が推定されます。

ブートストラップ法

ブートストラップ法とは、ある標本集団から母集団の性質を推定するための方法です。

復元抽出法により標本集団から標本集団と同じ数だけランダムに値を再抽出し、新しいデータセットを取得し統計値を計算します。

それを例えば1000回繰り返したら、1000個分のデータの平均値と標準偏差、信頼区間を算出することが出来ます。イメージとして以下のようになります。

BootStrap法は、推定値の信頼性評価を目的としています。再抽出データから得られたデータのばらつきを評価することでオリジナルデータから得られる統計値の性質を導く手法です。

信頼区間

母集団の平均・分散値があると推定される区間のことです。

95%信頼区間、99%信頼区間などを用いますがそれは、母集団の値がこの区間にある確率が95%、99%であることを表しています。

いいかえれば「100回サンプリングしたら、95回(99回)はこの範囲内に値があてはまる(という確率)」といえます。信頼度95%の信頼区間は次の式にて計算します。

〔例〕
正規分布の場合

m:標本集団の平均値

s:標本集団の標準偏差

n:標本集団の大きさ

m-1.96 (s / √n)<= X <= m+1.96 (s / √n)

この値が95%信頼区間となります。図で示すと以下の通りです。

両側検定を行う場合は、分布の両側、2.5%区間と97.5%区間を除いた部分が、100回サンプリングを行って、95回分のデータが収まる区間です。

AIC

赤池情報量基準(Akaike’s information criterion)の略。

AIC = -2×(最大対数尤度)+2×(説明変数の数)

で定義されます。

最大対数尤度とは、あるモデルに含まれるパラメータにおける真の値の最尤推定値対数です。説明変数はモデルに含まれる未知のパラメータです。

この式の右辺の第2項は推定するパラメータの個数が多くなるとその分不安定になるので、真のモデルから遠くなるという“悪さ”が評価されています。

つまり、AICの値が小さければ、それだけ良いモデルということになります。

クロス分割表の場合のAIC算出方法

AICはより小さい方がモデルとして適性である、と定義されています。

ここでは2×2分割表を例にして、独立モデル(IM)と従属モデル(DM)の求め方を説明します。まず、fig.1のような2×2分割表を定義します。

	B	B’	計
A	a	b	h
A’	c	d	n-h
計	k	n-k	n

独立モデル(IM)

以下の条件にて、fig.2のような分割表のモデルを作成します。

• A群とB群の間に関連が存在しない(事象A群と事象B群は独立)
• Aである確率= p、A’である確率= 1-p、Bである確率= q、B’である確率 =1-q
• AとA’、BとB’は排反事象 (p(A)+p(A’)=1、p(B)+p(B’)=1)

モデルを独立と仮定すると、A∩Bの起こる確率はそれぞれの起こる確率を掛け合わせたものになります。
よって以下のようなモデルを考えることが出来、モデルのパラメータ数= 2になります。

	B	B’	計
A	pq	p(1-q)	p
A’	(1-p)q	(1-p)(1-q)	(1-p)
計	q	1-q	1

この独立モデルの確率関数をPi(a,b,c,d)とすると、それぞれのマス目の要素同士が排反事象であるので、このモデルを多項分布に当てはめて以下の確率関数を算出出来ます。

Pi(a,b,c,d) = (n!/a!b!c!d!)×{pq^a×p(1-q)^b×(1-p)q^c×(1-p)(1-q)^d}
= (n!/a!b!c!d!)×{p^hqk(1-p)^n-h(1-q)^n-k}

pとqに対する対数尤度を求めると、対数尤度を取る

LPi(a,b,c,d) = h log p+k log q+(n-h)log(1-p)+(n-k)log(1-q)

p,qについて偏微分して、最尤推定すると

p=h/n q=k/n

これをLPi(a,b,c,d)に代入して最大対数尤度を求めると

MLPi(a,b,c,d) = h log h+k log k+(n-h)log(n-h)+(n-k)log(n-k)-2nlogn
AIC(IM) = -2×MLPi(a,b,c,d) + 2×2

により計算します。

従属モデル

以下の条件にて、fig.3のような分割表のモデルを作成します。

• A群とB群の間に関連が存在する(事象A群と事象B群は独立ではない)
• Aである確率= p、A’である確率= 1-p、Bである確率= q、B’である確率 =1-q
• AとA’、BとB’は排反事象 (p(A)+p(A’)=1、p(B)+p(B’)=1)

モデルを従属と仮定すると、A∩Bの起こる確率はp×qではなくなります。従って、p(A∩B) = p11、p(A∩B’) = p12、p(A’∩B) = p21、p(A∩B’) = p22、 とそれぞれの確率を定義します。この時、p11+p12+p21+p22 = 1という条件が存在するので、モデルのパラメータ数= 3になります。

	B	B’	計
A	p11	p12	p
A’	p21	p22	(1-p)
計	q	1-q	1

この従属モデルの確率関数をPd (a,b,c,d)とすると、それぞれのマス目同士は排反事象であることより、このモデルを多項分布に当てはめて、以下の確率関数を算出出来ます。

Pd(a,b,c,d) = (n!/a!b!c!d!)×(p11^ap12^bp21^cp22^d)

p11、p12、p21、p22に対する対数尤度を求めると、対数を取る

LLd = a log p11+b log p12+c log p21+d log p22

偏微分して最尤推定すると、

p11 = a/n p12 = b/n p21 = c/n p22 = d/n

これをLPd (a,b,c,d)に代入して最大対数尤度を求めると

MLLd = a log a+b log b+c log c+d log d-n log n
AIC(DM) = -2×MLLd + 2×3

により計算します。

AICはより小さい方がモデルとして適性である、と定義されているので、IM-DM > 0の場合にこの分割表には従属モデルを採用する、という判断が出来ます。より適正なモデルが導かれるための基準量を計算していることになります。用いる数値が分割表の実測値のみであるのも利点の一つです。

〔例〕
新薬Aと従来薬Bの薬効比較

	治癒	効果無し	計
新薬A	353	166	519
従来薬B	304	104	408
計	657	270	927

〔独立モデル:AIC(IM)〕
MLPi = 519log519 + 657log657 + 408log408 + 270log270 – 2・927log927 = -1195.124
AIC(IM) = (-2)× (-1195.124) + 2×2 = 2394.25

〔従属モデル:AIC(DM)〕
MLPd = 353log353 + 166log166 + 304log304 + 104log104 – 927log927 = -1192.775
AIC(DM) = (-2)×(-1192.775) + 2×3 = 2391.55

AIC(IM) – AIC(DM) = 2394.25 – 2391.55 = 2.7 {AIC(IM) > AIC(DM)}

従って、このデータには従属モデルを当てはめるべきであることが分かり、新薬Aは従来薬Bよりも多少は治癒効果が高いと判断出来ます。

フィッシャーの正確確率

2×2分割表において、いずれかのセルの期待値が5以下である場合、χ²検定を適用するのは不適当です。

その場合、フィシャーの正確確率を用いて分割表の直接生起確率を求めて有意確率(p-value)を算出する方法を使用します。例えばこのような2×2分割表があったとして、

	B	B’	計
A	13(a)	4(b)	17
A’	6(c)	14(d)	20
計	19	18	37

この分割表の存在確率は₁₇C₁₃×₂₀C₆/₃₇C₁₉ = 0.00522 となります。

周辺度数を固定した時、2×2分割表の自由度は1となるため、A∩Bの数を決定すれば他のセルに入る値も自動的に決定されます。フィッシャーの正確確率を求めるためにはA∩Bが0～17までの値を取った場合の分割表の存在確率を全て求めます。

また、2要因の関連の強さの指標として ad-bc を定義します。すると、以下のような結果が得られます。

ad-bcの値(もしくは絶対値)が測定値A∩B = 13の場合のad-bcより大きい場合を「極端な生起確率」とします。この値を合計したものが、有意確率となります。

片側検定の場合は(ad-bc > 158)、両側検定の場合は(|ad-bc| > 158)の生起確率を合計します。合計した確率が < 0.05 の場合、帰無仮説を棄却することが出来ます。

片側検定・両側検定

つまり、調べたい検定の帰無仮説・対立仮説の立て方によってどこに片側検定にするか両側検定にするかが決まります。 実際、片側検定と両側検定ではどう違うかというと、両側検定におけるP値0.05は、片側0.025になりますが、片側検定では片側のみで0.05なので、 当然片側検定の方が有意差が出やすくなります。通常は両側検定をしておけば間違いないと言えます。

〔例〕
フィッシャーの正確確率検定の場合

(1) 新しく開発された新薬Aが従来薬Bよりも有効であることを示す場合
この場合、「新薬 > 従来薬」という結果があるか否か、という検定であるので、片側検定を行います。

	治癒	効果無し	計
新薬A	13	4	17
従来薬B	6	14	20
計	19	18	37

p-value = 0.0058548

で新薬Aは従来薬Bよりも効果が高い、という結果になります。(仮想例)

(2) 男女間で読書時間に違いがあるかどうかを示す場合
この場合、違いがある場合の時間差は不明なので両側検定を行います。

	6時間以上	6時間未満	計
男性	13	4	17
女性	6	14	20
計	19	18	37

p-value = 0.0081381

で男性の方が女性よりも睡眠時間が長い、という結果になります。(仮想例)

このように片側検定を使うか、両側検定を使うかは、検定を行う条件により決まります。
データが同じでも片側か両側かによって計算されるp-valueが異なります。

オッズ比

オッズとは、ある状態の形成に寄与する因子が2つあった時に、その因子の寄与がどれほど大きいかを表す数値です。 オッズはある事象が発生する確率をpとしたときp / (1-p)で表されます。

	Case	Control
SNP(+)	a	b
SNP(-)	c	d
	a+c	b+d

上の分割表のCase群について考えると

「Case群におけるSNPのオッズ」=(a/(a+c))/(c/(a+c))=a/c

Control群について考えると

「Control群におけるSNPのオッズ」=(b/(b+d))/(d/(b+d))=b/d

さらにオッズ比という考え方では、Case群に対するControl群のSNPオッズを考えます。

「オッズ比」=Case群におけるSNPのオッズ/Control群におけるSNPのオッズ
=(a/c)/(b/d)=ad/bc

SNPが無い場合(SNP-)に比べ、ある場合(SNP+)ではad/bc倍 Case状態を形成しやすいとなります。

〔例〕
SNPがある場合の罹患オッズ比

	Case	Control	計
SNP(+)	13	4	17
SNP(-)	6	14	20
計	19	18	37

上の分割表では、

Case群の形成にSNP+が寄与するオッズ = 2.17、SNP‐が寄与するオッズ = 0.46

となります。つまり、SNPがある場合、Caseになる可能性はSNPが無い場合に比べて2.17倍高いということになります。

さらに上の例では、

オッズ比 = (13×14) / (4×6) = 7.58

つまり、SNPが無い場合(SNP-)に比べ、ある場合(SNP+)では 7.58倍 Case状態を形成しやすいということが分かります。

ディプロタイプ型について

ハプロタイプの組合せをディプロタイプと呼びます。表現型に寄与する遺伝子発現にはディプロタイプの型が大きく関係しています。

例えば、ある薬剤Xに対する薬剤耐性が遺伝子座Aによって支配されていると仮定します。遺伝子座Aの対立遺伝子Bは薬剤Xに対する耐性が高く、bは低いとします。この時、

遺伝子座Aの遺伝子型

B/B→耐性高
B/b→耐性中
b/b→耐性低

であることが予想されます。遺伝子座をハプロタイプに置き換えても同じことが言えるため、ディプロタイプ型を考えることはSNP解析にとって重要となります。

最尤推定法(尤度)

「尤度」:
ある事象Aがすでに起きている状態(実測データが得られている)からモデルを考えます。
モデルの仮説を立て、その仮説が正しいとした条件下で先ほど得られた実測値が得られる確率を尤度といいます。

「確率」:
ある事象Aがまだ起こっていない状態(実測データが得られていない)において、事象Aのモデルから確率を考えます。

〔例〕
事象A:コインを3回投げて2回表が出る (コインの表が出る確率 = p)

「確率」
コインに歪みなどが無いとして、

p = 0.5 (これは経験より考え出された可能性)

これより

₃C₂×(1/2)²×(1/2) = 0.375

確率を求める場合、予め表の出る確率が判明しているので、全体モデルの確率をそのまま計算出来ます。

---

「尤度」
コインを3回投げて2回表が出るというデータが得られました。しかしコインの表が出る確率 = pは不明とします。まずモデル全体を考えます。

₃C₂×p×p×(1-p) = L (likelihood)・・・(1)

尤度 = Lはこのモデルの確率です。
では3回中2回が表である、というデータの場合に最も尤度が大きくなるのは?

₃C₂×p×p×(1-p) = L
L = 3p²-3p³

両辺を微分して

dL = 6p-9p²
dL = 0 とおいて(最大値を求めるため)
p(6-9p) = 0
p = 6/9 = 0.667

つまり、3回中2回が表のデータ下においては、
表の出る確率 p = 0.667の時にモデル(1)は最も尤もらしい値をとります。
実際にp = 0.5, p = 0.6, p = 0.667, p = 0.7を(1)に当てはめて値を算出すると

p = 0.5 → 0.375(確率で考えた時の結果なので当たり前)
p = 0.6 → 0.432
p = 0.667 → 0.444
p = 0.7 → 0.384

というようにモデル(1)の確率を最も大きくするのは表の出る確率 p = 0.667の場合と分かります。
3回中2回表が出ているのだから、2/3 = 0.667 という直感的な計算結果と一致します。
これが最尤推定法です。

Permutation test

Permutation testは、2群の差を検定する際によく用いられる手法です。

2群それぞれのデータ数を固定したまま、データをランダムに入れ替えます。その結果から統計値を算出し、オリジナルの統計値と比較することで、2つの集団に差があるかどうかを検定します。2群に差があるのならば、オリジナル統計値よりも高い値出る可能性が低いはずです。

計測した分布がどの程度有意であるかの判断に用いることも出来ます。ブートストラップ法との違いは、データを入れ替えるだけの非復元抽出を行う点と、2群間の差に適用する点です。下のようなイメージになります。

〔例〕
2群の平均値の差

グループA(SNP+)	60kg	54kg	75kg	48kg	55kg
グループB(SNP-)	49kg	76kg	65kg	58kg	62kg

この2群に属する人の体重の平均値に差はあるかを調べます。

平均値は

平均値

グループA(SNP+)	60kg 54kg 75kg 48kg 55kg	→平均値:58.4kg
グループB(SNP-)	49kg 76kg 65kg 58kg 62kg	→平均値:62.0kg

平均値の差 = 3.6kg

この平均値の差が、AとBの母集団に違いがあるために得られた差であるのか偶然によるものなのかを判断するために、permutation testを行います。

まず、グループAとグループBの母集団に違いが無いと仮定すると、先ほどの平均値の差は、偶然得られた差であるということになります。

ここで、A=5、B=5のデータ数は固定して、ランダムにデータの並べ替えを行います。

並べ替え1回目

グループA(SNP+)	75kg 49kg 58kg 48kg 76kg	→平均値:61.2kg
グループB(SNP-)	55kg 62kg 65kg 60kg 54kg	→平均値:59.2kg

平均値の差 = 2.0kg

並べ替え2回目

グループA(SNP+)	58kg 65kg 49kg 76kg 55kg	→平均値:60.6kg
グループB(SNP-)	62kg 75kg 48kg 60kg 54kg	→平均値:59.8kg

平均値の差 = 0.8kg

・
・
・

この組合せは₁₀C₅ = 252通り存在します。

観測値である3.6kgという値が偶然生じた差であるならば、252通りの組合せの中に3.6kgよりも大きな値が頻繁に現れるはずです。

252通り全ての並べ替えを行い、観測値3.6kgよりも大きい値を取ったのは9通りであったとします。すると、偶然に3.6kg以上の差が生じてしまう可能性は9/252 = 0.0357 = 3.57%となります。

有意確率5%よりも小さな確率なので、グループAとグループBの母集団に違いがある可能性がある、と判断出来ます。

この0.0357をpermutation p-valueと言います。