最尤推定法に基づく正規方程式の導出

被説明変数 $\boldsymbol{y}$ を、説明変数 $\boldsymbol{X}$ と回帰係数 $\boldsymbol{\beta}$ の線形結合と誤差項 $\boldsymbol{\varepsilon}$ によって表現する線形回帰モデル

\newcommand{\b}[1]{\boldsymbol{#1}} \b{y} = \b{X} \b{\beta} + \b{\varepsilon}

(1)

を考える。また、以下を仮定する（古典的正規回帰モデル（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)

(2)

ただし、 $\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

(3)

対数尤度は

Undefined control sequence: \b at position 21: …n{align}
\ln L(\̲b̲{y}| \b{X}, \b{…

\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}

導出メモ

まず尤度を軽く整理する

\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}

(5)

対数をとると

\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}

(6)

（参考）対数の公式

\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}

(7)

回帰係数の推定¶

対数尤度を\b{\beta}について微分した勾配を0と置いて解くと

Undefined control sequence: \b at position 28: …}
\nabla \ln L(\̲b̲{y}| \b{X}, \b{…

\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}

となり、最尤推定量と最小二乗推定量が同じ方程式になることがわかる

導出メモ

対数尤度を\b{\beta}について微分した勾配を整理すると

Undefined control sequence: \b at position 28: …}
\nabla \ln L(\̲b̲{y}| \b{X}, \b{…

\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}

（\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とおいて

Undefined control sequence: \b at position 25: …{1}{\sigma^2} (\̲b̲{X}^\top \b{y})…

0 = \frac{1}{\sigma^2} (\b{X}^\top \b{y})^\top
  - \b{\beta}^\top \frac{1}{\sigma^2} \b{X}^\top \b{X}

として整理すると

Undefined control sequence: \b at position 1: \̲b̲{\beta}^\top \f…

\b{\beta}^\top \frac{1}{\sigma^2} \b{X}^\top \b{X}
  = \frac{1}{\sigma^2} (\b{X}^\top \b{y})^\top

Undefined control sequence: \b at position 20: …in{align}
&\to \̲b̲{\beta}^\top \f…

\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}

（(\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

(13)

ベクトル $\boldsymbol{b} = (b_1, b_2, \dots, b_n)^\top$ とベクトル $\boldsymbol{a}_i = (a_{i1}, a_{i2}, \cdots, a_{im})$ があるとき、スカラーの $b_i$ との積の和は

\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}

(14)

ベクトル\b{a}_iを集めた行列 $A$ を考えたとき、

Undefined control sequence: \b at position 1: \̲b̲{A}^\top \b{b}
…

\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}

ゆえに

\sum_{i=1}^n b_i \boldsymbol{a}_i^\top = (A^\top b)^\top

(16)

導出メモ

スカラーの場合、合成関数の微分より

\frac{d (y - \hat{y})^2 }{d \hat{y}} = 2 (y - \hat{y}) \times (-1)

(17)

の形になるがベクトルなので

Undefined control sequence: \b at position 23: …partial (y_i - \̲b̲{x}_i^\top \b{\…

\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$ の導出¶

Undefined control sequence: \b at position 36: …\partial \ln L(\̲b̲{y}| \b{X}, \b{…

\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}

これをゼロにする $\sigma^2$ は

Undefined control sequence: \b at position 76: …N_{i=1} (y_i - \̲b̲{x}_i^\top \b{\…

\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}

となる

\b{\beta}の最尤推定量 \hat{\b{\beta}}_{ML} = (\b{X}^\top \b{X})^{-1} \b{X}^\top \b{y}を代入すると

Undefined control sequence: \b at position 70: …N_{i=1} (y_i - \̲b̲{x}_i^\top \hat…

\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}

回帰係数の推定¶

標準偏差σ\sigmaσの導出¶

標準偏差 $\sigma$ の導出¶