非負値行列因子分解

非負値行列因子分解（ひふちぎょうれついんしぶんかい、英: Non-negative matrix factorization、略称：NMFまたはNNMF）、あるいは非負値行列近似^[1][2]は、多変量解析および線形代数における一群のアルゴリズムであり、行列 V を、通常2つの行列 W と H に分解するものである。このとき、3つの行列すべての要素が負でない（非負）という制約がある。

近似的非負値行列因子分解の図：行列 V は、2つの小さな行列 W と H に分解される。これらを乗算することで、V をおおよそ再構成できる。

この非負性制約により、得られた行列は直感的に解釈しやすくなる。また、音声スペクトログラムや筋活動などのデータでは、そもそも非負性が前提となっている。一般には問題を厳密に解けないため、通常は数値的に近似解を求める。

NMFは、天文学^[3][4]、コンピュータビジョン、文書クラスタリング^[1]、欠損データ補完^[5]、計量化学、音声信号処理、レコメンダシステム^[6]、およびバイオインフォマティクス^[7]など、さまざまな分野に応用されている。

歴史

計量化学において、非負値行列因子分解（NMF）は「自己モデリング曲線分解（Self Modeling Curve Resolution）」という名称で古くから知られていた^[8]。 この枠組みにおいては、右側の行列のベクトルは離散ベクトルではなく連続的な曲線として表現される。

また、1990年代にはフィンランドの研究グループによって「正値行列因子分解（Positive Matrix Factorization）」という名称で初期の研究が進められた^[9][10][11]。

「非負値行列因子分解（Non-negative Matrix Factorization）」という名称は、ダニエル・リーとセバスチャン・スンがアルゴリズムの性質を研究し、2種類の因子分解手法に対する単純かつ有用なアルゴリズムを発表した1999年以降、広く知られるようになった^[12][13]。

背景

行列 V を W と H の積として表すとする：

${\displaystyle \mathbf {V} =\mathbf {W} \mathbf {H} \,.}$

行列の積は、V の各列ベクトルを、W の列ベクトルの線形結合として計算することで実装できる。このとき、係数は H の列によって与えられる。すなわち、V の各列は以下のように計算できる：

${\displaystyle \mathbf {v} _{i}=\mathbf {W} \mathbf {h} _{i}\,,}$

ここで、v_i は V の第 i 列ベクトル、h_i は H の第 i 列ベクトルである。

行列の積では、因子行列の次元が積の行列（V）よりも大幅に小さくできる場合があり、この特性がNMFの基礎となっている。NMFは、元の行列と比べて大幅に次元が削減された因子を生成する。たとえば、V が m × n の行列、W が m × p の行列、H が p × n の行列であれば、p は m および n よりもかなり小さくすることができる。

以下にテキストマイニングに基づく例を示す：

• 入力行列 V は 10000 行 × 500 列の行列であり、行に単語、列に文書が対応する。すなわち、10000語によってインデックス付けされた500の文書がある。
• アルゴリズムに対して、10個の特徴を見つけるように指定したとすると、「特徴行列」W は 10000 行 × 10 列、「係数行列」H は 10 行 × 500 列の行列となる。
• W と H の積は 10000 行 × 500 列の行列となり、これは元の入力行列 V と同じサイズであり、うまく因子分解ができていれば、V の妥当な近似になっている。
• 先述の行列積の説明から、WH の各列は、W の10個の列ベクトルの線形結合として構成されることが分かる。係数は H の列によって与えられる。

この最後の点がNMFの本質である。元の文書（入力行列の各列ベクトル）を、少数の隠れた特徴（W の列ベクトル）によって構成されたものと考えることができる。NMFはこれらの特徴を自動的に抽出する。

特徴行列 W の各特徴（列ベクトル）は、単語の集合から成る「文書の原型（アーキタイプ）」と考えることができ、各単語のセル値は、その単語がその特徴においてどれだけ重要か（ランク）を示す。係数行列 H の各列は元の文書を表し、各セル値はその文書が各特徴に対してどれだけ関連しているか（ランク）を表す。元の文書は、W の各特徴ベクトルを、H の各セル値で重み付けした線形結合として再構成できる。

クラスタリング特性

NMF（非負値行列因子分解）には、内在的にクラスタリングの特性があるとされる^[14]。 すなわち、入力データの列ベクトル ${\displaystyle \mathbf {V} =(v_{1},\dots ,v_{n})}$ を自動的にクラスタリングする性質を持つ。

より具体的には、${\displaystyle \mathbf {V} \simeq \mathbf {W} \mathbf {H} }$ によって ${\displaystyle \mathbf {V} }$ を近似するために、${\displaystyle \mathbf {W} }$ および ${\displaystyle \mathbf {H} }$ を誤差関数（フロベニウスノルム）の最小化によって求めると

${\displaystyle \left\|V-WH\right\|_{F}}$

ここで、制約条件 ${\displaystyle W\geq 0,H\geq 0}$ を満たす。

さらに、${\displaystyle \mathbf {H} \mathbf {H} ^{T}=I}$（すなわち ${\displaystyle \mathbf {H} }$ の直交性）という制約を追加すると、この最小化問題は数学的に K-meansクラスタリング の最小化問題と等価になる^[14]。

このとき、求められた ${\displaystyle \mathbf {H} }$ はクラスタの所属情報を示し、たとえば ${\displaystyle \mathbf {H} _{kj}>\mathbf {H} _{ij}}$ がすべての i ≠ k に対して成り立つ場合、データ ${\displaystyle v_{j}}$ は第 ${\displaystyle k}$ クラスタに属すると考えられる。また、${\displaystyle \mathbf {W} }$ はクラスタの重心（セントロイド）を示し、第 ${\displaystyle k}$ 列が第 ${\displaystyle k}$ クラスタの重心を与える。このクラスタ重心の表現は、凸NMF（Convex NMF）によってさらに明確になる。

この直交制約 ${\displaystyle \mathbf {H} \mathbf {H} ^{T}=I}$ を明示的に課さない場合でも、${\displaystyle \mathbf {H} }$ はある程度の直交性を自然に持ち、クラスタリング特性は依然として維持される。

また、誤差関数として カルバック・ライブラー情報量を使用した場合、NMFは確率的潜在意味解析（PLSA）と等価になる。PLSA は、文書クラスタリングにおいて広く用いられている手法である^[15]。

種類

近似的非負値行列因子分解

通常、NMF における W の列数と H の行数は、WH が V に近似するように選ばれる。完全な分解は、W と H、および残差行列 U によって表され、

${\displaystyle \mathbf {V} =\mathbf {WH} +\mathbf {U} }$

となる。残差行列 U の要素は負の値または正の値をとる可能性がある。

W および H のサイズが V より小さい場合、それらの格納と操作が容易になる。また、V の各要素をより少ないデータで表現したい場合には、潜在的な構造を推定する必要があり、そのために小さい行列への因子分解が有効である。

凸非負値行列因子分解 (Convex NMF)

標準的な NMF では、因子 W ∈ R₊^{m × k} は任意の非負値をとる。凸NMF^[16] では、W の列が入力データベクトル ${\displaystyle (v_{1},\dots ,v_{n})}$ の凸結合となるよう制限する。これにより、W によるデータ表現の質が大幅に向上し、H はより疎で直交的になる傾向がある。

非負ランク因子分解

もし V の非負ランクがその通常のランクと等しい場合、V = WH は非負ランク因子分解（NRF）と呼ばれる^[17][18][19]。

V の NRF（存在する場合）を求める問題は、NP困難であることが知られている^[20]。

異なるコスト関数と正則化

非負値行列因子分解にはさまざまなタイプがあり、それらはV と WH の間の差異（損失関数）を測る方法や、W および H に対する正則化により区別される^[1]。

リーとスンによって研究された2つの単純な誤差関数は、二乗誤差（フロベニウスノルム）と、正の行列に対するカルバック・ライブラー情報量の拡張である。各誤差関数は、異なる反復更新規則に基づく異なる NMF アルゴリズムを導く。

二乗誤差を用いた NMF の場合、最小化する目的関数は以下のように表される：

${\displaystyle F(\mathbf {W} ,\mathbf {H} )=\left\|\mathbf {V} -\mathbf {WH} \right\|_{F}^{2}}$

画像データに対する別の種類の NMF は、全変動ノルムに基づいている^[21]。

また、L1正則化（ラッソ回帰に類似）を追加した NMF は、「非負値スパースコーディング」とも呼ばれ、スパースコーディング問題に類似している^[22][23]。 ただし、これも依然として NMF と呼ばれることがある^[24]。

オンラインNMF

多くの標準的なNMFアルゴリズムは、全データを一括で処理する前提で設計されており、すなわち、全体の行列が初めから利用可能である必要がある。しかし、データが膨大でメモリに収まりきらない場合や、データがデータストリームとして逐次的に提供される場合、このような手法は不都合である。例えば、レコメンダシステムにおける協調フィルタリングでは、ユーザー数やアイテム数が非常に多く、1ユーザーまたは1アイテムが追加されるたびに再計算するのは非効率である。これらのケースでは、標準のNMFとは異なるアルゴリズムが必要とされる^[25][26]。

畳み込みNMF

もし V の列が時間や空間といった次元に沿ってサンプリングされたデータ（例：音声信号、画像、動画）を表す場合、これらの次元に対する平行移動不変性を持つ特徴を学習するために「畳み込みNMF」が利用される。この場合、W はスパースであり、その各列はVの空間・時間次元に沿った平行移動によって共有される畳み込みカーネルを表す。Hの空間的・時間的プーリングを行い、その結果を再び入力として畳み込みNMFに使用することで、深い特徴階層を学習できる^[27]。

アルゴリズム

W と H を求める方法にはいくつかの選択肢がある。LeeとSeungによる乗法更新則^[13] は、実装の容易さから人気のある手法である。このアルゴリズムは以下のようになる：

初期化：W および H を非負値で初期化する。

次に、反復ごとに以下を計算して n を反復のインデックスとする：

${\displaystyle \mathbf {H} _{[i,j]}^{n+1}\leftarrow \mathbf {H} _{[i,j]}^{n}{\frac {((\mathbf {W} ^{n})^{T}\mathbf {V} )_{[i,j]}}{((\mathbf {W} ^{n})^{T}\mathbf {W} ^{n}\mathbf {H} ^{n})_{[i,j]}}}}$

および

${\displaystyle \mathbf {W} _{[i,j]}^{n+1}\leftarrow \mathbf {W} _{[i,j]}^{n}{\frac {(\mathbf {V} (\mathbf {H} ^{n+1})^{T})_{[i,j]}}{(\mathbf {W} ^{n}\mathbf {H} ^{n+1}(\mathbf {H} ^{n+1})^{T})_{[i,j]}}}}$

W と H が収束するまで繰り返す。

この更新は行列全体ではなく要素単位で実施される。

乗法更新における W および H の更新因子、すなわち ${\displaystyle {\frac {\mathbf {W} ^{\mathsf {T}}\mathbf {V} }{\mathbf {W} ^{\mathsf {T}}\mathbf {W} \mathbf {H} }}}$ および ${\displaystyle {\textstyle {\frac {\mathbf {V} \mathbf {H} ^{\mathsf {T}}}{\mathbf {W} \mathbf {H} \mathbf {H} ^{\mathsf {T}}}}}}$ は、${\displaystyle \mathbf {V} =\mathbf {W} \mathbf {H} }$ となるときはすべて1の行列になる。

最近では他のアルゴリズムも提案されている。いくつかの手法は、非負最小二乗法に基づく交互最適化に依拠しており、1ステップごとに H を固定して W を非負最小二乗法で解き、次に W を固定して H を解く。このとき W と H の最適化方法は同じでも異なっていてもよい。他の手法としては、射影最急降下法^[28]、有効制約法^[29]、最適勾配法、ブロック主ピボット法^[30] などが挙げられる。

逐次NMF

NMFの構成要素（W および H）を逐次構築する手法は、天文学における主成分分析（principal component analysis; PCA）との関係を示すために初めて使用された^[31]。PCAでは、構成要素の寄与は対応する固有値の大きさによってランク付けされる。一方、NMFでは構成要素が1つずつ順次に構築されるため、経験的に寄与をランク付けすることができる。

逐次NMFの構成要素の寄与は、カルフーネン・ロエヴェ展開（PCAの応用）と比較でき、これは固有値のプロットによって示される。PCAでは、構成要素の選択において典型的な方法は「エルボー」に基づいており、固有値プロットにおいて平坦な部分が存在すると、PCAがデータを効率的に捉えられていないことを示す。最後に急激な減少が生じる部分では、ランダムノイズを捉えてしまい、過学習の領域に入ることを示している^[32][33]。

逐次NMFの場合、固有値の代わりに残差分散の比率（Fractional Residual Variance, FRV）カーブが用いられ、これが連続的に減少し、PCAよりも高いレベルで収束することが示されている^[4]。これは、逐次NMFが過学習を抑える特性を持つことを示している。

厳密NMF

NMFの変種において、追加の制約が行列 V に対して成り立つ場合、（多項式時間で）厳密解を得ることができる。行列 V がそのランクと等しいランクを持つ一般化置換部分行列を含む場合、非負ランク因子分解を解くための多項式時間アルゴリズムが Campbell と Poole によって1981年に提案された^[34]。

Kalofolias および Gallopoulos（2012）^[35] はこの問題の対称なバージョン（V が対称行列であり、ランク r の主対角部分行列を含む）を解いた。彼らのアルゴリズムは V が密な場合に O(rm²) の時間で動作する。

Arora, Ge, Halpern, Mimno, Moitra, Sontag, Wu, & Zhu（2013）^[36]は、因子行列 W の一方が可分性条件を満たす場合に、厳密なNMFを実行するための多項式時間アルゴリズムを提案している。

他の手法との関係

1999年の論文「Learning the parts of objects by non-negative matrix factorization」において、Lee と Seung は主に画像のパーツベース分解のためにNMFを提案した^[37]。この研究では、NMFをベクトル量子化や主成分分析（PCA）と比較し、3つの手法すべてが因子分解として書き表せるにもかかわらず、課す制約が異なるため、異なる結果を生むことを示している。

NMF を確率的グラフィカルモデルとして表現。観測変数（V）は隠れ変数（H）とWによって接続され、Vは平均値https://wikimedia.org/api/rest_v1/media/math/render/svg/592990f983f4c227d8bd44456f45b6e22fc64c59を持つ確率分布から生成されるとみなせる。^[12]:5

一部のNMFは、より一般的な確率モデル「多項分布主成分分析（multinomial PCA）」の一形態と見なすことができる^[38]。NMFがカルバック・ライブラー情報量（KLダイバージェンス）を最小化する場合、これは確率的潜在意味解析（PLSA）と数学的に同値である。これは、最尤推定により学習される多項分布PCAの別の実装とみなされる^[39]。

NMFにおける最小二乗法目的関数を用いた場合、これはK-meansクラスタリングの緩和版と等価である。ここでは、因子行列 W がクラスタの重心を表し、H がクラスタの所属指標を示す^[40]。ただし、k-meansでは重心の非負制約がないため、最も近い類似手法は「semi-NMF」である^[41]。

NMFは、ベイジアンネットワークの2層有向グラフィカルモデルともみなすことができる。1層が観測変数、もう1層が隠れ変数である^[42]。

NMFはテンソルにも拡張できる。これは、たとえばPARAFACモデルに対応する非負制約付きバージョンとみなされる^[43][44]。

さらに、複数のデータ行列やテンソルを同時に因子分解し、一部の因子を共有する「共通因子分解」も存在する。これはセンサーフュージョンや関係学習に有用である^[45]。

NMFは非負二次計画法の一例でもあり、これはサポートベクターマシン（support-vector machine; SVM）とも共通する性質である。さらに、両者の関係性はアルゴリズムレベルでも見られ、一方の手法において提案された解法が他方にも応用することができる^[46]。

一意性

NMFの因子分解は一意ではない。行列とその逆行列を使って、2つの因子行列を次のように変換することができる^[47]。

${\displaystyle \mathbf {WH} =\mathbf {WBB} ^{-1}\mathbf {H} }$

ここで、新たに定義された行列 ${\displaystyle \mathbf {\tilde {W}} =\mathbf {WB} }$ および ${\displaystyle \mathbf {\tilde {H}} =\mathbf {B} ^{-1}\mathbf {H} }$ が非負行列であれば、これらもまた元の因子分解の別の表現となる。

${\displaystyle \mathbf {\tilde {W}} }$ と ${\displaystyle \mathbf {\tilde {H}} }$ の非負性は、${\displaystyle \mathbf {B} }$ が非負な一般化置換行列である場合には必ず満たされる。この単純なケースでは、スケーリングや置換に相当する変換となる。

NMFの非一意性に対してより強い制御を与える手法として、スパース性制約を加える方法がある^[48]。

応用

天文学

天文学において、NMF（非負値行列因子分解）は、天体物理学的信号が非負であるという特性から、次元削減手法として有望である。NMFはスペクトル観測^[49][3] や直接撮像観測^[4] に適用されており、天体の共通特性の解析や観測データの後処理に用いられている。

Blanton と Roweis（2007）によるスペクトル観測への応用では、天文観測に伴う不確かさを考慮した。Zhu（2016）^[31]はこの手法を拡張し、欠損データや並列計算への対応も導入した。これらの手法はRenら（2018）により系外惑星や原始惑星系円盤の直接撮像観測に応用されている。

Renら（2018）は、NMFの構成要素（コンポーネント）を逐次的に構築する際の安定性を証明し、NMFモデリングの線形性を保証した。この線形性は、恒星光と惑星や円盤からの散乱光を分離する上で重要である。

直接撮像観測において、非常に明るい恒星光の中から、暗い系外惑星や原始惑星系円盤を明らかにするために、統計的手法が多数利用されている^[50][51][32]。しかし、これらの手法ではしばしば対象天体（惑星や円盤）からの光も過剰に除去されてしまう（過学習／オーバーフィッティング）ため、真のフラックスを推定するにはフォワードモデリング（順方向モデリング）が必要となる^[52][33]。

このフォワードモデリングは点状光源（惑星など）に対しては最適化されているが、原始惑星系円盤のような不規則な形状の拡張構造には適していない。そのため、NMFの非負性およびスパース性（疎性）といった性質による過学習の抑制効果が注目されており、少数のスケーリング係数を用いたフォワードモデリングが可能になる^[4]。

データ補完

統計における欠損データの補完において、NMFはコスト関数を最小化する際に欠損データを無視できるため、欠損値を0として扱うことなく処理することが可能である^[5]。この性質により、NMFはデータ補完の数学的に正当な手法として利用されている^[5]。

まず、NMFコンポーネントが既知の場合、Renら（2020）は、データ補完中に欠損データが与える影響は2次の効果（セカンドオーダー効果）であることを証明した。次に、NMFコンポーネントが未知の場合においても、欠損データが与える影響は1次または2次の効果であることが示された。

NMFコンポーネントの取得方法によっては、上記の2つのステップは相互に独立、または依存的である可能性がある。また、NMFコンポーネント数を増やすことで補完品質が向上することが確認されており、その例はRenら（2020）の図4に示されている^[5]。

テキストマイニング

NMFはテキストマイニングへの応用が可能である。このプロセスでは、さまざまな単語の重み（通常は重み付けされた単語頻度情報）から構成された文書-単語行列が作成される。この行列は、「単語-特徴量行列」と「特徴量-文書行列」に分解される。

特徴量は文書の内容から導き出され、特徴量-文書行列は関連する文書のクラスタリングを表す。

具体的な応用例としては、PubMedの科学論文要旨の小規模サブセットに対して階層的NMFを使用した研究がある^[53]。

別の研究グループは、エンロンの電子メールデータセット^[54]の一部、65,033通のメッセージと91,133語の語彙を50のクラスタに分類した^[55]。

また、NMFは引用データにも応用されており、その一例として、英語版ウィキペディア記事と学術雑誌を、ウィキペディア内の外部科学論文への引用に基づいてクラスタリングする研究がある^[56]。

Aroraら（2013）は、NMFを用いたトピックモデルの学習に関して多項式時間アルゴリズムを提示している。このアルゴリズムでは、トピック行列が「可分性」条件を満たすことを仮定しており、この条件は実際の応用においてしばしば成立する^[36]。

Hassani、Iranmanesh、Mansouri（2019）は、NMFを利用した特徴凝集法（feature agglomeration）を提案しており、これは文書-単語行列をテキストクラスタリングに適したより小さな行列に変換する手法である^[57]。

スペクトルデータ解析

NMFはスペクトルデータの解析にも使用されており、その一例として宇宙物体や宇宙ごみの分類が挙げられる^[58]。

スケーラブルなインターネット距離予測

NMFはスケーラブルなインターネット距離（ラウンドトリップタイム）の予測にも応用されている。N個のホストからなるネットワークでは、NMFの助けを借りることで、すべての N² 個のエンドツーエンドリンクの距離を、O(N) の測定だけで予測できる。この種の手法は、最初に Internet Distance Estimation Service (IDES) で導入された^[59]。

非定常音声のノイズ除去

音声のノイズ除去は、音声信号処理における長年の課題である。ノイズが定常的（時間的に変化しない）である場合には、多くのアルゴリズムが存在する。たとえば、ウィーナーフィルタは加法的なガウス雑音に対して適している。

しかし、ノイズが非定常的（時間とともに変化する）である場合には、古典的なノイズ除去アルゴリズムは性能が劣る傾向がある。これは、非定常ノイズの統計情報を正確に推定するのが困難であるためである。Schmidtら^[60]は、非負スパースコーディング（NMFの一種）を用いて非定常ノイズ下の音声を除去する手法を提案した。このアプローチは、従来の統計的手法とはまったく異なる。

この手法の鍵となる考え方は、「クリーンな音声信号は音声辞書によってスパースに表現可能であるが、非定常ノイズはそれができない」というものである。同様に、非定常ノイズもノイズ辞書によってスパースに表現できるが、音声信号はそうではない。

NMFによるノイズ除去のアルゴリズムは以下のように進められる：

1. 音声辞書とノイズ辞書の2つを事前にオフラインで学習しておく。
2. ノイズのある音声信号が与えられたら、その短時間フーリエ変換の振幅スペクトルを計算する。
3. このスペクトルを、NMFにより2つの部分に分解する。一方は音声辞書によってスパースに表現でき、もう一方はノイズ辞書によってスパースに表現できる。
4. 最終的に、音声辞書で表現された部分がクリーンな音声信号として推定される。

集団遺伝学

スパースNMFは、集団遺伝学において、個体の混血係数（admixture coefficients）の推定や、集団サンプル中の個体の遺伝的クラスタリング、ゲノムの遺伝子混合評価に利用されている。

ヒトの遺伝的クラスタリングにおいては、NMFアルゴリズムによって得られる推定結果は、コンピュータプログラムSTRUCTUREと類似しているが、NMFのほうが計算効率に優れており、大規模なゲノムデータセットの解析が可能である^[61]。

バイオインフォマティクス

NMFは、バイオインフォマティクスにおいて、遺伝子発現やDNAメチル化データのクラスタリング、クラスタを最もよく表す遺伝子の抽出などに成功裏に応用されている^{[23][62][63][64]}。

がんの変異解析では、NMFを用いて多くのがんに共通する変異パターンを同定し、それらが異なる原因によって引き起こされている可能性を示している^[65]。 NMF技術は、細胞型、疾患のサブタイプ、集団構造、組織構成、腫瘍のクローン性などの変動要因の識別にも役立つ^[66]。

NMFの特定の変種である非負値三因子分解（Non-Negative Matrix Tri-Factorization, NMTF）は、ドラッグリポジショニング（既存薬の新たな適応探索）のタスクにおいて、承認薬に対する新たなタンパク質標的および治療適応を予測するために用いられている^[67]。

また、NMTFは抗がん剤の組み合わせ効果（薬剤の相乗効果）の推定にも利用されており、がん治療における有効な薬剤ペアの発見を支援している^[68]。

核医学イメージング

非負値行列因子分解（NMF）は、この分野では「因子分析（Factor Analysis）」とも呼ばれ、1980年代からSPECTおよびPETといった動的核医学イメージングにおける画像列の解析に使用されてきた^[69]。

NMFの非一意性（解が一つに定まらないこと）は、スパース性制約を導入することによって対処されている^[70]。

例えば、Boutchkoらは、脳PET画像において、クラスタリングによる初期化を用いた「クラスタリング初期化因子分析（Clustering Initiated Factor Analysis: CIFA）」を提案し、脳血流に関連する組織の分類に応用している^[71]。

さらに、SPECT画像における不整合な投影データから4次元動的画像を再構成するために、「スプライン初期化FADSアルゴリズム（SIFADS）」が提案されている^[72]。