2007-09-18

ケースコントロール関連以外でパーミュテーションテスト

こちらのレビュー(Crit Care. 2004; 8(3): 196–199.)が簡潔にまとまっている。
Statistics review 10: Further nonparametric methods
ちなみにこのレビューシリーズは、全14回(継続中？)
このシリーズ１−１４回へのリンクはこちら。

データ
- SNPのジェノタイプは３つ
- それぞれのジェノタイプについて、フェノタイプが「あり　か　なし　か」ではなくて、「量」。
- 「量」というとき、いわゆるQuantitative traitのようにさまざまな値をとる場合と、順序量の場合とがある。また、さまざまな値をとる場合にも、その値の分布が正規分布やその他の、「きれいな」分布をとる場合と、まったくそうなっていない場合とがある。今、「あり　か　なし　か」で扱えないデータをひとまず一括して扱ってしまいたい、というのが要請であるから、分布を仮定しないで、順序量にも使ってしまう方法を採用する
- ちなみに、「正規分布を仮定」できるようなQuantitative traitについて、複数のカテゴリ(この場合はジェノタイプ)の間の違いを調べ、そのときに、ジェノタイプ3種類には序列を持ち込まないとしたら、それはANOVA。線形の関係を持ち込むなら、linear regression
- SNPのジェノタイプに序列を持ち込まずに、各ジェノタイプに観測される、「きれいでない」量的データ(順序量を含む)に対して行う検定は、Kruskal-Wallis。
- SNPのジェノタイプに序列を持ち込むと、Jonckheere-Terpstra。
ケース・コントロール分割表への(お試し)適用
- 今、あるアレルに着目し、その本数によって、３ジェノタイプの量的情報を、「０，１，２」とする
- 今、ある質的形質に着目し、そのありなし、を「０，１」とする。
- ２つの見方ができる。
  - 場合１
    - ケース群とコントロール群について、観測したところ、「０，１，２」という量が観測された。今、ケース群・コントロール群という群分けと「０，１，２」という量の分布とに関係があるかどうかを考える。
    - ケース群・コントロール群という２群なので、群に対して順序を問題にする必要はないが、拡張することを考慮して、順序があるものとする。
      - ケースのAA,Aa,aaがn1=30,n2=40,n3=20人。コントロールがm1=40,m2=30,m3=20人だとする。
        
        このデータは２カラム、１カラム目に90行、２カラム目に90行で、各レコードは「０，１，２」で記録される。
    - 「０，１，２」に順位をつけ、Kruskal-Wallis検定が可能。Rでの実行方法と結果は、次の通り。エクセルでは、D22,E22のカラム。自由度は１。
- -

n1<-30
n2<-40
n3<-20
m1<-40
m2<-30
m3<-20
case<-c(rep(2,n1),rep(1,n2),rep(0,n3))
cont<-c(rep(2,m1),rep(1,m2),rep(0,m3))
kruskal.test(list(case,cont))

- - - その結果は

        Kruskal-Wallis rank sum test

data:  list(case, cont) 
Kruskal-Wallis chi-squared = 1.1506, df = 1, p-value = 0.2834

- - - このときの、分割表検定の結果(エクセル)は

2x3khai square	2.8571 	0.239651036
Chi-trend	0.952380952	0.329113996
Cochran-Armitage trend	0.947089947	0.330461144
Kruskal-Wallis	1.150596878	0.283424463

- - 場合２
    - AA群、Aa群、aa群の３群について、「ケース＝１」「コントロール＝０」の観測がなされる。ジェノタイプ３群という群分けと形質の「０，１」という量の分布とに関係があるかどうかを考える。
    - ジェノタイプ３群なので、群に対して順序を問題にしないか、問題にするかの選択が必要としてもよい。。
      - ケースのAA,Aa,aaがn1=30,n2=40,n3=20人。コントロールがm1=40,m2=30,m3=20人だとする。
        
        このデータは３カラム、１カラム目に70行、２カラム目に70行、３カラム目に40行で、各レコードは「０，１」で記録される。
    - 「０，１」に順位をつけ、Kruskal-Wallis検定が可能。Rでの実行方法と結果は、次の通り。エクセルでは、D23,D23のカラム。自由度は２。

AA<-c(rep(n1,1),rep(m1,0))
Aa<-c(rep(n2,1),rep(m2,0))
aa<-c(rep(n3,1),rep(m3,0))

kruskal.test(list(AA,Aa,aa))

        Kruskal-Wallis rank sum test

data:  list(AA, Aa, aa) 
Kruskal-Wallis chi-squared = 2.8413, df = 2, p-value = 0.2416

- - - Kruskal-Wallisは当然のことながら、ジェノタイプ３タイプについて、順序を考慮していないので、２ｘ３分割表のカイ自乗検定(自由度２)と同じ挙動をする。
    - このときの、分割表検定の結果(エクセル)は

2x3khai square	2.8571 	0.239651036
χsq (1/2)=	1.1429 	0.285049646
χsq (11/12+22)=	2.3377 	0.126278968
χsq (11+12/22)=	0.0000 	1
Chi-trend	0.952380952	0.329113996
Cochran-Armitage trend	0.947089947	0.330461144
Kruskal-Wallis	2.841269841	0.241560596

- - ２ｘ３分割表の各種検定とともにKruskal-Wallisの枠組みでの統計量とそのｐ値を算出するエクセルはこちら。(アップロード予定)
  - Kruskal-Wallis検定を別途行うエクセルはこちら。

- より面倒くさいデータにあてはめられるように、制約を課した手法を単純なデータにあてはめているので、すこしずつ保守的に寄ってくる。

- Kruskal-Wallisの検算はRで。
  - RでのKruskal-Wallisはこちらを。