Voice to Health

概要:音声から健康情報の抽出

モチベーション

音声認識ライブラリで音から文字に変換した場合、音声に含まれている情報の内、文章情報だけしか利用できない、音声に含まれているが今までノイズ(ゆらぎ)として捨てられていた情報の活用方法を探りたい。

調査内容

今回は健康情報の抽出をテーマとする。実際に風を引いた時の音声を使うことがベストだが、すぐに入手できないため、演技で元気の良い声と疲れている声を出し、それを録音した。録音した音声に対して特徴量の抽出と判別を行った。

結論

以下の方法で元気そうな声と疲れている声の判別に成功した。

  • 特徴量:
    • 振幅の平均値
    • ピッチのピーク
  • 判別方法:
    • SVM

set graph inline on ipython notebook

In [1]:
%matplotlib inline

import librarys

In [2]:
import glob
import numpy as np
from scipy.io import wavfile
import matplotlib.pyplot as plt
from sklearn.svm import SVC

Webカメラ付属のマイクで自分の声を録音し、wav形式で保存した。1ファイル1単語であり。単語はすべて『おはよう』で統一した。strongとファイル名に付いているファイルは元気よく発話した音声であり、weakと付いているものは元気がなさそうに発話した音声である。

In [3]:
filelist = glob.glob('./ohayou_*.wav')
filelist
Out[3]:
['.\\ohayou_strong1.wav',
 '.\\ohayou_strong2.wav',
 '.\\ohayou_strong3.wav',
 '.\\ohayou_strong4.wav',
 '.\\ohayou_strong5.wav',
 '.\\ohayou_strong6.wav',
 '.\\ohayou_weak1.wav',
 '.\\ohayou_weak2.wav',
 '.\\ohayou_weak3.wav',
 '.\\ohayou_weak4.wav',
 '.\\ohayou_weak5.wav',
 '.\\ohayou_weak6.wav']

音声ファイルを読み込み。連想配列の配列に保存。ステレオ録音のため一列目(おそらく右)の波形だけ使用。

In [4]:
snd_array = []
for fname in filelist:
    fs, snd = wavfile.read(fname)
    snd_array.append({'data':snd[:, 1], 'sample_rate':fs, 'fname':fname})

一つ目のファイルをプロット。

In [6]:
plot_data = snd_array[0]
plt.plot(np.arange(0.0, np.size( plot_data['data']), 1)/plot_data['sample_rate'], plot_data['data'])
plt.xlim([0,1])
plt.ylim([-30000,30000])
plt.title(plot_data['fname'])
plt.xlabel('[second]')
Out[6]:
<matplotlib.text.Text at 0x98ec330>

元気がなさそうに発音した時の音声をプロット。後ろの小さな山はため息(元気がなさそうに話すと無意識の内に出てしまった)。

In [7]:
plot_data = snd_array[7]
plt.plot(np.arange(0.0, np.size( plot_data['data']), 1)/plot_data['sample_rate'], plot_data['data'])
plt.xlim([0,1])
plt.ylim([-30000,30000])
plt.title(plot_data['fname'])
plt.xlabel('[ms]')
Out[7]:
<matplotlib.text.Text at 0x9912ad0>

今後何度も使うので。時間のインデックスを返す関数を定義しておく。

In [8]:
def time_index(data):
    import numpy as np
    return np.arange(0.0, np.size(data['data']), 1)/data['sample_rate']
In [9]:
time_index(snd_array[0])
Out[9]:
array([  0.00000000e+00,   2.26757370e-05,   4.53514739e-05, ...,
         5.49863946e-01,   5.49886621e-01,   5.49909297e-01])

すべての音声をプロットする。赤色は元気な声、青色は元気のない声(以降、特に断りのない限り同じ色分けを使う)。

In [10]:
for snd_data in snd_array:
    if 'weak' in snd_data['fname']:
        col = 'b'
    elif 'strong' in snd_data['fname']:
        col = 'r'
    plt.figure()
    plt.plot(time_index(snd_data), snd_data['data'], color=col)
    plt.xlim([0,1])
    plt.ylim([-2.**15, 2.**15])
    plt.title(plot_data['fname'])
    plt.xlabel('[second]')

振幅情報を抽出し、プロットする。振幅は元の音声データの絶対値を取り、500点平均したものとする。先頭に0を詰め数を合わせる。

In [11]:
for snd_data in snd_array:
    n_window = 500
    org_data = np.r_[np.zeros(n_window), snd_data['data']]
    scan_range = range(0, snd_data['data'].shape[0])
    amplitude_mean = np.zeros(len(scan_range))
    for i_time in scan_range:
        extract_data = org_data[i_time:i_time + n_window]
        amplitude_mean[i_time] = np.mean(np.abs(extract_data))
    snd_data['amp'] = amplitude_mean
    if 'weak' in snd_data['fname']:
        col = 'b'
    elif 'strong' in snd_data['fname']:
        col = 'r'
    plt.figure()
    plt.plot(time_index(snd_data), snd_data['amp'], color=col)

振幅が0.01以上かつ、0.1秒以上連続している区間を発話区間とする。

In [12]:
def extract_speak(snd, th_amp=1000, th_time=0.1):
    import numpy as np
    index_snd_buff = th_amp < snd['amp']
    index_snd = np.bool_(np.zeros(index_snd_buff.shape[0]))
    count = 0
    for i, v in enumerate(index_snd_buff):
        if v == False:
            if th_time * snd['sample_rate'] < count:
                index_snd[i-count:i] = True
            count = 0
        else:
            count += 1
    return index_snd

発話区間だけ抽出してプロット

In [13]:
for snd_data in snd_array:
    plot_data = snd_data['amp'][extract_speak(snd_data)]
    if 'weak' in snd_data['fname']:
        col = 'b'
    elif 'strong' in snd_data['fname']:
        col = 'r'
    plt.figure()
    plt.plot(time_index(snd_data)[:plot_data.shape[0]], plot_data, color=col)

自己相互相関を求める関数

In [14]:
def autocorr(dt):
    import scipy.signal as sig
    cor = sig.correlate(dt,dt,mode="full")
    return cor[cor.size/2:]

sin関数で挙動をテスト

In [15]:
sin_wave = np.sin(np.linspace(0,20,100))
plt.plot(sin_wave)
Out[15]:
[<matplotlib.lines.Line2D at 0x14e5b950>]
In [16]:
plt.plot(autocorr(sin_wave))
Out[16]:
[<matplotlib.lines.Line2D at 0x110c8d90>]

音声データの周波数を確認する(実行に30秒ほどかかる)

In [17]:
snd_data = snd_array[0]
snd_index = extract_speak(snd_data)
pitch = autocorr(np.float_(snd_data['data'][snd_index]))
In [18]:
plt.plot(np.arange(0.0, sum(snd_index), 1)/snd_data['sample_rate'], pitch)
Out[18]:
[<matplotlib.lines.Line2D at 0x14d5db70>]

0以外で最初のピークの位置を探したいがノイズが有るため抽出が難しい。(以下のグラフでは0.005[ms]をピークとして見つけたい)

In [19]:
plt.plot(np.arange(0.0, sum(snd_index), 1)/snd_data['sample_rate'],pitch)
plt.xlim([0,0.02])
Out[19]:
(0, 0.02)

ピッチ抽出の方針

人の声の限界として100[Hz]以下は出せない。そのため。自己相互相関は0.01まで見えればいい。高速化のため、音声の一部を切り出して自己相互相関を取り、その中で0.002~0.01の間の最大値をピークとする。

発話区間を抽出し、その中央部分を抜き取ってくる関数

In [20]:
def extract_speak_mid(snd, extract_length=0.05):
    snd_index = extract_speak(snd)
    extract_start = sum(snd_index)/2
    extract_end = sum(snd_index)/2+extract_length*snd_data['sample_rate']
    return snd_data['data'][snd_index][extract_start:extract_end]
In [21]:
snd_data = snd_array[0]
corr_data = extract_speak_mid(snd_data)
pitch = autocorr(np.float_(corr_data))

抽出区間

In [22]:
plt.plot(corr_data)
Out[22]:
[<matplotlib.lines.Line2D at 0x1524d730>]

自己相互相関。音声全体を使った場合と形が変わっているが、最初のピークは0.005の位置に残っている。

In [23]:
plt.plot(np.arange(0.0, pitch.shape[0], 1)/snd_data['sample_rate'], pitch)
plt.xlim([0,0.02])
Out[23]:
(0, 0.02)

0.002~0.01の間の最大値を取得する関数

In [24]:
def extract_pitch(corr, fs, min_pitch=0.002, max_pitch=0.01):
    import numpy as np
    return (np.argmax(corr[min_pitch*fs:max_pitch*fs]) + min_pitch*fs) / fs
In [25]:
extract_pitch(pitch, snd_data['sample_rate'])
Out[25]:
0.0048117913832199542

すべての音声に振幅の平均値と、ピッチの情報を追加

In [26]:
for snd_data in snd_array:
    use_index = extract_speak(snd_data)
    snd_data['amp_mean'] = np.mean(snd_data['amp'][use_index])
    pitch = autocorr(np.float_(extract_speak_mid(snd_data)))
    snd_data['pitch'] = extract_pitch(pitch, snd_data['sample_rate'])

振幅平均とピッチで音声データを散布図表示。元気のいい声は振幅が大きくピッチが小さい(高音)。元気のない声はそれと逆の性質を示している。綺麗に分かれているため、分離は容易と考えられる。

In [27]:
plt.hold(True)
for snd_data in snd_array:
    if 'weak' in snd_data['fname']:
        col = 'b'
    elif 'strong' in snd_data['fname']:
        col = 'r'
    plt.plot(snd_data['amp_mean'], snd_data['pitch'], '.', color=col)
plt.xlabel('amp_mean')
plt.ylabel('pitch')
plt.hold(False)

scikit-learnのSVMに合わせて学習データを作成。

In [83]:
vol_array = np.zeros(len(snd_array))
peak_array = np.zeros(len(snd_array))
label = np.zeros(len(snd_array))
for i, snd_data in enumerate(snd_array):
    vol_array[i] = snd_data['amp_mean']
    peak_array[i] = snd_data['pitch']
    if 'weak' in snd_data['fname']:
        label[i] = 1

data_training = [[x/1000, y*1000] for (x, y) in zip(vol_array, peak_array)] # 桁あわせ
label_training = [int(x) for x in label]
In [84]:
data_training
Out[84]:
[[6.4325938705373611, 4.8117913832199539],
 [11.985085977119212, 4.2448979591836737],
 [11.261728768038482, 4.2222222222222214],
 [4.2273977415806705, 4.8344671201814062],
 [4.7016150086609203, 4.0408163265306118],
 [6.1255228149128467, 4.2222222222222214],
 [4.4952506407786501, 6.2857142857142847],
 [4.6162763278823267, 5.9909297052154189],
 [5.4627958332919038, 5.7868480725623579],
 [2.4536728651327486, 6.5351473922902494],
 [2.0822530068408001, 6.3764172335600904],
 [5.0163243388526091, 6.2403628117913827]]
In [85]:
label_training
Out[85]:
[0, 0, 0, 0, 0, 0, 1, 1, 1, 1, 1, 1]

学習

In [86]:
estimator = SVC(C=10.0, kernel='linear')
estimator.fit(data_training, label_training)
Out[86]:
SVC(C=10.0, cache_size=200, class_weight=None, coef0=0.0, degree=3, gamma=0.0,
  kernel='linear', max_iter=-1, probability=False, random_state=None,
  shrinking=True, tol=0.001, verbose=False)

予測

In [87]:
label_prediction = estimator.predict(data_training)
print(label_prediction)
[0 0 0 0 0 0 1 1 1 1 1 1]

SVMの学習内容を表示。主にピッチで分離されている。上(ピッチが広い)ほど低周波のため、高周波は元気な声、低周波は元気のない声と分類できる。

In [91]:
X = np.array(data_training)
h = .02
x_min, x_max = X[:, 0].min()-1, X[:, 0].max()+1
y_min, y_max = X[:, 1].min()-1, X[:, 1].max()+1
xx, yy = np.meshgrid(np.arange(x_min, x_max, h),
                     np.arange(y_min, y_max, h))
Z = estimator.predict(np.c_[xx.ravel(), yy.ravel()])
Z = Z.reshape(xx.shape)
plt.contourf(xx, yy, Z, cmap=plt.cm.Paired, alpha=0.8)

plot_data = np.array(data_training)
label = np.array(label_training)
plt.hold(True)
plt.plot(plot_data[label==0, 0], plot_data[label==0, 1], '.', color='r')
plt.plot(plot_data[label==1, 0], plot_data[label==1, 1], '.', color='b')
plt.hold(False)
plt.xlabel('amp')
plt.ylabel('pitch')
Out[91]:
<matplotlib.text.Text at 0x11f5e310>