CAR に関する説明ですが,最初に読んだときは深澤君の説明がなんだかヘンな のでは……と思ったのですが,川渡でご教示いただいた GeoBUGS user manual 読んでみて,ようやくこのあたり理解できました.「ヘン」なのはむしろ car.normal() で使ってる intrinsic Gaussian CAR model のほうだったようで すね (便利なモデルだとは思いますが).
intrinsic Gaussian CAR model は隣接する場所間の local な共分散構造だけ きめて,調査地全体の global なばらつき (random effects 全体のばらつき) は決めない,というちょっとふうがわりなモデルだと今さらながら理解しまし た.この単純化のおかげで,ばらつきパラメーターは一個指定するだけで MCMC 計算できているわけです.
しかしながら空間構造を考えない GLMM なんかとの対応をきちんと説明するの は難しそうですね.私はこの点を誤解していて,普通の GLMM と car.normal() 使ったモデルが簡単に対応づけられる,と考えていました.
しかしながら実際は簡単ではないようで,となっているように思います.「対応がつかない」と言ってるのは, car.normal() の weight や tau をどうひねくりまわしても「場所差あり,た だし地点ごと独立同分布の GLMM」にはならない,ということです.
- 普通の場所差 GLMM: 調査地全体の random effects のばらつきは考慮している (あてはまりの良さと複雑さのバランスがとられる) が,どの二地点間も無相関であると仮定している (モデルで生成される確率場の variogram は水平線になる)
- car.normal() を使ったモデル: 近傍地点間の相関は推定されるが,調査地全体の場所差のばらつきをあらわすパラメーターはない (モデルで生成される確率場の variogram の最大値 (sill) が不定 (tau が定数であっても) …… のはず)
では,もう少し GLMM との対応がわかりやすいモデルはないかと考えると,お そらくこれが WinBUGS (GeoBUGS) では spatial.exp() を使ったモデルで phi が大きい状況 (距離とともに急速に無相関になる) というのが空間構造なしの GLMM に相当するかと思います.モデルで生成される確率場の variogram は原 点 (0, 0) をとおって上に凸の関数で最大値 (sill) に漸近する (これは spatial.exp() においては tau から決まるはず),つまり空間統計学の教科書 とかでみかけるもの (ただし nugget はゼロ) になります.
ベイズ企画集会においてこのあたりをごちゃごちゃ説明する必要はまったくな いと思いますが,今回つかっている空間統計モデルが一般的な CAR というより, (やや簡便版ということになるのでしょうか) intrinsic Gaussian CAR model ですよ,というのはスライドのはしっこにでも記載しておけばよいのかもしれ ません.
上のように考えてみると,このあたりは意外とめんどうなような気も してきました (他にも,たとえば空間相関のある random effects とは「移動 分散制約の効果なのか,環境不均質性の効果なのか……時間変化がわからない データではその区別がつかない」問題とか).
With Editor
状態のまま止まっている.
最初の revision のときは督促メイルだしたら対応したんだけど,
今回はまったく動く気配がないなぁ.
dpois()
がポアソン乱数 (整数)
を生成してないんでないか
(実数になる)
疑惑があって,
ですね.
model { Tau.noninformative <- 1.0E-2 for (i in 1:N) { Y[i] ~ dpois(lambda) y[i] ~ dpois(lambda) } log(lambda) <- beta beta ~ dnorm(0, Tau.noninformative) }この例では
Y[i]
は観測データなんだけど,
y[i]
は (データにしばられない) 単なるポアソン乱数になるはず
……
そうなっていると確認できた.
model { Tau.noninformative <- 1.0E-2 for (i in 1:N) { Y[i] ~ dpois(lambda[Group[i]]) } for (j in 1:N.group) { lambda[j] ~ dpois(mean) } log(mean) <- beta beta ~ dnorm(0.0, Tau.noninformative) }これは奇妙なモデルでポアソン分布の平均値である
lambda[j]
がポアソン分布にしたがう,
といった状況だ
……
で,この場合はやはりというべきか
> post.mcmc[, "lambda[1]"] [1] 14.65 16.96 15.10 17.87 15.35 15.08 15.56 16.87 13.78 12.98 12.89 13.72 [13] 14.21 14.55 13.61 15.93 16.36 15.92 14.20 14.46 16.02 15.52 13.64 15.63 [25] 14.14 13.24 16.72 11.70 14.97 17.59 15.24 15.64 15.82 13.35 15.93 15.97 [37] 16.46 17.37 16.46 13.53 14.96 14.85 15.29 15.38 14.56 13.79 15.47 15.34 [49] 15.07 14.41 13.14 15.10 17.02 14.03 14.25 13.95 13.06 14.41 14.50 14.29 ....と Gibbs sampling しているねえ. すなわち, WinBUGS 使い大学院生の (良くも悪くも) 独創的というほかない BUGS coding が WinBUGS 作者たちの想定した状況を超越してしまった, というコトのよーで (ひょっとしたらこれは WinBUGS の仕様どおりの動作, という可能性もゼロではないけど).
model { Tau.noninformative <- 1.0E-2 for (i in 1:N) { Y[i] ~ dpois(r.lambda[Group[i]]) } for (j in 1:N.group) { r.lambda[j] <- 1.0 * lambda[j] lambda[j] ~ dpois(mean) } log(mean) <- beta beta ~ dnorm(0.0, Tau.noninformative) }
dcat(p[])
を使ってみると
……
WinBUGS は「正しく」というべきか incompatibility
を指摘しているらしいナゾの trap 窓
(おそらく dpois()
のパラメーターに整数いれるな,
といった趣旨なんだろう)
をだして計算停止,
と.
なるほどね.