naoya_t@hatenablog

いわゆるチラシノウラであります

Metropolis-Hastingsアルゴリズム/PRML図11.9の再現

(このエントリはQiitaで12/1より開催中のMachine Learning Advent Calendar 2012の1日目の記事「Metropolis-Hastingsアルゴリズム」の転載です)

ご挨拶

今日から始まりました Machine Learning Advent Calendar 2012 幹事の @naoya_t です。

このアドベント・カレンダーの記事内容は、パターン認識・機械学習自然言語処理データマイニング等、データサイエンスに関する事でしたら何でもOKです。テーマに沿っていれば分量は問いません。(PRMLの読んだ箇所のまとめ、実装してみた、論文紹介、数式展開、etc.)

執筆する皆さんも読むだけの皆さんも共に楽しみましょう!

さて本題

いまPRML11章を読んでるので11章から何か。

明日のyag_aysさんがMCMCあたりと宣言されているので内容がもしかしたら被ってしまうかもしれないけれど、僕はOctaveだしyagさんは多分Rなので良しとしましょう!

Metropolisアルゴリズム

マルコフ連鎖モンテカルロ(MCMC)において、現在の状態から新しい状態をサンプリングする際、提案分布の対称性を仮定(q({\bf z}_A|{\bf z}_B)=q({\bf z}_B|{\bf z}_A))し、
\displaystyle A_k({\bf z}^*,{\bf z}^{(\tau)})=min\bigl(1,\frac{\tilde{p}({\bf z}^*)}{\tilde{p}({\bf z}^{(\tau)})}\bigr)
(11.33) を新しいサンプルを受理するか否かの基準として用いるアルゴリズム。

Metropolisアルゴリズムでは、サンプルが棄却された場合、現在の状態を新しい状態として引き続き用いる。(サンプルが重複することになる)

図11.9の再現

図11.9(PRML下巻p.253)では、Metropolisアルゴリズムを用いて2次元ガウス分布からサンプリングする例を図示しています。


一般に、等高線は contour() でも描けるのですが、図11.9のように標準偏差の線を1本だけ描きたい場合の使い方は分からなかったので、与えられたμとΣから楕円を1つ描く関数 plot_gauss2_contour() を用意しました。

Metropolis-Hastingsアルゴリズム

Metropolis-Hastingsアルゴリズム(§11.2.2)は、Metropolisアルゴリズムを提案分布が引数に対して対称な関数でない場合まで一般化したもの。

受理確率 A_k({\bf z}^*,{\bf z}^{(\tau)}) の計算が、
式 (11.33)
\displaystyle A_k({\bf z}^*,{\bf z}^{(\tau)})=min\bigl(1,\frac{\tilde{p}({\bf z}^*)}{\tilde{p}({\bf z}^{(\tau)})}\bigr)
から式 (11.44)
\displaystyle A_k({\bf z}^*,{\bf z}^{(\tau)})=min\bigl(1,\frac{\tilde{p}({\bf z}^*)q_k({\bf z}^{(\tau)}|{\bf z}^*)}{\tilde{p}({\bf z}^{(\tau)})q_k({\bf z}^*|{\bf z}^{(\tau)})}\bigr)
に差し替えられている。

Metropolis-Hastingアルゴリズムにおいて、提案分布の選択がアルゴリズムの性能に与える影響

等方ガウス分布を提案分布として用いる場合、標準偏差ρを(サンプリングしたい分布のσに対し)どの程度にしたら良いのかという話。(PRML下巻p.256, 図11.10参照)

1. ρ « σ_min
受理される遷移の割合は高い(9割前後)が、状態空間での前進はゆっくりとしたランダムウォークの形を取り、系列は長い時間にわたって相関を持つことになる。

2. ρ = σ_max
分散パラメータを大きく取れば、状態空間での前進は大きくなるが、同時に棄却率も高くなる。(この例では受理率7%前後)

3. ρ = σ_min
この位が良いらしい。(受理率7割前後)
ステップ数のオーダーは (\frac{\sigma_{max}}{\sigma_{min}})^2 になる。

4. (おまけ)10ステップ毎にサンプルを採用した場合の散布

終わりに

告知をいくつか。
PRML復々習レーン#7
12/15 14:00-20:00 @株式会社VOYAGE GROUP(渋谷・神泉)

PRMLカラオケクラスタひみつ集会Vol.1
12/15 20:30-22:30 @カラオケの鉄人 渋谷道玄坂

もう1つ。
先月CodeIQで簡単な機械学習(というか線形分類)の問題を出しました。
金貨の真贋を見分けよう
問題は12/20 10amまで公開していますので、ぜひ挑戦(あるいは秒殺)してみてください!

というわけで明日の担当は @yag_ays さんです。皆さん乞うご期待!