Theano: DeepLearning : ロジスティック回帰による MNIST 分類

DeepLearning 0.1 文書: ロジスティック回帰による MNIST 分類（翻訳/要約）

* DeepLearning 0.1 documentation: Classifying MNIST digits using Logistic Regression の簡単な要約です。

本文

[Note] このセクションは次の Theano コンセプトを理解していることを仮定しています: 共有変数 , 基本算術 ops , T.grad , floatX。GPU 上でのコードの実行を意図するならば、GPU も読んでください。

[Note] このセクションのためのコードはここからダウンロード可能です。GitHub も。

このセクションでは、もっとも基本的な分類器: ロジスティック回帰を実装するのに Theano がどのように使用されるかを示します。refresher としてもまた表記のよりどころとしても機能する、モデルの迅速な入門から始めます、そして数式がどのように Theano グラフ上にマップされるかを示します。

モデル

ロジスティック回帰は確率的な、線形分類器です。それは重み行列 W とバイアス・ベクトル b でパラメータ化されます。分類は入力ベクトルを、それぞれクラスに相当する、超平面の集合の上に射影されることでなされます。入力から超平面への距離は入力が該当するクラスのメンバーである確率を反映します。

数学的には、入力ベクトル x がクラス i のメンバーである確率、確率変数 Y の値は次のように書けます :

$P(Y=i|x, W,b) &= softmax_i(W x + b) \\ &= \frac {e^{W_i x + b_i}} {\sum_j e^{W_j x + b_j}}$

モデルの予測 $y_{pred}$ は確率が最大となるクラスです、明示的に :

$y_{pred} = {\rm argmax}_i P(Y=i|x,W,b)$

Theano でこれを行なうコードは次のようなものです :

        # 重み W を shape  (n_in, n_out) の行列として 0 で初期化する
        self.W = theano.shared(
            value=numpy.zeros(
                (n_in, n_out),
                dtype=theano.config.floatX
            ),
            name='W',
            borrow=True
        )
        # バイアス b を n_out 0 のベクトルとして初期化する
        self.b = theano.shared(
            value=numpy.zeros(
                (n_out,),
                dtype=theano.config.floatX
            ),
            name='b',
            borrow=True
        )

        # class-membership 確率の行列を計算するためのシンボリック式
        # ここで:
        # W は行列で、カラム-k はクラス-k のための分離超平面を表す。
        # x は行列で、行-j は入力訓練サンプル-j を表します。
        # b はベクトルで、要素-k は超平面-k の自由パラメータを表します。
        self.p_y_given_x = T.nnet.softmax(T.dot(input, self.W) + self.b)

        # 確率が最大となるクラスとしての予測をどのように計算するかのシンボリックな表現。
        self.y_pred = T.argmax(self.p_y_given_x, axis=1)

モデルのパラメータは訓練を通して永続的な状態を維持しなければならないので、 $W,b$ のために共有変数を割り当てます。これはそれら両者をシンボリックな Theano 変数として宣言しますが、その内容も初期化します。そして dot と softmax オペレータはベクトル $P(Y|x, W,b)$ を計算するために使用されます。結果 p_y_given_x はベクトル型のシンボリック変数です。

実際のモデル予測を得るためには、T.argmax オペレータを使用できます、これは p_y_given_x が最大となるインデックスを返します（i.e. 最大確率のクラス）。

さてもちろん、ここまでに定義したモデルはまだ有用なことは何もしていません、そのパラメータが依然として初期状態にありますから。次のセクションでは最適なパラメータをどのように学習するかをカバーします。

[Note] Theano ops の完全なリストは、list of ops 参照。

損失関数を定義する

最適なモデルパラメータんの学習は損失関数を最小化することも伴います。多クラス・ロジスティック回帰の場合には、損失として負の対数尤度を用いることが普通です。

これはモデルパラメータ $\theta$ のもとでデータセット $\cal{D}$ の尤度を最大化することと等価です。最初に尤度 $\cal{L}$ と損失 $\ell$ を定義することから始めましょう :

$\mathcal{L} (\theta=\{W,b\}, \mathcal{D}) = \sum_{i=0}^{|\mathcal{D}|} \log(P(Y=y^{(i)}|x^{(i)}, W,b)) \\ \ell (\theta=\{W,b\}, \mathcal{D}) = - \mathcal{L} (\theta=\{W,b\}, \mathcal{D})$

書籍全体が最小化のトピックに捧げられていますが、勾配降下はこれまでのところ任意の非線形関数を最小化するためのもっとも単純な方法です。このチュートリアルではミニバッチによる確率的勾配降下 (MSGD) の方法を用います。

次の Theano コードは与えられたミニバッチのための（シンボリック）損失を定義します :

        # y.shape[0] は y の（シンボリックな）行数です、i.e. ミニバッチのサンプル数 (n) です。
        # T.arange(y.shape[0]) はシンボリックなベクトルで、[0,1,2,... n-1] を含みます。
        # T.log(self.p_y_given_x) は サンプル毎の１行とクラス毎の１カラムを持つ Log-確率 (LP) の行列です。
        # LP[T.arange(y.shape[0]),y] はベクトル v で [LP[0,y[0]], LP[1,y[1]], LP[2,y[2]], ..., LP[n-1,y[n-1]]] を含みます。
        # そして T.mean(LP[T.arange(y.shape[0]),y]) は v の要素の（ミニバッチ・サンプルに渡る）平均です。
        # i.e., ミニバッチに渡る対数尤度の平均です。
        return -T.mean(T.log(self.p_y_given_x)[T.arange(y.shape[0]), y])

[Note] 損失はフォーマルには個々のエラー項の、データセットに渡る、合計 (sum) として定義されますが、実際には、コードでは平均 (T.mean) を使用します。これは学習率の選択についてミニバッチ・サイズへの依存を小さくすることを可能にします。

LogisticRegression クラスを作成する

さて、ロジスティック回帰の基本的な挙動をカプセル化する、LogisticRegression クラスを定義するために必要なツールの全てが揃いました。コードはここまでにカバーしてきたことと非常に類似していて、自明であるべきです。

class LogisticRegression(object):
    """多クラス ロジスティック回帰クラス

    ロジスティック回帰は重み行列: math:`W` とバイアス・ベクトル:math:`b` で完全に記述されます。
    分類はデータポイントの超平面の集合への射影で行なわれます、
    それへの距離はクラスに属する確率を決定するために使われます。
    """

    def __init__(self, input, n_in, n_out):
        """ ロジスティック回帰のパラメータを初期化します。

        :type input: theano.tensor.TensorType
        :param input: シンボリック変数、アーキテクチャ（一つのミニバッチ）の入力を記述します。

        :type n_in: int
        :param n_in: 入力ユニットの数、データポイントが存在するスペースの次元

        :type n_out: int
        :param n_out: 出力ユニットの数、ラベルが存在するスペースの次元

        """
        # start-snippet-1
        # 重み W を shape  (n_in, n_out) の行列として 0 で初期化する
        self.W = theano.shared(
            value=numpy.zeros(
                (n_in, n_out),
                dtype=theano.config.floatX
            ),
            name='W',
            borrow=True
        )
        # バイアス b を n_out 0 のベクトルとして初期化する
        self.b = theano.shared(
            value=numpy.zeros(
                (n_out,),
                dtype=theano.config.floatX
            ),
            name='b',
            borrow=True
        )

        # クラスに属する確率の行列を計算するためのシンボリック式
        # ここで:
        # W は行列で、カラム-k はクラス-k のための分離超平面を表す。
        # x は行列で、行-j は入力訓練サンプル-j を表します。
        # b はベクトルで、要素-k は超平面-k の自由パラメータを表します。
        self.p_y_given_x = T.nnet.softmax(T.dot(input, self.W) + self.b)

        # 確率が最大となるクラスとしての予測をどのように計算するかのシンボリックな表現。
        self.y_pred = T.argmax(self.p_y_given_x, axis=1)
        # end-snippet-1

        # モデルのパラメータ
        self.params = [self.W, self.b]

        # モデル入力を追跡
        self.input = input

    def negative_log_likelihood(self, y):
        """与えられたターゲット分布に基づく、このモデルの予測の負の対数尤度の平均を返します。

        .. math::

            \frac{1}{|\mathcal{D}|} \mathcal{L} (\theta=\{W,b\}, \mathcal{D}) =
            \frac{1}{|\mathcal{D}|} \sum_{i=0}^{|\mathcal{D}|}
                \log(P(Y=y^{(i)}|x^{(i)}, W,b)) \\
            \ell (\theta=\{W,b\}, \mathcal{D})

        :type y: theano.tensor.TensorType
        :param y: corresponds to a vector that gives for each example the
                  correct label

        Note: 合計 (sum) の代わりに平均 (sum) を使用します、
              その結果、学習率のバッチサイズへの依存が小さくなります。
        """
        # start-snippet-2
        # y.shape[0] は y の（シンボリックな）行数です、i.e. ミニバッチのサンプル数 (n) です。
        # T.arange(y.shape[0]) はシンボリックなベクトルで、[0,1,2,... n-1] を含みます。
        # T.log(self.p_y_given_x) は サンプル毎の１行とクラス毎の１カラムを持つ Log-確率 (LP) の行列です。
        # LP[T.arange(y.shape[0]),y] はベクトル v で [LP[0,y[0]], LP[1,y[1]], LP[2,y[2]], ..., LP[n-1,y[n-1]]] を含みます。
        # そして T.mean(LP[T.arange(y.shape[0]),y]) は v の要素の（ミニバッチ・サンプルに渡る）平均です。
        # i.e., ミニバッチに渡る対数尤度の平均です。
        return -T.mean(T.log(self.p_y_given_x)[T.arange(y.shape[0]), y])
        # end-snippet-2

    def errors(self, y):
        """Return a float representing the number of errors in the minibatch
        over the total number of examples of the minibatch ; zero one
        loss over the size of the minibatch

        :type y: theano.tensor.TensorType
        :param y: corresponds to a vector that gives for each example the
                  correct label
        """

        # check if y has same dimension of y_pred
        if y.ndim != self.y_pred.ndim:
            raise TypeError(
                'y should have the same shape as self.y_pred',
                ('y', y.type, 'y_pred', self.y_pred.type)
            )
        # check if y is of the correct datatype
        if y.dtype.startswith('int'):
            # the T.neq operator returns a vector of 0s and 1s, where 1
            # represents a mistake in prediction
            return T.mean(T.neq(self.y_pred, y))
        else:
            raise NotImplementedError()

このクラスを次のようにインスタンス化します :

    # 入力のためのシンボリック変数を生成する（x と y はミニバッチを表します。）
    # generate symbolic variables for input (x and y represent a
    # minibatch)
    x = T.matrix('x')  # データ、ラスター化された画像として表されます
    y = T.ivector('y')  # ラベル、[int] ラベルの 1D ベクトルとして表されます

    # ロジスティック回帰クラスを構築する
    # 各 MNIST 画像はサイズ 28*28 を持ちます。
    classifier = LogisticRegression(input=x, n_in=28 * 28, n_out=10)

訓練入力 $x$ と相当するクラス $y$ のためにシンボリック変数を割り当てることから始めます。x と y は LogisticRegression オブジェクトのスコープの外で定義されていることに注意してください。クラスはそのグラフを構築するために入力を必要としますので、それは __init__ 関数のパラメータとして渡されます。深層ネットワークを形作るためにそのようなクラス群のインスタンスを結合したい時には有用です。一つの層の出力は上の層の入力として渡すことができます。（このチュートリアルは多層ネットワークを構築しませんが、このコードはそれを行なう先々のチュートリアルで再利用されます。）

最後に、classifier.negative_log_likelihood インスタンス・メソッドを使用して、最小化する（シンボリックな）コスト変数を定義します。

    # 訓練中に最小化するコストは、シンボリックなフォーマットのモデルの負の対数尤度
    cost = classifier.negative_log_likelihood(y)

x は、コストの定義への暗黙的なシンボリックな入力である点に注意してください。何故なら、分類器のシンボリック変数は初期時に x の観点から定義されたからです。

モデルを学習する (Learning the Model)

多くのプログラミング言語 (C/C++, Matlab, Python) で MSGD を実装するためには、パラメータに関する損失の勾配のために手動で式を導出することから始めるでしょう: この場合は $\partial{\ell}/\partial{W}$ と $\partial{\ell}/\partial{b}$ ですが、これは複雑なモデルのためにはかなり技巧的になります、何故なら $\partial{\ell}/\partial{\theta}$ のための式が非常に複雑になるからで、特に数値的安定性を考慮に入れた場合に時に。

Theano では、この作業は大きく単純化されます。自動微分を行ない、そして数値的安定性を改善するために幾つかの数学的な変換を適用します。Theano で $\partial{\ell}/\partial{W}$ と $\partial{\ell}/\partial{b}$ の勾配を得るためには、単に次を行ないます :

    g_W = T.grad(cost=cost, wrt=classifier.W)
    g_b = T.grad(cost=cost, wrt=classifier.b)

g_W と g_b はシンボリック変数で、計算グラフの一部として使用されます。そして勾配降下の 1 ステップを実行する関数 train_model は次のように定義できます :

    #  (変数, 更新式) ペアのリストとして
    # モデルのパラメータをどのように更新するかを指定
    updates = [(classifier.W, classifier.W - learning_rate * g_W),
               (classifier.b, classifier.b - learning_rate * g_b)]

    # cost を返す Theano 関数 `train_model` をコンパイルする、
    # 同時にモデルのパラメータを `updates` で定義されたルールを基に更新する
    train_model = theano.function(
        inputs=[index],
        outputs=cost,
        updates=updates,
        givens={
            x: train_set_x[index * batch_size: (index + 1) * batch_size],
            y: train_set_y[index * batch_size: (index + 1) * batch_size]
        }
    )

ペアのリストを更新します。最初の要素はステップで更新されるシンボリック変数で、２番目の要素はその新しい値を計算するためのシンボリック関数です。同様に、辞書が与えられた時には、そのキーはシンボリック変数でその値はステップの間の置き換えを指定します。

そして関数 train_model は次のように定義されます :

入力はミニバッチ・インデックス index で、バッチサイズ（= これは固定されているので入力ではありません）と一緒に x を該当ラベル y とともに定義します。
戻り値は index で定義される x, y と関連する cost/loss です。
全ての関数呼び出しにおいて、最初に x と y を index で指定された訓練セットからのスライスで置き換えます。そして、ミニバッチに関連する cost を評価して updates リストで定義された演算を適用します。

train_model(index) が呼び出されるたびに、MSGD のステップを実行しながら、ミニバッチの cost を計算して返します。全体の学習アルゴリズムはこうしてデータセットの全てのサンプルに渡りループして、一度に一つのミニバッチのサンプル全てを考慮し、train_model 関数を繰り返し呼び出すことによって成ります。

モデルをテストする

モデルをテストする時には分類し損なった (misclassified) サンプルの数に興味があります（そして尤度だけでなく）。それ故に LogisticRegression クラスは余分のインスタンス・メソッドを持ち、これは各ミニバッチで分類し損なったサンプルの数を取得するためのシンボリック・グラフを構築します。

コードは次のようなものです :

    def errors(self, y):
        """ミニバッチ・サンプルの総数に渡ってミニバッチのエラー数を表す float を返す;
        ミニバッチ・サイズに渡る 0-1 損失

        :type y: theano.tensor.TensorType
        :param y: 各サンプルに対して正しいラベルを与えるベクトルに相当します。
        """

        # y が y_pred と同じ次元を持つかチェックする
        if y.ndim != self.y_pred.ndim:
            raise TypeError(
                'y should have the same shape as self.y_pred',
                ('y', y.type, 'y_pred', self.y_pred.type)
            )
        # y が正しい datatype かチェックする
        if y.dtype.startswith('int'):
            # T.neq オペレータは 0 と 1 のベクトルを返します。
            # ここで 1 は予測の誤りを表します。
            return T.mean(T.neq(self.y_pred, y))
        else:
            raise NotImplementedError()

そして関数 test_model と関数 validate_model を作成します、これはこの値を取得するために呼ぶことができます。すぐに見るように、validate_model は early-stopping 実装のキーです。（Early-Stopping 参照）これらの関数はミニバッチ index を取り、ミニバッチのサンプルに対して、モデルで誤分類した数を計算します。それらの唯一の違いは test_model はテストセットからのミニバッチを引き出し、validate_model は検証セットから引き出すことです。

    # Theano 関数をコンパイルします、これはミニバッチ上のモデによる誤りを計算します。
    test_model = theano.function(
        inputs=[index],
        outputs=classifier.errors(y),
        givens={
            x: test_set_x[index * batch_size: (index + 1) * batch_size],
            y: test_set_y[index * batch_size: (index + 1) * batch_size]
        }
    )

    validate_model = theano.function(
        inputs=[index],
        outputs=classifier.errors(y),
        givens={
            x: valid_set_x[index * batch_size: (index + 1) * batch_size],
            y: valid_set_y[index * batch_size: (index + 1) * batch_size]
        }
    )

Putting it All Together

最終的な成果は以下です。

（訳注：原文では全コードを引用。代わりに github code を参照。）

MNIST 数字を SGD ロジスティック回帰で分類することを学習するためには、DeepLearningTutorials フォルダの中からタイピングします :

python code/logistic_sgd.py

期待できる出力は次のような形です :

...
epoch 72, minibatch 83/83, validation error 7.510417 %
     epoch 72, minibatch 83/83, test error of best model 7.510417 %
epoch 73, minibatch 83/83, validation error 7.500000 %
     epoch 73, minibatch 83/83, test error of best model 7.489583 %
Optimization complete with best validation score of 7.500000 %,with test performance 7.489583 %
The code run for 74 epochs, with 1.936983 epochs/sec

Intel(R) Core(TM)2 Duo CPU E8400 @ 3.00 Ghz 上、コードはおよそ 1.936 epochs/sec で実行され、7.489% のテスト・エラーに 75 エポックで達します。GPU 上ではコードは殆ど 10.0 epochs/sec で動作します。このインスタンスのためにはバッチサイズ 600 を使用しました。

訓練モデルを使用した予測

sgd_optimization_mnist は、新しい検証エラーの最小値に達するたびに、モデルをシリアライズして pickle（訳注: pickle = python library の一つ）します。このモデルを reload して新しいデータのラベルを予測できます。予測関数はこれがどのようになされるかの例を示します。

def predict():
    """
    訓練されたモデルをどのようにロードしてラベルを予測するために使用するかの例です。
    """

    # 保存されたモデルをロードします
    classifier = pickle.load(open('best_model.pkl'))

    # predictor 関数をコンパイルします。
    predict_model = theano.function(
        inputs=[classifier.input],
        outputs=classifier.y_pred)

    # We can test it on some examples from test test
    dataset='mnist.pkl.gz'
    datasets = load_data(dataset)
    test_set_x, test_set_y = datasets[2]
    test_set_x = test_set_x.get_value()

    predicted_values = predict_model(test_set_x[:10])
    print("Predicted values for the first 10 examples in test set:")
    print(predicted_values)

[Footnotes]
[1] より小さなデータセットと単純なモデルのためには、より洗練された降下アルゴリズムがより効率的です。サンプルコード logistic_cg.py は、ロジスティック回帰タスクで SciPy の conjugate gradient solver を Theano と一緒にどのように使用するかを示します。

以上

月	火	水	木	金	土	日
	1	2	3	4	5	6
7	8	9	10	11	12	13
14	15	16	17	18	19	20
21	22	23	24	25	26	27
28	29	30	31