勉強 数学 機械学習/AI

ベイジアンニューラルネットワーク

 

 

ニューラルネットワークにベイズ推定の考え方を導入し応用したものに、ベイジアンニューラルネットワークがある。

通常のニューラルネットワークとの違いやその理論についてイメージをつかむべく整理していく。

 

ベイジアンニューラルネットワーク（BNN）とは

ベイジアンニューラルネットワーク(Bayesian Neural Network)は、近年注目を集めている一つの研究テーマである。この技術は、従来のニューラルネットワークにベイジアン理論を組み合わせることで、予測の不確実性を量化できるという特長がある。そのため、自動運転車、医療診断、金融取引など、確率的な要素が含まれる多くの問題に対して有用である。

普通のニューラルネットワークの弱点

NNのイメージを伝えるために、便宜的に収束したパラメータを定数で表している。

 

通常のニューラルネットワークは各層の重みパラメータの値は学習によって収束した定数として扱われ、そしてその予測値も定数として出力される。

さらには、一般にニューラルネットワークモデルを使った予測はどうやっても100パーセントの精度を達成できない。これはどうしても予測結果に対して不確実な要素が含まれてしまうからで、その不確実なものにしてしまう要素というのは例えば、データのノイズやそもそもの予測能力の不足などがある。

ニューラルネットワークの中身はブラックボックスとも言われており、どのような判断基準で予測しているのかわからないし、その予測結果をどの程度信頼できるのかもわからない。

なので、その出力が自信満々に普通に間違えている可能性がある。

 

通常のニューラルネットワークの場合この不確実な部分を分解して考えることができない。

ベイジアンニューラルネットワーク（BNN）の利点

 

ベイジアンニューラルネットワークでは、重みパラメータをある確率分布に基づいて生成される確率変数と考える。

そして、推論時においてモデルに入力が与えられるたびに、確率分布\(P(w)\)にサンプリングされた値が演算されて重みパラメータが更新されていく。

 

ベイジアンニューラルネットワークは、重みパラメータを確立変数とみなすことで重みパラメータの生成に関する確率分布を定義でき、予測値以外の統計値情報を出力できたり、重みパラメータ自体に情報量を定義することが可能となったりする。これにより、ニューラルネットワークの出力に対する不確実性を計測する様々な情報を手に入れることができるようになる。

 

さらに、ベイジアンな手法を採用することによって事前分布を使用していることから、データが足りない場合でもうまく動いたり、過学習を予測できたりする。

主なメリット、さらには実用例をまとめると

ベイジアンニューラルネットワークの主なメリット

  • 過学習の防止策としてデータが不十分な場合でも問題解決に役立つ（データ取得が高コストで、困難な実験研究を要する分野などでも使える）
  • BNNは通常のNNよりも普遍的に使用可能
  • 不確実性を定量的に扱える
  • 逐次的に学習できる

実用例・使いどころ

  • 意思決定プロセス内で意図しない行動を防ぐカギとなる可能性がある
  • 医療の画像診断
  • 分子生物学の画像認識
  • 車の自動運転
  • 時系列データの異常検知
  • 能動学習
  • 転移学習
  • 強化学習

などなど…

 

ベイズニューラルネットワーク（BNN)の問題点

利点ばかりではなく、問題点も勿論あって

ベイズ線形回帰の場合は、事後分布がガウス分布になるので厳密に評価できて、予測分布も閉形式で表すことができる。

しかし、多層ニューラルネットワークの場合では、ネットワーク関数がパラメータに極度に非線形に依存しているので、事後分布を厳密に計算するのは難しい。つまりベイズの枠組みで扱うことはできない。実際、事後分布の対数は非凸であり、誤差関数には複数の局所的極小点が存在する。

 

ではどうするのかというと事後分布を近似するしかないのである。

変分推論法とラプラス近似

変分推論法とラプラス近似はともに事後分布を近似する手法である。

変分推論法は、因子分解されたガウス近似を事後分布に用い、かつ完全な共分散を持つガウス分布を用いることで、ベイズニューラルネットワークに適用する方法である。

ラプラス近似は、最も確率が高い点（モード）を中心にガウス分布で事後分布を近似する。

ここからはラプラス近似の具体的な手法について解説していく。

 

ベイジアンニューラルネットワーク（回帰）における事後分布の計算

まず、ここから基本として考える回帰問題は、ある入力ベクトル\(\mathbf{x}\)から一次元の目標変数\(t\)を予測する問題だ。(簡単に多次元目標に拡張可能）ここでは、\(t\)がどのように\(\mathbf{x}\)から生成されるかという確率モデルを考える。具体的には、この条件付き確率\(p(t \mid \mathbf{x}, \mathbf{w}, \beta)\)はガウス分布で表され、その平均はニューラルネットワークの出力\(y(\mathbf{x}, \mathbf{w})\)で、\(\beta\)は精度（分散の逆数）

$$
p(t \mid \mathbf{x}, \mathbf{w}, \beta)=\mathcal{N}\left(t \mid y(\mathbf{x}, \mathbf{w}), \beta^{-1}\right)
$$

となる。

次に、このモデルのパラメータ（重み\(\mathbf{w}\)）に関する事前分布を考える。これもガウス分布を用い、平均を0、精度（分散の逆数）を\(\alpha\)とする。

$$
p(\mathbf{w} \mid \alpha)=\mathcal{N}\left(\mathbf{w} \mid \mathbf{0}, \alpha^{-1} \mathbf{I}\right)
$$

さて、このような設定のもとで、複数のデータ点\(\mathbf{x}_1, \ldots, \mathbf{x}_N\)とその目標値\(t_1, \ldots, t_N\)が与えられた場合、尤度関数は次のようになる。

$$
p(\mathcal{D} \mid \mathbf{w}, \beta)=\prod_{n=1}^N \mathcal{N}\left(t_n \mid y\left(\mathbf{x}_n, \mathbf{w}\right), \beta^{-1}\right)
$$

最後に尤度関数と事前分布が定義できたので、事後分布を計算すると

$$
p(\mathbf{w} \mid \mathcal{D}, \alpha, \beta) \propto p(\mathbf{w} \mid \alpha) p(\mathcal{D} \mid \mathbf{w}, \boldsymbol{\beta})
$$

ここまではいつも通りのベイズの流れ

先ほどの問題点としても言ったが、\(y(\mathbf{x}, \mathbf{w})\) が \(\mathbf{w}\)に非線形に依存しているからこの事後分布はガウス分布にならない。これは問題で、事前分布がガウス分布であるのであるならば、事後分布もガウス分布でないとベイズの枠組みで使いづらい。ゆえに、ラプラス近似を使用する必要があるのであった。

 

ラプラス近似を使ってパラメータの事後分布を計算

いきなり計算を始める前にまずは大まかな流れを確認しておこう

ラプラス近似の流れ

  1. 事後分布の最大値を探す: まず、事後確率が最も高くなるようなパラメータ（重みとバイアス）を見つける。これがガウス分布の「中心」になる。
  2. ヘッセ行列の近似:二乗和誤差関数部分の\(\mathbf{w}\)の各要素によるヘッセ行列をアルゴリズムで計算し近似
  3. ガウス分布で近似: 事後分布をこのガウス分布で近似する。
  4. 予測分布の計算: この近似を用いて、新しいデータに対する予測を行う。

確率分布をガウス分布で近似するラプラス近似のイメージ

ここからは一旦\(\alpha\),\(\beta\)は固定する

事後分布の最大値を探す

先ほど求めた事後分布は\(\mathcal{D}\)が与えられた後の、パラメータ\(\mathbf{w}\)の確率分布を表している。この事後分布が最も高くなる\(\mathbf{w}\)を求めるのが第一のステップであり、この\(\mathbf{w}\)を\(\mathbf{W}_{\text {MAP }}\)と呼ぶ。具体的には、事後分布の対数を取った式

$$
\ln p(\mathbf{w} \mid \mathcal{D})=-\frac{\alpha}{2} \mathbf{w}^{\mathrm{T}} \mathbf{w}-\frac{\beta}{2} \sum_{n=1}^N\left\{y\left(\mathbf{x}_n, \mathbf{w}\right)-t_n\right\}^2+\text { const }
$$

この式を反復的数値最適法で最大値を求めるのが便利であり、例えば「共役勾配法」というアルゴリズムを使う。ここで必要な微分の計算には誤差逆伝播が使われる。

勾配最適化（共役勾配法）について

ニューラルネットワーク１

                  ニューラルネットワークとは ニューラルネットワークは脳の神経細胞（ニューロ ...

続きを見る

対数事後分布を最大化して得られる\(\mathbf{W}_{\text {MAP }}\)は\(\mathbf{w}\)の初期値に依存する。

 

ヘッセ行列の近似

モード\(\mathbf{W}_{\text {MAP }}\)を見つけたら、事後分布の負の対数尤度の二階微分の行列を評価して、局所的にガウス分布を近似できる。これは先ほどの事後分布を直接計算して

$$
\mathbf{A}=-\nabla \nabla \ln p(\mathbf{w} \mid \mathcal{D}, \alpha, \beta)=\alpha \mathbf{I}+\beta \mathbf{H}
$$

で与えられる。ここでの\(\mathbf{H}\)は二乗和誤差関数の\(\mathbf{w}\)の各要素（要素数\(M\))による二階微分からなるヘッセ行列（\(M × M\)）を計算し近似するアルゴリズムはここでは割愛する。

 

ガウス分布で近似

ここがラプラス近似の核になる部分だが、ヘッセ行列が近似できたならガウス分布に対する標準的な結果を使って

$$
q(\mathbf{w} \mid \mathcal{D})=\frac{|\mathbf{A}|^{1 / 2}}{(2 \pi)^{M / 2}} \exp \left\{-\frac{1}{2}\left(\mathbf{w}-\mathbf{w}_{\text {MAP }}\right)^{\mathrm{T}} \mathbf{A}\left(\mathbf{w}-\mathbf{w}_{\text {MAP }}\right)\right\}=\mathcal{N}\left(\mathbf{w} \mid \mathbf{w}_{\text {MAP }}, \mathbf{A}^{-1}\right)
$$

つまり

$$
p(\mathbf{w} \mid \mathcal{D}, \alpha, \beta) \simeq q(\mathbf{w} \mid \mathcal{D})=\mathcal{N}\left(\mathbf{w} \mid \mathbf{w}_{\mathrm{MAP}}, \mathbf{A}^{-1}\right)
$$

事後分布を近似したガウス分布を得ることができる。

予測分布の計算

最後に、この近似した事後分布を使って、新しいデータ点\(\mathbf{x}\)に対する目標変数\(t\)の予測分布を計算する。具体的には、

$$
p(t \mid \mathbf{x}, \mathcal{D})=\int p(t \mid \mathbf{x}, \mathbf{w}) q(\mathbf{w} \mid \mathcal{D}) \mathrm{d} \mathbf{w}
$$

という形で求められる。

 

ネットワーク関数の線形近似

しかしながら、事後分布をガウス分布で近似したとしても、ネットワーク関数\(y(\mathbf{x}, \mathbf{w})\)が\(\mathbf{w}\)の非線形関数なので、結局先ほどの予測分布は解析的には扱いづらいのだ。

というわけで、次にモデルの出力関数（ここでは\(y(\mathbf{x}, \mathbf{w})\)）も線形近似する必要がある。そのためにテイラー展開を用いる。

ネットワーク関数のテイラー展開

具体的にはネットワーク関数の\(\mathbf{w}_{\text {MAP }}\)のまわりでのテイラー展開を考えてあげることで

$$
y(\mathbf{x}, \mathbf{w}) \simeq y\left(\mathbf{x}, \mathbf{w}_{\mathrm{MAP}}\right)+\mathbf{g}^{\mathbf{T}}\left(\mathbf{w}-\mathbf{w}_{\mathrm{MAP}}\right)
$$

となる。ここで\(\mathbf{w}_{\mathrm{MAP}}\)は事後分布の最大値を取る\(\mathbf{w}\)の値であり、\(\mathbf{g}\)は勾配（傾き）である。つまり、\(\mathbf{g}=\left.\nabla_{\mathbf{w}} y(\mathbf{x}, \mathbf{w})\right|_{\mathbf{w}=\mathbf{w}_{\mathrm{MAP}}}\)

 

これで条件付き確率

$$
p(t \mid \mathbf{x}, \mathbf{w}, \beta)=\mathcal{N}\left(t \mid y(\mathbf{x}, \mathbf{w}), \beta^{-1}\right)
$$

であったものが

$$
p(t \mid \mathbf{x}, \mathbf{w}, \beta) \simeq \mathcal{N}\left(t \mid y\left(\mathbf{x}, \mathbf{w}_{\mathrm{MAP}}\right)+\mathbf{g}^{\mathbf{T}}\left(\mathbf{w}-\mathbf{w}_{\mathrm{MAP}}\right), \beta^{-1}\right)
$$

と近似して表すことができた。

 

パーツがそろったので周辺化して予測分布を得よう

整理すると

得られたパーツ

重みパラメータ\(\mathbf{w}\)に関する事後分布（近似済み）

$$
p(\mathbf{w} \mid \mathcal{D}, \alpha, \beta) \simeq q(\mathbf{w} \mid \mathcal{D})=\mathcal{N}\left(\mathbf{w} \mid \mathbf{w}_{\mathrm{MAP}}, \mathbf{A}^{-1}\right) \tag{1}
$$

条件付き確率\(p(t \mid \mathbf{x}, \mathbf{w}, \beta) \)（近似済み）

$$
p(t \mid \mathbf{x}, \mathbf{w}, \beta) \simeq \mathcal{N}\left(t \mid y\left(\mathbf{x},　 \mathbf{w}_{\mathrm{MAP}}\right)+\mathbf{g}^{\mathbf{T}}\left(\mathbf{w}-\mathbf{w}_{\mathrm{MAP}}\right), \beta^{-1}\right) \tag{2}
$$

ここから周辺分布

$$
p(t \mid \mathbf{x}, \mathcal{D}, \alpha, \beta)=\mathcal{N}\left(t \mid y\left(\mathbf{x}, \mathbf{w}_{\mathrm{MAP}}\right), \sigma^2(\mathbf{x})\right)
$$

を得ることができる。ここで

$$
\sigma^2(\mathbf{x})=\beta^{-1}+\mathbf{g}^{\mathrm{T}} \mathbf{A}^{-1} \mathbf{g} .
$$

証明

この一般に成り立つ結果を使用する

ガウス分布の周辺分布と条件付き分布

\(\mathbf{x}\)の周辺ガウス分布と、\(\mathbf{x}\)が与えられたときの\(\mathbf{y}\)の条件付きガウス分布が次式で与えられたとする。

$$
\begin{aligned}
p(\mathbf{x}) & =\mathcal{N}\left(\mathbf{x} \mid \boldsymbol{\mu}, \mathbf{\Lambda}^{-1}\right)  \\
p(\mathbf{y} \mid \mathbf{x}) & =\mathcal{N}\left(\mathbf{y} \mid \mathbf{A x}+\mathbf{b}, \mathbf{L}^{-1}\right)
\end{aligned}
$$

(1)の式が一つ目の式(2)の式が二つ目の式に相当すると考える。

\(\mathbf{y}\)の周辺分布と、\(\mathbf{y}\)が与えられたときの\(\mathbf{x}\)の条件付き分布は

$$
\begin{aligned}
p(\mathbf{y}) & =\mathcal{N}\left(\mathbf{y} \mid \mathbf{A} \boldsymbol{\mu}+\mathbf{b}, \mathbf{L}^{-1}+\mathbf{A} \mathbf{\Lambda}^{-1} \mathbf{A}^{\mathrm{T}}\right)  \\
p(\mathbf{x} \mid \mathbf{y}) & =\mathcal{N}\left(\mathbf{x} \mid \mathbf{\Sigma}\left\{\mathbf{A}^{\mathrm{T}} \mathbf{L}(\mathbf{y}-\mathbf{b})+\boldsymbol{\Lambda} \boldsymbol{\mu}\right\}, \mathbf{\Sigma}\right)
\end{aligned}
$$

で与えられる。ただし、

$$
\boldsymbol{\Sigma}=\left(\boldsymbol{\Lambda}+\mathbf{A}^{\mathrm{T}} \mathbf{L} \mathbf{A}\right)^{-1}
$$

とする。

これを使って、変数の対応に注意しながら周辺化を行う。

$$
\begin{gathered}
p(t \mid \mathbf{x}, \mathcal{D}, \alpha, \beta)=\mathcal{N}\left(t \mid y\left(\mathbf{x}, \mathbf{w}_{\mathrm{MAP}}\right), \sigma^2(\mathbf{x})\right) \\
\sigma^2(\mathbf{x})=\beta^{-1}+\mathbf{g}^{\mathrm{T}} \mathbf{A}^{-1} \mathbf{g}
\end{gathered}
$$

示せた□

 

これによって、予測分布\(p(t \mid \mathbf{x}, \mathcal{D})\)は平均がネットワーク関数\(y\left(\mathbf{x}, \mathbf{w}_{\mathrm{MAP}}\right)\)で与えられるガウス分布になることがわかる。

 

ハイパーパラメータの最適化

先ほどまでは、ハイパーパラメータ\(\alpha\)と\(\beta\)は固定され既知として考えてきた。

ここからは、得られた予測分布にさらにハイパーパラメータの不確かさも考慮に入れたものを考えよう。

 

こういう時は、線形回帰モデルの時にやった経験ベイズのエビデンス理論を使用すれば、実用的に手続き的にハイパーパラメータの値を選ぶことができるのであった。

線形回帰モデルの経験ベイズ

経験ベイズってなに？線形回帰モデルをベイズ的に扱うための近似とエビデンス関数

  エビデンス近似 僕たちが統計的な問題を解くとき、大抵は何かを予測したり、何かのパターンを見つけたりしたいわけで、これを達成するために、僕たちは通常、データに適合するようなモデルを構築する ...

続きを見る

さっとわかる経験ベイズの流れ

  1. 周辺尤度関数\(p(\mathcal{D} \mid \alpha, \beta)\)（エビデンス関数）を同時分布を重みパラメータ\(\mathbf{w}\)に関して積分して得る
  2. 対数エビデンス関数を得る。
  3. 固有値の計算
  4. \(\alpha\)に関する対数エビデンス関数の最大化
  5. \(\beta\)に関する対数エビデンス関数の最大化

エビデンス関数を得る

ここでのエビデンス関数はどうなるのかというと

$$
p(\mathcal{D} \mid \alpha, \beta)=\int p(\mathcal{D} \mid \mathbf{w}, \beta) p(\mathbf{w} \mid \alpha) \mathrm{d} \mathbf{w}
$$

対数エビデンス関数を得る

となる。これを経験ベイズの手続きに従って対数化すると

対数エビデンス関数

$$
\ln p(\mathcal{D} \mid \alpha, \beta) \simeq-E\left(\mathbf{w}_{\mathrm{MAP}}\right)-\frac{1}{2} \ln |\mathbf{A}|+\frac{W}{2} \ln \alpha+\frac{N}{2} \ln \beta-\frac{N}{2} \ln (2 \pi)
$$

となる。ここで\(W\)は\(\mathbf{w}\)に含まれるパラメータの総数であり、正則化誤差関数は

$$
E\left(\mathbf{w}_{\mathrm{MAP}}\right)=\frac{\beta}{2} \sum_{n=1}^N\left\{y\left(\mathbf{x}_n, \mathbf{w}_{\mathrm{MAP}}\right)-t_n\right\}^2+\frac{\alpha}{2} \mathbf{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathbf{w}_{\mathrm{MAP}} .
$$

と定義されている。

証明

$$
p(\mathcal{D} \mid \alpha, \beta)=\int p(\mathcal{D} \mid \mathbf{w}, \beta) p(\mathbf{w} \mid \alpha) \mathrm{d} \mathbf{w}
$$

という分布をラプラス近似するには、以下を使用する

正規化係数のラプラス近似

$$
\begin{aligned}
Z & =\int f(\mathbf{z}) \mathrm{d} \mathbf{z} \\
& \simeq f\left(\mathbf{z}_0\right) \int \exp \left\{-\frac{1}{2}\left(\mathbf{z}-\mathbf{z}_0\right)^{\mathrm{T}} \mathbf{A}\left(\mathbf{z}-\mathbf{z}_0\right)\right\} \mathrm{d} \mathbf{z} \\
& =f\left(\mathbf{z}_0\right) \frac{(2 \pi)^{M / 2}}{|\mathbf{A}|^{1 / 2}}
\end{aligned}
$$

これは、ある分布\(f(\mathbf{z})\)とそのモード\(\mathbf{z}_0\)、そしてそのモード周りのヘッセ行列\(\mathbf{A}\)を使用して、その正規化係数\(Z\)を近似するもの。

つまり、\(\mathrm{f}(\mathrm{w})=\mathrm{p}(\mathrm{D} \mid \mathrm{w}, \beta) \mathrm{p}(\mathrm{w} \mid \alpha), \mathrm{Z}=\mathrm{p}(\mathrm{D} \mid \alpha, \beta)\)として、上の近似を適用すると、（\(\mathbf{W}\)は\(\mathbf{w}\)の次元数）

$$
\begin{array}{r}
\mathrm{p}(\mathrm{w} \mid \alpha)=\mathrm{N}\left(\mathrm{w} \mid 0, \alpha^{-1} \mathrm{I}\right) \\
\mathrm{p}(\mathrm{D} \mid \mathrm{w}, \beta)=\prod_{\mathrm{n}=1}^{\mathrm{N}} \mathrm{N}\left(\mathrm{t}_{\mathrm{n}} \mid \mathrm{y}\left(\mathrm{x}_{\mathrm{n}}, \mathrm{w}\right), \beta^{-1}\right)
\end{array}
$$

だったので、

$$
\begin{aligned}
\mathrm{Z} \simeq \mathrm{f}\left(\mathrm{w}_{\mathrm{MAP}}\right) & =\mathrm{p}\left(\mathrm{D} \mid \mathrm{w}_{\mathrm{MAP}}, \beta\right) \mathrm{p}\left(\mathrm{w}_{\mathrm{MAP}} \mid \alpha\right) \\
& =\prod_{\mathrm{n}=1}^{\mathrm{NN}}\left(\mathrm{t}_{\mathrm{n}} \mid \mathrm{y}\left(\mathrm{x}_{\mathrm{n}}, \mathrm{w}_{\mathrm{MAP}}\right), \beta^{-1}\right) \mathrm{N}\left(\mathrm{w}_{\mathrm{MAP}} \mid 0, \alpha^{-1} \mathrm{I}\right) \\
& =\prod_{\mathrm{n}=1}^{\mathrm{N}} \frac{1}{(2 \pi)^{1 / 2}} \frac{1}{\left(\beta^{-1}\right)^{1 / 2}} \exp \left[-\frac{1}{2 \beta^{-1}}\left\{\mathrm{t}_{\mathrm{n}}-\mathrm{y}\left(\mathrm{x}_{\mathrm{n}}, \mathrm{w}_{\mathrm{MAP}}\right)\right\}^2\right] \\
& \qquad \qquad ×\frac{1}{(2 \pi)^{\mathrm{w} / 2}} \frac{1}{\left|\alpha^{-1} \mathrm{I}\right|^{1 / 2}} \exp \{-\frac{1}{2} \mathrm{w}_{\mathrm{MAP}}^{\mathrm{T}}\left(\alpha^{-1} \mathrm{I}\right)^{-1} \mathrm{w}_{\mathrm{MAP}}\} \\
& =\prod_{\mathrm{n}=1}^{\mathrm{N}}\left(\frac{\beta}{2 \pi}\right)^{1 / 2} \exp \left[-\frac{\beta}{2}\left\{\mathrm{t}_{\mathrm{n}}-\mathrm{y}\left(\mathrm{x}_{\mathrm{n}}, \mathrm{w}_{\mathrm{MAP}}\right)\right\}^2\right]\left(\frac{\alpha}{2 \pi}\right) \exp \left(-\frac{\alpha}{2} \mathrm{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathrm{w}_{\mathrm{MAP}}\right)
\end{aligned}
$$

よって、

$$
\begin{aligned}
\ln (\mathrm{D} \mid \alpha, \beta) & \simeq \operatorname{lnf}\left(\mathrm{w}_{\mathrm{MAP}}\right)+\frac{\mathrm{W}}{2} \ln (2 \pi)-\frac{1}{2} \ln |\mathrm{A}| \\
& =\sum_{\mathrm{n}=1}^{\mathrm{N}}\left[\frac{1}{2}\{\ln \beta-\ln (2 \pi)\}-\frac{\beta}{2}\left\{\mathrm{t}_{\mathrm{n}}-\mathrm{y}\left(\mathrm{x}_{\mathrm{n}}, \mathrm{w}_{\mathrm{MAP}}\right)\right\}^2\right]\\
& \qquad \qquad +\frac{\mathrm{W}}{2}\{\ln \alpha-\ln (2 \pi)\}-\frac{\alpha}{2} \mathrm{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathrm{w}_{\mathrm{MAP}} + \frac{1}{2}\ln|\mathrm{A}| \\
& =-\left[\frac{\beta}{2} \sum_{\mathrm{n}=1}^{\mathrm{N}}\left\{\mathrm{t}_{\mathrm{n}}-\mathrm{y}\left(\mathrm{x}_{\mathrm{n}}, \mathrm{w}_{\mathrm{MAP}}\right)\right\}^2+\frac{\alpha}{2} \mathrm{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathrm{w}_{\mathrm{MAP}}\right]\\
& \qquad \qquad-\frac{1}{2} \ln |\mathrm{A}|+\frac{\mathrm{N}}{2} \ln \beta-\frac{\mathrm{N}}{2} \ln (2 \\
& =-\mathrm{E}\left(\mathrm{w}_{\mathrm{MAP}}\right)-\frac{1}{2} \ln |\mathrm{A}|+\frac{\mathrm{N}}{2} \ln \beta-\frac{\mathrm{N}}{2} \ln (2 \pi)+\frac{\mathrm{W}}{2} \ln \alpha
\end{aligned}
$$

よって示せた□

ヘッセ行列の固有値を計算する

次の流れはヘッセ行列の固有値を計算することである。まず固有方程式

$$
\beta \mathbf{H} \mathbf{u}_i=\lambda_i \mathbf{u}_i
$$

を考える。

そして、有効パラメータ\(\gamma\)を求める

$$
\gamma=\sum_{i=1}^W \frac{\lambda_i}{\alpha+\lambda_i}
$$

最適な\(\alpha\)を求める

線形回帰の通りの計算で\(\alpha\)の再推定値

$$
\alpha=\frac{\gamma}{\mathbf{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathbf{w}_{\mathrm{MAP}}}
$$

を得ることができる。

ここに注意

この結果は線形回帰の時には厳密に求めることができるが、このニューラルネットワーク場合ではそうではない。

なぜならば、\(\alpha\)の変化によってヘッセ行列\(\mathbf{H}\)が変化し、それによりさらに固有値が変化するので値が変化してしまう。ここではそれを無視している。

最適な\(\beta\)を求める

これまた線形回帰の通りの計算で\(\beta\)をエビデンス対数関数について最大化すると、再推定した値は

$$
\frac{1}{\beta}=\frac{1}{N-\gamma} \sum_{n=1}^N\left\{y\left(\mathbf{x}_n, \mathbf{w}_{\mathrm{MAP}}\right)-t_n\right\}^2
$$

となる。

この手順はハイパーパラメータ\(\alpha\)と\(\beta\)を交互に更新を行わなければならない。

 

最適化を振り返る

 

最適化の流れ

最適化の流れ,NNは(Neural Network)

 

 

ベイジアンニューラルネットワークをクラス分類に拡張

これまでは、回帰モデルを扱ってきたが、最後にこれをクラス分類に拡張できることを確認して終わろう。

 

ここでは、１つのロジスティックシグモイド出力を持つネットワークにようる２クラス分類問題を考える。（簡単にソフトマックス出力の多クラス分類に拡張可能）

このときの、対数尤度関数は

クラス分類時の対数尤度関数

$$
\ln p(\mathcal{D} \mid \mathbf{w})=\sum_{n=1}^N\left\{t_n \ln y_n+\left(1-t_n\right) \ln \left(1-y_n\right)\right\}
$$

ここで、\(t_n \in \{0,1\}\)は目標値であり、\(y_n \equiv y\left(\mathbf{x}_n, \mathbf{w}\right)\)である。データ点のラベリングは正しいとしているので、ノイズのハイパーパラメータ\(\beta\)は現れない。事前分布は等方ガウス分布

パラメータの事後分布を計算

これが決まったのであれば、回帰の場合と同じく対数事後分布を最大化して、その時の重みパラメータ\(\mathbf{W}_{\text {MAP }}\)を見つける。

この時は負の対数尤度関数に、正則化項を加えた正則化誤差関数

$$
E(\mathbf{w})=-\ln p(\mathcal{D} \mid \mathbf{w})+\frac{\alpha}{2} \mathbf{w}^{\mathrm{T}} \mathbf{w}
$$

を最小化することと等価である。同じく最適化アルゴリズムと誤差逆伝播を組み合わせて見つけることができる。

 

\(\mathbf{W}_{\text {MAP }}\)を見つけることができたら、負の対数尤度関数の二回微分からなるヘッセ行列\(\mathbf{H}\)を評価してあげて、

$$
\mathbf{A}=\alpha \mathbf{I} + \mathbf{H}
$$

として、結局回帰の時と同じ形で、ガウス分布で近似された事後分布は

$$
q(\mathbf{w} \mid \mathcal{D})=\mathcal{N}\left(\mathbf{w} \mid \mathbf{w}_{\mathrm{MAP}}, \mathbf{A}^{-1}\right)
$$

という形で書ける。

要するに、回帰の時とおんなじ。

\(\alpha\)を最適化

ハイパーパラメータ\(\alpha\)を最適化するには、同じく対数エビデンス関数を最大化するという流れになり、その式は

クラス分類の場合の対数エビデンス関数

$$
\ln p(\mathcal{D} \mid \alpha) \simeq-E\left(\mathbf{w}_{\mathrm{MAP}}\right)-\frac{1}{2} \ln |\mathbf{A}|+\frac{W}{2} \ln \alpha+\text { const }
$$

同じく\(W\)は\(\mathbf{w}\)に含まれるパラメータの総数であり、正則化誤差関数は

$$
E\left(\mathbf{w}_{\mathrm{MAP}}\right)=-\sum_{n=1}^N\left\{t_n \ln y_n+\left(1-t_n\right) \ln \left(1-y_n\right)\right\}+\frac{\alpha}{2} \mathbf{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathbf{w}_{\mathrm{MAP}}
$$

である。

証明

$$
p(\mathcal{D} \mid \alpha, \beta)=\int p(\mathcal{D} \mid \mathbf{w}, \beta) p(\mathbf{w} \mid \alpha) \mathrm{d} \mathbf{w}
$$

回帰の時と同じく、という分布をラプラス近似するには、以下を使用する。

正規化係数のラプラス近似

$$
\begin{aligned}
Z & =\int f(\mathbf{z}) \mathrm{d} \mathbf{z} \\
& \simeq f\left(\mathbf{z}_0\right) \int \exp \left\{-\frac{1}{2}\left(\mathbf{z}-\mathbf{z}_0\right)^{\mathrm{T}} \mathbf{A}\left(\mathbf{z}-\mathbf{z}_0\right)\right\} \mathrm{d} \mathbf{z} \\
& =f\left(\mathbf{z}_0\right) \frac{(2 \pi)^{M / 2}}{|\mathbf{A}|^{1 / 2}}
\end{aligned}
$$

これは、ある分布\(f(\mathbf{z})\)とそのモード\(\mathbf{z}_0\)、そしてそのモード周りのヘッセ行列\(\mathbf{A}\)を使用して、その正規化係数\(Z\)を近似するもの。

ここで\(f(\mathbf{w})=p(\mathcal{D} \mid \mathbf{w}, \beta) p(\mathbf{w} \mid \alpha), Z=p(\mathcal{D} \mid \alpha, \beta)\)として、上の式に代入して。

$$
\begin{aligned}
p(\mathcal{D} \mid \alpha, \beta) & \simeq p\left(\mathcal{D} \mid \mathbf{w}_{\mathrm{MAP}}, \beta\right) p\left(\mathbf{w}_{\mathrm{MAP}} \mid \alpha\right) \int \exp \left\{-\frac{1}{2}\left(\mathbf{w}-\mathbf{w}_{\mathrm{MAP}}\right)^{\mathrm{T}} \mathbf{A}\left(\mathbf{w}-\mathbf{w}_{\mathrm{MAP}}\right)\right\}d\mathbf{w} \\
& =f\left(\mathbf{w}_{\mathrm{MAP}}\right) \frac{(2 \pi)^{W / 2}}{|\mathbf{A}|^{1 / 2}}
\end{aligned}
$$

\(f\left(\mathbf{w}_{\mathrm{MAP}}\right)\)について展開すると

$$
\begin{aligned}
f\left(\mathbf{w}_{\mathrm{MAP}}\right)= & p\left(\mathcal{D} \mid \mathbf{w}_{\mathrm{MAP}}, \beta\right) p\left(\mathbf{w}_{\mathrm{MAP}} \mid \alpha\right) \\
= & \prod_{n=1}^N \mathcal{N}\left(t_n \mid y\left(\mathbf{x}_n, \mathbf{w}_{\mathrm{MAP}}\right), \beta^{-1}\right) \mathcal{N}\left(\mathbf{w}_{\mathrm{MAP}} \mid \mathbf{0}, \alpha^{-1} \mathbf{I}\right) \\
= & \prod_{n=1}^N\left(\frac{\beta}{2 \pi}\right)^{1 / 2} \exp \left[-\frac{\beta}{2}\left\{t_n-y\left(\mathbf{x}_n, \mathbf{w}_{\mathrm{MAP}}\right)\right\}^2\right] \\
& \frac{1}{(2 \pi)^{W / 2}} \frac{1}{\left|\alpha^{-1} \mathbf{I}\right|^{1 / 2}} \exp \left\{-\frac{1}{2} \mathbf{w}_{\mathrm{MAP}}^{\mathrm{T}}\left(\alpha^{-1} \mathbf{I}\right)^{-1} \mathbf{w}_{\mathrm{MAP}}\right\} \\
= & \prod_{n=1}^N\left(\frac{\beta}{2 \pi}\right)^{1 / 2} \exp \left[-\frac{\beta}{2}\left\{t_n-y\left(\mathbf{x}_n, \mathbf{w}_{\mathrm{MAP}}\right)\right\}^2\right]\left(\frac{\alpha}{2 \pi}\right)^{W / 2} \exp \left(-\frac{\alpha}{2} \mathbf{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathbf{w}_{\mathrm{MAP}}\right).
\end{aligned}
$$

これの対数を取って

$$
\begin{aligned}
\ln p(\mathcal{D} \mid \alpha, \beta) & \simeq \ln f\left(\mathbf{w}_{\mathrm{MAP}}\right)+\frac{W}{2} \ln (2 \pi)-\frac{1}{2} \ln |\mathbf{A}| \\
& =\sum_{n=1}^N\left[\frac{1}{2}\{\ln \beta-\ln (2 \pi)\}-\frac{\beta}{2}\left\{t_n-y\left(\mathbf{x}_n, \mathbf{w}_{\mathrm{MAP}}\right)\right\}^2\right] \\
& \qquad+\frac{W}{2}\{\ln \alpha-\ln (2 \pi)\}-\frac{\alpha}{2} \mathbf{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathbf{w}_{\mathrm{MAP}}+\frac{W}{2} \ln (2 \pi)-\frac{1}{2} \ln |\mathbf{A}| \\
& =-\left[\frac{\beta}{2} \sum_{n=1}^N\left\{t_n-y\left(\mathbf{x}_n, \mathbf{w}_{\mathrm{MAP}}\right)\right\}^2+\frac{\alpha}{2} \mathbf{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathbf{w}_{\mathrm{MAP}}\right]\\
& \qquad-\frac{1}{2} \ln |\mathbf{A}|+\frac{N}{2} \ln \beta-\frac{N}{2}\ln(2\pi)+\frac{W}{2}\ln(\alpha) \\
& =-E\left(\mathbf{w}_{\mathrm{MAP}}\right)-\frac{1}{2} \ln |\mathbf{A}|+\frac{W}{2} \ln \alpha+\frac{N}{2} \ln \beta-\frac{N}{2} \ln (2 \pi)
\end{aligned}
$$

よって示せた。□

この式を同様に計算すると、結局回帰の時の結果と同じなる。つまり、再推定方程式は

$$
\alpha=\frac{\gamma}{\mathbf{w}_{\mathrm{MAP}}^{\mathrm{T}} \mathbf{w}_{\mathrm{MAP}}}
$$

となる。

エビデンス理論（経験ベイズの流れ）を人工的２クラスデータ集合に用いた例、緑の曲線は最適な決定境界を示し、黒い曲線は８個の隠れユニットを持つ二層ネットワークに最尤推定を用いたフィッティングを示している。赤い曲線は正則化項を含む場合を示していて、\(\alpha\)は\(\alpha= 0\)を初期値としてエビデンスの手続きで最適化されている。

赤い線が今回のベイジアンニューラルネットワークに相当するが、黒い線の最尤推定の方法に比べて過学習が大幅に抑制されていることがわかると思う。

予測分布を考える

最後に予測分布を求めてみよう。

一つの方法はある意味分散を無視して

$$
p(t \mid \mathbf{x}, \mathcal{D}) \simeq p\left(t \mid \mathbf{x}, \mathbf{w}_{\mathrm{MAP}}\right)
$$

 

とすることで、もう一つの方法は、分散を考慮する。

 

２クラス分類では出力ユニットがロジスティックシグモイド活性化関数を使用していることから、出力が(0,1)に制限されているので、回帰の時のようにネットワークの出力そのものを線形近似するのは不適切である。なので活性化関数を通す前の出力ユニットの活性を線形近似することを考える。

$$
a(\mathbf{x}, \mathbf{w}) \simeq a_{\mathrm{MAP}}(\mathbf{x})+\mathbf{b}^{\mathrm{T}}\left(\mathbf{w}-\mathbf{w}_{\mathrm{MAP}}\right)
$$

ここで、\(a_{\mathrm{MAP}}(\mathbf{x})=a\left(\mathbf{x}, \mathbf{w}_{\mathrm{MAP}}\right)\)であり、ベクトル\(\mathbf{b} \equiv \nabla a\left(\mathbf{x}, \mathbf{w}_{\mathrm{MAP}}\right)\)は逆伝播によって求めることができるものである。

\(\mathbf{w}\)の事後分布のガウス近似とさらに\(\mathbf{w}\)の線形関数である\(a\)のモデルが得られたので、ネットワークの重みの分布から導かれる出力ユニットの活性の分布は

$$
p(a \mid \mathbf{x}, \mathcal{D})=\int \delta\left(a-a_{\mathrm{MAP}}(\mathbf{x})-\mathbf{b}^{\mathrm{T}}(\mathbf{x})\left(\mathbf{w}-\mathbf{w}_{\mathrm{MAP}}\right)\right) q(\mathbf{w} \mid \mathcal{D}) \mathrm{d} \mathbf{w}
$$

とあらわすことができる。ここで\(q(\mathbf{w} \mid \mathcal{D})=\mathcal{N}\left(\mathbf{w} \mid \mathbf{w}_{\mathrm{MAP}}, \mathbf{A}^{-1}\right)\)で与えられる事後分布のガウス近似である。

この分布は

平均：\(a_{\mathrm{MAP}} \equiv a\left(\mathbf{x}, \mathbf{w}_{\mathrm{MAP}}\right)\)
分散：\(\sigma_a^2(\mathbf{x})=\mathbf{b}^{\mathrm{T}}(\mathbf{x}) \mathbf{A}^{-1} \mathbf{b}(\mathbf{x})\)

のガウス分布であることがわかる。なので、結局のところ予測分布は

$$
p(t=1 \mid \mathbf{x}, \mathcal{D})=\int \sigma(a) p(a \mid \mathbf{x}, \mathcal{D}) \mathrm{d} a
$$

を用いて\(a\)についての周辺化をしなければならず、この計算は困難なので以下の近似を用いる

$$
\int \sigma(a) \mathcal{N}\left(a \mid \mu, \sigma^2\right) \mathrm{d} a \simeq \sigma\left(\kappa\left(\sigma^2\right) \mu\right)
$$

これより、予測分布は

２クラス分類の予測分布

$$
p(t=1 \mid \mathbf{x}, \mathcal{D})=\sigma\left(\kappa\left(\sigma_a^2\right) \mathbf{b}^{\mathrm{T}} \mathbf{w}_{\mathrm{MAP}}\right)
$$

ただし、

$$
\kappa\left(\sigma^2\right)=\left(1+\pi \sigma^2 / 8\right)^{-1 / 2}
$$

で定義される。

活性化関数が\(tanh\)である隠れユニット８個とロジスティックシグモイド出力ユニットを１個持つベイジアンニューラルネットワークにラプラス近似を用いた例である。重みパラメータの導出にはスケール化共役勾配法を用い、ハイパーパラメータ\(\alpha\)はエビデンス理論（経験ベイズ）を用いて最適化している。左の図はパラメータの点推定\(\mathbf{w}_{MAP}\)に基づく単純な近似（分散を無視）を用いた結果である。緑の曲線は\(y=0.5\)の決定境界を示していて、ほかの等高線はそれぞれの出力の確率が\(y = 0.1,0.3,0.7,0.9\)に対応している。右の図は分散を考慮したパターンを用いた結果

 

 

 

 

-勉強, 数学, 機械学習/AI