jiku log

データサイエンスの核心を掴む : 学びと発見の記録

「確率的機械学習 入門編I」を読む ~第4章 統計学 ④ベイズ統計学・その3~

はじめに

持橋大地・鈴木大慈 監訳「確率的機械学習 入門編I」は,世界的に評価の高いK.P.Murphy著 "Probabilistic Machine Learning (Book1)" の和訳であり,確率モデルに基づく機械学習,深層学習といった基礎が丁寧に整理されている。私は統計検定1級として数理統計の基礎は学んできたが,機械学習・深層学習は実務に応じて場当たり的に学んできた。実務での応用に向けて,機械学習・深層学習の基礎を体系的に学び,チームの技術力を底上げしたいと考えている。また読書メモに自身の理解をまとめることで,製造業に携わる若いエンジニアにとっても有益な知識を還元できればと考えている。

※なおボリュームが多い本なので,知っているところは端折りながら読み進めたい。


本記事は,「第4章 統計学」のベイズ統計学における,信用区間ベイズ機械学習,事後分布の近似に関する読書メモである。

4.6 ベイズ統計学

4.6.6 信用区間

中央区

本項では,推定の不確実性を定量化するための考え方である信用区間(credible interval)を説明している。

 100(1 - \alpha)%信用区間とは,連続な区間 C = (l, u)(下限と上限)であって,事後分布の質量のうち 1 - \alphaの割合を含むものである。


 \begin{align}
C_{\alpha}(\mathcal{D}) = (l, u) : P(l \leq \theta \leq u \mid \mathcal{D}) = 1 - \alpha \\ \\
\end{align}

これを満たす区間は多数存在しうるので,通常は両側の裾に (1-\alpha)/2の割合で質量をもつものを選び,これは中央区間と呼ばれる。

事後分布の関数形が既知であり,累積分布関数の逆関数 F^{-1}gが分かっている場合, l = F^{-1}(\alpha/2)および u = F^{-1}(1- \alpha/2)である。
ただしこの逆関数は一般には求めるのが難しい場合がある。この場合,モンテカルロ近似などで近似する方法がある。


中央区間の問題点は,分布が左右非対称のときに,確率密度が大きい領域が含まれなくなることがある点である。
下図は,中央区間を表している。この分布は右側の裾が長いため,左端の確率密度が高い領域が中央区間に含まれなくなってしまう,すなわち確率密度の値が大きいのに信用区間から外れてしまう領域が出現することである。

中央区

最大事後密度領域

この問題に対応するために,最大事後密度(highest posterior density, HPD)領域を考える。これは,事後確率がある閾値 p^*以上となる点を集めた領域で,その領域が指定された確率質量( 1 - \alpha)を最もコンパクトに含む領域のことを指す。

下図において,青線よりも上の領域に,確率密度関数の質量の 1 - \alphaが集中しているとする。

最大事後密度領域

すなわち青線 p^*の位置は,


 \begin{align}
1 - \alpha = \int_{\theta : p(\theta \mid \mathcal{D}) \gt p^*} p(\theta \mid \mathcal{D}) d\theta \\ \\
\end{align}
を満たす。この青線よりも内側の領域

 \begin{align}
C_{\alpha}(\mathcal{D}) = \{ \theta : p(\theta \mid \mathcal{D} \geq p^*) \} \\ \\
\end{align}
が最大事後密度領域である。

事後分布が多峰性を持つとき,最大事後密度領域は右図のように2つの領域に分かれることがある。

多峰性がある分布における中央区間(左)と最大事後密度領域(右)

4.6.7 ベイズ機械学習

事後分布 p(\boldsymbol{\theta} \mid \mathcal{D})を通じてモデルパラメータ \boldsymbol{\theta}ベイズ的に扱う枠組みを,ベイズ機械学習(Bayesian machine learning)と呼ぶ。

プラグイン近似

パラメータについて事後分布を計算した後,未知パラメータを周辺化することによって,入力が与えられた下での出力の事後予測分布が計算できる。


 \begin{align}
p(\boldsymbol{y} \mid \boldsymbol{x}, \mathcal{D}) = \int p(\boldsymbol{y} \mid \boldsymbol{x}, \boldsymbol{\theta}) p(\boldsymbol{\theta} \mid \mathcal{D}) d\boldsymbol{\theta} \\ \\
\end{align}

しかし,この積分計算は困難であることが多い。かわりに,MLEなどで点推定したパラメータ \hat{\boldsymbol{\theta}}を用いたデルタ関数によって,


 \begin{align}
 p(\boldsymbol{\theta} \mid \mathcal{D})  = \delta(\boldsymbol{\theta} - \hat{\boldsymbol{\theta}}) \\ \\
\end{align}
で,1点に集中したような事後分布を考えると

 \begin{align}
p(\boldsymbol{y} \mid \boldsymbol{x}, \mathcal{D}) = \int p(\boldsymbol{y} \mid \boldsymbol{x}, \boldsymbol{\theta}) p(\boldsymbol{\theta} \mid \mathcal{D}) d\boldsymbol{\theta} 
= p(\boldsymbol{y} \mid \boldsymbol{x}, \hat{\boldsymbol{\theta}}) \\ \\
\end{align}
のように,単に点推定の値を代入(プラグイン)したような事後予測分布が得られる。

プラグイン近似は過剰適合のおそれがある。そのため,いくつかのパラメータをモンテカルロ法で出力し,これらを用いて事後分布を近似する方法がある。

4.6.8 計算上の問題

事後分布の計算は,一般に困難なことが多い。本項では,近似事後分布推論(approximate posterior inference)の方法を説明している。

格子点近似

事後分布推論を近似するための最も単純なアプローチは,未知量 \boldsymbol{\theta}が取り得る値の空間を有限の候補 \boldsymbol{\theta}_1, ..., \boldsymbol{\theta}_Kで集合に分割し,総当たりで事後分布を近似する方法である。
しかしこの方法は,パラメータの次元数に対して指数的に候補が増加してしまうため,2,3次元以上になる場合は使うことができない。

二次(ラプラス)近似

ラプラス近似は,モード \hat{\boldsymbol{\theta}}周りでテイラー展開して2次近似する方法であり,多変量のガウス分布で近似する方法である。

事後分布を以下のように変形する。


 \begin{align}
p(\boldsymbol{\theta} \mid \mathcal{D}) &= \frac{p(\boldsymbol{\theta} , \mathcal{D}) }{ p(\mathcal{D} )}  = \frac{1}{Z} e^{-\mathcal{E}(\boldsymbol{\theta}) }   \\ \\
\mathcal{E}(\boldsymbol{\theta}) &= -\log p(\boldsymbol{\theta}, \mathcal{D}), \quad Z = p(\mathcal{D}) \\ \\
\end{align}
ここで, \mathcal{E}(\boldsymbol{\theta})はエネルギー関数と呼ばれる。

モード \hat{\boldsymbol{\theta}}(すなわち,エネルギーが最も低い状態)の周りでテイラー展開すると,


 \begin{align}
\mathcal{E}(\boldsymbol{\theta}) \approx \mathcal{E}(\hat{\boldsymbol{\theta}} ) + (\boldsymbol{\theta}  - \hat{\boldsymbol{\theta}} )^T \boldsymbol{g}
 + \frac{1}{2} (\boldsymbol{\theta} - \hat{\boldsymbol{\theta}} )^T \boldsymbol{H} (\boldsymbol{\theta} - \hat{\boldsymbol{\theta}} )  \\ \\
\end{align}
ただし,  \boldsymbol{g},  \boldsymbol{H}はそれぞれエネルギー関数の勾配とヘッセ行列である。

  \hat{\boldsymbol{\theta}}はモードであるため,1次の項は0になる。よって,


 \begin{align}
\hat{p}(\boldsymbol{\theta} \mid \mathcal{D}) = 
\frac{1}{Z} e^{-\mathcal{E}(\boldsymbol{\theta}) } \left[ - \frac{1}{2} (\boldsymbol{\theta} - \hat{\boldsymbol{\theta}} )^T \boldsymbol{H} (\boldsymbol{\theta} - \hat{\boldsymbol{\theta}} ) \right] \\ \\
\end{align}
となる。

この式には, \exp()の中に2次形式が出てくるので,多変量ガウス分布の形になる。よって,


 \begin{align}
\hat{p}(\boldsymbol{\theta} \mid \mathcal{D}) = \mathcal{N}(\boldsymbol{\theta} \mid \hat{\boldsymbol{\theta}}, \boldsymbol{H}) \\ \\
\end{align}
が得られる。

まとめと感想

今回は,「4 統計学」のベイズ統計学における,信用区間ベイズ機械学習,事後分布の近似についてまとめた。

信用区間は,推定されたパラメータの不確実性を評価するための方法であるが,中央区間・最大事後密度領域の特徴について把握することができた。
ベイズ統計を用いると,推定されたパラメータの分布を評価することができるが,実務上は得られる事後分布が対称な形になっている保証がないので,最大事後密度領域を用いる方が妥当だと感じた。

プラグイン近似は,実務ではよく用いる方法であるが,過剰適合の可能性があるとのことである。その対策としてモンテカルロ近似が紹介されていたが,計算自体も比較的容易なので,実務でも試してみたい。

ラプラス近似は,複雑な事後分布を多変量ガウス分布で近似するという方法であるが,多変量ガウス分布であれば扱いやすいうえに理解もしやすいので,実務上も有用だと感じた。


本記事を最後まで読んでくださり,どうもありがとうございました。