(57)【要約】【課題】残響や背景雑音が存在し、観測される特徴量が明確にクラスタリングできず、重なり合うような場合に対しても、音源定位／音源分離性能を向上させることができる音源分離定位方法を提供する。【解決手段】ブラインド音源分離（ＢＳＳ）の問題に対し、ＥＭアルゴリズムを適用した新しいアルゴリズムを提案する。具体的には、音源信号のスパース性に基づき、各時間周波数成分に寄与する音源は１個であるという観測モデルの下、各時間周波数成分に寄与している音源のインデックスを隠れ変数とみなし、最大尤度を与える音源方向と、各時間周波数成分への各音源の寄与率をＥＭアルゴリズムによって推定する。その結果、残響時間が多い場でも音源の分離定位を行うことが可能となった。【選択図】図１
2008145610.tif 000002
【特許請求の範囲】
【請求項１】
複数チャネルの観測信号に基づき、複数の音源からの信号を分離する音源分離定位方法において、
前記複数チャネルの観測信号のそれぞれを、時間周波数領域の信号に変換するステップと、
音源定位の初期値を決定する初期音源定位ステップと、
前記音源定位に基づき、前記時間周波数領域上において、各音源の信号を分離する音源分離ステップと、
前記分離結果に基づき、前記時間周波数領域上で各成分が属する音源の尤度を求め、この尤度の期待値である補助関数を最大値にする音源定位を求めて、音源の定位を確率的に行う繰り返し音源定位ステップと、
前記音源分離ステップと、前記繰り返し音源定位ステップと、を繰り返し実行する繰り返しステップと、
を含むことを特徴とする音源分離定位方法。
【請求項２】
請求項１記載の音源分離定位方法において、
前記音源定位は、前記音源の位置、前記音源の方向、前記音源からの音波の前記複数チャネル間の時間差、のいずれか１種又は２種以上を含むことを特徴とする音源分離定位方法。
【請求項３】
請求項１記載の音源分離定位方法において、
前記初期音源定位ステップは、
各観測信号の間の特徴量に基づきクラスタリングを行うステップと、
クラスタリングの結果に基づき各時間周波数成分ごとにいずれの音源に属するかを決定するステップと、
上記決定に基づきそれぞれの音源の方向を求めるステップと、
を含むことを特徴とする音源分離定位方法。
【請求項４】
請求項１記載の音源分離定位方法において、
前記音源分離ステップは、各時間周波数成分がどの音源に属するかを表す尤度に対する寄与を表す分配関数を求めることによって、音源分離を行うことを特徴とする音源分離定位方法。
【請求項５】
請求項４記載の音源分離定位方法において、
前記繰り返し音源定位ステップは、
前記分配関数を、対応する各時間周波数成分の尤度に乗算し、その結果を、全時間・全周波数にわたって加算する補助関数を計算するステップと、
前記分配関数のパラメータを所定の数値範囲で走査し、前記補助関数の値が最大となるようなパラメータを求めるステップと、
を含み、前記新たなパラメータに基づき、前記音源分離ステップにおいて尤度が求められることを特徴とする音源分離定位方法。
【請求項６】
請求項５記載の音源分離定位方法において、
前記パラメータ前記複数の音源方向の組であることを特徴とする音源分離定位方法。
【請求項７】
請求項５記載の音源分離定位方法において、
前記補助関数は、時間周波数領域上の関数であることを特徴とする音源分離定位方法。
【請求項８】
請求項１記載の音源分離定位方法において、
前記繰り返しステップは、定位する音源の方向の変化が所定量未満になった場合に、繰り返しを終了することを特徴とする音源分離定位方法。
【請求項９】
請求項５記載の音源分離定位方法において、
前記繰り返しステップは、前記補助関数の値の増加が所定量未満になった場合に、繰り返しを終了することを特徴とする音源分離定位方法。
【請求項１０】
Ｎチャネルの観測信号に基づき、Ｍ個の音源からの信号を分離する音源分離定位方法において、
前記Ｎチャネルの観測信号のそれぞれを、時間周波数領域の信号に変換するステップと、
音源定位の初期値を決定する初期音源定位ステップと、
前記音源定位に基づき、前記時間周波数領域上において、各音源の信号を分離する音源分離ステップと、
前記分離結果に基づき、前記時間周波数領域上で各成分が属する音源の音源方向の尤度を求め、この尤度の期待値である補助関数を最大値にする前記音源間の音源方向の値を確率的に求めることによって音源定位を行う繰り返し音源定位ステップと、
前記音源分離ステップと、前記繰り返し音源定位ステップと、を繰り返し実行する繰り返しステップと、
を含むことを特徴とする音源分離定位方法。ここで、前記Ｎ及びＭは２以上の整数である。
【請求項１１】
請求項１０記載の音源分離定位方法において、
前記音源定位は、前記各音源の方位、又は、前記記Ｎチャネル間の観測信号の時間差、強度比、誤差分散のいずれか１種又は２種以上を含むことを特徴とする音源分離定位方法。
【請求項１２】
請求項１０記載の音源分離定位方法において、
前記初期音源定位ステップは、
各観測信号の間の時間差又は強度比に基づきクラスタリングを行うステップと、
クラスタリングの結果に基づき各時間周波数成分ごとにいずれの音源に属するかを決定するステップと、
上記決定に基づきそれぞれの音源の方位として、前記時間差又は強度比を求めるステップと、
を含むことを特徴とする音源分離定位方法。
【請求項１３】
請求項１０記載の音源分離定位方法において、
前記音源分離ステップは、各時間周波数成分が属する音源の前記時間差又は強度比を表す尤度に対する寄与を表す分配関数を求めることによって、音源分離を行うことを特徴とする音源分離定位方法。
【請求項１４】
請求項１３記載の音源分離定位方法において、
前記繰り返し音源定位ステップは、
前記分配関数を、対応する各時間周波数成分の尤度に乗算し、その結果を、所定の数値範囲にわたって加算する補助関数を計算するステップと、
前記分配関数のパラメータである時間差又は強度比を所定の数値範囲で走査し、前記補助関数の値が最大となるような時間差又は強度比を求めるステップと、
を含み、前記新たな時間差又は強度比に基づき、前記音源分離ステップにおいて尤度が求められることを特徴とする音源分離定位方法。
【請求項１５】
請求項１４記載の音源分離定位方法において、
前記補助関数は、時間周波数領域上の関数であることを特徴とする音源分離定位方法。
【請求項１６】
請求項１０記載の音源分離定位方法において、
前記繰り返しステップは、定位する音源の方向の変化が所定量未満になった場合に、繰り返しを終了することを特徴とする音源分離定位方法。
【請求項１７】
請求項１４記載の音源分離定位方法において、
前記繰り返しステップは、前記補助関数の値の増加が所定量未満になった場合に、繰り返しを終了することを特徴とする音源分離定位方法。
【請求項１８】
コンピュータに、複数チャネルの観測信号に基づき、複数の音源からの信号を分離する音源分離定位方法を実行させるためのプログラムにおいて、
前記複数チャネルの観測信号のそれぞれを、時間周波数領域の信号に変換する手順と、
音源定位の初期値を決定する初期音源定位手順と、
前記音源定位に基づき、前記時間周波数領域上において、各音源の信号を分離する音源分離手順と、
前記分離結果に基づき、前記時間周波数領域上で各成分が属する音源の尤度を求め、この尤度の期待値である補助関数を最大値にする音源定位を求めて、音源の定位を確率的に行う繰り返し音源定位手順と、
前記音源分離ステップと、前記繰り返し音源定位ステップと、を繰り返し実行する繰り返し手順と、
を、前記コンピュータに実行させることを特徴とするプログラム。
【請求項１９】
コンピュータに、Ｎチャネルの観測信号に基づき、Ｍ個の音源からの信号を分離する音源分離定位方法を実行させるプログラムにおいて、
前記Ｎチャネルの観測信号のそれぞれを、時間周波数領域の信号に変換する手順と、
音源定位の初期値を決定する初期音源定位手順と、
前記音源定位に基づき、前記時間周波数領域上において、各音源の信号を分離する音源分離手順と、
前記分離結果に基づき、前記時間周波数領域上で各成分が属する音源の音源方向の尤度を求め、この尤度の期待値である補助関数を最大値にする前記音源間の音源方向の値を確率的に求めることによって音源定位を行う繰り返し音源定位手順と、
前記音源分離ステップと、前記繰り返し音源定位ステップと、を繰り返し実行する繰り返し手順と、
を、前記コンピュータに、実行させることを特徴とするプログラム。ここで、前記Ｎ及びＭは２以上の整数である。
【請求項２０】
請求項１又は１０記載の音源分離定位方法において、
前記繰り返し音源定位ステップは、さらに、
前記補助関数をσ２で微分した更新式を計算することによって雑音の分散σ２を求めるステップ、
を含むことを特徴とする音源分離定位方法。
【請求項２１】
請求項１又は１０記載の音源分離定位方法において、
前記繰り返しステップを終了した際に得られた分配関数に基づき、最も大きい寄与の音源のみを通過させる時間周波数マスキングによって、最終的な音源分離を行なう分離ステップを含むことを特徴とする音源分離定位方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、いわゆるブラインドの音源定位、音源分離技術に関する。
【背景技術】
【０００２】
１。序論
近年、実環境での音声認識やロボット聴覚などへの応用を目的として、複数の音源信号が混合した観測信号から、元の音源信号を分離する技術が着目されている。コンピュータによる音声認識や音環境の理解のために目的とする音声だけを分離する技術が強く望まれている。
【０００３】
特に、音源から観測点への伝達関数が未知の場合でも分離を実現するブラインド音源分離（Ｂｌｉｎｄ  Ｓｏｕｒｃｅ  Ｓｅｐａｒａｔｉｏｎ、  ＢＳＳ）の研究が活発に行えわれている。ブラインドとは、簡単に言えば、マイクロフォン等のセンサー側から話者の情報がないことを言う。もし、話者の方向がわかっている場合は、ガンマイクロフォン等の指向性の鋭いマイクロフォンを用いれば比較的簡易に特定の話者の音声のみを取得することができる。
【０００４】
ＢＳＳの有効な手法の一つとして独立性分分析が知られているが、この分析方法は、音源信号数が観測信号数以下の場合にしか適用できない。簡単に言えば、音源数よりマイクロフォンの数が多い場合にはこの独立成分分析が適用できる。
【０００５】
一方、観測信号より多くの音源信号を分離する枠組みとしては、音源信号のスパース性を利用した手法が提案されている。
【０００６】
この手法は、下記非特許文献１、非特許文献３、非特許文献６、非特許文献７、非特許文献９、非特許文献１１等に詳しい。
【０００７】
１．１  スパース性
対象とする信号のエネルギーがある領域（多くの場合、時間周波数領域）で一部の領域に集中し、その他の多くの領域でほぼ０であるような性質がある場合、それを信号のスパース性と呼ぶ。特に音声の場合は、フォルマント構造や調波構造、有声／無声の時間的変動などの要因から、複数の音声が混合している場合でも、エネルギーが集中している領域が重なり合うことは少ないことが知られている（下記非特許文献２参照）。
【０００８】
音声信号の場合、簡単に言えば、スパース性とは、時間周波数領域で音声の有意なエネルギーがまばらにしか存在しない性質を言う。
【０００９】
１．２  時間周波数マスキング
この音声のスパース性を利用した音源分離の代表的な先行研究として、時間周波数マスキングが挙げられる。これは目的の音源信号成分が支配的である時間周波数成分を、マスキング処理によって抜き出すことによって分離する手法である。
【００１０】
音声信号を、時間−周波数領域において解析が行えわれる場合、時間分解能は、およそ６０ｍｓｅｃ程度に選ばれることが多い。一方周波数分解能は、６０ｍｓｅｃの逆数、すなわちおよそ１６６Ｈｚ程度の分解能となる。
【００１１】
この時間分解能や周波数分解能は種々の値を選ぶことが可能であり、選び方に関する様々な研究もなされている。
【００１２】
本文では、一般的な６０ｍｓｅｃ、１６６Ｈｚ程度の分解能を例として説明を行うが、他の分解能を利用してもかまわない。
【００１３】
このように、時間−周波数領域を解析する場合、時間軸で約６０ｍｓｅｃ単位で、周波数軸で約１６６Ｈｚ単位で区分して解析を行う。この各区分（時間周波数成分）に対して、どの音源からの信号かを判断する。これは音声のスパース性を前提としているからである。そして、各区分（の信号）がどの音源の信号かを全て判断できれば、所望の音源からの信号が含まれる区分のみを取り出せば、所望の音源からの信号のみを取り出す（音源分離する）ことができる。
【００１４】
このような所望の信号のみを取り出すためのマスクは、各時間周波数成分ごとに「通過させる（取り出す）」か、「通過させない（除去する）」かを決めたマスクが用いられる。このマスクは「通過させる・させない」という２値のマスクであるので、バイナリマスクと呼ばれる。このようなマスキングの手法が時間周波数マスキングである。
【００１５】
実際には、各時間周波数成分（区分）において、複数の音源の信号がほぼ同じ比率で混在しており、いずれかの音源の音であるかを断定できない場合もあり得るが、そのような混在している時間周波数成分は一般に除いてしまっても音源分離した信号の品質（音声が聞き取れるか否か）にはあまり影響を及ぼさないことが知られている。
【００１６】
いずれの音源の音であるかを決定できた時間周波数成分についてのみ取り出すようにすれば、聴覚上、その音源からの音をほぼ聞き取れることが知られている。
【００１７】
さて、各時間周波数成分において、目的音源が支配的であるかどうかの判定には、２ｃｈの観測信号間の時間差δ（τ、ω）及び強度比ａ（τ、ω）が多く用いられる。各時間周波数成分における時間差は信号の到来方向と対応する特徴量なので、この特徴量に基づきクラスタリングすれば同じ音源から到来した成分だけを抜き出し、時間周波数平面のマスクを作ることができる。この技術が下記非特許文献１、非特許文献３、非特許文献８に記載されている。
【００１８】
時間周波数マスキングに関しては既に多くの研究がなされており、リアルタイムでの分離がＢａｅｃｋらにより研究されている（下記非特許文献１２）。また、その性能の上限や、マスキングの特徴であるミュージカルノイズを低減する時間周波数分解について荒木らが論じている（下記非特許文献４、非特許文献５）。
【００１９】
また、単純な時間周波数マスキング以上の性能を示す分離手法として、各時間周波数成分に寄与する音源の数が複数（観測信号数以下）であると仮定し、混合行列を推定して分離する手法（下記非特許文献７）が知られている。また、音源信号に事前分布を仮定して確率的に信号成分を推定する手法も提案されている（下記非特許文献９、非特許文献１１）。
【００２０】
しかしこれら従来の多くの手法では、残響や背景雑音が多く存在する環境下では、各音源から生じる強度比や時間差がばらついて互いに重なり合うため、特徴空間でのクラスタリングが困難になることが大きな問題であった。
【００２１】
１．３  混合ガウス分布の推定
一方、音声認識をはじめとした幅広い分野で用いられている手法として、混合ガウス分布モデルの推定手法がある。
【００２２】
このモデルを利用した推定問題においては、複数の分布が混在している状況で、分布の各パラメータの推定が行えわれる。このモデルは、主に、分布同士が互いに重なり合うような状況において用いられる。
【００２３】
このような場合においても、混合ガウス分布モデルによれば分布のパラメータの推定が精度良く行えわれることが知られている。
【００２４】
特に、混合ガウス分布モデルの特徴の一つとして、帰属度が連続値であることが挙げられる。上で述べた通常のクラスタリング手法では、各データをいずれかのクラスに０、１（離散値）で帰属させている。これに対して、混合ガウス分布モデルではこの連続値である帰属度を採用し、モデルフィッティングにより、全体にとして尤度が最大になるような解を求めている。
【００２５】
このような特徴が、上述した重なり合う分布の状況下においても効率の良い推定を可能にしている一つの鍵となっている。さらにまたこのモデルに関しては、ＥＭアルゴリズムと呼ばれる効率的な解法が存在することも、このモデルが良く採用される要因の一つと考えられる。
【００２６】
１．４  本発明の起点
本願発明者らは、この考え方（混合ガウス分布モデルの推定手法）を上記スパース性に基づくＢＳＳに適用することを鋭意検討し、本願発明をなすに至った。本発明によれば、残響環境下のように従来のＢＳＳでは特徴空間でのクラスタリングが困難になる状況でも、適切に時間周波数マスクを設計することが可能となった。本発明では、このような新しい手法を提案する。
【００２７】
なお、本発明では、従来は別々の処理であったクラスタリングによる音源定位と、マスキングによる音源分離と、を共通の目的関数を最大化させるという統一的な処理で行うことができることも利点の一つと考えられる。
【００２８】
以下、本発明の基礎となるスパース性を利用したＢＳＳを簡単に説明する。
【００２９】
２．スパース性に基づくＢＳＳの概要
２．１  時間周波数マスキング
スパース性に基づくＢＳＳで代表的な手法は時間周波数マスキングである（非特許文献１、非特許文献３参照）。今、Ｎ個の音源信号Ｓｎ（τ、ω）（ｎ＝１、。。。、Ｎ）が混合した観測信号Ｍ（τ、ω）の各時間周波数（τ、ω）において、最もエネルギーの大きい音源信号のインデックスをｋ（τ、ω）とする。ここで、Ｓｎは、短時間複素フーリエ変換（ＳＴＦＴ）の値である。τは時間であり、ωは角周波数である。
【００３０】
さて、時間周波数マスキングとは、理想的には
【数１】
2008145610.tif 000003
と表されるようなマスクを設計し、
【数２】
2008145610.tif 000004
のように観測信号に乗じることで、目的信号Ｓｊ（τ、ω）を推定するという手法である。
【００３１】
時間周波数マスキングは原理的に、個々の時間周波数成分を「通過させる」、「阻止する」、のいずれかしかできないため、理想的なマスクが設計できた場合の分離性能は、対象とする信号のスパース性の度合いに依存する。
【００３２】
これは、単なる「通過させる」「阻止する」というバイナリマスクではなく、部分的な通過を許す連続値マスクを用いたとしても同様である。
【００３３】
また実際には、得られた観測信号からいかにしてこのようなマスクを設計するかが重要な問題となる。
【００３４】
できるだけ信号をスパースに表現するための時間周波数分解法としては、これまでは短時間Ｆｏｕｒｉｅｒ変換を前提とし、音声を対象としたときの最適な窓幅が主に論じられてきた。例えば１６ｋＨｚサンプリングの場合には、１０２４点の窓幅が最も良いとの報告がある（下記非特許文献１参照）。また、短時間Ｆｏｕｒｉｅｒ変換に代えて、フィルタバンク分析を行えった場合には、聴覚的なフィルタバンクを用いたほうが良いとの報告もある（下記非特許文献１０参照）。
【００３５】
２．２  時間周波数マスクの従来の設計法
時間周波数マスクの設計法として、従来から用いられている設計法は、観測信号間の時間差・強度比のクラスタリングである。例えば２ｃｈの場合には、観測信号の各時間周波数成分の比ＭＲ（τ、ω）／ＭＬ（τ、ω）から、音源位置に依存した特徴量（強度比ａ、時間差δ）が抽出できる。２ｃｈであるので、右と左の１対のマイクロフォンを用いて２個の観測信号を得る。右側の観測信号をＭＲ（τ、ω）と表し、右側の観測信号をＭＬ（τ、ω）と表している。
【００３６】
これらの特徴量は、音源が移動しない限り、同一音源から到来した成分は時間、周波数にかかわらず、ほぼ同じ値をとると期待できるので、これらのクラスタリングにより時間周波数マスクを設計することができる。
【００３７】
Ｙｉｌｍａｚらは、観測信号の各時間周波数成分からａ、δを求め、その分布全体から複数の音源位置に対応する強度比ａｉ、時間差δｉを決定し、次に各時間周波数成分に対して尤度最大になる音源を選ぶクラスタリングを行えい、時間周波数マスクを設計する手法を提案している（下記非特許文献１参照）。ここで、ｉは、音源のインデックスである。
【００３８】
２．３  残響環境下での問題点
Ｙｉｌｍａｚらは、無響環境であれば時間周波数マスキングにより十分に分離が可能であるが、残響時間が長くなると、ａ、δは大きな分散を生じ、音源定位自体が困難になると報告している（下記非特許文献１参照）。残響時間が異なる環境で、３つの音源が存在する際、２個のマイクロフォンで観測される時間差の散布図が、図４、図５、図６に示されている。図４は残響時間０ｍｓｅｃであり、図５は残響時間５０ｍｓｅｃであり、図６は残響時間１７０ｍｓｅｃの場合を示している。これらは全時間周波数成分において位相差から時間差を検出し、縦軸にその成分のパワーをとってプロットしたものである。図４、図５、図６からわかる通り、残響時間が０の場合は、概ね３個のクラスタリングが可能である。しかし、残響時間が長くなるにつれて、多重反射の影響で時間周波数ごとに様々な時間差が生じ、クラスタリングが困難になっていく様子が図４、図５、図６に明確に示されている。
【００３９】
３．  先行特許文献
下記特許文献１には、マイクロフォン数以上の音源を扱うことができる音源定位と音源分離とを実現する技術が開示されている。この文献では、周波数と位相差の２次元データから、予め定められた図形を検出し、各図形に基づいて音源候補に対する情報を得るとされている。
【００４０】
下記特許文献２には、２本のマイクロフォンを用いて複数の音源を扱うことができる音源定位と音源分離とを実現する技術が開示されている。この文献では、上記特許文献１と同様に、周波数と位相差の２次元データから、予め定められた図形を検出し、各図形に基づいて、音源の数、音源の存在する範囲、音源の音の存在する時間、音源の成分、音源ごとの分離音声、等を得るとされている。
【００４１】
下記特許文献３には、ノイズの影響のある環境で、ＢＳＳを用いたバイナリマスキング処理と、バイノーラル信号分離処理と、を組み合わせて、音源分離を行う装置が開示されている。
【００４２】
下記特許文献４には、複数の信号が混合されている信号から、時間領域ＢＳＳを用いて音源分離し、音源分離した信号をサブバンド合成して原信号に対する信号を得る手法が開示されている。
【００４３】
【特許文献１】特開２００６−２６７４４４号公報
【特許文献２】特開２００６−２５４２２６号公報
【特許文献３】特開２００６−１５４３１４号公報
【特許文献４】特開２００３−２７１１６８号公報
【非特許文献１】O. Yilmaz and S. Rickard: "Blind Separation of Speech Mixtures via Time-Frequency Masking," IEEE Transaction on Signal Processing, Vol. 52, No. 7, pp 1830-1847, (2004)
【非特許文献２】S. Rickard and O. Yilmaz: "On the Approximate W-disjoint Orthogonality of Speech, " Proc. ICASSP, Vol. I, pp. 529-532, (2002)
【非特許文献３】S. Araki, H. Sawada, R. Mukai, S. Makino: "DOA Estimation for Multiple sparse sources with normalized observation vector clustering," ICASSP, Vol. V, pp 33-36 (2006)
【非特許文献４】S. Araki, S. Makino, H. Sawada, and R. Mukai: "Reducing Musical Noise by a Fine-Shift  Overlap-add Method Applied to Source Separation using a Time-Frequency Mask," Proc. ICASSP, vol. III, pp. 81-84, (2005)
【非特許文献５】S.Araki, R. Mukai, S. Makino, T. Nishikawa and H.Saruwatari,"The Fundamental Limitation of Frequency Domain Blind Source Separation for Convolutive Mixtures of Speech," IEEE Trans. on Speech Audio Processing, Vol. 11, No. 2, pp. 109-116 (2003)
【非特許文献６】L. Vielva, D. Erdogmus, C. Pantaleon, I. Santamaria, J. C.Principe: "Underdetermined Blind Source Separation in a Time-Varing environment," Proc. ICASSP Vol. III, pp3049-3052, (2002)
【非特許文献７】A. Blin, S. Araki and S. Makino: "A Sparseness-Mixing Matrix Estimation (SMME) Solving the Underdetermined BSS for Convolutive Mixtures," Proc. ICASSP, Vol. IV, pp85-88, (2004)
【非特許文献８】S. Winter, H. Sawada, S. Araki and S. Makino: "Overcomplete BSS for Convolutive Mixtures Based on Hierarchical Clustering," Proc. SAPA2004, S1.3, (2004)
【非特許文献９】S. Winter, H. Sawada, S. Makino: "On Real and Complex Valued L1-norm Minimization for Overcomplete Blind Source Separation," Proc. WASPAA2005, pp. 86-89, (2005)
【非特許文献１０】小野, 和泉, 嵯峨山, "音声のスパース性を最大化するフィルタバンクの検討," 日本音響学会2006 年春季研究発表会講演論文集, 1-5-23, pp. 551-552, (2006)
【非特許文献１１】C. Fevotte and S. J. Godsill, "A Bayesian Approach for Blind Separation of Sparse Sources," IEEE Trans. on Speech and Audio Processing, Vol. 14,
【非特許文献１２】M. Baeck and U. Zolzer, "Real-time Implementation of a source separation algorithm," Proc. of the 6th Int. Conference on Digital Audio Effects(DAFx-03), (2003)
【発明の開示】
【発明が解決しようとする課題】
【００４４】
本願発明は、このような背景に鑑みなされたものであり、ノイズの影響下であっても、音源定位・音源分離を行える手法を実現することである。特に、音源数がマイクロフォン数より多くても適用可能なブラインド音源定位・分離の実現を目的とする。
【課題を解決するための手段】
【００４５】
上述した課題に鑑み、本願発明者は、一意のクラスタリングではなく、確率による定式化を検討した（ポイント１）。さらに、同じ領域での繰り返し処理を行うことにした（ＥＭアルゴリズムの採用）（ポイント２）。以下、詳細に説明する。
【００４６】
４．    本発明の概要
４．１  本発明の特徴
今、ある時間周波数（τ、ω）で観測信号
【数３】
2008145610.tif 000005
が得られたとき、これがある音源方向θから到来した信号である確率を
【数４】
2008145610.tif 000006
と表す。なお、この観測信号は、時間周波数領域上の複素ベクトルである。本特許において太字はベクトルを表す。
【００４７】
以下、上記確率を単一方向尤度と呼ぶ。方向尤度を定めることができたならば、音源が１個の場合には、全観測データに対する対数尤度の和である下記（３）式
【数５】
2008145610.tif 000007
を最大化することによって、最尤音源方向θＭＬを求めることができる。
【００４８】
さて、本４章では、音源の数をＮと表す。便宜上他の章では異なる表記をする場合もある。
【００４９】
音源がこのようにＮ個存在する場合においても、音源信号がスパースであり、各時間周波数成分に寄与する音源がたかだか１個であるとみなせるならば、ｎ（ｎは１からＮまでの整数）番目の音源方向θｎは、ｎ番目の音源が寄与する時間周波数（τ、ω）の集合Ωｎに対する対数尤度の和である下記（４）式
【数６】
2008145610.tif 000008
を最大化することで推定される。ただし、Ωｎを求めること自体が時間周波数マスキングによる音源分離そのものであることに注意する。つまり、各時間周波数成分に寄与する音源がたかだか１個であるという仮定の下では、以下のことが言える。
【００５０】
１）Ωｎが求まれば（音源分離できれば）θｎが求まる（音源定位できる）
２）θｎ（ｎ＝１、。。。、Ｎ）が求まれば（音源定位できれば）、Ωｎ（ｎ＝１、。。。、Ｎ）が求まる（音源分離できる）。
【００５１】
すなわち、これらはみな相互に関連した関係にある。従来の枠組みでは、強度比・時間差などの特徴量検出後、その特徴空間における投票法やｋ−ｍｅａｎｓ法などのクラスタリングにより音源定位が先に行えわれ、その後、時間周波数マスクが設計されるという２段階の処理が行えわれることが多かった。
【００５２】
この問題はそもそも、時間周波数領域を個々の音源に帰属する成分に分けるクラスタリングの問題ととらえることができるが、クラスタリングの分野で良く扱われる類似の問題として、混合ガウス分布モデル（ＧＭＭ）の推定問題が知られている。この問題は、各データは複数のガウス分布のいずれかから生成されるが、どのガウス分布から生成されたものかはわからない、という条件下で、各ガウス分布の平均・分散を推定する問題である。この混合ガウス分布モデルの推定問題は、
１）データを
【数７】
2008145610.tif 000009
２）複数のガウス分布を各音源に対応する方向尤度分布
３）ガウス分布の平均・分散を音源方向やその他方向尤度を決定するパラメータ
と、それぞれ読み替えれば、本件の音源定位・音源分離の問題と全く同種の問題であることがわかる。
【００５３】
従来の通常のクラスタリング手法では、与えられたデータが個々のクラスに属しているか属していないかを０（属さない）、１（属する）で決定するのに対し、混合ガウス分布のモデルは、帰属率を連続値の確率として扱うため、分布同士が重なり合っているような場合でも、ロバストに推定を行うことができる。この考え方をスパース性に基づくＢＳＳに導入することによって、残響や背景雑音が存在し、観測される特徴量が明確にクラスタリングできず、重なり合うような場合に対して、音源定位／音源分離性能を向上させることができると、本願発明者らは考え、本発明をなすに至った。
【００５４】
すなわち、残響時間が多い場面においては、従来の単純なクラスタリングで音源を分離することは極めて困難であった。このような状況下では、上で述べたように、ある時間周波数成分に対して一意に音源を決めることはできないとして扱う枠組みが必要であると本願発明者らは考えたのである（着眼点１）。
【００５５】
また、このような混合分布の推定問題に対し最尤解を求めるための効率的な手法がＥＭアルゴリズム（Expectation Maximization Algorithm）として知られている。以下では各時間周波数成分に寄与する音源がたかだか１個であるというモデルの下で、ＢＳＳの問題を、各時間周波数成分に寄与する音源がどれであるかを隠れ変数とした最尤問題として定式化する。
【００５６】
４．２  ＥＭアルゴリズムによる定式化
本特許でまず扱う問題は、下記（５）式
【数８】
2008145610.tif 000010
を最大化する音源方向の組
【数９】
2008145610.tif 000011
を求めることである。これは方向の組であり、ベクトルである。本特許においては太字の記号はベクトルを表す。この上記式（５）（数８参照）が、ＥＭアルゴリズムにおける「目的関数」である。式（５）のｐ（尤度）は、下記式（６）（下記数１３参照）のように表される。また、音源方向である上記数９が推定したいパラメータであり、各時間周波数での音源のインデックスｋが隠れ変数となる。
【００５７】
ここで
【数１０】
2008145610.tif 000012
は、音源が
【数１１】
2008145610.tif 000013
方向に存在するときに、
【数１２】
2008145610.tif 000014
が観測される尤度である。これは、各時間周波数成分に寄与する音源が１個であるというモデルの下では、下記式（６）
【数１３】
2008145610.tif 000015
のように、周辺化して表すことができる。ここでｋ（τ、ω）は、（τ、ω）成分に寄与する音源のインデックスであり、実際には観測することができない隠れ変数である。混合ガウス分布の場合と比較すると、あるデータがいずれのガウス分布からの出力であるか、を示すガウス分布の番号にちょうど対応している。
【００５８】
もしｋ（τ、ω）が既知であったとすると、上記式（５）の尤度はｋ番目の音源方向にしか依存しないため、式（７）
【数１４】
2008145610.tif 000016
となり、前小節（３．１節）で議論した単一方向尤度と一致する。
【００５９】
ＥＭアルゴリズムでは、このような隠れ変数を含んだ最尤問題を、仮のパラメータ（ここでは、仮の音源方向）
【数１５】
2008145610.tif 000017
を使って定義されるＱ関数と呼ばれる補助関数を導入し、次のようなＥステップ・Ｍステップ
【数１６】
2008145610.tif 000018
という２個のステップの反復、すなわち式（８）
【数１７】
2008145610.tif 000019
の算出を繰り返すことで、パラメータの逐次推定を行う。ここで、本問題におけるＱ関数は、下記式（９）
【数１８】
2008145610.tif 000020
のように与えられる。ただし、この式（９）においては下記式（１０）（１１）
【数１９】
2008145610.tif 000021
のように与えられる。
【００６０】
ここで、「仮の」パラメータとは、音源方向を「仮に」決めたという意味である。このように音源方向を仮に決めてそれを用いて音源位置を定位し、音源位置に基づき音源方向を又定めるということを繰り返す。繰り返しの最初の初期値は、真の値ではなく、推定値であるので、「仮パラメータ」と呼んでいる。
【００６１】
上記式（９）（数１８）のように、Ｑ関数が個々のθｋのみに依存する関数の和に分解されるということは、本問題においては、複数音源定位が複数の単一音源定位問題に分解されることを意味している。Ｑ関数（補助関数）は、方向の尤度の対数（対数尤度）のいわば期待値である。
【００６２】
Ｅステップで計算される
【数２０】
2008145610.tif 000022
は分配関数とも呼ばれ、ある時間周波数成分
【数２１】
2008145610.tif 000023
の尤度への寄与を確率的に分配する。この結果、従来のクラスタリングでは一意に帰属音源を決められないデータであっても、この手法によればその曖昧さを確率として含んで取り扱う枠組みを実現可能である。
【００６３】
また、音源方向
【数２２】
2008145610.tif 000024
が決まれば、ｊ番目の音源を分離するバイナリマスクは、下記式（１２）
【数２３】
2008145610.tif 000025
のように設計することができる。
【００６４】
５．手段
本発明は、具体的には以下のような手段を採用する。
【００６５】
（１）上記課題を解決するために、本発明は、複数チャネルの観測信号に基づき、複数の音源からの信号を分離する音源分離定位方法において、前記複数チャネルの観測信号のそれぞれを、時間周波数領域の信号に変換するステップと、音源定位の初期値を決定する初期音源定位ステップと、前記音源定位に基づき、前記時間周波数領域上において、各音源の信号を分離する音源分離ステップと、前記分離結果に基づき、前記時間周波数領域上で各成分が属する音源の尤度を求め、この尤度の期待値である補助関数を最大値にする音源定位を求めて、音源の定位を確率的に行う繰り返し音源定位ステップと、前記音源分離ステップと、前記繰り返し音源定位ステップと、を繰り返し実行する繰り返しステップと、を含むことを特徴とする音源分離定位方法である。
【００６６】
（２）また、本発明は、上記（１）記載の音源分離定位方法において、前記音源定位は、前記音源の位置、前記音源の方向、前記音源からの音波の前記複数チャネル間の時間差、のいずれか１種又は２種以上を含むことを特徴とする音源分離定位方法である。
【００６７】
（３）また、本発明は、上記（１）記載の音源分離定位方法において、前記初期音源定位ステップは、各観測信号の間の特徴量に基づきクラスタリングを行うステップと、クラスタリングの結果に基づき各時間周波数成分ごとにいずれの音源に属するかを決定するステップと、上記決定に基づきそれぞれの音源の方向を求めるステップと、を含むことを特徴とする音源分離定位方法である。
【００６８】
（４）また、本発明は、上記（１）記載の音源分離定位方法において、前記音源分離ステップは、各時間周波数成分がどの音源に属するかを表す尤度に対する寄与を表す分配関数を求めることによって、音源分離を行うことを特徴とする音源分離定位方法である。
【００６９】
（５）また、本発明は、上記（４）記載の音源分離定位方法において、前記繰り返し音源定位ステップは、前記分配関数を、対応する各時間周波数成分の尤度に乗算し、その結果を、全時間・全周波数にわたって加算する補助関数を計算するステップと、前記分配関数のパラメータを所定の数値範囲で走査し、前記補助関数の値が最大となるようなパラメータを求めるステップと、を含み、前記新たなパラメータに基づき、前記音源分離ステップにおいて尤度が求められることを特徴とする音源分離定位方法である。
【００７０】
（６）また、本発明は、上記（５）記載の音源分離定位方法において、前記パラメータ前記複数の音源方向の組であることを特徴とする音源分離定位方法である。
【００７１】
（７）また、本発明は、上記（５）記載の音源分離定位方法において、前記補助関数は、時間周波数領域上の関数であることを特徴とする音源分離定位方法である。
【００７２】
（８）また、本発明は、上記（１）記載の音源分離定位方法において、前記繰り返しステップは、定位する音源の方向の変化が所定量未満になった場合に、繰り返しを終了することを特徴とする音源分離定位方法である。
【００７３】
（９）また、本発明は、上記（５）記載の音源分離定位方法において、前記繰り返しステップは、前記補助関数の値の増加が所定量未満になった場合に、繰り返しを終了することを特徴とする音源分離定位方法である。
【００７４】
（１０）上記課題を解決するために、本発明は、Ｎチャネルの観測信号に基づき、Ｍ個の音源からの信号を分離する音源分離定位方法において、前記Ｎチャネルの観測信号のそれぞれを、時間周波数領域の信号に変換するステップと、音源定位の初期値を決定する初期音源定位ステップと、前記音源定位に基づき、前記時間周波数領域上において、各音源の信号を分離する音源分離ステップと、前記分離結果に基づき、前記時間周波数領域上で各成分が属する音源の音源方向の尤度を求め、この尤度の期待値である補助関数を最大値にする前記音源間の音源方向の値を確率的に求めることによって音源定位を行う繰り返し音源定位ステップと、前記音源分離ステップと、前記繰り返し音源定位ステップと、を繰り返し実行する繰り返しステップと、を含むことを特徴とする音源分離定位方法。ここで、前記Ｎ及びＭは２以上の整数である。
【００７５】
ここで、ＭとＮの大小関係は問わない。本発明では、ＭはＮより大きくてもかまわない。
【００７６】
（１１）また、本発明は、上記（１０）記載の音源分離定位方法において、前記音源定位は、前記各音源の方位、又は、前記記Ｎチャネル間の観測信号の時間差、強度比、誤差分散のいずれか１種又は２種以上を含むことを特徴とする音源分離定位方法である。
【００７７】
（１２）また、本発明は、上記（１０）記載の音源分離定位方法において、前記初期音源定位ステップは、各観測信号の間の時間差又は強度比に基づきクラスタリングを行うステップと、クラスタリングの結果に基づき各時間周波数成分ごとにいずれの音源に属するかを決定するステップと、上記決定に基づきそれぞれの音源の方位として、前記時間差又は強度比を求めるステップと、を含むことを特徴とする音源分離定位方法である。
【００７８】
（１３）また、本発明は、上記（１０）記載の音源分離定位方法において、前記音源分離ステップは、各時間周波数成分が属する音源の前記時間差又は強度比を表す尤度に対する寄与を表す分配関数を求めることによって、音源分離を行うことを特徴とする音源分離定位方法である。
【００７９】
（１４）また、本発明は、上記（１３）記載の音源分離定位方法において、前記繰り返し音源定位ステップは、前記分配関数を、対応する各時間周波数成分の尤度に乗算し、その結果を、所定の数値範囲にわたって加算する補助関数を計算するステップと、前記分配関数のパラメータである時間差又は強度比を所定の数値範囲で走査し、前記補助関数の値が最大となるような時間差又は強度比を求めるステップと、を含み、前記新たな時間差又は強度比に基づき、前記音源分離ステップにおいて尤度が求められることを特徴とする音源分離定位方法である。
【００８０】
（１５）また、本発明は、上記（１４）記載の音源分離定位方法において、前記補助関数は、時間周波数領域上の関数であることを特徴とする音源分離定位方法である。
【００８１】
（１６）また、本発明は、上記（１０）記載の音源分離定位方法において、前記繰り返しステップは、定位する音源の方向の変化が所定量未満になった場合に、繰り返しを終了することを特徴とする音源分離定位方法である。
【００８２】
（１７）また、本発明は、上記（１４）記載の記載の音源分離定位方法において、前記繰り返しステップは、前記補助関数の値の増加が所定量未満になった場合に、繰り返しを終了することを特徴とする音源分離定位方法である。
【００８３】
（１８）また、本発明は、上記課題を解決するために、コンピュータに、複数チャネルの観測信号に基づき、複数の音源からの信号を分離する音源分離定位方法を実行させるためのプログラムにおいて、前記複数チャネルの観測信号のそれぞれを、時間周波数領域の信号に変換する手順と、前記時間周波数領域上で音源定位を行う初期音源定位手順と、前記音源定位に基づき、前記時間周波数領域上において、各音源の信号を分離する音源分離手順と、前記分離結果に基づき、前記時間周波数領域上で各成分が属する音源の尤度を求め、この尤度の期待値である補助関数を最大値にする音源定位を求めて、音源の定位を確率的に行う繰り返し音源定位手順と、前記音源分離ステップと、前記繰り返し音源定位ステップと、を繰り返し実行する繰り返し手順と、を、前記コンピュータに実行させることを特徴とするプログラムである。
【００８４】
（１９）また、本発明は、上記課題を解決するために、コンピュータに、Ｎチャネルの観測信号に基づき、Ｍ個の音源からの信号を分離する音源分離定位方法を実行させるプログラムにおいて、前記Ｎチャネルの観測信号のそれぞれを、時間周波数領域の信号に変換する手順と、前記時間周波数領域上で、観測信号間の相違量に基づき、音源方向を求めることによって音源定位を行う初期音源定位手順と、前記音源定位に基づき、前記時間周波数領域上において、各音源の信号を分離する音源分離手順と、前記分離結果に基づき、前記時間周波数領域上で各成分が属する音源の音源方向の尤度を求め、この尤度の期待値である補助関数を最大値にする前記音源間の音源方向の値を確率的に求めることによって音源定位を行う繰り返し音源定位手順と、前記音源分離ステップと、前記繰り返し音源定位ステップと、を繰り返し実行する繰り返し手順と、を、前記コンピュータに、実行させることを特徴とするプログラムである。ここで、前記Ｎは２以上の整数であり、前記ＭはＮより大きい整数である。
【００８５】
（２０）また、本発明は、上記（１）又は（１０）記載の音源分離定位方法において、前記繰り返し音源定位ステップは、さらに、前記補助関数をσ２で微分した更新式を計算することによって雑音の分散σ２を求めるステップ、を含むことを特徴とする音源分離定位方法である。
【００８６】
（２１）また、本発明は、上記（１）又は（１０）記載の音源分離定位方法において、前記繰り返しステップを終了した際に得られた分配関数に基づき、最も大きい寄与の音源のみを通過させる時間周波数マスキングによって、最終的な音源分離を行なう分離ステップと、を含むことを特徴とする音源分離定位方法である。
【００８７】
さて、上記（１）や（１０）の方法で、各時間周波数成分ごとに寄与率（分配関数）が最終的に求められた後、当然に、その求められた分配関数で分配して音源分離を行うことが考え方としては一般的であろう。
【００８８】
しかし、実際の分離においては、最終的な分離のためのマスクはバイナリマスクを用いた方が処理が簡単で現実的である。さらに、バイナリマスクでも聴感上は遜色ないことが経験的に知られている。
【００８９】
そこで、（２１）においては、各時間周波数成分ごとに、最も大きい寄与の音源のみを通過させるようなバイナリの時間周波数マスキングを行うことが好ましい。このようなマスクによって、各成分は、寄与が最大の音源に対してのみ「１：通過」、他の音源に対しては「０：不通過」となる。
【発明の効果】
【００９０】
以上述べたように、本発明によれば、音源定位を確率的に取り扱う枠組みを提供することによって、従来の音源定位・音源分離技術よりより精度の向上した音源分離定位方法が得られる。したがって、従来では音源分離できなかった残響の多い環境や、雑音の多い環境下でも音源分離、音源定位可能である。
【００９１】
また、本発明によれば、同じ時間周波数領域上の繰り返し処理によって、より精度の高い音源分離・音源定位を実現することができる。
【００９２】
また、この繰り返しは、ＥＭアルゴリズムと呼ばれており、本発明では、そのＥ−ｓｔｅｐが音源分離、そのＭ−ｓｔｅｐが音源定位にそれぞれ相当する。特に、ＥＭアルゴリズムの採用によって定位と分離で共通の目的関数（尤度）の最大化を行ええばよいという枠組みを提供することができた。
【００９３】
特に、本発明で採用するＥＭアルゴリズムでは、目的関数が尤度であるため、情報量規準などにより音源数推定が可能になる枠組みを実現できた。
【００９４】
さらに、本発明によれば、時間差・強度比などの特徴量空間ではなく、元の信号空間（時間周波数空間）での観測モデルにより尤度を定義しており、残響環境に対する拡散音場モデルなど、物理的な観測モデルの導入が可能である。
【００９５】
また、本発明によれば、繰り返しの際に雑音の分散も更新したので、観測信号中の雑音の分散も求めることができる。
【発明を実施するための最良の形態】
【００９６】
以下、図面に基づき、本発明を実施するための最良の形態を説明する。特に、具体的なアルゴリズムの導出を２ｃｈを例にして説明する。
【００９７】
６．本実施の形態の具体的アルゴリズムの導出
６．１  単一方向尤度の導出
上で述べた本発明の枠組み自体は、明らかに、一般論として、ＮｃｈのＢＳＳに適用可能である。ここで、Ｎとは自然数であり、マイクロフォンをＮ個用いて得たＮ個の音声信号（観測信号）に関するＢＳＳに適用できる。本６章では、チャネル数、すなわち観測信号の数としてこのＮを用いる。他の章では便宜上他の表記をする場合もある。
【００９８】
しかし、以下では具体的なアルゴリズムをわかりやすく説明するため、Ｎ＝２、すんわち２ｃｈのＢＳＳに話を絞って説明を進める。
【００９９】
具体的なアルゴリズムを決めるため必要なことは、単一方向尤度の与え方と、未知パラメータの選択である。
【０１００】
以下、単一方向尤度の算出プロセスを説明する。また、未知パラメータとしては、時間差δ、強度比ａ、誤差分散σが考えられるが、これらの内、どれを未知として扱い、どれを既知の値とするかは、種々の考え方があり、設計思想によっても異なってくる。本実施の形態では、後述するように強度比ａは全て１であるとし、時間差δを未知パラメータとして扱う。しかし、それ以外の選び方でもかまわない。
【０１０１】
一つの音源から発せられた信号が空間中を球面波として音源から伝播すると仮定すると、２ｃｈの観測信号
【数２４】
2008145610.tif 000026
の間には、音源位置に依存した時間差δｋ、強度比ａｋが生じる。以下では簡単のため、強度比に関してはａｋ＝１として平面波伝播を仮定し、音源の位置情報として方向θｋの代わりに、これとほぼ等価な情報である時間差δｋを用いる。
【０１０２】
本実施の形態では、時間差δｋを用いたが、強度比ａｋや、又は、時間差δｋと強度比ａｋの双方を用いることも好ましい。
【０１０３】
なお、マイクロフォンが検出する観測信号は、時間によって変化する振幅で表されるが、これを短時間フーリエ変換によって時間周波数領域に変換しておく。本特許では、観測信号は全て時間周波数領域に変換されたものとして取り扱っている。この様子が図１の（１）に示されている。特に図１では、横軸が時間、縦軸が周波数を表すグラフが２枚示され、２チャネル分の観測信号があることが示されている。図１には、本実施の形態の音源分離定位方法の流れを示す概念図が示されている。
【０１０４】
さて、２ｃｈ間の信号の間の時間差δｋを導入するためにステアリングベクトルを下記のように
【数２５】
2008145610.tif 000027
とベクトル表示すると、観測モデルは、下記式（１３）
【数２６】
2008145610.tif 000028
と表せる。ステアリングベクトルとは、音源方向に依存する複素ベクトルである。ただし、
【数２７】
2008145610.tif 000029
は、それぞれの観測信号に含まれる、残響、背景音を含む誤差であり、この誤差も観測信号等と同様の複素ベクトルであり太字で記されている。ここでＮＬ（τ、ω）、ＮＲ（τ、ω）は音源信号Ｓｋ（τ、ω）とは独立であると仮定する。
【０１０５】
なお、以下、音源信号や観測信号、誤差の記号から（τ、ω）を省略し、
【数２８】
2008145610.tif 000030
等と表記する。
【０１０６】
【数２９】
2008145610.tif 000031
のガウス雑音と仮定すると、その尤度は、
【数３０】
2008145610.tif 000032
と表される。上記式（１４）（数３０）には、まだ未知の変数Ｓｋが含まれているので、これをどのように扱うかにはいくつかの選択肢が考えられる。
【０１０７】
本実施の形態では、そのような選択肢の一つとして、Ｓｋとして、最尤値を採用する例を示す。すなわち、まず、下記式（１５）
【数３１】
2008145610.tif 000033
を上記式（１４）に代入し、さらに、下記式（１６）のように、
【数３２】
2008145610.tif 000034
と仮定すれば、単一方向尤度の具体的な形は、下記式（１７）
【数３３】
2008145610.tif 000035
となる。
【０１０８】
６．１．１  時間差δｋの初期値
さて、次節からｋ番目の音源に対応する時間差δｋを逐次更新していくことになるが、それには初期値が必要となる。本実施の形態では、初期値を従来のクラスタリングで取得している。これは既に説明したように、時間周波数領域における各成分に対して対応する音源を決定し、決定に基づき各音源の信号を取り出し、取り出した信号から、各音源の時間差δｋを求めている。図１の例では、例えば音源が３個の場合の例が図１（ｂ）に示されている。この場合、２本の観測信号間の時間差がδ１、δ２、δ３として求められる。これらはいわば初期値として扱われる。
【０１０９】
初期値の求め方は種々考えられる。もちろん、真の値に近い値の方が収束するための時間が短くなることは言うまでもない。
【０１１０】
６．２  Ｑ関数の導出
さて、ｋ番目の音源に対応する時間差δｋを求めるためのＱ関数（補助関数）は、上記式（１０）（数１９参照）より、下記式（１８）
【数３４】
2008145610.tif 000036
のように表される。ただし、δｋを含まない項は簡単のため、定数Ｃとして表した。この表現を見ると、分配関数の平方根
【数３５】
2008145610.tif 000037
は観測信号ＭＬ、ＭＲへ乗じられ、いわば、連続値のマスクに相当する働きを持つことがわかる。すなわち、Ｅステップで分配関数を計算しＱ関数を求める操作は、観測信号ＭＬ、ＭＲに連続値の分配関数を乗じてソフトに（なだらかに）音源分離を行えっているとみなすことができる。
【０１１１】
このＥステップの動作の概念が、図１（ｃ）に示されている。このように、バイナリマスクではなく、０〜１の連続的な透過率を持ったマスクによる音源の分離と考えられる。
【０１１２】
本実施の形態において特徴的なことは、このような連続的な値を有するマスクを用いて音源分離を行ったことである。この結果、音源の推定精度を向上させることが可能である。図１（ｃ）の例では音源が３個の例を示しており、Ｓ１、Ｓ２、Ｓ３に観測信号が分けられている。上述したようにこの分離は、連続的な値を有するマスクで行っている。この信号Ｓ１、Ｓ２、Ｓ３は、短時間複素フーリエ変換（ＳＴＦＴ）で表されている。図１（ｃ）に示されているように２チャネル分のグラフ（横軸は時間、縦軸は周波数）で２種の観測信号があることを示している。
【０１１３】
この分離は、全体の方向尤度（式（５）：（数８参照））を最大にする分解になっているが、音源信号の分離という観点から最適になっている保障はないため、本発明ではこの連続値マスクによる分離は音源定位のためにのみ行えい、分離は、式（１３）（数２６参照）のように最尤となる音源に基づくバイナリマスキングにより行うことにしている。
【０１１４】
Ｑ関数は、対数尤度の期待値である。本文の条件では、隠れ変数ｋが不明（ｋ：音源の数が不明）であるので、前のステップまでに推定されたパラメータの値を用いて、対数尤度の期待値を求めた結果がＱ関数である。
【０１１５】
６．３  パラメータの更新式
音源位置に対応した時間差δｋを更新するためには、下記の式（１９）を最大とするδｋを求める必要があるが、これは解析的には求まらないので、まず、適当に離散化した方向全てに対して
【数３６】
2008145610.tif 000038
を数値的に求め、以下の式（１９）のように更新することとした。
【０１１６】
【数３７】
2008145610.tif 000039
この結果、パラメータであるδｋが更新される。これはＥＭアルゴリズムのＭステップに相当する。この様子が図１（ｄ）に示されている。この操作は、抽出した各観測信号（時間周波数領域上の信号）から新しい時間差を求める作業となる。この新しい時間差δｋを用いて、再び音源分離（Ｅステップ）が実行される（図１（ｃ）参照））。
【０１１７】
一方、雑音の分散σ2 も、未知パラメータとして扱い、データから学習することができる。全Ｑ関数をσ2 で微分して０とおくことにより、下記の更新式（２０）を得る。
【０１１８】
【数３８】
2008145610.tif 000040
ただし、ここでのＮは全時間周波数成分の個数である。なお、上記の更新式では、全ての時間周波数成分に対し、雑音分散が等しいと仮定したが、これを周波数ごとに異なる値として求めることも可能である。
【０１１９】
各時間周波数成分を各音源に対し、どのように確率的に分配するかは、この分散の大きさに依存する。上記のように、分散を逐次的に推定することは本実施の形態では大きなポイントである。
【０１２０】
おおまかにいえば、非常に誤差分散が大きく雑音が大きい環境では、観測された音源方向が、雑音の影響である音源の方向とずれることも十分にありえる。そのため、誤差分散が大きい環境下では、どの音源にもある程度ずつ、エネルギーが分配されることになる。
【０１２１】
これに対し、誤差分散が小さく雑音が小さい環境では、観測された音源方向が正しい音源方向から大きくずれることはあまりないので、必然的に、観測に最も近い音源に、多くのエネルギーが分配されることになる。
【０１２２】
したがって、最終的な音源定位等からその「ばらつき」として分散を推定するのではなく、逐次的にこの分散を更新していくことが、残響環境下で分離と定位を良好に働かせるポイントの１つである。
【実施例】
【０１２３】
７．  シミュレーション実験による検証
２ｃｈＢＳＳに対しＥＭアルゴリズムを適用した提案手法を実装し、シミュレーションによりその分離性能を確認した。図２のように３つの音源及び２つのマイクロフォンを配置し、球面波伝播と残響を鏡像法のシミュレーションによって行った。分離性能の評価には、分離の前後での元音声に対するＳ／Ｎ比の改善値を用いた。音声データは研究用連続音声データベース（著作者：板橋秀一［日本音響学会／編］１９９１Ｖｏｌ。１−３）を使用した。（＝１０２４Ｐ）
また、Ｙｉｌｍａｚら（非特許文献１）の議論をもとに、サンプリング周期１６ｋＨｚ、フレーム長Ｔは２１０（＝１０２４）、シフトはＴ／２（５１２）、窓関数をＨａｍｍｉｎｇ窓として、観測信号を短時間Ｆｏｕｒｉｅｒ変換して時間周波数表現を得た。マイクロフォン間の距離は４ｃｍとした。ＥＭアルゴリズムの反復は、Ｑ関数の増加がある閾値以下になったら終了とした。
【０１２４】
比較対象とした従来法は、Ｙｉｌｍａｚらの手法（非特許文献１）に基づいた。まず、パワーで重みづけしたａ、δの２次元ヒストグラムを作成し、これに矩形関数を畳み込むことでスムージングを行えい、そのピークの頂点をａｉ、δｉとして推定した。次に各時間周波数成分に対する尤度が最大になる音源を通過させるマスクを作成し、分離信号を得た。なお、ここでｉは音源のインデックスである。
【０１２５】
音源定位結果を図３（１）、分離結果を図３（２）に示す。
【０１２６】
まず、残響がない場合には、本実施例の手法も、従来の手法もどちらもある程度の精度で定位を行えい分離が可能である。図３（１）に示すように、音源Ｓ１の真の時間差１０．４μｓに対して、従来手法では１０．３μｓ、本実施例の手法では９．８μｓと、双方ともかなり近い値を示している。音源Ｓ２に関しても、音源Ｓ２の真の時間差０．０μｓに対して、従来手法では０．０μｓ、本実施例の手法でも０．０μｓと、双方とも真の値と等しくなった。音源Ｓ３に関しても、音源Ｓ３の真の時間差−７．３μｓに対して、従来手法では−６．７μｓ、本実施例の手法でも−６．７μｓと、双方とも真の値と非常に近い値を示している。
【０１２７】
しかしながら、３７０μｓの残響環境下においては、従来手法の場合ではクラスタリングができず、したがって次のステップである音源分離が破綻している。これに対して、本実施例で提案する手法では、音源定位ができ、分離も行えていることが確認できよう。
【０１２８】
図３（１）に示すように、音源Ｓ１の真の時間差１０．４μｓに対して、残響影響下では従来手法では１．０μｓ、本実施例の手法では１０．３μｓとなり、従来手法では音源定位が破綻しているが、本実施例では真の値に近い値を示している。音源Ｓ２に関しても、音源Ｓ２の真の時間差０．０μｓに対して、従来手法では−４．２μｓ、本実施例の手法では０．０μｓとなり、従来手法では音源定位が破綻しているが、本実施例では真の値に近い値を示している。音源Ｓ３に関しても、音源Ｓ３の真の時間差−７．３μｓに対して、従来手法では−５．１μｓ、本実施例の手法では−８．８μｓとなり、本実施例の手法の方がより真の値と近い値を示している。
【０１２９】
また、本実施の形態において提案する手法における、σ2 の推定値と残響時間との関係が図３（３）に示されている。残響時間が長くなるにつれて、σ2 の推定値が大きくなっており、環境に応じて観測誤差の大きさを推定することが可能であることが理解されよう。ただし、σは観測信号の振幅と同じ単位を持つが特に明記していない。
【０１３０】
このように本実施の形態によれば、観測モデルに含まれる雑音項の大きさ（分散）も観測信号から推定可能であり、パラメータチューニング等が不要であるという効果を奏する。
【０１３１】
８．  まとめ
以上、本実施の形態では、ＮｃｈのＢＳＳ（Ｎは２以上の整数）に対し、各時間周波数成分の帰属音源を隠れ変数としてＥＭアルゴリズムを適用して音源分離を行う手法を提案した。提案した手法のＥＭアルゴリズムにおけるＥステップが音源分離を、Ｍステップが音源定位を行っていると見ることができる。
【０１３２】
従来手法では、音源定位・音源分離を２つの別個の処理によって分離を行っていたのに対し、本実施の形態で提案する手法は、音源定位と音源分離について共通の目的関数を導入し、この共通の目的関数を扱う統一した処理を実行する枠組みを採用している。
【０１３３】
また、上記実施例では、２ｃｈのＢＢＳについて、シミュレーション実験によって従来手法と、本実施の形態で提案する手法とを比較する実験を行えった。この実験結果から、従来手法では分離が困難になるような残響環境下においても、本特許で提案する手法によれば、分離を行うことができることが示された。
【０１３４】
９．  プログラム
本実施の形態で述べた音源分離・定位方法は、種々の装置、各種のＬＳＩ、種々のハードウェア・ソフトウェアで実行することが可能である。上で述べたシミュレーションもコンピュータ上で全てディジタル信号として音声を取り扱い、シミュレーションを実行している。
【０１３５】
好ましい一例としては、コンピュータ上で、上記の各種式を計算する処理はそのようなプログラムをコンピュータに実行させることによって実現することが挙げられよう。数式の計算は、コンピュータの一般的な動作であるため、そのようなプログラムを記述することは当業者にとって容易である。また、パラメータを所定の範囲で動かして最大値を求める等の繰り返し処理もコンピュータで一般に行える処理であるため、そのような繰り返しを伴うプログラムを記述することは当業者にとって容易である。
【０１３６】
また、そのようなプログラムは、コンピュータのハードディスク等の記録媒体に格納しておくことが好ましい。記録媒体は種々の光ディスクや磁気ディスク等の持ち運び可能な記録媒体に格納しておくことも好ましい。
【０１３７】
１０．  本件発明の特徴
本件発明では、スパース性に基づくＢＳＳの問題に対し、ＥＭアルゴリズムを適用した新しいアルゴリズムを提案した。具体的には、音源信号のスパース性に基づき、各時間周波数成分に寄与する音源は１個であるという観測モデルの下、各時間周波数成分に寄与している音源のインデックスを隠れ変数とみなし、最大尤度を与える音源方向と、各時間周波数成分への各音源の寄与率をＥＭアルゴリズムによって推定する。
【０１３８】
本件発明の方法は、
（１）Ｅ−ｓｔｅｐが音源分離、Ｍ−ｓｔｅｐが音源定位に相当し、定位と分離で共通の目的関数（尤度）の最大化が行えわれること、
（２）目的関数が尤度であるため、情報量規準などにより音源数推定が可能になる枠組みであること、
（３）時間差・強度比などの特徴量空間ではなく、元の信号空間（時間周波数空間）での観測モデルにより尤度を定義しており、残響環境に対する拡散音場モデルなど、物理的な観測モデルの導入が可能になること、
等の特長を有している。
【０１３９】
全体の枠組みは、一般にＮｃｈのＢＳＳに適用できる。上記実施例では特に２ｃｈのＢＳＳに議論を絞って具体的なアルゴリズムを導出し、シミュレーションによる本件手法の音源分離実験結果を示した。
【図面の簡単な説明】
【０１４０】
【図１】本実施の形態の音源分離定位方法の流れを示す概念図が示されている。
【図２】本実施例におけるシミュレーション実験の各構成の配置の様子を示す説明図である。
【図３】本実施例におけるシミュレーション実験結果を示す表の図である。
【図４】残響時間０ｍｓｅｃである場合の時間差の散布図である。
【図５】残響時間５０ｍｓｅｃである場合の時間差の散布図である。
【図６】残響時間１７０ｍｓｅｃである場合の時間差の散布図である。
【符号の説明】
【０１４１】
δ  時間差
Ｓ  観測信号（短時間複素フーリエ変換済）
ＭＬ  左側観測信号
ＭＲ  右側観測信号
【図１】
2008145610.tif 000041
【図２】
2008145610.tif 000042
【図３】
2008145610.tif 000043
【図４】
2008145610.tif 000044
【図５】
2008145610.tif 000045
【図６】
2008145610.tif 000046