データ科学2.pptx

Size: px

Start display at page:

Download "データ科学2.pptx"

ひろとちづ
2 years ago
Views:

1 データ科学多重検定 2 mul%ple test False Discovery Rate 藤博幸

2 前回の復習 1 多くの検定を繰り返す時には単純に個々の検定を繰り返すだけでは不十分 5% 有意水準ということは, 1000 回検定を繰り返すと, 50 回くらいは帰無仮説が正しいのに間違って棄却されてすまうじちがあるということ ex) 1 万個の遺伝子について正常細胞とガン細胞でそれぞれの遺伝子の発現に差があるかどうかを検定

3 RNA-seq で発現量を 5 回計測正常細胞ガン細胞遺伝子遺伝子正常細胞ガン細胞 1.5 遺伝子遺伝子正常細胞ガン細胞. 1.5 遺伝子遺伝子正常細胞遺伝子ガン細胞遺伝子遺伝子遺伝子正常細胞ガン細胞遺伝子遺伝子遺伝子遺伝子遺伝子

4 前回の復習 2 Family Wise Error Rate で制御する 1) Bonferroni 法 2) Holm 法

5 FWER: FP/(FP + TN) を制御 TN (True Nega*ve) 帰無仮説が正しく棄却されない FP (False Posi%ve) 帰無仮説が正しいのに棄却される

6 FWER は厳しすぎる FP ( 帰無仮説が正しいのみあやまって帰無仮説を棄却すること ; 正常細胞とガン細胞で遺伝子 A の発現に差がないのに棄却してしまうこと ) を抑えようとしているがそのため TP ( 帰無仮説が間違っていて帰無仮説が正しく棄却されること ; 正常細胞とガン細胞で遺伝子 A の発現に差があり帰無仮説が棄却されること ) であっても棄却されにくい ----à 発現に差のある遺伝子を検出しにくい ( お金と時間がかかっているのに )

7 FDR (False Discovery Rate) FP が混じっても良いので TP を増やすただし FP がどれだけ含めれているか (FDR) を推定してそれを新たな基準とする TP をどの程度許可するかの基準を有意水準の p- 値と区別するため q- 値とよぶ p- 値の分布に仮定を設けて FDR を計算

8 FDR = FP / (FP + TP) TP = true posi%ve: 対立仮説が正しく帰無仮説が棄却 FP = false posi%ve: 帰無仮説が正しいのに帰無仮説が棄却

9 Benjamini-Hochberg 法 1 m 回の多重検定の場合 (1) p-value を昇順に並べる (p 1 < p 2 < < p m ) (2) i = m とする (3) p i <= α i / m を満たす時 p 1 ~ p i を有意とするそうでなければ i を i 1 にして, 上の条件を確認する i = 1 になっても条件を満たさない場合は有意なものはないとする FWER の Holm 法が p 1 から順次 p-value が大きい方に向かうのに対し BH 法は p m から出発して, p-value が小さい方に向かうことに注意

10 Benjamini-Hochberg 法 2 p 1 < p 2 < < p i < p m とする i 番目 (p i ) を検討する α を false posi%ve を含む割合とする. m α は m 回の検定中 false posi%ve の期待回数 i / m は, m 回の検定中 i 個の占める割合すると, m α i / m = α i は i 回の検定の中で false posi%ve を起こす回数 (= 棄却されたのに本当は帰無仮説が正しいもの ) の期待値 FP の上限として α i を考える

11 Benjamini-Hochberg 法 3 今 p i 以下の p-value で棄却するすると posi%ve となる検定が i 個あるこの i 個には true posi%ve も false posi%ve も含まれているとすると FP + TP = i FDR = FP/( FP + TP) = FP/i < α i / i = α このようにして FDR を α 以下になるよう制御できる

12 Benjamini-Hochberg 法 4 p i <= α i / m を満たす時 p 1 ~ p i を有意とするこの式を書き直すと p i m / i <= α を満たす時 p 1 ~ p i を有意とする q i = p i m / i を p-value に対して q-value とよぶ

13 Benjamini-Hochberg 法 5 実際に行ってみよう p-value Gene Gene Gene Gene Gene p- 値のソート p-value Gene Gene Gene Gene Gene h_ps:// より

14 Benjamini-Hochberg 法 6 p-value q-value q- 値への変換 Gene Gene Gene Gene Gene x (5/1) x (5/2) x (5/3) x (5/4) x 5/5 0.21

15 Benjamini-Hochberg 法 7 FDR の閾値を 0.05 として p-value の大きいものから順番に検討 p-value q-value Gene Gene Gene Gene Gene x (5/1) x (5/2) x (5/3) x (5/4) x 5/ < 0.05 < 0.05 > 0.05 > 0.05 > 0.05 Gene 2 と Gene 5 が有意

16 R でやってみよう 1 Gene1~Gene5 の p-value をベクトルとして表現 > pv <- c(0.21, 0.001, 0.1, 0.06, 0.005) pv を昇順にソート > spv <- sort(pv) > spv [1] q-value を記憶させる空ベクトルを作成 > qv <- c() for 文で q-value に変換 > for (i in 1:length(spv)) qv <- c(qv, spv[i]*length(spv)/i) > qv [1] 前ページの表と同じ結果が得られたことを確認

17 R でやってみよう 2 p.adjust 関数を使ってみる > p.adjust(pv, method="bh", n=length(pv)) [1] > p.adjust(pv, method="bh", n=length(pv)) < 0.05 [1] FALSE TRUE FALSE FALSE TRUE 遺伝子の並び順のままで検定を行える

18 Benjamini-Hochberg 法の問題点 1 Benjamini-Hochberg 法のキモは q-value の計算で q i = p i m / i を p-value に対して q-value とよぶ p i を基準とした時その時の false posi%ve の数が p i m と仮定している点にある p i が一様分布していればこの仮定は正しい帰無仮説が全て正しい場合 (ex 個全ての遺伝子で正常細胞とガン細胞における発現量に差はない ) p-value は一様分布する

19 帰無仮説が全て正しい場合における p-value の一様分布の確認前回の多重検定での false posi%ve 生成と同じシミュレーションを実施 > N < # サンプル数を > p.values <- double(n) # p-value を格納するベクトル > for(i in 1:N) { + x <- rnorm(10, mean=0, sd=1) + y <- rnorm(10, mean=0, sd=1) + p.values[i] <- t.test(x,y)$p.value + } # 同じ正規分布から 10 こずつサンプルして平均値の差の検定 # 同じ分布からのサンプルなので差はないはず > hist(p.values) # p-value のヒストグラム作成

20 p-value がほぼ一様分布に従っていることがわかる

21 Benjamini-Hochberg 法の問題点 2 h_ps://sites.google.com/site/scriptocioinforma%cs/maikuroarei-guan-xi/fdr-zhi-yu-r

22 Benjamini-Hochberg 法の問題点 3 実際には差があるものが含まれるので p-value は 0 に近い方に偏った分布となる > N < > rm1 <- runif(n)*5 > rm2 <- runif(n)*5 > rs1 <- runif(n)*2 > rs2 <- runif(n)*2 > p.values <- double(n*2) > # 5000 個は同じ分布からサンプリング > for (i in 1:N) { + x <- rnorm(10, mean=0, sd=1) + y <- rnorm(10, mean=0, sd=1) + p.values[i] <- t.test(x,y)$p.value + } > # 残り 5000 個はランダムに生成した平均と標準偏差を持つ分布からサンプリング > for (i in 1:N) { + x <- rnorm(10, mean=rm1[i], sd=rs1[i]) + y <- rnorm(10, mean=rm2[i], sd=rs2[i]) + p.values[n+i] <- t.test(x,y)$p.value + } > hist(p.values)

23 p-value の一様性は成立していないこの点については Q-value 法や Local FDR 法などの改良があるが今回は説明を省く

24 帰無仮説の分布と対立仮説の分布の混合分布帰無仮説と対立仮説の密度比が各々 π 0 :1 - π 0 の比だったとする h_ps://sites.google.com/site/scriptocioinforma%cs/maikuroarei-guan-xi/fdr-zhi-yu-r h_ps:// の記述だともう少し複雑

25 ヒストグラムのどの場所で高さをとるかという問題を考えたのが Q-value Q-value では高さをとる場所を γ γ から推測される π 0 の値を π 0 (γ) としてプロットしこのデータを自然スプライン関数で回帰このスプライン関数を γ=0 まで外推した値を π 0 の予測値とする p-value のアナロジーとして使われている小文字の q-value とは違ってこれは Q value という手法名である点に注意 Q value はこの π 0 を BH 法の q-value にかけるそれ以外の計算手順は BH 法と同じなので先にこの π 0 の予測値をかけておいた p -value を BH 法に適用すれば Q-vakue を求める事ができる R のパッケージ qvalue をインストールすると Q-value を計算できる

26 h_p://strimmerlab.org/notes/fdr.html

27 h_ps:// h_ps:// h_ps://sites.google.com/site/scriptocioinforma%cs/maikuroarei-guan-xi/fdr-zhi-yu-r h_p://d.hatena.ne.jp/hoxo_m/ /p1 h_p://

EBNと疫学

EBNと疫学推定と検定 57 ( 復習 ) 記述統計と推測統計統計解析は大きく 2 つに分けられる記述統計推測統計記述統計観察集団の特性を示すもの代表値 ( 平均値や中央値 ) やばらつきの指標 ( 標準偏差など ) 図表を効果的に使う推測統計観察集団のデータから母集団の特性を推定する平均 / 分散 / 係数値などの推定 ( 点推定 ) 点推定値のばらつきを調べる ( 区間推定 ) 検定統計量を用いた検定