線形回帰におけるパラメータの求め方の整理

# はじめに
久々に統計の勉強をしていて、最小二乗法・最尤法・線形回帰・勾配降下法あたりで頭ぐちゃぐちゃになっていたので、整理するためにこの記事を書きます。

線形回帰は、回帰のなかでも最も基本的な考え方ですが、そのパラメータの推定は、複数の方法が存在します。
ここでは、以下の2章立てでまとめました。
- パラメータの解析的導出: 最小二乗法と最尤法
- 計算機における値の計算: 行列計算と勾配降下法


# 準備
説明変数を $\boldsymbol{x} = (x_1, x_2, ..., x_m)^T$、目的変数を $y$ とした、以下のような線形回帰を考えます。
``` math
y = \boldsymbol{w}^T \boldsymbol{x} + b \\
```
定数項が出ていると邪魔なので、1項目と合わせて、以下のように書き直します (変数名は面倒なので一緒にします)
``` math
y = \boldsymbol{w}^T \boldsymbol{x} \\
```

データが、$(\boldsymbol{x_1}, y_1), (\boldsymbol{x_2}, y_2), ...$ と存在することとします。


# パラメータの解析的導出
以下の2つの観点から、パラメータを推定します。
- 最小二乗法
- 最尤法

## 最小二乗法による求め方
最小二乗法では、誤差の2乗が最小となるように、パラメータ $\hat{w}$ を求めます。すなわち、
``` math
\begin{align}
\epsilon & = \frac{1}{2}\sum_i (y_i - w^T\boldsymbol{x_i})^2 \\
& = \frac{1}{2}|\boldsymbol{y} - X\boldsymbol{w}| ^2 \\
& = \frac{1}{2}(\boldsymbol{y} - X\boldsymbol{w})^T(\boldsymbol{y} - X\boldsymbol{w})
\end{align}
```
を最小とするような、パラメータ $\hat{w}$ を求めることになります。ここで、$X, \boldsymbol{y}$は、$\boldsymbol{x_i}, y_i$ を並べた、行列、縦ベクトルです。
この $\hat{w}$ を求めるには、微分して値が0になるような点を求めればよいので、以下の通りとなります。

``` math
\begin{align}
&\frac{\partial \epsilon}{\partial \boldsymbol{w}} = 0 \\
\iff &X^T(\boldsymbol{y} - X\boldsymbol{w}) = 0 \\
\iff &X^T\boldsymbol{y} - X^TX\boldsymbol{w} = 0 \\
\iff &\boldsymbol{w} = (X^TX)^{-1}X^T\boldsymbol{y} 
\end{align}
```

## 最尤法による求め方
線形回帰の考え方の根本には、「誤差は平均0の正規分布に従う」という仮定があります。つまり、

``` math
\epsilon \sim N(0, \sigma^2)
```

です。確率密度関数 $f$ は、標準偏差 $\sigma$ をパラメータとして、以下のように表現できます。

``` math
\begin{align}
f(\epsilon;\sigma) &= \frac{1}{\sqrt{2\pi\sigma^2}}exp(-\frac{\epsilon^2}{2\sigma^2}) \\
f(y;\sigma, \boldsymbol{w}) &= \frac{1}{\sqrt{2\pi\sigma^2}}exp(-\frac{(y-\boldsymbol{w}^T\boldsymbol{x})^2}{2\sigma^2})
\end{align}
```

我々は、事後確率を最大化するようなパラメータを求める必要がありますから、最尤法の観点から
以下のような尤度関数を最大化することになります。

``` math 
L(w, \sigma) = \prod_i \frac{1}{\sqrt{2\pi\sigma^2}}exp(-\frac{(y_i-\boldsymbol{w}^T\boldsymbol{x_i})^2}{2\sigma^2})
```

今回の関心事は $w$ですから、$w$だけに注目します ($\sigma$ は定数扱いします)。
上記の式の $L$ を最大化することは、以下の式で表される「負の対数尤度」を最小化することと同値です。

``` math 
\begin{align}
-log L(w) &= \frac{1}{2\sigma^2}\sum_i (y_i-\boldsymbol{w}^T\boldsymbol{x_i})^2 + const \\
& = \frac{1}{2\sigma^2}|\boldsymbol{y} - X\boldsymbol{w}|^2 + const
\end{align}
```
定数を除けば、$\frac{1}{2}|\boldsymbol{y} - X\boldsymbol{w}|^2$ を最小化することになるので、最小二乗法と最終的に同じ結論となります。




# 計算機における値の計算
ここまでは、解析的にパラメータを求めてきました。
ここからは、計算機でどのように実際の値を求めるか考えます。

## 行列計算
ここまでの議論で、以下のようなパラメータが求まりました。

``` math
\boldsymbol{\hat{w}} = (X^TX)^{-1}X^T\boldsymbol{y} 
```
あとは、この計算をするだけで、答えが求まります。


## 勾配降下法
前項には1点欠点があります。サンプル数が大きすぎたりすると、逆行列の計算などが難しく、計算機が悲鳴を上げてしまう
場合があります。こうした場合は、数理最適化の手法を用いて値を計算することになります。
数理最適化の手法は複数ありますが、一番基礎的な勾配降下法の観点から考えたいと思います。

勾配降下法の基本から、

``` math
\begin{align}
\boldsymbol{w_{next}} &= \boldsymbol{w} - \alpha\Delta\boldsymbol{w} \\
 &= \boldsymbol{w} - \alpha\frac{df(\boldsymbol{w})}{d\boldsymbol{w}} \\
\end{align}
```

を計算することになります。ここで、$f$ は最小化したい対象となる関数ですから、ここでは最小二乗誤差をターゲットとしましょう。すなわち、

``` math
\begin{align}
\boldsymbol{w_{next}} &= \boldsymbol{w} - \alpha\frac{d}{d\boldsymbol{w}}\frac{1}{2}| X\boldsymbol{w}-\boldsymbol{y} | ^2 \\
&= \boldsymbol{w} - \alpha X^T(X\boldsymbol{w} - \boldsymbol{y}) \\
\end{align}
```
を計算することになります。
(微分した際の符号の問題があり、これまでとは、$y$ と $X\boldsymbol{w}$ の順序を入れ替えています)

なお、ここではバッチ学習を前提として式を書きましたが、オンライン学習の場合や、確率的勾配降下法の場合も、ほぼ同様かと思います。

# 参考文献
- https://www.math.uwaterloo.ca/~hwolkowi/matrixcookbook.pdf
 - https://qiita.com/sugiyama34/items/9ad73435f547817283a3#%E6%9C%80%E5%B0%8F2%E4%B9%97%E3%81%AE%E5%A0%B4%E5%90%88
- https://www.data.math.ryukoku.ac.jp/course/probstat1_2021/w12.pdf