ハミモン、ゲットだぜ！
Hamilton Monte-Carlo法(Hybrid Monte-Carlo法)はMCMCによる分布関数のサンプリングを高速化させる手法の一つであり、近年StanやTheanoなど統計的解析を行うためのプログラミング言語に実装させており、注目を集めています。
今回は混合ガウス分布に対するその分布の軌跡をplotしました。
HMCの基本
基本的には通常のモンテカルロ法で扱う変数の組みxに対して同じ数の運動量変数pを用意し、ハミルトニアンHを
\( H=\frac{ p^2}{2}+V(x) \)
\( V(x) = -\log(p(x)) + C \)
と定義します。
- 運動量の変数(p)のランダムな更新
- 一定期間分ハミルトン方程式
\( \frac{dp}{dt}=-\frac{ \partial H}{\partial x} \)
\( \frac{dx}{dt}=\frac{ \partial H}{\partial p} \)
を積分して解く(leapflog法による更新の繰り返し)
- Metropolis Hasting法で更新されたxの値を採択するか否かを決定する。
を繰り返すという手順です。
p(x)が小さい(V(x)が大きい)領域では運動量に相当するpの値が大きくなり、大きな幅で遷移が行われ、効率的な分布が作成できると期待されます。
ただし積分がある分計算量的には1遷移あたりの計算量が通常のMCMCより大きくなる点に注意が必要です。
code
計算結果
codingとplot方法に関しては@teramonagiさんの[レプリカ交換モンテカルロ法（パラレル・テンパリング）による混合ガウス分布に従う乱数の生成 - My Life as a Mock Quant
を参考にさせていただきました。
通常のMCMC
遷移のヒストグラム acceptされなかった遷移は0に集中している。
Hamilton Monte-Carlo
遷移のヒストグラム 通常のMCMCよりも広い範囲に広がっている(ように見える)。
Reference
ほぼ同じことをかなり以前に紹介されています。参考文献の論文も有用です。
- TheanoのHMC sample codeへの日本語注釈
http://nbviewer.ipython.org/gist/xiangze/c2719235434bee796288
TheanoではTheano.T.grad()を用いてハミルトニアンの偏微分を数式として計算することが可能です。これによって確率分布関数の代入と高速な計算を両立させることが出来ます(ただし式のコンパイルに時間がかかります)。
主にTheanoの関数の説明になってしまっています。サンプルコードでは並列化が行われていないためにGPUの利点が発揮できてはいないことに注意が必要です。
- 作者: C.M.ビショップ,元田浩,栗田多喜夫,樋口知之,松本裕治,村田昇
- 出版社/メーカー: 丸善出版
- 発売日: 2012/02/29
- メディア: 単行本
- 購入: 6人 クリック: 14回
- この商品を含むブログを見る
PRML読み会での発表資料だそうです。
http://www.slideshare.net/wk77/bishop-prml-115116wk771006061152
