情報基礎A 「Cプログラミング」(ステップ7・最小二乗法・理論編・補足)

このページでは、最小二乗法で得た回帰直線のパラメータの「誤差」について考える。

理想と現実

ある量 $X$ に対する変量 $Y$ の関係が、分布が$N(0,\sigma^2)$であるような確率的な変動$\Xi$を含んだ線形の関係 $$ Y = a X + b + \Xi $$ にあるとしよう。 ただしこれは「天上」での話しであって、地上の我々は、$n$ 回のサンプリングを行って、$(x_i, y_i)$ の $n$ 組の値を得ることができるだけである。 このような状況で、最小二乗法によって推定した $a,b$ には、どれくらいのあいまいさが残されているのだろうか。

上記のモデルの「真の」パラメータ $a,b$ は $$ a = \frac{E(XY) - E(X) E(Y)}{E(XX) - E(X)^2} ,\\ b = \frac{E(XX) E(Y) - E(XY) E(X)}{E(XX) - E(X)^2} $$ で与えられるが、これらの式中に現れている平均値($E(\cdot)$)が、有限のサンプルデータから得た値と、どの程度食い違うかを考えなければならない。

サンプルから得たパラメータの推定値

ここで、$y_i=a x_i + b + \xi_i$とおいて、サンプルについての平均を使って、$a$の推定値(真の値と区別するため、ここでは$\tilde{a}$と表記する)を計算すると、 $$ \begin{align} \tilde{a} & = \frac{\frac{1}{n} \sum_i x_i (a x_i + b + \xi_i) - \left( \frac{1}{n} \sum_i x_i \right) \left(\frac{1}{n} \sum (a x_i + b + \xi_i) \right)}{\frac{1}{n} \sum_i {x_i}^2 - \left(\frac{1}{n} \sum_i {x_i}\right)^2} \\ & = a + \frac{\frac{1}{n} \sum_i \left( x_i - \frac{1}{n} \sum_j x_j \right) \xi_i}{\frac{1}{n} \sum_i {x_i}^2 - \left(\frac{1}{n} \sum_i {x_i}\right)^2} \end{align} $$ となる。 $\tilde{a}$の「揺らぎ」は、$x_i$の平均値からのずれに応じた重みを付けて、揺らぎ(誤差)$\xi_i$を加算した格好になっている。

さて、正規分布$N(0,\sigma^2)$に従うような独立な確率変数 $Z_i$ を「重み」$w_i$ を付けて加えたとき、その和 $\sum_i w_i Z_i$ が従う統計性を考えよう。 それぞれの項の分散は $(w_i \sigma)^2$ になる、すなわち、 $$ w_i Z_i \sim N(0, {(w_i \sigma)}^2) $$ であるから、総和全体では $$ \sum_i w_i Z_i \sim N\left(0, \sum_i (w_i \sigma)^2\right) $$ である。 このことを使って、現実と理想の差 $\tilde{a} - a$ の分散 ${\tilde{\sigma}_a}^2$ を計算する。 上の「重み」にあたるのが $$ w_i = \frac{\frac{1}{n} \left( x_i - \frac{1}{n} \sum_j x_j \right)}{\frac{1}{n} \sum_j {x_j}^2 - \left(\frac{1}{n} \sum_j {x_j}\right)^2} $$ であることに注意すれば、直ちに、 $$ {\tilde{\sigma}_a}^2 = \frac{ \frac{1}{n} \sigma^2 }{\frac{1}{n} \sum_i {x_i}^2 - \left(\frac{1}{n} \sum_i {x_i}\right)^2} $$ が得られる。 ここで、$\Xi$の分散$\sigma^2$として、回帰直線からの残差の二乗平均 $$ {\tilde{\sigma}}^2 = \frac{1}{n} \sum_i \left( y_i - a x_i - b \right)^2 $$ を用いて得られる推定値、すなわち、不偏分散 $$ \sigma^2 = \frac{n}{n-1} {\tilde{\sigma}}^2 $$ を用いる。

ここまでをまとめると、データから推定されたフィッティングパラメータ $\tilde{a}$ は、「理想」の値 $a$ から、分散 $$ {\tilde{\sigma}_a}^2 = \frac{ 1 } { \frac{1}{n} \sum_i {x_i}^2 - \left(\frac{1}{n} \sum_i {x_i}\right)^2 } \frac{{\tilde{\sigma}}^2}{n-1} $$ だけ「ばらついて」いる可能性があり、これがパラメータ自体の曖昧さを表すひとつの指標になっている。 その分母は、$x_i$の分散なので、データが横軸方向に大きく散布していれば、曖昧さは減じられる。 一方、縦軸方向に回帰直線から大きく分散していれば、当然ながら、曖昧さは増加する。

全く同様に考えることで、y切片 $\tilde{b}$についても、その分散 ${\tilde{\sigma}_b}^2$ を見積もることができて、結果は $$ {\tilde{\sigma}_b}^2 = \frac{\frac{1}{n} \sum_i {x_i}^2 } { \frac{1}{n} \sum_i {x_i}^2 - \left(\frac{1}{n} \sum_i {x_i}\right)^2 } \frac{{\tilde{\sigma}}^2}{n-1} $$ となる。 こちらは、$x=0$付近にデータ点が沢山あるほど、精度が高まる、という形の式になっている。