【個人的備忘録】Jetson Nanoでポテチの音声分類に挑戦した話

### 動機

　とある学生プロジェクトで，ポテチの咀嚼音から人がポテチを食べているかをリアルタイムで検出するというタスクが発生したため，Kerasを用いた音声分類に挑戦してみた．
　@cvusk さんの記事「ディープラーニングで音声分類」（https://qiita.com/cvusk/items/61cdbce80785eaf28349）を大いに参考にさせていただきました．

###　学習

####　教師データの準備

・ポテチを食べるときの特徴的な「パリッ」という音に着目した．0.5秒もあればこの音を検出できると判断し，時間窓の間隔は0.5秒とした．

・また，2種類の自然音（環境音）を用意した．**時間窓は同様に0.5秒**とした．
**※0.5秒でなくてもよいが，すべてのデータの時間はそろえること．**

・それぞれに対して40サンプル（合計120サンプル）用意した．
**※40サンプルじゃなくてもいい．任意．**

・音声の分割には**Audacity**を用いた．形式は**WAV**．サンプリングレートは44100[Hz]，量子化bit数は32bit

・ポテチは0，自然音1は1，自然音2は2としてラベル付けを行った．ラベル付け用に書いたコードは以下．

```python
#labeling.py

import csv
import random
import pandas as pd

l = []
for i in range(40):
    l.append(('j' + str(i + 1) + '.wav', 0, 'potechi'))
    l.append(('n' + str(i + 1) + '.wav', 1, 'natural'))
    l.append(('o' + str(i + 1) + '.wav', 2, 'natural2'))

random.shuffle(l)
print(l)

df = pd.DataFrame(l, columns = ['filename', 'target', 'category'])
print(df)

df.to_csv('eatSound.csv')
```

ラベリング結果

eatSound.csv

| No   | filename | target | category |
| ---- | -------- | ------ | -------- |
| 0    | o17.wav  | 2      | natural2 |
| 1    | o11.wav  | 2      | natural2 |
| 2    | o31.wav  | 2      | natural2 |
| 3    | j36.wav  | 0      | potechi  |
| 4    | j39.wav  | 0      | potechi  |
| 5    | o4.wav   | 2      | natural2 |
| 6    | n3.wav   | 1      | natural  |
| 7    | j24.wav  | 0      | potechi  |
| 8    | j27.wav  | 0      | potechi  |
| 9    | j17.wav  | 0      | potechi  |
| 10   | n27.wav  | 1      | natural  |
| 11   | o1.wav   | 2      | natural2 |
| 12   | o2.wav   | 2      | natural2 |
| 13   | n9.wav   | 1      | natural  |
| 14   | o38.wav  | 2      | natural2 |
| 15   | o37.wav  | 2      | natural2 |
| 16   | o27.wav  | 2      | natural2 |
| 17   | o8.wav   | 2      | natural2 |
| 18   | n17.wav  | 1      | natural  |
| 19   | j2.wav   | 0      | potechi  |
| 20   | j38.wav  | 0      | potechi  |
| 21   | o19.wav  | 2      | natural2 |
| 22   | n40.wav  | 1      | natural  |
| 23   | j4.wav   | 0      | potechi  |
| 24   | o3.wav   | 2      | natural2 |
| 25   | o26.wav  | 2      | natural2 |
| 26   | o18.wav  | 2      | natural2 |
| 27   | n35.wav  | 1      | natural  |
| 28   | n30.wav  | 1      | natural  |
| 29   | j15.wav  | 0      | potechi  |
| 30   | o36.wav  | 2      | natural2 |
| 31   | j40.wav  | 0      | potechi  |
| 32   | o21.wav  | 2      | natural2 |
| 33   | j29.wav  | 0      | potechi  |
| 34   | n2.wav   | 1      | natural  |
| 35   | j6.wav   | 0      | potechi  |
| 36   | o6.wav   | 2      | natural2 |
| 37   | j20.wav  | 0      | potechi  |
| 38   | n39.wav  | 1      | natural  |
| 39   | j31.wav  | 0      | potechi  |
| 40   | n13.wav  | 1      | natural  |
| 41   | j21.wav  | 0      | potechi  |
| 42   | n7.wav   | 1      | natural  |
| 43   | n21.wav  | 1      | natural  |
| 44   | n24.wav  | 1      | natural  |
| 45   | n10.wav  | 1      | natural  |
| 46   | n15.wav  | 1      | natural  |
| 47   | n22.wav  | 1      | natural  |
| 48   | j14.wav  | 0      | potechi  |
| 49   | j32.wav  | 0      | potechi  |
| 50   | j22.wav  | 0      | potechi  |
| 51   | o15.wav  | 2      | natural2 |
| 52   | j12.wav  | 0      | potechi  |
| 53   | j25.wav  | 0      | potechi  |
| 54   | n34.wav  | 1      | natural  |
| 55   | o13.wav  | 2      | natural2 |
| 56   | n1.wav   | 1      | natural  |
| 57   | o22.wav  | 2      | natural2 |
| 58   | n26.wav  | 1      | natural  |
| 59   | o5.wav   | 2      | natural2 |
| 60   | j9.wav   | 0      | potechi  |
| 61   | j35.wav  | 0      | potechi  |
| 62   | j5.wav   | 0      | potechi  |
| 63   | n19.wav  | 1      | natural  |
| 64   | n36.wav  | 1      | natural  |
| 65   | n32.wav  | 1      | natural  |
| 66   | o35.wav  | 2      | natural2 |
| 67   | n6.wav   | 1      | natural  |
| 68   | j18.wav  | 0      | potechi  |
| 69   | j37.wav  | 0      | potechi  |
| 70   | j34.wav  | 0      | potechi  |
| 71   | o39.wav  | 2      | natural2 |
| 72   | n33.wav  | 1      | natural  |
| 73   | j10.wav  | 0      | potechi  |
| 74   | n29.wav  | 1      | natural  |
| 75   | j13.wav  | 0      | potechi  |
| 76   | n37.wav  | 1      | natural  |
| 77   | n11.wav  | 1      | natural  |
| 78   | n20.wav  | 1      | natural  |
| 79   | o9.wav   | 2      | natural2 |
| 80   | o40.wav  | 2      | natural2 |
| 81   | j7.wav   | 0      | potechi  |
| 82   | n23.wav  | 1      | natural  |
| 83   | j28.wav  | 0      | potechi  |
| 84   | j19.wav  | 0      | potechi  |
| 85   | n31.wav  | 1      | natural  |
| 86   | o7.wav   | 2      | natural2 |
| 87   | o20.wav  | 2      | natural2 |
| 88   | o23.wav  | 2      | natural2 |
| 89   | n38.wav  | 1      | natural  |
| 90   | j3.wav   | 0      | potechi  |
| 91   | j23.wav  | 0      | potechi  |
| 92   | n8.wav   | 1      | natural  |
| 93   | n25.wav  | 1      | natural  |
| 94   | o28.wav  | 2      | natural2 |
| 95   | n28.wav  | 1      | natural  |
| 96   | n14.wav  | 1      | natural  |
| 97   | o16.wav  | 2      | natural2 |
| 98   | j33.wav  | 0      | potechi  |
| 99   | n16.wav  | 1      | natural  |
| 100  | o25.wav  | 2      | natural2 |
| 101  | o32.wav  | 2      | natural2 |
| 102  | o33.wav  | 2      | natural2 |
| 103  | j8.wav   | 0      | potechi  |
| 104  | o10.wav  | 2      | natural2 |
| 105  | n18.wav  | 1      | natural  |
| 106  | o34.wav  | 2      | natural2 |
| 107  | o14.wav  | 2      | natural2 |
| 108  | o30.wav  | 2      | natural2 |
| 109  | j26.wav  | 0      | potechi  |
| 110  | n4.wav   | 1      | natural  |
| 111  | j16.wav  | 0      | potechi  |
| 112  | o29.wav  | 2      | natural2 |
| 113  | j30.wav  | 0      | potechi  |
| 114  | o12.wav  | 2      | natural2 |
| 115  | n5.wav   | 1      | natural  |
| 116  | o24.wav  | 2      | natural2 |
| 117  | j11.wav  | 0      | potechi  |
| 118  | j1.wav   | 0      | potechi  |
| 119  | n12.wav  | 1      | natural  |

### 教師データの学習

　原理とかモデルとかは「ディープラーニングで音声分類」https://qiita.com/cvusk/items/61cdbce80785eaf28349 を参照．

#### 環境

・OS：Ubuntu 18.04. LTS
・GPU：Geforce GTX 1080
・Python 3.7.3
	--tensorflow-gpu	1.13.1+nv19.4
	--keras	2.3.1
	--numpy	1.18.1
	--pandas	1.0.0
	--librosa	0.7.2
	--matplotlib	3.1.3
	--seaborn 0.10.0
	--scikit-learn	0.22.1
	--IPython	7.13.0



#### 学習コード

・変数timeは最初に表示されるmelsp sizeの右の値にする．
今回はmelsp size：（128, 173）なのでtime = 173

・25%をテストデータとするので，
	train_num = 90
	test_num = 30　とする．

・batch sizeは2，epoch数は100とした．

```python
#learning.py

import random
import numpy as np
import pandas as pd
import librosa
import librosa.display
import matplotlib.pyplot as plt
import seaborn as sn
from sklearn import model_selection
from sklearn import preprocessing
import IPython.display as ipd
import os

# define directories
base_dir = "./"
meta_file = os.path.join(base_dir, "audio/eatSound.csv")
audio_dir = os.path.join(base_dir, "audio/audiofile/")

# load metadata
meta_data = pd.read_csv(meta_file)
meta_data

# get data size
data_size = meta_data.shape
print(data_size)

# arrange target label and its name
class_dict = {}
for i in range(data_size[0]):
    if meta_data.loc[i,"target"] not in class_dict.keys():
        class_dict[meta_data.loc[i,"target"]] = meta_data.loc[i,"category"]
class_pd = pd.DataFrame(list(class_dict.items()), columns=["labels","classes"])
class_pd

# 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()

# example data
x, fs = load_wave_data(audio_dir, meta_data.loc[0,"filename"])
melsp = calculate_melsp(x)
print("wave size:{0}\nmelsp size:{1}\nsamping rate:{2}".format(x.shape, melsp.shape, fs))
show_wave(x)
show_melsp(melsp, fs)

ipd.Audio(x, rate=fs)

# data augmentation: add white noise
def add_white_noise(x, rate=0.002):
    return x + rate*np.random.randn(len(x))

x_wn = add_white_noise(x)
melsp = calculate_melsp(x_wn)
print("wave size:{0}\nmelsp size:{1}\nsamping rate:{2}".format(x_wn.shape, melsp.shape, fs))
show_wave(x_wn)
show_melsp(melsp, fs)

ipd.Audio(x_wn, rate=fs)

# data augmentation: shift sound in timeframe
def shift_sound(x, rate=2):
    return np.roll(x, int(len(x)//rate))

x_ss = shift_sound(x)
melsp = calculate_melsp(x_ss)
print("wave size:{0}\nmelsp size:{1}\nsamping rate:{2}".format(x_ss.shape, melsp.shape, fs))
show_wave(x_ss)
show_melsp(melsp, fs)


ipd.Audio(x_ss, rate=fs)

# data augmentation: stretch sound
def stretch_sound(x, rate=1.1):
    input_length = len(x)
    x = librosa.effects.time_stretch(x, rate)
    if len(x)>input_length:
        return x[:input_length]
    else:
        return np.pad(x, (0, max(0, input_length - len(x))), "constant")

x_st = stretch_sound(x)
melsp = calculate_melsp(x_st)
print("wave size:{0}\nmelsp size:{1}\nsamping rate:{2}".format(x_st.shape, melsp.shape, fs))
show_wave(x_st)
show_melsp(melsp, fs)

ipd.Audio(x_st, rate=fs)

# get training dataset and target dataset
x = list(meta_data.loc[:,"filename"])
y = list(meta_data.loc[:, "target"])

x_train, x_test, y_train, y_test = model_selection.train_test_split(x, y, test_size=0.25, stratify=y)
print("x train:{0}\ny train:{1}\nx test:{2}\ny test:{3}".format(len(x_train),
                                                                len(y_train), 
                                                                len(x_test), 
                                                                len(y_test)))


a = np.zeros(50)
for c in y_test:
    a[c] += 1
print(a)

freq = 128
time = 173 #1723
# save wave data in npz, with augmentation
def save_np_data(filename, x, y, aug=None, rates=None):
    np_data = np.zeros(freq*time*len(x)).reshape(len(x), freq, time)
    np_targets = np.zeros(len(y))
    for i in range(len(y)):
        _x, fs = load_wave_data(audio_dir, x[i])
        if aug is not None:
            _x = aug(x=_x, rate=rates[i])
        _x = calculate_melsp(_x)
        np_data[i] = _x
        np_targets[i] = y[i]
    np.savez(filename, x=np_data, y=np_targets)

# save test dataset
if not os.path.exists("esc_melsp_test.npz"):
    save_np_data("esc_melsp_test.npz", x_test,  y_test)

# save raw training dataset
if not os.path.exists("esc_melsp_train_raw.npz"):
    save_np_data("esc_melsp_train_raw.npz", x_train,  y_train)

# save training dataset with white noise
if not os.path.exists("esc_melsp_train_wn.npz"):
    rates = np.random.randint(1,50,len(x_train))/10000
    save_np_data("esc_melsp_train_wn.npz", x_train,  y_train, aug=add_white_noise, rates=rates)

# save training dataset with sound shift
if not os.path.exists("esc_melsp_train_ss.npz"):
    rates = np.random.choice(np.arange(2,6),len(y_train))
    save_np_data("esc_melsp_train_ss.npz", x_train,  y_train, aug=shift_sound, rates=rates)

# save training dataset with stretch
if not os.path.exists("esc_melsp_train_st.npz"):
    rates = np.random.choice(np.arange(80,120),len(y_train))/100
    save_np_data("esc_melsp_train_st.npz", x_train,  y_train, aug=stretch_sound, rates=rates)

# save training dataset with combination of white noise and shift or stretch
if not os.path.exists("esc_melsp_train_com.npz"):
    np_data = np.zeros(freq*time*len(x_train)).reshape(len(x_train), freq, time)
    np_targets = np.zeros(len(y_train))
    for i in range(len(y_train)):
        x, fs = load_wave_data(audio_dir, x_train[i])
        x = add_white_noise(x=x, rate=np.random.randint(1,50)/1000)
        if np.random.choice((True,False)):
            x = shift_sound(x=x, rate=np.random.choice(np.arange(2,6)))
        else:
            x = stretch_sound(x=x, rate=np.random.choice(np.arange(80,120))/100)
        x = calculate_melsp(x)
        np_data[i] = x
        np_targets[i] = y_train[i]
    np.savez("esc_melsp_train_com.npz", x=np_data, y=np_targets)

import keras
from keras.models import Model
from keras.layers import Input, Dense, Dropout, Activation
from keras.layers import Conv2D, GlobalAveragePooling2D
from keras.layers import BatchNormalization, Add
from keras.callbacks import EarlyStopping, ModelCheckpoint
from keras.models import load_model

# dataset files
train_files = ["esc_melsp_train_raw.npz", 
               "esc_melsp_train_ss.npz",
               "esc_melsp_train_st.npz", 
               "esc_melsp_train_wn.npz",
               "esc_melsp_train_com.npz"]
test_file = "esc_melsp_test.npz"
train_num = 90#1500
test_num = 30#500

# define dataset placeholders
x_train = np.zeros(freq*time*train_num*len(train_files)).reshape(train_num*len(train_files), freq, time)
y_train = np.zeros(train_num*len(train_files))

# load dataset
for i in range(len(train_files)):
    data = np.load(train_files[i])
    x_train[i*train_num:(i+1)*train_num] = data["x"]
    y_train[i*train_num:(i+1)*train_num] = data["y"]
    print('clear')

# load test dataset
test_data = np.load(test_file)
x_test = test_data["x"]
y_test = test_data["y"]

# redefine target data into one hot vector
classes = 3
y_train = keras.utils.to_categorical(y_train, classes)
y_test = keras.utils.to_categorical(y_test, classes)

# reshape training dataset
x_train = x_train.reshape(train_num*5, freq, time, 1)
x_test = x_test.reshape(test_num, freq, time, 1)

#classes = 50
#y_test = keras.utils.to_categorical(y_test, classes)
#x_test = x_test.reshape(test_num, freq, time, 1)

print("x train:{0}\ny train:{1}\nx test:{2}\ny test:{3}".format(x_train.shape,
                                                                y_train.shape, 
                                                                x_test.shape, 
                                                                y_test.shape))

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 = 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 = 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 = 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_4 = cba(inputs, filters=32, kernel_size=(1,64), strides=(1,2))
x_4 = cba(x_4, filters=32, kernel_size=(64,1), strides=(2,1))
x_4 = cba(x_4, filters=64, kernel_size=(1,64), strides=(1,2))
x_4 = cba(x_4, filters=64, kernel_size=(64,1), strides=(2,1))

x = Add()([x_1, x_2, x_3, x_4])

x = cba(x, filters=128, kernel_size=(1,16), strides=(1,2))
x = cba(x, filters=128, kernel_size=(16,1), strides=(2,1))

x = GlobalAveragePooling2D()(x)
x = Dense(classes)(x)
x = Activation("softmax")(x)

model = Model(inputs, x)
model.summary()

# initiate Adam optimizer
opt = keras.optimizers.Adam(lr=0.00001, decay=1e-6, amsgrad=True)

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

# directory for model checkpoints
model_dir = "./models"
if not os.path.exists(model_dir):
    os.mkdir(model_dir)

# early stopping and model checkpoint# early  
es_cb = EarlyStopping(monitor='val_loss', patience=10, verbose=1, mode='auto')
#chkpt = os.path.join(model_dir, 'model_.{epoch:02d}_{val_loss:.4f}_{val_acc:.4f}.hdf5')
chkpt = os.path.join(model_dir, 'model_.{epoch:02d}_{val_loss:.4f}.hdf5')
cp_cb = ModelCheckpoint(filepath = chkpt, monitor='val_loss', verbose=1, save_best_only=True, mode='auto')


# between class data generator
class MixupGenerator():
    def __init__(self, x_train, y_train, batch_size=16, alpha=0.2, shuffle=True):
        self.x_train = x_train
        self.y_train = y_train
        self.batch_size = batch_size
        self.alpha = alpha
        self.shuffle = shuffle
        self.sample_num = len(x_train)

    def __call__(self):
        while True:
            indexes = self.__get_exploration_order()
            itr_num = int(len(indexes) // (self.batch_size * 2))

            for i in range(itr_num):
                batch_ids = indexes[i * self.batch_size * 2:(i + 1) * self.batch_size * 2]
                x, y = self.__data_generation(batch_ids)

                yield x, y

    def __get_exploration_order(self):
        indexes = np.arange(self.sample_num)

        if self.shuffle:
            np.random.shuffle(indexes)

        return indexes

    def __data_generation(self, batch_ids):
        _, h, w, c = self.x_train.shape
        _, class_num = self.y_train.shape
        x1 = self.x_train[batch_ids[:self.batch_size]]
        x2 = self.x_train[batch_ids[self.batch_size:]]
        y1 = self.y_train[batch_ids[:self.batch_size]]
        y2 = self.y_train[batch_ids[self.batch_size:]]
        l = np.random.beta(self.alpha, self.alpha, self.batch_size)
        x_l = l.reshape(self.batch_size, 1, 1, 1)
        y_l = l.reshape(self.batch_size, 1)

        x = x1 * x_l + x2 * (1 - x_l)
        y = y1 * y_l + y2 * (1 - y_l)

        return x, y

# train model
batch_size = 2
epochs = 100

training_generator = MixupGenerator(x_train, y_train)()
model.fit_generator(generator=training_generator,
                    steps_per_epoch=x_train.shape[0] // batch_size,
                    validation_data=(x_test, y_test),
                    epochs=epochs, 
                    verbose=1,
                    shuffle=True,
                    callbacks=[es_cb, cp_cb])

print('Finish learning.')
```



###　リアルタイム推論

　result配列に推論結果の確率が格納される．

```python
import time
import numpy as np
from keras.models import load_model
import librosa
import pyaudio
import matplotlib as plt

# メルスペクトログラム変換処理
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

# モデルの読み込み
print('start loading model...')
model = load_model('models/model.hdf5')#モデルのパス
print(model.summary())
print('finish loading.')
#　オーディオ設定

t = 0.5           # 計測時間[s]
samplerate = 44100  # サンプリングレート
fs = int(samplerate * t)          # フレームサイズ

#録音する関数
def record(samplerate, fs, t):
    pa = pyaudio.PyAudio()
    data = []
    dt = 1 / samplerate
    
    # ストリームの開始
    stream = pa.open(format=pyaudio.paInt16, channels=1, rate=samplerate,
                     input=True, frames_per_buffer=fs)
    
    # フレームサイズ毎に音声を録音していくループ
    for i in range(int(((t / dt) / fs))):
        frame = stream.read(fs)
        data.append(frame)
 
    # ストリームの終了
    stream.stop_stream()
    stream.close()
    pa.terminate()
    
    # データをまとめる処理
    data = b"".join(data)
    
    # データをNumpy配列に変換
    data = np.frombuffer(data, dtype='int16') / float(2**15)

    return data, i

while True:
    #推論
    x, i = record(samplerate, fs, time)
    xm = calculate_melsp(x)
    xm = xm.reshape(1, 128, 173, 1)
    result = model.predict(xm, batch_size=2)
    print('ポテチ：', result[0], '自然音１：', result[1],'自然音2：', result[2],)
```



##　最後に

・ラベル付けもうちょっと楽にしたい．

・説明が雑なので今後も記事をアップデートする．