①線形分類問題（２値分類） - 社畜エンジニア発掘戦線

第一問

設問１．入力x1、x2、正解ラベルtが (x1, x2, t) = (2, 2, 0)、(-2, -2, 1) で、それぞれにガウシアンノイズが付加されたデータセットを作成せよ
設問２．入力(x1, x2)に対する重みを(w1, w2)としてシグモイド関数を活性化関数として出力yを求めよ
設問３．クロスエントロピーを損失関数として損失量を求めよ
設問４．損失関数のW=(w1, w2)による微分式を求めよ
設問６．勾配降下法を用いてパラメータを更新し、イテレーションを30回繰り返した後の出力yを求めよ
設問７．損失量の推移をグラフ化せよ
まとめ

ディープラーニングでトライできる問題として、回帰問題に続き、「分類問題」にチャレンジしていきたいと思います。2次平面上にプロットされ、「0と1」にラベリングされている集合データをディープラーニングで分類してみます。

今回もPYTHONISTA3でコードを書きながら解いていってみましょう。

設問１．入力x1、x2、正解ラベルtが (x1, x2, t) = (2, 2, 0)、(-2, -2, 1) で、それぞれにガウシアンノイズが付加されたデータセットを作成せよ

さて、今回は分類問題にチャレンジです。設問を見ると、入力がx1とx2と増えています。これは２次元平面上にx座標、y座標としてプロットできるということですね。(2, 2)の座標にあたるデータの正解ラベルが「0」、(-2, -2)の方は「1」ということですね。イメージはこんな感じ。
f:id:sutokun:20190204093515j:plain

ではコードで書いていってみます。いつも通りnumpyとmatplotlibをインポートするんですが、今回は2次平面なのでmatplotlibのカラーマップも合わせてインポートします。

import numpy as np
import matplotlib.pyplot as plt
import matplotlib.cm as cm

次にサンプル数（num_of_sam）とガウシアンノイズの広がり具合（std_dv）を定義しておきます。データ数は40、広がり係数は1.8としましょう。

num_of_sam = 40
std_dv = 1.8

では実際に座標を設定します。

group1 = np.array([2,2])
group2 = np.array([-2,-2])

このままではただの点なので、ガウシアンノイズをのせたサンプル数を足します。

group1 = np.array([2,2])+np.random.randn(num_of_sam, 2)*std_dv
group2 = np.array([-2,-2])+np.random.randn(num_of_sam, 2)*std_dv

さらに、group1、group2をひとつのデータセットになるように結合します。このXの1列目がx1、2列目がx2にあたります。

X = np.vstack((group1, group2))

グラフ化してみます、ただノイズの広がりが1.8だと大きすぎて分からないのでグラフは「std_dv=0.5」で描いています。「X[:,0]」と「X[:,1]」がそれぞれx1、x2、グラフのx座標とy座標です。

plt.plot(X[:,0], X[:,1], 'o')
plt.show()

f:id:sutokun:20190204094518p:plain
うまくできています。しかし、今の状態では正解ラベルがないので全て同じ点（同じ色）の集合になっています。それぞれのグループに「0」と「1」を設定します。

t_group1 = np.zeros((num_of_sum, 1))
t_group1 = np.ones((num_of_sum, 1))

Xと同じように結合します。

t = np.vstack((t_group1, t_group2))

データセットの完成です。カラーマップでグラフ化してみましょう。ノイズの広がりは「std_dv=1.8」に戻します。「vmax、vmin」は色を表す値の最大、最小です。正解ラベルの0と1に対応しています。「cmap=cm.bwr」は色の種類を示しています。今回は0が青、1が赤の設定です。

plt.scatter(X[:,0], X[:,1], vmin=0, vmax=1, c=t[:,0], cmap=cm.bwr, marker='o', s=50)
plt.show()

f:id:sutokun:20190206113643p:plain
いい感じです。この問題は「線形分類」なので、直線でこの2グループを分類します。とにかく進めていきましょう。

設問２．入力(x1, x2)に対する重みを(w1, w2)としてシグモイド関数を活性化関数として出力yを求めよ

ディープラーニングの問題解決フローは同じなので、ご察しの通り次はニューラルネットワークの設計です。入力はx1とx2、出力は0か1かをラベリングするyひとつなので、こんな感じ。
f:id:sutokun:20190203061819j:plain

出力には非線形回帰で用いたシグモイド関数を使います。ここで使う関数を活性化関数と呼びます。正解ラベルは0か1なので、この間で推移するシグモイド関数が効果的です。今回は入力が2つになっていますが、非線形回帰問題とかに比べると結構シンプルですね。

まずはシグモイド関数が使えるようにPythonの関数として定義しておきます。

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

次に重みの初期値を設定します。w1とw2をまとめて行列Wで書きましょう。

W = np.random.randn(2,1)

準備が整ったので出力を求めましょう。XとWは線形結合（内積計算）でシグモイド関数を通します。

y = sigmoid(np.dot(X, W))

念のため、print(Y.shape)で出力の型を見ておきましょうか。正解ラベルtと型が一致しているはずです。

print(t) #(80, 1)
print(y) #(80, 1)

これで設問2は終了です。

設問３．クロスエントロピーを損失関数として損失量を求めよ

流れ通り、次は損失関数の設定ですね。回帰問題でスタンダードだった二乗和誤差と異なり、今回はクロスエントロピーを用います。これも場面によっていくつか定義式があるんですが、今回は次式の定義を用いて進めていきたいと思います。

クロスエントロピー： ${\displaystyle E = -\sum_{i} \Bigl( t_{i} \log y_{i} +(1-t_{i}) \log (1-y_{i}) \Bigr) }$

前回の二乗和誤差関数よりも複雑なので、こちらもシグモイド同様、Pythonの関数で書いておきましょう。

def loss(y, t):
    return -np.sum(np.multiply(t, np.log(y)) + np.multiply((1 - t), np.log(1 - y)))

では設問2で求めたYとデータセットの正解ラベルtからクロスエントロピーを用いて損失量を求めてみます。

E = loss(y, t)

print(E)してみるとちゃんとスカラー値が出てきます、これで設問2は終了です。

設問４．損失関数のW=(w1, w2)による微分式を求めよ

ここの流れももうおなじみです、連鎖率を用いて記述してみます。

${\displaystyle \frac{\partial E}{\partial W}\ = \frac{\partial y }{\partial W}\ \frac{\partial E }{\partial y}\ }$

ではクロスエントロピーのy微分から計算してみます。

${\displaystyle \frac{\partial E }{\partial y}\ = \frac{\partial }{\partial y}\ \Bigl(-\sum_{i} \Bigl( t_{i} \log y_{i} +(1-t_{i}) \log (1-y_{i}) \Bigr) \Bigr) = \sum_{i} \Bigl( -\frac{t_{i}}{y_{i}} + \frac{1 - t_{i}}{1 - y_{i}} \Bigr) }$

次にyのW微分です。yは入力XとパラメータWの内積をシグモイド関数に通したものです。シグモイド関数の微分を忘れないように。

${\displaystyle \frac{\partial y_{i} }{\partial W}\ = \frac{\partial }{\partial W}\ \Bigl( \sigma (X_{i} \cdot W) \Bigr) = \sum_{i} \Bigl( X_{i}y_{i}(1-y_{i}) \Bigr) }$

以上より、

${\displaystyle \frac{\partial E}{\partial W}\ = \sum_{i} \Bigl( X_{i}y_{i}(1-y_{i}) \Bigr) \sum_{i} \Bigl( -\frac{t_{i}}{y_{i}} + \frac{1 - t_{i}}{1 - y_{i}} \Bigr) = \sum_{i} \Bigl( X_{i}(y_{i}-t_{i}) \Bigr) }$

気持ちいい鬼約分ですね、かなりシンプルになりました。

ではこの内容をコードに落とします。シグマの方向に注意（axis=0）。

dW = np.sum(X*(y-t),axis=0)

以上で設問4は終了です。

設問６．勾配降下法を用いてパラメータを更新し、イテレーションを30回繰り返した後の出力yを求めよ

ここは今までどおり勾配降下法なので、コードに落とすだけです。とりあえず学習率は0.002としましょう。いつも通り損失量を保存するコードを追記しておきます。

勾配降下法： ${\displaystyle P_{new} = P_{old} - \mu \frac{\partial E}{\partial P}\ }$

learning_rate = 0.002
E_save = []

num_of_itr = 1
for i in range(num_of_itr):
    y = sigmoid(np.dot(X, W))
    E = loss(y, t)
    E_save = np.append(E_save, E)

    dW = np.sum(X*(y-t),axis=0)
    W = W - learning_rate*np.reshape(dW,(2,1))

dWは要素が2つの行列ですが、型が(1,2)に反転しているので、reshapeでもとのWと同じ型にそろえています。これでパラメータWが更新され、出力を求める準備が整いました。

さて、これで求められる出力をどのように表示するかが問題です。普通に数字だけ出しても「それで？」って感じです。実際に2次元グラフ上に描画してみましょう。色々なやり方があると思いますが、いちばん直感的かな〜と思うのが、2次元グラフの格子点を使って出力を表示してみる方法です。入力データセットXの(x1, x2)をそれぞれ格子点として用意し、各格子点に対して更新された重みWを使って出力yを求めます。
f:id:sutokun:20190203055111j:plain

ちょっとなんのこっちゃという感じですが、とりあえず進めます。

まずは格子点を用意します。どの範囲で表示するかをgrid_range、その範囲に何点の格子点を突っ込むかをresolutionで定義します。-5から5までの範囲を60分割するように設定します。

grid_range = 5
resolution = 60
x1_grid = x2_grid = np.linespace(-grid_range, grid_range, resolution)

これで格子点の外側の目盛りを作成できました。次に中の格子点を作っていきます。numpyのmeshgridを使うと楽に作成できます。

xx, yy = np.meshgrid(x1_grid, x2_grid)
X_grid = np.c_[xx.ravel(), yy.ravel()]

これで格子点の作成が完了しました。一度グラフにプロットしてみましょう。

plt.plot(X_grid[:,0], X_grid[:,1], 'o')
plt.show()

f:id:sutokun:20190203061659j:plain
目がチカチカしますが、ちゃんと格子点がグラフを埋めていますね。

ではこの格子点におけるそれぞれの出力を求めてみます。

Y_grid = sigmoid(np.dot(X_grid, W))

出力はsigmoid関数なので、0から1までの値を取ります。正解ラベルは0か1なのでこのどちらかに入るよう、四捨五入しましょう。

Y_predict = np.around(Y_grid)

ではこの出力をカラーマップで描画してみましょう。同じ色調で描くと色がきついので格子点はalphaで透かしています。

plt.scatter(X[:,0], X[:,1], vmin=0, vmax=1, c=T[:,0], cmap=cm.bwr, marker='o', s=50)
plt.scatter(X_grid[:,0], X_grid[:,1], vmin=0, vmax=1, c=Y_predict[:,0], cmap=cm.bwr, marker='o', s=50,alpha=0.3)
plt.show()

f:id:sutokun:20190206113733p:plain
いい具合に青と赤で分類できました。

設問７．損失量の推移をグラフ化せよ

【解決済】
ここで問題発生です。今回の問題では損失関数にクロスエントロピーを使っていますが、logの計算に不具合が生じます（0で割り算しちゃう系）。これはディープラーニングの本質と関係しているのかなんとも分からず、詳細を調べるのはめんどくさ…非効率なので、今回この設問はパスすることにします。回避方法がひらめいたら修正します（すんません）。【追記】少し調べたらイテレーション5回あたりからエラーが発生します。やはり学習が進んでyの値が更新されてくるとlog計算がだんだん難しくなってくる感じですかね。print(y)してみるとすごい綺麗な「1」が出力されていました。クロスエントロピーの式は「log(1-y)」を含むのでマイナス∞に発散したようです。じゃあなんでこんな綺麗な「1」が出力されるんだ？やっぱり問題が簡単すぎるのか…？【追記】嗚呼ぁぁああああ！！！！凡ミスぅうううう！！！パラメータを更新するときに学習率（learning_rate）をかけ忘れていました。つまり、学習率0.002かけるとか言っといて実は学習率1の天才ニューラルネットワークさんだったわけですね。そりゃ綺麗な「1」も出力されるよ。コレは恥ずかしい。このページのコードは修正しますがこの追記は戒めとして残しておきます。

（いろいろトラブルもありましたが）損失量をグラフ化してみます。