kerasでvgg16とGrad-CAMの実装による異常検出および異常箇所の可視化

#概要
画像認識による異常検出は最近工場現場の品質検査工程などで実用化が進められていますが、ディープラーニングの仕組みがブラックボックスとなっているため、AIが提示してくれた結果に不信感が生じやすいです。そのような課題を解消するために、Grad-CAMなどの判断根拠可視化手法が近年提案されました。
本稿では、kerasでvgg16モデルをファインチューニングし、DAGMデータセットの異常検知を試してみました。それから、Grad-CAMによる異常箇所の可視化も実装してみました。わりと良い結果が得られましたので、その手順と注意点をまとめて公開します。
特にGrad-CAMの実装にあたって、かなりエラーで苦しんでいたので、その注意点とコードをもしご参考になれれば嬉しいです。

#開発環境
Google Colaboratoryを使って、その便利さにめっちゃ感動されました。jupyter, python, kerasなどのライブラリーを一から導入する必要が無いため、ライブラリー間のバージョンによる互換性を全く意識せずにコードを書けて、機械学習の初心者にはとても優しい開発環境です。また、高速なGPUを無料で使えることが何より幸せです！具体的にどのGPUを使っているかは調べていませんが、筆者の経験だとNvidia Tesla K80よりも速いです。

使い方や詳しい説明はすでにいっぱい記事が書かれたので、本稿では割愛します。
以下のリンクよりご参照ください⬇️⬇️
[【秒速で無料GPUを使う】深層学習実践Tips on Colaboratory](https://qiita.com/tomo_makes/items/b3c60b10f7b25a0a5935)

#データセット
DAGM 2007というドイツで開催されたコンペで使われたデータセット

素材の模様に人為的につけられた欠陥の検出を目的とし、異常検出によく使われるデータセットです。
正常画像1000枚と異常画像150枚を1セットで、合計５セットあります。
![3.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/259387/66dc0400-1cc4-367c-eee4-66921dc06ebb.png)

データは以下のサイトより入手できます。
https://resources.mpi-inf.mpg.de/conference/dagm/2007/prizes.html
#実装
さてさて、ここからは実装のコードと注意点のご紹介です。
##Import
```ruby:vgg16_grad-cam.ipynb
from __future__ import print_function
import keras
from keras.applications import VGG16
from keras.models import Sequential, load_model, model_from_json
from keras import models, optimizers, layers
from keras.optimizers import SGD
from keras.layers import Dense, Dropout, Activation, Flatten
from sklearn.model_selection import train_test_split  
from PIL import Image 
from keras.preprocessing import image as images
from keras.preprocessing.image import array_to_img, img_to_array, load_img
from keras import backend as K 
import os
import numpy as np  
import glob  
import pandas as pd
import cv2
```
##Google Driverにマウント
```ruby:vgg16_grad-cam.ipynb
from google.colab import drive
drive.mount('/content/gdrive')
%cd ./gdrive/'My Drive'/"Colab Notebooks"
```
##データの前処理
DAGMのClass1データセットを使います。正常異常データを150枚ずつ取得し、事前にGoogle Driverの'/Colab Notebooks/DAGM/'にあるClass1とClass1_defの中にアップロードしておきます。

```ruby:vgg16_grad-cam.ipynb
num_classes = 2
folder = ["Class1","Class1_def"]                                     
image_size = 224
x = []
y = []
    
for index, name in enumerate(folder):
    dir = "./DAGM/" + name
    files = glob.glob(dir + "/*.png")    
    for i, file in enumerate(files):    
        image = Image.open(file)                       
        image = image.convert("RGB")    
        image = image.resize((image_size, image_size))
        data = np.asarray(image)        
        x.append(data)  
        y.append(index) 
        
x = np.array(x)   
y = np.array(y)  
 
x_train, x_test, y_train, y_test = train_test_split(x, y, test_size=0.2, random_state=111)
x_train = x_train.astype('float32')
x_test = x_test.astype('float32')
x_train /= 255
x_test /= 255
 
# ｙ　ラベルをワンホット表現に
y_train = keras.utils.to_categorical(y_train, num_classes)
y_test = keras.utils.to_categorical(y_test, num_classes)

print(x_train.shape[0], 'train samples')
print(x_test.shape[0], 'test samples')
```
##保存した重みの読み込み
事前に保存した学習済みモデルのパラーメーターの重みを読みこみます。初回は実行しないでください。

```ruby:vgg16_grad-cam.ipynb
model.load_weights('grad1_vgg16_weight_DAGM_C1.h5')
```
##vgg16モデル構築
vgg16モデルの紹介およびファインチューニングに関する記事が多く書かれているため、説明は省略します。
[Keras VGG16学習済みモデルでファインチューニングをやってみる
](http://cedro3.com/ai/keras-vgg16-fine-tuning/)[ディープラーニング実践入門 〜 Kerasライブラリで画像認識をはじめよう！](https://employment.en-japan.com/engineerhub/entry/2017/04/28/110000)

**ただし、Grad-CAMを実装するために、ここで１つ注意点があります！
モデルの構築にはSequentialモデルを使うと、Grad-CAMを実装するときはエラーが出まくってしまいます。本当の理由は未だにわかっていないのですが、筆者が考えている原因は以下の通りです。
・Grad-CAMの実装には、K.gradients()で各レイヤーのパラメーターの勾配を取得する必要があります。sequentialモデルを使うと、vgg16自体がモデルの１つのレイヤーになってしまうため（model.summary()でレイヤーの形を確認できる）、vgg16の畳み込み層のパラメーターの勾配を取得しようとすると、戻り値がnoneとなってしまいます。**
（もし間違っていれば、ぜひご指摘いただけるとありがたいです）



```ruby:vgg16_grad-cam.ipynb
vgg_conv = VGG16(weights='imagenet', include_top=False, input_shape=(image_size, image_size, 3))
last = vgg_conv.output

mod = Flatten()(last)
mod = Dense(1024, activation='relu')(mod)
mod = Dropout(0.5)(mod)
preds = Dense(2, activation='sigmoid')(mod)

model = models.Model(vgg_conv.input, preds)
model.summary()

epochs = 100
batch_size = 48
  
model.compile(loss='binary_crossentropy',
              optimizer=optimizers.SGD(lr=1e-4, momentum=0.9),
              metrics=['accuracy'])
```
```ruby:モデルのサマリー
_________________________________________________________________
Layer (type)                 Output Shape              Param #   
=================================================================
input_2 (InputLayer)         (None, 224, 224, 3)       0         
_________________________________________________________________
block1_conv1 (Conv2D)        (None, 224, 224, 64)      1792      
_________________________________________________________________
block1_conv2 (Conv2D)        (None, 224, 224, 64)      36928     
_________________________________________________________________
block1_pool (MaxPooling2D)   (None, 112, 112, 64)      0         
_________________________________________________________________
block2_conv1 (Conv2D)        (None, 112, 112, 128)     73856     
_________________________________________________________________
block2_conv2 (Conv2D)        (None, 112, 112, 128)     147584    
_________________________________________________________________
block2_pool (MaxPooling2D)   (None, 56, 56, 128)       0         
_________________________________________________________________
block3_conv1 (Conv2D)        (None, 56, 56, 256)       295168    
_________________________________________________________________
block3_conv2 (Conv2D)        (None, 56, 56, 256)       590080    
_________________________________________________________________
block3_conv3 (Conv2D)        (None, 56, 56, 256)       590080    
_________________________________________________________________
block3_pool (MaxPooling2D)   (None, 28, 28, 256)       0         
_________________________________________________________________
block4_conv1 (Conv2D)        (None, 28, 28, 512)       1180160   
_________________________________________________________________
block4_conv2 (Conv2D)        (None, 28, 28, 512)       2359808   
_________________________________________________________________
block4_conv3 (Conv2D)        (None, 28, 28, 512)       2359808   
_________________________________________________________________
block4_pool (MaxPooling2D)   (None, 14, 14, 512)       0         
_________________________________________________________________
block5_conv1 (Conv2D)        (None, 14, 14, 512)       2359808   
_________________________________________________________________
block5_conv2 (Conv2D)        (None, 14, 14, 512)       2359808   
_________________________________________________________________
block5_conv3 (Conv2D)        (None, 14, 14, 512)       2359808   
_________________________________________________________________
block5_pool (MaxPooling2D)   (None, 7, 7, 512)         0         
_________________________________________________________________
flatten_2 (Flatten)          (None, 25088)             0         
_________________________________________________________________
dense_3 (Dense)              (None, 1024)              25691136  
_________________________________________________________________
dropout_2 (Dropout)          (None, 1024)              0         
_________________________________________________________________
dense_4 (Dense)              (None, 2)                 2050      
=================================================================
Total params: 40,407,874
Trainable params: 40,407,874
Non-trainable params: 0
_________________________________________________________________
```
最後のプーリング層の１つ前のレイヤーが「block5_conv3」であることにご注目を、後ほどのGrad-CAMの実装で使われます。

##学習
```ruby:vgg16_grad-cam.ipynb
history = model.fit(x_train, y_train,
                    batch_size=batch_size,
                    epochs=epochs,
                    validation_data=(x_test, y_test),
                    shuffle=True)
```

##学習済みモデルの評価とaccuracy&lossの推移
```ruby:vgg16_grad-cam.ipynb
scores = model.evaluate(x_test, y_test, verbose=1)
print('Test loss:', scores[0])
print('Test accuracy:', scores[1])
 
### Plot accuracy & loss
import matplotlib.pyplot as plt 

acc = history.history["acc"]
val_acc = history.history["val_acc"]
loss = history.history["loss"]
val_loss = history.history["val_loss"]
epochs = range(1, len(acc) + 1)
 
#plot accuracy
plt.plot(epochs, acc, label = "Training acc" )
plt.plot(epochs, val_acc, label = "Validation acc")
plt.title("Training and Validation accuracy")
plt.legend()
plt.show()

#plot loss
plt.plot(epochs, loss,  label = "Training loss" )
plt.plot(epochs, val_loss, label = "Validation loss")
plt.title("Training and Validation loss")
plt.legend()
plt.show()
```
![スクリーンショット 2019-05-27 14.49.31.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/259387/d24d06ce-9c55-06f8-4643-3bc82722c45b.png)
100epoch回したら精度が98.3%となり、良い感じの結果が得られましたので、こちらの学習済みモデルで、Grad-CAMを実装します。



##学習済みモデルの重みの保存
```ruby:vgg16_grad-cam.ipynb
model.save_weights('grad_vgg16_weight_DAGM_C1.h5')

```
##Grad-CAMの実装
コード自体は下のリンクを参照していますので、詳細の説明は譲りますが、ここで何をやっているかを簡単に説明します。
**出力層に最も近いレイヤーを抽出し（vgg16の場合は「block5_conv3」）、このレイヤーのパラメーターの勾配を元に、画像の各領域が最終の出力に与える影響を計算し、その影響度の高さをヒートマップで表現します。**
[kerasでGrad-CAM 自分で作ったモデルで](https://qiita.com/haru1977/items/45269d790a0ad62604b3)
[ディープラーニングの注視領域の可視化
](https://qiita.com/bele_m/items/a7bb15313e2a52d68865)

```ruby:vgg16_grad-cam.ipynb
K.set_learning_phase(1) #set learning phase

def Grad_Cam(input_model, pic_array, layer_name):

    # 前処理
    pic = np.expand_dims(pic_array, axis=0)
    pic = pic.astype('float32')
    preprocessed_input = pic / 255.0

    # 予測クラスの算出
    predictions = input_model.predict(preprocessed_input)
    class_idx = np.argmax(predictions[0])
    class_output = input_model.output[:, class_idx]

    #  勾配を取得
    conv_output = input_model.get_layer(layer_name).output   # layer_nameのレイヤーのアウトプット
    grads = K.gradients(class_output, conv_output)[0]  # gradients(loss, variables) で、variablesのlossに関しての勾配を返す
    gradient_function = K.function([input_model.input], [conv_output, grads])  # input_model.inputを入力すると、conv_outputとgradsを出力する関数

    output, grads_val = gradient_function([preprocessed_input])
    output, grads_val = output[0], grads_val[0]

    # 重みを平均化して、レイヤーのアウトプットに乗じる
    weights = np.mean(grads_val, axis=(0, 1))
    cam = np.dot(output, weights)

    # 画像化してヒートマップにして合成
    cam = cv2.resize(cam, (224, 224), cv2.INTER_LINEAR) 
    cam = np.maximum(cam, 0) 
    cam = cam / cam.max()

    jetcam = cv2.applyColorMap(np.uint8(255 * cam), cv2.COLORMAP_JET)  # モノクロ画像に疑似的に色をつける
    jetcam = cv2.cvtColor(jetcam, cv2.COLOR_BGR2RGB)  # 色をRGBに変換
    jetcam = (np.float32(jetcam) + pic / 2)   # もとの画像に合成
    return jetcam
```
##テスト画像を指定
```ruby:vgg16_grad-cam.ipynb
pic_array = img_to_array(load_img('DAGM/Class1_def/12.png', target_size=(224, 224)))
pic = pic_array.reshape((1,) + pic_array.shape)
array_to_img(pic_array)
```
![12.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/259387/3b008dfb-d006-8421-2f2c-b0b812c91c01.png)
まずは異常系画像を試してみます。
若干見辛いのですが、画像の上部にシミがついていることがわかります。

##異常箇所の可視化
```ruby:vgg16_grad-cam.ipynb
picture = Grad_Cam(model, pic_array, 'block5_conv3')
picture = picture[0,:,:,]
array_to_img(picture)
```
![12_.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/259387/7b909e8f-df39-e684-ba86-d616cce8ef26.png)
完璧に異常箇所を示していることを確認できます。
ヒートマップの赤のところは、今回訓練したモデルがこの画像を異常画像に分類している根拠であると示しています。
##正常画像の判別
|  元画像  |  判別根拠  |
| --- | --- | 
| ![seijyou.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/259387/234b9fac-bb5d-166c-2b05-20421df23d1d.png)| ![sei.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/259387/e5bd64c0-c823-2cfd-ab7a-abd6ebedd298.png)| 
上記の出力結果からわかるように、正常画像を判定する場合は、画像全体を見渡して、異常箇所や欠陥があるかを確認していることが明らかになっています。
#終わりに
以上でvgg16とGrad-CAMの実装手順をまとめました。
ディープラーニング入門してからまだ1ヶ月程度なので、もし間違いがあればぜひご指摘ください！また、ご意見やご質問などありましたら、ぜひコメントで気軽に教えてくれると嬉しいです！
#参考
[【秒速で無料GPUを使う】深層学習実践Tips on Colaboratory](https://qiita.com/tomo_makes/items/b3c60b10f7b25a0a5935)
[Keras VGG16学習済みモデルでファインチューニングをやってみる
](http://cedro3.com/ai/keras-vgg16-fine-tuning/)[ディープラーニング実践入門 〜 Kerasライブラリで画像認識をはじめよう！](https://employment.en-japan.com/engineerhub/entry/2017/04/28/110000)
[kerasでGrad-CAM 自分で作ったモデルで](https://qiita.com/haru1977/items/45269d790a0ad62604b3)
[ディープラーニングの注視領域の可視化
](https://qiita.com/bele_m/items/a7bb15313e2a52d68865)