AutoEncoderでシンプルな音声の異常検知(Keras/TensorFlow)

#はじめに
AutoEncoderで音のデータでの異常検知をすることを目標に、AutoEncoderを実装してみました。ピー音（規制音）のようなシンプルな音を学習させて、違う音との差分を検出することで異常判定したいと思います。勉強中・試行錯誤しながらなので、色々指摘いただけたら嬉しいです。

こんなユーザが対象です。

* MNISTではないデータ（音声）を扱って、異常検知したい
* AutoEncoderでの異常検知を実装しながら理解を深めたい

##環境

* Anaconda3(Python3.6)
* Keras/Tensorflow
* logfilterbank

## AutoEncoder(オートエンコーダー)とは
オートエンコーダーは、機械学習において、ニューラルネットワークを使用した次元圧縮のためのアルゴリズムで、入力層と出力層を同じデータで学習させることで、入力されたデータを出力層で再現できる生成モデルです。

[オートエンコーダ(Wikipedia)](https://ja.wikipedia.org/wiki/%E3%82%AA%E3%83%BC%E3%83%88%E3%82%A8%E3%83%B3%E3%82%B3%E3%83%BC%E3%83%80)

今回は、その特性を利用して入力されたデータが再現できるのであれば正常、そうでないなら異常であると判定することで異常検知をしたいと思います。

## 簡単な流れ
1. 音声データを準備します（正常、異常）
2. フィルタバンクで音声から特徴次元抽出
3. AutoEncoderネットワークを実装
4. 学習と評価

## 1. 音声WAVデータを準備します（正常、異常）

まず、音声WAVデータを準備します。できるだけシンプルな、何か反復するような音でまずは試行したいと考え、今回は、ピー音（規制音）を準備しました。規制音は1kHzの正弦波の音源(10分音源)です。異常音データとしては、上記のピー音の最初の5分を2kHzにしたもの(10分音源)を作りました。下図を見てもらえれば分かるように波形やスペクトルが正常と異常で異なることが分かると思います。

* 正常音データ信号
![正常音データ波形](https://storage.googleapis.com/permanent_koizumi/20181229_signal(regularity).png "正常音データ波形")

* 正常音データ周波数スペクトル（FFT）
![正常音データ周波数スペクトル](https://storage.googleapis.com/permanent_koizumi/20181229_1000Hz_regurality.png "正常音データ周波数スペクトル")

* 異常音データ信号
![異常音データ波形](https://storage.googleapis.com/permanent_koizumi/20181229_signal(anomaly).png "異常音データ波形")

* 異常音データ周波数スペクトル（FFT）
![異常音データ周波数スペクトル](https://storage.googleapis.com/permanent_koizumi/20181229_2000Hz_anomaly.png "異常音データ周波数スペクトル")


## 2.フィルタバンクで音声から複数の特徴データ抽出
以下のサイトのpypthon-speech-featuresライブラリのフィルタバンクを利用することで、簡単に入力音声を複数の特徴データに分割することができます。
今回は、ウィンドウ幅を10msecとして、かつ、40次元に分割をしました。

[https://python-speech-features.readthedocs.io/en/latest/](https://python-speech-features.readthedocs.io/en/latest/)

```Python:
from python_speech_features import logfbank
import scipy.io.wavfile as wav

# 音声データをpython speech featuresで配列変換
(rate,sig) = wav.read("20181229_train_sound(regularity).wav")
fbank_feat = logfbank(sig,rate,winlen=0.01,nfilt=40)
```



## 3.AutoEncoderモデルを実装

まず、学習データとテストデータは60000データくらいあったため、前半30000でトレーニング、後半30000程度でテストとしておきました。

```Python:
n_traintest_sep = 30000

X_trn = fbank_feat[:n_traintest_sep][:]
X_tst = fbank_feat[n_traintest_sep:][:]
```

次に、データを正規化します。

```Python:
# 値の範囲を[0,1]に正規化
scaler = MinMaxScaler()
X_trn = scaler.fit_transform(X_trn)
X_tst = scaler.transform(X_tst)
```

AutoEncoderモデルを実装します。
[本家のサイトのDeepAutoEncoderサンプルを参考にしました。](https://blog.keras.io/building-autoencoders-in-keras.html)
Sequentialモデルを利用し、レイヤは36→18→9→18→36と設定しました。
活性化関数は、relu/sigmoidを利用しています。keras実装上、他の選択肢としては、[本家のサイト](https://keras.io/ja/activations/#_2)に記載されており、tanhや、LeakyReLUなども選択可能です。n_dimは、2項で分割した特徴次元数=40です。

```Python:
# AutoEncoderの構築

ae = Sequential()
ae.add(Dense(36, input_dim = n_dim, activation='relu'))
ae.add(Dense(18, activation='relu'))
ae.add(Dense(9, activation='relu', name = 'encoder'))
ae.add(Dense(18, activation='relu'))
ae.add(Dense(36, activation='relu'))
ae.add(Dense(n_dim, activation='sigmoid'))
```

モデルをcompileします。ロス関数は、MSE(Mean Squared Error)としています。keras実装上、他の選択肢としては多数あり、[本家のサイト](https://keras.io/ja/losses/)に記載されています。(mean_squared_errorや、mean_absolute_error、クラス分類時に利用するcategorical_crossentropyなど)
今回は、回帰問題のため、MSEを利用しています。

```Python:
ae.compile(loss = 'mse', optimizer ='adam')
```

最後にモデルをfitで訓練させます。訓練データのNumpy配列と、ターゲットのNumpy配列を指定します。AutoEncoderの場合は、訓練もターゲットも同一のX_trnとなります。バッチサイズやエポック数は、いくつか試して、収束が早く、十分な結果が得られたので、この値にしています。validation_data(評価データ)として、テストデータとテストのターゲットデータをタプルで設定します。今回は、(X_tst, X_tst)となります。

```Python:
records_ae = ae.fit(X_trn, X_trn,
                    epochs=100,
                    batch_size=128,
                    shuffle=True,
                    validation_data=(X_tst, X_tst))
```

ネットワークの概要表示は以下です。
想定通りに設定がされていることが分かります。

```Python:
# ネットワークの概要
ae.summary()
```

![AutoEncoderネットワーク概要](https://storage.googleapis.com/permanent_koizumi/20181229_autoencoder_summary2.png "AutoEncoderネットワーク概要")




## 4. 学習と評価

上記で学習をさせると、以下のようなロスヒストリーとなりました。
シンプルなデータ(ピー音（規制音）)だったため、かなりの精度になっていますが、今回は実装方法の勉強がメインなので特に気にせず進みます。
各種パラメータ（層の数や、活性化関数、ロス関数など）が異なる場合や、今回のようなシンプルな音データでない場合での評価などは後日実施予定です。

![ロスヒストリ](https://storage.googleapis.com/permanent_koizumi/20181229_loss_history.png "ロスヒストリ")

評価をする際は、正常(regular)と異常データ(anomaly)を、当該モデルに入力(predict)して、出力された結果データと入力データとの、MSEを計算することで、どれだけ異なっているかを判断します。学習で利用していない音声データを利用します。


```Python:
# 別の正常音データを入力してみて、MSEを見てみる
(rate_reg,sig_reg) = wav.read("20181229_test_sound(regularity).wav")
fbank_feat_reg = logfbank(sig_reg,rate_reg,winstep=0.01,nfilt=40)
X_reg = scaler.transform(fbank_feat_reg)

output_reg = ae.predict(X_reg)
mse_reg = mean_squared_error(X_reg, output_reg )
print("RegularityData mse:",mse_reg)

# 異常音データを入力してみて、MSEを見てみる
(rate_ano,sig_ano) = wav.read("20181229_test_sound(anomaly).wav")
fbank_feat_ano = logfbank(sig_ano,rate_ano,winstep=0.01,nfilt=40)
X_ano = scaler.transform(fbank_feat_ano)

output_ano = ae.predict(X_ano)
mse_ano = mean_squared_error(X_ano, output_ano)
print("AnomalyData mse:",mse_ano)
```

結果は、以下のようなものとなりました。
異常データを当該モデルに入力すると、正常データと比較してMSEの値が大きくなり、異常であると判定できる情報が得られました。

| 項目        | MSE                |
|:------------|:-------------------|
| 正常データ  | 0.1663619225413806  |
| 異常データ  | 31.446604284564717  |


以下に全コードを載せておきます。

```Python:全コード
from keras.utils import np_utils
from keras.models import Sequential, Model
from keras.layers import Activation, Dense, Dropout, Input
from keras.optimizers import Adam
import matplotlib.pyplot as plt
import os
import numpy as np
%matplotlib inline

#from python_speech_features import mfcc
from python_speech_features import logfbank
import scipy.io.wavfile as wav
from sklearn.preprocessing import MinMaxScaler
from sklearn.metrics import mean_squared_error

# =====================================================================
# 音声データをpython speech featuresで配列変換
(rate,sig) = wav.read("20181229_train_sound(regularity).wav")

print("sig:",sig,",len:",len(sig),",type:",type(sig))
print("rate:",rate)

fbank_feat = logfbank(sig,rate,winlen=0.01,nfilt=40)

# =====================================================================
# 学習データとテストデータは60000データくらいのうちの前半30000でトレーニング、後半30000でトレーニングとしておく。適当。
n_traintest_sep = 30000

X_trn = fbank_feat[:n_traintest_sep][:]
X_tst = fbank_feat[n_traintest_sep:][:]

# 値の範囲を[0,1]に変換
scaler = MinMaxScaler()
X_trn = scaler.fit_transform(X_trn)
X_tst = scaler.transform(X_tst)

# 入力データの次元数(=40)を取得
n_dim = X_trn.shape[1]


# =====================================================================
# 学習履歴をプロットする関数
   
# 損失関数値の履歴のプロット
def plot_history_loss(rec):
    plt.plot(rec.history['loss'],"o-",label="train",)
    plt.plot(rec.history['val_loss'],"o-",label="test")
    plt.title('loss history')
    plt.xlabel('epochs')
    plt.ylabel('loss')
    plt.legend(loc='upper right')
    plt.show()
    
    
    
# =====================================================================
# AutoEncoderの構築

ae = Sequential()
ae.add(Dense(36, input_dim = n_dim, activation='relu'))
ae.add(Dense(18, activation='relu'))
ae.add(Dense(9, activation='relu', name = 'encoder'))
ae.add(Dense(18, activation='relu'))
ae.add(Dense(36, activation='relu'))
ae.add(Dense(n_dim, activation='sigmoid'))

ae.compile(loss = 'mse', optimizer ='adam')
records_ae = ae.fit(X_trn, X_trn,
                    epochs=100,
                    batch_size=128,
                    shuffle=True,
                    validation_data=(X_tst, X_tst))

# 学習済み重みの保存
#ae.save_weights('autoencoder.h5')
# ネットワークの概要
ae.summary()
# 損失関数値の履歴のプロット
plot_history_loss(records_ae)

#testデータを入力してみて、MSEを見てみる
output_tst = ae.predict(X_tst)
mse_tst = mean_squared_error(X_tst, output_tst )
print("TestData mse:",mse_tst)

# 別の正常音データを入力してみて、MSEを見てみる
(rate_reg,sig_reg) = wav.read("20181229_test_sound(regularity).wav")
fbank_feat_reg = logfbank(sig_reg,rate_reg,winstep=0.01,nfilt=40)
X_reg = scaler.transform(fbank_feat_reg)

output_reg = ae.predict(X_reg)
mse_reg = mean_squared_error(X_reg, output_reg )
print("RegularityData mse:",mse_reg)

# 異常音データを入力してみて、MSEを見てみる
(rate_ano,sig_ano) = wav.read("20181229_test_sound(anomaly).wav")
fbank_feat_ano = logfbank(sig_ano,rate_ano,winstep=0.01,nfilt=40)
X_ano = scaler.transform(fbank_feat_ano)

output_ano = ae.predict(X_ano)
mse_ano = mean_squared_error(X_ano, output_ano)
print("AnomalyData mse:",mse_ano)
```


## おわりに

シンプルなDeep AutoEncoderで、異常検知を判定してみました。
今回は理想的なデータ(ピー音（規制音）)を作って異常検知してみましたが、より実際のデータなどを利用したうえでの異常判定ができるかや、AutoEncoderの中でも他の手法でも試行してみたいと思います。

## 参考

* [Kerasで学ぶAutoencoder - Elix Tech Blog](https://elix-tech.github.io/ja/2016/07/17/autoencoder.html)
* [DeepAutoEncoder - blog.keras.io](https://blog.keras.io/building-autoencoders-in-keras.html)
* [python-speech-features](https://python-speech-features.readthedocs.io/en/latest/)
* [Pythonで音声信号処理 - 人工知能に関する断創録](http://aidiary.hatenablog.com/entry/20110514/1305377659)