ニューラルネットワークを実装する [Part 3 隠れ層とバックプロパゲーション(誤差逆伝播法)編]

はじめに

前回の記事の続きです。Peters's NoteのPart 3を追っていきます。
非常にシンプルなモデルを使って、隠れ層が果たす役割や、バックプロパゲーション(誤差逆伝播法)の考え方を見つつ、実装していきます。

はじめに
目次
今回扱うニューラルネットワークのモデル
使用するデータの用意
非線形活性化関数
バックプロパゲーションで重みの最適化
- バックプロパゲーションの概要
  - 順伝播
  - 逆伝播
実際にバックプロパゲーションで重みを更新してみる
学習したモデルの可視化
GitHub
関連記事

今回扱うニューラルネットワークのモデル

今回は入力が1次元の2クラス分類問題を扱います。
ニューラルネットワークのモデルは、次の図のような隠れ層を1つだけ持ったシンプルなもの扱います。

f:id:kanohk:20171022174344p:plain — 出典: https://peterroelants.github.io/posts/neural_network_implementation_part03/#Hidden-layer

使用するデータの用意

今までと同様に、まずは使用するデータを用意します。入力 $x$ は1次元で、 $x$ はクラス $t \in \{0, 1\}$ に属します。 $t = 1$ だと青クラスで、 $t = 0$ だと赤クラスです。
次の画像のように、0の周辺には青クラスが、0から離れたところには赤クラスが存在しているようなデータを用意してあげます。
f:id:kanohk:20171022174948p:plain

このデータは前回と似たような形で、次のように用意できます。

# 入力サンプルを用意する
num_of_samples_per_class = 20
blue_mean = [0]
red_left_mean = [-2]
red_right_mean = [2]
std_dev = 0.5

x_blue = np.random.randn(num_of_samples_per_class, 1) * std_dev + blue_mean
x_red_left = np.random.randn(num_of_samples_per_class // 2, 1) * std_dev + red_left_mean
x_red_right = np.random.randn(num_of_samples_per_class // 2, 1) * std_dev + red_right_mean

# xに入力のベクトル、tにクラスのベクトルが入る
x = np.vstack((x_blue, x_red_left, x_red_right))
t = np.vstack((np.ones((x_blue.shape[0], 1)),
               np.zeros((x_red_left.shape[0], 1)),
               np.zeros((x_red_right.shape[0], 1))))

非線形活性化関数

隠れ層の活性化関数には非線形関数を使用します。通常シグモイドなどの関数が使われるようですが、今回はガウシアンRBFという関数を使用します。
$\displaystyle \phi(z) = e^{-z^2}$

この関数を使用する理由は、関数をプロットすると見えてきます。ガウシアンRBFは0の周辺で値が1に近くなり、0から離れると値が0に小さくなります。この関数を活性化関数と使用することで、0の周辺に集まっている青クラスと、0から離れている赤クラスをうまく分離できることを期待しているのです。
f:id:kanohk:20171022175408p:plain

というわけで、ガウシアンRBFを定義しておきます。

# RBFの定義
def rbf(z):
    return np.exp(-z**2)

後ほどの計算で使用しますが、RBFの微分は簡単に計算でき、
$\displaystyle \frac{d \phi(z)}{dz} = -2ze^{-z^2} = -2 z \phi(z)$
となります。

バックプロパゲーションで重みの最適化

バックプロパゲーション(誤差逆伝播法)を使用して、ネットワークの重みの更新をしていきます。
個人的にはこのバックプロパゲーションは非常にややこしく、何度聞いてもなかなか直感的に理解できない場所です。
Qiitaなどにも解説記事はたくさんあるので、私の説明でいまいち理解できない場合は、他の記事もいくつか目を通してみると良いかもしれません。
誤差逆伝播法のノート
 誤差逆伝播法をはじめからていねいに

バックプロパゲーションの概要

バックプロパゲーションには大きく分けて2つのステップがあります。

[順伝播] モデルへの入力をもとに、各層の出力を順方向に計算していき、最終的に出力層の値を計算する
[逆伝播] 出力層でのコスト(=コスト関数の値)をもとに、出力層から入力層の向きに(逆方向に)各層の重みを更新していく

順伝播

順伝播では、入力 $x$ をもとに、隠れ層の値、出力層の値を順に計算していきます。
隠れ層の値 $h$ は、先に定義したガウシアンRBFを用いて次のように書くことができます。
$\displaystyle h = \phi (x w_h) = e^{-(xw_h)^2}$

これをコードに落とすとこうなります。

# 隠れ層の値
def hidden_activations(x, wh):
    return rbf(x * wh)

この隠れ層の値をもとに、出力層の値 $y$ を計算します。出力層の活性化関数には、前回と同じロジスティック関数 $\sigma(z) = \frac{1}{1 + e^{-z}}$ を使用します。すると、
$\displaystyle y = \sigma (h w_o - 1) = \frac{1}{1 + e^{-h w_o - 1}}$
と計算できます。

ここで、隠れ層から出力層へ進むときに $h w_o - 1$ とバイアス -1 を加えている点に注意してください。本来はこのバイアスも学習すべきものなのですが、このチュートリアルでは問題を簡単にするために-1という定数をバイアスとして使用しているようです。

以上の数式をコードに落とすと次のようになります。

# ロジスティック関数
def logistic(z):
    return 1 / (1 + np.exp(-z))

# 出力層の値
def output_activations(h, wo):
    return logistic(h * wo - 1)

以上の関数をもとに、今回のニューラルネットワーク全体と、ニューラルネットワークを使って入力からクラスを予測する関数を定義します。

# モデルを定義
def nn(x, wh, wo):
    return output_activations(hidden_activations(x, wh), wo)

# モデルを使用してクラスを予測する関数
def nn_predict(x, wh, wo):
    return np.around(nn(x, wh, wo))

逆伝播

逆伝播のステップで、最急降下法を使用して重みを出力層側から順に更新します。
重みの更新方法は以前と変わらず $w(k+1) = w(k) - \Delta w(k)$ を使用します。 $\Delta w(k) = \mu \cdot \frac{\partial \xi}{\partial w}$ でした。

コスト関数にはクロスエントロピーを使用します。 $i$ 番目のサンプルのコストは次のようになります。
$\displaystyle \xi(t_i, y_i) = - \{ t_i \log(y_i) + (1 - t_i) \log (1 - y_i) \}$

コード上でもコスト関数を定義しておきます。

# コスト関数(クロスエントロピー)
def cost(y, t):
    return -np.sum(np.multiply(t, np.log(y)) + np.multiply((1 - t), np.log(1 - y)))

def cost_for_param(w, wh, wo, t):
    return cost(nn(w, wh, wo), t)

このコスト関数を3次元空間上にプロットすると次のようになります。
f:id:kanohk:20171022185059p:plain

いよいよ出力層の重み $w_o$ を更新します。 $w_o$ の更新には $\frac{\partial \xi_i}{\partial w_o}$ が必要ですが、前回この値はすでに使用していて、
$\displaystyle \frac{\partial \xi}{\partial w} = x (y - t)$
でした。出力層への入力は $h_i$ なので、 $\frac{\partial \xi_i}{\partial z_{oi}} = \delta_{oi}$ としたとき、
$\displaystyle \frac{\partial \xi_i}{\partial w_o} = h_i \frac{\partial \xi_i}{\partial z_{oi}} = h_i \delta_{oi}$
となります。

$\delta_{oi}$ をgradient_output, $\frac{\partial \xi_i}{\partial w_o}$ をgradient_weight_outputとしてコードで定義すると次のようになります。

# エラーを出力層の入力で微分したもの
def gradient_output(y, t):
    return y - t

# エラーを出力層の重みで微分したもの
def gradient_weight_out(h, grad_output):
    return h * grad_output

これで $w_o$ は更新できるのですが、問題は $w_h$ です。 $w_h$ の更新には $\frac{\partial \xi_i}{\partial w_h}$ が必要です。

連鎖律より、
$\displaystyle \frac{\partial \xi_i}{\partial w_h} = \frac{\partial z_{hi}}{\partial w_h} \frac{\partial \xi_i}{\partial z_{hi}}$
となります。
$\frac{\partial z_{hi}}{\partial w_h}$ は $z_{hi} = x w_h$ (隠れ層への入力)であることから $\frac{\partial z_{hi}}{\partial w_h} = x_i$ とすぐにわかります。

$\frac{\partial \xi_i}{\partial z_{hi}}$ は先ほどの $\delta$ を使うと $\delta_{hi}$ と書けます。つまり、
$\displaystyle \frac{\partial \xi_i}{\partial w_h} = x_i \delta_{hi}$
となります。

ここで $\delta_{hi}$ を詳しく見てみると、
$\displaystyle \delta_{hi} = \frac{\partial \xi_i}{\partial z_{hi}} = \frac{\partial h_i}{\partial z_{hi}} \frac{\partial z_{oi}}{\partial h_i} \frac{\partial \xi_i}{\partial z_{oi}}$
と連鎖律を使って変形できます。3つの微分をそれぞれ見てあげると、1つ目はRBFの微分を思い出すと
$\displaystyle \frac{\partial h_i}{\partial z_{hi}} = -2 z_{hi} h_i$
と書けます。2つ目は $z_{oi} = h_i w_o - 1$ なので
$\displaystyle \frac{\partial z_{oi}}{\partial h_i} = w_o$
となります。3つ目は先ほど $\delta_{oi}$ と置いたので、
$\displaystyle \frac{\partial \xi_i}{\partial z_{oi}} = \delta_{oi}$
となります。

つまり、
$\displaystyle \delta_{hi} = -2 z_{hi} h_i w_o \delta_{oi}$
となり、出力層のエラー $\delta_{oi}$ を使用して、その一つ手前の隠れ層のエラー $\delta_{hi}$ が求められることがわかります。

このように逆向きにエラーを計算していき重みを更新していく手順が逆伝播になります。

$\delta_h$ をgradient_hidden、 $\frac{\partial \xi}{\partial w_h}$ をgradient_weight_hiddenとしてPythonに落としたコードが次になります。

# エラーを隠れ層への入力で微分したもの
def gradient_hidden(wo, grad_output):
    return wo * grad_output

# エラーを隠れ層の重みで微分したもの
def gradient_weight_hidden(x, zh, h, grad_hidden):
    return x * -2 * zh * h * grad_hidden

実際にバックプロパゲーションで重みを更新してみる

では既に定義した関数を使って実際にバックプロパゲーションを実行してみます。
まずはバックプロパゲーションのイテレーション一回分の更新を行う関数を定義します。
必要な微分した関数などはすでにすべて定義してあるので、基本的にはそれらを呼び出して $w(k + 1) = w(k) - \Delta w(k)$ と重みを更新してあげるだけです。

# バックプロパゲーションでの重みの更新
def backprop_update(x, t, wh, wo, learning_rate):
    # まずは順伝播
    zh = x * wh # 隠れ層への入力
    h = rbf(zh) # 隠れ層の出力 (=出力層への入力)
    y = output_activations(h, wo)  # 出力層の出力
    
    # 次に逆伝播
    # 出力層の重みwoの更新に必要なΔwoを求める
    grad_output = gradient_output(y, t)
    d_wo = learning_rate * gradient_weight_out(h, grad_output)
    
    # 同様に、隠れ層の重みwhの更新に必要なΔwhを求める
    grad_hidden = gradient_hidden(wo, grad_output)
    d_wh = leraning_rate * gradient_weight_hidden(x, zh, h, grad_hidden)

    # 重みを両方とも更新する
    return (wh - d_wh.sum(), wo - d_wo.sum())

最後に、このbackprop_update関数を繰り返し呼び出して重みを更新していきます。
学習率は最初は大きめに、イテレーションが進むにつれて小さくなるようにしてあげています。
元記事ではイテレーションを50回走らせていましたが、手元の環境だとうまく学習できなかったのでイテレーションを30回にしてあります。

# バックプロパゲーションを実行する

# 初期重み
wh = 2
wo = -5

# 学習率
learning_rate = 0.2

# 学習率は初期値から徐々に小さくしていく
num_of_iterations = 30
lr_update = learning_rate / num_of_iterations

# 各イテレーションでの重みを保存しておく
w_cost_iter = [(wh, wo, cost_for_param(x, wh, wo, t))] 

# イテレーションをまわす
for i in range(num_of_iterations):
    learning_rate -= lr_update
    wh, wo = backprop_update(x, t, wh, wo, learning_rate)
    w_cost_iter.append((wh, wo, cost_for_param(x, wh, wo, t)))
    
print('最終的なコストは{:.2f}で、その時の重みは wh: {:.2f} and wo: {:.2f} です。'.format(cost_for_param(x, wh, wo, t), wh, wo))
# 最終的なコストは13.53で、その時の重みは wh: 0.76 and wo: 2.92 です。

この時の最適化の過程を可視化すると、このようになります。初期地点から色の濃い部分に徐々に近づいて行っている様子が確認できます。
f:id:kanohk:20171022193314p:plain

学習したモデルの可視化

前回と同様、学習したモデルを可視化してみます。学習したモデルで青クラスと分類される値は青色で、赤クラスと分類される値は赤色で塗りつぶした図が次になります。
f:id:kanohk:20171022193432p:plain

確かに、モデルがうまく入力を分類できていることが確認できます。なぜこのようにうまく分類できるかは、途中の隠れ層の出力を確認してみると理解できます。
用意したデータと学習した重みを使用してhidden_activations関数を読んであげると隠れ層の出力が取得でき、それを次のコードで可視化します。

# 隠れ層の出力を可視化する
plt.figure(figsize=(8,0.5))
plt.xlim(-0.01,1)
plt.ylim(-1,1)

# 入力にhidden_activation関数を適用して隠れ層の出力を計算して、プロット
plt.plot(hidden_activations(x_blue, wh), np.zeros_like(x_blue), 'b|', ms = 30) 
plt.plot(hidden_activations(x_red_left, wh), np.zeros_like(x_red_left), 'r|', ms = 30) 
plt.plot(hidden_activations(x_red_right, wh), np.zeros_like(x_red_right), 'r|', ms = 30) 

plt.gca().axes.get_yaxis().set_visible(False)
plt.title('Projection of the input samples by the hidden layer.')
plt.xlabel('h')
plt.show()

するとこのような画像が描画できます。赤クラスと青クラスが0.4あたりを境目にきれいに左右に分かれているようすが確認できます。
隠れ層の出力がこのような形になっているので、出力層ではこの左右をうまく分けるような重みを学習すればよいわけです。

f:id:kanohk:20171022193701p:plain

GitHub

いつも通りJupyter Notebookはこちらにおいておきます。
notebooks/Part3.ipynb at master · kanoh-k/notebooks · GitHub

次回記事: ニューラルネットワークを実装する [Part 4 バックプロパゲーションのベクトル化] - 今日も窓辺でプログラム

これまでの記事