4 最小二乗法

最小二乗法を感覚的に理解するために，最初に一次関数の最小二乗法を示す．その後，一 般的な線形最小二乗法について説明する．非線形最小二乗法は難しいので，範囲外とする．

4.1 簡単な例

最小二乗法というのは，データをある関数で最良近似する方法である．例えば，

$$(1.2,\,2.2)\quad(2.1,\,3.8)\quad(3.3,\,5.6)\quad(4.1,\,7.1)\quad(5,8.8)$$ (20)

の$(x,y)$ の実験データがあるとする．これを直線で近似$y=ax+b$ したい．どうすればよ いか?--という問題である．誤差の二乗が最小になる直線が最良近似とすることができる． これを最小二乗法(least squares method)と言う．式で表すと，誤差の二乗の和$E(a,b)$ は，

$$E(a,b)=\sum_{i=1}^n (y_i-ax_i-b)^2$$ (21)

となる． $(x_i,\,y_i)$ は，$i$ 番目のデータで，$N$ はデータの個数である．この誤差が 最小になる$a$ と$b$ を捜す．式(21)は$a$ にも$b$ にも2次式でその係 数は正の値なので最小値がある．誤差$E$ の最小値は，それぞれ偏微分した値がゼロとな るときに得ることができる．

$$\if 11 \frac{\partial E}{\partial a} \else \frac{\partial^{1} E}{\partial a^{1}}\fi =-\sum_{i=1}^N 2(y_i-ax_i-b)x_i=0 \if 11 \frac{\partial E}{\partial b} \else \frac{\partial^{1} E}{\partial b^{1}}\fi =-\sum_{i=1}^N 2(y_i-ax_i-b)=0$$ (22)

これは，$a$ と$b$ の連立方程式である．すなわち，

$$\left\{\quad \begin{aligned}a\sum_{i=1}^Nx_i^2+b\sum_{i=1}^N x_... ...y_i\\ a\sum_{i=1}^Nx_i+nb&=\sum_{i=1}^N y_i \end{aligned} \right.$$

である．これを解くと

$$a=\cfrac {N\sum_{i=1}^N x_iy_i-\sum_{i=1}^N x_i\sum_{i=1}^N y_i} {N\sum_{i=1}^N x_i^2-(\sum_{i=1}^N x_i)^2} b=\cfrac {\sum_{i=1}^Nx_i\sum_{i=1}^Ny_i-\sum_{i=1}^Nx_iy_i\sum_{i=1}^Nx_i} {N\sum_{i=1}^Nx_i^2-(\sum_{i=1}^Nx_i)^2}$$ (24)

となる．

最初に示したのデータについて計算してみると， $a=1.452119,\,b=0.708006$ となる．ゆ えに，最小二乗法による1次関数は

$$y=1.452119 x + 0.708006$$ (25)

となる．データをこの間数をプロットすると，図6のようになる．

グラフ作成ソフトウェアーは，最小二乗法によるデータのフィッティングをサポートして いるものが多い．Excelでも可能なはずである．また，私のWEBページでは web_gnuplot と称して，最小二乗法がweb上で使えるようになっている．実験データの整理に使うと良 い．というか，実験データの整理に最小二乗法を使わないことはあり得ない．

図 6: データを最小二乗法でフィット

ここでは，偏微分により最小二乗法の式を導いたが，線形代数の部分空間への射影を考え る方が簡単である．これについては，参考文献 [2]に詳しく書い てある．これは良い教科書なので，一読を勧める．

4.2 線形最小二乗法

一次関数の例を拡張して，線形最小二乗法の説明を行う．ただ，ここの線形最小二乗法の 説明は，結構いい加減なところもある．また，線形最小二乗法に適した連立方程式の解き 方も述べていない．最小二乗法について詳細に知りたければ，文献 [3]が大いに 参考になる．本格的な最小二乗法のプログラムを作成するときには，この文献を見ること を勧める．

データをフィットする線形最小二乗法の関数は，

$$y(x)=\sum_{k=1}^{M}a_kf_k(x)$$ (26)

である．ここで，$f_k(x)$ は基底関数(base function)と呼ばれるもので，$M$ はその個数 である．基底関数は，どのような関数でも良い．例えば，

$$y(x)=a_1+a_2 x+a_3x^3+a_4\sin(x)+a_4\cos^2(x^2)$$ (27)

とすることができる．

ここでは，得られたデータ点 $(x_1,\,y_1),\,(x_2,\,y_2),\,(x_2,\,y_2),\cdots,(x_N,\,y_N)$ に最も近い係数 $a_k$ を求めることが問題となる．近いというのは，一次関数の例と同じように，誤差の 二乗和が最小になることを言う．式(26)の場合，誤差の二乗は，

$$E(a_1,\,a_2,\,a_3,\cdots,a_M)=\sum_{i=1}^N \left[y_i-\sum_{k=1}^Ma_kf_k(x_i)\right]^2$$ (28)

となる．

誤差の二乗和 $E(a_1,\,a_2,\,a_3,\cdots,a_M)$ が最小になるように決めれられた係数 $a_1,\,a_2,\,a_3,\cdots,a_M$ が最良近似となる．これは，それぞれの係数で偏微分した 値がゼロになるようにすれば良い．一次関数の時と全く同じだ! 式で表すと次のようにな る．

$$\left\{ \begin{aligned}& \if 11 \frac{\partial E}{\partial a_1}... ...)\left[y_i-\sum_{k=1}^Ma_kf_k(x_i)\right]=0 \end{aligned} \right.$$

これをよく見ると，次のような連立方程式に書き直すことができる．

$$\left\{ \begin{aligned}&\sum_{i=1}^N f_1(x_i)f_1(x_i)a_1+ \sum_... ...M(x_i)f_M(x_i)a_M =\sum_{i=1}^N f_M(x_i)y_i \end{aligned} \right.$$

この連立方程式を解いて， $a_1,\,a_2,\,a_3,\cdots,a_M$ を求めれば，データをフィッ トする関数の式(26)の最適係数が求まったことになる．

連立方程式(30)の解法については，文献 [3]が詳しい．と りあえずガウス・ジョルダン法で解いてみて，ダメ(特異に近い)ならば，この文献を参考 にしてプログラムを作成せよ．特異値分解(singular value decomposition, SVD)が良い ようである．

ホームページ: Yamamoto's laboratory
著者: 山本昌志