🎯

【統計学】メトロポリス-ヘイスティング法(Metropolis-Hastings Algorithm)で詳細釣り合い条件が成立することの証明

2023/08/13に公開

はじめに

"統計科学のフロンティア12 計算統計II"を読んでいたところ, メトロポリス-ヘイスティング法における詳細釣り合い条件が成り立つことの確認が書いてあった. その証明はメトロポリス法と同様にして出来るとのことであった. 本記事ではその確認を実際に計算し確認する.

議論の背景

規約性と非周期性を持つ状態と遷移確率関数, 定常分布の存在を仮定すると, 下記の定理が成り立つことが知られている. しかし, この定理が一般的なため, これを満たすマルコフ連鎖を構成することが難しい. そこで詳細釣り合い条件を制約に追加することでこれを満たすマルコフ連鎖を構成する. 以下補足.

定義

状態\mathbf{x}がベクトルであることに注意しておく.

\mathbf{x}から\mathbf{x}'への遷移確率関数を以下のように定義する.

\begin{align} \pi(\mathbf{x} \to \mathbf{x'}) \end{align}

定常分布を以下のように定義する.

\begin{align} &P(\mathbf{x}), \text{Where}\\ &\forall \mathbf{x}, P(\mathbf{x'}) = \sum_{\mathbf{x} \in \mathcal{X}} P(\mathbf{x}) \pi(\mathbf{x} \to \mathbf{x'}) \end{align}

定理

状態と遷移確率関数について規約性と非周期性が成り立っているとする. すなわち, 以下が成り立つとする.

\begin{align} &あるMが存在して, 任意の\mathbf{x}, \mathbf{x}' に対して, \\ &\pi による遷移を丁度M回行うと, 有限の確率で\mathbf{x}から\mathbf{x}'へ到達できる. \end{align}

定常分布P(\mathbf{x})が存在すると仮定する.

このとき, 以下が成り立つ.

\begin{align} m \rightarrow \inftyでP_m(\mathbf{x})は初期状態\mathbf{x}_0に依存せず, P(\mathbf{x})に収束する. \end{align}

証明は参考文献を参照すること.

詳細釣り合い条件(detailed balance)

詳細釣り合い条件とは以下のような条件である.

\begin{align} \pi(\mathbf{x} \to \mathbf{x'}) P(\mathbf{x}) = \pi(\mathbf{x'} \to \mathbf{x}) P(\mathbf{x'}) \end{align}

制約がより強くなっていることに注意する. すなわち, 上記を満たせば前述の定理から定常分布への収束が保証される.

メトロポリス法(Metropolis Algorithm)とメトロポリス-ヘイスティング法(Metropolis-Hasting Algorithm)

定義

新たに提案分布Q(\mathbf{x} \to \mathbf{x'})を導入する. この分布は\mathbf{x}に対して\mathbf{x}'を遷移先の候補とする確率である

メトロポリス法のアルゴリズム

メトロポリス法では以下の制約をおく.

\begin{align} Q(\mathbf{x} , \mathbf{x'}) = Q(\mathbf{x'} , \mathbf{x}) \end{align}
  1. 候補先\mathbf{x}'を提案分布Q(\mathbf{x}, \mathbf{x'})からサンプリングする.
  2. \mathbf{x}'を受理する確率を以下のように定義する.
\begin{align} r = \min \left( 1, \frac{P(\mathbf{x'})}{P(\mathbf{x})}\right) \end{align}
  1. rに従って\mathbf{x}'を受理するか否かを決定する.
  2. \mathbf{x}'を受理した場合は\mathbf{x}'を新たな状態として採用する. そうでない場合は\mathbf{x}を新たな状態として採用する.

メトロポリス-ヘイスティング法のアルゴリズム

メトロポリス法では提案分布に仮定をおいたが, メトロポリス法はこの仮定を外したものである. 詳細釣り合い条件を満たすために, 以下のようにアルゴリズムを変更する.

  1. 候補先\mathbf{x}'を提案分布Q(\mathbf{x}, \mathbf{x'})からサンプリングする.
  2. \mathbf{x}'を受理する確率を以下のように定義する.
\begin{align} r = \min \left( 1, \frac{P(\mathbf{x'}) Q(\mathbf{x'} , \mathbf{x})}{P(\mathbf{x}) Q(\mathbf{x} ,\mathbf{x'})} \right) \end{align}
  1. rに従って\mathbf{x}'を受理するか否かを決定する.
  2. \mathbf{x}'を受理した場合は\mathbf{x}'を新たな状態として採用する. そうでない場合は\mathbf{x}を新たな状態として採用する.

メトロポリス法(Metropolis Algorithm)における詳細釣り合い条件が成り立つことの確認

P(\mathbf{x'})>P(\mathbf{x})が成り立つかどうかで場合分けする.なぜならば, この不等式が真であればrが必ず1になるためである.

P(\mathbf{x'})>P(\mathbf{x})の場合

\begin{align} &\pi(\mathbf{x} \to \mathbf{x}') = Q(\mathbf{x} , \mathbf{x'}) \\ &\pi(\mathbf{x}' \to \mathbf{x}) = \frac{Q(\mathbf{x}', \mathbf{x})P(\mathbf{x})}{P(\mathbf{x}')} \\ \end{align}
  • (12)式を変形して(11)式を代入する.
    すると,
\begin{align} \pi(\mathbf{x}' \to \mathbf{x}) &= Q(\mathbf{x} , \mathbf{x}') P(\mathbf{x}) / P(\mathbf{x}') \\ \pi(\mathbf{x}' \to \mathbf{x})&= \pi(\mathbf{x} \to \mathbf{x}') P(\mathbf{x})/P(\mathbf{x}') \\ \pi(\mathbf{x}' \to \mathbf{x})P(\mathbf{x}') &= \pi(\mathbf{x} \to \mathbf{x}') P(\mathbf{x}) \end{align}

よって, 詳細釣り合い条件が成り立つ.

ただし, 以下に注意する. (12)式では遷移元が\mathbf{x}'で遷移先が\mathbf{x}である.

\begin{align*} r &= \min \left( 1, \frac{P(\mathbf{\text{遷移先}})}{P(\text{遷移元})}\right) \\ &= \min \left( 1, \frac{P(\mathbf{x})}{P(\mathbf{x}')}\right) \end{align*}

P(\mathbf{x'})<P(\mathbf{x})の場合

\begin{align} &\pi(\mathbf{x} \to \mathbf{x}') = Q(\mathbf{x} , \mathbf{x}') \frac{P(\mathbf{x}')}{P(\mathbf{x})} \\ &\pi(\mathbf{x}' \to \mathbf{x}) = Q(\mathbf{x}', \mathbf{x}) \\ \end{align}
  • (16)式を変形して(17)式を代入する.
    すると,
\begin{align} \pi(\mathbf{x} \to \mathbf{x}') &= Q(\mathbf{x}', \mathbf{x}) \frac{P(\mathbf{x}')}{P(\mathbf{x})} \\ \pi(\mathbf{x} \to \mathbf{x}')&= \pi(\mathbf{x}' \to \mathbf{x}) P(\mathbf{x}')/P(\mathbf{x}) \\ \pi(\mathbf{x} \to \mathbf{x}')P(\mathbf{x}') &= \pi(\mathbf{x}' \to \mathbf{x}) P(\mathbf{x}) \end{align}

よって, 詳細釣り合い条件が成り立つ.

メトロポリス-ヘイスティング法(Metropolis-Hasting Algorithm)における詳細釣り合い条件

P(\mathbf{x'})Q(\mathbf{x}', \mathbf{x}) > P(\mathbf{x})Q(\mathbf{x}, \mathbf{x}')が成り立つかどうかで場合分けする.なぜならば, この不等式が真であればrが必ず1になるからである.

P(\mathbf{x'})Q(\mathbf{x}', \mathbf{x}) > P(\mathbf{x})Q(\mathbf{x}, \mathbf{x}')の場合

\begin{align} &\pi(\mathbf{x} \to \mathbf{x}') = Q(\mathbf{x} , \mathbf{x'}) \\ &\pi(\mathbf{x}' \to \mathbf{x}) = Q(\mathbf{x}', \mathbf{x}) \frac{P(\mathbf{x})Q(\mathbf{x}, \mathbf{x}')}{P(\mathbf{x}')Q(\mathbf{x}', \mathbf{x})} \\ \end{align}
  • (21)式を整理して(22)式を代入する.
    すると,
\begin{align} &\pi(\mathbf{x}' \to \mathbf{x}) = \frac{P(\mathbf{x})Q(\mathbf{x}, \mathbf{x}')}{P(\mathbf{x}')}\\ &\pi(\mathbf{x}' \to \mathbf{x})P(\mathbf{x}') = Q(\mathbf{x}', \mathbf{x})P(\mathbf{x}) \\ &\pi(\mathbf{x}' \to \mathbf{x})P(\mathbf{x}') = \pi(\mathbf{x}' \to \mathbf{x})P(\mathbf{x}) \end{align}

よって, 詳細釣り合い条件が成り立つ.

P(\mathbf{x'})Q(\mathbf{x}', \mathbf{x}) < P(\mathbf{x})Q(\mathbf{x}, \mathbf{x}')の場合

\begin{align} &\pi(\mathbf{x} \to \mathbf{x}') = Q(\mathbf{x} , \mathbf{x}') \frac{P(\mathbf{x}')Q(\mathbf{x}, \mathbf{x}')}{P(\mathbf{x})Q(\mathbf{x}, \mathbf{x}')} \\ &\pi(\mathbf{x}' \to \mathbf{x}) = Q(\mathbf{x}', \mathbf{x}) \\ \end{align}
  • (26)式を整理して(27)式を代入する.
    すると,
\begin{align} &\pi(\mathbf{x} \to \mathbf{x}') = Q(\mathbf{x}' , \mathbf{x})P(\mathbf{x}')/P(\mathbf{x}) \\ &\pi(\mathbf{x} \to \mathbf{x}')P(\mathbf{x}) = \pi(\mathbf{x}' \to \mathbf{x})P(\mathbf{x}) \end{align}

よって, 詳細釣り合い条件が成り立つ.

まとめ

メトロポリス法とメトロポリス-ヘイスティング法における詳細釣り合い条件が成り立つことを確認した. 今後はメトロポリス法とメトロポリス-ヘイスティング法の実装に取り組みたい.

参考文献

Discussion