前半は独立成分分析について解説し、後半でscikit-learnを用いて音声の分離を行います。

独立成分分析とは

独立成分分析とは複数の独立な変数がある比率で混合されているものを独立な変数に分離する分析手法です。
数式としては、次のように定式化されます。

x_i=\sum^q_{i=1}a_{ij} s_j ~~~ (i = 1,2,...,q)
s_jはそれぞれ、独立な非ガウス分布に従うとします。
これを行列を用いて表すと
x = As
となります。ここでAを混合行列と呼びます。

独立成分分析では混合行列Aを推定します。

（！注意！)
Aの列の交換、尺度変化については識別できません.。
例1) 列の交換
\left(\begin{array}{c} x_1 \\ x_2 \end{array}\right) = \left(\begin{array}{cc}a_{11} & a_{12} \\a_{21} & a_{22}\end{array}\right)\left(\begin{array}{c}s_1 \\ s_2 \end{array}\right), \left(\begin{array}{c} x_1 \\ x_2\end{array}\right) = \left(\begin{array}{cc}a_{12} & a_{11} \\a_{22} & a_{21}\end{array}\right)\left(\begin{array}{c} s_2 \\ s_1 \end{array}\right)

例2) 尺度
\left(\begin{array}{c} x_1 \\ x_2 \end{array}\right) = \left(\begin{array}{cc}a_{11} & a_{12} \\a_{21} & a_{22}\end{array}\right)\left(\begin{array}{c} s_1 \\ s_2 \end{array}\right),\left(\begin{array}{c} x_1 \\ x_2 \end{array}\right) = \left(\begin{array}{cc}c*a_{11} & a_{12} \\c*a_{21} & a_{22}\end{array}\right)\left(\begin{array}{c} s_1/c \\ s_2 \end{array}\right)
上記の例のように、データからは列の交換、尺度の変化については区別できず、推測することはできません。

混合行列の推定法

y = Wx \ y = s
となるような復元行列Sを推定します。
ここで、W = A^{-1}の時 y = xとなります。
sはそれぞれ独立であったので、yが独立になるようにWを推定します。

そこで、独立性の指標として相互情報量を用います。
相互情報量は以下で定義されます
I(y) = \left(\sum_{i = 1}^q H(yi)\right) - H(y) \ H(y) = E[-\log p(y)]
相互情報量はyの各成分が独立の時のみ0となります。

上記の期待値をサンプルの平均により近似し、次のように相互情報量を推定します。
\hat{I}(y) = \left(\sum_{i = 1}^q \frac{1}{n} \sum_{m=1}^n -\log p(y_i^m)\right) - \frac{1}{n} \sum_{m=1}^n -\log p(y^m)
y = Wx
上記を最小化するようにWを推定します。
推定方法としては不動点法などが知られています。

実際に独立成分分析を行う

では実際に音声を分離していきましょう

import numpy as np
import matplotlib.pyplot as plt
import IPython as IP
from sklearn.decomposition import FastICA

まずは音声を生成します。
高い音と低い音をそれぞれ生成し、混合していきます。

fs = 5000 #1秒に何回聞くか
duration = 1 # 何秒聞くか

#音声生成
f1 = 500 #振動数
f2 = 200 #振動数
t = np.linspace(0., duration, int(fs * duration))  # ( start, stop, num of data )
x1 = np.sin(f1 * (2. * np.pi) * t)
x2 = np.sin(f2 * (2. * np.pi) * t)

#波形を見てみる。
plt.figure(figsize=(20,5))
plt.plot(x1[:100])
plt.plot(x2[:100])

# 音を聞いてみる
# 音量注意！
IP.display.display(IP.display.Audio(x1, rate=fs, autoplay=True))

# 音を聞いてみる
# 音量注意！
IP.display.display(IP.display.Audio(x2, rate=fs, autoplay=True))

次に音声を混ぜてみます。

# 音を混ぜてみる
x_mix1 = 2*x1 + x2
x_mix2 = x1 + x2*5
plt.figure(figsize=(20,5))
plt.subplot(2,1,1)
plt.plot(x_mix1[:100])
plt.subplot(2,1,2)
plt.plot(x_mix2[:100])

# 音を聞いてみる
# 音量注意！
IP.display.display(IP.display.Audio(x_mix1, rate=fs, autoplay=True))

# 音を聞いてみる
# 音量注意！
IP.display.display(IP.display.Audio(x_mix2, rate=fs, autoplay=True))

上記で、音声の混合ができたので実際に独立成分分析にかけてみましょう

# 実際に独立成分分析を行う
X = np.concatenate((x_mix1.reshape(-1,1),x_mix2.reshape(-1,1)),1)
transformer = FastICA(n_components=2,random_state=0)
X_transformed = transformer.fit_transform(X)

たった3行で出来てしまいました。
分離された音声をそれぞれ見てみます。

# 復元された音声を見てみる
plt.figure(figsize=(20,5))
plt.subplot(2,1,1)
plt.plot(X_transformed[:,0][:100])
plt.subplot(2,1,2)
plt.plot(X_transformed[:,1][:100])

# 音を聞いてみる
# 音量注意！
IP.display.display(IP.display.Audio(X_transformed[:,0], rate=fs, autoplay=True))

# 音を聞いてみる
# 音量注意！
IP.display.display(IP.display.Audio(X_transformed[:,1], rate=fs, autoplay=True))

尺度のせいか、音量の変化が著しいですが、音声を聞いたみた感じはちゃんと分離されていると思います。

次回は独立成分分析を用いて行われる統計的因果探索の手法であるLinGAMについて解説します。

参考文献

清水 昌平 , 統計的因果探索(機械学習プロフェッショナルシリーズ), 講談社,2017