#0. 概要
多項式を用いてデータセットを近似する方法について説明する。今回は最急降下法等を使わず、最尤推定による近似を行う。
#1.0. 理論
理論について説明する。
##1.1. 推論
多項式基底による推論方法について説明する。
まず基底関数を以下のように多項式とする。
{\phi _{i}\left( x\right) =x^{i} (i=0,\ldots ,M-1)}
次に多項式係数$w$との線形結合によって、あるデータセットを$f(x)$で表現できるとする。
{f(x)=w^{T}\phi(x)=\sum ^{n-1}_{i=0}W_{i}\phi _{i}(x)}
ここに、上記式だけではデータセットを表現しきれないエラーも含めて式を展開すると以下になる。
{f(x)=w^{T}\phi(x)+\epsilon=\sum ^{n-1}_{i=0}W_{i}\phi _{i}(x) + \epsilon}\\
={f(x)=w_0+ w_1x^1 +\cdots + w_ix^i+ \epsilon}
$i$が次元数である。基本的に高次になればなるほど、複雑な関数を表現できるようになる。
しかし、データセットの複雑性に対して不必要に高い次数を採用した場合、オーバフィッティングする可能性が生じてくる。
#1.2. 多項式係数の推定
次に、多項式をあるデータセットにフィッティングさせることを考える。その場合、データセットと多項式による近似がどれくらい離れているかの評価関数(もとい損失関数)が必要である。一般的に損失関数は最小二乗法を用いる。
{E(w)=\dfrac {1}{2}\sum ^{N}_{n=1}\left( f\left( x_n\right) -t_{n}\right) ^{2}}
ここで$E(w)$があるデータセットと多項式による近似の距離(誤差)である。$t$が教師データである。$n$はデータセットのデータ数を示している。
すなわち、この誤差が全て0であれば、あるデータセットに対して多項式がフィッティングできたということになる。
{\frac{\partial E(w)}{\partial w_m}=\sum^N_{n=1}\phi_m(x_n)\left(\sum^M_{j=1}w_j\phi_j(x_n)-t_n\right)=0 \hspace{1.0em}(m=1,\ldots ,M-1)
}
偏微分の結果、以下のようになる
\frac{\partial E(w)}{\partial w_0}=\sum^N_{i=1}\left\{(w_0 + w_1x^1_{i}+ w_2x^2_{i} +\cdots + w_mx^m_{i})-t_i\right\}\cdot 1=0\\
\frac{\partial E(w)}{\partial w_1}=\sum^N_{i=1}\left\{(w_0 + w_1x^1_{i}+ w_2x^2_{i} +\cdots + w_mx^m_{i})-t_i\right\}x^1_{i}=0\\
\frac{\partial E(w)}{\partial w_2}=\sum^N_{i=1}\left\{(w_0 + w_1x^1_{i}+ w_2x^2_{i} +\cdots + w_mx^m_{i})-t_i\right\}x^2_{i}=0\\
\frac{\partial E(w)}{\partial w_m}=\sum^N_{i=1}\left\{(w_0 + w_1x^1_{i}+ w_2x^2_{i} +\cdots + w_mx^m_{i})-t_i\right\}x^m_{i}=0
各$w$は$E×x^m$で表現できることが確認できた。
次に各式を$t$周りに展開する。
まず上述の$\partial w_0$を出発式とする。偏微分の結果が全て0であればフィッティングできているということになる。
\frac{\partial E(w)}{\partial w_0}=\sum^N_{i=1}\left\{(w_0 + w_1x^1_{i}+ w_2x^2_{i} +\cdots + w_mx^m_{i})-t_i\right\}\cdot 1=0
外側の$1$を$\sum$に掛けてから展開し、$t$を移項する。
\sum^N_{i=1}1{\cdot}w_0+\sum^N_{i=1}1{\cdot}w_1x^1_{i}+\sum^N_{i=1}1{\cdot}w_2x^2_{i}+\cdots + \sum^N_{i=1}1{\cdot}w_mx^m_{i}=\sum^N_{i=1}t_i
$\sum^N_{i=1}1$は展開すると$1+1+1 \cdots$とN回続くので、$N$となる。最後に式を整理して以下が得られる。
N{\cdot}w_0+w_1\sum^N_{i=1}x^1_{i}+w_2\sum^N_{i=1}x^2_{i}+\cdots + w_m\sum^N_{i=1}x^m_{i}=\sum^N_{i=1}t_i
同様に、$\partial w_0$から$\partial w_2$まで展開すると、以下が得られる
N{\cdot}w_0+w_1\sum^N_{i=1}x^1_{i}+w_2\sum^N_{i=1}x^2_{i}+\cdots + w_m\sum^N_{i=1}x^m_{i}=\sum^N_{i=1}t_i\\
w_0\sum^N_{i=1}x^1_i+w_1\sum^N_{i=1}x^2_{i}+w_2\sum^N_{i=1}x^3_{i}+\cdots + w_m\sum^N_{i=1}x^{m+1}_{i}=\sum^N_{i=1}x_{i}t_i\\
w_0\sum^N_{i=1}x^2_i+w_1\sum^N_{i=1}x^3_{i}+w_2\sum^N_{i=1}x^4_{i}+\cdots + w_m\sum^N_{i=1}x^{m+2}_{i}=\sum^N_{i=1}x^2_{i}t_i
これを行列式として整理する。
{
\begin{bmatrix}
N & \sum^N_{i=1}x^1_{i} & \cdots & \sum^N_{i=1}x^m_{i}\\
\sum^N_{i=1}x^1_i & \sum^N_{i=1}x^2_{i} & \cdots & \sum^N_{i=1}x^{m+1}_{i} \\
\vdots & \vdots & \ddots & \vdots \\
\sum^N_{i=1}x^{m}_i & \sum^N_{i=1}x^{m+1}_{i} & \cdots & \sum^N_{i=1}x^{2m}_{i}
\end{bmatrix}
\begin{bmatrix}
w_0\\
w_1 \\
\vdots \\
w_m
\end{bmatrix}
=
\begin{bmatrix}
\sum^N_{i=1}t_i\\
\sum^N_{i=1}x_{i}t_i \\
\vdots \\
\sum^N_{i=1}x_{m}t_i
\end{bmatrix}
}
これで、勾配が0になることでフィッティングが出来ることが理解できた。後は$w$を解くだけである。
これから先は未知数$w$を推定する2つの解法を紹介する。
ヴァンデルモンドの行列式による解法
偏微分の行列式では、1乗の係数で全てのデータセットを評価し、2乗の係数で全てのデータセットを評価し・・・としていたが、究極的には勾配が0、すなわち多項式基底が$t$と一致すれば言い訳なので、以下の式でも同様である。
t_0 = (w_0 + w_1x^1_{0}+ w_2x^2_{0} +\cdots + w_mx^m_{0})\\
t_1 = (w_0 + w_1x^1_{1}+ w_2x^2_{1} +\cdots + w_mx^m_{1})\\
t_2 = (w_0 + w_1x^1_{2}+ w_2x^2_{2} +\cdots + w_mx^m_{2})\\
\vdots\\
t_i = (w_0 + w_1x^1_{i}+ w_2x^2_{i} +\cdots + w_mx^m_{i})\\
次に、この式を行列式に整理すると以下のような形になり、これはヴァンデルモンドの行列式と一致する。
また、2乗誤差を採用した偏微分方程式とも一致している。
無論、これまで計算してきた偏微分の結果（Sumが全てについた行列式）で計算しても結果は変わらないはずである。
最後に割り算によって答えを標準化するため。ここでは、簡単のため、Sumを取り除いて計算を行う。
{
\begin{bmatrix}
1 & x^1_{1} & x^2_{1} & \cdots & x^m_{i}\\
1 & x^1_{2} & x^2_{2} & \cdots & x^{m}_{i} \\
\vdots & \vdots & \vdots & \ddots & \vdots \\
1 & x^{1}_{i} & x^{2}_{i} & \cdots & x^{m}_{i}
\end{bmatrix}
\begin{bmatrix}
w_0\\
w_1 \\
w_2 \\
\vdots \\
w_m
\end{bmatrix}
=
\begin{bmatrix}
t_0\\
t_1 \\
t_2 \\
\vdots \\
t_i
\end{bmatrix}
}
すなわち、ベクトルで整理しなおすと
xw=y
次に$w=$の形に展開するため$x$を両方に掛ける。
x^Txw=x^Ty
最後に$x$を移項して、以下が得られる。
w=(x^Tx)^{-1}x^{T}y
Cramer's Ruleを用いた解法
まず、$w_m$は以下の行列式によって得られるものとする
w_m=\dfrac{det(M_m)}{det(M)}
このとき、クラメルの公式は以下によって与えられる。
すなわち、$x$は$w$の行列であるとして、以下のように整理できる。
無論、こちらもSumありでもなしでも答えは同じになるはずである。こちらはSumありで記述する。
b=\begin{bmatrix}
\sum^N_{i=1}t_i\\
\sum^N_{i=1}x_{i}t_i \\
\vdots \\
\sum^N_{i=1}x_{m}t_i
\end{bmatrix}\\
M={
\begin{bmatrix}
N & \sum^N_{i=1}x^1_{i} & \cdots & \sum^N_{i=1}x^m_{i}\\
\sum^N_{i=1}x^1_i & \sum^N_{i=1}x^2_{i} & \cdots & \sum^N_{i=1}x^{m+1}_{i} \\
\vdots & \vdots & \ddots & \vdots \\
\sum^N_{i=1}x^{m}_i & \sum^N_{i=1}x^{m+1}_{i} & \cdots & \sum^N_{i=1}x^{2m}_{i}
\end{bmatrix}
}\\
M0={
\begin{bmatrix}
\sum^N_{i=1}t_i & \sum^N_{i=1}x^1_{i} & \cdots & \sum^N_{i=1}x^m_{i}\\
\sum^N_{i=1}x_{i}t_i & \sum^N_{i=1}x^2_{i} & \cdots & \sum^N_{i=1}x^{m+1}_{i} \\
\vdots & \vdots & \ddots & \vdots \\
\sum^N_{i=1}x_{m}t_i & \sum^N_{i=1}x^{m+1}_{i} & \cdots & \sum^N_{i=1}x^{2m}_{i}
\end{bmatrix}
}\\
M1={
\begin{bmatrix}
N & \sum^N_{i=1}t_i & \cdots & \sum^N_{i=1}x^m_{i}\\
\sum^N_{i=1}x^1_i & \sum^N_{i=1}x_{i}t_i & \cdots & \sum^N_{i=1}x^{m+1}_{i} \\
\vdots & \vdots & \ddots & \vdots \\
\sum^N_{i=1}x^{m}_i & \sum^N_{i=1}x^{m+1}_{i} & \cdots & \sum^N_{i=1}x^{2m}_{i}
\end{bmatrix}
}\\
\vdots \\
M_{m}
この得られた各Mに対して行列式を解いていけば$w$が求まる。
なお、このMの求め方は非常に単純である。以下のような配列操作を行っているだけである。
このように$M$行列の列方向に$t$をスライドさせていくだけで$M_0, M_1, ...$が出来上がる。
#2.0. 例題を用いた推定
せっかくなのでクラメルの公式を用いて、解いてみる。以下をデータセットする。
x=\left\{-3, -2, -1, -0.2, 1, 3 \right\}
t=\left\{0.9, 0.8, 0.4, 0.2, 0.1, 0 \right\}
先ほどの行列式を基に、データセットから3次項までの係数$w_0, w_1, w_2$を求める。
それに先立ち、M、M0、M1、M2を計算する。
M0={
\begin{bmatrix}
\sum^N_{i=1}t_i & \sum^N_{i=1}x^1_{i} & \cdots & \sum^N_{i=1}x^m_{i}\\
\sum^N_{i=1}x_{i}t_i & \sum^N_{i=1}x^2_{i} & \cdots & \sum^N_{i=1}x^{m+1}_{i} \\
\vdots & \vdots & \ddots & \vdots \\
\sum^N_{i=1}x_{m}t_i & \sum^N_{i=1}x^{m+1}_{i} & \cdots & \sum^N_{i=1}x^{2m}_{i}
\end{bmatrix}
\begin{bmatrix}
w_0\\
w_1 \\
\vdots \\
w_m
\end{bmatrix}
=
\begin{bmatrix}
2.44 & -2.2 & 24.04\\
-4.64 & 24.04 & -8.008\\
11.808 & -8.008 & 180.0016
\end{bmatrix}
\begin{bmatrix}
w_0\\
w_1 \\
w_2
\end{bmatrix}}
M1={
\begin{bmatrix}
N & \sum^N_{i=1}t_i & \cdots & \sum^N_{i=1}x^m_{i}\\
\sum^N_{i=1}x^1_i & \sum^N_{i=1}x_{i}t_i & \cdots & \sum^N_{i=1}x^{m+1}_{i} \\
\vdots & \vdots & \ddots & \vdots \\
\sum^N_{i=1}x^{m}_i & \sum^N_{i=1}x^{m+1}_{i} & \cdots & \sum^N_{i=1}x^{2m}_{i}
\end{bmatrix}
\begin{bmatrix}
w_0\\
w_1 \\
\vdots \\
w_m
\end{bmatrix}\\
=
\begin{bmatrix}
6& 2.44 & 24.04\\
-2.2 & -4.64 & -8.008\\
24.04 & 11.808 & 180.0016
\end{bmatrix}
\begin{bmatrix}
w_0\\
w_1 \\
w_2
\end{bmatrix}
}
$M0$及び$M1$を求めてしまえば、クラメルの公式から行列を入れ替えるだけで、$M2$及び$M$を構築できる。
M2={
\begin{bmatrix}
6 &-2.2& 2.44 \\
-2.2 &24.04& -4.64\\
24.04 &-8.008& 11.808
\end{bmatrix}
\begin{bmatrix}
w_0\\
w_1 \\
w_2
\end{bmatrix}
}
M={
\begin{bmatrix}
6 &-2.2& 24.04 \\
-2.2 &24.04& -8.008\\
24.04 &-8.008& 180.0016
\end{bmatrix}
\begin{bmatrix}
w_0\\
w_1 \\
w_2
\end{bmatrix}
}
最後に求めた$M$から$w$を計算する。
w_0=\dfrac{det(M_0)}{det(M)}=\dfrac{2671.20}{11661.27}=0.2291\\
w_1=\dfrac{det(M_1)}{det(M)}=\dfrac{-1898.46}{11661.27}=-0.1628\\
w_2=\dfrac{det(M_2)}{det(M)}=\dfrac{323.76}{11661.27}=0.0278\\
よって、多項式基底のフィッティング結果は以下となる。
y=0.0278x^2-0.1628x+0.2991