ラビットチャレンジ 深層学習Day3 その2

Section2 LSTM

RNNの課題

  • 勾配消失

 時間を遡るほど勾配が消失してしまうため、長い時系列の学習が困難だった。
 勾配消失は、誤差逆伝搬法で深い層のニューラルネットワークを遡る時に見られる現象で、1より小さい値になる微分値が多数乗算されることにより引き起こされる。例えば、活性化関数として使用されるシグモイド関数は、微分値は最大でも0.25にしかならないため、2つ以上乗算されると、非常に小さな値となり、勾配はほぼ消失してしまう。

  • 確認テスト

シグモイド関数微分した時、入力値が0の時に最大値をとる。その値として正しいものを選択肢から選べ。
(1)0.15
(2)0.25
(3)0.35
(4)0.45
[解答]
(2)0.25

  • 勾配爆発

 勾配消失とは逆に、勾配が層を逆伝搬するごとに指数関数的に大きくなってしまう問題。活性化関数に恒関数を用いたりすると、この問題が起こりやすい。

  • 演習チャレンジ

RNNや深いモデルでは勾配の消失または勾配爆発が起こる傾向がある。勾配爆発を防ぐために勾配のクリッピングを行うという手法がある。具体的には勾配のノルムが閾値を超えたら、勾配のノルムを閾値に正規化するというものである。以下は勾配のクリッピングを行う関数である。
(さ)にあてはまるコードをかけ
f:id:tibet:20211129184855p:plain
[解答]
gradient*rate

LSTM(Long Short Term Memory)は、RNNの勾配消失や勾配爆発という問題を解決するために考案された構造である。
次に、LSTMの全体像を示す。
f:id:tibet:20211129195005p:plain
入力から出力の間がRNNの中間層に相当する。中間層で自己回帰(青点線)していることがわかる。
中間層内の構造の工夫がLSTMである。

1-2 CEC

 LSTMの中心はCEC(Constant Error carousel:定誤差カルーセル)である。このセルは、勾配を記憶する役割を果たしている。
 具体的には、CECの中では、勾配消失及び勾配爆発の解決法として、勾配を1としている。
 δ^{t-z-1}=δ^{t-z}(Wf'(u^{t-z-1}))=1
 \dfrac{\partial E}{\partial c^{t-1}}=\dfrac{\partial E}{\partial c^t}\dfrac{\partial}{\partial c^{t-1}}(a^t-c^{t-1})=\dfrac{\partial E}{\partial c^t}

CECの課題は、入力データについて時間依存度の関係なく重みが一律のため、学習が出来ないということである。
 つまり、学習と記憶の役割を中間層内で分けている。
 学習については、入力層の重荷については入力ゲート、出力層への重みについては出力ゲートでそれぞれ担っている。

2-2入力ゲートと出力ゲート

入力ゲート

f:id:tibet:20211129210828p:plain
入力ゲートは、CECに記憶する内容を調整する機能を持つ。
具体的には入力側の変数x(t)と、出力からの一つ前の数値の回帰h(t-1)にそれぞれ重みWiとUiを掛け合わせて演算した値をCECの入力との内積演算ユニットに送る。

出力ゲート

f:id:tibet:20211129211301p:plain
出力ゲートは、ECEからの出力を調整して出力層に送る機能を持つ。
具体的には入力側の変数x(t)と、出力からの一つ前の数値の回帰h(t-1)にそれぞれ重みWoとUoを掛け合わせて演算した値をCECの出力との内積演算ユニットに送る。

このように、入力・出力ゲートを追加することで、それぞれのゲートへの入力値の重みを重み行列W,Uで可変可能とすることで、CECの課題を解決している。

2-3忘却ゲート

CECは過去のすべての情報が保管されているが、逆に不要な場合でも削除が出来ない。
そこで、過去の上方がいらなくなった場合、そのタイミングで情報を忘却するための機能として、忘却ゲートが設けられた。
f:id:tibet:20211130095950p:plain
上記の赤枠が忘却ゲートである。
忘却ゲートは、入出力にそれぞれWf, Ufの重みを乗算した関数忘却関数f(t)を作る。次にCECの一つ前の時間の出力c(t-1)と内積を取り、CEC内で入力層からの信号i(t)a(t)と線形結合をとり、下記c(t)を出力する。
 c(t)=i(t)\cdot a(t)+f(t)\cdot c(t-1)
忘却関数f(t)と内積をとることで、CECの記憶の忘却度が調整される。

  • 確認テスト

以下の文章をLSTMに入力し、空欄に当てはまる単語を予測したいとする。文中の「とても」という言葉は空欄の予測において、無くなっても影響を及ぼさないと考えられる。
このような場合、どのゲートが作用すると考えられるか。
「映画おもしろかったね。ところで、とてもおなかが空いたからなにか____。」

[解答]
文意に影響を及ぼさないため、忘却ゲートが作用していると考えられる。

  • 演習チャレンジ

以下のプログラムはLSTMの順伝搬を行うプログラムである。ただし、_sigmoid関数は、要素ごとにシグモイドを作用させる関数である。
(け)にあてはまるコードを示せ。
f:id:tibet:20211130135130p:plain
[解答]
input_gate*a+forget_gate*c

覗き穴結合

CECの保存されている過去の情報を任意のタイミングで他のノードに伝搬させたり、任意のタイミングで忘却させたりしたいが、CEC自身の値は、ゲート制御に影響を与えていないという課題がある。
そこで、覗き穴結合というパスを作り、CEC自身の値に重み行列を介して伝搬可能にした構造を設けた。
f:id:tibet:20211130135346p:plain
上図中の赤枠で囲んだパスが覗き穴結合に当たり、入力ゲート、忘却ゲート、出力ゲートそれぞれにCECに出力にVi, Vf,Voという重みを乗算して伝達している。

Section 3 GRU

 従来のLSTMでは、パラメータ数が多いため、計算負荷が大きかった。そこで、パラメータを大幅に削って計算負荷を小さくし、かつ精度は同等以上を見込める構造として、GRUが考案された。
f:id:tibet:20211130171202p:plain
GRUはリセットゲート、更新ゲート、活性化関数を持つ内部ノードからなり、CECは持っていない。
数式化のために、下記のように信号と重みの定義をする。
x(t):入力信号
W_h:入力重み
W_r:入力からリセットゲートへの重み
U_r:出力からリセットゲートへの重み
b_h(t):リセットゲートのバイアス
r(t):リセットゲートの出力信号
 U_h:出力信号の回帰への重み
 f(x):内部ノードの活性化関数
W_z:入力から更新ゲートへの重み
U_z:出力から更新ゲートへの重み
b_z(t):更新ゲートのバイアス
z(t):更新ゲートの出力信号
h(t):出力信号
まず、リセットゲートの出力は、下記のように表せる。
 r(t)=Wrx(t)+Ur\cdot h(t-1)+b_h(t)
次段の内部演算ユニットの出力は下記のようになる。
 U_h\cdot (r(t)\cdot h(t-1))
よって、活性化関数f(x)を持つ内部ノードの出力は下記のようになる。
 f(W_hx(t)+U_h\cdot (r(t)\cdot h(t-1))
また、更新ゲートの出力は下記のようになる。
 z(t)=W_zx(t)+U_z\cdot h(t-1)+b_z(t)
よって、出力層の出力h(t)は下記のように表せる。
 h(t)=z(t)\cdot h(t-1)+(1-z(t))\cdot f(W_hx(t)+U_h\cdot (r(t)\cdot h(t-1))

  • 確認テスト

LSTMとCECが抱える問題について、それぞれ簡潔に述べよ
[解答]
LSTMはパラメータが多く、計算負荷が高いため、計算に時間がかかる。CECは勾配1で記憶するが、情報処理が出来ない。そこで、周りに情報処理のためのゲートを設けるため、構造が複雑になる。

実装演習

RNNで単語予測をするネットワークをtensorflowで記述した。

import tensorflow as tf
import numpy as np
import re
import glob
import collections
import random
import pickle
import time
import datetime
import os

# logging levelを変更
tf.logging.set_verbosity(tf.logging.ERROR)

class Corpus:
    def __init__(self):
        self.unknown_word_symbol = "<???>" # 出現回数の少ない単語は未知語として定義しておく
        self.unknown_word_threshold = 3 # 未知語と定義する単語の出現回数の閾値
        self.corpus_file = "./corpus/**/*.txt"
        self.corpus_encoding = "utf-8"
        self.dictionary_filename = "./data_for_predict/word_dict.dic"
        self.chunk_size = 5
        self.load_dict()

        words = []
        for filename in glob.glob(self.corpus_file, recursive=True):
            with open(filename, "r", encoding=self.corpus_encoding) as f:

                # word breaking
                text = f.read()
                # 全ての文字を小文字に統一し、改行をスペースに変換
                text = text.lower().replace("\n", " ")
                # 特定の文字以外の文字を空文字に置換する
                text = re.sub(r"[^a-z '\-]", "", text)
                # 複数のスペースはスペース一文字に変換
                text = re.sub(r"[ ]+", " ", text)

                # 前処理: '-' で始まる単語は無視する
                words = [ word for word in text.split() if not word.startswith("-")]


        self.data_n = len(words) - self.chunk_size
        self.data = self.seq_to_matrix(words)

    def prepare_data(self):
        """
        訓練データとテストデータを準備する。
        data_n = ( text データの総単語数 ) - chunk_size
        input: (data_n, chunk_size, vocabulary_size)
        output: (data_n, vocabulary_size)
        """

        # 入力と出力の次元テンソルを準備
        all_input = np.zeros([self.chunk_size, self.vocabulary_size, self.data_n])
        all_output = np.zeros([self.vocabulary_size, self.data_n])

        # 準備したテンソルに、コーパスの one-hot 表現(self.data) のデータを埋めていく
        # i 番目から ( i + chunk_size - 1 ) 番目までの単語が1組の入力となる
        # このときの出力は ( i + chunk_size ) 番目の単語
        for i in range(self.data_n):
            all_output[:, i] = self.data[:, i + self.chunk_size] # (i + chunk_size) 番目の単語の one-hot ベクトル
            for j in range(self.chunk_size):
                all_input[j, :, i] = self.data[:, i + self.chunk_size - j - 1]

        # 後に使うデータ形式に合わせるために転置を取る
        all_input = all_input.transpose([2, 0, 1])
        all_output = all_output.transpose()

        # 訓練データ:テストデータを 4 : 1 に分割する
        training_num = ( self.data_n * 4 ) // 5
        return all_input[:training_num], all_output[:training_num], all_input[training_num:], all_output[training_num:]


    def build_dict(self):
        # コーパス全体を見て、単語の出現回数をカウントする
        counter = collections.Counter()
        for filename in glob.glob(self.corpus_file, recursive=True):
            with open(filename, "r", encoding=self.corpus_encoding) as f:

                # word breaking
                text = f.read()
                # 全ての文字を小文字に統一し、改行をスペースに変換
                text = text.lower().replace("\n", " ")
                # 特定の文字以外の文字を空文字に置換する
                text = re.sub(r"[^a-z '\-]", "", text)
                # 複数のスペースはスペース一文字に変換
                text = re.sub(r"[ ]+", " ", text)

                # 前処理: '-' で始まる単語は無視する
                words = [word for word in text.split() if not word.startswith("-")]

                counter.update(words)

        # 出現頻度の低い単語を一つの記号にまとめる
        word_id = 0
        dictionary = {}
        for word, count in counter.items():
            if count <= self.unknown_word_threshold:
                continue

            dictionary[word] = word_id
            word_id += 1
        dictionary[self.unknown_word_symbol] = word_id

        print("総単語数:", len(dictionary))

        # 辞書を pickle を使って保存しておく
        with open(self.dictionary_filename, "wb") as f:
            pickle.dump(dictionary, f)
            print("Dictionary is saved to", self.dictionary_filename)

        self.dictionary = dictionary

        print(self.dictionary)

    def load_dict(self):
        with open(self.dictionary_filename, "rb") as f:
            self.dictionary = pickle.load(f)
            self.vocabulary_size = len(self.dictionary)
            self.input_layer_size = len(self.dictionary)
            self.output_layer_size = len(self.dictionary)
            print("総単語数: ", self.input_layer_size)

    def get_word_id(self, word):
        # print(word)
        # print(self.dictionary)
        # print(self.unknown_word_symbol)
        # print(self.dictionary[self.unknown_word_symbol])
        # print(self.dictionary.get(word, self.dictionary[self.unknown_word_symbol]))
        return self.dictionary.get(word, self.dictionary[self.unknown_word_symbol])

    # 入力された単語を one-hot ベクトルにする
    def to_one_hot(self, word):
        index = self.get_word_id(word)
        data = np.zeros(self.vocabulary_size)
        data[index] = 1
        return data

    def seq_to_matrix(self, seq):
        print(seq)
        data = np.array([self.to_one_hot(word) for word in seq]) # (data_n, vocabulary_size)
        return data.transpose() # (vocabulary_size, data_n)

class Language:
    """
    input layer: self.vocabulary_size
    hidden layer: rnn_size = 30
    output layer: self.vocabulary_size
    """

    def __init__(self):
        self.corpus = Corpus()
        self.dictionary = self.corpus.dictionary
        self.vocabulary_size = len(self.dictionary) # 単語数
        self.input_layer_size = self.vocabulary_size # 入力層の数
        self.hidden_layer_size = 30 # 隠れ層の RNN ユニットの数
        self.output_layer_size = self.vocabulary_size # 出力層の数
        self.batch_size = 128 # バッチサイズ
        self.chunk_size = 5 # 展開するシーケンスの数。c_0, c_1, ..., c_(chunk_size - 1) を入力し、c_(chunk_size) 番目の単語の確率が出力される。
        self.learning_rate = 0.005 # 学習率
        self.epochs = 1000 # 学習するエポック数
        self.forget_bias = 1.0 # LSTM における忘却ゲートのバイアス
        self.model_filename = "./data_for_predict/predict_model.ckpt"
        self.unknown_word_symbol = self.corpus.unknown_word_symbol

    def inference(self, input_data, initial_state):
        """
        :param input_data: (batch_size, chunk_size, vocabulary_size) 次元のテンソル
        :param initial_state: (batch_size, hidden_layer_size) 次元の行列
        :return:
        """
        # 重みとバイアスの初期化
        hidden_w = tf.Variable(tf.truncated_normal([self.input_layer_size, self.hidden_layer_size], stddev=0.01))
        hidden_b = tf.Variable(tf.ones([self.hidden_layer_size]))
        output_w = tf.Variable(tf.truncated_normal([self.hidden_layer_size, self.output_layer_size], stddev=0.01))
        output_b = tf.Variable(tf.ones([self.output_layer_size]))

        # BasicLSTMCell, BasicRNNCell は (batch_size, hidden_layer_size) が chunk_size 数ぶんつながったリストを入力とする。
        # 現時点での入力データは (batch_size, chunk_size, input_layer_size) という3次元のテンソルなので
        # tf.transpose や tf.reshape などを駆使してテンソルのサイズを調整する。

        input_data = tf.transpose(input_data, [1, 0, 2]) # 転置。(chunk_size, batch_size, vocabulary_size)
        input_data = tf.reshape(input_data, [-1, self.input_layer_size]) # 変形。(chunk_size * batch_size, input_layer_size)
        input_data = tf.matmul(input_data, hidden_w) + hidden_b # 重みWとバイアスBを適用。 (chunk_size, batch_size, hidden_layer_size)
        input_data = tf.split(input_data, self.chunk_size, 0) # リストに分割。chunk_size * (batch_size, hidden_layer_size)

        # RNN のセルを定義する。RNN Cell の他に LSTM のセルや GRU のセルなどが利用できる。
        cell = tf.nn.rnn_cell.BasicRNNCell(self.hidden_layer_size)
        outputs, states = tf.nn.static_rnn(cell, input_data, initial_state=initial_state)
        
        # 最後に隠れ層から出力層につながる重みとバイアスを処理する
        # 最終的に softmax 関数で処理し、確率として解釈される。
        # softmax 関数はこの関数の外で定義する。
        output = tf.matmul(outputs[-1], output_w) + output_b

        return output

    def loss(self, logits, labels):
        cost = tf.reduce_mean(tf.nn.softmax_cross_entropy_with_logits(logits=logits, labels=labels))

        return cost

    def training(self, cost):
        # 今回は最適化手法として Adam を選択する。
        # ここの AdamOptimizer の部分を変えることで、Adagrad, Adadelta などの他の最適化手法を選択することができる
        optimizer = tf.train.AdamOptimizer(learning_rate=self.learning_rate).minimize(cost)

        return optimizer

    def train(self):
        # 変数などの用意
        input_data = tf.placeholder("float", [None, self.chunk_size, self.input_layer_size])
        actual_labels = tf.placeholder("float", [None, self.output_layer_size])
        initial_state = tf.placeholder("float", [None, self.hidden_layer_size])

        prediction = self.inference(input_data, initial_state)
        cost = self.loss(prediction, actual_labels)
        optimizer = self.training(cost)
        correct = tf.equal(tf.argmax(prediction, 1), tf.argmax(actual_labels, 1))
        accuracy = tf.reduce_mean(tf.cast(correct, tf.float32))

        # TensorBoard で可視化するため、クロスエントロピーをサマリーに追加
        tf.summary.scalar("Cross entropy: ", cost)
        summary = tf.summary.merge_all()

        # 訓練・テストデータの用意
        # corpus = Corpus()
        trX, trY, teX, teY = self.corpus.prepare_data()
        training_num = trX.shape[0]

        # ログを保存するためのディレクトリ
        timestamp = time.time()
        dirname = datetime.datetime.fromtimestamp(timestamp).strftime("%Y%m%d%H%M%S")

        # ここから実際に学習を走らせる
        with tf.Session() as sess:
            sess.run(tf.global_variables_initializer())
            summary_writer = tf.summary.FileWriter("./log/" + dirname, sess.graph)

            # エポックを回す
            for epoch in range(self.epochs):
                step = 0
                epoch_loss = 0
                epoch_acc = 0

                # 訓練データをバッチサイズごとに分けて学習させる (= optimizer を走らせる)
                # エポックごとの損失関数の合計値や(訓練データに対する)精度も計算しておく
                while (step + 1) * self.batch_size < training_num:
                    start_idx = step * self.batch_size
                    end_idx = (step + 1) * self.batch_size

                    batch_xs = trX[start_idx:end_idx, :, :]
                    batch_ys = trY[start_idx:end_idx, :]

                    _, c, a = sess.run([optimizer, cost, accuracy],
                                       feed_dict={input_data: batch_xs,
                                                  actual_labels: batch_ys,
                                                  initial_state: np.zeros([self.batch_size, self.hidden_layer_size])
                                                  }
                                       )
                    epoch_loss += c
                    epoch_acc += a
                    step += 1

                # コンソールに損失関数の値や精度を出力しておく
                print("Epoch", epoch, "completed ouf of", self.epochs, "-- loss:", epoch_loss, " -- accuracy:",
                      epoch_acc / step)

                # Epochが終わるごとにTensorBoard用に値を保存
                summary_str = sess.run(summary, feed_dict={input_data: trX,
                                                           actual_labels: trY,
                                                           initial_state: np.zeros(
                                                               [trX.shape[0],
                                                                self.hidden_layer_size]
                                                           )
                                                           }
                                       )
                summary_writer.add_summary(summary_str, epoch)
                summary_writer.flush()

            # 学習したモデルも保存しておく
            saver = tf.train.Saver()
            saver.save(sess, self.model_filename)

            # 最後にテストデータでの精度を計算して表示する
            a = sess.run(accuracy, feed_dict={input_data: teX, actual_labels: teY,
                                              initial_state: np.zeros([teX.shape[0], self.hidden_layer_size])})
            print("Accuracy on test:", a)


    def predict(self, seq):
        """
        文章を入力したときに次に来る単語を予測する
        :param seq: 予測したい単語の直前の文字列。chunk_size 以上の単語数が必要。
        :return:
        """

        # 最初に復元したい変数をすべて定義してしまいます
        tf.reset_default_graph()
        input_data = tf.placeholder("float", [None, self.chunk_size, self.input_layer_size])
        initial_state = tf.placeholder("float", [None, self.hidden_layer_size])
        prediction = tf.nn.softmax(self.inference(input_data, initial_state))
        predicted_labels = tf.argmax(prediction, 1)

        # 入力データの作成
        # seq を one-hot 表現に変換する。
        words = [word for word in seq.split() if not word.startswith("-")]
        x = np.zeros([1, self.chunk_size, self.input_layer_size])
        for i in range(self.chunk_size):
            word = seq[len(words) - self.chunk_size + i]
            index = self.dictionary.get(word, self.dictionary[self.unknown_word_symbol])
            x[0][i][index] = 1
        feed_dict = {
            input_data: x, # (1, chunk_size, vocabulary_size)
            initial_state: np.zeros([1, self.hidden_layer_size])
        }

        # tf.Session()を用意
        with tf.Session() as sess:
            # 保存したモデルをロードする。ロード前にすべての変数を用意しておく必要がある。
            saver = tf.train.Saver()
            saver.restore(sess, self.model_filename)

            # ロードしたモデルを使って予測結果を計算
            u, v = sess.run([prediction, predicted_labels], feed_dict=feed_dict)

            keys = list(self.dictionary.keys())


            # コンソールに文字ごとの確率を表示
            for i in range(self.vocabulary_size):
                c = self.unknown_word_symbol if i == (self.vocabulary_size - 1) else keys[i]
                print(c, ":", u[0][i])

            print("Prediction:", seq + " " + ("<???>" if v[0] == (self.vocabulary_size - 1) else keys[v[0]]))

        return u[0]


def build_dict():
    cp = Corpus()
    cp.build_dict()

if __name__ == "__main__":
    #build_dict()

    ln = Language()

    # 学習するときに呼び出す
    #ln.train()

    # 保存したモデルを使って単語の予測をする
    ln.predict("some of them looks like")


結果の一部

lipoprotein : 1.3668956e-14
ldl : 1.4424414e-14
statin : 1.378644e-14
a--z : 1.3740787e-14
simvastatin : 1.413377e-14
bmi : 3.1720919e-07
covariates : 2.8344898e-06
yhl : 2.2330354e-07
yol : 9.992968e-11
obesity : 1.3501551e-09
evgfp : 6.1829963e-09
unintended : 4.6785074e-09
sizes : 2.569963e-07
obese : 1.9368042e-07
<???> : 2.9195742e-05
Prediction: some of them looks like et
  • 演習チャレンジ

GROもLSTMと同様にRNNの一種であり、単純な
RNNにおいて問題となる勾配消失問題を解決し、長期的な依存関係を学習することが出来る。LSTMに比べ変数の数やゲートの数が少なく、より単純なモデルであるが、タスクによってはLSTMより良い性能を発揮すr。以下のプログラムはGRUである。ただし、_sigmoid関数は要素ごとにシグモイド関数を作用させる関数である。
(こ)にあてはまるコードを書け。
f:id:tibet:20211130223210p:plain
[解答]
GRUの更新ゲートの出力なので、下記になる。
(1-z)*h+z*h_bar

  • 確認テスト

LSTMとGRUの違いを述べよ
[解答]
LSTMはCEC、入力ゲート、忘却ゲート、出力ゲートを備え、構成が複雑でパラメータが多い。
GRUは、リセットゲート、更新ゲートで更新され、LSTMより構成が単純でパラメータが少ない。

Section4 双方向RNN

 過去の上方だけでなく、未来の情報を加味することで、精度を向上させるモデル。文章の推敲や機械翻訳などに使用される。
f:id:tibet:20211130223758p:plain

  • 演習チャレンジ

以下は双方向RNNの順伝搬を行うプログラムである。順方向については、入力から中間層への重みW_f、1ステップ前の中間層出力から中間層への重みをU_f、逆方向に関しては同様にパラメータW_b、U_bを持ち、両者の中間表現を合わせた特徴から出力層への重みはVである。_rnn関数はRNNの順伝搬を表し中間層の系列を返す関数であるとする。(か)にあてはまるコードを示せ。
f:id:tibet:20211130223912p:plain
[解答]
hsは、順方向と逆方向の中間表現の特徴量を合わせたものなので、下記になる。

np.concatenate([h_f,h_b[::,-1]],axis=1)

ラビットチャレンジレポート 深層学習Day3 その1

Section1 RNN

1-1 RNN全体像

1-1-1 RNNとは

 RNNとは、時系列データに対応可能なニューラルネットワークである。

1-1-2 時系列データ

 時系列データとは、時間的順序を追って一定間隔ごとに観察され、
 しかも相互に統計的依存関係が認められるようなデータの系列

    • 音声データ
    • テキストデータ etc

時系列データは物理時間に従って現れるものだけでなく、テキストデータのようにある順序で系列として処理するデータも含まれる。
 CNNでは、このような系列処理は難しいため、RNNの技術が登場するようになった。これにより、自然言語処理技術が大きく飛躍することになった。

1-1-3 RNNとは

RNNは、以前の中間層の状態も利用して学習する構造になっている。
下図に、RNNの概念図を示す。
f:id:tibet:20211125172959p:plain
上図の左部が、RNNのネットワークの模式図で、入力層x、中間層z、出力層で構成されている。
RNNの特徴は、中間層zが自己回帰していることである。
分かりやすく展開したものが右部になる。
中間層の初期状態がz0とし、
時系列で[x1,x2,x3,x4....]というデータが入力されるとすると、
y1の処理には、中間層状態z1と前の中間状態に重みをかけたW×z0の線形結合を使用し、
y2の処理には、中間層状態z2と前の中間状態に重みをかけたW×z1の線形結合を使用し、
y3の処理には、中間層状態z3と前の中間状態に重みをかけたW×z2の線形結合を使用し、
.....
という形になる。
数式としては、下記のように書ける。
 u^t=W_{(in)}x^t+Wz^{t-1}+b
 z^t=f(W_{(in)}x^t+Wz^{t-1}+b)
 y^t=g(W_{(out)}z^t+c)
f(x):活性化関数
g(x):活性化関数
コードとしては、下記のように表せる。

u[:,t+1]=np.dot(X,W_in)+np.dot(z[:,t].reshape(1,-1),W)
z[:,t+1]=functions.sigmoid(u[:,t+1])
np.dot(z[:,t+1].reshape(1,-1),W_out)
y[:,t]=functions.sigmoid(np.dot(z[:,t+1].reshape(1,-1), W_out))
  • 確認テスト

RNNのネットワークには大きく分けて3つの重みがある。一つは入力から現在の中間層を定義する際にかけられる重み、一つは中間層から出力を定義する際にかけられる重みである。
残り1つの重みについて説明せよ。
[回答]
中間層から中間層に回帰する際の重み。

  • 演習チャレンジ

以下は再帰ニューラルネットワークにおいて、構文木を入力として、再帰的に文全体の表現ベクトルを得るプログラムである。ただし、ニューラルネットワークの重みパラメータはグローバル変数として定義してあるものとし、activation関数は何らかの活性化関数であるとする。木構造再帰的な辞書で定義してあり、rootが最も外側の辞書であると仮定する。
(く)にあてはまるコードは何か
f:id:tibet:20211126154417p:plain
[回答]
(2)
隣接単語から表現ベクトルを作るという処理は、隣接している表現leftとrightを合わせたものを特徴量としてそこに重みをかけることで実現する。

1-2 BPTT

1-2-1 BPTTとは

 Back Propagation Through Timeで、RNNにおける逆誤差伝搬法である。計算結果から微分を逆算することで、不要な再帰的計算を避けることが出来る。

1-2-2 BPTTの数学的記述

RNNのモデル式を再掲する。
 u^t=W_{(in)}x^t+Wz^{t-1}+b
 z~t=f(W_{(in)}x^t+Wz^{t-1}+b)
 v^t=W_{(out)}z~t+c
y^t=g(W_{(out)}z^t+c)

u[:,t+1]=np.dot(X,W_in)+np.dot(z[:,t].reshape(1,-1),W)
z[:,t+1]=functions.sigmoid(u[:,t+1])
np.dot(z[:,t+1].reshape(1,-1),W_out)
y[:,t]=functions.sigmoid(np.dot(z[:,t+1].reshape(1,-1), W_out))

BPTTの微分の導出と、対応するコードを示す。
[]は、すべてのtに渡る計算をすることを示す。
まず、損失関数Eの入力重みWinによる微分
 \dfrac{\partial E}{\partial W_{(in)}}=\dfrac{\partial E}{\partial u^t}[\dfrac{u^t}{W_{(in)}}]^T=δ^t[x^t]^T

np.dot(X.T,delta[:t].reshape(1,-1))

Eの出力重みWoutによる微分
 \dfrac{\partial E}{\partial W_{(out)}}{\partial W}=\dfrac{\partial E}{\partial v^t}[\dfrac{v^t}{W_{(out)}}]^T=δ^{out,t}[x^t]^T

np.dot(z[:,t+1].reshape(-1,1),delta_out[:,t].reshape(-1,1))

Eの中間層重みによる微分
 \dfrac{\partial E}{\partial W}=\dfrac{\partial E}{\partial u^t}[\dfrac{\partial u^t}{\partial W}]^T=δ^t[z^{t-1}]^T

np.dot(z[:,t].reshape(-1,1),delta[:,t].reshape(1,-1))

Eのバイアスbによる微分
 \dfrac{\partial E}{\partial b}=\dfrac{\partial E}{\partial u^t}\dfrac{\partial u^t}{\partial b}^T=δ^t
Eのバイアスcによる微分
 \dfrac{\partial E}{\partial b}=\dfrac{\partial E}{\partial u^t}\dfrac{\partial u^t}{\partial c}^T=δ^{out,t}

  • 確認テスト

下図のy1をx,s0,s1,w,woutを用いて数式で表せ。
*バイアスは任意の文字で定義せよ。
*また中間層の出力にシグモイド関数g(x)を作用させよ。
f:id:tibet:20211127225100p:plain
[回答]
 y_1=g(W_{(out)}s_1+c)
 s_1=f(W_{(in)}x_1+Ws_0+b)

ここで、時間的に遡っての微分δt-z-1について計算してみる。
 \dfrac{\partial E}{\partial u^t}=\dfrac{\partial E}{\partial v^t}\dfrac{\partial v^t}{\partial u^t}=\dfrac{\partial E}{\partial v^t}\dfrac{\partial(W_{(out)}f(u^t)+c)}{\partial u^t}=δ^{out,t}W_{(out)}^Tf'(u^t)=δ^t
 δ^{t-1}=\dfrac{\partial E}{\partial u^{t-1}}=\dfrac{\partial E}{\partial u^t}\dfrac{\partial u^t}{\partial u^{t-1}}=δ^t\left(\dfrac{\partial u^t}{\partial z^{t-1}}\dfrac{\partial z^{t-1}}{\partial u^{t-1}}\right)=δ^t{Wf'(u^{t-1})}
 δ^{t-z-1}=δ^{t-z}(Wf'(u^{t-z-1})
第一式のコードは下記のように書ける

delta[:,t]=(np.dot(delta[:,t+1].t,W.T)+np.dot(delta_out[:,t].T,W_out.T))*functions.d_sigmoid(u[:,t+1])

最終的に各パラメータの更新式は以下のようになる。
 W^{t+1}_{(in)}=W^t_{(in)}-ε\dfrac{\partial E}{\partial W_{(in)}}=W^t_{(in)}-ε\displaystyle\sum^{T_t}_{z=0}δ^{t-z}[x^{t-z}]^T
 W^{t+1}_{(out)}=W^t_{(out)}-ε\dfrac{\partial E}{\partial W_{(out)}}=W^t_{(in)}-εδ^{out,t}[z^{t}]^T
 W^{t+1}=W^t-ε\dfrac{\partial E}{\partial W}=W^t-ε\displaystyle\sum^{T_t}_{z=0}δ^{t-z}[z^{t-z-1}]^T
 b^{t+1}=b^t-ε\dfrac{\partial E}{\partial b}=b^t-ε\displaystyle\sum^{T_t}_{z=0}δ^{t-z}
 c^{t+1}=c^t-ε\dfrac{\partial E}{\partial c}=c^t-εδ^{out,t}
上部3つの式に関するコードは下記のように書ける

W_in-=learning_rate*W_in_grad
W_out=learnign_rate*W_out_grad
W-=learning_rate*W_grad
BPTTの全体像

最後に、損失関数の中で、時間にわたる計算がどのようにあらわされるかを示す。
損失関数をE=loss(x1,x2)、教師データをdとすると、RNNの損失関数は下記のように表せる。
 E^t=loss(y^t,d^t)
             =loss(g(W_{(out)}z^t+c),t^t)
             =loss(g(W_{(out)}f(W_{(in)}x^t+Wz_{t-1}+b);c),dt)
上記のように、yを展開すると、 z^{t-1}の項が出てくる。さらに、 z^{t-1}を展開すると、 z^{t-2}の項が出てくる。このように、時間をさかのぼって計算することが出来る。

RNN,BPTTの実装演習

バイナリを加算して桁上げする操作を、RNNで実行する。

import numpy as np
from common import functions
import matplotlib.pyplot as plt

def d_tanh(x):
    return 1/(np.cosh(x) ** 2)

# データを用意
# 2進数の桁数
binary_dim = 8
# 最大値 + 1
largest_number = pow(2, binary_dim)
# largest_numberまで2進数を用意
binary = np.unpackbits(np.array([range(largest_number)],dtype=np.uint8).T,axis=1)

input_layer_size = 2
hidden_layer_size = 16
output_layer_size = 1

weight_init_std = 1
learning_rate = 0.1

iters_num = 10000
plot_interval = 100

# ウェイト初期化 (バイアスは簡単のため省略)
W_in = weight_init_std * np.random.randn(input_layer_size, hidden_layer_size)
W_out = weight_init_std * np.random.randn(hidden_layer_size, output_layer_size)
W = weight_init_std * np.random.randn(hidden_layer_size, hidden_layer_size)
# Xavier
# W_in = np.random.randn(input_layer_size, hidden_layer_size) / (np.sqrt(input_layer_size))
# W_out = np.random.randn(hidden_layer_size, output_layer_size) / (np.sqrt(hidden_layer_size))
# W = np.random.randn(hidden_layer_size, hidden_layer_size) / (np.sqrt(hidden_layer_size))
# He
# W_in = np.random.randn(input_layer_size, hidden_layer_size) / (np.sqrt(input_layer_size)) * np.sqrt(2)
# W_out = np.random.randn(hidden_layer_size, output_layer_size) / (np.sqrt(hidden_layer_size)) * np.sqrt(2)
# W = np.random.randn(hidden_layer_size, hidden_layer_size) / (np.sqrt(hidden_layer_size)) * np.sqrt(2)


# 勾配
W_in_grad = np.zeros_like(W_in)
W_out_grad = np.zeros_like(W_out)
W_grad = np.zeros_like(W)

u = np.zeros((hidden_layer_size, binary_dim + 1))
z = np.zeros((hidden_layer_size, binary_dim + 1))
y = np.zeros((output_layer_size, binary_dim))

delta_out = np.zeros((output_layer_size, binary_dim))
delta = np.zeros((hidden_layer_size, binary_dim + 1))

all_losses = []

for i in range(iters_num):
    
    # A, B初期化 (a + b = d)
    a_int = np.random.randint(largest_number/2)
    a_bin = binary[a_int] # binary encoding
    b_int = np.random.randint(largest_number/2)
    b_bin = binary[b_int] # binary encoding
    
    # 正解データ
    d_int = a_int + b_int
    d_bin = binary[d_int]
    
    # 出力バイナリ
    out_bin = np.zeros_like(d_bin)
    
    # 時系列全体の誤差
    all_loss = 0    
    
    # 時系列ループ
    for t in range(binary_dim):
        # 入力値
        X = np.array([a_bin[ - t - 1], b_bin[ - t - 1]]).reshape(1, -1)
        # 時刻tにおける正解データ
        dd = np.array([d_bin[binary_dim - t - 1]])
        
        u[:,t+1] = np.dot(X, W_in) + np.dot(z[:,t].reshape(1, -1), W)
        z[:,t+1] = functions.sigmoid(u[:,t+1])
#         z[:,t+1] = functions.relu(u[:,t+1])
#         z[:,t+1] = np.tanh(u[:,t+1])    
        y[:,t] = functions.sigmoid(np.dot(z[:,t+1].reshape(1, -1), W_out))


        #誤差
        loss = functions.mean_squared_error(dd, y[:,t])
        
        delta_out[:,t] = functions.d_mean_squared_error(dd, y[:,t]) * functions.d_sigmoid(y[:,t])        
        
        all_loss += loss

        out_bin[binary_dim - t - 1] = np.round(y[:,t])
    
    
    for t in range(binary_dim)[::-1]:
        X = np.array([a_bin[-t-1],b_bin[-t-1]]).reshape(1, -1)        

        delta[:,t] = (np.dot(delta[:,t+1].T, W.T) + np.dot(delta_out[:,t].T, W_out.T)) * functions.d_sigmoid(u[:,t+1])
#         delta[:,t] = (np.dot(delta[:,t+1].T, W.T) + np.dot(delta_out[:,t].T, W_out.T)) * functions.d_relu(u[:,t+1])
#         delta[:,t] = (np.dot(delta[:,t+1].T, W.T) + np.dot(delta_out[:,t].T, W_out.T)) * d_tanh(u[:,t+1])    

        # 勾配更新
        W_out_grad += np.dot(z[:,t+1].reshape(-1,1), delta_out[:,t].reshape(-1,1))
        W_grad += np.dot(z[:,t].reshape(-1,1), delta[:,t].reshape(1,-1))
        W_in_grad += np.dot(X.T, delta[:,t].reshape(1,-1))
    
    # 勾配適用
    W_in -= learning_rate * W_in_grad
    W_out -= learning_rate * W_out_grad
    W -= learning_rate * W_grad
    
    W_in_grad *= 0
    W_out_grad *= 0
    W_grad *= 0
    

    if(i % plot_interval == 0):
        all_losses.append(all_loss)        
        print("iters:" + str(i))
        print("Loss:" + str(all_loss))
        print("Pred:" + str(out_bin))
        print("True:" + str(d_bin))
        out_int = 0
        for index,x in enumerate(reversed(out_bin)):
            out_int += x * pow(2, index)
        print(str(a_int) + " + " + str(b_int) + " = " + str(out_int))
        print("------------")

lists = range(0, iters_num, plot_interval)
plt.plot(lists, all_losses, label="loss")
plt.show()

結果の一部表示

iters:0
Loss:1.245858363786796
Pred:[0 0 0 0 0 0 0 0]
True:[1 0 1 0 0 0 0 1]
82 + 79 = 0
------------
iters:100
Loss:0.9697292981990132
Pred:[0 1 1 1 1 0 1 1]
True:[1 0 0 1 1 0 1 1]
55 + 100 = 123
------------
.....(途中一部省略)
.------------
iters:9800
Loss:0.0009243998471599557
Pred:[0 1 1 0 0 1 0 1]
True:[0 1 1 0 0 1 0 1]
68 + 33 = 101
------------
iters:9900
Loss:0.0025480327477791424
Pred:[0 0 1 0 0 0 0 1]
True:[0 0 1 0 0 0 0 1]
19 + 14 = 33
------------

f:id:tibet:20211126154215p:plain
学習回数の増加に従って、PredとTrueの値が一致してきていることも分かる。
また、繰り返して学習することで、損失関数が低減していることがわかる。

  • コード演習問題

下図は、BPTTを行うプログラムである。なお簡単化のため活性化関数は恒等関数であるとする。
また、calculate_dout関数は損失関数を出力に関して偏微分した値を返す関数であるとする。
(お)にあてはまる値はどれか
f:id:tibet:20211128170542p:plain
[解答]
中間出力 h_{t}は過去の中間層出力 h_{t-1},h_[t-2},,,h_1に依存する。
 \dfrac{dh_{t}}{dh_{t-1}}=U
であることに注意すると、過去にさかのぼるたびにUがかけられる。つまり、(お)は、delta_t=delta_t.dot(U)となる。

ラビットチャレンジレポート 深層学習Day2 その3

Section4 畳み込みニューラルネットワークの概念

畳み込みニューラルネットワーク(CNN)は各次元で連続性のあるデータを扱うのに適したネットワークである。
CNNの構造の代表例は下記である。
入力層

畳み込み層

畳み込み層

プーリング層

畳み込み層

畳み込み層

プーリング層

全結合層

出力層

基本的な構造として有名なものは、LeNet(入力データ行列32×32)がある。

4-1畳み込み層

 畳み込み層では、画像の場合、建て、横、チャンネルの3次元のデータをそのまま学習し、次に伝えることが出来る。
 具体的には、入力画像Xに、フィルター行列Wを左端から順番にアダマール積をとって行列を出力する。パディング0、ストライド1の場合、入力データがm×m行列、フィルター行列が(m-t)×(m-t)ならば、出力行列は(t+1)×(t+1)行列になる。

4-1-1. 畳み込みの演算概念:バイアス

 一般に畳み込み演算では、入力データとフィルタの演算後、バイアスを加える。バイアスの設定によって収束の速さが変化する。

4-1-1.畳み込みの演算概念:パディング

 一般に畳み込み演算を行うと、入力データに対して、出力データは小さくなる。例えば、入力データが4×4行列、フィルタが3×3行列の場合、出力データは2×2行列になる。
 このため、何度も畳み込みを繰り返すと出力データが小さくなりすぎたり、出力データが計算機処理には適さないサイズ(2^mが適する)になる場合がある。
 これを防ぐためにパディング処理がある。
 パディングは、元の入力画像の周りにダミーにデータを配置する手法である。
 元の入力データ 
 f:id:tibet:20211121174014p:plain
 パディングデータ
 f:id:tibet:20211121174510p:plain
 ダミーデータは0のほか、隣接データを使用するなどの手法がある。

4-1-2.畳み込みの演算概念:ストライド

 畳み込み時に入力データ上をフィルタが動くステップ間隔をストライドと呼ぶ。
 ストライドは1のほか、2以上をとることもある。ストライドが大きいほど出力データは小さくなる。

4-1-3.畳み込みの演算概念:チャネル

 入力データやフィルタの行列の枚数のことをチャネル言う。

全結合層で画像を学習した場合の欠点

 画像の場合、縦、横、チャネルの各次元のデータの位置関係が重要だが、全結合層の場合、1次元のデータとして処理され、各次元間の関連性が反映されない。
 そのため、画像を全結合層で学習するのは難しい。

4-2 プーリング層

 プーリング層は、特徴マップを空間的な局所ごとに代表値にまとめて空間解像度を下げる(ダウンサンプリング)層である。これにより、入力データのフィルタ形状による位置ずれを吸収する役割がある。
 プーリング層の具体的な操作は、畳み込みフィルタと同じく、ある大きさのウィンドウで左上からストライド分ずつ移動しながら代表値をとる。
 代表値として、ウインドウ内の最大値をとる場合をMax Pooling、ウィンドウ内の平均値をとる場合をAverage Poolingとそれぞれ呼ぶ。

  • 確認テスト

サイズ6×6の入力画像を、サイズ2×2のフィルタで畳み込んだ時の出力画像のサイズを答えよ。なおストライドとパディングは1とする。
[解答]
出力画像の高さOH、幅HWは下記の公式で導出できる
 O_H=\dfrac{入力画像の高さ+2×パディング高さ-フィルタ高さ}{ストライド}+1
 O_W=\dfrac{入力画像の幅+2×パディング幅-フィルタ幅}{ストライド}+1
これらを用いると、出力画像のサイズは7×7である。 

畳み込みの実装演習

画像データを2次元配列に変換する関数、畳み込み処理の関数、マックスプーリングの関数を定義し、CNNネットワークを実装する。

画像データを2次元配列に変換する(im2col)

import pickle
import numpy as np
from collections import OrderedDict
from common import layers
from common import optimizer
from data.mnist import load_mnist
import matplotlib.pyplot as plt

# 画像データを2次元配列に変換
'''
input_data: 入力値
filter_h: フィルターの高さ
filter_w: フィルターの横幅
stride: ストライド
pad: パディング
'''
def im2col(input_data, filter_h, filter_w, stride=1, pad=0):
    # N: number, C: channel, H: height, W: width
    N, C, H, W = input_data.shape
    # 切り捨て除算
    out_h = (H + 2 * pad - filter_h)//stride + 1
    out_w = (W + 2 * pad - filter_w)//stride + 1

    img = np.pad(input_data, [(0,0), (0,0), (pad, pad), (pad, pad)], 'constant')
    col = np.zeros((N, C, filter_h, filter_w, out_h, out_w))

    for y in range(filter_h):
        y_max = y + stride * out_h
        for x in range(filter_w):
            x_max = x + stride * out_w
            col[:, :, y, x, :, :] = img[:, :, y:y_max:stride, x:x_max:stride]

    col = col.transpose(0, 4, 5, 1, 2, 3) # (N, C, filter_h, filter_w, out_h, out_w) -> (N, filter_w, out_h, out_w, C, filter_h)    
    
    col = col.reshape(N * out_h * out_w, -1)
    return col

畳み込みのクラスを作成する

class Convolution:
    # W: フィルター, b: バイアス
    def __init__(self, W, b, stride=1, pad=0):
        self.W = W
        self.b = b
        self.stride = stride
        self.pad = pad
        
        # 中間データ(backward時に使用)
        self.x = None   
        self.col = None
        self.col_W = None
        
        # フィルター・バイアスパラメータの勾配
        self.dW = None
        self.db = None

    def forward(self, x):
        # FN: filter_number, C: channel, FH: filter_height, FW: filter_width
        FN, C, FH, FW = self.W.shape
        N, C, H, W = x.shape
        # 出力値のheight, width
        out_h = 1 + int((H + 2 * self.pad - FH) / self.stride)
        out_w = 1 + int((W + 2 * self.pad - FW) / self.stride)
        
        # xを行列に変換
        col = im2col(x, FH, FW, self.stride, self.pad)
        # フィルターをxに合わせた行列に変換
        col_W = self.W.reshape(FN, -1).T

        out = np.dot(col, col_W) + self.b
        # 計算のために変えた形式を戻す
        out = out.reshape(N, out_h, out_w, -1).transpose(0, 3, 1, 2)

        self.x = x
        self.col = col
        self.col_W = col_W

        return out

    def backward(self, dout):
        FN, C, FH, FW = self.W.shape
        dout = dout.transpose(0, 2, 3, 1).reshape(-1, FN)

        self.db = np.sum(dout, axis=0)
        self.dW = np.dot(self.col.T, dout)
        self.dW = self.dW.transpose(1, 0).reshape(FN, C, FH, FW)

        dcol = np.dot(dout, self.col_W.T)
        # dcolを画像データに変換
        dx = col2im(dcol, self.x.shape, FH, FW, self.stride, self.pad)

        return dx
||>
プーリングのクラス
>|python|
class Pooling:
    def __init__(self, pool_h, pool_w, stride=1, pad=0):
        self.pool_h = pool_h
        self.pool_w = pool_w
        self.stride = stride
        self.pad = pad
        
        self.x = None
        self.arg_max = None

    def forward(self, x):
        N, C, H, W = x.shape
        out_h = int(1 + (H - self.pool_h) / self.stride)
        out_w = int(1 + (W - self.pool_w) / self.stride)
        
        # xを行列に変換
        col = im2col(x, self.pool_h, self.pool_w, self.stride, self.pad)
        # プーリングのサイズに合わせてリサイズ
        col = col.reshape(-1, self.pool_h*self.pool_w)
        
        #maxプーリング
        arg_max = np.argmax(col, axis=1)
        out = np.max(col, axis=1)
        out = out.reshape(N, out_h, out_w, C).transpose(0, 3, 1, 2)

        self.x = x
        self.arg_max = arg_max

        return out

    def backward(self, dout):
        dout = dout.transpose(0, 2, 3, 1)
        
        pool_size = self.pool_h * self.pool_w
        dmax = np.zeros((dout.size, pool_size))
        dmax[np.arange(self.arg_max.size), self.arg_max.flatten()] = dout.flatten()
        dmax = dmax.reshape(dout.shape + (pool_size,)) 
        
        dcol = dmax.reshape(dmax.shape[0] * dmax.shape[1] * dmax.shape[2], -1)
        dx = col2im(dcol, self.x.shape, self.pool_h, self.pool_w, self.stride, self.pad)
        
        return dx

CNNのクラスを定義する

class SimpleConvNet:
    # conv - relu - pool - affine - relu - affine - softmax
    def __init__(self, input_dim=(1, 28, 28), conv_param={'filter_num':30, 'filter_size':5, 'pad':0, 'stride':1},
                 hidden_size=100, output_size=10, weight_init_std=0.01):
        filter_num = conv_param['filter_num']        
        filter_size = conv_param['filter_size']
        filter_pad = conv_param['pad']
        filter_stride = conv_param['stride']
        input_size = input_dim[1]
        conv_output_size = (input_size - filter_size + 2 * filter_pad) / filter_stride + 1
        pool_output_size = int(filter_num * (conv_output_size / 2) * (conv_output_size / 2))

        # 重みの初期化
        self.params = {}
        self.params['W1'] = weight_init_std * np.random.randn(filter_num, input_dim[0], filter_size, filter_size)
        self.params['b1'] = np.zeros(filter_num)
        self.params['W2'] = weight_init_std * np.random.randn(pool_output_size, hidden_size)
        self.params['b2'] = np.zeros(hidden_size)
        self.params['W3'] = weight_init_std * np.random.randn(hidden_size, output_size)
        self.params['b3'] = np.zeros(output_size)

        # レイヤの生成
        self.layers = OrderedDict()
        self.layers['Conv1'] = layers.Convolution(self.params['W1'], self.params['b1'], conv_param['stride'], conv_param['pad'])
        self.layers['Relu1'] = layers.Relu()
        self.layers['Pool1'] = layers.Pooling(pool_h=2, pool_w=2, stride=2)
        self.layers['Affine1'] = layers.Affine(self.params['W2'], self.params['b2'])
        self.layers['Relu2'] = layers.Relu()
        self.layers['Affine2'] = layers.Affine(self.params['W3'], self.params['b3'])

        self.last_layer = layers.SoftmaxWithLoss()

    def predict(self, x):
        for key in self.layers.keys():
            x = self.layers[key].forward(x)
        return x
        
    def loss(self, x, d):
        y = self.predict(x)
        return self.last_layer.forward(y, d)

    def accuracy(self, x, d, batch_size=100):
        if d.ndim != 1 : d = np.argmax(d, axis=1)
        
        acc = 0.0
        
        for i in range(int(x.shape[0] / batch_size)):
            tx = x[i*batch_size:(i+1)*batch_size]
            td = d[i*batch_size:(i+1)*batch_size]
            y = self.predict(tx)
            y = np.argmax(y, axis=1)
            acc += np.sum(y == td) 
        
        return acc / x.shape[0]

    def gradient(self, x, d):
        # forward
        self.loss(x, d)
        
        # backward
        dout = 1
        dout = self.last_layer.backward(dout)
        layers = list(self.layers.values())
        
        layers.reverse()
        for layer in layers:
            dout = layer.backward(dout)

        # 設定
        grad = {}
        grad['W1'], grad['b1'] = self.layers['Conv1'].dW, self.layers['Conv1'].db
        grad['W2'], grad['b2'] = self.layers['Affine1'].dW, self.layers['Affine1'].db
        grad['W3'], grad['b3'] = self.layers['Affine2'].dW, self.layers['Affine2'].db

        return grad

mnistデータセットを用いて、CNNの学習と評価を行う。

from common import optimizer

# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(flatten=False)

print("データ読み込み完了")

# 処理に時間のかかる場合はデータを削減 
x_train, d_train = x_train[:5000], d_train[:5000]
x_test, d_test = x_test[:1000], d_test[:1000]


network = SimpleConvNet(input_dim=(1,28,28), conv_param = {'filter_num': 30, 'filter_size': 5, 'pad': 0, 'stride': 1},
                        hidden_size=100, output_size=10, weight_init_std=0.01)

optimizer = optimizer.Adam()

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10



for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]
    
    grad = network.gradient(x_batch, d_batch)
    optimizer.update(network.params, grad)

    loss = network.loss(x_batch, d_batch)
    train_loss_list.append(loss)

    if (i+1) % plot_interval == 0:
        accr_train = network.accuracy(x_train, d_train)
        accr_test = network.accuracy(x_test, d_test)
        accuracies_train.append(accr_train)
        accuracies_test.append(accr_test)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))               

lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果の図示
f:id:tibet:20211121185204p:plain
訓練制度、テスト精度とも1.0に近く、高い精度で学習と汎化が出来ていることが確認できる。

Section5 最新のCNN

ここでは、Alexnetを取り上げる。

Alexnetとは

 Alex Krizhevskyが博士指導教官であるGeoffrey Hintonらとともに設計したネットワーク。
 2012年の国際的な画像認識コンテストであるILSVRC2012に参加し、エラー率が次点より10.8%も低いという圧倒的なスコアで優勝した。
 Alexの優勝によって、特徴量設計を人間が行うのではなく、たくさんのデータの中から機械が見つけ出す深層学習の特徴によって大きな優位性が実現できることが示され、深層学習ブームの端緒となった。
 Alexnetの論文はコンピュータビジョンで発表された最も影響力のある論文と考えられている。

Alexnetの構造

 8層の構造になっており、3つの畳み込み層、2つのプーリング層、及び3つの全結合層によって構成されている。活性化関数にはReLuを使用している。
f:id:tibet:20211121191115p:plain

全結合層の処理
  • Flatten:前層のデータを、1次元に並べる処理。13×13×256行列なら、4096×1の1次元行列に変換する。
  • Global max pooling: 1チャネルずつMax Poolingを行う処理。13×13×256行列なら、13×13行列にMax Poolingを適用し、256×1の1次元行列に変換する。
  • Global average pooling: 1チャネルずつAverage Poolingを行う処理。13×13×256行列なら、13×13行列にAverage Poolingを適用し、256×1の1次元行列に変換する。

Flattenより、Global Max PoolingやGlobal Average Poolingのほうが特徴量をよくとらえ、より少ない要素数で、より高い性能を発揮する場合が多い。

実装演習

ここでは、畳み込み層6つ、プーリング層3つの構成で、出力の全結合層のドロップアウトを適用したDeepConvNetを作成し、学習評価する。
DeepConvNetの作成

import pickle
import numpy as np
from collections import OrderedDict
from common import layers
from data.mnist import load_mnist
import matplotlib.pyplot as plt
from common import optimizer

class DeepConvNet:
    '''
    認識率99%以上の高精度なConvNet

    conv - relu - conv- relu - pool -
    conv - relu - conv- relu - pool -
    conv - relu - conv- relu - pool -
    affine - relu - dropout - affine - dropout - softmax
    '''
    def __init__(self, input_dim=(1, 28, 28),
                 conv_param_1 = {'filter_num':16, 'filter_size':3, 'pad':1, 'stride':1},
                 conv_param_2 = {'filter_num':16, 'filter_size':3, 'pad':1, 'stride':1},
                 conv_param_3 = {'filter_num':32, 'filter_size':3, 'pad':1, 'stride':1},
                 conv_param_4 = {'filter_num':32, 'filter_size':3, 'pad':2, 'stride':1},
                 conv_param_5 = {'filter_num':64, 'filter_size':3, 'pad':1, 'stride':1},
                 conv_param_6 = {'filter_num':64, 'filter_size':3, 'pad':1, 'stride':1},
                 hidden_size=50, output_size=10):
        # 重みの初期化===========
        # 各層のニューロンひとつあたりが、前層のニューロンといくつのつながりがあるか
        pre_node_nums = np.array([1*3*3, 16*3*3, 16*3*3, 32*3*3, 32*3*3, 64*3*3, 64*4*4, hidden_size])
        wight_init_scales = np.sqrt(2.0 / pre_node_nums)  # Heの初期値
        
        self.params = {}
        pre_channel_num = input_dim[0]
        for idx, conv_param in enumerate([conv_param_1, conv_param_2, conv_param_3, conv_param_4, conv_param_5, conv_param_6]):
            self.params['W' + str(idx+1)] = wight_init_scales[idx] * np.random.randn(conv_param['filter_num'], pre_channel_num, conv_param['filter_size'], conv_param['filter_size'])
            self.params['b' + str(idx+1)] = np.zeros(conv_param['filter_num'])
            pre_channel_num = conv_param['filter_num']
        self.params['W7'] = wight_init_scales[6] * np.random.randn(pre_node_nums[6], hidden_size)
        print(self.params['W7'].shape)
        self.params['b7'] = np.zeros(hidden_size)
        self.params['W8'] = wight_init_scales[7] * np.random.randn(pre_node_nums[7], output_size)
        self.params['b8'] = np.zeros(output_size)

        # レイヤの生成===========
        self.layers = []
        self.layers.append(layers.Convolution(self.params['W1'], self.params['b1'], 
                           conv_param_1['stride'], conv_param_1['pad']))
        self.layers.append(layers.Relu())
        self.layers.append(layers.Convolution(self.params['W2'], self.params['b2'], 
                           conv_param_2['stride'], conv_param_2['pad']))
        self.layers.append(layers.Relu())
        self.layers.append(layers.Pooling(pool_h=2, pool_w=2, stride=2))
        self.layers.append(layers.Convolution(self.params['W3'], self.params['b3'], 
                           conv_param_3['stride'], conv_param_3['pad']))
        self.layers.append(layers.Relu())
        self.layers.append(layers.Convolution(self.params['W4'], self.params['b4'],
                           conv_param_4['stride'], conv_param_4['pad']))
        self.layers.append(layers.Relu())
        self.layers.append(layers.Pooling(pool_h=2, pool_w=2, stride=2))
        self.layers.append(layers.Convolution(self.params['W5'], self.params['b5'],
                           conv_param_5['stride'], conv_param_5['pad']))
        self.layers.append(layers.Relu())
        self.layers.append(layers.Convolution(self.params['W6'], self.params['b6'],
                           conv_param_6['stride'], conv_param_6['pad']))
        self.layers.append(layers.Relu())
        self.layers.append(layers.Pooling(pool_h=2, pool_w=2, stride=2))
        self.layers.append(layers.Affine(self.params['W7'], self.params['b7']))
        self.layers.append(layers.Relu())
        self.layers.append(layers.Dropout(0.5))
        self.layers.append(layers.Affine(self.params['W8'], self.params['b8']))
        self.layers.append(layers.Dropout(0.5))
        
        self.last_layer = layers.SoftmaxWithLoss()

    def predict(self, x, train_flg=False):
        for layer in self.layers:
            if isinstance(layer, layers.Dropout):
                x = layer.forward(x, train_flg)
            else:
                x = layer.forward(x)
        return x

    def loss(self, x, d):
        y = self.predict(x, train_flg=True)
        return self.last_layer.forward(y, d)

    def accuracy(self, x, d, batch_size=100):
        if d.ndim != 1 : d = np.argmax(d, axis=1)

        acc = 0.0

        for i in range(int(x.shape[0] / batch_size)):
            tx = x[i*batch_size:(i+1)*batch_size]
            td = d[i*batch_size:(i+1)*batch_size]
            y = self.predict(tx, train_flg=False)
            y = np.argmax(y, axis=1)
            acc += np.sum(y == td)

        return acc / x.shape[0]

    def gradient(self, x, d):
        # forward
        self.loss(x, d)

        # backward
        dout = 1
        dout = self.last_layer.backward(dout)

        tmp_layers = self.layers.copy()
        tmp_layers.reverse()
        for layer in tmp_layers:
            dout = layer.backward(dout)

        # 設定
        grads = {}
        for i, layer_idx in enumerate((0, 2, 5, 7, 10, 12, 15, 18)):
            grads['W' + str(i+1)] = self.layers[layer_idx].dW
            grads['b' + str(i+1)] = self.layers[layer_idx].db

        return grads

mnistデータによる学習と評価を行う。

(x_train, d_train), (x_test, d_test) = load_mnist(flatten=False)

# 処理に時間のかかる場合はデータを削減 
x_train, d_train = x_train[:5000], d_train[:5000]
x_test, d_test = x_test[:1000], d_test[:1000]

print("データ読み込み完了")

network = DeepConvNet()  
optimizer = optimizer.Adam()

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10


for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]
    
    grad = network.gradient(x_batch, d_batch)
    optimizer.update(network.params, grad)

    loss = network.loss(x_batch, d_batch)
    train_loss_list.append(loss)

    if (i+1) % plot_interval == 0:
        accr_train = network.accuracy(x_train, d_train)
        accr_test = network.accuracy(x_test, d_test)
        accuracies_train.append(accr_train)
        accuracies_test.append(accr_test)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))               

lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果の図示
f:id:tibet:20211121194151p:plain
訓練データ、テストデータ双方に98%以上の高い正解率を示した。

ラビットチャレンジレポート 深層学習DAY2 その2

Section2 学習率最適化手法について

学習率とは

学習を通して教師データとの誤差を最小にするネットワークを作成するために取られる手法を勾配降下法という。勾配降下法では、パラメータ更新する際の変化を規定する係数εを学習率という。
· w^{(t+1)}=w^{(t)}-ε∇E

  • 学習率の値が大きい場合

最適値にいつまでもたどりつかない

  • 学習率の値が小さい場合

収束まで時間がかかる。
局所極小解に収束する場合がある。

学習率の決め方

  • 初期の学習率設定指針
    • 初期の学習率を大きく設定し、徐々に小さくしていく。
    • パラメータごとに学習率を可変させる

学習率最適化手法を利用する 

学習率最適化手法:モメンタム

 誤差をパラメータで微分したものと学習率の席を減算した後、現在の重みに前回の重みを減算した値と慣性の積を加算する。誤差関数の勾配のほかに、重みの一回当たりの変化(≒速度)の要素を加えるので、モメンタム(慣性)と呼ばれる。
 w^{(t+1)}=w^{(t)}+V_t
 V_t=μV_{(t-1)}-ε∇E
慣性:μ

  • モメンタムの特徴
    • 局所最適解にならず、大域的最適解にたどり着きやすい
    • 谷間についてから最も低い位置(最適値)に行くまでの時間が早い
    • 傾きが急だと慣性が大きすぎて最適解を通り過ぎてしまったり、中々収束しないという課題がある。

学習率最適化手法:Adagrad

Adagradは、学習率を学習の進み方に応じて再定義する手法である。具体的には、誤差関数の微分の2乗を加算するパラメータhを導入し、これを用いて学習率を調整する。
 w^{(t+1)}=w^{(t)}-ε\dfrac{1}{\sqrt{h_t}+θ}∇E
 h_0=θ
 h_t=h_{(t-1)}+(∇E)^2

  • Adagradの特徴
    • 勾配の緩やかな斜面に対して最適値に近づける
    • 学習率が徐々に小さくなるので鞍点問題を引き起こすことがある。

学習率最適化手法:RMSProp

 自己減衰係数αを導入して、Adagradの鞍点問題を起こしやすいという欠点を改良したもの。
 w^{(t+1)}=w^{(t)}-ε\dfrac{1}{\sqrt({h_t}+θ}∇E
 h_t=αh_t+(1-α)(∇E)^2

  • RMSPropの特徴
    • 局所最適解にならず、大域的最適解にたどり着きやすい
    • ハイパーパラメータの調整が必要な場面が少ない

  • 確認テスト

モメンタム・AdaGrad・RMSPropの特徴をそれぞれ簡潔に説明せよ
[解答]

  • モメンタム:勾配降下法に慣性項Vtを導入して大域的最適解にたどり着きやすくしたもの
  • AdaGrad:勾配降下法の学習率が∇Eに応じて変化するよう調整したもので、緩やかな勾配に強い。一方で、鞍点問題を起こしやすいという課題がある。
  • RMSProp:Adagradの鞍点問題を引き起こしやすいという欠点を改良したもの。

学習率最適化手法:Adam

 モメンタムとRMSPropの長所を取り入れたアルゴリズムである。比較的スムーズに裁定買いにたどり着くことが出来る。
 w_t=w_{(t-1)}-η\dfrac{\hat{m}_t}{\sqrt{\hat{ν}_t}+ε}
 \hat{m}_t=\dfrac{m_t}{1-β_1^t}
 \hat{ν}_t=\dfrac{ν_t}{1-β_2^t}
 m_t=β_1m_{t-1}+(1-·β_1)∇E
 ν_t=β_2ν_{t-1}+(1-β_2)∇E^2

実装演習

SGD、モメンタム、Adagrad、RMSProp、Adamの学習推移を比較し、Optimizerの違いが学習にどのような影響があるかを確認する。

  • SDG
import numpy as np
from collections import OrderedDict
from common import layers
from data.mnist import load_mnist
import matplotlib.pyplot as plt
from multi_layer_net import MultiLayerNet


# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True, one_hot_label=True)

print("データ読み込み完了")

# batch_normalizationの設定 =======================
# use_batchnorm = True
use_batchnorm = False
# ====================================================


network = MultiLayerNet(input_size=784, hidden_size_list=[40, 20], output_size=10, activation='sigmoid', weight_init_std=0.01,
                       use_batchnorm=use_batchnorm)

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100
learning_rate = 0.01

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    # 勾配
    grad = network.gradient(x_batch, d_batch)
    
    for key in ('W1', 'W2', 'W3', 'b1', 'b2', 'b3'):
        network.params[key] -= learning_rate * grad[key]
        
        loss = network.loss(x_batch, d_batch)
        train_loss_list.append(loss)
    
    
    if (i + 1) % plot_interval == 0:
        accr_test = network.accuracy(x_test, d_test)
        accuracies_test.append(accr_test)        
        accr_train = network.accuracy(x_batch, d_batch)
        accuracies_train.append(accr_train)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))

        
lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果の図示
f:id:tibet:20211120145058p:plain
Accuracyが向上せず、学習が進んでいない。

  • モメンタム
# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True, one_hot_label=True)

print("データ読み込み完了")

# batch_normalizationの設定 =======================
# use_batchnorm = True
use_batchnorm = False
# ====================================================

network = MultiLayerNet(input_size=784, hidden_size_list=[40, 20], output_size=10, activation='sigmoid', weight_init_std=0.01,
                       use_batchnorm=use_batchnorm)

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100
learning_rate = 0.3
# 慣性
momentum = 0.9

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    # 勾配
    grad = network.gradient(x_batch, d_batch)
    if i == 0:
        v = {}
    for key in ('W1', 'W2', 'W3', 'b1', 'b2', 'b3'):
        if i == 0:
            v[key] = np.zeros_like(network.params[key])
        v[key] = momentum * v[key] - learning_rate * grad[key]
        network.params[key] += v[key]

        loss = network.loss(x_batch, d_batch)
        train_loss_list.append(loss)
        
    if (i + 1) % plot_interval == 0:
        accr_test = network.accuracy(x_test, d_test)
        accuracies_test.append(accr_test)        
        accr_train = network.accuracy(x_batch, d_batch)
        accuracies_train.append(accr_train)

        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))
        
        
lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果の図示
f:id:tibet:20211120145247p:plain
学習が進み、600回程度でほぼ収束していることがわかる。

  • Adagrad
# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True, one_hot_label=True)

print("データ読み込み完了")

# batch_normalizationの設定 =======================
# use_batchnorm = True
use_batchnorm = False
# ====================================================

network = MultiLayerNet(input_size=784, hidden_size_list=[40, 20], output_size=10, activation='sigmoid', weight_init_std=0.01,
                       use_batchnorm=use_batchnorm)

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100
learning_rate = 0.1

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    # 勾配
    grad = network.gradient(x_batch, d_batch)
    if i == 0:
        h = {}
    for key in ('W1', 'W2', 'W3', 'b1', 'b2', 'b3'):
        if i == 0:
            h[key] = np.full_like(network.params[key], 1e-4)
        else:
            h[key] += np.square(grad[key])
        network.params[key] -= learning_rate * grad[key] / (np.sqrt(h[key]))

        loss = network.loss(x_batch, d_batch)
        train_loss_list.append(loss)        
        
        
    if (i + 1) % plot_interval == 0:
        accr_test = network.accuracy(x_test, d_test)
        accuracies_test.append(accr_test)        
        accr_train = network.accuracy(x_batch, d_batch)
        accuracies_train.append(accr_train)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))
        
        
lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果の図示
f:id:tibet:20211120145457p:plain
最終的に学習は収束しているが、モメンタムより学習に時間がかかっている。

  • RMSProp
# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True, one_hot_label=True)

print("データ読み込み完了")

# batch_normalizationの設定 =======================
# use_batchnorm = True
use_batchnorm = False
# ====================================================

network = MultiLayerNet(input_size=784, hidden_size_list=[40, 20], output_size=10, activation='sigmoid', weight_init_std=0.01,
                       use_batchnorm=use_batchnorm)

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100
learning_rate = 0.01
decay_rate = 0.99

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    # 勾配
    grad = network.gradient(x_batch, d_batch)
    if i == 0:
        h = {}
    for key in ('W1', 'W2', 'W3', 'b1', 'b2', 'b3'):
        if i == 0:
            h[key] = np.zeros_like(network.params[key])
        h[key] *= decay_rate
        h[key] += (1 - decay_rate) * np.square(grad[key])
        network.params[key] -= learning_rate * grad[key] / (np.sqrt(h[key]) + 1e-7)

        loss = network.loss(x_batch, d_batch)
        train_loss_list.append(loss)                
        
    if (i + 1) % plot_interval == 0:
        accr_test = network.accuracy(x_test, d_test)
        accuracies_test.append(accr_test)        
        accr_train = network.accuracy(x_batch, d_batch)
        accuracies_train.append(accr_train)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))
        
        
lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果の図示
f:id:tibet:20211120145634p:plain
収束は、Adagradより早い

  • Adam
# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True, one_hot_label=True)

print("データ読み込み完了")

# batch_normalizationの設定 =======================
# use_batchnorm = True
use_batchnorm = False
# ====================================================

network = MultiLayerNet(input_size=784, hidden_size_list=[40, 20], output_size=10, activation='sigmoid', weight_init_std=0.01,
                       use_batchnorm=use_batchnorm)

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100
learning_rate = 0.01
beta1 = 0.9
beta2 = 0.999

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    # 勾配
    grad = network.gradient(x_batch, d_batch)
    if i == 0:
        m = {}
        v = {}
    learning_rate_t  = learning_rate * np.sqrt(1.0 - beta2 ** (i + 1)) / (1.0 - beta1 ** (i + 1))    
    for key in ('W1', 'W2', 'W3', 'b1', 'b2', 'b3'):
        if i == 0:
            m[key] = np.zeros_like(network.params[key])
            v[key] = np.zeros_like(network.params[key])
            
        m[key] += (1 - beta1) * (grad[key] - m[key])
        v[key] += (1 - beta2) * (grad[key] ** 2 - v[key])            
        network.params[key] -= learning_rate_t * m[key] / (np.sqrt(v[key]) + 1e-7)                
        
        loss = network.loss(x_batch, d_batch)
        train_loss_list.append(loss)        
        
    if (i + 1) % plot_interval == 0:
        accr_test = network.accuracy(x_test, d_test)
        accuracies_test.append(accr_test)        
        accr_train = network.accuracy(x_batch, d_batch)
        accuracies_train.append(accr_train)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))
                

lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果の図示
f:id:tibet:20211120145737p:plain
RMSPropと同程度に学習の収束が早く、かつ他のOptimizerより滑らかに精度が向上していることが確認できる。

Section3 過学習

過学習とは

訓練誤差は非常に小さくなるが、テスト誤差はそれほど下がらず、両者で学習曲線が乖離すること。
原因としては、一般にニューラルネットワークのパラメータは非常に多く、自由度が高いため、特定の訓練データに特化して学習が進みやすいことにある。
過学習を抑制する手法として正則化がある。

正則化

正則化とは、ネットワークの自由度(総数、ノード数、パラメータ数etc)を制約することである。

  • 確認テスト

機械学習で使われる線形モデル(線形回帰、主成分分析..)の正則化は、モデルの重みを制限することで可能となる。
前述の線形モデルの正則化手法の中にリッジ回帰という手法があり、その特徴として正しいものを選択しなさい
(a)ハイパーパラメータを大きな値に設定するとすべての重みが限りなく0に近づく
(b)ハイパーパラメータ0に設定すると、非線形回帰となる
(c)バイアス項についても、正則化される
(d)リッジ回帰の場合、隠れ層に対して正則化項を加えることになる
[解答]
(a)
正則化項の影響が大きくなりすぎるため。
ハイパーパラメータは、適切な値を選ぶ必要がある。

Weight decay(荷重減衰)

過学習の原因として、重みが大きな値をとることがあげられる。大きい値をとる重みとは学習において重要な値とみなされる。
 過学習が起こりそうな重みの大きさ以下で重みをコントロールし、かつ重みの大きさにばらつきを出せば過学習が起こりにくいと考えられる。そこで対策として、誤差に対して正則化項を加算して重みを抑制する手法が考案された。

L1,L2正則化

誤差関数に、pノルムを正則化項として加算する。
 E_n(w)+\dfrac{1}{p}λ||x||_p
 ||x||_p=(|x|^p+..+|x_n|^p)^{(1/p)}

p=1の場合、L1正則化と呼ぶ。
p=2の場合、L2正則化と呼ぶ。
L1ノルムは、要素間のマンハッタン距離を計算することに相当する。L2ノルムは、要素間のユークリッド距離を計算することに相当する。

 微分不能な点があるため、微分を用いた最適化手法が使用出来ないが、最適点である重みをゼロにするスパース化が出来るため、計算効率を高くすることが出来る。

 スパース化は難しいが、全点微分可能なため、微分を用いた最適化手法を使用できる。

  • 確認テスト

下図について、L1正則化を表しているグラフはどちらか答えよ。
f:id:tibet:20211120194521p:plain
[解答]
 右のLasso回帰を表すグラフである。

  • 例題チャレンジ1

f:id:tibet:20211120200447p:plain
[解答]
(4)
Ridge正則項を加えた誤差関数は下記のように表せる。
 E_t=E_t+λ||w||^2
両辺を w_i微分すると、
 grad=\dfrac{\partial E_t}{\partial w_i}=\dfrac{\partial E_{t-1}}{\partial w_i}+λ\cdot 2|w_i|
 rate=2λと置けば、後半の項はrate*paramと表せる。

  • 例題チャレンジ2

f:id:tibet:20211120201107p:plain
[解答]
(3)
Lasso正則項を加えた誤差関数は下記のように表せる。
 E_t=E_t+λ|w|
 w_iで両辺を微分した場合、正則化項はλ*( w_iの符号に応じて1,0,-1)となる。
よって、sign(param)である。

ドロップアウト

ドロップアウトとは、ランダムにノードを削除して学習させる手法。ノード数を減少させることで過学習を防ぐ。
メリットとしては、データ量を変化させずに異なるモデルを学習させると解釈できることである。

実装演習

まず、対策なしの過学習の状態を確認する。その後、L2正則化、L1正則化ドロップアウトの効果をそれぞれ確認し、最後にL1正則化ドロップアウトを複合した場合の効果を確認する。

対策なしで訓練データとテストデータの精度の差異を確認する。

import numpy as np
from collections import OrderedDict
from common import layers
from data.mnist import load_mnist
import matplotlib.pyplot as plt
from multi_layer_net import MultiLayerNet
from common import optimizer


(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True)

print("データ読み込み完了")

# 過学習を再現するために、学習データを削減
x_train = x_train[:300]
d_train = d_train[:300]

network = MultiLayerNet(input_size=784, hidden_size_list=[100, 100, 100, 100, 100, 100], output_size=10)
optimizer = optimizer.SGD(learning_rate=0.01)

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10


for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    grad = network.gradient(x_batch, d_batch)
    optimizer.update(network.params, grad)

    loss = network.loss(x_batch, d_batch)
    train_loss_list.append(loss)
        
    if (i+1) % plot_interval == 0:
        accr_train = network.accuracy(x_train, d_train)
        accr_test = network.accuracy(x_test, d_test)
        accuracies_train.append(accr_train)
        accuracies_test.append(accr_test)

        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))        

lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果の表示
f:id:tibet:20211120205640p:plain
訓練データに対しては、accuracy=1.0と非常に精度高く予測で来ているが、テストデータに対しては、0.77程度と訓練データに対して精度が低く、過学習が起きていることが確認できる。

次にL2正則化を導入する。

from common import optimizer

(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True)

print("データ読み込み完了")

# 過学習を再現するために、学習データを削減
x_train = x_train[:300]
d_train = d_train[:300]


network = MultiLayerNet(input_size=784, hidden_size_list=[100, 100, 100, 100, 100, 100], output_size=10)


iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100
learning_rate=0.01

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10
hidden_layer_num = network.hidden_layer_num

# 正則化強度設定 ======================================
weight_decay_lambda = 0.1
# =================================================

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    grad = network.gradient(x_batch, d_batch)
    weight_decay = 0
    
    for idx in range(1, hidden_layer_num+1):
        grad['W' + str(idx)] = network.layers['Affine' + str(idx)].dW + weight_decay_lambda * network.params['W' + str(idx)]
        grad['b' + str(idx)] = network.layers['Affine' + str(idx)].db
        network.params['W' + str(idx)] -= learning_rate * grad['W' + str(idx)]
        network.params['b' + str(idx)] -= learning_rate * grad['b' + str(idx)]        
        weight_decay += 0.5 * weight_decay_lambda * np.sqrt(np.sum(network.params['W' + str(idx)] ** 2))

    loss = network.loss(x_batch, d_batch) + weight_decay
    train_loss_list.append(loss)        
        
    if (i+1) % plot_interval == 0:
        accr_train = network.accuracy(x_train, d_train)
        accr_test = network.accuracy(x_test, d_test)
        accuracies_train.append(accr_train)
        accuracies_test.append(accr_test)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))               


lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

f:id:tibet:20211120205824p:plain
訓練データのaccuracyは少し下がっている(0.94程度)がテストデータに対しては0.74程度と対策なしに比較して精度の差が若干小さくなっている。

L1正則化を導入する。

(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True)

print("データ読み込み完了")

# 過学習を再現するために、学習データを削減
x_train = x_train[:300]
d_train = d_train[:300]

network = MultiLayerNet(input_size=784, hidden_size_list=[100, 100, 100, 100, 100, 100], output_size=10)


iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100
learning_rate=0.1

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10
hidden_layer_num = network.hidden_layer_num

# 正則化強度設定 ======================================
weight_decay_lambda = 0.005
# =================================================

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    grad = network.gradient(x_batch, d_batch)
    weight_decay = 0
    
    for idx in range(1, hidden_layer_num+1):
        grad['W' + str(idx)] = network.layers['Affine' + str(idx)].dW + weight_decay_lambda * np.sign(network.params['W' + str(idx)])
        grad['b' + str(idx)] = network.layers['Affine' + str(idx)].db
        network.params['W' + str(idx)] -= learning_rate * grad['W' + str(idx)]
        network.params['b' + str(idx)] -= learning_rate * grad['b' + str(idx)]        
        weight_decay += weight_decay_lambda * np.sum(np.abs(network.params['W' + str(idx)]))

    loss = network.loss(x_batch, d_batch) + weight_decay
    train_loss_list.append(loss)        
        
    if (i+1) % plot_interval == 0:
        accr_train = network.accuracy(x_train, d_train)
        accr_test = network.accuracy(x_test, d_test)
        accuracies_train.append(accr_train)
        accuracies_test.append(accr_test)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))               
                
lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

f:id:tibet:20211120210207p:plain
一部の重みが0になっているため、学習の進み方が特徴的になっている。

ドロップアウトを導入する。

class Dropout:
    def __init__(self, dropout_ratio=0.5):
        self.dropout_ratio = dropout_ratio
        self.mask = None

    def forward(self, x, train_flg=True):
        if train_flg:
            self.mask = np.random.rand(*x.shape) > self.dropout_ratio
            return x * self.mask
        else:
            return x * (1.0 - self.dropout_ratio)

    def backward(self, dout):
        return dout * self.mask
from common import optimizer
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True)

print("データ読み込み完了")

# 過学習を再現するために、学習データを削減
x_train = x_train[:300]
d_train = d_train[:300]

# ドロップアウト設定 ======================================
use_dropout = True
dropout_ratio = 0.15
# ====================================================

network = MultiLayerNet(input_size=784, hidden_size_list=[100, 100, 100, 100, 100, 100], output_size=10,
                        weight_decay_lambda=weight_decay_lambda, use_dropout = use_dropout, dropout_ratio = dropout_ratio)
optimizer = optimizer.SGD(learning_rate=0.01)
# optimizer = optimizer.Momentum(learning_rate=0.01, momentum=0.9)
# optimizer = optimizer.AdaGrad(learning_rate=0.01)
# optimizer = optimizer.Adam()

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10


for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    grad = network.gradient(x_batch, d_batch)
    optimizer.update(network.params, grad)

    loss = network.loss(x_batch, d_batch)
    train_loss_list.append(loss)    
    
    if (i+1) % plot_interval == 0:
        accr_train = network.accuracy(x_train, d_train)
        accr_test = network.accuracy(x_test, d_test)
        accuracies_train.append(accr_train)
        accuracies_test.append(accr_test)

        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))        
        
lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

f:id:tibet:20211120210420p:plain
ドロップアウトをしたことで、学習の進み方は遅くなっている。

最後にドロップアウトとL1正則化の組み合わせを評価する。

from common import optimizer
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True)

print("データ読み込み完了")

# 過学習を再現するために、学習データを削減
x_train = x_train[:300]
d_train = d_train[:300]

# ドロップアウト設定 ======================================
use_dropout = True
dropout_ratio = 0.08
# ====================================================

network = MultiLayerNet(input_size=784, hidden_size_list=[100, 100, 100, 100, 100, 100], output_size=10,
                        use_dropout = use_dropout, dropout_ratio = dropout_ratio)

iters_num = 1000
train_size = x_train.shape[0]
batch_size = 100
learning_rate=0.01

train_loss_list = []
accuracies_train = []
accuracies_test = []
hidden_layer_num = network.hidden_layer_num

plot_interval=10

# 正則化強度設定 ======================================
weight_decay_lambda=0.004
# =================================================

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    grad = network.gradient(x_batch, d_batch)
    weight_decay = 0
    
    for idx in range(1, hidden_layer_num+1):
        grad['W' + str(idx)] = network.layers['Affine' + str(idx)].dW + weight_decay_lambda * np.sign(network.params['W' + str(idx)])
        grad['b' + str(idx)] = network.layers['Affine' + str(idx)].db
        network.params['W' + str(idx)] -= learning_rate * grad['W' + str(idx)]
        network.params['b' + str(idx)] -= learning_rate * grad['b' + str(idx)]        
        weight_decay += weight_decay_lambda * np.sum(np.abs(network.params['W' + str(idx)]))

    loss = network.loss(x_batch, d_batch) + weight_decay
    train_loss_list.append(loss)        
        
    if (i+1) % plot_interval == 0:
        accr_train = network.accuracy(x_train, d_train)
        accr_test = network.accuracy(x_test, d_test)
        accuracies_train.append(accr_train)
        accuracies_test.append(accr_test)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))               
        
lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

f:id:tibet:20211120210747p:plain
ドロップアウト単体に対して学習の進み方は早くなり、テストデータに対しても精度の飽和が見られず、過学習しにくいことが確認できる。

ラビットチャレンジレポート 深層学習DAY2その1

Section1 勾配消失問題

誤差逆伝搬法

 誤差逆伝搬法とは、計算結果から微分を逆算することで、不要な再帰的計算を避けて微分を算出する手法である。

  • 確認テスト

 連鎖律の原理を使い、下記の式のdz/dxを求めよ
 z=t^2
 t=x+y
[解答]
 \dfrac{dz}{dt}=2t
 \dfrac{dt}{dx}=1
 \dfrac{dz}{dx}=\dfrac{dz}{dt}\dfrac{dt}{dx}=2t\cdot 1=2t=2)x+y)
 

勾配消失問題とは

誤差逆伝搬法が下位層に進んでいくに連れて、勾配がどんどん緩やかになっていく問題。
そのため、勾配降下法による更新では、階層のパラメータはほとんど変わらず、訓練は最適値に収束しなくなる。
 シグモイド関数を使用した場合、値の絶対値が大きいと勾配がほぼ0になるため、勾配が消失する。
 また、微分の連鎖律に起因する問題もある。
例として、出力から3層下位の重みwの連鎖律の式を下記に示す。
 \dfrac{\partial E}{\partial w^{(3)}}=\dfrac{\partial E}{\partial y}\dfrac{\partial y}{\partial u^{(3)}}\dfrac{\partial u^{(3)}}{\partial z^{(2)}}\dfrac{\partial z^{(2)}}{\partial u^{(2)}}\dfrac{\partial u^{(2)}}{\partial z^{(2)}}\dfrac{\partial z^{(2)}}{\partial u^{(1)}}\dfrac{\partial u^{(1)}}{\partial w^{(1)}}
この式を見ると、活性化関数zの入力uに対する微分 \dfrac{\partial z}{\partial u}が2回出てくる。活性化関数としてシグモイド関数を使用した場合、シグモイド関数微分値は最大で0.25程度である。これを二回乗算すると、0.0625となってしまい、勾配がほぼ消失してしまうように見える。
これが勾配消失問題の例である。

  • 確認テスト

シグモイド関数微分した時、入力値が0の時に最大値をとる。その値として正しいものを選択肢から選べ。
(1)0.15
(2)0.25
(3)0.35
(4)0.45

[解答]
(2)0.25

勾配消失問題の解決法として知られているものは次の3つである。

  • 活性化関数の選択
  • 重みの初期値設定
  • バッチ正規化

活性化関数

勾配消失問題回避に効果のある活性化関数として広く用いられているものに、ReLu関数がある。
 f(x)= {x(x>0), 0(x \leq 0)}
ReLu関数は、x>0では微分すると1のため、何度乗算しても減衰しないという特徴がある。
また、x≦0では0のため、活性ノードが疎になる(スパース化)ことで、計算量の軽減や精度向上につながるという利点がある。

重みの初期設定

ニューラルネットワークにおける重みの初期値は、学習の精度や速度に大きく影響する非常に重要な問題である。
 通常の標準正規分布から生成した初期値では、0や1に偏ってしまったり、標準偏差を非常に狭くしたものでは、0.5付近に偏ってしまったりと、偏りや勾配消失が起こりやすくなってしまう。
 そこで、初期値生成の手法として、「Xavierの初期値」と「Heの初期値」が提案された。

Xavierの初期値

 前層のノード数がn個の時、平均0、標準偏差 \dfrac{1}{\sqrt n} である正規分布から初期値を生成する。主にシグモイド関数を活性化関数としたときに効果を発揮する初期値生成法である。

Heの初期値

 前層のノード数がn個の時、平均0、標準偏差  \dfrac{2}{\sqrt n} である正規分布から初期値を生成する手法。主にReLu関数を活性化関数としたときに効果を発揮する初期値生成法である。

  • 確認テスト

 重みの初期値を0に設定すると、どのような問題が発生するか。簡潔に説明せよ。
[解答]
 学習が出来なくなる。重みを均質にすると、学習時に同じように更新されてしまい、特徴の表現が出来なくなる。

バッチ正規化

 勾配降下法の一つで、データを小さな塊に分けて学習をするミニバッチ学習がある。このミニバッチ学習をそのままのデータで行うと、バッチごとにデータの偏りが出て、勾配消失問題が起こる等学習が上手く出来ない場合がある。それを抑制するために各バッチのデータの分布を正規化することをバッチ正規化という。

  • 確認問題

 一般に考えられるバッチ正規化の効果を2点あげよ。
[解答]
-学習が安定化する
-パラメータのスケールや初期値の影響が小さくなるので、高い学習率を設定することが可能になり、学習スピードが速くなる。

  • バッチ正規化の手順

1 データの平均値をとる
 μ_t=\dfrac{1}{N_t}\displaystyle \sum_{i=1}^{N_t}x_ni
2 データの分散をとる
 σ_t^2=\dfrac{1}{N_t}\displaystyle \sum_{i=1}^{N_t}(x_{ni}-μ_t)^2
3 平均ゼロにし、標準偏差にθを加えて分布をコントロール
 \hat{x}_{ni}=\dfrac{x_{ni}-μ_t}{\sqrt{σ_t^2+θ}}
4 変倍とバイアスをつけて変換
 y_ni=ɤx_ni+β

実装演習

シグモイド関数の活性化関数による勾配消失と、ReLu関数や、重み初期値の設定方法による勾配消失問題の回避効果をニューラルネットワークモデルを用いて検証する。
まず、ニューラルネットワークを定義する。

import numpy as np
from common import layers
from collections import OrderedDict
from common import functions
from data.mnist import load_mnist
import matplotlib.pyplot as plt

class MultiLayerNet:
    '''
    input_size: 入力層のノード数
    hidden_size_list: 隠れ層のノード数のリスト
    output_size: 出力層のノード数
    activation: 活性化関数
    weight_init_std: 重みの初期化方法
    '''
    def __init__(self, input_size, hidden_size_list, output_size, activation='relu', weight_init_std='relu'):
        self.input_size = input_size
        self.output_size = output_size
        self.hidden_size_list = hidden_size_list
        self.hidden_layer_num = len(hidden_size_list)
        self.params = {}

        # 重みの初期化
        self.__init_weight(weight_init_std)

        # レイヤの生成, sigmoidとreluのみ扱う
        activation_layer = {'sigmoid': layers.Sigmoid, 'relu': layers.Relu}
        self.layers = OrderedDict() # 追加した順番に格納
        for idx in range(1, self.hidden_layer_num+1):
            self.layers['Affine' + str(idx)] = layers.Affine(self.params['W' + str(idx)], self.params['b' + str(idx)])
            self.layers['Activation_function' + str(idx)] = activation_layer[activation]()

        idx = self.hidden_layer_num + 1
        self.layers['Affine' + str(idx)] = layers.Affine(self.params['W' + str(idx)], self.params['b' + str(idx)])

        self.last_layer = layers.SoftmaxWithLoss()

    def __init_weight(self, weight_init_std):
        all_size_list = [self.input_size] + self.hidden_size_list + [self.output_size]
        for idx in range(1, len(all_size_list)):
            scale = weight_init_std
            if str(weight_init_std).lower() in ('relu', 'he'):
                scale = np.sqrt(2.0 / all_size_list[idx - 1])
            elif str(weight_init_std).lower() in ('sigmoid', 'xavier'):
                scale = np.sqrt(1.0 / all_size_list[idx - 1])

            self.params['W' + str(idx)] = scale * np.random.randn(all_size_list[idx-1], all_size_list[idx])
            self.params['b' + str(idx)] = np.zeros(all_size_list[idx])

    def predict(self, x):
        for layer in self.layers.values():
            x = layer.forward(x)

        return x

    def loss(self, x, d):
        y = self.predict(x)

        weight_decay = 0
        for idx in range(1, self.hidden_layer_num + 2):
            W = self.params['W' + str(idx)]

        return self.last_layer.forward(y, d) + weight_decay

    def accuracy(self, x, d):
        y = self.predict(x)
        y = np.argmax(y, axis=1)
        if d.ndim != 1 : d = np.argmax(d, axis=1)

        accuracy = np.sum(y == d) / float(x.shape[0])
        return accuracy

    def gradient(self, x, d):
        # forward
        self.loss(x, d)

        # backward
        dout = 1
        dout = self.last_layer.backward(dout)

        layers = list(self.layers.values())
        layers.reverse()
        for layer in layers:
            dout = layer.backward(dout)

        # 設定
        grad = {}
        for idx in range(1, self.hidden_layer_num+2):
            grad['W' + str(idx)] = self.layers['Affine' + str(idx)].dW
            grad['b' + str(idx)] = self.layers['Affine' + str(idx)].db

        return grad

このモデルを用いて、まず活性化関数にシグモイドを用い、重みの初期値を、標準偏差0.01、平均0の正規分布で発生させたものを用いて学習効果を確認する。

# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True, one_hot_label=True)

print("データ読み込み完了")

network = MultiLayerNet(input_size=784, hidden_size_list=[40, 20], output_size=10, activation='sigmoid', weight_init_std=0.01)

iters_num = 2000
train_size = x_train.shape[0]
batch_size = 100
learning_rate = 0.1

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    # 勾配
    grad = network.gradient(x_batch, d_batch)
    
    for key in ('W1', 'W2', 'W3', 'b1', 'b2', 'b3'):
        network.params[key] -= learning_rate * grad[key]
    
    loss = network.loss(x_batch, d_batch)
    train_loss_list.append(loss)
    
    if (i + 1) % plot_interval == 0:
        accr_test = network.accuracy(x_test, d_test)
        accuracies_test.append(accr_test)        
        accr_train = network.accuracy(x_batch, d_batch)
        accuracies_train.append(accr_train)

        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))
        

lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果のグラフ表示
f:id:tibet:20211119175631p:plain
Accuracy(予測精度)が学習の回数を重ねても向上せず、学習が進んでいないことがわかる。

次に、活性化関数にReLu関数を使用し、重みの初期値は前回と同じく標準偏差が0.01、平均0の正規分布に従って設定する。

# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True, one_hot_label=True)

print("データ読み込み完了")

network = MultiLayerNet(input_size=784, hidden_size_list=[40, 20], output_size=10, activation='relu', weight_init_std=0.01)

iters_num = 2000
train_size = x_train.shape[0]
batch_size = 100
learning_rate = 0.1

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    # 勾配
    grad = network.gradient(x_batch, d_batch)
    
    for key in ('W1', 'W2', 'W3', 'b1', 'b2', 'b3'):
        network.params[key] -= learning_rate * grad[key]
    
    loss = network.loss(x_batch, d_batch)
    train_loss_list.append(loss)
    
    if (i + 1) % plot_interval == 0:
        accr_test = network.accuracy(x_test, d_test)
        accuracies_test.append(accr_test)        
        accr_train = network.accuracy(x_batch, d_batch)
        accuracies_train.append(accr_train)

        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))
        
        
lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

結果の図示
f:id:tibet:20211119194844p:plain
500回前後から精度が向上し、1000回程度まで急速に学習が進み、その後緩やかに向上している。
勾配消失問題が緩和していることがわかる。

次に、活性化関数をシグモイド関数に戻し、重みの初期値をXavierを用いて設定する。

# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True, one_hot_label=True)

print("データ読み込み完了")

network = MultiLayerNet(input_size=784, hidden_size_list=[40, 20], output_size=10, activation='sigmoid', weight_init_std='Xavier')

iters_num = 2000
train_size = x_train.shape[0]
batch_size = 100
learning_rate = 0.1

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    # 勾配
    grad = network.gradient(x_batch, d_batch)
    
    for key in ('W1', 'W2', 'W3', 'b1', 'b2', 'b3'):
        network.params[key] -= learning_rate * grad[key]
    
    loss = network.loss(x_batch, d_batch)
    train_loss_list.append(loss)
    
    if (i + 1) % plot_interval == 0:
        accr_test = network.accuracy(x_test, d_test)
        accuracies_test.append(accr_test)        
        accr_train = network.accuracy(x_batch, d_batch)
        accuracies_train.append(accr_train)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))
        

lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

f:id:tibet:20211119195909p:plain
シグモイド関数を使っているにも関わらず、学習回数に応じて精度は向上しており、勾配消失が緩和していることがわかる。

最後に、活性化関数をReLu関数にして、重みの初期値をHeで設定する。

# データの読み込み
(x_train, d_train), (x_test, d_test) = load_mnist(normalize=True, one_hot_label=True)

print("データ読み込み完了")

network = MultiLayerNet(input_size=784, hidden_size_list=[40, 20], output_size=10, activation='relu', weight_init_std='He')

iters_num = 2000
train_size = x_train.shape[0]
batch_size = 100
learning_rate = 0.1

train_loss_list = []
accuracies_train = []
accuracies_test = []

plot_interval=10

for i in range(iters_num):
    batch_mask = np.random.choice(train_size, batch_size)
    x_batch = x_train[batch_mask]
    d_batch = d_train[batch_mask]

    # 勾配
    grad = network.gradient(x_batch, d_batch)
    
    for key in ('W1', 'W2', 'W3', 'b1', 'b2', 'b3'):
        network.params[key] -= learning_rate * grad[key]
    
    loss = network.loss(x_batch, d_batch)
    train_loss_list.append(loss)
    
    if (i + 1) % plot_interval == 0:
        accr_test = network.accuracy(x_test, d_test)
        accuracies_test.append(accr_test)        
        accr_train = network.accuracy(x_batch, d_batch)
        accuracies_train.append(accr_train)
        
        print('Generation: ' + str(i+1) + '. 正答率(トレーニング) = ' + str(accr_train))
        print('                : ' + str(i+1) + '. 正答率(テスト) = ' + str(accr_test))
        

lists = range(0, iters_num, plot_interval)
plt.plot(lists, accuracies_train, label="training set")
plt.plot(lists, accuracies_test,  label="test set")
plt.legend(loc="lower right")
plt.title("accuracy")
plt.xlabel("count")
plt.ylabel("accuracy")
plt.ylim(0, 1.0)
# グラフの表示
plt.show()

f:id:tibet:20211119200237p:plain
学習の初期からaccuracyが向上しており、Heを使わない場合に対して、学習のスピードが著しく向上していることがわかる。

  • 例題チャレンジ

以下は特徴データdata_x、ラベルデータdata_tに対してミニバッチ学習を行うプログラムである。
下記の(き)に当てはまるコードを答えよ。
f:id:tibet:20211119200438p:plain
[解答]

data_x[i:i_end],data_t[i:i_end]

ラビットチャレンジレポート 深層学習DAY1

Section 1 入力層~活性層

  • 確認テスト

f:id:tibet:20211118121514p:plain
[答え]
ディープラーニングは明示的なプログラムの代わりに、ニューラルネットワークを用いて各種パラメータを最適化することによって、入力値から目的とする出力地に変換する数学モデルを構築すること
最適化の最終目的:③重み[W]、④バイアス[b]


f:id:tibet:20211118124810p:plain
入力層は、図のように入力  x_i及びバイアスbが入力される層である。各入力層から中間層に至るネットワークには、各入力に重みw_iが乗算される。中間層ではそれらを線形結合して総入力uが作られる。
中間層の出力zは、活性化関数f()を用いて、下記のように書ける。
 z=f(u)
zは、次の中間層の入力となる。
一般には、このように中間層(隠れ層ともいう)を複数積層することで、ニューラルネットワークが形作られる。
数式化すると、重みベクトルW, 入力ベクトルxは次のように書ける。
 W=[w_1....w_I^T]
 x=[x_1....x_I^T]
 u=w_1x_1+w_2x_2+..+b=Wx+b

  • 確認テスト:

この図式に動物分類の実例を入れてみよう
[答え]
 x_1=10 体長
 x_2=300 体重
 x_3=300 髭の本数
 x_4=15 毛の平均長
 x_5=50 耳の大きさ

  • 確認テスト:

この式をパイソンで書け
 u=w_1x_1+w_2x_2+..+b=Wx+b
[答え]

u=np.dot(x,W)+b
  • 確認テスト:

1-1のファイルから中間層の出力を定義しているソースを抜き出せ
[答え]

    # 2層の総出力
    z2 = functions.relu(u2)


実装演習

numpyと関数のインポート及び表示用関数定義

import numpy as np
import functions

def print_vec(text, vec):
    print("*** " + text + " ***")
    print(vec)
    print("")

単相、単ユニットの順伝搬による入力から中間層出力までのモデル

# 重み
W = np.array([[0.1], [0.2]])
print_vec("重み", W)

# バイアス
b = np.array([0.1, 0.2, 0.3])
print_vec("バイアス", b)

# 入力値
x = np.array([1.0, 5.0, 2.0, -1.0])
print_vec("入力", x)

#  総入力
u = np.dot(x, W) + b
print_vec("総入力", u)

# 中間層出力
z = functions.sigmoid(u)
print_vec("中間層出力", z)

出力結果

*** 重み ***
[[0.1 0.2 0.3]
 [0.2 0.3 0.4]
 [0.3 0.4 0.5]
 [0.4 0.5 0.6]]

*** バイアス ***
[0.1 0.2 0.3]

*** 入力 ***
[ 1.  5.  2. -1.]

*** 総入力 ***
[1.4 2.2 3. ]

*** 中間層出力 ***
[0.80218389 0.90024951 0.95257413]

Section2 活性化関数

ニューラルネットにおいて、次の層への出力の大きさを決める非線形の関数
入力値の値によって、次の層への信号のON/OFFや強弱を定める働きを持つ

  • 確認テスト

線形と非線形の違いを図に書いて簡易に説明せよ。
[答え]
f:id:tibet:20211118133641p:plain
線形関数は、左図のように入出力の関係が直線のイメージで表され、加法性と斉次性を持つもの。

  • 加法性:  f(x+y)=f(x)+f(y)
  • 斉次性:  f(kx)=kf(x)

非線形関数は、加法性と斉次性を持たない関数。左図のように入出力の関係が直線にならないイメージ。

活性化関数の種類

中間層用の活性化関数としては、次がよく知られている。

出力層用の活性化関数としては、下記がよく知られている。

活性化関数:ステップ関数

数式
 f(x)= 1(x \geq 0) or 0 (x<0)
閾値を超えたら発火する関数であり、出力は常に0か1。パーセプトロンで利用された関数。
0-1の間を表現できず、線形分離可能なモノしか学習できなかった。

活性化関数:シグモイド関数

数式:
 \dfrac{1}{1+e^{-u}}
0-1の間を緩やかに変化する関数で、ステップ関数に比べて信号の強弱を伝えられるようになり、ニューラルネットワーク普及のきっかけとなった。
一方で、逆誤差伝搬法を用いると、微分値が小さいために微分値を何度も掛け合わせると消失してしまい、深いニューラルネットワークでは勾配消失問題を引き起こすという課題があった。

活性化関数:ReLu関数

数式:
 f(x)= x (x>0) or 0 (x \leq 0)
現在最も使われている関数。x>0では、微分値が常に1となるため、勾配消失問題が起きにくい。
また、 x\leq 0では常に0のため、層の活性化がスパースとなり、計算量が少なく、精度も向上しやすい。

  • 確認テスト

配布されたソースコードより、該当する場所を抜き出せ。
f:id:tibet:20211118141401p:plain
[答え]

z1=functions.sigmoid(u)

実装演習

3層複数ユニットのネットワーク作成

# ウェイトとバイアスを設定
# ネートワークを作成
def init_network():
    print("##### ネットワークの初期化 #####")
    network = {}
    
    network['W1'] = np.array([
        [0.1, 0.3, 0.5],
        [0.2, 0.4, 0.6]
    ])
    network['W2'] = np.array([
        [0.1, 0.4],
        [0.2, 0.5],
        [0.3, 0.6]
    ])
    network['W3'] = np.array([
        [0.1, 0.3],
        [0.2, 0.4]
    ])
    network['b1'] = np.array([0.1, 0.2, 0.3])
    network['b2'] = np.array([0.1, 0.2])
    network['b3'] = np.array([1, 2])

    print_vec("重み1", network['W1'] )
    print_vec("重み2", network['W2'] )
    print_vec("重み3", network['W3'] )
    print_vec("バイアス1", network['b1'] )
    print_vec("バイアス2", network['b2'] )
    print_vec("バイアス3", network['b3'] )

    return network

# プロセスを作成
# x:入力値
def forward(network, x):
    
    print("##### 順伝播開始 #####")

    W1, W2, W3 = network['W1'], network['W2'], network['W3']
    b1, b2, b3 = network['b1'], network['b2'], network['b3']
    
    # 1層の総入力
    u1 = np.dot(x, W1) + b1
    
    # 1層の総出力
    z1 = functions.relu(u1)
    
    # 2層の総入力
    u2 = np.dot(z1, W2) + b2
    
    # 2層の総出力
    z2 = functions.relu(u2)

    # 出力層の総入力
    u3 = np.dot(z2, W3) + b3
    
    # 出力層の総出力
    y = u3
    
    print_vec("総入力1", u1)
    print_vec("中間層出力1", z1)
    print_vec("総入力2", u2)
    print_vec("出力1", z1)
    print("出力合計: " + str(np.sum(z1)))

    return y, z1, z2

# 入力値
x = np.array([1., 2.])
print_vec("入力", x)

# ネットワークの初期化
network =  init_network()

y, z1, z2 = forward(network, x)

結果

*** 入力 ***
[1. 2.]

##### ネットワークの初期化 #####
*** 重み1 ***
[[0.1 0.3 0.5]
 [0.2 0.4 0.6]]

*** 重み2 ***
[[0.1 0.4]
 [0.2 0.5]
 [0.3 0.6]]

*** 重み3 ***
[[0.1 0.3]
 [0.2 0.4]]

*** バイアス1 ***
[0.1 0.2 0.3]

*** バイアス2 ***
[0.1 0.2]

*** バイアス3 ***
[1 2]

##### 順伝播開始 #####
*** 総入力1 ***
[0.6 1.3 2. ]

*** 中間層出力1 ***
[0.6 1.3 2. ]

*** 総入力2 ***
[1.02 2.29]

*** 出力1 ***
[0.6 1.3 2. ]

出力合計: 3.9

Section3 出力層

誤差関数

 誤差関数は、ネットワークの出力値と、あらかじめ訓練データとして用意しておいた正解との誤差を評価する関数である。
特に知られているものは、MSE(平均二乗誤差)である。
 E_n(w)=\dfrac{1}{2}\displaystyle\sum_{J=1}^{l}(y_j-d_j)^2=\dfrac{1}{2}||(y-d)||^2
他にも、平均絶対誤差(MAE)や、平均二乗誤差の平方根(RMSE)、平均二乗対数誤差などが知られている。

-確認テスト
・なぜ引き算でなく二乗するかを述べよ。
[答え]
誤差の和をとって評価するため、正の値と負の値が存在すると打ち消し合って正確に評価が出来ない。二乗すれば常に正だから。
・MSEの1/2はどういう意味を持つか述べよ
[答え]
誤差逆伝搬法の計算時に微分をとるので、計算がしやすいから。

出力層の種類

  • 回帰
    • 恒等関数  f(u)=u

入力がそのまま出力になる関数。誤差関数としては二乗誤差を使う

入力が0~1の間をとる出力に変換される関数。誤差関数としては交差エントロピーを用いる

  • 多クラス分類
    • ソフトマックス関数  f(i,u)=\dfrac{e^{ui_}}{\sum_{K=1}^{K}e^{u_k}

多値の確率が出力される関数。誤差関数は交差エントロピーを用いる。

  • 確認テスト

 f:id:tibet:20211118193407p:plain
[答え]

def softmax(x):
#この関数全体が①に当たる。入力は
 if x.dim ==2 # 入力が2次元の場合
  x=x.T #入力を転置する
  x=x-np.max(x,axis=0) #オーバーフロー対策
  y=np.exp(x)/np.sum(np.exp(x), axis=0) #分子が②、分母が③
 return y.T # yを転置して出力
# 入力が2次元以上の場合
 x=x-np.max(x) #オーバーフロー対策
 return np.exp(x)/np.sum(np.exp(x)) #分子が②、分母が③

  • 確認テスト

交差エントロピーの下記の式で①~②に該当するソースコードを示し、1行ずつ処理の説明をせよ
f:id:tibet:20211118202618p:plain
[答え]

def cross_entropy_error(d, y): #全体の関数は①
 if y.ndim==1: #yが1次元の場合
   d=d.reshape(1, d.size)
   y=y.reshape(1, y.size)
 #教師データがone-hot-vectorの場合、正解ラベルのインデックスに変換
   if d.size==y.size: #dとyのサイズが等しい場合
   d=d.argmax(axis=1) #dにdの最大値のインデックスを代入
   batch_size=y.shape[0] #バッチサイズをyの行数にする(通常1)
 return -np.sum(np.log(y[np.arange(batch_size),d]+1e-7))/batch_size #こちらの数式が②の部分

実装演習

  • 多クラス分類(2-3-4ネットワーク)

入力層サイズ:3
隠れ層サイズ:50
出力層サイズ:6
出力層活性化関数:softmax

def init_network():
    print("##### ネットワークの初期化 #####")

    network = {}
    
    input_layer_size = 3
    hidden_layer_size=50
    output_layer_size = 6
   
    network['W1'] = np.random.rand(input_layer_size, hidden_layer_size)
    network['W2'] = np.random.rand(hidden_layer_size,output_layer_size)

    network['b1'] =  np.random.rand(hidden_layer_size)
    network['b2'] =  np.random.rand(output_layer_size)
    
    print_vec("重み1", network['W1'] )
    print_vec("重み2", network['W2'] )
    print_vec("バイアス1", network['b1'] )
    print_vec("バイアス2", network['b2'] )

    return network

# プロセスを作成
# x:入力値
def forward(network, x):
    
    print("##### 順伝播開始 #####")
    W1, W2 = network['W1'], network['W2']
    b1, b2 = network['b1'], network['b2']
    
    # 1層の総入力
    u1 = np.dot(x, W1) + b1

    # 1層の総出力
    z1 = functions.relu(u1)

    # 2層の総入力
    u2 = np.dot(z1, W2) + b2
    
    # 出力値
    y = functions.softmax(u2)
    
    print_vec("総入力1", u1)
    print_vec("中間層出力1", z1)
    print_vec("総入力2", u2)
    print_vec("出力1", y)
    print("出力合計: " + str(np.sum(y)))
        
    return y, z1

## 事前データ
# 入力値
x = np.array([1., 2.,  3.])

# 目標出力
d = np.array([0, 0, 0, 1, 0, 0])

# ネットワークの初期化
network =  init_network()

# 出力
y, z1 = forward(network, x)

# 誤差
loss = functions.cross_entropy_error(d, y)

## 表示
print("\n##### 結果表示 #####")
print_vec("出力", y)
print_vec("訓練データ", d)
print_vec("交差エントロピー誤差",  loss)

結果(一部)

##### 結果表示 #####
*** 出力 ***
[6.10163245e-08 1.58053510e-08 9.99981036e-01 8.10503745e-16
 7.56282983e-06 1.13244940e-05]
shape: (6,)

*** 訓練データ ***
[0 0 0 1 0 0]
shape: (6,)

*** 交差エントロピー誤差 ***
16.118095642853284
shape: ()

Section4 勾配降下法

勾配降下法

 誤差関数E(w)を最小化するパラメータwを発見的に探索するために活用される。
 誤差関数のパラメータによる偏微分に学習率をかけたものを重みから差し引き、新たな重みとして更新する。
 w^{(t+1)}=w~{(t)}-ε∇E
 ∇E=\dfrac{\partial E}{\partial w}=\displaystyle \left( \dfrac{\partial E}{\partial w_1}, ... , \dfrac{\partial E}{\partial w_M}\right)

  • 確認問題

下記に該当するソースコードを探してみる。
 w^{(t+1)}=w~{(t)}-ε∇E
 ∇E=\dfrac{\partial E}{\partial w}=\displaystyle \left( \dfrac{\partial E}{\partial w_1}, ... , \dfrac{\partial E}{\partial w_M}\right)
[回答]

network[key]-=learning_rate*grad[key]
grad=backward(x,d,z1,y)

勾配降下法は、学習率の値によって学習率の値が大きく異なる

  • 学習率が小さすぎる場合:収束までの時間がかかり、かつ局所解に陥る場合がある。
  • 学習率が大きすぎる場合:解が収束しない場合がある。

 勾配降下法の学習率の決定、収束性向上のためのアルゴリズムは様々な手法が考え出されている。

確率的勾配降下法(SGD)

全サンプルの平均誤差を使用する勾配降下法に対し、確率的勾配降下法は、ランダムに抽出したサンプルの誤差を使用する。
メリット

  • データが冗長な場合の計算コストの低減
  • 局所極小解に収束するリスクの軽減
  • オンライン学習が出来る

  • 確認テスト

オンライン学習とは何か2行でまとめよ
[回答]
学習データが入ってくるたびにその都度、新たなデータのみを使用して学習し、パラメータを随時更新する学習法。

ミニバッチ勾配降下法

勾配降下法と確率的勾配降下法の中間のような手法で、ランダムに分割したデータの集合(ミニバッチ)に属するサンプルの平均誤差を使用して学習する。
メリット
確率的勾配降下法のメリットを損なわず、計算機の計算資源を有効利用できる。
 →CPUを利用したスレッド並列化やGPUを利用したSIMD並列化が可能

  • 確認テスト

 w^{(t+1)}=w^{(t)}-ε∇E_t
この数式の意味を図に書いて説明せよ
[解答]
f:id:tibet:20211118214200p:plain
図のようにエポックが1増えるごとに、-ε∇E分だけ、重み w_tが更新されることを意味している。

実装演習

確率的勾配降下法の演習
3層のニューラルネットワークで、epoch数が増えるごとに損失関数の値がどのように変化するかを図示する。

# サンプルとする関数
#yの値を予想するAI

def f(x):
    y = 3 * x[0] + 2 * x[1]
    return y

# 初期設定
def init_network():
    network = {}
    nodesNum = 10
    network['W1'] = np.random.randn(2, nodesNum)
    network['W2'] = np.random.randn(nodesNum)
    network['b1'] = np.random.randn(nodesNum)
    network['b2'] = np.random.randn()
    return network

# 順伝播
def forward(network, x):
    W1, W2 = network['W1'], network['W2']
    b1, b2 = network['b1'], network['b2']
    u1 = np.dot(x, W1) + b1
    z1 = functions.relu(u1)
    
    u2 = np.dot(z1, W2) + b2
    y = u2    
    return z1, y

# 誤差逆伝播
def backward(x, d, z1, y):
    # print("\n##### 誤差逆伝播開始 #####")    

    grad = {}
    
    W1, W2 = network['W1'], network['W2']
    b1, b2 = network['b1'], network['b2']

    # 出力層でのデルタ
    delta2 = functions.d_mean_squared_error(d, y)
    # b2の勾配
    grad['b2'] = np.sum(delta2, axis=0)
    # W2の勾配
    grad['W2'] = np.dot(z1.T, delta2)

    delta1 = np.dot(delta2, W2.T) * functions.d_sigmoid(z1)

    delta1 = delta1[np.newaxis, :]
    # b1の勾配
    grad['b1'] = np.sum(delta1, axis=0)
    x = x[np.newaxis, :]
    # W1の勾配
    grad['W1'] = np.dot(x.T, delta1)
    
    return grad

# サンプルデータを作成
data_sets_size = 100000
data_sets = [0 for i in range(data_sets_size)]

for i in range(data_sets_size):
    data_sets[i] = {}
    # ランダムな値を設定
    data_sets[i]['x'] = np.random.rand(2)
    
    ## 試してみよう_入力値の設定
    # data_sets[i]['x'] = np.random.rand(2) * 10 -5 # -5〜5のランダム数値
    
    # 目標出力を設定
    data_sets[i]['d'] = f(data_sets[i]['x'])
    
losses = []
# 学習率
learning_rate = 0.07

# 抽出数
epoch = 1000

# パラメータの初期化
network = init_network()
# データのランダム抽出
random_datasets = np.random.choice(data_sets, epoch)

# 勾配降下の繰り返し
for dataset in random_datasets:
    x, d = dataset['x'], dataset['d']
    z1, y = forward(network, x)
    grad = backward(x, d, z1, y)
    # パラメータに勾配適用
    for key in ('W1', 'W2', 'b1', 'b2'):
        network[key]  -= learning_rate * grad[key]

    # 誤差
    loss = functions.mean_squared_error(d, y)
    losses.append(loss)

print("##### 結果表示 #####")    
lists = range(epoch)


plt.plot(lists, losses, '.')
# グラフの表示
plt.show()

結果
f:id:tibet:20211118214813p:plain
Epochが増加すると、損失関数の値が低くなって収束していく様子がわかる。

Section5 誤差逆伝搬法

ニューラルネットワークの誤差勾配を計算する時、算出された誤差を、出力層から順に微分し、前の層へと伝搬させることで、最小限の計算で各パラメータでの微分値を解析的に計算する手法。

-確認テスト
誤差逆伝搬法では不要な再帰的処理を避けることが出来る。すでに行った計算結果を保持しているソースコードを抽出せよ。

    #  出力層でのデルタ
    delta2 = functions.d_sigmoid_with_loss(d, y)
    #  中間層でのデルタ
    delta1 = np.dot(delta2, W2.T) * functions.d_relu(z1)

誤差勾配の計算には、順方向の計算結果と微分の連鎖律を利用する。
 E(y)=\dfrac{1}{2}\displaystyle \sum_{j=1}^{J}(y_j-d_j)^2=\dfrac{1}{2}||y-d||^2 :誤差関数

 y=u^{(L)} :出力層の活性化関数
 u^{(l)}=w^{(i)}z^{(l-1)}+b^{(l)} :総入力の計算

連鎖律
 \dfrac{\partial E}{\partial w_{ji}^{(2)}}= \dfrac{\partial E}{\partial y}\dfrac{\partial y}{\partial u}\dfrac{\partial u}{\partial w_{ji}^{(2)}}
第一項の計算
 \dfrac{\partial E(y)}{\partial y}=\dfrac{\partial}{\partial y}\dfrac{1}{2}||y-d||^2=y-d
第二項の計算
 \dfrac{\partial y(u)}{\partial u}=\dfrac{\partial u}{\partial u}=1
第三項の計算
 \dfrac{\partial u(w)}{\partial w_{ji}}=\dfrac{\partial}{\partial w_{ji}}(w^{(l)}Z^{(l-1)}z^{(l-1)}+b^{(l)})=(0....z_j...0)^T

これらの結果より、誤差関数の重みによる微分は下記のように計算できる。
 \dfrac{\partial E}{\partial y} \dfrac{\partial y}{\partial u} \dfrac{\partial u}{\partial w_{ji}^{(2)}}=(y_j-d_j)z_i

-確認テスト
2つの空欄に該当するソースコードを探せ
[解答]
 \dfrac{\partial E}{\partial y}

delta2=functions.d_mean_squared_error(d,y)

 \dfrac{\partial E}{\partial y}\dfrac{\partial y}{\partial u}

delta2=functions.d_mean_squared_error(d,y)

 \dfrac{\partial E}{\partial y} \dfrac{\partial y}{\partial u} \dfrac{\partial u}{\partial w_{ji}^{(2)}}

grad['W2']=np.dot(z1.T,delta2)

実装演習

3層のニューラルネットワークによる誤差逆伝搬法のデモ

def init_network():
    print("##### ネットワークの初期化 #####")

    network = {}
    network['W1'] = np.array([
        [0.1, 0.3, 0.5],
        [0.2, 0.4, 0.6]
    ])

    network['W2'] = np.array([
        [0.1, 0.4],
        [0.2, 0.5],
        [0.3, 0.6]
    ])

    network['b1'] = np.array([0.1, 0.2, 0.3])
    network['b2'] = np.array([0.1, 0.2])
    
    print_vec("重み1", network['W1'])
    print_vec("重み2", network['W2'])
    print_vec("バイアス1", network['b1'])
    print_vec("バイアス2", network['b2'])

    return network

# 順伝播
def forward(network, x):
    print("##### 順伝播開始 #####")

    W1, W2 = network['W1'], network['W2']
    b1, b2 = network['b1'], network['b2']
    
    u1 = np.dot(x, W1) + b1
    z1 = functions.relu(u1)
    u2 = np.dot(z1, W2) + b2
    y = functions.softmax(u2)
    
    print_vec("総入力1", u1)
    print_vec("中間層出力1", z1)
    print_vec("総入力2", u2)
    print_vec("出力1", y)
    print("出力合計: " + str(np.sum(y)))

    return y, z1

# 誤差逆伝播
def backward(x, d, z1, y):
    print("\n##### 誤差逆伝播開始 #####")

    grad = {}

    W1, W2 = network['W1'], network['W2']
    b1, b2 = network['b1'], network['b2']
    #  出力層でのデルタ
    delta2 = functions.d_sigmoid_with_loss(d, y)
    #  b2の勾配
    grad['b2'] = np.sum(delta2, axis=0)
    #  W2の勾配
    grad['W2'] = np.dot(z1.T, delta2)
    #  中間層でのデルタ
    delta1 = np.dot(delta2, W2.T) * functions.d_relu(z1)
    # b1の勾配
    grad['b1'] = np.sum(delta1, axis=0)
    #  W1の勾配
    grad['W1'] = np.dot(x.T, delta1)
        
    print_vec("偏微分_dE/du2", delta2)
    print_vec("偏微分_dE/du2", delta1)

    print_vec("偏微分_重み1", grad["W1"])
    print_vec("偏微分_重み2", grad["W2"])
    print_vec("偏微分_バイアス1", grad["b1"])
    print_vec("偏微分_バイアス2", grad["b2"])

    return grad
    
# 訓練データ
x = np.array([[1.0, 5.0]])
# 目標出力
d = np.array([[0, 1]])
#  学習率
learning_rate = 0.01
network =  init_network()
y, z1 = forward(network, x)

# 誤差
loss = functions.cross_entropy_error(d, y)

grad = backward(x, d, z1, y)
for key in ('W1', 'W2', 'b1', 'b2'):
    network[key]  -= learning_rate * grad[key]

print("##### 結果表示 #####")    


print("##### 更新後パラメータ #####") 
print_vec("重み1", network['W1'])
print_vec("重み2", network['W2'])
print_vec("バイアス1", network['b1'])
print_vec("バイアス2", network['b2'])

結果

##### ネットワークの初期化 #####
*** 重み1 ***
[[0.1 0.3 0.5]
 [0.2 0.4 0.6]]

*** 重み2 ***
[[0.1 0.4]
 [0.2 0.5]
 [0.3 0.6]]

*** バイアス1 ***
[0.1 0.2 0.3]

*** バイアス2 ***
[0.1 0.2]

##### 順伝播開始 #####
*** 総入力1 ***
[[1.2 2.5 3.8]]

*** 中間層出力1 ***
[[1.2 2.5 3.8]]

*** 総入力2 ***
[[1.86 4.21]]

*** 出力1 ***
[[0.08706577 0.91293423]]

出力合計: 1.0

##### 誤差逆伝播開始 #####
*** 偏微分_dE/du2 ***
[[ 0.08706577 -0.08706577]]

*** 偏微分_dE/du2 ***
[[-0.02611973 -0.02611973 -0.02611973]]

*** 偏微分_重み1 ***
[[-0.02611973 -0.02611973 -0.02611973]
 [-0.13059866 -0.13059866 -0.13059866]]

*** 偏微分_重み2 ***
[[ 0.10447893 -0.10447893]
 [ 0.21766443 -0.21766443]
 [ 0.33084994 -0.33084994]]

*** 偏微分_バイアス1 ***
[-0.02611973 -0.02611973 -0.02611973]

*** 偏微分_バイアス2 ***
[ 0.08706577 -0.08706577]

##### 結果表示 #####
##### 更新後パラメータ #####
*** 重み1 ***
[[0.1002612  0.3002612  0.5002612 ]
 [0.20130599 0.40130599 0.60130599]]

*** 重み2 ***
[[0.09895521 0.40104479]
 [0.19782336 0.50217664]
 [0.2966915  0.6033085 ]]

*** バイアス1 ***
[0.1002612 0.2002612 0.3002612]

*** バイアス2 ***
[0.09912934 0.20087066]

ラビットチャレンジレポート 機械学習後半

非線形回帰モデル

複雑な非線形構造を内在する現象に関しては、非線形モデルを適用する。

  • 規定展開法
    • 回帰関数として基底関数と呼ばれる非線形関数とパラメータベクトルの線形結合を使用する
    • よく使われる既定関数は、多項式関数やガウス型既定関数、スプライン関数

多項式型関数:  Φ_j=x^j
ガウス型基底:  Φ_j(x)=exp \left(\dfrac{(x-μ_j)^2}{2h_j}\right)

    • 規定展開法も、最小二乗法や最尤法など線形回帰と同じフレームワークで推定可能

説明変数  x_i=(x_{i1},x_{i2},...,x_{im})\in R^m
非線形関数ベクトル  Φ(x_i)=(Φ_1(x_i),Φ_2(x_i),.....,Φ_k(x_i))^T \in R^m
非線形関数の計画行列  Φ^{(train)}=(Φ_1(x_i),Φ_2(x_i),.....,Φ_k(x_n))^T \in R^{n×k}
最尤法による予測値  \hat{y}=Φ(Φ^{(train)T}Φ^{(train)})^{-1}Φ^{(train)T}y^{(train)}

未学習(underfitting)と過学習(overfitting)

  • 学習データに対して、十分小さな誤差が得られてないモデル→未学習
  • 対策:表現力の高いモデルを利用する
  • 学習データに対しては非常に小さな誤差を得られるが、テストデータに対する誤差が大きいモデル→過学習
    • 対策1:学習データの数を増やす
    • 対策2:不要な既定関数を削除
    • 対策3:正則化法を利用して表現力を抑止
  • 不要な基底関数の削除
    • 解きたい問題に対して多すぎる基底関数を用意してしまうと過学習の問題が起こりやすいため、適切な数の基底関数まで削除する(オッカムのカミソリ)
    • ただし、どの程度の基底関数が適切かの判断は難しい。

(参考)赤池情報基準(AIC)
統計モデルの良さを評価する基準。AICはモデルの複雑さとデータとの適合度のバランスをとるための指標として使われる。AIC最小のモデルを選択すれば、多くの場合良いモデルが作成できるとされる。
 AIC=-2lnL+2k
Lは最大尤度、kは自由パラメータの数

  • 正則化法(罰則化法)
    • 「モデルの複雑さに伴って、その値が大きくなる正則化項」を導入した関数を最小化

 S_ɤ=(y-Φw)^T(y-Φw)+ɤR(w)
ɤR(w):正則化

汎化性能の評価法

    • 学習に使用した入力だけでなく、新たな入力に対する予測性能
  • ホールドアウト法
    • 有限のデータを学習用とテスト用の2つに分割し、予測精度や誤り率を推定。
    • 手軽に出来る点が長所だが、手元に大量にデータが大量にない場合は、学習データと検証データで偏りが起きやすく、良い性能評価を与えないという欠点がある。
  • クロスバリデーション(交差検証)

 学習データをn分割し、1番目を検証用に割り当て、2番目以降を学習データに割り当てて精度評価を行い、次に2番目を検証用に割り当てて精度検証し、次に3番目を検証に割り当て..とn回繰り返し、精度の平均をとる方法。この平均精度をCV値と呼び、データの偏りが起きにくいため、汎化性能指標として広く用いられる。

パラメータチューニング

  • グリッドサーチ

すべてのチューニングパラメータの組み合わせで評価値を算出し、最も良い評価値を持つチューニングパラメータを持つ組み合わせを、「いいモデルのパラメータ」として採用する

非線形回帰モデルのハンズオン

まず環境設定をする。

import numpy as np
import matplotlib.pyplot as plt
import seaborn as sns
%matplotlib inline
#seaborn設定
sns.set()
#背景変更
sns.set_style("darkgrid", {'grid.linestyle': '--'})
#大きさ(スケール変更)
sns.set_context("paper")
n=100
def true_func(x):
    z = 1-48*x+218*x**2-315*x**3+145*x**4
    return z 
def linear_func(x):
    z = x
    return z 

データ生成

 # 真の関数からデータ生成
data = np.random.rand(n).astype(np.float32)
data = np.sort(data)
target = true_func(data)

#  ノイズを加える
noise = 0.5 * np.random.randn(n) 
target = target  + noise

# ノイズ付きデータを描画

plt.scatter(data, target)

plt.title('NonLinear Regression')
plt.legend(loc=2)

f:id:tibet:20211113225328p:plain

まず、線形回帰モデルを使用してみる。

from sklearn.linear_model import LinearRegression
clf = LinearRegression()
data = data.reshape(-1,1)
target = target.reshape(-1,1)
clf.fit(data, target)

p_lin = clf.predict(data)
plt.scatter(data, target, label='data')
plt.plot(data, p_lin, color='darkorange', marker='', linestyle='-', linewidth=1, markersize=6, label='linear regression')
plt.legend()
print(clf.score(data, target))

f:id:tibet:20211113225658p:plain
データに対して、うまく学習出来ていない様子がわかる。

次に、L2制約(Ridge推定)が付いたガウス型基底関数を用いた非線形回帰モデル(KernelRidge)を用いて推定する。

from sklearn.kernel_ridge import KernelRidge
clf = KernelRidge(alpha=0.0002, kernel='rbf')
clf.fit(data, target)
p_kridge = clf.predict(data)
plt.scatter(data, target, color='blue', label='data')
plt.plot(data, p_kridge, color='orange', linestyle='-', linewidth=3, markersize=6, label='kernel ridge')
plt.legend()
#plt.plot(data, p, color='orange', marker='o', linestyle='-', linewidth=1, markersize=6)

f:id:tibet:20211113230325p:plain
上手くフィッティングしている様子がわかる。

ロジスティック回帰

  • 分類問題(クラス分類)
    • ある入力からクラスに分類する問題
    • 2値分類問題(0or1)

説明変数:  x=(x_1,x_2,..,x_m)^T\in R^m
目的変数:  y\in {0,1}

 σ(x)=\dfrac{1}{1+exp(-a)}

 シグモイド関数微分シグモイド関数で表現できる。
 \dfrac{\partial σ(x)}{\partial x}=aσ(x)(1-σ(x))

ロジスティック回帰モデル

 P(Y=1|x)=σ(w_0+\displaystyle  \sum_{k=1}^n (w_kx_k))

最尤推定

  • 尤度関数:データを固定し、パラメータを変化させる
  • 尤度関数を最大化するようなパラメータを選ぶ推定方法を最尤推定という
  • ベルヌーイ分布:数学において確率pで1、確率1-pで0をとる2値の離散分布。ロジスティック回帰のベースの確率分布。

 P(y)=p^y(1-p)^{1-y}

  • n回のベルヌーイ試行の同時確率による尤度関数

 P(y_1,y_2,...,y_n;p)=\prod_{i=1}^np^{y_i}(1-p)^{1-y_i}

ロジスティック回帰モデルの最尤推定

 P(y_1,y_2,...,y_n|w_0,w_1,...,w_m)=\prod_{i=1}^np^{y_i}(1-p)_{1-y_i}
                                                          =\prod_{i=1}^nσ(w^Tx_i)^{y_i}(1-σ(w^Tx_i))^{1-y_i}
                                                          =L(w)

  • ロジスティック回帰モデルの尤度関数は、通常対数をとる

 理由としては、積を和に変換できる、確率の掛け算は多重に行うと極小で桁落ちになり実務上取り扱いが出来ないなど。
  E(w)=-logL(w)=-\displaystyle \sum_{i=1}^{n}{y_ilogp_i+(1-y_i)log(1-p_i)}

勾配降下法

  • 反復学習によりパラメータを逐次的に更新するアプローチ
  • 対数尤度関数をパラメータで微分してゼロになる値を求める必要があるが、解析的にこの値を求めることは不可能。

 w(k+1)=w^k-η\dfrac{\partial E(w)}{\partial w}

  • ロジスティック回帰モデルでは、下記のように書ける。

 w(k+1)=w^k-η\displaystyle \sum_{i=1}^{n}(y_i-p_i)x_i]

 勾配降下法では、パラメータを更新するのにN個すべてのデータを求める必要がある。
 そのため、メモリの問題や、計算時間の問題が発生する。
 確率的勾配降下法では、データをランダムに選んでパラメータを更新する。
 勾配降下法でパラメータを一回更新更新するのと同じ計算量でパラメータをn回更新できる。

機械学習の精度指標

  • 混合行列(confusion matrix)
検証用データの結果
positive negative
モデルの予測結果 positive 真陽性 偽陰性 ]
negative 偽陽性 真陰性
  • 精度(accuracy)

 \dfrac {TP+TF}{TP+FN+FP+TN}

  • 再現率(Recall)

 全陽性群の中で、真陽性のデータの割合
 適用例:病変診断
 \dfrac{TP}{TP+FN}

  • 適合率(Precision)

モデルが陽性と予測したものの中で、真陽性のデータの割合
 適用例:スパムフィルタ
 \dfrac{TP}{TP+FP} 

 PresitionとRecallの調和平均
  \dfrac{2Recall\cdot Precision}{Recall+Precision}

ロジスティック回帰ハンズオン

  • テーマ:タイタニック
  • 課題:年齢が30歳で男の乗客は生き残れるか

環境設定、データ読み込み

import pandas as pd
import numpy as np
train=pd.read_csv('train.csv')
train_df=pd.DataFrame(train)

学習に使うデータのみ選定し、nullを中央値で補完

train_df=train_df.drop(["PassengerId","Name","Ticket","Fare","Cabin"], axis=1,inplace=True)
train_df['Age']=train_df['Age'].fillna(train_df['Age'].median())

性別で、maleを1、Femaleを0に数値化

train_df['Sex2']=train_df['Sex'].apply(lambda x:1 if x=='male' else 0 )

ターゲットデータとラベルを作成

data=train_df.loc[:,['Sex2','Age']]
label=train_df.loc[:,['Survived']]

ロジスティック回帰モデルを作成

from sklearn.linear_model import LogisticRegression as lr
model=lr()
model.fit(data,label)

30歳男性(Sex2=1, Age=30)を予測

model.predict([[1,30]])

array([0])
結果は、0つまり死亡の予測だった。
probaで(死亡(0)の予測確率、生存(1)の予測確率)を算出する。

model.predict_proba([[1,30]])

array(0.80695278, 0.19304722)
それぞれ80.7%、19.3%の予測結果となった。

主成分分析(PCA)

  • 多変量データの持つ構造をより少数個の指標に圧縮する次元削減の手法
    • 情報の量を分散の大きさととらえる
    • 線形返還後の変数の分散が最大となる射影軸を探索する
  • 以下の制約付き最適化問題を解く

目的関数:  arg max_{a\in{R^m}} a_j^TVar(\bar{X})a_j
制約条件 :  a_j^Ta_j=1

ラグランジュ関数を最大にする係数ベクトルを探索する(微分して0になる点を探索)
ラグランジュ関数
 E(a_j)=a_j^TVar(bar{X}a_j-λ(a_j^Ta_j-1)

    • 実務的な解法

1.分散共分散行列を計算
2.固有値問題を解く
3.mこの固有値固有ベクトルのペアを算出
4.寄与率から削減する次元数を決定

      • 寄与率

第k主成分の分散の全文さんに対する割合

      • 累積寄与率

第1-k主成分まで圧縮した時の情報損失量の割合

主成分分析のハンズオン

  • 設定:乳がん検査データ
  • 課題:32次元のデータを2次元上に次元圧縮した際に、うまく判別できるか

環境設定

import pandas as pd
from sklearn.model_selection import train_test_split
from sklearn.preprocessing import StandardScaler
from sklearn.linear_model import LogisticRegressionCV
from sklearn.metrics import confusion_matrix
from sklearn.decomposition import PCA
import matplotlib.pyplot as plt
%matplotlib inline

cancer_df=pd.read_csv(cancer.csv)

目的変数をdiagnosis、説明変数を3列目以降として、ロジスティック回帰で分類
>|python}
cancer_df.drop('Unnamed: 32', axis=1, inplace=True)
y = cancer_df.diagnosis.apply(lambda d: 1 if d == 'M' else 0)
X = cancer_df.loc[:, 'radius_mean':]
# 学習用とテスト用でデータを分離
X_train, X_test, y_train, y_test = train_test_split(X, y, random_state=0)

# 標準化
scaler = StandardScaler()
X_train_scaled = scaler.fit_transform(X_train)
X_test_scaled = scaler.transform(X_test)

# ロジスティック回帰で学習
logistic = LogisticRegressionCV(cv=10, random_state=0)
logistic.fit(X_train_scaled, y_train)

# 検証
print('Train score: {:.3f}'.format(logistic.score(X_train_scaled, y_train)))
print('Test score: {:.3f}'.format(logistic.score(X_test_scaled, y_test)))
print('Confustion matrix:\n{}'.format(confusion_matrix(y_true=y_test, y_pred=logistic.predict(X_test_scaled))))
|

Train score: 0.988
Test score: 0.972
検証スコア97%で分類できることを確認
目的変数と不要な変数1つを除いた30次元のデータについて、寄与率を観察

pca = PCA(n_components=30)
pca.fit(X_train_scaled)
plt.bar([n for n in range(1, len(pca.explained_variance_ratio_)+1)], pca.explained_variance_ratio_)

f:id:tibet:20211114212728p:plain
上位2変数の寄与度が非常に大きいため、2次元に圧縮しても分類情報が保たれる可能性は高い。
2次元まで圧縮する。

pca = PCA(n_components=2)
X_train_pca = pca.fit_transform(X_train_scaled)
print('X_train_pca shape: {}'.format(X_train_pca.shape))
# X_train_pca shape: (426, 2)

# 寄与率
print('explained variance ratio: {}'.format(pca.explained_variance_ratio_))
# explained variance ratio: [ 0.43315126  0.19586506]

# 散布図にプロット
temp = pd.DataFrame(X_train_pca)
temp['Outcome'] = y_train.values
b = temp[temp['Outcome'] == 0]
m = temp[temp['Outcome'] == 1]
plt.scatter(x=b[0], y=b[1], marker='o') # 良性は○でマーク
plt.scatter(x=m[0], y=m[1], marker='^') # 悪性は△でマーク
plt.xlabel('PC 1') # 第1主成分をx軸
plt.ylabel('PC 2') # 第2主成分をy軸

f:id:tibet:20211114212923p:plain

k近傍法

  • 分類問題のための機械学習手法
  • 最近棒のデータをk個とってきて、それらが最も多く所属するクラスに識別する
  • kを変化させると結果も変わる

k-means

1.各クラスタ中心の初期値を設定する。
2.各データ点に対して、各クラスタ中心との最も距離の近いクラスタを割り当てる
3.各クラスタの重心を計算する
4.収束するまで2,3の処理を繰り返す
中心の初期値を変えると、クラスタリングの結果も変わりえる

SVM

サポートベクトルマシン(SVM)は1990年代に提案され、2000年代前半に急速に発展した手法で、データ分析の場において近年最も注目を集めているモデルである。
オリジナルのSVMは、2クラス分類問題を対象としている。その後、回帰問題や教師無し問題などへも応用されている。

  • 決定境界と分類境界

2クラス分類問題では、特徴ベクトルxがどちらのクラスに属するか判定するために、次の決定関数が使われる。
 f(x)=w^Tx+b
データを分類するには、ある入力データxに対して決定関数f(x)を計算し、その符号によって2つのクラスに分類する。
一般に、x1-x2-..-xm平面とf(x)の平面の交線f(x)=0が特徴ベクトルxを2つのクラスに分ける境界線になっている。一般にこのような境界を分類境界と呼ぶ。

  • 線形サポートベクトル分類(ハードマージン)

 一般に訓練データを分離できる分類境界は複数存在する。SVMでは、それぞれのクラスのデータが分類境界からなるべく離れるようにして「最適な」分類境界を決定する。分類境界を挟んで2つのクラスがどのくらい離れているかをマージンと呼ぶ。なるべく大きなマージンを持つ分類境界を探す考え方をマージン最大化と呼ぶ。
訓練データを完全に分類できる分離境界があることを分離可能と呼ぶ。分離可能性を仮定したSV分類のことを、一般にハードマージンと呼ぶ。
分類境界の決定にかかわるのは分類境界に最も近いデータxiのみのため、xiをサポートベクトルと呼ぶ。

  • 線形サポートベクトル分類(ソフトマージン)

 ハードマージンでは訓練データを完璧に分類できる決定関数f(x)が存在すると仮定をおいたが、現実の多くの問題にとってこの仮定は強すぎる。そこで、分離可能でない問題までSV分類を拡張したものを、ソフトマージンと呼ぶ。マージン内に入るデータやご分類されたデータに対する誤差を表す変数をスラック変数と呼ぶ。

 線形分離ではうまく分類できないデータに対して、特徴ベクトルを非線形変換してその空間での線形分類を行う「カーネルトリック」と呼ばれる手法で非線形分離を行うことも可能である。
 非線形分離をするためには、写像Φ(x)で高次元空間へ変換する必要がある。しかし、最適化問題をそのまま解くことは不可能である。なぜならΦ(x)の内積部分の計算量が莫大になるからである。
そのため、この部分の計算を簡略化するためのテクニックはカーネルトリックと呼ばれ、Φの内積部分をと呼ばれる関数で置き換える方法である。
 K(x_i,X_j)=Φ(x_i)^TΦ(x_j)
決定関数f(x)はカーネル関数を用いて
[tex: f(x)=w^TΦ(x)+b=\displaystyle \sum_{i=1}^{n}α_iy_iK(x_i,x)+b
カーネル関数さえ決めれば、最適化や決定係数の計算にはΦ(x)は一切必要なくなる。
代表的なカーネル関数は次の3つ
多項式カーネル
ガウスカーネル
・シグモイドカーネル

SVMハンズオン

環境設定

%matplotlib inline
import numpy as np
import matplotlib.pyplot as plt
  • 線形分離可能な訓練データ生成
def gen_data():
    x0 = np.random.normal(size=50).reshape(-1, 2) - 2.
    x1 = np.random.normal(size=50).reshape(-1, 2) + 2.
    X_train = np.concatenate([x0, x1])
    ys_train = np.concatenate([np.zeros(25), np.ones(25)]).astype(np.int)
    return X_train, ys_train
X_train, ys_train = gen_data()
plt.scatter(X_train[:, 0], X_train[:, 1], c=ys_train)

f:id:tibet:20211114223222p:plain
学習
>|python}
t = np.where(ys_train == 1.0, 1.0, -1.0)

n_samples = len(X_train)
# 線形カーネル
K = X_train.dot(X_train.T)

eta1 = 0.01
eta2 = 0.001
n_iter = 500

H = np.outer(t, t) * K

a = np.ones(n_samples)
for _ in range(n_iter):
grad = 1 - H.dot(a)
a += eta1 * grad
a -= eta2 * a.dot(t) * t
a = np.where(a > 0, a, 0)
|

予測

index = a > 1e-6
support_vectors = X_train[index]
support_vector_t = t[index]
support_vector_a = a[index]

term2 = K[index][:, index].dot(support_vector_a * support_vector_t)
b = (support_vector_t - term2).mean()
xx0, xx1 = np.meshgrid(np.linspace(-5, 5, 100), np.linspace(-5, 5, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T

X_test = xx
y_project = np.ones(len(X_test)) * b
for i in range(len(X_test)):
    for a, sv_t, sv in zip(support_vector_a, support_vector_t, support_vectors):
        y_project[i] += a * sv_t * sv.dot(X_test[i])
y_pred = np.sign(y_project)
# 訓練データを可視化
plt.scatter(X_train[:, 0], X_train[:, 1], c=ys_train)
# サポートベクトルを可視化
plt.scatter(support_vectors[:, 0], support_vectors[:, 1],
                    s=100, facecolors='none', edgecolors='k')
# 領域を可視化
#plt.contourf(xx0, xx1, y_pred.reshape(100, 100), alpha=0.2, levels=np.linspace(0, 1, 3))
# マージンと決定境界を可視化
plt.contour(xx0, xx1, y_project.reshape(100, 100), colors='k',
                     levels=[-1, 0, 1], alpha=0.5, linestyles=['--', '-', '--'])

# マージンと決定境界を可視化
plt.quiver(0, 0, 0.1, 0.35, width=0.01, scale=1, color='pink')

f:id:tibet:20211114223511p:plain

  • 線形分離不可能な訓練データ生成
factor = .2
n_samples = 50
linspace = np.linspace(0, 2 * np.pi, n_samples // 2 + 1)[:-1]
outer_circ_x = np.cos(linspace)
outer_circ_y = np.sin(linspace)
inner_circ_x = outer_circ_x * factor
inner_circ_y = outer_circ_y * factor

X = np.vstack((np.append(outer_circ_x, inner_circ_x),
               np.append(outer_circ_y, inner_circ_y))).T
y = np.hstack([np.zeros(n_samples // 2, dtype=np.intp),
               np.ones(n_samples // 2, dtype=np.intp)])
X += np.random.normal(scale=0.15, size=X.shape)
x_train = X
y_train = y
plt.scatter(x_train[:,0], x_train[:,1], c=y_train)

f:id:tibet:20211114223656p:plain

学習
カーネルとしてRBFカーネルを利用する

def rbf(u, v):
        sigma = 0.8
        return np.exp(-0.5 * ((u - v)**2).sum() / sigma**2)
    
X_train = x_train
t = np.where(y_train == 1.0, 1.0, -1.0)

n_samples = len(X_train)
# RBFカーネル
K = np.zeros((n_samples, n_samples))
for i in range(n_samples):
    for j in range(n_samples):
        K[i, j] = rbf(X_train[i], X_train[j])

eta1 = 0.01
eta2 = 0.001
n_iter = 5000

H = np.outer(t, t) * K

a = np.ones(n_samples)
for _ in range(n_iter):
    grad = 1 - H.dot(a)
    a += eta1 * grad
    a -= eta2 * a.dot(t) * t
    a = np.where(a > 0, a, 0)

予測

index = a > 1e-6
support_vectors = X_train[index]
support_vector_t = t[index]
support_vector_a = a[index]

term2 = K[index][:, index].dot(support_vector_a * support_vector_t)
b = (support_vector_t - term2).mean()
xx0, xx1 = np.meshgrid(np.linspace(-1.5, 1.5, 100), np.linspace(-1.5, 1.5, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T

X_test = xx
y_project = np.ones(len(X_test)) * b
for i in range(len(X_test)):
    for a, sv_t, sv in zip(support_vector_a, support_vector_t, support_vectors):
        y_project[i] += a * sv_t * rbf(X_test[i], sv)
y_pred = np.sign(y_project)
# 訓練データを可視化
plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train)
# サポートベクトルを可視化
plt.scatter(support_vectors[:, 0], support_vectors[:, 1],
                    s=100, facecolors='none', edgecolors='k')
# 領域を可視化
plt.contourf(xx0, xx1, y_pred.reshape(100, 100), alpha=0.2, levels=np.linspace(0, 1, 3))
# マージンと決定境界を可視化
plt.contour(xx0, xx1, y_project.reshape(100, 100), colors='k',
                     levels=[-1, 0, 1], alpha=0.5, linestyles=['--', '-', '--'])

f:id:tibet:20211114223835p:plain

  • ソフトマージンSVM

重なりのある訓練データの生成

x0 = np.random.normal(size=50).reshape(-1, 2) - 1.
x1 = np.random.normal(size=50).reshape(-1, 2) + 1.
x_train = np.concatenate([x0, x1])
y_train = np.concatenate([np.zeros(25), np.ones(25)]).astype(np.int)
plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train)

f:id:tibet:20211114224014p:plain
学習
スラック変数を導入する

X_train = x_train
t = np.where(y_train == 1.0, 1.0, -1.0)

n_samples = len(X_train)
# 線形カーネル
K = X_train.dot(X_train.T)

C = 1
eta1 = 0.01
eta2 = 0.001
n_iter = 1000

H = np.outer(t, t) * K

a = np.ones(n_samples)
for _ in range(n_iter):
    grad = 1 - H.dot(a)
    a += eta1 * grad
    a -= eta2 * a.dot(t) * t
    a = np.clip(a, 0, C)

予測

index = a > 1e-8
support_vectors = X_train[index]
support_vector_t = t[index]
support_vector_a = a[index]

term2 = K[index][:, index].dot(support_vector_a * support_vector_t)
b = (support_vector_t - term2).mean()
xx0, xx1 = np.meshgrid(np.linspace(-4, 4, 100), np.linspace(-4, 4, 100))
xx = np.array([xx0, xx1]).reshape(2, -1).T

X_test = xx
y_project = np.ones(len(X_test)) * b
for i in range(len(X_test)):
    for a, sv_t, sv in zip(support_vector_a, support_vector_t, support_vectors):
        y_project[i] += a * sv_t * sv.dot(X_test[i])
y_pred = np.sign(y_project)
# 訓練データを可視化
plt.scatter(x_train[:, 0], x_train[:, 1], c=y_train)
# サポートベクトルを可視化
plt.scatter(support_vectors[:, 0], support_vectors[:, 1],
                    s=100, facecolors='none', edgecolors='k')
# 領域を可視化
plt.contourf(xx0, xx1, y_pred.reshape(100, 100), alpha=0.2, levels=np.linspace(0, 1, 3))
# マージンと決定境界を可視化
plt.contour(xx0, xx1, y_project.reshape(100, 100), colors='k',
                     levels=[-1, 0, 1], alpha=0.5, linestyles=['--', '-', '--'])

f:id:tibet:20211114224202p:plain