【Python/PyRoomAcoustics】ILRMAでブラインド音源分離

ブラインド音源分離の手法の一つであるILRMA（独立低ランク行列分析に基づくブラインド音源分離）がPythonモジュールのPyRoomAcousticsで使えたので試してみました。

今回は話者分離（スピーカーダイアライゼーション）での使用です。

収録する際のマイクロホンや音源の種類や位置，部屋の形状等の情報（混合系）が未知の状況でも音源を分離することができる統計的信号処理アルゴリズムを**ブラインド音源分離（blind audio source separation: BASS）** や、**ブラインド信号源分離（blind source separation: BSS）** といいます。

# ILRMAとは
ILRMAについての詳細はこちらから。

- 論文　http://d-kitamura.net/pdf/paper/Kitamura2017EA10.pdf
- 北村さんの音源分離デモ http://d-kitamura.net/demo-ILRMA.html
- 北村さんのオリジナルソースコード（MATLAB）　https://github.com/d-kitamura/ILRMA

複数のマイクロホンを使うことでブラインド音源分離ができます。優決定系（欲しい音源数より収録マイクの方が多い）で使えます。欲しい音源数より収録マイクの方が少ない劣決定系では**MNMF**などを使ってください。

# PythonでILRMA
理論は論文で確認するとして、PyRoomAcousticsで簡単に実装できるので紹介します。

- PyRoomAcousticsのIRLMA https://pyroomacoustics.readthedocs.io/en/pypi-release/pyroomacoustics.bss.ilrma.html

```pip install pyroomacoustics``` でインストール

ColaboratoryなりJupyterLabなりで実行してみてください

追記：ありがたいことに、北村先生ご本人からアドバイスを得られたので、少し書き直しました。

```ilrma.py
#cal
import pyroomacoustics as pr
import scipy.signal as sig
from sklearn.decomposition import PCA
#system
import glob
import numpy as np
#audio
from IPython.display import display, Audio
import scipy.io.wavfile as wav

#パラメータ
data_dir = "content" #音源wavの入ってるディレクトリ
N = 2 #音源数
RT = 0.3 #室内残響時間

#データ読み込み
data_list = [] #ORIGINAL
for f in glob.glob(data_dir+'/*.wav'):
    RATE,data = wav.read(f)
    data_list.append(data)
t_data = np.stack(data_list, axis=1)


# 主成分分析（PCA）で次元下げ
pca = PCA(n_components=N, whiten=True)
H = pca.fit_transform(t_data)

#STFT()
#STFT窓長を室内残響時間より長くします
n = RT*RATE
i = 1
while ( i * 2 <= n):
    i *= 2
seg = i * 2
stft_list = [] #STFT
for i in range(N):
    _,_,Z = sig.stft(H[:,i],nperseg=seg)
    stft_list.append(Z.T)
f_data = np.stack(stft_list, axis=2)

# ILRMA計算
Array = pr.bss.ilrma(f_data, n_src=None, n_iter=100, proj_back=True, W0=None, n_components=2, return_filters=False, callback=None)

# 出力
sep = [] #ここに音源データが入るので煮るなり焼くなり
for i in range(N):
    x=sig.istft(Array[:,:, -(i+1)].T, nperseg=seg)

    sep.append(x[1])
    display(Audio(x[1], rate=RATE))
```

##出力結果
２人での会話を４chで収録し計算してみました。２人の音声が分離されてますね。![After.png](https://qiita-image-store.s3.ap-northeast-1.amazonaws.com/0/543519/485f05b9-6ad3-fa31-5387-469866910ecf.png)


#簡単な説明

## STEP1 事前処理
まず分離したい音源以上の数のマイクロホンを用いて収録します。
マイクがない場合、シミュレーションで作るのもありです。
[【Python/PyRoomAcoustics】Pythonで室内音響シミュレーション](https://qiita.com/oozzZZZZ/items/8dd90f17aad100b8ec43)

ILRMAは決定的な系で解くので、いきなりPyRoomAcousticsのIRLMAにぶち込むのではなく主成分分析（PCA）で音源数次元へデータを減らしてやる必要があります。PCAはscikit-learnに簡単関数があります。

```pca.py
# 主成分分析（PCA）で次元下げ
pca = PCA(n_components=N, whiten=False)
H = pca.fit_transform(t_data)
```

## STEP2 STFT
短時間フーリエ変換（STFT）をしてデータを時間フレーム×周波数ビン×次元数の３次元配列に。
STFT窓長については理論上室内の残響時間より長くしてやる必要があるそうです。
（これを知らずにデフォルトで計算したところうまくいきませんでした）

```STFT.py
#STFT()
#STFT窓長は室内残響時間より長くします
n = RT*RATE
i = 1
while ( i * 2 <= n):
    i *= 2
seg = i * 2
stft_list = [] #STFT
for i in range(N):
    _,_,Z = sig.stft(H[:,i],nperseg=seg)
    stft_list.append(Z.T)
f_data = np.stack(stft_list, axis=2)
```

## STEP3 ILRMA計算
``` n_iter``` はイテレーション（反復計算数）です。デフォルトは20回。100回くらいやれば大抵収束するらしいです。他のパラメータは[PyRoomAcousticsのドキュメント](https://pyroomacoustics.readthedocs.io/en/pypi-release/pyroomacoustics.bss.ilrma.html)をみてください。


```ilrma.py
# ILRMA計算
Array = pr.bss.ilrma(f_data, n_src=None, n_iter=100, proj_back=True, W0=None, n_components=2, return_filters=False, callback=None)
```
## STEP4 ISTFT
逆STFTして出力！　おしまい！！

```ISTFT.py
# 出力
sep = [] #ここに音源データが入るので煮るなり焼くなり
for i in range(N):
    x=sig.istft(Array[:,:, -(i+1)].T, nperseg=seg)
    sep.append(x[1])
    display(Audio(x[1], rate=RATE))
```

# 最後に

音声を聞いてみると残響だけ入っていましたがこの辺りは残響処理でなんとかなりそうです。
北村先生アドバイスをいただきありがとうございました！