目的
データ$(x_i,y_i)$が複数（$i=1 \cdots n$）与えられたとき、次のモデルのパラメータ$(a,b,c)$を線形最小二乗法により最適化する。
\[ y=a\cos(x + b) + c \]方法
- 加法定理と$a' = a\cos b$、$b' = a\sin b$で変形し、3変数連立一次方程式の形にする。
- 平均自乗誤差を最小化する問題として線形最小二乗法の問題として解く。
- 求まった最適な($a',b',c$)から$a = \sqrt{a'^2 + b'^2}$、$b = \tan^{-1}(b'/a')$として元の変数を求める。
ソースコード
Jupyter NotebookのipynbファイルをGitHubで公開している．
STEP 1: 加法定理と変数変換でパラメータの線形和に
実際には、データ$(x,y)$には観測誤差が含まれて定式の等号は成り立たないので観測誤差をまとめて$\epsilon_i$と表すと、$i$番目のデータ$(x_i,y_i)$の関係は以下のように表現できる。
\[ \epsilon_i = a\cos(x_i + b) + c - y_i \]ここでの最適化とはこの全てのデータに対して$\epsilon_i$が小さくなる$(a,b,c)$を探すことなので、その指標として平均自乗誤差$E = \sum_i^n \epsilon_i^2/n$を最小化する。加法定理で変形し、$a' = a\cos b$、$b' = a\sin b$と置くと以下のように変形できる。
\[ \begin{array}{ll} \epsilon_i &= a\cos(x_i + b) + c - y_i\\ &= a\{\cos(x_i)\cos(b) - \sin(x_i)\sin(b)\} + c - y_i\\ &= a'\cos(x_i) - b'\sin(x_i)\ + c - y_i \end{array} \]STEP 2: 誤差関数を微分する
今、平均自乗誤差$E$を最小化しようとしているので、最適な値$(a',b',c)$では$E$が極値になることが必要条件である。
\[ \frac{\partial E}{\partial a'} = 0,\ \frac{\partial E}{\partial b'} = 0,\ \frac{\partial E}{\partial c} = 0 \]以下、微分に関係ない定数$n$を左辺に移項して上式を実際に計算する。
\[ \begin{array}{ll} n\frac{\partial E}{\partial a'} &= \frac{\partial}{\partial a'}\sum_i^n \epsilon_i^2\\ &= \sum_i^n \frac{\partial}{\partial a'} \epsilon_i^2\ = \sum_i^n 2\epsilon_i \frac{\partial}{\partial a'}\epsilon_i\\ &= \sum_i^n 2\{a'\cos(x_i) - b'\sin(x_i)\ + c - y_i \}\cos(x_i)\\ &= 2\{a'\sum_i^n\cos^2(x_i) - b'\sum_i^n \sin(x_i)\cos(x_i) + c\sum_i^n \cos(x_i) - \sum_i^n y_i\cos(x_i) \} = 0\\ % n\frac{\partial E}{\partial b'} &= \sum_i^n - 2\{a'\cos(x_i) - b'\sin(x_i)\ + c - y_i \}\sin(x_i) = 0\\ &= -2\{a'\sum_i^n\cos(x_i)\sin(x_i) - b'\sum_i^n\sin^2(x_i) + c\sum_i^n \sin(x_i) - \sum_i^n y_i\sin(x_i) \} = 0\\ % n\frac{\partial E}{\partial c} &= \sum_i^n 2\{a'\cos(x_i) - b'\sin(x_i)\ + c - y_i \} = 0\\ &= 2\{a'\sum_i^n\cos(x_i) - b'\sum_i^n\sin(x_i) + cn - \sum_i^n y_i \} = 0\ \end{array} \]上記3つの方程式は、パラメータ$(a',b',c)$のみに注目すれば3変数の連立一次方程式である。大学生の線形代数の知識を使って、以下のように行列の式にすれば簡単に解ける。
\[ \left[ \begin{array}{ccc} \sum_i^n\cos^2(x_i)& -\sum_i^n \sin(x_i)\cos(x_i)& \sum_i^n \cos(x_i)\\ \sum_i^n\cos(x_i)\sin(x_i)& -\sum_i^n\sin^2(x_i)& \sum_i^n \sin(x_i)\\ \sum_i^n\cos(x_i)& -\sum_i^n\sin(x_i)& n \end{array} \right] \left[ \begin{array}{c} a'\\ b'\\ c \end{array} \right] = \left[ \begin{array}{c} \sum_i^n y_i\cos(x_i)\\ \sum_i^n y_i\sin(x_i)\\ \sum_i^n y_i \end{array} \right] \]STEP 3: 変数変換を戻して答えを求める
求まった最適な($a',b',c$)から$a = \sqrt{a'^2 + b'^2}$、$b = \tan^{-1}(b'/a')$として元の変数を求める。プログラム上では関数atan2を使うと$a'=0$の場合分けの必要がない。使い方は以下のとおり（ただし、a_=$a'$、b_=$b'$）。
b = atan2(a_, b_)
