Dropoutによる近似ベイズ推論でDeep Learningの不確かさを評価する

Deep learningの推定結果の不確かさってどうやって評価するのか疑問を持っていました。
Dropoutを使ったサンプリングをすることで不確かさ評価をできるということなので、kerasで実装して検証してみました。

以下の検証に関するコードはgithubにあげてあります。
github.com

1. 本記事の概要

記事の目的
Dropoutによるベイズ推論と不確かさ評価の概要
不確かさ評価の検証
まとめ

2. 記事の目的

Dropoutを使ったDLの不確かさ評価の概要をつかむ。
kerasでの実装方法の記録。

3. Dropoutによるベイズ推論と不確かさ評価の概要

Dropoutによるベイズ推論

Deep Learning（以降DL）を使うことで回帰問題や分類問題の答えを精度よく推定することができます。しかし、精度が良いと言っても線形回帰等と同じく、その推定結果はかならず不確かさを含みます。
ベイズ推論を使えば推定結果がどのようにばらつくのかがわかるので不確かさを評価できるのですが、DLの式は複雑なのでベイズ推論を適用することは難しかったとか。（詳しくはよく理解していないです笑）

以下の論文では、DLの過学習を防ぐ手法であるDropoutを使って推論(推定)を行うことが近似的にベイズ推論になっていることが証明されています。
https://arxiv.org/pdf/1506.02142.pdf
（こういう機械学習の手法が統計的・数学的にこんな意味をもっているんだぜーってやつ、めっちゃわくわくしますよね。）

詳しくは以下の日本語ブログでも詳細されています。
ニューラルネットへのベイズ推定 - Bayesian Neural Network - nykergoto’s blog
Dropoutによる近似ベイズ推論について – 戦略コンサルで働くデータサイエンティストのブログ

僕が理解した範囲で、できるだけ平易な言葉でイメージを表現しておきます。

ベイズ推論の枠組みをDLに適用する場合、DLの各重みは一定値でなく分布を持つと考えます。これにより推定結果も分布（不確かさ）を持ちます。
各重みの分布から各重みを何回もサンプリングし、何回も推定を行うことで推定結果の分布を求めることができます。しかし、そもそも各重みの分布を求めるのはすごくしんどい。。
Dropoutは各重みをランダムに0にする効果があります。通常はDropoutは学習時にしか使いませんが、Dropoutを推論時にも使うことで各重みの分布からサンプリングしていることになります。
そのため、Dropoutを使って推定を何回も繰り返すことで、推定結果の分布（不確かさ）を求めることができます。

不確かさ評価手順

Dropoutを使うことで、以下の手順で不確かさを評価できます。

Dropout、L2正則化を使ってDLを学習させます。
同じ入力xに対して複数回DLの推定を行います。このとき、学習時と同じ率のDropoutを使います。
推定した結果をヒストグラムにする、分散をもとめる、エントロピーを求めるなどして不確かさを評価します。ヒストグラムの裾が広いほど、分散が大きいほど、エントロピーが大きいほど不確かさが大きいと言えます。

原論文では上記のようにDropoutを使いながらサンプリングすることをモンテカルロドロップアウト(MC dropout)と呼んでいます。
分散は、「普通の分散」+「ドロップアウト率、L2正則化係数、サンプリング数等に関連した量」で求められます。詳しくは原論文に書かれています。

この不確かさ評価は、通常のDLの学習・推論とやっていることは大差ありません。異なるところは、推論時にもDropoutを使う、同じ入力に対して複数回推定を行う、というところぐらいです。そのため、keras等の既存のDLライブラリを使うことができるのでお手軽です。
ただし、複数回推定を行わないといけないので、推定に時間がかかります。

4. 不確かさ評価の検証

ここでは、DropoutによるDLの不確かさ評価の検証をしてみます。
ゴールは、不確かさ評価のイメージをつかむことです。

検証内容

検証の設定は以下のとおりです。二値分類問題で不確かさがどのようなふるまいをするか確認します。

対象は犬猫画像の二値分類問題。
- 画像データとしてcifar10を使用。
- モデルにはCNNを使用。構造は、入力画像→畳み込み層→…→全結合層→Dropout→全結合層→Dropout→...→出力(sigmoid)
- 犬を1、猫を0とラベルをつけてCNNを通常どおり学習させます。できあがったCNNは犬である確率を出力します。
学習済みCNNで推論時にもDropoutを使い、推定値の標準偏差(原論文で提案された分散ではなく、普通の方法で算出)を計算します。また、推定値の平均も計算します。
- 1つの入力に対して100回モンテカルロドロップアウトサンプリングを実施。
- 学習に使っていないラベル（飛行機、自動車、鳥、鹿、カエル、馬、船、トラック）についても評価を実施。

f:id:st1990:20190731002809p:plain — モデル構造

実装

実装にはpython, kerasを使用しました。

学習は通常どおり行います。

import keras
from keras.preprocessing.image import ImageDataGenerator
from keras.models import Sequential
from keras.layers import Dense, Dropout, Activation, Flatten, BatchNormalization
from keras.layers import Conv2D, MaxPooling2D
from keras import regularizers
import numpy as np

class ClassifierCnn:
    """
    classify dog and cat in cifar10
    """
    def __init__(self):
        return

    def built_model(self, dropout_rate, l2, input_shape=None):
        # constants
        self.DROPOUT_RATE = dropout_rate //0.5
        self.L2 = l2 //0.0001

        if input_shape is None:
            # assume cifar10 image
            input_shape = (32, 32, 3)

        # model structure
        model = Sequential()
        model.add(Conv2D(32, (3, 3), padding='same', 
                         kernel_regularizer=regularizers.l2(self.L2), input_shape=input_shape))
        model.add(Activation('relu'))
        model.add(Conv2D(32, (3, 3), kernel_regularizer=regularizers.l2(self.L2)))
        model.add(Activation('relu'))
        model.add(MaxPooling2D(pool_size=(2, 2)))

        model.add(Conv2D(64, (3, 3), padding='same', kernel_regularizer=regularizers.l2(self.L2)))
        model.add(Activation('relu'))
        model.add(Conv2D(64, (3, 3), kernel_regularizer=regularizers.l2(self.L2)))
        model.add(Activation('relu'))
        model.add(MaxPooling2D(pool_size=(2, 2)))

        model.add(Flatten())
        model.add(Dense(512, kernel_regularizer=regularizers.l2(self.L2)))
        model.add(Activation('relu'))
        model.add(Dropout(self.DROPOUT_RATE))
        model.add(Dense(64, kernel_regularizer=regularizers.l2(self.L2)))
        model.add(Activation('relu'))
        model.add(Dropout(self.DROPOUT_RATE))
        model.add(Dense(1, kernel_regularizer=regularizers.l2(self.L2)))
        model.add(Activation('sigmoid'))

        self.model = model
        self.model.summary()
        return

    def train_model(self, x_train, y_train, x_test, y_test, epochs, batch_size):
        # compile
        self.model.compile(loss='binary_crossentropy', optimizer='adam', metrics=['accuracy'])

        # datagen
        datagen = ImageDataGenerator(
            featurewise_center=False,  # set input mean to 0 over the dataset
            samplewise_center=False,  # set each sample mean to 0
            featurewise_std_normalization=False,  # divide inputs by std of the dataset
            samplewise_std_normalization=False,  # divide each input by its std
            zca_whitening=False,  # apply ZCA whitening
            zca_epsilon=1e-06,  # epsilon for ZCA whitening
            rotation_range=20,  # randomly rotate images in the range (degrees, 0 to 180)
            # randomly shift images horizontally (fraction of total width)
            width_shift_range=0.1,
            # randomly shift images vertically (fraction of total height)
            height_shift_range=0.1,
            shear_range=0.,  # set range for random shear
            zoom_range=0.1,  # set range for random zoom
            channel_shift_range=0.,  # set range for random channel shifts
            # set mode for filling points outside the input boundaries
            fill_mode='nearest',
            cval=0.,  # value used for fill_mode = "constant"
            horizontal_flip=True,  # randomly flip images
            vertical_flip=False,  # randomly flip images
            # set rescaling factor (applied before any other transformation)
            rescale=None,
            # set function that will be applied on each input
            preprocessing_function=None,
            # image data format, either "channels_first" or "channels_last"
            data_format=None,
            # fraction of images reserved for validation (strictly between 0 and 1)
            validation_split=0.0)
        datagen.fit(x_train)

        # fit
        self.model.fit_generator(datagen.flow(x_train, y_train, batch_size=batch_size),
                            epochs=epochs,
                            steps_per_epoch=int(x_train.shape[0] / batch_size),
                            validation_data=(x_test, y_test),
                            )

推論時にもDropoutを使うため、読み込んだモデルのDropout層をDropout()(training=True)と入れ替えています。
※参考：【Keras】 Kerasで訓練時だけじゃなく予測時にもドロップアウトを使う - Qiita
こうすることで推論時にもDropoutが使われるため、複数回model.predictを実行することでモンテカルロドロップアウトサンプリングができます。

from keras.layers import Dropout
from keras.models import Model

class MontecarloDropout:
    def __init__(self):
        return

    def build_model(self, model_file_path):
        """
        keras modelからMontecarloDropoutに対応したモデルを作成
        build monte carlo dropout model base on keras_model.
        """

        model = self.__load_model(model_file_path)

        # change dropout layer to dropout layer that can use dropout in inference.
        # ドロップアウト層を推論時にもドロップアウトできるドロップアウト層に変更する。
        for ily, layer in enumerate(model.layers):
            # input layer
            if ily == 0:
                input = layer.input
                h = input
            # is dropout layer ?
            if 'dropout' in layer.name:
                # change dropout layer
                h = Dropout(layer.rate)(h, training=True)
            else:
                h = layer(h)

        self.model = Model(input, h)
        return

    def md_predict(self, xs, sampling_num):
        """
        predict with using monte carlo dropout sampling.
        return prediction average, std

        xs : input sample array. xs = x0, x1, x2, ...
        """
        pre_ys = []
        for ismp in range(sampling_num):
            pre_y = self.model.predict(xs)
            pre_ys.append(pre_y)
        pre_ys = np.array(pre_ys)

        # calculate ave, std
        pre_ave = np.average(pre_ys, axis=0)
        pre_std = np.std(pre_ys, axis=0)

        return pre_ave, pre_std

他の部分はgithubを参照してください。
github.com

検証結果

作成したモデルの精度

モデルの精度(犬猫分類の正解率)は以下のとおりになりました。MC Dropoutはアンサンブルのような効果があるので精度が良くなるかと思っていましたが、大差ないですね。

CNN: 　　　学習データ 92.3%, テストデータ 86.3%
MC Dropout: 学習データ 92.2%, テストデータ 86.3%
※MC Dropoutの推定値は、入力1つに対して100回モンテカルロドロップアウトして得られた100個の推定値を平均したものです。

MC Dropoutでの推定値 vs CNNの推定値（学習データ）

MC Dropoutでの推定値 vs CNNの推定値（テストデータ）

犬猫画像の分類結果と不確かさ

学習に使った犬猫画像の分類結果とMC Dropoutで求めた不確かさ(標準偏差)の関係をみていきます。
推定値の標準偏差に対する推定値のグラフ、推定値の標準偏差のヒストグラムを以下に示します。
推定値が0.5に近いほど標準偏差が大きく、0や1では不確かさがほぼないという結果になっています。推定値=0.5は「犬か犬以外の確率が半々」ということを表しているので、どっちかよくわからんという状況ということで不確かさが大きくなるのでしょうか。

f:id:st1990:20190731000343p:plain — 推定値の標準偏差 vs 推定値（学習データ）

f:id:st1990:20190731000419p:plain — 推定値の標準偏差 vs 推定値（テストデータ）

f:id:st1990:20190731002322p:plain — 推定値の標準偏差 vs 推定値（犬のみ）（学習データ）

f:id:st1990:20190731002353p:plain — 推定値の標準偏差 vs 推定値（猫のみ）（学習データ）

f:id:st1990:20190731000437p:plain — 推定値の標準偏差のヒストグラム

テストデータのROC曲線を下図に示します。
「threshold=y」ラベルは、犬か猫か判定する閾値を単純に推定値としてROC曲線を描いたものです。推定値>閾値であれば犬と判定されます。
「threshold=0.5+astd」ラベルは、閾値を0.5+astd (-∞<a<∞)としてROC曲線を描いたものです。推定値>0.5+a*stdであれば犬と判定されます。推定値が0.5に対して不確かさを含めてどの程度余裕があるかを表しています。
結果を見ると…ほぼ同じですね。

f:id:st1990:20190816021850p:plain — ROC曲線(テストデータ)

同程度の推定値で、標準偏差が小さい場合と大きい場合の画像を以下に示します。標準偏差の大きい画像は人が見ても判別が難しいものであることを期待していましたが、傾向がよくわからん。。笑

f:id:st1990:20190731001405p:plain — 推定値0.48～0.52で標準偏差が小さい方から25位までの画像

f:id:st1990:20190731001542p:plain — 推定値0.48～0.52で標準偏差が大きい方から25位までの画像

f:id:st1990:20190731001608p:plain — 推定値0.68～0.72で標準偏差が小さい方から25位までの画像

f:id:st1990:20190731001629p:plain — 推定値0.68～0.72で標準偏差が大きい方から25位までの画像

f:id:st1990:20190731001651p:plain — 推定値0.78～0.82で標準偏差が小さい方から25位までの画像

f:id:st1990:20190731001708p:plain — 推定値0.78～0.82で標準偏差が大きい方から25位までの画像

学習に使っていないラベルでの不確かさ評価結果

学習に使っていないラベル（飛行機、自動車、鳥、鹿、カエル、馬、船、トラック）についても不確かさ評価を実施しました。
推定値の標準偏差のヒストグラム、推定値の標準偏差に対する推定値のグラフを以下に示します。
0.5付近の推定値が増えるので、ヒストグラムのピークとなる標準偏差の値が大きくなったことがわかります。しかし、標準偏差vs推定値の分布は犬猫画像と変わりなく見えます。学習に使っていないラベルは学習データのデータ分布外なので、推定値が0や1に近くても標準偏差が大きくなるかと思っていたのですが。。