Hazy Ideas

日々の勉強の気づきを書き出しています

R: GLMのベイズモデル化

*個人勉強用メモです。

 

MCMCによるパラメータ推定とベイズモデリングを、一般化線形モデルに適応する例を考える。

 

簡単な例として、ランダム項なしのモデル 平均 λ = exp(β1 + β2X)を考える。

尤度関数 L(β1, β2) = Π (y | λ) = Π (y | β1, β2, X)

ベイズの事後分布は尤度と事前分布だったので、観測データYが与えられた時の

p(β1, β2 |Y) ∝ p(Y | β1, β2) p(β1) p(β2)

 

β1, β2の事前分布を考える。

無情報事前分布:無情報のように見える分布、よく平均0、標準偏差が広い分布が使われる。

乱数によってYを生み出して、メトロポリス法で実装する。

 

MCMCサンプリングの数について。

初期の捨てるステップをburninで指定する。

サンプリングの反復数も指定する。

収束診断には、収束の指数であるR指数を用いる。

収束しない場合は、モデリングが不適切、コード間違い、データ間違いなどが考えられる。

 

得られた周辺事後分布が収束しているかの判断の例。

事後分布の95%信用区間で推定したい場合、両端2.5%に独立なMCMCサンプルが入るようにサンプリング回数を決める。(ちょっと理解不足)

 

ギブスサンプリング:効率と汎用性が高いサンプリング方法

新しい値の確率分布を作り、その確率分布のランダムサンプルを次の新しい値とする。

最初にβ1, β2の初期値を設定し、

p(β1 | Y, β2)に従う乱数を発生させて、β1とする

p(β2 | Y, β1)に従う乱数を発生させて、β2とする

これを繰り返す。

性質がよくわかっている確率分布のときは、それに適した乱数発生方法がある。