JP6448567B2 - 音響信号解析装置、音響信号解析方法、及びプログラム - Google Patents

Description

本発明は、音響信号解析装置、音響信号解析方法、及びプログラムに係り、特に、複数のマイクロホンで取得した音響信号を用いて、残響除去と音源分離を行うための音響信号解析装置、方法、及びプログラムに関する。

複数のマイクロホンで取得した多チャネル信号を処理し、音源の空間情報を手がかりにして音源分離などを行う枠組をマイクロホンアレー信号処理という。

近年、マイクロホンアレー信号処理に関して、ボイスレコーダ、ノートパソコン、スマートフォン、及びビデオカメラ等の身の回りにある様々な録音機器による多チャンネル録音を用いた、アドホックマイクロホンアレーの研究が盛んに行われている。アドホックマイクロホンアレーを用いた場合、特殊な装置や配線を要する従来のマイクロホンアレーを用いた場合に比べて手軽かつ安価にマイクロホンアレーシステムを構築できるため、注目を集めている。

従来のマイクロホンアレーで商用化されているものの多くは、各マイクロホンが小規模に集中配置されているため、録音チャンネル間の音声信号の時間差が音源分離のための手がかりとなるのに対し、アドホックマイクロホンアレーを用いた場合、従来のマイクロホンアレーに比べてマイクロホンを広範囲に分散して配置することが容易となるため、録音チャンネル間の時間差に加えて、更に音声信号の強度比も音源分離のための手がかりとなる。

一般に、音声信号に残響及び雑音が重畳され、観測信号が得られるプロセスを順問題と捉えると、アドホックマイクロホンアレーで集音した観測信号から目的音声のみを分離抽出する問題は逆問題とみなすことができる。雑音、或いは室内伝達系の情報が未知の場合で、且つ、マイクロホン数より音源数が多いという劣決定条件の場合、この逆問題には解が無数に存在しうるため、解を絞り込むための何らかの仮定が必要となる。

一方、近年、劣決定条件における音源分離手法の一例として、非負値行列因子分解(Non-negative Matrix Factorization: NMF) の多チャンネル拡張を用いたアプローチが注目されている（非特許文献１、２）。

NMF とは、非負値行列を二つの非負値行列（基底行列と係数行列）の積に分解することをいい、スペクトログラムを非負値行列と見なしてNMF を適用することはスペクトログラムを低ランクの非負値行列で近似することに相当し、各時刻のスペクトルを基底行列の列数分のスペクトルテンプレートの非負結合で説明しようとしていることを意味する。そして、NMF の多チャンネル拡張は、各音源のパワースペクトログラムにこの構造を仮定した多チャンネル音源分離手法の一例である。なお、音源数よりマイクロホン数が多いという過剰決定条件でのNMF の多チャンネル拡張の例も提案されている（非特許文献３）。

A. Ozerov and C. F_evotte, "Multichannel nonnegative matrix factorization in convolutive mixtures for audio source separation," IEEE Transactions on Audio, Speech and Language Processing, vol.18, no. 3, pp. 550-563, 2010. A. Ozerov, C. F_evotte, R. Blouet and J.L. Durrieu, "Multichannel nonnegative tensor factorization with structured constraints for user-guided audio source separation," in Proc. 2011 IEEE International Conference on Acoustics, Speech and Signal Processing (ICASSP2011), pp. 257-260, 2011. Hirokazu Kameoka, Takuya Yoshioka, Mariko Hamamura, Jonathan Le Roux, Kunio Kashino,"Statistical model of speech signals based on composite autoregressive system with application to blind source separation," in Proc. 9th International Conference on Latent Variable Analysis and Signal Separation (LVA/ICA 2010), LNCS 6365, pp. 245-253, Sep. 2010. T. Yoshioka, Tomohiro Nakatani, Masato Miyoshi, and Hiroshi G. Okuno, "Blind separation and dereverberation of speech mixtures by joint optimization," IEEE Transactions on Audio, Speech and Language Processing, vol. 19, no. 1, pp. 69-84, 2011.

上述した従来のアプローチでは、室内伝達系に時不変性などの制約が置かれ、その条件の下で逆問題が定式化されるが、アドホックマイクロホンアレーの枠組では手軽にアレーシステムが構築できる利点がある一方で、各マイクロホンの位置は固定されないため音源及びマイクロホンの相対位置関係が録音中に変化しやすいという脆弱性を有している。このように録音中に音源及びマイクロホンの相対位置関係が変化する場合、室内伝達系に対する上述の仮定が成立しなくなり、当該仮定の下で設計されたアルゴリズムは、音源分離に関して高い性能を発揮できなくなるという問題があった。

本発明は、上記の事情を鑑みてなされたものであり、音源及びマイクロホンの相対位置関係が変化する時変残響環境下であっても、複数の音源の音声が重畳した観測信号から、各音源信号を精度よく分離することができる音響信号解析装置、方法、及びプログラムを提供することを目的とする。

上記の目的を達成するために本発明に係る音響信号解析装置は、Ｊ本の各マイクロホンｊで集音した音響信号の時系列データを入力として、各時刻ｌにおける各周波数ｋの観測信号時間周波数成分ｙ_j,k,lを出力する時間周波数解析部と、音源ｉから前記マイクロホンｊまでの時刻ｎだけ遅延して集音される伝達特性を表す時変ステアリングベクトルの振幅成分Ａ_j,i,k,n、複数フレームのスペクトルを連結したスペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底スペクトルＷ_i,k,m,τ、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底オンセットＨ_i,m,τの各々に初期値を設定する初期値設定部と、（j,k,l）の全ての組み合わせにおける、前記観測信号時間周波数成分ｙ_j,k,lと、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底スペクトルＷ_i,k,m,τ、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底オンセットＨ_i,m,τ、各音源ｉ及び各時刻ｎの前記振幅成分Ａ_j,i,k,nに基づいて算出される前記マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lとの距離が小さくなるように、前記基底スペクトルＷ_i,k,m,τと、前記基底オンセットＨ_i,m,τと、前記振幅成分Ａ_j,i,k,nとを更新するパラメータ更新部と、予め定められた終了条件を満たすまで、前記パラメータ更新部による更新を繰り返し行う終了判定部と、を含んで構成されている。

本発明に係る音響信号解析方法は、時間周波数解析部と、初期値設定部と、パラメータ更新部と、終了判定部とを含む音響信号解析装置における音響信号解析方法であって、前記時間周波数解析部が、Ｊ本の各マイクロホンｊで集音した音響信号の時系列データを入力として、各時刻ｌにおける各周波数ｋの観測信号時間周波数成分ｙ_j,k,lを出力し、前記初期値設定部が、音源ｉから前記マイクロホンｊまでの時刻ｎだけ遅延して集音される伝達特性を表す時変ステアリングベクトルの振幅成分Ａ_j,i,k,n、複数フレームのスペクトルを連結したスペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底スペクトルＷ_i,k,m,τ、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底オンセットＨ_i,m,τの各々に初期値を設定し、前記パラメータ更新部が、（j,k,l）の全ての組み合わせにおける、前記観測信号時間周波数成分ｙ_j,k,lと、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底スペクトルＷ_i,k,m,τ、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底オンセットＨ_i,m,τ、各音源ｉ及び各時刻ｎの前記振幅成分Ａ_j,i,k,nに基づいて算出される前記マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lとの距離が小さくなるように、前記基底スペクトルＷ_i,k,m,τと、前記基底オンセットＨ_i,m,τと、前記振幅成分Ａ_j,i,k,nとを更新し、前記終了判定部が、予め定められた終了条件を満たすまで、前記パラメータ更新部による更新を繰り返し行う。

本発明に係るプログラムは、上記の音響信号解析装置の各部としてコンピュータを機能させるためのプログラムである。

以上説明したように、本発明の音響信号解析装置、方法、及びプログラムによれば、観測信号時間周波数成分ｙ_j,k,lと、スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底スペクトルＷ_i,k,m,τ、スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底オンセットＨ_i,m,τ、各音源ｉ及び各時刻ｎにおける時変ステアリングベクトルの振幅成分Ａ_j,i,k,nに基づいて算出されるマイクロホンｊのパワースペクトログラムモデルＸ_j,k,lとの距離が小さくなるように、基底スペクトルＷ_i,k,m,τと、基底オンセットＨ_i,m,τと、振幅成分Ａ_j,i,k,nとを更新することを繰り返すことにより、複数の音源の音声が重畳した観測信号から、各音源信号を精度よく分離することができる、という効果が得られる。

本発明の実施の形態に係る音響信号解析装置の構成を示す概略図である。 本発明の実施の形態に係る音響信号解析装置における音響信号解析処理ルーチンの内容を示すフローチャートである。 音響信号解析装置の評価実験を実施する部屋の環境を示す図である。 音響信号解析装置の評価実験における、壁の反射係数に対するSource to Distortion Ratio(SDR)の変化の一例を示すグラフである。 音響信号解析装置の評価実験における、伝達系への外乱に対するSDRの変化の一例を示すグラフである。

以下、図面を参照して本発明の実施の形態を詳細に説明する。

（発明の原理）

まず、本発明の提案モデルについて説明する。

（時間周波数領域における畳み込み混合モデル）

音源からマイクロホンアレーへの伝達系が線形時不変であり、また、残響成分が時間周波数解析の窓長内に収まっていると仮定できる場合、マイクロホンアレーで得られる信号は音源の瞬時混合で記述することができる。

一方、窓長を超える残響成分が無視できない状況下では、マイクロホンで観測される観測信号は、例えば、非特許文献４に示されているように、時間周波数領域の畳み込み混合モデルを用いて（１）式のように表される。

ここで、変数ｉ、ｋは、それぞれ音源と周波数のインデックスを表し、ｌ、ｎは時間フレームのインデックスを表す。ｙ^_k,l∈Ｃ^Ｊは、マイクロホンアレーで観測される観測信号ベクトルであり、ｊはマイクロホンのインデックスを表し、Ｊはマイクロホン数を表す。ａ^_i,k,nは、各音源からマイクロホンへのステアリングベクトルであり、時間がｎフレーム遅れて到来する反響成分に対応するものとする。また、ｓ_i,k,lは、各音源の時間周波数領域における複素スペクトログラムである。以降、ベクトル、行列又は確率変数を表す変数には“^”を付すものとする。

今、音源或いはマイクロホンの位置が時間変化すると仮定した場合、ステアリングベクトルａ^_i,k,nは時刻ｌに依存し、（１）式の混合過程は、（２）式のような時変ステアリングベクトルａ^_i,k,n,lを用いて表される。

ここで、各音源ｉの複素スペクトログラムが複素ガウス分布に従う、すなわち、ｓ_i,k,l〜Ｎ_Ｃ(０,Ｐ_i,k,l)と仮定すれば、マイクロホンでの観測信号ベクトルｙ^_k,lは、（３）式に示す分布に従う。なお、Ｐ_i,k,lは、音源ｉの時刻ｌにおけるパワースペクトログラムを表し、ａ^^H _i,k,n,lは、時変ステアリングベクトルａ^_i,k,n,lのエルミート転置を表す。

ここで、時変ステアリングベクトルａ^ _i,k,n,lを絶対値と偏角の要素に分解して表した式を（４）式に示す。

今、音源或いはマイクロホンの微小移動等、音響信号解析環境の軽微な変化については、時変ステアリングベクトルａ^ _i,k,n,lの振幅成分を時不変、位相成分を時変と仮定する特殊な時変系を設定し、このような混合過程を「半時変形」と呼ぶことにする。すなわち、|ａ^_j,i,k,n,l|は時刻ｌに依存しないことになる。したがって、|ａ^_j,i,k,n,l|＝Ａ_j,i,k,nとすれば、（４）式は（５）式のように表すことができる。

（５）式を（３）式に代入すれば、（６）式が得られる。

（非負値テンソル二重畳み込みモデル）

アドホックマイクロホンアレーでは、アレー素子が非同期であることによる音声信号のサンプリング周波数の僅かなずれ、並びに音源或いはマイクロホンの僅かな位置の変化が、通常のマイクロホンアレーに比べて起こりやすいため、半時変系の混合過程として取り扱う必要がある。

半時変系の混合過程を取り扱う際の一つの解決策は、時変ステアリングベクトルａ^_{ｉ,ｋ,ｎ,ｌ}の時間変化量をオンラインで推定して補償した後に、公知のアレー信号処理を適用することである。

一方、時変ステアリングベクトルａ^_i,k,n,lの位相成分を確率的に変動する確率変数として扱う方法も考えられるが、ここでは、後者の考え方にしたがって、時変ステアリングベクトルａ^_i,k,n,lの位相成分に、次の２つの条件を設定する。

（条件１）φ_j,i,k,n,lとφ_j',i,k,n,l'(j≠ｊ^'またはｌ≠ｌ^')は互いに独立である

（条件２）φ_j,i,k,n,lは区間[０,２π)で一様分布に従う。

（条件１）及び（条件２）の下で位相成分φ_j,i,k,n,lに関して周辺化を行うと、Ｅ[ψ^_i,k,n,l ψ^^H _i,k,n,l]は単位行列となるため、（６）式は（７）式のように表される。

次に、音源ｉのパワースペクトログラムＰ_i,k,lについて検討する。

周波数のインデックスｋの最大値をＫ、及び時間フレームのインデックスｌの最大値をＬとすれば、NMFの多チャンネル拡張では、Ｋ行Ｌ列で表されるパワースペクトログラムベクトルＰ^_i＝(Ｐ_i,k,l)_K×Lを２つの非負値行列の積で表される。これは、音源ｉの時刻ｌにおけるパワースペクトログラムＰ_i,k,lが、限られた数のスペクトルテンプレートの非負結合で表されるという仮定に基づくものである。

しかし、音声には単一の時間フレームｌにおけるスペクトルのみならず、そのダイナミクス、すなわち、局所的な時間変化パターンに大きな特徴があるものと考えられる。したがって、各時間フレームｌにおけるスペクトルを音声を構成する要素と考えるよりも、複数の時間フレームｌに亘るスペクトルを連結したものを、音声を構成する要素単位とみなした方が、音声をより特徴付けた好ましい表現と考えることができる。

そこで、本実施の形態では、複数の時間フレームｌに亘るスペクトルを連結したスペクトログラム素片のテンプレートとアクティベーション系列とを畳み込む混合モデルによって、パワースペクトログラムＰ_i,k,lのモデル化を行う。具体的には、公知の非負値行列因子逆畳み込み(Nonnegative Matrix Factor Deconvolution:NMFD)の考え方を、パワースペクトログラムＰ_i,k,lのモデル化に適用する。この場合、音源ｉの時刻ｌにおけるパワースペクトログラムＰ_i,k,lは（８）式で表される。

ここで、Ｗ_i,k,m,τは音源ｉのパワースペクトログラムＰ_i,k,lの基底スペクトルを表し、Ｈ_i,m,lは基底オンセットを表す。また、変数ｍは基底のインデックスを表し、変数τは時間フレームのインデックスを表す。なお、τ＝｛０｝（｛ｘ｝は、“ｘ”が集合の要素に含まれていることを示す）の場合は、NMFを用いた混合モデルに一致する。

（パラメータ推定）

次に、最尤推定を用いたパラメータ推定について説明する。まず、（７）式の対数尤度関数の負値をＣ_ＭＬとして（９）式のように定義する。なお、Ｙ^は、観測信号ｙ_j,k,lに対応する確率変数、Ｐ^は、パワースペクトログラムＰ_i,k,lに対応する確率変数、及びＡ^は、時変ステアリングベクトルａ^_i,k,n,lの絶対値に対応する確率変数である。また、観測信号ｙ_j,k,lは、観測信号時間周波数成分の一例である。

（９）式を用いてＣ_MLを具体的に計算すれば（１０）式が得られる。

ここで、記号＝^cは、定数項以外の他の項は等しいことを表す。

すなわち、最尤推定を用いたパラメータ推定は、マイクロホンｊの観測パワースペクトログラムＹ_j,k,lと、音源ｉのパワースペクトログラムＰ_i,k,lとして（８）式を適用した場合における、マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lとの板倉斎藤距離の最小化問題に帰着し、目的関数である（１１）式をＡ_i,j,k,n≧０、Ｗ_i,k,m,τ≧０、及びＨ_i,m,τ≧０の制約の下で最適化する問題になる。

そこで、（１１）式を公知の補助関数法の原理に基づいて最適化する。なお、ここでは、板倉斎藤距離を内包した、より一般的な乖離度規準であるβダイバージェンスを規準として最適化アルゴリズムを導出する。

最適化アルゴリズムの詳細な導出過程の説明は省略するが、（１１）式で表される目的関数に対して、Jensenの不等式と接線不等式とを用いて補助関数を設計することで、（１３）〜（１５）式に示す乗法更新式を得ることができる。

ただし、ρ(β)は、βの値に応じて設定される値であり、（１６）式で定義される。

（システム構成）

次に、アドホックマイクロホンアレーで取得した、伝達系への外乱及び残響成分が重畳された複数の音源の観測信号を解析して、複数の音源の各々の音源信号に分離する音響信号解析装置に本発明を適用した場合を例にして、本発明の実施の形態を説明する。

図１に示すように、本発明の実施の形態に係る音響信号解析装置は、ＣＰＵと、ＲＡＭと、後述する音響信号解析処理ルーチンを実行するためのプログラムを記憶したＲＯＭとを備えたコンピュータで構成され、機能的には次に示すように構成されている。

音響信号解析装置１００は、入力部１０と、演算部２０と、記憶部３０と、出力部４０とを備えている。

入力部１０により、伝達系への外乱及び残響成分が重畳された複数の音源ｉを含む観測信号ｙ_j[ｌ]の時系列データが入力される。記憶部３０は、入力部１０により入力された観測信号ｙ_j[ｌ]の時系列データを記憶する。また、記憶部３０は、後述する各処理での結果を記憶すると共に、本処理ルーチンで用いる各パラメータの初期値を記憶している。

演算部２０は、時間周波数解析部２１と、初期設定部２２と、パラメータ更新部２３と、終了判定部２４と、信号変換部２５とを備えている。

時間周波数解析部２１は、例えばマイクロホンｊの時系列信号としての観測された観測信号ｙ_j[ｌ]を入力として、マイクロホンｊの観測パワースペクトログラムＹ_j,k,lを計算する。また、計算したマイクロホンｊの観測パワースペクトログラムＹ_j,k,lを、記憶部３０に記憶しておく。より詳細には、時間周波数解析部２１は、例えばマイクロホンｊで観測された観測信号の時系列データを入力として、短時間フーリエ変換（Short-Time Fourier Transform:STFT）を用いて時間周波数解析を行うことにより、マイクロホンｊの観測パワースペクトログラムＹ_j,k,lを計算する。

初期設定部２２は、後述する処理で用いる各パラメータＡ_j,i,k,n、Ｗ_i,k,m,τ、及びＨ_i,m,τの各初期値を設定する。なお、各パラメータの初期値は、例えば乱数を用いて適当な値に設定すればよい。この場合、Ａ_j,i,k,n、Ｗ_i,k,m,τ、及びＨ_i,m,τの各パラメータの初期値は非負値となるように設定する。

パラメータ更新部２３は、（ｉ、ｋ、ｌ）の全ての組み合わせの各々について、記憶部３０に記憶されているＷ_i,k,m,τ及びＨ_i,m,τに基づいて、上記（８）式に従って、音源ｉのパワースペクトログラムＰ_i,k,lを計算し、記憶部３０に格納する。

また、パラメータ更新部２３は、（ｊ、ｋ、ｌ）の全ての組み合わせの各々について、記憶部３０に記憶されているＰ_i,k,l及びＡ_j,i,k,nに基づいて、上記（１２）式に従って、マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lを計算し、記憶部３０に格納する。

また、パラメータ更新部２３は、（ｊ,i,ｋ,ｎ）の全ての組み合わせの各々について、記憶部３０に記憶されているＹ_j,k,l、Ｘ_j,k,l、Ａ_j,i,k,n、及びＰ_i,k,lに基づいて、上記（１１）式の目的関数を小さくするように、（１３）式に従って、時変ステアリングベクトルの振幅成分Ａ_j,i,k,nを更新し、記憶部３０に格納する。この際、パラメータ更新部２３は、更新したＡ_j,i,k,nを用いて、上記（１２）式に従って、マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lを更新し、記憶部３０に格納する。

また、パラメータ更新部２３は、（ｉ,ｋ,ｍ,τ）の全ての組み合わせの各々について、記憶部３０に記憶されているＹ_j,k,l、Ｘ_j,k,l、Ａ_j,i,k,n、Ｗ_i,k,m,τ、及びＨ_i,m,τに基づいて、上記（１１）の目的関数を小さくするように、（１４）式に従って、音源ｉのパワースペクトログラムＰ_i,k,lの基底スペクトルＷ_i,k,m,τを更新し、記憶部３０に格納する。この際、パラメータ更新部２３は、更新したＷ_i,k,m,τを用いて、上記（８）式及び（１２）式に従って、音源ｉのパワースペクトログラムＰ_i,k,l及びマイクロホンｊのパワースペクトログラムモデルＸ_j,k,lを更新し、記憶部３０に格納する。

更に、パラメータ更新部２３は、（ｉ,ｍ,τ）の全ての組み合わせの各々について、記憶部３０に記憶されているＹ_j,k,l、Ｘ_j,k,l、Ａ_j,i,k,n、Ｈ_i,m,τ、及びＷ_i,k,m,τに基づいて、上記（１１）の目的関数を小さくするように、（１５）式に従って、音源ｉのパワースペクトログラムＰ_i,k,lの基底オンセットＨ_i,m,τを更新し、記憶部３０に格納する。この際、パラメータ更新部２３は、更新したＨ_i,m,τを用いて、上記（８）式及び（１２）式に従って、音源ｉのパワースペクトログラムＰ_i,k,l及びマイクロホンｊのパワースペクトログラムモデルＸ_j,k,lを更新し、記憶部３０に格納する。

終了判定部２４は、予め定められた終了条件を満足するか否かを判定し、終了条件を満足していない場合には、パラメータ更新部２３の各処理を繰り返す。終了判定部２４は、終了条件を満足したと判定した場合には、信号変換部２５による処理に移行する。

信号変換部２５は、記憶部３０に記憶されている音源ｉのパワースペクトログラムＰ_i,k,lの基底スペクトルＷ_i,k,m,τ、及び音源ｉのパワースペクトログラムＰ_i,k,lの基底オンセットＨ_i,m,τに基づいて、複数の音源ｉの各々について、当該音源ｉの音源信号を生成して出力部４０に出力する。出力部４０は、複数の音源ｉの各々の音源信号を出力する。

なお、終了条件としては、繰り返し回数がＬ-1回目の目的関数（１１）式の値と、繰り返し回数がＬ回目の目的関数（１１）式の値との差が、予め定めた閾値よりも小さくなったことを用いればよい。あるいは、終了条件として、繰り返し回数が、予め定められた上限回数に到達したことを用いてもよい。

（音響信号解析装置の作用）

次に、本実施の形態に係る音響信号解析装置１００の作用について説明する。Ｊ本のマイクロホンｊからなるアドホックマイクロホンアレーで取得された、伝達系への外乱及び残響成分が重畳された複数の音源ｉを含む観測信号ｙ_j[ｌ]の時系列データが、音響信号解析装置１００に入力され、記憶部３０に格納される。そして、音響信号解析装置１００において、図２に示す音響信号解析処理ルーチンが実行される。

まず、ステップＳ１００において、記憶部３０から、観測信号ｙ_j[ｌ]を読み込み、当該観測信号ｙ_j[ｌ]に対して、短時間フーリエ変換を用いた時間周波数分析を行い、マイクロホンｊの観測パワースペクトログラムＹ_j,k,lを算出すると共に、得られたマイクロホンｊの観測パワースペクトログラムＹ_j,k,lを記憶部３０に記憶する。

そして、ステップＳ１０２において、乱数を用いて、Ａ_j,i,k,n、Ｗ_i,k,m,τ、及びＨ_i,m,τの各初期値を設定して、記憶部３０に記憶する。

次に、ステップＳ１０４では、ステップＳ１０２で設定されたＷ_i,k,m,τ及びＨ_i,m,τに基づいて、上記（８）式に従って、音源ｉのパワースペクトログラムＰ_i,k,lを各（ｉ、ｋ、ｌ）の組み合わせについて算出して、記憶部３０に格納する。更に、ステップＳ１０４では、ステップＳ１０２で設定されたＡ_j,i,k,nと、本ステップで算出されたＰ_i,k,lに基づいて、上記（１２）式に従って、マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lを各（ｊ、ｋ、ｌ）の組み合わせについて算出して、記憶部３０に格納する。

ステップＳ１０６では、ステップＳ１００で算出されたＹ_j,k,lと、ステップＳ１０２で設定されたＡ_j,i,k,nと、ステップＳ１０４で算出されたＰ_i,k,l及びＸ_j,k,lに基づいて、上記（１２）式に従って、時変ステアリングベクトルの振幅成分Ａ_j,i,k,nを各（ｊ、ｉ、ｋ、ｎ）の組み合わせについて更新して、記憶部３０に格納する。

ステップＳ１０８では、ステップＳ１０６で更新されたＡ_j,i,k,nと、ステップＳ１０４で算出したＰ_i,k,lに基づいて、（１２）式に従って、マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lを各（ｊ、ｋ、ｌ）の組み合わせについて更新して、記憶部３０に格納する。

ステップＳ１１０では、ステップＳ１００で算出されたＹ_j,k,lと、ステップＳ１０８で更新されたＸ_j,k,lと、ステップＳ１０６で更新されたＡ_j,i,k,nと、ステップＳ１０２で設定されたＨ_i,m,τ及びＷ_i,k,m,τに基づいて、（１４）式に従って、音源ｉのパワースペクトログラムＰ_i,k,lの基底スペクトルＷ_i,k,m,τを各（ｉ、ｋ、ｍ、τ）の組み合わせについて更新して、記憶部３０に格納する。

ステップＳ１１２では、ステップＳ１１０で更新されたＷ_i,k,m,τと、ステップＳ１０２で設定されたＨ_i,m,τに基づいて、（８）式に従って、音源ｉのパワースペクトログラムＰ_i,k,lを各（ｉ、ｋ、ｌ）の組み合わせについて更新して、記憶部３０に格納する。また、本ステップで更新されたＰ_i,k,lと、ステップＳ１０６で更新されたＡ_j,i,k,nに基づいて、（１２）式に従って、マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lを各（ｊ、ｋ、ｌ）の組み合わせについて更新して、記憶部３０に格納する。

ステップＳ１１４では、ステップＳ１００で算出されたＹ_j,k,lと、ステップＳ１０２で設定されたＨ_i,m,τと、ステップＳ１１２で更新されたＸ_j,k,lと、ステップＳ１０６で更新されたＡ_j,i,k,nと、ステップＳ１１０で更新されたＷ_i,k,m,τに基づいて、（１５）式に従って、音源ｉのパワースペクトログラムＰ_i,k,lの基底オンセットＨ_i,m,τを各（ｉ、ｍ、τ）の組み合わせについて更新して、記憶部３０に格納する。

ステップＳ１１６では、ステップＳ１１０で更新されたＷ_i,k,m,τと、ステップＳ１１４で更新されたＨ_i,m,τに基づいて、（８）式に従って、音源ｉのパワースペクトログラムＰ_i,k,lを各（ｉ、ｋ、ｌ）の組み合わせについて更新して、記憶部３０に格納する。また、本ステップで更新されたＰ_i,k,lと、ステップＳ１０６で更新されたＡ_j,i,k,nに基づいて、（１２）式に従って、マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lを各（ｊ、ｋ、ｌ）の組み合わせについて更新して、記憶部３０に格納する。

次のステップＳ１１８では、ステップＳ１００で算出したＹ_j,k,lと、ステップＳ１１６で更新されたＸ_j,k,lに基づいて、（１１）式に従って、目的関数の値を算出して、記憶部３０に記憶する。そして、前回のステップＳ１１８で算出した目的関数の値を記憶部３０から読み込み、今回のステップＳ１１８で算出した目的関数の値と、前回のステップＳ１０８で算出した目的関数の値との差分が、予め記憶部３０に記憶されている予め定められた閾値よりも小さいか否かを判定し、差分が予め定められた閾値以上の場合には、終了条件を満足していないと判断して、上記ステップＳ１０６へ戻り、上記ステップＳ１０６〜ステップＳ１１８の処理を繰り返す。

一方、差分が予め定められた閾値未満の場合には、終了条件を満足したと判断して、ステップＳ１２０で、記憶部３０に記憶されている音源ｉのパワースペクトログラムＰ_i,k,lの基底スペクトルＷ_i,k,m,τ、及び音源ｉのパワースペクトログラムＰ_i,k,lの基底オンセットＨ_i,m,τに基づいて、複数の音源ｉの各々について、当該音源ｉの音源信号を生成して出力部４０から出力して、音響信号解析処理ルーチンを終了する。

（実施結果）

次に、本実施の形態に係る手法の有効性を示す目的で、２つの残響環境下での劣決定条件における音源分離実験を行った。

まず、残響環境下での頑健さの度合いを確認するため、異なる残響の強さを持つ環境での音源分離実験を行った。次に、伝達系への外乱に対する頑健さの度合いを確認するため、マイクロホンの位置が観測信号取得中に変化した場合における音源分離実験を行った。

なお、劣決定条件として、部屋の中に音源数を３つ、マイクロホン数を２つ設置し、部屋の形状と、音源の位置と、マイクロホンの位置から鏡像法を用いてインパルス応答を生成した。

図３に、音源分離実験を行う部屋（以降、単に「部屋」と称す）の形状と、各音源の位置及び各マイクロホンの位置を示す。なお、図３に示すように、部屋は、５ｍ×１０ｍの長方形状の大きさを有し、Ｓ１〜Ｓ３で示される“×”は音源の位置を、Ｍ１及びＭ２で示される“●”はマイクロホンの位置をそれぞれ示している。

音源分離実験で使用する音源としては、ＡＴＲ音声対話データベースの３話者１５発話を用いた。３話者の音源１〜３のうち、音源１、２の話者は女性であり、音源３の話者は男性である。

また、残響の強さは、部屋の壁の反射係数を調整することで変化させた。具体的には、部屋の壁の反射係数を0.5にした場合、残響時間は60msとなり、部屋の壁の反射係数を0.8にした場合、残響時間は210msとなる。

音源分離実験では、公知の多チャンネルNMF（Multichannel extensions of Non-negative Matrix Factorization:MNMF）を用いた手法と、本実施の形態に係る提案手法を用いた場合とについて比較した。MNMFの場合、瞬時混合モデルが仮定されているため、STFTのフレーム外に残響成分が存在する場合は、性能が低下することが予想される。

また、用意した各話者の１５発話のうち、１つの発話を分離用の信号とし、残りの１４個の発話を教師データとして事前の学習に用いた。なお、本実施の形態に係る提案手法の事前学習にはNMFD、MNMFにはNMFを用いて、基底スペクトルの学習を行い、各音源１〜３に対して、それぞれ４０個と２０個の基底を学習した。また、距離尺度としては、一般化KLダイバージェンスを用い、評価指標としては、Source-to-distortion ratio(SDR)を用いた。この際、STFTのフレーム長を32msとし、シフト長は16msとした。

図４に、部屋の壁の反射係数の変化に対する各比較手法の性能の変化の一例を示すグラフを示す。図４に示すグラフの横軸は部屋の壁の反射係数を表すと共に、縦軸はSDRを表し、グラフ５０がMNMFを用いた場合、グラフ５１が本実施の形態に係る提案手法を用いた場合のグラフを示す。なお、各マイクロホンにおける各音源のSDRを平均した値を、対応する部屋の壁の反射係数に対するSDRとして表している。

図４に示すように、反射係数の低い、瞬時混合モデルが成立しているとみなせる範囲においては、MNMFは本実施の形態に係る提案手法に比べて良い性能を示している。しかし、反射係数が高くなるにしたがって、本実施の形態に係る提案手法では、スペクトログラム素片のテンプレートとアクティベーション系列とを畳み込む混合モデルを採用していることから、残響成分を精度よく推定し、MNMFより良い性能を示す傾向が見られる。

次に、伝達系への外乱に対する頑健さの度合いを確認するため、マイクロホンの位置が観測信号取得中に変化した場合における音源分離実験について説明する。

本実験では、図３に示すマイクロホンＭ１、Ｍ２で取得した観測信号と、マイクロホンＭ１、Ｍ２の位置からそれぞれΔｘｍずらした位置で取得した観測信号とを繋げることにより、本音源分離実験で用いる観測信号を生成した。この場合、伝達系の振幅成分に比べて、位相成分には大きな外乱が生じることになる。なお、本音源分離実験における部屋の壁の反射係数は0.8に設定した。

図５に、観測信号取得位置のずれ量Δｘに対する各比較手法の性能の変化の一例を示すグラフを示す。図５に示すグラフの横軸は観測信号取得位置のずれ量Δｘを表すと共に、縦軸はSDRを表し、グラフ５２がMNMFを用いた場合、グラフ５３が本実施の形態に係る提案手法を用いた場合のグラフを示す。

図５に示すように、MNMFにおける性能を示すグラフ５２は、Δｘが大きくなるにしたがって、すなわち、外乱が大きくなるにしたがって、SDRの値が大きく落ち込んでいることがわかる。具体的には、Δｘ＝０の場合のSDRに対して、Δｘ＝０．２の場合のSDRは約８６％の値まで減少している。

一方、提案手法の場合、Δｘ＝０の場合のSDRに対して、Δｘ＝０．２の場合のSDRは約９３％の値までしか減少していない。したがって、本実施の形態に係る提案手法は、MNMFを用いる場合に比べて、伝達系への外乱に対しての音源分離性能の低下が少なく、伝達系への外乱に対する頑健さの度合いが高いということができる。

このように、本発明に係る提案手法では、音源或いはマイクロホンの微小移動等、音響信号解析環境の軽微な変化に対して半時変形モデルを設定し、複数の時間フレームに亘るスペクトルを連結したスペクトログラム素片のテンプレートとアクティベーション系列とを畳み込む混合モデルを適用することで、音源及びマイクロホンの相対位置関係が変化する時変残響環境下であっても、複数の音源の音声が重畳した観測信号から、各音源の音源信号を精度よく分離することができる。

なお、本発明は、上述した実施形態に限定されるものではなく、この発明の要旨を逸脱しない範囲内で様々な変形や応用が可能である。

例えば、上述の音響信号解析装置は、内部にコンピュータシステムを有しているが、「コンピュータシステム」は、ＷＷＷシステムを利用している場合であれば、ホームページ提供環境（あるいは表示環境）も含むものとする。

また、本願明細書中において、プログラムが予めインストールされている実施形態として説明したが、当該プログラムを、コンピュータ読み取り可能な記録媒体に格納して提供することも可能である。

１０ 入力部
２０ 演算部
２１ 時間周波数解析部
２２ 初期設定部
２３ パラメータ更新部
２４ 終了判定部
２５ 信号変換部
３０ 記憶部
４０ 出力部
１００ 音響信号解析装置

Claims (5)

1. Ｊ本の各マイクロホンｊで集音した音響信号の時系列データを入力として、各時刻ｌにおける各周波数ｋの観測信号時間周波数成分ｙ_j,k,lを出力する時間周波数解析部と、
   音源ｉから前記マイクロホンｊまでの時刻ｎだけ遅延して集音される伝達特性を表す時変ステアリングベクトルの振幅成分Ａ_j,i,k,n、複数フレームのスペクトルを連結したスペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底スペクトルＷ_i,k,m,τ、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底オンセットＨ_i,m,τの各々に初期値を設定する初期値設定部と、
   （j,k,l）の全ての組み合わせにおける、前記観測信号時間周波数成分ｙ_j,k,lと、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底スペクトルＷ_i,k,m,τ、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底オンセットＨ_i,m,τ、各音源ｉ及び各時刻ｎの前記振幅成分Ａ_j,i,k,nに基づいて算出される前記マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lとの距離が小さくなるように、前記基底スペクトルＷ_i,k,m,τと、前記基底オンセットＨ_i,m,τと、前記振幅成分Ａ_j,i,k,nとを更新するパラメータ更新部と、
   予め定められた終了条件を満たすまで、前記パラメータ更新部による更新を繰り返し行う終了判定部と、
   を含む音響信号解析装置。
2. 前記観測信号時間周波数成分ｙ_j,k,lと、前記マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lとの距離の尺度として板倉斎藤距離を用いた
   請求項１記載の音響信号解析装置。
3. 時間周波数解析部と、初期値設定部と、パラメータ更新部と、終了判定部とを含む音響信号解析装置における音響信号解析方法であって、
   前記時間周波数解析部が、Ｊ本の各マイクロホンｊで集音した音響信号の時系列データを入力として、各時刻ｌにおける各周波数ｋの観測信号時間周波数成分ｙ_j,k,lを出力し、
   前記初期値設定部が、音源ｉから前記マイクロホンｊまでの時刻ｎだけ遅延して集音される伝達特性を表す時変ステアリングベクトルの振幅成分Ａ_j,i,k,n、複数フレームのスペクトルを連結したスペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底スペクトルＷ_i,k,m,τ、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底オンセットＨ_i,m,τの各々に初期値を設定し、
   前記パラメータ更新部が、（j,k,l）の全ての組み合わせにおける、前記観測信号時間周波数成分ｙ_j,k,lと、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底スペクトルＷ_i,k,m,τ、前記スペクトル素片に対応する各時刻τにおける、各基底ｍ及び各周波数ｋの非負値の基底オンセットＨ_i,m,τ、各音源ｉ及び各時刻ｎの前記振幅成分Ａ_j,i,k,nに基づいて算出される前記マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lとの距離が小さくなるように、前記基底スペクトルＷ_i,k,m,τと、前記基底オンセットＨ_i,m,τと、前記振幅成分Ａ_j,i,k,nとを更新し、
   前記終了判定部が、予め定められた終了条件を満たすまで、前記パラメータ更新部による更新を繰り返し行う
   音響信号解析方法。
4. 前記観測信号時間周波数成分ｙ_j,k,lと、前記マイクロホンｊのパワースペクトログラムモデルＸ_j,k,lとの距離の尺度として板倉斎藤距離を用いた
   請求項３記載の音響信号解析方法。
5. コンピュータを、請求項１又は請求項２記載の音響信号解析装置の各部として機能させるためのプログラム。