機械学習基礎理論独習

誤りがあればご指摘いただけると幸いです。数式が整うまで少し時間かかります。リンクフリーです。

勉強ログです。リンクフリーです
目次へ戻る

エビデンス近似

MAP推定とリッジ回帰

MAP推定とリッジ回帰の関係という記事で
事前分布を平均は{\bf 0}、共分散行列は等方的な分布\alpha^{-1}{\bf I}_Mとしたときに
\bf wについてMAP推定することと、
二乗和誤差関数と二次正則化項の和を最小化することが等価であることを示しました。

これを利用して、リッジ回帰(L2正則化)のパラメータ\lambdaを推定してみます。

MAP推定の記事より、事後分布は

\begin{eqnarray}
&&p({\bf w}|{\bf X},{\bf t})=\mathcal{N}({\bf w}|{\bf m}_N,{\bf S}_N)\tag{1}\\
&&{\bf m}_N={\bf S}_N({\bf S}_0^{-1}{\bf m}_0+\beta{\boldsymbol\Phi}^\top{\bf t})\tag{2}\\
&&{\bf S}_N^{-1}={\bf S}_0^{-1}+\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}\tag{3}\\
\end{eqnarray}

でした。{\bf m}_0={\bf 0},{\bf S}_0^{-1}=\alpha{\bf I}_Mとおくと、

\begin{eqnarray}
&&{\bf m}_N=\beta{\bf S}_N{\boldsymbol\Phi}^\top{\bf t}\tag{4}\\
&&{\bf S}_N^{-1}=\alpha{\bf I}_M+\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}\tag{5}\\
\end{eqnarray}

となります。

事前分布は

\begin{eqnarray}
&&p({\bf w}|\alpha)=\mathcal{N}({\bf w}|{\bf 0},\alpha^{-1}{\bf I}_M)\tag{6}\\
\end{eqnarray}

です。

尤度は

\begin{eqnarray}
p({\bf t}|{\bf w},\beta)&=&\prod_{n=1}^N\mathcal{N}(t_n|{\bf w}^\top{\boldsymbol\phi}({\bf x}_n),\beta^{-1})\\
&=&\prod_{n=1}^N\left(\frac{\beta}{2\pi}\right)^{\frac{1}{2}}\exp\left(-\frac{\beta}{2}\left(t_n-{\boldsymbol\phi}({\bf x}_n)^\top{\bf w}\right)^2\right)\\
&=&\frac{1}{(2\pi\beta^{-1})^{\frac{N}{2}}}\exp\left(-\frac{\beta}{2}\sum_{n=1}^N \left(t_n-{\boldsymbol\phi}({\bf x}_n)^\top{\bf w}\right)^2\right)\\
&=&\frac{1}{(2\pi)^{\frac{N}{2}}|\beta^{-1}{\bf I}_N|^{\frac{1}{2}}}\exp\left(-\frac{\beta}{2}({\bf t}-{\boldsymbol\Phi}{\bf w})^\top({\bf t}-{\boldsymbol\Phi}{\bf w})\right)\\
&=&\mathcal{N}({\bf t}|{\boldsymbol\Phi}{\bf w},\beta^{-1}{\bf I}_N)\tag{7}\\
\end{eqnarray}

です。

周辺尤度

周辺尤度を計算します。

\begin{eqnarray}
p({\bf t}|\alpha,\beta)&=&\int p({\bf t},{\bf w}|\alpha,\beta){\rm d}{\bf w}\\
&=&\int p({\bf t}|{\bf w},\beta)p({\bf w}|\alpha){\rm d}{\bf w}\\
&=&\int \mathcal{N}({\bf t}|{\boldsymbol\Phi}{\bf w},\beta^{-1}{\bf I}_N)\mathcal{N}({\bf w}|{\bf 0},\alpha^{-1}{\bf I}_M){\rm d}{\bf w}\tag{8}\\
\end{eqnarray}

(8)は以下の公式を使えば、積分計算せずに求まります。

\begin{eqnarray}
&&p(\bf x) &=& \mathcal{N}(\mathbf x | \boldsymbol\mu, \mathbf\Lambda^{-1})\tag{9}\\
&&p(\bf y | \bf x) &=& \mathcal{N}(\mathbf y | \mathbf A \mathbf x + \mathbf b, \mathbf{L}^{-1}) \tag{10}\\
\end{eqnarray}

のとき、

\begin{eqnarray}
p(\mathbf y) = \mathcal{N}(\mathbf y | \mathbf A {\boldsymbol\mu} + \mathbf b , \mathbf{L}^{-1} + \mathbf A \mathbf \Lambda^{-1} \mathbf A^{\top}) \tag{11}
\end{eqnarray}

が成り立ちます。

この公式に{\bf x}={\bf w},{\boldsymbol\mu}={\bf 0},{\boldsymbol\Lambda}^{-1}=\alpha^{-1}{\bf I}_M,{\bf y}={\bf t},{\bf A}={\boldsymbol\Phi},{\bf L}^{-1}=\beta^{-1}{\bf I}_N,{\bf b}={\bf 0}のように当てはめると

\begin{eqnarray}
&&p({\bf t}|\alpha,\beta)=\mathcal{N}({\bf t}|{\bf 0},\beta^{-1}{\bf I}_N+\alpha^{-1}{\boldsymbol\Phi}{\boldsymbol\Phi}^\top)\tag{12}\\
\end{eqnarray}

となります。
対数を取ります。

\begin{eqnarray}
\ln p({\bf t}|\alpha,\beta)&=&\ln \mathcal{N}({\bf t}|{\bf 0},\beta^{-1}{\bf I}_N+\alpha^{-1}{\boldsymbol\Phi}{\boldsymbol\Phi}^\top)\\
&=&-\frac{N}{2}\ln 2\pi-\frac{1}{2}\ln|\beta^{-1}{\bf I}_N+\alpha^{-1}{\boldsymbol\Phi}{\boldsymbol\Phi}^\top|-\frac{1}{2}{\bf t}^\top(\beta^{-1}{\bf I}_N+\alpha^{-1}{\boldsymbol\Phi}{\boldsymbol\Phi}^\top)^{-1}{\bf t}\tag{13}\\
\end{eqnarray}

ここで公式|a{\bf A}|=a^N|{\bf A}|,|{\bf I}_N+{\bf AB}^\top|=|{\bf I}_M+{\bf A}^\top{\bf B}|より、(13)の第2項は

\begin{eqnarray}
 |\beta^{-1}{\bf I}_N+\alpha^{-1}{\boldsymbol\Phi}{\boldsymbol\Phi}^\top|&=&\beta^{-N}|{\bf I}_N+\beta\alpha^{-1}{\boldsymbol\Phi}{\boldsymbol\Phi}^\top|\\
&=&\beta^{-N}|{\bf I}_M+\beta\alpha^{-1}{\boldsymbol\Phi}^\top{\boldsymbol\Phi}|\\
&=&\beta^{-N}\alpha^{-M}|\alpha{\bf I}_M+\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}|\\
&=&\beta^{-N}\alpha^{-M}|{\bf A}|\tag{14}
\end{eqnarray}

(14){\bf A}=\alpha{\bf I}_M+\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}とおきました。

ウッドベリーの公式

\begin{eqnarray}
({\bf A}+{\bf B}{\bf D}^{-1}{\bf C})^{-1}={\bf A}^{-1}-{\bf A}^{-1}{\bf B}({\bf D}+{\bf C}{\bf A}^{-1}{\bf B})^{-1}{\bf C}{\bf A}^{-1}\tag{15}
\end{eqnarray}

{\bf A}=\beta^{-1}{\bf I}_N,{\bf B}=\alpha^{-1}{\boldsymbol\Phi},{\bf C}={\boldsymbol\Phi}^\top,{\bf D}^{-1}={\bf I}_Mとして、(13)の第3項に適用します。

\begin{eqnarray}
 -\frac{1}{2}{\bf t}^\top(\beta^{-1}{\bf I}_N+\alpha^{-1}{\boldsymbol\Phi}{\boldsymbol\Phi}^\top)^{-1}{\bf t}&=&-\frac{1}{2}{\bf t}^\top(\beta{\bf I}_N-\beta\alpha^{-1}{\boldsymbol\Phi}({\bf I}_M+{\boldsymbol\Phi}^\top\beta\alpha^{-1}{\boldsymbol\Phi})^{-1}{\boldsymbol\Phi}^\top\beta){\bf t}\\
&=&-\frac{1}{2}{\bf t}^\top(\beta{\bf I}_N-\beta^2\alpha^{-1}{\boldsymbol\Phi}(\alpha^{-1}(\alpha{\bf I}_M+\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}))^{-1}{\boldsymbol\Phi}^\top){\bf t}\\
&=&-\frac{1}{2}{\bf t}^\top(\beta{\bf I}_N-\beta^2\alpha^{-1}{\boldsymbol\Phi}(\alpha(\alpha{\bf I}_M+\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi})^{-1}){\boldsymbol\Phi}^\top){\bf t}\\
&=&-\frac{1}{2}{\bf t}^\top(\beta{\bf I}_N-\beta^2{\boldsymbol\Phi}(\alpha{\bf I}_M+\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi})^{-1}{\boldsymbol\Phi}^\top){\bf t}\\
&=&-\frac{1}{2}{\bf t}^\top(\beta{\bf I}_N-\beta^2{\boldsymbol\Phi}{\bf A}^{-1}{\boldsymbol\Phi}^\top){\bf t}\\
&=&-\frac{1}{2}(\beta{\bf t}^\top{\bf t}-{\beta^2}{\bf t}^\top{\boldsymbol\Phi}{\bf A}^{-1}{\boldsymbol\Phi}^\top{\bf t})\\
&=&-\frac{1}{2}(\beta{\bf t}^\top{\bf t}-{\beta^2}{\bf t}^\top{\boldsymbol\Phi}{\bf A}^{-1}{\bf A}{\bf A}^{-1}{\boldsymbol\Phi}^\top{\bf t})\\
&=&-\frac{1}{2}(\beta{\bf t}^\top{\bf t}-({\beta}{\bf t}^\top{\boldsymbol\Phi}{\bf A}^{-1}){\bf A}(\beta{\bf A}^{-1}{\boldsymbol\Phi}^\top{\bf t}))\\
&=&-\frac{1}{2}(\beta{\bf t}^\top{\bf t}-(\beta{\bf A}^{-1}{\boldsymbol\Phi}^\top{\bf t})^\top{\bf A}(\beta{\bf A}^{-1}{\boldsymbol\Phi}^\top{\bf t}))\\
&=&-\frac{1}{2}(\beta{\bf t}^\top{\bf t}-{\bf m}_N^\top{\bf A}{\bf m}_N)\\
&=&-\frac{1}{2} \left ( \beta \mathbf t^\top \mathbf t-2 \mathbf m^\top_N \mathbf A \mathbf m_N +  \mathbf m^\top_N \mathbf A \mathbf m_N \right ) \\
&=& - \frac { 1 } { 2 } \left ( \beta \mathbf { t } ^\top \mathbf { t } - 2  \mathbf { m } ^ \top_ { N } \mathbf { A } \left ( \beta \mathbf { A } ^ { - 1 } \mathbf { \Phi } ^\top \mathbf { t } \right ) +  \mathbf { m } ^\top_ { N }  \left ( \alpha \mathbf { I } _ { M } + \beta \mathbf { \Phi } ^\top \mathbf { \Phi } \right ) \mathbf { m } _ { N } \right ) \\
&=&-\frac{1}{2} \left ( \beta \mathbf t^\top \mathbf t-2 \mathbf m^\top_N  \boldsymbol \Phi  ^\top \mathbf t \beta + \beta  \mathbf m^\top_N  \boldsymbol \Phi  ^\top\boldsymbol \Phi\mathbf m_N+\alpha \mathbf m^\top_N \mathbf m_N \right ) \\
&=&-\frac{1}{2} \left ( \beta  \left ( \mathbf t-\boldsymbol \Phi \mathbf m_N \right ) ^\top \left ( \mathbf t-\boldsymbol \Phi \mathbf m_N \right )  +\alpha \mathbf m_N^\top \mathbf m_N \right ) \\
&=& -\frac{\beta}{2}  \left | \left | \mathbf t-\boldsymbol \Phi \mathbf m_N \right | \right | ^ 2-\frac{\alpha}{2}\mathbf m_N^\top \mathbf m_N\tag{16}
\end{eqnarray}

(14),(16)(13)に代入します。

\begin{eqnarray}
\ln p({\bf t}|\alpha,\beta)&=&-\frac{N}{2}\ln 2\pi-\frac{1}{2}\ln\beta^{-N}\alpha^{-M}|{\bf A}|-\frac{\beta}{2}  \left | \left | \mathbf t-\boldsymbol \Phi \mathbf m_N \right | \right | ^ 2-\frac{\alpha}{2}\mathbf m_N^\top \mathbf m_N\\
&=&\frac{M}{2}\ln\alpha+\frac{N}{2}\ln\beta-E({\bf m}_N)-\frac{1}{2}\ln|{\bf A}|-\frac{N}{2}\ln(2\pi)\tag{18}\\
\end{eqnarray}

(18)E({\bf m}_N)は以下のようにおきました。

\begin{eqnarray}
E({\bf m}_N)=\frac{\beta}{2}  \left | \left | \mathbf t-\boldsymbol \Phi \mathbf m_N \right | \right | ^ 2+\frac{\alpha}{2}\mathbf m_N^\top \mathbf m_N\tag{19}
\end{eqnarray}

周辺尤度の最大化

対数周辺尤度を最大化する\alpha,\betaを求めるのですが、解析的に求めるのが難しい為、
\alpha,\betaの更新式を導出し、繰り返し法で解を求めます。

行列\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}を行列\bf Pで対角化します。

\begin{eqnarray}
{\bf P}^{-1}(\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}){\bf P}={\rm diag}(\lambda_1,\ldots,\lambda_M)\tag{20}
\end{eqnarray}

すると、行列{\bf A}にも対角化ができます。

\begin{eqnarray}
{\bf P}^{-1}{\bf A}{\bf P}={\rm diag}(\alpha+\lambda_1,\ldots,\alpha+\lambda_M)\tag{21}
\end{eqnarray}

(21)の導出については、後で説明します。
よって、

\begin{eqnarray}
 |{\bf A}|=\prod_{i=1}^M(\lambda_i+\alpha)\tag{22}
\end{eqnarray}

となります。

\ln|{\bf A}|\alpha微分します。

\begin{eqnarray}
\frac{d}{d\alpha}\ln|{\bf A}|&=&\frac{d}{d\alpha}\ln\prod_{i=1}^M(\lambda_i+\alpha)\\
&=&\frac{d}{d\alpha}\sum_{i=1}^M\ln(\lambda_i+\alpha)\\
&=&\sum_{i=1}^M\frac{1}{\lambda_i+\alpha}\tag{23}\\
\end{eqnarray}

対数周辺尤度(18)\alpha微分して=0とおきます。

\begin{eqnarray}
&&\frac{\partial}{\partial\alpha}\ln p({\bf t}|\alpha,\beta)=0\\
&&\Leftrightarrow \frac{M}{2\alpha}-\frac{1}{2}{\bf m}_N^\top{\bf m}_N-\frac{1}{2}\sum_{i=1}^M\frac{1}{\lambda_i+\alpha}=0\\
&&\Leftrightarrow\alpha{\bf m}_N^\top{\bf m}_N=M-\sum_{i=1}^M\frac{\alpha}{\lambda_i+\alpha}\\
&&\Leftrightarrow\alpha{\bf m}_N^\top{\bf m}_N=\sum_{i=1}^M\frac{\lambda_i}{\lambda_i+\alpha}\\
&&\Leftrightarrow\alpha=\frac{\gamma}{{\bf m}_N^\top{\bf m}_N}\tag{24}\\
\end{eqnarray}

(24)\gammaを以下のようにおきました。

\begin{eqnarray}
\gamma=\sum_{i=1}^M\frac{\lambda_i}{\lambda_i+\alpha}\tag{25}
\end{eqnarray}

(24)\alphaの更新式です。
(24)の右辺の\gamma\alphaを含むので、(24)\alphaに関する陰関数になっていることに注意してください。

固有値\lambda_i\betaに比例するので、\lambda_i=\beta a_iとおいて、\beta微分すると、

\begin{eqnarray}
\frac{d}{d\beta}\lambda_i&=&\frac{d}{d\beta}\beta a_i\\
&=&a_i\\
&=&\frac{\lambda_i}{\beta}\tag{26}
\end{eqnarray}

となります。
\ln|{\bf A}|\beta微分します。

\begin{eqnarray}
\frac{d}{d\beta}\ln|{\bf A}|&=&\frac{d}{d\beta}\ln\prod_{i=1}^M(\lambda_i+\alpha)\\
&=&\frac{d}{d\beta}\sum_{i=1}^M\ln(\lambda_i+\alpha)\\
&=&\sum_{i=1}^M\frac{d}{d(\lambda_i+\alpha)}\ln(\lambda_i+\alpha)\frac{d}{d\beta}(\lambda_i+\alpha)\\
&=&\frac{1}{\beta}\sum_{i=1}^M\frac{\lambda_i}{\lambda_i+\alpha}\\
&=&\frac{\gamma}{\beta}\tag{27}
\end{eqnarray}

対数周辺尤度(18)\beta微分して=0とおきます。

\begin{eqnarray}
&&\frac{\partial}{\partial\beta}\ln p({\bf t}|\alpha,\beta)=0\\
&&\Leftrightarrow\frac{N}{2\beta}-\frac{1}{2}||{\bf t}-{\boldsymbol\Phi}{\bf m}_N||^2-\frac{\gamma}{2\beta}=0\\
&&\Leftrightarrow\frac{1}{\beta}=\frac{1}{N-\gamma}||{\bf t}-{\boldsymbol\Phi}{\bf m}_N||^2\tag{28}\\
\end{eqnarray}

(27)\betaの更新式です。
(27)の右辺の\gamma\lambda_iを含んでおり、\lambda_i\betaに依存するので、
(27)\betaに関する陰関数になっています。

\alpha,\betaに関する更新式が求まりましたので、収束するまで繰り返せば、近似解\hat{\alpha},\hat{\beta}が求まります。
MAP推定とリッジ回帰の関係の記事より、リッジ回帰の\lambda\alpha / \betaなので、

\begin{eqnarray}
\hat{\lambda}=\frac{\hat{\alpha}}{\hat{\beta}}\tag{29}
\end{eqnarray}

となります。

(21)の導出

\alpha{\bf I}固有値は全て\alpha(重解)なので、固有ベクトル\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}固有ベクトルに合わせることができます。
よって、対角化する直交行列は\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}を対角化する直交行列{\bf P}と同じものとなります。

以下で、式(21)を導出します。

\begin{eqnarray}
&&{\bf P}^{-1}(\alpha{\bf I}){\bf P}+{\bf P}^{-1}(\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}){\bf P}={\rm diag}(\alpha,\ldots,\alpha)+{\rm diag}(\lambda_1,\ldots,\lambda_M)\\
&&\Leftrightarrow{\bf P}^{-1}(\alpha{\bf I}+\beta{\boldsymbol\Phi}^\top{\boldsymbol\Phi}){\bf P}={\rm diag}(\alpha+\lambda_1,\ldots,\alpha+\lambda_M)\\
&&\Leftrightarrow{\bf P}^{-1}{\bf A}{\bf P}={\rm diag}(\alpha+\lambda_1,\ldots,\alpha+\lambda_M)\tag{30}
\end{eqnarray}

偉人の名言

f:id:olj611:20210331184521p:plain:h300
できると思えばできる、できないと思えばできない。
これは、無情なゆるぎない法則である。
パブロ・ピカソ

参考文献

パターン認識機械学習 上巻 p164-p168
パターン認識機械学習の学習 p24-p25

動画

目次へ戻る