【MCMC法シリーズその2】ハミルトニアンモンテカルロ法と物理学の関連性について

こんにちは。スキルアップAI編集部です。

機械学習の手法の中には、物理学の知識と非常に密接に関わっているものが多くあります。 今回はマルコフ連鎖モンテカルロ (Markov Chain Monte Carlo; MCMC) 法の1つである、ハミルトニアンモンテカルロ (Hamiltonian Monte Carlo; HMC) 法について触れていきます。MCMC法は、パラメータの事後分布を推定する方法として、機械学習の様々な場面で利用されます。 HMC法は、物理学で発展した概念をMCMC法に応用した重要な手法です。
本記事では、HMC法と物理学の関連性について、厳密性は無視して直感的な説明を試みます。

1.マルコフ連鎖モンテカルロ (MCMC) 法とは

MCMC法は『マルコフ連鎖』と『モンテカルロ法』の2つの用語を組み合わせた単語です。
『マルコフ連鎖』とは、今の状態は1個前の状態だけによって決まるという考え方です。例えば、天気の統計モデルを考える際、ある日の天気を決めるのはその前日の天気だけと割り切り、2日前の天気の影響は無視します。これは統計モデルを考えやすくするために重要な考え方です。

『モンテカルロ法』とは、乱数を振って数値的にパラメータの値を推定する方法です。有名な例として、円周率の近似値をモンテカルロ法で推定する方法が挙げられます。図1に、モンテカルロ法で円周率を求める方法のイメージを示します。モンテカルロ法を用いて円周率を求めるには、まず、大量の赤点をランダムに発生させます。その後、青色領域と白色領域に入った赤点の数をそれぞれ数え、その比を計算します。その比を用いると面積が求まるので、その面積を用いて円周率を算出します。

図1. モンテカルロ法で円周率を求める方法
(参考文献[1]より引用)

以上の2語を併せたのがMCMC法です。MCMC法では大まかに以下の手順でパラメータ θ の事後分布を推定します。

MCMC法の手順の大まかな流れ
1. パラメータを θ とし、その初期値 θ(0) を定めます。
2. ある規則に従って、パラメータの候補値 θ(a) を得ます。
3. 現在のパラメータを更新するか (θ(t) ＝ θ(a)) 、もしくは更新しないかを確率的に決定します。
4. 一定回数 2 と 3 の操作を繰り返します。

この手順のうち、2 と 3 の違いによって様々なMCMC法が提唱されています。
今回注目しているHMC法は、この 2 と 3 の操作にハミルトニアンという物理学の道具を応用した手法です。

2.ハミルトニアンとは

ハミルトニアンは、解析力学や量子力学といった物理学の諸分野でよく登場する概念です。厳密性を無視して言えば、ハミルトニアンは、ある物体 (質点) の力学的エネルギーと言えます。もう少し詳細に言うと、ハミルトニアン H(p,h) は、運動量 p を変数に持つ運動エネルギー K(p) と位置 (高さ) h を変数に持つ位置エネルギー (ポテンシャルエネルギー) U(h) を足した物理量です。

【式1】

図2. 摩擦を無視した坂を運動する質点の概念図
※黒い丸は質点を表す

単純な例を図2に示します。ここでは、摩擦が無視できる坂を転がる質点を想定します。最初頂上から初速度0で運動を始めると、次第に速さを増していき、最下点に到達すると運動エネルギーが最大になります。しかし、どの地点を見ても運動エネルギーと位置エネルギーの和、すなわち力学的エネルギーは一定のままです。

【式2】

ただし、 m は質点の質量、 v は速度、 g は重力加速度、そして h は物体の位置 (高さ) を意味しています。

この系で言えば、この力学的エネルギーを運動量 p = mv と位置 (高さ) h で表現したものがハミルトニアンと言えます。運動量は質量と速度をかけた値として定義されます。

この運動量 p を用いて上記の式を変形すると、

【式3】

となります。これがハミルトニアンです。

証明は省きますが、ハミルトニアンが一定という条件から以下のハミルトンの運動方程式を導出することができます。

【式4】

これは重要な方程式であり、HMC法では、この式に従ってパラメータの更新を行います。

また、位相空間という概念も重要です。位相空間とは、運動量 p と位置 h の軸を持つ空間です。位相空間には、ハミルトニアンの等高線を描くことができます。

図2. 位相空間の概念図
※横軸が運動量、縦軸が位置 (高さ)を表し、
曲線はハミルトニアンが一定の場所を表す

3.ハミルトニアンモンテカルロ(HMC)法について

【式5】

標準正規分布とは平均0、分散1の特別な正規分布です。この分布について、以下の変形を施します。

ここで、以下のような置き換えをします。

【式6】

すると、

【式7】

となります。【式5】とおくのは、対数確率の高さを物理空間の低さとみなすとことに相当します。

このようにハミルトニアンを用いて確率分布を表せるという結果が導けます。ただし、物理学の文脈とHMC法の文脈では変数の意味が変わっています。物理学では、変数は運動量 p と位置 (高さ) h だったのに対し、HMC法では乱数 p とサンプリングしたい量 θ となっています。

あとは、【式6】における p と θ で構成されるハミルトニアンについて、ハミルトンの運動方程式【式3】を用いて p と θ を更新します。その更新には、計算精度を考慮してリープフロッグ法などがよく用いられますが、今回は割愛します。この部分の詳細は文献[2]などが詳しいです。

4.HMC法と物理学の関連性

物理学でのハミルトニアンとHMC法でのハミルトニアンには、以下の対応関係を見ることができます。

・物理学での運動量 p とHMC法での乱数 p
・物理学での位置(高さ) h とHMC法でのパラメータ θ
・物理学での位置エネルギーとHMC法でのパラメータ θ の事後分布

これより、図3のように、HMC法におけるサンプリングの流れを物理学的なイメージに繋げることができます。

図3. 物理学とHMC 法の位相空間上でのアナロジー

(左) 物理学の観点。横軸を運動量、縦軸を位置 (高さ) h とみる。(1) 質点に運動量 p を加えると別のハミルトニアンに移る。 (2)質点はハミルトニアン一定曲線上を運動する。
(右) HMC法の観点。横軸を乱数 p の値で縦軸をサンプリングするパラメータ θ とみる。(1)パラメータ θ の現在点はある確率に基づいて現在のパラメータ θ を別のハミルトニアンへ移る。(2)パラメータ θ の現在点は、ハミルトニアン一定の線上を移動する。この軌跡を記録していく。この軌跡がパラメータ θ のサンプリング結果になる。

最後に、HMC法の大まかな流れをまとめます。

HMC法の計算手順の大まかな流れ
1. パラメータを θ とし、その初期値 θ(0) を定めます。
2. 正規分布から乱数 p をサンプルします。
3. ハミルトンの運動方程式に従って、 θ(a)、 p(a) を計算します。 計算にはリープフロッグ法などを利用します。
4. 確率 min(1,r) でパラメータを更新 (θ(t) = θ(a)) するか、もしくは更新しないか (θ(t) = θ(t-1)) を決定します。 ただし、 r = exp(H(θ(t),p(t)) – H(θ(a),p(a))) と定義されます。
5. 一定回数だけ 2 から 4 の操作を繰り返します。

5.おわりに

本記事では、HMC法と物理学の関連性について紹介しました。機械学習の手法の中には、物理学と関連したものがHMC法以外にもあります。ご興味のある方はぜひ調べてみてください。

6.参考文献

[1] 学びTimes : 高校数学の美しい物語 モンテカルロ法と円周率の近似計算， https://manabitimes.jp/math/1182 (2021年)．
[2] 豊田秀樹: 基礎からのベイズ統計学　ハミルトニアンモンテカルロ法による実践的入門， 朝倉書店 (2015年)．

☆☆☆
スキルアップAIのメールマガジンでは会社のお知らせや講座に関するお得な情報を配信しています。
配信を希望される方はこちら

また、SNSでも様々なコンテンツをお届けしています。興味を持った方は是非チェックしてください♪
Twitterはこちら
Facebookはこちら
LinkedInはこちら
スキルアップAI公式YouTube AIビジネスチャンネルはこちら