はじめてのMCMC (ギブス・サンプリング)
こんにちはtatsyです。
前回のメトロポリス・ヘイスティングス法に引き続きギブス・サンプリングについて解説したいと思います。どうでも良い話ですが、「ギブズ」サンプリングではなく「ギブス」サンプリングなんですね。いや、本当どうでもいい話です。
ちなみに前回の記事はこちらです。
ギブス・サンプリングとは？
M-H法ではガウス関数などのサンプリングが容易な任意の提案分布を用いました。M-H法は適用先の不変分布が評価できる (この時、規格化定数は未知でも良い)ことだけが必要な条件で、適用範囲が広いのが特徴です。
その一方で、M-H法は採択確率によって、新しいサンプルを採択するのか棄却するのかを決定するため、新しい未知のサンプルを多く得づらいという問題があります。
この問題を解決した手法の一つがギブス・サンプリングで、不変分布(=サンプリングを行いたい確率密度関数)が特殊な形をしている場合に、新しく提案された値を常に採択しながら効率的にサンプリングを進めることができます。
提案分布が次の形で求めたい不変分布と、以下の関係を満たすと仮定するのが、ギブス・サンプリングのポイントです。
$$ q(x | x^t) = \begin{cases} p(x_i | x_{-i}^t) & \text{if } x_{-i} = x_{-i}^t \\ 0 & \text{otherwise} \end{cases} $$
なお、この式では、$n$次元空間の点$x \in \mathbb{R}^n$について、特定の次元$i$を除いた$\{ x_1, \ldots, x_{i-1}, x_{i+1}, x_n \}$を$x_{-i}$と略記しています。
この式は、提案分布により新しく提案される$x$は、特定の次元$x_i$のみが変化する形で、事前の状態$x^t$から遷移が起こることを表しています。
面白いことに、上記の条件の元では、1ステップにつき次元1つだけを更新するようなマルコフ連鎖を考えると、採択確率が常に1になります。
具体的に、とある次元$i$が変化するとして、詳細つり合い条件を考えてみると、採択確率$a(x^{t+1} | x^t)$は次のように書けます。
$$ \begin{array}{rcl} a(x^{t+1} | x^t) &=& \min \left[ 1, \frac{p(x^{t+1}) q(x^t | x^{t+1})}{p(x^t) q(x^{t+1} | x^t)} \right] \\ &=& \min \left[ 1, \frac{p(x^{t+1}) p(x_i^t | x_{-i}^{t+1})}{p(x^t) p(x_i^{t+1} | x_{-i}^t)} \right] \\ &=& \min \left[1, \frac{p(x_i^{t+1}|p(x_{-i}^{t+1}) p(x_{-i}^{t+1}) p(x^t | x^{t+1})}{p(x_i^t | x_{-i}^t) p(x_{-i}^t) p(x^{t+1} | x^t)} \right] \\ &=& \min \left[ 1, \frac{p(x_{-i}^{t+1})}{p(x_{-i}^t)} \right] = 1 \end{array} $$
上の式では3番目の等式においてベイズの定理を使い、最後の等式においては、特定の次元$i$以外が変化していないことを利用しています。
具体例
それでは具体的な不変分布として $$ p(x, y) = \frac{4 \pi}{\sqrt{3}} \exp \left( -\frac{1}{2} (x^2 -xy + y^2) \right) $$ を考えてみます。
この式は2変数の関数なので、ギブス・サンプリングに用いる提案分布として$x$と$y$のそれぞれを更新するものを考えます。この場合、$x$の提案分布は$y$を定数として不変分布を整理すれば良いので、以下のような提案分布が得られます。 $$ q(x | y) \propto \exp\left( -\frac{1}{2} (x^2 -xy + \frac{1}{4}y^2) \right) \cdot \exp \left(-\frac{3}{8}y^2 \right) \propto \exp \left( -\frac{1}{2} \left( x - \frac{y}{2} \right)^2 \right) $$ 同様にして$y$に対する提案分布は以下のように書けます。 $$ q(y | x) \propto \exp\left( -\frac{1}{2} \left( y - \frac{x}{2} \right)^2 \right) $$
これを利用したギブス・サンプリングのコードが以下になります (結果はコードをスクロールしてください)。
結果から分かる通り、対象の二変数確率分布がおよそ正しくサンプリングできていることが分かります。
まとめ
という感じでギブス・サンプリングも解説おしまいです。M-H法よりも適応範囲は狭いですが、使えるときには力を発揮する方法です。
またこの記事のつづきでハイブリッド・モンテカルロについても解説しておりますので、もしよろしければそちらもご覧ください。
今回の記事を書くにあたり、いくつか参考にした文献を挙げておきます(書籍はamazonに飛びます)。
- パターン認識と機械学習 下 (著: Bishop)
- マルコフ連鎖モンテカルロ法 (著: 豊田秀樹)
- An Introduction of MCMC for Machine Learning (Andrieu et al. 2003, Machine Learning)
今回はここまでです。最後までお読みいただきありがとうございました。