ブログ

Adam(Adaptive moment)の導出と直感的理解 〜momentum、AdaGrad〜

DeepLearningの学習でよく用いられるAdam(Adaptive moment)ですが、式が少々複雑なので直感的な理解が難しいです。当記事ではAdamを構成する主な考え方のmomentumとAdaGradなどを元にAdamを表し、その直感的な理解に関して取りまとめを行いました。

深層学習」の$3.1$節の「確率的勾配降下法(SGD)」と$3.5$節の「SGDの改良」を参考に作成を行いました。

Adamを構成する考え方

確率的勾配降下法(SGD)

$(x_1,y_1), …, (x_n,y_n)$が観測された際に$f(x_i,w)$によって$y_i$を予測することを考える。DeepLearningを考えるにあたって$f(x_i,w)$は多層パーセプトロンやCNNなどが用いられるが、$f(x_i,w) = x_i^{\mathrm{T}}w$のようにシンプルな回帰を考えてもこれは同様である。この学習にあたっては、クロスエントロピー誤差関数などを用いてパラメータ$w$の損失関数$E(w)$を定義し、$E(w)$の最小化を行う。

このときの目的は$E(w)$を最小にする$w$を計算することであり、この計算にあたって勾配降下法のような最適化のアルゴリズムが用いられる。線形回帰では正規方程式の解を計算することができたが、GLMやDeepLearningでは解を直接解くことが難しいので勾配法などが用いられる。

ここで損失関数$E(w)$はサンプル$(x_i,y_i)$ごとに計算を行うことができるので、$i$番目のサンプルの損失関数を$E_i(w)$、サンプルセット全体の損失関数を$E(w)$のように表すことを考える。

このとき$E(w)$は$E_i(w)$を用いて下記のように表すことができる。
$$
\large
\begin{align}
E(w) = \sum_{i=1}^{n} E_i(w)
\end{align}
$$

上記に対し、$E(w)$を用いて勾配法を実行することをDeepLearningの分野ではバッチ学習(batch learning)という。バッチ学習では下記のような漸化式を用いて$w$を繰り返し演算によって求める。
$$
\large
\begin{align}
w_{t+1} &= w_{t} – \alpha \nabla E(w)
\end{align}
$$

バッチ学習のような手法は局所解から抜け出せない場合もありうるので、確率的に選んだサンプルに対して学習を行うことで大域最適解を探すにあたって確率的勾配降下法(SGD)の考え方が重要になる。

確率的勾配降下法(SGD; Stochastic Gradient Descent)では全体の損失関数の$E(w)$ではなく、サンプル$(x_i,y_i)$に対する損失関数の$E_i(w)$を用いることで下記のように$w$の更新を行う。
$$
\large
\begin{align}
w_{t+1} = w_{t} – \alpha \nabla E_i(w)
\end{align}
$$

一方で確率的勾配降下法ではサンプル$1$つずつに関して計算することで、計算の並列化が難しくパラメータの収束がなかなか進まない。この解決にあたって、サンプル$1$つずつではなく数十〜数百ほどのサンプルを同時に取り扱うことで効率化が行われる。

このような複数サンプルはミニバッチ(minibatch)といわれ、現行の多くのDeepLearningの学習においてはミニバッチに基づいてパラメータの推定が行われる。

$w_{t}$の更新に用いるミニバッチを$\mathcal{D}_t$、$\mathcal{D}_t$に対応する勾配を$E_t(w)$のようにおくとき、$E_t(w)$は下記のように表される。
$$
\large
\begin{align}
E_j(w) = \sum_{i \in \mathcal{D}_t} E_i(w)
\end{align}
$$

また、ミニバッチの損失関数$E_t(w)$を元に勾配法を考えると下記のような数式で表すことができる。
$$
\large
\begin{align}
w_{t+1} &= w_{t} – \alpha \nabla E_t(w) \quad (1)
\end{align}
$$

なお、ミニバッチを用いた学習のことも確率的勾配降下法ということは抑えておくとよい。

・注意事項
参照元の「深層学習」の表記では学習率は$\epsilon$で表されるが、AdaGrad以降で$\varepsilon$が出てきて判別しにくいので当記事では$\epsilon$から$\alpha$に置き換えを行った。

momentum

確率的勾配法(SGD)を用いるにあたってはパラメータの更新毎に異なるサンプルセットを用いることで勾配がばらつくということが起こる。これにより学習が不安定になりやすいが、この対応にあたって、前のステップの修正量を用いるモメンタム(momentum)という手法がよく用いられる。

momentumではミニバッチ$\mathcal{D}_{t-1}$に対する修正量を$v_{t} = w_{t} – w_{t-1}$のように定義し、下記のような式に基づいて$w$に関する計算を行う。
$$
\large
\begin{align}
w_{t+1} = w_{t} – \alpha \nabla E_{t}(w) + \mu v_t
\end{align}
$$

上記の計算式の$w_{t+1} = w_{t} – \alpha \nabla E_{t}(w)$はミニバッチを用いるSGDと同じ式だが、momentumでは$\mu v_t = \mu(w_{t} – w_{t-1})$がさらに加算される。

ここで$\mu$は加算の割合を制御するハイパーパラメータであり、通常は$0.5$から$0.9$の範囲から選ばれる。

また、$(2)$式は$v_{t} = w_{t} – w_{t-1}$であることに基づいて下記のように表すこともできる。
$$
\large
\begin{align}
w_{t+1} &= w_{t} – \alpha \nabla E_{t}(w) + \mu v_t \\
w_{t+1} – w_{t} &= \mu v_t – \alpha \nabla E_{t}(w) \\
v_{t+1} &= \mu v_t – \alpha \nabla E_{t}(w) \quad (3)
\end{align}
$$

$\Delta w$を用いた式表記

シンプルな表記を行うにあたって$w$の$k$番目のパラメータに関して$\displaystyle \Delta w_{t,k} = w_{t+1,k}-w_{t,k}, g_{t,k} = \frac{\partial}{\partial w_k} E_t(w)$の表記を導入することを考える。ここで$\Delta w_{t,j}$は「momentum」で導入した$v_t$と同様の式であることに注意する。

確率的勾配降下法を表す$(1)$式は$\Delta w_t$を用いて下記のように書き直すことができる。
$$
\large
\begin{align}
w_{t+1,k} &= w_{t,k} – \alpha \frac{\partial}{\partial w_k} E_t(w) \\
w_{t+1,k} – w_{t,k} &= – \alpha \frac{\partial}{\partial w_k} E_t(w) \\
\Delta w_{t,k} &= – \alpha g_{t,k}
\end{align}
$$

適応的調整(AdaGrad)

AdaGradでは前項の表記に基づいて下記のような$\Delta w_{t,k}$を元に$w$の繰り返し演算を行う手法である。
$$
\large
\begin{align}
\Delta w_{t,k} = – \frac{\alpha}{\sqrt{\sum_{t’=1}^{t} g_{t’,k}^2+\varepsilon}} g_{t,k}
\end{align}
$$

上記の式は、「$k$方向の累積の修正量が多い場合は$k$方向の修正量を減らす」というように解釈することができる。また$k$方向の累積の修正量の$2$乗の和は$0$になることもあり得ることから$0$で割ることを避けるにあたって$\varepsilon$が導入されている。

このような$w$の修正量の調整は「適応的(Adaptive)」と表現され、AdaGradは”Adaptive Gradient”の略だと考えることができる。

移動平均(RMS Prop, Adadelta)

前項の「AdaGrad」では学習開始時からの$g_{t,k}^2$の総和を考えたが、これでは学習が進むにつれて修正量が$0$になるので、これを避けるにあたって総和ではなく移動平均を取る手法が考案されている。

移動平均を取る主な手法には「RMSProp」と「Adadelta」があるので、以下それぞれの移動平均の取り方に関して取りまとめ、その後に移動平均の基本式の理解に関して確認を行う。

・RMSPropで用いられる移動平均
RMSPropでは勾配の二乗和の移動平均を$\langle g_k^2 \rangle_t = \gamma \langle g_k^2 \rangle_{t-1} + (1-\gamma) g_{t,k}^2$のように定め、下記のように重み$w$の修正量の計算を行う。
$$
\large
\begin{align}
\Delta w_{t,k} = – \frac{\alpha}{\sqrt{\langle g_k^2 \rangle_t+\varepsilon}} g_{t,k}
\end{align}
$$

・Adadeltaで用いられる移動平均

・移動平均の式の理解
移動平均の式の理解にあたっては平均の漸化式的な計算を行う式に基づいて理解すればよい。移動平均に関して考える前に$a_1, a_2, …, a_n$の平均を$\bar{a}_n$とおき、$\bar{a}_{n}$を$\bar{a}_{n}$と$a_{n}$を用いて表すことを考える。
$$
\large
\begin{align}
\bar{a}_{n} &= \frac{1}{n} \sum_{i=1}^{n} a_i \\
&= \frac{1}{n} \sum_{i=1}^{n-1} a_i + \frac{1}{n} a_n \\
&= \frac{1}{n} \times \frac{n-1}{n-1} \times \sum_{i=1}^{n-1} a_i + \frac{1}{n} a_n \\
&= \frac{n-1}{n} \times \frac{1}{n-1} \sum_{i=1}^{n-1} a_i + \frac{1}{n} a_n \\
&= \frac{n-1}{n} \bar{a}_{n-1} + \frac{1}{n} a_n
\end{align}
$$

上記に対し、移動平均を$\tilde{a}_{n} = \gamma \tilde{a}_{n-1} + (1-\gamma) a_n = (1-\gamma) a_n + \gamma \tilde{a}_{n-1}$のように表すと考える。この式を元に$\tilde{a}_{n}$は下記のように導出できる。
$$
\large
\begin{align}
\tilde{a}_{n} &= (1-\gamma) a_{n} + \gamma \tilde{a}_{n-1} \quad (4) \\
&= (1-\gamma) a_{n} + \gamma ( (1-\gamma) a_{n-1} + \gamma \tilde{a}_{n-2}) \\
&= (1-\gamma) a_{n} + \gamma (1-\gamma) a_{n-1} + \gamma^2 ( (1-\gamma) a_{n-2} + \gamma \tilde{a}_{n-3}) \\
&= … \\
&=(1-\gamma) ( a_{n} + \gamma a_{n-1} + \gamma^2 a_{n-2} … )
\end{align}
$$

ここで$0 < \gamma < 1$の$\gamma$を考えると、上記は等比数列の和であることが確認でき、移動平均の定義を表すことが確認できる。

Adam

Adamの式表記

Adamでは勾配の$1$次と$2$次のモーメントをそれぞれ$m_{t,k}, v_{t,k}$のように定義し、それぞれを下記の移動平均で計算を行う。
$$
\large
\begin{align}
m_{t,k} &= \beta_{1} m_{t-1,k} + (1-\beta_{1}) g_{t,k} \quad (5) \\
v_{t,k} &= \beta_{2} v_{t-1,k} + (1-\beta_{2}) g_{t,k}^2 \quad (6)
\end{align}
$$

ここで上記は偏差を含むので、下記のように$m_{t,k}, v_{t,k}$の補正を行った$\hat{m}_{t,k}, \hat{v}_{t,k}$を定義する。
$$
\large
\begin{align}
\hat{m}_{t,k} &= \frac{m_{t,k}}{(1-\beta_{1}^t)} \\
\hat{v}_{t,k} &= \frac{v_{t,k}}{(1-\beta_{2}^t)}
\end{align}
$$

上記を元に$\Delta w_{t,k}$を下記のように考える。
$$
\large
\begin{align}
\Delta w_{t,k} = – \frac{\alpha}{\sqrt{\hat{v}_{t,k}}+\varepsilon} \hat{m}_{t,k}
\end{align}
$$

Adamでは上記を用いてパラメータ$w$の推定を行う。

Adamの式の導出

$(5)$式の導出にあたっては、$(3)$式を$(4)$式に基づいて再構成を行ったと考えることができる。$(3)$式を成分$k$に関して考える場合、$m_{t,k}$を用いて下記のように表すことができる。
$$
\large
\begin{align}
v_{t+1,k} &= \mu v_{t,k} – \alpha \frac{\partial}{\partial w_k} E_{t}(w) \quad (3)’ \\
m_{t,k} &= \mu m_{t-1,k} – \alpha g_{t,k}
\end{align}
$$

このとき上記のパラメータ$\mu, \alpha$を$(4)$式のような移動平均の漸化式で表すことを考える。
$$
\large
\begin{align}
m_{t,k} &= \mu m_{t-1,k} – \alpha g_{t,k} \\
m_{t,k} &= \beta_{1} m_{t-1,k} + (1-\beta_{1}) g_{t,k} \quad (5)
\end{align}
$$

上記のように$(5)$式が導出できる。

・注意事項
当項の導出は「深層学習」の表記から筆者が推測したものに過ぎないので、要出典であることに注意が必要である。

偏差の補正の詳細

「Adamの式表記」では$(5)$式と$(6)$式の補正を行ったが、当項では補正の詳細に関して式変形の確認を行う。Adamの論文では$(6)$式を元に解説が行われているので、以下同様に$(6)$式の導出を確認する。

$v_t$に関して$v_0=0$とすると、$(6)$式に基づいて$v_t$は下記のように変形を行うことができる。
$$
\large
\begin{align}
v_{t,k} &= (1-\beta_{2}) g_{t,k}^2 + \beta_{2} v_{t-1,k} \quad (6) \\
&= (1-\beta_{2}) g_{t,k}^2 + \beta_{2} ((1-\beta_{2}) g_{t-1,k}^2 + \beta_{2} v_{t-2,k}) \\
&= (1-\beta_{2}) g_{t,k}^2 + \beta_{2} (1-\beta_{2}) g_{t-1,k}^2 + \beta_{2}^2 ((1-\beta_{2}) g_{t-2,k}^2 + \beta_{2} v_{t-3,k}) \\
&= \cdots \\
&= (1-\beta_{2}) \sum_{i=1}^{t} \beta_{2}^{t-i} g_{i}^2 + \beta_{2}^n v_{0,k} \\
&= (1-\beta_{2}) \sum_{i=1}^{t} \beta_{2}^{t-i} g_{i}^2 \quad (7)
\end{align}
$$

ここで上記の$(7)$式の両辺の期待値を考える。
$$
\large
\begin{align}
\mathbb{E}[v_{t,k}] &= \mathbb{E} \left[ (1-\beta_{2}) \sum_{i=1}^{t} \beta_{2}^{t-i} g_{i,k}^2 \right] \quad (7) \\
&= (1-\beta_{2}) \sum_{i=1}^{t} \beta_{2}^{t-i} \mathbb{E} \left[ g_{i,k}^2 \right] \\
&= (1-\beta_{2}) \sum_{i=1}^{t} \beta_{2}^{t-i} \mathbb{E} \left[ g_{t,k}^2 \right] \\
&= \mathbb{E} \left[ g_{t,k}^2 \right] (1-\beta_{2}) \sum_{j=1}^{t} \beta_{2}^{j-1} \\
&= \mathbb{E} \left[ g_{t,k}^2 \right] (1-\beta_{2}) \frac{1-\beta_{2}^t}{1-\beta_{2}} \\
&= \mathbb{E} \left[ g_{t,k}^2 \right] (1-\beta_{2}^t) \quad (8)
\end{align}
$$

ここで$v_{t,k}$を用いて$g_{t,k}^2$の推定を行うと考えるとき、上記の$(8)$式に基づいて不偏性を考える場合、$v_{t,k}/(1-\beta_{2}^t)$を元に推定を行うと良いことがわかる。この変形が$\displaystyle \hat{v}_{t,k} = \frac{v_{t,k}}{(1-\beta_{2}^t)}$に対応する。

またここで行った変形を$m_{t,k}, \hat{m}_{t,k}$に同様に適用することができることも合わせて抑えておくと良い。

・参考
Adam論文 Section.$3$ Initialization Bias Correction

Adamの解釈・直感的理解

Adamの改良

参考

・深層学習

・Adam論文

多項分布と最尤法を用いたクロスエントロピー(cross entropy)誤差の導出

DeepLearningなどの手法で近年よく用いられるクロスエントロピー(cross entropy)誤差ですが、実は多項分布に対して最尤法を用いることで導出を行うことができます。当記事では導出の大枠と関連で出てくるソフトマックス関数に関して取りまとめを行いました。

「深層学習」の$2.4$節の「問題の定式化:出力層と損失関数の設計」を元に再構成や追記を行うことで作成を行いました。

前提の理解

ベルヌーイ分布

クロスエントロピー(cross entropy)の式の基本形はベルヌーイ分布の確率関数の積に対し、対数を取ることで導出することができる。よって、ベルヌーイ分布の確率関数の理解が重要になる。一方でベルヌーイ分布の式の形は一見少々難しく見えるので注意が必要である。

確率変数$X$が確率を表すパラメータ$p \in [0,1]$を用いたベルヌーイ分布$Bin(n=1,p)$に従うとき、$X=0, X=1$の確率$P(X=0|p), P(X=1|p)$はそれぞれ下記のように表すことができる。
$$
\large
\begin{align}
P(X=0|p) &= 1 – p \\
P(X=1|p) &= p
\end{align}
$$

上記の確率関数は$X=k, k=0,1$のように考えることで$1$つの式にまとめることができる。
$$
\large
\begin{align}
P(X=k|p) = p^{k} (1-p)^{1-k}, \quad k=0, 1
\end{align}
$$

上記の右辺に$k=1$を代入すると$p^{1} (1-p)^{1-1}=p$、$k=0$を代入すると$p^{0} (1-p)^{1-0}=1-p$がそれぞれ成立する。数式は一見難しく見えるが、$k$は$0$か$1$以外の値を取らないことに基づいて$p^{0}=1, (1-p)^{0}=1$が成立することを考えればシンプルに考えることができる。

多項分布

ベルヌーイ分布は確率変数が$2$値の場合を取り扱う確率分布だが、多項分布では確率変数が$K>2$の値を取りうる場合を取り扱う。一般的な多項分布は$2$項分布と同様に複数回の試行に基づいた確率を考えるが、ここではベルヌーイ分布に対応して$1$回の試行のみに関して考える。

多項分布$\mathrm{Mul}(1,p_1,p_2,…,p_K)$の確率関数はベルヌーイ分布と同様に考えることで下記のように表すことができる。
$$
\large
\begin{align}
P(X=k|p_1,p_2,…,p_K) = \prod_{i=1}^{K} p_{i}^{d_i}, \quad k = 1, 2, …, K
\end{align}
$$

ここで上記の$d_i$は$i=k$のときに$d_i=1$、$i \neq k$のときに$d_i=0$となる変数である。ベルヌーイ分布と同様に$p_{i}^{0}=1, p_{i}^{1}=p_{i}$のように考えられることから、上記を元に$P(X=2|p_1,p_2,…,p_K) = p_{2}, P(X=5|p_1,p_2,…,p_K) = p_{5}$のようにそれぞれ考えることができる。$\displaystyle \prod$が出てくる分さらに難しく見えるが、ほとんどが$p_{i}^{0}=1$のようになることさえ把握しておけばシンプルに理解できる。

最尤法

下記などで詳しく取り扱った。
https://www.amazon.co.jp/dp/B08FYMTYBW/

・参考
問題演習:基本的な最尤法の流れ

導出の詳細

クロスエントロピー(cross entropy)の導出

$n$対のサンプル$(x_1,y_1), …, (x_n,y_n), \quad y_i \in \{1,2,…,K\}$に対し、$y_i=k$のとき$d_{ik}=1$、$y_i \neq k$のとき$d_{ik}=0$となるように$d_{ik}$を定める。

このとき同時確率関数$P(y_1,…,y_K|p_1,p_2,…,p_K)$は下記のように表すことができる。
$$
\large
\begin{align}
P(y_1,…,y_K|p_1,p_2,…,p_K) &= \prod_{i=1}^{n} P(y_i|p_1,p_2,…,p_K) \\
&= \prod_{i=1}^{n} \prod_{k=1}^{K} p_{k}^{d_{ik}}
\end{align}
$$

上記の同時確率関数はパラメータ$p_1,p_2,…,p_K$に関する尤度$L(p_1,p_2,…,p_K)$であると考えることができるので、$L(p_1,p_2,…,p_K)$は下記のように表される。
$$
\large
\begin{align}
L(p_1,p_2,…,p_K) &= P(y_1,…,y_K|p_1,p_2,…,p_K) \\
&= \prod_{i=1}^{n} \prod_{k=1}^{K} p_{k}^{d_{ik}}
\end{align}
$$

ここで上記において、GLMなどと同様に「$p_k$は$f_k(x_i|w)$のように$x_i$にパラメータを作用させて得られる」と考えると下記のように書き直すことができる。
$$
\large
\begin{align}
L(w) &= P(y_1,…,y_K|x_1,…,x_n,w) \\
&= \prod_{i=1}^{n} \prod_{k=1}^{K} f_k(x_i|w)^{d_{ik}}
\end{align}
$$

クロスエントロピー誤差$E(w)$は上記の対数を取り、「『尤度最大化』から『誤差関数の最小化』に変えるにあたってマイナスをかけることで導出できる。
$$
\large
\begin{align}
E(w) &= – \log{L(w)} \\
&= – \log \left[ \prod_{i=1}^{n} \prod_{k=1}^{K} f_k(x_i|w)^{d_{ik}} \right] \\
&= – \sum_{i=1}^{n} \sum_{k=1}^{K} d_{ik} \log{f_k(x_i|w)}
\end{align}
$$

ここで上記の$d_{ik}$は実測値の$y_1,…,y_n$から作成し、$f_k(x_i|w)$は$x_1,…,x_n$を用いて作成した予測値に対応することは抑えておくと良い。

ソフトマックス関数

前項で定義した$f_k(x_i|w)$は確率の予測値であるので、$0$から$1$の値を取るのが望ましい。この際にロジスティック回帰などのGLMではリンク関数を設定するが、DeepLearningでも同様のことを考える。

この際によく用いられるのがソフトマックス関数(softmax function)であり、各出力層の値を$u_k(x_i|w)$とするとき、$f_k(x_i|w)$は$u_k(x_i|w)$を用いて下記のように定義される。
$$
\large
\begin{align}
f_k(x_i|w) = \frac{\exp(u_k(x_i|w))}{\sum_{j=1}^{K} \exp(u_j(x_i|w))}
\end{align}
$$

参考

・深層学習

正規分布(Normal Distribution)のパラメータと確率密度関数の形状の変化

正規分布(Normal Distribution)の平均$\mu$や分散$\sigma^2$はそれぞれ平均とばらつき度合いを表すパラメータとされますが、それぞれの値によって確率密度関数の形状がどのように変わるかに関して考察されている場合が少ないと思われます。

そこで当記事では正規分布のパラメータが変わると確率密度関数の形状がどのように変化するかに関して取りまとめを行いました。具体的には「①平均を表す$\mu$が変わると確率密度関数$f(x)$が$x$方向に平行移動する」、「②$\sigma$が変化すると確率密度関数の最大値が変わる」、「③正規分布の分散の$1/2$乗の$\sigma$が確率密度関数の変曲点に一致する」などを確率密度関数の数式を元に確認を行いました。

確率密度関数と$\mu$の変化による平行移動

正規分布$N(\mu,\sigma^2)$の確率密度関数

確率変数$X$が$X \sim N(\mu,\sigma^2)$のように正規分布$N(\mu,\sigma^2)$に従う際の確率密度関数を$f(x)$とおくと、$f(x)$は下記のように表すことができる。
$$
\large
\begin{align}
f(x) = \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left( – \frac{(x-\mu)^2}{2 \sigma^2} \right)
\end{align}
$$

$\mu$の変化と確率密度関数の平行移動

・$N(0,\sigma^2)$
$$
\large
\begin{align}
f(x) = \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left( – \frac{x^2}{2 \sigma^2} \right)
\end{align}
$$

・$N(\mu,\sigma^2)$
$$
\large
\begin{align}
f(x) = \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left( – \frac{(x-\mu)^2}{2 \sigma^2} \right)
\end{align}
$$

上記で表した$N(0,\sigma^2), N(\mu,\sigma^2)$に関する確率密度関数の式より、$N(\mu,\sigma^2)$の確率密度関数は$N(0,\sigma^2)$の確率密度関数を$x$方向に$\mu$だけ平行移動させた関数であると考えることができる。

$\sigma$の値の変化による確率密度関数の形状の変化

確率密度関数の最大値

指数関数$f(x) = \exp(x)$は単調増加の関数であり、$x$が増加するとき$f(x)$も$x$に合わせて増加する。よって、確率密度関数の$\displaystyle \frac{1}{\sqrt{2 \pi \sigma^2}} \exp \left( – \frac{x^2}{2 \sigma^2} \right)$は$\displaystyle \frac{x^2}{2 \sigma^2} \geq 0$であることより、$\displaystyle \frac{x^2}{2 \sigma^2} = 0$であるときに最大値$\displaystyle \frac{1}{\sqrt{2 \pi \sigma^2}} \exp(0) = \frac{1}{\sqrt{2 \pi \sigma^2}}$をとることがわかる。

ここで$\displaystyle \frac{1}{\sqrt{2 \pi \sigma^2}} = \frac{1}{\sqrt{2 \pi} \sigma}$を$\sigma$に関して考えると、確率密度関数の最大値が$\sigma$の逆数に比例することがわかる。

これより、「分散が大きくなるにつれ確率密度関数の最大値が小さく」なり、「分散が小さくなるにつれて確率密度関数の最大値が大きく」なることがわかる。これは分散がばらつきの度合いを表すということと直感的に反しない結果であると考えることができる。

確率密度関数の変曲点

正規分布の確率密度関数の形状は、確率密度関数の$x$に関する$1$階微分、$2$階微分の式を元に増減表を描くことで数学的に把握することができる。

以下、$\displaystyle g(x) = \exp \left( – \frac{x^2}{2 \sigma^2} \right)$に対し、$1$階微分の$g'(x)$と$2$階微分の$g^{”}(x)$を計算する。
$$
\large
\begin{align}
g'(x) &= \exp \left( – \frac{x^2}{2 \sigma^2} \right) \times \left( – \frac{x^2}{2 \sigma^2} \right)’ \\
&= – \frac{x}{\sigma^2} \exp \left( – \frac{x^2}{2 \sigma^2} \right) \\
g^{”}(x) &= – \frac{1}{\sigma^2} \exp \left( – \frac{x^2}{2 \sigma^2} \right) – \frac{x}{\sigma^2} \exp \left( – \frac{x^2}{2 \sigma^2} \right) \times \left( – \frac{x^2}{2 \sigma^2} \right)’ \\
&= – \frac{1}{\sigma^2} \exp \left( – \frac{x^2}{2 \sigma^2} \right) – \frac{x}{\sigma^2} \exp \left( – \frac{x^2}{2 \sigma^2} \right) \times \left( – \frac{x}{\sigma^2} \right) \\
&= – \frac{1}{\sigma^2} \exp \left( – \frac{x^2}{2 \sigma^2} \right) + \frac{x^2}{\sigma^4} \exp \left( – \frac{x^2}{2 \sigma^2} \right) \\
&= \frac{1}{\sigma^4} (x^2-\sigma^2) \exp \left( – \frac{x^2}{2 \sigma^2} \right) \\
&= \frac{1}{\sigma^4} (x + \sigma)(x – \sigma) \exp \left( – \frac{x^2}{2 \sigma^2} \right)
\end{align}
$$

上記より、$g(x)$に関して下記のような増減表を作成することができる。
$$
\large
\begin{array}{|c|*9{c|}}\hline x & -\infty & \cdots & -\sigma & \cdots & 0 & \cdots & \sigma & \cdots & \infty \\
\hline g'(x)& 0 & + & + & + & 0 & – & – & – & 0 \\
\hline g^{”}(x)& + & + & 0 & – & – & – & 0 & + & + \\
\hline g(x)& 0 & \nearrow & \nearrow & \nearrow & 1 & \searrow & \searrow & \searrow & 0 \\
\hline
\end{array}
$$

上記より、$g(x)$は$x=0$で最大値をとり、$x = \pm \sigma$で変曲点をとることがわかる。

また、$f(x)$を$x$に着目して考えると定数$C$を用いて$f(x)=Cg(x)$のように表せるので、$f(x)$の変曲点は$g(x)$の変曲点に一致すると考えることができる。

したがって、正規分布$N(0,\sigma^2)$の確率密度関数の変曲点は$x = \pm \sigma$であり、これを$x$方向に平行移動させた$N(\mu,\sigma^2)$の確率密度関数の変曲点は$x = \mu \pm \sigma$であると考えることができる。

参考

・連続型確率分布
https://www.hello-statisticians.com/explain-terms-cat/probdist2.html

正規方程式(normal equation)の導出と直感的理解

正規方程式(normal equation)は回帰式に対して最小二乗法(least square method)などを用いてパラメータの推定を行う際に出てくる数式です。$y = \beta_0 + \beta_1 x$のような単回帰の式の場合はシンプルな一方で、重回帰などを考えるにあたって行列表記を行うと難しそうに見えるので、当記事では行列に基づく正規方程式の導出の詳細を確認し、その直感的な理解に関して取りまとめを行いました。

内容の作成にあたっては「現代数理統計学(学術図書出版社)」の$11.1$節の「回帰モデル」と$11.2$節の「回帰モデルの推定」を主に参考にしました。

基本事項の確認

計画行列と回帰式

$$
\large
\begin{align}
y_i &= \beta_0 + \beta_1 x_i + \varepsilon_i \\
\varepsilon_i & \sim N(0,\sigma^2), \quad i.i.d.
\end{align}
$$

単回帰などを行うにあたっては、$x \in \mathbb{R}, y \in \mathbb{R}$のようなスカラー値を元に上記のように式を表す。同様に$k$次元の重回帰の式を考えると、下記のように表すことができる。
$$
\large
\begin{align}
y_i &= \beta_0 + \beta_1 x_{i1} + … \beta_k x_{ik} + \varepsilon_i \\
\varepsilon_i & \sim N(0,\sigma^2), \quad i.i.d.
\end{align}
$$

ここで上記の$2$式はどちらも$i$番目のサンプルに対応しており、さらに重回帰の式では複数変数を取り扱っており、要素を一つ一つ表すとなかなか表記が大変になる。

ここで表記をシンプルに行うにあたっては$n$対のサンプルに関する計画行列(design matrix)という行列を定義すると以下のように式の整理が行いやすい。
$$
\large
\begin{align}
y &= X \beta + \varepsilon \\
y &= \left(\begin{array}{c} y_{1} \\ \vdots \\ y_{n} \end{array} \right) \\
X &= \left(\begin{array}{cccc} 1 & x_{11} & \cdots & x_{1k} \\ \vdots & \vdots & \ddots & \vdots \\ 1 & x_{n1} & \cdots & x_{nk} \end{array} \right) \\
\beta &= \left(\begin{array}{c} \beta_{0} \\ \beta_{1} \\ \vdots \\ \beta_{k} \end{array} \right) \\
\varepsilon & \sim N_{n}(0,\sigma^2I_{n})
\end{align}
$$

ここで上記の$X$を計画行列(design matrix)、$\beta$を回帰係数ベクトルという。

残差ベクトル$e$と誤差の二乗和

計画行列と回帰式」で表したような実際に得られた$y$に対して$X \beta$を用いて計算した予測値との差を考え、$e = y – X \beta$のように定義するとき、$e$は残差ベクトルという。

また、下記で表したように$e^{\mathrm{T}}e$は誤差の二乗和に一致する。
$$
\large
\begin{align}
e^{\mathrm{T}}e &= (y – X \beta)^{\mathrm{T}} (y – X \beta) \\
&= \sum_{i=1}^{n} (y_i – (\beta_{0} + \beta_{1}x_{i1} + … + \beta_{1}x_{ik}))^{2}
\end{align}
$$

正規分布の最尤推定=最小二乗法

最小二乗法は「残差ベクトル$e$と誤差の二乗和」で確認した誤差二乗和の$e^{\mathrm{T}}e$が最小になるような$\beta$をパラメータの推定量と考える手法である。ここで「最小二乗法」によって導出される結果が「正規分布の最尤推定」の結果に一致することを以下に示す。

観測値を元に作成する計画行列$X \in \mathbb{R}^{n \times (k+1)}$と$X$に対応する$y \in \mathbb{R}^{n}$に対し、パラメータのベクトル$\beta$に関する尤度関数$L(\beta)$は$y$に関する同時確率密度関数を元に下記のように定義される。
$$
\large
\begin{align}
L(\beta) &= f(y|X,\beta,\sigma^2) \\
&= \frac{1}{(2 \pi \sigma^2)^{n/2}} \exp \left( -\frac{1}{2 \sigma^2}(y – X \beta)^{\mathrm{T}} (y – X \beta) \right) \\
&= \exp \left( -\frac{1}{2 \sigma^2}(y – X \beta)^{\mathrm{T}} (y – X \beta) – \frac{n}{2}\log{(2 \pi \sigma^2)} \right)
\end{align}
$$

上記に対して「$L(\beta)$が$\beta$に関して最大 $\iff$ $\log{L(\beta)}$が$\beta$に関して最大」である。このとき$\log{L(\beta)}$は下記のように表すことができる。
$$
\large
\begin{align}
\log{L(\beta)} &= \log \left[ \exp \left( -\frac{1}{2 \sigma^2}(y – X \beta)^{\mathrm{T}} (y – X \beta) – \frac{n}{2}\log{(2 \pi \sigma^2)} \right) \right] \\
&= -\frac{1}{2 \sigma^2}(y – X \beta)^{\mathrm{T}} (y – X \beta) – \frac{n}{2}\log{(2 \pi \sigma^2)}
\end{align}
$$

上記の$\beta$に関する最大化を考えることは$e^{\mathrm{T}}e = (y – X \beta)^{\mathrm{T}} (y – X \beta)$の$\beta$に関する最小化を考えることに一致することが式より確認できる。

よって、「最小二乗法」によって導出される結果が「正規分布の最尤推定」の結果に一致することが示される。

ベクトル・行列に関する微分の公式

・$a^{\mathrm{T}}b$の微分
$$
\large
\begin{align}
\left(\begin{array}{c} \frac{\partial}{\partial a_{1}} \\ \vdots \\ \frac{\partial}{\partial a_{q}} \end{array} \right) a^{\mathrm{T}}b &= \left(\begin{array}{c} \frac{\partial}{\partial a_{1}} \\ \vdots \\ \frac{\partial}{\partial a_{q}} \end{array} \right) \sum_{k=1}^{p} a_{k} b_{k} \\
&= \left(\begin{array}{c} b_{1} \\ \vdots \\ b_{q} \end{array} \right) = b
\end{align}
$$

・対称行列$C$に関して$a^{\mathrm{T}}Ca$の微分
$$
\large
\begin{align}
\left(\begin{array}{c} \frac{\partial}{\partial a_{1}} \\ \vdots \\ \frac{\partial}{\partial a_{q}} \end{array} \right) a^{\mathrm{T}}Ca &= 2 \left(\begin{array}{ccc} c_{11} & \cdots & c_{1q} \\ \vdots & \ddots & \vdots \\ c_{q1} & \cdots & c_{qq} \end{array} \right) \left(\begin{array}{c} a_{1} \\ \vdots \\ a_{q} \end{array} \right) \\
&= 2Ca
\end{align}
$$

詳しい導出は「現代数理統計学(学術図書出版社)」の章末課題$11.2$で取り扱った。

対称行列$X^{\mathrm{T}}X$の半正定値性

「特異値分解(SVD)」の「半正定値行列の固有値」で詳しく取り扱ったが、対称行列$X^{\mathrm{T}}X$は半正定値性を持つ。

正規方程式の導出

導出の詳細

正規分布の最尤推定=最小二乗法」の議論で、$\beta$に関する対数尤度$\log{L(\beta)}$は下記のように導出した。
$$
\large
\begin{align}
\log{L(\beta)} = -\frac{1}{2 \sigma^2}(y – X \beta)^{\mathrm{T}} (y – X \beta) – \frac{n}{2}\log{(2 \pi \sigma^2)}
\end{align}
$$

上記は$\beta^{\mathrm{T}} X^{\mathrm{T}} y, y^{\mathrm{T}} X \beta$がどちらもスカラーであり$\beta^{\mathrm{T}} X^{\mathrm{T}} y = y^{\mathrm{T}} X \beta$が成立することに基づいて下記のように変形できる。
$$
\large
\begin{align}
\log{L(\beta)} &= -\frac{1}{2 \sigma^2}(y – X \beta)^{\mathrm{T}} (y – X \beta) – \frac{n}{2}\log{(2 \pi \sigma^2)} \\
&= -\frac{1}{2 \sigma^2}(y^{\mathrm{T}}y – y^{\mathrm{T}} X \beta – \beta^{\mathrm{T}}X^{\mathrm{T}}y + \beta^{\mathrm{T}} X^{\mathrm{T}} X \beta) – \frac{n}{2}\log{(2 \pi \sigma^2)} \\
&= -\frac{1}{2 \sigma^2}(y^{\mathrm{T}}y – 2 \beta^{\mathrm{T}}X^{\mathrm{T}}y + \beta^{\mathrm{T}} X^{\mathrm{T}} X \beta) – \frac{n}{2}\log{(2 \pi \sigma^2)}
\end{align}
$$

上記に対してベクトル$\displaystyle \beta = \left(\begin{array}{c} \beta_{0} \\ \vdots \\ \beta_{p} \end{array} \right)$での微分を考える。
$$
\large
\begin{align}
\left(\begin{array}{c} \frac{\partial}{\partial \beta_{0}} \\ \vdots \\ \frac{\partial}{\partial \beta_{p}} \end{array} \right) \log{L(\beta)} = 2 X^{\mathrm{T}} X \beta – 2 X^{\mathrm{T}}y
\end{align}
$$

上記が零ベクトルに一致する$\beta$を$\hat{\beta}$とおくと、$\hat{\beta}$は下記のように表すことができる。
$$
\large
\begin{align}
2 X^{\mathrm{T}} X \hat{\beta} – 2 X^{\mathrm{T}}y &= 0 \\
X^{\mathrm{T}} X \hat{\beta} &= X^{\mathrm{T}}y \\
\hat{\beta} &= (X^{\mathrm{T}} X)^{-1} X^{\mathrm{T}}y
\end{align}
$$

ここで$\hat{\beta}$の導出に用いた$X^{\mathrm{T}} X \beta = X^{\mathrm{T}}y$を正規方程式(normal equation)とよぶ。また、上記の$\hat{\beta}$は最尤推定量に一致する。

・参考
「現代数理統計学(学術図書出版社)」 章末課題$11.3$

正規方程式の直感的理解

正規方程式の直感的な理解を行うにあたって、$y_i = \beta_{0} + \beta_{1} x_i$の式に対し、そのまま偏微分を行うパターンと$y = X \beta$を用いるパターンの$2$通りで導出を行い、それぞれの式の確認を行う。

・$y_i = \beta_{0} + \beta_{1} x_i$の微分
$y_i = \beta_{0} + \beta_{1} x_i$に関しては誤差の関数$E(\beta_{0},\beta_{1})$は下記のように表される。
$$
\large
\begin{align}
E(\beta_{0},\beta_{1}) = \sum_{i=1}^{n} (y_{i} – (\beta_{0} + \beta_{1} x_i))^2
\end{align}
$$

ここで上記の$E(\beta_{0},\beta_{1})$に対し、$\beta_{0},\beta_{1}$で偏微分を行うとそれぞれ下記のように計算できる。
$$
\large
\begin{align}
\frac{\partial E(\beta_{0},\beta_{1})}{\partial \beta_{0}} &= – 2 \sum_{i=1}^{n} (y_{i} – (\beta_{0} + \beta_{1} x_i)) \\
\frac{\partial E(\beta_{0},\beta_{1})}{\partial \beta_{1}} &= – 2 \sum_{i=1}^{n} (y_{i} – (\beta_{0} + \beta_{1} x_i)) x_i
\end{align}
$$

上記が$0$に一致すると考えるとき、それぞれの式は下記のように整理できる。
$$
\large
\begin{align}
\frac{\partial E(\beta_{0},\beta_{1})}{\partial \beta_{0}} &= 0 \\
– 2 \sum_{i=1}^{n} (y_{i} – (\beta_{0} + \beta_{1} x_i)) &= 0 \\
\sum_{i=1}^{n} (\beta_{0} + \beta_{1} x_i) &= \sum_{i=1}^{n} y_{i} \quad (1)
\end{align}
$$

$$
\large
\begin{align}
\frac{\partial E(\beta_{0},\beta_{1})}{\partial \beta_{1}} &= 0 \\
– 2 \sum_{i=1}^{n} (y_{i} – (\beta_{0} + \beta_{1} x_i)) x_i &= 0 \\
\sum_{i=1}^{n} (\beta_{0} x_i + \beta_{1} x_i^2) &= \sum_{i=1}^{n} x_{i}y_{i} \quad (2)
\end{align}
$$

・$y = X \beta$の微分
$y = X \beta$に関しては誤差の関数$E(\beta)$は下記のように表される。
$$
\large
\begin{align}
E(\beta) &= (y – X \beta)^{\mathrm{T}} (y – X \beta) \\
&= y^{\mathrm{T}}y – 2 \beta^{\mathrm{T}} X^{\mathrm{T}} y + \beta^{\mathrm{T}} X^{\mathrm{T}} X \beta
\end{align}
$$

ここで上記の$E(\beta)$に対し、ベクトル$\beta$で偏微分を行うと下記のように計算できる。
$$
\large
\begin{align}
\frac{\partial}{\partial \beta} E(\beta) = 2 X^{\mathrm{T}} X \beta – 2 X^{\mathrm{T}} y
\end{align}
$$

上記が$0$に一致すると考えるとき、それぞれの式は下記のように整理できる。
$$
\large
\begin{align}
\frac{\partial}{\partial \beta} E(\beta) &= 0 \\
2 X^{\mathrm{T}} X \beta – 2 X^{\mathrm{T}} y &= 0 \\
X^{\mathrm{T}} X \beta &= X^{\mathrm{T}} y \quad (3)
\end{align}
$$

・式の対応
ここで$(1),(2)$式を$(3)$式の形式で表すことを考える。
$$
\large
\begin{align}
\left(\begin{array}{c} \displaystyle \sum_{i=1}^{n} (\beta_{0} + \beta_{1} x_i) \\ \displaystyle \sum_{i=1}^{n} (\beta_{0} x_i + \beta_{1} x_i^2) \end{array} \right) &= \left(\begin{array}{c} \displaystyle \sum_{i=1}^{n} y_{i} \\ \displaystyle \sum_{i=1}^{n} x_{i}y_{i} \end{array} \right) \\
\left(\begin{array}{c} \displaystyle n \beta_{0} + \beta_{1} \sum_{i=1}^{n} x_i \\ \displaystyle \beta_{0} \sum_{i=1}^{n} x_i + \beta_{1} \sum_{i=1}^{n} x_i^2 \end{array} \right) &= \left(\begin{array}{c} \displaystyle \sum_{i=1}^{n} y_{i} \\ \displaystyle \sum_{i=1}^{n} x_{i}y_{i} \end{array} \right) \\
\left(\begin{array}{c} n & \displaystyle \sum_{i=1}^{n} x_i \\ \displaystyle \sum_{i=1}^{n} x_i & \displaystyle \sum_{i=1}^{n} x_i^2 \end{array} \right) \left(\begin{array}{c} \beta_{0} \\ \beta_{1} \end{array} \right) &= \left(\begin{array}{ccc} 1 & \cdots & 1 \\ x_1 & \cdots & x_n \end{array} \right) \left(\begin{array}{c} y_{1} \\ \vdots \\ y_{n} \end{array} \right) \\
\left(\begin{array}{ccc} 1 & \cdots & 1 \\ x_1 & \cdots & x_n \end{array} \right) \left(\begin{array}{cc} 1 & x_1 \\ \vdots & \vdots \\ 1 & x_n \end{array} \right) \left(\begin{array}{c} \beta_{0} \\ \beta_{1} \end{array} \right) &= \left(\begin{array}{ccc} 1 & \cdots & 1 \\ x_1 & \cdots & x_n \end{array} \right) \left(\begin{array}{c} y_{1} \\ \vdots \\ y_{n} \end{array} \right) \\
X^{\mathrm{T}} X \beta &= X^{\mathrm{T}} y
\end{align}
$$

上記のように式の対応を考えることができる。よって、$X^{\mathrm{T}} X \beta = X^{\mathrm{T}} y$の形で表される正規方程式は直感的には左辺が予測に基づく式、右辺が実測値に基づく式であるとそれぞれ解釈しておくと良い。

・参考
最小二乗法・相関係数・決定係数
https://www.hello-statisticians.com/explain-terms-cat/regression1.html

Ch.11 「線形モデル」の章末問題の解答例 〜現代数理統計学(学術図書出版社)〜

当記事は「現代数理統計学(学術図書出版社)」の読解サポートにあたってChapter.11の「線形モデル」の章末問題の解説について行います。

基本的には書籍の購入者向けの解説なので、まだ入手されていない方は購入の上ご確認ください。また、解説はあくまでサイト運営者が独自に作成したものであり、書籍の公式ページではないことにご注意ください。(そのため著者の意図とは異なる解説となる可能性はあります)

↓下記が公式の解答なので、正確にはこちらを参照ください。
https://www.gakujutsu.co.jp/text/isbn978-4-7806-0860-1/

章末の演習問題について

問題11.1の解答例

$(11.14)$式を元に$(11.13)$式は下記のように変形することができる。
$$
\large
\begin{align}
f(y) &= \frac{1}{(2 \pi \sigma^2)^{n/2}} \exp \left( -\frac{1}{2 \sigma^2}(y – X \beta)^{\mathrm{T}}(y – X \beta) \right) \\
&= \exp \left( -\frac{1}{2 \sigma^2} y^{\mathrm{T}}y + \frac{1}{\sigma^2} \beta^{\mathrm{T}} X^{\mathrm{T}}y – \frac{\beta^{\mathrm{T}} X^{\mathrm{T}} X \beta}{\sigma^2} – \frac{n}{2} \log{(2 \pi \sigma^2)} \right) \quad (1)
\end{align}
$$

上記の$(1)$式に対して、$\displaystyle \psi_{1} = -\frac{1}{2 \sigma^2}, \psi_{2} = \frac{\beta^{\mathrm{T}}}{\sigma^2}$とおき、$\displaystyle T_{1} = y^{\mathrm{T}}y, T_{2} = X^{\mathrm{T}}y, c(\psi) = \frac{n}{2} \log{(2 \pi \sigma^2)}$とおくと、$(1)$式は$(6.17)$式の指数型分布族の式に一致する。よって$(y^{\mathrm{T}}y, X^{\mathrm{T}}y)$が完備十分統計量であると考えることができる。

・参考 $(6.17)$式による指数型分布族の定義
$$
\large
\begin{align}
f(x,\psi) = h(x) \exp \left( \sum_{j=1}^{k} T_{j}(x) \psi_{j} – c(\psi) \right)
\end{align}
$$

問題11.2の解答例

・$a^{\mathrm{T}}b$に関する導出
$a^{\mathrm{T}}b$は下記のように成分を用いて表すことができる。
$$
\large
\begin{align}
a^{\mathrm{T}}b &= \left(\begin{array}{ccc} a_{1} & \cdots & a_{q} \end{array} \right) \left(\begin{array}{c} b_{1} \\ \vdots \\ b_{q} \end{array} \right) \\
&= \sum_{k=1}^{p} a_{k} b_{k}
\end{align}
$$

よって上記の$a_{i}$での偏微分は$\displaystyle \frac{\partial}{\partial a_{i}} a^{\mathrm{T}}b = b_{i}$のように考えることができる。また、これをベクトル表記することで下記が成立することも示すことができる。
$$
\large
\begin{align}
\left(\begin{array}{c} \frac{\partial}{\partial a_{1}} \\ \vdots \\ \frac{\partial}{\partial a_{q}} \end{array} \right) a^{\mathrm{T}}b &= \left(\begin{array}{c} \frac{\partial}{\partial a_{1}} \\ \vdots \\ \frac{\partial}{\partial a_{q}} \end{array} \right) \sum_{k=1}^{p} a_{k} b_{k} \\
&= \left(\begin{array}{c} b_{1} \\ \vdots \\ b_{q} \end{array} \right) = b
\end{align}
$$

・$a^{\mathrm{T}}Ca$に関する導出
$a^{\mathrm{T}}Ca$は下記のように成分を用いて表すことができる。
$$
\large
\begin{align}
a^{\mathrm{T}}Ca &= \left(\begin{array}{ccc} a_{1} & \cdots & a_{q} \end{array} \right) \left(\begin{array}{ccc} c_{11} & \cdots & c_{1q} \\ \vdots & \ddots & \vdots \\ c_{q1} & \cdots & c_{qq} \end{array} \right) \left(\begin{array}{c} a_{1} \\ \vdots \\ a_{q} \end{array} \right) \\
&= \left(\begin{array}{ccc} a_{1} & \cdots & a_{q} \end{array} \right) \left(\begin{array}{c} \sum_{l=1}^{q} c_{1j}a_{j} \\ \vdots \\ \sum_{j=1}^{q} c_{ql}a_{l} \end{array} \right) \\
&= a_{1} \sum_{l=1}^{q} c_{ql}a_{l} + … + a_{q} \sum_{l=1}^{q} c_{ql}a_{l} \\
&= \sum_{k=1}^{q} a_{k} \sum_{l=1}^{q} c_{kl}a_{l} \\
&= \sum_{k=1}^{q}\sum_{j=1}^{q} c_{kl}a_{k}a_{l} \\
&= \sum_{k=1}^{q} c_{kk}a_{k}^2 + 2 \sum_{i<j}^{q} c_{kl}a_{k}a_{l}
\end{align}
$$

上記より、偏微分$\displaystyle \frac{\partial}{\partial a_{i}} a^{\mathrm{T}}Ca$は下記のように計算できる。
$$
\large
\begin{align}
\frac{\partial}{\partial a_{i}} a^{\mathrm{T}}Ca &= \frac{\partial}{\partial a_{i}} \left( \sum_{k=1}^{q} c_{kk}a_{k}^2 + 2 \sum_{k<l}^{q} c_{kl}a_{k}a_{l} \right) \\
&= 2c_{ii}a_{i}^2 + 2 \sum_{j \neq i}^{q} c_{ij}a_{j} \\
&= 2 \sum_{j=1}^{q} c_{ij}a_{j}
\end{align}
$$

また、上記をベクトル表記することで下記が成立することも示すことができる。
$$
\large
\begin{align}
\left(\begin{array}{c} \frac{\partial}{\partial a_{1}} \\ \vdots \\ \frac{\partial}{\partial a_{q}} \end{array} \right) a^{\mathrm{T}}Ca &= \left(\begin{array}{c} \frac{\partial}{\partial a_{1}} \\ \vdots \\ \frac{\partial}{\partial a_{q}} \end{array} \right) \left( \sum_{k=1}^{q} c_{kk}a_{k}^2 + 2 \sum_{k<l}^{q} c_{kl}a_{k}a_{l} \right) \\
&= 2 \left(\begin{array}{c} \sum_{j=1}^{q} c_{1j}a_{j} \\ \vdots \\ \sum_{j=1}^{q} c_{qj}a_{j} \end{array} \right) \\
&= 2 \left(\begin{array}{ccc} c_{11} & \cdots & c_{1q} \\ \vdots & \ddots & \vdots \\ c_{q1} & \cdots & c_{qq} \end{array} \right) \left(\begin{array}{c} a_{1} \\ \vdots \\ a_{q} \end{array} \right) = 2Ca
\end{align}
$$

問題11.3の解答例

$$
\large
\begin{align}
\left(\begin{array}{c} \frac{\partial}{\partial \beta_{0}} \\ \vdots \\ \frac{\partial}{\partial \beta_{p}} \end{array} \right) Q(\beta) = 2 X^{\mathrm{T}} X \beta – 2 X^{\mathrm{T}}y
\end{align}
$$
$Q(\beta) = (y – X \beta)^{\mathrm{T}}(y – X \beta) = y^{\mathrm{T}}y – 2 \beta^{\mathrm{T}}X^{\mathrm{T}}y + \beta^{\mathrm{T}} X^{\mathrm{T}} X \beta$に対して、上記が成立することを示す。

$X^{\mathrm{T}} X$は対称行列であるので、$11.2$の導出結果より、下記がそれぞれ成立する。
$$
\large
\begin{align}
\left(\begin{array}{c} \frac{\partial}{\partial \beta_{0}} \\ \vdots \\ \frac{\partial}{\partial \beta_{p}} \end{array} \right) \beta^{\mathrm{T}}X^{\mathrm{T}}y &= X^{\mathrm{T}}y \\
\beta^{\mathrm{T}} X^{\mathrm{T}} X \beta &= 2 X^{\mathrm{T}} X \beta
\end{align}
$$

よって、下記のように考えることができる。
$$
\large
\begin{align}
\left(\begin{array}{c} \frac{\partial}{\partial \beta_{0}} \\ \vdots \\ \frac{\partial}{\partial \beta_{p}} \end{array} \right) Q(\beta) &= \left(\begin{array}{c} \frac{\partial}{\partial \beta_{0}} \\ \vdots \\ \frac{\partial}{\partial \beta_{p}} \end{array} \right) (y^{\mathrm{T}}y – 2 \beta^{\mathrm{T}}X^{\mathrm{T}}y + \beta^{\mathrm{T}} X^{\mathrm{T}} X \beta) \\
&= – 2 X^{\mathrm{T}}y + 2 X^{\mathrm{T}} X \beta \\
&= 2 X^{\mathrm{T}} X \beta – 2 X^{\mathrm{T}}y
\end{align}
$$

問題11.4の解答例

$X$の列数と同じ要素の数の任意のベクトル$a$に関して$a^{\mathrm{T}} X^{\mathrm{T}} X a = (Xa)^{\mathrm{T}} (Xa) \geq 0$が成立するので、$X^{\mathrm{T}} X$は半正定値行列であることがわかる。

以下、問題文を元に背理法を用いて示す。$X^{\mathrm{T}} X$が半正定値行列であることより、正定値行列でない場合は$a^{\mathrm{T}} X^{\mathrm{T}} X a = (Xa)^{\mathrm{T}} (Xa) = 0$が成立すると考えられ、この時に「一次独立となる場合がある」と仮定する。

このとき$(Xa)^{\mathrm{T}} (Xa) = 0$より、$Xa$は要素が$0$である零ベクトルに一致する。これは$X$の列が一次従属であることを意味する。逆に一次従属であるときは$(Xa)^{\mathrm{T}} (Xa) = 0$が成立する。

上記の議論より、「$X$の列が一次従属 $\iff$ $a^{\mathrm{T}} X^{\mathrm{T}} X a = (Xa)^{\mathrm{T}} (Xa) = 0$」が成立する。これは「$X^{\mathrm{T}} X$が正定値行列でない場合に一次独立となる場合がある」という仮定に反する。

よって「$X$の列が一次独立 $\iff$ $a^{\mathrm{T}} X^{\mathrm{T}} X a = (Xa)^{\mathrm{T}} (Xa) > 0$」が示される。

・考察
問題文にあるような背理法を用いることで逆にわかりにくいように思われた。「一次独立=$Xa$が零ベクトルでない」と考え、この際に$a^{\mathrm{T}} X^{\mathrm{T}} X a = (Xa)^{\mathrm{T}} (Xa) > 0$であることから示す方がシンプルであると思われる。

問題11.5の解答例

式の簡易化にあたって$\tau = \sigma^2$とおき、集約尤度関数(concentrated likelihood function)を表す$(11.13)$式に代入すると下記が得られる。
$$
\large
\begin{align}
f(y) &= \frac{1}{(2 \pi \tau)^{n/2}} \exp \left( -\frac{1}{2 \tau}(y – X \beta)^{\mathrm{T}}(y – X \beta) \right) \\
&= \exp \left( -\frac{1}{2 \tau}(y – X \beta)^{\mathrm{T}}(y – X \beta) – \frac{n}{2} \log{(2 \pi \tau)} \right)
\end{align}
$$

以下、$\tau$に関する偏微分を考えることで、$\log{f(y)}$を最大にする$\tau$の導出を行う。
$$
\large
\begin{align}
\frac{\partial \log{f(y)}}{\partial \tau} &= \frac{\partial}{\partial \tau} \left( -\frac{1}{2 \tau}(y – X \beta)^{\mathrm{T}}(y – X \beta) – \frac{n}{2} \log{(2 \pi \tau)} \right) \\
&= \frac{1}{2 \tau^2}(y – X \beta)^{\mathrm{T}}(y – X \beta) – \frac{n}{2} \times \frac{2 \pi}{2 \pi \tau} \\
&= \frac{1}{2 \tau^2} \left( (y – X \beta)^{\mathrm{T}}(y – X \beta) – n \tau \right)
\end{align}
$$

上記は$\tau$に関して単調減少であるので、$\displaystyle \tau = \frac{(y – X \beta)^{\mathrm{T}}(y – X \beta)}{n} = \frac{e^{\mathrm{T}}e}{n}$のときに$\log{f(y)}$と$f(y)$は最大値を取ることがわかる。

ここまでの導出により$(11.23)$式の$\displaystyle \hat{\sigma}_{ML}^{2} = \frac{e^{\mathrm{T}}e}{n}$が成立することが確認できる。

問題11.6の解答例

問題11.7の解答例

$$
\large
\begin{align}
X = \left(\begin{array}{cccc} 1 & 0 & \cdots & 0 \\ \vdots & \vdots & \vdots & \vdots \\ 1 & 0 & \cdots & 0 \\ 0 & 1 & \cdots & 0 \\ \vdots & \vdots & \vdots & \vdots \\ 0 & 1 & \cdots & 0 \\ 0 & 0 & \cdots & 0 \\ \vdots & \vdots & \vdots & \vdots \\ 0 & 0 & \cdots & 0 \\ 0 & 0 & \cdots & 1 \\ \vdots & \vdots & \vdots & \vdots \\ 0 & 0 & \cdots & 1 \end{array} \right)
\end{align}
$$

上記のように$X$を定めるとき、$X^{\mathrm{T}}X$は下記のように考えることができる。
$$
\large
\begin{align}
X^{\mathrm{T}}X &= \left(\begin{array}{cccccccccccc} 1 & \cdots & 1 & 0 & \cdots & 0 & 0 & \cdots & 0 & 0 & \cdots & 0 \\ 0 & \cdots & 0 & 1 & \cdots & 1 & 0 & \cdots & 0 & 0 & \cdots & 0 \\ \vdots & \cdots & \vdots & \vdots & \cdots & \vdots & \vdots & \cdots & \vdots & \vdots & \cdots & \vdots \\ 0 & \cdots & 0 & 0 & \cdots & 0 & 0 & \cdots & 0 & 1 & \cdots & 1 \end{array} \right) \left(\begin{array}{cccc} 1 & 0 & \cdots & 0 \\ \vdots & \vdots & \vdots & \vdots \\ 1 & 0 & \cdots & 0 \\ 0 & 1 & \cdots & 0 \\ \vdots & \vdots & \vdots & \vdots \\ 0 & 1 & \cdots & 0 \\ 0 & 0 & \cdots & 0 \\ \vdots & \vdots & \vdots & \vdots \\ 0 & 0 & \cdots & 0 \\ 0 & 0 & \cdots & 1 \\ \vdots & \vdots & \vdots & \vdots \\ 0 & 0 & \cdots & 1 \end{array} \right) \\
&= \left(\begin{array}{cccc} n_1 & 0 & \cdots & 0 \\ 0 & n_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & n_k \end{array} \right)
\end{align}
$$

よって$(X^{\mathrm{T}}X)^{-1}$は下記のように表すことができる。
$$
\large
\begin{align}
(X^{\mathrm{T}}X)^{-1} &= \left(\begin{array}{cccc} n_1 & 0 & \cdots & 0 \\ 0 & n_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & n_k \end{array} \right)^{-1} \\
&= \left(\begin{array}{cccc} 1/n_1 & 0 & \cdots & 0 \\ 0 & 1/n_2 & \cdots & 0 \\ \vdots & \vdots & \ddots & \vdots \\ 0 & 0 & \cdots & 1/n_k \end{array} \right)
\end{align}
$$

また、$X^{\mathrm{T}}Y$は下記のように考えることができる。
$$
\large
\begin{align}
X^{\mathrm{T}}Y &= \left(\begin{array}{cccccccccccc} 1 & \cdots & 1 & 0 & \cdots & 0 & 0 & \cdots & 0 & 0 & \cdots & 0 \\ 0 & \cdots & 0 & 1 & \cdots & 1 & 0 & \cdots & 0 & 0 & \cdots & 0 \\ \vdots & \cdots & \vdots & \vdots & \cdots & \vdots & \vdots & \cdots & \vdots & \vdots & \cdots & \vdots \\ 0 & \cdots & 0 & 0 & \cdots & 0 & 0 & \cdots & 0 & 1 & \cdots & 1 \end{array} \right) \left(\begin{array}{c} Y_{11} \\ Y_{12} \\ \vdots \\ Y_{1n_{1}} \\ Y_{21} \\ \vdots \\ Y_{kn_{k}} \end{array} \right) \\
&= \left(\begin{array}{c} \sum_{j=1}^{n_1} Y_{1j} \\ \sum_{j=1}^{n_2} Y_{2j} \\ \vdots \\ \sum_{j=1}^{n_k} Y_{kj} \end{array} \right)
\end{align}
$$

よって、$\hat{\beta} = (X^{\mathrm{T}}X)^{-1}X^{\mathrm{T}}Y$は下記のように表せる。
$$
\large
\begin{align}
\hat{\beta} &= (X^{\mathrm{T}}X)^{-1}X^{\mathrm{T}}Y \\
&= \left(\begin{array}{c} \bar{Y}_{1} \\ \vdots \\ \bar{Y}_{k} \end{array} \right) \\
\bar{Y}_{i} &= \frac{1}{n_i} \sum_{j=1}^{n_i} Y_{ij}
\end{align}
$$

上記は$(11.33)$式に一致する。

問題11.8の解答例

・$(11.36)$式
$$
\large
\begin{align}
\alpha_{1} + \cdots + \alpha_{k} = 0
\end{align}
$$

上記で表した$(11.36)$式の制約は、下記のような$(11.37)$式の形式で表すことができる。
$$
\large
\begin{align}
\mu &= \frac{\mu_1 + \cdots + \mu_k}{k} = \bar{\mu} \\
\alpha_i &= \mu_i – \bar{\mu}
\end{align}
$$

上記より、$\mu_1, \cdots , \mu_k$から$\mu, \alpha_1, \cdots , \alpha_{k}$が定まると考えることができる。また、逆方向も$\mu_i = \mu + \alpha_i$で定まると考えることができる。これより$1$対$1$対応が確認できる。

・$(11.39)$式
$$
\large
\begin{align}
n_1 \alpha_{1} + \cdots + n_k \alpha_{k} = 0
\end{align}
$$

上記で表した$(11.39)$式に関しても同様に下記のような$(11.40)$式の形式で表すことができる。
$$
\large
\begin{align}
\mu &= \frac{n_1 \mu_1 + \cdots + n_k \mu_k}{n} = \bar{\mu} \\
\alpha_i &= \mu_i – \bar{\mu} \\
n &= n_1 + \cdots + n_k
\end{align}
$$

上記より、$\mu_1, \cdots , \mu_k$から$\mu, \alpha_1, \cdots , \alpha_{k}$が定まると考えることができる。また、逆方向も$\mu_i = \mu + \alpha_i$で定まると考えることができる。これより$1$対$1$対応が確認できる。

問題11.9の解答例

問題11.10の解答例

問題11.11の解答例

問題11.12の解答例

問題11.13の解答例

問題11.14の解答例

問題11.15の解答例

問題11.16の解答例

問題11.17の解答例

問題11.18の解答例

問題11.19の解答例

問題11.20の解答例

問題11.21の解答例

行列の定義まとめ 〜正方行列、正則行列、正定値行列 etc〜|線形代数入門【1】

行列・線形代数を学ぶにあたって様々な行列の名称が出てきますが、正方行列、正則行列、正規行列、正定値行列、直交行列、対称行列、エルミート行列など数多くの名称が出てくるので抑えておくのがなかなか大変です。そこで当記事では行列の定義に関して取りまとめを行いました。

基本的な行列定義

正方行列(square matrix)

行と列の要素の数が同じ行列を正方行列(square matrix)という。$n$行$n$列の正方行列$A$は下記のように表すことができる。
$$
\large
\begin{align}
A &= (a_{ij}) \\
&= \left(\begin{array}{ccc} a_{11} & \cdots & a_{1n} \\ \vdots & \ddots & \vdots \\ a_{n1} & \cdots & a_{nn} \end{array} \right)
\end{align}
$$

転置行列(transposed matrix)

$m$行$n$列の行列$A$に対して$A$の$(i,j)$要素と$(j,i)$要素を入れ替えてできる$n$行$m$列の行列を転置行列(transposed matrix)という。転置行列を$A^{\mathrm{T}}$とおくと、$A^{\mathrm{T}}$は下記のように表すことができる。
$$
\large
\begin{align}
A^{\mathrm{T}} &= (a_{ji}) \\
&= \left(\begin{array}{ccc} a_{11} & \cdots & a_{1n} \\ \vdots & \ddots & \vdots \\ a_{m1} & \cdots & a_{mn} \end{array} \right)^{\mathrm{T}} \\
&= \left(\begin{array}{ccc} a_{11} & \cdots & a_{1m} \\ \vdots & \ddots & \vdots \\ a_{n1} & \cdots & a_{nm} \end{array} \right)
\end{align}
$$

随伴行列(adjoint matrix)

随伴行列は$A^{\dagger}$などで表す。行列$A$の随伴行列は$A$を転置し、その複素共役を取ったものが該当する。

行列の対称性に基づく行列の定義

対称行列(symmetric matrix)

エルミート行列

$A^{\dagger}$を随伴行列とするとき、$A=A^{\dagger}$が成立する$A$をエルミート行列という。
実対称行列の複素数に対する拡張版の概念である。

分散共分散行列

下記などで詳しく取り扱った。
https://www.hello-statisticians.com/explain-terms-cat/pca1.html

行列の直交性に基づく行列の定義

直交行列(orthogonal matrix)

ユニタリ行列(unitary matrix)

行列分解(Matrix factorization)まとめ① 〜特異値分解(SVD)〜

特異値分解(SVD; singular value decomposition)は行列を分解する手法の一つで、$n \times n$の行列にとどまらず$m \times n$の行列に用いることができます。当記事では特異値分解や関連する行列分解の手法に関して取りまとめを行いました。

特異値分解

特異値分解の概要

$$
\large
\begin{align}
\mathbf{A} = \mathbf{U} \mathbf{\Sigma} \mathbf{V}^{*}
\end{align}
$$

$m$行$n$列の行列$\mathbf{A}$に対して$m$行$m$列の$\mathbf{U}$、$m$行$n$列の$\mathbf{\Sigma}$、$m$行$n$列の$\mathbf{V}^{*}$を用いて上記のような分解を考えることができる。これを行列$\mathbf{A}$の特異値分解という。

ここで$\mathbf{\Sigma}$は$i$行$i$列の対角成分のみが値を持ち、他の値が$0$の行列である。また、$\mathbf{V}^{*}$は$\mathbf{V}$の随伴行列であるが、ここでは実数のみを考えることから以下転置行列の$\mathbf{V}^{T}$で表す。

$n \times n$対称行列の固有値分解

$n \times n$対称行列(square matrix)の$\mathbf{B}$に関して固有値分解を行うことを考える。まず、$n \times n$の単位行列を$\mathbf{I}_n$とおき、$\det(\mathbf{B}-\lambda \mathbf{E})$の相異なる解を値の大きい順に並べたものを$\lambda_1, \lambda_2, …, \lambda_n$のようにおく。このとき、$\lambda_i$は行列$\mathbf{B}$の固有値であり、$\lambda_i$に対応する固有ベクトルを$\mathbf{u}_i$のようにおくと固有値・固有ベクトルの定義より下記の式が成立する。
$$
\large
\begin{align}
\mathbf{B} \mathbf{u}_i = \lambda_i \mathbf{u}_i
\end{align}
$$

ここで固有値$\lambda_i$と固有ベクトル$\mathbf{u}_i$を用いて下記のように行列$\mathbf{\Lambda}, \mathbf{U}$を定義する。
$$
\large
\begin{align}
\mathbf{\Lambda} &= \left(\begin{array}{ccc} \lambda_{1} & \cdots & 0 \\ \vdots & \ddots & \vdots \\ 0 & \cdots & \lambda_{n} \end{array} \right) \\
\mathbf{U} &= \left(\begin{array}{ccc} \mathbf{u}_{1} & … & \mathbf{u}_{n} \end{array} \right)
\end{align}
$$

上記の$\mathbf{U}$は直交行列であり、$\mathbf{U}^{\mathrm{T}} \mathbf{U} = \mathbf{I}_{n}$が成立する。$\mathbf{U}$が直交行列となることは下記で詳しく取り扱った。

このとき、$\mathbf{B}, \mathbf{\Lambda}, \mathbf{U}$に関して下記が成立する。
$$
\large
\begin{align}
\mathbf{B} \mathbf{U} &= \mathbf{U} \mathbf{\Lambda} \\
\mathbf{B} &= \mathbf{U} \mathbf{\Lambda} \mathbf{U}^{-1} \\
&= \mathbf{U} \mathbf{\Lambda} \mathbf{U}^{\mathrm{T}}
\end{align}
$$

$n \times n$対称行列$\mathbf{B}$の固有値分解は上記のように行うことができる。本稿の主題である特異値分解はこの固有値分解の拡張で理解することができる。

$\mathbf{A} \mathbf{A}^{\mathrm{T}}$・$\mathbf{A}^{\mathrm{T}} \mathbf{A}$の計算

$m$行$n$列の行列$\mathbf{A}$に対して$\mathbf{A} \mathbf{A}^{\mathrm{T}}$を考えるとこれは$m$行$m$列の対称行列となる。また、$\mathbf{A}^{\mathrm{T}} \mathbf{A}$を考えるとこれは$n$行$n$列の対称行列となる。

$$
\large
\begin{align}
\mathbf{A} = \left(\begin{array}{ccc} a_{11} & a_{12} & a_{13} \\ a_{21} & a_{22} & a_{23} \end{array} \right)
\end{align}
$$
以下、上記のように$3$行$2$列で定義した$\mathbf{A}$を元に$\mathbf{A} \mathbf{A}^{\mathrm{T}}$と$\mathbf{A}^{\mathrm{T}} \mathbf{A}$の行列をそれぞれ確認する。

・$\mathbf{A} \mathbf{A}^{\mathrm{T}}$
$$
\large
\begin{align}
\mathbf{A} \mathbf{A}^{\mathrm{T}} &= \left(\begin{array}{ccc} a_{11} & a_{12} & a_{13} \\ a_{21} & a_{22} & a_{23} \end{array} \right) \left(\begin{array}{cc} a_{11} & a_{21} \\ a_{12} & a_{22} \\ a_{13} & a_{23} \end{array} \right) \\
&= \left(\begin{array}{cc} a_{11}^2+a_{12}^2+a_{13}^2 & a_{11}a_{21}+a_{12}a_{22}+a_{13}a_{23} \\ a_{21}a_{11}+a_{22}a_{12}+a_{23}a_{13} & a_{11}^2+a_{12}^2+a_{13}^2 \end{array} \right) \\
&= \left(\begin{array}{cc} a_{11}^2+a_{12}^2+a_{13}^2 & a_{11}a_{21}+a_{12}a_{22}+a_{13}a_{23} \\ a_{11}a_{21}+a_{12}a_{22}+a_{13}a_{23} & a_{21}^2+a_{22}^2+a_{23}^2 \end{array} \right)
\end{align}
$$

・$\mathbf{A}^{\mathrm{T}} \mathbf{A}$
$$
\large
\begin{align}
\mathbf{A}^{\mathrm{T}} \mathbf{A} &= \left(\begin{array}{cc} a_{11} & a_{21} \\ a_{12} & a_{22} \\ a_{13} & a_{23} \end{array} \right) \left(\begin{array}{ccc} a_{11} & a_{12} & a_{13} \\ a_{21} & a_{22} & a_{23} \end{array} \right) \\
&= \left(\begin{array}{ccc} a_{11}^2+a_{21}^2 & a_{11}a_{12}+a_{21}a_{22} & a_{11}a_{13}+a_{21}a_{23} \\ a_{12}a_{11}+a_{22}a_{21} & a_{12}^2+a_{22}^2 & a_{12}a_{13}+a_{22}a_{23} \\ a_{13}a_{11}+a_{23}a_{21} & a_{13}a_{12}+a_{23}a_{22} & a_{13}^2+a_{23}^2 \end{array} \right) \\
&= \left(\begin{array}{ccc} a_{11}^2+a_{21}^2 & a_{11}a_{12}+a_{21}a_{22} & a_{11}a_{13}+a_{21}a_{23} \\ a_{11}a_{12}+a_{21}a_{22} & a_{12}^2+a_{22}^2 & a_{12}a_{13}+a_{22}a_{23} \\ a_{11}a_{13}+a_{21}a_{23} & a_{12}a_{13}+a_{22}a_{23} & a_{13}^2+a_{23}^2 \end{array} \right)
\end{align}
$$

上記の例では$\mathbf{A}$に対して$\mathbf{A} \mathbf{A}^{\mathrm{T}}, \mathbf{A}^{\mathrm{T}} \mathbf{A}$がそれぞれ対称行列となることが確認できる。

また、ここで考えた$m$行$n$列の行列$\mathbf{A}$に対する$\mathbf{A} \mathbf{A}^{\mathrm{T}}, \mathbf{A}^{\mathrm{T}} \mathbf{A}$はどちらも半正定値行列である。このことは任意の$m$次元ベクトル$\displaystyle \mathbf{x} = \left(\begin{array}{c} x_{1} \\ \vdots \\ x_{m} \end{array} \right)$と$n$次元ベクトル$\displaystyle \mathbf{y} = \left(\begin{array}{c} y_{1} \\ \vdots \\ y_{n} \end{array} \right)$に対して下記が成立することよりそれぞれ確認できる。
$$
\large
\begin{align}
\mathbf{x}^{\mathrm{T}} \mathbf{A} \mathbf{A}^{\mathrm{T}} \mathbf{x} &= (\mathbf{A}^{\mathrm{T}}\mathbf{x})^{\mathrm{T}} (\mathbf{A}^{\mathrm{T}} \mathbf{x}) \geq 0 \\
\mathbf{y}^{\mathrm{T}} \mathbf{A}^{\mathrm{T}} \mathbf{A} \mathbf{y} &= (\mathbf{A}\mathbf{y})^{\mathrm{T}} (\mathbf{A}\mathbf{y}) \geq 0
\end{align}
$$

半正定値行列の固有値

半正定値行列の固有値は非負である。

前項で示した$\mathbf{A} \mathbf{A}^{\mathrm{T}}, \mathbf{A}^{\mathrm{T}} \mathbf{A}$が半正定値行列であることより、$\mathbf{A} \mathbf{A}^{\mathrm{T}}, \mathbf{A}^{\mathrm{T}} \mathbf{A}$の固有値は非負であることを示すことができる。

また、分散共分散行列に関しても$\mathbf{x}^{\mathrm{T}} \mathbf{X}^{\mathrm{T}} \mathbf{X} \mathbf{x} = (\mathbf{X} \mathbf{x})^{\mathrm{T}} (\mathbf{X} \mathbf{x}) \geq 0$のような表し方ができることから、半正定値行列であることが確認できる。

分散共分散行列の数式表現の詳細はPCAの解説の際に関連で取り扱ったので、合わせてご確認ください。

以下、半正定値行列の固有値が非負であることを示す。

$n \times n$対称行列の固有値分解の内容に基づいて$\mathbf{B}=\mathbf{U} \mathbf{\Lambda} \mathbf{U}^{\mathrm{T}}$と表せると考えるとき、$\mathbf{x}^{\mathrm{T}}\mathbf{B}\mathbf{x}$は下記のように変形することができる。
$$
\large
\begin{align}
\mathbf{x}^{\mathrm{T}}\mathbf{B}\mathbf{x} &= \mathbf{x}^{\mathrm{T}} \mathbf{U} \mathbf{\Lambda} \mathbf{U}^{\mathrm{T}} \mathbf{x} \\
&= (\mathbf{U}^{\mathrm{T}} \mathbf{x})^{\mathrm{T}} \mathbf{\Lambda} (\mathbf{U}^{\mathrm{T}} \mathbf{x})
\end{align}
$$

ここで上記の式における$\mathbf{U}^{\mathrm{T}} \mathbf{x}$を$\mathbf{y} = \mathbf{U}^{\mathrm{T}} \mathbf{x}$のようにおくことを考える。このとき$\mathbf{y}^{\mathrm{T}}\mathbf{\Lambda}\mathbf{y} = \mathbf{x}^{\mathrm{T}}\mathbf{B}\mathbf{x} \geq 0$が成立するが、$\mathbf{\Lambda}$が対角行列であることから$\mathbf{\Lambda}$の全ての値が$0$以上でなければ$\mathbf{y}^{\mathrm{T}}\mathbf{\Lambda}\mathbf{y} \geq 0$は成立しない。

よって半正定値行列の固有値は非負であると考えることができる。

特異値分解の導出

$$
\large
\begin{align}
\mathbf{A} = \mathbf{U} \mathbf{\Sigma} \mathbf{V}^{\mathrm{T}}
\end{align}
$$

行列$\mathbf{A}$が上記のように表すことができるとき、直交行列$\mathbf{U}$と$\mathbf{V}$はそれぞれ$\mathbf{A} \mathbf{A}^{\mathrm{T}}$、$\mathbf{A}^{\mathrm{T}} \mathbf{A}$の固有ベクトルに基づく行列であることを以下に示す。

$$
\large
\begin{align}
\mathbf{A} \mathbf{A}^{\mathrm{T}} &= (\mathbf{U} \mathbf{\Sigma} \mathbf{V}^{\mathrm{T}}) (\mathbf{U} \mathbf{\Sigma} \mathbf{V}^{\mathrm{T}})^{\mathrm{T}} \\
&= \mathbf{U} \mathbf{\Sigma} \mathbf{V}^{\mathrm{T}} \mathbf{V} \mathbf{\Sigma}^{\mathrm{T}} \mathbf{U}^{\mathrm{T}} \\
&= \mathbf{U} \mathbf{\Sigma} \mathbf{\Sigma}^{\mathrm{T}} \mathbf{U}^{\mathrm{T}}
\end{align}
$$

$$
\large
\begin{align}
\mathbf{A}^{\mathrm{T}} \mathbf{A} &= (\mathbf{U} \mathbf{\Sigma} \mathbf{V}^{\mathrm{T}})^{\mathrm{T}} (\mathbf{U} \mathbf{\Sigma} \mathbf{V}^{\mathrm{T}}) \\
&= \mathbf{V} \mathbf{\Sigma}^{\mathrm{T}} \mathbf{U}^{\mathrm{T}} \mathbf{U} \mathbf{\Sigma} \mathbf{V}^{\mathrm{T}} \\
&= \mathbf{V} \mathbf{\Sigma}^{\mathrm{T}} \mathbf{\Sigma} \mathbf{V}^{\mathrm{T}}
\end{align}
$$

上記を$n \times n$対称行列の固有値分解の式と同様に考えることで$\mathbf{U}$が$\mathbf{A} \mathbf{A}^{\mathrm{T}}$の固有ベクトルに基づく行列であり、$\mathbf{V}$が$\mathbf{A}^{\mathrm{T}} \mathbf{A}$の固有ベクトルに基づく行列であることがわかる。

SVDと行列のノルム

上記などを参考に行列$\mathbf{A} \in \mathbb{R}^{m \times n}$のフロベニウスノルムを$||\mathbf{A}||_{F}$、$2$-normを$||\mathbf{A}||_{2}$とおく。このときSVDにおける対角行列(diagonal matrix)$\mathbf{\Sigma} = \mathrm{diag}(\sigma_{1}, \cdots \sigma_{p}), \, p = \min{(m, n)}$の最大成分を$\sigma_{\max}$とおくと、$||\mathbf{A}||_{F}$と$||A||_{2}$は下記のように表される。
$$
\large
\begin{align}
|\mathbf{A}||_{F} &= \sqrt{\sigma_{1}^{2} + \cdots \sigma_{p}^{2}} \quad (1) \\
||\mathbf{A}||_{2} &= \sigma_{\max} \quad (2)
\end{align}
$$

フロベニウスノルムについては次節でも取り扱った。

$(1)$式の導出

$\mathbf{A}=\mathbf{U} \mathbf{\Sigma} \mathbf{V}^{T}$の$\mathbf{U}$と$\mathbf{V}$が直交行列であることから下記が成立する。
$$
\large
\begin{align}
||\mathbf{A}||_{F} &= ||\mathbf{U} \mathbf{\Sigma} \mathbf{V}^{T}||_{F} \\
&= ||\mathbf{\Sigma}||_{F}
\end{align}
$$

ここで$\mathbf{\Sigma}=\mathrm{diag}(\sigma_{1}, \cdots \sigma_{p})$より下記が成立する。
$$
\large
\begin{align}
||\mathbf{A}||_{F} &= ||\mathbf{\Sigma}||_{F} \\
&= \sqrt{\sigma_{1}^{2} + \cdots \sigma_{p}^{2}} \quad (1)
\end{align}
$$

$(2)$式の導出

$\mathbf{A}=\mathbf{U} \mathbf{\Sigma} \mathbf{V}^{T}$の$\mathbf{U}$と$\mathbf{V}$が直交行列であることから下記が成立する。
$$
\large
\begin{align}
||\mathbf{A}||_{2} &= ||\mathbf{U} \mathbf{\Sigma} \mathbf{V}^{T}||_{2} \\
&= ||\mathbf{\Sigma}||_{2}
\end{align}
$$

ここで$2$-normの定義より下記が成立する。
$$
\large
\begin{align}
||\mathbf{A}||_{2} &= ||\mathbf{\Sigma}||_{2} \\
&= \sigma_{\max} \quad (2)
\end{align}
$$

特異値分解の近似

前節では固有値分解の考え方に基づいて特異値分解の導出を確認したが、行列の行数や列数が大きい場合は固有値や固有ベクトルの計算が複雑になり、同様の計算を行うのは難しい。したがって実際に行列の分解を行うにあたっては近似的な手法が基本的に用いられる。

近似を考える際は$\mathbf{A}$の分解した行列を元に$\mathbf{A}$をどのくらい復元できるかと考えることで近似を評価すると考えれば良い。以下ではフロベニウスノルムによって近似を評価し、これに対して勾配法を適用することによって分解を行うベーシックな手順について取りまとめる。

フロベニウスノルム(Frobenius norm)

特異値分解の近似を考えるにあたっては、$\mathbf{A} = \mathbf{U} \mathbf{\Sigma} \mathbf{V}^{\mathrm{T}} = \mathbf{U}’ \mathbf{V}’$のような分解に対し、行列ノルムの$||\mathbf{A}-\mathbf{U}’ \mathbf{V}’||$を考え、与えられた$\mathbf{A}$に対してノルムを最小にするように$\mathbf{U}’, \mathbf{V}’$を繰り返し演算によって計算すれば良い。

このように行列のノルムを考える際のシンプルな考え方がフロベニウスノルム(Frobenius norm)である。$m$行$n$列の行列$\mathbf{A} = (a_{ij})$に関するフロベニウスノルム$|\mathbf{A}|_{F}$は下記のように定義できる。
$$
\large
\begin{align}
||\mathbf{A}||_{F} = \sqrt{ \sum_{i=1}^{m} \sum_{j=1}^{n} a_{ij}^2 }
\end{align}
$$

フロベニウスノルムは行列の各要素の二乗の和を計算し、計算した二乗和に対して平方根を考えた値に一致する。フロベニウスノルムはヒルベルト=シュミットノルム(Hilbert–Schmidt norm)ともいわれることも抑えておくと良い。

上記のフロベニウスノルムの定義から二つの行列$\mathbf{X}, \mathbf{Y}$の差の大小を考えるにあたって、$||\mathbf{X}-\mathbf{Y}||_{F}^{2}$を用いて考えることは妥当であると思われる。

勾配降下法による近似

前項までの内容を元に考えると、$||\mathbf{A}-\mathbf{\tilde{U}}
\mathbf{\tilde{V}}||_{F}^{2}$を最小にする$\mathbf{\tilde{U}}, \mathbf{\tilde{V}}$を計算することで特異値分解を行うことができることがわかる。

$\mathbf{\tilde{U}}$全体に関する勾配を$\nabla_{\mathbf{\tilde{U}}} ||\mathbf{A}-\mathbf{\tilde{U}}\mathbf{\tilde{V}}||_{F}^{2}$とおくと、$\nabla_{\mathbf{U}’} ||\mathbf{A}-\mathbf{\tilde{U}} \mathbf{\tilde{V}}||_{F}^{2}$は下記のように計算できる。
$$
\large
\begin{align} \nabla_{\mathbf{\tilde{U}}} ||\mathbf{A}-\mathbf{\tilde{U}}
\mathbf{\tilde{V}}||_{F}^{2} = -(\mathbf{A}-\mathbf{\tilde{U}}
\mathbf{\tilde{V}}) \mathbf{\tilde{V}}^{\mathrm{T}}
\end{align}
$$

同様に$\mathbf{\tilde{V}}$全体に関する勾配を$\nabla_{\mathbf{\tilde{V}}} ||\mathbf{A}-\mathbf{\tilde{U}}\mathbf{\tilde{V}}||_{F}^{2}$とおくと、$\nabla_{\mathbf{\tilde{V}}} ||\mathbf{A}-\mathbf{\tilde{U}} \mathbf{\tilde{V}}||_{F}^{2}$は下記のように計算できる。
$$
\large
\begin{align} \nabla_{\mathbf{\tilde{V}}} ||\mathbf{A}-\mathbf{\tilde{U}}
\mathbf{\tilde{V}}||_{F}^{2} = -\mathbf{A}^{\mathrm{T}} \mathbf{\tilde{U}} + \mathbf{\tilde{V}} \mathbf{\tilde{U}}^{\mathrm{T}} \mathbf{\tilde{U}}
\end{align}
$$

上記の導出は「関係データ学習」の$4.4.1$節の「$1$次交互勾配降下法」で取り扱われているので、詳しくはそちらを参照ください。

$l^{2}$正則化と勾配降下法

参考

・実$2$次形式と正定値行列
https://www.chaos.cs.tsukuba.ac.jp/ila/chapter5.pdf
・特異値分解とその応用
https://www.chaos.cs.tsukuba.ac.jp/ila/chapter7.pdf
・関係データ学習

超幾何分布(Hypergeometric distribution)の定義と期待値・分散

超幾何分布(Hypergeometric distribution)は$2$種からなる有限の集団から無作為非復元サンプリングを行う際に考える分布である。当記事では超幾何分布の定義や期待値・分散の導出、有限母集団修正の修正項の導出に関して取り扱った。
「基礎統計学Ⅰ 統計学入門(東京大学出版会)」の$6.1$節の「超幾何分布」や「統計学実践ワークブック」の$5$章の「離散型分布」を参考に作成を行なった。

・離散確率分布まとめ
https://www.hello-statisticians.com/explain-terms-cat/probdist1.html

超幾何分布の定義とその理解

超幾何分布の定義

$2$種の$A,B$から構成される$N$個の集団があり、$A,B$はそれぞれ$M, N-M$個であると考えるとき、この集団から$n$個取り出したとき、$A$の個数を確率変数$Y$で表すと考える。このとき、$n$個を取り出す際に復元抽出(sampling with replacement)を行う場合は$Y$は二項分布$\displaystyle \mathrm{Bin} \left( n,\frac{M}{N} \right)$に従う一方で、非復元抽出(sampling without replacement)を行う場合は$Y$は超幾何分布$\mathrm{HyperGeo}(N,M,n)$に従う。

以下、非復元抽出の場合を考えるにあたって、超幾何分布の確率関数を$P(Y=y|N,M,n)$とおくと、$P(Y=y|N,M,n)$は下記のように表される。
$$
\large
\begin{align}
P(Y=y|N,&M,n) = \frac{{}_M C_{y} \times {}_{N-M} C_{n-y}}{{}_{N} C_{n}} \\
{}_a C_b &= 0, \quad if \quad a<b \quad or \quad b<0
\end{align}
$$

上記の式の解釈にあたっては、「集団全体の$N$から$n$個抽出する組み合わせ」と「($A$に関して$M$から$y$個抽出する組み合わせ) × ($B$に関して$N-M$から$n-y$個抽出する組み合わせ)」の割合であると考えると良い。またここでの問題の定義より、${}_a C_b$に対して$a < b$と$b < 0$が成立するときはないので${}_a C_b = 0$と定義した。元の定義上は$a, b$の定義域を定める方が直接的だが、$a < b$と$b < 0$のときは${}_a C_b = 0$と定める方が取り扱いがシンプルになるのでこの式の形で抑えておく方がわかりやすいと思われる。

超幾何分布に基づく具体的な確率の確認

以下、「基礎統計学Ⅰ 統計学入門(東京大学出版会)」の例に基づいて確認を行う。$N=1000, M=200, n=5$のとき$Y=0, Y=1$のときの確率関数$P(Y=0|N=1000,M=200,n=5), P(Y=1|N=1000,M=200,n=5)$の計算をそれぞれ行う。

・$P(Y=0|N=1000,M=200,n=5)$
$$
\large
\begin{align}
P(Y=0|&N=1000,M=200,n=5) = \frac{{}_{200} C_{0} \times {}_{1000-200} C_{5-0}}{{}_{1000} C_{5}} \\
&= \frac{800 \cdot 799 \cdot 798 \cdot 797 \cdot 796}{5!} \times \frac{5!}{1000 \cdot 999 \cdot 998 \cdot 997 \cdot 996} \\
&= 0.326859…
\end{align}
$$

・$P(Y=1|N=1000,M=200,n=5)$
$$
\large
\begin{align}
P(Y=1|&N=1000,M=200,n=5) = \frac{{}_{200} C_{1} \times {}_{1000-200} C_{5-1}}{{}_{1000} C_{5}} \\
&= \frac{200}{1!} \times \frac{800 \cdot 799 \cdot 798 \cdot 797}{4!} \times \frac{5!}{1000 \cdot 999 \cdot 998 \cdot 997 \cdot 996} \\
&= 0.410626…
\end{align}
$$

超幾何分布を用いた資源調査

超幾何分布の期待値・分散の計算と有限母集団修正

$i$回目の抽出でグループ$A$を引いた場合に$X_i=1$、グループ$B$を引いた場合に$X_i=0$のように確率変数$X_i$を定義する。このとき確率変数$Y$は$Y = X_1 + X_2 + \cdots + X_n$のように表せることに基づいて超幾何分布の期待値$E[Y]$や分散$V[Y]$の導出を行える。

超幾何分布の期待値$E[Y]$の計算

$E[X_i]$は次のように表せる。
$$
\large
\begin{align}
E[X_i] &= 0 \cdot P(X_i=0) + 1 \cdot P(X_i=1) \\
&= P(X_i=1) \\
&= \frac{M}{N}
\end{align}
$$

よって、超幾何分布の期待値$E[Y]$は下記のように考えることができる。
$$
\large
\begin{align}
E[Y] &= E[X_1 + X_2 + \cdots + X_n] \\
&= nE[X_i] \\
&= \frac{nM}{N}
\end{align}
$$

超幾何分布の分散$V[Y]$の計算

$E[X_i^2], E[X_iX_j], \, i \neq j$は次のように表せる。
$$
\large
\begin{align}
E[X_i^2] &= 0^2 \cdot P(X_i=0) + 1^2 \cdot P(X_i=1) \\
&= P(X_i=1) \\
&= \frac{M}{N} \\
E[X_iX_j] &= 0 \cdot (P(X_i=0,X_j=0)+P(X_i=1,X_j=0)+P(X_i=0,X_j=1)) + 1 \cdot P(X_i=1,X_j=1) \\
&= P(X_i=1,X_j=1) \\
&= \frac{M(M-1)}{N(N-1)}
\end{align}
$$

このとき、確率変数$X_i, X_j$に関して下記のように$V[X_i], \mathrm{Cov}(X_i,X_j)$が得られる。
$$
\large
\begin{align}
V[X_i] &= E[X_i^2] – E[X_i]^2 \\
&= \frac{M}{N} – \left( \frac{M}{N} \right)^2 \\
&= \frac{M(N-M)}{N^2} \\
\mathrm{Cov}(X_i,X_j) &= E[X_iX_j] – E[X_i]E[X_j] \\
&= \frac{M(M-1)}{N(N-1)} – \frac{M^2}{N^2} \\
&= \frac{MN(M-1)}{N^2(N-1)} – \frac{M^2(N-1)}{N^2(N-1)} \\
&= \frac{\cancel{M^2N} – MN – \cancel{M^2N} + M^2}{N^2(N-1)} \\
&= -\frac{M(N-M)}{N^2(N-1)}
\end{align}
$$

よって、超幾何分布の分散$V[Y]$は下記のように考えることができる。
$$
\large
\begin{align}
V[Y] &= V[X_1 + X_2 + \cdots + X_n] \\
&= nV[X_i] + 2 \cdot {}_{n} C_{2} \mathrm{Cov}(X_i,X_j) \\
&= nV[X_i] + n(n-1) \mathrm{Cov}(X_i,X_j) \\
&= \frac{nM(N-M)}{N^2} – \frac{n(n-1)M(N-M)}{N^2(N-1)} \\
&= \frac{nM(N-M)}{N^2(N-1)} [(N-\cancel{1})-(n-\cancel{1})] \\
&= \frac{nM(N-M)(N-n)}{N^2(N-1)} \\
&= n \cdot \frac{M}{N} \left( 1-\frac{M}{N} \right) \times \frac{N-n}{N-1}
\end{align}
$$

有限母集団修正

「復元抽出」の場合は二項分布を元に考えることができる。$Z \sim \mathrm{Bin}(n,p)$の分散を$V[Z]$とおくと$V[Z]=np(1-p)$のように表せる。ここで$p=M/N$のように表すとき、$V[Z]$は下記のように表せる。
$$
\large
\begin{align}
V[Z] &= np(1-p) \\
&= n \cdot \frac{M}{N} \left( 1-\frac{M}{N} \right)
\end{align}
$$

このとき$V[Y]/V[Z]$は下記のように計算できる。
$$
\large
\begin{align}
\frac{V[Y]}{V[Z]} = \frac{N-n}{N-1}
\end{align}
$$

上記が有限母集団修正にあたっての修正項に一致する。

参考

・統計学実践ワークブック$5$章: 離散型分布
・統計学実践ワークブック$5$章: 問題$5.3$

教師あり学習・教師なし学習・強化学習はそれぞれどのように理解すべきか

機械学習(machine leanring)の大別にあたって「教師あり学習(Supervised Learning)」、「教師なし学習(Unsupervised Learning)」、「強化学習(Reinforcement Learning)」の三つが挙げられることが多いですが、入門コンテンツの多くはミスリードになるような解説が多く、特に「教師なし学習」に関して論文などで調査を行う際などにわからなくなる場合があります。当記事ではその解決にあたって、それぞれの定義に関して詳しく確認を行います。

基本的な考え方

「教師あり学習(Supervised Learning)」、「教師なし学習(Unsupervised Learning)」、「強化学習(Reinforcement Learning)」の三つが機械学習の大別に挙げられることが多いですが、「教師なし学習」を「クラスタリング」などのアルゴリズムと対応させて解説を行うことでミスリードになりやすいので注意が必要です。このような解説がなされることで「教師あり/教師なし」に対応してそれぞれアルゴリズムがあると考えがちですが、この考え方は適切ではありません。

「教師あり/教師なし」に対応するのは「アルゴリズム・学習の仕組み」ではなく、「前処理も含んだ処理全体」です。たとえばDeep Learningを考えた際に、下記のようにそれぞれ「教師あり学習」、「教師なし学習」、「強化学習」の三つの考え方が対応します。

・教師あり学習
-> ImageNetの分類 etc
・教師なし学習
-> Auto Encoder、GAN、BERT etc
・強化学習
-> Deep Q-Network etc

上記のようにDeep Learningを考えた際に「教師あり学習」、「教師なし学習」、「強化学習」の三パターンを考えられることから、「ニューラルネットワーク」のような「学習の仕組み」に基づいて三つの大別を行うことができません。

ここでDeepLearningなどを考える際によく出てくるのが「人間がアノテーションを与えるかどうか」という視点です。たとえばImageNetの分類では画像$X$に対し、正解ラベル$y$を人間が作成し学習を行いますが、Auto EncoderやGANは画像$X$があれば学習を行うことができます。

教師あり学習と教師なし学習(Auto Encoder)

上図で「通常の教師あり学習」と「教師なし学習に対応するAuto Encoder」を元に具体的な図示を行いましたが、ニューラルネットワークの入力はどちらも$X$である一方で、出力は人間がラベル付けを行なった(human made)$y$と、入力と同様な$X$でそれぞれ異なります。

ここで着目すべきは手法自体はどちらもニューラルネットワークである一方で、出力にラベルを与えたかどうかで「教師あり/教師なし」を判断する場合も多いことです。したがって「ニューラルネットワーク」のような「学習のアルゴリズム」で判断するよりはfit(X,y)と書くかfit(X)と書くかのように関数のように理解する方が良いと思います。fit(X)と書く場合は「教師なし学習」と考える一方で、内部で教師に該当するyを生成しfit(X,y)を実行する場合もあり、この場合はfit(X)の視点で見れば「教師なし学習」、fit(X,y)の視点で見れば「教師あり学習」のように考えることができると思います。

次節で「教師なし学習」を詳しく確認するにあたってBERTの事前学習(pre-training)を取り扱いますが、BERTのMLM(Masked Language Modeling)を題材にすることで「教師なし学習」がより理解できると思います。

ここまで「教師あり/教師なし」に関して確認を行いましたが、「強化学習」も同様に「強化学習に対応するアルゴリズムがある」というよりは、$X$の代わりにエージェントと呼ばれるプレイヤーを定義して学習を行う場合を「強化学習」と呼ぶという理解の方が良いと思います。

強化学習に関しても詳しく確認すると難しくなるので、詳細は次節で確認を行います。

詳細

教師あり学習

教師なし学習

前節では$X$をニューラルネットワークの入力と出力に用いるAuto Encoderを元に「教師なし学習」を確認しましたが、当項ではBERTの事前学習(pre-training)を例に「教師を乱数を用いて”機械的に”作る」場合の確認を行います。

BERTの事前学習であるMLM(Masked Language Modeling)では文章に対応する系列を取り扱うDeepLearningの主流の手法であるTransformerに対し、一定確率でランダムに$X$にあたる入力を消去し、これを出力であてるというタスク設定を行います。

より具体的には$X = (x_1,x_2,x_3,x_4,…) = ($私, は, 本屋, に, 行った, 。, 教科書, を, 買った, 。$)$という文があった際に入力を$($私, は, MASK, に, 行った, 。, 教科書, を, 買った, 。$)$に変化させ、出力で「本屋」を予測するようなタスクをMLM(Masked Language Modeling)では取り扱うと理解すれば良いです。少々分かりにくいので図で表すと下記のように表すことができます。

BERTにおけるMasked Language Modelingの学習の全容

ここで上記で表したMasked Language Modelingは入力$X$に対して$y$を予測している問題のようにも見えるので「教師あり学習」であると考えがちですが、出力にあたる「本屋」は乱数などを用いて作成したものであり人手でアノテーションを行なったわけではありません。

このような例もあるなど、DeepLearningを考える際には「教師あり/教師なし」という用語は単に「アノテーションを作成するコスト」に着目することが多いです。この解釈にあたっては、「ニューラルネットワーク」単体で見れば常にfit(X,y)であるのに対して、BERTのように前処理まで含めばfit(X)であり教師なし学習という見方ができるとわかりやすいのではないかと思います。

強化学習

強化学習(Reinforcement Learning)に関しては直感的にわかりやすい定義が少ないですが、意思決定主体のエージェント(Agent)を定義し「エージェントの最適行動を考える問題」と理解すると良いと思います。

もう少し詳しく考えると「エージェントが何らかの状態(state)を観測した際に意思決定(action)を行いその結果新しい状態と報酬を観測する場合、一連の意思決定に基づく報酬を最大化する意思決定を考える問題」と考えると良いです。

これに関連して上図のような図がよく用いられますが、この図を元に理解するにあたっては下記に注意すると良いです。

1) 強化学習は逐次的意思決定問題(Sequential Decision Making Problem)を取り扱う
-> 強化学習は系列処理の問題に意思決定の概念を導入し、意思決定によって発生する報酬を元に意思決定方針の最適化を考えると理解しやすいです。図の理解にあたっては1つの時点における意思決定に関連する報酬のみを表示しており、実際は逐次処理の一部であることには注意が必要です。
より具体的には「テレビゲームを行う際に、画面を見てコントローラーでコマンドを入力し、状況が遷移する際に報酬が発生すれば獲得する」を1フレームと見た際に、この1フレームの概略が図に対応すると考えれば良いと思います。

2) 大まかに理解する際は環境(Environment)ではなくエージェント(Agent)ありきで考えるべき
-> 「強化学習には環境とエージェントがある」という解説がされやすいですが、方針の最適化を行うのは環境ではなくエージェントであることは注意しておくと良いです。環境は迷路や囲碁などのようにはっきり定まる場合もあればAtariゲームのように何らかの定義が必要な場合もあります。
したがって、大まかに把握する際はエージェントが自身を取り巻く環境に対し何らかの意思決定を行い報酬などの反応を通して方針を学習すると考えると良いです。

3) DeepLearningなどを用いる際はエージェントの意思決定をDeepLearningを用いて近似する
-> 1)、2)で取り扱ったように強化学習ではエージェントの最適行動方針を学習します。この考えに基づいてDeepLearningを導入するにあたっては、エージェントの最適な意思決定をDeepLearningでいかに近似するかがメインの課題になると理解しておくと良いと思います。

強化学習に関しては下記などで別途取りまとめを行なっているので下記も参考にしてみてください。
https://www.hello-statisticians.com/explain-terms-cat/reinforcement_learning1.html

また、下記のChapter$1$などに強化学習の定義に関してまとめられているので、専門的・学術的なリファレンスに用いる場合は下記を参照すると良いと思います。

その他

半教師あり学習(semi-supervised learning)

深層学習[第2版]の11.4節の内容を元に以下取りまとめます。半教師あり学習(semi-supervised learning)はラベル付きの訓練データ$\mathcal{D}_L$とラベルのない入力のみのデータ$\mathcal{D}_{UL}$がある際に両者を効果的に用いて目的タスクを学習する方法です。

半教師ありの手法に関しては、主に一致性正則化(consistency regularization)、擬似ラベル(pseudo label)、エントロピー最小化(entropy minimization)の$3$つの手法があるので以下それぞれに関して取りまとめを行います。

・一致性正則化(consistency regularization)
一致性正則化(consistency regularization)はラベルのないサンプルの$\mathbf{x}$を中身が変わらない範囲で$\mathbf{x} \to \mathbf{x}’$のように変動させるとき、予測結果が変わらないであろうことに基づく手法です。

関数$f$を元に出力の近似を行う際に出力の変化を測る尺度を$\delta(f(\mathbf{x}),f(\mathbf{x}’))$のようにおき、これに基づいて正則化項を構成し誤差関数に加えることで$\mathcal{D}_{UL}$に基づいて正則化を行うことができると考えることができます。

・擬似ラベル(pseudo label)
擬似ラベル(pseudo label)は先に$\mathcal{D}_L$を用いて教師ネットワーク$T$を学習させたのちに、$\mathcal{D}_{UL}$を教師ラベルを用いて分類し、生徒ネットワーク$S$を教師ネットワークの出力に基づいて学習させるという手法です。

・エントロピー最小化(entropy minimization)
エントロピー最小化(entropy minimization)はネットワークの出力結果の各クラスのスコアのエントロピー$\displaystyle – \sum_{k=1}^{K} f(\mathbf{x})_k \log{f(\mathbf{x})_k}$を最小にする手法です。予測結果の事後確率の分布が1つのクラスのみが大きくなるときエントロピーは小さくなるので、$\mathcal{D}_{UL}$の予測のエントロピーが小さくなるように学習させると考えると良いです。

まとめ

参考

・機械学習中級(筆者作成)
https://www.amazon.co.jp/gp/product/B08CRSY7XN/

・深層学習[第2版]

強化学習まとめ② 〜動的計画法 (Dynamic Programming)〜

強化学習(Reinforcement Learning)は数式が多く一見難しそうに見えますが、統計学の確率分布の表記に慣れていれば実はそれほど難しくありません。そこで何回かにわけて強化学習の基本トピックに関して取り扱いを行います。
第$2$回は最適ベルマン方程式の導出にあたって用いる「動的計画法(Dynamic Programming)」を取り扱います。Richard S. Suttonの”Reinforcement Learning, second edition: An Introduction“のChapter$4$を参考に作成を行いました。

基本的な理論の確認

概要

動的計画法(Dynamic Programming)はマルコフ決定過程(MDP; Markov Decision Process)の完全な環境モデル(perfect model of environment)が与えられるとき最適方策(optimal policy)を計算するのに用いられる手法である。

ここで「完全な環境モデル(perfect model of environment)」は、第$1$回で取り扱った状態遷移確率(state-transition probabilities)の$p(s’,r|s,a) = P(S_{t+1}=s’,R_{t+1}=r|S_t=s,A_{t}=a)$が対応すると考えておけば良い。要するに状態$S_t$とエージェントの行動$A_t$が定まった際に、次の状態$S_{t+1}$と同時に発生する報酬$R_{t+1}$の確率分布が既知である場合をこのように考える。

この環境モデル(environment model)の考え方は「モデルフリー(model free)強化学習」と「モデルベース(model based)強化学習」にも関連する。「モデルベース(model based)強化学習」は近年ではalpha Goやalpha Zeroなどのボードゲームにおける強化学習でよく研究されており、環境モデルから最適方策を導出する方法をプランニング(planning)と表すことも抑えておくと良い。

動的計画法を考える際の主要な考え方は価値関数(value function)を「良い方策(good policy)」を探索する際に用いるということである。ここで用いる価値関数は第$1$回で確認したベルマン方程式(Bellman equation)の数式の形を用いる。
$$
\large
\begin{align}
v_{\pi}(s) &= \mathbb{E}_{\pi} [G_{t}|S_{t}(s)] \\
&= \mathbb{E}_{\pi} [R_{t+1} + \gamma v_{\pi}(S_{t+1})|S_{t}(s)] \\
&= \sum_{a} \pi(a|s) \sum_{s’,r} p(s’,r|s,a) \left[ r + \gamma v_{\pi}(s’) \right] \quad (1)
\end{align}
$$

動的計画法のアルゴリズムは上記で表したベルマン最適方程式の式を用いて表現を行う。式の中に出てきた$\pi(a|s)$と$p(s’,r|s,a)$は期待値を考える際の確率分布と同様に理解しておくと良い。

より具体的には「動的計画法」に基づいて「方策評価」を行い、その「方策評価」に基づいて「方策反復法(Policy Iteration)」や「価値反復法(Value Iteration)」を構築し、繰り返し計算によって「最適方策」の計算を行う。

動的計画法の再確認

前項では「強化学習」における「動的計画法」の活用にあたっての概要を確認を行なったが、「強化学習」に関連して「動的計画法(Dynamic Programming)」の解説がなされることは少ないように思われる。

「動的計画法(Dynamic Programming)」に関してはWikipediaを元に考えるが、「対象となる問題を複数の部分問題に分割し、部分問題の計算結果を記録しながら解いていく手法を総称」と記載されている。

この具体例に「フィボナッチ数列」が紹介されているが、$a_5 = a_4 + a_3 = (a_3+a_2) + (a_2+a_1) = …$のように考えるのではなく、$a_2 = a_1 + a_0, a_3 = a_2 + a_1$のように考えることを動的計画法を用いたプログラムに挙げられている。

要するに、漸化式を元に$a_0, a_1, a_2, a_3, …$のように順に値を計算していくことを動的計画法(Dynamic Programming)と総称していると理解すると良いように思われる。

ここで次項で扱う「反復方策評価」は漸化式を元に値を順に計算する手法であるので、動的計画法に基づくと考えられる。

一方で「動的計画法」は多義的な表現であり、もう少し適切な表現の方があるようにも思われる。これに関しては、「ベルマン方程式」を考案したベルマン(Richard Ernest Bellman)が「動的計画法」も同時に発表していたことに起因すると思われる。

動的計画法に関しては「トップダウン方式(分割統治法)」と「ボトムアップ方式」の二つがあり、「分割統治法」をメインで紹介されることも多いのでこの辺はミスリードになる場合もあるように思われた。「ボトムアップ方式」は漸化式に関連する問題の解法を考えるにあたって直感的にそのまま用いる解法であるので、「動的計画法」に関して学ぶ際には既知である場合も多いと考えられるが、強化学習の動的計画法は単にシンプルな「ボトムアップ方式」を用いていると考えれば良いと思われる。

当項は単に筆者の考察であり、出典が存在しないのでご注意ください。また、関連する文献をご存知の方がおられたらご指摘ください。

・参考
動的計画法
リチャード・E・ベルマン
ベルマン方程式

反復方策評価(Iterative Policy Evaluation)

反復方策評価(Iterative Policy Evaluation)は$1.1$節の$(1)$式に対して動的計画法の考え方に基づいて繰り返し演算を適用する手法である。
$$
\large
\begin{align}
v_{k+1}(s) & \equiv \mathbb{E}_{\pi} [R_{t+1} + \gamma v_{k}(S_{t+1})|S_{t}(s)] \\
&= \sum_{a} \pi(a|s) \sum_{s’,r} p(s’,r|s,a) \left[ r + \gamma v_{k}(s’) \right] \quad (2)
\end{align}
$$

反復方策評価は具体的には上記で表した繰り返し演算に基づいて価値関数の$v_{\pi}(s)$の近似を行う。ここで$S_{t}=s \to S_{t+1}=s’$が状態の変化を表す一方で、$k \to k+1$は価値関数の変化を表すことに注意が必要である。

方策改善(Policy Improvement)

反復方策評価によって価値関数$v_{\pi}$が収束した際に、価値関数に基づいて方策を改善することが必要になる。この一連の流れを方策改善(Policy Improvement)という。

方策改善は方策改善定理(Policy Improvement Theorem)に基づいており、価値関数に基づいて方策を改善することは既存の方策が最適でない限りは既存のものよりもより良い方策が得られるという考え方である。

方策反復法(Policy Iteration)

方策反復法は方策評価(Policy Evaluation)と方策改善(Policy Improvement)を繰り返す手法である。アルゴリズムは下記のように表される。

上記の一連の処理を方策反復法(Policy Iteration)という。一連の流れは下記のようにも概略を表すことができる。
$$
\large
\begin{align}
\pi_0 \to v_{\pi_0} \to \pi_1 \to v_{\pi_1} \to … \to \pi_{*} \to v_{\pi_{*}}
\end{align}
$$

価値反復法(Value Iteration)

まとめ

具体例の確認