しまてく

学んだ技術を書きためるブログ

ITエンジニアのための機械学習理論入門 6

第6章 k平均法: 教師なしモデルの基礎

  • k平均法は教師なし学習モデル。
  • 与えられたトレーニングセットをクラスタリングする方法。

  • いままでの教師あり学習では、トレーニングセットは目的変数$t_n$を持っていた。

  • そうしたデータを分析し、未知のデータに対して推測をしてきた。

  • 一方教師なし学習では目的変数を持たない。

  • そのため、分類をするための何らかの基準を決める必要がある。
  • この基準を「二重歪み」と呼び、これをなるべく最小化するようにする。

k平均法によるクラスタリング

  • トレーニングセット $ {(x_n, y_n)}_{n=1}^N $ が与えられたとする。
  • 感覚的に、このデータは二つに分類できそうな気がする。
  • そこで、k平均法を用いて、これらのデータを二つのクラスターに分類する。
    • 分類する数は事前に決める必要がある。

実際の分類の手順

  1. トレーニングセットに含まれるデータは ${\bf x}_n = (x_n, y_n)^T$と表記する。
  2. 適当な代表点 ${\mu_k}_{k=1}^2$ を決める。
  3. 各点と代表点との距離 $|{\bf x}_n - \mu_k|$ を計算して、距離が短い方の代表点に属するものとする。
  4. ここで、どちらの代表点に所属するかを表す変数$r_{nk}$を定義しておく。

$$ r_{nk} = \left\{ \begin{array}{l} 1 \ \ \ \ {\bf x}_n がk番目の代表点に属する場合 \\ 0 \ \ \ \ それ以外の場合 \end{array} \right . \tag{6.1} $$

  1. 次に、(3)で分類したクラスターの「重心」を新たな代表点とする。
  2. 重心の公式を使って式を立てると以下のようになる。
    • $N$は$k$番目の代表点に属する点の個数。
    • 分子の$\sum$は$k$番目の代表点に属する点についてのみ加える。

$$ \mu_k = \dfrac{\sum{\bf x}_n}{N_k} (k = 1,2) \tag{6.2} $$

  1. (6)を(6.1)を用いてスマートに表現すると、

$$ \mu_k = \dfrac{\sum_{n=1}^Nr_{nk}{\bf x}_n}{\sum_{n=1}^Nr_{nk}} \tag{6.3} $$

  1. (5)から(7)をクラスターに所属する点に変化がなくなるまで繰り返す。
  2. 最終的に得られた ${\mu_k}_{k=1}^2$ がそれぞれのクラスターの代表点となる。

画像データへの応用

  • ある画像から代表色を抽出するという問題を解く。
  • 最も簡単なのは、RGBの各値を3次元空間の点とみなす方法。
  • 画像のすべてのピクセルを3次元空間にマッピングすると色の分布を可視化できる。
  • すべてのデータに対してクラスタリングをすることで代表色を決められる。
  • ピクセルを代表色で置き換えることで、画像の減色処理が実現できる。

k平均法の数学的根拠

  • k平均法では「二乗歪み」というグループ分けの基準を用いている。

$$ J = \sum_{n=1}^N\sum_{k=1}^Kr_{nk}|{\bf x}_n - \mu_k|^2 \tag{6.4} $$

  • $r_{nk}$の定義(6.1)を考えると、これは『「自分の属するクラスターの代表点からの距離の2乗」の合計』になっている。
  • つまり、$J$を小さくするということは、「代表点のなるべく近くにデータが集まるように分類する」ということ。
  • 以下は、k平均法の手続きによって、この「二乗歪み」が最終的に極小値に達することを証明する。

証明

  • トレーニングセットは、任意の特定次元ベクトルの集合 ${{\bf x}_n}_{n=1}^N$とする。
  • K個のクラスターに分類するとして、代表点を${\mu_k}_{k=1}^K$とする。
  • 各データが属するクラスタ-を以下のように表す。 $$ r_{nk} = \left\{ \begin{array}{l} 1 \ \ \ {\bf x}_n がk番目の代表点に属する場合 \\ 0 \ \ \ それ以外の場合 \end{array} \right. \tag{6.5} $$

  • 現在の分類における「二乗歪み」は以下の式で表せる。

$$ J = \sum_{n=1}^N\sum_{k=1}^Kr_{nk}|{\bf x}_n - \mu_k|^2 \tag{6.6} $$

以下では、k平均法の手続きに従って $r_{nk}$ と $\mu_k$ を修正していくと$J$の値は減少し、極小値に達することを示す。

  • まず、各データが属するクラスターを選択しなおす。
  • この時、各データ${\bf x}_n$について、代表点からの距離$|{\bf x}_n - \mu_k|$が最も小さいクラスターを選択する。
    • Jの意味を考えれば、この操作でJの値が大きくなることはない。
    • この操作は、$r_{nk}$を以下の条件で再定義していることになる。

$$ r_{nk} = \left\{ \begin{array}{l} 1 \ \ \ k = \mathop{\rm arg\,min}\limits_{k'}|{\bf x}_n - \mu_{k'}| \\ 0 \ \ \ それ以外の場合 \end{array} \right. \tag{6.7} $$

  • ($\mathop{\rm arg\,min}\limits_{x}f_x$は、$f_x$を最小にする$x$の値)
  • 次に、現在の分類状態において、各クラスターの代表点を取り直す。
  • (6.6)を最小にするという条件で、 $\mu_k$を選択する。
  • (6.6)を$\mu_k$について見ると、下に凸な二次関数となっている。
    • これは偏微分係数が0になるという条件で最小化するということ。
  • $J$を成分表記すると以下のようになる。($[x_n]_i$という記号はベクトル${\bf x}_n$の$i$成分を表す)

$$ J = \sum_{n=1}^N \sum_{k=1}^K \left\{r_{nk}\sum_i([{\bf x}_n]i - [\mu_k]i)^2 \right\} \tag{6.8} $$

  • これを特定の成分$[\mu_k]_i$で偏微分する。
  • $k$と$i$が特定の成分なので、$k = k'$, $i = i'$として、$J$を以下のように変形する。

$$ \begin{eqnarray} J' &=& \sum_{n=1}^Nr_{nk'}([{\bf x}_n]i' - [\mu_{k'}]i')^2 \nonumber \\ &=& \sum_{n=1}^Nr_{nk'}([{\bf x}_n]i'^2 -2[{\bf x}_n]i' [\mu_{k'}]i' + [\mu_{k'}]i'^2) \nonumber \end{eqnarray} $$

  • これを$[\mu_{k'}]_{i'}$で偏微分する。

$$ \dfrac{\partial J}{\partial [\mu_{k'}]_{i'}} = -2 \sum_{n=1}^Nr_{nk'}([{\bf x}_n]_{i'} - [\mu_{k'}]_{i'}) $$

  • 元々は特定の成分$[\mu_k]_i$なので、書き直すと以下になる。

$$ \dfrac{\partial J}{\partial [\mu_k]_i} = -2 \sum_{n=1}^Nr_{nk}([{\bf x}_n]_i - [\mu_k]_i) \tag{6.9} $$

  • これが0になるという条件から、$[\mu_k]_i$は以下のように決まる。

$$ \begin{eqnarray} -2 \sum_{n=1}^Nr_{nk}([{\bf x}_n]_i - [\mu_k]_i) &=& 0 \nonumber \\ -2 \sum_{n=1}^Nr_{nk}[{\bf x}_n]_i + 2 \sum_{n=1}^Nr_{nk}[\mu_k]_i &=& 0 \nonumber \\ \sum_{n=1}^Nr_{nk}[\mu_k]_i &=& \sum_{n=1}^Nr_{nk}[{\bf x}_n]_i \nonumber \\ [\mu_k]_i &=& \dfrac{\sum_{n=1}^Nr_{nk}[{\bf x}_n]_i}{\sum_{n=1}^Nr_{nk}} \tag{6.10} \end{eqnarray} $$

  • 成分表記からベクトル表記に戻すと、以下の結果が得られる。

$$ \begin{eqnarray} \mu_k &=& \dfrac{\sum_{n=1}^Nr_{nk}{\bf x}_n}{\sum_{n=1}^Nr_{nk}} \tag{6.11} \end{eqnarray} $$

  • これは各クラスターの重心を新たな代表点に取る(6.3)と同じものになっている。
  • 従って(6.3)の手続きによっても$J$が大きくなることはない。
  • 以上により、k平均法の操作を繰り返すと、Jの値は必ず小さくなるか、それ以上変化しない極小値に達する。

ただし、二乗歪み$J$は有限回の操作で極小値に達するが、完全に値が変化しなくなるまでには計算に時間がかかるため、$J$の減少分が0.1%未満になるなどで計算を打ち切る。

  • この例では3次元空間の点$x$だったので、通常のユークリッド距離$|{\bf x}_n - \mu_k|$で計算をした。
  • ただし、これ以外の距離を用いてもk平均法の手続きを実施することはできる。
  • 例えば、「文書」のカテゴリー分けで、「文書間の距離」を定義すれば利用可能。
  • 何らかの方法で類似度を計算し、類似度が高いほど距離が近くなるようにする。
    • 単語の出現頻度などで判断する。

怠惰学習モデルとしてk近傍法

  • k近傍法は教師あり学習に分類される。
  • しかし一般的な目的変数とはことなり、近傍データによって自分の目的変数を推測する。
  • 近い方からK個分のデータを見て推測するとして、
    • K=1の場合はトレーニングセットのデータにオーバーフィッティングする。
    • K=3などの場合はよりなめらかに分類できる。
  • しかし、k近傍法は2つの点で問題がある。
    1. 新たなデータの分類に時間がかかる。
      • データが与えられた場合、すべてのデータについて参照して計算しなおす必要がある。
    2. 分析のモデルが明確ではない。
      • 与えられたデータからたまたまそのように分類できた、というだけ。
      • 仮説と実証ではなく、単純に事実からの判断でしかないため将来に生かせない。

ITエンジニアのための機械学習理論入門 5-1

5章 ロジスティック回帰とROC曲線: 学習モデルの評価方法

  • ロジスティック回帰は、パーセプトロンと同じ分類アルゴリズムの一つ。
  • 確率を用いた最尤推定法でパラメーターを決定する点が異なる。
  • 確率を用いる事で、未知のデータの属性を推定する際に違いがある。
    • 「このデータは$t=1$である」という単純なものではなく
    • 「このデータが$t=1$である確率は70%」というような感じになる。
  • さらにROC曲線を用いて、学習モデルを評価する方法を解説する。

5.1 分類問題への最尤推定法の適用

  • 最尤推定法では、適当に決めた「あるデータが得られる確率」から、 トレーニングセットとして与えられる確率(尤度関数)を計算する。
  • 尤度関数が最大になるという条件から、設定した確率の式のパラメーターを求める。

5.1.1 データの発生確率の設定

  • まずはパーセプトロンと同様に2種類のデータを分類する直線を表す線形関数$f(x, y)$を用意する。

$$ f(x, y) = \omega_0 + \omega_1x + \omega_2y \tag{5.1} $$

  • 図5.1(本文参照)のように $f(x, y) = 0$で分割線が決まった場合、
  • 分割線に直行する方向に移動すると $- \infty < f(x, y) < \infty $の範囲で$f(x, y)$の値が変化する。
  • 次に、$(x, y)$平面上の任意の点で、得られたデータが$t=1$である確率を考える。
  • 図5.1(本文参照)において、分割線から右上に離れるほど$t=1$である確率は高いと考える。
  • 逆に左下に離れれば低くなると考える。
  • 分割線上では$t=1$と$t=0$の確率は同じ($\dfrac{1}{2}$)。
  • また、今回データは2つしか無いので、$t=1$となる確率を$P$とすると$t=0$となる確率は$1-P$となる。
  • 得られた全てのデータについて、$t=1$である確率を対応させると、図5.2の下のグラフになる。
  • これをロジスティック関数のグラフと言う。

  • また、ロジスティック関数は以下の式で表せる。

$$ σ(a) = \dfrac{1}{1 + e^{-a}} \tag{5.2} $$

  • この式において、$a$の値を$-\infty$から$\infty$まで変化させると $ σ{(a)} $は0から1に向かって滑らかに変化する。
  • $a=0$ではちょうど$\dfrac{1}{2}$になる。
  • この関数の$a$として$f(x, y)$を代入すると、図5.2の対応関係が得られる。

以上をまとめると、

  • 点$(x, y)$で得られたデータが$t=1$である確率は次式で表せる。

$$ P(x, y) = σ{(\omega_0 + \omega_1x + \omega_2y)} \tag{5.3} $$

  • 反対に、$t=0$となる確率は$1 - P(x, y)$となる。

この確率から、トレーニングセットとして与えられたデータ${(x_n, y_n, t_n)}_{n=1}^{N}$ が得られる確率を考える。

  • あるデータ$(x_n, y_n, t_n)$が得られる確率は、$t=1$と$t=0$で次のように場合分けされる。

$$ t=1 の場合 : P(x_n, y_n) \tag{5.4} $$

$$ t=0 の場合 : 1 - P(x_n, y_n) \tag{5.5} $$

  • 数学的な技巧を使うと以下のようにまとめて書くことができる。

$$ P_n = P(x_n, y_n)^{t_n}{1 - P(x_n, y_n)}^{1-t_n} \tag{5.6} $$

  • ここで(5.6)に(5.3)を代入すると次のようになる。

$$ P_n = σ(\omega_0 + \omega_1x_n + \omega_2y_n)^{t_n}{1 - σ(\omega_0 + \omega_1x_n + \omega_2y_n)}^{1-t_n} $$

  • $σ(\omega_0 + \omega_1x_n + \omega_2y_n)$は、「n番目のデータが$t=1$である確率」を表す。
  • これを$z_n$と定義し、さらに行列を用いて表すと以下のようになる。

$$ P_n = z_n^{t_n}(1 - z_n)^{1-t_n} \tag{5.9} $$

  • $z_n$の定義は以下になる。 $$ z_n = σ(w^T\phi_n) \tag{5.10} $$

  • $w$と$\phi$はパーセプトロンの計算で用いた(4.12)、(4.13)と同じ。

トレーニングセットに含まれる全てのデータをまとめて考えると、これらが得られる確率$P$は、各データが得られる確率($P_n$)(5.9)の積になる。

$$ P = Π_{n=1}^{N}P_n = Π_{n=1}^{N}z_n^{t_n}(1-z_n)^{1-t_n} \tag{5.13} $$

  • 確率$P$は、(5.10)を通して、求めるべき係数$w$の関数になっている。

5.1.2 最尤推定法によるパラメーターの決定

$$ w_{new} = w_{old} -(Φ^TRΦ)^{-1}Φ^T(z-t) \tag{5.14} $$

  • $t$はトレーニングセットの各データの属性値$t_n$を並べたベクトル
  • $\Phi$は各データの座標を表すベクトル$\Phi_n$を横ベクトルにして並べたNx3行列。
  • $z$は(5.10)の$z_n$を並べたベクトル。
  • $R$は$z_n(1-z_n)$を対角成分とする対角行列。

$$ t = \left( \begin{array}{c} t_1 \\ \vdots \\ t_n \end{array} \right) \tag{5.15} $$

$$ \Phi = \left( \begin{array}{c} 1 & x_1 & y_1 \\ 1 & x_2 & y_2 \\ \vdots & \vdots & \vdots \\ 1 & x_n & y_n \end{array} \right) \tag{5.16} $$

$$ z = \left( \begin{array}{c} z_1 \\ \vdots \\ z_n \end{array} \right) \tag{5.17} $$

$$ R = \left( \begin{array}{c} z_1(1-z_1) & \cdots & \cdots & 0 \\ \vdots & z_2(1-z_2) & \cdots & \cdots \\ \vdots & \vdots & \ddots & \cdots \\ \vdots & \vdots & \vdots & z_n(1-z_n) \end{array} \right) \tag{5.18} $$

  • ここで、$t$と$\Phi$はトレーニングセットのデータから決まる「定数」のベクトル/行列。
  • 一方$z$と$R$に含まれる$z_n$は(5.10)を通してパラメーター$w$に依存している。
  • つまり、パラメーター $ w_{old} $ が与えられた際に、これを用いて$z$と$R$を計算しておき、それを(5.14)に代入することで修正された新しいパラメーター $ w_{new} $ が得られる。
  • この $ w_{new} $ を次の $ w_{old} $として繰り返し計算していく。
  • 最終的に(5.13)の確立$P$が最大値に達する。
  • ※具体的な導出は「5.3 付録 - IRLS法の導出」を参照。
  • この方法は「IRLS (Iteratively Reweighted Least Squares)法、あるいは反復再重み付け最小二乗法」という。
  • 後述のニュートン法との類似性から、(5.14)の計算を繰り返し$P$の値が最大値に近づくに つれて、パラメーター$w$の変化の割合は小さくなる。
  • サンプルでは以下の条件で計算を打ち切る。

$$ \dfrac{{|| w_{new} - w_{old} ||}^2}{|| w_{old} ||} < 0.001 \tag{5.19} $$

  • この条件は、トレーニングセットに対してオーバーフィッティングを発生させないためのガード。

5.1.3 サンプルコードによる確認

「05-logistic_vs_perceptron.py」の実行結果は以下。 f:id:cimadai:20160523100500p:plain

  • 実線はロジスティック回帰の結果、破線がパーセプトロンの結果。
  • 4角トレーニングセットはランダムで、それぞれがどのように分類されるか。

  • パーセプトロンでは、「$n=1$〜$N$について、正しく分類されていなければパラメーター$w$を修正する」という処理を30回繰り返し、終了する。
  • ロジスティック回帰では、(5.14)でパラメーター$w$を修正した際に(5.19)の条件が成立した場合に計算を終了する。
  • (5.14)の修正を30回繰り返しても(5.19)が成り立たない場合はその時点で終了。
  • 図を見ると、ロジスティック回帰ではよりデータの中央で分割していることがわかる。
  • これはパーセプトロン確率的勾配降下法では、一度全てのデータが正しく分類されると、そこでパラメーターの変化が停滞するため。
  • ロジスティック回帰の場合は、トレーニングセットのデータが得られる全体的な確率を最大化しようとするため、よりもっともらしい分割線が選択される。

ITエンジニアのための機械学習理論入門 4-1

第4章 パーセプトロン : 分類アルゴリズムの基礎

  • 複数の観測値があるときに、$ t=\pm 1 $ に分類するような直線を発見する問題。
  • 最小二乗法に類似の「誤差関数」を使って解く。
  • ただし「紙と鉛筆による計算」では解けない。
  • 数値計算をしてパラメータの修正を繰り返す「確立的勾配降下法」を利用する。

4.1 確率的勾配降下法アルゴリズム

パラメトリックモデルの3つのステップ(再掲)

  1. パラメーターを含むモデル(数式)を設定する
  2. パラメーターを評価する基準を定める
  3. 最良の評価を与えるパラメーターを決定する

4.1.1 平面を分割する直線の方程式

  • まずステップ1として「パラメーターを含むモデル」を設定する。
  • 今回は($x$, $y$)平面上のデータを分類する直線を求めることが目的。
  • 直線の式と言えば $ y = ax + b$ だが、ここでは$x$と$y$を対等に扱うため、以下の関数を設定する。

$$ f(x, y) = \omega_0 + \omega_1x + \omega_2y \tag{4.1} $$

  • この時、($x$, $y$)平面を分割する直線は以下の式で表せる。

$$ f(x, y) = 0 \tag{4.2} $$

  • (4.2)の直線(境界線)から離れれば離れるほど$f$($x$, $y$)の絶対値は大きくなる。

  • ($x$, $y$)平面をこのように分割するのは、$ t=\pm 1 $ の2種類の属性を持つデータを分類すること。

  • ここでは次のルールでデータを分類するものとする。

$$ f(x, y) > 0 \Rightarrow t = +1 \tag{4.3} $$

$$ f(x, y) < 0 \Rightarrow t = -1 \tag{4.4} $$

  • この時、トレーニングデータとして与えられたデータ $ {( x_n, y_n, t_n )}_{n=1}^N $ について、正しく分類されているかは以下で判定する。

$$ f(x_n, y_n) \times t_n > 0 : 正解 \tag{4.5} $$

$$ f(x_n, y_n) \times t_n \le 0 : 不正解 \tag{4.6} $$

  • $ t=\pm 1 $ のどちらの場合においても、同じルールで正誤判定ができるのがポイント。
  • すべての $( x_n, y_n, t_n )$ について、(4.5)が成り立つような直線、すなわち(4.1)の係数($\omega_0, \omega_1, \omega_2$)を見つけることがゴール。
  • ということで、ステップ2としてパラメーター($\omega_0, \omega_1, \omega_2$)の評価基準を定める。

4.1.2 誤差関数による分類結果の評価

  • パラメータの誤差として、正しく分類できない点、すなわち(4.6)が成り立つ点があった場合に、それを誤差として計算する。
  • ここで言う誤差は、誤判定に対するペナルティと考えてもいい。
  • この誤差が小さいほど正しい分類に近いと言える。
  • 具体的な誤差の値としては、「境界線から離れるほど $f(x, y)$ の絶対値が大きくなる」という特徴を利用して次の値を採用する。

$$ ある点の誤差 E_n = \mid f(x_n, y_n) \mid \tag{4.7} $$

  • (4.7)は、正しく分類できなかった点についてのみ計算することに注意。
  • 誤って分類された点の誤差を合計したものを、分類の誤差 $E$ とする。

$$ E = \sum_{n}{E_n} = \sum_{n} \mid f(x_n, y_n) \mid \tag {4.8} $$

  • (4.8)は誤分類された点の誤差の合計を表している。
  • このような点は(4.6)を満たすので、次の関係式が成り立つ。

$$ \mid f(x_n, y_n) \mid = - f(x_n, y_n) \times t_n \tag {4.9} $$

  • (4.9)および $ f(x, y) $ の定義 (4.1)を用いると、(4.8)は次のように表せる。

$$ E = - \sum_{n}{(\omega_0 + \omega_1x_n + \omega_2y_n ) t_n} \tag {4.10} $$

  • あるいは、ベクトルを用いて次のように表すこともできる。

$$ E = - \sum_{n}{t_n w^T \phi_n} \tag {4.11} $$

  • ここで、$w$と$\phi_n$は次で定義されるベクトルになる。

$$ w = \left( \begin{array}{r} \omega_0 \\ \omega_1 \\ \omega_2 \end{array} \right) \tag{4.12} \ $$

$$ \phi_n = \left( \begin{array}{r} 1 \ \\ x_n \\ y_n \end{array} \right) \begin{array}{r} \longleftarrow バイアス項 \\ \\ \ \end{array} \tag{4.13} $$

  • $w$は求めるべき係数をならべたベクトル。
  • $\phi_n$はトレーニングセットに含まれるデータの座標$(x_n, y_n)$を並べたもの。
  • ただし、係数$w_0$に対応する成分として、第1成分に定数1を入れてある。
  • これはやや技巧的ではあるが、(4.10)をベクトル形式で書き直すために導入した。
  • このような形で追加する定数項を「バイアス項」と呼ぶ。
  • (4.11)で計算される誤差$E$が小さくなるほど、トレーニングセットは適切に分類されていると言える。
  • すべてのデータが完全に正しく分類できれば、$E=0$になる。
  • 最後のステップ3は、(4.11)の誤差$E$を最小にするパラメーター$w$を求めること。
  • ただこれが難しく、本節のテーマである「確率的勾配降下法」を用いて解く。

4.1.3 勾配ベクトルによるパラメーターの修正

  • 最小二乗法ではパラメーターによる偏微分係数が0になるという条件から、二乗誤差$E_D$を最小にする係数$w$を決定することができた。
  • これと同様に、(4.10)の誤差$E$の偏微分係数を0と置いてみる。

$$ \dfrac{\partial \ E}{\partial \ \omega_m} = 0 \ \ \ (m = 0, 1, 2) \tag {4.14} $$

  • あるいは、ベクトル形式で、勾配ベクトルが0になるとしても同じ。

$$ \nabla E(w) = - \sum_{n}{t_n \phi_n} = 0 \tag {4.15} $$

  • 一般に、勾配ベクトルは次式で定義されるベクトルになる。

$$ \nabla E(w) = \left( \begin{array}{r} \dfrac{\partial \ E}{\partial \ \omega_0} \\ \dfrac{\partial \ E}{\partial \ \omega_1} \\ \dfrac{\partial \ E}{\partial \ \omega_2} \end{array} \right) \tag{4.16} $$

  • しかし、(4.14)や(4.15)を変形して係数$w$を表す式を得ることは不可能。
  • (4.15)を見るとわかるように、この式の中にはそもそも$w$が含まれていない。
  • そこで、純粋な式変形で$w$を求めることはやめ、数値計算を用いて$w$の値を修正しながら、誤差$E$がなるべく小さくなるものを求めることにする。
  • ここでヒントになるのが、勾配ベクトルの幾何学的な性質になる。
  • 例えば、$(x, y)$平面上に、原点(0, 0)を谷底とするようなすり鉢状の形状をした2変数関数$h(x, y)$があったとする。

$$ 例: h(x, y) = \dfrac{3}{4}(x^2 + y^2) \tag{4.17} $$

  • この場合、勾配ベクトルは次のように計算される。

$$ \nabla h = \left( \begin{array}{r} \dfrac{\partial \ h}{\partial \ x} \\ \dfrac{\partial \ h}{\partial \ y} \end{array} \right) = \left( \begin{array}{r} \dfrac{3}{2}x \\ \dfrac{3}{2}y \end{array} \right) \tag{4.18} $$

  • この時、任意の点$(x, y)$において、どちらかの方向に少し移動すると、hの値が増減するが、勾配ベクトル$\nabla h$は、$h$の値が最も増加する方向、すなわち「斜面をまっすぐに這い上がる方向」を表す。
  • また、勾配ベクトルの大きさは、その点における斜面(接平面)の傾きを表す。
  • つまり、各点における勾配ベクトルの方向に移動していけば$h(x, y)$の値はどんどん大きくなる。
  • 逆に、勾配ベクトルの反対方向に進めば、$h(x, y)$の値は小さくなる。
  • 現在点を$x{old}$とし、次の点$x{new}$を次の式で決定するアルゴリズムとして表現できる。

$$ x_{new} = x_{old} - \nabla h \tag{4.19} $$

  • 座標$x$を(4.19)に従って何度も更新していくと、次第に原点に近づき、原点にたどり着いた(勾配ベクトルが0になった)時に停止する。
  • 途中で$\nabla h$が大きすぎて原点を超えることがあるが、何度も行ったり来たりしながら次第に収束する。

$$ - \nabla E(w) = \sum_{n}{t_n \phi_n} \tag{4.20} $$

  • 確率的勾配降下法」という名前の「勾配降下」というのは、勾配ベクトルの反対方向にパラメーターを修正して、「誤差の谷」を降りていく、という考えから。
  • では「確率的」というのはなにか?
    • (4.20)を見ると、右辺は「正しく分類されなかった点」についての和になる。
    • 点が100個ある場合は100個分の ${t_n \phi_n}$ を合計して、その値を$w$に加えるという計算が必要。
    • しかしデータ数が膨大になると、事前に${t_n \phi_n}$を計算するのは難しいので、サンプリングを行う。
    • サンプリングは、正しく分類されていない点をどれか一つ選んで、とりあえずその分だけのパラメーターを修正する。

$$ w_{new} = w_{old} + t_n\phi_n \tag{4.21} $$

  • さらに、修正された新しい$w$の下で、正しく分類されていない点を一つ選び、繰り返し(4.21)の修正を行う。
  • このように「正しく分類されていない点」をランダムに選びながら、パラメーターを修正していくことこそが「確率的勾配降下法」となる。
  • ただし、今回の例では「ランダムに選ぶ」のも面倒なので、かたっぱしから計算していく。
  • n = 1~Nについて、処理が終わったらさらにもう一度n=1~Nについて同じ処理を行う。
  • これを何度も繰り返し、すべての点を正しく分類する直線にたどり着くと、「$(x_n, y_n)$が正しく分類されていない」という条件を満たさなくなるので、これ以上パラメーター$w$は変化しなくなるためアルゴリズムは終了。
  • すべての点を正しく分類できるかはトレーニングセット次第。
  • 仮にすべての点を正しく分類できる直線が存在する場合は、この処理でいつかその直線にたどり着くことができる。(Novikovの定理として数学的に証明されている)
  • しかしそのような直線が存在しない場合は、いつまでも$w$の値は変化し続ける。
    • 従って一定回数で打ち切るように処理をする必要がある。

おしまい

今回は計算式すくなめ。 ナブラ($\nabla$)というベクトル微分演算記号を覚えた!

ITエンジニアのための機械学習理論入門 3-2

3.2 単純化した例による解説

いままでやっていた例題をもう少し単純化する。

  • 今までは複数の観測点 $ \{x_n\}_{n=1}^N $ における観測値を予測することが目標だった。
  • 今回は観測点をある点に固定し、繰り返し観測値$t$を得るものとする。
  • ある値を中心に散らばったデータ群 $ \{t_n\}_{n=1}^N $ が得られる。

平均μ、標準偏差σの正規分布である特定のデータ $t=t_n$ が得られる確率 (ほぼ再掲)

$$ N(t_n\ |\ μ, σ^2) = \dfrac{1}{\sqrt{2πσ^2}}e^{-\dfrac{1}{2σ^2}\{t_n-μ\}^2} \tag{3.23} $$

すべての観測点 $\{t_n\}_{n=1}^N$ について合わせて考えると、 一連のデータ群が得られる確率Pは、それぞれの確率の積になる。

$$ \begin{eqnarray} P &=& N(t_1\ |\ μ, σ^2) × N(t_1\ |\ μ, σ^2) × \cdot\cdot\cdot × N(t_N\ |\ μ, σ^2) \nonumber \\ &=& \Pi_{n=1}^NN(t_n\ |\ μ, σ^2) \tag{3.24} \end{eqnarray} $$

この確率はμとσの2つのパラメータに依存しており、μとσを変数とする尤度関数Pが得られたことになる。

このPを最大化するようなμとσを求め、平均および標準偏差の推定値として採用する。

計算

まずは (3.23) を (3.24) に代入して整理する。

$$ \begin{eqnarray} P &=& \Pi_{n=1}^NN(t_n\ |\ μ, σ^2) \nonumber \\ &=& \Pi_{n=1}^N\left(\dfrac{1}{\sqrt{2πσ^2}}e^{-\dfrac{1}{2σ^2}\{t_n-μ\}^2}\right) \nonumber \\ &=& \Pi_{n=1}^N\left(\dfrac{1}{\sqrt{2πσ^2}}\right)\Pi_{n=1}^N\left(e^{-\dfrac{1}{2σ^2}\{t_n-μ\}^2}\right) \nonumber \\ &=& \Pi_{n=1}^N\left\{\left(\dfrac{1}{2πσ^2}\right)^\dfrac{1}{2}\right\}\Pi_{n=1}^N\left\{exp\left({-\dfrac{1}{2σ^2}\{t_n-μ\}^2}\right)\right\} \nonumber \\ &=& \left(\dfrac{1}{2πσ^2}\right)^{\dfrac{N}{2}}exp\left({-\dfrac{1}{2σ^2}\sum_{n=1}^{N}\{t_n-μ\}^2}\right) \tag{3.25} \\ \end{eqnarray} $$

(3.25)にはパラメータσが $\dfrac{1}{σ^2}$ という形のみで含まれるので、 $$ β = \dfrac{1}{σ^2} \tag{3.26} $$ とおく。

計算を単純にするために対数をとり、対数尤度関数lnPを最大化する。

$$ \begin{eqnarray} P &=& \left(\dfrac{1}{2πσ^2}\right)^{\dfrac{N}{2}}exp\left({-\dfrac{1}{2σ^2}\sum_{n=1}^{N}(t_n-μ)^2}\right) \nonumber \\ &=& \left(\dfrac{β}{2π}\right)^{\dfrac{N}{2}}exp\left({-\dfrac{β}{2}\sum_{n=1}^{N}(t_n-μ)^2}\right) \nonumber \\ lnP &=& ln\left\{\left(\dfrac{β}{2π}\right)^{\dfrac{N}{2}}exp\left({-\dfrac{β}{2}\sum_{n=1}^{N}(t_n-μ)^2}\right)\right\} \nonumber \\ &=& ln\left(\dfrac{β}{2π}\right)^{\dfrac{N}{2}} + ln\ exp\left({-\dfrac{β}{2}\sum_{n=1}^{N}(t_n-μ)^2}\right) \nonumber \\ &=& \dfrac{N}{2}ln\left(\dfrac{β}{2π}\right) + ln\ exp\left({-\dfrac{β}{2}\sum_{n=1}^{N}(t_n-μ)^2}\right) \nonumber \\ &=& \dfrac{N}{2}lnβ - \dfrac{N}{2}ln2π -\dfrac{β}{2}\sum_{n=1}^{N}(t_n-μ)^2 \tag{3.27} \\ \end{eqnarray}$$

対数は単調増加なのでlnPが最大になることと、Pが最大になることは同値。 対数なので傾き0で最大となるので、対数尤度関数を最大にする $(μ, β)$ は次の条件で決まる。

$$ \dfrac{∂(lnP)}{∂μ} = 0 \tag{3.28} $$

$$ \dfrac{∂(lnP)}{∂β} = 0 \tag{3.29} $$

(3.28)は(3.27)を代入して、

$$ \begin{eqnarray} \dfrac{∂\left(\dfrac{N}{2}lnβ - \dfrac{N}{2}ln2π -\dfrac{β}{2}\sum_{n=1}^{N}(t_n-μ)^2 \right)}{∂μ} &=& 0 \nonumber \\ \dfrac{∂\left(-\dfrac{β}{2}\sum_{n=1}^{N}(t_n-μ)^2 \right)}{∂μ} &=& 0 \nonumber \\ \dfrac{∂\left(-\dfrac{β}{2}\sum_{n=1}^{N}(t_n^2 - 2t_nμ + μ^2) \right)}{∂μ} &=& 0 \nonumber \\ \dfrac{∂\left(\dfrac{β}{2}\sum_{n=1}^{N}t_n^2 + β\sum_{n=1}^{N}t_nμ + -\dfrac{β}{2}\sum_{n=1}^{N}μ^2 \right)}{∂μ} &=& 0 \nonumber \\ \dfrac{∂\left(β\sum_{n=1}^{N}t_nμ -\dfrac{β}{2}\sum_{n=1}^{N}μ^2 \right)}{∂μ} &=& 0 \nonumber \\ β\sum_{n=1}^{N}t_n -β\sum_{n=1}^{N}μ &=& 0 \nonumber \\ β\sum_{n=1}^{N}t_n -βNμ &=& 0 \nonumber \\ β(\sum_{n=1}^{N}t_n -Nμ) &=& 0 \tag{3.30} \\ \sum_{n=1}^{N}t_n &=& Nμ \nonumber \\ μ &=& \dfrac{1}{N}\sum_{n=1}^{N}t_n \tag{3.31} \\ \end{eqnarray}$$

この時、(3.31)の右辺は観測をN回試行した際の平均値(標本平均)を表す。

一方、(3.29)の左辺は(3.27)を用いると、

$$ \begin{eqnarray} \dfrac{∂\left(\dfrac{N}{2}lnβ - \dfrac{N}{2}ln2π -\dfrac{β}{2}\sum_{n=1}^{N}(t_n-μ)^2 \right)}{∂β} &=& 0 \nonumber \\ \dfrac{∂\left(\dfrac{N}{2}lnβ -\dfrac{β}{2}\sum_{n=1}^{N}(t_n-μ)^2 \right)}{∂β} &=& 0 \nonumber \\ \dfrac{N}{2β} -\dfrac{1}{2}\sum_{n=1}^{N}(t_n-μ)^2 &=& 0 \nonumber \\ \dfrac{N}{2β} &=& \dfrac{1}{2}\sum_{n=1}^{N}(t_n-μ)^2 \nonumber \\ \dfrac{1}{β} &=& \dfrac{1}{N}\sum_{n=1}^{N}(t_n-μ)^2 \tag{3.33} \\ \end{eqnarray} $$

(3.26)より、$ β = \dfrac{1}{σ^2} $ なので、$ σ^2 = \dfrac{1}{β} $ と変形でき、(3.33)は $ σ^2 = \dfrac{1}{N}\sum_{n=1}^{N}(t_n-μ)^2 $ となる。

つまり、観測データ $ \{t_n\}_{n=1}^N $ の分散(標本分散)になっており、分散$σ^2$の推定値として、観測データの標本分散を採用する事を意味する。

3.2.3 推定量の評価方法 (一致性と普遍性)

  • 最尤推定法は必ずしも正解(真の母数)を与えるものではない。
    • これはどんな推測アプローチも同様。
  • データサイエンスの目的は「過去のデータから未来を予測する」
    • 言い換えれば「有限個のデータからその背後の一般的な事実を推測する」
  • 様々な推測方法があり、それぞれにメリット/デメリットがある。
    • 機械学習の結果を鵜呑みにするのではなく、モデルの汎化能力を評価しよう。
    • テストセットによる検証やクロスバリデーションなど。
  • (図3.11)からは「標準偏差$σ$(分散$σ^2$)の推定値は実際の値よりも小さくなりがち」という事が観測できる
    • ならば推定値を少し大きくすればいい?
      • どの程度大きくすればいいか?
      • 「推定量の評価」をする。

推定量

  • 何らかの陸地に基づいて、推定値を計算する方法が得られた場合、その計算方法を「推定量」と呼ぶ。
  • 具体的なデータから計算された値(推定値)と区別するために、こういう呼び方になっている。
  • 様々な推定量の中でも「一致性」と「不偏性」を持つものは性質のよい推定量。

一致性

  • 観測するデータ数Nを増やすと、$μ$と$σ^2$の推定値は真の母数に近づく。
  • このようにデータ数を大きくすることで真の値に近づくことを「一致性」と呼ぶ。

  • 一般に一致性を持つ推定量を「一致推定量」と呼ぶ。

不偏性

  • 偏りが無い、という性質。
  • 集める推定値の数を増やしていくと「μの推定値の平均」は真の母数である0に近づく。
  • 「何度も推定を繰り返した際に、推定値の平均が真の母数に近づく」という性質。
  • 大きい方にも外れるが、同じだけ小さい方にも外れる。→平均すると変わらない。

今週はここまで。

ITエンジニアのための機械学習理論入門 3-1

途中の計算式がわからないので補完しながら予習メモ。(といいつつ結構書いてるけど大丈夫なんだろうか...)

3.1 確率モデルの利用

  1. 最尤推定
  2. 「あるデータが得られる確率」を設定して、そこから最良のパラメータを決定する。
  3. パラメトリックモデルの3つのステップ
    1. パラメータを含むモデル(数式)を設定する。
    2. パラメーターを評価する基準を定める。
    3. 最良の評価を与えるパラメータを決定する。

3.1.1 「データ発生確率」の設定

  • 一般に回帰分析では、データの背後にある関数関係を推測する。
  • しかし最小二乗法で見たように「すべての点を正確に通る関数」は未来予測には役立たない。
    • オーバーフィッティングのため。
  • どの程度の誤差を持つのか?を含めて分析する必要がある。
  • 最小二乗法は、誤差による散らばりの中心部分を表す。
  • ここにデータが持つ誤差(散らばり具合)を含めて推測するのが確率モデル。
    • ビジネスでは「どの程度の範囲で予想が外れそうか」が大事。
  • このデータの背後にはM次多項式の関数があり、さらに標準偏差σの誤差が含まれている。と仮定する。
    • M次多項式は最小二乗法のところと同じ。
    • 誤差についての仮定が新しい。

まず、最小二乗法と同じように、特徴変数$x$と目的変数$t$の間にはM次多項式の関係があるとする。

$$ \begin{eqnarray} f(x) &=& ω_0 + ω_1x + ω_2x^2 + \cdot\cdot\cdot + ω_Mx^M \nonumber \\ &=& \sum_{m=0}^{M}ω_mx^m \tag{3.1} \end{eqnarray} $$

その上で、「観測点$x_n$における観測値$t$は、$f(x_n)$を中心として、およそ$f(x_n)\pmσ$の範囲に散らばる」と考える。

分布

$μ$を中心として、およそ$μ\pmσ$の範囲に散らばる乱数は、平均$μ$、分散$σ^2$の正規分布で表現できる。

正規分布は釣鐘型の確率で散らばる乱数なので、この値は次の関数で与えられる。

$$ N(x\ |\ μ,σ^2) = \dfrac{1}{\sqrt{2πσ^2}}e^{-\dfrac{1}{2σ^2}(x-μ)^2} \tag{3.2} $$

ここで、乱数で散らばる値$t$と、散らばりの中心$f(x_n)$を代入して考えると、

$$ N(t\ |\ f(x_n), σ^2) = \dfrac{1}{\sqrt{2πσ^2}}e^{-\dfrac{1}{2σ^2}{t-f(x_n)}^2} \tag{3.3} $$

となる。

ここで言っている観測値$t$は「次に観測される値」の事を言っていて、 トレーニングセットとして与えられている$t_n$(すでに観測された特定の値)とは違うことに注意。

$t_0$を具体的な値として、$t = t_0$が得られるような確率は次の式で求められる。

$$ N(t_0\ |\ f(x_n), σ^2) = \dfrac{1}{\sqrt{2πσ^2}}e^{-\dfrac{1}{2σ^2}{t_0-f(x_n)}^2} \tag{3.4} $$

ここまででパラメトリックメソッドの(1)ステップが終わり。


このモデルに含まれるパラメータ

  1. (3.1)に含まれる係数 $ {ω_m}_{m=0}^M $
  2. (3.3)に含まれる標準偏差$σ$ の2つ。

次のステップでは、これらのパラメーター値を評価する基準を設定して、最良のパラメーター値を決定する。

「データに含まれる誤差」は何故正規分布か?

  • 実は、「正規分布」でも「t分布」でも「ロジスティック分布」でも「一様分布」でもなんでもいい。
  • でも計算が楽ですぐに評価できるという観点から、「正規分布」を採用する。
  • 仮説・検証を高速に繰り返せることはとてもアドバンテージ。
  • うまく行かなかったら、「なぜその仮説ではうまくいかないのか」を解明することで先に進める。

尤度関数によるパラメータの評価

(3.1)と(3.3)を使って、「トレーニングセットに含まえるデータ ${(x_n,t_n)}_{n=1}^N$ が得られる確率」を計算してみる。

ある特定の観測点$x_n$について考えると、そこで$t_n$が得られる確率は以下の式で表される。

$$ N(t_n\ |\ f(x_n), σ^2) = \dfrac{1}{\sqrt{2πσ^2}}e^{-\dfrac{1}{2σ^2}{t_n-f(x_n)}^2} \tag{3.7} $$

すべての観測点 $ {x_n}_{n=1}^N $ について合わせて考えると、全体としてトレーニングセット ${(x_n,t_n)}_{n=1}^N$ が得られる確率Pは、それぞれの確率の積になる。

$$ \begin{eqnarray} P &=& N(t_1\ |\ f(x_1), σ^2) x N(t_1\ |\ f(x_1), σ^2) x \cdot\cdot\cdot x N(t_N\ |\ f(x_N), σ^2) \nonumber \\ &=& Π_{n=1}^NN(t_n\ |\ f(x_n), σ^2) \tag{3.8} \end{eqnarray} $$

この確率はパラメーター ($ {ω_m}_{m=0}^M $ と $σ$) によって値がかわるので、これらのパラメーターの関数と考えられる。

ここまででパラメトリックメソッドの(2)ステップが終わり。


このように、「トレーニングセットのデータが得られる確率」をパラメーターの関数とみなしたものを 尤度関数と呼ぶ。

そして

「観測されたデータ(トレーニングセット)は、最も発生確率が高いデータに違いない」

という仮説を立てる。

これは適当な仮説で、この仮説が正しいという保証はないが、自分は「トレーニングセットでレアなデータを得られるほど運はよくないはずだ」と考えているようなもの。

とにかくこの仮説が正しいとして、(3.8)の確率Pが最大になるようにパラメータを決定する手法を「最尤推定法」と呼ぶ。

計算

(3.8)のPを最大化するパラメーターを求めるために、まず(3.7)を(3.8)に代入する。

$$ \begin{eqnarray} P &=& Π_{n=1}^N{\dfrac{1}{\sqrt{2πσ^2}}}e^{-\dfrac{1}{2σ^2}{t_n-f(x_n)}^2} \nonumber \\ &=& \left({\dfrac{1}{2πσ^2}}\right)^{\dfrac{N}{2}}Π_{n=1}^Ne^{-\dfrac{1}{2σ^2}{t_n-f(x_n)}^2} \nonumber \\ &=& \left({\dfrac{1}{2πσ^2}}\right)^{\dfrac{N}{2}}e^{-\dfrac{1}{2σ^2}\sum_{n=1}^{N}{t_n-f(x_n)}^2} \nonumber \\ &=& \left({\dfrac{1}{2πσ^2}}\right)^{\dfrac{N}{2}}exp \left[ {-\dfrac{1}{2σ^2}\sum_{n=1}^{N}{t_n-f(x_n)}^2} \right] \tag{3.9} \end{eqnarray} $$

ここで、(3.9)の指数関数の中を見ると最小二乗法で仕様した二乗誤差$E_D$と同じものがふくまれていることがわかる。

二乗誤差の定義は改めて以下の通り。

$$ E_D = \dfrac{1}{2}\sum_{n=1}^{N}{t_n-f(x_n)}^2 \tag{3.10} $$

これを用いると、尤度関数は次のように表現できる。

$$ \begin{eqnarray} P &=& \left( \dfrac{1}{2πσ^2} \right)^\dfrac{N}{2} exp \left( -\dfrac{1}{σ^2}E_D \right) \nonumber \\ &=& \left( \dfrac{1}{2πσ^2}\right)^\dfrac{N}{2} e^{-\dfrac{1}{σ^2}E_D} \tag{3.11} \end{eqnarray} $$

ここでパラメーターに対する依存性を確認する。 (3.11)にはパラメーターσが、$\dfrac{1}{σ^2}$という形でしか存在しない。 そこで、後の計算を楽にするために、

$$ β = \dfrac{1}{σ^2} \tag{3.12} $$

とおいて、σの代わりにβをパラメーターとみなして計算を進める。 また、二乗誤差$E_D$は、多項式の係数$w=(w_0,\cdot\cdot\cdot,w_M)^T$に依存している。 従って、パラメーター $(β,w)$ に対する依存性を明示すると、以下のようになる。

$$ P(β,w) = \left( \dfrac{β}{2π}\right)^\dfrac{N}{2} e^{-βE_D(w)} \tag{3.13} $$

これを最大化する $(β,w)$ を求めれば良い。ここでさらに計算を簡単にするために、 Pの対数lnP (自然対数$log_ex = lnx$) を最大化することにする。

$$ \begin{eqnarray} lnP(β,w) &=& ln \left( \left( \dfrac{β}{2π}\right)^\dfrac{N}{2} e^{-βE_D(w)} \right) \nonumber \\ &=& ln \left( \dfrac{β}{2π}\right)^\dfrac{N}{2} + lne^{-βE_D(w)} \nonumber \\ &=& \dfrac{N}{2}ln \dfrac{β}{2π} -βE_D(w) \nonumber \\ &=& \dfrac{N}{2}(lnβ - ln2π) -βE_D(w) \nonumber \\ &=& \dfrac{N}{2}lnβ - \dfrac{N}{2}ln2π -βE_D(w) \tag{3.14} \end{eqnarray} $$

対数は単調増加関数なので、lnPが最大になることと、Pが最大になることは同値。 一般に $lnP$ を対数尤度関数と呼ぶ。

対数尤度関数を最大にする $(β,w)$ は、次の条件で決まる。

$$ \begin{eqnarray} & \dfrac{∂(lnP)}{∂ω_m}& = 0   (m = 0, \cdot\cdot\cdot, M) \tag{3.15} \\ \nonumber \\ & \dfrac{∂(lnP)}{∂β}& = 0 \tag{3.16} \end{eqnarray} $$

まず(3.14)を(3.15)に代入すると、次が得られる。

$$ \begin{eqnarray} \dfrac{∂\left(\dfrac{N}{2}lnβ - \dfrac{N}{2}ln2π -βE_D(w)\right)}{∂ω_m} &=& 0 \nonumber \\ \dfrac{∂}{∂ω_m}\left(\dfrac{N}{2}lnβ \right) - \dfrac{∂}{∂ω_m}\left(\dfrac{N}{2}ln2π \right) - \dfrac{∂}{∂ω_m}\left( βE_D(w) \right) &=& 0 \nonumber \\ \dfrac{∂}{∂ω_m}\left( βE_D(w) \right) &=& 0 \nonumber \\ β\dfrac{∂E_D}{∂ω_m}(w) &=& 0 \nonumber \\ \dfrac{∂E_D}{∂ω_m} &=& 0   (m=0,\cdot\cdot\cdot,M) \tag{3.17} \ \end{eqnarray} $$

これは二乗誤差を最小にする条件と全く同じなので、2.1.3と同じ計算を行って同じ結論が得られる。 つまり多項式の係数 $ {ω_m}_{m=0}^M $ は最小二乗法と同じ値に決まる。

一方、(3.14)を(3.16)に代入すると、

$$ \begin{eqnarray} \dfrac{∂\left(\dfrac{N}{2}lnβ - \dfrac{N}{2}ln2π -βE_D(w)\right)}{∂β} &=& 0 \nonumber \\ \dfrac{N}{2β} -E_D &=& 0 \nonumber \\ \dfrac{1}{β} &=& \dfrac{2E_D}{N} \tag{3.18} \ \end{eqnarray} $$

これに(3.12)を代入すると、標準偏差$σ$を決定する式が得られる。

$$ \begin{eqnarray} (β = \dfrac{1}{σ^2})    \dfrac{1}{\dfrac{1}{σ^2}} &=& \dfrac{2E_D}{N} \nonumber \\ σ^2 &=& \dfrac{2E_D}{N} \nonumber \\ σ &=& \sqrt{\dfrac{2E_D}{N}} = E_{R M S} \tag{3.19} \ \end{eqnarray} $$

この $ E_{R M S} $ は(2.20)で定義した平方根二乗誤差となる。

つまり、(3.19)は、トレーニングセットに含まれるデータの 「多項式で推定される値$f(x_n)$に対する平均的な誤差」を標準偏差$σ$の推定値として採用することを意味する。

対数関数の微分

$$ \dfrac{d}{dx}lnx = \dfrac{1}{x} $$


今週はとりあえずここまで。

ITエンジニアのための機械学習理論入門 2-1

毎週少しずつ読んでいるんですが、なかなか理解が遅いので記事としてまとめてみます。

完全に自分用メモ。

2.1.3 数学徒の小部屋

$$ 誤差 E_D = \dfrac{1}{2} \sum_{n=1}^{N} (\sum_{m'=0}^{M}ω_{m'}x_{n}^{m'} - t_{n})^2 \tag{2.4} $$

ここで(2.4)を最小にする $ {ω_m}_{m=0}^M $ を決定する。

$$ \dfrac{∂E_{D}}{∂ω_{m}} = 0     (m = 0, …, M)\tag{2.5} $$

に、(2.4)を代入すると以下の式になる。

$$ \dfrac{∂{ \dfrac{1}{2} \sum_{n=1}^{N} (\sum_{m'=0}^{M}ω_{m'}x_{n}^{m'} - t_{n})^2 }}{∂ω_m} = 0 $$

ここで、$ E_D $に代入するmは$ \sum_{m=0}^{M} $の、冪級数を表すmなので、 偏微分したいωmのmとは違う。

よって、冪級数を表すmをm'として置く。 そこで、

$$ \dfrac{1}{2} \sum_{n=1}^{N} (\sum_{m'=0}^{M}ω_{m'}x_{n}^{m'} - t_{n})^2 $$

を$ ω_m $で偏微分する。

実際の偏微分は、つい2乗を展開したくなるのが これは間違いで、合成関数の偏微分を利用する。


【合成関数の偏微分

$$ \dfrac{∂f(g(x, y))}{∂x} = f'(g(x, y)) \cdot \dfrac{∂g(x, y)}{∂x} $$


ここで、

$$ \begin{eqnarray} g(x) &=& \sum_{m'=0}^{M}ω_{m'}x_{n}^{m'} - tn \nonumber \\ f(x) &=& g(x)^2 \nonumber \end{eqnarray} $$

と置くと、

$$ \begin{eqnarray} f'(x) &=& 2g(x) \nonumber \\ &=& 2(\sum_{m'=0}^{M}ω_{m'}x_n^{m'} - tn) \nonumber \end{eqnarray} $$

となる。

また、$ \dfrac{∂g(x, y)}{∂x} $ は $ \dfrac{∂(\sum_{m'=0}^{M}ω_{m'}x_n^{m'} - tn)}{∂ω_m} = x_n^m $ となる。

従って $ \dfrac{∂{ \dfrac{1}{2} \sum_{n=1}^{N} ( \sum_{m'=0}^{M} ω_{m'} x_n^{m'} - t_n)^2 }}{∂ω_m} $ は

$$ \begin{eqnarray} \dfrac{1}{2} \sum_{n=1}^{N} \cdot 2( \sum_{m'=0}^{M}ω_{m'}x_n^{m'} - t_n) \cdot x_n^m &=& 0 \nonumber \\ \sum_{n=1}^{N}(\sum_{m'=0}^{M}ω_{m'}x_n^{m'} - t_n) \cdot x_n^m &=& 0 \nonumber \tag{2.7} \end{eqnarray} $$

となる。

やや作為的だが、これを次のように変形する。

$$ \sum_{m'=0}^{M}ω_{m'} \sum_{n=1}^{N}x_n^{m'}x_n^m - \sum_{n=1}^{N}t_nx_n^m = 0 \tag{2.8} $$

ここで、$ x_n^m $ を (n, m)成分とするNx(M+1)行列φを用いると、これは行列形式で書き直せる。

$$ w^Tφ^Tφ-t^Tφ = 0 \tag{2.9} $$

(2.8)から(2.9)の変換の正しさを逆方向の式変換で考える。

(2.9)の $ w^Tφ^Tφ $ に注目する。

$$ w = \left( \begin{array}{c} ω_0 \\ ω_1 \\ ︙ \\ ω_M \end{array} \right),      φ = \left( \begin{array}{c} x_1^0 & x_1^1 & \cdots & x_1^M \\ x_2^0 & x_2^1 & \cdots & x_2^M \\ \vdots & \vdots & \ddots & \vdots \\ x_N^0 & x_N^1 & \cdots & x_N^M \end{array} \right),      φ^T = \left( \begin{array}{c} x_1^0 & x_2^0 & \cdots & x_N^0 \\ x_1^1 & x_2^1 & \cdots & x_N^1 \\ \vdots & \vdots & \ddots & \vdots \\ x_1^M & x_2^M & \cdots & x_N^M \end{array} \right) \tag{2.10} $$

具体的な数値で、N = 2, M = 2と置くと、

$$ w = \left( \begin{array}{c} ω_0 \\ ω_1 \\ ω_2 \end{array} \right),      φ = \left( \begin{array}{c} x_1^0 & x_1^1 & x_1^2 \\ x_2^0 & x_2^1 & x_2^2 \end{array} \right),      φ^T = \left( \begin{array}{c} x_1^0 & x_2^0 \\ x_1^1 & x_2^1 \\ x_1^2 & x_2^2 \end{array} \right) $$

ここで、$ φ^Tφ $ を計算すると、

$$ \begin{eqnarray} φ^Tφ &=& \left( \begin{array}{c} x_1^0 & x_2^0 \\ x_1^1 & x_2^1 \\ x_1^2 & x_2^2 \end{array} \right) \left( \begin{array}{c} x_1^0 & x_1^1 & x_1^2 \\ x_2^0 & x_2^1 & x_2^2 \end{array} \right) \nonumber\\ &=& \left( \begin{array}{c} x_1^0x_1^0 + x_2^0x_2^0 & x_1^0x_1^1 + x_2^0x_2^1 & x_1^0x_1^2 + x_2^0x_2^2 \\ x_1^1x_1^0 + x_2^1x_2^0 & x_1^1x_1^1 + x_2^1x_2^1 & x_1^1x_1^2 + x_2^1x_2^2 \\ x_1^2x_1^0 + x_2^2x_2^0 & x_1^2x_1^1 + x_2^2x_2^1 & x_1^2x_1^2 + x_2^2x_2^2 \end{array} \right) \nonumber \\ &=& \left( \begin{array}{c} \sum_{n=1}^{N}x_n^0x_n^0 & \sum_{n=1}^{N}x_n^0x_n^1 & \sum_{n=1}^{N}x_n^0x_n^2 \\ \sum_{n=1}^{N}x_n^1x_n^0 & \sum_{n=1}^{N}x_n^1x_n^1 & \sum_{n=1}^{N}x_n^1x_n^2 \\ \sum_{n=1}^{N}x_n^2x_n^0 & \sum_{n=1}^{N}x_n^2x_n^1 & \sum_{n=1}^{N}x_n^2x_n^2 \end{array} \right) \nonumber \\ \end{eqnarray} $$

次に $ w^Tφ^Tφ $ を計算すると、

$$ \begin{eqnarray} w^Tφ^Tφ &=& \left( \begin{array}{c} ω_0 & ω_1 & ω_2 \end{array} \right) \left( \begin{array}{c} \sum_{n=1}^{N}x_n^0x_n^0 & \sum_{n=1}^{N}x_n^0x_n^1 & \sum_{n=1}^{N}x_n^0x_n^2 \\ \sum_{n=1}^{N}x_n^1x_n^0 & \sum_{n=1}^{N}x_n^1x_n^1 & \sum_{n=1}^{N}x_n^1x_n^2 \\ \sum_{n=1}^{N}x_n^2x_n^0 & \sum_{n=1}^{N}x_n^2x_n^1 & \sum_{n=1}^{N}x_n^2x_n^2 \end{array} \right) \nonumber \\ &=& \left( \begin{array}{c} ω_0(\sum_{n=1}^{N}x_n^0x_n^0) + ω_1(\sum_{n=1}^{N}x_n^1x_n^0) + ω_2(\sum_{n=1}^{N}x_n^2x_n^0) & ω_0(\sum_{n=1}^{N}x_n^0x_n^1) + ω_1(\sum_{n=1}^{N}x_n^1x_n^1) + ω_2(\sum_{n=1}^{N}x_n^2x_n^1) & ω_0(\sum_{n=1}^{N}x_n^0x_n^2) + ω_1(\sum_{n=1}^{N}x_n^1x_n^2) + ω_2(\sum_{n=1}^{N}x_n^2x_n^2) \end{array} \right) \nonumber \\ &=& \left( \begin{array}{c} \sum_{m'=0}^{N}ω_{m'}\sum_{n=1}^{N}x_n^{m'}x_n^0 & \sum_{m'=0}^{N}ω_{m'}\sum_{n=1}^{N}x_n^{m'}x_n^1 & \sum_{m'=0}^{N}ω_{m'}\sum_{n=1}^{N}x_n^{m'}x_n^2 \end{array} \right) \nonumber \\ &=& \sum_{m'=0}^{N}ω_{m'}\sum_{n=1}^{N}x_n^{m'}x_n^m     (m=0,...,M) \nonumber \\ \end{eqnarray} $$

よって、(2.8)は(2.9)に式変換できる。

話は戻って、(2.9) $ w^Tφ^Tφ-t^Tφ = 0 $を転置を取ってwについて式変換すると、

$$ \begin{eqnarray} w &=& \dfrac{φ^Tt}{φ^Tφ} \nonumber \\ &=& {(φ^Tφ)}^{-1}φ^Tt \tag{2.11} \end{eqnarray} $$

φとtはそれぞれトレーニングセットに含まれる観測データから決まるものなので、 (2.11)は、与えられたトレーニングセットを用いて多項式の係数wを決定する式になっている。

(2.11)の $ φ^Tφ $ は逆行列を持つのか?という確認。

$E_D$の2階偏微分係数を表すヘッセ行列を用いて説明する。


二階偏微分

$$ \dfrac{∂^2f}{∂x∂y} = \dfrac{∂}{∂y}\left(\dfrac{∂f}{∂x}\right) = \dfrac{∂}{∂x}\left(\dfrac{∂f}{∂y}\right) = \dfrac{∂^2f}{∂y∂x} $$


ヘッセ行列Hは、次の成分を持つ(M+1)x(M+1)の正方行列となっている。

$$ H_{mm'} = \dfrac{∂^2E_D}{∂ω_m∂ω_{m'}}       (m, m' = 0, ..., M) \tag{2.12} $$

$$ \begin{eqnarray} H_{mm'} &=& \dfrac{∂^2E_D}{∂ω_m∂ω_{m'}} \nonumber \\ &=& \dfrac{∂}{∂ω_{m'}}\left(\dfrac{∂E_D}{∂ω_{m}}\right) \nonumber \\ &=& \dfrac{∂}{∂ω_{m'}}\left(\sum_{m'=0}^{M}ω_{m'} \sum_{n=1}^{N}x_n^{m'}x_n^m - \sum_{n=1}^{N}t_nx_n^m\right) \nonumber \\ &=& \sum_{n=1}^{N}x_n^{m'}x_n^m \nonumber \\ \end{eqnarray} $$ $※ \dfrac{∂E_D}{∂ω_m} はすでに行っている計算なので結果を代入している。 $

(2.10)を用いると、(2.11)で逆行列を取っている部分の行列がヘッセ行列に一致することがわかる。

$$ H = φ^Tφ \tag{2.14} $$

この時、$M+1 \le N $、 すなわち係数の個数 M+1がトレーニングセットのデータ数N以下であれば、 ヘッセ行列は正定値であることがわかる。

正定値というのは、任意のベクトルu ($u \ne 0$)に対して、$u^THu \gt 0 $が成立することを言う。

今の例では、以下の式になる。 $$ u^TH_u = u^Tφ^Tφu = {|| φu ||}^2 > 0 $$

この不等号が成立するのは、 $φu \ne 0$の場合に限るが、$φ$の定義(2.10)を思い出すと、$φu = 0$は 要素数がM+1のベクトルuに対するN本の斉次な連立一次方程式になるので、$M+1 \le N$の場合、 自明でない解$u \ne 0$を見つけることはできない。

したがって、$u \ne 0$は必ず成立して、ヘッセ行列$φ^Tφ$は正定値となる。

さらに、正定値な行列は逆行列をもつことが証明できるので、逆行列$(φ^Tφ)^{-1}$が確かに存在し、 停留点は(2.11)で一意に決まる

そしてヘッセ行列が正定値であることから、この停留点は$E_D$の極小値を与えることが示される。 これで(2.11)は$E_D$を最小にする唯一のwを与える事が示された。

一方、$ M+1 \gt N $、すなわち、係数の個数がトレーニングセットのデータを超える場合は、 ヘッセ行列は半正定値($u^THu \ge 0$)となるため、$E_D$を最小にするwは複数存在し、一意に決定されなくなる。