アカデミック

【超初心者向け】wavファイルからpythonでMCMS(時間/周波数分解能の異なるメル周波数スペクトログラム)を抽出する

時間周波数分解能の異なるメル」周波数スペクトログラムを抽出したいゾ。

時間周波数分解能ってどうやて変えるんだろう?

今回は,基本的な音響特徴量であるログメルスペクトログラムとMFCCをPythonで抽出する方法をお伝えしていこうと思います。本記事はpython実践講座シリーズの内容になります。その他の記事は,こちらの「Python入門講座/実践講座まとめ」をご覧ください。

【超初心者向け】python入門講座/実践講座まとめ目次入門講座 1.実行環境 2.文字の出力 3.データ型 4.変数 5.更新と変換 6.比較演算子 7.論理演算子 8.条件分岐 9.リスト...

今回は他の記事とは違ってライブラリに頼ります。コーディングに関して未熟な部分がたくさんあると思いますので,もし何かお気づきの方は教えていただけると幸いです。また,誤りについてもご指摘していただけると非常に助かります。

結論

●時間分解能=STFTの窓幅 or ホップサイズで制御
●周波数分解能=STFTの窓幅で制御
※窓幅でトレードオフが決まります。

STFTの概要

STFTからスペクトログラムを作る操作の概要。
窓幅とホップサイズが時間分解能と周波数分解能に影響を与えることが分かる。

メル周波数スペクトログラムの概要

スペクトログラムをメル周波数スペクトログラムに変換する操作の概要。メルフィルタの窓の数が変換後の周波数方向のビンの数に相当することが分かる。

メル周波数スぺクトログラムとは,簡単にまとめると以下のような特徴量です。

●wavファイル(生の音)にSTFT(短時間フーリエ変換)を施して
●メルフィルタバンクを適用した特徴量

短時間フーリエ変換では,音の周波数に関する時間変化を表すスペクトログラムという特徴量を得ることができます。そのスペクトログラムを,人間の聴覚特性にフィットするような形に整形(=メルフィルタバンクを適用する)したものがメル周波数対数スぺクトログラムです。

おそらく一般的には「メル周波数スペクトログラム」と呼ばれていますが,この記事では対数変換(デシベルに変換)を施すことを強調するために「メル周波数対数スペクトログラム」と記述することにします。

ちなみに,MFCC(メル周波数ケプストラム係数)は,メル周波数対数スペクトログラムをさらに離散コサイン変換してケプストラム領域に飛ばし,そのうち低次元の係数をとってきた特徴量のことを指しています。一昔前までは,音響特徴量といえばMFCCでしたが,深層学習が利用されるようになってからはメル周波数対数スペクトログラムもよく利用されています。

時間分解能

時間分解能は「STFTの窓幅」によって決められます。つまり,短時間フーリエ変換を行う際の窓関数の幅を小さくすれば,それだけ時間方向の解像度が上がります。これは逆を考えると分かりやすく,窓幅を大きくすればそれだけ長い時間情報を1つの情報に凝縮してしまっているため,それだけ時間方向に荒くスキャンしているということになります。

また,時間分解能はSTFTの窓関数のホップサイズにも依存します。つまり,ホップサイズが小さければ小さいほど時間方向に細かくスキャンすることになるため,時間分解能は高くなります。一方で,ホップサイズを大きくすればそれだけ時間方向に荒くスキャンすることになるため,時間分解能は低くなります。

一般的に,ホップサイズは「窓幅の半分」「窓幅の1/4」などのように窓幅を基準にして考えられることが多いです。そのため,時間分解能は窓幅に依存すると言われることが多いです。

一方で,周波数分解能も「STFTの窓幅」によって決められます。これは非常に単純で,窓幅を大きくすればするほどたくさんの周波数帯を含めることができるため,周波数方向に細かくスキャンしていることになります。ちなみに,ホップサイズを変えてもスキャンできる周波数帯は変わらないため,周波数分解能は変わりません。

以上を踏まえれば,時間分解能と周波数分解能はどちらも窓幅に依存し,トレードオフの関係にあることが分かります。ちなみに,これを「不確定性原理」と呼びます。

ログメル周波数を抽出する関数

def extract_logmel(wav, sr, n_mels, n_fft, hop_length, fmin, fmax): #(logmel_dim, time_frame)
    audio, _ = librosa.load(wav, sr=sr)
    logmel = librosa.power_to_db(librosa.feature.melspectrogram(y=audio, sr=sr, n_mels=n_mels, n_fft=n_fft, hop_length=hop_length, fmin=fmin, fmax=fmax))
    return logmel

実際に走らせるコード

x = "hoge.wav"

fs = 44100
n_fft1 = int(0.023 * fs) # 窓幅は23ms
n_fft2 = int(0.046 * fs) # 窓幅は46ms
n_fft2 = int(0.093 * fs) # 窓幅は93ms
hop_length = int(0.01* fs) #ホップサイズは10ms
fmin = 27.5 # メル周波数の下限は27.5Hz
fmax = 16000 # メル周波数の上限は16kHz

logmel1 = extract_logmel(wav=x, sr=fs, n_mels=80, n_fft=n_fft1, hop_length=hop_length, fmin=fmin, fmax=fmax)
logmel2 = extract_logmel(wav=x, sr=fs, n_mels=80, n_fft=n_fft2, hop_length=hop_length, fmin=fmin, fmax=fmax)
logmel3 = extract_logmel(wav=x, sr=fs, n_mels=80, n_fft=n_fft3, hop_length=hop_length, fmin=fmin, fmax=fmax)

窓幅を変えることで時間周波数分解能の異なるスペクトログラムを生成しました。n_mels=80としているのは80個の窓を用意しているということです。そうすることで,周波数方向に80個のビンが出現します。また,fminとfmaxという引数を指定することでめるフィルタの下限と上限を設定することができます。詳しくはlibrosaのドキュメントをご覧ください。

【librosaドキュメント】
librosa.filters.mel
librosa.feature.melspectrogram

まとめ

時間周波数分解能の異なるスペクトログラムを抽出する方法をお伝えしました。

COMMENT

メールアドレスが公開されることはありません。 * が付いている欄は必須項目です