TadaoYamaokaの開発日記

個人開発しているスマホアプリや将棋AIの開発ネタを中心に書いていきます。

Pythonで3層パーセプトロンの誤差逆伝播を実装してみる

DeepLearning Python 機械学習

入力層、隠れ層、出力層で構成される単純なニューラルネットワークで、
誤差逆伝播を計算します。

f:id:TadaoYamaoka:20160410125104p:plain

隠れ層の活性化関数はsigmoid、
出力層の活性化関数はsoftmaxとします。

誤差関数(損失関数)は、交差エントロピーを使用します。

それぞれの式は以下の通りです。

sigmoid:\displaystyle \frac{1}{1+\mathrm{e}^{-{\bf u}}}

softmax:\displaystyle \frac{\mathrm{e}^{u_j}}{\sum_{k} \mathrm{e}^{u_k}}

softmaxの交差エントロピー:-\sum_{k} d_k \log(y_k)

ここで、{\bf u}は各層の入力ベクトル、{\bf d}は教師データのベクトル、{\bf y}は出力層の出力ベクトルとします。

ネットワークのパラメータの行列は、
1層目を{\bf W^1}=\begin{pmatrix} w_{11}^1 & w_{12}^1 \\ w_{21}^1 & w_{22}^1 \end{pmatrix}、
2層目を{\bf W^2}=\begin{pmatrix} w_{11}^2 & w_{12}^2 \\ w_{21}^2 & w_{22}^2 \end{pmatrix}で表します。

行列の成分w_{ji}^nの添え字は、jが出力側のベクトル成分、iが入力側のベクトル成分を表します。
上付きの数字nは、何番目の層かを表します。

順伝播

順伝播は、入力{\bf x}から出力{\bf y}を以下の式で計算します。

\displaystyle {\bf u^1} = {\bf W^1} {\bf x} = \begin{pmatrix} w_{11}^1 & w_{12}^1 \\ w_{21}^1 & w_{22}^1 \end{pmatrix} \begin{pmatrix} x_1 \\ x_2 \end{pmatrix} = \begin{pmatrix} w_{11}^1 x_1 + w_{12}^1 x_2 \\ w_{21}^1 x_1 + w_{22}^1 x_2 \end{pmatrix}

\displaystyle {\bf z^1} = sigmoid({\bf u^1})

\displaystyle {\bf u^2} = {\bf W^2} {\bf z^1} = \begin{pmatrix} w_{11}^2 & w_{12}^2 \\ w_{21}^2 & w_{22}^2 \end{pmatrix} \begin{pmatrix} z_1^1 \\ z_2^1 \end{pmatrix} = \begin{pmatrix} w_{11}^2 z_1^1 + w_{12}^2 z_2^1 \\ w_{21}^2 z_1^1 + w_{22}^2 z_2^1 \end{pmatrix}

\displaystyle {\bf z^2} = sofmax({\bf u^2})

\displaystyle {\bf y} = {\bf z^2}

出力{\bf y}は、出力層の活性化関数にsoftmaxを使用しているので、ベクトルの各成分がクラスを表し、値がそのクラスに属する確率を示します。
入力{\bf x}が与えられたときにクラスを分類する識別器として機能します。


Pythonでは以下の通り記述します。

import numpy as np

def sigmoid(u):
    return 1 / (1 + np.exp(-u))

def softmax(u):
    e = np.exp(u)
    return e / np.sum(e)

def forward(x):
    global W1
    global W2
    u1 = x.dot(W1)
    z1 = sigmoid(u1)
    u2 = z1.dot(W2)
    y = softmax(u2)
    return y, z1

xはnumpyのarrayオブジェクトで、入力ベクトルの成分をx1,x2とすると、

x = np.array([[x1, x2]])

となります。

W1、W2はネットワークのパラメータの行列で、

W1 = np.array([[w1_11, w1_21], [w1_12, w1_22]])

のようになります。W2も同様です。
上の方で説明した{\bf W^1}を転置した形になっていることに注意してください。
これは入力ベクトルを横ベクトルで表しているため計算の都合上、転置しています。


x.dot(W1)は、xを横ベクトル、W1を転置した形で表しているので、
\displaystyle \begin{pmatrix} x_1 & x_2 \end{pmatrix} \cdot \begin{pmatrix} w^1_{11} & w^1_{21} \\ w^1_{12} & w^1_{22} \end{pmatrix} = \begin{pmatrix} x_1 w^1_{11} + x_2 w^1_{12} & x_2 w^1_{21} + x_2 w^1_{22} \end{pmatrix}
を計算しています。.dotは行列の積(ドット積)を表します。

同様に、z1.dot(W2)は、z1の成分をz^1_1,z^1_2とすると、
\displaystyle \begin{pmatrix} z^1_1 & z^1_2 \end{pmatrix} \cdot \begin{pmatrix} w^2_{11} & w^2_{21} \\ w^2_{12} & w^2_{22} \end{pmatrix} = \begin{pmatrix} z^1_1 w^2_{11} + z^1_2 w^2_{12} & z^1_2 w^2_{21} + z^1_2 w^2_{22} \end{pmatrix}
を計算しています。


順伝播を使って、具体的な値で計算すると以下のようになります。

W1 = np.array([[0.1, 0.3], [0.2, 0.4]])
W2 = np.array([[0.1, 0.3], [0.2, 0.4]])

x = np.array([[1, 0.5]])
y, z1 = forward(x)

print z1
print y

結果は以下の通りになります。

>>> print z1
[[ 0.549834    0.62245933]]
>>> print y
[[ 0.44165237  0.55834763]]

誤差逆伝播

次に誤差逆伝播を行い、ネットワークのパラメータを教師データに近づくように調整します。

教師データと順伝播の出力結果の差分を誤差関数(損失関数)で計ります。
誤差関数(損失関数)が0に近づくように調整します。

誤差関数(損失関数)は、誤差逆伝播を行うには微分が行いやすいことが条件となりますので、交差エントロピーを使います。

softmaxの交差エントロピーは、微分が引き算で計算できるという特性があります。

2層目のネットワークのパラメータの勾配

誤差関数Eをネットワークのパラメータで微分する式は、微分の連鎖則を使うと、以下の通りとなります。

\displaystyle \frac{\partial E}{\partial w^2_{ji}} = \frac{\partial E}{\partial u^2_{j}} \frac{\partial u^2_{j}}{\partial w^2_{ji}} \tag{1}


式(1)の右辺の第1項は、softmaxの交差エントロピーは以下の式なので、

\displaystyle E = -\sum_{k} d_k \log(y_k)

これを出力層の入力で微分すると、

\displaystyle \frac{\partial E}{\partial u^2_{k}} = y_k - d_k

となります。(途中の計算式は省略。詳しくは専門書を参照するなどしてください。)

式(1)の右辺の第2項は、
\displaystyle u^2_j = \sum_i w^2_{ji} z^1_i
なので、
\displaystyle \frac{\partial u^2_{ji}}{\partial w^2_{ji}} = z^1_i
となります。(右辺でijが左辺と一致しない項は微分すると0になる。)

したがって、式(1)は、
\displaystyle \frac{\partial E}{\partial w^2_{ji}} = (y_j - d_j) z^1_i \tag{2}
となります。

行列で書き下すと以下の通りになります。

\displaystyle \begin{pmatrix} \frac{\partial E}{\partial w^2_{11}} & \frac{\partial E}{\partial w^2_{12}} \\ \frac{\partial E}{\partial w^2_{21}} & \frac{\partial E}{\partial w^2_{22}} \end{pmatrix} = \begin{pmatrix} (y_1 - d_1)z^1_1 & (y_1 - d_1)z^1_2 \\ (y_2 - d_2)z^1_1 & (y_2 - d_2)z^1_2 \end{pmatrix}

これで、2層目のネットワークのパラメータの勾配が計算できたことになります。

1層目のネットワークのパラメータの勾配

続いて、1層目のネットワークのパラメータの勾配を計算します。

先ほどの微分の連鎖則を1層目のネットワークのパラメータにも適用すると、

\displaystyle \frac{\partial E}{\partial w^1_{ji}} = \frac{\partial E}{\partial u^1_{j}} \frac{\partial u^1_{j}}{\partial w^1_{ji}} \tag{3}

となります。

右辺の1項目は、u^1がEに与える影響は、u^2の変化で表せることから、微分の連鎖則を使うと、

\displaystyle \frac{\partial E}{\partial u^1_{j}} = \sum_k \frac{\partial E}{\partial u^2_k} \frac{\partial u^2_k}{\partial u^1_j}
で表せます。

右辺の2項目は、\displaystyle u^2_k = \sum_k w^2_{kj} z^1_j = \sum_k w^2_{kj} f(u^1_j)であることから、
\displaystyle \frac{\partial u^2_k}{\partial u^1_j} = w^2_{kj} f'(u^1_j)
と表せます。(活性化関数をfで表している。)

ここで、
\displaystyle \frac{\partial E}{\partial u^l_k} = \delta^l_j
と置き換えると、

\displaystyle \frac{\partial E}{\partial u^1_j} = \delta^1_j = \sum_k \delta^2_k (w^2_{kj} f'(u^1_j)) \tag{4}
と表せます。

式(4)を一般化すると、

\displaystyle \delta^l_j = \sum_k \delta^{l+1}_k (w^{l+1}_{kj} f'(u^l_j)) \tag{5}

となります。


式(3)の右辺の第2項は、\displaystyle u^1_j = \sum_i w^1_{ji} x_iであることから、
\displaystyle \frac{\delta u^1_j}{\delta w^1_{ji}} = x_i
となります。(右辺のjiが左辺のjiと一致しない項は微分すると0になる。)

したがって、式(3)は、

\displaystyle \frac{\partial E}{\partial w^1_{ji}} = \delta^1_j x_i \tag{6}

となり、式(4)から\delta^1_jは、
\displaystyle \delta^2_k = \frac{\partial E}{\partial u^2_k} = y_k - d_k \tag{7}
を使って、

\displaystyle \delta^1_j = \sum_k (y_k - d_k) (w^2_{kj} f'(u^1_j)) \tag{8}

となります。

隠れ層の活性化関数はsigmoidで、

\displaystyle sigmoid({\bf u^1}) = {\bf z^1} \tag{9}

式(9)の微分は、

\displaystyle sigmoid'({\bf u^1}) = {\bf z^1} \circ (1 - {\bf z^1}) \tag{10}

となります。

式(8)のf'を式(10)で置き換えると、式(8)は、

\displaystyle \delta^1_j = \sum_k \delta^2_k (w^2_{kj} z^1_j (1 - z^1_j)) \tag{11}

となります。


したがって、式(6)は、

\displaystyle \frac{\partial E}{\partial w^1_{ji}} = \sum_k (y_k - d_k) (w^2_{kj} z^1_j (1 - z^1_j)) x_i \tag{12}

となり、1層目のネットワークのパラメータの勾配が計算できます。


式(5)のように、前の層の\delta^{l+1}_kから、勾配が計算できることから、誤差逆伝播と呼ばれています。

以上の通り計算した勾配に定数をかけた値で、ネットワークのパラメータを更新します。


数式で書くとややこしいですが、イメージとしては、教師データとの差分が大きければ、影響の大きかった前の層のパラメータをより多く調整するということを行っています。


Pythonでは以下の通り記述します。

def back_propagation(x, z1, y, d):
    global W1
    global W2
    delta2 = y - d
    grad_W2 = z1.T.dot(delta2)

    sigmoid_dash = z1 * (1 - z1)
    delta1 = delta2.dot(W2.T) * sigmoid_dash
    grad_W1 = x.T.dot(delta1)

    W2 -= learning_rate * grad_W2
    W1 -= learning_rate * grad_W1
各行の解説
    delta2 = y - d

は、上記の式(6)を計算しています。

    grad_W2 = z1.T.dot(delta2)

は、上記の式(2)を計算しています。
z1、delta2はともに横ベクトルになっているので、z1を転置して、
\displaystyle \begin{pmatrix} \frac{\partial E}{\partial w^2_{11}} & \frac{\partial E}{\partial w^2_{21}} \\ \frac{\partial E}{\partial w^2_{12}} & \frac{\partial E}{\partial w^2_{22}} \end{pmatrix} = \begin{pmatrix} z^1_1 \\ z^1_2 \end{pmatrix} \cdot \begin{pmatrix} \delta^2_1 & \delta^2_2 \end{pmatrix} = \begin{pmatrix} z^1_1 \delta^2_1 & z^1_1 \delta^2_2 \\ z^1_2 \delta^2_1 & z^1_2 \delta^2_2 \end{pmatrix}
を計算しています。
これで、2層目のネットワークのパラメータの勾配grad_W2が計算できました。

    sigmoid_dash = z1 * (1 - z1)

は、式(10)を計算しています。

    delta1 = delta2.dot(W2.T) * sigmoid_dash

は、式(11)を計算しています。
行列で書き下すと、
\displaystyle \begin{pmatrix} \delta^1_1 & \delta^1_2 \end{pmatrix} = \begin{pmatrix} \delta^2_1 & \delta^2_2 \end{pmatrix} \cdot \begin{pmatrix} w^2_{11} & w^2_{12} \\ w^2_{21} & w^2_{22} \end{pmatrix} * \begin{pmatrix} sigmoid'_1 & sigmoid'_2 \end{pmatrix} \\ = \begin{pmatrix} \delta^2_1 w^2_{11} + \delta^2_2 w^2_{21} & \delta^2_1 w^2_{12} + \delta^2_2 w^2_{22} \end{pmatrix} * \begin{pmatrix} sigmoid'_1 & sigmoid'_2 \end{pmatrix} \\ = \begin{pmatrix} (\delta^2_1 w^2_{11} + \delta^2_2 w^2_{21}) sigmoid'_1 & (\delta^2_1 w^2_{12} + \delta^2_2 w^2_{22}) sigmoid'_2 \end{pmatrix}
となります。(*は要素同士の積を計算します。)

    grad_W1 = x.T.dot(delta1)

は、式(6)を計算しています。
xは横ベクトルなので転置を行っています。行列で書き下すと、
\displaystyle \begin{pmatrix} \frac{\partial E}{\partial w^1_{11}} & \frac{\partial E}{\partial w^1_{21}} \\ \frac{\partial E}{\partial w^1_{12}} & \frac{\partial E}{\partial w^1_{22}} \end{pmatrix} = \begin{pmatrix} x_1 \\ x_2 \end{pmatrix} \cdot \begin{pmatrix} \delta^1_1 & \delta^1_2 \end{pmatrix} = \begin{pmatrix} x_1 \delta^1_1 & x_1 \delta^1_2 \\ x_2 \delta^1_1 & x_2 \delta^1_2 \end{pmatrix}
となります。
これで、1層目のネットワークのパラメータの勾配grad_W1が計算できました。

最後に、求めた勾配を使って、ネットワークのパラメータW2、W1を更新します。

    W2 -= learning_rate * grad_W2
    W1 -= learning_rate * grad_W1

learning_rateは、学習係数で通常は非常に小さい値(0.00001など)を使用します。

行列の演算のみで誤差逆伝播が計算できたことに注目してください。
これは演算にGPUを利用する際に重要になります。


誤差逆伝播を使って、具体的な値で計算すると以下のようになります。

W1 = np.array([[0.1, 0.3], [0.2, 0.4]])
W2 = np.array([[0.1, 0.3], [0.2, 0.4]])
learning_rate = 0.005

# 順伝播
x = np.array([[1, 0.5]])
y, z1 = forward(x)

# 誤差逆伝播
d = np.array([[1, 0]]) # 教師データ
back_propagation(x, z1, y, d)

print W1
print W2

結果は以下の通りになります。

>>> print W1
[[ 0.0998618   0.29986879]
 [ 0.1999309   0.39993439]]
>>> print W2
[[ 0.10153499  0.29846501]
 [ 0.20173774  0.39826226]]


以上で、説明は終わりです。

本来は、ネットワークのパラメータにバイアス項を加えますが、説明を簡単にするため省略しています。

数式の説明は、この本を参考にしました。

Pythonのコードは、この本を参考にしました。