ベイズ統計の理論と方法：ノート5.2節その1（平均場近似）

以下の本を読みます。キャラクターは架空のものです。解釈の誤りは筆者に帰属します。お気付きの点がありましたらコメント等でご指摘いただけますと幸いです。

ベイズ統計の理論と方法

作者:渡辺澄夫
出版社/メーカー: コロナ社
発売日: 2012/03/01
メディア: 単行本

書籍サポートページ

ベイズ統計の理論と方法：ノート1：主要な量の定義。
ベイズ統計の理論と方法：ノート2：指数型分布の場合の具体的な事後分布とか。
ベイズ統計の理論と方法：ノート2.5：コイン投げの例とか。
ベイズ統計の理論と方法：ノート3： 1章の章末問題と、真の分布と確率モデルの関係。
ベイズ統計の理論と方法：ノート4章その0.1：開集合とか。
ベイズ統計の理論と方法：ノート4：キュムラント母関数。
ベイズ統計の理論と方法：ノート4章その1：多様体の導入だけ。
ベイズ統計の理論と方法： 4章の章末問題【2】（実対数閾値とその多重度）

f:id:cookie-box:20200101101603p:plain:w60

…データから混合指数型分布 $p(x|w) = \sum_{k=1}^K a_k v(x) \exp \bigr( f (b_k) \cdot g(x) \bigl)$ を推定するのに、

「 $k$ の目が出る確率が $a_k v(x) \exp \bigr( f (b_k) \cdot g(x) \bigl)$ であるようなサイコロ」

を仮想的に考えて、

「 $n$ 個のデータ $x^n$ を受け取ったとき、パラメータ $a_{1:K}, b_{1:K}$ は何であったろうか」という問題を
「 $n$ 個のデータ $x^n$ を受け取ったとき、パラメータ $a_{1:K}, b_{1:K}$ と各データのサイコロの目 $y^n$ は何であったろうか」という問題に

読み換えたのですよね？そうすることで事後分布 $P(y^n, a_{1:K}, b_{1:K} | x^n)$ を、独立な分布の積である試行分布 $q(y^n) r(a_{1:K}, b_{1:K})$ で平均場近似することが実現できるというのはわかります。わかりますが、そもそもそんなサイコロを考えなかったらどうなるんです？

f:id:cookie-box:20200101101614p:plain:w60

…まずサイコロを導入した場合はどうなったの？

えっと、そもそも、 $p(w_1, w_2)= \exp \bigl( -\beta H(w_1, w_2) \bigr) / Z(\beta)$ という形で表される目標分布を、独立な分布の積である試行分布 $q(w_1, w_2) = q_1(w_1)q_2(w_2)$ で（カルバック・ライブラー情報量 $D(q || p)$ の意味で）最もよく近似したとき、 $q_1(w_1)q_2(w_2)$ を $p(w_1, w_2)$ の平均場近似とよぶのですね。

それではどのような $q_1(w_1)q_2(w_2)$ が $p(w_1, w_2)$ の平均場近似なのか考えると、これらのカルバック・ライブラー情報量は具体的に、

$D(q || p) = \int q_1(w_1)q_2(w_2) \log \bigl( q_1(w_1)q_2(w_2)/p(w_1, w_2) \bigr) dw$
$\quad = \int q_1(w_1)q_2(w_2) \log q_1(w_1)q_2(w_2) dw - \int q_1(w_1)q_2(w_2) \log p(w_1, w_2) dw$
$\quad = \int q_1(w_1)q_2(w_2) \log q_1(w_1)q_2(w_2) dw + \int q_1(w_1)q_2(w_2) \bigl( \beta H(w_1, w_2) + \log Z(\beta) \bigr) dw$
$\quad = \int q_1(w_1)q_2(w_2) \log q_1(w_1)q_2(w_2) dw + \beta \int q_1(w_1)q_2(w_2) H(w_1, w_2) dw - \beta F(\beta)$
$\quad = \int q_1(w_1) \log q_1(w_1) dw_1 + \int q_2(w_2) \log q_2(w_2) dw_2 + \beta \int q_1(w_1)q_2(w_2) H(w_1, w_2) dw - \beta F(\beta)$
$\quad \geqq 0$

となりますから（ただし、 $F(\beta) = - \beta^{-1} \log Z(\beta)$ は自由エネルギーです）、これを最小化する $q_1(w_1)q_2(w_2)$ が $p(w_1, w_2)$ の平均場近似ということになります。等号成立は $p(w_1, w_2) = q_1(w_1)q_2(w_2)$ のときですが、いま、試行確率分布は成分どうしが独立という制約がありますし、この等号が達成できる保証はありません。ところで、上式の下から2行目の、試行分布に依存しない $- \beta F(\beta)$ を除いた部分を、平均場自由エネルギー $\beta \hat{F}(\beta)$ とよぶとのことです。もし試行分布が目標分布を実現できれば真の自由エネルギーに等しくなる量ですが、一般的には真の自由エネルギー以上になる量ですね。最小化すべき目的関数の別表現といった感じがします。それで、もし $\beta \hat{F}(\beta)$ が最小値になっているならば、この最小化問題のラグランジュ関数

$\int q_1(w_1) \log q_1(w_1) dw_1 + \int q_2(w_2) \log q_2(w_2) dw_2 + \beta \int q_1(w_1)q_2(w_2) H(w_1, w_2) dw_1 dw_2$
$\quad + \lambda_1 ( \int q_1(w_1)dw_1 - 1 ) + \lambda_2 ( \int q(w_2)dw_2 - 1 )$

の $q_1, \, q_2$ に関する1次の変化分が恒等的に $0$ になっているはずです。このことから、

$\quad q_1(w_1) \propto \exp \bigl( - \beta \int q_2(w_2) H(w_1, w_2) dw_2 \bigr) \propto \exp \bigl(\int q_2(w_2) \log p(w_1, w_2) dw_2 \bigr)$
$\quad q_2(w_2) \propto \exp \bigl( - \beta \int q_2(w_1) H(w_1, w_2) dw_1 \bigr) \propto \exp \bigl( \int q_2(w_1) \log p(w_1, w_2) dw_1 \bigr)$

が要請され、これを自己無矛盾条件とよぶと。それで自己無矛盾条件を満たす試行分布をみつけよということですが…一般にラグランジュ関数の停留点は最小点である保証はないのであくまで最小点であることの必要条件になると思うんですが、この場合も極大点であることはないんでしょうか。そもそも必ず最小点があるんでしょうか。テキストの雰囲気的には「自己無矛盾条件を満たす試行分布の中には必ず平均場近似はある」といった感じですが。

あと、平均場自由エネルギー $\beta \hat{F}(\beta)$ も何がエネルギーなのかわかりませんね。いうほど通常の自由エネルギーも何のことかわかりませんが。

自由エネルギー $\beta F(\beta) = -\log Z(\beta)$ って、確率分布が $p(w_1, w_2)= \exp \bigl( -\beta H(w_1, w_2) \bigr) / Z(\beta)$ の形でかかれていることを前提としていて、このとき点 $(w_1, w_2)$ の選択情報量（本当は $p(w_1, w_2)$ は確率じゃなくて確率密度だから選択情報量とはいわないけど）は $\log \bigl( 1/p(w_1, w_2) \bigr)$ だけど、「もし正規化を忘れて選択情報量を $\log \left( 1/\exp \bigl( -\beta H(w_1, w_2) \bigr) \right)$ と考えてしまうと $\beta F(\beta)$ だけ底上げされているよ！ $\beta F(\beta)$ だけ差し引く必要があるよ！」という量なのかな。 $w_1, w_2$ によらずいつも $\log \left( 1/\exp \bigl( -\beta H(w_1, w_2) \bigr) \right)$ から $\beta F(\beta)$ が差し引かれているよね。

$\displaystyle \log \left( \frac{1}{p(w_1, w_2)} \right) = \log \left( \frac{1}{\exp \bigl( -\beta H(w_1, w_2) \bigr) } \right) + \log Z(\beta)= \log \left( \frac{1}{\exp \bigl( -\beta H(w_1, w_2) \bigr) } \right) - \beta F(\beta)$

選択情報量を「その出来事を知らされたときの驚きの量」みたいに捉えるなら、「必ず起きることを知らされたときの驚きの量をゼロにするために差し引くべきオフセット」が自由エネルギーだといえるかもしれない。じゃあ確率分布 $q_1(w_1)q_2(w_2)$ によるエンコードで同じ点 $(w_1, w_2)$ の選択情報量を測ったものはというと、 $\log \left( 1/\exp \bigl( -\beta H(w_1, w_2) \bigr) \right)$ に比べて $\log \left( 1/\exp \bigl( -\beta H(w_1, w_2) \bigr) \right) + \log q_1(w_1)q_2(w_2)$ だけ差し引かれたものになっているよね。トートロジカルだけど。

$\displaystyle \log \left( \frac{1}{q_1(w_1)q_2(w_2)} \right) = \log \left( \frac{1}{\exp \bigl( -\beta H(w_1, w_2) \bigr) } \right) - \log \left( \frac{1}{\exp \bigl( -\beta H(w_1, w_2) \bigr) } \right) - \log q_1(w_1)q_2(w_2)$

この差し引くオフセットは $w_1, w_2$ に依存しているけど、確率分布 $q_1(w_1)q_2(w_2)$ で平均すると平均場自由エネルギー $\beta \hat{F}(\beta)$ になるね。だからやっぱりこれも「驚きの量から差し引くべきオフセットの平均」になっている。このときの驚きの量のエンコーディングは最適ではないかもしれないけど。
$\beta \hat{F}(\beta) = \int q_1(w_1)q_2(w_2) \log \Bigl( q_1(w_1)q_2(w_2)/ \exp \bigl( -\beta H(w_1, w_2) \bigr) \Bigr) dw_1 dw_2$
$\quad \displaystyle =\int q_1(w_1)q_2(w_2) \left( \log \left( \frac{ 1}{ \exp \bigl( -\beta H(w_1, w_2) \bigr)} \right) + \log q_1(w_1)q_2(w_2) \right) dw_1 dw_2$

ただ「平均場自由エネルギー」を「驚きの量から差し引くべき平均的なオフセット」と捉えても、「なぜ平均場自由エネルギーを最小化しなければならないのか」に示唆は得られないな…。いまは「カルバック・ライブラー情報量」つまり「真の驚きの量からの平均的なムダ」を最小化したいから、それと定数 $\beta F(\beta)$ ずれているだけの「（仮の）驚きの量から差し引くべきオフセットの平均」を最小化したいはずだ。エンコードがずれているほど差し引くべきオフセットは平均的にかさむからね。

だからその説明で誰かに伝わりますかね…？まあそれで冒頭の推定をするんですが、つまり、パラメータ $a_{1:K}, \, b_{1:K}$ の分布 $\varphi(a_{1:K}, b_{1:K})$ を更新していくわけです。

パラメータ $a_{1:K}$ の分布はディリクレ分布 ${\rm Dir}(a_{1:K}|\phi)$ ですね。ハイパーパラメータ $\phi \in \mathbb{R}_{+}^K$ を更新していくことになります。ディリクレ分布は $K$ 面あるサイコロの各目が出る割合を生成する確率分布ですね。ある割合における確率密度は「各目が $\phi_k - 1$ 回ずつ出たときにその割合であった確率密度」といえます。
パラメータ $b_{1:K}$ の分布は $\prod_{k=1}^K \bigl(1 / z(\eta_k) \bigr) \exp \bigl( \eta_k \cdot f(b_k) \bigr)$ です。12ページでやった指数型分布の共役な事前分布ですね。ハイパーパラメータ $\eta \in \mathbb{R}^{J \times K}$ を更新していくことになると思います。

12ページの指数型分布のときにはもう事後分布を直接計算しにいきましたが、ここでは謎の確率変数 $Y$ を導入します。これは実質的には $1, \cdots, K$ のいずれかの値を取る確率変数で、元々の確率変数 $X$ と $Y$ の同時分布が混合分布の $k$ 番目の山のみになるようなものです。言い換えると、 $Y$ は $K$ 面サイコロではあるのですが、 $\mathbb{R}^N$ 空間内の点によってサイコロの各面の出る割合が変動します。その割合は、その点における混合分布の $1, \cdots, K$ 番目の山の密度の比になります。 $a_{1:K}, \, b_{1:K}$ から山たちが決まり、山たちから「各点における密度」「各点における各山の密度の比」が決まるので、 $x^n, y^n, a_{1:K}, b_{1:K}$ の同時分布はこうです。 $\displaystyle P(x^n, y^n, a_{1:K}, b_{1:K}) = \varphi(a_{1:K}, b_{1:K}) \prod_{i=1}^n p(x_i, y_i |a_{1:K}, b_{1:K})$ これをさらに観測データで条件付けた $P(y^n, a_{1:K}, b_{1:K} | x^n)$ を $q(y^n) r(a_{1:K}, b_{1:K})$ で近似したいのですが、 $P(y^n, a_{1:K}, b_{1:K} | x^n) \propto P(x^n, y^n, a_{1:K}, b_{1:K})$ なのでこのまま自己無矛盾条件に代入して差し支えありません。そうすると結局、 $q(y^n)$ は、各データの $K$ 面サイコロの目があるパターンのときだけ $1$ をとる確率分布になります（これは誤り）。 $r(a_{1:K}, b_{1:K})$ の方はハイパーパラメータを更新した形になります。それで「各データのサイコロの目のパターン」と「ハイパーパラメータ」を交互に更新していけばよかろうという算段のようなんです。

…この手続きで、「謎の確率変数 $Y$ 」と「交互の更新」がしれっと導入されましたが、そもそもこんなことをしてよいのでしょうか？「謎の確率変数 $Y$ 」は、いわば各データの実家が何番目の山かというようなものになっていますが、本当はどの山が実家とかないでしょう？そこには山の重なりのような分布があるだけで、各データがどの山からきたのかなんて人間の妄想に過ぎないはずなんです。それに、その「各データがどの山からきたのか」と「ハイパーパラメータ」を交互に更新するというのも、直感的には解に向かっていきそうですが、何の根拠があってそんなことをしているのか…。

2020-01-04 追記

「こんなことをしてもよいのか」どうかは、「何をするためにそうしたのか」によるだろう。今は何をしたいんだっけ。

何をしたいんだっけって、データから混合分布を推定したいんですよ。

推定するって何？

推定するって何って…この本では、4ページで定義した事後分布によって確率モデルを平均したものこそが推定分布でしたよね。だからこの定義に沿って事後分布を得ようとするのが筋のはずです。

$\displaystyle \begin{split} p(a_{1:K}, b_{1:K}|X^n) &= \frac{1}{Z_n(1)} \varphi(a_{1:K}, b_{1:K}) \prod_{i=1}^n p(X_i |a_{1:K}, b_{1:K} ) \\ &= \frac{1}{Z_n(1)} \varphi(a_{1:K}, b_{1:K}) \prod_{i=1}^n \sum_{k=1}^K a_k v(X_i) \exp \bigr(f (b_k) \cdot g(X_i) \bigl) \end{split}$ …山が1つのみのときは共役な事前分布がありましたが、山が複数あると山どうしのクロスタームが出てきてしまいますね…少なくとも山の個数だけ共役な事前分布を用意すればよいということにはなりません…。

まあ解析的に求まる例ではないんだよね。だからモンテカルロ法で解くとか、別の方法で解っぽい何かを探すことになる。5.2節の「平均場近似」はその「解っぽい何かを探す別の方法」だから、「4ページで定義した事後分布で確率モデルを得ること」はもはや放棄されている。「こんなことをしてもよいのか」というよりは「こうすることにした」んだね。…といっても、「だからどうしてそうすることにしたのか」の回答にはなっていないな。いま解析的には解けないけど「パラメータ $a_{1:K}, b_{1:K}$ の分布をよりよさそうなものに更新する」ことを達成したい。だから次善の策として平均場近似では別の分布でこう近似することにした。

本当の事後分布とのカルバック・ライブラー距離の極小点を目指すことにする。
極小点が満たすべき条件を更新式として利用できるようにするために、パラメータを2グループ（以上）に分けて、こちらのグループの分布を固定してあちらのグループの分布を更新する → 今度はさっき固定していたグループの分布を更新する → …という繰り返しができるようにする。この目的のため、パラメータの分布を独立な分布の積とする。

…実際に平均場近似が編み出された歴史的経緯はさておき、目的からブレークダウンするなら平均場近似というフレームワークは上のように解釈できるだろう。

なるほど、だから「謎の確率変数」が導入された…説明になっていないんですが。「交互の更新」の方にはそれで説明が与えられましたが。

うん、上のフレームワーク自体は「謎の確率変数の導入」を要請しないからね。ただ今回は、157ページの $\log P$ は「謎の確率変数」を導入しなくては「あるパラメータたちの分布を固定して別のパラメータたちの分布のハイパーパラメータを更新する」という形に持ち込めないはず。「謎の確率変数」を導入しないと結局山どうしのクロスタームが出ちゃうしね。

「謎の確率変数」の導入はそういう背景が…しかし、「じゃあこういう確率変数を導入しよう」というのも突飛な気がしますが…。

そんな突飛かな？「山どうしのクロスターム消えてほしい」→「各データの実家の山は一つまで」って思えば自然に導入できるんじゃない？

いうほど自然に導入できますかね？？それに、クロスタームを消してしまうということはそれだけ正しいやり方を損ねているでしょう？「各データの実家の山は一つまで」というのも結構大胆な仮定にみえるんですが、それが心配というか…。

もう少しつづく