「ベイズ因子による検定」の編集履歴(バックアップ)一覧はこちら
「ベイズ因子による検定」(2007/07/13 (金) 22:21:58) の最新版変更点
追加された行は緑色になります。
削除された行は赤色になります。
=ベイズ因子とは=
*古典的な仮説検定の代わりとなるベイジアンな方法
*データベクトルxのもとで、二つの仮説M1とM2のいずれかのモデルを選択する際のベイズ因子Kは、
::<math>K=\frac{p(x|M_1)}{p(x|M_2)}</math>
:である。ここで<math>p(x|M_i)</math>はモデルiの周辺尤度と呼ばれる。
*尤度比検定と似て(最尤法とは異なり)、ベイジアンではこれをパラメータで平均を取る。一般にモデルM1、M2はパラメータベクトル<math>\theta_1,\theta_2</math>でパラメータ化され、Kは
::<math>K=\frac{p(x|M_1)}{p(x|M_2)}=\frac{\int p(\theta_1|M_1)p(x|\theta_1,M_1)d\theta_1}{\int p(\theta_2|M_2)p(x|\theta_2,M_2)d\theta_2}</math>
:と変形される。
*Kの対数は、xを与えたときのM1のM2に対するweight of evidenceと呼ばれることがある。
K > 1は、モデルM1がM2よりも強く支持されるということを示唆する。ハロルド・ジェフリーのKの解釈スケールはこうである
:{|
! K !! dB !! Strength of evidence
|-
! < 1:1
| <center> < 0 </center>
| Negative (supports M<sub>2</sub>)
|-
!1:1 to 3:1
| <center>0 to 5</center>
| Barely worth mentioning
|-
!3:1 to 10:1
| <center>5 to 10</center>
| Substantial
|-
!10:1 to 30:1
| <center> 10 to 15
| Strong
|-
!30:1 to 100:1
| <center>15 to 20</center>
| Very strong
|-
!>100:1
| <center>>20</center>
| Decisive
|}
=具体例、WTCCC論文のメソッドより=
*N個体からなるセット(ケースN<sub>1</sub>、コントロールN<sub>2</sub>)のSNPタイピングを行った。
*Y<sub>i</sub>は表現型で、個体iがケースのとき1、コントロールのとき2である。
*Z<sub>i</sub>は個体iにおけるアレル1の個数である。
{| border=1
|-
| || 0 || 1 || 2
|-
| Case || s<sub>0</sub> || s<sub>1</sub> || s<sub>2</sub>
|-
| Control || r<sub>0</sub> || r<sub>1</sub> || r<sub>2</sub>
|}
*モデルは3通り
**M<sub>0</sub> 関連なし
**M<sub>1</sub> additive effect on the log-odds scale
**M<sub>2</sub> general 3 parameter model of association
*M1とM0の間のベイズ因子は次のように定義する。
::<math>BF_1=\frac{P(D|M_1)}{P(D|M_0)}=\frac{\int P(D|\theta_1,M_1)P(\theta_1|M_1) d\theta_1}{\int P(D|\theta_0,M_0)P(\theta_0|M_0)d\theta_0}</math>
:ここでDは観察データであり、<math>\theta_1</math>と<math>\theta_0</math>はそれぞれモデルM<sub>1</sub>、M<sub>0</sub>のパラメータである。
*頻度論者のように尤度を最大化する代わりに、事前分布を与えてパラメータを重み付けし積分することができる。
*いずれのモデルにおいても尤度の[[ロジスティック回帰]]モデルを用いる。
::<math>P(D|\theta)=\prod_{i=1}^N p_i^{Y_i}(1-p_i)^{1-Y_i}</math>
:モデルM<sub>1</sub>においては
::<math>\theta_1=(\mu,\gamma) </math> <math>log \frac{p_i}{1-p_i}=\mu+\gamma Z_i</math>
:モデルM<sub>0</sub>では
::<math>\theta_0=(\mu)</math> <math>log\frac{p_i}{1-p_i}=\mu</math>
*事前分布<math>P(\theta_1|M_1)=P(\mu,\gamma|M_1)</math>をはっきりさせる必要がある。
**<math>\mu</math>はベースラインのオッズである。ケースとコントロールの数に影響されるが、症例対照研究では人為的にケースの数が大きくなっている。そのため[[確率分布#正規分布|正規分布]]<math>N(\alpha_1,\beta_1)</math>を<math>\mu</math>の事前分布として使用する。実際には<math>\mu \sim N(0,1)</math>とした。
*<math>\gamma</math>はリスクアレルの数に応じた対数オッズの増大を表しており、<math>e^{\gamma}</math>はオッズ比の加算効果を示す。よい事前情報がある。
**common diseaseの原因となるgenetic variantsはリスクアレルのオッズ比が1-2であり、特に1-1.5であろうと信じられている。そんなわけで事前分布<math>\gamma \sim N(\alpha_2,\beta_2)</math>をとる。
*まとめると、事前分布はこういう形になる。
::<math>P(\theta_1|M_1) \propt \frac{1}{\beta_1}e^{-\frac{(\mu-\alpha_1)^2}{2\beta^2_1}}\frac{1}{\beta_2}e^{-\frac{(\gamma-\alpha_2)^2}{2\beta_2^2}}</math>
*事前分布はベイズ因子に大きな影響を与える。<math>\mu</math>については両モデルに共通のdiffuseな分布を用い、<math>\gamma</math>に焦点を当てた比較を行うこととした。
*周辺尤度を評価する。
::<math>P(D|M_1)=\int P(D|\theta_1,M_1)P(\theta_1|M_1) d\theta_1</math>
:に対して[[公式#ラプラス近似|ラプラス近似]]を行う。
=ベイズ因子とは=
*古典的な仮説検定の代わりとなるベイジアンな方法
*データベクトルxのもとで、二つの仮説M1とM2のいずれかのモデルを選択する際のベイズ因子Kは、
::<math>K=\frac{p(x|M_1)}{p(x|M_2)}</math>
:である。ここで<math>p(x|M_i)</math>はモデルiの周辺尤度と呼ばれる。
*尤度比検定と似て(最尤法とは異なり)、ベイジアンではこれをパラメータで平均を取る。一般にモデルM1、M2はパラメータベクトル<math>\theta_1,\theta_2</math>でパラメータ化され、Kは
::<math>K=\frac{p(x|M_1)}{p(x|M_2)}=\frac{\int p(\theta_1|M_1)p(x|\theta_1,M_1)d\theta_1}{\int p(\theta_2|M_2)p(x|\theta_2,M_2)d\theta_2}</math>
:と変形される。
*Kの対数は、xを与えたときのM1のM2に対するweight of evidenceと呼ばれることがある。
K > 1は、モデルM1がM2よりも強く支持されるということを示唆する。ハロルド・ジェフリーのKの解釈スケールはこうである
:{|
! K !! dB !! Strength of evidence
|-
! < 1:1
| <center> < 0 </center>
| Negative (supports M<sub>2</sub>)
|-
!1:1 to 3:1
| <center>0 to 5</center>
| Barely worth mentioning
|-
!3:1 to 10:1
| <center>5 to 10</center>
| Substantial
|-
!10:1 to 30:1
| <center> 10 to 15
| Strong
|-
!30:1 to 100:1
| <center>15 to 20</center>
| Very strong
|-
!>100:1
| <center>>20</center>
| Decisive
|}
=具体例、WTCCC論文のメソッドより=
*N個体からなるセット(ケースN<sub>1</sub>、コントロールN<sub>2</sub>)のSNPタイピングを行った。
*Y<sub>i</sub>は表現型で、個体iがケースのとき1、コントロールのとき2である。
*Z<sub>i</sub>は個体iにおけるアレル1の個数である。
{| border=1
|-
| || 0 || 1 || 2
|-
| Case || s<sub>0</sub> || s<sub>1</sub> || s<sub>2</sub>
|-
| Control || r<sub>0</sub> || r<sub>1</sub> || r<sub>2</sub>
|}
*モデルは3通り
**M<sub>0</sub> 関連なし
**M<sub>1</sub> additive effect on the log-odds scale
**M<sub>2</sub> general 3 parameter model of association
*M1とM0の間のベイズ因子は次のように定義する。
::<math>BF_1=\frac{P(D|M_1)}{P(D|M_0)}=\frac{\int P(D|\theta_1,M_1)P(\theta_1|M_1) d\theta_1}{\int P(D|\theta_0,M_0)P(\theta_0|M_0)d\theta_0}</math>
:ここでDは観察データであり、<math>\theta_1</math>と<math>\theta_0</math>はそれぞれモデルM<sub>1</sub>、M<sub>0</sub>のパラメータである。
*頻度論者のように尤度を最大化する代わりに、事前分布を与えてパラメータを重み付けし積分することができる。
*いずれのモデルにおいても尤度の[[ロジスティック回帰]]モデルを用いる。
::<math>P(D|\theta)=\prod_{i=1}^N p_i^{Y_i}(1-p_i)^{1-Y_i}</math>
:モデルM<sub>1</sub>においては
::<math>\theta_1=(\mu,\gamma) </math> <math>log \frac{p_i}{1-p_i}=\mu+\gamma Z_i</math>
:モデルM<sub>0</sub>では
::<math>\theta_0=(\mu)</math> <math>log\frac{p_i}{1-p_i}=\mu</math>
*事前分布<math>P(\theta_1|M_1)=P(\mu,\gamma|M_1)</math>をはっきりさせる必要がある。
**<math>\mu</math>はベースラインのオッズである。ケースとコントロールの数に影響されるが、症例対照研究では人為的にケースの数が大きくなっている。そのため[[確率分布#正規分布|正規分布]]<math>N(\alpha_1,\beta_1)</math>を<math>\mu</math>の事前分布として使用する。実際には<math>\mu \sim N(0,1)</math>とした。
*<math>\gamma</math>はリスクアレルの数に応じた対数オッズの増大を表しており、<math>e^{\gamma}</math>はオッズ比の加算効果を示す。よい事前情報がある。
**common diseaseの原因となるgenetic variantsはリスクアレルのオッズ比が1-2であり、特に1-1.5であろうと信じられている。そんなわけで事前分布<math>\gamma \sim N(\alpha_2,\beta_2)</math>をとる。
*まとめると、事前分布はこういう形になる。
::<math>P(\theta_1|M_1) \propt \frac{1}{\beta_1}e^{-\frac{(\mu-\alpha_1)^2}{2\beta^2_1}}\frac{1}{\beta_2}e^{-\frac{(\gamma-\alpha_2)^2}{2\beta_2^2}}</math>
*事前分布はベイズ因子に大きな影響を与える。<math>\mu</math>については両モデルに共通のdiffuseな分布を用い、<math>\gamma</math>に焦点を当てた比較を行うこととした。
*周辺尤度を評価する。
::<math>P(D|M_1)=\int P(D|\theta_1,M_1)P(\theta_1|M_1) d\theta_1</math>
:に対して[[公式#ラプラス近似|ラプラス近似]]を行う。
::<math>\log P(D|M_1) \approx \log P(D|\hat{\theta}_1,M_1)+log P(\hat{\theta}_1|M_1)+\frac{d}{2} \log(2\pi)-\frac{1}{2}\log |A|</math>
:ここで<math>\hat{\theta}_1</math>は<math>P(D|\theta_1,M_1)P(\theta_1|M_1)</math>を最大化するような<math>\theta_1</math>で、<i>maximum a posteriori</i>(MAP)推定量として知られる。
*<math>A</math>は<math>\hat{\theta}_1</math>は<math>P(D|\theta_1,M_1)P(\theta_1|M_1)</math>のnegative Hessianで、dは<math>\theta_1</math>のdimensionである。
*<math>\hat{\theta}_1</math>を求めるためNewton-Raphson法を用いたが、それで収束しないときにはline-search法を用いた。