Wizard Notes

Python, JavaScript を使った音楽信号分析の技術録、作曲活動に関する雑記

PythonでSSM(自己類似度行列)を使って楽曲構成(Aメロ,サビなど)を分析

計算機を使った音楽分析技術の一つとして、楽曲構成(Aメロ、Bメロ、サビ)を分析する方法があります。

しかし、日本語で技術の具体的な実現方法や実装について説明のある文献はあまり多くありません。また、Pythonの音楽分析ライブラリ LibROSA にはまだ実装されていません*1

そこで、PythonでLibROSAを使いつつ、楽曲構成(Aメロ、Bメロ、サビ)の分析を実装してみました。

以下の文献および実装を元に書きましたので、技術の細かい所はこちらをご参照ください。

概要

今回の実装では、処理の流れは以下のようになっています。

  1. リサンプリング(ダウンサンプリング)
  2. 非定常成分(打楽器音)の除去
  3. クロマグラムの算出
  4. 自己類似度行列 (SSM: Self-Similarity-Matirx)の算出
  5. SSM のスムージング
  6. SSMをTime-lag表現に変換
  7. Novelty Function算出

f:id:Kurene:20210203082840p:plain

各処理の詳細

1. リサンプリング(ダウンサンプリング)

ダウンサンプリングは必ずしも必要な処理ではありませんが、以下の理由でリサンプリングしています。

  • クロマグラム算出の計算量を減らす
  • 高い周波数成分はクロマグラムではノイズとなることがある

今回の実装では、元のサンプリング周波数の10分の1にダウンサンプリングしています。

sr = sr_orig // 10
x = resampy.resample(x, sr_orig, sr)

2. 非定常成分(打楽器音)の除去

これもマストではないですが、打楽器成分はクロマグラムを使った分析ではノイズとなりやすいので入れました。

手法としては、HPSS(調波打楽器音分離)を使っています。

HPSSの詳細については以下の記事をご参照ください。

実装としては、LibROSAを使えば以下のように簡単に書くことができます。

n_fft      = 512
hop_length = n_fft // 2
x_stft = librosa.stft(x, n_fft=n_fft, hop_length=hop_length)
x_stft_harm, _ = librosa.decompose.hpss(x_stft, margin=1.2)
y = librosa.util.fix_length(librosa.istft(x_stft_harm, hop_length=hop_length), len(x))

3. クロマグラムの算出

クロマグラムについては、以下の記事をご覧ください。

なお、楽曲構成の分析では、クロマグラムのほかにも

を使う方法があります。

今回は、楽曲構成における各セクションはコード進行が異なるはずであると仮定し、クロマグラムを利用しました。

注意として、4 の自己類似度行列(SSM)は、時間フレーム数×時間フレーム数の行列となります。

従って、楽曲の時間が長いと計算量が非常に大きくなります。

計算量を減らす手段としては、まずはクロマグラム等の特徴表現行列の時間フレーム数が少なくなるようにhop_lengthなどを調整するのがよいと思います。

y_chroma = librosa.feature.chroma_cqt(y=y, sr=sr, hop_length=hop_length, n_octaves=6)
y_chroma = y_chroma / np.max(y_chroma)
np.nan_to_num(y_chroma, nan = 0.0)

4. 自己類似度行列 (SSM: Self-Similarity-Matirx)の算出

楽曲構造の分析に利用するSSMでは、2つの時間フレーム間の特徴表現ベクトルの類似度が算出されています。その結果、SSMの内部には対角方向へのパスと、ブロック構造が表れます。これが、楽曲のコードやテンポがどのように変化しているかを分析する手がかりとなります。

f:id:Kurene:20210203093707p:plain
https://www.audiolabs-erlangen.de/resources/MIR/FMP/C4/C4S2_SSM.htmlより引用

SSMの実際の例としては、ページ最下部の分析結果をご覧ください。

なお、類似度としては、

などの利用が考えられます。

無音区間でのゼロ除算などによるNaNの出現に注意しつつ算出してください。また、算出される行列は、時間フレーム×時間フレームの正方行列となります。

本実装では、コサイン類似度を利用しました。

def cossim(mx):
    norm = (mx * mx).sum(0, keepdims=True) ** .5 + 1e-16
    return mx.T @ mx / norm / norm.T
ssm = cossim(y_chroma**2)

5. SSM のスムージング

SSM をスムージング(平滑化処理)することで、7で算出するNovelty functionを滑らかにし、楽曲構造を分析しやすい結果にすることができます。

スムージング方法としては、

  • Medianフィルタ
  • ガウシアンフィルタ
  • 平均

など、画像処理などでよく使われている平滑化フィルタを採用するとよいと思います。

本実装では、参考文献 [1]に紹介されている平均を取る方法を採用しました。

@jit(nopython=True, nogil=True)
def mx_smoothing(mx, N, L):
    forward  = np.zeros((N+L, N+L))
    forward[0:N, 0:N]      = mx
    for k in range(1, L):
        mx += forward[k:k+N, k:k+N]
    mx /= L

len_smooth = 50
N = ssm.shape[0]
ssm_smooth = ssm.copy()
mx_smoothing(ssm_smooth, N, len_smooth)

6. SSMをTime-lag表現に変換

SSMの対角方向のパスが水平方向になるように変換することで、水平軸を時間方向としてどのように楽曲構造が変わっているか分かりやすい形にすることができます、また、分析のための計算効率を向上させることができます。

f:id:Kurene:20210203094044p:plain
https://www.audiolabs-erlangen.de/resources/MIR/FMP/C4/C4S4_StructureFeature.htmlより引用

参考文献 [1] によると、循環型のTime-lag表現はnp.rollを使うことで簡単に実装できます。

N = ssm.shape[0]
struct = np.zeros((N, N))
for k in range(N):
    struct[:, k] = np.roll(ssm_smooth[:, k], -k)

7. Novelty Function算出

6 のtime-lag表現行列より、水平方向に眺めてみると、楽曲構造の変わり目では行列の列ベクトルの構造が大きく変わることに気が付きます。

従って、特徴表現の変化点を検出するためのNovelty Functionは直前の時間フレームとの誤差として定義することができます。そのため、Novelty Functionの値が大きくなっている時間フレームでは、楽曲構造に変化があると考えられます。

for k in range(1, N):
    novelty[k] = np.linalg.norm(struct[:, k] - struct[:, k-1])

デモ

Perfume の無限未来を分析してみました。

分析結果(Novelty Function)

Novelty Function の値が大きくなっているところで楽曲構成におけるセクションが変わっていることが分かります。また、それ以外にも、繰り返しの頭で値が大きくなっていることが聞き比べることで分かります。これは、同じセクション内で同じコード進行の繰り返しがあるためです。

f:id:Kurene:20210203092509p:plain

分析結果(SSM, Time-lag 表現)

f:id:Kurene:20210203091841p:plain

なお、タイムラインのプロットは以下をご参考ください。

matplotlibでタイムラインチャートをプロット - Wizard Notes

*1:2021年2月現在