So-net無料ブログ作成
  • ブログをはじめる
  • ログイン

R: NIMBLEでWAICを計算してみた

きのうにひきつづき、NIMBLEをつかってみるテストです。WAICを計算させてみました。

ZIPモデルのWAICを計算してみた」で使ったデータを使います。

以下がコードです。

はじめのモデルでは潜在変数zをつかってコーディングしています。25行目で、zが既知の部分には値(1)をあたえています。

結果です。

> zip1.out$summary
$chain1
            Mean    Median    St.Dev. 95%CI_low 95%CI_upp
lambda 2.4278981 2.4050451 0.41310146 1.6811911 3.3384034
omega  0.4784749 0.4726144 0.08675762 0.3246998 0.6461032

$chain2
            Mean    Median   St.Dev. 95%CI_low 95%CI_upp
lambda 2.4633608 2.4342948 0.4136240 1.7018781 3.3253671
omega  0.4754656 0.4760111 0.0857376 0.3188948 0.6431468

$chain3
            Mean    Median    St.Dev. 95%CI_low 95%CI_upp
lambda 2.4569056 2.4495721 0.41610300 1.7153744 3.2865996
omega  0.4754655 0.4713965 0.08927651 0.3116644 0.6607087

$all.chains
            Mean    Median    St.Dev. 95%CI_low 95%CI_upp
lambda 2.4493882 2.4282182 0.41442723 1.6949683 3.3245902
omega  0.4764687 0.4733892 0.08725235 0.3168532 0.6559968


> zip1.out$WAIC
[1] 180.881

推定値はあっているのですが、計算されたWAICの値がかなりおおきくなっています。データ中のzの一部がNAのせいかもしれませんが、原因は不明です。また、データにzをあたえないと、WAICは計算されません。

つぎのモデルでは、ZIPの分布を定義しています。NIMBLEの例題のコードを参考にしました。この定義では、ゼロになる確率を定義しているので、BUGSコードではY[i] ~ dZIP(lambda, 1 - omega)としています。

結果です。

> zip2.out$summary
$chain1
            Mean    Median    St.Dev. 95%CI_low 95%CI_upp
lambda 2.4506730 2.4037902 0.42448499 1.7255186 3.3503215
omega  0.4763018 0.4759791 0.08501409 0.3177539 0.6516389

$chain2
            Mean    Median    St.Dev. 95%CI_low 95%CI_upp
lambda 2.4499844 2.4263678 0.40543207 1.7266570 3.2979968
omega  0.4716324 0.4675331 0.08836584 0.3077814 0.6524603

$chain3
            Mean    Median    St.Dev. 95%CI_low 95%CI_upp
lambda 2.4465500 2.4206025 0.41498684 1.6929659 3.2913297
omega  0.4777543 0.4771137 0.08818478 0.3178834 0.6529767

$all.chains
            Mean    Median    St.Dev. 95%CI_low 95%CI_upp
lambda 2.4490691 2.4202375 0.41490637 1.7116558 3.3098488
omega  0.4752295 0.4734133 0.08721187 0.3153975 0.6522954


> zip2.out$WAIC
[1] 109.9841

WAICの値も、Stanとlooで計算させたものと近くなっています。それでもややちがうのですが、その理由はまだ追求していません。


nice!(2)  コメント(0) 
共通テーマ:日記・雑感

nice! 2

コメント 0

コメントを書く

お名前:
URL:
コメント:
画像認証:
下の画像に表示されている文字を入力してください。

Facebook コメント