卒業研究論文
論 文 題 目

正則化非負値行列因子分解による非負低ランク行列補完

提出年月日
学

科

氏

名

令和 6年 2月 22日

電 気 情 報 工 学 科
和気 佑弥

印

指導教員(主査)

北村 大地 講師

印副査

柿元 健 准教授

印長

漆原 史朗 教授

印学科

香川高等専門学校

Nonnegative low-rank matrix completion using
regularized nonnegative matrix factorization
Yuya Wake
Department of Electrical and Computer Engineering
National Institute of Technology, Kagawa College
Abstract
Real-world data often include missing values, resulting from collection or subsequent processing, influenced byvarious factors such as human errors or technical limitations. As the dataset size increases, missing values become
more prevalent. Directly handling missing values in data analysis is often impractical, leading to reduced usable
data size. Imputation becomes crucial for effective utilization of the dataset. In particular, matrix completion
methods are essential and have widely been investigated. Nonnegative matrix factorization (NMF) is a matrix
decomposition algorithm that focuses on the negativity of real-world data. NMF is an unsupervised learning andshows high matrix completion accuracy resulting in a cost-reduced application compared to deep-learning-basedor other supervised approaches. The matrix decomposition property in NMF can be controlled by introducing
regularization. However, the suitable regularization type in matrix completion has not extensively been explored.
In this thesis, I investigate the performance of NMF-based matrix completion using two regularization types;
cosine similarity regularization and 𝐿1-norm regularization. The former emphasizes orthogonality between basis
vectors, while the latter induces sparsity of the coefficient matrix in NMF. Experiments are conducted using
artificial random matrices and audio spectrograms with missing values, and the performance of the matrix
completion is compared. Results indicate that the proposed methods achieve a significant accuracy improvement
for the matrix completion. Both regularizations provide stable estimation of missing values even when the numberof basis vectors in NMF and the oracle rank of the input matrix are mismatched. In the audio spectrogram
experiment, artificial distortions due to the missing values are mitigated, which is substantial enhancementcompared to simple NMF without regularization.
Keywords: nonnegative matrix factorization, matrix completion, missing value, regularization

（和訳）
実世界で取り扱われるデータには，収集段階，またはその後のデータ処理により，データ内に欠損値が発生する．欠損値の発生要因は人為的なミスや技術的な限界など多様であり，データサイズが大きくなるほどその数は増す．多くの場合，データ解析において欠損値をそのまま扱うことはできず，結果として使用可能なデータサイズが減少してしまう．このため，情報資源を有効活用する上で欠損値に対し本来あるべき適切な値を推定して置き換える補完操作，特に対象データが 2 次元データであれば行列補完が重要であり，様々な手法が提案されている．行列補完手法の一つとして，現実において取り扱うべきデータがもつ非負性に着目した非負値行列因子分解（NMF）と呼ばれる行列分解アルゴリズムがある．NMF は事前学習不要のアルゴリズムでありながら高精度な補完性能を持ち，コスト面で深層学習やその他教師あり学習を必要とする手法に勝る．NMF は正則化によって分解行列の特性が変化する
ことが知られており，様々な正則化が提案されている．一方で，これらの正則化の行列補完への有効性はあまり検討されていない．本研究では，基底行列間の直交化を目的としたコサイン類似度正則化及び係数行列要素のスパース化を目的とした𝐿1ノルム正則化を行列補完 NMF に適用することで，補完性能の変化を調査する．人工欠損乱数行列と音声スペクトログラムを用いた補完実験を行い，正則化を用いない NMF による補完手法と比較した．実験により，人工欠損乱数行列の補完において，正則化項を加えることで従来法よりコスト値が大きく改善することが示された．また，従来法では性能が安定しない
基底数の設定が実際のデータと乖離した条件において，いずれの正則化も安定した結果をもたらすことを確認した．音声スペクトログラムによる実験では，僅かな歪みの改善が見られ，正則化を用いない
NMF による補完手法に対し大きな改善を得た．

目次
第1章

緒言

1

1.1

本論文の背景 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

1

1.2

本論文の目的 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

3

1.3

本論文の構成 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

4

NMF を用いた行列補完

5

2.1

はじめに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5

2.2

NMF . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5

2.2.1

NMF の概要 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

5

2.2.2

NMF の定式化 . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

6

NMF の行列補完への応用 . . . . . . . . . . . . . . . . . . . . . . . . . . .

10

2.3.1

バイナリマスク行列を用いた NMF に基づく行列補完 . . . . . . . .

10

2.3.2

バイナリマスク行列を用いた NMF に基づく行列補完の課題 . . . . .

11

本章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

12

NMF の正則化項付与

13

3.1

はじめに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

13

3.2

CSRSNMF . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

13

3.3

L1 -RNMF . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

15

3.4

本章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

16

正則化付き NMF に基づく行列補完

17

4.1

はじめに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

17

4.2

動機 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

17

4.3

CSRNMF に基づく行列補完 . . . . . . . . . . . . . . . . . . . . . . . . . .

18

4.3.1

コサイン類似度正則化の NMF への適用 . . . . . . . . . . . . . . .

18

4.3.2

補助関数法に基づく上限関数の導出 . . . . . . . . . . . . . . . . . .

18

4.3.3

CSRNMF の更新式の導出 . . . . . . . . . . . . . . . . . . . . . . .

19

4.3.4

CSRNMF の行列補完への適用 . . . . . . . . . . . . . . . . . . . .

21

L1 -RNMF に基づく行列補完 . . . . . . . . . . . . . . . . . . . . . . . . .

22

L1 -RNMF に基づく行列補完の目的関数 . . . . . . . . . . . . . . .

22

第2章

2.3

2.4
第3章

第4章

4.4

4.4.1

4.4.2

補助関数法に基づく上限関数の導出 . . . . . . . . . . . . . . . . . .

23

4.4.3

L1 -RNMF に基づく行列補完の更新式の導出 . . . . . . . . . . . . .

23

4.5

正則化項による NMF の挙動変化 . . . . . . . . . . . . . . . . . . . . . . .

25

4.6

本章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

28

人工欠損行列の補完実験

29

5.1

はじめに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

29

5.2

実験方法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

29

5.3

実験結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

30

5.4

本章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

36

欠損音声スペクトログラムの補完実験

38

6.1

はじめに . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

38

6.2

実験背景 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

38

6.2.1

STFT . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

38

6.2.2

時間周波数マスク . . . . . . . . . . . . . . . . . . . . . . . . . . .

40

6.3

実験方法 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

41

6.4

実験結果 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

46

6.5

本章のまとめ . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

49

結言

50

第5章

第6章

第7章
謝辞

51

参考文献

51

付録 A

補助関数法で利用される不等式

54

Jensen の不等式 . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . . .

54

付録 B

人工欠損行列の補完実験結果

55

付録 C

欠損音声スペクトログラムの補完実験結果

65

A.1

1

第1章

緒言
1.1 本論文の背景
今日のデータ解析では情報資源の肥大化によりビッグデータを扱う機会が増加している．取り扱われるデータは音声データ，画像データ，顧客の購買データなど，その種類は多様であり，それぞれの収集背景を持つ．このようなデータには往々にして欠損値が発生する．
欠損値とは，データの一部が破損して本来の値がわからない状態を言い，データの収集背景などに左右されて様々な要因により生じる．例として，Fig. 1.1 に音声信号，画像データ及び顧客の購買データにおける欠損の生成過程を示す．音声信号の場合，入力される音声が録音機
器のダイナミックレンジを超えると，それらの音は最大値，または最小値で表されることにな
るため，本来の大きさが表現できず欠損となる．画像データの場合，撮影段階において被写体が明るすぎる，または暗すぎると白飛びや黒潰れが発生し，これらが発生した部分の本来の画素値を取得できず欠損となる．顧客の購買データの場合，会員情報収集時において記入が不十分であると，一部の情報が不明になることで欠損が発生する．その他のデータにおいても，使用する機器の特性や不備，人為的なミスなど，あらゆる場面で欠損が発生することがある．これらの欠損は，その後に続く処理の中で望ましくない影響を及ぼすことや，そもそも欠損を含むデータを活用することができなくなる等の問題を引き起こす．
以上のように，データにおける欠損値の発生は少なくなく，データサイズが大きくなるほどにその数は増加する．一方で，データ分析に用いる手法には欠損値を含むデータを扱えないものも多くあり，分析前に何らかの方法によって欠損値を削除する必要がある．欠損値の削除方
法としては，欠損値を含む要素を切り捨てる方法と欠損部分に相当する値を予測し挿入する方
法の 2 通りがあり，後者は欠損代入（missing value imputation），欠損補間（missing valueinterpolation），欠損補完（missing value completion）等と呼ばれる．欠損データの切り捨てを行う場合，当然データ量が削減されるため，欠損の割合によっては解析結果やデータ自体が使い物にならない場合が生じる．このように，欠損値の対処において元データの性質を損なわないことは最も重要であり，補完が検討される場合も多い．
欠損データの補完を検討する場合，その手法には様々なものが提案されている．中でも，最も手軽な方法は単一代入法に分類される手法である．単一代入法は複数の欠損値部分に共通し

2

第 1 章 緒言Fig. 1.1. Examples of missing values: (a) audio recording, (b) image capturing, and (c) surveying
customer data.

て一つの値を代入する方法であり，データ全体の平均値を挿入する平均値補完や中央値を挿入する中央値補完などの種類がある．このような手法は手軽に適用できる一方，局所的な変化を持つ変数の一部が欠損した場合や特定の要因に従って欠損が発生する場合においては，一定値の挿入によってデータの特性を変質させる恐れがある．
そこで，より本来のデータに近い値で補完するために，データの本質的構造に基づく補完手法も多く提案されている [1]．近年，データの本質的構造に即した補完には，深層学習を代表とする学習データを用いた教師あり手法が多く提案されている [2]．これらの手法は，補完対象のデータを適切に学習することにより非常に高精度な補完を行うことができる．一方で，一定の精度を得るには莫大なコストがかかるため，小規模なデータ解析やサンプルデータを多く用意できない場合においては適用が困難となる．このため，コストを抑えつつデータの本質的構造を捉えた補完が可能な手法が求められている．

1.2 本論文の目的

3

1.2 本論文の目的
本論文では，非負値のみを含む 2 次元データの欠損補完問題を取り扱う．この問題は欠損を含む非負行列の行列補完問題として解釈できる．非負値のみを取り扱うことは，本論文が適用できる対象を限定することになるが，この制約が実際のデータにおいて大きな問題となることは比較的少ない．その理由として，実際に観測される様々なデータは非負値であることが多い点が挙げられる．
顧客データを例にあげると，5 人の顧客が 3 種類の商品の内，何をいくつ購入したかというデータは Fig. 1.2(a) のように 3 × 5 の非負値行列で表される．この時，−1 回購入したという負のデータは通常あり得ない．また，8 ビット白黒画像は Fig. 1.2(b) のように 0 から 255 の輝度を持つ非負値行列として表される．音響信号も例外ではなく，後述する短時間フーリエ変換（short-time Fourier transform: STFT）を用いて 1 次元データの時間信号を 2 次元データの時間周波数信号に変換することで複素行列が得られ，その各要素の振幅値を取った振幅行列を定義すれば，Fig. 1.2(c) のように非負値行列が現れる．この行列は多くの音響信号処理で加工の対象となる特徴量である．このように，実世界において非負値行列で表されるものは多く，このような非負値行列の構造抽出において有効なアルゴリズムとして，非負値行列因子分解（nonnegative matrix factorization: NMF） [3, 4] が知られている．

NMF は対象とするデータに非負制約を課し，低ランク性を仮定する行列分解アルゴリズム
であり，補完手法にも応用される．NMF に基づく行列補完は，学習不要で高精度な補完手法
であり，データの本質的構造を捉えながらコストを抑えることができる．上述の通り，NMFの非負制約は実世界において比較的問題になることはなく，また，NMF の仮定する低ランク性についても同様である．変数間に相関を持ったデータである場合，データサイズが巨大であればあるほど，相対的に少数の構造やパターンによって構成されることが多い．即ち，低ラン
ク性を持つデータとなる．例えば，1 万人の購買データがあるとき，その中から「食料品をよく買うグループ」や「娯楽品をよく買うグループ」といった特徴的なグループに分類することができる．このように，1 万人の顧客データをカテゴライズしていくことで，実際には 100 種
類程度の傾向に分類でき，この傾向の組み合わせによって顧客のほとんどの特性を説明できる
場合，行列においては低ランクであると言える．このようなパターンや傾向を見つけることができないデータとは，即ち，完全にランダムな要素を持つデータであり，実問題で扱うことは少ない．
以上より，実世界には非負かつ低ランクなデータは多く，NMF は幅広いデータに対して適用することができる．また，これらの制約により，NMF は学習データが不要でありながら高精度な補完を可能としている．一方で，高精度な補完を行うには，低ランク性が保証され，かつそのランクが判明しているという理想的な状態である必要がある．実際には，観測される
データはすでに欠損していることや観測中に付加されるノイズの影響によって，このような条件を満たすことは困難である．そして，本質的に低ランクなデータがその低ランク性を損なう

4

第 1 章 緒言Fig. 1.2. Nonnegative matrices appearing in real situations: (a) customer purchase data, (b)monochrome image data, (c) time-frequency amplitude of audio signal.

場合や元データのランクが不明な場合には補完精度が低下してしまう．そこで，本論文ではNMF に基づく行列補完について改善を加えることで，より現実的な条件においても高精度な補完品質を安定させることを目的とする．

1.3 本論文の構成
2 章では，本論文で主に扱う NMF の定式化と行列補完への適用手法について説明し，従来手法の持つ問題点について述べる．3 章では，NMF に適用することで補完性能の向上が期待できる正則化項について説明する．ここでは，従来手法とともにコサイン類似度正則化項とL1 ノルム正則化項について説明し，その効果について述べる．4 章では，コサイン類似度正則
化付き NMF 及び L1 ノルム正則化付き NMF の定式化を行い，その特性について確認する．

5 章では，人工欠損行列の補完実験によって，コサイン類似度正則化付き NMF 及び L1 ノル
ム正則化付き NMF の補完精度を従来手法と比較し，その性能を評価する．6 章では，音響信号の補完実験を行い，より実用的な条件下における補完性能を評価する．7 章で本論文をまと
める．

5

第2章

NMF を用いた行列補完
2.1 はじめに本章では，2.2 節で NMF の基礎理論について説明し，定式化を行う．2.3 節では NMF を用いた行列補完の手法について説明し，またその課題についても述べる．2.4 節で，本章についてまとめる．

2.2 NMF
2.2.1 NMF の概要
NMF とは非負の観測行列を 2 つの非負行列の行列積で近似的に分解する数理アルゴリズムである．このとき，2 つの非負行列の行列積による近似行列を低ランクに制約することで，観測行列の低ランク性を仮定した分解が可能になる．NMF は，この非負制約及び低ランク近似の両性質により，観測行列に頻出する少数の潜在的なパターンを抽出することが可能な教師なし学習である．実世界では，音の振幅時間周波数行列（振幅スペクトログラム）や画像等の非負行列で表されるデータが多く存在する．そのような非負値のデータを，非負値の構成要素に分解するために NMF を適用する．得られる非負行列は解釈が容易であるため，商品の購買データパターン抽出 [5] や音源分離 [6] 等のパターン認識を活用する分野で広く用いられている．

NMF による行列の近似分解は次の式で表される．
X ≈ WH(2.1)

I×K
ここで，X ∈ RI×J
≥0 は全要素が非負の観測行列であり，W = [w1 w2 · · · wK ] ∈ R≥0 及び

H = [h1 h2 · · · hK ]T ∈ RK×Jは NMF で推定すべき非負変数行列で，W は基底行列，H≥0は係数行列と呼ばれる．I 及び J はそれぞれ観測行列の行数と列数である．また，W の列ベ
クトルである wk は基底ベクトルと呼ばれ，K は基底ベクトルの個数を示す．k = 1, 2, · · · , K
は基底ベクトルのインデクスを示す．NMF は基底ベクトルの数 K を K ≪ min(I, J) とな

6

第 2 章 NMF を用いた行列補完

Fig. 2.1. Decompose observation matrix using NMF.

るように設定することで，観測行列 X を式 (2.1) のように低ランク近似する．Fig. 2.1 は，NMF による行列の分解近似の例である．このように，観測行列を非負値の基底行列と係数行
c を式 (2.1) のように列で表現することで，ランクを基底ベクトルの数 K に抑えた近似行列 X得ることができる．

NMF を適用するにあたって注意すべき点として，NMF で得られる W と H には原理的に次の 2 つの任意性が存在する．

• 順序の任意性（permutation ambiguity）
W に含まれる K 本の基底ベクトル w1 , w2 , · · · , wK の順序を入れ替えても，H に含まれる係数ベクトル h1 , h2 , · · · , hK を同様の順序に入れ替えれば，W H は一致する．


例： w1

 
 hT

1
w2
T = w2
h2

 
 hT
2
w1
hT
1

• 大きさの任意性（scale ambiguity）
W に 含 ま れ る K 本 の 基 底 ベ ク ト ル w1 , w2 , · · · , wK の そ れ ぞ れ に 対 し て
C1 , C 2 , · · · , C K と い う ス カ ラ ー を 乗 じ て も ， H に 含 ま れ る 係 数 ベ ク ト ル−1h1 , h2 , · · · , hK の そ れ ぞ れ に 対 し て C1−1 , C2−1 , · · · , CK
を 乗 じ れ ば ，W H は

一致する．



例： w1

w2

 
 hT

1
T = C 1 w1
h2

C2 w2



 C1−1 hT
1
C2−1 hT
2

また，W 及び H の推定問題は解析的に解けるわけではなく，次項で説明する最適化問題を解く必要がある．さらに，局所解の存在により常に同じ結果を得られるとは限らず．与える初期値や用いる最適化アルゴリズムによって異なる分解となる．

2.2.2 NMF の定式化
NMF では，次の最適化問題を解くことで変数行列 W 及び H を推定する．
Minimize D(X|W H) s.t. wik , hkj ≥ 0 ∀i, j, k
W ,H

(2.2)

2.2 NMF

7

ここで，xij ，wik ，及び hkj はそれぞれ X ，W ，及び H の要素であり，i = 1, 2, · · · , I
及び j = 1, 2, · · · , J はそれぞれ観測行列の行インデクス及び列インデクスを示す．また，D(·|·) は 2 つの入力行列間の乖離度を測る関数である．本論文では，次式で表される一般化
Kullback–Leibler（KL）ダイバージェンス [7] DKL (·|·) を用いる．

XX
aij
DKL (A|B) =aij log− (aij − bij )
bij
i
j

(2.3)

ここで，aij 及び bij はそれぞれ A ∈ RI×J 及び B ∈ RI×J の（i, j ）要素を表す．
一般化 KL ダイバージェンスに基づく NMF（KLNMF）の反復最適化更新式を求める．また，
本論文では NMF の反復更新式の導出に補助関数法（majorization-minimization algorithm）[8] と呼ばれる最適化アルゴリズムを用いる．まず，KLNMF の最適化問題は式 (2.3) より，
次式のようになる．

Minimize J (W , H) s.t. wik , hkj ≥ 0 ∀i, j, k
W ,H
"
!#
XX
X
xij
J (W , H) =xij log P− xij −
wik hkj
k wik hkj
i
j

(2.4)
(2.5)

k

この目的関数を変形し次式を得る．

J (W , H) =XXixij log xij − xij log

Xjwik hkj − xij +kX!
wik hkj(2.6)

k式 (2.6) の目的関数の wik 及び hkj に関する偏微分を計算する場合，第 2 項の k に関する総和を含む対数関数の偏微分が原因で W や H の停留点を求めることができない．そこで，負対数関数は凸関数であることから，付録 (A.1) に示す Jensen の不等式を適用することで，式(2.6) を上から抑える補助関数を設計する．式 (2.6) の第 2 項に対して，次の不等式を考える．
− log

X

wik hkj = − log

k

≤−

X
k

ここで，δijk ≥ 0 は補助変数であり

X

δijk

wik hkj
δijk

δijk log

wik hkj
δijk

k

(2.7)

P

k δijk = 1 を満たす．式 (2.7) は，Jensen の不等式において，凸関数 f : R → R を f (·) = − log (·) として適用した不等式である．この不等式を式(2.6) の第 2 項に適用することで，補助関数を次式のように定義できる．
J (W , H) ≤ J + (W , H, ∆)
=

XX
i

j

xij log xij − xijXkXwik hkj− xij +
wik hkjδijk logδijkk

ここで，∆ ∈ RI×J×K
は補助変数 δijk で構成される 3 階のテンソルである．
≥0

!
(2.8)

8

第 2 章 NMF を用いた行列補完

次に，補助関数 J + (W , H, ∆) に対する補助変数 ∆ の最小化について考える．いま，補助
変数には等式制約条件

Pk δijk = 1 が課せられているため，次の等式制約条件付き最適化問題

を考える．

Minimize J + (W , H, ∆) s.t.

X∆δijk = 1 ∀i, j

(2.9)

k

最適化問題 (2.9) の Lagrange 関数は次式で与えられる．

X

L(W , H, ∆, λ) = J + − λ!
δijk − 1(2.10)k

従って，δijk 及び λ の偏微分をそれぞれ 0 とおくと次式が得られる．

wik hkj∂L= −xij log−1 −λ=0
∂δijk
δijk
X
∂L
=−
δijk + 1 = 0∂λ(2.11)(2.12)k

式 (2.12) は等式制約条件

PK

k=1 δijk = 1 そのものである．また，式 (2.11)，δijk > 0，及びwik , hkj ≥ 0 より次式となる．
δijk = wik hkj exp

λ
−1
xij


(2.13)式 (2.13) の両辺に対して k について総和をとると次式を得る．

X

よって，等式制約条件

Pδijk =kX

wik hkj exp

k

λ
−1
xij


(2.14)

k δijk = 1 より次式を得る，

X
wik hkj exp

k


λ
−1 =1
xij

(2.15)

未定乗数 λ は k に依らない定数であるため，結局，次式となる．



λexp 1 −
xij


=

X

wik hkj(2.16)k

これを式 (2.13) に代入することで，補助変数の最小解が次式を得る．

wik hkj∀i, j, k′′k′ wik hk jδijk = P(2.17)

従って，式 (2.17) が補助関数 (2.9) の補助変数に関する最小解を与える．今，変数の更新の反復回数を t とし，各変数の上付き文字として表示すると，次のようになる．
(t) (t)wik hkj(t+1)
δijk = P(t) (t)k′ wik′ hk′ j∀i, j, k(2.18)

2.2 NMF

9

続いて，設計した補助関数 J + (W , H, ∆) を本来の変数 W 及び H について最小化する．
初めに，補助関数 (2.8) を wik で偏微分して 0 とおき，次式を得る．


X
∂J +
δijk hkj=−xij δijk+ hkj = 0
∂wik
w
h
δ
ik
kj
ijk
j

(2.19)式 (2.19) を整理すると次式のようになる．

−

X
j

xij

δijk X+
hkj = 0
wik
j
X
X
wik
hkj =
xij δijkjjPj xij δijk
P
j hkj

wik =

(2.20)式 (2.20) が補助関数 (2.8) を最小化する wik となる．反復回数 t を記述し，次式を得る．

P

(t+1)
j xij δijk
P (t)
j hkj

(t+1)
wik
=

(2.21)

同様に，補助関数 (2.8) を hkj で偏微分して 0 とおくと，hkj の最小解も次のように得られる．


X
∂J +
δijk wik=−xij δijk+ wik = 0∂hkj
wik hkj δijk
iX
δijk X+
wik = 0
−
xij
hkj
i
i
X
X
hkj
wik =
xij δijk
i

i=1

hkj =

P
i xij δijk
P
i wik

(2.22)

反復回数 t を考慮して記述すると，次のようになる．

P
(t+1)
hkj =

(t+1)
i xij δijk
P (t+1)i wik(2.23)

従って，式 (2.7) で補助変数を更新し，式 (2.21) 及び式 (2.23) で本来の変数を更新する．この 2 つの手順を繰り返すことで，目的関数 (2.6) の値を最小化する W 及び H を求めることができる．
なお，補助変数の反復更新式 (2.18) を式 (2.21) 及び式 (2.23) にそれぞれ代入することで，

10

第 2 章 NMF を用いた行列補完

次式のような統合された反復更新式に書き換えることができる．

P
(t+1)

wik

P

=

P
(t+1)

hkj

wik′ hk′ j

(t)

h

(t)

xij
(t)

k′

P

(t)

wik′ hk′ j

hkj(2.24)

(t)
j hkj

w

(t+1)

h

(t)

P ik (t+1)kj (t)

P
P
(t)

k′

P

j

i xij

=

(t)

(t)
j hkj

P
(t)
= wik

w

P ik (t)kj (t)

j xij

i

k′

wik′

hk′ j

(t+1)
i wik
xij

P

= hkj

(t+1)

(t+1)

k

′ wik′

P

(t)

hk′ j

wik

(2.25)

(t+1)

i wik

式 (2.24) 及び式 (2.25) は，次のように行列形式で表現することができる．

W ←W ⊙X
T
WH H
1(I×J) H T(2.26)H←H⊙

W T WXH
W T 1(I×J)(2.27)

ここで，見やすさのために反復回数 t の表記は省略し，変数更新を表す演算子 ← を用いている．また，1(I×J) 全ての要素が 1 の I × J の行列であり，演算子 ⊙ は要素毎の積（Hadamard

product），行列の分数は要素毎の商を表す．以上より，KLNMF の反復更新式が導出された．

2.3 NMF の行列補完への応用
2.3.1 バイナリマスク行列を用いた NMF に基づく行列補完
NMF に基づく行列補完には，ナイーブな手法として欠損要素のインデクスを行列とした，
バイナリマスク行列を用いる方法がある．この手法の概要は Fig. 2.2 の通りである．NMF に基づく行列補完では，初めに，観測行列の欠損要素を 0，その他を 1 とするバイナリマスク行
列 M ∈ {0, 1}I×J を生成する．生成されたバイナリマスク行列を観測行列及び近似行列のそれぞれに要素積することで，欠損要素を無視した行列補完が可能となる．即ち，観測可能な要素についてのみ，観測行列と近似行列の目的関数値を最小化することで，観測不能な要素についても行列全体の代数構造に基づいて何らかの非負値を補完する手法である．

NMF に基づく行列補完は，次式の最適化問題で表される．
Minimize J (W , H) s.t. wik , hkj ≥ 0 ∀i, j, k
W ,H
!#
"
X
XX
xij
− xij −
wik hkj
J (W , H) =mij xij log P
k wik hkjijk(2.28)(2.29)

2.3 NMF の行列補完への応用

11

Fig. 2.2. NMF-based matrix completion using binary masks. Missing values are ignored inbcalculation of KL divergence during estimation of approximated matrix X.

ここで，mij はバイナリマスク行列 M の要素である．バイナリマスク行列は定数行列である
ため，2.2.2 項と同様の手順で次式の更新式を導出できる．

P
(t+1)

wik

j

(t)

k′

P

= wik

P
(t+1)
(t)
hkj = hkj(t)xij
P mij(t)(t) hkjiwik′ hk′ j(t)(2.30)

j mij hkj

(t+1)
xij
P mij(t+1)(t) wik
k′

P

wik′

hk′ j

(t+1)
i mij wik(2.31)式 (2.30) 及び式 (2.31) は行列形式で表現すると次式となる．

W ←W ⊙(M ⊙X)T
WH H
M HT(2.32)⊙X)
W T (M
WH
H←H⊙(2.33)
W TM
以上で，バイナリマスク行列を用いた NMF に基づく行列補完の更新式が導出された．

2.3.2 バイナリマスク行列を用いた NMF に基づく行列補完の課題
バイナリマスク行列を用いた NMF に基づく行列補完の適用例を Fig. 2.3 に示す．Fig. 2.3
は，ランク 5 の 25 行 25 列の正解行列の内，50% の要素を欠損させた観測行列に対し，バイ
ナリマスク行列を用いた NMF に基づく行列補完を適用したときの近似行列である．このように，NMF に基づく行列補完は補完精度が高く，Fig. 2.3 においても，目視では正解行列と近似行列の差異を確認できない精度である．しかしながら，これは理想的な状態における欠損補完の例である．即ち，欠損前の観測行列のランクが 5 であることを既知とし，NMF の基底ベクトル数を K = 5 に設定した結果を示している．さらに，観測行列にノイズも付与されておらず，完全な最適化が達成された場合には式 (2.29) の目的関数の値を 0 にすることができるケースである．

NMF において，欠損前の観測行列（正解行列）のランクの推測精度（即ち K の設定値）がc の近似精度に与える影響は大きい．NMF が低ランク近似するうえで，基推定結果である X

12

第 2 章 NMF を用いた行列補完

Fig. 2.3. Example of performing NMF-based matrix completion using binary mask. This method
can estimate highly accurate completion results when number of bases is set to ideal
value.

底行列の各基底ベクトルの役割は大きく，その本数となるランクは非常に重要な設定値となる．このため，NMF に用いるランク K が実際の正解行列のランク（欠損前の観測行列のランク）と 1 異なるだけでも，行列補完の精度が大きく低下し，Fig. 2.3 のような結果は得られな
くなる．また，Fig 2.3 の例では，正解行列と観測行列の間に欠損は発生しているが，観測行列にノイズは与えられていない．実際にはいかなる観測系においてもノイズは発生してしまう
ため，多くの場合で観測行列は正解行列のランクによらずフルランクとなる．よって，NMFの仮定する低ランク性は成立せず，結果として NMF の近似精度に影響を与える．
バイナリマスク行列を用いた NMF に基づく行列補完は，低ランクかつ正解行列の基底数が既知の時，高精度な補完を可能とする．一方で，実際の観測行列は観測時のノイズによってフルランクになり，またノイズと欠損が生じる前の行列のランクの推定も困難であるため，その性能を十分に発揮することが難しい．そのため，バイナリマスク行列を用いた NMF に基づく
行列補完は，観測行列に含まれるノイズやランクの推定誤差（NMF で用意する基底ベクトルの個数に関するモデル誤差）に対して頑健になるような工夫が求められる．

2.4 本章のまとめ
本章では，NMF の基礎理論とバイナリマスク行列を用いた NMF に基づく行列補完の手法について説明した．この手法は，本論文における従来手法であり，この手法で達成することのできる行列補完の性能がベースラインとなる．また，バイナリマスク行列を用いた NMF に基づく行列補完が抱える問題について述べた．次章では，観測行列に含まれるノイズや用意する
基底ベクトルの個数への頑健性を持つことが期待できる正則化について説明する．

13

第3章

NMF の正則化項付与
3.1 はじめに本章では，NMF に対してこれまで提案されている正則化項について説明する．3.2 節では，音源分離の分野で用いられる半教師あり NMF（semi-supervised NMF : SNMF）について示し，SNMF の改良型であるコサイン類似度正則化条件付き SNMF（cosine-similarity-

regularized SNMF : CSRSNMF）に用いられているコサイン類似度正則化項について述べる．3.3 節では，L1 ノルム正則化付き NMF（L1 -norm-regularized NMF : L1 -RNMF）について説明し，L1 ノルム正則化項がもたらす効果について述べる．3.4 節では本章についてまとめる．

3.2 CSRSNMFNMF には，音源分離の分野等で用いられる SNMF [9, 10] という手法がある．これは，混合信号の振幅スペクトログラムを，目的音源及び非目的音源の振幅スペクトログラムの和に近似的に分解する手法である．目的音源と非目的音源の振幅スペクトログラムはいずれも基底行
列と係数行列の行列積で表される．より具体的に，SNMF では，混合信号の振幅スペクトロ
グラム Y ∈ RI×J
≥0 を次の最適化により分離する．

Minimize D(Y |F G + HU ) s.t. gkj , hil , ulj ≥ 0 ∀ i, j, k, l
G,H,U

(3.1)

ここで，非目的音源の基底行列 H のランクは L であり，l = 1, 2, · · · , L は非目的音源の基底
行列 H に含まれる基底ベクトルのインデクスを示す．F ∈ RI×J
≥0 は事前に分離対象である目
的音源の学習データから特徴量を学習した教師基底行列であり，G ∈ RK×J
は目的音源の係≥0
L×J
数行列である．また，H ∈ RI×L
≥0 及び U ∈ R≥0 はそれぞれ非目的音源を表すための基底行

列と係数行列であり，gkl ，hil 及び ulj はそれぞれ G，H 及び U の非負要素である．

SNMF は混合信号中に目的音源と非目的音源が類似する要素を含むときに分離精度が劣化する問題がある [10]．これは，類似する要素が教師基底行列 F または非目的音源の基底行列

H のいずれを用いても表現できてしまうためである．即ち，類似成分を F と H が取り合う14

第 3 章 NMF の正則化項付与

Cosine similarity

Fig. 3.1. Two vectors and their cosine similarity. Cosine similarity becomes zero when two
vectors are orthogonal.

ことで，目的音源と非目的音源が完全に分離できず混合してしまうことが原因である．
この問題に対し，非目的音源の基底行列 H をできるだけ教師基底行列 F と直交になるように推定する改良を加えた CSRSNMF [11] が提案されている．これは，F と H の各基底ベク
トルの総組み合わせのコサイン類似度の総和を正則化項として目的関数に導入することで，基底行列 F 及び H 間ができるだけ直交するように誘導する手法である．基底行列 F 及び H間の直交性を促すことにより，教師基底行列 F が含む要素を非目的音源の基底行列 H が含まないようになり，混合信号中の類似した成分を取り合う現象を抑制することができる．
直交の度合いを計る指標として用いているコサイン類似度に基づく正則化項とは次式で表される．

Rcos (F , H) =XX

cos(fk , hl )

(3.2)

P
fik hklcos(fk , hl ) = P 2 i 1 P 2 1( i fik ) 2 ( i hil ) 2

(3.3)

k

l

ここで，fk 及び hl はそれぞれ F 及び H の k 番目及び l 番目の基底ベクトル（列ベクトル）である．式 (3.3) のコサイン類似度は，Fig. 3.1 に示すようにベクトル間の類似性を，そのベクトル間のコサイン値を基に表現する尺度である．2 つのベクトルの組が一次従属の関係にあるとき，コサイン類似度は 1 または −1 となり，直交関係にあるとき，コサイン類似度は 0 となる．なお NMF の文脈では，基底ベクトルは非負の要素から成ることが保証されているため，式 (3.3) のコサイン類似度は必ず 0 から 1 の値を取る．

CSRSNMF の最適化問題と目的関数は次式によって定義される．
Minimize JCSRSNMF s.t. gkj , hil , ulj ≥ 0 ∀ i, j, k, l

(3.4)

JCSRSNMF = DKL (X|F G + HU ) + µRcos (F , H)

(3.5)

G,H,U

ここで，µ > 0 は正則化項の重み係数である．式 (3.5) より，CSRSNMF は通常の KLNMFにコサイン類似度を最小化する正則化項を付与する形となっていることがわかる．前述の通3.3 L1 -RNMF

15

り，NMF におけるコサイン類似度の最小値は 0 となるため，式 (3.5) における正則化項は変数 H の推定において，できるだけ F と直交することを誘導する効果がある．
一般に，最適化問題における目的関数をデータ近似項（data ﬁdelity term）と呼ぶ．また，変数の推定において何らかの性質や構造を誘導するために追加する項を正則化項（reguralizationterm）と呼ぶ．式 (3.5) は，第一項がデータ近似項であり，第二項が正則化項に該当する．正則化項の強さは重み係数 µ で制御され，この値が大きいほど最適化に正則化が強く作用することになる．

3.3 L1 -RNMF
NMF の文脈では，係数行列の要素をスパース（sparse，疎）に誘導することで，観測行列の特徴や潜在的な構造をより強く反映した基底行列を頑健に推定する手法が提案されている．この手法はスパース NMF と呼ばれ，独自のスパース性基準を定義して正則化に用いる手法 [12]や係数行列の L1 ノルムを正則化に用いる手法 [13, 14] 等が提案されている．特に L1 ノルムを用いた手法は，有名な least absolute shrinkage and selection operator（lasso）回帰 [15]の正則化を NMF に導入したものであり，一般的な正則化といえる．

L1 ノルムは || · ||1 と表され，行列 X において，次式のように定義される．
||X||1 =XXi|xij |(3.6)

j

これは即ち，行列の全要素の絶対値総和を表す．この値が小さくなることは，行列 X 中の非零要素の個数が少ないことを意味する．但し，本来は L0 擬ノルムが行列の非零要素の個数を表したものであり，L1 ノルムは L0 擬ノルムの最良な凸緩和であるため，最適化の正則化として L0 擬ノルムの代わりに用いられるものである．

NMF 等のパラメータの更新によって目的関数を最小化するアルゴリズムにおいて，L1 ノルム正則化項はそのパラメータをスパースにする効果を持つことが知られている [16]．このような L1 ノルム正則化項を NMF の係数行列に対して付与したものが L1 -RNMF である．

L1 -RNMF は次式によって定義される．
Minimize JL1 −RNMF s.t. wik , hkj ≥ 0 ∀i, j, k

(3.7)

JL1 −RNMF = D(X|W H) + ν||H||1

(3.8)

W ,H

ここで，ν は L1 ノルム正則化項の重み係数である．

L1 -RNMF は NMF の持つ低ランク性に加え，係数行列にスパース性を持つことになる．この結果，観測行列を特徴的な値によって少ないパラメータで表そうとすることにより，より強力に特徴量の抽出が可能となる．

16

第 3 章 NMF の正則化項付与

3.4 本章のまとめ
本章ではコサイン類似度正則化と L1 ノルム正則化の 2 つの正則化について説明した．いずれの正則化も NMF の特徴抽出をノイズの発生や基底ベクトルの個数のモデル誤差に対して頑健にする効果が期待できる手法である．次章では，本論文の提案手法として，この 2 つの正則化をを付与した NMF に基づく行列補完の定式化を行う．

17

第4章

正則化付き NMF に基づく行列補完
4.1 はじめに本章では，コサイン類似度正則化及び L1 ノルム正則化の 2 つの正則化を，NMF に基づく
行列補完へと適用する．4.2 節では，NMF に基づく行列補完に正則化を加える動機について述べる．4.3 節では，SNMF ではなく通常の NMF におけるコサイン類似度正則化付き NMF（cosine-similarity-regularized NMF: CSRNMF）の定式化と行列補完への適用について説明する．4.4 節では，L1 -RNMF に基づく行列補完の定式化について説明する．4.5 節では，正則化項が NMF にもたらす影響について確認し，4.6 節で本章についてまとめる．

4.2 動機
2.3 節に記載の通り，通常の NMF に基づく行列補完は，欠損やノイズが発生する前の観測行列である正解行列のランクが既知で，かつ低ランクである場合において高精度な補完を可能とする．一方，実問題においてこの条件を満たす場合は少ない．実際には正解行列が低ランクである場合においても，ノイズが発生することによりフルランク化することで低ランク性を失ってしまう．さらに，欠損が生じて観測行列となった後に，正解行列のランクを正確に推定することは難しい．そのため NMF の実用上は基底ベクトル数 K を適当な値や経験的な値で定めることが多く，当然正解行列のランクと適当に定めた K の間にはモデル誤差が生じてしまう．観測行列の低ランクな近似行列を得ようとする NMF において，基底行列の各基底ベク
トルは非常に重要な役割を果たすにもかかわらず，正解行列のランクと基底ベクトル数 K の設定値が異なるならば，近似精度に多大な影響を与える問題がある．以上から，NMF に基づ
く行列補完の精度向上には，観測行列に含まれるノイズや正解行列のランクと基底ベクトル数のモデル誤差に対し，より頑健な性能を持つ必要があると考えられる．
本論文では，コサイン類似度及び L1 ノルムの正則化を用いることで，観測行列に含まれるノイズの影響を受けにくくなるように改良された NMF に基づく行列補完を提案する．これらの正則化はいずれも，K 個の基底ベクトルの個々の独立性を高める誘導を行うものであり，ノ

18

第4章

正則化付き NMF に基づく行列補完

イズのような低ランク性に寄与しない要素を無視する効果が期待できる．そのため，観測信号
中のノイズに対する頑健性が向上することが期待される．また，前述の効果は正解行列のランクと基底ベクトル数 K 間のモデル誤差に対する頑健性ももたらす可能性がある．

3.2 節において説明したコサイン類似度正則化は，CSRSNMF においては，目的音源と非目的音源の基底行列間を直交させるために活用される．本論文の提案手法である行列補完の用途に応用するためには，2 つの基底行列間ではなく，1 つの基底行列内に含まれる各基底ベクト
ルを相互に直交させることを考える．

4.3 CSRNMF に基づく行列補完
4.3.1 コサイン類似度正則化の NMF への適用
コサイン類似度による正則化項を KLNMF に適用する．目的関数は次式となる．

Minimize JCSRNMF (W , H) s.t. wik , hkj ≥ 0 ∀i, j, k

(4.1)

W ,H

JCSRNMF (W , H) =DKL (X|W H) + µ

K
K
X
X

Rcos (wk , wk′ )(4.2)

k=1 k′ =k+1

=

XX
j

X
xij
xij log P− xij +
wik hkj
k wik hkj

K
X

K
X

i

+µ

!

k

P

wik wik′
P 2i 1 P 2 1( i wik ) 2 ( i wik′ ) 2
k=1 k′ =k+1

(4.3)

ここで，k ′ = k + 1, k + 2, · · · , K は基底インデクスであり，常に k < k ′ の関係にある．

4.3.2 補助関数法に基づく上限関数の導出
前項で示した式 (4.3) から，JCSRNMF の上限関数を求める．煩雑になるのを避けるため，
式 (4.2) の第一項と第二項を分けて考える．
式 (4.2) の第一項は，2.2.2 項と同様に，Jensen の不等式を利用して次のように変形できる．

XXj

xij logxij − xij

X

i

k

ここで，δijk > 0 は補助変数であり，

X
wik hkjδijk log− xij +
wik hkj
δijk

!
(4.4)

k

P

k δijk = 1 を満たす．等式条件は次式となる．

δijk = P

wik hkj′′
k′ wik hk j

(4.5)

4.3 CSRNMF に基づく行列補完

19

式 (4.2) の第二項は次のように変形できる．

µ

K
K
X
X

P

K
K
X
X

wik wik′
P 2i 1 P 2 1( i wik ) 2 ( i wik′ ) 2
k=1 k′ =k+1
P
K
K
X
X
ik wik′
i w
P
=µ
2
P 2 1
wik
( i wik′)2
k=1 k′ =k+1
i ϵik ϵik


!− 12
K
K
X
X
X
XX
3
′
w
ik 

≤µ
wi2′ k′(4.6)
wi′ k′ ϵi2′′ k
′′ k
w
i
′
′
′′
′

Rcos (wk , wk′ ) = µ

k=1 k′ =k+1i

k=1 k =k+1ii

式 (4.6) は Jensen の不等式を利用して変形している．ここで，ϵik > 0 はPi ϵik = 1 を満たす補助変数であり，式 (4.6) の等号成立条件は次式で与えられる．

w2
ϵik = P ik 2
i′ w i′ k(4.7)

以上より，式 (4.2) の第一項及び第二項を式 (4.4) 及び式 (4.6) に置き換えた JCSRNMF の+上限関数 JCSRNMF
は次式となる．

+JCSRNMF ≤ JCSRNMF=!
X
wik hkjj xij logxij − xijδijk log− xij +
wik hkj
δijk
i
k
k


!− 12
K
K
XX
X X
X
3
wi′ k wi′ k′ ϵi2′′ k(4.8)+µ
wi2′ k′
w
i′′ k
′
′′
′
′

XX

k=1 k =k+1

X

i

i

i

これにより，式 (4.8) を最小化することで最適化を行う．次項では基底行列及び係数行列の更新式を導出する．

4.3.3 CSRNMF の更新式の導出
基底行列 W の更新式は，式 (4.8) を wik で偏微分し 0 と置いた次式から得られる．
+
X xij δijk
∂JCSRNMF
=−
+ hkj
∂wik
wik
j


!− 12 
X
X
X 3 1
3 X
1
2

−

wi2′ k′+µwi′ k′ wi′ k + wik′
ϵi2′ k
2 ϵik
w
wi′ k
ik
′
′
′
′
k =k+1

=0

i

i ̸=i

i ̸=i

(4.9)

20

第4章

正則化付き NMF に基づく行列補完

2式 (4.9) の両辺に wik
を乗じると，次式となる．

−wik

X

2
xij δijk + wik

X

j

X

X

k′ =k+1i′hkj − µjXXk′ =k+1

i′

+µ

!− 12
wi2′ k′

3
2
ϵik

X

!− 12

2
wik
wik′

wi2′ k′wi′ k′ wi′ ki′ ̸=iX 3 1
ϵi2′ k=0
wi′ k′(4.10)i ̸=i

式 (4.10) を wik の次数ごとにまとめると，次式が得られる．


X
X
2 
wik
hkj + µwik′j

+ wik −k′ =k+1

X



X


X 3 1

ϵi2′ k
wi′ k
′

!− 21
wi2′ k′

i′



i ̸=i

X

X

k′ =k+1i′32xij δijk  + −µϵikj!− 12 
wi2′ k′Xwi′ k′ wi′ k  = 0 (4.11)

i′ ̸=i2
ここで，次式を用いて式 (4.11) の wik，wik 及び第三項の係数をそれぞれ置き換える．

aik =

X
j

bik = −

X

hkj + µ

X

X

wik′


X 3 1
2


ϵik
wi′ k
′

wi2′ k′i′k′ =k+1

3
2
i′ k(4.12)i ̸=i

xij δijk(4.13)

jcik = −µϵ!− 12 XX

k′ =k+1i′!− 12 wi2′ k′X
wi′ k′ wi′ k (4.14)i′ ̸=i式 (4.12)，式 (4.13)，及び式 (4.14) を式 (4.11) に適用すると，次式が得られる．
2
aik wik+ bik wik + cik = 0(4.15)+
式 (4.15) は wik の二次方程式となるため，上限関数 JCSRNMFの wik に関する停留点は二次方程式の解の公式で次のように与えられる．

wik =−bik ±

pb2ik − 4aik cik
2aik(4.16)

ここで，式 (4.16) の ± は wik が非負になるように決定する．wik の更新則は等号成立条件の式 (4.5) 及び式 (4.7) を式 (4.16) に代入することで求めることができる．しかし，非常に複雑

4.3 CSRNMF に基づく行列補完になるため，代わりに aik ，bik 及び cik の更新式をそれぞれ示す．

aik ←

X

hkj + µ

wik′

X

X

xij P

j


cik ← −µw2
P ik!− 12 2i′ wi′ kXwi2′ kX!− 23 wi2′′ k

(4.17)i′′i′ ̸=iwik hkj′′
k′ wik hk j 32

wi2′ k′i′k′ =k+1

jbik ← −

X21

(4.18)

XX

k′ =k+1i′!− 12 
wi2′ k′Xwi′ k′ wi′ k (4.19)i′ ̸=i以上より，更新後の wik は式 (4.17)，式 (4.18)，及び式 (4.19) によって更新した aik ，bik ，及び cik を式 (4.16) の計算に用いることで得られる．なお，µ = 0 のとき，式 (4.16)，式 (4.17)，式 (4.18)，及び式 (4.19) は 2.2.2 項で示した KLNMF の基底行列 W の更新式 (2.24) と一致
する．
係数行列 H の更新式は，式 (4.8) を hkj で偏微分し 0 と置いた次式から得られる．
+
X xij δijk
∂JCSRNMF=−+ wik = 0
∂hkj
hkj
i

(4.20)

式 (4.20) を hkj について解き，等式条件の式 (4.5) を代入すると次式が得られる．

1 X−xij δijk + hkj wik = 0
hkj i
X
X
hkj
wik =
xij δijkiiPi xij δijk
hkj = Pi wik
P
P wik hkj
i xij
wik′ hk′ j
P k′
=
wik
P i xijPw
iwik′ hk′ j ik
k′
P
= hkj
i wik

(4.21)

これより，hkj の更新則は次式となる．

P
hkj ← hkj

i

P xij
w
wik′ hk′ j ik
k′
P
i wik

(4.22)

従って，CSRNMF の係数行列 H は式 (4.22) によって更新される．これは 2.2.2 項で示した
通常の KLNMF における係数行列 H の更新式 (2.25) と一致する．

4.3.4 CSRNMF の行列補完への適用
CSRNMF に基づく行列補完は，2.3 節で述べた手法と同様のアイデアによって実現する．
従って，観測行列の欠損要素に対応する要素を 0，その他の要素を 1 とするバイナリマスク行

22

第4章

正則化付き NMF に基づく行列補完

列 M を定義することで，CSRNMF に基づく行列補完は次式で表される．

Minimize JCSRNMF (W , H) s.t. wik , hkj ≥ 0 ∀i, j, k
(4.23)
W ,H
X X
XX
X
Rcos (wk , wk′ )JCSRNMF (W , H) =mij DKL (xij |wik hkj ) + µij

k=1 k′ =k+1

k(4.24)
ここで，バイナリマスク行列は W や H に依らない定数行列であるため，基底行列及び係数行列の更新式の導出は 4.3.3 項と同様に行うことができる．よって，基底行列の更新式は式(4.16) 同様に aik ，bik 及び cik を用いることで次式となる．
−bik ±wik =p
b2ik − 4aik cik
2aik(4.25)

このとき，aik ，bik 及び cik はそれぞれ次式となる．

aik ←

X

X

mij hkj + µ

wik′

X

cik ← −µ

wi2′ k′X

X

wi2′ kwi2′′ k



 32

w2
P ik

2
i′ wi′ k(4.26)

i′′

i′ =
̸ i

wik hkj′′
k′ wik hk j!− 23 mij xij Pj!− 12 i′k′ =k+1

jbik ← −

X

(4.27)

XX

k′ =k+1i′!− 12 wi2′ k′X
wi′ k′ wi′ k (4.28)i′ ̸=iまた，係数行列の更新式は次式となる．

P
hkj ← hkji

P mij xij wik
wik′ hk′ j
Pk′i mij wik

(4.29)

以上で，基底行列及び係数行列の更新式が示された．

4.4 L1 -RNMF に基づく行列補完
4.4.1 L1 -RNMF に基づく行列補完の目的関数
L1 -RNMF に基づく行列補完は，観測行列の欠損要素に対応する要素を 0，その他の要素を1 とするバイナリマスク行列 M を定義することで，次式で表される．
Minimize JL1 −RNMF (W , H) s.t. wik , hkj ≥ 0 ∀i, j, k
W ,H
XX
X
JL1 −RNMF (W , H) =mij DKL (xij |wik hkj ) + ν||H||1
i

=

j

XX
i

j

k

mij

X
xij
− xij +
wik hkj
xij log Pw hk ik kjk(4.30)(4.31)
!
+ν

XX
k

hkj

j

(4.32)

4.4 L1 -RNMF に基づく行列補完

23

4.4.2 補助関数法に基づく上限関数の導出
前項で示した式 (4.32) から，JL1 −RNMF の上限関数を求める．式 (4.31) の第一項と第二項を分けて考えると，第二項については偏微分は容易であるため，第一項のみについて導出すればよい．
式 (4.31) の第一項は，2.2.2 項及び 4.3.2 項と同様に，Jensen の不等式を利用して次のように変形できる．

XXixij logxij − xij

mij

j

X
k

ここで，δijk は補助変数であり，

X
wik hkj− xij +
wik hkjδijk log
δijk!
(4.33)kP

δ = 1 を満たす．等式条件は次式となる．
k ijk
δijk = P

wik hkj′′
k′ wik hk j

(4.34)

よって，式 (4.31) の第一項を式 (4.33) に置き換えた J2 の上限関数 J2+ は次式となる．

JL1 −RNMF ≤JL+1 −RNMF
=

XX
i

+ν

mij

xij logxij − xij

j

XX
k

X
k

X
wik hkjδijk log− xij +
wik hkj
δijk

!

k

hkj

(4.35)j

これにより，最小化すべき上限関数が求まったため，次項では基底行列及び係数行列の更新式を導出する．

4.4.3 L1 -RNMF に基づく行列補完の更新式の導出
基底行列 W の更新式は式 (4.35) を wik で偏微分し 0 と置いた次式から得られる．

X∂JL+1 −RNMF
=−
∂wik
j



mij xij δijk
+ mij hkj
wik


=0

(4.36)

24

第4章

正則化付き NMF に基づく行列補完式 (4.36) を wik について解き，等式条件の式 (4.34) を代入すると次式が得られる．

−

X
1 X
mij xij δijk +mij hkj = 0
wik j
j

P
m h
1
j ij kj
=P
wik
m x δ
j ij ij ijkPm x δ
j ij ij ijk
wik = Pm hj ij kjPw hkj
m x P ikj ij ij
w ′ hk ′ j
k′ ik
wik =
P
m h
j ij kj
P P mij xij
hkjj′ wik′ hk′ j
wik = wikPkm hj ij kj(4.37)

これより，wik の更新式は次式となる．

P
wik ← wikj

P mij xij hkj′ wik′ hk′ jPkm hj ij kj(4.38)

従って，L1 -RNMF に基づく行列補完の基底行列 W は式 (4.38) によって更新され，これは2.3.1 項で示した通常の NMF に基づく行列補完における基底行列 W の更新式 (2.32) と一致
する．
係数行列 H の更新式は式 (4.35) を hik で偏微分し 0 と置いた次式から得られる．

X  mij xij δijk
∂JL+1 −RNMF=−+ mij wik + ν = 0
∂hkj
hkj
i

(4.39)

式 (4.39) を hkj について解くと，次式を得る．

X
1 X−mij xij δijk +mij wik + ν = 0
hkj i
i
P
1 X
mij xij δijk = i mij wik + ν
hkj iP
mij xij δijk
hkj = P i
m w +ν
i ij ik(4.40)

補助変数 δijk に等式条件の式 (4.34) を代入すると次式が得られる．

P

kj
mij xij P ik
w ′ hk′ j
k′ ik
hkj =
P
m w +ν
i ij ikP P
mij xij wik
i′ wik′ hk′ j= hkj P k
m w +ν
i ij ikwhi(4.41)

4.5 正則化項による NMF の挙動変化

25

これより，hkj の更新式は次式となる．

P P
mij xij wikihkj = hkj Pik′wik′ hk′ j

mij wik + ν(4.42)

従って，L1 -RNMF に基づく行列補完の基底行列 H は式 (4.42) によって更新される．なお，

ν = 0 のとき，式 (4.42) は 2.3.1 項で示した通常の NMF に基づく行列補完における係数行列
H の更新式 (2.33) と一致する．
注意すべき点として L1 -RNMF において，2.1 節で述べたスケールの任意性は大きな問題となる．これは，最小化問題の目的関数の式 (4.31) において，正則化項 ∥H∥1 は係数行列 H のスケールに依存して小さくなるためである．例えば，限りなく 0 に近い値 ε > 0 を H に乗じて εH と更新だけで ∥H∥1 を限りなく 0 に近づけることができる．このとき，基底行列 Wは ε−1 W と更新することで W H の値が不変となり，式 (4.31) の第 1 項の値は変わらない．
従って，最小化問題 (4.30) は厳密には最小解の無い問題であり解くことはできない [17]．この問題に対して，基底ベクトルを正規化した最適化問題を考えて対処する手法 [14, 17] も提案されているが，本論文では NMF の文脈でよく用いられてきたヒューリスティックな対処法である更新毎の W 及び H の正規化を用いた L1 -RNMF を取り扱う．具体的には，式 (4.38)
及び式 (4.42) で W 及び H を更新した後に，正規化係数 C1 , C2 , . . . , CK を用いて次式の正規化処理を適用する．

Ck =

X

wik

(4.43)

i

wik
∀i, k
Ck
hkj = Ck hkj ∀k, j

wik =(4.44)(4.45)

この正規化式の前後で目的関数の式 (4.31) の値は不変であるため，この処理により間接的に

L1 ノルムの正則化を変数の最適化に影響させ続けることができる．

4.5 正則化項による NMF の挙動変化
4.3 節及び 4.4 節によって，正則化付き NMF の反復更新式を導出した．ここでは，コサイン類似度正則化項及び L1 ノルム正則化項が NMF に与える影響について確認する．

Fig. 4.1 に区間 (0, 1) の一様分布乱数を用いて生成した観測行列を示す．また，この観測行
列を KLNMF，CSRNMF，及び L1 -RNMF によって基底行列と係数行列に分解した結果をFigs. 4.2–4.4 にそれぞれ示す．Figs. 4.2–4.4 は，基底ベクトル数 K を 5 としてそれぞれのNMF によって分解した結果である．ここで，CSRNMF の重み係数 µ は 100000，L1 -RNMFの重み係数 ν は 100000 としている．

Figs. 4.2–4.4 より，それぞれの NMF における近似行列を比較すると，通常の NMF が最も正確に観測行列の細かな値の変化まで近似できていることがわかる．一方で，正則化付き

26

第4章

正則化付き NMF に基づく行列補完

Fig. 4.1. Observed matrix of size 100 × 50 and rank 5.

Fig. 4.2. Example of matrix decomposition by KLNMF when size of observed matrix is 100 × 50
and its rank is 5.

NMF による近似行列は，いずれも観測行列の各行，各列の値を平滑化したような近似となっており，観測行列のパターン，特徴を強調するような結果となっていることがわかる．

Fig. 4.3 より，コサイン類似度正則項付き NMF の基底行列に注目すると，0 に近い値を持つ要素が増加し，各行において値を持つ要素がほとんど 1 つになっていることがわかる．これより，基底行列の列ベクトルを基底ベクトルとして，各基底ベクトル同士のコサイン類似度を考えると，内積が 0 になることから，それぞれの基底ベクトル同士において直交していると考えられる．よって，コサイン類似度正則項によって基底行列の基底ベクトル間を直交できる
ことが確認でき，これにより補完後の推定行列はパターンが著しく強調された行列となっていることが見て取れる．実際には欠損していない要素は観測値をそのまま利用できる為，補完結
果は Fig. 4.3 の推定行列ほど改変されたものにはならないが，それでもこの結果は CSRNMF

4.5 正則化項による NMF の挙動変化

27

Fig. 4.3. Example of matrix decomposition by CSRNMF when size of observed matrix 100 × 50,
its rank is 5, and regularization coeﬃcient is set to µ = 100000.

Fig. 4.4. Example of matrix decomposition by L1 -RNMF when size of observed matrix 100 × 50,
its rank is 5, and regularization coeﬃcient is set to ν = 10000.

が強力な特徴抽出（基底ベクトルの推定）をもたらすことを示唆している．

Fig. 4.4 より，L1 -RNMF の基底行列に注目すると，係数行列の値を持つ要素が Fig. 4.2 に比べ少なくなっていることがわかる．コサイン類似度正則項付き NMF が基底行列要素が減少していたのに対し，L1 -RNMF においては係数行列の要素が減少しており，正解行列の各要素を少ない要素で表現できるように分解していることがわかる．これより，L1 ノルム正則化によって，値を持つ要素を抑制できていることが確認でき，観測行列の持つパターンを鮮明にした近似行列を求められていることがわかる．
以上の結果は，理想的な場合において KLNMF が最良であることも表している．正則化の

28

第4章

正則化付き NMF に基づく行列補完

恩恵は，4.2 節で示したように，観測行列にノイズが含まれている場合や，正解行列のランクと基底ベクトル数 K 間で誤差が生じた場合にもたらされることを期待している．次章では，欠損補完というタスクにおいて，そのような条件で各正則化項がどの程度性能改善に寄与するかについて調査する．

4.6 本章のまとめ
本章では，CSRNMF 及び L1 -RNMF に基づく行列補完の定式化をそれぞれ示した．4.5 節
で述べた通り，それぞれの正則化項はおおむね目的の通りの挙動を示すことを確認できた．次章では，人工的に欠損行列を生成し，欠損行列に対してそれぞれの NMF に基づく行列補完を実験することで，性能を評価する．

29

第5章

人工欠損行列の補完実験
5.1 はじめに本章では，4 章で示した行列補完正則化付き NMF の性能について，人工的に作成した欠損
乱数データを用いた実験により検証する．5.2 節では，補完実験の方法について説明する．5.3
節では，補完実験の結果について解説し，各正則化付き NMF の性能を比較する．5.4 節で本章についてまとめる．

5.2 実験方法
正則化項の付与による行列補完の性能変化を確認するために，通常の KLNMF，CSRNMF，
及び L1 -RNMF の 3 手法を対象とした人工欠損行列の補完実験を行う．本実験では，乱数によって生成した正解行列 T に対しノイズを加算した後に欠損を加えた観測行列 X を用意し，調査対象とする NMF によって行列補完を行い，その補完性能を調査する．
観測行列の生成過程を，Fig. 5.1 に示す．初めに，ランク Ktrue の正解行列 T ∈ RI×J を
サイズが I × Ktrue 及び Ktrue × J の乱数行列の行列積により生成した．乱数には区間 (0, 1)の一様分布乱数を用いた．次に，正解行列 T に上乗せするノイズ行列 N ∈ RI×J
≥0 を作成した．ノイズ行列の生成には区間 (0, 1) の一様分布及び母数 σ = 1.0 のレイリー分布に従う乱数の 2 種類を用いた．レイリー分布とは，時間領域において発生した正規分布に従うノイズ（白色雑音）が，STFT と絶対値処理により振幅スペクトログラムに変換されたときに従う確率分布であり，次章で述べる音響信号の振幅スペクトログラムにおける欠損補完を応用先に考慮した実験条件である．次に，正解行列 T とノイズ行列 N の大きさが所望の信号雑音比（signal-to-noise ratio: SNR）となるような係数 d を次式で求めた．

d = 10

SNR−r
20

SNR = 10 log10



∥T ∥2
∥N ∥2



(5.1)
(5.2)

30

第5章

人工欠損行列の補完実験

Fig. 5.1. Process ﬂow of producing observed matrix.

ここで，r は所望の SNR dB，∥ · ∥2 は次式で定義される L2 ノルムである．

∥X∥2 =

sX X
i

x2ij

(5.3)

j

従って，所望の SNR r を持つ行列は T + dN となる．このようにノイズを正解行列 T に付加した後，全ての要素を確率 α ∈ [0, 1] で欠損させた．
以上の手順により，NMF の補完対象となる観測行列 X を作成した．観測行列 X をもとに，どれだけ正解行列 T に近い補完ができるかを評価することで，各手法の性能を比較する．
正解行列 T にどれだけ近づけることができたかという評価基準には欠損を与えた全要素ごとの一般化 KL ダイバージェンスの総和を用いた．即ち，値が小さいほど，より補完性能が高いと評価できる．
実験手順と各処理で使用するパラメータを Fig. 5.2 に示す．Fig. 5.2 に示す通り，行列の補完処理においては NMF に用いるランクの推定値と正則化項の重み係数を調整する．本実験
では，正解行列のサイズを I = 100 かつ J = 200，ランク Ktrue = 10 に固定し，Fig. 5.2 の各パラメータについて網羅的に条件を変更して実験を行った．詳細な実験条件は Table 5.1 に示す．

5.3 実験結果
Figs. 5.3–5.8 は，ノイズにレイリー分布を用いたときの，SNR r を 0.88 dB，5.99 dB，及び 21.54 dB と設定した場合の各正則化付き NMF による行列補完の実験結果をそれぞれ示したヒートマップ図である．また，その他の条件における結果を付録 B に示す．ヒートマップの縦軸は欠損率 α，横軸は NMF の基底ベクトル数 K （即ち，正解行列 T のランクの予測値であり K = 10 で誤差の無い理想的な条件）を示している．各セルは対応する α 及び K の条
件下において，正解行列及び補完行列間の誤差に関する評価値を，赤色から青色のグラデーションで示している．ヒートマップ上の各色が示す具体的な値は，ヒートマップ右側のカラー
バーによって確認することができる．このカラーバーの値は，正則化項の無い NMF に基づく

5.3 実験結果

31

Fig. 5.2. Experimental conditions changed in evaluation and performance criterion of matrix
completion.

Table 5.1. Experimental conditions

Noise type

Uniform dist./Rayleigh dist.

Number of basis vectors K

6, 7, · · · , 20

SNR r [dB]10−2 to 103 divided by 40
on a logarithmic scale

Missing rate α

0 to 0.8 divided by 20
on a linear scale

Regularizationcoeﬃcient µ and ν

10−3 to 105 divided by 40
on a logarithmic scale

Initialization of
W and H

Uniformly distributed random values
in the range (0,1)

Number of iterations
for parameter updates

500

Number of trials

10 with diﬀerent random seeds

行列補完（従来手法）の結果の値との差分を表している．従って，0 の色は提案手法と従来手法が全く同じ性能であることを表し，負値になるほど提案手法の方が高精度な補完を達成できていることを表している．即ち，ヒートマップにおいてはより赤いほど性能が向上している．

Fig. 5.3 はノイズにレイリー分布を用いた時の，SNR r = 0.88 dB の条件下の各欠損率
α 及び基底ベクトル数 K における CSRNMF の補完実験結果である．この時，重み係数µ = 61584.82 である．Fig. 5.3 の結果より，欠損率 α がおよそ 50% より大きい範囲について，CSRNMF の補完精度が KLNMF を上回っている．一方で，欠損率 α が 50% を下回る範囲のほとんどにおいては改善が見られなかった．

Fig. 5.4 はノイズにレイリー分布を用いた時の，SNR r = 5.99 dB の条件下の結果である．

32

第5章

人工欠損行列の補完実験

Fig. 5.3. Improvement value obtained by CSRNMF: noise type is Rayleigh dist., r = 0.88 dB,
and µ = 61584.82.

この時，重み係数 µ = 0.55 である．Fig. 5.4 より，基底ベクトル数 K が正解行列のランクである 10 に対して小さく，かつ欠損率 α が小さい範囲で改善を示している．補完精度が改善した範囲について，SNR r = 0.88 dB の条件下における改善量と比べると小さくなっている．
一方で，その他の広い範囲については目立つ改善は得られなかった．

Fig. 5.5 はノイズにレイリー分布を用いた時の，SNR = 21.54 dB の条件下の各欠損率 α 及び基底ベクトル数 K における CSRNMF の補完実験結果である．この時，重み係数 µ = 0.55
である．Fig. 5.5 より，基底ベクトル数 K が正解行列のランクである 10 に対して小さく，かつ欠損率 α が小さい範囲で改善を示していることがわかる．一方で，この他の範囲については
目立つ改善は得られなかった．
以上より，CSRNMF はノイズが大きい時，欠損率 α が大きい範囲において基底ベクトル数

K によらず同程度の改善を示し，ノイズが小さい時，欠損率 α 及び基底ベクトル数 K が正解
行列のランクより小さい時において改善を示す．ノイズ及び欠損率 α が大きい条件とは，正解行列の特徴抽出が最も困難な状態であると言える．また，ノイズが小さく欠損率 α 及び基底ベ
クトル数 K が正解行列のランクより小さい条件というのは，正解行列に近い値を取る要素の多くが既知の中で，正解行列を構成する特徴の中でより重要な特徴を抽出する必要がある条件であると言える．このような条件において改善を示したことより，コサイン類似度正則化項が観測行列から正解行列の構造を抽出する効果を NMF に与えていると考えられる．つまり，コサイン類似度正則化項により，KLNMF よりも強力な特徴抽出が可能であると言える．

5.3 実験結果

33

Fig. 5.4. Improvement value obtained by CSRNMF: noise type is Rayleigh dist., r = 5.99 dB,
and µ = 0.55.

Fig. 5.6 はノイズにレイリー分布を用いた時の，SNR r = 0.88 dB の条件下の各欠損率 α 及び基底ベクトル数 K における L1 -RNMF の補完実験結果である．この時，重み係数 ν = 1.44
である．Fig. 5.6 より，実験を行った全ての範囲について，L1 -RNMF の補完精度が KLNMF
を上回っていることがわかる．総じて，補完精度は欠損率 α が低いほど改善量が大きく，同じ欠損率 α について基底ベクトル数 K の値によらず同程度の性能であることがわかる．

Fig. 5.7 はノイズにレイリー分布を用いた時の，SNR r = 5.99 dB の条件下の各欠損率
α 及び基底ベクトル数 K における L1 -RNMF の補完実験結果である．この時，重み係数
ν = 0.89 である．Fig. 5.7 より，実験を行ったほぼ全ての範囲について，L1 -RNMF の補完精度が KLNMF を上回っていることがわかる．総じて SNR r = 0.88 dB の条件時と同様に，補完精度は欠損率 α が低いほど改善量が大きく，同じ欠損率 α について基底ベクトル数 K の値によらず同程度の性能であることがわかる．但し，欠損率 α = 80% に近づくと，改善がほとんど無いか，KLNMF の方が補完精度が良いところが現れ，SNR r = 0.88 dB の条件時よりは改善が減少している．

Fig. 5.8 はノイズにレイリー分布を用いた時の，SNR r = 21.54 dB の条件下の各欠損率 α 及び基底ベクトル数 K における L1 -RNMF の補完実験結果である．この時，重み係数 ν = 0.08
である．Fig. 5.8 より，欠損率 α がおよそ 50% を下回る範囲について，L1 -RNMF の補完精度
が KLNMF を上回っていることがわかる．総じて SNR r = 0.88 dB 及び SNR r = 5.99 dBの条件下同様，補完精度は欠損率 α が低いほど改善量が大きく，同じ欠損率 α について基底

34

第5章

人工欠損行列の補完実験

Fig. 5.5. Improvement value obtained by CSRNMF: noise type is Rayleigh dist., r = 21.54 dB,
and µ = 0.55.

Fig. 5.6. Improvement value obtained by L1 -RNMF: noise type is Rayleigh dist., r = 0.88 dB,
and ν = 1.44.

5.3 実験結果

35

Fig. 5.7. Improvement value obtained by L1 -RNMF: noise type is Rayleigh dist., r = 5.99 dB,
and ν = 0.89.

ベクトル数 K の値によらず同程度の性能であることがわかる．一方で，欠損率 α が 50% を超えると KLNMF の補完精度を下回りはじめ，SNR r = 5.99 dB の条件下より広い範囲で改善が得られなかった．
以上より，L1 -RNMF はノイズの大小に関係なく，欠損率 α が小さい時に基底ベクトル数

K によらず大きな改善を示し，欠損率 α が大きくなるにつれ改善量が小さくなっていた．即ち，ノイズの影響を受けずに補完精度を向上させており，基底ベクトル数 K によって改善量が変化しないことから基底ベクトル数 K の誤差に対する特性は KLNMF による行列補完と同様であると言える．これより，L1 ノルム正則化項はノイズに対する頑健性を向上させており，
行列補完の精度を底上げする正則化項であると言える．

CSRNMF と L1 -RNMF の補完精度を比較すると，いずれの条件においても L1 -RNMF の補完精度が上回ることがわかる．但し，CSRNMF は改善量は小さい一方で，SNR r = 0.88 dBのようにノイズが特に大きい時において欠損率 α が大きい範囲で良い性能を発揮しており，欠損率 α が小さい範囲を得意とする L1 -RNMF とは異なる特徴を持つことがわかる．それぞれの正則化付き NMF の実験結果の特徴をまとめる．
いずれの正則化付き NMF に共通することとして，ノイズが大きい時ほど改善量が大きく，
ノイズが小さくなると改善量も減少していた．これは，Fig. 5.9 に示す結果より説明できる．

Fig. 5.9 は各 SNR の条件における，KLNMF による補完行列と正解行列との一般化 KL ダイ
バージェンスにおける距離の推移である．ノイズにはレイリー分布を用いている．Fig. 5.9 よ

36

第5章

人工欠損行列の補完実験

Fig. 5.8. Improvement value obtained by L1 -RNMF: noise type is Rayleigh dist., r = 21.54 dB,
and ν = 0.08.

り，ノイズが小さくなるにつれ，KLNMF による行列補完の精度が増す．つまり，ノイズが小さくなると行列補完の精度が増すことで一般化 KL ダイバージェンスの値が小さくなり，改善の余地が減少する．このため，改善があった場合もその値が小さくなっている．また，いずれの正則化付き NMF においても，ノイズが大きい時には重み係数 µ 及び ν が大きいほどよい性能を示し，ノイズが小さい時には重み係数の値が小さい時によい性能を示している．これは，正則化の影響を大きくすることで，Fig. 4.3 及び Fig. 4.4 の基底行列または係数行列のような推定する傾向が強くなる．このような基底行列または係数行列は，観測行列の中に含まれる「低ランク性に寄与しない成分」を無視するような推定となっている．そのため，観測行列に含まれるノイズの影響をあまり受けずに正解行列の低ランク構造を抽出することができ，結果としてノイズに対する頑健性が向上しているといえる．

5.4 本章のまとめ
本章では，人工欠損行列を用いた補完実験により，それぞれの正則化付き NMF の行列補完の精度について検証した．いずれの正則化付き NMF においても改善が得られた一方で，補完対象がノイズ行列であるため，実用性が検証可能な実験条件であるとは言えない．そこで，次章では音響信号の振幅スペクトログラムに欠損を与え，これの補完実験によって補完性能の検証を行う．

5.4 本章のまとめ

37

Fig. 5.9. Evaluation criterion for each SNR conditions: : noise type is Rayleigh dist. and K = 10.

38

第6章

欠損音声スペクトログラムの補完
実験
6.1 はじめに本章では，音源分離等の処理を通して音響信号に発生しうる欠損値の補完について，正則化
付き NMF の補完精度の評価を行う．6.2 節では，NMF によって音声信号を扱う上で欠かせない STFT について説明する．また，本実験における欠損原因となる時間周波数マスクについて説明する．6.3 節では実験の詳細と各種条件について説明する．6.4 節では実験結果を示し，各正則化付き NMF の補完性能と特徴について述べる．6.5 節で本章についてまとめる．

6.2 実験背景
本実験では，音響信号の振幅スペクトログラムに発生する欠損について補完実験を行う．本節では，実験において，音声信号を非負値行列（振幅スペクトログラム）に変換するために用いる STFT と本実験が仮定する欠損発生メカニズムである時間周波数マスクについて説明する．

6.2.1 STFT
観測された音響信号を非負行列として表現する手法として，STFT がある．STFT とは

Fig. 6.1 に示すように，一次元の時間信号を二次元の時間周波数信号に変換する処理である．
STFT の分析窓関数の長さ及びシフト長をそれぞれ Q 及び τ としたとき，時間領域の信号
z(p)（p = 1, 2, · · · , P は離散時間サンプルのインデクスを表す）の j 番目の短時間区間（時間フレーム）の信号は次式で表される．
Tz (j) = [z ((j − 1)τ + 1) , z ((j − 1)τ + 2) , · · · , z ((j − 1)τ + Q)]hiT= z (j) (1), z (j) (2), · · · , z (j) (q), · · · , z (j) (Q)
∈ RQ

(6.1)

6.2 実験背景

39

Fig. 6.1. Process of STFT.

ここで，j = 1, 2, · · · , J 及び q = 1, 2, · · · , Q は，それぞれ時間フレーム及び時間フレーム内のサンプルのインデクスを示す．また，時間フレーム数 J は次式によって与えられる．

J=

P
τ

(6.2)

但し，信号長 L はセグメント数 J が整数となるように各時間フレームの信号の両端にゼロを挿入する処理（ゼロパディング）が施されている．式 (6.1) で定義される時間フレームの信号

を全ての j についてまとめた全時間フレームの信号を z = z (1) z (2) · · · z (J) ∈ RQ×J と表記すると，STFT の処理は次式のように表される．

Z = STFTω (z) ∈ CI×J

(6.3)

ここで，ω = [ω(1), ω(2), · · · , ω(Q)]T ∈ RQ は STFT で用いる窓関数である．スペクトログ
ラム Z の (i, j) 番目の要素は次式で表される．

zij =

X
q


ω(q)z

(j)

(q)exp−ι2π(q − 1)(i − 1)
F


(6.4)

ここで F は ⌊ F2 ⌋ + 1 = I を満たす整数（⌊·⌋ は床関数），i = 1, 2, · · · , I は周波数ビンのイン
デクス，ι は虚数単位を示している．このように，時間領域の信号に対して一定幅の短時間ごとに分析窓関数を乗じて離散フーリエ変換を行うことで，横軸が時間，縦軸が周波数のスペク
トログラムと呼ばれる複素行列 Z で表すことができる．
音響信号に NMF を適用する場合，時間信号に STFT を適用して得られるスペクトログラ
ム（複素行列）の振幅値（振幅スペクトログラム）やパワー値（パワースペクトログラム）を非負値観測行列 X ∈ RI×J
≥0 として扱うことが一般的である．振幅スぺクトログラム X を NMFで分解した様子が Fig. 6.2 である．振幅スペクトログラム X は NMF を用いることで基底行

40

第6章

欠損音声スペクトログラムの補完実験

Fig. 6.2. NMF for audio signals, where K = 2.

列 W と係数行列 H の行列積で近似されている．Fig. 6.2 において，基底行列 W は 2 つの基底ベクトルから成る（K = 2 に設定されている）．即ち，観測行列 X のスペクトログラムを 2 本のスペクトルパターン（基底ベクトル）で近似的に表現している．同様に，係数行列 H
は，W 中の 2 本の基底ベクトルがどの時間にどの程度の強さで生起するかを表す時間的強度変化を行ベクトルとして含んでいる．
このように，音響信号であっても，STFT によってスペクトログラムを得ることで NMF を適用することができる．NMF で抽出される非負構成要素は解釈が容易であり，音響信号中のスペクトルパターンを教師なし学習できるため，音楽信号解析 [18] や音源分離 [7, 10] 等に頻繁に適用される．本実験では，単に欠損のある音響信号を NMF に基づく行列補完の観測行列として適用するため，STFT によってパワースペクトログラムに変換する．

6.2.2 時間周波数マスク
次に，時間周波数マスクについて説明する．時間周波数マスクとは，観測信号のスペクトロ
グラムのある時間周波数要素に対して，目的の分離信号の成分が存在しているかどうかを表す
二次元行列である．ソフトマスクであれば，時間周波数マスクは 0 から 1 までの値で構成され，バイナリマスク行列であれば，0 又は 1 の値で構成される．この概要を Fig. 6.3 に示す．
図のような赤，緑，青の音源から成る混合信号から赤の音源成分のみを取り出したい状況を仮定する．この時，赤の音源の部分を 1 それ以外を 0 とするようなバイナリの時間周波数マスクを構成し，混合信号のスペクトログラムと時間周波数マスクを要素毎に乗算することで，赤の音源とその他の音源を分離することが可能である．
このような時間周波数マスクはそのまま音源分離に適用され，いかに目的音源のみを通すマ
スクを作成できるかが検討されている（例えば，[19, 20, 21] 等）．適切な時間周波数マスクを生成できれば，簡単に音源分離可能な一方で，特にバイナリ時間周波数マスクによって分離された音源は，スペクトログラムの中に 0 を含むことになる．即ち，元の目的音源に対し欠損した状態であり，これは，人間の聴覚には歪みとして現れる．
本実験では，バイナリ時間周波数マスクを適用して音源分離された音響信号に対し，歪みを

6.3 実験方法

41

Time-frequency maskfor extracting red source

Frequency [Hz]Frequency [Hz]Spectrogram of mixture signals

Time [s]Time [s]Black and white show zero and ones, respectively.

(a)Frequency [Hz]Extracted red source

Time [s](b)Fig. 6.3. Overview of time-frequency mask: (a) mixture signal and time-frequency binary mask
for red source and (b) extracted red source by time-frequency binary masking.

軽減することを目的として提案手法で欠損補完を施すことを想定する．バイナリ時間周波数マ
スクは欠損のある音源を作成するための，完全に分離された（混合前の）音源信号のパワース
ペクトログラムから生成する．これは即ち，前段のバイナリ時間周波数マスクに基づく音源分
離アルゴリズムが理想的な性能を達成した場合において，提案手法による欠損補完が分離音の品質をどの程度改善させることができるか，ということを比較・検証することに相当する．

6.3 実験方法
本実験では，音源に対しバイナリ時間周波数マスクをかけることで欠損を生成し，これを観測行列とした補完実験を行う．音源には SiSEC2011 [22] データセットから楽曲 bearlin-

roads のアコースティックギター音源（Gt.）及びピアノ音源（Pf.）の 2 つの音源を用いる．
Table 6.1 に使用音源を示す．SiSEC2011 はプロ奏者の音楽信号が収録されたデータセットであり，bearlin-roads には音源分離への活用を想定してトレーニングセットとテストセットの2 つのデータセットが用意されている．本実験において学習は不要のため，トレーニングセッ
トのみを用いる．バイナリ時間周波数マスクは，各音源のパワースペクトログラムの各要素の値から生成する．
今，アコースティックギター音源のパワースペクトログラム Z (Gt) ∈ RI×J
≥0 及びピアノ音源

42

第6章

欠損音声スペクトログラムの補完実験

Table 6.1. Source signals

Dataset nameInstrument typedev1 bearlin-roads snip 85 99 acoustic guit main.wav

Acoustic guitar

dev1 bearlin-roads snip 85 99 piano.wav

Piano

Table 6.2. STFT conditions

Audio length

14.0 sSampling rate

44.1 kHzWindow function in STFT

Blackman window

Window length in STFT

92.9 ms

Shift length in STFT

46.4 msのパワースペクトログラム Z (Pf) ∈ RI×J
≥0 が与えられた時，ピアノ音源のバイナリ時間周波数
マスク M (Gt) ∈ {0, 1}I×J 及びギター音源のバイナリ時間周波数マスク X (Pf) ∈ {1, 0}I×J
は次式によって得る．

(
1
0
(
1
(Pf)
mij =
0

(Gt)

(Gt)

ここで，zij

(Pf)

，zij

(Gt)

，mij

∀i, j

(6.5)

(Gt)

if mij = 0∀i, jotherwise(Pf)，及び mij(Pf)

if zij ≥ zij
otherwise

(Gt)
mij =

(6.6)

はそれぞれ Z (Gt) ，Z (Pf) ，M (Gt) ，及び M (Pf) の要素である．即ち，式 (6.5) 及び式 (6.6) は，2 つの観測信号のパワースペクトログラムの各要素について大小比較し，大きい要素を持つ観測信号のバイナリ時間周波数マスクを 1，もう一方を 0 として生成する．完全に分離された信号のスペクトログラム Z (Gt) 及び Z (Pf) を利用しているので，実際の音源分離では作ることができない理想的なバイナリ時間周波数マスクである．
欠損を含む観測行列は，観測信号のパワースペクトログラムにバイナリ時間周波数を掛け合
I×J
わせることにより生成する．アコースティックギター音源の観測行列を X (Gt) ∈ R≥0，ピア

ノ音源の観測行列を X (Pf) ∈ RI×J
≥0 とすると，次式より得られる．

X (Gt) = M (Gt) ⊙ Z (Gt)
X

(Pf)

=M

(Pf)

⊙Z

(Pf)

(6.7)
(6.8)

以上より，観測行列 X (Gt) 及び X (Pf) を得る．Fig. 6.4 に使用音源の時間信号とパワー
スペクトログラムを，Fig. 6.5 にバイナリ時間周波数マスクと観測行列をそれぞれ示す．また，この時の STFT の条件を Table 6.2 に示す．本実験では，観測行列 X (Gt) 及び X (Pf)

6.3 実験方法

43

Fig. 6.4. Souce signals and their power spectrograms Z (Gt) and Z (Pf) used in experiments.

を KLNMF 及び 2 つの正則化付き NMF でそれぞれ補完を行い，性能を評価する．バイナリ
時間周波数マスクによって得られる分離音源は波長の一部に 0 を含むことで歪みを有する．
行列補完により適切に欠損値が埋められた場合，このような歪みは減少すると考えられる．
また，補完時に目的音源のパターンを強く抽出することができれば，非目的音源の干渉を抑えることができる．従って，音源分離分野においてよく用いられる指標である信号対歪み比（source-to-distortion ratio: SDR），信号対干渉比（source-to-interference ratio: SIR），及び信号群対歪み比（sources-to-artifacts ratio: SAR）[23] の改善量によって性能を評価する．
これらは分離信号の品質を測る指標であり，SIR は値が高いほど非目的音源の干渉が小さいことを，SAR は値が高いほど信号に残留する人工歪みが小さいことをそれぞれ示す．また，

SDR は SIR と SAR の両指標を加味した総合的な音源分離性能を評価する指標である．
今，目的音源 s(p) における推定信号 ŝ(p) は以下の式に示すように分解することができる．

ŝ(p) = starget (p) + sinterf (p) + eartif (p)(6.9)

ここで，starget (p)，sinterf (p)，及び eartif (p) はそれぞれ推定信号 ŝ(p) 中の目的音源成分，残留した非目的音源成分，及び音源分離によって生じた人工的な歪み成分を表す．SDR は次式で表される．

P
SAR = 10 log10 P

p |starget (p)|22
p |sinterf (p) + eartif (p)|[dB]

(6.10)

従って，高い SDR 値を達成するには sinterf (p) 及び eartif (p) が少なく，starget (p) が高精度に

44

第6章

欠損音声スペクトログラムの補完実験

Fig. 6.5. Binary time-frequency masks and separated power spectrograms, which are used as
observed matrix in matrix completion experiment. Black elements in observed matrix
show missing values.

推定されている必要がある．SIR は，starget (p) 及び sinterf (p) を用いて次式で表される．

P

2
p |starget (p)|SIR = 10 log10 P
2
p |sinterf (p)|[dB]

(6.11)

従って，高い SIR 値を達成するには sinterf (p) をより小さくする必要がある．SAR は，starget (p)，sinterf (p)，及び eartif (p) を用いて次式で表される．
P
SAR = 10 log102
p |starget (p) + sinterf (p)|P2
l |eartif (p)|[dB]

(6.12)

従って，高い SAR 値を達成するには目的音源及び非目的音源によって構成される分離音源の大きさに対し，歪み成分の大きさをより小さくする必要がある．以上より，これら 3 つの指標を確認することで補完性能について評価することができる．
実験手順と各処理で使用するパラメータを Fig. 6.6 にまとめる．本実験では，使用音源をパワースペクトログラムに変換した後，バイナリ時間周波数マスクによって観測行列 X (Gt) 及6.3 実験方法

45

Fig. 6.6. Experimental conditions changed in evaluation.

Fig. 6.7. Cumulative singular values of power spectrogram of acoustic guitar and piano sources.

び X (Pf) を生成する．観測行列に対し各設定値の NMF を適用し行列補完を行い，時間領域に変換した補完信号と観測信号を比較し，補完精度を評価する．

NMF の基底ベクトル数 K の設定値候補を算出するために，Z (Gt) 及び Z (Pf) がどの程度のランクで良く近似できるかを確認する．具体的には，累積特異値が 80% 程度に到達するのに必要な特異ベクトルの個数を確認し，その周辺の値を NMF の基底ベクトル数 K の設定値候補とする．実際に Z (Gt) 及び Z (Pf) の特異値分解により累積特異値を求めた結果を Fig. 6.7に示す．これより，いずれの音源においてもランク 50 程度で累積特異値の割合が全体の 80%を超えることがわかる．また，観測行列は 2049 行 302 列の行列であり，NMF は低ランク近似であることから，基底ベクトル数 K の設定値は 4 から 100 までとした．正則化付き NMFの重み係数 µ 及び ν は 10−3 から 105 まで対数スケールで 40 点推移させた．各実験条件において NMF の基底行列及び係数行列の初期化に用いる乱数シードを 10 種類用意し，その平均を評価値とした．

46

第6章

欠損音声スペクトログラムの補完実験

6.4 実験結果
Figs. 6.8–6.10 に，基底ベクトル数 K 毎の各 NMF による補完信号の SDR 改善量 (SDRi)，
SIR 改善量 (SIRi)，及び SAR 改善量 (SDRi) の推移をそれぞれ示す．これらの改善量は，バ
イナリ時間周波数マスクによって欠損したパワースペクトログラムをそのまま（補完せずに）時間領域へ戻した信号の SDR 値，SIR 値及び SAR 値との差であり，値が高いほど行列補完
処理によって品質が改善されたことを示す．

Fig. 6.8 より，両音源において CSRNMF 及び KLNMF は SDR が無補完信号に対し低下しており，分離信号として劣化したことがわかる．KLNMF と CSRNMF を比較すると，基底ベクトル数 K 値によって前後するものの同程度の値をとっているといえる．一方で，L1 -RNMF
は SDRi が常にほぼ 0 となっており，僅かな改善もみられる．また，KLNMF や CSRNMFと異なり，基底ベクトル数 K の変化に対して一定の値をとっている．

Fig. 6.9 より，アコースティックギター音源において SIR は，KLNMF 及び CSRNMF 共に無補完信号に対して劣化はあるものの，KLNMF が平均的に 25 dB 程度低下したのに対し

CSRNMF は 8 dB 程度の低下に抑えられている．一方で，ピアノ音源については CSRNMFが無補完信号に対して広い基底ベクトル数 K 値において改善を示しており，最大で 10 dB 程度まで改善した．また，L1 -RNMF は SDR 同様ほぼ 0 の値を推移しており，ピアノ音源において極僅かな改善量が見られる．

Fig. 6.10 より，アコースティックギター音源において SAR は，KLNMF 及び CSRNMF
共に改善が見られず，基底ベクトル数 K 値が大きくなるにつれていずれの両手法とも減少していることがわかる．L1 -RNMF は基底ベクトル数 K が低い時に僅かな改善を示し，基底ベ
クトル数 K が大きくなるにつれ，KLNMF 及び CSRNMF 同様に SAR が低下していることがわかる．一方で，ピアノ音源においては，KLNMF 及び CSRNMF が低下したのに対しL1 -RNMF は広い基底ベクトル数 K 値の範囲の多くの点で 1 dB 程度の改善を示した．
以上より，L1 -RNMF は SDR，SIR 及び SAR において，音源によらず 0 に近い値をとっており，一部において僅かな改善が見られた．L1 -RNMF において特徴的なのは，ピアノ音源
における SAR が 1 dB 程度上昇したことである，ピアノ音源においてのみこのような結果を示したことから，本実験におけるアコースティックギター音源とピアノ音源の特徴の違いが，

L1 -RNMF の性能に影響したと考えられる．Fig. 6.4 のスペクトログラムにおいて特に強い振幅を持つことを示す黄色い要素が，Fig. 6.5 においてアコースティックギター音源は多く欠損したのに対し，ピアノ音源は保持していることがわかる．これより，ピアノ音源のスペクトロ
グラムは重要な特徴を残した行列であるといえる．このような行列のため，欠損値をピアノ音源に矛盾なく補完することができ，結果として SAR が改善したと考えられる．CSRNMF は

SDR 及び SAR が，多くの条件において KLNMF と同等か低い結果となった．一方で，SIR
においては KLNMF に対し明らかに性能が向上しており，ピアノ音源においては無補完信号に対し 10 dB 程度の改善を示した．このような結果の理由として，L1 -RNMF 同様ピアノ音

6.4 実験結果

47

Fig. 6.8. Average SDRi obtained by matrix completion based on KLNMF, CSRNMF, and L1 RNMF when µ = 0.001, and ν = 42.813.

Fig. 6.9. Average SIRi obtained by matrix completion based on KLNMF, CSRNMF, and L1 RNMF when µ = 61585, and ν = 100000.

源のスペクトログラムが重要な特徴を保持しているためであると考えられる．無補完信号における干渉信号とはバイナリ時間周波数マスクによる 0 の信号である．CSRNMF によってこれを補完するにあたり，より元の信号の特徴を含んだピアノ音源において性能が向上したと考えられる．また，SDR が全体的に低下している要因として，バイナリ時間周波数マスクによる音源分離が理想的状態が大きい．欠損値を補完する値が元音源と異なるとき，その誤差が人工的歪みや非目的音源の干渉となってしまい，音源分離の総合的指標となる SDR は改善しにくいと言える．
最後に，各手法について，SDR，SIR 及び SAR から総合的に評価する．Fig. 6.11 に，基底

48

第6章

欠損音声スペクトログラムの補完実験

Fig. 6.10. Average SARi obtained by matrix completion based on KLNMF, CSRNMF, and L1 RNMF when µ = 100000, and ν = 3.793.

Fig. 6.11. Average SDR, SIR, and SAR values when K = 30, µ = 6.159, and ν = 297.635.

ベクトル数 K が 30 の時における無補完信号を含めた各手法の最大値を示す．また，その他の基底ベクトル数 K における各手法の最大値を付録 C に示す．Fig. 6.11 において，CSRNMF
及び L1 -RNMF の重み係数はそれぞれ µ = 6.159，ν = 297.635 である．また，この時の各値
を Table 6.3 に示す．これより，KLNMF は全ての指標において無補完信号から低下する一方，CSRNMF は SIR が，L1 -RNMF は SDR 及び SAR の改善が見られ，補完手法として性能向上が確認できる．総合的な音源補完の性能としては，L1 -RNMF が最も安定していると言え，SIR を最も重視する場合には CSRNMF が有効であると考えられる．いずれの場合においても，正則化項による効果が確認でき，補完性能の向上が得られた．

6.5 本章のまとめ

49

Table 6.3. Average SDR, SAR, and SIR values when K = 30, µ = 6.159, and ν = 297.635.

Acourstic guitar

Piano

MethodSDR [dB]

SIR [dB]

SAR [dB]

SDR [dB]

SIR [dB]

SAR [dB]

Without
completion

1.612

22.9001.667

13.78417.22416.485KLNMF-14.645-0.305-9.508-12.20212.534-11.901

CSRNMF-11.366

26.762-11.349-20.223

16.670-20.098

L1 -RNMF

1.682

21.5591.757

13.79617.05216.658

6.5 本章のまとめ
本章では，実際の音源を用いた補完実験を行い，それぞれの手法における性能を評価した．

SCRNMF による行列補完は SDR が改善しなかったものの，SIR は大きく改善する場合が確認できた．L1 -RNMF による行列補完は SDR，SIR 及び SAR において僅かな改善が見られ，

KLNMF 及び CSRNMF に対し最も高い補完性能が得られた．これより，正則化項が働くことで補完性能に一定の改善を得られることが確認できた．

50

第7章

結言
本論文では，あらゆるデータにおいて発生する欠損値を低コストかつ高精度に補完するため，2 つの正則化付き NMF を検討した．コサイン類似度正則化は NMF に付与することにより，基底行列に含まれる各基底ベクトルを相互に直交させることが確認できた．L1 ノルム
正則化は NMF に付与することにより，係数行列の要素がスパースになるよう促すことが確認できた．人工欠損行列を用いた実験において，CSRNMF が特にノイズの大きい条件について従来の KLNMF を用いた手法より高精度な補完性能を発揮した．また，L1 -RNMF は多くの条件において，KLNMF に対し補完性能が向上しており，基底ベクトル数 K のモデル
誤差及びノイズに対する頑健性が確認できた．データセットの音源を用いた実験においては，

L1 -RNMF が安定した補完性能を見せたものの劇的な音源の改善は見られなかった．また，
CSRNMF は特定の条件において SIR のみ 10 dB 程度の大きな改善が確認できたが，これに伴い SDR 及び SAR は低下した．
本論文では，正則化を NMF に適用することで行列補完性能を向上させることができることを確認した．一方で，正則化項の影響度を決める重み係数によってその性能は大きく変化し，各実験において条件ごとに最適な値を示す重み係数は異なっていた．現状，最適な値は網羅的に重み係数値を変更して調べる必要があり，多大な時間を要する．この問題を解決するため，効率的な最適な重み係数のチューニング方法の検討が必要となる．

51

謝辞
本論文は，香川高等専門学校電気情報工学科北村研究室にて行われた研究に基づくものです．
初めに，本研究を進めるにあたり，ご多忙のところ多大なご指導を賜った，北村大地講師に深く感謝いたします．北村大地講師にはこの 1 年間で，本研究分野に関する基礎的な知識をはじめ，資料作成，実験結果の考察など細部に至るまで熱心にご指導下さいました．また，自分の進路についても非常に丁寧にサポート頂き，本当にありがとうございました．
本論文の副査を担当頂いた柿元健准教授には，論文の構成や記述に関して大変有益な助言を頂き，大変お世話になりました．ここに厚く御礼申し上げます．
北村研究室の先輩である川口翔也氏，溝渕悠朔氏，村田佳斗氏，綾野翔馬氏にはソフトウェ
アの使い方やプログラミングの記法，実験方法の助言など，研究に必要な多くのことを支援頂きました．研究について高い技術と熱意を持ちながら，研究室内をユーモアあふれる空間にして下さったことで，臆することなく様々な学びを得ることができました．また，綾野翔馬氏には本論文執筆において様々な助言を頂き，自身の研究についてより理解を深めることができました．心より感謝致します．
北村研究室の同期である加藤大輝氏，鈴木慶氏，松本愛花氏には，日頃のディスカッション
をはじめ，研究室生活を楽しいものにして頂きました．特に，加藤大輝氏とは行動を共にすることも多く，その中で多くのインスピレーションを与えていただきました，モティベーションの高い同期に囲まれたこの環境だからこと，研究を進めることができました．
最後に，実家を離れ香川高専へ進学することを快く後押ししてくれた家族には感謝の念に堪えません．本科 5 年間においても，安心して寮生活を送れるよう暖かく見守って下さり，本当にありがとうございました．重ねて，本研究に協力して下さった皆様には心より感謝申し上げ
ます．

52

参考文献
[1] L. Ren, T. Wang. A. S. Seklouli, H. Zhang, and A. Bouras, “A review on missing
values for main challenges and methods,” Information Systems, vol. 119, 2023.
[2] T. Emmanuel, T. Maupong, D. Mpoeleng, T. Semong, B. Mphago, and O. Tabona,“A survey on missing data in machine learning,” Journal of Big Data, vol. 8, 2021.
[3] D. D. Lee and H. S. Seung, “Learning the parts of objects by non-negative matrix
factorization,” Nature, vol. 401, no. 6755, pp. 788–791, 1999.
[4] D. D. Lee and H. S. Seung, “Algorithms for non-negative matrix factorization,” inProc. Neural Info. Process. Syst., pp. 556–562, 2000.
[5] T. Kujiraoka, S. Fumiaki, and I. Syohei, “Extraction of customer satisfaction topics
regarding product delivery using non-negative matrix factorization,” in Proc. IEEE
International Conference on Industrial Engineering and Engineering Management,
pp. 225–229, 2017.
[6] D. Kitamura, H. Saruwatari, H. Kameoka, Y. Takahashi, K. Kondo, and S. Nakamura, “Multichannel signal separation combining directional clustering and nonnegative matrix factorization with spectrogram restoration,” IEEE/ACM Trans. Audio,
Speech, Lang. Process., vol. 23, no. 4, pp. 654–669, 2015.
[7] D. FitzGerald, M. Cranitch, and E. Coyle, “On the use of the beta divergence formusical source separation,” in Proc. Irish Signals Syst. Conf., 2009.
[8] D. R. Hunter and K. Lange, “Quantile regression via an MM algorithm,” J. Comput.
Graph. Stat., vol. 9, no. 1, pp. 60–77, 2000.
[9] P. Smaragdis, B. Raj, and M. Shashanka, “Supervised and semi-supervised separationof sounds from single-channel mixtures,” in Proc. LVA/ICA, pp. 414–421, 2007.
[10] D. Kitamura, H. Saruwatari, K. Yagi, K. Shikano, Y. Takahashi, and K. Kondo,“Music signal separation based on supervised nonnegative matrix factorization withorthogonality and maximum-divergence penalties,” IEICE Trans. Fundamentals,
vol. E97–A, no. 5, pp. 1113–1118, 2014.
[11] Y. Iwase and D. Kitamura, “Supervised audio source separation based on nonnegative
matrix factorization with cosine similarity penalty,” IEICE Trans. Fundamentals,
vol. E105–A, no. 6, pp. 906–913, 2022.

参考文献

53

[12] P. O. Hoyer, “Non-negative matrix factorization with sparseness constraints,” The
Journal of Machine Learning Research 5, pp. 1457–1469, 2004.
[13] W. Liu, N. Zheng, and X. Lu, “Non-negative matrix factorization for visual coding,”
in Proc. ICASSP, pp. III–293, 2003.
[14] J. Le Roux, F. Weninger and J. R. Hershey, “Sparse NMF - half-baked or well done?,”Mitsubishi Electric Research Lab. Technical Report, TR2015–023, 2015.
[15] R. Tibshirani, “Regression shrinkage and selection via the Lasso,” J. Royal Statist.
Soc. Ser. B 58, pp. 267–288, 1996.
[16] M. A. Rasmussen, and R. Bro, “A tutorial on the Lasso approach to sparse modeling,”
Chemometrics and Intelligent Laboratory Systems, vol. 119, pp. 21–31, 2012.
[17] A. Marmin, J. Goulart, and C. Févotte, “Majorization-minimization for Sparse Nonnegative Matrix Factorization with the β-divergence.” IEEE Trans. Signal Process.,
pp. 1435–1447, 2023.
[18] C. Févotte, N. Bertin, and J.-L. Durrieu, “Nonnegative matrix factorization with
the Itakura-Saito divergence. With application to music analysis,” Neural Comput.,vol. 21, no. 3, pp. 793–830, 2009.
[19] H. Sawada, S. Araki, and S. Makino, “Underdetermined convolutive blind source separation via frequency bin-wise clustering and permutation alignment,” IEEE Trans.
Audio, Speech, Lang. Process., vol. 19, no. 3, pp. 516–527, 2011.
[20] Y. Luo and N. Mesgarani, “Conv-TasNet: Surpassing ideal time–frequency magnitude
masking for speech separation,” IEEE/ACM Trans. Audio, Speech, Lang. Process.,
vol. 27, no. 8, pp. 1256–1266, 2019.
[21] S. Oyabu, D. Kitamura, and K. Yatabe, “Linear multichannel blind source separation
based on time-frequency mask obtained by harmonic/percussive sound separation,”in Proc. ICASSP, pp. 201–205, 2021.
[22] S. Araki, F. Nesta, E. Vincent, Z. Koldovský, G. Nolte, A. Ziehe, and A. Benichoux,“The 2011 signal separation evaluation campaign (SiSEC2011) - audio source separation -,” in Proc. LVA/ICA, pp. 414–422, 2012.
[23] E. Vincent, R. Gribonval, and C. Févotte, “Performance measurement in blind audio source separation,” IEEE Trans. Audio, Speech, Lang. Process., vol. 14, no. 4,
pp. 1462–1469, 2006.

54

付録 A

補助関数法で利用される不等式
A.1 Jensen の不等式
補題 A.1. （Jensen の不等式）αi > 0 を，Pi αi = 1 を満たす補助変数とする．関数 f (x)が凸関数であるとき，xi (i = 1, . . . , I) に対して以下の不等式が成立する．

f

I
X
i=1

!
α i xi≤IX

αi f (xi )

(A.1)

i=1

f (x) が狭義凸関数であるとき，不等式中の等号が成立するための条件は x1 = · · · = xi =· · · = xI である．

55

付録 B

人工欠損行列の補完実験結果
本付録では，5 章で行った人工欠損行列の補完実験の内，本文で示した条件以外の実験結果を掲載する．Figs. B.1–B.11 にノイズに一様分布を用いた時の各 SNR 値における実験結果を示す．また，Figs. B.12–B.19 にノイズにレイリー分布を用いた時の各 SNR 値における実験結果を示す．

Fig. B.1. Improvement value obtained by CSRNMF: noise type is uniform dist., r = 0 dB,
and µ = 1274.28. The smaller the value and the redder the color, the better the
performance.

56

付録 B 人工欠損行列の補完実験結果

Fig. B.2. Improvement value obtained by CSRNMF: noise type is uniform dist., r = 3.16 dB,
and µ = 5455.59.

Fig. B.3. Improvement value obtained by CSRNMF: noise type is uniform dist., r = 5.99 dB,
and µ = 6.16.

57

Fig. B.4. Improvement value obtained by CSRNMF: noise type is uniform dist., r = 21.54 dB,
and µ = 0.88.

Fig. B.5. Improvement value obtained by CSRNMF: noise type is uniform dist., r = 40.84 dB,
and µ = 0.55.

58

付録 B 人工欠損行列の補完実験結果

Fig. B.6. Improvement value obtained by CSRNMF: noise type is uniform dist., r = 77.43 dB,
and µ = 0.55.

Fig. B.7. Improvement value obtained by L1RNMF: noise type is uniform dist., r = 0 dB, andν = 1.44.

59

Fig. B.8. Improvement value obtained by L1 RNMF: noise type is uniform dist., r = 3.16 dB,
and ν = 0.89.

Fig. B.9. Improvement value obtained by L1 -RNMF: noise type is uniform dist., r = 5.99 dB,
and ν = 0.55.

60

付録 B 人工欠損行列の補完実験結果

Fig. B.10. Improvement value obtained by L1 -RNMF: noise type is uniform dist., r = 21.54 dB,
and ν = 0.08.

Fig. B.11. Improvement value obtained by L1 -RNMF: noise type is uniform dist., r = 40.84 dB,
and ν = 0.01.

61

Fig. B.12. Improvement value obtained by CSRNMF: noise type is Rayleigh dist., r = 0 dB,
and µ = 61584.82.

Fig. B.13. Improvement value obtained by CSRNMF: noise type is Rayleigh dist., r = 3.16 dB,
and µ = 784.76.

62

付録 B 人工欠損行列の補完実験結果

Fig. B.14. Improvement value obtained by CSRNMF: noise type is Rayleigh dist., r = 40.84 dB,
and µ = 0.21.

Fig. B.15. Improvement value obtained by CSRNMF: noise type is Rayleigh dist., r = 77.43 dB,
and µ = 0.21.

63

Fig. B.16. Improvement value obtained by L1 -RNMF: noise type is Rayleigh dist., r = 0 dB,
and ν = 2.34.

Fig. B.17. Improvement value obtained by L1 -RNMF: noise type is Rayleigh dist., r = 3.16 dB,
and ν = 0.89.

64

付録 B 人工欠損行列の補完実験結果

Fig. B.18. Improvement value obtained by L1 -RNMF: noise type is Rayleigh dist., r = 40.84 dB,
and ν = 0.01.

Fig. B.19. Improvement value obtained by L1 -RNMF: noise type is Rayleigh dist., r = 77.43 dB,
and ν = 0.01.

65

付録 C

欠損音声スペクトログラムの補完実
験結果
本付録では，6 章で行った欠損音声スペクトログラムの補完実験の内，本文で示した条件以外の実験結果を掲載する．Figs. C.1–C.17 に，各基底ベクトル数 K 値における各手法の最大値を示す．

Fig. C.1. Average SDR, SIR, and SAR values when K = 4, µ = 0.001, and ν = 23357.214.

66

付録 C

欠損音声スペクトログラムの補完実験結果

Fig. C.2. Average SDR, SIR, and SAR values when K = 10, µ = 0.004, and ν = 8858.668.

Fig. C.3. Average SDR, SIR, and SAR values when K = 16, µ = 0.001, and ν = 483.293.

67

Fig. C.4. Average SDR, SIR, and SAR values when K = 22, µ = 10, and ν = 784.760.

Fig. C.5. Average SDR, SIR, and SAR values when K = 28, µ = 3.793, and ν = 784.760.

68

付録 C

欠損音声スペクトログラムの補完実験結果

Fig. C.6. Average SDR, SIR, and SAR values when K = 28, µ = 2.335, and ν = 297.635.

Fig. C.7. Average SDR, SIR, and SAR values when K = 40, µ = 3.703, and ν = 784.760.

69

Fig. C.8. Average SDR, SIR, and SAR values when K = 46, µ = 3.793, and ν = 483.293.

Fig. C.9. Average SDR, SIR, and SAR values when K = 52, µ = 0.001, and ν = 297.635.

70

付録 C

欠損音声スペクトログラムの補完実験結果

Fig. C.10. Average SDR, SIR, and SAR values when K = 58, µ = 1.438, and ν = 297.635.

Fig. C.11. Average SDR, SIR, and SAR values when K = 64, µ = 2.336, and ν = 183.298.

71

Fig. C.12. Average SDR, SIR, and SAR values when K = 70, µ = 0.336, and ν = 297.635.

Fig. C.13. Average SDR, SIR, and SAR values when K = 76, µ = 0.207, and ν = 483.293.

72

付録 C

欠損音声スペクトログラムの補完実験結果

Fig. C.14. Average SDR, SIR, and SAR values when K = 82, µ = 0.546, and ν = 183.298.

Fig. C.15. Average SDR, SIR, and SAR values when K = 88, µ = 0.546, and ν = 183.298.

73

Fig. C.16. Average SDR, SIR, and SAR values when K = 94, µ = 0.207, and ν = 183.298.

Fig. C.17. Average SDR, SIR, and SAR values when K = 100, µ = 0.127, and ν = 183.298.