次のPRML勉強会が10章の変分ベイズなので予習。

10章概要

まずは節タイトルをを眺めていく。

  • 冒頭: EM法での厳密解の式
  • 10.1 factorized approximation
  • 10.2 混合ガウス分布
  • 10.3 変分線形回帰
  • 10.4 指数族分布
  • 10.5 ローカル変分法
  • 10.6 変分ロジスティック回帰
  • 10.7 Expectation Propgation

EM法で求める時にqなどに適当に仮定を置く、というのが10章のテーマの模様。 その考えだと10.1、10.2、10.4の3つはタイトルからだいたい内容が予想出来る。

逆に10.3と10.5は何をやるのか良く分からないな。 ちょっと軽く10.3を見てみよう。

10.3の最初の所を見ると、線形回帰で、アルファとベータが確率分布するfully bayesianなアプローチの時に、このハイパーパラメータによる積分が計算不能なので近似を使う、という話らしい。

10.91のようにfactorize出来る、という仮定での議論なのかな。

10.6の変分ロジスティック回帰は4.5でラプラス近似使ってた所をローカル変分法というのに変えてやってみる話っぽい。 4.5って何やる話だっけ？

ちらっと見直すと、こちらもロジスティック回帰のパラメータのベイズ的な取扱いの模様。

軽く見直すとパラメータの事後分布が4.142の形になって、そこにラプラス近似を使う、みたいな話がp 218にある。

これになんか違う近似を入れる、という話なのかなぁ。

まぁ概要的にはこんな物でしょう。 では実際に読み進めて行きます。

10.1 factorized distribution

なんとなく読んでいったら、10.1.2で何やってるか良く分からなくなってきたのでちゃんとここにメモを書いて読み直す。

10.1.1では分布が潜在変数Zの各要素でfactorize出来るという仮定の近似を扱っている（10.5式）

これをEM法のLに入れて変形すると10.6となり、これを各 \(z_i\) について最大化していく事を考えている。

すると答えは10.9になるとの事。同時分布のj以外での期待値にするのが、qjとしてはもっとも良い物となる、という結論。

10.1.2は何をやってるか

やってる事は二変数のガウス分布が真の分布の時に、こいつをfactorization近似を適用してみて、近似が何を捨ててるのかを実際に見てみよう、という事。

さて、factorizationの近似をした以上、各変数についてqを求めていける訳だ。 それは先程の、\(z_1\)を求める時は同時分布の\(z_2\)での期待値を求める、という奴だよな。

なんでここで分からなくなったのかを冷静に考えよう。pとqの区別が良く分からなくなったのだよな。

EM法のqとはなんだったか？

EM法としては、qというのはなんだったのか。

p450の9.4くらいに戻ってみよう。

ある観測変数Xと、潜在変数Zがあって、さらにパラメータのシータがある。

で、以下を最大にするシータを求めたい。

左辺だけで考えると、これは観測値をもとにパラメータを求める、最尤推定の式となっている。

右辺はコンプリートなデータの時の同時分布を潜在変数についてマージナライズした物となっている。

本来は左辺だけで片付く話だが、これが難しいケースの時に、けれど何らかのZを導入すると右辺については扱いやすい、という時にEM法を用いてこのシータを求める、という話だった。

さて、qというのはここでシータを求める時に、便宜的に導入されたものだ。

q自体には大した意味は無く、なんらかの関数形でこいつとシータを交互に変形していくとシータが求められる、という物に過ぎない。

ただ、qは毎回以下のKLダイバージェンスを最小化していく。

ここでpとかを適当に流すから良く分からなくなるのだな。ちゃんと注目しよう。

毎回最小化する時に、現在のシータの値での

の分布と一致させる。で、シータはだんだん真の値に近づく、という事で、この分布の近似値となっているのだな。

シータは省くと、観測値が与えられた時の潜在変数の分布か。

10.1.2のpとqとはなんぞや？

さて、10.1.2に戻ろう。

qを求める、というのは、観測変数が与えられた時の隠れ変数の分布を求める事と同じだ。

で、factorizeされる、というのは、これらの隠れ変数が独立だ、という仮定と等価となる。

pはなんぞや？というと、本来はxとzの同時分布なのだが、ここにはXの記述は無い。

だが、もともとXは計算には関係ない。Xには、ある観測値を代入してその尤度を求めている訳だから、実際の定数が代入されている。

という訳でこのXの事は単純な定数と思っておけば良いはず。

という事でpは、同時分布にXに観測値を代入したものか。

それが二変数のガウシアンだ、と。 で、その時にqを求める。 これはXの元でのZの条件付き確率分布。

分かってきた。同時分布にXの観測値を代入した値は、そのXの実現確率に引っ張られた小さい値のはずだ。

で、qは条件付き確率だからそうでは無い。

結局、10.1.2は何をやっているのか？

これはEM法の時に求めた結果のパラメータの所に、パラメータの期待値を入れた物に似た形となっている。

パラメータを点推定した物と、分布を仮定した場合にその期待値で点推定の部分を置き換えた物とがだいたい似ているのは自然。

また、逆に言えば差分が最尤推定でバイアスされている部分と言える。

さて、このrはハイパーパラメータに依存していて、ハイパーパラメータはrに依存しているのだから、この2つを交互にもとめていく必要がありそう。

具体的な手順は

  1. 現時点でのハイパーパラメータを元にrを計算する。この為にはハイパーパラメータを元にパラメータの期待値やモーメントを計算する事になる。
  2. 計算したrを元にパラメータの分布を求める。この時には更新されたハイパーパラメータも手に入る。

この2つを交互に進める事で、数値解が求まる。 これがp479の真ん中あたりに書いてある内容。

10.2.2以降の粗筋

ここからは軽そうなのでまとめて。

10.2.2は、lower boundを直接計算してみる、という話。 これは各イテレーションで低下しなくちゃいけない値として計算時のチェックに使える。 また、関数形を知っていればパラメータに対してこのlower boundを最大化するパラメータを出そうとすると、10.2.1で計算したのと同じ結果になるらしい。

10.2.3は特定の新しい観測値が得られる確率密度の話。 新しい観測値に関連した潜在変数zの値が得られて、それを用いてこの確率密度は10.78が得られて、zについて足し合わせて10.79が出る。
こいつは計算不能な事が多いので、パラメータの分布を変分近似したqに置き換えて計算していくと、10.81のようになるらしい。

10.2.4はKの値の選択。同じXの分布を生むパラメータの組み合わせが複数あるのでその分のペナルティ項を加えて観測データの尤度を比較する事でKを決められるとか。
いまいちなんで等価なパラメータの組み合わせがあると問題なのかは良く分からなかった。

Old Faithfulはイエローストーンの間欠泉の名前とか。

また、このやり方だとK個のモデルをトレーニングしないといけないので、もっと簡易的な方法として\(\pi\)の点推定をするというのがあるとか。

普通はこれも分布しているとして扱っていたのをここだけ点推定に変えるとゼロになる項が出てきてその分の潜在変数を落とせるのかな。

10.2.5は近似として導入したfactorize以外にも、厳密解の所で存在する独立性から、追加のfactorizeが発生している事と、その重要性について。
factorize出来ていると数値計算的に記憶していないといけない要素が激減するので、そいつを有効利用すべき、という話。

この追加で発生するfactorizeが成立する為の条件は条件付き独立で表現出来て、これはグラフィカルモデルのd-sepから分かる、とか。

10.3 線形回帰

ここまで読み進めた今、別段難しい事は無い。線形回帰でハイパーパラメータのアルファが確率分布する場合の話だ。
内容もストレートフォワード。

ガンマ分布は以下の形。

なお、ガンマ関数は以下。

Wikpedia:ガンマ関数

10.4 exponential family

これも大した話では無く、これまで具体的に混合ガウス分布などてやってた結果を指数関数族に拡張した場合の計算。

パラメータはデータ点の数に依存せず決まった数で、潜在変数はデータ点と同じだけある、という仮定のもとに計算をする。

10.4.1は指数関数族あんま関係ない気もするが、factorizeの近似の一般化として、ベイジアンネットワークで記述されるfactorizeのケースの時の計算について触れられていて、10.123の形で近似を入れてメッセージパッシングの形で扱うとかなり一般的に扱える、とか。 これはお話だけ。
Infer.NETとかこの辺の話なのかしら？

10.5 local variational methods

ただ読んでいったら10.130のあたりで分からなくなったので、ひたすら計算を追う。

手を動かしてみれば別段難しい事は無いので、式だけ。

10.6 ロジスティック回帰の変分ベイズ

ロジスティック関数に対して、lower boundで抑えて結果を計算して、このlower boundをパラメータで最大化する事で近似解とする、という話。

さらにハイパーパラメータが確率分布するように取り扱うには、wとアルファがfactorizeできるという近似を入れた変分ベイズで出せる、という話だが…

個人的に全然テンション上がらないのであんまり真面目にノートを追う気も起こらない。 やってる事はこれまでの集大成の感じはあるのだが、結果が二値分類がちょっと正確に出来ます、って…

計算はごついけれど、ストーリー自体はこれまでと同じで追いやすいので、まぁいいか、という事でサボり。

10.7 expectation propagation

2.4の復習

10.186の所で2.226が出て来るのと、この周辺はexponential familyの話が多いので、簡単に復習しておく。

2.226の導出周辺

定義は以下。

積分して1になる事から、以下。

以下、イェータの満たすべき条件を考える。 両辺をイェータで微分してみて以下。

ガウス分布との対応

Exponential familyの定義と見比べて、以下。

KLダイバージェンスの定義式

reverseが出てきて混乱したので書いておく。

左側で積分。

10.7.1 計算メモ

ゆとりなので全部追う気は無いが、結果を読み解ける程度のメモはしていこうと思う。

10.216の導出。

この第一項は、定数係数を除くと以下みたいな式となっている。

この計算は昔やった記憶があるがめんどかった気がする。
まぁ10.216になりそうな気もするからいいか。

3つの統計量を9章と比べて考える

さて、10.51から10.53は、9章のMステップと似ているので併記しておこう。

9.26に対応する物が無いが、これは確率分布するからそういう物かもしれない。

基本的にある観測値が、幾つかのzにまたがっていて、それぞれの寄与でxが実現されている、と考える。 例えば、k=1に0.3だけ所属していて、k=2に0.7だけ所属している、 みたいな部分的に所属している、と考えた時の実効的な個数などを問題とする。

10.51は観測値のうちkに所属している実効的な数。

10.52はxのうちkに所属している物だけを集めて平均をとった物。

10.53はその分散を意味する。

これが9.24や9.25に対応するのは明らかに見える。 むしろこれらが違う変数を割り当てている方が注意が必要で、これはパラメータが確率分布するので、これらの値がそのままパラメータになる訳じゃなく、事後確率を更新する事で求まる、という事か。

感覚的には事前確率とこれらの値の平均というか、間くらいの所に落ち着くのだろう。

パラメータの事後分布を求める

Zの事後分布が出たので、次はパラメータの事後分布を求めていく。

10.54からがその話となる。

まず同時分布の式を出し、そこからおのおのの分布を変分ベイズで出す。 同時分布はほぼ定義より10.54となる。 あとは\(\pi, \mu, \Lambda\) を別々に出せば良い。 といっても\(\mu, \Lambda\)はお互い依存しあっているので、ちょっと注意は必要そう。

で、\(\pi\)の分布を求めると、10.57となる。 10.58を見ると、予想通り事前分布と事後分布の間くらいの所に落ち着いた分布になっていそう。

10.59の導出はゆとりなのでgive up。 ただ、結果を解釈する事はやっておこう。

事前分布をGaussian-Wishartと仮定していたので、事後分布もパラメータを少し更新したこの分布となっているのだろう。

気分的には中心が先程の\(\bar{x}\) と\(m_0\)の間くらいに\(\mu\)の中心が、\(W_0\)と \(S\)の間くらいに\(\Lambda\)の中心が来るのだろう。

10.61と10.62がまさにそういう感じの式となっている。

その他細かい結果もだいたいハイパーパラメータが事前分布と標本から求めた類似の統計量との間くらいになるようになっているので、結果自体は納得出来る。

パラメータとrの式から数値解を求める

10.64からは、実際に数値解を求める話となる。

rを求めるにはパラメータの期待値が必要で、パラメータの分布はrに依存している。 EM法と似た構成となっている。

さて、パラメータの期待値はハイパーパラメータを用いて表す事が出来て、それが10.64から10.66となる。

これらの期待値をrの定義に入れてやるとrが出る。 それが10.67。

これはEM法の時に求めた結果のパラメータの所に、パラメータの期待値を入れた物に似た形となっている。

パラメータを点推定した物と、分布を仮定した場合にその期待値で点推定の部分を置き換えた物とがだいたい似ているのは自然。

また、逆に言えば差分が最尤推定でバイアスされている部分と言える。

さて、このrはハイパーパラメータに依存していて、ハイパーパラメータはrに依存しているのだから、この2つを交互にもとめていく必要がありそう。

具体的な手順は

  1. 現時点でのハイパーパラメータを元にrを計算する。この為にはハイパーパラメータを元にパラメータの期待値やモーメントを計算する事になる。
  2. 計算したrを元にパラメータの分布を求める。この時には更新されたハイパーパラメータも手に入る。

この2つを交互に進める事で、数値解が求まる。 これがp479の真ん中あたりに書いてある内容。