畳みこみニューラルネットを0から実装する (第2回)
2022年1月に加筆を行って記事を書き直しました！
こんにちはtatsyです。
畳みこみニューラルネットを0から実装する記事の第2回です。前回の第1回ではMNISTから画像とラベルのデータを読み込みました。第1回の記事のご確認は以下のリンクからお願いします。
第2回の記事の内容
第2回となる今回はシグモイド関数を使った通常のニューラルネットを対象に、学習パラメータの更新に用いられる誤差逆伝搬法 (back propagation)について解説します。
ニューラルネットの概要
詳しい解説はその他の専門書にゆずるとして、ニューラルネットとはベクトル値ベクトル関数 (ベクトルをとってベクトルを返す関数) の合成だと考えることができます。
例えばMNISTの例ならニューラルネットの入力は画像で、出力はラベルです。
入力の画像は画素値を要素として画素数と同じだけの要素数を持つベクトルとみなすことができます。出力のラベルは(学習の都合上ですが)十次元のベクトルとみなせます。正解は0から9までの数字ですから、正解をxとするとx番号目の要素だけが1でそれ以外が0の十次元ベクトルです。
関数の合成だ、といったのはニューラルネットが複数のレイヤーから構成されるためです。ここでいうレイヤーからレイヤーへの遷移が上で述べたベクトル値ベクトル関数というわけです。
入力層、隠れ層、出力層からなる三層のニューラルネットを例に取れば、ニューラルネット全体が、「入力層から隠れ層への関数」と「隠れ層から出力層への関数」という2つの関数の合成だというわけです。
ニューラルネットの学習
ここまで、ニューラルネットはベクトルを入力としてベクトルを出力とする関数(=ベクトル値ベクトル関数)とみなせる、という話をしてきました。以下では、ニューラルネット全体が関数$f$であると考えて話を勧めます。
ニューラルネットは多くの問題に適用可能ですが、ここでは入出力のデータの組$(\mathbf{y}, \mathbf{x})$が与えられたときに$\mathbf{y} = f(\mathbf{x})$となるような$f$を決める、という回帰問題を考えてみましょう。
ここで、隠れ層を1つもつ三層のニューラルネットを考えます。まず、入力層から隠れ層への関数を$f_1$とします。すると、隠れ層で受け取る出力は、 $$ \mathbf{h} = f_1(\mathbf{x}) $$ です。同じように隠れ層から出力層への関数を$f_2$とすると出力層で受け取るベクトルは、 $$ \mathbf{y} = f_2(\mathbf{h}) = f_2(f_1(\mathbf{x})) $$ となります。
したがって、訓練データを元にしてなるべく$\mathbf{y}$と$f_2(f_1(\mathbf{x}))$が近くなるように$f_1, f_2$を定めるというのが学習の目的になります。
例えばデータの出力ベクトル$\mathbf{y}$の近さを図る指標としてユークリッドノルムを採用すると、 $$ \begin{equation} L(\mathbf{x}, \mathbf{y}) = \frac{1}{2} \| \mathbf{y} - f_2(f_1(\mathbf{x})) \|^2 \label{eq:l2-loss} \end{equation} $$ を最小化するような$f_1, f_2$を求められれば良いことになります。
全結合層の概要
それでは$f_1, f_2$にもう少し具体的な関数を入れて考えてみます。ニューラルネットの層は、一般的に「線形の演算」と「非線形の演算」の組み合わせによりできています。
この際、非線形の計算を、より単純な形で扱うために、線形の演算部分には学習により決定されるパラメータをもたせ、非線形の演算はシグモイド関数に代表されるような「活性化関数」によって行います (例外としてPReLUのようなパラメータを持つ活性化関数もある)。
すなわち、関数$f_i$は一般的に以下のような形で書き表せます。 $$ f_i(\mathbf{x}) = \sigma(\mathbf{W}_i \mathbf{x} + \mathbf{b}_i) $$ この式で$\sigma$は非線形の活性化関数で、その引数として行列の掛け算とベクトルの足し算という線形の演算結果が入っています。
これを用いて、改めてニューラルネットの全体の関数を書き直してみると、 $$ f(\mathbf{x}) = \sigma(\mathbf{W}_2 \sigma(\mathbf{W}_1 \mathbf{x} + \mathbf{b}_1) + \mathbf{b}_2) $$ のように書けます。
ニューラルネットを用いた機械学習は$f$に含まれるパラメータである$(\mathbf{W}_1, \mathbf{b}_1, \mathbf{W}_2, \mathbf{b}_2)$を\eqref{eq:l2-loss}のような目的関数を最小化するように最適化することでした。
最適化には最急降下法のような勾配に基づく最適化法を用いるので、目的関数$L$の$\mathbf{W}_i, \mathbf{b}_i$に関する勾配を求めます。なお、後ほど損失関数を変更してご説明する都合上、以下では誤差関数$L$を具体的に書き下すことはせずに勾配を計算してみます。
まずは、より出力に近い$\mathbf{W}_2, \mathbf{b}_2$に関する勾配を考えてみます。$L$の引数になっている$\mathbf{y}$は$\mathbf{W}_2, \mathbf{b}_2$の関数なので、合成関数の微分を考えると、 $$ \begin{aligned} \frac{\partial L}{\partial \mathbf{W}_2} &= \frac{\partial L}{\partial \mathbf{y}} \frac{\partial \mathbf{y}}{\partial \mathbf{W}_2} = \frac{\partial L}{\partial \mathbf{y}} \sigma’(\mathbf{W}_2 \mathbf{h} + \mathbf{b}_2) \mathbf{h}^\top \\ \frac{\partial L}{\partial \mathbf{b}_2} &= \frac{\partial L}{\partial \mathbf{y}} \frac{\partial \mathbf{y}}{\partial \mathbf{b}_2} = \frac{\partial L}{\partial \mathbf{y}} \sigma’(\mathbf{W}_2 \mathbf{h} + \mathbf{b}_2) \end{aligned} $$
この様な形で、一番上のレイヤーに関する微分は比較的簡単に求まるのですが、これを下へ下へと微分していくのは一見かなり骨が折れそうです。
ですが、実は下のレイヤーのパラメータでの偏微分は上のレイヤーでのパラメータの偏微分とチェインルールと呼ばれる関係を持っており、この関係を用いて各レイヤーでの偏微分を簡単に計算しようというのが誤差逆伝播です。
チェインルールとは、高校までの数学の言葉で言えば、要するに「合成関数」の微分で、複数の関数が「鎖状」に合成されている時、その合成関数の微分が、合成されている各関数の「微分の積」で表せる、ということです。
実際に、$f$を$\mathbf{W}_1, \mathbf{b}_1$で微分しようとすると、以下のような式が得られることが分かります。 $$ \begin{aligned} \frac{\partial L}{\partial \mathbf{W}_1} &= \frac{\partial L}{\partial \mathbf{y}} \frac{\partial \mathbf{y}}{\partial \mathbf{h}} \frac{\partial \mathbf{h}}{\partial \mathbf{W}_1} = \frac{\partial L}{\partial \mathbf{h}} \frac{\partial f_1(\mathbf{x})}{\partial \mathbf{W}_1} \\ \frac{\partial L}{\partial \mathbf{b}_1} &= \frac{\partial L}{\partial \mathbf{y}} \frac{\partial \mathbf{y}}{\partial \mathbf{h}} \frac{\partial \mathbf{h}}{\partial \mathbf{b}_1} = \frac{\partial L}{\partial \mathbf{h}} \frac{\partial f_1(\mathbf{x})}{\partial \mathbf{b}_1} \end{aligned} $$ この式を注意深く見てみると、ネットワーク中に合わられる関数$f_i$のパラメータに関する微分というのは、以下の2つの成分の積であることが分かります。
- 関数$f_i$の出力となる変数で誤差関数を微分したもの
- 関数$f_i$を目的のパラメータで微分したもの
この性質を用いると、誤差関数に近い方にある関数から順に微分を計算していくことで、効率的にパラメータに関する微分を計算することができます。
パラメータの更新
今回のコードではパラメータの更新に単純な最急降下法を用います。最急降下法はパラメータの勾配、例えば$\frac{\partial L}{\partial \mathbf{W}}$を計算後、ステップ幅$\gamma$をかけてパラメータを更新します。 $$ \begin{aligned} \Delta \mathbf{W}^{t} &= \frac{\partial L}{\partial \mathbf{W}^{t}} \\ \mathbf{W}^{t+1} &= \mathbf{W}^t - \gamma \Delta \mathbf{W}^{t} \end{aligned} $$ 上記の式で上付きの添字は$t$回更新された、ということを意味しています。
基本的には上記の式である程度は最適化可能なのですが、現在のニューラルネットの学習は計算コストと学習効率の観点から、データセットの一部だけを用いて勾配を計算する確率的最急降下法 (stochastic gradient descent)を用いるのが一般的です。
学習に使うデータセットの一部は、特にバッチと呼ばれていて、通常、数個から数十個の学習データのペアだけを用いて勾配を計算します。ですが、計算される勾配がどれほどデータ全体から計算される勾配に近いかはバッチの性質にもよります。したがって、バッチの選び方がたまたま一部のデータに偏っていた場合、計算される勾配の性質が悪く、最適化に失敗してしまう可能性があります。
この問題を防ぐためには、いろいろなやり方が提案されていますが、今回は単純な慣性 (momentum)を勾配の安定化手法を用います。ここで言う慣性とは、勾配方向が急激に変化するのを防ぐ成分で、パラメータとしては慣性の強さを表す$\mu$を使います。慣性を加えたパラメータの更新式は以下のように書けます。 $$ \begin{aligned} \Delta \mathbf{W}^{t} &= \mu \mathbf{W}^{t-1} + \gamma \frac{\partial L}{\partial \mathbf{W}^{t}} \\ \mathbf{W}^{t+1} &= \mathbf{W}^t - \Delta \mathbf{W}^{t} \end{aligned} $$ この式に示すように、慣性を用いた確率的最急降下法では、過去の勾配方向をある程度維持しつつパラメータを更新することで、バッチ中のデータの偏りによる不適切なパラメータ更新を防いでいます。
なお、Yan LeCunによるオリジナルの畳みこみニューラルネットの論文[LeCun98]ではルーベンバーグ・マルカート法 (Levenberg Marquardt method) という最急降下法とニュートン法を組み合わせたアルゴリズムが用いられていますが、2022年現在も、基本的には最急降下法に基づく最適化アルゴリズムが主流で、二階微分まで考えるルーベンバーグ・マルカート法のようなアルゴリズムはあまり使われていないように思います。
パラメータの初期化
多くの最適化問題で最適化される変数に対して適当な初期値が与えられるように、ニューラルネット中のパラメータ (全結合層であれば$\mathbf{W}$, $\mathbf{b}$)も、最初に適当な初期値で初期化しておいて、そこから学習を開始します。
よくある最適化問題では最適化される変数について何らかの事前知識が得られていて、それにしたがって初期値を設定することも多いのですが、ニューラルネットの場合にはパラメータも多く、定式化についてもモデル化される対象に即したものとは言いづらい曖昧なものなので、パラメータはランダムネスを持った値で初期化するのが一般的です。
その際、よく使われるテクニックにXavier Glorot uniformあるいはXavier Glorot normalというものがあります。これは、Xavier Glorot氏とYoshua Bengio氏によって2010年に発表された以下の論文で提案されたパラメータの初期化手法で、前向きの関数の評価と誤差逆伝播の最中に各レイヤーを評価する時、入出力の変数が持つ値の分散が一定に保たれることを仮定した重みを導入します。
- X. Glorot and Y. Bengio, “Understanding the difficulty of training deep feedforward neural networks,” AISTATS 2010. [link]
より具体的には、入出力の変数の数をそれぞれ$n_{in}, n_{out}$であるレイヤーの重み変数$\mathbf{W}$を初期化する場合、正規分布で初期化するなら、その標準偏差を $$ stddev = \frac{\sqrt{2}}{\sqrt{n_{in} + n_{out}}} $$ とし、一様分布で初期化するならば、一様分布の下端ならびに上端を、それぞれ $$ lower = -\frac{\sqrt{6}}{\sqrt{n_{in} + n_{out}}}, upper = \frac{\sqrt{6}}{\sqrt{n_{in} + n_{out}}} $$ とするというものです。
実装例: 全結合層
全体コード: https://github.com/tatsy/educnn/blob/master/sources/fully_connected_layer.h
以下の実装例では、順方向の関数の評価をforwardで、誤差逆伝播の関数をbackwardで書いています。パラメータはXavier Glorotの重みを用いて正規分布でパラメータを初期化しています。
forward関数の中では、関数の出力を計算するとともに、誤差逆伝播に使用するために、input_、output_という変数に関数の入出力を保存しています。
backward関数では、先ほどのチェインルールの式に従って、各パラメータ $\mathbf{W}$, $\mathbf{b}$の勾配を計算しています。計算された勾配は、慣性付き確率的最急降下法の式に基づいて、パラメータの更新に用います。
class FullyConnectedLayer : public AbstractLayer {
public:
// Public methods
FullyConnectedLayer()
: AbstractLayer() {
}
FullyConnectedLayer(int input_size, int output_size)
: AbstractLayer()
, input_size_(input_size)
, output_size_(output_size)
, W(output_size, input_size)
, b(1, output_size)
, dW(output_size, input_size)
, db(1, output_size) {
// X. Glorot's standard deviation for parameter initialization
// X. Glorotによるパラメータ初期化のための標準偏差
const double xg_stddev = sqrt(2.0 / (input_size_ + output_size_));
// Parameter initialization
// パラメータの初期化
Random &rng = Random::getInstance();
for (int i = 0; i < output_size; i++) {
for (int j = 0; j < input_size; j++) {
W(i, j) = rng.normal() * xg_stddev;
dW(i, j) = 0.0;
}
b(0, i) = 0.0;
db(0, i) = 0.0;
}
}
virtual ~FullyConnectedLayer() {
}
const Matrix &forward(const Matrix &input) override {
// Simple linear operation (y = Wx + b)
// 単純な線形演算 (y = Wx + b)
const int batchsize = (int)input.rows();
input_ = input;
output_ = input * W.transpose() + b.replicate(batchsize, 1);
return output_;
}
Matrix backward(const Matrix &dLdy, double lr = 0.1, double momentum = 0.5) override {
// Assum x and y are input and output of this layer, hence back-prop transforms dLdy to dLdx.
// xとyがこのレイヤーの入出力だと仮定. 誤差逆伝播のためにdLdyをdLdxに変換する
const int batchsize = (int)dLdy.rows();
const Matrix dLdx = dLdy * W;
// Momentum SGD
// 慣性つき確率的最急降下法
const Matrix current_dW = dLdy.transpose() * input_;
const Matrix current_db = dLdy.colwise().sum();
dW = momentum * dW + lr * current_dW;
db = momentum * db + lr * current_db;
// Update parameters
// パラメータの更新
W -= dW;
b -= db;
return dLdx;
}
private:
// Private parameters
int input_size_ = 0;
int output_size_ = 0;
Matrix W = {};
Matrix b = {};
Matrix dW = {};
Matrix db = {};
}; // class FullyConnectedLayer
まとめ
今回の記事では、全結合層を例に誤差逆伝播の方法について紹介しました。合わせて確率的最急降下法とXavier Glorotの重みを用いたパラメータの初期化手法についても紹介しています。
次回は誤差関数の定義とバッチを用いた学習全体の流れについてご紹介する予定です。
今回も最後までお読みいただき、ありがとうございました。