Lecture notes The lecture notes will be handed out at the lectures. As photocopying the notes is much cheaper than printing them, please do not print out the notes.
12月くらいからMCMCの勉強しだして、いくつか代表的なアルゴリズムによるサンプリングをやったのでまとめておく。 Example of Rejection Sampling - yasuhisa's blog Example of importance sampling - yasuhisa's blog Example of Metropolis Hastings Algorithm - yasuhisa's blog Metropolis Hastings Algorithmの続き - yasuhisa's blog Gibbs Sampler Algorithmによって多変量正規分布からのサンプル抽出を行なう - yasuhisa's blog あとはモデルによって色々変わるけど、根幹となるアルゴリズムはできたからまあよいか。 これで一応自分で作れるという感じにはなったので、MCMC
もちろんRで。 2変数正規分布 2変数の正規分布からのサンプリングはここのをそのまままねすると B <- 2000 x1 <- rep(NA,B+1) x2 <- rep(NA,B+1) x1[1] <- -2 x2[1] <- 1 for (i in 1:B){ x1[i+1] <- rnorm(1,1+0.7*(x2[i]-2),sqrt(1-0.7^2)) x2[i+1] <- rnorm(1,2+0.7*(x1[i+1]-1),sqrt(1-0.7^2)) } 相関係数もこんな感じになって、 > cor(x1,x2) [1] 0.712666 散布図もこんな感じ。うまくいってますね。 面倒なのでbuild-inとか考えてない(ry。 3変数正規分布 Gibbs Sampler Algorithmを使っているので、full conditioningした分布が分かってないといけないです
今さら感はあるが、MCMC (Markov Chain Monte Carlo; マルコフ連鎖モンテカルロ)を使えるようになろうと、まずは簡単な例から試してみた。 手始めに、正規乱数から生成した標本の平均と標準偏差を推定してみる。 やはりRを使用。MCMCpackパッケージを あらかじめインストールしておいて、呼び出す。MCMCpack中のMCMCmetrop1R()関数を利用して、メトロポリス法によるMCMC推定をおこなう。 library(MCMCpack) 乱数系列を初期化。 set.seed(1) 平均10、標準偏差3の乱数を1000個生成して、xに入れる。 m <- 10 s <- 3 x <- rnorm(1000, m, s) MCMC推定に使用する関数を用意する。betaは要素数2のベクトル。beta[1]が平均、beta[2]が標準偏差で、betaを推定する。関数の返り値
リリース、障害情報などのサービスのお知らせ
最新の人気エントリーの配信
処理を実行中です
j次のブックマーク
k前のブックマーク
lあとで読む
eコメント一覧を開く
oページを開く