130

日本音響学会誌 75 巻 3 号（2019），pp. 130–138

小特集—位相情報を考慮した音声音響信号処理—

複素生成モデルに基づく
非負値行列因子分解と音源分離への応用*
北 村 大 地（香川高等専門学校）∗∗
43.60.Cg, Uv

も十分小さな値に設定されるため，式 (1) は観測

1. は じ め に

行列 X を限られた数（K 個）の基底で表現する

非負値行列因子分解（nonnegative matrix fac-torization: NMF）[1, 2] は，非負要素から成る 2

低ランク近似と解釈でき，少数の有意な非負特徴量が基底 wk 及び係数 hk として得られる。

NMF の変数は次式の最適化問題で推定される。

次元データ，すなわち非負行列から，少数の有意な非負特徴量を抽出する数理アルゴリズムである。

min D(X|W H) s.t. wik , hkj ≥ 0 ∀i, j, kW ,H

非負行列を観測データとするパターン抽出に効果(4)的であり，その解釈の容易さから様々な分野で活用されてきたが，とりわけ音響分野において高度に発展してきた歴史がある。
今，解析対象となる非負観測行列を X ∈ RI×J
≥0

とおく。ここで，I と J はそれぞれ X の行数と列数である。このとき，NMF は次式で示される
近似分解を行う。

X ≈ X̂ = W H

wk h T
=
k

(1)
(2)

k

ここで，W = (w1 · · · wK ) ∈ RI×Kは基底行
≥0列と呼ばれ，X 中の有意な非負の特徴量ベクトル（基底）wk ∈ RI×1
≥0 を列に持つ行列である。また，

H = (h1 · · · hK )T

∈ RK×Jは係数行列やアク≥0

ティベーション行列と呼ばれ，基底 wk の係数ベ

クトル hk ∈ RJ×1
≥0 を行に持つ行列である。式 (1)を要素ごとに記述すると次式となる。

xij ≈ x̂ij =

wik hkj(3)

k

ここで，xij ，x̂ij ，wik ，及び hkj はそれぞれ X ，X̂ ，W ，及び H の非負要素である。添え字の i,j, 及び k はそれぞれ X の行，列，及び基底のインデクスを表す。一般に基底数 K は I や J より∗

Nonnegative matrix factorization based on complex
generative model and its application to audio source
separation.
∗∗
Daichi Kitamura (National Institute of Technology,
Kagawa Collage, Takamatsu, 761–8058)

ここで D(·|·) は二つの行列の類似度を測る任意の関数であり，二乗 Euclid 距離 [2]，一般化 Kulback–Leibler（KL）ダイバージェンス [2]，Itakura–Saito（IS）ダイバージェンス [3] 等がよく用いられる。
従って式 (4) は，非負観測行列 X を良く近似する
ような非負低ランクモデル行列 W H を推定する問題となる。式 (4) の解は閉形式では与えられな
いため，W と H を何等かの方法で初期化（例えば [4] 等）したうえで，反復計算により解を求めるアルゴリズムが提案されている [5, 6]。
音響信号処理においては，信号を短時間フーリエ変換（short-time Fourier transform: STFT）して得られる複素スペクトログラム C ∈ CI×J（行列C の複素要素を cij と定義する）を非負化したもの
を観測行列 X とおくことが一般的である。このとき，I は周波数ビン数，J は時間フレーム数に対応する。特に，振幅スペクトログラム (X = |C|.1 )やパワースペクトログラム (X = |C|.2 ) に NMFを適用する例が多い1 。Fig. 1 に，パワースペクト
ログラムを NMF で分解した例を示す。観測行列には音高の異なる二つのピアノ音が含まれているが，基底数を K = 2 と設定して分解することで，各音の調波構造（スペクトル）が基底 w1 及び w2に，各音の時間強度が係数 h1 及び h2 にそれぞれ現れている。このように，NMF は音響信号中の1行列に対するドット付き指数と絶対値記号を組み合わせた演算子 | · |.p を，要素ごとの絶対値の p 乗を要素に持つ行列と定義する。すなわち，[|C|.p ]ij = |cij |p である。

131

複素生成モデルに基づく非負値行列因子分解と音源分離への応用

Fig. 1 NMF decomposition for audio signals.

頻出スペクトル及びその時間的な強度変化を基底及び係数として抽出できる教師無し学習である。

Fig. 2 Wrong mixture model assumed in NMF for audio signals.

NMF は音響信号から直感的な特徴量を簡便に抽出できるため，音源分離 [7–11]，自動採譜 [12–14]，
音響イベント検出 [15]，音超解像 [16] 等，様々なペクトルに関して鈍感ではあるが，インパルス信タスクに応用されてきた。本稿では，非負スペク

知覚に大いに寄与する場合もある。また，NMF

トログラム X ではなく，複素スペクトログラムで推定した非負スペクトログラム成分 wk hT
k を逆

C を NMF でモデル化する手法と，複素生成モデ

号と白色雑音のように，位相スペクトルが人間のSTFT で時間領域に戻す際にも，wk hT
k に何等かルによる統計的解釈及びその拡張理論について紹の位相を付与する必要があり，例えば観測複素ス

介し，音源分離問題への応用例を示す。

ペクトログラム C の位相をそのまま付与する方

2. NMF に基づく音響モデリングの問題点

法，C に対して Wiener フィルタを適用する方法，
位相復元法 [17–19] を適用する方法等がある。

NMF において非負行列 X は，式 (2) に示されるように，ランク 1 の非負行列 wk hT
k の和として手法が提案されている。本稿では，特に，NMF

近似表現される。この分解モデルは，本来非負での分解成分 wik hkj を複素数に拡張した複素 NMF

あるような観測データ（白黒画像の輝度値や顧客（complex NMF: CNMF）[20–22] 及び複素 Gaussの商品購買数等）に対しては正当性があるが，音分布生成モデルに基づく NMF [3] とその理論拡響信号の分解に対しては次のような問題が生じる。

張手法 [23–25] を取り上げて紹介する。

Fig. 2 に示すように，本来音響信号の混合は時間領域の信号の和であり，時間周波数領域では各この問題を解決するために，これまでに様々な

3. 位相スペクトルを考慮した CNMF

音響信号の複素スペクトログラムの和に対応する。

従来の NMF では，位相スペクトルの影響を無一方で，音響信号に NMF を適用することは，「複視して少数の非負スペクトログラム成分 wk hT
k を

素スペクトログラムの和」を「非負のランク 1 スペ

推定しており，現実の物理現象との乖離がある。そクトログラム成分 wk hT
k の和」で近似表現するこ

こで，少数の音響成分を非負値から複素数値に拡とに相当する。一般に二つの複素数値 c1 及び c2張した CNMF が提案されている [20]。CNMF でに対して非負値の加法性 |c1 + c2 | = |c1 | + |c2 |は，複素スペクトル cij を次式のように分解する。

p

p

p

は p = 0 で成り立たないため，NMF で仮定される非負スペクトログラム成分の加法性は物理的にcij ≈ ĉij ==音響信号が混合する際の位相ずれによる打ち消し

れは NMF を用いた音響信号処理特有の問題であり，NMF が複素行列ではなく非負行列を対象としていることが原因である。人間の聴覚は位相ス(5)wik hkj ejφij,k(6)k合いが取り扱われないため，NMF によるモデル
まま基底及び係数を推定していることになる。こĉij,kk

間違ったモデルである。より具体的には，二つの化では，位相に起因する物理現象が考慮されない√−1 である。また，ĉij,k ∈ C は複素スペクトル成分であり，|ĉij,k | = wik hkj 及びarg(ĉij,k ) = φij,k である。すなわち，ランク 1 の振幅スペクトログラム成分 wk hT
k と位相スペクト
I×J
ログラム成分 Φk ∈ R[0,2π]（要素は φij,k ）から成
ここで，j =

132

日本音響学会誌 75 巻 3 号（2019）

る複素スペクトログラム成分 Ĉk ∈ CI×J（要素はĉij,k ）を定義し，観測複素スペクトログラム C をC ≈ Ĉ =(7)Ĉkkと近似するモデルである。

CNMF で推定すべき変数は wik ，hkj ，及び φij,kであり，次式の生成モデルに基づく最尤推定として定式化される。

Fig. 3 Isotropic complex Gaussian distribution.

cij = ĉij + εij(8)εij ∼ NC (0, σ 2 )(9)数に加えて，新たに位相を推定しなければならな1|c − μ|2
NC (c; μ, σ ) = 2 exp −
(10)
πσ
σ2

いため，最適化の困難性に起因する極端な初期値

ここで，NC (c; μ, σ 2 ) は平均 μ，分散 σ 2 の点対称

4. 複素生成モデルに基づく NMF（等方性）複素 Gauss 分布であり，μ = 0 の場合は

CNMF は，音響信号の混合現象を少数の複素ス

2Fig. 3 に示すような複素平面上の原点対称分布である。これは，式 (10) が μ = 0 において複素確
率変数のノルム |c| にのみ依存する確率密度関数であることを示している。CNMF は式 (7) のように少数の複素スペクトログラム成分 Ĉk の和で観測複素スペクトログラム C を近似するが，その近似誤差 εij が時間周波数ごとに独立な Fig. 3 の分
布に従うと仮定すると，変数 wik ，hkj ，及び φij,kの最尤推定として次の最小化問題を得る。


2

 


min
ĉij,k cij −

wik ,hkj ,φij,ki,jks.t. wik , hkj ≥ 0 ∀i, j, k,
wik = 1 ∀k

依存の問題が指摘されている [26]。

ペクトログラム成分の和でモデル化することで，2 章に述べた問題を解決している。一方，一部の類似度関数に基づく NMF の最適化問題を，複素スペクトログラムの生成モデルに基づく非負パラメータの最尤推定問題として解釈できることが明らかにされており，非負スペクトログラム成分の加法性を統計的な意味で正当化できる方法が発見されている。本章では，NMF の複素生成モデルに基づく解釈について解説し，その理論的な拡張について紹介する。

4.1 IS ダイバージェンスに基づく NMF の統
計的解釈

IS ダイバージェンスは次式で定義される。

i

D(c|σ) =これは二乗 Euclid 距離に基づく NMF の複素数版と解釈でき，補助関数法を用いることで通常のNMF と同様に最適化が可能である [20]。更に近年では，CNMF の距離規範を一般化 KL ダイバージェンスや β ダイバージェンス（二乗 Euclid 距
離，一般化 KL ダイバージェンス，及び IS ダイ
バージェンスを含むより一般的な類似度関数）へと拡張する手法が提案されている [21, 22]。
以上のように，CNMF は振幅スペクトログラム

|c|2
|c|2
−
log
−1
σ2
σ2

(11)

式 (4) における類似度関数 D を式 (11) で定義し，

|c|2 = xij 及び σ 2 = k wik hkj と置いた NMF（ISNMF）の最適化問題は，変数 wik 及び hkj に無関係な項を省略して記述すると次式のようになる。

min

W ,H


i,j




x
 ij
+ log
wik hkj
k wik hkjks.t. wik , hkj ≥ 0 ∀i, j, k

(12)の低ランク構造を仮定しつつ複素スペクトログラ

この最小化問題において，非負観測データの定義

ム成分の和で観測信号を近似する物理的に正しいを xij = |cij |2 （パワースペクトル）とした場合，モデルとなっており，NMF の直接的な複素数拡

以下に示すような統計的解釈が可能となる。

張と解釈される。しかしながら，通常の NMF 変

今，観測された複素スペクトログラム C の各時

133

複素生成モデルに基づく非負値行列因子分解と音源分離への応用

間周波数における成分 cij が cij =k cij,k として少数の複素スペクトル成分 cij,k の和で構成されていると仮定する。また，各複素スペクトル成分
は，次式のように Fig. 3 の原点対称複素 Gauss 分布から生成されたものと仮定する。
2cij,k ∼ NC (0, σij,k
)

(13)

2
ここで，分散 σij,k> 0 は周波数 i と時間 j の両方に依存して変化する非負パラメータである。このとき，複素 Gauss 分布の安定性（あるいは再生性）によ
2り，cij,k の和である cij もまた分散 σij=

2
k σij,kの複素 Gauss 分布から生成される。
2cij ∼ NC (0, σij) 2= NC (0, k σij,k )Fig. 4 Local Gaussian model assumed in ISNMF.

(14)

式 (14) は観測スペクトルの生成モデルとして，時間周波数ごとに分散の異なる原点対称複素 Gauss2
分布 p(cij ) = NC (0, σij) を仮定している。この生成モデルを Fig. 4 に示す。パワーの大きな時間周

波数グリッドでは大きな分散値を持つため振幅の大きい複素数値 cij を生成し易く，逆にパワーの
小さな時間周波数グリッドではほとんど零付近の複素数値 cij しか生成しないような分布となっている。また，いずれの時間周波数グリッドも原点対称な分布であることから，位相 arg(cij ) に関し

この生成モデルは多チャネル信号（多変量分布）にも拡張され，local Gaussian model（LGM）[27]と呼ばれている。
次に，式 (16) の複素生成モデルに基づいて，分2
散 σij,kの最尤推定問題を考える。観測データ Cの尤度関数 L は次式で表される。

L=

NC (0,i,j

2k σij,k )

|c |2
 2 exp −  ij 2=
π k σij,kk σij,ki,j1ては無情報（一様）な分布である。各時間周波数2
グリッドにおける分散 σij
は複素スペクトル成分

cij のパワーの期待値であるため

k

=

(17)
更に，負対数尤度関数は− log L|cij |22
 2 +log
σij,k+log π
=
σ
k
ij,k
i,j

2
σij
= E[|cij |2 ]E[|cij,k |2 ]
=





k

2
σij,k

(15)

k

が成立する。ここで，E[·] は観測データに対する期待値を表す。式 (14) の生成モデルが周波数 i と
時間 j に対して統計的に独立であると仮定すると，
観測複素スペクトログラム C の（時間周波数全体
における）生成モデルは次式となる。

(18)となる。従って，分散の最尤値を得るには，負対2を求めればよい。
数尤度 (18) を最小化する σij,kここで，式 (18) と ISNMF の最小化問題 (12) と2
比較すると，非負観測データ xij 及び分散 σij,kを2= wik hkj と置く
それぞれ xij = |cij |2 及び σij,kとき，両式は定数項の違いを除いて一致する。

C ∼ p(c11 , c12 , · · · , cIJ )以上より，ISNMF をパワースペクトログラム

p(cij )

X = |C|.2 に対して適用することは，式 (16) の2複素生成モデル（LGM）を仮定した分散 σij,k
の

=
i,j

NC (0,=i,j

2k σij,k )(16)

最尤推定問題と等価であることが分かる。また，
この複素生成モデルでは，非負成分 wik hkj の和134

日本音響学会誌 75 巻 3 号（2019）2
が，式 (15) に示す分散 σij,k= E[|cij,k |2 ] の和にトル成分 cij,k が式 (20) から生成されるとき，時対応していることが分かる。この事実は，複素ス

間周波数ごとに定義される尺度母数は複素スペク

ペクトル成分の混合（cij =k cij,k ）を非負パ
 22
ラメータである分散値の和（σij =k σij,k ）でトル成分の振幅の期待値 γij,k = E[|cij,k |] に対応表現できることを意味し，音響信号に NMF を適

測複素スペクトル cij =用する際の「複素数を非負化して少数の成分に分また解」という一見場当たり的にみえる手続きを期待値の意味で正当化できる。従って，パワースペク

する。複素 Cauchy 分布は安定性を持つため，観cij ∼ CC (0,k cij,k の生成モデルも

k γij,k )(21)

トログラムを ISNMF で少数の非負スペクトログとなり，振幅スペクトル成分 γij,k = wik hkj の加ラム成分 wk hT
k の和に近似分解する場合，Fig. 2

法性が期待値の意味で正当化される。

に示すような時間領域での信号の混合現象を，時この生成モデルに基づく NMF（Cauchy NMF）

ラムの加法性を期待値の意味で正当化するために[23] が提案され，パワースペクトログラム |C|.2 及
び振幅スペクトログラム |C|.1 の加法性を正当化できる NMF がそれぞれ ISNMF [3] 及び Cauchy
NMF [23] として与えられた。また，複素 Gauss
分布と複素 Cauchy 分布の両方を含む一般化である複素 Student’s t 分布に基づく NMF（tNMF）[24] も提案されており，自由度パラメータ ν を変は，安定性を持つ分布，すなわち安定分布を生成モ化させることで両分布の中間的な生成モデルを仮デルとして仮定する必要がある。安定性とは，同定することが可能となっている。ただし，安定性

じ分布族から独立に生成された二つの確率変数 v1

は Gauss 分布（ν → ∞）と Cauchy 分布（ν = 1）と v2 に対して，その線形結合 av1 + bv2 と，別のの場合にのみ成立するため，非負スペクトログラ

確率変数 dv + e（ここで a > 0，b > 0，d > 0，

ムの加法性は両分布でのみ正当化される。tNMF

及び e は定数）がいずれも同一の分布族から生成も LGM と同様に多チャネル音源分離へと応用さされる性質を指す [32]。このような安定分布を生れ [33, 34]，より高精度な分離を達成している。な成モデルとして仮定した場合，確率変数の和をそお，各 NMF の最適化アルゴリズム（反復更新式）の分布のパラメータの和としてモデル化することについては，それぞれの文献を参照されたい。

が可能となる。事実，原点対称複素 Gauss 分布は5. 複素一般化 Gauss 分布に基づく NMF とスパース雑音除去への応用

間周波数領域の非負値の和でモデル化することの正当性が保障される。この LGM に基づく定式化
は，多チャネル信号の音源分離へと応用されている [28–31]。

4.2 LGM の理論拡張に基づく NMF
ISNMF のように音響信号の非負スペクトログ

安定分布の特殊形であり，次式のように複素数値の和を分散の和で表現できる。

本章では，LGM の新たな一般化として，複素c1 ∼ NC (0, σ12 ), c2 ∼ NC (0, σ22 )一般化 Gauss 分布（generalized Gaussian distri-c1 + c2 ∼ NC (0, σ12 + σ22 )(19)

同様に，複素数値の和を 1 次の期待値（尺度母数）γ の和で表現できる分布として，次式がある [32]。
2

CC (c; 0, γ) =
2π

−1/2

|c|2 +

γ2
2−1/2 γ32(20)

ここで，CC (c; 0, γ) は最頻値 0 かつ尺度母数 γ の原点対称複素 Cauchy 分布である2 。複素スペク

文献 [23] では，式 (20) の尺度母数 γ の係数 2−1/2 はない。これは文献 [23] が内部で参照している文献 [32]（pp. 80–81）の誤植に由来する。係数 2−1/2 を付与する
ことで，|c| = γ のときに尤度が最大となる。
2

bution: GGD）[25] に基づく NMF（GGDNMF）
について解説し，その性質を活用したスパース雑音除去の例を紹介する。

5.1 GGDNMF の生成モデル
GGDNMF では，ISNMF で仮定されていた
LGM を複素 GGD へと拡張している。すなわ
ち，平均値 0，形状母数 ρ > 0，時間周波数ごとに異なる尺度母数 σij を持つ原点対称複素
GGD GC (cij ; 0, ρ, σij ) を，観測複素スペクトロ
グラム C の生成モデルとして仮定する。
C∼

GC (cij ; 0, ρ, σij )
i,j

135

複素生成モデルに基づく非負値行列因子分解と音源分離への応用

Fig. 5 Isotropic complex GGD with σij = 1.5.

ρ 

2 |cij |
=
exp −
2
1− ρ
2 Γ(2/ρ)
ρ σij
2
πσ
i,j
2

ρ1− ρ



ij

(22)
p
=
σij



wik hkj(23)

k

ここで，Γ(·) はガンマ関数であり，p は低ランク近似の領域に対応するドメインパラメータである。
すなわち，p = 1 の場合は振幅スペクトログラム|C|.1 ，p = 2 の場合はパワースペクトログラム|C|.2 を W と H で低ランク近似することに対応する。複素 GGD GC (cij ; 0, ρ, σij ) は Fig. 5 に示
すとおり，ρ = 2 及び ρ = 1 のとき，それぞれ原点対称な複素 Gauss 分布及び複素 Laplace 分布に一致する。また，ρ > 2 では劣 Gauss（platykurtic），
ρ < 2 では優 Gauss（leptokurtic）な分布となる。
5.2 複素 GGD から導かれるダイバージェンス
式 (22) の生成モデルと NMF における類似度関数の関係を明らかにするために，対数尤度差（deviance）から複素 GGD に基づくダイバージェンスを導出する。複素 GGD の対数尤度は次式となる。
log L=log GC (c; 0, ρ, σ)
=log

ρ
2

2
1− ρ

21−ρ πΓ(2/ρ)

−2 log σ−

2 |c|
ρ σ

ρ

のダイバージェンスは，IS ダイバージェンスの式(11) と比較すると，形状母数 ρ に関する一般化となっていることが分かる。また，より一般化された類似度関数である α-β ダイバージェンス [35] の一部となっていることも明らかにされている [25]。

5.3 最適化アルゴリズム
GGDNMF の最小化すべき関数は次式である。
|cij |ρD(cij |σij ) =ρ( k wik hkj ) p
i,j
i,j


ρ
+ log
wik hkj (25)pkただし，定数項を除いて示している。式 (25) を最小化する wik 及び hkj は補助関数法により導出でき，次式の更新式が得られる [25]。

⎡

⎤ ρ+p
zij
j (  w  h  )2 hkjp

wik ← wik ⎣ ⎡zij = |cij |

ρ
p

ik

kj

1

hkj
j
k wik hkji (hkj ← hkj ⎣ 

k



i

1− ρ
σij p

zij
k wik hkj



)

2

wik

1
wik
k wik hkj

⎦

p
⎤ ρ+p

⎦

p

5.4 スパース雑音除去への応用
複素 GGD における ρ < 2 のようにヘビーテイ
ルな分布に基づく NMF の応用例として，信号に∂log L/∂σ=0 より，σ の最尤値は σML =|c| である。対数尤度差 D= log L(σML )− log L(σ)≥0 は
 ρ
2 |c|
2
D(c σ)=−2 log |c| − +2 log σ+ρρ σ ρ 
 ρ
2
|c|
|c|
=
−log
−1 (24)
ρσ
σ

クトル減算法等の非線形な信号処理を適用した際となる。式 (24) は非負かつ σ = |c| のとき唯一 0に生じる人工的な成分であり，ミュージカルノイとなるため，ダイバージェンスの公理を満たす。こズと呼ばれている。この観測信号（Fig. 6(b)）に対足しあわされたスパース雑音（時間周波数領域でスパースに分布する成分）の除去が挙げられる。
実験に用いた原信号と雑音を付与した観測信号をFig. 6 に示す。このようなスパース雑音は，スペ136

日本音響学会誌 75 巻 3 号（2019）

Fig. 6 Power spectrograms of (a) original and (b) observed noisy signals.

Fig. 8 Average SDR improvements of (a) tNMF and(b) GGDNMF for various ν, ρ, and p.

に対して行い，tNMF と GGDNMF に関して平
Fig. 7 Examples of estimated power spectrogram by(a) ISNMF (SDR: −13.52 dB), (b) Cauchy NMF(SDR: 3.77 dB), (c) tNMF (SDR: 7.26 dB), and (d)GGDNMF (SDR: 7.38 dB).

均 SDR 改善量を示したグラフである。これをみると，仮定すべき最適な分布の裾の重さ（各線における最高性能の横軸位置）が観測信号のドメイン pに応じて変化していく様子が分かる。これは，p の増加に伴って値の大きなスパース雑音成分が強調

して，パワースペクトログラムに対する ISNMF，され，外れ値としての影響が強くなるためである。

振幅スペクトログラムに対する Cauchy NMF，

SDR 改善量の観点では，tNMF と GGDNMF は.p

任意の信号ドメイン（|C| ）に対する tNMF 及同程度の雑音除去性能を示している。

び GGDNMF を適用した結果を Fig. 7 に示す。

6. お わ り に

ここで，NMF の基底数はすべて K = 30 としている。また，tNMF は ν = 2，p = 0.5 とし，

本稿では，音響信号に NMF を適用する際の問GGDNMF は ρ = p = 0.1 に設定している。評価
値には source-to-distortion ratio（SDR）[36] の改善量を用いた。Fig. 7 を見ると，ISNMF 以外

題点として，位相スペクトログラムの影響が無視さの複素数拡張である CNMF を紹介した。また，パの手法はスパース雑音が比較的高精度に除去されワースペクトログラムを ISNMF で分解することていることが確認できる。これは，LGM をよりが LGM に基づく最尤推定と等価である解釈を示

ヘビーテイルな分布へと拡張した Cauchy NMF，し，NMF で仮定される非負スペクトログラムの分

tNMF，及び GGDNMF が，スパース雑音を外れ解の正当化が可能であることも紹介した。最後に，値として扱いながら低ランク近似した結果である。
すなわち，ヘビーテイルな分布を仮定する NMF

LGM の理論拡張である Cauchy NMF，tNMF，
及び GGDNMF について，スパース雑音除去へのでは，低ランク性に寄与しないパワーの大きな成適用例を示した。

れる問題について説明し，その解決法として NMF

分を無視しながら基底 wk と係数 hk を推定できるため，Fig. 6(b) のような信号からも頑健に低ランク構造を抽出できる。

Fig. 8 は，Fig. 7 と同様の実験を 6 種類の信号

謝

辞

本研究の一部は SECOM 科学技術支援財団，ヤマハ株式会社，及び JSPS 科研費 17H06572 の助成を受けた。

複素生成モデルに基づく非負値行列因子分解と音源分離への応用
文

献[ 1 ] D. D. Lee and H. S. Seung, “Learning the partsof objects by non-negative matrix factorization,” Nature, 401(6755), 788–791 (1999).
[ 2 ] D. D. Lee and H. S. Seung, “Algorithms for nonnegative matrix factorization” Proc. NIPS, pp. 556–
562 (2000).
[ 3 ] C. Févotte, N. Bertin and J.-L. Durrieu, “Nonnegative matrix factorization with the Itakura-Saito
divergence. With application to music analysis,” Neural Comput., 21, 793–830 (2009).
[ 4 ] D. Kitamura and N. Ono, “Eﬃcient initialization for nonnegative matrix factorization based on
nonnegative independent component analysis,” Proc.
IWAENC (2016).
[ 5 ] M. Nakano, H. Kameoka, J. Le Roux, Y. Kitano,
N. Ono and S. Sagayama, “Convergence-guaranteed
multiplicative algorithms for nonnegative matrix factorization with β-divergence,” Proc. MLSP, pp. 283–
288 (2010).
[ 6 ] C. Févotte and J. Idier, “Algorithms for nonnegative matrix factorization with the β-divergence,”
Neural Comput., 23, 2421–2456 (2011).
[ 7 ] T. Virtanen, “Monaural sound source separation
by nonnegative matrix factorization with temporalcontinuity and sparseness criteria,” IEEE Trans. Audio Speech Lang. Process., 15, 1066–1074 (2007).
[ 8 ] P. Smaragdis, B. Raj and M. Shashanka, “Supervised and semi-supervised separation of sounds
from single-channel mixtures,” Proc. ICA, pp. 414–421 (2007).
[ 9 ] H. Kameoka, M. Nakano, K. Ochiai, Y. Imoto,
K. Kashino and S. Sagayama, “Constrained and regularized variants of non-negative matrix factorization incorporating music-speciﬁc constraints,” Proc.
ICASSP, pp. 5365–5368 (2012).
[10] D. Kitamura, H. Saruwatari, K. Yagi, K. Shikano,Y. Takahashi and K. Kondo, “Music signal separation based on supervised nonnegative matrix factorization with orthogonality and maximum-divergencepenalties,” IEICE Trans. Fundam., E97-A, 1113–
1118 (2014).
[11] D. Kitamura, H. Saruwatari, H. Kameoka, Y.
Takahashi, K. Kondo and S. Nakamura, “Multichannel signal separation combining directional clusteringand nonnegative matrix factorization with spectrogram restoration,” IEEE/ACM Trans. Audio SpeechLang. Process., 23, 654–669 (2015).
[12] P. Smaragdis and J. C. Brown, “Non-negative
matrix factorization for polyphonic music transcription,” Proc. WASPAA, pp. 177–180 (2003).
[13] S. A. Raczyński, N. Ono and S. Sagayama, “Multipitch analysis with harmonic nonnegative matrix
approximation,” Proc. ISMIR, pp. 381–386 (2007).
[14] N. Bertin, R. Badeau and E. Vincent, “Enforcing harmonicity and smoothness in Bayesian nonnegative matrix factorization applied to polyphonicmusic transcription,” IEEE Trans. Audio SpeechLang. Process., 18, 538–549 (2010),[15] T. Heittola, A. Mesaros, T. Virtanen and A.
Eronen, “Sound event detection in multisource environments using source separation” Proc. CHiME,
pp. 36–40 (2011).
[16] D. Bansal, B. Raj and P. Smaragdis, “Bandwidth
expansion of narrowband speech using non-negative
matrix factorization,” Proc. Interspeech, pp. 1505–

137

1508 (2005).
[17] D. W. Griﬃn and J. S. Lim, “Signal estimationfrom modiﬁed short-time Fourier transform,” IEEE
Trans. Acoust. Speech Signal Process., 32, 236–243(1984).
[18] J. Le Roux, H. Kameoka, N. Ono and S.
Sagayama, “Fast signal reconstruction from magnitude STFT spectrogram based on spectrogram consistency,” Proc. DAFx, pp. 397–403 (2010).
[19] 矢田部浩平, 升山義紀, 草野 翼, 及川靖広, “位相変換による複素スペクトログラムの表現,” 音響学会誌, 75,
147–155 (2019).
[20] H. Kameoka, N. Ono, K. Kashino and S.
Sagayama, “Complex NMF: A new sparse representation for acoustic signals,” Proc. ICASSP, pp. 3437–
3440 (2009).
[21] H. Kameoka, H. Kagami and M. Yukawa, “Complex NMF with the generalized Kullback-Leibler divergence,” Proc. ICASSP, pp. 56–60 (2017).
[22] P. Magron and T. Virtanen, “Towards complex nonnegative matrix factorization with the betadivergence,” Proc. IWAENC, pp. 156–160 (2018).
[23] A. Liutkus, D. Fitzgerald and R. Badeau,“Cauchy nonnegative matrix factorization,” Proc.
WASPAA (2015).
[24] K. Yoshii, K. Itoyama and M. Goto, “Student’s t nonnegative matrix factorization and positive
semideﬁnite tensor factorization for single-channel
audio source separation,” Proc. ICASSP, pp. 51–55(2016).
[25] 北村大地, 高宗典玄, 最上伸一, 三井祥幹, 猿渡 洋,
高橋 祐, 近藤多伸, “ヘビーテイルな分布に基づく非負値行列因子分解を用いたスパース雑音除去,” 音講論集,
pp. 441–444 (2018).
[26] P. Magron, R. Badeau and B. David, “Phase
recovery in NMF for audio source separation: An
insightful benchmark,” Proc. ICASSP, pp. 81–85(2015).
[27] A. Ozerov, E. Vincent and F. Bimbot, “A general
ﬂexible framework for the handling of prior information in audio source separation,” IEEE Trans. Audio
Speech Lang. Process., 20, 1118–1133 (2012).
[28] A. Ozerov and C. Févotte, “Multichannel nonnegative matrix factorization in convolutive mixtures foraudio source separation,” IEEE Trans. Audio SpeechLang. Process., 18, 550–563 (2010).
[29] H. Sawada, H. Kameoka, S. Araki and N. Ueda,“Multichannel extensions of non-negative matrix factorization with complex-valued data,” IEEE Trans.
Audio Speech Lang. Process., 21, 971–982 (2013).
[30] D. Kitamura, N. Ono, H. Sawada, H. Kameoka
and H. Saruwatari, “Determined blind source separation unifying independent vector analysis and nonnegative matrix factorization,” IEEE/ACM Trans.
Audio Speech Lang. Process., 24, 1626–1641 (2016).
[31] D. Kitamura, N. Ono, H. Sawada, H. Kameoka
and H. Saruwatari, “Determined blind source separation with independent low-rank matrix analysis,”
Audio Source Separation, S. Makino, Ed. (Springer,
Cham, 2018), pp. 125–155.
[32] G. Samorodnitsky and M. S. Taqqu, Stable NonGaussian Random Processes: Stochastic Models with
Infinite Variance (Chapman & Hall/CRC Press,
Boca Raton, 1994).
[33] K. Kitamura, Y. Bando, K. Itoyama and K.
Yoshii, “Student’s t multichannel nonnegative matrix factorization for blind source separation,” Proc.

138
IWAENC (2016).
[34] S. Mogami, D. Kitamura, Y. Mitsui, N. Takamune, H. Saruwatari and N. Ono, “Independent lowrank matrix analysis based on complex Student’s tdistribution for blind audio source separation,” Proc.
MLSP (2017).
[35] A. Cichocki, S. Cruces and S. Amari, “Generalized alpha-beta divergences and their application torobust nonnegative matrix factorization,” Entropy,
13, 134–170 (2011).
[36] E. Vincent, R. Gribonval and C. Fevotte, “Performance measurement in blind audio source separation,” IEEE Trans. Audio Speech Lang. Process., 14,1462–1469 (2006).

日本音響学会誌 75 巻 3 号（2019）北村 大地
香川高等専門学校専攻科修了 (2012)。
奈良先端科学技術大学院大学博士前期課
程修了 (2014)。総合研究大学院大学博士後期課程修了 (2017)。博士（情報学）。現在，香川高等専門学校助教。主として統計的信号処理，アレイ信号処理，機械学
習，音源分離の研究に従事。日本音響学会
学生・若手フォーラム幹事会員（2014–現在）。日本音響学会粟屋潔学術奨励賞 (2015)，日本学術振興会育志賞 (2017)，IEEE SPS Japan Best Paper Award(2017)，電気通信普及財団テレコムシステム技術賞奨励賞(2018)，日本音響学会独創研究奨励賞板倉記念 (2018) 等受賞。