(57)【要約】【課題】各構成音が混合した観測信号から、各構成音を精度よく分離することができる。【解決手段】時間周波数展開部２４が、各構成音が混合された観測信号の時系列データを入力として、観測複素スペクトログラムを出力する。パラメータ推定部３６が、各構成音及び各周波数の基底スペクトル、各構成音の各時刻における音量を表すアクティベーションパラメータ、各構成音及び各周波数の各時刻における位相スペクトル、及び観測複素スペクトログラムから求められる、観測複素スペクトログラムを、各構成音に分配した各構成音の複素スペクトルの絶対値と、基底スペクトル及びアクティベーションパラメータの積との誤差の大きさを表すＩダイバージェンス規準を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルと、各構成音及び各周波数の基底スペクトルと、各構成音の各時刻におけるアクティベーションパラメータとを推定する。【選択図】図１
2018036409.tif 000002
【特許請求の範囲】
【請求項１】
各構成音が混合された観測信号の時系列データを入力として、各時刻及び各周波数の複素スペクトルを表す観測複素スペクトログラムを出力する時間周波数展開部と、
前記時間周波数展開部により出力された前記観測複素スペクトログラム、及び予め学習された音響信号の各構成音及び各周波数における振幅スペクトルを表す基底スペクトルに基づいて、各構成音及び各周波数の基底スペクトル、各構成音の各時刻における音量を表すアクティベーションパラメータ、各構成音及び各周波数の各時刻における位相スペクトル、及び前記観測複素スペクトログラムから求められる、前記観測複素スペクトログラムの各時刻及び各周波数の複素スペクトルを、各構成音に分配した各構成音の複素スペクトルの絶対値と、前記基底スペクトル及び前記アクティベーションパラメータの積との誤差の大きさを表すＩダイバージェンス規準を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルと、各構成音及び各周波数の基底スペクトルと、各構成音の各時刻におけるアクティベーションパラメータとを推定するパラメータ推定部と、
を含む信号解析装置。
【請求項２】
前記Ｉダイバージェンス規準は、以下の式で表わされる請求項１記載の信号解析装置。
【数１】
2018036409.tif 000087
2018036409.tif 000088
ただし、Ｘl,k,mは、各時刻ｍ、各周波数ｋ、及び各構成音ｌの複素スペクトルを表し、Ｈk,lは、各周波数ｋ及び各構成音ｌの基底スペクトルを表し、Ｕl,mは、各構成音ｌの各時刻ｍにおけるアクティベーションパラメータを表し、Ｙk,mは、各時刻ｍ及び各周波数ｋの観測複素スペクトルを表し、
2018036409.tif 000089
は、アクティベーションパラメータ
2018036409.tif 000090
のスパース正則化項である。
【請求項３】
前記パラメータ推定部は、
前記Ｉダイバージェンス規準の上界関数である補助関数を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルと、各構成音及び各周波数の基底スペクトルと、各構成音の各時刻におけるアクティベーションパラメータとを更新するパラメータ更新部と、
予め定められた収束条件を満たすまで、前記パラメータ更新部による更新を繰り返させる収束判定部と、
を含む請求項１又は２記載の信号解析装置。
【請求項４】
時間周波数展開部と、パラメータ推定部とを含む信号解析装置における信号解析方法であって、
前記時間周波数展開部が、各構成音が混合された観測信号の時系列データを入力として、各時刻及び各周波数の複素スペクトルを表す観測複素スペクトログラムを出力し、
前記パラメータ推定部が、前記時間周波数展開部により出力された前記観測複素スペクトログラム、及び予め学習された音響信号の各構成音及び各周波数における振幅スペクトルを表す基底スペクトルに基づいて、各構成音及び各周波数の基底スペクトル、各構成音の各時刻における音量を表すアクティベーションパラメータ、各構成音及び各周波数の各時刻における位相スペクトル、及び前記観測複素スペクトログラムから求められる、前記観測複素スペクトログラムの各時刻及び各周波数の複素スペクトルを、各構成音に分配した各構成音の複素スペクトルの絶対値と、前記基底スペクトル及び前記アクティベーションパラメータの積との誤差の大きさを表すＩダイバージェンス規準を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルと、各構成音及び各周波数の基底スペクトルと、各構成音の各時刻におけるアクティベーションパラメータとを推定する、
信号解析方法。
【請求項５】
前記Ｉダイバージェンス規準は、以下の式で表わされる請求項４記載の信号解析方法。
【数２】
2018036409.tif 000091
2018036409.tif 000092
ただし、Ｘl,k,mは、各時刻ｍ、各周波数ｋ、及び各構成音ｌの複素スペクトルを表し、Ｈk,lは、各周波数ｋ及び各構成音ｌの基底スペクトルを表し、Ｕl,mは、各構成音ｌの各時刻ｍにおけるアクティベーションパラメータを表し、Ｙk,mは、各時刻ｍ及び各周波数ｋの観測複素スペクトルを表し、
2018036409.tif 000093
は、アクティベーションパラメータ
2018036409.tif 000094
のスパース正則化項である。
【請求項６】
前記パラメータ推定部が推定することでは、
パラメータ更新部が、前記Ｉダイバージェンス規準の上界関数である補助関数を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルと、各構成音及び各周波数の基底スペクトルと、各構成音の各時刻におけるアクティベーションパラメータとを更新し、
収束判定部が、予め定められた収束条件を満たすまで、前記パラメータ更新部による更新を繰り返させることを含む請求項４又は５記載の信号解析方法。
【請求項７】
請求項１〜請求項３の何れか１項に記載の信号解析装置の各部としてコンピュータを機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、信号解析装置、方法、及びプログラムに係り、特に、パラメータを推定する信号解析装置、方法、及びプログラムに関する。
【背景技術】
【０００２】
近年、非負値行列因子分解(Non-negative Matrix Factorization; NMF)に基づく音源分離手法が注目されている。観測信号の振幅（あるいはパワー）スペクトログラムを非負値行列と見なし、これを二つの非負値行列の積で近似することは、各時刻における観測スペクトルを複数のスペクトルテンプレートの重みつき和でフィッティングしていることに相当する。従って、NMFにより求まる二つの非負値行列を用いることで観測スペクトログラムを各スペクトルテンプレートに対応する音源成分に分解することができる。このアプローチではスペクトルの加法性を仮定しているが、振幅スペクトルやパワースペクトルは信号を非線形変換したものであるため、スペクトルは実際には加法的ではない。このため、NMFにより観測スペクトログラムを加法的な成分に分解したとしても必ずしも各音源に対応したスペクトログラムが得られるとは限らない。NMF法には以上の点に解決すべき課題があった。
【０００３】
この課題を解決するため、NMF法における各スペクトルテンプレートに位相スペクトルを時変パラメータとして付加することで、NMF法と類似した考え方により複素スペクトル領域での信号分解を可能にする複素NMF（非特許文献１）と呼ぶ枠組が提案されている。従来、複素NMFにおけるパラメータ推定問題（NMF法における行列分解に相当）は、観測複素スペクトログラムと複素スペクトログラムモデルとのEuclid距離（差の絶対値の二乗）を規準とした最適化問題として定式化されており、効率的なパラメータ推定アルゴリズムが提案されている。
【先行技術文献】
【非特許文献】
【０００４】
【非特許文献１】亀岡弘和, 小野順貴, 柏野邦夫, 嵯峨山茂樹, "複素NMF: 新しいスパース信号分解表現と基底系学習アルゴリズム," 日本音響学会2008 年秋季研究発表会講演論文集, 2-8-13, pp. 657-660, Sep. 2008.
【発明の概要】
【発明が解決しようとする課題】
【０００５】
Euclid距離以外の規準でのパラメータ推定アルゴリズムはまだ報告されていない。NMFではEuclid距離、Iダイバージェンス（一般化Kullback-Leibler ダイバージェンス）、板倉齋藤距離、ダイバージェンスを規準とした場合のパラメータ推定アルゴリズムがそれぞれ提案されているが、経験的にIダイバージェンスを規準とした場合に高い音源分離性能が得られることが知られている。このことから、複素NMFにおいてもIダイバージェンスなどのEuclid距離以外の規準を用いた場合のパラメータ推定アルゴリズムを導出することができればより高い音源分離性能が得られる可能性がある。
【０００６】
本発明では、上記事情を鑑みて成されたものであり、各構成音が混合した観測信号から、各構成音を精度よく分離することができる信号解析装置、方法、及びプログラムを提供することを目的とする。
【課題を解決するための手段】
【０００７】
上記目的を達成するために、本発明に係る信号解析装置は、各構成音が混合された観測信号の時系列データを入力として、各時刻及び各周波数の複素スペクトルを表す観測複素スペクトログラムを出力する時間周波数展開部と、前記時間周波数展開部により出力された前記観測複素スペクトログラム、及び予め学習された音響信号の各構成音及び各周波数における振幅スペクトルを表す基底スペクトルに基づいて、各構成音及び各周波数の基底スペクトル、各構成音の各時刻における音量を表すアクティベーションパラメータ、各構成音及び各周波数の各時刻における位相スペクトル、及び前記観測複素スペクトログラムから求められる、前記観測複素スペクトログラムの各時刻及び各周波数の複素スペクトルを、各構成音に分配した各構成音の複素スペクトルの絶対値と、前記基底スペクトル及び前記アクティベーションパラメータの積との誤差の大きさを表すＩダイバージェンス規準を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルと、各構成音及び各周波数の基底スペクトルと、各構成音の各時刻におけるアクティベーションパラメータとを推定するパラメータ推定部と、を含んで構成されている。
【０００８】
本発明に係る信号解析方法は、時間周波数展開部と、パラメータ推定部とを含む信号解析装置における信号解析方法であって、前記時間周波数展開部が、各構成音が混合された観測信号の時系列データを入力として、各時刻及び各周波数の複素スペクトルを表す観測複素スペクトログラムを出力し、前記パラメータ推定部が、前記時間周波数展開部により出力された前記観測複素スペクトログラム、及び予め学習された音響信号の各構成音及び各周波数における振幅スペクトルを表す基底スペクトルに基づいて、各構成音及び各周波数の基底スペクトル、各構成音の各時刻における音量を表すアクティベーションパラメータ、各構成音及び各周波数の各時刻における位相スペクトル、及び前記観測複素スペクトログラムから求められる、前記観測複素スペクトログラムの各時刻及び各周波数の複素スペクトルを、各構成音に分配した各構成音の複素スペクトルの絶対値と、前記基底スペクトル及び前記アクティベーションパラメータの積との誤差の大きさを表すＩダイバージェンス規準を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルと、各構成音及び各周波数の基底スペクトルと、各構成音の各時刻におけるアクティベーションパラメータとを推定する。
【０００９】
また、本発明のプログラムは、コンピュータを、上記の信号解析装置を構成する各部として機能させるためのプログラムである。
【発明の効果】
【００１０】
以上説明したように、本発明の信号解析装置、方法、及びプログラムによれば、各構成音及び各周波数の基底スペクトル、各構成音の各時刻における音量を表すアクティベーションパラメータ、各構成音及び各周波数の各時刻における位相スペクトル、及び前記観測複素スペクトログラムから求められる、前記観測複素スペクトログラムの各時刻及び各周波数の複素スペクトルを、各構成音に分配した各構成音の複素スペクトルの絶対値と、前記基底スペクトル及び前記アクティベーションパラメータの積との誤差の大きさを表すＩダイバージェンス規準を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルと、各構成音及び各周波数の基底スペクトルと、各構成音の各時刻におけるアクティベーションパラメータとを推定することにより、各構成音が混合した観測信号から、各構成音を精度よく分離することができる。
【図面の簡単な説明】
【００１１】
【図１】本発明の実施の形態に係る信号解析装置の機能的構成を示すブロック図である。
【図２】本発明の実施の形態に係る信号解析装置における学習処理ルーチンを示すフローチャート図である。
【図３】本発明の実施の形態に係る信号解析装置におけるパラメータ推定処理ルーチンを示すフローチャート図である。
【図４】実験結果を示す図である。
【発明を実施するための形態】
【００１２】
以下、図面を参照して本発明の実施の形態を詳細に説明する。
【００１３】
＜既存手法＞
＜NMF法＞
周波数と時刻に対応するインデックスをk=1,…,K、m=1,…,Mとし、観測信号の複素スペクトログラムを
2018036409.tif 000003
、Yk,mの絶対値またはその絶対値二乗を
2018036409.tif 000004
とする。
【００１４】
NMF法は観測振幅（またはパワー）スペクトログラム
2018036409.tif 000005
を非負行列
2018036409.tif 000006
の積で近似する手法で、
2018036409.tif 000007
の各要素を
【００１５】
【数１】
2018036409.tif 000008
【００１６】
と近似することが目的となる。
【００１７】
2018036409.tif 000009
の各列ベクトルをスペクトルテンプレートと解釈すると、式(1)は、各時刻の観測スペクトルをスペクトルテンプレートの重みつき和で近似していることに相当する。ただし、lはスペクトルテンプレートに対応するインデックスである。式(1)の分解表現が得られれば、例えばWienerフィルタ
【００１８】
【数２】
2018036409.tif 000010
【００１９】
により各構成音の複素スペクトログラム
2018036409.tif 000011
を得ることができる。
【００２０】
式(1)の分解表現を得る問題は
2018036409.tif 000012
の誤差の規準を用いた最適化問題として定式化することができる。誤差規準としては二乗誤差
2018036409.tif 000013
、Iダイバージェンス
2018036409.tif 000014
【００２１】
【数３】
2018036409.tif 000015
【００２２】
などが用いられる。以上の目的関数を局所最小化する
【００２３】
2018036409.tif 000016
は次節で述べる補助関数法に基づく反復アルゴリズムにより探索することができる。ただし、(・)k,mは行列のk行m列目の要素を表す。
【００２４】
＜補助関数法＞
パラメータθの目的関数を
2018036409.tif 000017
とするとき
【００２５】
【数４】
2018036409.tif 000018
【００２６】
を満たす
2018036409.tif 000019
を
2018036409.tif 000020
の補助関数、
2018036409.tif 000021
を補助変数と呼ぶ。式(5)を満たすような補助関数を設計できれば、パラメータθと補助変数
2018036409.tif 000022
を
【００２７】
【数５】
2018036409.tif 000023
【００２８】
により反復的に更新することで目的関数
2018036409.tif 000024
を局所最小化するθを得ることができる。従来のNMFの反復アルゴリズムは上述の原理に基づいて導かれる。
【００２９】
＜複素NMF＞
式(1)では、振幅スペクトルまたはパワースペクトルが加法的であるという仮定の下で混合音スペクトルがモデル化されているが、振幅スペクトルもパワースペクトルも実際には加法的ではない。そこで、加法性が成り立つ複素スペクトル領域でNMF法と同様の考え方による信号分解を可能にする枠組として提案されたのが複素NMF（非特許文献１）である。
【００３０】
前述のNMFのモデルにおける各スペクトルテンプレートに時変の位相スペクトル
2018036409.tif 000025
をパラメータとして付加することで
【００３１】
【数６】
2018036409.tif 000026
【００３２】
のような形の複素スペクトログラムのモデルを立てることができる。観測信号の複素スペクトログラムを
2018036409.tif 000027
とすると、複素NMFはこのモデルと観測複素スペクトログラムとの間のEuclid距離
【００３３】
【数７】
2018036409.tif 000028
【００３４】
を最小化する問題として定式化される。ただし、上式の第二項は
2018036409.tif 000029
の要素がスパースになるように誘導する目的の正則化項であり、0＜p≦2、λ＞0とする。この目的関数を最小にする
2018036409.tif 000030
を解析的に得ることはできないが、補助関数法により局所最適解を得るための反復アルゴリズムを導くことができる。まず第一項に関し、二次関数は凸関数であるため、Jensenの不等式
【００３５】
【数８】
2018036409.tif 000031
【００３６】
が成り立つ。ただし、
2018036409.tif 000032
である。次に第二項に関し、0＜p≦2より、
【００３７】
【数９】
2018036409.tif 000033
【００３８】
が成り立つ。ただし、Vl,m≧0である。以上より、
2018036409.tif 000034
を補助変数にもつ
【００３９】
【数１０】
2018036409.tif 000035
【００４０】
のような補助関数を設計することができる。また、βl,k,mは、βl,k,m＞0、Σlβl,k,m＝1を満たす任意の定数である。以上の補助関数を最小にする補助変数Xl,k,m、Vl,mは、
【００４１】
【数１１】
2018036409.tif 000036
【００４２】
となる。式(13)は、各ステップで推定されたl番目の構成音に対応する複素スペクトログラムモデル
2018036409.tif 000037
に全体のモデル化誤差Yk,m−Fk,mのβl,k,m倍を足したものをXl,k,mと置いていることを意味する。このことから、Xl,k,mはΣlXl,k,m=Yk,mとなるように分配されたl番目の構成音に対応する複素スペクトログラムの成分と解釈することができる。同様に、以上の補助関数を最小にするパラメータ
2018036409.tif 000038
も解析的に
【００４３】
【数１２】
2018036409.tif 000039
【００４４】
と求まる。なお、βl,k,m＞0、Σlβl,k,m=1を満たす範囲であれば、βl,k,mを各ステップで更新したとしてもアルゴリズムの収束性は保証されており、上記非特許文献１では
2018036409.tif 000040
のような更新式が採用されている。以上より、複素NMFアルゴリズムは下記のようにまとめることができる。
【００４５】
Step0)  初期値設定
【００４６】
Step1)  補助変数の更新
【００４７】
【数１３】
2018036409.tif 000041
【００４８】
Step2)  モデルパラメータの更新
【００４９】
【数１４】
2018036409.tif 000042
【００５０】
＜問題の定式化＞
【００５１】
従来の複素NMFでは各構成音の位相スペクトログラム
2018036409.tif 000043
はモデルパラメータであり、Xl,k,mは各構成音の複素スペクトログラムの成分を意味する補助変数であった。ここで、
2018036409.tif 000044
とXl,k,mに対する扱いを逆転させ、
2018036409.tif 000045
を補助変数、Xl,k,mをモデルパラメータと扱った場合の、上述の複素NMFアルゴリズムと等価な手続きを行う反復アルゴリズムについて考える。理由は後述するが、実はこの反復アルゴリズムは次の最適化問題(P0)
【００５２】
【数１５】
2018036409.tif 000046
【００５３】
の局所最適解を得るためのアルゴリズムとなっている。このことを以下で確認する。
【００５４】
表記の簡単化のため
2018036409.tif 000047
とすると、複素NMFの補助関数
2018036409.tif 000048
の第一項は、
【００５５】
【数１６】
2018036409.tif 000049
【００５６】
と書ける。ここでg(z)=−Re(c*z)(|c|=1)という関数とf(z)=-|z|という関数をRe(z)とIm(z)の2変数の関数として見たとき、g(z)がf(z)にcの方向で接する接平面になっていることに注意すると、
【００５７】
【数１７】
2018036409.tif 000050
【００５８】
が成り立つ。従って、
【００５９】
【数１８】
2018036409.tif 000051
【００６０】
2018036409.tif 000052
【００６１】
が言え、cl,k,mを補助変数、Xl,k,mをモデルパラメータと見なせば、
2018036409.tif 000053
は最適化問題(P0)の目的関数の補助関数にもなっていることが分かる。従って、(P0)の局所最適解を探索するアルゴリズムも
2018036409.tif 000054
を用いて導くことができ、上述した複素NMFアルゴリズムにおいて式(18)と式(22)を入れ替えたものとなる。以上の関係より、(P0)を複素NMFの双対形式と呼ぶことにする。
【００６２】
以上の複素NMFの双対形式では、|Xl,k,m|とHk,lUl,mの二乗誤差が目的関数となっているが、|Xl,k,m|とHk,lUl,mはいずれも非負値のため、二乗誤差の代わりにIダイバージェンス規準を用いて両者の誤差を測ることもできる。そこで、
【００６３】
【数１９】
2018036409.tif 000055
【００６４】
のような最適化問題を考える。この問題をIダイバージェンス規準の複素NMFと定義する。ただし、
2018036409.tif 000056
は
2018036409.tif 000057
に対するスパース正則化項であり、
【００６５】
【数２０】
2018036409.tif 000058
【００６６】
とする。
【００６７】
＜パラメータ推定アルゴリズムの導出＞
＜補助関数の設計＞
2018036409.tif 000059
,
はXl,k,mに関して微分不可能な点を含んでおり、
2018036409.tif 000060
を最小化するXl,k,mは解析的に求まらない。しかし、下記に示す補助関数を用いることで
2018036409.tif 000061
を局所最小化するHk,l、Ul,m、Xl,k,mを得るためのアルゴリズムを導くことができる。まず、対数関数は凹関数より、対数関数の接線が対数関数を下回らないことを利用すると、
【００６８】
【数２１】
2018036409.tif 000062
【００６９】
が成立する。明らかに等号はξl,k,m = |Xl,k,m|のとき成立する。式(35)の第二項には依然として微分不可能な点が含まれる。ここで、|Xl,k,m|の係数を
2018036409.tif 000063
と置くと、Dl,k,mはξl,k,m、Hk,l、Ul,m の値に応じて正にも負にもなりうる。そこで、
【００７０】
【数２２】
2018036409.tif 000064
および
2018036409.tif 000065
【００７１】
が成り立つ（等号はbl,k,m= |Xl,k,m|、
2018036409.tif 000066
のとき成立する）ことを利用し、Dl,k,mの符号に合わせて上記のいずれかの不等式を適用することで、
【００７２】
【数２３】
2018036409.tif 000067
【００７３】
という不等式の形にまとめることができる。ただし、
【００７４】
【数２４】
2018036409.tif 000068
【００７５】
である。
【００７６】
また、正則化項
2018036409.tif 000069
については、
【００７７】
【数２５】
2018036409.tif 000070
【００７８】
が成り立つ。
【００７９】
以上より、
【００８０】
【数２６】
2018036409.tif 000071
【００８１】
を
2018036409.tif 000072
【００８２】
の補助関数として用いることができる。
【００８３】
＜更新式＞
Xl,k,mにはΣlXl,k,m=Yk,mという制約があるため、ラグランジアン
【００８４】
【数２７】
2018036409.tif 000073
2018036409.tif 000074
【００８５】
を用いて補助関数
2018036409.tif 000075
を最小にするXl,k,mの更新式を求めることができる。
【００８６】
2018036409.tif 000076
についての偏微分を0と置くことで
【００８７】
【数２８】
2018036409.tif 000077
【００８８】
を得る。これをΣlXl,k,m=Yk,mに代入することで、未定乗数λ*k,mを
【００８９】
2018036409.tif 000078
【００９０】
より得ることができる。式(48)を式(46)に代入することで，Xl,k,mの更新式
【００９１】
【数２９】
2018036409.tif 000079
【００９２】
が得られる。次に、目的関数
2018036409.tif 000080
を最小にするHk,l、Ul,mの更新式は、
2018036409.tif 000081
をそれぞれの変数に関して偏微分したものを0と置くことで
【００９３】
【数３０】
2018036409.tif 000082
【００９４】
のように得られる。
【００９５】
以上をまとめると次のようなアルゴリズムを得る。
【００９６】
Step0)初期値設定
【００９７】
Step1)補助変数の更新
【００９８】
【数３１】
2018036409.tif 000083
2018036409.tif 000084
【００９９】
Step2)モデルパラメータの更新
【０１００】
【数３２】
2018036409.tif 000085
【０１０１】
ただし、ｃl,k,mは、各構成音ｌ及び各周波数ｋの各時刻ｍにおける位相スペクトルを表す。
【０１０２】
＜本発明の実施の形態に係る信号解析装置の構成＞
次に、本発明の実施の形態に係る信号解析装置の構成について説明する。図１に示すように、本発明の実施の形態に係る信号解析装置１００は、ＣＰＵと、ＲＡＭと、後述する学習処理ルーチン及びパラメータ推定処理ルーチンを実行するためのプログラムや各種データを記憶したＲＯＭと、を含むコンピュータで構成することができる。この信号解析装置１００は、機能的には図１に示すように入力部１０と、演算部２０と、出力部９０と、を含んで構成されている。
【０１０３】
入力部１０は、各構成音の構成音信号の時系列データを受け付ける。また、入力部１０は、複数の構成音が混じっている音響信号（以後、観測信号）の時系列データを受け付ける。
【０１０４】
演算部２０は、時間周波数展開部２４と、基底スペクトル学習部２８と、基底スペクトル記憶部３０と、パラメータ推定部３６と、構成音信号生成部３８と、を含んで構成されている。
【０１０５】
時間周波数展開部２４は、各構成音信号の時系列データに基づいて、各時刻における各周波数の複素スペクトルを表す振幅スペクトログラム又はパワースペクトログラムを計算する。なお、第１の実施の形態においては、短時間フーリエ変換やウェーブレット変換などの時間周波数展開を行う。
【０１０６】
また、時間周波数展開部２４は、観測信号の時系列データに基づいて、各時刻ｍにおける各周波数ｋの観測複素スペクトルＹk,mを表す振幅スペクトログラム又はパワースペクトログラムである観測複素スペクトログラムＹを計算する。
【０１０７】
基底スペクトル学習部２８は、時間周波数展開部２４によって計算された、各構成音信号の各時刻における各周波数の複素スペクトルに基づいて、従来技術である複素ＮＭＦを用いて、各構成音ｌ及び各周波数kにおける振幅スペクトルを表す基底スペクトルＨk,lを推定する。
【０１０８】
基底スペクトル記憶部３０は、基底スペクトル学習部２８によって推定された、各構成音ｌ及び各周波数kにおける振幅スペクトルを表す基底スペクトルＨk,lを記憶している。
【０１０９】
パラメータ推定部３６は、時間周波数展開部２４により出力された観測複素スペクトログラムＹ、及び基底スペクトル記憶部３０に記憶された各構成音及び各周波数における基底スペクトルに基づいて、各構成音及び各周波数の基底スペクトルＨk.l、各構成音の各時刻における音量を表すアクティベーションパラメータＵl,m、各構成音及び各周波数の各時刻における位相スペクトル
2018036409.tif 000086
、及び観測複素スペクトログラムＹから求められる、観測複素スペクトログラムの各時刻及び各周波数の複素スペクトルを、各構成音に分配した各構成音の複素スペクトルＸl,k,mの絶対値と、基底スペクトルＨk.l及びアクティベーションパラメータＵl,mの積との誤差の大きさを表す上記（３０）式のＩダイバージェンス規準を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mと、各構成音及び各周波数の基底スペクトルＨk.lと、各構成音の各時刻におけるアクティベーションパラメータＵl,mとを推定する。
【０１１０】
具体的には、パラメータ推定部３６は、初期値設定部４０、補助変数更新部４２、パラメータ更新部４４、及び収束判定部４６を備えている。
【０１１１】
初期値設定部４０は、各構成音信号の基底スペクトルＨの初期値として、基底スペクトル記憶部３０に記憶された各構成音信号の各周波数における基底スペクトルＨk,lを設定する。また、初期値設定部４０は、各構成音信号のアクティベーションパラメータＵl,mと、各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mとに初期値を設定する。
【０１１２】
補助変数更新部４２は、初期値である、又は前回更新した、各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mと、各構成音信号のアクティベーションパラメータＵl,mとに基づいて、上記（５４）式〜（５７）式に従って、各構成音ｌ、各周波数k、及び各時刻ｍに対するζl,k,m、ｂl,k,m、ｃl,k,m、各構成音ｌ及び各時刻ｍに対するＶl,mを更新する。
【０１１３】
パラメータ更新部４４は、時間周波数展開部２４により出力された観測複素スペクトログラムＹと、補助変数更新部４２によって更新された各構成音ｌ、各周波数k、及び各時刻ｍに対するζl,k,m、ｂl,k,m、ｃl,k,m、各構成音ｌ及び各時刻ｍに対するＶl,mと、初期値である、又は前回更新した、各構成音信号の各周波数における基底スペクトルＨk,lと、各構成音信号のアクティベーションパラメータＵl,mと、各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mとに基づいて、上記（５８）式〜（６３）式に従って、各構成音信号の各周波数における基底スペクトルＨk,lと、各構成音信号のアクティベーションパラメータＵl,mと、各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mとを推定する。
【０１１４】
収束判定部４６は、収束条件を満たすか否かを判定し、収束条件を満たすまで、補助変数更新部４２における更新処理と、パラメータ更新部４４における更新処理とを繰り返させる。
【０１１５】
収束条件としては、例えば、繰り返し回数が、上限回数に到達したことを用いることができる。あるいは、収束条件として、上記（３０）式の規準の値と前回の規準の値との差分が、予め定められた閾値以下であることを用いることができる。
【０１１６】
構成音信号生成部３８は、パラメータ推定部３６において取得した各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mとに基づいて、各構成音の構成音信号を生成し、出力部９０から出力する。
【０１１７】
＜本発明の実施の形態に係る信号解析装置の作用＞
次に、本発明の実施の形態に係る信号解析装置１００の作用について説明する。まず、入力部１０において各構成音の構成音信号の時系列データを受け付けると、信号解析装置１００は、図２に示す学習処理ルーチンを実行する。
【０１１８】
まず、ステップＳ１００では、入力部１０において受け付けた各構成音の構成音信号の時系列データに基づいて、各構成音信号の各時刻における各周波数の複素スペクトルを計算する。
【０１１９】
次に、ステップＳ１０４では、上記ステップＳ１００において取得した各構成音信号の各時刻における各周波数の複素スペクトルに基づいて、従来技術である複素ＮＭＦにより、各構成音ｌ及び各周波数kにおける振幅スペクトルを表す基底スペクトルＨk,lを推定し、基底スペクトル記憶部３０に格納して、学習処理ルーチンを終了する。
【０１２０】
次に、入力部１０において、各構成音が混在した観測信号の時系列データを受け付けると、信号解析装置１００は、図３に示すパラメータ推定処理ルーチンを実行する。
【０１２１】
まず、ステップＳ１２０では、入力部１０において受け付けた観測信号の時系列データに基づいて、観測複素スペクトログラムＹを計算する。
【０１２２】
ステップＳ１２２では、各構成音ｌの基底スペクトルＨk,lの初期値として、基底スペクトル記憶部３０に記憶された各構成音及び各周波数における基底スペクトルを設定する。また、各構成音のアクティベーションパラメータＵl,mと、各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mとに初期値を設定する。
【０１２３】
ステップＳ１２４では、初期値である、又は後述するステップＳ１２５で前回更新した、各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mと、各構成音信号のアクティベーションパラメータＵl,mとに基づいて、上記（５４）式〜（５７）式に従って、各構成音ｌ、各周波数k、及び各時刻ｍに対するζl,k,m、ｂl,k,m、ｃl,k,m、各構成音ｌ及び各時刻ｍに対するＶl,mを更新する。
【０１２４】
次に、ステップＳ１２５では、上記ステップＳ１２０で得られた観測複素スペクトログラムＹと、上記ステップＳ１２４で更新された、各構成音ｌ、各周波数k、及び各時刻ｍに対するζl,k,m、ｂl,k,m、ｃl,k,m、各構成音ｌ及び各時刻ｍに対するＶl,mと、初期値である、又は前回更新した、各構成音信号の各周波数における基底スペクトルＨk,lと、各構成音信号のアクティベーションパラメータＵl,mと、各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mとに基づいて、上記（５８）式〜（６３）式に従って、各構成音信号の各周波数における基底スペクトルＨk,lと、各構成音信号のアクティベーションパラメータＵl,mと、各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mとを推定する。
【０１２５】
次に、ステップＳ１２８では、収束条件を満たすか否かを判定する。収束条件を満たした場合には、ステップＳ１３０へ移行し、収束条件を満たしていない場合には、ステップＳ１２４へ移行し、ステップＳ１２４〜ステップＳ１２５の処理を繰り返す。
【０１２６】
ステップＳ１３０では、上記ステップＳ１２５で最終的に更新された各時刻、各周波数、及び各構成音の複素スペクトルＸl,k,mに基づいて、各構成音信号を生成し、出力部９０から出力して、パラメータ推定処理ルーチンを終了する。
【０１２７】
＜実験例＞
本実施の形態の手法の有効性を確認するため、SiSEC 2013 Database の"professionally produced music recordings"（5 音源からなるポピュラー音楽の音響信号）を用いて教師あり音源分離実験を行い、二乗誤差規準NMF(EU-NMF)、Iダイバージェンス規準NMF(KL-NMF)、板倉齋藤距離規準NMF(IS-NMF)、二乗誤差規準複素NMF(EU-CNMF)，Iダイバージェンス規準複素NMF(KL-CNMF)による各音源のSignal to Distortion Ratio (SDR) 改善値を比較した。図４にその結果を示す。５音源のうち４音源においてKL-CNMF が他手法に比べて高い改善値を得られていることが確認できた。
【０１２８】
以上説明したように、本発明の実施の形態に係る信号解析装置によれば、各構成音及び各周波数の基底スペクトル、各構成音の各時刻における音量を表すアクティベーションパラメータ、各構成音及び各周波数の各時刻における位相スペクトル、及び観測複素スペクトログラムから求められる、観測複素スペクトログラムの各時刻及び各周波数の複素スペクトルを、各構成音に分配した各構成音の複素スペクトルの絶対値と、基底スペクトル及び前記アクティベーションパラメータの積との誤差の大きさを表すＩダイバージェンス規準を小さくするように、各時刻、各周波数、及び各構成音の複素スペクトルと、各構成音及び各周波数の基底スペクトルと、各構成音の各時刻におけるアクティベーションパラメータとを推定することにより、各構成音が混合した観測信号から、各構成音を精度よく分離することができる。
【０１２９】
なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。
【０１３０】
例えば、各構成音信号の基底スペクトルを学習する処理と、観測信号からパラメータを推定するパラメータ推定とを別々の装置で行うように構成してもよい。
【０１３１】
また、更新するパラメータの順番には任意性があるため、上記の実施の形態の順番に限定されない。
【０１３２】
また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能であるし、ネットワークを介して提供することも可能である。
【符号の説明】
【０１３３】
１０   入力部
２０   演算部
２４   時間周波数展開部
２８   基底スペクトル学習部
３０   基底スペクトル記憶部
３６   パラメータ推定部
３８   構成音信号生成部
４０   初期値設定部
４２   補助変数更新部
４４   パラメータ更新部
４６   収束判定部
９０   出力部
１００ 信号解析装置
【図１】
2018036409.tif 000095
【図２】
2018036409.tif 000096
【図３】
2018036409.tif 000097
【図４】
2018036409.tif 000098