機械学習における線形回帰モデルの最大事後確率推定
定理
データ$\mathbf{x}_{i} \in \mathbb{R}^{n}$とそのラベル$y_{i} \in \mathbb{R}$の関係が以下の線形モデルであると仮定する。
$$ y_{i} = \mathbf{w}^{\mathsf{T}} \mathbf{x}_{i} + \epsilon_{i}, \qquad i = 1, \ldots, K \tag{1} $$
事後確率が最大となるパラメータ$\mathbf{w}_{\text{MAP}}$は次の通りである。$\mathbf{y} = \begin{bmatrix} y_{1} & \cdots & y_{K} \end{bmatrix}^{\mathsf{T}}$と$\mathbf{X} = \begin{bmatrix} \mathbf{x}_{1} & \cdots & \mathbf{x}_{K} \end{bmatrix}^{T} \in \mathbb{R}^{n \times K}$について、
事前分布が正規分布の場合: $$ \mathbf{w}_{\text{MAP}} = \left(\dfrac{1}{\sigma^{2}} \mathbf{X}^{T}\mathbf{X} + \boldsymbol{\Sigma}^{-1} \right)^{-1} \left(\dfrac{1}{\sigma^{2}} \mathbf{X}^{T} \mathbf{y} + \boldsymbol{\Sigma}^{-1} \boldsymbol{\mu} \right) $$
このとき、$\boldsymbol{\mu}$と$\boldsymbol{\Sigma}$はそれぞれ$\mathbf{w}$の事前分布の平均ベクトルと共分散行列である。
事前分布がラプラス分布の場合:
最適解の明示的な形は存在しない。
説明
以下の$(3)$を見れば、$\mathbf{w}$の事前分布として標準正規分布$N(\mathbf{0}, \mathbf{I})$を仮定すると、リッジ回帰と同じ問題であることが分かる。
$$ \argmin_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} + \dfrac{1}{2} \| \mathbf{w} \|_{2}^{2} \right] $$
一方、事前分布としてラプラス分布$\operatorname{Laplace}(0, b)$を仮定すると、ラッソ回帰と同じ問題である。
$$ \argmin_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} + \dfrac{1}{b} \| \mathbf{w} \|_{1} \right] $$
事前分布=正規分布
$(1)$で$\mathbf{w} \in \mathbb{R}^{n}$は母数パラメータであり、$\epsilon_{i} \sim N(0, \sigma^{2})$をガウシアンノイズと仮定する。$\epsilon_{i}$が$N(0, \sigma^{2})$に従うと仮定したので、$y_{i} = \mathbf{w}^{\mathsf{T}} \mathbf{x}_{i} + \epsilon_{i}$は$N(\mathbf{w}^{\mathsf{T}} \mathbf{x}_{i}, \sigma^{2})$に従う。
$$ y_{i} \sim N(\mathbf{w}^{\mathsf{T}} \mathbf{x}_{i}, \sigma^{2}) $$
最大事後確率推定は次を満足する$\mathbf{w}_{\text{MAP}}$を求めるものである。
$$ \mathbf{w}_{\text{MAP}} = \argmax_{\mathbf{w}} p(\mathbf{y} | \mathbf{w}, \mathbf{X}) p(\mathbf{w}) \tag{2} $$
$p(\mathbf{y} | \mathbf{w}, \mathbf{X})$は尤度であり、$p(\mathbf{w})$は事前確率である。尤度関数は次の通りである。
$$ p(\mathbf{y} | \mathbf{w}, \mathbf{X}) = \dfrac{1}{(2\pi \sigma^{2})^{K/2}} \exp \left[ -\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} \right] $$
また、$\mathbf{w}$の事前分布が以下のような多変量正規分布に従うと仮定する。
$$ \mathbf{w} \sim N(\boldsymbol{\mu}, \boldsymbol{\Sigma}), \qquad p(\mathbf{w}) = \dfrac{1}{\sqrt{(2\pi)^{n} \det \boldsymbol{\Sigma}}} \exp \left[ -\dfrac{1}{2} (\mathbf{w} - \boldsymbol{\mu})^{\mathsf{T}} \boldsymbol{\Sigma}^{-1} (\mathbf{w} - \boldsymbol{\mu}) \right] $$
事後確率が指数関数で表現されるため、対数尤度を考慮することが計算上便利である。
$$ \begin{align*} \mathbf{w}_{\text{MAP}} &= \argmax_{\mathbf{w}} \log (p(\mathbf{y} | \mathbf{w}, \mathbf{X}) p(\mathbf{w})) \\ &= \argmax_{\mathbf{w}} \left[-\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} -\dfrac{1}{2} (\mathbf{w} - \boldsymbol{\mu})^{\mathsf{T}} \boldsymbol{\Sigma}^{-1} (\mathbf{w} - \boldsymbol{\mu}) \right] \\ &= \argmin_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} + \dfrac{1}{2} (\mathbf{w} - \boldsymbol{\mu})^{\mathsf{T}} \boldsymbol{\Sigma}^{-1} (\mathbf{w} - \boldsymbol{\mu}) \right] \tag{3}\\ \end{align*} $$
したがって、上式を微分して、$\mathbf{0}$となるような$\mathbf{w}$が$\mathbf{w}_{\text{MAP}}$である。勾配を計算してみよう。
$$ \begin{align*} & \nabla_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} +\dfrac{1}{2} (\mathbf{w} - \boldsymbol{\mu})^{\mathsf{T}} \boldsymbol{\Sigma}^{-1} (\mathbf{w} - \boldsymbol{\mu}) \right] \\ &= \nabla_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} (\mathbf{y} - \mathbf{X}\mathbf{w})^{T}(\mathbf{y} - \mathbf{X}\mathbf{w}) +\dfrac{1}{2} (\mathbf{w} - \boldsymbol{\mu})^{\mathsf{T}} \boldsymbol{\Sigma}^{-1} (\mathbf{w} - \boldsymbol{\mu}) \right] \\ &= \nabla_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} (\mathbf{y} - \mathbf{X}\mathbf{w})^{T}(\mathbf{y} - \mathbf{X}\mathbf{w})\right] + \nabla_{\mathbf{w}} \left[\dfrac{1}{2} (\mathbf{w} - \boldsymbol{\mu})^{\mathsf{T}} \boldsymbol{\Sigma}^{-1} (\mathbf{w} - \boldsymbol{\mu}) \right] \\ \end{align*} $$
微分の規則は以下を参照しよう。
内積 $$ \frac{ \partial f(\mathbf{x})}{ \partial \mathbf{x} } = \frac{ \partial (\mathbf{w}^{T}\mathbf{x})}{ \partial \mathbf{x} } = \frac{ \partial (\mathbf{x}^{T}\mathbf{w})}{ \partial \mathbf{x} } = \mathbf{w} $$
ノルム $$ \nabla f(\mathbf{x}) = \dfrac{\partial \left\| \mathbf{x} \right\|^{2}}{\partial \mathbf{x}} = \dfrac{\partial (\mathbf{x}^{T}\mathbf{x})}{\partial \mathbf{x}} = 2\mathbf{x} $$
二次形式 $$ \dfrac{\partial f(\mathbf{x})}{\partial \mathbf{x}} = \dfrac{\partial (\mathbf{x}\mathbf{R}\mathbf{x})}{\partial \mathbf{x}} = (\mathbf{R} + \mathbf{R}^{T})\mathbf{x} $$ $\mathbf{R}$が対称行列であれば、 $$ \dfrac{\partial f(\mathbf{x})}{\partial \mathbf{x}} = 2\mathbf{R}\mathbf{x} $$
最初の項から微分を計算すると次のようになる。
$$ \begin{align*} &\nabla_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} (\mathbf{y} - \mathbf{X}\mathbf{w})^{T}(\mathbf{y} - \mathbf{X}\mathbf{w})\right] \\ &=\nabla_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} \left( \mathbf{y}^{T}\mathbf{y} - 2\mathbf{y}^{T}\mathbf{X}\mathbf{w} + \mathbf{w}^{T}\mathbf{X}^{T}\mathbf{X}\mathbf{w} \right) \right] \\ &= \dfrac{1}{2\sigma^{2}}\left[- 2\mathbf{X}^{T}\mathbf{y} + 2\mathbf{X}^{T}\mathbf{X}\mathbf{w} \right] \\ &= \dfrac{1}{\sigma^{2}}\left[-\mathbf{X}^{T}\mathbf{y} + \mathbf{X}^{T}\mathbf{X}\mathbf{w} \right] \\ &= -\dfrac{1}{\sigma^{2}} \mathbf{X}^{T} \left(\mathbf{y} - \mathbf{X}\mathbf{w} \right) \\ \end{align*} $$
次に第二項を計算すると以下のようになる。共分散行列は対称行列であり、対称行列の逆行列も対称行列であることから
$$ \nabla_{\mathbf{w}} \left[\dfrac{1}{2} (\mathbf{w} - \boldsymbol{\mu})^{\mathsf{T}} \boldsymbol{\Sigma}^{-1} (\mathbf{w} - \boldsymbol{\mu}) \right] \\ = \boldsymbol{\Sigma}^{-1}(\mathbf{w} - \boldsymbol{\mu}) $$
したがって次を得る。
$$ -\dfrac{1}{\sigma^{2}} \mathbf{X}^{T} \left(\mathbf{y} - \mathbf{X}\mathbf{w}_{\text{MAP}} \right) + \boldsymbol{\Sigma}^{-1}(\mathbf{w}_{\text{MAP}} - \boldsymbol{\mu}) = \mathbf{0} $$
$\mathbf{w}_{\text{MAP}}$を求めるために上式を解くと次のようになる。
$$ \begin{align*} && - \dfrac{1}{\sigma^{2}} \mathbf{X}^{T} \mathbf{y} + \dfrac{1}{\sigma^{2}} \mathbf{X}^{T}\mathbf{X}\mathbf{w}_{\text{MAP}} + \boldsymbol{\Sigma}^{-1} \mathbf{w}_{\text{MAP}} - \boldsymbol{\Sigma}^{-1} \boldsymbol{\mu} = \mathbf{0} \\ \implies && \left(\dfrac{1}{\sigma^{2}} \mathbf{X}^{T}\mathbf{X} + \boldsymbol{\Sigma}^{-1} \right) \mathbf{w}_{\text{MAP}} = \dfrac{1}{\sigma^{2}} \mathbf{X}^{T} \mathbf{y} + \boldsymbol{\Sigma}^{-1} \boldsymbol{\mu} \\ \implies && \mathbf{w}_{\text{MAP}} = \left(\dfrac{1}{\sigma^{2}} \mathbf{X}^{T}\mathbf{X} + \boldsymbol{\Sigma}^{-1} \right)^{-1} \left(\dfrac{1}{\sigma^{2}} \mathbf{X}^{T} \mathbf{y} + \boldsymbol{\Sigma}^{-1} \boldsymbol{\mu} \right) \end{align*} $$
事前分布=ラプラス分布
いくつかの仮定や計算は、事前分布がラプラス分布であることを除けば、上と同じである。最大事後確率推定は次を満足する$\mathbf{w}_{\text{MAP}}$を求めるものである。
$$ \mathbf{w}_{\text{MAP}} = \argmax_{\mathbf{w}} p(\mathbf{y} | \mathbf{w}, \mathbf{X}) p(\mathbf{w}) $$
$p(\mathbf{y} | \mathbf{w}, \mathbf{X})$は尤度であり、$p(\mathbf{w})$は事前確率である。尤度関数は次の通りである。
$$ p(\mathbf{y} | \mathbf{w}, \mathbf{X}) = \dfrac{1}{(2\pi \sigma^{2})^{K/2}} \exp \left[ -\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} \right] $$
また、$\mathbf{w}$の事前分布に関して、それぞれの$w_{i}$が独立にラプラス分布$\operatorname{Laplace}(\mu, b)$に従うと仮定しよう。
$$ p(\mathbf{w}) = \prod_{i=1}^{n} \dfrac{1}{2b} \exp \left( -\dfrac{|w_{i} - \mu|}{b} \right) $$ $$ \implies \log p(\mathbf{w}) = -n\log 2b - \sum_{i=1}^{n} \left[ \dfrac{|w_{i} - \mu|}{b} \right] $$
事後確率が指数関数で表現されるため、対数尤度を考慮することが計算上便利である。
$$ \begin{align*} \mathbf{w}_{\text{MAP}} &= \argmax_{\mathbf{w}} \log (p(\mathbf{y} | \mathbf{w}, \mathbf{X}) p(\mathbf{w})) \\ &= \argmax_{\mathbf{w}} \left[-\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} - \sum_{i=1}^{n} \dfrac{|w_{i} - \mu|}{b} \right] \\ &= \argmin_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} + \dfrac{1}{b} \| \mathbf{w} - \mu\mathbf{1} \|_{1} \right] \\ \end{align*} $$
ここで、$\mathbf{1} \in \mathbb{R}^{n}$はすべての成分が1のベクトルである。ここで$\mu = 0$とすれば、ラッソ回帰と同じ形態になる。
$$ \argmin_{\mathbf{w}} \left[\dfrac{1}{2\sigma^{2}} \| \mathbf{y} - \mathbf{X}\mathbf{w} \|_{2}^{2} + \dfrac{1}{b} \| \mathbf{w} \|_{1} \right] $$
残念ながら、この場合、最適解の閉じた形closed formがないことが知られている。