最尤推定法に基づく正規方程式の導出
被説明変数\(\boldsymbol{y}\)を、説明変数\(\boldsymbol{X}\)と回帰係数\(\boldsymbol{\beta}\)の線形結合と誤差項\(\boldsymbol{\varepsilon}\)によって表現する線形回帰モデル
\[
\newcommand{\b}[1]{\boldsymbol{#1}}
\b{y} = \b{X} \b{\beta} + \b{\varepsilon}
\]
を考える。また、以下を仮定する（古典的正規回帰モデル（Classical Normal Regression Model: CNRM） の仮定）
説明変数\(\b{X}\)は非確率的である
\(E[Y] = \b{X \beta}\)であり、したがって誤差項の期待値はゼロ：\(E(\b{\varepsilon}) = 0\)
誤差項の分散\(\sigma^2\)は一定（均一分散）であり、共分散はゼロ（独立性）：\(Var(\b{\varepsilon}) = \sigma^2 \b{I}\)
\(\b{X}\)の階数は\(k\)：\(rank(\b{X}) = k\)（\(\b{X^\top X}\)に逆行列が存在することの仮定）
\(Y\)は正規分布に従う（正規性）：\(Y \sim \mathcal{N}(\b{X\beta}, \sigma^2\b{I}), \quad \b{\varepsilon} \sim \mathcal{N}(0, \sigma^2\b{I})\)
特に3（標本がi.i.d.）と5（正規性）は最尤推定のために必要で、誤差項\(\varepsilon_i = y_i - \boldsymbol{x}_i^\top \boldsymbol{\beta}\)が従う分布型を仮定することで最尤推定が可能になる。
尤度としては、平均\(\boldsymbol{x}_i^\top \boldsymbol{\beta}\)の正規分布で観測値\(y_i\)が得られる確率\(\mathcal{N}(y_i| \boldsymbol{x}_i^\top \boldsymbol{\beta}, \sigma^2I)\)を用いて
\[
\newcommand{\b}[1]{\boldsymbol{#1}}
L(\b{y}| \b{X}, \b{\beta}, \sigma) = \prod^N_{i=1} \mathcal{N}(y_i| \b{x}_i^\top \b{\beta}, \sigma^2I)
\]
ただし、\(\mathcal{N}\)は正規分布
\[
\mathcal{N}(x|\mu, \sigma^2)
= \frac{1}{\sqrt{2\pi}\sigma}
\exp \left\{ - \frac{(x-\mu)^2}{2\sigma^2} \right\}
, \hspace{1em} -\infty < x < \infty
\]
対数尤度は
\[\begin{split}
\begin{align}
\ln L(\b{y}| \b{X}, \b{\beta}, \sigma)
&= \sum^N_{i=1} \ln \mathcal{N}(y_i| \b{x}_i^\top \b{\beta}, \sigma^2)
\\
&= - \frac{N}{2} \ln (2\pi)
- \frac{N}{2} \ln \sigma^2
- \frac{1}{2\sigma^2} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})^2
\end{align}
\end{split}\]
導出メモ
まず尤度を軽く整理する
\[\begin{split}
\begin{align}
L(\mu, \sigma^2)
&= \prod_{i=1}^n \frac{1}{\sqrt{2\pi\sigma^2}}
\exp \left \{ - \frac{(x_i - \mu)^2}{2\sigma^2} \right \}
\\
&= \left( \frac{1}{\sqrt{2\pi\sigma^2}} \right)^n
\prod_{i=1}^n
\exp \left \{ - \frac{(x_i - \mu)^2}{2\sigma^2} \right \}
\end{align}
\end{split}\]
対数をとると
\[\begin{split}
\begin{align}
\ln L(\mu, \sigma^2)
&= n \ln \left( \frac{1}{\sqrt{2\pi\sigma^2}} \right)
+ \sum^n_{i=1} \ln \left(
\exp \left\{ - \frac{(x_i - \mu)^2}{2\sigma^2} \right\}
\right)
\\
&= n \ln ( 2\pi\sigma^2 )^{-1/2}
+ \sum^n_{i=1} \left(
- \frac{(x_i - \mu)^2}{2\sigma^2}
\ln e
\right)
\\
&= - \frac{n}{2} \ln ( 2\pi\sigma^2 )
- \sum^n_{i=1} \frac{(x_i - \mu)^2}{2\sigma^2}
\\
&= - \frac{n}{2} \ln (2\pi)
- \frac{n}{2} \ln \sigma^2
- \frac{1}{2\sigma^2} \sum^n_{i=1} (x_i - \mu)^2
\end{align}
\end{split}\]
（参考）対数の公式
\[\begin{split}
\begin{align}
\log 1 &= 0\\
\log_e e &= 1\\
\log x^a &= a\log x\\
\log \sqrt{a} &= \log a^{1/2} = \textstyle \frac{1}{2} \log a\\
\log ab &= \log a + \log b\\
\log \frac{a}{b} &= \log a - \log b
\end{align}
\end{split}\]
回帰係数の推定
対数尤度を\(\b{\beta}\)について微分した勾配を0と置いて解くと
\[\begin{split}
\begin{align}
\nabla \ln L(\b{y}| \b{X}, \b{\beta}, \sigma)
&= \frac{1}{\sigma^2} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})\b{x}_i^\top\\
&= \frac{1}{\sigma^2} (\sum^N_{i=1} y_i \b{x}_i^\top)
- \b{\beta}^\top \frac{1}{\sigma^2} (\sum^N_{i=1} \b{x}_i \b{x}_i^\top)\\
&= \frac{1}{\sigma^2} (\b{X}^\top \b{y})^\top
- \b{\beta}^\top \frac{1}{\sigma^2} \b{X}^\top \b{X}\\
&= 0\\
\to \hat{\b{\beta}} &= (\b{X}^\top \b{X})^{-1} \b{X}^\top \b{y}\\
\end{align}
\end{split}\]
となり、最尤推定量と最小二乗推定量が同じ方程式になることがわかる
導出メモ
対数尤度を\(\b{\beta}\)について微分した勾配を整理すると
\[\begin{split}
\begin{align}
\nabla \ln L(\b{y}| \b{X}, \b{\beta}, \sigma)
&= \frac{1}{\sigma^2} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})\b{x}_i^\top\\
&= \frac{1}{\sigma^2} \sum^N_{i=1} y_i \b{x}_i^\top
- \frac{1}{\sigma^2} \sum^N_{i=1} \b{x}_i^\top \b{\beta} \b{x}_i^\top\\
&= \frac{1}{\sigma^2} \sum^N_{i=1} y_i \b{x}_i^\top
- \frac{1}{\sigma^2} \sum^N_{i=1} \b{\beta}^\top \b{x}_i \b{x}_i^\top\\
&= \frac{1}{\sigma^2} (\sum^N_{i=1} y_i \b{x}_i^\top)
- \b{\beta}^\top \frac{1}{\sigma^2} (\sum^N_{i=1} \b{x}_i \b{x}_i^\top)\\
&= \frac{1}{\sigma^2} (\b{X}^\top \b{y})^\top
- \b{\beta}^\top \frac{1}{\sigma^2} \b{X}^\top \b{X}
\end{align}
\end{split}\]
（\(\sum^N_{i=1} y_i \b{x}_i^\top = (\b{X}^\top \b{y})^\top\)は行ベクトルで\(\b{\beta}^\top \b{X}^\top \b{X}\)も行ベクトルなので演算できる）
これを0とおいて
\[
0 = \frac{1}{\sigma^2} (\b{X}^\top \b{y})^\top
- \b{\beta}^\top \frac{1}{\sigma^2} \b{X}^\top \b{X}
\]
として整理すると
\[
\b{\beta}^\top \frac{1}{\sigma^2} \b{X}^\top \b{X}
= \frac{1}{\sigma^2} (\b{X}^\top \b{y})^\top
\]
\[\begin{split}
\begin{align}
&\to \b{\beta}^\top \frac{1}{\sigma^2} \b{X}^\top \b{X} = \frac{1}{\sigma^2} (\b{X}^\top \b{y})^\top\\
&\to \b{\beta}^\top (\b{X}^\top \b{X}) = (\b{X}^\top \b{y})^\top\\
&\to \b{\beta}^\top = (\b{X}^\top \b{X})^{-1} (\b{X}^\top \b{y})^\top\\
&\to \b{\beta} = (\b{X}^\top \b{X})^{-1} \b{X}^\top \b{y}\\
\end{align}
\end{split}\]
（\((\b{X}^\top \b{X})\)は対角行列なので\((\b{X}^\top \b{X})^\top = (\b{X}^\top \b{X})\)）
※内積の結果はスカラーになるので、\(\b{x}^\top \b{\beta} = \b{\beta}^\top \b{x}\)である
2つの行列\(A, B\)があるとき、その行ベクトル\(a_i, b_i\)の外積\(\boldsymbol{a}_i \boldsymbol{b}_i^\top\)と行列の積は一致する
\[
A^\top B = \sum^n_{i=1} \boldsymbol{a}_i \boldsymbol{b}_i^\top
\]
ベクトル\(\boldsymbol{b} = (b_1, b_2, \dots, b_n)^\top\)とベクトル\(\boldsymbol{a}_i = (a_{i1}, a_{i2}, \cdots, a_{im})\)があるとき、スカラーの\(b_i\)との積の和は
\[\begin{split}
\begin{align}
&\sum_{i=1}^n b_i \boldsymbol{a}_i^\top\\
&= \sum_{i=1}^n (b_i a_{i1}, \cdots, b_i a_{im})\\
&= \begin{pmatrix}
b_1 a_{11} & \cdots & b_1 a_{1m}\\
+ & & +\\
\vdots & \cdots & \vdots\\
+ & & +\\
b_n a_{n1} & \cdots & b_n a_{nm}\\
\end{pmatrix}\\
\end{align}
\end{split}\]
ベクトル\(\b{a}_i\)を集めた行列\(A\)を考えたとき、
\[\begin{split}
\b{A}^\top \b{b}
= \begin{pmatrix}
b_1 a_{11} + \cdots + b_n a_{n1}\\
\vdots\\
b_1 a_{m1} + \cdots + b_n a_{nm}\\
\end{pmatrix}
\end{split}\]
ゆえに
\[
\sum_{i=1}^n b_i \boldsymbol{a}_i^\top
= (A^\top b)^\top
\]
導出メモ
スカラーの場合、合成関数の微分より
\[
\frac{d (y - \hat{y})^2 }{d \hat{y}} = 2 (y - \hat{y}) \times (-1)
\]
の形になるがベクトルなので
\[
\frac{\partial (y_i - \b{x}_i^\top \b{\beta})^2 }{\partial \beta_j}
= - 2 (y_i - \b{x}_i^\top \b{\beta}) \b{x}_i^\top
\]
これを各\(j\)について計算してベクトル（勾配）とし、それを\(i\)ごとに和する
標準偏差\(\sigma\)の導出
\[\begin{split}
\begin{align}
\frac{\partial \ln L(\b{y}| \b{X}, \b{\beta}, \sigma) }{ \partial \sigma }
&= \frac{\partial}{ \partial \sigma }
\left(- \frac{N}{2} \ln (2\pi)
- \frac{N}{2} \ln \sigma^2
- \frac{1}{2\sigma^2} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})^2 \right)\\
&= - \frac{N}{2} \frac{1}{\sigma^2} 2\sigma
+ \frac{4\sigma}{4\sigma^4} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})^2\\
&= - \frac{N}{\sigma}
+ \frac{1}{\sigma^3} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})^2\\
\end{align}
\end{split}\]
これをゼロにする\(\sigma^2\)は
\[\begin{split}
\begin{align}
&- \frac{N}{\sigma} + \frac{1}{\sigma^3} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})^2 = 0\\
&\to \frac{N}{\sigma} = \frac{1}{\sigma^3} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})^2\\
&\to \frac{\sigma^3}{\sigma} = \frac{\sigma^3}{N \sigma^3} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})^2\\
&\to \sigma^2 = \frac{1}{N} \sum^N_{i=1} (y_i - \b{x}_i^\top \b{\beta})^2\\
\end{align}
\end{split}\]
となる
\(\b{\beta}\)の最尤推定量 \(\hat{\b{\beta}}_{ML} = (\b{X}^\top \b{X})^{-1} \b{X}^\top \b{y}\)を代入すると
\[\begin{split}
\begin{align}
\hat{\sigma}^2_{ML}
&= \frac{1}{N} \sum^N_{i=1} (y_i - \b{x}_i^\top \hat{\b{\beta}}_{ML})^2\\
&= \frac{ \b{\varepsilon}^\top \b{\varepsilon} }{N}\\
\end{align}
\end{split}\]