ベイズ推論により混合分布のパラメータ推論をやってみる　〜混合数の推論１〜

【概要】

混合分布（混合モデル）はモデルを潜在変数でスイッチする構造を持ったモデルであり、実用的な観点でも面白いです
弊ブログでは数回にわたって、混合分布を使って遊んでみています（これが4記事目）
第4弾、第5弾では、これまで既知としていた混合数（クラスタ数）を含めて推論してみます（記事が長くなったので分割しました）
まだまだ勉強中なので、間違いなど指摘いただけると助かります

【目次】

はじめに
混合モデル
- 混合モデルにおける課題
- 課題に対するアプローチ
ディリクレ過程の実装
- 中華料理店過程(CRP; Chinese Restaurant Process)
- 棒折り過程(SBP; Stick Breaking Process)
まとめ
参考文献

はじめに

機械学習や統計の問題では、手元にあるデータを解釈して応用しようとしますね。この時、明に暗になんらかの「確率モデル」を仮定しているはずです。

確率モデルの中でも、混合分布（混合モデル）は、複数の確率モデルの組み合わせとして定義されており、複雑なデータ構造を表現できます。応用としても、クラス分類や多クラスの回帰など面白い実用例があります。

ということで当ブログではここまで3回にわたって、混合モデルを使って遊んできました。（第3回からここまでだいぶ時間が空いてしまいましたが。。。）
第4弾となる本記事では、これまで既知としていた混合数（クラスタ数）もデータに合わせて推論してみたいと思います。ディリクレ過程混合モデル(DPMM; Dirichlet Process Mixture Model)と呼ばれるモデルです。ほんとは一つの記事でまとめたかったのですが、長くなってしまったので、二つに分割します。DPMM自体の実装については次回で扱います。

【ベイズ推論により混合分布のパラメータ推論をやってみる記事一覧】

第1弾： learning-with-machine.hatenablog.com
第2弾： learning-with-machine.hatenablog.com
第3弾： learning-with-machine.hatenablog.com
第4弾：当記事
第5弾： learning-with-machine.hatenablog.com

【トップに戻る】

混合モデル

混合モデルについては、これまでの当ブログの記事や参考文献3、参考文献4などを参照ください。複数の確率モデルの組み合わせで定義されており、潜在変数 $s$ でそれらの確率モデルをスイッチする構造を持った確率モデルです。複数の確率モデルを組み合わせることができるため、複雑な構造のデータをモデリングするのに利用できます。

グラフィカルモデルで表現すると以下の通りです。

f:id:hippy-hikky:20200308143035p:plain

通常の混合モデルでは、データ $x$ がどのコンポーネント(クラスタ)から生成されたかを割り当てる変数が $s$ です。つまり $s$ は、クラスタのIDを示す離散変数なので、カテゴリ分布から生成されるとします。そして、カテゴリ分布のパラメータであるそれぞれのコンポーネントの割合を示す変数 $\pi$ はディリクレ分布から生成されるとするとします。

f:id:hippy-hikky:20200825235330p:plain:h70

データ $x_n$ は、割り当てられたコンポーネントのパラメータ $\theta_k$ に基づいて生成されるとします。

f:id:hippy-hikky:20200825235748p:plain:h60

混合モデルにおける課題

上記の混合モデルでは、コンポーネントの数は $s$ の種類、つまり、 $\pi$ の次元数で決まります。

現実の問題では、予めクラスタ数が決まっている問題しか扱えないのであれば柔軟性がありません。既知のクラスタ数で学習したモデルに対して、未知のクラスのデータを入力した場合、通常は既知クラスのいずれかのクラスタに分類してしまいます（e.g. Open Set Recognition）。

そこで、クラスタ数も含めてデータから構造を学習(推論)したいと言う欲求が現れます。

課題に対するアプローチ

上記の欲求に対して、コンポーネント数を無限として無限クラスの混合モデルを考えるのがディリクレ過程混合モデルです。他にも、混合数を複数パターン推論してAIC等の基準でモデル選択するとの考え方もありますが、複数パターンのモデルの推論が必要なため計算効率が悪いです。

無限次元の混合モデルを扱うためのアプローチとして以下の二つのアプローチがあります。

$\pi$ を周辺化除去して、 $s_n$ について無限クラスの分割方法に関する事前分布を考える
無限次元の $\pi$ を表現する事前分布を考える

なおここで、「無限次元」と表現していますが、データは有限です。なので、実際に観測されるクラスタ数は有限になります（なので計算ができる）。予め明示的にクラスタ数を設定しているのではなく、クラスタ数も確率的に変動するという点で「無限次元」の混合モデルと表現されています（と理解しています）。

【トップに戻る】

ディリクレ過程の実装

正直に言うと、ディリクレ過程について説明できるほど数学的な理解はできていません（いずれ理解をまとめて記事書こうと思います）。表面的な理解のままとりあえず実装を進めて行きたいと思います。

上記の二つのアプローチを実装するための代表的な二つの確率過程があります。 1のアプローチでのディリクレ過程の実現例として「中華料理店過程」と呼ばれる確率過程があります。また、2のアプローチの実現例として「棒折り過程」と呼ばれる確率過程があります。

以下ではそれぞれの確率過程を実装して無限次元のクラスタの表現というのがどのようなイメージなのかを見ていきたいと思います。

中華料理店過程(CRP; Chinese Restaurant Process)

CRPは、N個のデータ $X (=[x_1, x_2, \cdots, x_N ])$ を分割（クラスタリング）することを考えた際に、分割数を予め定めずどのような分割の仕方があるのかについての確率過程です。分割の仕方、つまり、クラスタのIDを示す潜在変数 $s$ の種類数や割り当ての確率を含めての事前分布になります（ $\pi$ が現れてこない→周辺化除去されている）。

CRPは以下の流れで分割数の事前分布を構成します（詳しくは参考資料2を参照してください）。

最初の客（データ）は任意のテーブル（クラスタ）に着席する（割り当てられる）
２番以降の客（データ）は、すでに $n_i(>0)$ 人着席しているテーブルiに確率 $\frac{n_i}{(n-1+\alpha)}$ で着席し、誰も座っていないテーブルに確率 $\frac{\alpha}{(n-1+\alpha)}$ で着席する

テーブルがクラスタを表し、お客さんがデータを示しています。このようにしてN個のデータを予め定めないK個のクラスタに分割する仕方を表すのがCRPです。パラメータは $\alpha$ で、大きいとそれだけ新規のテーブルに座る確率が上がるため、多数のクラスタができやすくなります。

CRPを実装し、1,000個のデータがどのように分割されるのかシミュレーションしてみたのが以下の結果です。

f:id:hippy-hikky:20200826004446p:plain — 2種類のパラメータ(2, 10)でのCRPの結果．横軸がデータ数を示し，縦軸がテーブル数(クラスタ数)を示す．パラメータが大きいほど多数のテーブルにお客さんが座ることがわかる．

それぞれのクラスタにどのくらいのデータが属するのかを示したのが以下の図です。

f:id:hippy-hikky:20200826004724p:plain — 各テーブルに何人のお客さんが座っているのかを示した棒グラフ．横軸はテーブル番号，座っている縦軸はお客さんの数．

$\alpha$ が大きいと、多数のクラスタができやすくなるということが図からわかります。なお、これはあくまでも「事前分布」なので、データがこのようにクラスタリングされるということを言っているわけではないことに注意です（これは私も理解に時間がかかりました。DPMMの実装例を確認すると気持ちがわかりやすいと思います。）。

CRPの実装例は以下のnotebookを参照してください。

棒折り過程(SBP; Stick Breaking Process)

CRPは、上述の通り、有限のデータがどのように分割されるのかを示すものでした。つまり、明示的に無限次元を扱ってはいなかったのです。これに対して、無限次元の混合比( $\pi = [\pi_1, \pi_2, \cdots$ ])を明示的に扱うのが棒折り過程（SBP）です。

$\pi_k$ はk番のクラスタの割合を示す変数であり、 $\sum^K_{k=1}\pi_k = 1.0$ です。つまり、1.0を無限個に分割するという考えで無限次元の混合比を扱います。そのアルゴリズムは以下の通りです。

長さ1.0の棒を $v_1 : (1-v_1)$ の比で折る（ $\pi_1 = v_1$ ）
残りの棒(長さ $(1-v_1)$ )を $v_2 : (1-v_2)$ の比で折る（ $\pi_2 = v_2$ ）
以下同様に繰り返し

ここで、 $v_i$ は割合であり、ベータ分布( $\mathrm{Beta}(1, \alpha)$ )から生成します。

f:id:hippy-hikky:20200826010348p:plain:h70

パラメータ $\alpha$ はベータ分布のパラメータであり、大きい場合には小さい値がより生成されやすくなります。そのため、 $v_i$ が小さい多数のクラスタが生成されるようになります。
ベータ分布( $\mathrm{Beta}(1, \alpha)$ )の確率密度は以下の通りです。

f:id:hippy-hikky:20200826010607p:plain — ベータ分布の確率密度．パラメータが異なる3つの確率密度をプロット．パラメータが大きいほど低い値にピークが立つため，生成される値は小さくなる．

SBPを実装し、混合比がどのように割り当てられるのかをシミュレートしてみました。なお実装にあたって、無限次元のベクトルを生成することはできないので、十分大きな次元を予め設定します（上限を設定しない仕組みも提案されているらしいですが論文読んで無いです。参考文献3を参考にしてください）。

f:id:hippy-hikky:20200826011044p:plain — 2種類のパラメータ(2, 10)でのSBPの結果．縦軸が混合比の累積確率を示し，横軸が分割番号(クラスタ番号)を示す．

この図から、 $\alpha=2$ では累積確率のカーブが急峻であり、 $\pi_i$ の値が大きい少数のクラスタで累積確率がほぼ1に達していることがわかります。これがCRPでのテーブル数が少ないケースに対応しています。

次に、それぞれの $\pi_i$ の分布を確認してみます。

f:id:hippy-hikky:20200826011434p:plain — 分割した棒の大きさ(各クラスタの混合比)の分布．

$\alpha=2$ のケースでは、大きい値の少数の棒がありますが、他はほぼ棒が見えないです。つまり、クラスタ数が少なくなりがちになることを示しています。

ここでSBPでは理論的には無限次元の $\pi$ を考えるのですが、上図のように有限な数のクラスタ $\pi_k$ が意味のある値になり、以降のクラスタの $\pi_k$ はほぼ0になります。棒を折っていくと原理的には無限回折ることはできますが、ある程度折れば微小になると言うのはイメージつくかなと思います。これが無限の混合数のイメージで、意味のあるクラスタは無限のクラスタの中の一部ということだと理解しています。

以上の結果を生成する、SBPの実装例は以下のnotebookを参照してください。

【トップに戻る】

まとめ

今回は、混合モデルにおける混合数の推論まで行うDPMMの準備として、ディリクレ過程を計算機で実装するための二つの確率過程を確認しました。

これらは共に混合モデルの確率変数の事前分布となります。私自身がそうだったのですが、これだけでは何が嬉しいのかわかりにくいんじゃないかと思います。実際にDPMMとしてディリクレ過程を利用すると、なんとなく理解ができそうな気がします。DPMMの実装については次回扱いますので、よければ参考にしていただけたらと思います。

【トップに戻る】

参考文献

[1] Dirichlet process mixtures for density estimation, PyMC3 Examples, 2020/08/20アクセス
[2] 石井, 上田, 続・わかりやすいパターン認識, オーム社, 2014(初版)
- これは実際わかりやすかった。扱ってる内容がハードなので難しいんだけど、それでもわかりやすかったです。これは必読です。