(57)【要約】【課題】様々な遅延量を持つ混合音信号について遅延量を考慮した高精度な分離処理を可能とする構成を提供する。【解決手段】入力音信号を時間周波数領域に変換して生成した観測信号スペクトログラムを時間周波数領域において畳み込み混合された観測信号として解釈し、畳み込み混合を解く独立成分分析により信号分離を行なう、または観測信号スペクトログラムに対する短時間フーリエ変換（ＳＴＦＴ）により生成したモジュレーション・スペクトログラムを瞬時混合として解釈し、瞬時混合を解く独立成分分析によって信号分離結果を生成する。本構成により直接波、反射波など様々な遅延量を持つ混合音信号について遅延量を考慮した高精度な分離処理が実現される。【選択図】図１５
2011215649.tif 000002
【特許請求の範囲】
【請求項１】
複数の信号が混合した信号を入力して個別の信号に分離する信号分離装置であり、
入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換手段と、
前記信号変換手段の生成した観測信号スペクトログラムから信号分離結果を生成する信号分離手段を有し、
前記信号分離手段は、
前記観測信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により、信号分離結果を生成する構成であることを特徴とする信号分離装置。
【請求項２】
前記瞬時混合ＩＣＡを適用した処理は、時間周波数領域の観測信号と分離行列から時間周波数領域の分離信号を生成し、生成した時間周波数領域の分離信号と、多次元確率密度関数から導出される多次元スコア関数によって計算される分離行列とがほぼ収束するまで分離行列を修正し、修正した分離行列を適用して時間周波数領域の分離信号を生成する処理であることを特徴とする請求項１に記載の信号分離装置。
【請求項３】
前記信号分離手段は、
複数の信号入力源の観測信号各々に対応して生成される複数の観測信号スペクトログラムシフトセットを積み重ねた複数チャンネル対応の観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡを適用して信号分離結果を生成することを特徴とする請求項１に記載の信号分離装置。
【請求項４】
前記信号分離手段は、
前記シフトの際に生じた隙間をゼロまたはゼロに近い値、または前記観測信号スペクトログラムの両端の値をコピーして設定して、前記観測信号スペクトログラムシフトセットを生成することを特徴とする請求項１に記載の信号分離装置。
【請求項５】
前記信号分離手段は、
前記シフトをシフトではみ出した一端のデータを他端にコピーする巡回シフト処理を実行することを特徴とする請求項１に記載の信号分離装置。
【請求項６】
前記信号分離手段は、
最小シフト量を０、最大シフト量を観測信号から分離結果を生成する際のフレームタップ数［Ｌ'］として設定した複数のシフトデータを生成し、生成した異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成することを特徴とする請求項１に記載の信号分離装置。
【請求項７】
前記信号分離手段は、
周波数に応じて前記フレームタップ数［Ｌ'］を変更して前記観測信号スペクトログラムシフトセットを生成することを特徴とする請求項１に記載の信号分離装置。
【請求項８】
前記信号分離手段は、
前記観測信号スペクトログラムに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により第１の信号分離結果を生成し、該第１の信号分離結果から、どの音源にも対応しないと判定される不要チャンネル除去処理を実行し、該除去処理後に残存する観測信号スペクトログラムをフレーム方向へシフトさせて観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡを適用して信号分離結果を生成する構成であることを特徴とする請求項１に記載の信号分離装置。
【請求項９】
複数の信号が混合した信号を入力して個別の信号に分離する信号分離装置であり、
入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換手段と、
前記信号変換手段の生成した観測信号スペクトログラムから信号分離結果を生成する信号分離手段を有し、
前記信号分離手段は、
前記観測信号スペクトログラムに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により信号分離結果Ｙ１〜Ｙｎを生成し、
信号分離結果Ｙ１〜Ｙｎの各々に対応する信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により残響除去処理を実行し、残響除去済みスペクトログラムの統合処理によって、残響を除去した信号分離結果を生成する構成であることを特徴とする信号分離装置。
【請求項１０】
前記瞬時混合ＩＣＡを適用した処理は、時間周波数領域の観測信号と分離行列から時間周波数領域の分離信号を生成し、生成した時間周波数領域の分離信号と、多次元確率密度関数から導出される多次元スコア関数によって計算される分離行列とがほぼ収束するまで分離行列を修正し、修正した分離行列を適用して時間周波数領域の分離信号を生成する処理であることを特徴とする請求項９に記載の信号分離装置。
【請求項１１】
信号分離装置において、複数の信号が混合した信号を入力して個別の信号に分離する信号分離方法であり、
信号変換手段が、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換ステップと、
信号分離手段が、前記信号変換ステップにおいて生成した観測信号スペクトログラムから信号分離結果を生成する信号分離ステップを有し、
前記信号分離ステップは、
前記観測信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により、信号分離結果を生成するステップであることを特徴とする信号分離方法。
【請求項１２】
信号分離装置において、複数の信号が混合した信号を入力して個別の信号に分離する信号分離方法であり、
信号変換手段が、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換ステップと、
信号分離手段が、前記信号変換ステップにおいて生成した観測信号スペクトログラムから信号分離結果を生成する信号分離ステップを有し、
前記信号分離ステップは、
前記観測信号スペクトログラムに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により信号分離結果Ｙ１〜Ｙｎを生成し、
信号分離結果Ｙ１〜Ｙｎの各々に対応する信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により残響除去処理を実行し、残響除去済みスペクトログラムの統合処理によって、残響を除去した信号分離結果を生成するステップであることを特徴とする信号分離方法。
【請求項１３】
信号分離装置において、複数の信号が混合した信号を入力して個別の信号に分離する信号分離処理を実行させるコンピュータ・プログラムであり、
信号変換手段に、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成させる信号変換ステップと、
信号分離手段に、前記信号変換ステップにおいて生成した観測信号スペクトログラムから信号分離結果を生成させる信号分離ステップを有し、
前記信号分離ステップは、
前記観測信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により、信号分離結果を生成させるステップであることを特徴とするコンピュータ・プログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、信号分離装置、および信号分離方法、並びにコンピュータ・プログラムに関する。さらに、詳細には、本発明は、複数の信号が混合された信号を独立成分分析（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ；ＩＣＡ）を用いて信号毎に分離する信号分離装置、および信号分離方法、並びにコンピュータ・プログラムに関する。
【背景技術】
【０００２】
複数の原信号が未知の係数によって線形に混合されているときに、統計的独立性のみを用いて原信号を分離・復元するという独立成分分析（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ；ＩＣＡ）の手法が信号処理の分野で注目されている。この独立成分分析を応用することで、例えば話者とマイクロホンとが離れた場所にあり、マイクロホンで話者の音声以外の音を拾ってしまうような状況でも、音声信号を分離・復元することが可能となる。
【０００３】
ＩＣＡとは、多変量分析の一種であり、信号の統計的な性質を利用して多次元信号を分離する手法のことである。ＩＣＡ自体の詳細については、例えば非特許文献１（「入門・独立成分分析」（村田昇著、東京電機大学出版局））などを参照されたい。
【０００４】
まず、時間周波数領域の独立成分分析を用いて、複数の信号（特に音信号）が混合された信号を時間周波数領域で分離する方法について説明し、次にその方法が持つ問題点について述べる。図１に示すように、Ｎ個の音源（信号源）から異なる音が鳴っていて、それらをｎ個のマイク（センサー）で観測するという状況を考える。複数の音源が発した音（原信号）がマイクに届く場合、マイクの取得する音は直接波、反射波が含まれ、各音源との距離に基づく時間遅れなどがあるため、ある１つのマイクｊ（ただし１≦ｊ≦ｎ）で観測される信号（観測信号）は以下に示す式［１．１］のように、原信号と伝達関数との畳み込み演算を全音源について総和した式として表わせる（以下では「畳み込み混合」と呼ぶ）。さらに、全てのマイク１〜ｎについての観測信号を一つの式で表わすと、式［１．２］のように表わせる。ただし、ｘ（ｔ），ｓ（ｔ）はそれぞれｘｋ（ｔ），ｓｋ（ｔ）を要素とする列ベクトルであり、Ａ［ｌ］はａｋｊを要素とするｎ×Ｎの行列である。（以降では、ｎ＝Ｎとする。）
【数１】
2011215649.tif 000003
【０００５】
このような畳み混合を解く手法として、以下の２つの方法が知られている。
（１）時間領域で畳み込み混合を直接解く。（時間領域逆畳み込み）
（２）観測信号を時間周波数領域に変換し、瞬時混合問題として解く。
以下では、それぞれの方法について説明する。
【０００６】
（１）時間領域で畳み込み混合を直接解く（時間領域逆畳み込み）手法について
上記した式［１．２］の畳み込みを解くために、以下に示す式［２．１］のような、観測信号の畳み込み混合の式を用意する。
【数２】
2011215649.tif 000004
【０００７】
上記の式［２．１］のような、観測信号の畳み込み混合の式を用意し、そして、分離結果ｙ（ｔ）の成分であるｙ１（ｔ）〜ｙｎ（ｔ）が全てのｔについて最も独立となるように、分離行列Ｗ［０］〜Ｗ［Ｌ'］を決める（以降、Ｗ［０］〜Ｗ［Ｌ'］をまとめて分離フィルターと呼ぶ）。そのためには、式［２．１］〜式［２．４］を、分離行列および分離結果が収束するまで繰り返す（以降、このような繰り返しを「学習」と呼ぶ。また、分離行列を更新する式やΔＷを算出する式などを「学習規則」と呼ぶ）。なお、式［２．３］のＥｔ［］はｔについての平均を表わす。同式のφはスコア関数または活性化関数と呼ばれる関数である。なお、時間領域での畳み込み混合を解く式の詳細については、例えば、非特許文献２（「詳解独立成分分析」（Ａａｐｏ  Ｈｙｖａｒｉｎｅｎｎ  ほか著。東京電機大学出版局）１９．２畳み込み混合の暗中分離，１９．２．３  自然勾配法）を参照されたい。
【０００８】
（２）観測信号を時間周波数領域に変換し、瞬時混合問題として解く手法について
時間領域の畳み込み混合は、時間周波数領域では瞬時混合で表わされることが知られており、その特徴を利用したのが時間周波数領域のＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）である。時間周波数領域ＩＣＡ自体については、上述の非特許文献２（「詳解独立成分分析」（Ａａｐｏ  Ｈｙｖａｒｉｎｅｎｎ  ほか著。東京電機大学出版局「１９．２．４フーリエ変換法」）や、特許文献１（特開２００６−２３８４０９「音声信号分離装置・雑音除去装置および方法」）などを参照されたい。
【０００９】
時間周波数領域の独立成分分析では、上記式［１．２］のｘ（ｔ）からＡおよびｓ（ｔ）を直接推定するのではなく、ｘ（ｔ）を時間周波数領域の信号に変換し、Ａおよびｓ（ｔ）に対応する信号を時間周波数領域で推定する。以下では、主に本発明と関係ある点を説明する。上記式［１．２］の両辺を短時間フーリエ変換すると、以下に示す式［３．１］が近似的に得られる。すなわち、信号ベクトルｘ（ｔ）、ｓ（ｔ）を長さＬの窓で短時間フーリエ変換したものをそれぞれＸ（ω，ｔ），Ｓ（ω，ｔ）とし、行列Ａ（ｔ）を同様に短時間フーリエ変換したものをＡ（ω）とすると、時間領域の上記式［１．２］は時間周波数領域の下記式［３．１］で表すことができる。但し、ωは周波数ビンの番号を示し（１≦ω≦Ｍ）、ｔはフレーム番号を示す（１≦ｔ≦Ｔ）。時間周波数領域の独立成分分析では、式［３．１］のＳ（ω，ｔ）、Ａ（ω）を時間周波数領域で推定することになる。
【数３】
2011215649.tif 000005
【００１０】
上記式［３．１］において、ωは周波数ビンの番号、ｔはフレームの番号である。ωを固定すると、この式は瞬時混合と見なせる。そこで、観測信号を分離するには、式［３．５］のような式を用意し、Ｙ（ω，ｔ）の各成分が最も独立になるように分離行列Ｗ（ω）を決める。
【００１１】
なお、周波数ビンの個数は、本来は窓の長さＬと同一であり、各周波数ビンは、−Ｒ／２からＲ／２まで（Ｒはサンプリング周波数）をＬ等分したそれぞれの周波数成分を表す。但し、負の周波数成分は正の周波数成分の共役複素数であり、Ｘ（−ω）＝ｃｏｎｊ（Ｘ（ω））（ｃｏｎｊ（・）は共役複素数）として求めることができる。
時間周波数領域でＳ（ω，ｔ）、Ａ（ω）を推定するには、先ず、下記式（４）のような式を考える。式［３．５］において、Ｙ（ω，ｔ）はｙｋ（ｔ）を長さＬの窓で短時間フーリエ変換したＹｋ（ω，ｔ）を要素とする列ベクトルを表し、Ｗ（ω）はｗｉｊ（ω）を要素とするｎ行ｎ列の行列（分離行列）を表す。
【００１２】
従来の時間周波数領域ＩＣＡでは、パーミュテーション問題と呼ばれる、「どの成分がどのチャンネルに分離されるか」が周波数ビンごとに異なるという問題が発生していたが、この問題については、発明者自身の特許出願である上述の特許文献１（特開２００６−２３８４０９「音声信号分離装置・雑音除去装置および方法」）によって、ほぼ解決した。
【００１３】
本発明は、この公開されている特許出願である特開２００６−２３８４０９の理論的発展版であるため、以下では特開２００６−２３８４０９の特徴についても説明する。
【００１４】
従来、すなわち、上述の特許文献１（特開２００６−２３８４０９）で示される手法が開示される以前は、時間周波数領域の分離の式として、周波数ビンごとの式である［３．５］を用い、さらに周波数ビンごとに独立性を最大にするような分離行列Ｗ（ω）求めていた。
【００１５】
すなわち、ωを固定してｔを変化させたときにＹ１（ω，ｔ）〜Ｙｎ（ω，ｔ）が統計的に独立となる（実際には、独立性が最大となる）ようなＷ（ω）を求める。なお、後述するが時間周波数領域の独立成分分析ではパーミュテーション（ｐｅｒｍｕｔａｔｉｏｎ）及びスケーリングの不定性があるため、Ｗ（ω）＝Ａ（ω）−１以外にも解が存在する。統計的に独立となるＹ１（ω，ｔ）〜Ｙｎ（ω，ｔ）が全てのωについて得られたら、それらを逆フーリエ変換することで、時間領域の分離信号ｙ（ｔ）を得ることができる。
【００１６】
時間周波数領域における従来の独立成分分析の概略について説明する。ｎ個の音源が発するお互いに独立な原信号をｓ１〜ｓｎとし、それらを要素とするベクトルをｓとする。マイクロホンで観測される観測信号ｘは、原信号ｓに上記式［１．２］の畳み込み・混合演算を施したものである。次に、観測信号ｘに対して短時間フーリエ変換を施し、時間周波数領域の信号Ｘを得る。Ｘの要素をＸｋ（ω，ｔ）とすると、Ｘｋ（ω，ｔ）は複素数値をとる。Ｘｋ（ω，ｔ）の絶対値である｜Ｘｋ（ω，ｔ）｜を色の濃淡で表現した図をスペクトログラムという。スペクトログラムは、例えば、横軸をｔ（フレーム番号）、縦軸をω（周波数ビン番号）としてＸｋ（ω，ｔ）の絶対値である｜Ｘｋ（ω，ｔ）｜を色の濃淡で表現した図である。続いて、信号Ｘの各周波数ビンにＷ（ω）を乗算することで分離信号Ｙを得る。そして、分離信号Ｙを逆フーリエ変換することで時間領域の分離信号ｙを得る。
【００１７】
しかし、上述した時間周波数領域の独立成分分析では、信号の分離処理を周波数ビン毎に行っており、周波数ビンの間の関係は考慮していない。そのため、分離自体は成功しても、周波数ビンの間でスケーリング及び分離先の不統一が発生する可能性がある。このうち、スケーリングの不統一については、音源毎に観測信号を推定する方法により解決できる。一方、分離先の不統一とは、例えばω＝１ではＹ１にＳ１由来の信号が現れるのに対してω＝２ではＹ１にＳ２由来の信号が現れる、というような現象のことであり、パーミュテーション（置換）の問題と呼ばれている。
【００１８】
これに対し、特許文献１（特開２００６−２３８４０９）では、スペクトログラム全体での分離を表わす式である以下に示す数式［４．４］を用い、スペクトログラム全体での独立性を最大にする分離行列Ｗを求めるという手法を採用している。
【数４】
2011215649.tif 000006
【００１９】
具体的には、スペクトログラム全体での独立性として、式［４．５］で表わされるＫｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量Ｉ（Ｙ）を導入し、Ｉ（Ｙ）を最小にする分離行列Ｗを求めている。独立成分分析において、独立性をどのような尺度で表現するか、また、どのようなアルゴリズムで独立性を最大化するかについては、種々のバリエーションが存在する。その１つの手法として、独立性をＫｕｌｌｂａｃｋ−Ｌｅｉｂｌｅｒ情報量（ＫＬ情報量）がある。Ｋｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量Ｉ（Ｙ）は、スペクトログラムごとのエントロピーの総和から、スペクトログラム全体の同時エントロピーを引いたものであり、全てのスペクトログラムがお互いに独立となった場合に最小（理想的には０）となる。
【００２０】
ＫＬ情報量Ｉ（Ｙ）は前述したように上記式［４．５］のように定義される。この式［４．５］において、Ｈ（Ｙｋ）は各チャンネルについてのスペクトログラム１枚分のエントロピーを表し、Ｈ（Ｙ）は全チャンネルについてのスペクトログラム１枚分の同時エントロピーを表す。ｎ＝２のときのＨ（Ｙｋ）とＨ（Ｙ）との関係を図２に示す。図２において、Ｐ（Ｙｋ（ｔ））は、Ｙｋ（ｔ）の確率密度関数であり、Ｈ（Ｙｋ）は各チャンネルについてのスペクトログラム１枚分のエントロピーである。Ｋｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量Ｉ（Ｙ）は、スペクトログラムごとのエントロピー１１，１２の総和から、スペクトログラム全体の同時エントロピー１３を引いたものであり、全てのスペクトログラムがお互いに独立となった場合に最小（理想的には０）となる。
【００２１】
スペクトログラム全体でのＫＬ情報量Ｉ（Ｙ）を最小にするためには、式［５．１］〜［５．３］をＷおよびＹが収束するまで繰り返す。
【数５】
2011215649.tif 000007
【００２２】
なお、式［５．３］に出てくるΔＷ（ω），Ｗ（ω），Ｙ（ω，ｔ）は、それぞれΔＷ，Ｗ，Ｙ（ｔ）からω番目の周波数ビンに対応する要素を抽出した部分行列である。こうすることで、パーミュテーション問題のない分離結果を得ることが可能となった。
【００２３】
しかし、上記の２つの畳み混合を解く手法、すなわち、
（１）時間領域で畳み込み混合を直接解く。（時間領域逆畳み込み）
（２）観測信号を時間周波数領域に変換し、瞬時混合問題として解く。
この２つの手法には、課題がある。すなわち、
（１）時間領域で畳み込み混合を直接解く。（時間領域逆畳み込み）
この方法については、収束が遅いという問題がある。遅収束の理由は、分離フィルターの係数が変化すると波形全体が変化することや、分離フィルターの更新式の計算量がタップ数Ｌ'の２乗に比例することなどが挙げられる。そのため、フィルターのタップ数Ｌ'が大きい場合は、分離フィルターの初期値として、収束値にできる限り近い値を事前に求めておかないと、実用的な時間で分離することは難しい。実環境の残響に対応させるためには、少なくとも数千のオーダーのタップ数が必要であるため、（１）の方法では数千の２乗の計算量が必要である。
【００２４】
一方、（２）観測信号を時間周波数領域に変換し、瞬時混合問題として解く。この方法は、短時間フーリエ変換（ＳＴＦＴ）の窓長と分離精度との間にトレードオフが存在することが課題である。観測信号が長い残響、すなわち大きなタップ数の畳み込み混同である場合、それを時間周波数領域の瞬時混合で表わすためには、ＳＴＦＴの窓の長さ（＝タップ数）も大きくする必要がある。（窓長＜残響長の場合は、残響が複数のフレームにまたがるため、瞬時混合では表現できない。）しかし、窓長を長くしすぎると、かえって分離精度が落ちることが知られている。なお、このトレードオフについては、例えば、以下の文献を参照されたい。
特許文献２（特開２００３−２７１１６８「信号抽出方法および信号抽出装置、信号抽出プログラムとそのプログラムを記録した記録媒体」）
非特許文献３（「サブバンド処理によるブラインド音源分離に関する検討」荒木章子・Ｒｏｂｅｒｔ  Ａｉｃｈｎｅｒ・牧野昭二・西川剛樹・猿渡洋  日本音響学会講演論文集２００２年３月  ｐｐ．６１９−−６２０）
非特許文献４（「帯域分割型ＩＣＡを用いたＢｌｉｎｄ  Ｓｏｕｒｃｅ  Ｓｅｐａｒａｔｉｏｎにおける帯域分割数の最適化」西川剛樹・荒木章子・牧野昭二・猿渡洋  日本音響学会講演論文集２００１年３月  ｐｐ．５６９−−５７０）
【００２５】
窓長を長くすると分離精度が落ちる原因は、窓長を長くする（＝タップ数を大きくする）ほど、生成されたスペクトログラムの時間方向の変化、すなわち時間方向エンベロープの変化がなだらかになることである。時間周波数領域ＩＣＡは、エンベロープ同士の独立性に着目して観測信号を分離するが、なだらかなエンベロープ同士の独立性は急激に変化するエンベロープ同士の独立性に比べて低めに算出される傾向がある。つまり、異なる音源に由来するエンベロープ同士であっても「相関がある」と判定される可能性があるため、結果として分離精度が悪くなる。
【００２６】
上述したように、（２）観測信号を時間周波数領域に変換し、瞬時混合問題として解く手法における問題点は、短時間フーリエ変換（ＳＴＦＴ）の窓長と分離精度との間にトレードオフが存在することである。以下、窓長と分離精度とのトレードオフに関して発明者本人が行なった実験結果を示す。図３は、ＳＴＦＴの窓長と時間周波数領域ＩＣＡの分離精度の関係をプロットしたグラフである。
【００２７】
図３において、横軸はＳＴＦＴの窓長（６４，１２８，２５６，５１２，１０２４，２０４８，４０９６）、縦軸は分離精度の尺度であるｓｉｇｎａｌ−ｉｎｔｅｒｆｅｒｅｎｃｅ−ｒａｔｉｏ（ＳＩＲ）であり、実線は、（２）観測信号を時間周波数領域に変換し、瞬時混合問題として解く手法として、特開２００６−２３８４０９の方法で分離した結果のＳＩＲである（実験の詳細は後述する）。また、図３上段のグラフは波形ベースのＳＩＲであり、下段は周波数ビンベースのＳＩＲである。また、図４は、横軸の窓長を実際の秒数で表わしたグラフである。どちらのグラフも、中間に分離精度のピークがあることが分かる。（波形ベースでは窓長＝１０２４がピーク。周波数ビンベースでは窓長＝５１２がピーク。）
【００２８】
すなわち、時間周波数領域のＩＣＡでは、長い残響に対応させようとしてＳＴＦＴの窓を長くしても、ある程度を超えるとかえって分離性能が落ちるという問題がある。
【００２９】
以上をまとめると、独立成分分析（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ；ＩＣＡ）の手法の方式である以下の２つの方式、すなわち、
（１）時間領域で畳み込み混合を直接解く。（時間領域逆畳み込み）
（２）観測信号を時間周波数領域に変換し、瞬時混合問題として解く。
この２つの方式のいずれも、タップ数の大きい畳み込み混合に対しては分離精度が不十分であるという問題が存在する。
【００３０】
なお、「残響長よりも短い窓を用いてＳＴＦＴした場合は、スペクトログラム上でも依然として畳み込みが残る」との想定に対応した処理を開示した技術して非特許文献５（Ｓｅｒｖｉｅｒｅ，Ｃ．  Ｓｅｐａｒａｔｉｏｎ  ｏｆ  ｓｐｅｅｃｈ  ｓｉｇｎａｌｓ  ｕｎｄｅｒ  ｒｅｖｅｒｂｅｒａｎｔ  ｃｏｎｄｉｔｉｏｎｓ．Ｉｎ  Ｐｒｏｃ．ＥＵＳＩＰＣＯ０４，ｐｐ．１６９３-１６９６（２００４））がある。
【００３１】
上記非特許文献５では、観測信号を時間周波数領域上の畳み込み混合であると考え、それを解く方法として、時間周波数領域での逆畳み込みのアルゴリズムを提案している。すなわち、「時間周波数領域において畳み込み混合を直接解く」方式に近い処理である。しかし、この非特許文献５において開示しているアルゴリズムは、２入力・２出力、すなわち、音声信号の出力音源が２つ、入力部としてのマイクを２つとした場合に限られたものである。また、本文献では、分離および逆畳み込みを周波数ビンごとに個別に行なう構成であり、パーミュテーション（ｐｅｒｍｕｔａｔｉｏｎ）問題と呼ばれる、「どの成分がどのチャンネルに分離されるか」が周波数ビンごとに異なるという問題が発生する。
【００３２】
上述したように、複数の信号が混合された音信号の分離処理について開示した従来技術はいくつか存在するが、独立成分分析（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ；ＩＣＡ）を用いて信号毎の高精度な分離処理を実現する信号分離処理においては、
（１）窓長（＝分析フレームの長さ）を超える残響に対する対処、
（２）パーミュテーション（ｐｅｒｍｕｔａｔｉｏｎ）問題に対する対処、
（３）２つの入出力を超える入出力構成に対する対処、
これらの様々な問題に対する十分な解決策が提示されていないのが現状である。
【先行技術文献】
【特許文献】
【００３３】
【特許文献１】特開２００６−２３８４０９
【特許文献２】特開２００３−２７１１６８
【非特許文献】
【００３４】
【非特許文献１】「入門・独立成分分析」（村田昇著、東京電機大学出版局）特開２００２−３４２１９８
【非特許文献２】「詳解独立成分分析」（Ａａｐｏ  Ｈｙｖａｒｉｎｅｎｎ  ほか著。東京電機大学出版局）
【非特許文献３】「サブバンド処理によるブラインド音源分離に関する検討」荒木章子・Ｒｏｂｅｒｔ  Ａｉｃｈｎｅｒ・牧野昭二・西川剛樹・猿渡洋  日本音響学会講演論文集２００２年３月  ｐｐ．６１９−−６２０
【非特許文献４】「帯域分割型ＩＣＡを用いたＢｌｉｎｄ  Ｓｏｕｒｃｅ  Ｓｅｐａｒａｔｉｏｎにおける帯域分割数の最適化」西川剛樹・荒木章子・牧野昭二・猿渡洋  日本音響学会講演論文集２００１年３月  ｐｐ．５６９−−５７０
【非特許文献５】Ｓｅｒｖｉｅｒｅ，Ｃ．  Ｓｅｐａｒａｔｉｏｎ  ｏｆ  ｓｐｅｅｃｈ  ｓｉｇｎａｌｓ  ｕｎｄｅｒ  ｒｅｖｅｒｂｅｒａｎｔ  ｃｏｎｄｉｔｉｏｎｓ．Ｉｎ  Ｐｒｏｃ．ＥＵＳＩＰＣＯ０４，ｐｐ．１６９３-１６９６（２００４）
【発明の開示】
【発明が解決しようとする課題】
【００３５】
本発明は、このような状況に鑑みてなされたものであり、複数の信号が混合された音信号を独立成分分析（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ；ＩＣＡ）を用いて信号毎の高精度な分離処理を実現する信号分離装置、および信号分離方法、並びにコンピュータ・プログラムを提供することを目的とするものであり、特に、タップ数の大きい畳み込み混合に対する分離精度を向上させた信号分離装置、および信号分離方法、並びにコンピュータ・プログラムを提供することを目的とする。
【課題を解決するための手段】
【００３６】
本発明の第１の側面は、
複数の信号が混合した信号を入力して個別の信号に分離する信号分離装置であり、
入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換手段と、
前記信号変換手段の生成した観測信号スペクトログラムから信号分離結果を生成する信号分離手段を有し、
前記信号分離手段は、
前記観測信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により、信号分離結果を生成する構成であることを特徴とする信号分離装置にある。
【００３７】
さらに、本発明の信号分離装置の一実施態様において、前記瞬時混合ＩＣＡを適用した処理は、時間周波数領域の観測信号と分離行列から時間周波数領域の分離信号を生成し、生成した時間周波数領域の分離信号と、多次元確率密度関数から導出される多次元スコア関数によって計算される分離行列とがほぼ収束するまで分離行列を修正し、修正した分離行列を適用して時間周波数領域の分離信号を生成する処理であることを特徴とする。
【００３８】
さらに、本発明の信号分離装置の一実施態様において、前記信号分離手段は、複数の信号入力源の観測信号各々に対応して生成される複数の観測信号スペクトログラムシフトセットを積み重ねた複数チャンネル対応の観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡを適用して信号分離結果を生成することを特徴とする。
【００３９】
さらに、本発明の信号分離装置の一実施態様において、前記信号分離手段は、前記シフトの際に生じた隙間をゼロまたはゼロに近い値、または前記観測信号スペクトログラムの両端の値をコピーして設定して、前記観測信号スペクトログラムシフトセットを生成することを特徴とする。
【００４０】
さらに、本発明の信号分離装置の一実施態様において、前記信号分離手段は、前記シフトをシフトではみ出した一端のデータを他端にコピーする巡回シフト処理を実行することを特徴とする。
【００４１】
さらに、本発明の信号分離装置の一実施態様において、前記信号分離手段は、最小シフト量を０、最大シフト量を観測信号から分離結果を生成する際のフレームタップ数［Ｌ'］として設定した複数のシフトデータを生成し、生成した異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成することを特徴とする。
【００４２】
さらに、本発明の信号分離装置の一実施態様において、前記信号分離手段は、周波数に応じて前記フレームタップ数［Ｌ'］を変更して前記観測信号スペクトログラムシフトセットを生成することを特徴とする。
【００４３】
さらに、本発明の信号分離装置の一実施態様において、前記信号分離手段は、前記観測信号スペクトログラムに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により第１の信号分離結果を生成し、該第１の信号分離結果から、どの音源にも対応しないと判定される不要チャンネル除去処理を実行し、該除去処理後に残存する観測信号スペクトログラムをフレーム方向へシフトさせて観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡを適用して信号分離結果を生成する構成であることを特徴とする。
【００４４】
さらに、本発明の第２の側面は、
複数の信号が混合した信号を入力して個別の信号に分離する信号分離装置であり、
入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換手段と、
前記信号変換手段の生成した観測信号スペクトログラムから信号分離結果を生成する信号分離手段を有し、
前記信号分離手段は、
前記観測信号スペクトログラムに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により信号分離結果Ｙ１〜Ｙｎを生成し、
信号分離結果Ｙ１〜Ｙｎの各々に対応する信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により残響除去処理を実行し、残響除去済みスペクトログラムの統合処理によって、残響を除去した信号分離結果を生成する構成であることを特徴とする信号分離装置にある。
【００４５】
さらに、本発明の信号分離装置の一実施態様において、前記瞬時混合ＩＣＡを適用した処理は、時間周波数領域の観測信号と分離行列から時間周波数領域の分離信号を生成し、生成した時間周波数領域の分離信号と、多次元確率密度関数から導出される多次元スコア関数によって計算される分離行列とがほぼ収束するまで分離行列を修正し、修正した分離行列を適用して時間周波数領域の分離信号を生成する処理であることを特徴とする。
【００４６】
さらに、本発明の第３の側面は、
信号分離装置において、複数の信号が混合した信号を入力して個別の信号に分離する信号分離方法であり、
信号変換手段が、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換ステップと、
信号分離手段が、前記信号変換ステップにおいて生成した観測信号スペクトログラムから信号分離結果を生成する信号分離ステップを有し、
前記信号分離ステップは、
前記観測信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により、信号分離結果を生成するステップであることを特徴とする信号分離方法にある。
【００４７】
さらに、本発明の第４の側面は、
信号分離装置において、複数の信号が混合した信号を入力して個別の信号に分離する信号分離方法であり、
信号変換手段が、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換ステップと、
信号分離手段が、前記信号変換ステップにおいて生成した観測信号スペクトログラムから信号分離結果を生成する信号分離ステップを有し、
前記信号分離ステップは、
前記観測信号スペクトログラムに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により信号分離結果Ｙ１〜Ｙｎを生成し、
信号分離結果Ｙ１〜Ｙｎの各々に対応する信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により残響除去処理を実行し、残響除去済みスペクトログラムの統合処理によって、残響を除去した信号分離結果を生成するステップであることを特徴とする信号分離方法にある。
【００４８】
さらに、本発明の第５の側面は、
信号分離装置において、複数の信号が混合した信号を入力して個別の信号に分離する信号分離処理を実行させるコンピュータ・プログラムであり、
信号変換手段に、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成させる信号変換ステップと、
信号分離手段に、前記信号変換ステップにおいて生成した観測信号スペクトログラムから信号分離結果を生成させる信号分離ステップを有し、
前記信号分離ステップは、
前記観測信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により、信号分離結果を生成させるステップであることを特徴とするコンピュータ・プログラムにある。
【００４９】
なお、本発明のコンピュータ・プログラムは、例えば、様々なプログラム・コードを実行可能なコンピュータ・システムに対して、コンピュータ可読な形式で提供する記憶媒体、通信媒体、例えば、ＣＤやＦＤ、ＭＯなどの記録媒体、あるいは、ネットワークなどの通信媒体によって提供可能なコンピュータ・プログラムである。このようなプログラムをコンピュータ可読な形式で提供することにより、コンピュータ・システム上でプログラムに応じた処理が実現される。
【００５０】
本発明のさらに他の目的、特徴や利点は、後述する本発明の実施例や添付する図面に基づくより詳細な説明によって明らかになるであろう。なお、本明細書においてシステムとは、複数の装置の論理的集合構成であり、各構成の装置が同一筐体内にあるものには限らない。
【発明の効果】
【００５１】
本発明の一実施例の構成によれば、複数の信号が混合した入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成し、観測信号スペクトログラムから信号分離結果を生成する信号分離処理において、観測信号スペクトログラムを時間周波数領域において畳み込み混合された観測信号として解釈し、時間周波数領域の畳み込み混合を解く処理により信号分離結果を生成する、あるいは、観測信号スペクトログラムに対する時間方向の短時間フーリエ変換（ＳＴＦＴ）によりモジュレーション・スペクトログラムを生成してモジュレーション・スペクトログラムを瞬時混合として解釈して信号分離結果を生成する構成としたので、直接波、反射波など様々な遅延量を持つ混合された音信号について、遅延量を考慮した高精度な分離処理が実現される。
【図面の簡単な説明】
【００５２】
【図１】複数の信号が混合された音信号の分離処理に適用する音情報の取得構成例について説明する図である。
【図２】各チャンネルについてのスペクトログラム１枚分のエントロピーＨ（Ｙｋ）と、全チャンネルについてのスペクトログラム１枚分の同時エントロピーＨ（Ｙ）の関係を示す図である。
【図３】ＳＴＦＴの窓長と時間周波数領域ＩＣＡの分離精度の関係をプロットしたグラフを示す図である。
【図４】横軸の窓長を実際の秒数として、時間周波数領域ＩＣＡの分離精度を表わしたグラフを示す図である。
【図５】複数の信号が混合された音信号の分離処理に適用する音情報の取得構成例について説明する図である。
【図６】時間領域の畳み混合が時間周波数領域において瞬時混合でなく畳み込み混合であると考える概念について説明する図である。
【図７】短時間フーリエ変換（ＳＴＦＴ）について説明する図である。
【図８】スペクトログラムＸから、時間方向に再び短時間フーリエ変換（ＳＴＦＴ）したＸ'（モジュレーション・スペクトログラム）への変換について説明する図である。
【図９】スペクトログラムＸから、時間方向に再び短時間フーリエ変換（ＳＴＦＴ）したＸ'（モジュレーション・スペクトログラム）への変換について説明する図である。
【図１０】エントロピーＨ（Ｙ'ｋ）の計算方法について説明する図である。
【図１１】観測信号スペクトログラムに対して、フレーム番号をずらしながら縦に積み重ねたベクトルを生成する処理について説明する図である。
【図１２】観測信号スペクトログラムＸについて、ｔ−ｌ（エル）番目からｔ−ｌ（エル）＋Ｌ'番目のフレームを畳み込むことで分離結果を生成する操作について説明する図である。
【図１３】シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理について説明する図である。
【図１４】シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理のシーケンスについて説明するフローチャートを示す図である。
【図１５】本発明の信号分離装置の構成例をについて説明する図である。
【図１６】本発明の信号分離装置の構成例をについて説明する図である。
【図１７】本発明の信号分離装置の処理シーケンスを説明するフローチャートを示す図である。
【図１８】本発明の信号分離装置の実行する分離処理の詳細シーケンスについて説明するフローチャートを示す図である。
【図１９】本発明の信号分離装置の実行する分離処理の詳細シーケンスについて説明するフローチャートを示す図である。
【図２０】周波数ごとにフレームタップ数［Ｌ'］の値を異ならせる処理について説明する図である。
【図２１】２段階分離によるチャンネル数削減を行なって信号分離を行なう処理についてのフローチャートを示す図である。
【図２２】残響除去処理について説明するフローチャートを示す図である。
【図２３】本発明の信号分離装置の効果の確認のための実験装置の構成について説明する図である。
【図２４】本発明の信号分離装置の効果の確認のための実験結果の評価データを示す図である。
【図２５】本発明の信号分離装置の効果の確認のための実験結果の評価データを示す図である。
【図２６】本発明の信号分離装置の効果の確認のための実験結果の評価データを示す図である。
【図２７】本発明の信号分離装置の効果の確認のための実験結果の評価データを示す図である。
【図２８】信号分離処理の評価実験を行なった環境について説明する図である。
【図２９】信号分離処理の評価実験に適用した音源について説明する図である。
【図３０】信号分離処理の評価実験に適用した音源の入出力パターンについて説明する図である。
【図３１】信号分離処理の評価実験における観測信号の例について説明する図である。
【図３２】信号分離処理の評価実験において得られたシフト＆積み重ね（図１１参照）を行なった結果について説明する図である。
【図３３】信号分離処理の評価実験において得られた分離結果およびＳＩＲについて説明する図である。
【図３４】信号分離処理の評価実験において得られた評価結果について説明する図である。
【図３５】信号分離処理の評価実験において得られた評価結果について説明する図である。
【発明を実施するための形態】
【００５３】
以下、図面を参照しながら本発明の信号分離装置、および信号分離方法、並びにコンピュータ・プログラムの詳細について説明する。本発明は、前述したように複数の原信号が混合されて取得された混合信号の信号解析によって原信号を分離・復元する処理を実行する信号分離処理を行なうものであり、独立成分分析（ＩＣＡ：Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）による信号分離処理を行なう構成である。
【００５４】
具体的には、図５に示すように、Ｎ個の音源１１１−１〜１１１−Ｎから異なる音が鳴っていて、それらをｎ個のマイク１２１−１〜１２１−ｎで観測するような状況において、マイク１２１−１〜１２１−ｎによって取得された混合信号に基づいて独立成分分析（ＩＣＡ）による信号分離処理を行なう。
【００５５】
先に説明したように、ある１つのマイクｊ（ただし１≦ｊ≦ｎ）で観測される信号（観測信号）は前述した式［１．１］のように、原信号と伝達関数との畳み込み演算を全音源について総和した式として表わせる（「畳み込み混合」）。さらに、全てのマイク１〜ｎについての観測信号を一つの式で表わすと、前述した式［１．２］のように表わせ、この畳み混合を解く手法として２つの方法、すなわち、
（１）時間領域で畳み込み混合を直接解く。（時間領域逆畳み込み）
（２）観測信号を時間周波数領域に変換し、瞬時混合問題として解く。
これらの手法があったが、
（２）観測信号を時間周波数領域に変換し瞬時混合問題として解くという処理を行なう前提として、従来の時間周波数領域ＩＣＡの枠組みでは、時間領域の畳み混合が時間周波数領域では瞬時混合で表わされると考えていた。それに対し本発明では、時間周波数領域でも依然として畳み込み混合であると考える。この概念について図６を参照して説明する。
【００５６】
図６に示す図６（ａ）は、原信号、すなわち、図５に示す各音源１１１−１〜１１１−Ｎの出力する原信号のスペクトログラムを縦に積み重ねたものである。それぞれの音源のスペクトログラムをＳ１，Ｓ２、両者を縦に積み上げたものをＳとする。なお、スペクトログラムは、前述したように、横軸をｔ（フレーム番号）、縦軸をω（周波数ビン番号）としてＸｋ（ω，ｔ）の絶対値である｜Ｘｋ（ω，ｔ）｜を色の濃淡で表現した図である。
【００５７】
図６（ａ）に示す原信号のスペクトログラムにおいて、ｔ番目のフレームの信号をベクトルで表現したものをＳ（ｔ）とおく。なお、スペクトログラムの１フレーム分をスペクトルと呼ぶ。
【００５８】
従来は、Ｓ（ｔ）がフレーム遅延なしでマイクに届くと考えていたが、本発明ではフレーム遅延があると考える。すなわち、図５を参照して説明すると、それぞれの音源１１１−１〜１１１−Ｎでスペクトルという名のベクトルが独立に発生し、それらが０以上の遅延を伴ってセンサーとしてのマイク１２１−１〜１２１−ｎに届く。これらには直接波と反射波が含まれる。
【００５９】
異なる音源からの直接波、また直接波と反射波、さらには単純な反射と複雑な反射など、様々な信号がマイクによって取得されることになり、その信号には様々な遅延量が存在すると推定される。ここで遅延の最大値をＬ＋１とすると、図６（ａ）に示す原信号のスペクトログラムにおけるｔ番目のフレーム信号のベクトル表現であるスペクトルＳ（ｔ）の影響は観測信号のｔ番目からｔ＋Ｌ番目のフレームに及ぶことになる。
【００６０】
図６（ｂ）は、観測信号のスペクトログラムであり、各マイク１２１−１〜１２１−ｎによって取得された観測信号について、短時間フーリエ変換（ＳＴＦＴ）を実行して生成した観測信号のスペクトログラムＸである。
【００６１】
短時間フーリエ変換（ＳＴＦＴ）について、図７を用いて説明する。例えば図５に示すような環境においてｋ番目のマイクによって収録された観測信号ｘｋを図７（ａ）に示す。この観測信号ｘｋから一定長を切り出した切り出しデータであるフレーム１７１〜１７３にハニング窓やサイン窓等の窓関数を作用させる。なお、切り出した単位をフレームと呼ぶ。切り出す長さ（サンプルポイント数）は、従来法の時間周波数領域ICAにおいて最も高精度の分離結果が得られる長さ（図３によれば、512ポイントまたは1024ポイント付近）と同じ値でよい。
１フレーム分のデータに対して、離散フーリエ変換（有限区間のフーリエ変換のこと。略称DFT）または高速フーリエ変換（FFT）を施すことにより、周波数領域のデータであるスペクトルＸｋ（ｔ）を得る（ｔはフレーム番号）。
【００６２】
切り出すフレームの間には、図に示すフレーム１７１〜１７３のように重複があってもよく、そうすることで連続するフレームのスペクトルＸｋ（ｔ−１）〜Ｘｋ（ｔ＋１）を滑らかに変化させることができる。また、スペクトルをフレーム番号に従って並べたものをスペクトログラムと呼ぶ。図７（ｂ）がスペクトログラムの例である。
【００６３】
なお、短時間フーリエ変換（ＳＴＦＴ）において切り出すフレーム間に重複がある場合は、逆フーリエ変換（ＦＴ）においてもフレームごとの逆変換結果（波形）波形を重複つきで重ね合わせる。これをオーバラップ加算（ｏｖｅｒｌａｐ  ａｄｄ）という。逆変換結果は、オーバラップ加算（ｏｖｅｒｌａｐ  ａｄｄ）の前にサイン窓等の窓関数を再び作用させても良く、これを、ｗｅｉｇｈｔｅｄ  ｏｖｅｒｌａｐ  ａｄｄ（ＷＯＬＡ）という。ＷＯＬＡにより、フレーム間の不連続性に由来するノイズを低減することができる。
【００６４】
図６（ｂ）は、図７を参照して処理によって得られる観測信号のスペクトログラムであり、スペクトログラムを縦に積み重ねたものである。それぞれのセンサー（マイク）のスペクトログラムをＸ１，Ｘ２、両者を縦に積み上げたものをＸとする。観測信号のスペクトログラムをＸとすると、Ｘ（ｔ）からＸ（ｔ＋Ｌ）のＬ＋１個のフレームは原信号スペクトルＳ（ｔ）の影響を受けている。逆にいうと、図６（ｂ）の観測信号中のｔ番目のフレームの観測信号Ｘ（ｔ）は、それより前のＬ＋１フレーム分の原信号の影響を受けている。
【００６５】
このように、観測信号中のｔ番目のフレームの観測信号Ｘ（ｔ）が、それより前のＬ＋１フレーム分の原信号の影響を受けていることを考慮すると、観測信号Ｘ（ｔ）は、以下に示す式［６．１］のような畳み込み混合で表わすことができる。
【数６】
2011215649.tif 000008
【００６６】
上記式［６．１］は、先に説明した式［１．２］と似ているが、式［６．１］は時間周波数領域の式であることに注意されたい。Ｌ＝０の場合は、従来の瞬時混合と等価となる。すなわち、観測信号Ｘ（ｔ）が原信号スペクトルＳ（ｔ）のみによって影響されている場合はＬ＝０となり、従来の瞬時混合と等価となる。
両者の畳み込みを区別するため、
式［１．２］のＬを「時間タップ数」、
式［６．１］のＬを「フレームタップ数」、
と定義する。
【００６７】
なお、上記式［６．１］は、ＳＴＦＴにおいてフレームのシフト幅を１とした場合には厳密に成立する。また、フレームのシフト幅を２以上に設定した場合でも、近似的に成立する。この点についての詳細は、発明者自身の論文である以下の文献を参照されたい。
［Ｈｉｒｏｅ，Ａ．Ｂｌｉｎｄ  Ｖｅｃｔｏｒ  Ｄｅｃｏｎｖｏｌｕｔｉｏｎ：Ｃｏｎｖｏｌｕｔｉｖｅ  Ｍｉｘｔｕｒｅ  Ｍｏｄｅｌｓ  ｉｎ  Ｓｈｏｒｔ−Ｔｉｍｅ  Ｆｏｕｒｉｅｒ  Ｔｒａｎｓｆｏｒｍ  Ｄｏｍａｉｎ．Ｉｎ  Ｍ．Ｅ．  Ｄａｖｉｅｓ  ｅｔ  ａｌ．（Ｅｄｓ．）：ＩＣＡ  ２００７，ＬＮＣＳ  ４６６６，ｐｐ．４７１-４７９，２００７．］
【００６８】
短時間フーリエ変換（ＳＴＦＴ）の窓長よりも残響の時間の方が長い場合、残響の影響は１フレームでは完結せず、複数のフレームに及ぶ。複数のフレームに跨がる残響は、時間周波数領域での畳み込みとして表現できるため、本発明で導入する「時間周波数領域での畳み込み混合」という考え方によって、ＳＴＦＴの窓長を越える残響も除去することが可能である。
【００６９】
すなわち、先に説明した図３のＳＴＦＴの窓長と時間周波数領域ＩＣＡの分離精度の関係をプロットしたグラフを例にして説明すると、長い窓（２０４８や４０９６など）の代わりに、短めの窓（５１２や１０２４）と複数のフレームタップ（１６や３２）との組み合わせが可能となり、長い窓のトレードオフを回避しつつ長い窓と同等のタイムスパン（時間タップ数とフレーム方向へのシフト幅とフレームタップ数とから算出される時間）を確保することが可能となる。
【００７０】
また、時間領域逆畳み込みと比べ、ずっと少ないタップ数の畳み込みで済む（数十のオーダー）ため、時間領域逆畳み込みの課題も回避できる。なお、以降では、
原信号から観測信号が生成される際のフレームタップ数をＬという文字で表わし、
一方、
観測信号から分離結果を生成する際のフレームタップ数をＬ'と表わす。
Ｌはその環境の残響時間およびＳＴＦＴの窓長とシフト幅から決まる値である。一方、Ｌ'はＬとは異なる値に設定することができる。（Ｌ'＝０とすると、従来法と等価になる。）
【００７１】
観測信号のフレームタップ数Ｌは、以下の式で計算することができる。
Ｌ＝Ｔｒ×Ｆｓ／Ｓ
ただし、
Ｔｒ：環境の残響時間
Ｆｓ：サンプリング周波数
Ｓ：ＳＴＦＴのシフト幅
である。
【００７２】
例えば、
残響時間Ｔｒ＝０．３秒、
サンプリング周波数Ｆｓ＝１６０００Ｈｚ、
シフト幅Ｓ＝２５６、
とすると、
原信号から観測信号が生成される際のフレームタップ数Ｌは、
Ｌ＝１８．７５である。
すなわち、残響の影響は１９フレームに及ぶ（端数切り上げ）ことが分かる。
【００７３】
観測信号Ｘから分離結果Ｙ、すなわち、図６（ｂ）の観測信号Ｘから図６（ｃ）の分離結果Ｙを生成するためのフレームタップ数Ｌ'については、Ｌが既知であれば（すなわち、残響時間が既知であれば）、
Ｌ'＝αＬ
とすればよい（αは適切な正の実数）。
Ｌが未知である場合、Ｌ'は例えば以下のいずれかの方法で決定することができる。
【００７４】
第一の方法は、Ｌ'＝６４やＬ'＝１００といった一定の値に決め打つことである。基本的に、Ｌ'が大きくなるほど計算量も増えるため、計算量と分離性能との兼ね合いからＬ'を決定しても良い。
第二の方法は、何らかの方法で残響時間を測定し、その残響時間から上記の式で求めたＬの値の定数倍をＬ'、すなわちＬ'＝αＬとする方法である。残響時間の測定方法としては、例えば装置自体に装備されたスピーカーからインパルス性の音を発し、その音が十分減衰するまでの時間を計測する。
第三の方法は、既知の原信号から生成された観測信号に対してさまざまなＬ'の下で分離を行ない、最も良い分離結果をもたらすＬ'の値を採用することである。そのためには、例えば装置の周辺にスピーカーを複数設置し、それぞれから既知の音を鳴らし、それらの音を複数のマイクで観測する。その観測結果に対して、異なるＬ'（例えばＬ'＝０〜１００のそれぞれの値）で分離結果を生成する。分離結果と原信号とから、後述のＳＩＲ（ｓｉｇｎａｌ−ｉｎｔｅｒｆｅｒｅｎｃｅ  ｒａｔｉｏ）という分離性能尺度を計算し、最高のＳＩＲをもたらすＬ'を採用する。環境が同じであれば、原信号が未知の場合でも、そのＬ'が最高の分離信号をもたらす可能性が高い。
【００７５】
例えば上記いずれかの方法によって、Ｌ'、すなわち、観測信号Ｘから分離結果Ｙを生成するためのフレームタップ数Ｌ'、具体的には、例えば、図６（ｂ）の観測信号Ｘから図６（ｃ）の分離結果Ｙを生成するためのフレームタップ数Ｌ'を決定し、このフレームタップ数Ｌ'を用いて、観測信号の複数の連続フレームから分離結果を生成する。
【００７６】
時間周波数領域において畳み込み混合された観測信号を分離する処理方式としては、たとえば以下の手法のいずれかを適用することが可能である。
（１）時間周波数領域において、畳み込み混合を直接解く。
（２）スペクトログラムを時間方向へもう一度短時間フーリエ変換（ＳＴＦＴ）し、瞬時混合として解く。
（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理によって解く。
以下では、それぞれの方法について説明する。
【００７７】
なお、上記の「（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理」は、「（１）時間周波数領域において、畳み込み混合を直接解く」と同等の分離処理を実現する方式であり、観測信号スペクトログラムをシフトしながら積み重ねた後、その結果に対して従来の時間周波数領域の瞬時混合ＩＣＡを適用する手法である。詳細については後段で説明する。
【００７８】
（１）時間周波数領域において、畳み込み混合を直接解く
まず、時間周波数領域において、畳み込み混合を直接解くことで、時間周波数領域において畳み込み混合された観測信号を分離する処理について説明する。
【００７９】
再び図６を参照して説明する。前述の通り、原信号スペクトログラムのｔ番目のフレームＳ（ｔ）は、観測信号のｔ番目からｔ＋Ｌ番目のフレームに影響を与える。従って、原信号の１フレーム分を推定するためには、観測信号がＬフレーム分かそれ以上必要である。その値をＬ'とする。
【００８０】
分離信号中のｔ番目のフレームを基準にした場合、例えば、図６（ｃ）の分離信号中のＹ（ｔ）を基準として考えると、Ｓ（ｔ）を推定するためには少なくとも以降のＬ＋１フレーム分のデータが必要である。そこで、原信号の推定結果（＝分離結果）であるＹ（ｔ）を、先に示した式［６．３］のように、観測信号Ｘ（ｔ）からＸ（ｔ＋Ｌ'）までの畳み込み混合として表わす。
【００８１】
一方、分離信号中のｔ＋Ｌ'番目のフレームを基準にした場合、例えば、図６（ｃ）の分離信号中のＹ（ｔ＋Ｌ'）を基準として考えると、Ｓ（ｔ）を推定するためには直前のＬ＋１フレーム分のデータが必要である。そこで、分離信号Ｙ（ｔ）は式［６．２］のようにＸ（ｔ−Ｌ'）からＸ（ｔ）までの畳み込み混合として表わす。
両者の式は、Ｓ（ｔ）とのフレームのずれが異なるが、本質的には等価であるため、以下では式［６．２］からＹ（ｔ）を推定する方法について説明する。
【００８２】
混合は同じ周波数ビンでのみ起こると仮定する（すなわち、伝播の途中で周波数の変調が起こることはないと仮定する）と、全周波数ビンの混合の式である式［６．１］は、周波数ビンごとの混合の式である式［６．４］のように書き表すことができる。その仮定の下では、式［６．２］の分離行列Ｗ［ｌ］は、式［６．５］のように、対角行列から構成される行列として表わせるため、Ｗ［ｌ］を推定するためには式［６．５］の非零の成分のみを推定すればよい。
【００８３】
式［６．２］から学習規則（ΔＷの式）を求める処理は以下のように行なう。スペクトログラム全体での独立性を表わす尺度として、式［４．５］で計算されるＫｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量Ｉ（Ｙ）を考える。なお、この手法は、特開２００６−２３８４０９に記載したと同様の処理である。
【００８４】
Ｙ（ｔ）の成分であるＹ１（ｔ）〜Ｙｎ（ｔ）を互いに独立にするためには、式［４．５］のＫｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量Ｉ（Ｙ）を最小にする分離行列Ｗ［０］〜Ｗ［Ｌ'］を求めれば良いわけである。なお、特開２００６−２３８４０９において示した方法は瞬時混合だったため、分離行列は１つだけ推定すればよかったが、本発明では、Ｌ'＋１個のフレームの畳み込み混合であるため、分離行列もＬ'＋１個推定する必要がある。
【００８５】
ここで「Ｙ１（ｔ）〜Ｙｎ（ｔ）は互いに独立」という仮定（チャンネル間の独立性）の他に、「Ｙｋ（ｔ−Ｌ'）〜Ｙｋ（ｔ）も互いに独立」という仮定（フレーム間の独立性）を設けると、最終的に、以下に示す式［７．１］の学習規則が導出される。
【数７】
2011215649.tif 000009
【００８６】
すなわち、分離行列Ｗ［０］〜Ｗ［Ｌ'］を求めるためには、式［６．２］，［７，１］，［７．８］をＷ［０］〜Ｗ［Ｌ'］が収束するまで（または一定回数）繰り返す。ただし、式［７，１］のΔＷ［ｌ］（ω），Ｗ［ｌ］（ω）は、それぞれΔＷ［ｌ］とＷ［ｌ］から周波数ビンωに対応する要素を抽出した部分行列（式［６．６］）であり、Ｒω［ｌ］は式［７．２］で計算されるクロス項である。式［７．２］のφω（Ｙ（ｔ））は、スコア関数からなるベクトルであり（式［７．４］）、これは本出願人の先の出願（特開２００６−２３８４０９）において示したスコア関数からなるベクトルと同一である。スコア関数は確率密度関数の対数微分として定義され（式［７．５］）、特開２００６−２３８４０９において開示したと同様に多変量のスコア関数を使用することで、パーミュテーションの発生を阻止できる。
【００８７】
スコア関数の具体例は特開２００６−２３８４０９において説明したと同一でよく、例えば式［７．６］を用いる。この式において、αｋ（ω），ｍ，γｋ（ω）は正の実数、βｋ（ω）は非負の実数である。簡単な例として、式［７．７］を適用してもよい。
【００８８】
式［７．８］において、ηは学習率と呼ばれる正の実数である。ηは例えば０．１といった定数でも良いが、式［７．９］のように適用的に算出しても良い。ただし、この式において‖Ｗ（ω）‖はＷ［０］（ω）〜Ｗ［Ｌ］（ω）の全要素の２乗和（式［７．１０］）、‖ΔＷ（ω）‖も同様にΔＷ［０］（ω）〜ΔＷ［Ｌ'］（ω）の全要素の２乗和、η０はηの上限値を表わす正の実数である。式［７．８］を用いると、学習の始めはηが比較的小さな値となり（‖ΔＷ（ω）‖が大きいため）、Ｗ（ω）がオーバーフローするのを回避できる。一方で学習の終わりではηが比較的大きな値となり（‖ΔＷ（ω）‖がゼロ行列に近いため）、Ｗ（ω）が目的の値に早く収束する。
【００８９】
なお、式［６．２］の代わりに式［６．３］を用いる場合は、学習において式［６．３］，［７，１］，［７．８］を繰り返す。ただし、式［７．１］のＲω［ｌ］として、式［７．２］の代わりに式［７．３］を用いる。
【００９０】
上記で式［７．１］を導出する際に、「Ｙｋ（ｔ−Ｌ'）〜Ｙｋ（ｔ）も互いに独立」という仮定を置いていたが、代わりに「Ｙｋ（ｔ−Ｌ'）〜Ｙｋ（ｔ）は互いに依存」という仮定を置くと、別の学習規則である以下に示す式［８．１］が得られる（式［７．１］は共通）。
【数８】
2011215649.tif 000010
【００９１】
式［７．２］と式［８．１］との違いはスコア関数の引数にあり、式［７．２］はＹ（ｔ）のみを引数としているのに対して式［８．１］はＹ（ｔ−Ｌ'）〜Ｙ（ｔ）を全て引数としている。このスコア関数は式［８．４］  で定義され、この式に現われるＰ（Ｙｋ（ｔ），...，Ｙｋ（ｔ−Ｌ'））は、隣接するＬ'＋１個のフレームのデータが同時に発生する確率を表わしている。そのため、式［８．１］を用いると、隣り合ったフレーム間での依存関係が分離行列に一層反映され得る。スコア関数の例としては、式［８．５］（式［８．６］はさらにその具体例）が挙げられる。
【００９２】
なお、式［８．１］は式［６．２］に対応した式である。式［６．２］の代わりに式［６．３］を用いる場合は、式［８．２］が対応する。
【００９３】
なお、上記の説明では、独立性の尺度としてKullback-Leiblar情報量を採用していたが、他の尺度を用いても良い。Kullback-Leiblar情報量以外で独立性を表す尺度としては、非正規性や尖度（kurtosis）などがあり、それらの量を最大または最小とするように分離行列を更新してもよい。
【００９４】
（２）スペクトログラムを時間方向へもう一度ＳＴＦＴし、瞬時混合として解く
次に、スペクトログラムを時間方向へもう一度短時間フーリエ変換（ＳＴＦＴ）し、瞬時混合として解くことで、時間周波数領域において畳み込み混合された観測信号を分離する処理について説明する。
【００９５】
畳み込みをタップ数よりも長い窓長で短時間フーリエ変換（ＳＴＦＴ）すると、畳み込みはただの積に変換される。これは、時間周波数領域の畳み込み混合についても同様である。すなわち、時間周波数領域の畳み込み混合である上述した式［６．４］を時間方向に再び短時間フーリエ変換（ＳＴＦＴ）すると、以下に示す式［９．１］が得られる。ただし、Ｘ'，Ａ'，Ｓ'は式［６．４］のＸ，Ａ，Ｓの各要素を短時間フーリエ変換（ＳＴＦＴ）した結果である。
【数９】
2011215649.tif 000011
【００９６】
上記式［９．１］は瞬時混合の式であり、観測信号を独立な成分へ分離するためには式［９．２］を考えればよい。
【００９７】
ここで、図８、図９を用いて、スペクトログラムＸから、時間方向に再び短時間フーリエ変換（ＳＴＦＴ）したＸ'（モジュレーション・スペクトログラム）への変換について説明する。比較のため、波形ｘからスペクトログラムＸへの変換についても説明する。
【００９８】
図８（ａ）は、観測信号の波形である（この図ではチャンネル数＝２としてあるが、
チャンネル数は任意である）。
図８（ｂ）は、観測信号の波形（図８（ａ））を短時間フーリエ変換（ＳＴＦＴ）することで生成されたスペクトログラムである（チャンネルごとにＳＴＦＴを行ない、それぞれの結果を縦に並べて表示してある）。窓長＝Ｎでフーリエ変換するとＮ個の周波数成分が得られるが、実数データの変換においては負の周波数成分は正の周波数成分の共役複素数の関係にある（共役対称と呼ぶ）ため、直流成分と正の周波数成分とのＮ／２＋１＝Ｍ本の周波数ビンだけ考慮すればよい。図８（ｂ）に示す周波数ビン２０１は、周波数ビンの１本を示している。なお、通常、スペクトログラムはＸの絶対値をプロットしたものを指すが、ここではＸ自体もスペクトログラムと呼ぶ。原信号Ｓや分離結果Ｙについても同様である。
【００９９】
ここでさらに、図８（ｂ）に示すスペクトログラムＸについて、周波数ビンごとに短時間フーリエ変換（ＳＴＦＴ）を行なう。スペクトログラムをもう一度ＳＴＦＴして生成されたデータをモジュレーション・スペクトログラムと呼ぶ。２度目の短時間フーリエ変換（ＳＴＦＴ）の窓長をＬ'とすると、１本の周波数ビン、例えば図８（ｂ）のビン２０１からＬ'本のビンが生成されるため、それを奥行き方向で表現する。これが、図９（ｃ）に示すビン２０２であり、これらを集積した結果が、図９（ｃ）に示すデータとなる。
【０１００】
すなわち、図９（ｃ）は、図８（ｂ）に示すスペクトログラムＸについて、周波数ビンごとに短時間フーリエ変換（ＳＴＦＴ）を行なって生成したモジュレーション・スペクトログラムであり、図９（ｃ）に示すような直方体の構造のモジュレーション・スペクトログラムＸ'で表わされる。奥行き方向も周波数成分であるが、波形の周波数成分ではなくてエンベロープの周波数成分である。２度目の短時間フーリエ変換（ＳＴＦＴ）では変換前データも複素数であるため、変換結果は共役対称にはならない。従って、Ｌ'本のビンは全て考慮しなければならない。
【０１０１】
新たに生成されたビンを、奥行き方向の代わりに縦方向に配置する。すなわち、図９（ｃ）に示すビン２０２を図９（ｄ）に示すビン２０３のように配置すると、モジュレーション・スペクトログラムは図９（ｄ）に示すように平面でも表現可能である。図９（ｄ）に示すモジュレーション・スペクトログラムＸ'は、図８（ｂ）に示すスペクトログラムＸと一見似ているが、周波数ビンの意味合いが異なることに注意されたい。（１チャンネル辺りのビンの本数は、図８（ｂ）に示すスペクトログラムＸがＬ'本、図９（ｄ）に示すモジュレーション・スペクトログラムＸ'がＭ×Ｌ'本である。
【０１０２】
ここで先に示した数式［９．ｎ］に戻ると、式［９．１］，［９．２］のＸ'に相当するのが、図９（ｃ）に示す立体版モジュレーション・スペクトログラムＸ'であり、ωが縦方向の周波数ビンを、ω２が奥行き方向のビンを表わしている。式［９．２］において、（ω，ω２）というペアをまとめてω'というインデックスで表わすと、式［９．３］が得られる。式［９．３］は、図９（ｄ）に示す平面版モジュレーション・スペクトログラムＸ'に対応している。
【０１０３】
式［９．２］または式［９．３］から学習規則（ΔＷの式）を求めるためには、以下のように考える。モジュレーション・スペクトログラム全体での独立性を表わす尺度として、式［９．５］で計算されるＫｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量を考える。この式は式［４．５］とほぼ同一だが、Ｈ（Ｙｋ'）は１チャンネル分のモジュレーション・スペクトログラムから算出されるエントロピー、Ｈ（Ｙ'）はモジュレーション・スペクトログラム全体から算出される同時エントロピーである。Ｈ（Ｙ'ｋ）の計算方法について、図１０を参照して説明する。
【０１０４】
図１０は、図９（ｃ）に示す立体版モジュレーション・スペクトログラムＸ'に相当する。すなわち、例えば、観測信号の波形（図８（ａ））を短時間フーリエ変換（ＳＴＦＴ）することで生成された図８（ｂ）に示すスペクトログラムＸに対して、さらに周波数ビンごとに短時間フーリエ変換（ＳＴＦＴ）を行なって生成したモジュレーション・スペクトログラムに相当する。この図１０に示す立体版モジュレーション・スペクトログラムＸ'において、例えば、１番目のチャンネルのエントロピー計算においては、図１０における１番目のフレームのモジュレーション・スペクトログラムＹ１'（ｔ）２２１は、平面を表わし、それを引数とする多変量確率密度関数Ｐ（Ｙ１'（ｔ））２２２にＹ１'（ｔ）を代入することで、エントロピーＨ（Ｙ１'）２２３を求める。
【０１０５】
式［９．３］は、変数名の違いを除くと式［３．５］と同一である。従って、学習規則を導出するためには、式［５．２］の変数名を付け替えればよく、結果として式［９．５］を得る。すなわち、式［９．３］，［９．５］，［９．６］をＷ'が収束するまで繰り返せば、Ｙ１'（ｔ）〜Ｙｎ'（ｔ）が互いに独立になる。
【０１０６】
互いに独立になったモジュレーション・スペクトログラムＹ１'〜Ｙｎ'それぞれに対して逆フーリエ変換とｏｖｅｒｌａｐ  ａｄｄを作用させると、互いに独立なスペクトログラムＹ１〜Ｙｎが得られる。
【０１０７】
上記の説明では、独立性の尺度としてKullback-Leiblar情報量を採用していたが、（１）と同様に他の尺度を用いても良い。
また、分離行列更新の式として自然勾配法に基づく式を導出したが、代わりに他のアルゴリズムを用いてもよい。他のアルゴリズムとしては、正規直交制約つき勾配法や不動点法やニュートン法などがあるが、この点は従来の瞬時混合ＩＣＡと同様である。
【０１０８】
（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理によって解く。
次に、シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理によって、時間周波数領域において畳み込み混合された観測信号を分離する処理について説明する。
【０１０９】
この３番目の処理方式は、先に説明した［（１）時間周波数領域において、畳み込み混合を直接解く方式］とほぼ同等の分離処理を実現する方式であり、本出願人の先の特許出願である特開２００６−２３８４０９において開示した瞬時混合ＩＣＡ処理を利用して実現される。
【０１１０】
本方式は、例えば観測信号スペクトログラムをシフトしながら積み重ねた後、その結果に対して時間周波数領域の瞬時混合ＩＣＡ、すなわち、本出願人の先の特許出願である特開２００６−２３８４０９において開示した瞬時混合ＩＣＡを適用することで実現される。この３番目の手法の適用により、パーミュテーション（置換）問題が解決され、かつ直接波、反射波など様々な遅延量を持つ混合された音信号について、遅延量を考慮した高精度な分離処理が実現される。
【０１１１】
まず、この３番目の手法について説明する前に、観測信号の分離処理において発生するパーミュテーション（置換）問題と、この問題を解決した本出願人の先の特許出願である特開２００６−２３８４０９の瞬時混合ＩＣＡの概要について、再度、簡潔に説明する。
【０１１２】
ｎ個の音源が発するお互いに独立な原信号をｓ１〜ｓｎとし、それらを要素とするベクトルをｓとしたとき、マイクロホンで観測される観測信号ｘは、原信号ｓに前述した式［１．２］の畳み込み・混合演算を施したものとなる。次に、観測信号ｘに対して短時間フーリエ変換を施し、時間周波数領域の信号Ｘを得る。Ｘの要素をＸｋ（ω，ｔ）とすると、Ｘｋ（ω，ｔ）は複素数値をとる。Ｘｋ（ω，ｔ）の絶対値である｜Ｘｋ（ω，ｔ）｜を色の濃淡で表現した図が例えば、図６（ｂ）に示す観測信号のスペクトログラムである。このスペクトログラムは、例えば図５に示すマイク１２１−１〜１２１−ｎによって取得された観測信号について、短時間フーリエ変換（ＳＴＦＴ）を実行して生成した観測信号のスペクトログラムＸである。
【０１１３】
スペクトログラムは、例えば、横軸をｔ（フレーム番号）、縦軸をω（周波数ビン番号）としてＸｋ（ω，ｔ）の絶対値である｜Ｘｋ（ω，ｔ）｜を色の濃淡で表現した図である。続いて、信号Ｘの各周波数ビンに、分離行列Ｗ（ω）を乗算することで分離信号Ｙを得る。そして、分離信号Ｙを逆フーリエ変換することで時間領域の分離信号ｙを得ることができる。
【０１１４】
しかし、前述したように、従来の時間周波数領域の独立成分分析では、信号の分離処理を周波数ビン毎に行っており、周波数ビンの間の関係は考慮していない。そのため、分離自体は成功しても、周波数ビンの間でスケーリング及び分離先の不統一が発生する可能性がある。スケーリングの不統一については、音源毎に観測信号を推定する方法により解決できるが、分離先の不統一、例えばω＝１ではＹ１にＳ１由来の信号が現れるのに対してω＝２ではＹ１にＳ２由来の信号が現れるといったパーミュテーション（置換）問題は解決できない。
【０１１５】
本出願人の先の特許出願である特開２００６−２３８４０９は、このパーミュテーション（置換）問題を解決する手法を開示した。すなわち、スペクトログラム全体での分離を表わす式については、先に説明した以下に示す数式［４．４］を用い、スペクトログラム全体での独立性を最大にする分離行列Ｗを求めるという手法を採用した。
【数１０】
2011215649.tif 000012
【０１１６】
具体的には、スペクトログラム全体での独立性として、式［４．５］で表わされるＫＬ（Ｋｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ）情報量Ｉ（Ｙ）を導入し、Ｉ（Ｙ）を最小にする分離行列Ｗを求める。ＫＬ情報量Ｉ（Ｙ）は、スペクトログラムごとのエントロピーの総和から、スペクトログラム全体の同時エントロピーを引いたものであり、全てのスペクトログラムがお互いに独立となった場合に最小（理想的には０）となる。
【０１１７】
ＫＬ情報量Ｉ（Ｙ）の定義式［４．５］において、Ｈ（Ｙｋ）は各チャンネルについてのスペクトログラム１枚分のエントロピーを表し、Ｈ（Ｙ）は全チャンネルについてのスペクトログラム１枚分の同時エントロピーを表す。
【０１１８】
例えば、ｎ＝２のときのＨ（Ｙｋ）とＨ（Ｙ）との関係は、先に図２を参照して説明したとおりである。図２において、Ｐ（Ｙｋ（ｔ））は、Ｙｋ（ｔ）の確率密度関数であり、Ｈ（Ｙｋ）は各チャンネルについてのスペクトログラム１枚分のエントロピーである。Ｋｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量Ｉ（Ｙ）は、スペクトログラムごとのエントロピー１１，１２の総和から、スペクトログラム全体の同時エントロピー１３を引いたものであり、全てのスペクトログラムがお互いに独立となった場合に最小（理想的には０）となる。
【０１１９】
スペクトログラム全体でのＫＬ情報量Ｉ（Ｙ）を最小にするためには、先に説明したように以下に示す式［５．１］〜［５．３］をＷおよびＹが収束するまで繰り返す。
【数１１】
2011215649.tif 000013
【０１２０】
なお、式［５．３］に出てくるΔＷ（ω），Ｗ（ω），Ｙ（ω，ｔ）は、それぞれΔＷ，Ｗ，Ｙ（ｔ）からω番目の周波数ビンに対応する要素を抽出した部分行列である。こうすることで、パーミュテーション問題のない分離結果を得ることが可能となった。
【０１２１】
３番目の処理方式は、この特開２００６−２３８４０９において開示した時間周波数領域の瞬時混合ＩＣＡを適用する手法である。特開２００６−２３８４０９において示されている時間周波数領域の瞬時混合ＩＣＡを適用した処理は、具体的には、信号分離処理として、時間周波数領域の観測信号と初期値が代入された分離行列から時間周波数領域の分離信号を生成し、生成した時間周波数領域の分離信号と、多次元確率密度関数から導出される多次元スコア関数によって計算される分離行列とがほぼ収束するまで分離行列の修正を行い、この修正された分離行列を適用して時間周波数領域の分離信号を生成する処理として実行される。詳細については、特開２００６−２３８４０９に開示された通りである。
【０１２２】
以下説明する３番目の処理方式、すなわち、［（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理によって時間周波数領域において畳み込み混合された観測信号を分離する処理では、この特開２００６−２３８４０９において開示した時間周波数領域の瞬時混合ＩＣＡを適用する。具体的には、例えば観測信号スペクトログラムをシフトしながら積み重ねた後、その結果に対して時間周波数領域の瞬時混合ＩＣＡを適用する。以下、この３番目の手法について説明する。
【０１２３】
本手法では、音声入力部である複数のマイクの各々の観測信号スペクトログラムに対して、フレーム番号をずらしながら縦に積み重ねたベクトルを生成する。例えばｋ番目のマイクに相当するｋチャンネル目の観測信号スペクトログラム、すなわち、上記の式［４．１］の観測信号スペクトログラムＸｋ（ｔ）に対して、フレーム番号をずらしながら縦に積み重ねたベクトルを考える。さらに、それを全チャンネル分積み重ねたベクトルを考える。以下に示す式［１１．１］のベクトルＸ''（ｔ）である。式［１１．１］のベクトルＸ''（ｔ）は、ｎチャンネル分のベクトルを含む。なお、チャンネルごとのベクトルをＸｋ''（ｔ）と示している。
【０１２４】
【数１２】
2011215649.tif 000014
【０１２５】
上記数式［１１．１］に示すベクトルＸ''（ｔ）の作成手順について、図１１以下を参照して説明する。図１１は、各マイクの入力信号に基づいて生成されるチャンネルごとの観測信号スペクトログラムＸｋ（ｔ）から、チャンネルごとのベクトルＸ''（ｔ）を生成する処理について説明する図である。図１１（ａ）に示すデータ３０１、すなわちＸｋは観測信号１チャンネル分のスペクトログラム、すなわちｋ番目のマイクに相当するｋチャンネル目の観測信号スペクトログラムＸｋであり、先に説明した図６（ｂ）のＸ１やＸ２に相当する。
【０１２６】
このＸｋを、左にｌ（小文字のエル）フレームずつシフトした結果をＸｋ［ｌ］とする。図１１（ｂ）には、観測信号スペクトログラムＸｋを、シフト量ｌをｌ＝０〜Ｌ'まで、順次、変化させながら縦に複数積み重ねた構造を示している。
データ３１１−０が、シフト量ｌ＝０、
データ３１１−１が、シフト量ｌ＝ｌ（エル）フレーム
：
データ３１１−Ｌ'が、シフト量ｌ＝Ｌ'フレーム
である。
なお、Ｌ'は前述したように、観測信号から分離結果を生成する際のフレームタップ数である。
【０１２７】
１つの観測信号スペクトログラムから、これらの複数の異なるフレーム方向へのシフト量を持つ観測信号スペクトログラムシフトセットを生成し、これを観測信号スペクトログラムシフトセット［Ｘ''］とする。観測信号スペクトログラムシフトセット［Ｘ''］から１フレーム分を切り出すと、図１１（ｂ）に示す数式３１２となる。この数式は、上述した式［１１．１］に含まれる１つのチャンネル対応のベクトル［Ｘｋ（ｔ）］に対応する。
【０１２８】
式［１１．１］は、先に説明したように、１つのチャンネル対応の観測信号スペクトログラムＸｋ（ｔ）を、フレーム番号をずらしながら縦に積み重ねて生成される観測信号スペクトログラムシフトセットからなるベクトルを、さらに全チャンネル分積み重ねて生成される複数チャンネル対応観測信号スペクトログラムシフトセットからなるベクトルである。
【０１２９】
なお、図１１（ｂ）に示すように、シフトの際に生じた隙間、すなわち、図１１（ｂ）に示すデータ３１１−０〜３１１−Ｌ'の斜線部分には、ゼロに近い値を代入しておくか、両端の値（Ｘ（１）やＸ（Ｔ）など）をコピーして設定する。または、後述のゼロ除算対策が施されている場合には、ゼロを代入しても構わない。または、両端の隙間を取り除き、中間のＴ−Ｌ'フレーム分のデータを用いるようにしても構わない。さらには、通常のシフト処理ではなく、長さＴの巡回シフト（シフトではみ出した左端のデータを右端にコピーする）を適用する設定としてもよい。以下において説明する処理例は、巡回シフトによって生成した観測信号スペクトログラムシフトセット［Ｘ''］を適用した処理例について説明する。
【０１３０】
図１１（ｂ）に示すような、シフト処理と積み重ね処理で生成された観測信号スペクトログラムシフトセット［Ｘ''］を元の観測信号スペクトログラム［Ｘ］と比較すると、
観測信号スペクトログラム［Ｘ］はｎチャンネル分のスペクトログラム、
であるのに対して、
観測信号スペクトログラムシフトセット［Ｘ''］は見かけ上、ｎ×（Ｌ'＋１）チャンネル分のスペクトログラムを含む。ｎがマイク数に相当するチャンネル数、（Ｌ'＋１）が１つのチャンネルに対応して設定されるシフトデータ数である。
【０１３１】
この観測信号スペクトログラムシフトセット［Ｘ''］を、ｎ×（Ｌ'＋１）チャンネルの観測信号スペクトログラムとして、本出願人の先の特許出願である特開２００６−２３８４０９において開示した瞬時混合ＩＣＡを適用した方法で分離処理を行なう。この処理によって、先に説明した「（１）時間周波数領域において、畳み込み混合を直接解く」方式と同等な分離を行なうことができる。以降では、その原理について説明する。
【０１３２】
観測信号スペクトログラムＸについて、ｔ−ｌ（エル）番目からｔ−ｌ（エル）＋Ｌ'番目のフレームを畳み込むことで分離結果を生成する操作について考察する。すなわち、図１２に示すように、Ｘ（ｔ−ｌ（エル））からＸ（ｔ−ｌ（エル）＋Ｌ'）までのＬ'＋１フレームから１フレーム分の分離結果を生成する操作であり、この処理によって分離信号Ｙ［ｌ］（ｔ）を得るための算出式は式［１１．２］で表わされる。
【０１３３】
分離結果をＹ［ｌ］（ｔ）とおく。分離信号Ｙ［ｌ］（ｔ）はＬ'＋１フレームの間の畳み込みであるため、係数の行列はＬ'＋１個必要となるが、さらに、シフトフレーム数［ｌ（エル）］によっても異なる値をとるため、分離行列［Ｗ］については、２種類の添字をつけて、
Ｗ［ｌ，０］〜Ｗ［ｌ，Ｌ'］
と表わす。すなわち、分離行列［Ｗ］は、シフトフレーム数［ｌ（エル）］と、各々のシフトスペクトログラムに応じて設定する。
【０１３４】
式［１１．３］および式［１１．４］は、式［１１．２］に現れる部分行列の詳細であり、式［１１．５］は、式［１１．４］に現れる部分行列の詳細を示している。
分離信号［Ｙ［ｌ］（ｔ）］と分離行列［Ｗ［ｌ，τ］］は、それぞれ、各チャンネルの成分に対応したベクトルや行列からなる。なお、Ｗに対する添字τはτ＝０〜Ｌ'である、
【０１３５】
ここで、式［１１．６］で示される、
分離結果：Ｙ［０］（ｔ）〜Ｙ［Ｌ'］（ｔ）
をすべて包含したベクトル
分離結果ベクトルＹ''（ｔ）と、
式［１１．７］で示される複数の分離行列、
Ｗ［０，０］〜Ｗ［Ｌ'＋１，Ｌ'］
をすべて包含した行列
Ｗ''
これらのベクトル［Ｙ''（ｔ）］と、行列［Ｗ''］とを用いると、分離処理を示す式は、単純に式［１１．８］、すなわち、
Ｙ''（ｔ）＝Ｗ''Ｘ''（ｔ）・・・［１１．８］
このように示すことができる。
【０１３６】
従来法として先に説明した、特許文献１（特開２００６−２３８４０９）では、スペクトログラム全体での分離を表わす式として、先に説明した数式［４．４］、すなわち、
Ｙ（ｔ）＝ＷＸ（ｔ）・・・［４．４］
上記式を用いた処理としているが、式［１１．８］と式［４．４］と比較すると、式［１１．８］は単にチャンネル数がｎからｎ×（Ｌ'＋１）に増えたものとして、式［４．４］を適用したとみなすことができる。
【０１３７】
すなわち、図１３に示すように、複数チャンネル分の観測信号スペクトログラムシフトセット［Ｘ''］は、Ｘ１''〜Ｘｎ''によって構成され、これらをｎ×（Ｌ'＋１）の個別のチャンネルに対応する観測信号スペクトログラムとして考えれば、式［１１．８］は単にチャンネル数がｎからｎ×（Ｌ'＋１）に増えたものとして、式［４．４］を適用したとみなすことができる。
【０１３８】
従って、ｎチャンネル分の観測信号スペクトログラムＸを、図１１を参照して説明した方法によってｎ×（Ｌ'＋１）チャンネルに拡張し、その結果である観測信号スペクトログラムシフトセット［Ｘ''］に対して、特開２００６−２３８４０９の学習式である式［５．１］〜［５．３］を繰り返し適用すると、分離結果であるＹ''と分離行列Ｗ''とが得られるのである。
【０１３９】
ただし、式［５．１］〜［５．３］の変数の詳細である式［５．４］〜［５．７］において、ｎはｎ×（Ｌ'＋１）に読み替えること。また、ｋは１≦ｋ≦ｎではなく、１≦ｋ≦ｎ×（Ｌ'＋１）を表わすインデックスとなる。
【０１４０】
分離結果であるＹ''はｎ×（Ｌ'＋１）チャンネル分のスペクトログラムを含んでいるが、所望のものはｎチャンネル分（またはｎ未満）であるため、必要に応じてスペクトログラムの選択を行なう。選択の方法としては、分離結果Ｙ''の中からＹ１［０］，Ｙ２［０］，…，Ｙｎ［０］のように特定のシフト量［ｌ（エル）］に該当する成分のみを残すといった方法が適用可能である。
【０１４１】
または、観測信号から分離結果を生成する際のフレームタップ数をＬ'の値の決定方法と同様に、既知の信号を用いて最適なフレーム方向へのシフト量［ｌ（エル）］を求めても良い。すなわち、既知の信号をスピーカー等から鳴らして本発明の手法で集音および分離を行なった後、分離結果のＹｋ［０］〜Ｙｋ［Ｌ'］それぞれについて分離精度の尺度であるＳＩＲ（ｓｉｇｎａｌ−ｉｎｔｅｒｆｅｒｅｎｃｅ−ｒａｔｉｏ）を計算する。そして最も高い分離精度（ＳＩＲ）をもたらすシフト数：ｌ（エル）に対応した分離結果［Ｙｋ［ｌ］］を選択するといった処理が可能である。
【０１４２】
この（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理によって、時間周波数領域において畳み込み混合された観測信号を分離する処理のシーケンスを説明するフローチャートを図１４に示す。図１４のフローにおける各ステップの処理について説明する。
【０１４３】
まず、ステップＳ１１において、観測信号スペクトログラムをシフトしながら積み重ねる。この処理は、図１１を参照して説明した処理であり、各マイクによって取得される観測信号から生成される観測信号スペクトログラムに対して、シフトフレーム（ｌ（エル））単位で順次シフトとして、フレームタップ数をＬ'に相当するシフト量になるまでシフトデータを生成して積み重ねて、観測信号スペクトログラムシフトセット［Ｘ''］を生成する。
【０１４４】
次にステップＳ１２において、瞬時混合ＩＣＡ（または、変更したスコア関数）を用いて分離結果Ｙ''を求める。すなわち、観測信号スペクトログラムシフトセット［Ｘ''］に対して、特開２００６−２３８４０９の学習式である式［５．１］〜［５．３］を繰り返し適用して、分離結果であるＹ''と分離行列Ｗ''とを算出する。ただし、式［５．１］〜［５．３］の変数の詳細である式［５．４］〜［５．７］において、ｎはｎ×（Ｌ'＋１）に読み替えること。また、ｋは１≦ｋ≦ｎではなく、１≦ｋ≦ｎ×（Ｌ'＋１）を表わすインデックスとなる。
【０１４５】
なお、スコア関数は確率密度関数の対数微分として定義され、式［５．７］によって定義される。先の［（１）時間周波数領域において、畳み込み混合を直接解く］方式に式［７．５］について説明したように、特開２００６−２３８４０９において開示したと同様、多変量のスコア関数を使用することで、パーミュテーションの発生を阻止できる。このスコア関数を用いた処理については後述する。
【０１４６】
次に、ステップＳ１３において、必要に応じて分離結果Ｙ''の中から所望のスペクトログラムを選択する。すなわち、上述したように、分離結果であるＹ''はｎ×（Ｌ'＋１）チャンネル分のスペクトログラムを含んでいるが、所望のものはｎチャンネル分（またはｎ未満）であるため、必要に応じてスペクトログラムの選択を行なう。
【０１４７】
選択方法としては、分離結果Ｙ''の中からＹ１［０］，Ｙ２［０］，…，Ｙｎ［０］のように特定のシフト量［ｌ（エル）］に該当する成分のみを残すといった方法が適用可能である。この際、最も高い分離精度（ＳＩＲ）をもたらすシフト数：ｌ（エル）に対応した分離結果［Ｙｋ［ｌ］］を選択するといった処理が可能である。
【０１４８】
上記で説明した方法は、先に、［（１）時間周波数領域において、畳み込み混合を直接解く］方式において説明した式［７．２］〜式［７．５］を用いた方法とほぼ同等の処理を実行していることに相当する。すなわち、本処理例では、ｎ×（Ｌ'＋１）チャンネルの信号をお互いに独立となるように分離する処理であり、例えば、図１３を参照して説明すると、複数チャンネル分の観測信号スペクトログラムシフトセット［Ｘ''］を適用した結果として得られる分離結果である信号Ｙ''３３１の中のスペクトログラム１枚分である信号Ｙ１［０］３４１は、他音源に由来する信号Ｙｎ［０］３４３やＹｎ［Ｌ'］３４４と独立となるだけでなく、同一音源に由来するはずのＹ１［Ｌ'］３４２とも独立になる。
【０１４９】
一方、上記で使用するスコア関数（式［５．７］ほか）を変更することで、先に、［（１）時間周波数領域において、畳み込み混合を直接解く］方式において説明した式［８．１］〜式［８．４］を用いた方式と同等の処理も行なうことが可能である。
【０１５０】
［（１）時間周波数領域において、畳み込み混合を直接解く］方式において説明した式［８．１］〜式［８．４］を用いた方式は、「Ｙｋ（ｔ−Ｌ'）〜Ｙｋ（ｔ）は互いに依存」という仮定に基づく処理である。本処理例［（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理］においても、分離結果の依存を考慮した処理が可能である。すなわち、図１３を参照して説明すると、Ｙ１［０］３４１は、Ｙｎ［０］３４３やＹｎ［Ｌ'］３４４とは独立で、Ｙ１［Ｌ'］３４２とは依存関係があるという分離を行なうことができる。以下では、その方法について説明する。
【０１５１】
同一の音源に由来する分離結果Ｙｋ［０］，…，Ｙｋ［Ｌ'］の間で依存関係を持たせるには、先に説明したΔＷ（ω）の算出式［５．２］の代わりに、以下に示す式［１２．１］を用いる。
【０１５２】
【数１３】
2011215649.tif 000015
2011215649.tif 000016
【０１５３】
ただし、式［１２．１］中のＹ''（ω，ｔ）およびＷ''（ω）は、それぞれＹ''およびＷ''からω番目の周波数ビンの成分を抽出したベクトルおよび行列であり、式［１２．３］および式［１２．４］のように表わされる。φω（Ｙ''（ｔ））は、式［１２．５］で表わされる通り、ｎ×（Ｌ'＋１）個のスコア関数を要素に持つベクトルである。（スコア関数の具体例は後述する。）
【０１５４】
式［１２．５］と式［５．６］との違いは、スコア関数の引数にある。すなわち、式［５．６］をｎ×（Ｌ'＋１）チャンネルに拡張した場合は、ｎ×（Ｌ'＋１）個のスコア関数が全て異なる引数をとるのに対し、式［１２．５］ではφｋω［０］（Ｙｋ''（ｔ））〜φｋω［Ｌ'］（Ｙｋ''（ｔ））は同一の引数Ｙｋ''（ｔ）をとるため、引数はｎ種類である。
【０１５５】
スコア関数φｋω［ｌ］（Ｙｋ''（ｔ））は、Ｙｋ''（ｔ）（すなわちＹｋ［０］，…，Ｙｋ［Ｌ'］）を引数とする多次元（多変量）確率密度関数の対数微分として定義される（式［１２．５］）。このように、一つの確率密度関数に複数の引数を含め、そこから導出されるスコア関数を用いてＩＣＡの学習を行なうと、引数となっている要素同士は依存関係を持つ（独立にはならない）ことが理論的に示されている。すなわち、再び図１３を参照して説明すると、Ｙ１［０］，…，Ｙ１［Ｌ'］の組である信号Ｙ１''３５１は、組の内部の各要素は依存関係を持つが、別の組の信号、例えばＹｎ''３５２とは独立になる。
【０１５６】
ここで、多次元確率密度関数とスコア関数の具体例について説明する。多次元確率密度関数の一種に、球状分布と呼ばれるものがある。これは、以下に示す式［１３．１］の通り、スカラーを引数にとる関数にベクトルのＬ２ノルムを代入することで生成される（「∝」は比例を表わす）。
【０１５７】
【数１４】
2011215649.tif 000017
【０１５８】
Ｌ２ノルムは、各要素の（絶対値の）２乗和の平方根であり、式［１３．２］のｍに２を代入することで得られる。球状分布の例として、式［１３．３］のような指数分布に基づくもの（γは正の実数）を用いると、対応するスコア関数として式［１３．４］が導出される。この式を［１２．５］に代入すればよい。
【０１５９】
なお、先に、［（１）時間周波数領域において、畳み込み混合を直接解く］方式において説明した式［７．６］と同様に、式［１３．４］に対しても変更を及ぼして構わない。その例を式［１３．５］に示す。変更の例としては、
１）ゼロ除算を防ぐために分母に正の値βｋ［ｌ］（ω）を加える。さらに、その値として、ｋやｌ（エル）やωごとに異なるものを用いる。
２）Ｌ２ノルムの代わりにＬ−ｍノルム（式［１３．２］）を用いる。
３）スコア関数の係数Ｋの代わりに、ｋやｌ（エル）やωごとに異なる正値γｋ［ｌ］（ω）を用いる
このようなことを行なう。
【０１６０】
式［１２．１］は自然勾配法に基づく更新式であるが、それ以外のアルゴリズムも使用可能である。例えば、「独立性による等分散適用的分離」（Ｅｑｕｉｖａｒｉａｎｔ  Ａｄａｐｔｉｖｅ  Ｓｅｐａｒａｔｉｏｎ  ｖｉａ  Ｉｎｄｅｐｅｎｄｅｎｃｅ：  ＥＡＳＩ）と呼ばれる、信号の無相関化と分離とを同時に行なうアルゴリズムに基づく更新式は、式［１２．２］の通りである。このアルゴリズムを用いると、自然勾配法に比べて少ないループ回数で学習を収束させることができる。
【０１６１】
なお、式［１２．１］および式［１２．２］において行列の要素の対称性に注目すると、計算量を削減することが可能である。その点について説明する。
【０１６２】
式［１２．１］のＥｔ［・］の内部は、式［１２．７］のような（Ｌ'＋１）ｎ×（Ｌ'＋１）ｎの行列に展開される（上線は共役複素数を表わす）。この式の各要素について平均を取る際に、各要素の１番目の項であるφｋω［α］（Ｙｋ''（ｔ））と２番目の要素であるＹｉ［β］（ω，ｔ）との間（α，βは０≦α，β≦Ｌ'を満たす整数）で相対的なシフト量が同じであれば、平均後の値はほぼ同じ値となる。すなわち、式［１２．８］の関係が成り立つ。特に、シフトとして前述の巡回シフトを用いた場合は、全く同一の値となる。
【０１６３】
この性質を用いると、式［１２．７］の｛（Ｌ'＋１）ｎ｝２個の要素のうち、実際に値を計算する必要のあるのは２（Ｌ'＋１）ｎ２個だけでよく、残りの要素については、式［１２．８］に従って値を再利用すればよい。
【０１６４】
同様に、式［１２．２］についても、計算量の削減が可能である。Ｅｔ［・］の内部の３項のうち、１番目は式［１２．１］と同様の計算が行なえる。また、２番目については１番目の項を求めた後、単にエルミート転置を計算すればよい（式［１２．９］）。３番目の項については、式［１２．１０］の変形を行なうことで、計算量削減が可能となる。ただし、式［１２．１０］のＸ''（ω，ｔ）は式［１１．１］からω番目の周波数ビンに対応した要素を抽出したベクトルであり、式［１２．１１］のように表わせる。
【０１６５】
１７あるＥｔ［Ｘ''（ω，ｔ）Ｘ''（ω，ｔ）Ｈ］は学習のあいだ常に一定である。従って、Ｅｔ［Ｘ''（ω，ｔ）Ｘ''（ω，ｔ）Ｈ］については、学習前に一度だけ計算しておけばよく、学習中に毎回平均操作を行なう必要はない。すなわち、式［１２．１０］の左辺よりも右辺の方が計算量を削減できるのである。
【０１６６】
さらに、Ｅｔ［Ｘ''（ω，ｔ）Ｘ''（ω，ｔ）Ｈ］の計算においては、式［１２．８］と同様の対称性である式［１２．１２］および、対角線に対する対称性である式［１２．１３］が成り立つ。そのため、｛（Ｌ'＋１）ｎ｝２個の要素のうち、実際に値を計算する必要があるのは（Ｌ'＋１）ｎ２個だけでよい。
【０１６７】
［具体的構成例および処理例］
本発明の信号分離装置の構成例を図１５および図１６に示す。図１５は時間周波数領域で畳み込み混合を解く方式を実行する信号分離装置、図１６はモジュレーション・スペクトログラムに変換してから瞬時混合を解く方式を実行する信号分離装置に対応する構成例である。
【０１６８】
（１）時間周波数領域で畳み込み混合を解く方式を実行する構成
先に、図１５に示す時間周波数領域で畳み込み混合を解く方式を実行する信号分離装置の構成および処理について説明する。なお、以下に説明する処理の統括的な制御は制御部４０９において実行される。制御部４０９は、例えば、装置の記憶部（図示せず）に予め記憶された以下に説明する処理を実行するプログラムに従って処理を制御する。以下、各構成部の処理について説明する。複数の音源が発する独立な音を複数のマイク４０１で観測し、ＡＤ変換部４０２において入力アナログ信号をデジタル信号に変換してデジタル観測信号を得る。
【０１６９】
デジタル観測信号は短時間フーリエ変換（ＳＴＦＴ）部４０３に入力され短時間フーリエ変換処理が行なわれ、観測信号のスペクトログラムを得る。ここまでの処理は、例えば、図６（ｂ）に示す観測信号のスペクトログラムＸを得る処理に相当する。
【０１７０】
信号分離部４０４は、短時間フーリエ変換（ＳＴＦＴ）部４０３で生成した観測信号のスペクトログラムＸを独立な成分に分離する。この図１５に示す信号分離装置は、時間周波数領域において畳み込み混合された観測信号の分離処理として時間周波数領域において畳み込み混合を直接解く方式を適用しており、式［６．２］および［７．２］・［７．１］・［７．８］の演算を分離行列および分離結果が十分収束するまで（または一定回数）繰り返し行なうことで観測信号の分離処理を実行する。この分離処理によって、図６（ｃ）に示す分離結果Ｙを得る。
【０１７１】
なお、畳み込み演算部４０８で行なう処理は、先に、図６を参照して説明した処理に従った処理である。すなわち、観測信号中のｔ番目のフレームの観測信号Ｘ（ｔ）が、遅延の最大値をＬ＋１として、前のＬ＋１フレーム分の原信号の影響を受けていることを考慮した処理である。すなわち、観測信号Ｘ（ｔ）を前述したフレームタップ数Ｌを適用した式［６．１］の畳み込み混合で表し、さらに、図６（ｃ）の分離信号中のＹ（ｔ＋Ｌ'）を基準とし、Ｓ（ｔ）を推定するため、直前のＬ＋１フレーム分のデータを考慮して、分離信号Ｙ（ｔ）を式［６．２］のようにＸ（ｔ−Ｌ'）からＸ（ｔ）までの畳み込み混合として表わして、式［６．２］と式［７．２］を適用して行なう畳み込み演算である。
【０１７２】
観測信号Ｘから分離結果Ｙ、すなわち、図６（ｂ）の観測信号Ｘから図６（ｃ）の分離結果Ｙを生成するためのフレームタップ数Ｌ'は、前述したように、Ｌが既知であれば（すなわち、残響時間が既知であれば）、L'＝αLとすればよく(αは適切な正の実数)、Ｌが未知である場合、Ｌ'は、例えば、以下のいずれかの方法で決定する。
（ａ）Ｌ'＝６４やＬ'＝１００といった一定の値に決め打つ。
（ｂ）残響時間を測定し、その残響時間から求めたＬの値をＬ'とする。
（ｃ）さまざまなＬ'の下で分離を行ない、最も良い分離結果をもたらすＬ'の値を採用する。例えばＳＩＲ（ｓｉｇｎａｌ−ｉｎｔｅｒｆｅｒｅｎｃｅ  ｒａｔｉｏ）という分離性能尺度を計算し、最高のＳＩＲをもたらすＬ'を採用する。
上記いずれかの方法によって、Ｌ'、すなわち、観測信号Ｘから分離結果Ｙを生成するためのフレームタップ数Ｌ'、具体的には、例えば、図６（ｂ）の観測信号Ｘから図６（ｃ）の分離結果Ｙを生成するためのフレームタップ数Ｌ'を決定し、このフレームタップ数Ｌ'を用いて、観測信号の複数の連続フレームから分離結果を生成する。
【０１７３】
リスケーリング部４０５では、分離信号の各周波数ビンに対してスケールを揃えるリスケーリング処理を行なう。リスケーリングとは、周波数ビンごとのスケールを調整する処理である。また、分離処理前に観測信号に正規化（平均や分散の調整）を行なっていた場合は、ここで元に戻す。
【０１７４】
逆フーリエ変換部４０６は、逆フーリエ変換によって分離信号のスペクトログラムを時間領域の信号へと変換する。変換された信号は、必要に応じて後段処理実行部４０７へ送られる。後段の処理とは、スピーカーからの再生や音声認識などである。なお、後段の処理によっては、逆フーリエ変換部を省略することも可能である。
【０１７５】
このように、図１５に示す信号分離装置は、複数の音信号が混合した信号を入力して個別の音信号に分離する信号分離装置であり、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換手段（ＳＴＦＴ部４０３）と、信号変換手段（ＳＴＦＴ部４０３）の生成した観測信号スペクトログラムから信号分離結果を生成する信号分離手段（信号分離部４０４）を有し、信号分離手段（信号分離部４０４）は、観測信号スペクトログラムを時間周波数領域において畳み込み混合された観測信号として解釈して、畳み込み演算部４０８における畳み込み演算の実行により信号分離結果を生成する。
【０１７６】
なお、信号変換手段（ＳＴＦＴ部４０３）は、入力信号に対して短時間フーリエ変換（ＳＴＦＴ）を実行して時間周波数領域に変換し観測信号スペクトログラムを生成する処理を実行する。
【０１７７】
また、信号分離手段（信号分離部４０４）は、フレーム番号（ｔ）の分離信号Ｙ（ｔ）を、観測信号Ｘ（ｔ−Ｌ'）〜Ｘ（ｔ）の畳み込み混合として設定し、分離信号Ｙ（ｔ）に含まれる個別の音声信号成分であるＹ１（ｔ）〜Ｙｎ（ｔ）各々の独立性を高める処理により信号分離結果を生成する。具体的には、分離信号Ｙ（ｔ）に含まれる個別の音信号成分であるＹ１（ｔ）〜Ｙｎ（ｔ）各々の独立性を高める処理として、独立性算出尺度であるＫｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量Ｉ（Ｙ）を適用し、Ｋｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量Ｉ（Ｙ）を最小にする分離行列の更新処理により信号分離結果を生成する。
【０１７８】
なお、（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理によって時間周波数領域において畳み込み混合された観測信号を分離する処理を実行する装置構成としては、例えば、図１５に示す構成から畳み込み演算部４０８を省いた構成が適用できる。なお、信号分離部において実行する処理は異なる。
【０１７９】
（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理を行なう装置では、ＳＴＦＴ部４０３が、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する信号変換手段として機能し、信号分離部４０４は、信号変換手段の生成した観測信号スペクトログラムから信号分離結果を生成する処理を行なう構成であり、信号分離部４０４では、先に図１１〜図１４を参照して説明したように、観測信号スペクトログラムをフレーム方向へシフトさせて、各々が異なるシフト量を持つデータを積み重ねた観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により、信号分離結果を生成する。なお、瞬時混合ＩＣＡを適用した処理は、特開２００６−２３８４０９に開示された手法、すなわち、時間周波数領域の観測信号と分離行列から時間周波数領域の分離信号を生成し、生成した時間周波数領域の分離信号と、多次元確率密度関数から導出される多次元スコア関数によって計算される分離行列とがほぼ収束するまで分離行列を修正し、修正した分離行列を適用して時間周波数領域の分離信号を生成する処理として実行される。
【０１８０】
（２）モジュレーション・スペクトログラムに変換してから瞬時混合を解く方式を実行する構成
次に、図１６に示すモジュレーション・スペクトログラムに変換してから瞬時混合を解く方式を実行する信号分離装置の構成、および処理について説明する。なお、以下に説明する処理の統括的な制御は制御部４６１において実行される。制御部４６１は、例えば、装置の記憶部（図示せず）に予め記憶された以下に説明する処理を実行するプログラムに従って処理を制御する。以下、各構成部の処理について説明する。複数の音源が発する独立な音を複数のマイク４５１で観測し、ＡＤ変換部４５２において入力アナログ信号をデジタル信号に変換してデジタル観測信号を得る。
【０１８１】
デジタル観測信号は、第１短時間フーリエ変換（ＳＴＦＴ）部４５３に入力され、短時間フーリエ変換（ＳＴＦＴ）処理が行なわれ、観測信号のスペクトログラムを得る。この段階で得られる信号は、例えば、図８（ｂ）に示すスペクトログラムＸである。さらに、第１段階の短時間フーリエ変換（ＳＴＦＴ）処理によって得られた観測信号のスペクトログラムを第２短時間フーリエ変換（ＳＴＦＴ）部４５４に入力して、周波数ビンごとに再び短時間フーリエ変換（ＳＴＦＴ）を実行し、モジュレーション・スペクトログラムを得る。
【０１８２】
この第２短時間フーリエ変換（ＳＴＦＴ）部４５４における短時間フーリエ変換（ＳＴＦＴ）によって得られるモジュレーション・スペクトログラムが例えば図９（ｃ），（ｄ）に示すモジュレーション・スペクトログラムＸ'である。
【０１８３】
信号分離部４５５は、モジュレーション・スペクトログラムＸ'を入力して、このモジュレーション・スペクトログラムＸ'を独立な成分に分離する。この分離処理は、先に図１０を参照して説明した処理である。すなわち、図１０は、図９（ｃ）に示す立体版モジュレーション・スペクトログラムＸ'に相当し、この図１０に示す立体版モジュレーション・スペクトログラムＸ'において、例えば、１番目のチャンネルのエントロピー計算においては、図１０における１番目のフレームのモジュレーション・スペクトログラムＹ１'（ｔ）２２１は、平面を表わし、それを引数とする多変量確率密度関数Ｐ（Ｙ１'（ｔ））２２２にＹ１'（ｔ）を代入することで、エントロピーＨ（Ｙ１'）２２３を求める。式［９．３］は、変数名の違いを除くと式［３．５］と同一である。従って、学習規則を導出するためには、式［５．２］の変数名を付け替えればよく、結果として式［９．５］を得る。すなわち、式［９．３］，［９．５］，［９．６］をＷ'が収束するまで繰り返せば、Ｙ１'（ｔ）〜Ｙｎ'（ｔ）が互いに独立になる。
【０１８４】
次に、第１リスケーリング部４５６でモジュレーション・スペクトログラムに対してリスケーリングを行なう。リスケーリングとは、周波数ビンごとのスケールを調整する処理である。さらに、第１逆フーリエ変換（ＦＴ）部４５７で、リスケーリングされたモジュレーション・スペクトログラムに対して逆フーリエ変換（ＦＴ）処理を実行して、モジュレーション・スペクトログラムをスペクトログラムへ変換する。その後、再び第２リスケーリング部４５８でリスケーリングを行なった後、第２逆フーリエ変換（ＦＴ）部４５９で、リスケーリングされたスペクトログラムに対して逆フーリエ変換（ＦＴ）処理を実行して、スペクトログラムを波形に変換する。波形に変換された信号は、必要に応じて後段処理実行部４６１へ送られ、必要に応じた後段処理を実行する。後段の処理とは、スピーカーからの再生や音声認識などである。
【０１８５】
このように、図１６に示す信号分離装置は、複数の信号が混合した信号を入力して個別の信号に分離する信号分離装置であり、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成する第１信号変換手段（第１ＳＴＦＴ部４５３）と、第１信号変換手段（第１ＳＴＦＴ部４５３）の生成した観測信号スペクトログラムに対するデータ変換を実行しモジュレーション・スペクトログラムを生成する第２信号変換手段（第２ＳＴＦＴ部４５４）と、第２信号変換手段（第２ＳＴＦＴ部４５４）の生成した前記モジュレーション・スペクトログラムから信号分離結果を生成する信号分離手段（信号分離部４５５）を有し、信号分離手段（信号分離部４５５）は、モジュレーション・スペクトログラムを瞬時混合として解釈し信号分離結果を生成する。
【０１８６】
第１信号変換手段（第１ＳＴＦＴ部４５３）は、入力信号に対して短時間フーリエ変換（ＳＴＦＴ）を実行して時間周波数領域に変換し観測信号スペクトログラムを生成する。さらに、第２信号変換手段（第２ＳＴＦＴ部４５４）は、観測信号スペクトログラムに対して時間方向の短時間フーリエ変換（ＳＴＦＴ）を実行しモジュレーション・スペクトログラムを生成する。
【０１８７】
信号分離手段（信号分離部４５５）は、モジュレーション・スペクトログラムに含まれる分離信号対応の信号成分Ｙ１'〜Ｙｎ'各々の独立性を高める処理により信号分離結果を生成する。具体的には、分離信号対応の信号成分Ｙ１'〜Ｙｎ'各々の独立性を高める処理として、独立性算出尺度であるＫｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量を適用し、Ｋｕｌｌｂａｃｋ−Ｌｅｉｂｌａｒ情報量を最小にする分離行列の更新処理により信号分離結果を生成する。
【０１８８】
なお、逆フーリエ変換手段（第１逆ＦＴ部４５７）は、信号分離手段（信号分離部４５５）において得られた分離信号対応の信号成分Ｙ１'〜Ｙｎ'各々に対して逆フーリエ変換を実行して分離信号対応のスペクトログラムＹ１〜Ｙｎを生成する。
【０１８９】
本発明の信号分離装置の実行する処理のシーケンスの一例について図１７に示すフローチャートを参照して説明する。ステップＳ１０１において、マイクで音を観測する。例えば先に図５を参照して説明したように、複数の音源から出力される音の混合信号をマイクで取得する。次に、ステップＳ１０２において、観測信号に対する短時間フーリエ変換（ＳＴＦＴ）処理を実行しスペクトログラムを得る。短時間フーリエ変換は、先に図７を参照して説明した処理であり、この処理によって、スペクトログラムを得る。このスペクトログラムは、例えば図６（ｂ）に示すスペクトログラムである。
【０１９０】
次に、ステップＳ１０３において、観測信号のスペクトログラムに対し、ＩＣＡによる分離処理を行なう。分離処理の処理シーケンスの詳細については後述する。分離結果に対し、ステップＳ１０４において、必要に応じて逆フーリエ変換（ＩＦＴ）を実行し、その後、必要に応じてステップＳ１０５において後段処理を実行する。
【０１９１】
ステップＳ１０３において実行する分離処理の詳細シーケンスについて、図１８と図１９に示すフローチャートを参照して説明する。
図１８、図１９に示す分離処理シーケンスは、それぞれ、先に、図１５、図１６を参照して説明した信号分離装置において実行する分離処理の具体的シーケンであり、
図１８は、図１５の信号分離装置の実行する時間周波数領域で畳み込み混合を解く方式における分離処理、
図１９は、図１６の信号分離装置の実行するモジュレーション・スペクトログラムに変換してから瞬時混合を解く方式における分離処理、
これらの詳細シーケンスである。
【０１９２】
まず、図１８を参照して、図１５の信号分離装置の実行する時間周波数領域で畳み込み混合を解く方式における分離処理、すなわち時間周波数領域での逆畳み込みを行なう分離処理シーケンスについて説明する。
【０１９３】
最初に、ステップＳ２０１において、観測信号スペクトログラムに対して正規化を行なう。本処理における正規化処理は、スペクトログラムの各周波数ビンに対して、平均を０に設定し分散を１とする処理、または以降の処理に都合のよい値に調整する処理である。次に、ステップＳ２０２において、分離行列の初期化処理、すなわち、分離行列Ｗ［τ］に初期値を代入する。初期値は、Ｗ［０］に対しては単位行列を、τ＞０の分離行列Ｗ［τ］に対してはゼロ行列を代入すればよい。または、前回の学習で求まった分離行列が存在している場合は、それを初期値として用いても良い。
【０１９４】
ステップＳ２０３〜Ｓ２１０は学習のループであり、分離行列および分離結果が収束するまでこのループを繰り返す。すなわち、
ステップＳ２０３：分離行列が収束したか否かの判定、
ステップＳ２０４：分離信号Ｙの計算、
ステップＳ２０５：周波数ビンループの開始（ω＝１，...，Ｍ）、
ステップＳ２０６：フレームタップループの開始（τ＝０，...，Ｌ）、
ステップＳ２０７：τ番目のフレームタップに対応する増分ΔＷ［τ］の計算、
ステップＳ２０８：フレームタップループの終了、
ステップＳ２０９：ΔＷ［０］（ω）〜Ｗ［Ｌ'］（ω）の更新、
ステップＳ２１０：周波数ビンループの終了、
これらのステップからなるループを繰り返し実行する。
【０１９５】
ステップＳ２０４の分離結果Ｙの計算には、先に説明した式［６．２］、または式［６．３］を用いる。（Ｙ＝［Ｙ（１），...，Ｙ（Ｔ）］とする。）ステップＳ２０５〜Ｓ２１０は周波数ビンについてのループであり、Ｍを周波数ビンの本数として、１≦ω≦Ｍを満たす各周波数（ω）について、ステップＳ２０６〜Ｓ２０９を繰り返す。なお、ループの代わりに、周波数ビンごとの並列処理を行なう構成としてもよい。なお、本出願と同一出願人の先の特許出願であり公開された特開２００６−２３８４０９において示した手法法では、推定する分離行列は１つ（または周波数ビンごとに１つ）だけであったが、本発明ではフレームタップの個数だけ分離行列を推定する必要がある。そこで、フレームタップの個数分だけループを回す（ステップＳ２０６〜Ｓ２０８）。
【０１９６】
ステップＳ２０７では、τ番目のフレームタップに対応する増分ΔＷ［τ］（ω）を求める。ΔＷ［τ］（ω）の計算には、式［７．１］を用いる。前述の通り、この式［７．１］のＲω［ｌ］は、分離結果Ｙの計算に式［６，２］、または式［６．３］のいずれの式を用いたかによって異なる。
【０１９７】
分離結果Ｙの計算に式［６．２］を用いた場合はＲω［ｌ］の計算に式［７．２］または式［８．１］を用い、分離結果Ｙの計算に式［６．３］を用いた場合はＲω［ｌ］の計算に式［７．３］または式［８．２］を用いる。
【０１９８】
ステップＳ２０６〜Ｓ２０８のフレームタップのループを抜けた後、ステップＳ２０９において、式［７．８］を用いて分離行列ΔＷ［０］（ω）〜Ｗ［Ｌ'］（ω）を更新する。なお、この処理は、ステップＳ２１０の後で全周波数ビンの分をまとめて行なっても構わない。（一方、フレームタップの内部には入れられないことに注意。）
【０１９９】
ステップＳ２０５〜Ｓ２１０の周波数ビンのループを抜けた後、再び、ステップＳ２０３の収束チェックに戻る。ステップＳ２０４において、で分離行列が収束した（または、所定の回数だけループした）と判定された場合は、分岐を右に進みステップＳ２１１に移行する。
【０２００】
なお、ステップＳ２０３における分離行列が収束したか否かの判定は、例えばΔＷのノルム‖ΔＷ‖（行列のノルムは、例えば前記した式［７．１０］で計算する）がある値を下回ったかどうか（または‖ΔＷ‖／‖Ｗ‖がある値を下回ったかどうか）で判断してもよいし、または、単純に一定回数のループ数を予め設定し、そのループ数を実行してもよい。
【０２０１】
ステップＳ２０３において、分離行列がまだ収束していないと判定された場合は、ステップＳ２０４〜Ｓ２１０の処理を繰り返し実行する。ステップＳ２０４において、で分離行列が収束した（または、所定の回数だけループした）と判定された場合は、分岐を右に進みステップＳ２１１に移行し、ステップＳ２１１において、リスケーリングを行なう。リスケーリングとは、周波数ビンごとのスケールを調整する処理である。また、正規化処理ステップ（Ｓ２０１）で周波数ビンの平均や分散を変更した場合は、ここで必要に応じて元に戻す。
【０２０２】
なお、ステップＳ２１１において実行するリスケーリングの係数は、以下のようにして求める。先に示した式［７．１１］を用いて、ある周波数ビンにおいて観測信号と分離結果との２乗誤差が最小となるようなスケールを求める（具体的には、最小二乗法などを用いる）。そして分離結果を、そのスケールを乗じた値に更新する（式［７．１２］）。また、必要に応じて、分離行列自体も同様に更新する（式［７．１３］）。
【０２０３】
または、以下のように行なっても良い。式［７．１４］を用いて、観測信号を分離結果と定数との線形和で表現する。スケールαｋ１（ω）〜αｋｎ（ω）および定数項βｋ（ω）は、式［７．１５］で求める（具体的には、最小二乗法などを用いる）。スケールが求まったら、式［７．１６］を用いて分離結果を更新する。（必要に応じて、分離行列も更新する。）
【０２０４】
なお、式［７．１４］に出てくるαｋｊ（ω）Ｙｊ（ω，ｔ）を全て出力すると、ｓｉｎｇｌｅ  ｉｎｐｕｔ  ｍｕｌｔｉｐｌｅ  ｏｕｔｐｕｔ  （ＳＩＭＯ）形式の出力が得られる。ＩＣＡのＳＩＭＯ出力とは、「観測信号を、それぞれの音源に由来する成分に分解する」ことであり、例えばαｋｊ（ω）Ｙｊ（ω，ｔ）は、Ｙｊをｉ番目の音源の推定結果とすると、「ｋ番目のマイクで観測される信号の内、ｉ番目の音源に由来する成分」を表わしている。以上で、時間周波数領域で畳み込み混合を解く場合についてのフローチャートの解説を終わる。
【０２０５】
次に、モジュレーション・スペクトログラム領域で瞬時混合を解く場合の処理について、図１９に示すフローチャートを参照して説明する。図１９は、図１６の信号分離装置の実行するモジュレーション・スペクトログラムに変換してから瞬時混合を解く方式における分離処理の詳細シーケンスである。
【０２０６】
ステップＳ３０１は、観測信号スペクトログラムに対して正規化を行なう。この処理は、図１８のフローにおけるステップＳ２０１の正規化処理と同様の処理であり、スペクトログラムの各周波数ビンに対して、平均を０に設定し分散を１とする処理、または以降の処理に都合のよい値に調整する処理である。ステップＳ３０２では、周波数ビンごとに短時間フーリエ変換（ＳＴＦＴ）を行ない、モジュレーション・スペクトログラム、すなわち、図９（ｃ），（ｄ）に示すモジュレーション・スペクトログラムＸ'を生成する。
【０２０７】
なお、このモジュレーション・スペクトログラムの生成には、先に図１６を参照して説明したように、デジタル観測信号に対する第１短時間フーリエ変換（ＳＴＦＴ）部４５３における短時間フーリエ変換（ＳＴＦＴ）処理を行なって観測信号のスペクトログラム（例えば、図８（ｂ）に示すスペクトログラムＸ）を得て、さらに、第１段階の短時間フーリエ変換（ＳＴＦＴ）処理によって得られた観測信号のスペクトログラムを第２短時間フーリエ変換（ＳＴＦＴ）部４５４に入力して、周波数ビンごとに再び短時間フーリエ変換（ＳＴＦＴ）を実行することが必要となる。この第２短時間フーリエ変換（ＳＴＦＴ）部４５４における短時間フーリエ変換（ＳＴＦＴ）によって得られるモジュレーション・スペクトログラムが例えば図９（ｃ），（ｄ）に示すモジュレーション・スペクトログラムＸ'である。
【０２０８】
モジュレーション・スペクトログラムは、図９（ｃ），（ｄ）に示すように、立方体形式（式［９．２］に相当）と、平面形式（式［９．３］に相当）とがあるが、以降の説明では平面形式を用いる。すなわち、図９（ｃ）に示す縦方向と奥行き方向の両方のビンをまとめてω'というインデックスで表現する。
【０２０９】
ステップＳ３０３では、モジュレーション・スペクトログラムの各ビンω'に対して、再び正規化を行なう。学習のループの前に、ステップＳ３０４では、分離行列Ｗ'に初期値を代入しておく。初期値は、単位行列でよいが、前回の学習で求まった分離行列でも良い。
【０２１０】
ステップＳ３０５〜ステップＳ３１０は学習のループであり、分離行列Ｗ'が収束するまで（または一定回数）繰り返す。ステップＳ３０５における収束性判定は、図１８を参照して説明したステップＳ２０３における処理と同様の判定であり、分離行列が収束したか否かの判定は、例えばΔＷ'のノルム‖ΔＷ'‖（行列のノルムは、例えば前記した式［７．１０］で計算する）がある値を下回ったかどうか（または‖ΔＷ'‖／‖Ｗ'‖がある値を下回ったかどうか）で判断してもよいし、または、単純に一定回数のループ数を予め設定し、そのループ数を実行してもよい。
【０２１１】
ステップＳ３０６で、分離結果モジュレーション・スペクトログラムであるＹ'を計算する。この計算は、式［９．３］を全てのω'とｔに対して行なえばよい。
【０２１２】
ステップＳ３０７〜Ｓ３１０は、図９（ｃ）に示すモジュレーション・スペクトログラムの各ビンω'、すなわち縦方向と奥行き方向の両方のビンω'についてのループである。なお、各ビンについての繰り返し処理とするループの代わりに、各ビンについての処理を並列処理として実行してもよい。ステップＳ３０８では、分離行列の増分を計算し（式［９．５］）、ステップＳ３０９においてで分離行列を更新する（式［９．６］）。
【０２１３】
ステップＳ３１０において、ループを抜けた後、再びステップＳ３０５の収束性判定に戻る。ステップＳ３０５において、分離行列が収束した（または所定の回数ループした）と判定された場合は、条件分岐を右に進む。ステップＳ３１１において、リスケーリングを行なう。リスケーリングは、各ビンのスケールを調整する処理である。分離結果のモジュレーション・スペクトログラムに対してリスケーリングを行なう。リスケーリングの方法は、先に図１８を参照して説明したステップＳ２１１の処理とほぼ同様であり、式［７．１１］〜［７．１６］の式のＹ，Ｘ，Ｗを適宜Ｙ'，  Ｘ'，Ｗ'に置き換えた式に基づいて行なう。また、必要に応じて、ステップＳ３０１の正規化を元に戻す処理も行なう。
【０２１４】
次に、ステップＳ３１２において、モジュレーション・スペクトログラムをスペクトログラムへ変換する逆フーリエ変換（ＦＴ）を実行する。その際、必要に応じてＷＯＬＡなどを行なう。すなわち、逆フーリエ変換（ＦＴ）においてもフレームごとの逆変換結果（波形）波形を重複つきで重ね合わせる。これをオーバラップ加算（ｏｖｅｒｌａｐ  ａｄｄ）という。逆変換結果は、オーバラップ加算（ｏｖｅｒｌａｐ  ａｄｄ）の前にサイン窓等の窓関数を再び作用させても良く、これを、ｗｅｉｇｈｔｅｄ  ｏｖｅｒｌａｐ  ａｄｄ（ＷＯＬＡ）という。ＷＯＬＡにより、フレーム間の不連続性に由来するノイズを低減することができる。
【０２１５】
さらに、ステップＳ３１３において、スペクトログラムに対してリスケーリングを行なう。これは、ステップＳ３１１のリスケーリングと同様の処理である。
【０２１６】
なお、図１７のフローにおけるステップＳ１０４と、図１９に示すフローにおけるステップＳ３１２において実行する逆フーリエ変換（ＦＴ）においては、分離結果のスペクトログラムやモジュレーション・スペクトログラムの他に、分離行列自体も、必要に応じて逆フーリエ変換（ＦＴ）を施す。
【０２１７】
［変形例］
次に、上述した実施例を変形した実施例について説明する。上記の実施例では、分離結果を生成する際に適用するフレームタップＬ'は、すなわち観測信号から分離結果を生成する際のフレームタップＬ'は、全ての周波数で一定値を用いていた。これを、周波数毎で一律とするのではなく周波数ごとにフレームタップＬ'の値を変えても良い。
【０２１８】
例えば、高い周波数の成分は、低い周波数の成分と比べて急激に減衰するため、残響時間は短い。そこで、高い周波数に対応する周波数ビンでは、フレームタップＬ'の値を低い周波数ビンよりも少なくしても良い。こうすることで、分離性能を保ったまま計算量を削減することができる。
【０２１９】
また、図１６の信号分離装置、および図１９に示すフローチャートを参照して説明した方法、すなわち、モジュレーション・スペクトログラムに変換してから瞬時混合を解く方式における分離処理においては、２度目の短時間フーリエ変換（ＳＴＦＴ）において、周波数ビンごとにフレームタップ数Ｌ'を異ならせる他に、シフト幅を異ならせることも可能である。ただし、周波数ビンごとにフレームタップ数やシフト幅を異ならせると、モジュレーション・スペクトログラムにおいて１フレームあたりの時間長が異なる可能性が出てくる。
【０２２０】
例えば、２度目の短時間フーリエ変換（ＳＴＦＴ）において、
低い周波数はタップ数＝３２・シフト幅＝１６を用い、
高い周波数はタップ数＝１６・シフト幅＝８を用いると、
変換後のモジュレーション・スペクトログラムにおける１フレームあたりの時間長は、低い周波数が高い周波数の倍となる。すなわち、低い周波数の方が単位時間当たりのフレーム数が少ない（半分である）。
【０２２１】
１フレーム当たりの時間長が一定である場合は、図１０に示すように、モジュレーション・スペクトログラムからＹｋ'（ｔ）２２１を切り出して独立性を計算することができるが、一定でない場合はそれが難しい。そういう場合は、以下に説明する方法（方法１〜３）のいずれかを用いることで、フレームの不一致に対処する。
【０２２２】
（方法１）フレームデータの間引き
生成されたモジュレーション・スペクトログラムにおいて、単位時間当たりのフレーム数が多い方のビンからデータを間引くことで、フレーム数が少ない方のビンとデータの個数を合わせる。上記の３２タップ・１６シフトと１６タップ・８シフトの例では、１６タップ・８シフトの短時間フーリエ変換（ＳＴＦＴ）を行なったビンに対して、一つ置きにデータを間引くと、両者で単位時間当たりのフレーム数が一致する（＝１フレーム当たりの時間が同じ）ようになる。
【０２２３】
（方法２）フレームデータの補間
上述の（方法１）とは逆に、個数が少ない方を多い方に合わせる方法。上記の３２タップ・１６シフトと１６タップ・８シフトの例では、３２タップ・１６シフトを行なったビンに対して、データの補間を行なう。例えば、前後のフレームデータの平均を取ることで、中間に新しいデータを挿入する。
【０２２４】
（方法３）フレームデータの重複
上述の（方法２）と同様に、個数が少ない方を多い方に合わせる方法。上記の３２タップ・１６シフトと１６タップ・８シフトの例では、３２タップ・１６シフトを行なったビンに対して、データを２回ずつ重複させることで、１６タップ・８シフトのビンとデータ数を合わせる。
【０２２５】
次に、図１１〜図１４を参照して説明した（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理による方法、すなわち、観測信号スペクトログラムをシフトしながら積み重ね、それを瞬時混合ＩＣＡ（例えば、特開２００６−２３８４０９に記載の手法）によって、時間周波数領域において畳み込み混合された観測信号を分離する処理において、
「周波数ごとに［Ｌ'］の値、すなわち、観測信号から分離結果を生成する際のフレームタップ数［Ｌ'］の値を異ならせる」という変形例について説明する。
【０２２６】
この変形例、すなわち、「周波数ごとにフレームタップ数［Ｌ'］の値を異ならせる」という変形例を実現するためには、以下のようにすればよい。周波数ビンごとに異なるＬ'をＬ'（ω）と表記する。図１１を参照して説明したシフト処理において、シフト量が［Ｌ'（ω）］を超えたら、その周波数ビンのデータを［０］に置き換える。その方法について、図２０を参照して説明する。
【０２２７】
周波数ビンごとに異なるフレームタップ数［Ｌ'（ω）］の値を、周波数ビン番号ωに応じて以下のように変更したいとする。（Ｍはスペクトログラム１枚あたりの周波数ビン数）
１≦ω＜Ｍ／４では、Ｌ'（ω）＝２
Ｍ／４≦ω＜Ｍ／２では、Ｌ'（ω）＝１
Ｍ／２≦ω＜Ｍでは、Ｌ'（ω）＝０
これを実現するためには、先に、図１１を参照して説明したシフト処理によって生成されるデータＸｋ［０］，Ｘｋ［１］，Ｘｋ［２］に対して以下の操作を行なう。
Ｘｋ［０］は、Ｘｋそのまま。（全ての周波数ビンで、シフト＝０は必要）、
Ｘｋ［１］は、Ｍ／２≦ωの周波数ビンを０でマスク。（Ｍ／２≦ωでは、１以上のシフトは不要）
Ｘｋ［２］は、Ｍ／４≦ωの周波数ビンを０でマスク。（Ｍ／４≦ωでは、２以上のシフトは不要）
【０２２８】
具体的には、図２０（ｂ）に示すように、黒く塗りつぶされているデータ部分５１１２が、０でマスクされる箇所である。なお、実際の処理においては、マスクされた部分のメモリを確保する必要はなく、スペクトログラムへのアクセスの際にマスク該当箇所をスキップすれば、処理時間やメモリ量の増加を防ぐことができる。
【０２２９】
なお、本発明の前処理として、従来の時間周波数領域の瞬時混合ＩＣＡ（例えば特開２００６−２３８４０９）を組み合わせると、処理時間の増加をある程度抑えることができる。以降では、両者の組み合わせについて説明する。以下の各処理例について、順次説明する。
（１）基本的な２段階分離
（２）チャンネル数の削減
（３）残響除去として利用
【０２３０】
（１）基本的な２段階分離
従来の時間周波数領域の瞬時混合ＩＣＡにおいて、残響よりも短い分析フレーム（または分析窓）を用いた場合、複数のフレームにまたがる妨害音は除去しきれない。その反面、本発明よりも（１回目のＳＴＦＴにおける分析フレーム長が同じであれば）計算量は少ない。そこで、最初に従来の時間周波数領域ＩＣＡで分離を行い、その結果のスペクトログラムを本発明の方法でさらに分離すれば、最初から本発明だけを用いる場合と比べて少ない時間で同等の分離精度を達成することができる。
【０２３１】
特に、本発明における「（１）時間周波数領域において、畳み込み混合を直接解く」方式を用いる場合は、従来法と本発明とをシームレスに動作させることが可能である。すなわち、式［７．２］および式［８．１］（または、式［７．３］および式［８．２］）においてＬ'＝０とすると従来法と等価になるという特徴を利用することができ、図１８に示したフローにおけるステップＳ２０３〜Ｓ２１０の学習ループにおいて、ループ回数が少ないうちはＬ'＝０、ループ回数がある値を超えたらＬ'を本来の値とすればよい。または、ループ回数の増加にともなって、Ｌ'を少しずつ増加させても構わない。
【０２３２】
（２）チャンネル数の削減
一般に、ＩＣＡの計算量は、チャンネル数の２乗に比例する。そのため、チャンネル数を削減することができれば、計算量を大幅に削減することができる。２段階分離を用いると、本発明のステップのチャンネル数を削減することもできる。その方法について説明する。
【０２３３】
時間周波数領域のＩＣＡにおいて、音源数よりもマイク数の方が多い場合、出力チャンネルのうちのいくつかは、どの音源にも対応しないと判定される信号が出力される。例えば、マイク数＝４・音源数＝３の場合、出力チャンネルのうち３つは音源に対応しているが、残りの１つはどの音源にも対応しない、背景雑音と残響音とが混ざったような信号が出力される。このような出力は、他のチャンネルと比べてパワーが極端に小さかったり、他のどのチャンネルとも相関があったりするため、容易に検出できる。
【０２３４】
そこで、２段階分離においては、図２１に示すフローチャートに従って、まず、ステップＳ５０１で時間周波数領域の瞬時混合ＩＣＡによる分離処理を行なう。この処理は、特開２００６−２３８４０９において開示した処理として実行可能である。その後、ステップＳ５０２において、例えば「どの音源にも対応しないと判定される出力」（不要チャンネル）を除去した後、上述した本発明に従った処理、すなわち、
（１）時間周波数領域において、畳み込み混合を直接解く。
（２）スペクトログラムを時間方向へもう一度短時間フーリエ変換（ＳＴＦＴ）し、瞬時混合として解く。
（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理によって解く。
これらのいずれかの処理によって、時間周波数領域において畳み込み混合された観測信号を分離する処理を実行すれば、分離処理における計算量を削減することができる。なお、入力チャンネル数＝音源数であれば分離は可能であるため、ステップＳ５０２においてチャンネル数を削減しても分離精度には影響しない。
【０２３５】
例えば、上記の（１）時間周波数領域において、畳み込み混合を直接解く方式に、この２段階処理を適用した場合は、信号分離手段が、観測信号スペクトログラムに対して、瞬時混合ＩＣＡ（Ｉｎｄｅｐｅｎｄｅｎｔ  Ｃｏｍｐｏｎｅｎｔ  Ａｎａｌｙｓｉｓ）を適用した処理により第１の信号分離結果を生成し、該第１の信号分離結果から、どの音源にも対応しないと判定される不要チャンネル除去処理を実行して、除去処理後に残存する観測信号スペクトログラムに対して時間周波数領域の畳み込み混合を解く処理を実行して信号分離結果を生成する。
【０２３６】
また、（２）スペクトログラムを時間方向へもう一度短時間フーリエ変換（ＳＴＦＴ）し、瞬時混合として解く方式にこの２段階処理を適用した場合は、
第１信号変換手段が、入力信号を時間周波数領域に変換し観測信号スペクトログラムを生成し、不要チャンネル除去手段が、第１信号変換手段の生成した観測信号スペクトログラムに対して、瞬時混合ＩＣＡを適用した処理により第１の信号分離結果を生成し、この信号分離結果から、どの音源にも対応しないと判定される不要チャンネル除去処理を実行し、さらに、第２信号変換手段が、不要チャンネルが除去された観測信号スペクトログラムに対してデータ変換を実行してモジュレーション・スペクトログラムを生成し、信号分離手段が、モジュレーション・スペクトログラムから信号分離結果を生成するといった処理となる。
【０２３７】
また、（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理にこの２段階処理を適用した場合は、
観測信号スペクトログラムに対して、瞬時混合ＩＣＡを適用した処理により第１の信号分離結果を生成し、生成した第１の信号分離結果から、どの音源にも対応しないと判定される不要チャンネル除去処理を実行し、除去処理後に残存する観測信号スペクトログラムをフレーム方向へシフトさせて観測信号スペクトログラムシフトセットを生成し、生成した観測信号スペクトログラムシフトセットに対して、再度、瞬時混合ＩＣＡを適用して信号分離結果を生成する構成となる。
【０２３８】
（３）残響除去として利用
本発明の分離処理、すなわち、
（１）時間周波数領域において、畳み込み混合を直接解く。
（２）スペクトログラムを時間方向へもう一度短時間フーリエ変換（ＳＴＦＴ）し、瞬時混合として解く。
（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理によって解く。
これらの処理のうち、３番目の「シフト積み重ね＋従来法」を用いた場合、分離自体は前処理の従来法で行ない、本発明は残響除去のみを行なうという役割分担も可能である。こうすることで、計算量はＯ（｛ｎ×（Ｌ'＋１）｝２）からＯ（ｎ×ｎ×（Ｌ'＋１））に削減される。以下はその方法について説明する。
【０２３９】
先に図１１を参照して説明した「スペクトログラムのシフト＋積み重ね」を行なうと、１チャンネル分のスペクトログラムが見かけ上はＬ'＋１チャンネルに拡張される。その結果に対して、従来の時間周波数領域の瞬時混合ＩＣＡを用いてＬ'＋１チャンネルの入力として処理すると、結果としてＬ'＋１チャンネル分のスペクトログラムが生成される。このような処理を行なっても、音源ごとの成分に分離されるわけではないが、複数のフレームにまたがった成分を取り除く効果、すなわち残響除去の効果はある。そこで、従来の時間周波数領域の瞬時混合ＩＣＡで分離を行い、ｎチャンネル分の分離結果スペクトログラムを生成した後、各チャンネルに対して前述の「残響除去」を行なうという組み合わせが考えられる。そのような処理について、図２２に示すフローチャートを参照して説明する。
【０２４０】
まず、ステップＳ６０１において、時間周波数領域の瞬時混合ＩＣＡによる分離処理を行なう。この処理は、特開２００６−２３８４０９において開示した処理として実行可能である。分離結果としてｎチャンネル分のスペクトログラムＹ１〜Ｙｎが生成される。以降の処理は、ｎチャンネル分のスペクトログラムＹ１〜Ｙｎに対して別個に行なう。第１チャンネル対応のスペクトログラムＹ１に対する処理がステップＳ６１１〜Ｓ６１３、第ｎチャンネル対応のスペクトログラムＹｎに対する処理がステップＳ６２１〜Ｓ６２３である。なお、ステップＳ６０１の瞬時混合ＩＣＡによる分離処理が終了した時点で、先に図２０を参照して説明したと同様、不要チャンネル（どの音源にも対応しないと判定される出力）を除去する処理を行なう構成としてもよい。
【０２４１】
ステップＳ６１１〜Ｓ６１３の処理は、先に説明した［（３）シフト積み重ねと瞬時混合ＩＣＡを組み合わせた処理によって解く］方式の処理シーケンスである図１４に示すフローのステップＳ１１〜Ｓ１３の処理に対応する処理である。ただし、図１４に示すフローのステップＳ１１の処理はｎチャンネル分のスペクトログラムをｎ×（Ｌ'＋１）に拡張する処理だったのに対し、図２１に示すフローのステップＳ６１１の処理は、１チャンネル分をＬ'＋１チャンネル分に拡張する処理である。また、ステップＳ６１２の残響除去処理は、処理自体は、図１４に示すフローにおけるステップＳ１２の処理と同一の処理となるが、前述の理由により、このステップＳ６１２の処理の効果は音源の分離ではなくて残響除去として実行される。ステップＳ６１３の処理は、Ｌ'＋１チャンネル分の残響除去済みスペクトログラムから所望の一つを選択する処理であり、図１４のステップＳ１３の処理と同様の処理である。
ステップＳ６２１〜Ｓ６２３の処理は、処理対象が異なるチャンネル対応の信号Ｙｎである点を除いてはステップＳ６１１〜Ｓ６１３の処理と同様である。
【０２４２】
全ての出力チャンネル（ただし不要チャンネルは削除して構わない）に対して残響除去と選択とが完了したら、ステップＳ６３１において、残ったスペクトログラムを統合する。例えば縦に積み重ねる処理を実行する。この処理によって、複数のフレームにまたがった成分を取り除く処理、すなわち残響除去処理が実現される。
【０２４３】
［本発明に従った信号分離処理における効果の検証］
上述した本発明の方法により、従来の時間周波数領域ＩＣＡを超える分離性能が出ることを実験で確かめた。以下、この実験結果に基づいて本発明に従った信号分離処理による効果について説明する。
【０２４４】
最初に、実験の条件について説明する。
音データの収録を、図２３に示す環境（オフィスの部屋）で行なった。
マイク数＝４（間隔＝７．５ｃｍ）、音源数＝３であり、音源として以下のＷｅｂページで公開されているものを用いた。
原信号：
ＩＣＡ'９９  ＳＹＮＴＨＥＴＩＣ  ＢＥＮＣＨＭＡＲＫＳ
ｈｔｔｐ：／／ｓｏｕｎｄ．ｍｅｄｉａ．ｍｉｔ．ｅｄｕ／ｉｃａ−ｂｅｎｃｈ／ｓｏｕｒｃｅｓ／
ｓｒｃ１：ｂｅｅｔ．ｗａｖ
ｓｒｃ２：ｂｅｅｔ９．ｗａｖ
ｓｒｃ３：  ｍｉｋｅ．ｗａｖ
なお、収録はそれぞれの音源を単独に鳴らした状態で行ない、後で計算機上で混合している。
【０２４５】
実験は以下の条件で行なった。
サンプリング周波数：１６ｋＨｚ
ＳＴＦＴの窓長：６４，１２８，２５６，５１２，１０２４，（２０４８，４０９６）
ＳＴＦＴのシフト幅：窓長の１／２
窓：短時間フーリエ変換（ＳＴＦＴ）時にサイン窓、逆フーリエ変換（ＦＴ）時に再びサイン窓
η０＝０．５（式）
ループ回数＝２００  ｏｒ  ４００
方式：
（方式１）式［５．２］      （従来法に相当）
（方式２）式［７．１］＆式［７．２］（以降「逆方向畳み込み」）
（方式３）式［９．５］（以降「再ＳＴＦＴ」）
スコア関数：式［７．７］を使用
スコア関数のγの値：
（方式１＆２）γ＝ｓｑｒｔ（Ｍ）      Ｍ：周波数ビンの本数
（方式３）    γ＝ｓｑｒｔ（Ｌ'Ｍ）
フレームタップ：
（方式２）Ｌ'＝４，５，８，１０，１５，１６，２０，２５，３０，３２
（方式３）Ｌ'＝４，８，１６，３２
【０２４６】
評価尺度として、波形ベースのｓｉｇｎａｌ−ｉｎｔｅｒｆｅｒｅｎｃｅ−ｒａｔｉｏ（ＳＩＲ）と周波数ビンベースのＳＩＲとを用いている。以下で、ＳＩＲの計算方法について説明する。
【０２４７】
ｋ番目のチャンネルに対応した分離結果（波形）をｙｋ（ｔ）とし、原信号ｓ１（ｔ）〜ｓＮ（ｔ）の線形結合でｙｋ（ｔ）を近似することを考える（以下に示す式［１０．１］）。
【数１５】
2011215649.tif 000018
【０２４８】
ｓ１（ｔ）〜ｓＮ（ｔ）の係数λ１〜λＮは、式［１０．２］の二乗誤差を最小にすることで求まる。
ｙｋ（ｔ）をｉ番目の音源ｓｉ（ｔ）の推定結果と見なした場合、ＳＩＲはｓｉ（ｔ）とそれ以外の音源とのパワー比として定義される（式［１０．３］）。
出力チャンネル数（＝マイク数）をｎとすると、１つの音源に対してＳＩＲはｎ通り計算されるが、その内の最大値を音源ｉのＳＩＲと定義する（式［１０．４］）。以降の実験結果では、３つの音源のからそれぞれ求めたＳＩＲを、さらに平均している。
周波数ビンベースのＳＩＲは、周波数ビンごとにＳＩＲを計算した後、全ての周波数ビンについて平均を取ることで計算する（式［１０．６］）。
【０２４９】
以下では、実験結果について説明する。以下、実験結果を表として示す。
各表において、
窓長：ＳＴＦＴの窓長、
ｆｒｍ−ｔａｐはフレームタップ数、
ＳＩＲ（ｗａｖｅ）は波形ベースのＳＩＲ、
ＳＩＲ（ｂｉｎ）は周波数ビンベースのＳＩＲ、
を表わす。
【０２５０】
以下、各表において、
（１）方式１（従来法）、２００回ループ
（２）方式２（式［６．１］，［７．１］，［７．２］）、２００回ループ
（３）方式３（式［９．２］，［９．５］）、２００回ループ
（４）方式１（従来法）、４００回ループ
（５）方式２（式［６．１］，［７．１］，［７．２］）、４００回ループ
（６）方式３（式［９．２］，［９．５］）、４００回ループ
これらの実験結果を示す。
【０２５１】
【表１】
2011215649.tif 000019
【０２５２】
【表２】
2011215649.tif 000020
【０２５３】
【表３】
2011215649.tif 000021
【０２５４】
【表４】
2011215649.tif 000022
【０２５５】
【表５】
2011215649.tif 000023
【０２５６】
【表６】
2011215649.tif 000024
【０２５７】
【表７】
2011215649.tif 000025
【０２５８】
【表８】
2011215649.tif 000026
【０２５９】
図２４は、以下の３方式による分離結果についての評価データである。
（１）方式１（従来法）、２００回ループ
（２）方式２（式［６．１］，［７．１］，［７．２］）、２００回ループ
（３）方式３（式［９．２］，［９．５］）、２００回ループ
これらの３方式を実行した場合の結果データに基づくＳＩＲデータをプロットしたものであり、
（ａ）波形ベースのＳＩＲ（ｓｉｇｎａｌ−ｉｎｔｅｒｆｅｒｅｎｃｅ−ｒａｔｉｏ）
（ｂ）周波数ビンベースのＳＩＲ
これらのＳＩＲデータをプロットしたものである。横軸がＳＴＦＴの窓長、縦軸がＳＩＲである。各グラフにおいて、
＊（実線）：方式１、
◆：が方式２、
＋：が方式３、
である。
いくつかの設定において、方式２と方式３は、従来法を上回っているのが確認できる。
【０２６０】
次に、横軸として以下の式で計算されるタイムスパンを用いてプロットした評価データを図２５に示す。
ｔｉｍｅ＿ｓｐａｎ＝｛（ｆｒａｍｅ＿ｔａｐ−１）×ｆｒａｍｅ＿ｓｈｉｆｔ＋ｗｉｎｄｏｗ＿ｌｅｎ｝／ｓｒａｔｅ
なお、
ｆｒａｍｅ＿ｔａｐ：フレームタップ数（＝Ｌ'）
ｗｉｎｄｏｗ＿ｌｅｎ：窓長（一度目のSTFTの切り出し区間長）
ｆｒａｍｅ＿ｓｈｉｆｔ：窓シフト幅（今回は窓長の１／２）
ｓｒａｔｅ：サンプリング周波数（１６ｋＨｚ）
【０２６１】
図２５も、以下の３方式による分離結果についての評価データである。
（１）方式１（従来法）、２００回ループ
（２）方式２（式［６．１］，［７．１］，［７．２］）、２００回ループ
（３）方式３（式［９．２］，［９．５］）、２００回ループ
これらの３方式を実行した場合の結果データに基づくＳＩＲデータをプロットしたものであり、
（ａ）波形ベースのＳＩＲ（ｓｉｇｎａｌ−ｉｎｔｅｒｆｅｒｅｎｃｅ−ｒａｔｉｏ）
（ｂ）周波数ビンベースのＳＩＲ
これらのＳＩＲデータをプロットしたものである。横軸が上述したタイムスパン（Ｔｉｍｅ＿ｓｐａｎ）の窓長、縦軸がＳＩＲである。各グラフにおいて、
＊（実線）：方式１、
◆：が方式２、
＋：が方式３、
である。
【０２６２】
従来は、長いタイムスパンをカバーするためには短時間フーリエ変換（ＳＴＦＴ）の窓長を長くするしかなく、それがＳＩＲの低下を招いていた。それに対し本発明では、短めの窓と複数のフレームタップという組み合わせを用いることで、ＳＩＲを低下させずに同等のタイムスパンをカバーすることができる。
【０２６３】
図２６は、以下の３方式による分離結果についての評価データである。
（４）方式１（従来法）、４００回ループ
（５）方式２（式［６．１］，［７．１］，［７．２］）、４００回ループ
（６）方式３（式［９．２］，［９．５］）、４００回ループ
これらの３方式を実行した場合の結果データに基づくＳＩＲデータをプロットしたものであり、
（ａ）波形ベースのＳＩＲ（ｓｉｇｎａｌ−ｉｎｔｅｒｆｅｒｅｎｃｅ−ｒａｔｉｏ）
（ｂ）周波数ビンベースのＳＩＲ
これらのＳＩＲデータをプロットしたものである。横軸がＳＴＦＴの窓長、縦軸がＳＩＲである。各グラフにおいて、
＊（実線）：方式１、
◆：が方式２、
＋：が方式３、
である。
【０２６４】
次に、分離処理のループ回数を400回に増やして同様の評価実験を行なった。
図２６に示すデータに対応するデータとして、横軸としてタイムスパンを用いてプロットした評価データを図２７に示す。以下の３方式による分離結果についての評価データである。
（４）方式１（従来法）、４００回ループ
（５）方式２（式［６．１］，［７．１］，［７．２］）、４００回ループ
（６）方式３（式［９．２］，［９．５］）、４００回ループ
これらの３方式を実行した場合の結果データに基づくＳＩＲデータをプロットしたものであり、
（ａ）波形ベースのＳＩＲ（ｓｉｇｎａｌ−ｉｎｔｅｒｆｅｒｅｎｃｅ−ｒａｔｉｏ）
（ｂ）周波数ビンベースのＳＩＲ
これらのＳＩＲデータをプロットしたものである。横軸が上述したタイムスパン（Ｔｉｍｅ＿ｓｐａｎ）の窓長、縦軸がＳＩＲである。各グラフにおいて、
＊（実線）：方式１、
◆：が方式２、
＋：が方式３、
である。
【０２６５】
図２６、図２７においても、方式２，方式３は、従来法を上回る設定が存在する。このように、本発明によって、従来の時間周波数領域ＩＣＡが持っていた「窓長と分離性能とのトレードオフ」という課題を回避することが可能である。
【０２６６】
次に、別種のデータについての評価実験について説明する。図２８は、収録環境であるオフィス環境の見取り図である。図に示すようにほぼ７５０ｃｍ×３７５ｃｍの長方形の室内で実験を行なった。なお、図に示すように完全な長方形ではなく、一方は高さ１５３ｃｍのパーティションによって区切られた空間である。部屋の残響時間は０．３秒よりやや短い値である。（以降では、０．２７５秒としてプロットしてある。）
【０２６７】
音源として、以下の３種類の音を用意した。（各信号のスペクトログラムを図２９に示す。）
音源１（ｓｒｃ１）：女性１名の発話（以降、女声またはＦ）
音源２（ｓｒｃ２）：男性１名の発話（以降、男声またはＭ）
音源３（ｓｒｃ３）：以下のＵＲＬで公開されているストリートノイズ（以降、雑踏またはＳ）ｈｔｔｐ：／／ｓｏｕｎｄ．ｍｅｄｉａ．ｍｉｔ．ｅｄｕ／ｉｃａ−ｂｅｎｃｈ／ｓｏｕｒｃｅｓ／ｓｔｒｅｅｔ．ｗａｖ
【０２６８】
図中のｓｐ１〜ｓｐ４の各スピーカーから上記の音をそれぞれ再生し、５ｃｍ間隔で並べた４本のマイク（ｍｉｃ１〜ｍｉｃ４）で収録した。次に、図３０に示す８通りの組み合わせで、スビーカｓｐ１〜ｓｐ４からの音声出力を行なって４本のマイク（ｍｉｃ１〜ｍｉｃ４）によって入力するデータの解析を行なった。女性１名の発話を［Ｆ］、男性１名の発話を［Ｍ］、ストリートノイズを［Ｓ］、音声出力なしを［０］として、
（１）ｓｐ１＝Ｓ，ｓｐ２＝０、ｓｐ３＝Ｆ、ｓｐ４＝Ｍ
（２）ｓｐ１＝Ｓ，ｓｐ２＝０、ｓｐ３＝Ｍ、ｓｐ４＝Ｆ
（３）ｓｐ１＝Ｆ，ｓｐ２＝Ｓ、ｓｐ３＝０、ｓｐ４＝Ｍ
（４）ｓｐ１＝Ｍ，ｓｐ２＝Ｓ、ｓｐ３＝０、ｓｐ４＝Ｍ
（５）ｓｐ１＝０，ｓｐ２＝０、ｓｐ３＝Ｆ、ｓｐ４＝Ｍ
（６）ｓｐ１＝０，ｓｐ２＝０、ｓｐ３＝Ｍ、ｓｐ４＝Ｆ
（７）ｓｐ１＝Ｆ，ｓｐ２＝０、ｓｐ３＝０、ｓｐ４＝Ｍ
（８）ｓｐ１＝Ｍ，ｓｐ２＝０、ｓｐ３＝０、ｓｐ４＝Ｍ
これらの８つのパターンである。
なお、実験では、（１）〜（８）の各パターンについて、観測信号が４秒の場合と８秒の場合とについて実験しているため、観測信号のバリエーションは合計で８×２＝１６通り存在する。
【０２６９】
観測信号の例を図３１に示す。これは、図３０に示すパターン中の［Ｔａｋｅ  Ｎｏ．＝３］に該当する。すなわち、以下の出力パターンである。
（３）ｓｐ１＝Ｆ，ｓｐ２＝Ｓ、ｓｐ３＝０、ｓｐ４＝Ｍ
図３１（ａ）に示すＸ１〜Ｘ４の４枚のスペクトログラムは、図２８に示す４本のマイク（ｍｉｃ１〜ｍｉｃ４）で観測された観測信号である。図３１（ｂ）は周波数ビンごとのＳＩＲである。４枚のスペクトログラムの間で、４つの音源の混ざり具合はほぼ同じであることが分かる。
【０２７０】
音源分離実験は以下の３つの方式について行なった。すなわち、前述の実験から（方式２）を省き、代わりに前述の「（３）シフト積み重ね＋瞬時混合ＩＣＡ」（の中の１番目の方法）を方式４として行なった。
（方式１）式［５．２］      （従来法に相当）
（方式３）式［９．５］（以降「再ＳＴＦＴ」）
（方式４）式［１１．１］＆式［５．２］（以降「シフト積み重ね」）
【０２７１】
実験の条件は以下とおりである。
共通：
サンプリング周波数＝１６ｋＨｚ
サンプルビット数＝１６
観測信号の長さ：４秒および８秒
【０２７２】
（方式１）
ＳＴＦＴの窓長：２５６，５１２，１０２４，２０４８，４０９６，８１９２
ＳＴＦＴのシフト幅：窓長の１／４
窓：短時間フーリエ変換（ＳＴＦＴ）時にハニング窓、逆フーリエ変換（ＦＴ）時には窓なし。
η０＝０．３
ループ回数＝４００
スコア関数のγの値：γ＝ｓｑｒｔ（Ｍ）    Ｍ：周波数ビンの本数
【０２７３】
ただし、観測信号＝４秒、ＳＴＦＴの窓長＝８１９２  の場合のみ、シフト幅は１／８の１０２４を使用した。（１／４シフトではフレーム数が少なくなりすぎるため。）
【０２７４】
（方式３）
ＳＴＦＴ（１回目）の窓長：５１２
ＳＴＦＴ（１回目）のシフト幅：窓長の１／４
窓（１回目）：短時間フーリエ変換（ＳＴＦＴ）時にハニング窓、逆フーリエ変換（ＦＴ）時には窓なし。
η０＝０．３
ループ回数＝４００
スコア関数のγの値：γ＝ｓｑｒｔ（Ｍ（Ｌ'＋１））      Ｍ：周波数ビンの本数
ＳＴＦＴ（２回目）の窓長：Ｌ'＋１＝４，８，１６，３２
ＳＴＦＴ（２回目）のシフト幅：窓長の１／８（端数切り上げ）
窓（２回目）：短時間フーリエ変換（ＳＴＦＴ）時にハミング窓、逆フーリエ変換（ＦＴ）時には窓なし。
【０２７５】
２回目のＳＴＦＴでハニング窓ではなくてハミング窓を使った理由は、タップ数が小さい場合でも両端のサンプルを有効に使用するためである。（ハニング窓は両端が０なので、有効なサンプル数が２個少なくなってしまう。）
【０２７６】
（方式４）
ＳＴＦＴ（１回目）の窓長：５１２
ＳＴＦＴ（１回目）のシフト幅：窓長の１／４
窓（１回目）：短時間フーリエ変換（ＳＴＦＴ）時にハニング窓、逆フーリエ変換（ＦＴ）時には窓なし。
η０＝０．３
ループ回数＝４００
スコア関数のγの値：γ＝ｓｑｒｔ（Ｍ（Ｌ'＋１））      Ｍ：周波数ビンの本数
フレームタップ：Ｌ'＋１＝２，４，８，１２
【０２７７】
図３２と図３３は、図３１の観測信号に対して方式４で処理をした結果である。Ｌ'＝１（すなわち２タップ分）でシフト＆積み重ね（図１１参照）を行なった結果が図３２であり、それを８チャンネルの観測信号として分離をした結果が図３３である。図３３（ａ）が分離結果のスペクトログラム、図３３（ｂ）が周波数ビンごとのＳＩＲである。図３３（ａ）の分離結果のスペクトログラムを見ると、原信号と以下のように対応している。
Ｙ２［１］，Ｙ４［０］：音源１
Ｙ３［０］，Ｙ３［１］：音源２
Ｙ１［０］，Ｙ１［１］：音源３
Ｙ２［０］，Ｙ４［１］：対応なし
【０２７８】
分離度合いを表わす尺度として、周波数ビンごとの改善ＳＩＲの平均を計算した。図３３を例にとると、図３３（ａ）に示す分離結果スペクトログラムの各チャンネル（Ｙ１［０］，〜Ｙ４［１］）について、最も強く現れている音源に対してＳＩＲを計算し、さらに全周波数ビンで平均をとった。例えばＹ２［１］では音源１が最も強く現れているので、音源１へのＳＩＲを計算する。Ｙ４［０］に対しても同様に計算し、両者で値が大きい方を音源１の分離度合いとする。音源２・３についても同様に計算してから３つの間で平均を取ることで、全体の分離度合いとする。この値から観測信号の平均ＳＩＲ、すなわち図３３（ｂ）に示す周波数ビンごとのＳＩＲのプロットを全周波数で平均したものを引くと、改善ＳＩＲが計算される。
【０２７９】
最後に、８回のテイクの間で平均を取ることで、１つの実験パラメーターについての分離度合いを計算する。なお、観測信号が４秒の場合と８秒の場合とは、別々に集計した。集計結果は、図３４と図３５に示す通りである。図３４が観測信号＝４秒の場合、図３５が８秒の場合であり、いずれも縦軸が改善ＳＩＲ、横軸がタイムスパン（対数表示）である。両図の縦の破線は部屋の残響時間であり、０．２７５秒としている。３本の折れ線は、それぞれ従来法（方式１）・再ＳＴＦＴ（方式３）・シフト積み重ね（方式４）に対応する。
【０２８０】
図３４、図３５に示すように、従来法では、ＳＴＦＴの窓長（分析フレーム長）を長くしていっても、ある値（観測信号＝４秒で１０２４、８秒で２０４８）で分離精度がピークに達し、それより長くするとかえって分離精度が悪化する。これは、窓を長くし過ぎると、ＳＴＦＴ結果の時間分解能が落ちるためである。時間分解能の低下は、観測信号が短い場合ほど強く影響するため、観測信号＝４秒の場合の方が８秒の場合よりも短い窓長でピークに達する。一方で、窓が短い場合は、時間分解能は高いものの、フレームをまたがった成分が多くなる（残響が１フレームに収まらなくなる）ため、十分な分離精度が出ない。
【０２８１】
それに対し、本発明である方法３・方法４では、短い窓（この実験では５１２）でＳＴＦＴした結果に対してさらに複数フレーム用いて分離を行なうため、時間分解能の低下を抑える一方で複数フレームにまたがった成分に対しても対応できる。そのため、従来法と同一のタイムスパンで比較した場合はいっそう高い分離精度を達成でき、また、ピークの分離精度で比較した場合は、より長いタイムスパンでより高い分離精度を達成できる。
【０２８２】
以上、特定の実施例を参照しながら、本発明について詳解してきた。しかしながら、本発明の要旨を逸脱しない範囲で当業者が該実施例の修正や代用を成し得ることは自明である。すなわち、例示という形態で本発明を開示してきたのであり、限定的に解釈されるべきではない。本発明の要旨を判断するためには、特許請求の範囲の欄を参酌すべきである。
【０２８３】
なお、明細書中において説明した一連の処理はハードウェア、またはソフトウェア、あるいは両者の複合構成によって実行することが可能である。ソフトウェアによる処理を実行する場合は、処理シーケンスを記録したプログラムを、専用のハードウェアに組み込まれたコンピュータ内のメモリにインストールして実行させるか、あるいは、各種処理が実行可能な汎用コンピュータにプログラムをインストールして実行させることが可能である。
【０２８４】
例えば、プログラムは記録媒体としてのハードディスクやＲＯＭ（Read Only Memory)に予め記録しておくことができる。あるいは、プログラムはフレキシブルディスク、ＣＤ−ＲＯＭ(Compact Disc Read Only Memory)，ＭＯ(Magneto optical)ディスク，ＤＶＤ(Digital Versatile Disc)、磁気ディスク、半導体メモリなどのリムーバブル記録媒体に、一時的あるいは永続的に格納（記録）しておくことができる。このようなリムーバブル記録媒体は、いわゆるパッケージソフトウエアとして提供することができる。
【０２８５】
なお、プログラムは、上述したようなリムーバブル記録媒体からコンピュータにインストールする他、ダウンロードサイトから、コンピュータに無線転送したり、ＬＡＮ(Local Area Network)、インターネットといったネットワークを介して、コンピュータに有線で転送し、コンピュータでは、そのようにして転送されてくるプログラムを受信し、内蔵するハードディスク等の記録媒体にインストールすることができる。
【０２８６】
なお、明細書に記載された各種の処理は、記載に従って時系列に実行されるのみならず、処理を実行する装置の処理能力あるいは必要に応じて並列的にあるいは個別に実行されてもよい。また、本明細書においてシステムとは、複数の装置の論理的集合構成であり、各構成の装置が同一筐体内にあるものには限らない。
【産業上の利用可能性】
【０２８７】
以上、説明したように、本発明の一実施例の構成によれば、複数の音信号が混合した入力信号を時間周波数領域に変換して観測信号スペクトログラムを生成し、観測信号スペクトログラムから信号分離結果を生成する信号分離処理において、観測信号スペクトログラムを時間周波数領域において畳み込み混合された観測信号として解釈し、畳み込み混合を解く独立性分析の実行により信号分離結果を生成する、あるいは、観測信号スペクトログラムに対する時間方向の短時間フーリエ変換（ＳＴＦＴ）によりモジュレーション・スペクトログラムを生成してモジュレーション・スペクトログラムを瞬時混合として解釈し、瞬時混合を解く独立性分析の実行により信号分離結果を生成する構成としている。そのため、直接波、反射波など様々な遅延量を持つ混合された音信号について、遅延量を考慮した高精度な分離処理が実現される。
【符号の説明】
【０２８８】
１１，１２  エントロピー
１３  同時エントロピー
１１１  音源
１２１  マイク
２０１  周波数ビン
２０２  ビン
２０３  ビン
２２１  モジュレーション・スペクトログラム
２２２  多変量確率密度関数
２２３，２２４  エントロピー
４０１  マイク
４０２  ＡＤ変換部
４０３  ＳＴＦＴ部
４０４  信号分離部
４０５  リスケーリング部
４０６  逆ＦＴ部
４０７  後段処理実行部
４０８  畳み込み演算部
４０９  制御部
４５１  マイク
４５２  ＡＤ変換部
４５３  第１ＳＴＦＴ部
４５４  第２ＳＴＦＴ部
４５５  信号分離部
４５６  第１リスケーリング部
４５７  第１逆ＦＴ部
４５８  第２リスケーリング部
４５９  第２逆ＦＴ部
４６０  後段処理実行部
４６１  制御部
【図３】
2011215649.tif 000027
【図４】
2011215649.tif 000028
【図１１】
2011215649.tif 000029
【図１３】
2011215649.tif 000030
【図１４】
2011215649.tif 000031
【図１５】
2011215649.tif 000032
【図１６】
2011215649.tif 000033
【図１７】
2011215649.tif 000034
【図１８】
2011215649.tif 000035
【図１９】
2011215649.tif 000036
【図２０】
2011215649.tif 000037
【図２１】
2011215649.tif 000038
【図２２】
2011215649.tif 000039
【図２３】
2011215649.tif 000040
【図２４】
2011215649.tif 000041
【図２５】
2011215649.tif 000042
【図２６】
2011215649.tif 000043
【図２７】
2011215649.tif 000044
【図３０】
2011215649.tif 000045
【図３４】
2011215649.tif 000046
【図３５】
2011215649.tif 000047
【図１】
2011215649.tif 000048
【図２】
2011215649.tif 000049
【図５】
2011215649.tif 000050
【図６】
2011215649.tif 000051
【図７】
2011215649.tif 000052
【図８】
2011215649.tif 000053
【図９】
2011215649.tif 000054
【図１０】
2011215649.tif 000055
【図１２】
2011215649.tif 000056
【図２８】
2011215649.tif 000057
【図２９】
2011215649.tif 000058
【図３１】
2011215649.tif 000059
【図３２】
2011215649.tif 000060
【図３３】
2011215649.tif 000061