Deep Learningの学習結果（重み）はブラックボックスで、隠れ層のユニット（特に深い層の！）が一体何を学習したのかがよくわからないと長年言われてきた。しかし、今回紹介する方法を使うとニューラルネットが何を学習したのか目で見える形で表現できる。

畳み込みニューラルネットで学習したフィルタの可視化というと以前やったように学習した第1層のフィルタの重みを直接画像として可視化する方法がある。

しかし、畳み込みフィルタのサイズは基本的に数ピクセル（MNISTの例では5x5ピクセル程度）のとても小さな画像なのでこれを直接可視化しても何が学習されたか把握するのはとても難しい。たとえば、MNISTを学習した畳み込みニューラルネット（2016/11/20）のフィルタを可視化しても各フィルタがどの方向に反応しやすいかがわかる程度だ。

f:id:aidiary:20161120200900p:plain:w300

各フィルタが何を学習したかを可視化する別のアプローチとして各フィルタの出力を最大化するような入力画像を生成する手法が提案された。この生成画像はニューラルネットの入力画像と同じサイズなのでフィルタが最終的に何を学習したのかより把握しやすいというメリットがある。

今回は畳み込みニューラルネットの一種であるVGG16を対象に学習したフィルタを可視化してみたい。あとでMNISTやCIFAR-10を学習したCNNや最新のResNetでもやってみるつもり。

基本的に本家のKeras Blogの How convolutional neural networks see the world を参考にした。しかし、この記事は

keras.applications.vgg16モジュールが導入される前に書かれている
正則化などの工夫が入っておらず生成される画像が美しくない

という問題がある。

そこでこの記事ではKerasのvgg16モジュールを使って書き換えるとともに、よりきれいな画像が生成されるようにいくつか工夫した。まあ少しは差別化しないとね（笑）

リポジトリ：dream

ナイーブな方法

基本的なアプローチは、畳み込みニューラルネットの指定したフィルタの出力を最大にする入力画像を勾配法を用いて更新することである。数式で書くと

$x \leftarrow x + \eta \frac{\partial a_i (x)}{\partial x}$

となる。ここで、 $x$ は入力画像、 $\eta$ は学習率、 $a_i (x)$ は画像 $x$ を入力したときのi番目のフィルタの出力（activation）だ。

ここで、偏微分がいわゆる勾配を表し、入力画像 $x$ をちょっと変えたときにフィルタの出力 $a_i (x)$ がどれくらい変化するかを表している。入力画像 $x$ をこの勾配方向に徐々に移動する勾配上昇法（gradient ascent）を用いてフィルタの出力を最大化する入力画像を求めている。

この式はニューラルネットの重みの更新式と比較するとわかりやすい。ニューラルネットの重みの更新式は

$w \leftarrow w - \eta \frac{\partial L(w)}{\partial w}$

これはよく目にする式だ。ここで、 $w$ は重み、 $\eta$ は学習率、 $L(w)$ は損失だ。

ここで、偏微分は勾配を表し、重み $w$ をちょっと変えたときに損失 $L(w)$ がどれくらい変化するかを表している。重み $w$ をこの勾配の負の方向に徐々に移動する勾配降下法（gradient descent）を用いて損失を最小化する重みを求めている。

こんな感じで比較してみると両者の基本原理はよく似ていることがわかる。

ニューラルネットの重み更新では重みの初期値としてランダムな値を使った。フィルタの可視化でも同様に入力画像の初期値として下記のようなランダムな画像を使う。一方、次回紹介するDeep Dreamはこの画像の初期値として任意の画像を使う。

f:id:aidiary:20170215230458p:plain:w300

Kerasで書くと下のようになる。

def visualize_filter(layer_name, filter_index, num_loops=200):
    """指定した層の指定したフィルタの出力を最大化する入力画像を勾配法で求める"""
    # 指定した層の指定したフィルタの出力の平均
    activation_weight = 1.0
    if layer_name == 'predictions':
        # 出力層だけは2Dテンソル (num_samples, num_classes)
        activation = activation_weight * K.mean(layer.output[:, filter_index])
    else:
        # 隠れ層は4Dテンソル (num_samples, row, col, channel)
        activation = activation_weight * K.mean(layer.output[:, :, :, filter_index])

    # 層の出力の入力画像に対する勾配を求める
    # 入力画像を微小量変化させたときの出力の変化量を意味する
    # 層の出力を最大化したいためこの勾配を画像に足し込む
    grads = K.gradients(activation, input_tensor)[0]

    # 正規化トリック
    # 画像に勾配を足し込んだときにちょうどよい値になる
    grads /= (K.sqrt(K.mean(K.square(grads))) + K.epsilon())

    # 画像を入力して層の出力と勾配を返す関数を定義
    iterate = K.function([input_tensor], [activation, grads])

    # ノイズを含んだ画像（4Dテンソル）から開始する
    x = np.random.random((1, img_height, img_width, 3))
    x = (x - 0.5) * 20 + 128

    # 勾配法で層の出力（activation_value）を最大化するように入力画像を更新する
    cache = None
    for i in range(num_loops):
        activation_value, grads_value = iterate([x])
        # activation_valueを大きくしたいので画像に勾配を加える
        step, cache = rmsprop(grads_value, cache)
        x += step
        print(i, activation_value)

    # 画像に戻す
    img = deprocess_image(x[0])

    return img

この関数の二つの引数は畳込みニューラルネットのどの層の何番目のフィルタを最大化したいかを指定する。畳み込みニューラルネット内の畳込み層または出力層を指定できる。たとえば、出力層の65番目のフィルタを最大化したいときは下のように指定する。

visualize_filter('predictions', 65)

層の名前はmodel.summary()で表示される。predictionsはVGG16の出力層に付けられた名前である。VGG16の出力層は1000ユニットから成るが65番目のユニットはsea snake（ウミヘビ）クラスを意味する。1個だけではつまらないので16個まとめて画像化すると下のようになる。左上にフィルタ番号に対応するクラス名を表示した。

f:id:aidiary:20170215230857p:plain

まあこの可視化画像を見てもはっきり言ってよくわからない。ウミヘビは何かウネウネしてるぞ、ゴリラに何か目があるぞとかそれくらい（笑）

でもこれらの画像をニューラルネットに入力するとウミヘビである確率が99%、ゴリラである確率が99%ととんでもない結果が出てくる・・・人間にはまったくウミヘビやゴリラと認識できない画像なのに、ニューラルネットは99%の確率でウミヘビ・ゴリラと判定してしまうのだ。ここら辺は参考文献で挙げた「ニューラルネットは簡単にだませる」という論文を参照。

このナイーブな方法ではあまりはっきりとした画像が生成できなかったが、いくつかの改良を加えるとよりはっきりとした鮮やかな画像が生成できることが知られている。そのために自然な画像らしさを表す事前知識（natural image priors）を正則化項として導入すればよい。

VGG16の出力をsoftmaxからlinearに

正則化項を導入する前に一つ重要な改善を入れよう。VGG16の出力層の出力は分類のためにsoftmax活性化関数を適用しており、1000クラスの合計が1.0になるように正規化されている。この場合、あるクラスの出力を最大化するためには (1) そのクラスの出力を上げる以外に、(2) 他のクラスの出力を下げるという方法もとれてしまう。この依存関係で最適化が惑わされてうまくいかないという指摘がある（Simonyan 2013）。そこで、最初の工夫としてsoftmaxをlinearにしてみた。

VGG16のsoftmax活性化関数をかける前の出力（つまり線形出力）が取れれば簡単なのだがKerasにそのようなプロパティは用意されていない。もしVGG16のsoftmax活性化関数が

model.add(Dense(1000))
model.add(Activation('softmax'))

のように独立した層として追加されていれば、一つ前のsoftmaxをかける前の出力は簡単に取れたのだが、VGG16は下のように層ではなく引数として実装されていた。

model.add(Dense(1000, activation='softmax'))

これだとsoftmaxをかける前の出力は取り出せない・・・仕方ないのでKerasのvgg16.pyをコピーして一部修正したVGG16クラスを用意した。121行目だけ下のように書き換えればOK。学習した重みはそのまま使える。

if include_top:
    # Classification block
    x = Flatten(name='flatten')(x)
    x = Dense(4096, activation='relu', name='fc1')(x)
    x = Dense(4096, activation='relu', name='fc2')(x)
    # ここの活性化関数をsoftmaxからlinearに置き換え
    x = Dense(1000, activation='linear', name='predictions')(x)

この工夫を入れてもう一度可視化すると下のようにコントラストが少し強くなる。だけどまだ何が何だかよくわからない。

f:id:aidiary:20170215231801p:plain

Lpノルム正則化の導入

次の工夫はLpノルム正則化である。Lpノルムの定義は

$||x||_p = (|x_1|^p + |x_2|^p + \cdots + |x_n|^p)^\frac{1}{p}$

である。pを $\infty$ に発散させると最大ノルムになる。

$||x||_\infty = \max (|x_1|, |x_2|, \cdots, |x_n|)$

Kerasで実装すると下のようになる。

def normalize(img, value):
    return value / np.prod(K.int_shape(img)[1:])

# Lpノルム正則化項
# 今回の設定ではactivationは大きい方がよいため正則化のペナルティ項は引く
p = 6.0
lpnorm_weight = 10.0
if np.isinf(p):
    lp = K.max(input_tensor)
else:
    lp = K.pow(K.sum(K.pow(K.abs(input_tensor), p)), 1.0 / p)
activation -= lpnorm_weight * normalize(input_tensor, lp)

ここで、Lpノルムの $x$ は生成画像であることに注意。またactivationは最大化を目指しているためペナルティである正則化項は差し引くことにも注意。noramlize()は画像のサイズによらないようにする正規化である。

この正則化項は画像のノルムを求めているのでピクセル値が極端に小さな値や大きな値になったときにペナルティをかけるような正則化だと考えられる。そのようなピクセル値が極端な値になる画像は自然ではないという事前知識を表している。

Lpノルムのpや正則化項の重み（lpnorm_weight）はハイパーパラメータなのだが参考文献を参考に設定した。いくつか異なる設定を試してみるとよいかも。Lpノルム正則化を入れると下の画像が生成された。あまり効果ない？

f:id:aidiary:20170215231957p:plain

Total Variation正則化の導入

もう一つTotal Variation正則化を導入しよう。これは効果絶大だった！画像 y のTotal Variationは

f:id:aidiary:20170216215852p:plain:w400

で定義される。つまり、Total Variationが大きい画像ほど隣通しのピクセルの差分が大きな画像と解釈できる。Kerasで書くと下のようになる（beta = 1のとき式と一致）。

# Total Variationによる正則化
beta = 2.0
tv_weight = 10.0
a = K.square(input_tensor[:, 1:, :-1, :] - input_tensor[:, :-1, :-1, :])
b = K.square(input_tensor[:, :-1, 1:, :] - input_tensor[:, :-1, :-1, :])
tv = K.sum(K.pow(a + b, beta / 2.0))
activation -= tv_weight * normalize(input_tensor, tv)

配列のスライスがややこしいが具体的な行列（5x5くらいの）で計算してみると隣通しのピクセルの差を取っているのがすぐにわかる。aは画像の行方向の差分でbは画像の列方向の差分を意味する。

ここでは、Total Variationが大きいほどペナルティがかかるようになっているためTotal Variationを小さく（=ピクセル間の差を小さく）する方向に正則化がはたらく。つまり、Total Variation正則化は隣通しのピクセルの輝度差が小さい滑らかな画像にする正則化であると考えられる。

この正則化を導入すると下の画像のようにきれいな画像が出てくる。ここまでくると元のクラスが何なのかうっすらと見えてくる！ナイーブな方法に比べると正則化の導入によってかなり改善しているのがわかる。

f:id:aidiary:20170215232009p:plain