MFCC
MFCCの手順を簡潔にまとめた。
実際に使用する際はlibrosaなどのライブラリを用いて１行で実装するのがいいと思う。
MFCCとは
音声認識で使用される特徴量抽出の方法。
流れ
- wavファイルから音声波形を読み取る
- 短い時間で切り出す
- 高周波数を強調する
- 窓関数を適応する
- フーリエ変換する
- メルフィルタバンクをかける
- 離散コサイン変換
1.wavファイルから音声波形を読み取る
def read_wave(file_name):
"""
パス先のwaveファイルを読み込む
"""
wf = wave.open(file_name, "r")
fs = wf.getframerate()
x = wf.readframes(wf.getnframes())
x = np.frombuffer(x, dtype="int16") / 32768.0
wf.close()
return x, float(fs)
wav, fs = read_wave(DATA_DIR + "a.wav")
2.短い時間で切り出す
def cut_wave(wave_data, cuttime=0.04):
center = len(wave_data) / 2 # 中心のサンプル番号
cuttime = 0.04 # 切り出す長さ [s]
wav = wave_data[int(center - cuttime/2*fs) : int(center + cuttime/2*fs)]
return wav
wave_data, fs = read_wave(DATA_DIR + "a.wav")
wav = cut_wave(wave_data)
※グラフを出力するコードは省略しています。
3.高周波数を強調する
def preEmphasis(wave, p=0.97):
"""プリエンファシスフィルタ"""
# 係数 (1.0, -p) のFIRフィルタを作成
return scipy.signal.lfilter([1.0, -p], 1, wave)
wave_data, fs = read_wave(DATA_DIR + "a.wav")
wav = cut_wave(wave_data)
wav = preEmphasis(wav)
4.窓関数を適応する
def window(wave):
"""
窓関数
"""
hanningWindow = np.hanning(len(wave))
wave = wave * hanningWindow
return wave
wave_data, fs = read_wave(DATA_DIR + "a.wav")
wav = cut_wave(wave_data)
wav = preEmphasis(wav)
wav = window(wav)
5.フーリエ変換する
def fft(wave):
"""
波形に対してフーリエ変換を行う関数
:param: wave_data
:return Adft: 振幅スペクトル
:return Pdft: パワースペクトル
"""
nfft = 2048 # FFTのサンプル数
# 離散フーリエ変換
dft = np.fft.fft(wave, n)
# 振幅スペクトル
Adft = np.abs(dft)[:int(nfft/2)]
# パワースペクトル
Pdft = np.abs(dft)[:int(nfft/2)] ** 2
return Adft, Pdft
wave_data, fs = read_wave(DATA_DIR + "a.wav")
wav = cut_wave(wave_data)
wav = preEmphasis(wav)
wav = window(wav)
Adft, Pdft = fft(wav)
6.メルフィルタバンクをかける
def hz2mel(f):
"""Hzをmelに変換"""
return 1127.01048 * np.log(f / 700.0 + 1.0)
def mel2hz(m):
"""melをhzに変換"""
return 700.0 * (np.exp(m / 1127.01048) - 1.0)
def melFilterBank(spec, fs=8820.0, nfft=2048, numChannels=20):
"""メルフィルタバンクを作成"""
# ナイキスト周波数（Hz）
fmax = fs / 2
# ナイキスト周波数（mel）
melmax = hz2mel(fmax)
# 周波数インデックスの最大数
nmax = nfft / 2
# 周波数解像度（周波数インデックス1あたりのHz幅）
df = fs / nfft
# メル尺度における各フィルタの中心周波数を求める
dmel = melmax / (numChannels + 1)
melcenters = np.arange(1, numChannels + 1) * dmel
# 各フィルタの中心周波数をHzに変換
fcenters = mel2hz(melcenters)
# 各フィルタの中心周波数を周波数インデックスに変換
indexcenter = np.round(fcenters / df)
# 各フィルタの開始位置のインデックス
indexstart = np.hstack(([0], indexcenter[0:numChannels - 1]))
# 各フィルタの終了位置のインデックス
indexstop = np.hstack((indexcenter[1:numChannels], [nmax]))
filterbank = np.zeros((int(numChannels), int(nmax)))
for c in np.arange(0, numChannels):
# 三角フィルタの左の直線の傾きから点を求める
increment= 1.0 / (indexcenter[c] - indexstart[c])
for i in np.arange(indexstart[c], indexcenter[c]):
i = int(i)
c = int(c)
filterbank[c, i] = (i - indexstart[c]) * increment
# 三角フィルタの右の直線の傾きから点を求める
decrement = 1.0 / (indexstop[c] - indexcenter[c])
for i in np.arange(indexcenter[c], indexstop[c]):
i = int(i)
filterbank[c, i] = 1.0 - ((i - indexcenter[c]) * decrement)
mspec = np.log10(np.dot(spec, filterbank.T))
return mspec, fcenters
wave_data, fs = read_wave(DATA_DIR + "a.wav")
wav = cut_wave(wave_data)
wav = preEmphasis(wav)
wav = window(wav)
Adft, Pdft = fft(wav)
Adft, fcenters = melFilterBank(Adft)
Pdft, fcenters = melFilterBank(Pdft)
7.離散コサイン変換
def dct(spec, nceps=12):
"""
離散コサイン変換
"""
spec = scipy.fftpack.realtransforms.dct(spec, type=2, norm="ortho", axis=-1)
return spec[:nceps]
wave_data, fs = read_wave(DATA_DIR + "a.wav")
wav = cut_wave(wave_data)
wav = preEmphasis(wav)
wav = window(wav)
A_spec, P_spec = fft(wav)
A_spec, fcenters = melFilterBank(A_spec)
P_spec, fcenters = melFilterBank(P_spec)
A_spec = dct(A_spec)
P_spec = dct(P_spec)