特許 6288561 ブラインド信号分離方法およびその装置

(19)【発行国】日本国特許庁(JP)

(12)【公報種別】特許公報(B2)

(11)【特許番号】6288561

(24)【登録日】2018年2月16日

(45)【発行日】2018年3月7日

(54)【発明の名称】ブラインド信号分離方法およびその装置

(51)【国際特許分類】

   G10L 21/0272 20130101AFI20180226BHJP

【ＦＩ】

   G10L21/0272 100Z

【請求項の数】4

【全頁数】19

(21)【出願番号】特願2014-101793(P2014-101793)

(22)【出願日】2014年4月24日

(65)【公開番号】特開2015-210512(P2015-210512A)

(43)【公開日】2015年11月24日

【審査請求日】2017年3月16日

(73)【特許権者】

【識別番号】509257961

【氏名又は名称】齋藤 晋哉

(73)【特許権者】

【識別番号】506178874

【氏名又は名称】大石 邦夫

(73)【特許権者】

【識別番号】513234798

【氏名又は名称】古川 利博

(72)【発明者】

【氏名】齋藤 晋哉

(72)【発明者】

【氏名】大石 邦夫

(72)【発明者】

【氏名】古川 利博

【審査官】 菊池 智紀

(56)【参考文献】

【文献】 特開２０１１−８１３７３（ＪＰ，Ａ）

【文献】 特開２００８−３３３０７（ＪＰ，Ａ）

【文献】 齋藤晋哉 他，"逆行列の補題を用いた適応ALSアルゴリズムにおけるブラインド音源分離"，電子情報通信学会技術研究報告，２００７年 ２月２７日，Vol.106，No.570，pp.63-67

【文献】 SAITO, Shinya et al.，"AN APPROACH TO CONVOLUTIVE BACKWARD-MODEL BLIND SOURCE SEPARATION BASED ON JOINT DIAGONALIZATION"，Proc. of EUSIPCO2012，２０１２年 ８月２７日，pp.579-583

(58)【調査した分野】（Int.Cl.，ＤＢ名）

Ｇ１０Ｌ ２１／０２０８−２１／０３６４，

１５／２０

ＩＥＥＥ Ｘｐｌｏｒｅ

(57)【特許請求の範囲】

【請求項1】

未知の畳み込み混合系により混在した互いに統計的に独立な未知信号源信号を、観測信号のみからブラインドで推定する方法であって、時間周波数領域においてエポック時刻毎に観測信号ベクトルから空間相関行列を求め、各エポック時刻において周波数ビンの中からフロベニウスノルムが最大となる空間相関行列を求め、そのノルムで全ての同エポック時刻の空間相関行列を正規化した後、正規化された空間相関行列から混合行列を推定するために、制約条件付きフォワードモデル型最小２乗型同時対角化問題をラグランジュの未定乗数法と反復法の組み合わせにより解法する手順と、手順の収束後、累乗法を１回用いて最小２乗型同時対角化問題の解となる混合行列を求め、混合行列のランクが落ちていた場合は、混合行列がフルランクになるように特異値分解を用いて基底を補い、次いで、ランク落ちした混合行列の特異値とその総和が変わらないように特異値を補い、最小２乗型一般化逆行列を用いて混合行列から分離行列を推定し、バックワード型最小２乗型同時対角化問題を最小２乗法で解いて対角行列を求め、同一信号源から発生した信号の周波数ビンの電力比に相関があることを利用してパーミュテーション問題を、多数決を利用して解法することによって信号分離精度の高い分離行列を求めることを特徴とするブラインド信号分離方法。

【請求項2】

請求項１記載の最小２乗型同時対角化問題の解法を用いたブラインド信号分離方法において、信号源とブラインド信号分離装置出力間の伝達関数を表した混合行列と分離行列の縦続接続モデルにラグランジュの未定乗数を導入した制約条件付きフォワードモデル型最小２乗型同時対角化問題とバックワード型最小２乗型同時対角化問題を導入し、更に分離行列に遅延を与えることによって信号源とブラインド信号分離装置出力間の因果的なモデルを推定することを特徴とするブラインド信号分離方法。

【請求項3】

請求項１乃至請求項２のいずれか１項に記載の最小２乗型同時対角化問題の解法を用いたブラインド信号分離方法において、基準周波数ビンを複数選択し、基準周波数ビン間において最も電力比の相関が大きいパーミュテーション行列を推定し、複数の基準周波数ビンから基準周波数ビンを１つ選択し、全ての周波数ビン間で電力比の相関が最も大きいパーミュテーション行列を推定する手順を、全ての基準周波数ビンが１度選択されるまで繰り返した後、各周波数ビンに複数割り当てられたパーミュテーション行列から多数決によって周波数ビンに割り当てられるパーミュテーション行列を決定し、多数決によってパーミュテーション行列が決定できない場合には、最も相関値が大きいパーミュテーション行列を採用することを特徴とするブラインド信号分離方法。

【請求項4】

請求項１乃至請求項３のいずれか１項に記載のブラインド信号分離方法を用いて信号源分離を行うように構成されていることを特徴とするブラインド信号分離方法を用いたブラインド信号分離装置。

【発明の詳細な説明】

【技術分野】

【0001】

本発明は、未知の畳み込み混合系により混在した互いに統計的に独立な未知信号源信号を、観測信号のみから推定するブラインド信号分離方法に係わり、特に、最小２乗型同時対角化問題の解法を用いて高い精度で信号を分離することができるブラインド信号分離装置に関する。

【背景技術】

【0002】

複数の未知信号源信号が未知の畳み込み混合系により混在されて観測されるとき、観測信号を分離して混在前の未知信号源信号を推定する処理をブラインド信号分離という。ブラインド信号分離では、未知信号源信号間の統計的独立性のみを条件として、観測信号から未知信号源信号を推定する方法であり、信号源の位置或いは観測信号の到来方向の推定を必ずしも必要としない方法である。

【0003】

勾配法に基づく周波数領域の信号分離方法が非特許文献１で提案されている。勾配法に基づき周波数ビン毎に分離行列を更新することによって分離信号を推定する。

【0004】

最小２乗型同時対角化問題の解法を用いたブラインド信号分離方法が非特許文献２と非特許文献３で提案されている。これらの方法は、最小２乗型同時対角化問題の解法による混合行列の推定、最小２乗型一般化逆行列を用いた混合行列からの分離行列の推定、スケーリング問題の解法、パーミュテーション問題の解法の４つの手順から成る。

【0005】

周波数ビン間の電力比の相関に基づきパーミュテーション行列を推定する方法が非特許文献４で提案されている。

【先行技術文献】

【0006】

【非特許文献1】 Ｌ．Ｐａｒｒａ ａｎｄ Ｃ．Ｓｐｅｎｃｅ，”Ｃｏｎｖｏｌｕｔｉｖｅ ｂｌｉｎｄ ｓｅｐａｒａｔｉｏｎ ｏｆ ｎｏｎ−ｓｔａｔｉｏｎａｒｙ Ｓｏｕｒｃｅｓ”，ＩＥＥＥ Ｔｒａｎｓａｃｔｉｏｎｓ ｏｎ Ｓｐｅｅｃｈ ａｎｄ Ａｕｄｉｏ Ｐｒｏｃｅｓｓｉｎｇ，ｖｏｌ．８，ｎｏ．３，ｐｐ．３２０−３２７，Ｍａｙ，２０００

【非特許文献2】 Ａ．Ｙｅｒｅｄｏｒ，”Ｎｏｎ−ｏｒｔｈｏｇｏｎａｌ ｊｏｉｎｔ ｄｉａｇｏｎａｌｉｚａｔｉｏｎ ｉｎ ｔｈｅ ｌｅａｓｔ−ｓｑｕａｒｅｓ ｓｅｎｓｅ ｗｉｔｈ ａｐｐｌｉｃａｔｉｏｎ ｉｎ ｂｌｉｎｄ ｓｏｕｒｃｅ ｓｅｐａｒａｔｉｏｎ”，ＩＥＥＥ Ｔｒａｎｓａｃｔｉｏｎｓ ｏｎ Ｓｉｇｎａｌ Ｐｒｏｃｅｓｓｉｎｇ，ｖｏｌ．５０，ｎｏ．７，ｐｐ．１５４５−１５５３，Ｊｕｌｙ ２００２

【非特許文献3】 Ｋ．Ｒａｈｂａｒ ａｎｄ Ｊ．Ｐ．Ｒｅｉｌｌｙ，”Ａ ｆｒｅｑｕｅｎｃｙ ｄｏｍａｉｎ ｍｅｔｈｏｄ ｆｏｒ ｂｌｉｎｄ ｓｏｕｒｃｅ ｓｅｐａｒａｔｉｏｎ ｏｆ ｃｏｎｖｏｌｕｔｉｖｅ ａｕｄｉｏ ｍｉｘｔｕｒｅｓ”，ＩＥＥＥ Ｔｒａｎｓａｃｔｉｏｎｓ ｏｎ Ｓｐｅｅｃｈ ａｎｄ Ａｕｄｉｏ Ｐｒｏｃｅｓｓｉｎｇ，ｖｏｌ．１３，ｎｏ．５，ｐｐ．８３２−８４４，Ｓｅｐｔ．２００５

【非特許文献4】 Ｈ．Ｓａｗａｄａ，Ｓ．Ａｒａｋｉ，ａｎｄ Ｓ．Ｍａｋｉｎｏ，”Ｍｅａｓｕｒｉｎｇ ｄｅｐｅｎｄｅｎｃｅ ｏｆ ｂｉｎ−ｗｉｓｅ ｓｅｐａｒａｔｅｄ ｓｉｇｎａｌｓ ｆｏｒ ｐｅｒｍｕｔａｔｉｏｎ ａｌｉｇｎｍｅｎｔ ｉｎ ｆｒｅｑｕｅｎｃｙ−ｄｏｍａｉｎ ＢＳＳ”，Ｐｒｏｃ．ＩＥＥＥ Ｉｎｔ．Ｓｙｍｐ．Ｃｉｒｃｕｉｔｓ Ｓｙｓｔ．，ｐｐ．３２４７−３２５０，Ｍａｙ ２００７

【発明の概要】

【発明が解決しようとする課題】

【0007】

従来のブラインド信号分離方法で使用されている勾配法は演算量が少ないが、収束が遅く、実環境下で十分な信号分離精度を得ることができない。

【0008】

従来のブラインド信号分離方法では、時間周波数領域で周波数ビンとエポック毎、即ち、ポイント毎に観測信号の空間相関行列をそのフロベニウスノルムで正規化しているため、未知信号源信号が音声の場合、音声の特徴を表す周波数特性が失われ、周波数ビンに対して振幅特性が一定になり、信号分離精度が劣化する。

【0009】

従来のブラインド信号分離方法では、時間周波数領域で観測信号の空間相関行列から最小２乗型同時対角化問題の解法を用いて混合行列を推定し、次いで最小２乗型一般化逆行列を用いて混合行列から分離行列を推定しているため、最小２乗型同時対角化問題に分離行列の推定が含まれず、仮に混合行列を精度良く推定できても分離行列が精度良く推定できるとは限らず、信号分離精度が劣化する場合がある。

【0010】

従来のブラインド信号分離方法では、時間周波数領域で観測信号の空間相関行列から制約条件を課すことなく最小２乗型同時対角化問題の解法を用いて混合行列を推定し、最小２乗型一般化逆行列を用いて混合行列から分離行列を推定する前に、塁乗法により制約条件を課しているので、手順が収束に要する反復回数は大幅に増加する。

【0011】

更に、従来のブラインド信号分離方法において混合行列の最小２乗型一般化逆行列を単に分離行列として取り扱うと、加法性雑音が存在する場合、雑音の影響が大きくなって信号分離性能の劣化につながる。

【0012】

同一信号源から発生した信号の隣接周波数ビンに相関があることを利用した従来のパーミュテーション問題の解法では、１度誤りが生じると、これ以降の解法が誤り続ける確率が非常に高くなる。

【0013】

本発明はこのような事情を鑑みてなされたものであり、信号源とブラインド信号分離装置出力間を表した混合行列と分離行列の縦続接続モデルを最小２乗型同時対角化問題に取り込み、更にモデルに因果性を与えることによって信号分離性能の向上と反復手順の収束高速化の両立、更に、同一信号源から発生した信号の周波数ビン間の電力比相関性と多数決を組み合わせによってパーミュテーション問題を正確に解法することを目的とする。

【課題を解決するための手段】

【0014】

このような目的に応えるために本発明（請求項１記載の発明）に係るブラインド信号分離方法は、未知の畳み込み混合系により混在した互いに統計的に独立な未知信号源信号を、観測信号のみからブラインドで推定する方法であって、時間周波数領域においてエポック時刻毎に観測信号ベクトルから空間相関行列を求め、各エポック時刻において周波数ビンの中からフロベニウスノルムが最大となる空間相関行列を求め、そのノルムで全ての同エポック時刻の空間相関行列を正規化した後、正規化された空間相関行列から混合行列を推定するために、制約条件付きフォワードモデル型最小２乗型同時対角化問題をラグランジュの未定乗数法と反復法の組み合わせにより解法する手順と、手順の収束後、累乗法を１回用いて最小２乗型同時対角化問題の解となる混合行列を求め、混合行列のランクが落ちていた場合は、混合行列がフルランクになるように特異値分解を用いて基底を補い、次いで、ランク落ちした混合行列の特異値とその総和が変わらないように特異値を補い、最小２乗型一般化逆行列を用いて混合行列から分離行列を推定し、バックワード型最小２乗型同時対角化問題を最小２乗法で解いて対角行列を求め、同一信号源から発生した信号の周波数ビンに相関があることを利用してパーミュテーション問題を、多数決を利用して解法することによって信号分離精度の高い分離行列を求めることを特徴とする。

【0015】

本発明（請求項２記載の発明）に係るブラインド信号分離方法は、信号源とブラインド信号分離装置出力間の伝達関数を表した混合行列と分離行列の縦続接続モデルにラグランジュの未定乗数を導入した制約条件付きフォワードモデル型最小２乗型同時対角化問題とバックワード型最小２乗型同時対角化問題を導入し、更に分離行列に遅延を与えることによって信号源とブラインド信号分離装置出力間の因果的なモデルを推定することを特徴とする。

【0016】

本発明（請求項３記載の発明）に係るブラインド信号分離方法は、基準周波数ビンを複数選択し、基準周波数ビン間において最も電力比の相関が大きいパーミュテーション行列を推定し、複数の基準周波数ビンから基準周波数ビンを１つ選択し、全ての周波数ビン間で電力比の相関が最も大きいパーミュテーション行列を推定する手順を、全ての基準周波数ビンが１度選択されるまで繰り返した後、各周波数ビンに複数割り当てられたパーミュテーション行列から多数決によって周波数ビンに割り当てられるパーミュテーション行列を決定し、多数決によってパーミュテーション行列が決定できない場合には、最も相関値が大きいパーミュテーション行列を採用することを特徴とする。

【0017】

すなわち、本発明によれば、観測信号の空間相関行列を制約条件付き最小２乗型同時対角化問題の対象に、ラグランジュの未定乗数法と反復法を組み合わせることによって最小２乗型同時対角化問題の近似解、即ち、混合行列を求めることができ、推定性能の向上と反復回数の低減を両立させる。

【0018】

また、本発明によれば、各周波数ビンに複数割り当てられたパーミュテーション行列から多数決によって当該周波数ビンに割り当てられるパーミュテーション行列を決定し、全体の信号分離性能を向上させる。

【発明の効果】

【0019】

本発明によれば、時間周波数領域においてエポック時刻毎に観測信号ベクトルから空間相関行列を求め、各エポック時刻において周波数ビンの中からフロベニウスノルムが最大となる空間相関行列を求め、そのノルムで全ての同エポック時刻の空間相関行列を正規化した後、これに制約条件付き最小２乗型同時対角化問題を適用する。制約条件付き最小２乗型同時対角化問題では、ラグランジュの未定乗数法と反復法を組み合わせることによって近似解、即ち、混合行列を求める。本発明に係るブラインド信号分離方法及びブラインド信号分離装置では、混合行列のランクが落ちていた場合は、混合行列がフルランクになるように特異値分解を用いて基底と特異値を補い、最小２乗型一般化逆行列を用いて混合行列から分離行列を推定し、バックワード型最小２乗型同時対角化問題を最小２乗法で解いて対角行列を求めるため、少ない反復回数で近似解に収束するという効果がある。

【0020】

また、本発明に係るブラインド信号分離方法及びブラインド信号分離装置では、信号源とブラインド信号分離装置出力間の伝達関数を表した混合行列と分離行列の縦続接続モデルにラグランジュの未定乗数を導入した制約条件付きフォワードモデル型最小２乗型同時対角化問題とバックワード型最小２乗型同時対角化問題を導入し、更に分離行列に遅延を与えることによって信号源とブラインド信号分離装置出力間の因果的なモデルを推定することによって雑音の影響を最小にする分離行列を推定できるという効果がある。

【0021】

更に、パーミュテーション問題の解法では、基準周波数ビンを複数選択し、基準周波数ビン間において最も電力比の相関が大きいパーミュテーション行列を推定し、複数の基準周波数ビンから基準周波数ビンを１つ選択し、全ての周波数ビン間で電力比の相関が最も大きいパーミュテーション行列を推定する手順を全ての基準周波数ビンが１度選択されるまで繰り返した後、各周波数ビンに複数割り当てられたパーミュテーション行列から多数決によって周波数ビンに割り当てられるパーミュテーション行列を決定し、多数決によってパーミュテーション行列が決定できない場合には、最も相関値が大きいパーミュテーション行列を採用することを特徴とすることによって本発明に係るブラインド信号分離方法及びブラインド信号分離装置には、信号分離性能を高めることができるという効果がある。

【図面の簡単な説明】

【0022】

【図1】 本発明に係るブラインド信号分離方法の実施の形態を示す図である。

【図2】 本発明に係るブラインド信号分離方法において最小２乗型同時対角化問題の対角化行列のラグランジェ未定乗数法の適用について説明するためのフローチャートである。

【図3】 本発明に係るブラインド信号分離方法において最小２乗型同時対角化問題の対角行列の解法について説明するためのフローチャートである。

【図4】 本発明に係るブラインド信号分離方法のパーミュテーション問題の解法について説明するためのフローチャートである。

【図5】 本発明に係る信号分離方法の実施例１における信号源（スピーカ）とマイクロホンの位置関係を表す平面図である。

【図6】 本発明に係る信号分離方法の実施例１における部屋の残響時間と信号分離性能の関係を示す図である。

【図7】 本発明に係る信号分離方法の実施例１におけるＳＮＲと信号分離性能の関係を示す図である。

【図8】 実施例２における従来のパーミュテーション問題の解法によるパーミュテーション行列の割り当て結果を示す図である。

【図9】 実施例２における本発明に係るブラインド信号分離方法によるパーミュテーション行列の割り当て結果を示す図である。

【図10】 分離行列の推定に本発明に係るブラインド信号分離方法を使用して、パーミュテーション問題の解法には非特許文献４を使用した場合と、本発明に係るブラインド信号分離方法の実施例２におけるパーミュテーション行列の割り当て正答率と信号分離性能を比較する図である。

【発明を実施するための形態】

【0023】

本発明に係るブラインド信号分離方法の実施の形態について図面を参照して説明する。

【0024】

１．畳み込み混合モデル
図１に示すように、時刻ｔにおいてＮ個の信号源１１、１２、…、１Ｎから発せられた信号源信号ｓ_ｊ（ｔ）が畳み込み混合されてｘ_ｉ（ｔ）として観測される。ｓ_ｊ（ｔ）は平均０で互いに統計的独立な非定常信号である。ｈ_ｉｊ（ｔ）は信号源１ｊからマイクロホン２ｉまでの経路の時不変なインパルス応答で、因巣的で非最小位相系である。また、ｎ_ｉ（ｔ）はマイクロホン２ｉに加わる平均０、分散σ^２のガウス性白色雑音で、ｓ_ｊ（ｔ）と統計的独立である。時刻ｔにおいてＪ個のマイクロホン２１、２２、…、２Ｊで観測される観測信号ｘ_ｉ（ｔ）は式（１）で表される。ここで、Ｊ≧Ｎ≧２とする。

【0025】

【数1】

ここで、＊は畳み込み演算を表す。

【0026】

観測信号ｘ_ｉ（ｔ）を３１で短時間フーリエ変換すると、フレーム時刻ｍにおける観測信号は式（２）により表される。式（２）において、ｗｉｎ（ｔ）は窓関数、Ｋは短時間フーリエ変換の点数、Ｔ_ｂは２つの重複窓間のシフトサイズ、ω_ｋ＝２πｋ／Ｋ、ｋ＝０，１，…，Ｋ−１、をそれぞれ表す。離散フーリエ変換の点数Ｋがｈ_ｉｊ（ｔ）のインパルス応答長より十分に大きいとき観測信号は式（３）により近似される。ここで、ｈ_ｉｊ（ｔ）のＫ点離散フーリエ変換をｈ_ｉｊ（ω_ｋ）、ｓ_ｊ（ｔ）に窓関数を乗算した後、Ｋ点短時間フーリエ変換で時間周波数領域に変換したフレーム時刻ｍの信号源信号をｓ_ｊ（ω_ｋ，ｍ）、同様に、ｎ_ｉ（ｔ）に窓関数を乗算した後、Ｋ点短時間フーリエ変換で時間周波数領域に変換したフレーム時刻ｍの雑音をｎ_ｉ（ω_ｋ，ｍ）とそれぞれ表記している。また、式（３）において、ｘ（ω_ｋ，ｍ）はフレーム時刻ｍに各マイクロホンでの観測信号ベクトル、ｓ（ω_ｋ，ｍ）はフレーム時刻ｍに各信号源信号ベクトル、混合行列Ｈ（ω_ｋ）はＮ個の信号源からＪ個のマイクロホンまでの混合行列、ｎ（ω_ｋ，ｍ）は雑音ベクトルでそれぞれ式（４）、（７）、（５）、（８）により定義される。入手可能なエポック時刻の総数をＭとすると、１≦ｍ≦Ｍとなる。信号源信号の共分散行列はＰ_ｓ（ω_ｋ，ｍ）＝Ｅ［ｓ（ω_ｋ，ｍ）ｓ（ω_ｋ，ｍ）^Ｈ］∈Ｒ^Ｎ×Ｎで、対角行列となる。Ｅ［・］と上付き添字^Ｈは期待値と複素共役転置をそれぞれ表す。また、上付き添字^ＴとＲ^Ｎ×Ｎは転置とＮ×Ｎの実数空間を表す。

【0027】

【数2】

【0028】

信号を分離するには５１、５２、…、５Ｋで周波数ビン毎に式（９）を満足する分離行列Ｗ（ω_ｋ）を推定し、６０で信号源の割り当てを定めるパーミュテーション行列Π（ω_ｋ）∈Ｒ^Ｎ×Ｎを決定する。周波数ビン毎に独立にΠ（ω_ｋ）を決定しても信号が完全に分離する保証はなく、同一信号源から発生した信号の隣接または近接周波数ビンに相関があることを利用してパーミュテーション行列Π（ω_ｋ）を決定する。

【0029】

【数3】

ここで、Ｄ（ω_ｋ）∈Ｃ^Ｎ×Ｎは周波数ビン毎に異なる任意の対角行列である。

【0030】

スケーリング問題とパーミュテーション問題を順に解法した後、７１、７２、…、７Ｋでｘ（ω_ｋ，ｍ）に左から分離行列Ｗ（ω_ｋ）を乗算すると、周波数ビンω_ｋにおける分離信号ｙ（ω_ｋ，ｍ）は式（１０）で表される。尚、スケーリング問題の解法については後述する。式（１０）を８０で短時間逆フーリエ変換と重複加算によって時間領域に変換すると分離信号ｙ_ｉ（ｔ）が求められる。雑音の分散σ^２が十分に小さいとき、ｙ_ｉ（ｔ）≒ｓ_ｉ（ｔ）になる。尚、分離信号ベクトルｙ（ω_ｋ，ｍ）は式（１１）により表される。

【0031】

【数4】

【0032】

本発明に関するブラインド信号分離方法について、図１乃至図４を参照して詳細に説明する。図２乃至図４は、図１の３１における短時間フーリエ変換後、ブラインド信号分離システム４０において本発明により周波数ビン毎に推定される分離行列の算出手順を示したものである。図４は、分離行列の算出後、ブラインド信号分離システム４０において本発明によりパーミュテーション行列の算出手順を示したものである。

【0033】

観測信号ｘ（ω_ｋ，ｍ）の共分散行列Ｐ_ｘ（ω_ｋ，ｍ）∈Ｃ^Ｊ×Ｊは式（１２）で与えられる。式（１３）の制約条件を課して式（１４）を満足する対角化行列Ｂ（ω_ｋ）と対角行列Λ（ω_ｋ，ｍ）を求めると、式（９）よりＢ（ω_ｋ）とＷ（ω_ｋ）の関係は式（１５）で与えられる。ただし、Ｃ^Ｊ×ＪはＪ×Ｊの複素空間を表す。

【0034】

【数5】

ただし、Ｉは単位行列である。

【0035】

２．最小２乗型同時対角化問題とその解法
観測信号ｘ（ω_ｋ，ｍ）の共分散行列Ｐ_ｘ（ω_ｋ，ｍ）の推定値Ｐ_ｘ（ω_ｋ，ｍ）を正規化して、式（１７）を最小にする対角化行列Ｂ（ω_ｋ）と対角行列Λ（ω_ｋ，ｍ）を求める。式（１７）は最小２乗型同時対角化問題の解法として知られている。

【数6】

【0036】

本発明では、制約条件付き最小２乗型同時対角化問題を解法することによって混合行列Ｂ（ω_ｋ）を推定した後、分離行列Ｗ（ω_ｋ）を求める。次いで、対角行列Λ（ω_ｋ，ｍ）を推定するために、最小２乗法を用いて分離行列Ｗ（ω_ｋ）を用いた評価量を最小化する。本発明では、混合行列と対角行列の推定を交互に繰り返す。音声は低域周波数帯にフォルマントと呼ばれる振幅スペクトルのピークを有している。この音声波形の特徴を失うことなく、式（１９）によって観測信号ｘ（ω_ｋ，ｍ）の共分散行列Ｐ_ｘ（ω_ｋ，ｍ）を正規化することが、本発明の特徴の一つである。

【0037】

２．１ 時間周波数領域における観測信号の共分散行列の正規化
時間領域観測信号は式（１８）の短時間フーリエ変換によって時間周波数領域に変換される。

【0038】

【数7】

式（１７）において、ｗ（ｔ）は窓関数、Ｋは短時間フーリエ変換の点数、Ｔ_ｓは２つの重複窓間のシフトサイズ、Ｔ_ｓはエポックサイズ、Ｎ_ｓ＋１は各エポックにおける総重複フレーム数でＫ＋Ｎ_ｓＴ_ｓ≦Ｔ_ｂ、ｌ＝０，１，…，Ｎ_ｓをそれぞれ表す。

【0039】

ステップＳ１０１においてエポック時刻ｍにおける時間周波数領域観測信号の共分散行列を式（１９）によって推定される。

【0040】

【数8】

【0041】

混合行列Ｂ（ω_ｋ）と分離行列Ｗ（ω_ｋ）を縦続に接続したとき、そのインパルス応答は式（１５）を最小にすることによって求められる。式（２１）の評価量ｅ（ω_ｋ）をＷ（ω_ｋ）によって微分すると、分離行列Ｗ（ω_ｋ）は式（２２）によって求められる。Ｂ（ω_ｋ）^ＨＢ（ω_ｋ）のランクがＮのときのみ、式（２１）の評価量ｅ（ω_ｋ）は零になる。一方、Ｂ（ω_ｋ）^ＨＢ（ω_ｋ）のランクがＮ未満のとき、ｅ（ω_ｋ）は零より大きくなる。そこで、制約条件‖ｂ_ｊ（ω_ｋ）‖_２＝１に制約条件ｒａｎｋ（Ｂ（ω_ｋ）^ＨＢ（ω_ｋ））＝Ｎを付け加え、混合行列Ｂ（ω_ｋ）の推定のための最小２乗型同時対角化問題を解法する。ここで、ｂ_ｊ（ω_ｋ）はＢ（ω_ｋ）のｊ番目の列ベクトル、‖・‖_２はユークリッドノルム、ｒａｎｋ（Ａ）は行列Ａのランクをそれぞれ表す。

【0042】

【数9】

ここで、‖・‖_Ｆはフロベニウスノルムを表す。

【0043】

ステップＳ１０２においてＰ_ｘ（ω_ｋ，ｍ）を式（２３）によって正規化する。

【0044】

【数10】

【0045】

２．２ 対角化行列の解法
制約条件‖ｂ_ｊ（ω_ｋ）‖_２＝１を課したフォワードモデル型最小２乗型同時対角化問題を周波数ビンω_ｋ毎に解くことによって、対角化行列Ｂ（ω_ｋ）、即ち、混合行列を求める。
評価量を式（２４）に示す。式（２４）はフォワードモデル型最小２乗型同時対角化問題として知られている。

【0046】

【数11】

ここで、γ_ｉはラグランジェの未定乗数を表す。

【0047】

ベクトル表現を用いると、式（２４）の評価関数は式（２６）のように表現することができる。ここで、ｒ_ｘ（ω_ｋ，ｍ）、Ｇ（ω_ｋ）、ｄ（ω_ｋ，ｍ）、Ｇ（ω_ｋ）ｄ（ω_ｋ，ｍ）はそれぞれ式（２８）〜（３１）により表される。ただし、ｖｅｃ｛Ａ｝は行列Ａの列を積み重ね

れぞれ表す。λ_ｉは対角行列Ａのｉ番目の要素を表す。

【0048】

【数12】

【0049】

反復法によって混合行列Ｂ（ω_ｋ）を求めるためにｚ_ｉ（ω_ｋ）とＴ（ω_ｋ）を式（３２）、（３３）によって定義され、ステップＳ１０３において作成される。

【0050】

【数13】

【0051】

ステップＳ１０４においてｇ_ｊ（ω_ｋ）を求める際、ｇ_ｊ（ω_ｋ）（ｊ≠ｉ）を式（３４）のように、ｚ_ｊ（ω_ｋ）を式（３４）のようにそれぞれ定数に設定して、Ｆ_ｊ（ω_ｋ）を計算する。ここで、ｇ_ｊ（ω_ｋ）はＧ（ω_ｋ）のｊ番目の列ベクトルである。

【0052】

【数14】

【0053】

式（２６）の制約条件付きフォワードモデル型最小２乗型同時対角化問題を式（３５）のように書き直すことができ、その近似解ｇ_ｊ（ω_ｋ）はステップＳ１０５において式（３６）のラグランジェの未定乗数法によって求められる。

【0054】

【数15】

ここで、ｕｎｖｅｃ｛Ａ｝は、Ｊ^２×１の列ベクトルＡをＪ×Ｊの行列に変換することを表す。

【0055】

ステップＳ１０７においてｊ＝１，２，…，Ｎについて誤差の限界がε_Ｇの近似解ｇ_ｊ（ω_ｋ）を反復法によって推定した後、ステップＳ１０８において累乗法を１回用いて式（３７）を最小にするｂ_ｊ（ω_ｋ）を算出する。次いで、ステップＳ１０９において式（３８）のようにＢ（ω_ｋ）を特異値分解する。ここで、ｔｒ［Ａ］は行列Ａのトレースを表す。

【0056】

【数16】

【0057】

ステップＳ１１０においてＢ（ω_ｋ）^ＨＢ（ω_ｋ）のランクがＮ未満のとき、式（２１）の評価量ｅ（ω_ｋ）を零にするために、Ｂ（ω_ｋ）を式（４４）の行列によって置き換える。ここで、ステップＳ１１１において正規直交基底ｖ_１（ω_ｋ），ｖ_２（ω_ｋ），…，ｖ_ｒ（ω_ｋ）によって張られる空間に直交する空間の正規直交基底ｖ_ｒ＋１（ω_ｋ），ｖ_ｒ＋２（ω_ｋ），…，ｖ_Ｎ（ω_ｋ）、同様に、ステップＳ１１２において正規直交基底ｕ_１（ω_ｋ），ｕ_２（ω_ｋ），…，ｕ_ｒ（ω_ｋ）によって張られる空間に直交する空間の正規直交基底ｕ_ｒ＋１（ω_ｋ），ｕ_ｒ＋２（ω_ｋ），…，ｕ_Ｎ（ω_ｋ）はそれぞれ求められる。

【0058】

【数17】

ここで、δ（ω_ｋ）＞０とする。ステップＳ１１３において追加される特異値δ（ω_ｋ）は、式（４４）の右辺の√Ｎ／（√Ｎ＋δ（ω_ｋ）Ｎ）によって条件ｔｒ［Σ（ω_ｋ）］＝√Ｎを満足するように設定される。

【0059】

２．３ 対角行列の解法
ステップＳ１１４においてＢ（ω_ｋ）から分離行列Ｗ（ω_ｋ）を式（２２）の最小２乗型一般化逆行列によって求める。式（１４）、（１５）よりΛ（ω_ｋ，ｍ）の左からＷ（ω_ｋ）Ｂ（ω_ｋ）、右からＢ（ω_ｋ）^ＨＷ（ω_ｋ）^Ｈをそれぞれ乗算すると、式（４７）を得る。

【数18】

【0060】

Ｐ_ｘ（ω_ｋ，ｍ）の推定値を使用して誤差Ψ（ω_ｋ，ｍ）を式（４８）で定義すると、式（４９）のバックワード型最小２乗型同時対角化問題に最小２乗法を適用すると、対角行列Λ（ω_ｋ，ｍ）はステップＳ１１５において式（５０）で推定される。式（５０）においてｄｉａｇ［Ａ］は行列Ａの対角行列を表す。

【0061】

【数19】

【0062】

誤差の限界がε_Ｃの近似解ｇ_ｊ（ω_ｋ）と近似解Λ（ω_ｋ，ｍ）を推定するまで、上記のアルゴリズムはステップＳ１１６において繰り返される。

【0063】

３．パーミュテーション問題の解法
ステップＳ１１７において基準周波数ビンを複数選択し、基準周波数ビン間において電力比の相関に基づきパーミュテーション行列を推定するためにΞ（ω_ｋ）を式（５３）によりステップＳ１１８で算出する。

【0064】

【数20】

ここで、Ｔｒ（・）は行列のトレースを表す。また、Ｑは行列の各行に１となる要素が１箇所、その他の要素は０で、１となる要素の位置が他の行と重複しない行列の集合である。

【0065】

ステップＳ１１９において基準周波数ビン間で電力比の相関が最も大きいパーミュテーション行列を式（５４）によって推定する。

【0066】

ステップＳ１２０において複数の基準周波数ビンから１つの基準周波数ビンを任意に選択し、ステップＳ１２１において全ての周波数ビン間で電力比の相関が最も大きいパーミュテーション行列を式（５４）によって推定する。

【0067】

ステップＳ１２２において選択された全ての基準周波数ビンが１度選択されるまで、上記のパーミュテーション行列の推定手順を繰り返す。この結果、基準周波数ビンを除き、各周波数ビンに複数のパーミュテーション行列が割り当てられることになる。

【0068】

ステップＳ１２３において多数決によって周波数ビンに割り当てられるパーミュテーション行列を決定する。ただし、多数決によってパーミュテーション行列が決定できない場合には、最も相関値が大きいパーミュテーション行列を採用する。式（５５）のように観測信号ｘ（ω_ｋ，ｍ）にΠ（ω_ｋ）Ｗ（ω_ｋ）を左から乗算して分離信号ｙ（ω_ｋ，ｍ）を得る。

【実施例1】

【0069】

４．１ 評価データ
図５のように４．４５×３．５５×２．５メートルの部屋に３個の信号源（スピーカ）１１、１２、１３を半径１．２メートルの円の円周上に、円の中心に位置する一辺が２０センチメートルの正三方形の頂点に３個のマイクロホン２１、２２、２３をそれぞれ配置した。尚、図５は信号源（スピーカ）とマイクロホンの位置関係を示す平面図である。部屋の残響時間は１００ミリ秒から９００ミリ秒に設定し、標本化周波数８ｋＨｚ、量子化ビット数１６ビットで信号源とマイクロホンの間のインパルス応答は人工的に発生させた。実験条件は、１０００秒の音声データ、Ｋ＝８１９２点の短時間フーリエ変換、エポック当たり重複率９９％の２３個のフレームの使用、窓関数にはハニング窓を用いた。ＳＮＲは５ｄＢ間隔で０〜３０ｄＢの範囲で変化させた。マイクロホン２１、２２、２３のＳＮＲの設定方法については４．２で説明する。本発明に係るブラインド信号分離方法では、ε_Ｇ＝ε_Ｃ＝１０^−６、δ（ω_ｋ）＝σ_ｒ（ω_ｋ）を用いている。スケーリング問題は周波数ビン毎に分離行列の行ベクトルを正規化することによって解法した。Ｃ言語で作成したプログラムをインテル製コアｉ７−２６００ ３．４ＧＨｚプロセッサを用いて実行した。信号源信号からマイクロホンまでの経路は時不変のインパルス応答で、因果的で非最小位相系であるので、因果的な分離行列を実現するために、Π（ω_ｋ）^−１Ｄ（ω_ｋ）^−１Ｗ（ω_ｋ）にｅ^−ｊπｋを乗算した後、逆離散フーリエ変換をして分離フィルタのインパルス応答を得た。

【0070】

４．２ 評価指標
ブラインド信号分離方法の信号分離性能を次の方法で評価した．式（５６）によって観測信号における所望信号源信号と干渉信号の電力の比、式（５７）によって出力信号における所望信号源信号と干渉信号の電力の比をそれぞれ計算し、ブラインド信号分離装置の各出力の信号分離性能を求める。各出力の平均を信号分離性能とした。γ_ｉｊ（ｔ）は式（５８）のΓ（ω_ｋ）のｉ行ｊ列の要素を、ｗ_ｉｊ（ｔ）はＷ（ω_ｋ）の要素をそれぞれ離散逆フーリエ変換したものである。また、分離行列の推定アルゴリズムにおいて収束に要した反復回数と計算時間も評価指標とする。ＳＮＲは、最適な分離行列ｅ^−ｊπｋＤ（ω_ｋ）^−１（Ｈ（ω_ｋ）^ＨＨ（ω_ｋ））^−１Ｈ（ω_ｋ）^Ｈとパーミュテーション行列Π_ｏｐｔ（ω_ｋ）を使用して観測信号から信号源信号を分離した後、分離信号ｙ_ｉ（ｔ）に含まれる雑音と干渉信号の電力と所望信号源信号の電力の比によって計算した。最適なパーミュテーション行列Π_ｏｐｔ（ω_ｋ）は式（５９）によって求めた。また、非ブラインド法は、受信信号を使用して分離行列を計算した後、混合行列が入手可能であるとして、式（６０）によってパーミュテーション行列を求めた。即ち、推定した分離行列に最適なパーミュテーション行列を求めることになり、ブラインド信号分離装置の性能の上限を与えることになる。

【0071】

【数22】

ここで、Ｃ_ｏｐｔ（ω_ｋ）＝ｅ^−ｊπｋＤ（ω_ｋ）^−１（Ｈ（ω_ｋ）^ＨＨ（ω_ｋ））^−１Ｈ（ω_ｋ）^ＨＨ（ω_ｋ）、Ｃ（ω_ｋ）＝Ｗ（ω_ｋ）Ｈ（ω_ｋ）である。

【0072】

４．３ 評価対象
勾配法を用いたバックワードモデル型ブラインド信号分離方法（非特許文献１）、最小２乗型同時対角化問題の解法を用いた２種類のフォワードモデル型ブラインド信号分離方法（非特許文献２、非特許文献３）を比較対象とする。従来のブラインド信号分離方法（非特許文献１、非特許文献２、非特許文献３）と本発明に係るブラインド信号分離方法における分離行列の推定精度を比較するため、パーミュテーション行列の推定法は共通の手法（非特許文献４）を使用した。尚、基準周波数ビンの番号には６１４を用いた。

【0076】

４．４ 評価結果
部屋の残響時間と信号分離性能の関係を図６に、ＳＮＲと信号分離性能の関係を図７にそれぞれ示す。尚、図６においてＳＮＲは２０ｄＢに設定している。両図において太字の数字が最も優れた性能を表している。両図の信号分離性能から明らかなように、計算時間では非特許文献１より劣るものの、本発明に係るブラインド信号分離方法が従来のブラインド信号分離方法よりも最も高い信号分離性能（高い出力ＳＩＲ）を最も少ない反復回数で得ることができた。この要因はラグランジェの未定乗数法を最小２乗型同時対角化問題に導入したこと、推定した混合行列がランク落ちしていた場合、補空間を補い分離行列を推定したことが高い信号分離性能の実現に貢献したと考えられる。また、非ブラインド法の出力ＳＩＲ、即ち、ブラインド信号分離装置の上限に近い値を、本発明に係るブラインド信号分離装置が実現できることが分かる。

【実施例2】

【0077】

５．１ 評価データ
図５のように４．４５×３．５５×２．５メートルの部屋に３個の信号源（スピーカ）１１、１２、１３を半径１．２メートルの円の円周上に、円の中心に位置する一辺が２０センチメートルの正三方形の頂点に３個のマイクロホン２１、２２、２３をそれぞれ配置した。部屋の残響時間は７００ミリ秒に設定し、標本化周波数８ｋＨｚ、量子化ビット数１６ビットで信号源とマイクロホンの間のインパルス応答は人工的に発生させた。実験条件は、１０００秒の音声データ、Ｋ＝８１９２点の短時間フーリエ変換、エポック当たり重複率８０％の２個のフレームの使用、窓関数にはハニング窓を用いた。ＳＮＲは２０ｄＢに設定した。本発明に係るブラインド信号分離方法では、ε_Ｇ＝ε_Ｃ＝１０^−６、δ（ω_ｋ）＝σ_ｒ（ω_ｋ）、基準周波数ビンの番号は６１６、６１７、６１８を用いている。スケーリング問題は周波数ビン毎に分離行列の行ベクトルを正規化することによって解法した。信号源信号からマイクロホンまでの経路は時不変のインパルス応答で、因果的で非最小位相系であるので、因果的な分離行列を実現するために、Π（ω_ｋ）^−１Ｄ（ω_ｋ）^−１Ｗ（ω_ｋ）にｅ^−ｊπｋ乗算した後、逆離散フーリエ変換をして分離フィルタのインパルス応答を得た。

【0078】

５．２ 評価指標
信号源とマイクロホンの個数が共に３である場合、式（６１）に示す６種類のパーミュテーション行列の何れか１つが各周波数ビンに割り当てられる。割り当てられたパーミュテーション行列が、任意のパーミュテーション行列に一致する割合と信号分離性能を計算する。

【0079】

【数23】

【0080】

５．３ 評価対象
同一信号源から発生した信号の周波数ビン間の電力比に相関があることを利用したパーミュテーション問題の解法（非特許文献４）を比較対象とする。基準周波数ビンの番号には６１４を用いた。従来のブラインド信号分離方法（非特許文献１、非特許文献２、非特許文献３）と本発明に係るブラインド信号分離方法におけるパーミュテーション行列の推定精度を比較するため、分離行列の推定法は共通の手法（本発明に係る分離行列推定法）を使用した。

【0081】

５．４ 評価結果
図８と図９に各周波数ビンに割り当てられたパーミュテーション行列の番号を×印で示す。パーミュテーション行列Π_ｉとパーミュテーション行列の番号ｉの関係を式（６１）に示している。図８と図９では、各周波数ビンでパーミュテーション行列の番号３に割り当てられると未知信号源への割り当てが揃うことになる。したがって、番号３を除く他の番号への割り当ては間違いになる。低周波数帯域（０〜２ｋＨｚ）と全周波数帯域におけるパーミュテーション行列の割り当て結果を図１０にまとめる。本発明に係るパーミュテーション行列の推定法が非特許文献４の方法に比べ正答率が向上していることが分かる。また、信号分離性能においても、本発明に係るパーミュテーション行列の推定法が高い出力ＳＩＲを達成することができた。

【符号の説明】

【0082】

１１〜１Ｎ…信号源、２１〜２Ｊ…マイクロホン、３１…短時間フーリエ変換、４０…ブラインド信号分離システム、５１、５２、…、５Ｋ…最小２乗型同時対角化問題の解法、６０…パーミュテーション問題の解法、７１、７２、…、７Ｋ…畳み込み演算、８０…離散逆フーリエ変換と重複加算

【図1】

【図2】

【図3】

【図4】

【図5】

【図6】

【図7】

【図8】

【図9】

【図10】