第3章-エビデンス近似をSageで試す

第1章-Sageを使って線形回帰を試してみるでは、 与えられたモデルパラメータM, α、βに対して解いていましたが、 実際には自分でベストなモデルパラメータM, α、βを求める必要があります。 パターン認識と機械学習 の３章ではエビデンス関数を使って最適なパラメータを求める方法が説明されています。

ここでは、Sageを使ってエビデンス関数の評価、最適なモデルパラメータの推定を試してみます。

エビデンス関数の評価

エビデンス関数の対数表現は、式(3.86) $$ \ln p(t|\alpha, \beta) = \frac{M}{2} \ln \alpha + \frac{N}{2} \ln \beta -E(m_N) - \frac{1}{2} \ln | A | + \frac{N}{2} \ln (2 \pi) $$ で与えられます。

ここで、$E(m_N)$は、式(3.82) $$ E(m_N) = \frac{\beta}{2}|| t = \Phi m_N ||^2 + \frac{\alpha}{2} m_N^T m_N $$ $m_N$は、式(3.84) $$ m_N = \beta A^{-1} \Phi^T t $$ Aは、式(3.81) $$ A = \alpha I + \beta \Phi^T \Phi $$

目指すは、図3.14です。

準備

第1章-Sageを使って線形回帰を試してみる で使ったデータと同じものを座標Xと目的値tにセットし、関数Φを定義します。

エビデンスの計算

まずは、α、βを固定値

• $\alpha = 5 \times 10^{-3}$
• $\beta = 11.1$

一番良いエビデンスは、M=4あたりになります。 理由は分かりませんが、図3.14のようなM=3以降の急激な下降は見られません。

β_MLを使ってエビデンスを計算

そこで、βを平均値の分散$\beta_{ML}$を使って計算してみました。 $$ \frac{1}{\beta_{ML}} = \frac{1}{N} \sum (y(x_n, W_ML) - t_n)^2 $$

期待に反し、M=0でも高い値となり、あまり良くありません。

どうして図3.14と合わないのか？

どうして図3.14のような結果にならないのかと調べたのですが、他にも同様の計算をした 方がいらして、私と同じ傾向になったとの記述がありました。 基底関数を色々変えて、線形回帰のエビデンスを計算してみた

私の推測では、ln |A|の計算をA.norm()と間違えたのではないかと思われます。以下にA.norm()としたときの図を示します。

最適なα、βを求める

α、βを固定にした場合に、エビデンスが最大はM=4だったので、M=4に対する 最適なα、βを求めてみます。

エビデンスを最大にするαは、式(3.92) $$ \alpha = \frac{\gamma}{m_N^T m_N} $$ ここで$\gamma$は、式(3.87) $$ (\beta \Phi^T \Phi) u_i = \lambda_i u_i $$ の固有値から、式(3.91) $$ \gamma = \sum_{i} \frac{\lambda_i}{\alpha + \lambda_i} $$ で求まります。

βについても、式(3.95) $$ \frac{1}{\beta} = \frac{1}{N - \gamma} \sum_{n=1}^{N} ( t_n - m_N^T \phi(x_n))^2 $$

$\gamma$は、$\alpha$に依存し、$m_N$も$\alpha$に依存するため、すぐに最適な値を得ることが できません。

1. 適当な $\alpha, \beta$を初期値とする
2. $m_N$, $\gamma$を求め
3. 式(3.95)から$\beta$を求める
4. 式(3.92)から$\alpha$を求め、2番目からを繰り返す

以下の例では、

• $\alpha = 5 \times 10^{-3}$
• $\beta = 11.1$

γと2αEw(mN)のグラフ

M=4の多項式回帰に対して、図3.16(a)

のように$\gammaと2\alpha E_w(m_N)$のグラフを$\ln \alpha$の関係を図化してみます。

形が図3.16とは異なり、曲線の交差する点もα=0.024だとln αの値は負になるのに、図3.16は正の部分で 交差しています。（何かちがう！）

残差と対数エビデンスのグラフ

同様に、M=4の多項式回帰に対して、図3.16(b)

のように残差と対数エビデンスを$\ln \alpha$の関係を図化してみます。

やはり形が、図3.16とは異なります。

M=4に対する最適な解

M=4の多項式回帰に対する最適なα、βを使って$sin(2 \pi x)$、観測データ、フィッティング曲線、分散をグラフに表示してみます。

M=9の最適解（図1.17）

どうしてエビデンス関数が合わないのか

つぎにどうして、図3.16と合わないのか、調べてみます。

図3.16の説明をよく読むと、基底関数が多項式ではなく、ガウス基底関数になっていました。

基底関数のチェック

ガウス基底関数がどのようなものなのか、図3.1を再現しながら、見てみましょう。

以下に、

• 多項式の基底関数
• ガウス基底関数
$$ \phi_j(x) = exp \left \{ - \frac{(x - \mu_j)^2}{2 s^2} \right \} $$
• シグモイド基底関数
$$ \phi_j(x) = \sigma \left( \frac{x - \mu_j}{s} \right) $$ $$ \sigma(a) = \frac{1}{1 + exp(-a)} $$

ガウス基底関数、シグモイド基底関数では、μの値として-1から1（問題によって0から1の場合もある）に とり、それを基底関数の数で等分した値をとるみたいです。 線形回帰モデルとかを参考にしました。

ガウス基底関数を使ったγと2αEw(mN)のグラフ

M=9個のガウス基底関数を使ったγと2αEw(mN)のグラフを以下に示します。

関数の形状と交差点の位置が図3.16と合います。

ガウス基底関数を使った残差と対数エビデンスのグラフ

同様に、ガウス基底関数を使った残差と対数エビデンスのグラフを示します。

残差の形状がやや違いますが、残差の最低値の付近に対数エビデンスの最大値が位置する 部分はほぼ合っています。

ガウス基底関数を回帰結果

おまけに、ガウス基底関数の回帰結果と分散の分布、オリジナルのsin曲線、データ点を プロットした図を示します。ガウス基底を使った回帰も良く合っていることが分かります。