CNNを使ったヒグラシの音声認識（モデル構築編）～夏休みの自由研究～

# はじめに
そろそろ夏休みに入った小中学校も多いかと思います。夏休みといえば必ずと言って出されていた夏休みの自由研究。そんな夏休みの自由研究に役立つのではと思いセミの音声認識を紹介します。（半分冗談ですが、CNNを使って自由研究を仕上げてくるつよつよ小中学生もいてもおかしくない時代なのかなと思います。）

音声の分類が可能であれば、ヒグラシが鳴いている時間帯の記録を自動で収集でき、例えば気象データと合わせてヒグラシがよく鳴く気象条件を導き出せたりするのではないかなと思います。とりあえず本記事での目標は分類器を作るところまでとしました。

# 実行環境
google colaboratory を使います。googleアカウントがあればだれでも使用でき、環境構築不要、無料でGPU（計算が速い環境）を使うことができるといったメリットがあります(詳しくはググってみてください)。今回はGPU環境を用いました。

# 実行ディレクトリの構成
少しわかりにくいかと思い、図にしました。
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/864771/9ce38a5a-2afb-a5cb-ba96-c9e758a3c099.png)
今回は図のようにフォルダおよびファイルを配置し
・class0にヒグラシの音声（ほとんど自力で集めました笑）
・class1にヒグラシが鳴いていない音声
・class2に雨の日の音声（おまけ）
を入れ検証をしました。

もちろん、分類class数を増やしてもいけると思いますし、セミの音声でなくても同じように分類できると思います。

# コードの紹介
作成に当たり、<a href="https://qiita.com/cvusk/items/61cdbce80785eaf28349">こちらの記事</a>を参考に少々工夫を加えてみました。
データを集め、Google drive内のディレクトリの図のように構成を行い、sound_class.ipynb内に以下のコードをコピペして実行していけば動くと思います。

###①Colab特有の呪文
マウント、ディレクトリ移動を行うところです。

```mount.py
from google.colab import drive
drive.mount('/content/drive/')
``` 

```move.py
%cd /content/drive//MyDrive/Colab Notebooks/
``` 

###②データの下処理
wavファイルを読み込み、メルスペクトログラムにする。というようなことをやっています。<a href="https://qiita.com/cvusk/items/61cdbce80785eaf28349">こちらの記事</a>
（音声データを波形の画像データにするようなイメージです）

まず、使うモジュールをインポートします。

```module.py
import os
import glob
import numpy as np
import librosa
import librosa.display
import matplotlib.pyplot as plt
``` 

次に、関数の定義（ここはほぼコピペです。すみません汗）

``` define.py
# load a wave data
def load_wave_data(audio_dir, file_name):
    file_path = os.path.join(audio_dir, file_name)
    x, fs = librosa.load(file_path, sr=44100)
    return x,fs

# change wave data to mel-stft
def calculate_melsp(x, n_fft=1024, hop_length=128):
    stft = np.abs(librosa.stft(x, n_fft=n_fft, hop_length=hop_length))**2
    log_stft = librosa.power_to_db(stft)
    melsp = librosa.feature.melspectrogram(S=log_stft,n_mels=128)
    return melsp

# display wave in plots
def show_wave(x):
    plt.plot(x)
    plt.show()

# display wave in heatmap
def show_melsp(melsp, fs):
    librosa.display.specshow(melsp, sr=fs)
    plt.colorbar()
    plt.show()

# wavfile division
def wav_div_nparr(fname):
    x, fs = load_wave_data('', fname)
    xls = []
    for i in range(0,len(x)-fs,fs):
        xls.append(np.copy(x[i:i+fs]))
    return np.array(xls)
``` 

データの変換部分(フォルダに入っている音声ファイルのデータをCNNに使える形にデータを変換します。）

```data_to_data.py
folder = 'sounddata/'
files = glob.glob(folder+'class*/*.wav')
Xls = []
yls = []
for file in files:
    label = file[file.find('class'):file.find('.wav')].split('/')[0][5:] # classの後の数字をとってくる
    x = wav_div_nparr(file)
    for i in range(x.shape[0]):
        melsp = calculate_melsp(x[i])
        Xls.append(melsp)
        yls.append(label)
X = np.array(Xls)
X = X.reshape(X.shape[0],X.shape[1],X.shape[2],1)
Y = np.array(yls).astype(int)
print(X.shape,Y.shape)
```

データの分割（学習に使うデータと検証に使うデータを分けます）
また、ラベルデータをカテゴリカル変数化しています。

``` data_split.py
from sklearn.model_selection import train_test_split
from tensorflow.keras.utils import to_categorical
x_train, x_test, y_train, y_test = train_test_split(
    X,
    Y,
    random_state = 0,
    test_size = 0.2
)
# y to categorical
classes = np.max(Y)+1 #今回は3クラス分類
y_train = to_categorical(y_train, classes)
y_test = to_categorical(y_test, classes)
```

メルスペクトログラムの可視化（本筋ではない確認のための部分）
例えばこんな感じです。

```melsp_display.py
fs = 44100
X = wav_div_nparr('sounddata/class2/rain.wav') # ファイル名は存在するものを指定
melsp = calculate_melsp(X[10])
show_wave(X[10])
show_melsp(melsp, fs)
```
オリジナルの波形
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/864771/4978c4f9-e5c4-22c6-6b67-bc87a8986668.png)
メルスペクトログラム
![image.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/864771/3086fd36-50d5-e2ba-7514-ce120e3dcde2.png)


###③モデルの定義
ここもほぼ参考記事のコピペになってしまいますが、工夫としてMaxPooling層を間に挟み、広い範囲の特徴量抽出をより軽い計算で回るようにしました。また、dropout層を挟むことで過学習に強いモデルとしました。

```model_define.py
from keras.optimizers import Adam
from keras.models import Model,Input
from keras.layers import Dense, Dropout, Activation
from keras.layers import Conv2D, GlobalAveragePooling2D,MaxPooling2D,Flatten
from keras.layers import BatchNormalization, Add

def cba(inputs, filters, kernel_size, strides):
    x = Conv2D(filters, kernel_size=kernel_size, strides=strides, padding='same')(inputs)
    x = BatchNormalization()(x)
    x = Activation("relu")(x)
    return x

# define CNN
inputs = Input(shape=(x_train.shape[1:]))

x_1 = cba(inputs, filters=32, kernel_size=(1,8), strides=(1,2))
x_1 = cba(x_1, filters=32, kernel_size=(8,1), strides=(2,1))
x_1 = MaxPooling2D(pool_size = (2, 2))(x_1)
x_1 = cba(x_1, filters=64, kernel_size=(1,8), strides=(1,2))
x_1 = cba(x_1, filters=64, kernel_size=(8,1), strides=(2,1))

x_2 = cba(inputs, filters=32, kernel_size=(1,16), strides=(1,2))
x_2 = cba(x_2, filters=32, kernel_size=(16,1), strides=(2,1))
x_2 = MaxPooling2D(pool_size = (2, 2))(x_2)
x_2 = cba(x_2, filters=64, kernel_size=(1,16), strides=(1,2))
x_2 = cba(x_2, filters=64, kernel_size=(16,1), strides=(2,1))

x_3 = cba(inputs, filters=32, kernel_size=(1,32), strides=(1,2))
x_3 = cba(x_3, filters=32, kernel_size=(32,1), strides=(2,1))
x_3 = MaxPooling2D(pool_size = (2, 2))(x_3)
x_3 = cba(x_3, filters=64, kernel_size=(1,32), strides=(1,2))
x_3 = cba(x_3, filters=64, kernel_size=(32,1), strides=(2,1))

x = Add()([x_1, x_2, x_3])
x = MaxPooling2D(pool_size = (2, 2))(x)
x = Dropout(0.25)(x)
x = cba(x, filters=64, kernel_size=(1,8), strides=(1,2))
x = cba(x, filters=64, kernel_size=(8,1), strides=(2,1))
x = GlobalAveragePooling2D()(x)

x = Dense(64)(x)
x = Activation("relu")(x)
x = Dense(classes)(x)
x = Activation("softmax")(x)

model = Model(inputs, x)

# initiate Adam optimizer
opt = Adam(learning_rate=0.001, decay=1e-6, amsgrad=True)

# Let's train the model using Adam with amsgrad
model.compile(loss='categorical_crossentropy',
              optimizer=opt,
              metrics=['accuracy'])
```

###④モデルの学習
一行です笑(epochsが全学習データを学習に使う回数で多いほど計算に時間がかかります)

```fit.py
model.fit(x_train,y_train,epochs=10,batch_size=30,validation_data=(x_test, y_test))
```
実行結果

```
Epoch 1/10
46/46 [==============================] - 12s 185ms/step - loss: 0.5362 - accuracy: 0.7706 - val_loss: 1.5023 - val_accuracy: 0.5565
Epoch 2/10
46/46 [==============================] - 7s 162ms/step - loss: 0.0579 - accuracy: 0.9784 - val_loss: 2.3608 - val_accuracy: 0.5159
Epoch 3/10
46/46 [==============================] - 8s 165ms/step - loss: 0.0367 - accuracy: 0.9854 - val_loss: 0.6912 - val_accuracy: 0.7304
Epoch 4/10
46/46 [==============================] - 8s 168ms/step - loss: 0.0406 - accuracy: 0.9869 - val_loss: 0.0366 - val_accuracy: 0.9884
Epoch 5/10
46/46 [==============================] - 8s 171ms/step - loss: 0.0181 - accuracy: 0.9935 - val_loss: 0.0274 - val_accuracy: 0.9942
Epoch 6/10
46/46 [==============================] - 8s 172ms/step - loss: 0.0109 - accuracy: 0.9994 - val_loss: 0.0203 - val_accuracy: 0.9971
Epoch 7/10
46/46 [==============================] - 8s 170ms/step - loss: 0.0047 - accuracy: 0.9996 - val_loss: 0.0138 - val_accuracy: 0.9942
Epoch 8/10
46/46 [==============================] - 8s 167ms/step - loss: 0.0050 - accuracy: 0.9979 - val_loss: 0.0013 - val_accuracy: 1.0000
Epoch 9/10
46/46 [==============================] - 8s 165ms/step - loss: 0.0025 - accuracy: 0.9998 - val_loss: 0.0019 - val_accuracy: 1.0000
Epoch 10/10
46/46 [==============================] - 8s 164ms/step - loss: 0.0162 - accuracy: 0.9944 - val_loss: 0.0426 - val_accuracy: 1.0000
<keras.callbacks.History at 0x7f59980f2690>
```

###⑤モデルの評価
学習データ(train)と検証データ(test)についてマトリックスと評価指標による評価を行います。

評価用モジュールのインポート

```metrics_import.py
from sklearn.metrics import confusion_matrix,classification_report
```
学習データに対しての評価

```train_result.py
pred_train = model.predict(x_train)
print(confusion_matrix(np.argmax(y_train,axis=1),np.argmax(pred_train,axis=1)))
print(classification_report(np.argmax(y_train,axis=1),np.argmax(pred_train,axis=1)))
```
実行結果

```
[[479   0   0]
 [  2 572   0]
 [  0   2 325]]
              precision    recall  f1-score   support

           0       1.00      1.00      1.00       479
           1       1.00      1.00      1.00       574
           2       1.00      0.99      1.00       327

    accuracy                           1.00      1380
   macro avg       1.00      1.00      1.00      1380
weighted avg       1.00      1.00      1.00      1380
```



検証用データに対しての評価

```test_result.py
pred_test = model.predict(x_test)
print(confusion_matrix(np.argmax(y_test,axis=1),np.argmax(pred_test,axis=-1)))
print(classification_report(np.argmax(y_test,axis=1),np.argmax(pred_test,axis=-1)))
```
実行結果

```
[[119   0   0]
 [  0 147   0]
 [  0   0  79]]
              precision    recall  f1-score   support

           0       1.00      1.00      1.00       119
           1       1.00      1.00      1.00       147
           2       1.00      1.00      1.00        79

    accuracy                           1.00       345
   macro avg       1.00      1.00      1.00       345
weighted avg       1.00      1.00      1.00       345
```
ほとんど正解してますね。（汎化性能には不安が残りますが、まずまずでしょう。）

最後にモデルの保存（今後別の日に取ったデータで検証することもあると思うので保存しておきます。）

```modelsave.py
model.save('modelA')
```

#おわりに
機械学習についてよく学んでいる方々からすると、この方法でのモデル評価では不十分だと思う方もいると思います。実際その通りだと思ってまして、今回のテーマに合わせては例えば別日でのデータを別の検証用データとして用意する必要があると考えています。この辺りは今後取り組んでいきたいと考えています。進捗あればまた記事にしようと思いますので、よろしくお願いいたします。