読了: Plummer, et al.(2006) RのCODAパッケージ

Plummer, M., Best, N., Cowles, K., Vines, K. (2006) CODA: Convergence Diagnosis and Output Analysis for MCMC. R News, 6(1), 7-11.
 読んだのはなんでも記録しておこうということで… MCMC診断用のRパッケージの古株、codaの解説記事。

 もともとcodaパッケージは昔のBAGSチームによってS Plus用に書かれたものだった。ちなみにboaパッケージもcodaから派生したもの。
 まずオブジェクトの作成から。mcmc(), mcmc.list()を用意している。後者は並列チェーンの出力を表す。read.coda(), read.jags(), read.openbags()とかでも作れる。
 mcmcクラスにはplotメソッドがあってトレース・プロットと密度プロットを描く。ほかにautocorr.plot()でacfが描ける。crosscorr.plot()というのもある。
 summaryメソッドで要約統計量が出せる。SEは時点ゼロのスペクトル密度を推定するspectrum0()で出せるんだけど、自己相関が高すぎるとクラッシュするので、代わりにspectrum0.ar()を使うように。effectiveSize()というものもある。
 診断手法としては、geweke.diag(), gelman.diag(), heidel.diag(), raftery.diag()がある。云々。

 … 以上は2006年段階の話なので、最新版(0.19-4)のヘルプをみて補完すると…
 プロット系ではcumuplot()ってのができてますね。横軸がイテレーションで縦軸が四分位点。
 要約系ではHPDinterval()というのができてる。HPD(事後密度が高い区間)を出せる。
 診断系では、geweke.plot(), gelman.plot()というのもある。gelman.diag()はmcmc.listしか食わないみたいだ(つまり単一チェーンのR-hatは依然として出せない)。

 多変量ESSはどうやらcodaじゃ無理みたいで、mcmcseパッケージを使うしかなさそうだ。
 それにしても、Rにはmcmc標本を扱うパッケージが乱立していて不安になる。いまざっと探したところでは、古株のcodaとboa, Stan公式が出しているbayesplot, みんな大好きなshinystanの他に、lattice系で描画するplotMCMC, なんだか知らんがMCMCvis, そしてtidyverse系のtidybayes(たぶんデータハンドリングのみ), ggplot2系のggmcmc, というのがあるらしい。いろいろあるなあ。