社団法人 電子情報通信学会
THE INSTITUTE OF ELECTRONICS,INFORMATION AND COMMUNICATION ENGINEERS

信学技報
TECHNICAL REPORT OF IEICE.

［ポスター講演］マルチチャネル音源分離のためのネスト型基底・音源混合モデルに基づく時間周波数クラスタリング
板倉 光佑†坂東 宜昭†

中村 栄太†

糸山 克寿†

吉井 和佳†

河原 達也†† 京都大学大学院 情報学研究科
E-mail: †{itakura,bando,enakamura,itoyama,yoshii,kawahara}@sap.ist.i.kyoto-u.ac.jp
あらまし

本稿では，時間周波数クラスタリングによるマルチチャネル音源分離のためのネスト型混合モデルについて述べる．時間周波数クラスタリングによる音源分離では，混合音は各音源の観測モデルの混合モデルに基づいて生成されるとする．この混合モデルを推定するための特徴量として各マイクでの音の位相と各音源のパワースペクトログラムを用いることができる．提案法ではこのパワースペクトログラムを基底の混合モデルによりモデル化する．これにより，提案法では混合音が音源の混合モデルで，各音源が基底の混合モデルでモデル化される．これをネスト型基底・音源混合モデルと呼ぶ．評価実験により，音声などの混合音に対して提案法により SDR と SIR が向上することを確認した．
キーワード

マルチチャネル音源分離，時間周波数クラスタリング，潜在的ディリクレ配分法

1. は じ め にデル化を行った．また，パワースペクトログラムの低ランク性を用いたモデル化も行われている [9]．この手法では，単チャネマイクロホンアレイを用いたマルチチャネル音源分離におル音源分離でよく用いられる非負値行列因子分解 (NMF) [10]いて，これまでに多くの手法が提案されてきた．広く用いられのようにパワースペクトログラムが基底スペクトルとアクティている手法のうちの一つに，独立成分分析 (ICA) [1] がある．

ベーションの積の和を用いて表現される．マルチチャネル音源

ICA は各音源の統計的な独立性を仮定することにより分離行分離においても，空間モデルだけではなく音源モデルも考慮す列を推定する．この ICA をもとに独立ベクトル分析 (IVA) [2]

ることにより分離性能の向上が期待される．

や FastICA [3] などのさまざまな手法が提案されているが，こ本稿では，このパワースペクトログラムをさらに混合モデル

れらの手法は共通してマイク数が音源数より少ない劣決定条件を用いてモデル化した２段階のネスト型混合モデルによる音源では分離できないという問題点がある．

分離法を提案する．提案法では，NMF のようにパワースペク

これに対し，劣決定条件でも分離が可能な手法として時間・周トログラムを基底スペクトルとアクティベーションに分解し，波数クラスタリングに基づく音源分離法が着目されている [4–9]．

それに加えて基底スペクトルのスパース性を仮定する．つまり，このアプローチでは，各音源スペクトログラムが時間・周波数パワースペクトログラムの各時間周波数ビンではいずれか一つ領域でスパースであると仮定することで，混合音スペクトログの基底の成分のみが観測されるとする．これにより，提案法でラムの各時間・周波数ビンにおける観測はそれぞれいずれか一は空間モデルが音源の混合モデル，音源モデルが基底の混合モつの音源成分が直接観測されたものであるとみなす．つまり，

デルを用いてモデル化される．この２つの混合モデルをネスト

この仮定では，混合音の観測モデルは各音源の観測モデルの混型基底・音源混合モデルと呼ぶ．提案法ではネスト型基底・音合モデルとして扱われる．この混合モデルを推定するため，マ源混合モデルに対して潜在的ディリクレ配分法 (LDA) の枠組

イク間の音の位相差とパワー差を特徴量として用いた混合ワみを用いてギブスサンプリングを行うことで二つの混合モデル

トソン分布のクラスタリング [4–7] や各マイクでの音の位相とを同時に推定する．図 1 に提案法における混合音の生成モデル

パワーを特徴量として用いた混合ガウス分布のクラスタリンを示す．

グ [8, 9] による分離法が提案されている．
特徴量にマイク間のパワー差ではなく各マイクでのパワーを2. 提 案 法

用いる利点として，音源の混合過程のような空間モデルだけで提案法では，ベイズモデルを用いてネスト型基底・音源混合なく各音源のパワースペクトログラムの性質に基づいた音源のモデルのモデル化を行い，LDA の枠組みを用いてそのモデルのモデル化も同時にできるという点が挙げられる．大塚ら [8] は推定を行う．本章ではそのモデル化と推定方法について述べる．

パワースペクトログラムのスパース性に基づいたモデル化を行

2. 1 モ デ ル 化

うため，このパワーに対しスパースとなる事前分布を用いたモ

ここでは提案法のモデルの定式化について述べる．提案法で—1—それ以外のときは 0 となる．このとき，観測 xtf は次のような分布にしたがって生成される．

xtf ∼!"#$%&'()*+,-K"

ztf kNC (0, λtf k G−1fk )

(6)

k=1

ここで，空間相関行列は音源の種類ではなく音源の方向に依存するため，空間相関行列 Gf k を音源ごとに独立な変数ではな!"#$%&'()*+

く，方向ごとに独立な変数 Gf d として考える．このとき，式,-$(a) 音源モデル(6) は音源 k の方向を示すベクトル sk = [sk1 , · · · , skD ]T を用いると次のように表される．

!"#$!"#$012"

K,Dxtf ∼%&'011()&*'+,-./%&'

ztf k skd
NC (0, λtf k G−1fd )

(7)

k,d=1

ただし，sk は 1 of D 表現のベクトルであり，音源 k が方向 dにあるときは skd = 1，それ以外のときは 0 となる．
次に音源 k のパワー λtf k のモデルについて考える．提案(b) 空間モデル

図 1 提案法の生成モデル．音源モデルでは各音源のパワースペクトログラムが基底・アクティベーション・マスクにより構成される．
マスクにより一つの基底のみが選択される (黒色の部分が選択さ法では NMF のように基底スペクトル wklf とアクティベー
ション hklt を用いて λtf k を表現する．ただし，NMF では

λtf k =!

lwklf hklt と全ての基底の和でパワーを表現するのにれた部分)．空間モデルでは混合音スペクトログラムが音源スペ対し，提案法では各時間周波数ビンごとに一つの基底 l′ を選択

クトログラム・伝達関数・マスクにより構成される．マスクによして λtf k = wkl′ f hkl′ t と表す．その基底を選択するための変り一つの音源のみが選択される．

数を utf k = [utf k1 , · · · , utf kL ]T とすると式 (7) は次のようには時間領域の信号に対して短時間フーリエ変換 (STFT) を行うことにより得られる時間周波数領域の信号に対してモデル化を行う．まず K 個の音源を M 個のマイクを用いて録音するとし，時刻 t，周波数 f での観測 xtf と音源信号 ytf を以下のように定義する．

xtf = [xtf 1 , · · · , xtf M ]T ∈ CMytf = [ytf 1 , · · · , ytf K ]T ∈ CK(1)(2)

このとき周波数領域での瞬時混合を仮定すると，観測は以下のように表される．

xtf =K!

af k · ytf k(3)k=1

ただし，af k は周波数 f での音源 k の伝達関数である．ここで，ytf k が次のような複素ガウス分布に従うとする．

ytf k ∼ NC (0, λtf k )(4)

2λtf k = E[ytfk ] は時刻 t，周波数 f での音源 k のパワーを表す．

このとき音源 k のみを観測した時の観測 xtf k は次のような複素ガウス分布に従う．

xtf k ∼ NC (0, λtf k G−1
fk )

(5)

表される．

"

K,D,L

xtf ∼

#NC 0, wklf hklt G−1
fdk,d,l=1$ztf k skd utf kl(8)

ただし，utf k は 1 of L 表現のベクトルであり，基底 l が用いられるときは utf kl = 1，それ以外のときは 0 となる．

2. 2 事前分布の設計
提案法では式 (8) のパラメータに対し，ぞれぞれ適切な事前分布を与えることで推論を行う．ここではその事前分布の与え方について述べる．まず，ztf ，sk ，utf k はクラスタリングによる推論を行うためにカテゴリカル分布から生成されるとする：

ztf | πt ∼ Categorical(πt )sk | φ ∼ Categorical(φ)

utf k | ψtk ∼ Categorical(ψtk )

(9)
(10)
(11)

ここで，ハイパーパラメータ πt ，φ，ψtk は観測に依存して変動するため推論を必要とする．したがって πt ，φ，ψtk はカテゴリカル分布と共役なディリクレ分布から生成されるとする：

πt ∼ Dirichlet(aπ0 1K )(12)

ψtk ∼ Dirichlet(aψ0 1L )(14)

φ ∼ Dirichlet(aφ0 1D )(13)

ここで，G−1f k は周波数 f での音源 k の空間相関行列であり，−1H
Gf k = af k aHはエルミート共役を示す．
f k である．ただし ∗ここで，1N は要素が全て 1 の N 次元のベクトルとし，a∗0 はち，各時間周波数ビンにおいて観測される音は高々一つである

前分布として次のような分布を与える．

ここで，音源スペクトログラムがスパースである，すなわとし，そのときの混合音の観測モデルについて考える．まず各時間周波数ビンにおいて観測される音源を示すための変数を

ztf = [ztf 1 , · · · , ztf K ]T とする．ただし，ztf は 1 of K 表現のベクトルであり，音源 k が観測されるときは ztf k = 1 となり，ハイパーパラメータとする．また，空間相関行列 Gf d ，基底wklf ，アクティベーション hklt は式 (8) と共役になるように事
Gf d ∼ WC (ν, G0f d )(15)hklt ∼ Gamma(ah0 , bh0 )(17)w
wklf ∼ Gamma(aw0 , b0 )(16)

—2—ここで，ν ，a∗0 ，b∗0 はハイパーパラメータであり，WC は複素0ウィシャート分布 (付録参照) とする．

2. 3 推論

提案法では，観測データ集合 X に対するすべてのパラメータ1100の事後分布 p(G, Z, S, U , π, φ, ψ, W , H|X) を最大とするパラ80
280メータを求めることを目標とする．ただし，これらのパラメー

428タを解析的に求めることは困難なので，提案法ではこれらのパラメータをギブスサンプリングにより求めることとする．ただ3し，事後分布に含まれるパラメータは G，Z ，S ，U ，W ，H ，

ψ ，π ，φ の 9 つがあるが，提案法ではこのうちの π ，φ，ψ は積分消去を行い，残りのパラメータ Θ = {G，Z ，S ，U ，W ，H} を求めることとする．

ギブスサンプリングではそれぞれのパラメータの事後分布を求め，それらの事後分布からサンプリングを繰り返すことにより推定を行う．それぞれの事後分布は事前分布と尤度関数の積により求めることができ，以下のようになる．

Gf d | X, Θ¬Gf d ∼ WC (νf′ d , G′f d )

(18)

図2

マイク配置

ここで，n¬tftk は時刻 t 周波数 f でのサンプルを除いて，時刻 tにおいて音源 k に割り当てられた時間周波数ビンの数を表し，c¬kd は音源 k を除いて方向 d に割り当てられた音源の数を表す．
また，ntkl (nf kl ) は時刻 t(周波数 f ) において音源 k ，基底 l にk割り当てられた時間周波数ビンの数を示し，n¬tf
は ntkl からtkl

時刻 t 周波数 f で音源 k に割り当てられた要素を除いたものである．

2. 4 分離音の生成(19)

提 案 法 で は ，時 間 周 波 数 領 域 で 音 源 方 向 ご と に マ ス

sk | X, Θ¬sk ∼ Categorical(φ′k )(20)

ク を 推 定 す る こ と に よ り 分 離 音 を 生 成 す る ．ギ ブ ス サww
wklf | X, Θ¬wklf ∼ GIG(γklf, ρw
klf , τklf )(22)′ztf | X, Θ¬ztf ∼ Categorical(πtf)′utf kl | X, Θ¬utf kl ∼ Categorical(ψtfk)

hhhklt | X, Θ¬hklt ∼ GIG(γklt, ρhklt , τklt)

(21)
(23)

ここで，Θ¬∗ は Θ から ∗ の要素のみを除いた集合とする．また，GIG は一般化逆ガウス分布 (付録参照) を示す．ここで，′ハイパーパラメータ νf′ d ，G′f d ，πtf，φ′k ，ψ ，γ ∗ ，ρ∗ ，τ ∗ は次のようになる．

!

T ,K

νf′ d = ν +(24)

ztf k skdt,k=1! xtf xtf H

T ,K,L0 −1
G′−1+
f d = (Gf d )t,k,l=1

wklf hklt

ztf k skd utf kl(25)&%&' H()skd utf kl
D,Lxtf Gf d xtf" && Gf d &&′πtf
exp−k =& wklf hklt &wklf hklt
d,l=1

× (aπ0 + n¬tftk )φ′kd ="T ,F,L

(26)&%&' H()ztf k utf kl& Gf d &&& exp − xtf Gf d xtf& wklf hklt &wklf hklt

t,f,l=1× (aφ0 + c¬kd )&%&' H()ztf k skd
xtf Gf d xtf" && Gf d &&′ψtf=exp−kl& wklf hklt &wklf hklt

(27)

Dwww
γklf
= aw
0 − M nf kl , ρklf = b0w
τklf=! xtf Gf d xtf

T ,Dt,d=1hklt

ztf k skd utf klhγklt= ah0 − M ntkl , ρhklt = bh0h
τklt=! xtf Gf d xtf

F,Df,d=1wklf

ztf k skd utf kl

W (i) , H (i) , G(i) , Z (i) , S (i) , U (i) とすると，時刻 t 周波数 f でdの音源方向 d に対するマスク Mtf
は次のようになる．
dMtf
=I !
K
1!
(i) (i)z sI i=1 k=1 tf k kd(33)d
この Mtfを用いて次の式により方向 d の分離音を生成することができる．
dxdtf = Mtf
xtf

(34)

3. 評 価 実 験
提案法の分離性能を評価するため，シミュレーションにより混合した音を用いた実験を行った．比較手法として，IVA [11]，マルチチャネル NMF (MNMF) [12]，音源モデルをスパースとし空間モデルに LDA を用いた分離法 (Sparse-LDA) [8]，音源モデルに NMF，空間モデルに LDA を用いた分離法 (NMF-LDA) [9] を用いた．これに対し，提案法では音源モデル・空間モデルともに LDA を用いているためここでは提案法を LDA-LDA と呼ぶ．また，Sparse-LDA では音源数の同時推定も行うが，条件を対等にするため音源数は既知とした．

3. 1 実 験 条 件

d=1¬tf k× (aψ0 + ntkl )ンプリングで i 回目の試行により得られるサンプルを

(28)
(29)

図 2 に音源とマイクの配置を示す．残響時間 400 ms のインパルス応答を用いて 3 音源を混合した音声を用いた．マイク数は 4 とした．混合音には，音声のみの混合音と音楽のみ(30)の混合音，音声と音楽の混合音をそれぞれ 10 個ずつ使用し(31)ス文 [14] から選択した．サンプリング周波数は 16 kHz とし，(32)た．用いる音楽と音声は SISEC [13] と JNAS の音素バラン

STFT では窓幅 512 のハミング窓をシフト幅 256 で使用した．
基底数 L = 20 とし，ハイパーパラメータは，ν = M + 1，whw
haπ0 = aφ0 = 10，aψ
0 = a0 = a0 = b0 = b0 = 1 とした．また，—3—

デルを用いた音源分離法を提案した．ネスト型基底・音源混合

表 1 音楽による評価

SDR

SIRSARIVA

0.3 dB

4.9 dB

5.7 dB

MNMF

1.0 dB

6.2 dB

6.7 dB

Sparse-LDA

0.7 dB

7.4 dB

4.1 dB

NMF-LDA1.2 dB

9.6 dB

3.5 dB

LDA-LDA

1.1 dB

8.8 dB

3.7 dB

表 2 音声による評価

SDR

SIR

SAR

IVA

3.4 dB

7.5 dB

7.1 dB

MNMF

4.8 dB

10.0 dB

7.7 dB

Sparse-LDA

5.5 dB

15.1 dB

6.3 dB

NMF-LDA

4.9 dB

15.0 dB

5.7 dB

LDA-LDA

5.9 dB

16.8 dB

6.4 dB

表 3 音楽 + 音声による評価

SDR

SIR

SAR

IVA

0.1 dB

5.3 dB

5.3 dB

MNMF

1.8 dB

8.6 dB

6.1 dB

Sparse-LDA

2.4 dB

11.5 dB

4.5 dB

NMF-LDA

1.3 dB

10.6 dB

3.9 dB

LDA-LDA2.6 dB

13.2 dB

4.2 dB−1
G0f d = (af d aHとし，af d には無響室で 5◦ 間隔
f d + 0.01 × I)で録音したインパルス応答を用いた．ギブスサンプリングの試行回数は 50 回とし，はじめの 30 回は burn-in として棄却した．
評価尺度として，Signal-to-distortion ratio (SDR)，Signal-to-inference ratio (SIR)，Signal-to-artificial ratio (SAR) を用いた．SDR は総合的な分離性能，SIR は目的音以外の音の除去性能，SAR は分離音の歪みの少なさを表す尺度である．

3. 2 実 験 結 果表 1，2，3 に実験結果を示す．それぞれの条件において最も数値が大きくなったものを太字で示した．SDR と SIR は音楽
においては NMF-LDA が最も大きく，音声と音楽+音声においては LDA-LDA が最も大きくなった．このことから，より低ランク性の強い音源である音楽に対しては低ランク近似を用いる

NMF-LDA の方が分離精度は高く，比較的低ランク性の弱い音源である音声などが含まれる混合音に対しては LDA-LDA の方が分離精度が高いことがわかった．また，Sparse-LDA と比較すると SDR と SIR はすべての条件において提案法の方が大きくなった．このことから提案法の音源モデルの有効性が確認できた．ただし，SAR は全ての条件において MNMF が最も大きくなった．また，Sparse-LDA，NMF-LDA，LDA-LDA では SAR にはそれほど大きな差はなかった．これは，LDA を用いた分離法では完全に排他的な割り当てが行われるのに対し，MNMF ではソフトな割り当てを行えることから自然なモデル化が可能となり分離音の歪みが小さくなったのではないかと推測される．したがって，SAR の改善を行うためには空間モデルの改善が必要であると考えられる．

4. お わ り に
本稿では，空間モデルを音源の混合モデルで記述し，音源モモデルを用いることで，位相情報だけでなく，音源のパワースペクトログラムの構造も考慮した音源分離を可能とした．実験の結果，音声などの比較的低ランクでない音源が含まれる混合音に対しては提案法により従来法よりも SDR や SIR が向上することを確認した．今後は，提案法のオンライン化や音源数・基底数の推定を行う．

付録

複素ウィシャート分布と一般化逆ガウス分布の確率密度関数は以下のとおりである．

WC (G|ν, G0 ) =GIG(y|γ, ρ, τ ) =|G|ν−M exp(−tr(G(G0 )−1 ))|G0 |ν π M (M −1)/2"M −1m=0

Γ(ν − m)exp{(γ − 1) log y − ρy − τ /y}√2τ γ/2 Kγ (2 ρτ )

(A·1)
(A·2)

ただし Kγ は第 2 種変形ベッセル関数である．

謝辞 本研究の一部は、JSPS 科研費 24220006, 15K12063 の支援を受けた．

文

献[1] A. Hyvärinen et al., Independent component analysis, John
Wiley & Sons, 2004.
[2] I. Lee et al., “Fast fixed-point independent vector analysis
algorithms for convolutive blind source separation,” Signal
Processing, vol.87, no.8, pp.1859–1871, 2007.
[3] A. Hyvarinen, “Fast and robust fixed-point algorithms forindependent component analysis,” IEEE Transactions on
Neural Networks, vol.10, no.3, pp.626–634, 1999.
[4] I. Jafari et al., “On the use of the watson mixture model
for clustering-based under-determined blind source separation,” INTERSPEECH, pp.988–992, 2014.
[5] N. Ito et al., “Permutation-free convolutive blind source
separation via full-band clustering based on frequencyindependent source presence priors,” ICASSP, pp.3238–3242, 2013.
[6] L. Drude et al., “Blind speech separation based on complex
spherical k-mode clustering,” ICASSP, pp.141–145, 2016.
[7] H. Sawada et al., “Underdetermined convolutive blindsource separation via frequency bin-wise clustering and permutation alignment,” TASLP, vol.19, no.3, pp.516–527,
2011.
[8] T. Otsuka et al., “Bayesian nonparametrics for microphonearray processing,” TASLP, pp.493–504, 2014.
[9] K. Itakura et al., “A unified bayesian model of timefrequency clustering and low-rank approximation for multichannel source separation,” EUSIPCO, to appear, 2016.
[10] P. Smaragdis et al., “Non-negative matrix factorization for
polyphonic music transcription,” WASPAA, pp.177–180,
2003.
[11] N. Ono, “Stable and fast update rules for independent vector analysis based on auxiliary function technique,” WASPAA, pp.189–192, 2011.
[12] H. Sawada et al., “Multichannel extensions of non-negative
matrix factorization with complex-valued data,” TASLP,
vol.21, no.5, pp.971–982, 2013.
[13] S. Araki et al., “The 2011 signal separation evaluation campaign (sisec2011):-audio source separation,” Latent Variable
Analysis and Signal Separation, pp.414–422, Springer, 2012.
[14] K. Itou et al., “The design of the newspaper-based japaneselarge vocabulary continuous speech recognition corpus,” ICSLP, pp.3261–3263, 1998.

デルを基底の混合モデルで記述したネスト型基底・音源混合モ—4—