ベイズ混合モデルにおける近似推論③ ～崩壊型ギブスサンプリング～

 さて、今回は混合モデルに対する効率的な学習アルゴリズムとして、崩壊型ギブスサンプリング（Collapsed Gibbs Sampling）*1を紹介します。ベースとなるアイデア自体はギブスサンプリングそのものなのですが、標準的な方法と比べて種々の計算テクニックが必要になってくるのと、実際に推論精度が向上することが多いのでここで紹介したいと思います。 

また、今回の実験で使用したソースコードはGithubにおいてあります。

GitHub - sammy-suyama/MLBlog: source codes used in MLBlog

 

[必要な知識]

下記をさらっとだけ確認しておくといいです。

• ポアソン混合モデルにおけるギブスサンプリング
• ベイズ推論
• カテゴリカル分布 

 

１、崩壊型ギブスサンプリングのアイデア

基本的なアイデアとしては、対象の確率モデル（同時分布）から一部の確率変数をあらかじめ周辺化除去（marginalise out）しておき、より数の少ない確率変数をサンプルすることによってサンプリングの効率性を高めようとするものです*2。

ここでは、崩壊型ギブスサンプリングの位置付けを明確するために、いくつかの異なるギブスサンプリングの手法を比較紹介します。

ある確率分布があり、この分布から,,をすべて同時にサンプルするのが不可能であるとしましょう。この様な問題に対して、一般的なギブスサンプリングは次のように確率変数を１つ１つ別々にサンプリングするという作戦を取ります。

これは「直接同時にサンプルするのが難しいのであれば、他の値を固定したより簡易な条件付き分布から順次サンプルすればいい」というアイデアなのでしたね。

次に、ブロッキングギブスサンプリング（blocking Gibbs sampling）というのもあり、この場合はいくつかの変数をまとめてサンプリングします。

「すべて同時にサンプルするのは難しいが、,の同時分布なら簡単に計算できる」というような場合ではこちらを使った方が断然性能が良いです。

さて、本題の崩壊型ギブスサンプリングでは上記2つの手法とは異なる作戦を取ります。この手法では、モデル中のいくつかの変数を周辺化除去することによって、より少ない数の変数だけをサンプルします。（Collapsed spaceでサンプルする、と呼ぶようです。）

最初のに対する積分除去が解析的に出来ることと、それによって得られた同時分布からの交互サンプリングが簡単に行えることが適用できる条件です。単純なギブスサンプリングと比べて、そもそもサンプルする変数の種類が少ないので、より速く目的の分布からのサンプリングが得られることが期待できます。

一般的にはブロッキングギブスサンプリングや崩壊型ギブスサンプリングの方が、単純なギブスサンプリングより性能が良い傾向にあります。しかし単純ギブスサンプリグの方が1回あたりの計算が軽いので、問題によってはより早く所望の結果が得られる場合があります。実用上では、まず最初の２つの手法を検討してみて、それでも精度がイマイチだった場合には崩壊型ギブスサンプリングを試してみるのがオススメです。*3

 

２、混合モデルからパラメータを周辺化除去

まずは標準的な混合モデルを考えてみましょう。ここではが個の観測データの集合、が各データに対する隠れ変数の集合、がの分布を制御するパラメータ、が隠れ変数の分布に対するパラメータです。

がパラメータに対する事前分布です。一般的に混合モデルでは、はカテゴリカル分布を使います。はの発生をモデル化するための確率分布であり、今回は例としてポアソン分布を使いますが、用途次第ではガウス分布でも多項分布でも構いません。

さて、この式に対するグラフィカルモデルを次の様に表現してみることにします。

混合モデルに対する崩壊型ギブスサンプリングでは、パラメータをモデルからあらかじめ周辺化除去するのが一般的なやり方です。

こんな感じに積分除去すると、対応するグラフィカルモデルは次のようになります。一応これをここでは周辺化モデルって呼ぶことにします。

当然、積分除去されたやはもはや周辺化モデルには現れません。しかし、そのおかげでモデルが簡単になったかというとそうでもないようですね。隠れ変数は元のモデルでは（が与えられたとき）独立性を保っていたもの、が積分除去された途端にすべての隠れ変数が依存関係を持つようになってしまいました。同様に、もが積分除去されたために、すべてのデータ点がお互いに依存するようになってしまいました*4 。このようなモデルからをサンプリングする計算式を導くのが目標になります。

 

３、周辺化モデルからのサンプリング

さて、周辺化モデルは非常にタイトな変数依存関係を持ってしまってしまっているので、を一度に同時にサンプルすることはできません。しかし、,・・・、のように、１つずつサンプルすることが可能である場合があります。あるを、以外のサンプルがすでに与えられた状態でサンプルできないか考えてみることにします。この条件付き分布は次のように計算することが出来ます。

1行目では、ただ単にベイズの定理の分母を無視したもので事後分布を書き直しただけです。2行目ではそれをバラバラに展開しました。真ん中の項はを除いたグラフィカルモデルを考えれば、独立性から、

であることが示せるため、結果今サンプルしたいとは無関係になるので3行目では除去されています。

さて、出来上がった最後の式をよく見てみると、目的の確率分布がとの２つの分布に分解されることがわかりますね。これらの確率分布はどういう意味を持っているのでしょうか。

まず始めにですが、次のようなに関するある種の予測分布であると解釈することができます。

まず、すでにサンプルされていると仮定しているを観測データかのように扱うことで、のある種の事後分布のようなものが計算できます。この事後分布は、次の様にベイズの定理を使って簡単に計算することができます。

今回ははカテゴリカル分布で、はディリクレ分布です。したがって、次のようなディリクレ分布事後分布が求まります。

さて、このに関する事後分布を使って予測分布を計算してみることにします。結局、有限シンボルの離散分布であることがわかるので、結果はカテゴリカル分布として表現することができます。

同様にして、残るの項も次の様にに関する（の条件付き）予測分布になることがわかります。

次のセクションで解説しますが、ポアソン・ガンマモデルではこの積分は解析的に実行できます。*5

ところでは実際にはデータとして与えられているので、より正確に言うとこの項は予測尤度とでも言った方が正しいかもしれません。肝心な点は、具体的なを与えることによってこの尤度が評価できるという点です。に対する種類の実現値１つ１つに対して式を評価し、正規化してあげれば、をサンプルするための離散分布が計算できることになります。

 

４、ポアソン・ガンマモデルにおけるベイズ推論

さて、ここまでで崩壊型ギブスサンプリングの説明はほぼ終わりなのですが、残された課題は次の条件付き予測尤度の計算ですね。

先ほどと同様、手順としてはまず右側の事後分布を計算し、その後で尤度関数と掛け合わせてからを積分除去します。

まず事後分布の計算からやってみましょう。今回は次元のポアソン・ガンマモデルを仮定しているのではポアソン分布のパラメータになります。こちらもベイズの定理を使うと、

として計算をすることができ、事後分布は当然ガンマ分布になります。

さて、この事後分布を使って予測分布を計算してみましょう。幸い事後分布がガンマ分布のままなので、このまま解析的にパラメータを積分除去することができ、得られる分布は次のような負の二項分布（Negative Binomial Distribution）になります。

ここでは、負の二項分布の定義として次のような表現を使いました。

混合モデルにおいて崩壊型ギブスサンプリングが適用できる条件としては、観測モデルの事後分布と予測分布がともに解析的に実行できることです。今回のようなポアソン・ガンマモデル以外でも、共役の関係にある分布を選べばこれは必ず実行できます。  

 

5、動作確認実験

さて、今回導いたポアソン混合モデルにおける崩壊型ギブスサンプリングアルゴリズムをJuliaで実装し、動作確認をしてみます。ソースコードはGithubにあります。

観測データ数N＝１０００、次元数D＝２、クラスタ数K＝８として混合モデルを学習させてみました。次の図は１００回までのサンプル結果をアニメーションにしています。 

左下の方のクラスタはかなり密集していますが、なんとかうまく分離できているように見えますね。右図でサンプル回数が５０回を超えてもまだパタパタしている箇所がありますが、それらのデータ点はクラスの所属が「あいまい」であることを示しています。

 

さて、せっかくですので前回導出した標準的なギブスサンプリングと簡単な性能比較をしてみます。下の図は横軸をサンプル回数、縦軸を最後のサンプルに対する尤度にしてグラフを描いたものです*6。GS（ギブスサンプリング）、Collapsed（崩壊型ギブスサンプリング）ともに５回ずつ実験しています。

傾向として、崩壊型ギブスサンプリングの方が少ない繰り返し回数で、より高い尤度に達しているようですね。最終的な到達点も崩壊型ギブスサンプリングの方が高い傾向があります。

ちなみに計算時間で比較すると、１００回のサンプリングでギブスサンプリングが平均7.03秒、崩壊型ギブスサンプリングが平均7.15秒と、あまり差がないようです。 一般的には、崩壊型ギブスサンプリングにおける予測尤度の計算が遅くなる可能性がありますが、今回の場合は特にボトルネックにはなっていないようです。

 

6、まとめ

さて、今回はベイズポアソン混合モデルを題材にして、崩壊型ギブスサンプリングのアルゴリズムを導いてみました。今回の簡易実験が示す通り、標準的なギブスサンプリングアルゴリズムよりも性能が高くなる傾向があるようです。オススメとしては、まず標準のギブスサンプリングを導いてみて、その性能で満足いかないようであれば崩壊型ギブスサンプリングを導出してみるのが良いかと思います。ただし、解析的にパラメータを積分除去するためには、今回のポアソン分布とガンマ分布の例のように、観測モデルが共役の関係性を持っている必要があります。

同様のアイデアは時系列モデルである隠れマルコフモデル（Hidden Markov Model, HMM）、トピックモデルにおけるLatent Dirichlet Allocation（LDA）に適用されているほか、中華料理店過程（Chinese Restaurant Process, CRP）やインド料理店過程（Indian Buffet Process, IBP）などのベイジアンノンパラメトリクスの手法にも応用されています。

 

*1:参考: Machine Learning: A Probabilistic Perspective by Kevin Murphy | The MIT Press

*2:このことから、崩壊型ギブスサンプリングは周辺化ギブスサンプリングとも呼ばれることもあるようです。

*3:ただし、モデルによっては解析的にうまくアルゴリズムを導けない場合や、解析的な周辺化のために積分近似を行う場合もあります。

*4:これは一般的にはであるためです。

*5:他にもガウス・ウィシャートモデルのように、共役関係が成り立っていれば解析計算できます。

*6:尤度で性能を比較するのはあまり望ましくないですが、今回は完全に同じモデルからアルゴリズムを導いているので、尤度でもほぼ「良いクラス割り当て」を反映しているとしてしまってもいいんじゃないかと思います。