最小二乗法

データセットを4次関数で最小二乗近似した例

最小二乗法(さいしょうにじょうほう,さいしょうじじょうほう；最小自乗法とも書く、英: least squares method)は、測定で得られた数値の組を、適当なモデルから想定される1次関数、対数曲線など特定の関数を用いて近似するときに、想定する関数が測定値に対してよい近似となるように、残差の二乗和を最小とするような係数を決定する方法、あるいはそのような方法によって近似を行うことである。

前提条件

測定データはモデル関数と誤差の和で表せる。モデル関数が正しければ誤差は偶然誤差だから、測定データとモデル関数の差の分散を最小にすることによって、モデル関数を決定することができる。物理現象の測定データには、誤差が含まれ、それは系統誤差と偶然誤差を含んでいる。この内、偶然誤差は、測定における信号経路の微視的現象に由来するならば、正規分布であると期待されることが多い。また、社会調査などの誤差理由の特定が困難な場合でも誤差が正規分布になると期待する考え方もある。

誤差が正規分布に従わない場合、最小二乗法によって得られたモデル関数は尤もらしくないことに注意する必要がある。偶然誤差が正規分布していない場合、系統誤差が無視できない位大きくそれをモデル関数に含めていない場合、測定データに正規分布から大きく外れた外れ値を含む場合などが該当する。

基礎的な考え方

話を簡単にするため、測定値は x, y の二次元の平面に分布するものとし、想定される分布が y = f(x) の形である場合を述べる。想定している関数 f は、既知の関数 g(x) の線型結合で表されていると仮定する。すなわち、

$f(x)=\sum_{k=1}^{m} a_k g_k(x)$

例えば、g_k(x)=x^k-1 は、多項式近似であり、特に m=2 の時は f(x) = a₁ + a₂x という直線による近似（線形回帰）になる。図は多項式近似で m=5 の例。

今、測定で得られた、次のような数値の組の集合があるとする。

$(x, y) = (x_1, y_1), (x_2, y_2), \ldots , (x_n, y_n)$

これら (x, y) の分布が、y = f(x) という関数に従うと仮定した時、想定される理論値は (x₁, f(x₁)), (x₂, f(x₂)), ..., (x_n, f(x_n)) ということになり、実際の測定値との残差は、各 i につき |y_i - f(x_i)| ということになる。 この残差の大きさは、xy-平面上での (x_i, y_i) と (x_i, f(x_i)) との距離でもある。

ここで、理論値からの誤差の分散の推定値は残差の平方和

$J = \sum_{i=1}^n (y_i - f(x_i))^2$

で与えられるから、J が最小になるように想定分布 f を(すなわち a_kを)、定めればよいということになる。

それには、上式は a_k を変数とする関数と見なすことができるので、J を a_k について偏微分したものをゼロと置く。こうして得られた m 個の連立方程式(正規方程式)を解き、a_k を決定すればよい。

一次方程式への近似

今、

$(x, y) = (x_1, y_1), (x_2, y_2), \ldots , (x_n, y_n)$

という測定結果が得られたとする。求めたい1次方程式の式を

$y=ax+b\,$

とおくと、aとbは次式で求められる。

$a=\frac{\displaystyle n\sum_{k=1}^n x_ky_k-\sum_{k=1}^n x_k\sum_{k=1}^n y_k}{\displaystyle n\sum_{k=1}^n x^2_k-\left( \sum_{k=1}^n x_k \right)^2}$

$b=\frac{\displaystyle \sum_{k=1}^n x^2_k\sum_{k=1}^n y_k-\sum_{k=1}^n x_ky_k\sum_{k=1}^n x_k}{\displaystyle n\sum_{k=1}^n x^2_k-\left( \sum_{k=1}^n x_k \right)^2}$

正規方程式による解法

当てはめたい関数 f は、

$f(x)= (g_1(x), g_2(x), \ldots, g_m(x))(a_1, a_2, \ldots, a_m)^T$

と、行列で表すことができる。ここに上付き添字 T は転置行列を表す。すると、最小にする関数 J は

J = (Ga − y)^T(Ga − y)

と、表される。ここに、G は、G_ij = g_j(x_i) なる成分を持つ行列、$a=(a_1, a_2, \ldots, a_m)^T$ で、$y=(y_1, y_2, \ldots, y_m)^T$。

前章で述べたように J を a のそれぞれの成分で偏微分してゼロと置いた m 個の式(正規方程式)は行列を用いて、

G^TGa = G^Ty

と、表される。行列 G^T G の逆行列が存在すれば、それを用いて a を求めることができる。数値的に解くには、LU分解や、コレスキー分解を用いることができる。

J の行列表現は、優決定の(=方程式の数が変数の数より大きい)連立一次方程式 Ga = y を近似的に解くという意味にも解釈できる。これは特異値分解(あるいは擬似逆行列)を用いて解ける。こうすれば、正規方程式の行列 G^T G が正則でない場合も解くことができる。

拡張

多次元

想定される分布が媒介変数 t を用いて (x, y) = (f(t), g(t)) の形(あるいは f, g は複数の媒介変数によって決まるとしても同様)であっても考察される。

すなわち、測定値 (x_i, y_i) がパラメータ t_i に対する (f(t_i), g(t_i)) を理論値として近似されているものと考えるのである。

この場合、各点の理論値 (f(t_i), g(t_i)) と測定値 (x_i, y_i) の間に生じる残差は

$\sqrt{(x_i - f(t_i))^2 + (y_i - g(t_i))^2}$

である。故に、残差平方和は

$\sum_{i=1}^{n}\left\{(x_i - f(t_i))^2 + (y_i -g(t_i))^2\right\}$

となるから、この値が最小であるように、f, g を決定するのである。

このように、n 組の (x, y) の測定値 (x_i, y_i) (i = 1, 2, ..., n) を n 組の (x₁, x₂, $\ldots$ x_m) の測定値 (x_1i, x_2i, ..., x_mi) (i = 1, 2, ..., n) に拡張したものも考察することができる。

測定の誤差が既知の場合

n 回の測定における誤差があらかじめ分かっている場合を考える。異なる測定方法で測定した複数のデータ列を結合する場合などでは、測定ごとに誤差が異なることはしばしばある。誤差が正規分布していると考え、その標準偏差 $\sigma_i (i=1,2,\ldots,n )$ で、誤差の大きさを表す。すると、誤差が大きい測定より、誤差が小さい測定の結果により重みをつけて近似関数を与えるべきであるから、

$J'= \sum_{i=1}^n \frac{(y_i - f(x_i))^2}{\sigma_i^2}$

を、最小にするように f を定める方がより正確な近似を与える。

毎回の測定が独立ならば、測定値の尤度は exp(-J') に比例する。そこで、上記の J' を最小にする f は、最尤推定値であるとも解釈できる。また、J' は自由度 n-m のカイ二乗分布に従うので、それを用いてモデル f の妥当性を検定することもできる。

毎回の測定誤差が同じ場合、J' を最小にするのは J を最小にするのと同じ意味になる。

非線型最小二乗法

もし、f が、a_k の線型結合で表されないときは、正規方程式を用いた解法は使えず、数値的に a_k の近似値を求める必要がある。例えば、ガウス・ニュートン法やLevenberg-Marquardt法が用いられる。とくにLevenberg-Marquardt法は多くの多次元非線形関数でパラメータを発散させずに効率よく収束させる(探索する)方法として知られている。

異常値の除去

前提条件の節で述べたように、測定データを最小二乗法によって近似する場合、異常値が含まれていると極端に近似の尤もらしさが低下することがある。また、様々な要因によって異常値を含む測定はしばしば得られるものである。

誤差が正規分布から極端に外れた異常値を取り除くための方法として修正トンプソン－τ法が用いられる。

関連項目

• 曲線あてはめ