JP6845373B2 - 信号分析装置、信号分析方法及び信号分析プログラム - Google Patents

Description

本発明は、信号分析装置、信号分析方法及び信号分析プログラムに関する。

従来、２個の源信号が混在する状況において、それぞれ異なる位置で取得された複数の観測信号から、各源信号に由来する成分である源信号成分を推定する方法が提案されている。

N.Q.K. Duong, E. Vincent and R. Gribonval, "Under-Determined Reverberant Audio Source Separation Using a Full-Rank Spatial Covariance Model", IEEE Transactions on Audio, Speech, and Language Processing, vol.18, no. 7, pp. 1830−1840, Sep. 2010.

ここで、図２４を用いて、従来の信号分析装置について説明する。図２４は、従来の信号分析装置の構成を示す図である。図２４に示すように、従来の信号分析装置１Ｐは、マルチチャネルウィーナーフィルタ部１０ａを有する。なお、ベクトル、行列又はスカラーであるＡに対し、“＾Ａ”と記載する場合は「“Ａ”の直上に“＾”が記された記号」と同じであるとする。

マルチチャネルウィーナーフィルタ部１０ａは、観測信号ｙ_ｉ（ｋ）（ｉは観測信号のインデックス、ｋはサンプル点のインデックス）と、源信号のサンプル点ごとの分散をモデル化するパラメータである分散パラメータｖ_ｊ（ｋ）（ｊは源信号のインデックス）と、源信号の空間的特性をモデル化するパラメータである空間共分散行列Ｒ_ｊと、に基づいて、源信号成分の推定値を計算する。

具体的には、マルチチャネルウィーナーフィルタ部１０ａは、すべての観測信号を要素とするＩ次元縦ベクトルである観測信号ベクトルｙ（ｋ）を次の（１）式により作成する。なお、以降において、上付きのＴは、転置を表す。

観測信号ベクトルｙ（ｋ）は、２個の源信号に由来する成分（源信号成分）ｘ_１（ｋ），ｘ_２（ｋ）の和として、次の（２）式によってモデル化される。

マルチチャネルウィーナーフィルタ部１０ａは、さらに、観測信号ベクトルｙ（ｋ）にＩ×Ｉの行列であるマルチチャネルウィーナーフィルタＷ_ｊ（ｋ）を適用することにより、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を次の（３）式により計算する。

しかしながら、従来の信号分析装置１Ｐでは、大きい計算量を必要とする逆行列演算をサンプル点ごとに行う必要があった。具体的には、従来の信号分析装置１Ｐでは、ｋごとに逆行列（ｖ_１（ｋ）Ｒ_１＋ｖ_２（ｋ）Ｒ_２）^−１を計算する必要があった。このため、従来の信号分析装置１Ｐでは、源信号に由来する成分である源信号成分の推定値を計算するために大きい計算量が必要であった。

本発明は、上記に鑑みてなされたものであって、源信号成分の推定値を少ない計算量で計算することができる信号分析装置、信号分析方法及び信号分析プログラムを提供することを目的とする。

上述した課題を解決し、目的を達成するために、本発明の信号分析装置は、Ｊ個（Ｊは２以上の整数）の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列になるような行列である同時無相関化行列Ｐまたは／及び、そのエルミート転置Ｐ^Ｈを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するためのパラメータとして得る信号パラメータ分析部を有することを特徴とする。

また、本発明の信号分析装置は、２個の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_１，Ｒ_２に対するｐ_ｉ ^ＨＲ_２ｐ_ｍ＝０（ｉ，ｍは１以上Ｉ以下の互いに異なる任意の整数）を満たす、それぞれ異なるＩ個（Ｉは２以上の整数）の一般化固有ベクトルｐ_１，・・・，ｐ_Ｉを列ベクトルとする行列である同時無相関化行列Ｐまたはそのエルミート転置Ｐ^Ｈを、Ｉ個の位置で取得される観測信号による観測信号ベクトルについて２個の源信号に対応する成分を無相関化するためのパラメータとして得る信号パラメータ分析部を有することを特徴とする。

また、本発明の信号分析装置は、Ｊ個（Ｊは２以上の整数）の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列になるような行列Ｐを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するための同時無相関化行列として、同時無相関化行列のエルミート転置を、観測信号ベクトルに掛けることにより、Ｊ個の源信号に対応する成分が無相関化された同時無相関化済観測信号ベクトルを得る同時無相関化部を有することを特徴とする。

また、本発明の信号分析装置は、Ｊ個（Ｊは２以上の整数）の混在する源信号から得た同時無相関化された共分散行列に基づいて、源信号の空間的特性をモデル化するパラメータである空間共分散行列を求める空間共分散行列計算部、または／及び同時無相関化された共分散行列に基づいて、源信号のサンプル点ごとの分散をモデル化するパラメータである分散パラメータを求める分散パラメータ計算部を有することを特徴とする。

本発明によれば、源信号成分の推定値を少ない計算量で計算することができる。

図１は、第１の実施形態に係る信号分析装置の構成の一例を示す図である。 図２は、第１の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。 図３は、第２の実施形態に係る信号分析装置の構成の一例を示す図である。 図４は、第２の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。 図５は、第３の実施形態に係る信号分析装置の構成の一例を示す図である。 図６は、第３の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。 図７は、従来の方法及び第３の実施形態を用いた音源分離実験の実験条件を示す図である。 図８は、従来の方法及び第３の実施形態を用いた音源分離実験結果（real time factor）を示す図である。 図９は、従来の方法及び第３の実施形態を用いた音源分離実験結果（音源分離性能）を示す図である。 図１０は、第４の実施形態に係る信号分析装置の構成の一例を示す図である。 図１１は、第４の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。 図１２は、第５の実施形態に係る信号分析装置の構成の一例を示す図である。 図１３は、第５の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。 図１４は、従来の方法及び第５の実施形態を用いた音源分離実験の実験条件を示す図である。 図１５は、従来の方法及び第５の実施形態を用いた音源分離実験結果（real time factor）を示す図である。 図１６は、従来の方法及び第５の実施形態を用いた音源分離実験結果（音源分離性能）を示す図である。 図１７は、第６の実施形態に係る信号分析装置の構成の一例を示す図である。 図１８は、第６の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。 図１９は、周波数領域信号分離部の構成の一例を示す図である。 図２０は、周波数領域信号分離処理の処理手順の一例を示すフローチャートである。 図２１は、第８の実施形態に係る信号分析装置の構成の一例を示す図である。 図２２は、第８の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。 図２３は、プログラムが実行されることにより、信号分析装置が実現されるコンピュータの一例を示す図である。 図２４は、従来の信号分析装置の構成を示す図である。

以下に、本願に係る信号分析装置、信号分析方法及び信号分析プログラムの実施形態を図面に基づいて詳細に説明する。また、本発明は、以下に説明する実施形態により限定されるものではない。なお、以下では、ベクトル、行列又はスカラーであるＡに対し、“＾Ａ”と記載する場合は「“Ａ”の直上に“＾”が記された記号」と同じであるとする。ベクトル、行列又はスカラーであるＡに対し、“~Ａ”と記載する場合は「“Ａ”の直上に“~”が記された記号」と同じであるとする。また、以下では、信号分離は、信号分析の下位概念とする。

［第１の実施形態］
まず、第１の実施形態に係る信号分析装置の構成、処理の流れ及び効果について説明する。なお、第１の実施形態においては、２個の源信号が混在する状況において、それぞれ異なる位置でマイクロホンなどのセンサにより取得されたＩ個（Ｉは２以上の整数）の観測信号ｙ_ｉ（ｋ）（ｉ＝１，・・・，Ｉは観測信号のインデックス、ｋはサンプル点のインデックス）と、源信号のサンプル点ごとの分散（パワー）をモデル化するパラメータである分散パラメータｖ_ｊ（ｋ）（ｊ＝１，２は源信号のインデックス）と、源信号の空間的特性をモデル化するパラメータである空間共分散行列Ｒ_ｊと、が信号分析装置に入力されるものとする。

なお、本第１の実施形態における「信号」の例としては、音、脳波、無線信号などが挙げられる。本第１の実施形態における「観測信号」とは、複数のセンサ（センサアレイ）により取得された信号を指す。センサアレイの例としては、マイクロホンアレイ、脳波計・脳磁計、アンテナアレイなどが挙げられる。また、本第１の実施形態における「源信号」とは、観測信号の構成要素である信号を指す。観測信号は２個の源信号の和としてモデル化される。

また、本第１の実施形態における「信号分離」の問題設定は、２個の目的信号（たとえば音声）が混在する状況においてこれらに対応する成分を推定する問題設定（たとえば音源分離）に加えて、目的信号（たとえば音声）とノイズ（たとえばテレビから流れる音楽や雑踏における喧噪）が混在する状況においてこれらに対応する成分（または目的信号に対応する成分のみ）を推定する問題設定（ノイズ除去とも呼ばれる）をも含むものとする。前者の場合における源信号は、目的信号（たとえば音声）である。後者の場合における源信号は、目的信号（たとえば音声）とノイズ（たとえばテレビから流れる音楽や雑踏における喧噪）とである。なお、「雑踏における喧噪」は、多数のノイズが混ざったものであるが、この場合まとめて１つの源信号とみなす。

［信号分析装置の構成及び処理］
まず、図１及び図２を用いて、第１の実施形態に係る信号分析装置の構成及び信号分析処理の処理手順について説明する。図１は、第１の実施形態に係る信号分析装置の構成の一例を示す図である。図２は、第１の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。第１の実施形態に係る信号分析装置１は、例えば、ＲＯＭ（Read Only Memory）、ＲＡＭ（Random Access Memory）、ＣＰＵ（Central Processing Unit）等を含むコンピュータ等に所定のプログラムが読み込まれて、ＣＰＵが所定のプログラムを実行することで実現される。

図１に示すように、信号分析装置１は、観測信号ベクトル作成部１０、信号パラメータ分析部２０、同時無相関化部３０、シングルチャネルウィーナーフィルタ部４０及び同時無相関化逆変換部５０を有する。

まず、信号分析装置１の各部の概要について説明する。観測信号ベクトル作成部１０は、入力された観測信号ｙ_ｉ（ｋ）を取得し（ステップＳ１０）、（４）式に示す、すべての観測信号からなるＩ次元縦ベクトルである観測信号ベクトルをサンプル点ごとに作成する（ステップＳ１１）。

信号パラメータ分析部２０は、入力された空間共分散行列Ｒ_ｊを取得し、一般化固有値問題を解くことにより、一般化固有値行列Λ及び一般化固有ベクトル行列Ｐを計算する（ステップＳ１２）。

同時無相関化部３０は、観測信号ベクトル作成部１０からの観測信号ベクトルｙ（ｋ）と、信号パラメータ分析部２０からの一般化固有ベクトル行列Ｐと、に基づいて、同時無相関化された観測信号ベクトルｙ´（ｋ）を計算する（ステップＳ１３）。具体的には、同時無相関化部３０は、（５）式に示すように、観測信号ベクトルｙ（ｋ）に、一般化固有ベクトル行列Ｐのエルミート転置Ｐ^Ｈを掛けることにより、同時無相関化された観測信号ベクトルｙ´（ｋ）を計算する。

ここで、「´（上付きのプライム）」は、同時無相関化を施された変数であることを表す。なお、同時無相関化とは、観測信号ベクトルｙ（ｋ）に含まれる、２個の源信号に対応する成分（源信号成分）を同時に無相関化する線型変換を指す。

シングルチャネルウィーナーフィルタ部４０は、入力された分散パラメータｖ_ｊ（ｋ）を取得し、信号パラメータ分析部２０からの一般化固有値行列Λと、同時無相関化部３０からの同時無相関化された観測信号ベクトルｙ´（ｋ）と、に基づいて、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を計算する（ステップＳ１４）。具体的には、シングルチャネルウィーナーフィルタ部４０は、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を（６），（７）式により計算する。

ただし、ベクトルαに対し、［α］_ｌはαの第ｌ要素を表し、行列Ａに対し、［Ａ］_ｌｍはＡの（ｌ，ｍ）要素を表す。また、変数の上の∧（ハット）は、推定値であることを表す。

同時無相関化逆変換部５０は、信号パラメータ分析部２０からの一般化固有ベクトル行列Ｐと、シングルチャネルウィーナーフィルタ部４０からの同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）と、に基づいて、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算し、出力する（ステップＳ１５）。具体的には、同時無相関化逆変換部５０は、（８）式に示すように、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）に対し、同時無相関化部３０で掛けた一般化固有ベクトル行列Ｐのエルミート転置Ｐ^Ｈの逆行列（Ｐ^Ｈ）^−１を掛けることにより、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算し、出力する。

［第１の実施形態の処理の背景］
次に、第１の実施形態の処理の背景について説明する。なお、以下では、観測信号が複素数で表される場合について説明する（観測信号が実数で表される場合の説明は、ほぼ同様だから省略する）。第１の実施形態においては、観測信号ベクトルｙ（ｋ）は、２個の源信号に対応する成分（源信号成分）ｘ_１（ｋ），ｘ_２（ｋ）の和として、次に示す（９）式でモデル化される。

源信号成分ｘ_ｊ（ｋ）は、平均が０、共分散行列がＳ_ｊ（ｋ）であるとし、２個の源信号成分ｘ_１（ｋ），ｘ_２（ｋ）は互いに無相関であるとする。２個の源信号成分ｘ_１（ｋ），ｘ_２（ｋ）は、ともに目的信号（たとえば音声）に対応する源信号成分であっても良いし、一方が目的信号（たとえば音声）に対応する源信号成分で他方がノイズ（たとえばテレビからの雑音や雑踏における喧噪）に対応する源信号成分であっても良い。

観測信号ベクトルｙ（ｋ）から源信号成分ｘ_ｊ（ｋ）を推定する方法としては、マルチチャネルウィーナーフィルタに基づく方法がある。この方法では、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を（１０）式により計算する。

すなわち、この方法では、観測信号ベクトルｙ（ｋ）に対し、Ｉ×Ｉの行列であるマルチチャネルウィーナーフィルタＳ_ｊ（ｋ）（Ｓ_１（ｋ）＋Ｓ_２（ｋ））^−１を掛けることにより、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算する。

［従来の信号分析方法］
ここで、従来の信号分析方法（例えば、非特許文献１に記載の信号分析方法）では、（１０）式に基づいて源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算する。その際、源信号成分ｘ_ｊ（ｋ）の共分散行列Ｓ_ｊ（ｋ）は、次の（１１）式によりモデル化される。なお、（１１）式におけるＥは、期待値演算である。

ここで、ｖ_ｊ（ｋ）は、源信号のサンプル点ごと（ｋごと）の分散をモデル化するパラメータである分散パラメータであり、Ｒ_ｊは、源信号の空間的特性をモデル化するパラメータである空間共分散行列である。従来の信号分析方法では、（１０）式に（１１）式を代入して得られる（１２）式を用いて、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算していた。

上述の従来の信号分析方法の問題点として、膨大な計算量を要するため、大規模なデータセットの処理（サンプル点ｋの数が膨大になる）や、リアルタイム処理（短時間での処理が要求される）、計算能力の低い機器などへの適用が困難だった。

実際、従来技術では、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算する際に、逆行列計算（具体的には、（ｖ_１（ｋ）Ｒ_１＋ｖ_２（ｋ）Ｒ_２）^−１の計算）をサンプル点毎に（つまり多数回）行う必要があった。逆行列計算は、Ｉ次行列に対しＯ（Ｉ^３）という大きな計算量を要する演算である。また、サンプル点の数は、典型的には数百〜数千万程にもなる。そのため、従来の信号分析方法では、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算する際に、膨大な計算量を必要としていた。

［第１の実施の形態に係る信号分析方法］
これに対し、第１の実施形態では、２個の源信号成分ｘ_１（ｋ），ｘ_２（ｋ）の同時無相関化に基づいて、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を少ない計算量で計算する方法を提供する。以下、第１の実施形態における考え方について説明する。

各サンプル点ｋに対し信号α（ｋ）の共分散行列Ｅ（α（ｋ）α（ｋ）^Ｈ）が対角行列であるとき、α（ｋ）は無相関であるという。これは、簡単に言えば、「ベクトルα（ｋ）の異なる要素同士は相関を持たない」ということである。いま仮に、（１１）式の源信号成分ｘ_ｊ（ｋ）の共分散行列ｖ_ｊ（ｋ）Ｒ_ｊ（ｊは１または２）がいずれも対角行列であるとすると、（１２）式における逆行列演算（ｖ_１（ｋ）Ｒ_１＋ｖ_２（ｋ）Ｒ_２）^−１は単なる対角要素ごとの逆数演算により実現できるため、少ない計算量で計算できる。しかしながら、通常、源信号成分ｘ_ｊ（ｋ）の異なる要素同士は互いに強い相関を持つため、その共分散行列ｖ_ｊ（ｋ）Ｒ_ｊは対角行列ではない。そこで、本第１の実施形態では、２個の源信号成分ｘ_１（ｋ），ｘ_２（ｋ）を同時無相関化することを考える。

２個の源信号成分ｘ_１（ｋ），ｘ_２（ｋ）の同時無相関化とは、（１３）式及び（１４）式に示すように、ｘ_１（ｋ），ｘ_２（ｋ）に正則行列Ｖを掛ける線型変換を施して、線型変換後の源信号成分ｘ´_１（ｋ），ｘ´_２（ｋ）が同時に無相関になるようにする（すなわち、ｘ´_１（ｋ），ｘ´_２（ｋ）の共分散行列Ｓ´_１（ｋ），Ｓ´_２（ｋ）が同時に対角行列になるようにする）ことを指すものとする。この同時無相関化により、源信号成分ｘ_ｊ（ｋ）の共分散行列ｖ_ｊ（ｋ）Ｒ_ｊはいずれも対角行列に変換される。したがって、同時無相関化により、（１２）式における逆行列演算（ｖ_１（ｋ）Ｒ_１＋ｖ_２（ｋ）Ｒ_２）^−１は、対角行列に対する逆行列演算に変わり、単なる対角要素ごとの逆数演算により実現できるようになるため、少ない計算量で計算できるようになる。

ｘ´_１（ｋ），ｘ´_２（ｋ）の共分散行列Ｓ´_１（ｋ），Ｓ´_２（ｋ）を計算すると、次の（１５）〜（１９）式のようになる。

同時無相関化の条件は、これらが同時に対角行列になること、すなわち以下の（２０）式及び（２１）式が同時に成り立つことである。なお、Λ_１，Λ_２は或る対角行列である。

すなわち、同時無相関化の目標は、空間共分散行列Ｒ_１，Ｒ_２に対し、左から正則行列Ｖを、右からＶのエルミート転置Ｖ^Ｈを、それぞれ掛けたときに、これらが同時に対角行列になるような正則行列Ｖを求めること（同時対角化）である。

このような行列Ｖは、一般化固有値問題に基づいて、以下のようにして求めることができる。（２２）式に示す、空間共分散行列Ｒ_１，Ｒ_２に対する一般化固有値問題を解く。

すると、Ｉ個の一般化固有値λ_１，・・・，λ_Ｉと、λ_１，・・・，λ_Ｉにそれぞれ対応するＩ個の一般化固有ベクトルｐ_１，・・・，ｐ_Ｉと、が得られる。ここで、一般化固有ベクトルｐ_１，・・・，ｐ_Ｉは、次の（２３）式を満たすように取ることができる。

これらの一般化固有ベクトルｐ_１，・・・，ｐ_Ｉを用いて、行列Ｐを次の（２４）式により定めると、これは、次の（２５）式及び（２６）式を満たす（例えば、参考文献１「G. Strang, “Linear Algebra and Its Applications”, Harcourt Brace Jovanovich, 1988.」参照）。

ただし、Ｉは単位行列、Λは一般化固有値λ_１，・・・，λ_Ｉからなる（２７）式の対角行列である。

したがって、行列Ｖを次の（２８）式により定めると、これは同時無相関化の条件である（２０）式及び（２１）式を満たす（Λ_１=Λ、Λ_２=Ｉである）。

以下、（２７）式における一般化固有値からなる対角行列Λを一般化固有値行列と呼び、（２４）式における一般化固有ベクトルからなる行列Ｐを一般化固有ベクトル行列と呼ぶことがある。また、上述のように、Ｐは同時無相関化を実現する行列でもあるから、その意味でＰを同時無相関化行列とも呼ぶ。

信号パラメータ分析部２０は、空間共分散行列Ｒ_１，Ｒ_２に対する一般化固有値λ_１，・・・，λ_Ｉと、λ_１，・・・，λ_Ｉにそれぞれ対応する一般化固有ベクトルｐ_１，・・・，ｐ_Ｉであって（２３）式を満たすものと、を計算し、一般化固有値λ_１，・・・，λ_Ｉを対角要素とする対角行列である一般化固有値行列Λと、一般化固有ベクトルｐ_１，・・・，ｐ_Ｉを列ベクトルとする正則行列である一般化固有ベクトル行列Ｐと、を、（２７）式と（２４）式とにより作成する。

上述のように、２個の源信号成分ｘ_１（ｋ），ｘ_２（ｋ）は、（２４）式で定義される行列Ｐのエルミート転置Ｐ^Ｈを掛ける線型変換により、同時無相関化される。観測信号ベクトルのモデルを示す（９）式の両辺に対し、この線形変換を施すと、次の（２９）式を得る。

ここで、（１８）式、（１９）式、（２５）式、（２６）式及び（２８）式より、ｘ´_１（ｋ），ｘ´_２（ｋ）の共分散行列Ｓ´_１（ｋ），Ｓ´_２（ｋ）は、次の（３０）式及び（３１）式のように、ともに対角行列である。

ただし、Λは（２７）式で定義される対角行列である。（２９）式にマルチチャネルウィーナーフィルタを適用する場合、ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）は、次の（３２）式で与えられる。

（３２）式における共分散行列Ｓ´_１（ｋ），Ｓ´_２（ｋ）は、（３０）式及び（３１）式に示すように対角行列である。したがって、（３２）式は、以下の（３６）式及び（４０）式に示すように非常に簡単になる。このように、同時無相関化により、（１２）式における逆行列演算（ｖ_１（ｋ）Ｒ_１＋ｖ_２（ｋ）Ｒ_２）^−１は、（３３）式と（３７）式のように対角行列に対する逆行列演算（ｖ_１（ｋ）Λ＋ｖ_２（ｋ）Ｉ）^−１に変わり、単なる対角要素ごとの逆数演算により実現できるようになるため、少ない計算量で計算できるようになる。

シングルチャネルウィーナーフィルタ部４０は、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を（３６）式及び（４０）式により計算する。

（３６）式及び（４０）式の処理は、シングルチャネルウィーナーフィルタに基づく処理と解釈することもできる。そのことを説明するために、まず、シングルチャネルウィーナーフィルタとは何であるかを説明する。マルチチャネルウィーナーフィルタが複数（Ｉ個）の観測信号を用いて推定を行うのに対し、シングルチャネルウィーナーフィルタは１個の観測信号のみを用いて推定を行う。シングルチャネルウィーナーフィルタでは、１個の観測信号（スカラー）ｙ（ｋ）を、２個の源信号成分（スカラー）ｘ_１（ｋ），ｘ_２（ｋ）の和として次の（４１）式でモデル化する。

ここで、観測信号が１個しかないことに対応して、ｙ（ｋ），ｘ_１（ｋ），ｘ_２（ｋ）はいずれもスカラーである。源信号成分ｘ_ｊ（ｋ）は、平均が０、分散がσ^２ _ｊ（ｋ）であるとし、２個の源信号成分ｘ_１（ｋ），ｘ_２（ｋ）は互いに無相関であるとする。シングルチャネルウィーナーフィルタに基づく方法では、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）は、次の（４２）式により計算される。

すなわち、この方法では、観測信号ｙ（ｋ）に対し、シングルチャネルウィーナーフィルタσ^２ _ｊ（ｋ）／（σ^２ _１（ｋ）＋σ^２ _２（ｋ））を掛けることにより、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算する。

ここで、（２９）式を要素ごとに見た場合、第ｉ要素より次の（４３）式を得る。

上式に基づくシングルチャネルウィーナーフィルタを用いて、［ｙ´（ｋ）］_ｉから［ｘ´_ｊ（ｋ）］_ｉを推定することを考える。一般に、確率ベクトルαの共分散行列がＡであるとき、αの第ｌ要素［α］_ｌの分散は、Ａの（ｌ，ｌ）要素［Ａ］_ｌｌであることから、［ｘ´_ｊ（ｋ）］_ｉの分散は［Ｓ´_ｊ（ｋ）］_ｉｉである。よって、（３０）式及び（３１）式より、［ｘ´_１（ｋ）］_ｉの分散はｖ_１（ｋ）λ_ｉであり、［ｘ´_２（ｋ）］_ｉの分散はｖ_２（ｋ）である。したがって、シングルチャネルウィーナーフィルタに基づく［ｘ´_ｊ（ｋ）］_ｉの推定値［＾ｘ´_ｊ（ｋ）］_ｉは、次の（４４）式及び（４５）式で与えられる。

これらはそれぞれ、（３６）式及び（４０）式の第ｉ要素に他ならない。

さて、上記のようにして得られた同時無相関化された源信号成分ｘ´_ｊ（ｋ）＝Ｐ^Ｈｘ_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）は、推定が上手く行っていれば、次の（４６）式のように真値に近い値を取ると期待される。

上式の両辺に行列（Ｐ^Ｈ）^−１を掛けると、次の（４７）式を得る。

そこで、次の（４８）式により、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算する。

このように、第１の実施形態によれば、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）の計算において、サンプル点ごとの逆行列演算を行う必要がないため、少ない計算量で信号分離を実現することができる。なお、上式は逆行列（Ｐ^Ｈ）^−１を含むが、これはｋに依存しないため、サンプル点ごとではなく１回計算するだけで良いことに注意する。

従来の信号分析方法における源信号成分ｘ_ｊ（ｋ）の推定値を＾ｘ_ｊ ^ｃｏｎｖ（ｋ）で表し、第１の実施形態における源信号成分ｘ_ｊ（ｋ）の推定値を＾ｘ_ｊ ^ｐｒｏｐ（ｋ）で表す。＾ｘ_ｊ ^ｃｏｎｖ（ｋ）と＾ｘ_ｊ ^ｐｒｏｐ（ｋ）とは理論上一致することが示される。ここで、＾ｘ_ｊ ^ｃｏｎｖ（ｋ）は、（４９）式によって与えられる。

また、（５）式と、（３３）式〜（４０）式と、（４８）式と、により、＾ｘ_ｊ ^ｐｒｏｐ（ｋ）は次の（５０）式及び（５１）式に等しい。

以下では、ｊ＝１の場合について説明する。なお、ｊ＝２の場合については、ｊ＝１の場合と同様であるから、説明を省略する。（２５）式及び（２６）式によって、以下の（５２）式及び（５３）式が導き出される。

これらを、（４９）式に代入すると、次の（５４）式〜（５７）式を得る。

上述のことから、従来の信号分析方法と第１の実施形態とで、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）の計算結果は理論上一致する。ただし、もちろんコンピュータによる計算誤差が生じることはある。

このように、本実施の形態１に係る信号分析装置１では、信号パラメータ分析部２０は、２個の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_１，Ｒ_２に対するｐ_ｉ ^ＨＲ_２ｐ_ｍ＝０（ｉ，ｍは１以上Ｉ以下の互いに異なる任意の整数）を満たす、それぞれ異なるＩ個（Ｉは２以上の整数）の一般化固有ベクトルｐ_１，・・・，ｐ_Ｉを列ベクトルとする行列である同時無相関化行列Ｐまたはそのエルミート転置Ｐ^Ｈを、Ｉ個の位置で取得される観測信号による観測信号ベクトルについて２個の源信号に対応する成分を無相関化するためのパラメータとして得る。信号分析装置１では、この無相関化するためのパラメータを用いて同時無相関化を行うことによって、フィルタリングでは、計算量の少ない要素ごとの掛け算を行うだけで足り、源信号成分の推定値の計算量を低減することができる。

同時無相関化部３０は、信号パラメータ分析部２０が得た同時無相関化行列Ｐのエルミート転置Ｐ^Ｈを、入力されたＩ個の観測信号による観測信号ベクトルに掛けることにより、２個の源信号に対応する成分が無相関化された同時無相関化済観測信号ベクトルを得る。

そして、シングルチャネルウィーナーフィルタ部４０は、同時無相関化部３０が得た同時無相関化済観測信号ベクトルの各要素に、同時無相関化行列Ｐに基づく行列Ｐ^ＨＲ_ｊＰ（ｊは１以上２以下の整数）の対角成分に基づくシングルチャネルウィーナーフィルタを適用して同時無相関化済源信号成分推定値によるベクトルを得る。

さらに、同時無相関化逆変換部５０は、シングルチャネルウィーナーフィルタ部４０が得た同時無相関化済源信号成分推定値によるベクトルに、同時無相関化行列Ｐのエルミート転置Ｐ^Ｈの逆行列（Ｐ^Ｈ）^−１を掛けることにより、源信号成分の推定値を得る。

［第１の実施形態の効果］
このような構成を有する第１の実施形態によれば、上記のように、従来の信号分析方法と比べて、信号分離性能を低下させることなく、より少ない計算量で信号分離を実現することができる。すなわち、第１の実施形態では、観測信号ベクトルに対し、同時無相関化を行ってから、フィルタリングを行い、最後に同時無相関化の逆変換を行うことにより、源信号成分の推定値を計算する。ここで、同時無相関化を実現する行列Ｐは、空間共分散行列に対する一般化固有値問題を解くことにより、求めることができる。

この結果、第１の実施形態では、従来技術と比べて少ない計算量で源信号成分の推定値を計算できる。これは主に、以下の理由による。

第一に、第１の実施の形態では、同時無相関化により、フィルタリングは要素ごとの掛け算という計算量の小さい処理になり、サンプル点ごとの逆行列演算という計算量の大きい処理が不要になる。第二に、第１の実施の形態では、サンプル点によらない空間共分散行列に対する一般化固有値問題を解いて得られる行列Ｐを用いて同時無相関化を行うため、一般化固有値問題はサンプル点ごとではなく１回だけ解けば良い。これは、源信号成分の共分散行列が、サンプル点によらない空間共分散行列のスカラー倍であるとモデル化したことにより可能になったことである。第三に、第１の実施の形態では、同時無相関化の逆変換を行う際には、同時無相関化の行列Ｐ^Ｈの逆行列演算を行う必要があるが、同じ理由により、これはサンプル点ごとではなく１回だけ行えば良い。

さらに、第１の実施の形態では、観測信号ベクトルに対して、同時無相関化を行ってから、フィルタリングを行い、最後に同時無相関化の逆変換を行う、第１の実施形態における処理と、観測信号ベクトルに対して直接フィルタリングを行う従来技術における処理とは、等価である。このため、第１の実施形態によれば、従来技術と比べて信号分離性能を劣化させることなく、従来技術より少ない計算量で源信号成分の推定値を計算することができる。

［第１の実施形態の変形例１］
次に、第１の実施形態の変形例１について説明する。第１の実施形態では、一般化固有ベクトルｐ_１，・・・，ｐ_Ｉは、（２３）式の制約条件を満たすものとした。本第１の実施形態の変形例１では、一般化固有ベクトルｐ_１，・・・，ｐ_Ｉが（２３）式以外の制約条件を満たすものとした場合の信号分析装置について説明する。本第１の実施形態の変形例１に係る信号分析装置における信号パラメータ分析部では、空間共分散行列Ｒ_１，Ｒ_２に対する（５８）式に示す一般化固有値問題を解く。

この結果、Ｉ個の一般化固有値λ_１，・・・，λ_Ｉとλ_１，・・・，λ_Ｉにそれぞれ対応するＩ個の一般化固有ベクトルｐ_１，・・・，ｐ_Ｉと、が得られる。ここで、上に述べた実施の形態１の例では、一般化固有ベクトルｐ_１，・・・，ｐ_Ｉは、次の（５９）式を満たすように取る。

これに対して、第１の実施形態の変形例１では、一般化固有ベクトルｐ_１，・・・，ｐ_Ｉを、次の（６０）式を満たすように取る（各σ_ｌは所定の正数）。

これらの一般化固有ベクトルｐ_１，・・・，ｐ_Ｉを用いて、行列Ｐを次の（６１）式により定めると、これは、（６２）式及び（６３）式を満たす。

ただし、Λ_１，Λ_２は、次の（６４）式及び（６５）式に示す対角行列である。

したがって、行列Ｖを次の（６６）式により定めると、これは、同時無相関化の条件（２０）式及び（２１）式を満たす。ただし、（６４）式及び（６５）式のΛ_１及びΛ_２が、（２０）式及び（２１）式におけるΛ_１及びΛ_２に相当する。

以下、（２７）式における一般化固有値からなる対角行列Λを一般化固有値行列、（６１）式における一般化固有ベクトルからなる行列Ｐを一般化固有ベクトル行列と呼ぶことがある。このようにして求めた行列Ｐを、同時無相関化行列として用いることができる。本第１の実施形態の変形例１に係る信号分析装置におけるシングルチャネルウィーナーフィルタ部では、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を、次の（６７）式及び（６８）式により計算する。

他の各部における処理は、上記の第１の実施形態と同様である。

このように、第１の実施形態の変形例１では、信号パラメータ分析部２０は、２個の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_１，Ｒ_２に対するｐ_ｉ ^ＨＲ_２ｐ_ｍ＝０（ｉ，ｍは１以上Ｉ以下の互いに異なる任意の整数）を満たす、それぞれ異なるＩ個（Ｉは２以上の整数）の一般化固有ベクトルｐ_１，・・・，ｐ_Ｉを列ベクトルとする行列である同時無相関化行列Ｐまたはそのエルミート転置Ｐ^Ｈを、Ｉ個の位置で取得される観測信号による観測信号ベクトルについて２個の源信号に対応する成分を無相関化するためのパラメータとして得る。

［第１の実施形態の変形例２］
次に、第１の実施形態の変形例２について説明する。第１の実施形態の変形例２では、信号パラメータ分析部における処理のもう一つの変形例について説明する。第１の実施形態では、信号パラメータ分析部２０は、固有値分析に基づいて、Ｐ^ＨＲ_１ＰとＰ^ＨＲ_２Ｐとが対角行列になる行列Ｐとそのときの対角行列Λ＝Ｐ^ＨＲ_１Ｐを求めていた。これに対し、本第１の実施形態の変形例２では、信号パラメータ分析部２０は、最適化に基づいて、Λ_１＝Ｐ^ＨＲ_１Ｐと、Λ_２＝Ｐ^ＨＲ_２Ｐとができるだけ、ある対角行列Σ_１，Σ_２に近づくような行列ＰとそのときのΛ_１＝Ｐ^ＨＲ_１Ｐ及びΛ_２＝Ｐ^ＨＲ_２Ｐを求める。

すなわち、第１の実施形態の変形例２では、信号パラメータ分析部２０は、次のような行列ＰとそのときのΛ_１＝Ｐ^ＨＲ_１Ｐ及びΛ_２＝Ｐ^ＨＲ_２Ｐを求める（（６９）式及び（７０）式参照）。

具体的には、次の（７１）式のコスト関数を行列Ｐと対角行列Σ_１とΣ_２に関して最小化することにより、行列Ｐを求めることができる。

ここで、||Ａ||_Ｆは、行列Ａのフロベニウスノルムを表す。例えば、上式のコスト関数の最小化は、参考文献２「A. Yeredor, “Non-Orthogonal Joint Diagonalization in the Least-Squares Sense with Application in Blind Source Separation”, IEEE Transactions on Signal Processing, vol. 50, no. 7, pp. 1545−1553, Jul. 2002.」に記載のアルゴリズムに従って行うことができる。このようにして求めた行列Ｐを、同時無相関化行列として用いることができる。

また、Λ_１＝Ｐ^ＨＲ_１Ｐ及びΛ_２＝Ｐ^ＨＲ_２Ｐに基づき、第１の実施形態の変形例１と同様にして、シングルチャネルウィーナーフィルタ部４０における同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を計算することができる。他の各部における処理は、上に述べた第１の実施形態と同様である。

このように、第１の実施形態の変形例２では、信号パラメータ分析部２０は、Ｊ個（この変形例ではＪは２）の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列に近づくように求めた行列である同時無相関化行列Ｐまたは／及びそのエルミート転置Ｐ^Ｈを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するためのパラメータとして得る。

［第１の実施形態の変形例３］
次に、第１の実施形態の変形例３について説明する。本第１の実施形態の変形例３では、第１の実施形態の変形例２の方法に基づき、源信号が２個の場合に限らず、源信号が一般にＪ個の場合に拡張する場合について説明する。

なお、第１の実施形態の変形例３においては、Ｊ個の源信号が混在する状況において、それぞれ異なる位置で取得されたＩ個の観測信号ｙ_ｉ（ｋ）と、源信号のサンプル点ごとの分散（パワー）をモデル化するパラメータである分散パラメータｖ_ｊ（ｋ）（ｊ＝１，・・・，Ｊは源信号のインデックス）と、源信号の空間的特性をモデル化するパラメータである空間共分散行列Ｒ_ｊと、が信号分析装置１に入力されるものとする。

第１の実施の形態の変形例３では、信号パラメータ分析部２０は、最適化に基づいて、Λ_１＝Ｐ^ＨＲ_１Ｐ，・・・，Λ_Ｊ＝Ｐ^ＨＲ_ＪＰが、できるだけ、ある対角行列Σ_１，・・・，Σ_Ｊに近づくような行列Ｐと、そのときのΛ_１＝Ｐ^ＨＲ_１Ｐ，・・・，Λ_Ｊ＝Ｐ^ＨＲ_ＪＰを求める。

すなわち、第１の実施形態の変形例３では、信号パラメータ分析部２０は、次の（７２）式〜（７３）式に示すような行列Ｐと、そのときのΛ_１＝Ｐ^ＨＲ_１Ｐ，・・・，Λ_Ｊ＝Ｐ^ＨＲ_ＪＰを求める。

具体的には、次の（７４）式のコスト関数を行列Ｐと各対角行列Σ_ｊ（ｊは１以上Ｊ以下の整数）に関して最小化することにより、行列Ｐを求めることができる。

例えば、上式のコスト関数の最小化は、参考文献２に記載のアルゴリズムに従って行うことができる。本第１の実施形態の変形例３では、このようにして求めたＰを、同時無相関化行列として用いることができる。また、Λ_１＝Ｐ^ＨＲ_１Ｐ，・・・，Λ_Ｊ＝Ｐ^ＨＲ_ＪＰに基づき、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を計算することができる。

次に、源信号がＪ個の場合の信号分析装置１の構成と処理について説明する。信号分析装置１は、観測信号ベクトル作成部１０、信号パラメータ分析部２０、同時無相関化部３０、シングルチャネルウィーナーフィルタ部４０及び同時無相関化逆変換部５０を有する。観測信号ベクトル作成部１０は、入力された観測信号ｙ_ｉ（ｋ）を取得し、すべての観測信号からなるＩ次元縦ベクトルである観測信号ベクトルｙ（ｋ）をサンプル点ごとに作成する。観測信号ベクトルｙ（ｋ）は、次の（７５）式で示される。

信号パラメータ分析部２０は、入力された空間共分散行列Ｒ_ｊを取得する。

そして、上述の（７６）〜（７８）式で定義される各行列Λ_ｊ（ｊは１以上Ｊ以下の整数）が対角行列になるような同時無相関化行列Ｐと、そのときのΛ_１，・・・，Λ_Ｊを計算する。

同時無相関化部３０は、観測信号ベクトル作成部１０からの観測信号ベクトルｙ（ｋ）と、信号パラメータ分析部２０からの同時無相関化行列Ｐと、に基づいて、同時無相関化された観測信号ベクトルｙ´（ｋ）を計算する。具体的には、（７９）式に示すように、同時無相関化部３０は、観測信号ベクトルｙ（ｋ）に同時無相関化行列Ｐのエルミート転置Ｐ^Ｈを掛けることにより、同時無相関化された観測信号ベクトルｙ´（ｋ）を計算する。

なお、同時無相関化とは、観測信号ベクトルｙ（ｋ）に含まれる、Ｊ個の源信号に対応する成分（源信号成分）を同時に無相関化する線型変換を指す。

シングルチャネルウィーナーフィルタ部４０は、入力された分散パラメータｖ_ｊ（ｋ）を取得し、信号パラメータ分析部２０からの行列Λ_ｊと、同時無相関化部３０からの同時無相関化された観測信号ベクトルｙ´（ｋ）と、に基づいて、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を計算する。具体的には、シングルチャネルウィーナーフィルタ部４０は、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を次の（８０）式により計算する。

同時無相関化逆変換部５０は、信号パラメータ分析部２０からの同時無相関化行列Ｐと、シングルチャネルウィーナーフィルタ部４０からの同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）と、に基づいて、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算し、出力する。具体的には、同時無相関化逆変換部５０は、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）に対し、（８１）式のように、同時無相関化部３０で掛けた行列Ｐ^Ｈの逆行列（Ｐ^Ｈ）^−１を掛けることにより、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算し、出力する。

このように、第１の実施形態の変形例３では、信号パラメータ分析部２０は、Ｊ個（Ｊは２以上の整数）の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列に近づくように求めた行列である同時無相関化行列Ｐまたは／及びそのエルミート転置Ｐ^Ｈを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するためのパラメータとして得る。

［第１の実施形態及び第１の実施形態の変形例１〜３の効果］
このように、本第１の実施形態及び本第１の実施形態の変形例１〜３によれば、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を少ない計算量で計算することができる。特に、第１の実施形態の固有値分析に基づく構成によれば、従来技術と比べて信号分離性能を劣化させることなく、従来技術より少ない計算量で信号分離を行うことができる。

なお、第１の実施形態の変形例２，３では、対角行列Σ_ｊに近づくように行列Λ_ｊを求める。ただし、行列Λ_ｊは実際には対角行列にならないことがある。この場合には、対角行列に近い行列Λ_ｊ＝行列Ｐ^ＨＲ_ｊＰ」の「対角成分」をλとして使ってシングルチャネルウィーナーフィルタの処理を行う。したがって、第１の実施形態の変形例２，３のシングルチャネルウィーナーフィルタの処理では、行列Λ_ｊの非対角成分は用いず、対角成分のみを用いるので、非対角成分は特定しない構成としても良い。

［第２の実施形態］
次に、第２の実施形態に係る信号分析装置の構成、処理の流れ及び効果について説明する。

第１の実施形態で説明したように、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊが既知の場合には、これらに基づいて、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算できる。一方、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊが未知の場合には、これらを推定する必要がある。

そこで、本第２の実施形態においては、２個の源信号が混在する状況において、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊと、源信号成分ｘ_ｊ（ｋ）と、を反復法により同時に推定する方法について説明する。第２の実施形態においては、それぞれ異なる位置で取得されたＩ個の観測信号ｙ_ｉ（ｋ）が信号分析装置に入力されるものとする。

［信号分析装置の構成及び処理］
図３及び図４を用いて、第２の実施形態に係る信号分析装置の構成及び信号分析処理の処理手順について説明する。図３は、第２の実施形態に係る信号分析装置の構成の一例を示す図である。図４は、第２の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。

図３に示すように、信号分析装置２０１は、初期化部（図示しない）、観測信号ベクトル作成部２１０、信号パラメータ分析部２２０、同時無相関化部２３０、シングルチャネルウィーナーフィルタ部２４０、事後共分散行列更新部２５０、分散パラメータ更新部２６０、空間共分散行列更新部２７０、収束判定部（図示しない）及び同時無相関化逆変換部２８０を有する。

まず、信号分析装置２０１の各部の概要について説明する。図示しない初期化部は、分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ）と、空間共分散行列Ｒ_ｊの初期推定値＾Ｒ_ｊ ^（０）と、を設定する（ステップＳ２０）。これらは、例えば乱数に基づいて設定すればよい。

観測信号ベクトル作成部２１０は、入力された観測信号ｙ_ｉ（ｋ）を取得し（ステップＳ２１）、すべての観測信号からなるＩ次元縦ベクトルである観測信号ベクトルｙ（ｋ）をサンプル点ごとに作成する（ステップＳ２２）。

信号パラメータ分析部２２０は、空間共分散行列更新部２７０からの同時無相関化された空間共分散行列Ｒ´_ｊの推定値＾Ｒ´_ｊ（ただし、例外として、信号パラメータ分析部２２０における最初の処理の際には、入力された空間共分散行列Ｒ_ｊの初期推定値＾Ｒ_ｊ ^（０））に基づいて、一般化固有値問題を解くことにより、一般化固有値行列Λ及び一般化固有ベクトル行列Ｐを更新する（ステップＳ２３）。

同時無相関化部２３０は、観測信号ベクトル作成部２１０からの観測信号ベクトルｙ（ｋ）と、信号パラメータ分析部２２０からの一般化固有ベクトル行列Ｐと、に基づいて、同時無相関化された観測信号ベクトルｙ´（ｋ）を更新する（ステップＳ２４）。

シングルチャネルウィーナーフィルタ部２４０は、信号パラメータ分析部２２０からの一般化固有値行列Λと、同時無相関化部２３０からの同時無相関化された観測信号ベクトルｙ´（ｋ）と、分散パラメータ更新部２６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）（ただし、例外として、シングルチャネルウィーナーフィルタ部２４０における最初の処理の際には、入力された分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ）と、に基づいて、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を更新する（ステップＳ２５）。

事後共分散行列更新部２５０は、信号パラメータ分析部２２０からの一般化固有値行列Λと、分散パラメータ更新部２６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）（ただし、例外として、事後共分散行列更新部２５０における最初の処理の際には、入力された分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ））と、に基づいて、同時無相関化された事後共分散行列Σ´（ｋ）を更新する（ステップＳ２６）。

分散パラメータ更新部２６０は、信号パラメータ分析部２２０からの一般化固有値行列Λと、シングルチャネルウィーナーフィルタ部２４０からの同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）と、事後共分散行列更新部２５０からの同時無相関化された事後共分散行列Σ´（ｋ）と、に基づいて、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）を更新する（ステップＳ２７）。

空間共分散行列更新部２７０は、シングルチャネルウィーナーフィルタ部２４０からの同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）と、事後共分散行列更新部２５０からの同時無相関化された事後共分散行列Σ´（ｋ）と、分散パラメータ更新部２６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、に基づいて、同時無相関化された空間共分散行列Ｒ´_ｊの推定値＾Ｒ´_ｊを更新する（ステップＳ２８）。

図示しない収束判定部は、収束したかどうかの判定を行う（ステップＳ２９）。収束判定部が収束していないと判定した場合（ステップＳ２９：Ｎｏ）、ステップＳ２３の信号パラメータ分析部２２０での処理に戻って、処理が継続される。

一方、収束判定部が収束したと判定した場合（ステップＳ２９：Ｙｅｓ）、同時無相関化逆変換部２８０での以下の処理が行われる。

具体的には、同時無相関化逆変換部２８０は、信号パラメータ分析部２２０からの一般化固有ベクトル行列Ｐと、シングルチャネルウィーナーフィルタ部２４０からの同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）と、に基づいて、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算し、出力する（ステップＳ３０）。

次に、信号分析装置２０１の各部の詳細について説明する。なお、初期化部については上述の通りである。

観測信号ベクトル作成部２１０は、入力された観測信号ｙ_ｉ（ｋ）を取得し、すべての観測信号からなるＩ次元縦ベクトルである観測信号ベクトルｙ（ｋ）を次の（８２）式により作成する。

信号パラメータ分析部２２０は、信号パラメータ分析部２２０における最初の処理の際に、初期化部からの空間共分散行列Ｒ_ｊの初期推定値＾Ｒ_ｊ ^（０）を受け取る。そして、信号パラメータ分析部２２０は、＾Ｒ_１ ^（０），＾Ｒ_２ ^（０）に対する一般化固有値問題を解くことにより、一般化固有値λ_１，・・・，λ_Ｉと、λ_１，・・・，λ_Ｉにそれぞれ対応する一般化固有ベクトルｐ_１，・・・，ｐ_Ｉであって、（８３）式を満たすものを計算する。

これによって、信号パラメータ分析部２２０は、一般化固有値行列Λ及び一般化固有ベクトル行列Ｐを次の（８４）式及び（８５）式により更新する。

信号パラメータ分析部２２０は、信号パラメータ分析部２２０における２回目以降の処理の際に、空間共分散行列更新部２７０からの同時無相関化された空間共分散行列Ｒ´_ｊの推定値＾Ｒ´_ｊを受け取って、＾Ｒ´_１，＾Ｒ´_２に対する一般化固有値問題を解くことにより、一般化固有値λ_１，・・・，λ_Ｉと、λ_１，・・・，λ_Ｉにそれぞれ対応する一般化固有ベクトルｑ_１，・・・，ｑ_Ｉであって、（８６）式を満たすものを計算する。

そして、信号パラメータ分析部２２０は、一般化固有値行列Λ及び行列Ｑを（８７）式及び（８８）式により更新する。

そして、信号パラメータ分析部２２０は、一般化固有ベクトル行列Ｐを次の（８９）式により更新する。

同時無相関化部２３０は、（９０）式に示すように、観測信号ベクトルｙ（ｋ）に一般化固有ベクトル行列Ｐのエルミート転置Ｐ^Ｈを掛けることにより、同時無相関化された観測信号ベクトルｙ´（ｋ）を更新する。

シングルチャネルウィーナーフィルタ部２４０は、同時無相関化された源信号成分ｘ´_ｊ（ｋ）の推定値＾ｘ´_ｊ（ｋ）を（９１）式及び（９２）式により更新する。

ただし、例外として、シングルチャネルウィーナーフィルタ部２４０における最初の処理の際に、上の２式において、分散パラメータ更新部２６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）の代わりに、初期化部からの分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ）を用いる。

事後共分散行列更新部２５０は、同時無相関化された事後共分散行列Σ´（ｋ）を、次の（９３）式により更新する。

ただし、例外として、事後共分散行列更新部２５０における最初の処理の際には、上式において、分散パラメータ更新部２６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）の代わりに、初期化部からの分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ）を用いる。

分散パラメータ更新部２６０は、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）を次の（９４）式により更新する。

空間共分散行列更新部２７０は、同時無相関化された空間共分散行列Ｒ´_ｊの推定値＾Ｒ´_ｊを次の（９５）式により更新する。（９５）式において、Ｋはサンプル点の数である。

同時無相関化逆変換部２８０は、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を次の（９６）式により計算し、出力する。

［第２の実施形態に係る信号分析装置の処理の背景］
［従来の信号分析方法］
次に、信号分析装置２０１の処理の背景について説明する。まず、従来技術について説明する。２個の源信号が混在する状況において、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊと、源信号成分ｘ_ｊ（ｋ）と、を反復法により同時に推定する方法が、非特許文献１に開示されている。

ここで、従来技術における観測信号ベクトルｙ（ｋ）のモデル化について説明する。（９）式における源信号成分ｘ_ｊ（ｋ）は、次の（９７）式のように平均０、共分散行列ｖ_ｊ（ｋ）Ｒ_ｊの複素ガウス分布に従うものとする。

ここで、Ｎ（ｚ；μ，Σ）は、平均がμ、共分散行列がΣの複素ガウス分布（（９８）式）を表す。

また、２個の源信号成分ｘ_１（ｋ），ｘ_２（ｋ）は互いに無相関であると仮定する。観測信号ベクトルｙ（ｋ）が与えられた条件下での源信号成分ｘ_ｊ（ｋ）の事後確率密度は、次の（９９）式で与えられる（例えば、非特許文献１を参照）。

ここで、平均＾ｘ_ｊ（ｋ）及び共分散行列Σ_ｊ（ｋ）は、次の（１００）式及び（１０１）式により与えられる。

（１００）式における平均＾ｘ_ｊ（ｋ）は、（１０）式におけるマルチチャネルウィーナーフィルタに基づく源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）に他ならない（両者に同じ記号を用いているのはそのためである）。したがって、（１０）式におけるマルチチャネルウィーナーフィルタに基づく源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）は、観測信号ベクトルｙ（ｋ）が与えられた条件下での源信号成分ｘ_ｊ（ｋ）の事後確率密度を最大化する点と解釈できることが分かった。（１０１）式はｊに依らないから、（１０２）式と書く。

従来技術では、下記第１の処理〜第４の処理により、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊと、源信号成分ｘ_ｊ（ｋ）と、を同時に推定する。

第１の処理では、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊと、を初期化する。これらは、例えば、乱数に基づいて初期化すればよい。

第２の処理では、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊと、に基づいて、源信号成分ｘ_ｊ（ｋ）を推定する。第１の実施形態では、源信号成分ｘ_ｊ（ｋ）の１つの推定値＾ｘ_ｊ（ｋ）を計算したのに対し、ここでは源信号成分ｘ_ｊ（ｋ）の事後確率密度（（９９）式参照）を更新する。具体的には、（９９）式の平均＾ｘ_ｊ（ｋ）と共分散行列Σ（ｋ）（（１０２）式の通りΣ_１（ｋ）及びΣ_２（ｋ）に等しい。（１０１）式も参照。）を次の（１０３）式及び（１０４）式により更新する。

第３の処理では、源信号成分ｘ_ｊ（ｋ）の事後確率密度の平均＾ｘ_ｊ（ｋ）と共分散行列Σ（ｋ）に基づいて、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊと、を次の（１０５）式及び（１０６）式により更新する。ここで、ｔｒはトレースである。

第４の処理では、第２及び第３の処理を、収束するまで交互に反復し、収束した場合に、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を出力する。

従来技術の問題点として、膨大な計算量を要するため、サンプル点ｋの数が膨大になる大規模なデータセットの処理、短時間での処理が要求されるリアルタイム処理、或いは、計算能力の低い機器などへの適用が困難だった。

実際、従来技術では、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）と事後共分散行列Σ（ｋ）を計算する際に、逆行列計算（具体的には（ｖ_１（ｋ）Ｒ_１＋ｖ_２（ｋ）Ｒ_２）^−１の計算）と行列乗算（具体的には（ｖ_１（ｋ）Ｒ_１）（ｖ_１（ｋ）Ｒ_１＋ｖ_２（ｋ）Ｒ_２）^−１（ｖ_２（ｋ）Ｒ_２）の計算）とを、サンプル点毎に（つまり多数回）行う必要があった。逆行列計算と行列乗算とはともに、Ｉ次行列に対し、Ｏ（Ｉ^３）という大きな計算量を要する演算である。また、サンプル点の数は、典型的には数百〜数千万程にもなる。

このため、従来の信号分析方法では、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）と事後共分散行列Σ（ｋ）を計算する際に、膨大な計算量を必要としていた。さらに、従来では、信号分析装置全体として、上記の第４の処理における反復を行う必要がある（反復回数は典型的には数十回にもなる）ため、計算量はさらに膨大となる。

［第２の実施の形態に係る信号分析方法］
これに対し、第２の実施形態では、２個の源信号成分ｘ_１（ｋ），ｘ_２（ｋ）の同時無相関化に基づいて、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）と事後共分散行列Σ（ｋ）と空間共分散行列Ｒと分散パラメータｖとを少ない計算量で計算する方法を提供する。この結果、信号分析装置２０１全体としても、少ない計算量で信号分析を実現することができる。以下、第２の実施形態における考え方について説明する。

従来技術（非特許文献１参照）の更新則（１０３）式〜（１０６）式において、何回目の反復であるかを明示すると、以下の（１０７）式〜（１１０）式のようになる。なお、ｎ＝１，・・・，Ｎは何回目の反復であるかを表すインデックス、Ｎは反復回数である。

空間共分散行列Ｒ_１，Ｒ_２の推定値＾Ｒ_１ ^{（ｎ−１）}，＾Ｒ_２ ^{（ｎ−１）}の一般化固有値問題に対して、第１の実施形態と同様にして定義される一般化固有値行列及び一般化固有ベクトル行列を、Λ^{（ｎ−１）}及びＰ^{（ｎ−１）}で表す。これらは、（１１１）式及び（１１２）式を満たす。

すなわち、Λ^{（ｎ−１）}及びＰ^{（ｎ−１）}は（１１３）式及び（１１４）式を満たす。Ｐ^{（ｎ−１）}は同時無相関化（同時対角化）を実現する行列であるから、その意味でＰ^{（ｎ−１）}を同時無相関化行列とも呼ぶ。

（１０７）式からサンプル点ごとの逆行列演算を消すために、（１０７）式に（１１３）式と（１１４）式とを代入して整理すると、第１の実施形態と同様に、次の（１１５）式、（１１６）式を得る。

同様に、（１０８）式からサンプル点ごとの逆行列演算を消すために、（１０８）式に、（１１３）式と（１１４）式とを代入すると、次の（１１８）式〜（１２１）式を得る。

これで、サンプル点ごとの逆行列演算を消すことができた。しかしながら、（１２１）式は、まだ、サンプル点ごとの行列の乗算を含んでいる。そこで、Σ^（ｎ）（ｋ）の代わりに、（１２２）式により定義される、Σ´^（ｎ）（ｋ）を更新することを考える。

（１２１）式より、Σ´^（ｎ）（ｋ）は、（１２３）式で表せる。

（１２３）式には、サンプル点ごとの行列の乗算が含まれない。Σ´^（ｎ）（ｋ）は同時無相関化の線型変換を施された事後共分散行列とみなせるから、以下、同時無相関化された事後共分散行列と呼ぶ。

そして、（１１５）式と（１１６）式とにおいても、＾ｘ_ｊ ^（ｎ）（ｋ）の代わりに、（１２４）式で定義される＾ｘ´_ｊ ^（ｎ）（ｋ）を更新することを考える。

（１１５）式及び（１１６）式より、（１２５）式及び（１２６）式となる。

＾ｘ´_ｊ ^（ｎ）（ｋ）は同時無相関化された源信号成分である。

次に、（１０９）式に着目すると、これは、＾ｘ_ｊ ^（ｎ）（ｋ）とΣ^（ｎ）（ｋ）に依存している。上述のように、第２の実施形態では、＾ｘ_ｊ ^（ｎ）（ｋ）とΣ^（ｎ）（ｋ）の代わりに＾ｘ´_ｊ ^（ｎ）（ｋ）とΣ´^（ｎ）（ｋ）とを更新する。そこで、（１０９）式も、＾ｘ´_ｊ ^（ｎ）（ｋ）とΣ´^（ｎ）（ｋ）とを用いて計算できる形に書き換える必要がある。（１２４）式及び（１２２）式より、＾ｘ_ｊ ^（ｎ）（ｋ）とΣ^（ｎ）（ｋ）は（１２７）式及び（１２８）式と表せる。

これらを（１０９）式に代入すると（１２９）式〜（１３３）式を経て（１３４）式が得られる。

（１１０）式についても同様に、＾ｘ´_ｊ ^（ｎ）（ｋ）とΣ´^（ｎ）（ｋ）とを用いて計算できる形に書き換えるために、（１２７）式と（１２８）式とを代入すると、次の（１３５）式及び（１３６）式のようになる。

ここでも、＾Ｒ_ｊ ^（ｎ）の代わりに、（１３７）式で定義される＾Ｒ´_ｊ ^（ｎ）を更新することにする。

これにより、逆行列演算と行列の乗算とをいくつか減らせるため、計算量を削減することができる。また、＾Ｒ´_ｊ ^（ｎ）は、（１３６）式より、（１３８）式のように表せる。

＾Ｒ´_ｊ ^（ｎ）は同時無相関化の線型変換を施された空間共分散行列とみなせるから、以下、同時無相関化された空間共分散行列と呼ぶ。

最後に、一般化固有値行列Λ及び一般化固有ベクトル行列Ｐの更新について説明する。最初の反復においては、空間共分散行列Ｒ_１，Ｒ_２の初期推定値＾Ｒ_１ ^（０），＾Ｒ_２ ^（０）の一般化固有値問題を解くことにより、Λ^（０）及びＰ^（０）を求めればよい。

２回目以降の反復においては、同時無相関化された空間共分散行列の推定値＾Ｒ´_１ ^（ｎ），＾Ｒ´_２ ^（ｎ）が得られている。これを同時無相関化の逆変換により＾Ｒ_１ ^（ｎ），＾Ｒ_２ ^（ｎ）に戻してから、一般化固有値問題を適用することも可能だが、やや無駄な計算が生じる。

そこで、以下のようにすれば、より効率的に計算できる。まず、＾Ｒ´_１ ^（ｎ），＾Ｒ´_２ ^（ｎ）の一般化固有値問題を解いて、一般化固有値を対角成分とする対角行列Λ^（ｎ）と、一般化固有ベクトルを列ベクトルとする行列Ｑ^（ｎ）と、を求める。これらは、（Ｑ^（ｎ））^Ｈ（＾Ｒ´_１ ^（ｎ））Ｑ^（ｎ）＝Λ^（ｎ）及び（Ｑ^（ｎ））^Ｈ（＾Ｒ´_２ ^（ｎ））Ｑ^（ｎ）＝Ｉを満たすものとする。

次に、Ｐ^（ｎ）を（１３９）式により更新する。

このとき、Λ^（ｎ）及びＰ^（ｎ）は、＾Ｒ_１ ^（ｎ），＾Ｒ_２ ^（ｎ）に対する一般化固有値行列及び一般化固有ベクトル行列になる。これは、次の（１４０）式〜（１４７）式のようにして確認できる。

［第２の実施形態の効果］
このように、第２の実施形態によれば、（同時無相関化された）源信号成分の推定値及び（同時無相関化された）事後共分散行列と空間共分散行列と分散パラメータとを従来技術と比べて少ない計算量で計算できるため、従来技術と比べて少ない計算量で信号分離を実現できる。

これは以下の理由による。第一に、第１の実施形態と同様に、同時無相関化に基づくことで、サンプル点ごとの逆行列演算が不要になる。第二に、事後共分散行列Σ^（ｎ）（ｋ）そのものではなく、同時無相関化された事後共分散行列Σ´^（ｎ）（ｋ）を更新するようにすることで、サンプル点ごとの行列の乗算が不要になる。

また、上記のように、従来の信号分析方法の処理と第２の実施形態の処理とは理論上等価であるため、従来の信号分析方法と第２の実施形態とで信号分析結果は理論上一致する（ただし、もちろんコンピュータによる計算誤差が生じることはある）。したがって、第２の実施形態によれば、従来の信号分析方法と比べて、信号分析性能を低下させることなく、より少ない計算量で信号分析を実現することができる。

また、第１の実施形態とは異なり、第２の実施形態によれば、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊが未知の場合でも、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊと、源信号成分ｘ_ｊ（ｋ）と、を同時に推定することができる。

［第３の実施形態］
次に、第３の実施形態に係る信号分析装置の構成、処理の流れ及び効果について説明する。

第３の実施形態においては、２個の音源信号が混在する状況において、それぞれ異なる位置でマイクロホンにより取得されたＩ個（Ｉは２以上の整数）の観測信号ｙ_ｉ（ｔ）が信号分析装置に入力されるものとする。観測信号ｙ_ｉ（ｔ）は、周波数分解され、周波数ごとに第２の実施形態に記載の方法で信号分離がなされ、最後に時間領域に変換される。

［信号分析装置の構成及び処理］
図５及び図６を用いて、第３の実施形態に係る信号分析装置の構成及び信号分析処理の処理手順について説明する。図５は、第３の実施形態に係る信号分析装置の構成の一例を示す図である。図６は、第３の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。

図５に示すように、信号分析装置３０１は、周波数領域変換部３１０、周波数領域信号分離部３２０及び時間領域変換部３３０を有する。

まず、信号分析装置３０１の各部の概要について説明する。周波数領域変換部３１０は、入力された観測信号ｙ_ｉ（ｔ）を受け取り（ステップＳ４０）、短時間フーリエ変換などにより時間周波数領域の観測信号ｙ_ｉ（ｋ，ｆ）を計算する（ステップＳ４１）。ここで、ｋは時間フレームのインデックスであり、ｆは周波数ビンのインデックスである。以下、ｙ_ｉ（ｋ，ｆ）をｙ_ｉ ^ｆ（ｋ）で表す。

周波数領域信号分離部３２０は、ｆごとに第２の実施形態の方法を適用して、周波数領域での信号分離を行う（ステップＳ４２）。すなわち、ｆ＝１に対しては、ｙ_ｉ ^１（ｋ）を入力として、第２の実施形態の方法により源信号成分の推定値＾ｘ_１ ^１（ｋ），＾ｘ_２ ^１（ｋ）を得る。

そして、ｆ＝２に対しては、ｙ_ｉ ^２（ｋ）を入力として、第２の実施形態の方法により、源信号成分の推定値＾ｘ_１ ^２（ｋ），＾ｘ_２ ^２（ｋ）を得る。各ｆについてこれらを繰り返す。

そして、ｆ＝Ｆに対しては、ｙ_ｉ ^Ｆ（ｋ）を入力として、第２の実施形態の方法により、源信号成分の推定値＾ｘ_１ ^Ｆ（ｋ），＾ｘ_２ ^Ｆ（ｋ）を得る。以下、＾ｘ_ｊ ^ｆ（ｋ）を、＾ｘ_ｊ（ｋ，ｆ）で表す。

時間領域変換部３３０は、時間周波数領域の源信号成分の推定値＾ｘ_ｊ（ｋ，ｆ）に基づいて、逆短時間フーリエ変換などにより、時間領域の源信号成分の推定値＾ｘ_ｊ（ｔ）を計算し、出力する（ステップＳ４３）。

このように、第３の実施形態では、時間周波数領域において、周波数ビンｆごとに第２の実施形態に係る方法を適用する。これに対し、非特許文献１には、時間周波数領域において、周波数ビンｆごとに第２の実施形態に記載した従来技術を適用する方法が開示されている。

この従来技術の問題点として、膨大な計算量を要するため、大規模なデータセットの処理や、リアルタイム処理（短時間での処理が要求される）、計算能力の低い機器などへの適用が困難だった。実際、従来技術では、逆行列計算と行列乗算とを、周波数ビンごと、フレームごと、反復ごとに（つまり多数回）行う必要があった。このため、従来の信号分析方法では、膨大な計算量を必要としていた。

例えば、後述の実験では、８秒の観測信号に対して、周波数ビンの数は２５６、フレームの数は２４９、反復数は１０であり、処理に７５０秒を要した。これが８００００秒（およそ２２時間）のデータセットになれば、処理に７５０００００秒（およそ８７日）掛かる計算になる。

［第３の実施形態の効果］
本発明の効果を確認するために、従来の方法及び第３の実施形態を用いた音源分離実験について説明する。

図７は、従来の方法及び第３の実施形態を用いた音源分離実験の実験条件を示す図である。図８は、従来の方法及び第３の実施形態を用いた音源分離実験結果（real time factor）を示す図である。図９は、従来の方法及び第３の実施形態を用いた音源分離実験結果（音源分離性能）を示す図である。

８秒間の英語音声に、図７に示す実験室で測定された残響インパルス応答を畳み込んで残響音声を生成し、複数人の残響音声を加算することで観測信号を生成した。このとき、観測信号に対して、各方法を用いて音源分離を行った場合の性能を図８及び図９に示す。図８は、real time factor（観測信号長を１としたときの、処理に要した時間）である。また、図９は、音源分離性能である信号対歪み比（値が大きいほど音源分離性能が高い）である。

図８及び図９に示すように、本第３の実施形態に係る方法によれば、従来の方法と比べて音源分離性能を低下させることなく、２５０倍以上の計算時間削減が実現できたことが分かる。なお、残響時間が１３０ミリ秒及び４４０ミリ秒の場合には、本第３の実施形態に係る方法による信号対歪み比は、従来の方法と比べてわずかに低下しているが、これはコンピュータによる計算の誤差によるものである。このように、本第３の実施形態に係る方法は従来の方法と理論的には等価だが、実際には計算誤差に起因して性能差が生じることもあるものの、それは図８のようにごく僅かな差に過ぎない。

このように、第３の実施形態によれば、従来の方法と比べて音源分離性能を低下させることなく、計算時間の大きな削減が実現できる。

［第４の実施形態］
次に、第４の実施形態に係る信号分析装置の構成、処理の流れ及び効果について説明する。

第４の実施形態は、第２の実施形態において、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊを、補助関数法を用いて推定する構成に変更したものである。

［信号分析装置の構成及び処理］
図１０及び図１１を用いて、第４の実施形態に係る信号分析装置の構成及び信号分析処理の処理手順について説明する。図１０は、第４の実施形態に係る信号分析装置の構成の一例を示す図である。図１１は、第４の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。

図１０に示すように、信号分析装置４０１は、初期化部（図示しない）、観測信号ベクトル作成部４１０、信号パラメータ分析部４２０、同時無相関化部４３０、分散パラメータ更新部４４０、空間共分散行列更新部４５０、収束判定部（図示しない）を有する。

まず、信号分析装置４０１の各部の概要について説明する。図示しない初期化部は、分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ）と、空間共分散行列Ｒ_ｊの初期推定値＾Ｒ_ｊ ^（０）と、を設定する（ステップＳ５１）。これらは、例えば乱数に基づいて設定すればよい。

観測信号ベクトル作成部４１０は、入力された観測信号ｙ_ｉ（ｋ）を取得し（ステップＳ５２）、すべての観測信号からなるＩ次元縦ベクトルである観測信号ベクトルｙ（ｋ）をサンプル点ごとに作成する（ステップＳ５３）。

信号パラメータ分析部４２０は、空間共分散行列更新部４５０からの空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊに基づいて、一般化固有値問題を解くことにより、一般化固有値行列Λ及び一般化固有ベクトル行列Ｐを更新する(ステップＳ５４)。具体的には、信号パラメータ分析部４２０は、空間共分散行列Ｒ_１，Ｒ_２の推定値＾Ｒ_１，＾Ｒ_２の一般化固有値問題に対して、第１の実施形態と同様にして定義される一般化固有値行列Λ及び一般化固有ベクトル行列Ｐを求める。これらは、（１４８）式及び（１４９）式を満たすものとする。

ただし、例外として、信号パラメータ分析部４２０における最初の処理の際には、信号パラメータ分析部４２０における上述の処理において、空間共分散行列更新部４５０からの空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊの代わりに、入力された空間共分散行列Ｒ_ｊの初期推定値^Ｒ_ｊ ^（０）を用いる。

同時無相関化部４３０は、観測信号ベクトル作成部４１０からの観測信号ベクトルｙ（ｋ）と、信号パラメータ分析部４２０からの一般化固有ベクトル行列Ｐと、に基づいて、同時無相関化された観測信号ベクトルｙ´（ｋ）を（９０）式により更新する（ステップＳ５５）。

分散パラメータ更新部４４０は、同時無相関化部４３０からの同時無相関化された観測信号ベクトルｙ´（ｋ）と、信号パラメータ分析部４２０からの一般化固有値行列Λと、に基づいて、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）を次の（１５０）式及び（１５１）式により更新する（ステップＳ５６）。

ここで、右辺の＾ｖ_１（ｋ）及び＾ｖ_２（ｋ）は、分散パラメータ更新部４４０における前回の処理で得られた分散パラメータの推定値（ただし、例外として、分散パラメータ更新部４４０における最初の処理の際には、初期化部からの分散パラメータの初期推定値）を表す。

空間共分散行列更新部４５０は、同時無相関化部４３０からの同時無相関化された観測信号ベクトルｙ´（ｋ）と、信号パラメータ分析部４２０からの一般化固有ベクトル行列Ｐ及び一般化固有値行列Λと、分散パラメータ更新部４４０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、に基づいて、空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊを更新する（ステップＳ５７）。具体的には、空間共分散行列更新部４５０は、まず次の（１５２）式により行列Ｃ_ｊを更新する。

次に、空間共分散行列更新部４５０は、次の（１５３）式によりベクトルｚ（ｋ）を更新する。

次に、空間共分散行列更新部４５０は、次の（１５４）式により行列Ｄ_ｊを更新する。

次に、空間共分散行列更新部４５０は、次の（１５５）式により空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊを更新する。

図示しない収束判定部は、収束したかどうかの判定を行う（ステップＳ５８）。収束判定部が収束していないと判定した場合（ステップＳ５８：Ｎｏ）、信号パラメータ分析部４２０での処理（ステップＳ５４）に戻って、処理が継続される。

一方、収束判定部が収束したと判定した場合（ステップＳ５８：Ｙｅｓ）、分散パラメータ更新部４４０が分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）を出力し、空間共分散行列更新部４５０が空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊを出力する（ステップＳ５９）。

[第４の実施形態の背景]
第４の実施形態の背景について説明する。補助関数法に基づけば、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊは、次のように推定できる。分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊの推定値はあらかじめ初期化され、以下に述べる反復処理により更新される。分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）は、次の（１５６）式により更新される。

空間共分散行列Ｒ_ｊの更新においては、まず（１５７）式及び（１５８）式により行列Ｃ_ｊ及び行列Ｄ_ｊが更新され、次に、（１５９）式により空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊが更新される。

ここで、半正定値エルミート行列Ａの平方根Ａ^１／２は、（１６０）式をＡの固有値分解とすると、（１６１）式により定義される。ここで、ＵはＡの固有ベクトルからなるユニタリ行列、ΛはＡの固有値からなる対角行列、Σ^１／２はΣの各対角要素[Σ]_ｉｉをその平方根（[Σ]_ｉｉ）^１／２で置き換えることにより得られる対角行列である。

上記の分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊの更新は、収束するまで反復される。

上記の分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊの更新則は、参考文献３「Kazuyoshi Yoshii，Ryota Tomioka，Daichi Mochihashi，and Masataka Goto，“Infnite Positive Semidefinite Tensor Factorization for Source Separation of Mixture Signals”，Proceedings of the 30th International Conference on Machine Learning，2013．」、参考文献４「Hiroshi Sawada，Hirokazu Kameoka，Shoko Araki，and Naonori Ueda，“Multichannel Extensions of Non−Negative Matrix Factorization With Complex−Valued Data”，IEEE TRANSACTIONS ON AUDIO， SPEECH， AND LANGUAGE PROCESSING，VOL．21，No.5，MAY 2013．」に開示されている方法と同様の方法により導出することができる。

上記の方法の問題点として、（１５６）式〜（１５８）式の各処理において、サンプル点ごとの逆行列演算を行う必要があるため、膨大な計算量を必要とする。

これに対し、本第４の実施形態では、同時無相関化（同時対角化）に基づいて、サンプル点ごとの逆行列演算を回避することにより、計算量を削減することができる。実際、（１４８）式及び（１４９）式を＾Ｒ_１、＾Ｒ_２について解いたものを（１５６）式〜（１５８）式に代入することにより、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊの推定値の更新則が以下のように得られる。

分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）の更新則は、次の（１６２）式及び（１６３）式のようになる。

ここで、ｙ´（ｋ）は、（１６４）式の同時無相関化済観測信号ベクトルであり、Λ_１は一般化固有値行列Λであり、Λ_２は単位行列Ｉである。

空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊの更新則は次の（１６５）式及び（１６６）式のようになる。まず、行列Ｃ_ｊに関しては、次式のようになる。

行列Ｄ_ｊに関しては次の（１６７）式及び（１６８）式のようになる。

ただし、ベクトルｚ（ｋ）は、次の（１６９）式〜（１７１）式により定義される。

［第４の実施形態の効果］
このように、同時無相関化（同時対角化）により、サンプル点ごとの逆行列演算を不要にし、計算量を削減することができる。具体的には、同時無相関化（同時対角化）に基づき、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）の更新式を（１６２）式のように対角行列Λ_ｊ＝Ｐ^ＨＲ_ｊＰを用いて表すことができ、その結果、（１６３）式のように逆行列演算を不要にすることができる。空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊの更新則についても同様である。

なお、（１５５）式の処理では、（１５９）式の処理と同様、逆行列演算を行う必要があるが、これはサンプル点ごとではなく各反復にて１回だけ行えば良いため、その計算量は比較的小さい。同様に、（１５５）式の処理では、（１５９）式の処理と同様、行列の平方根の計算のために固有値分解を行う必要があるが、これはサンプル点ごとではなく各反復にて２回だけ行えば良いため、その計算量は比較的小さい。また、（１６６）式、（１６８）式、及び（１５５）式の処理では、逆行列演算と同様に計算量の大きい演算である行列の乗算を含んでいるが、これもサンプル点ごとに行う必要はないため、その計算量は比較的小さい。

［第４の実施形態の変形例］
なお、本実施形態においては、源信号が２個の場合について説明したが、第１の実施形態の変形例３と同様の考え方に基づけば、源信号が２個の場合に限らず、源信号が一般にＪ個の場合に拡張することができる。

上記のように、第１〜第４の実施形態では、信号パラメータ分析部が、Ｊ個（Ｊは２以上の整数）の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列になるような行列である同時無相関化行列Ｐまたは／及び、そのエルミート転置Ｐ^Ｈを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するためのパラメータとして得る。信号分析装置１，２０１，３０１，４０１では、この無相関化するためのパラメータを用いて同時無相関化を行うことによって、フィルタリングでは、計算量の少ない要素ごとの掛け算を行うだけで足り、源信号成分の推定値の計算量を低減することができる。

［第５の実施形態］
次に、第５の実施形態に係る信号分析装置の構成、処理の流れ及び効果について説明する。第５の実施形態は、第１の実施形態の変形例３のバリエーションである。第１の実施形態の変形例３では、分散パラメータｖ_ｊ（ｋ）と空間共分散行列Ｒ_ｊが信号分析装置に入力されるとの前提条件がある。これに対し、第５の実施形態では、分散パラメータｖ_ｊ（ｋ）と空間共分散行列Ｒ_ｊとが信号分析装置に入力される必要はなく、分散パラメータｖ_ｊ（ｋ）と、空間共分散行列Ｒ_ｊをモデル化するパラメータである同時無相関化行列Ｐ及び同時無相関化された空間共分散行列Λ_ｊと、源信号成分ｘ_ｊ（ｋ）との同時推定を実現する。

本第５の実施形態では、一般にＪ個の源信号が混在し、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊが未知の状況で、分散パラメータｖ_ｊ（ｋ）と、空間共分散行列Ｒ_ｊをモデル化するパラメータである同時無相関化行列Ｐ及び同時無相関化された空間共分散行列Λ_ｊと、源信号成分ｘ_ｊ（ｋ）と、を同時に推定することにより、信号分離を実現する方法について説明する。本第５の実施形態では、それぞれ異なる位置で取得されたＩ個の観測信号ｙ_ｉ（ｋ）が信号分析装置に入力されるものとする。

［信号分析装置の構成及び処理］
図１２及び図１３を用いて、第５の実施形態に係る信号分析装置の構成及び信号分析処理の処理手順について説明する。図１２は、第５の実施形態に係る信号分析装置の構成の一例を示す図である。図１３は、第５の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。

図１２に示すように、信号分析装置５０１は、初期化部（図示しない）、観測信号ベクトル作成部５１０、同時無相関化部５３０、シングルチャネルウィーナーフィルタ部５４０、事後共分散行列更新部５５０、分散パラメータ更新部５６０、空間共分散行列更新部５７０、信号パラメータ分析部５２０、収束判定部（図示しない）及び同時無相関化逆変換部５８０を有する。

まず、信号分析装置５０１の各部の概要について説明する。図示しない初期化部は、分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ）と、同時無相関化行列Ｐの初期推定値＾Ｐ^（０）と、同時無相関化された空間共分散行列Λ_ｊの初期推定値＾Λ_ｊ ^（０）と、を設定する（ステップＳ６１）。これらは、例えば乱数に基づいて設定すればよい。

観測信号ベクトル作成部５１０は、入力されたそれぞれ異なる位置で取得されたＩ個の観測信号ｙ_ｉ（ｋ）を取得し（ステップＳ６２）、Ｉ次元縦ベクトルである観測信号ベクトルｙ（ｋ）をサンプル点ごとに（８２）式により作成する（ステップＳ６３）。

同時無相関化部５３０は、観測信号ベクトル作成部５１０からの観測信号ベクトルｙ（ｋ）と、信号パラメータ分析部５２０からの同時無相関化行列Ｐの推定値＾Ｐ（ただし例外として、同時無相関化部５３０における最初の処理の際には、初期化部からの同時無相関化行列Ｐの初期推定値＾Ｐ^（０））と、に基づいて、同時無相関化された観測信号ベクトルｙ´（ｋ）を次の（１７２）式により更新する（ステップＳ６４）。

ただし、同時無相関化部５３０における最初の処理の際には、（１７２）式の右辺において、＾Ｐを＾Ｐ^（０）で置き換える。

シングルチャネルウィーナーフィルタ部５４０は、同時無相関化部５３０からの同時無相関化された観測信号ベクトルｙ´（ｋ）と、分散パラメータ更新部５６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）（ただし、例外として、シングルチャネルウィーナーフィルタ部５４０における最初の処理の際には、初期化部からの分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ））と、空間共分散行列更新部５７０からの同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊ（ただし、例外として、シングルチャネルウィーナーフィルタ部５４０における最初の処理の際には、初期化部からの同時無相関化された空間共分散行列Λ_ｊの初期推定値＾Λ_ｊ ^（０））と、に基づいて、同時無相関化された源信号成分ｘ_ｊ´（ｋ）の推定値＾ｘ_ｊ´（ｋ）を次の（１７３）式により更新する（ステップＳ６５）。

ただし、シングルチャネルウィーナーフィルタ部５４０における最初の処理の際には、（１７３）式の右辺において、＾Λ_ｊ（ｊ＝１，・・・，Ｊ）を＾Λ_ｊ ^（０）（ｊ＝１，・・・，Ｊ）に、＾ｖ_ｊ（ｋ）（ｊ＝１，・・・，Ｊ；ｋ＝１，・・・，Ｋ）を＾ｖ_ｊ ^（０）（ｋ）（ｊ＝１，・・・，Ｊ；ｋ＝１，・・・，Ｋ）に、それぞれ置き換える。

事後共分散行列更新部５５０は、分散パラメータ更新部５６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）（ただし、例外として、事後共分散行列更新部５５０における最初の処理の際には、初期化部からの分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ））と、空間共分散行列更新部５７０からの同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊ（ただし、例外として、事後共分散行列更新部５５０における最初の処理の際には、初期化部からの同時無相関化された空間共分散行列Λ_ｊの初期推定値＾Λ_ｊ ^（０））と、に基づいて、同時無相関化された事後共分散行列Σ_ｊ´（ｋ）を次の（１７４）式により更新する（ステップＳ６６）。

ただし、事後共分散行列更新部５５０における最初の処理の際には、（１７４）式の右辺において、＾Λ_ｊを＾Λ_ｊ ^（０）に、＾ｖ_ｊ（ｋ）を＾ｖ_ｊ ^（０）（ｋ）に、それぞれ置き換える。

分散パラメータ更新部５６０は、シングルチャネルウィーナーフィルタ部５４０からの同時無相関化された源信号成分ｘ_ｊ´（ｋ）の推定値＾ｘ_ｊ´（ｋ）と、事後共分散行列更新部５５０からの同時無相関化された事後共分散行列Σ_ｊ´（ｋ）と、空間共分散行列更新部５７０からの同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊ（ただし、例外として、分散パラメータ更新部５６０における最初の処理の際には、初期化部からの同時無相関化された空間共分散行列Λ_ｊの初期推定値＾Λ_ｊ ^（０））と、に基づいて、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）を次の（１７５）式により更新する（ステップＳ６７）。

ただし、分散パラメータ更新部５６０における最初の処理の際には、（１７５）式の右辺において、＾Λ_ｊを＾Λ_ｊ ^（０）に置き換える。

空間共分散行列更新部５７０は、シングルチャネルウィーナーフィルタ部５４０からの同時無相関化された源信号成分ｘ_ｊ´（ｋ）の推定値＾ｘ_ｊ´（ｋ）と、事後共分散行列更新部５５０からの同時無相関化された事後共分散行列Σ_ｊ´（ｋ）と、分散パラメータ更新部５６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、に基づいて、同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊを次の（１７６）式により更新する（ステップＳ６８）。

信号パラメータ分析部５２０は、シングルチャネルウィーナーフィルタ部５４０からの同時無相関化された源信号成分ｘ_ｊ´（ｋ）の推定値＾ｘ_ｊ´（ｋ）と、事後共分散行列更新部５５０からの同時無相関化された事後共分散行列Σ_ｊ´（ｋ）と、分散パラメータ更新部５６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、空間共分散行列更新部５７０からの同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊと、（信号パラメータ分析部５２０における同時無相関化行列Ｐの推定値＾Ｐの最初の更新の際には、更に、初期化部からの同時無相関化行列Ｐの初期推定値＾Ｐ^（０）と、）を受け取って、次の（１７７）式により同時無相関化行列Ｐの推定値＾Ｐを更新する（ステップＳ６９）。ただし、ｅ_ｉは単位行列の第ｉ列を表し、行列Ａに対して［Ａ］_ｉは、Ａの第ｉ列を表す。

ただし、信号パラメータ分析部５２０における同時無相関化行列Ｐの推定値＾Ｐの最初の更新の際には、（１７７）式の右辺において、＾Ｐを＾Ｐ^（０）で置き換える。なお、信号パラメータ分析部５２０は、（１７７）式を複数回連続して適用することにより同時無相関化行列Ｐの推定値＾Ｐを更新してもよい。

図示しない収束判定部は、収束したかどうかの判定を行う（ステップＳ７０）。収束判定部が収束していないと判定した場合（ステップＳ７０：Ｎｏ）、同時無相関化部５３０によるステップＳ６４の処理に戻って、処理が継続される。

一方、収束判定部が収束したと判定した場合（ステップＳ７０：Ｙｅｓ）、同時無相関化部５３０における処理（ステップＳ７１）と、シングルチャネルウィーナーフィルタ部５４０における処理（ステップＳ７２）と、同時無相関化逆変換部５８０における処理（ステップＳ７３）と、が行われる。同時無相関化部５３０及びシングルチャネルウィーナーフィルタ部５４０における処理については上述の通りであるから、同時無相関化逆変換部５８０での処理について説明する。同時無相関化逆変換部５８０は、信号パラメータ分析部５２０からの同時無相関化行列Ｐの推定値＾Ｐと、シングルチャネルウィーナーフィルタ部５４０からの同時無相関化された源信号成分ｘ_ｊ´（ｋ）の推定値＾ｘ_ｊ´（ｋ）と、に基づいて、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を次の（１７８）式により計算し、出力する（ステップＳ７３）。

［信号分析装置の処理の背景］
次に、本第５の実施形態に係る信号分析装置５０１の処理の背景について説明する。本第５の実施形態においては、観測信号ベクトルｙ（ｋ）は、Ｊ個の源信号成分ｘ_１（ｋ），・・・，ｘ_Ｊ（ｋ）の和として、次の（１７９）式でモデル化される。

源信号成分ｘ_ｊ（ｋ）は、（９７）式のように平均０、共分散行列Ｓ_ｊ（ｋ）＝ｖ_ｊ（ｋ）Ｒ_ｊの複素ガウス分布に従うとする。ここで、ｖ_ｊ（ｋ）はｊ番目の源信号のサンプル点ごとの分散をモデル化する分散パラメータであり、Ｒ_ｊはｊ番目の源信号の空間的特性をモデル化するパラメータである空間共分散行列である。分散パラメータｖ_ｊ（ｋ）は正であり、空間共分散行列Ｒ_ｊは正定値エルミート行列であると仮定する。Ｊ個の源信号成分ｘ_１（ｋ），・・・，ｘ_Ｊ（ｋ）は互いに無相関であるとする。

観測信号ベクトルｙ（ｋ）から源信号成分ｘ_ｊ（ｋ）を推定する方法としては、マルチチャネルウィーナーフィルタに基づく方法がある。この方法では、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を次の（１８０）式により計算する。

本第５の実施形態のように、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊが未知の状況で（１８０）式により源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算する（すなわち信号分離を実現する）ためには、これらのパラメータを推定する必要がある。

［従来の信号分析方法］
分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊが未知の状況で、分散パラメータｖ_ｊ（ｋ）と、空間共分散行列Ｒ_ｊと、源信号成分ｘ_ｊ（ｋ）と、を同時に推定することにより、信号分離を実現する方法が、非特許文献１に開示されている。

従来技術では、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊを最尤推定する。この最尤推定は、下記の（１８１）式の最適化問題として定式化される。

ここで、Ψはパラメータの全体を表し、分散パラメータｖ_ｊ（ｋ）（ｊ＝１，・・・，Ｊ；ｋ＝１，・・・，Ｋ）及び空間共分散行列Ｒ_ｊ（ｊ＝１，・・・，Ｊ）からなる。Ｌ_１（Ψ）は、次の（１８２）式の対数尤度関数を表す。

下の（１８３）式は、行列Ａが正定値エルミート行列であることを表す。

従来技術では、ＥＭアルゴリズムに基づいて最尤推定を行う。このＥＭアルゴリズムでは、ＥステップとＭステップとが交互に反復される。

Ｅステップでは、パラメータΨの現在の推定値＾Ψに基づく、源信号成分ｘ_ｊ（ｋ）の事後確率ｐ（ｘ_ｊ（ｋ）｜ｙ（ｋ），＾Ψ）を計算する。この事後確率は、次の（１８４）式のように複素ガウス分布であることが示される。

ただし、平均＾ｘ_ｊ（ｋ）及び共分散行列Σ_ｊ（ｋ）は、パラメータの現在の推定値＾Ψを用いて次の（１８５）式及び（１８６）式のように表される（＾ｖ_ｊ（ｋ）は分散パラメータｖ_ｊ（ｋ）の現在の推定値であり、＾Ｒ_ｊは空間共分散行列Ｒ_ｊの現在の推定値である）。

したがって、Ｅステップでは、平均＾ｘ_ｊ（ｋ）及び共分散行列Σ_ｊ（ｋ）を、（１８５）式及び（１８６）式により更新する。なお、（１８５）式は、（１８０）式のマルチチャネルウィーナーフィルタによる源信号成分ｘ_ｊ（ｋ）の推定値に他ならない。以下、＾ｘ_ｊ（ｋ）を源信号成分ｘ_ｊ（ｋ）の推定値と呼び、Σ_ｊ（ｋ）を事後共分散行列と呼ぶ。

Ｍステップでは、Ｅステップで更新された源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）及び事後共分散行列Σ_ｊ（ｋ）に基づいて、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊを更新する。これらは、下記の更新則（（１８７）式及び（１８８）式）により更新される。

しかしながら、従来技術では、（１８５）式の源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）の更新と、（１８６）式の事後共分散行列Σ_ｊ（ｋ）の更新と、において、逆行列計算と行列乗算とをサンプル点ごとに行う必要があった。具体的には、（１８９）式に示す逆行列計算（反復ごとサンプル点ごとに１回）と、（１９０）式に示す行列乗算（反復ごとサンプル点ごとに２Ｊ回）とを行う必要があった。

そのため、従来技術では膨大な計算量を要していたため、リアルタイムでの信号分離や、センサネットワーク等における分散処理など、サンプル点の数Ｋが大きい（観測信号長が長い）場合や計算資源が限られている状況への適用が困難であった。

［第５の実施形態に係る信号分析方法］
これに対し、本第５の実施形態では、Ｊ個の源信号成分ｘ_１（ｋ），・・・，ｘ_Ｊ（ｋ）の同時無相関化に基づいて、Ｅステップにおける源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）及び事後共分散行列Σ_ｊ（ｋ）の更新を少ない計算量で実現することにより、少ない計算量での信号分離を可能にする。以下、本第５の実施形態における考え方について説明する。

いま仮に、（１８５）式及び（１８６）式における源信号成分ｘ_ｊ（ｋ）の共分散行列ｖ_ｊ（ｋ）Ｒ_ｊがいずれも対角行列であるとすると、（１８５）式及び（１８６）式における逆行列演算及び行列乗算は単なる対角要素ごとの逆数演算及び乗算に帰着するため、少ない計算量で計算できる。しかしながら、通常、ベクトルｘ_ｊ（ｋ）の異なる要素（つまり異なる位置で観測されたｊ番目の源信号）は互いに相関を持つため、その共分散行列ｖ_ｊ（ｋ）Ｒ_ｊの非対角要素は０ではない、すなわち共分散行列ｖ_ｊ（ｋ）Ｒ_ｊは対角行列ではない。

そこで本第５の実施形態では、Ｊ個の源信号成分ｘ_１（ｋ），・・・，ｘ_Ｊ（ｋ）を、正則行列Ｐのエルミート転置Ｐ^Ｈによる次の（１９１）式の変換により同時無相関化することを考える。

変換後のｘ_ｊ´（ｋ）の共分散行列は、次の（１９２）式及び（１９３）式で与えられる。

変換後のｘ_ｊ´（ｋ）（ｊ＝１，・・・，Ｊ）がすべて無相関になる条件、すなわち、（１９３）式で与えられるこれらの共分散行列がすべて対角行列になる条件は、次の（１９４）式〜（１９６）式である。

既に述べたように、信号数がＪ＝２の場合には、一般化固有値問題を解くことにより、（１９４）式〜（１９６）式を満たす正則行列Ｐ及び対角行列Λ_１，・・・，Λ_Ｊを得ることができる。これに対し、本第５の実施形態では、Ｊ＝２の場合に限らず一般の信号数Ｊの場合に適用可能なアプローチとして、未知のＪ個の空間共分散行列Ｒ_１，・・・，Ｒ_Ｊが、ある正則行列Ｐにより同時対角化されるという制約条件、すなわち、ある正則行列Ｐ及びあるＪ個の対角行列Λ_１，・・・，Λ_Ｊに対し、（１９４）式〜（１９６）式が成立するという制約条件の下、正則行列Ｐ及び対角行列Λ_１，・・・，Λ_Ｊを最尤推定する。上述の制約条件は、ある正則行列Ｐ及びあるＪ個の対角行列Λ_１，・・・，Λ_Ｊに対し、次の各式（（１９７）式〜（１９９）式）が成立することと同値である。

本第５の実施形態では、空間共分散行列Ｒ_１，・・・，Ｒ_Ｊは、パラメータである正則行列Ｐとパラメータである対角行列Λ_１，・・・，Λ_Ｊとを用いて、（１９７）式〜（１９９）式のように表される（モデル化される）と仮定する。この場合、空間共分散行列Ｒ_１，・・・，Ｒ_Ｊを推定することは、行列Ｐ（同時無相関化を実現する行列であるから、以下、同時無相関化行列と呼ぶ）及び対角行列Λ_１，・・・，Λ_Ｊ（同時無相関化行列Ｐを用いて（１９４）式〜（１９６）式のように空間共分散行列Ｒ_１，・・・，Ｒ_Ｊを変換して得られる行列であるから、以下、同時無相関化された空間共分散行列と呼ぶ）を推定することに帰着する。

そこで、本第５の実施形態では、同時無相関化行列Ｐと、同時無相関化された空間共分散行列Λ_ｊと、分散パラメータｖ_ｊ（ｋ）と、を最尤法により同時推定する。これにより、一般の信号数Ｊに対しても、源信号成分の同時無相関化（源信号成分の共分散行列の同時対角化）が実現し、少ない計算量での信号分離が可能になる。

上述の最尤推定は、次の最適化問題（（２００）式）として定式化される。

ここで、Θはすべてのパラメータの全体を表し、分散パラメータｖ_ｊ（ｋ）（ｊ＝１，・・・，Ｊ；ｋ＝１，・・・，Ｋ）と、同時無相関化行列Ｐと、同時無相関化された空間共分散行列Λ_ｊ（ｊ＝１，・・・，Ｊ）と、からなる。Ｌ_２（Θ）は、次の（２０１）式の対数尤度関数であり、（１８２）式に（１９７）式〜（１９９）式を代入することにより得られる。なお、ＧＬ（Ｉ，Ｃ）はすべての正則なＩ次複素行列からなる集合である。

本第５の実施形態では、ＥＭアルゴリズムにより、上述の最尤推定を実現する。このＥＭアルゴリズムはＥステップとＭステップとからなる。

Ｅステップでは、パラメータΘの現在の推定値＾Θに基づいて、源信号成分の事後確率を更新する。従来技術の場合と同様に、ｘ_ｊ（ｋ）の事後確率ｐ（ｘ_ｊ（ｋ）｜ｙ（ｋ），＾Θ）は、（１８４）式の右辺の複素ガウス分布であり、その平均は、（１８５）式で、その共分散行列は、（１８６）式で与えられることが示せる。ただし、従来技術の場合とは異なり、（１８５）式及び（１８６）式における空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊは、同時無相関化行列Ｐの推定値＾Ｐ及び同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊを用いて、次の（２０２）式〜（２０４）式のように表される。

（１８５）式及び（１８６）式に（２０２）式〜（２０４）式を代入することにより、次の（２０５）式〜（２０８）式を得る。

したがって、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を同時無相関化行列Ｐの推定値＾Ｐにより基底変換して得られる＾ｘ_ｊ´（ｋ）（以下、同時無相関化された源信号成分の推定値と呼ぶ）は、（１７３）式により更新でき、事後共分散行列Σ_ｊ（ｋ）を同時無相関化行列Ｐの推定値＾Ｐにより基底変換して得られるΣ_ｊ´（ｋ）（以下、同時無相関化された事後共分散行列と呼ぶ）は、（１７４）式により更新できる。その際、同時対角化により、逆行列演算及び行列乗算は、対角行列に対する逆行列演算及び行列乗算になっているため、その計算量は、Ｏ（Ｉ^３）ではなくＯ（Ｉ）であり少ない計算量で計算できる。

Ｍステップでは、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊと、同時無相関化行列Ｐの推定値＾Ｐと、を下記のＱ関数の最大化に基づいて更新する（（２０９）式及び（２１０）式参照）。ただし、（２０９）式及び（２１０）式ではΘに依らない定数項は無視してあり、したがって等号はその両辺がΘに依らない定数の差を除いて等しいことを表している。

（２１０）式のＱ関数を分散パラメータｖ_ｊ（ｋ）に関して偏微分して０とおくと、次の（２１１）式及び（２１２）式を得る。

これで、（１７５）式の分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）の更新式が導かれた。

また、（２１０）式のＱ関数を同時無相関化された空間共分散行列Λ_ｊに関して偏微分したものが零行列に等しいと置くと、次の（２１３）式及び（２１４）式を得る。

これで、（１７６）式の同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊの更新則が導かれた。なお、（２１３）式において、絶対値二乗｜・｜^２はベクトルの要素ごとに計算するものとし、ベクトルαに対してｄｉａｇ（α）は、αの要素を対角要素とする対角行列を表す。

（２０９）式のＱ関数を同時無相関化行列Ｐの複素共役Ｐ^＊に関して偏微分すると、次の（２１５）式を得る。

これが零行列に等しいと置いて両辺のｖｅｃを取ると、次の（２１６）式を得る。

ここで、行列Ａに対し、ｖｅｃ（Ａ）は、行列Ａの列ベクトルを縦に並べて得られる列ベクトルを表す。これに関し、次の（２１７）式の公式が成り立つ。

また、「×」を丸で囲んだ記号は、クロネッカ積を表す。ブロック対角行列の性質より、（２１６）式は次の（２１８）式〜（２２０）式と等価である。

そこで、不動点法に基づいて、（１７７）式により、同時無相関化行列Ｐの推定値＾Ｐを更新すれば良い。このように不動点法に基づくことで、勾配法や自然勾配法とは異なり学習率の調整が不要であるため、安定な最適化が可能である。

以下、事後確率の詳細な導出を述べる。隠れ変数ｘ_１（ｋ），・・・，ｘ_Ｊ−１（ｋ）の対数事後確率は、次の（２２１）式〜（２２３）式のように書き下せる。

ここで、等号「＝」の上に「ｃ」を付した記号は、両辺がｘ_１（ｋ），・・・，ｘ_Ｊ−１（ｋ）に依らない定数の差を除いて等しいことを表す。行列Ａは、Ｉ次の単位行列ＩをＪ−１個横に並べた次の（２２４）式の行列を表す。また、空間共分散行列Ｒ_ｊ（ｊ＝１，・・・，Ｊ）は（１９７）式〜（１９９）式で与えられる。なお、本明細書では、Ｉを観測信号の個数（スカラー）と単位行列（行列）との２通りの意味で用いているが、Ｉがスカラーであるか行列であるかは文脈より明らかであるから、これにより混乱が生じるおそれはない。

~ｘ（ｋ）は、ｘ_１（ｋ），・・・，ｘ_Ｊ−１（ｋ）を縦に並べた次の（２２５）式のベクトルを表す。

（２２３）式より、隠れ変数ｘ_１（ｋ），・・・，ｘ_Ｊ−１（ｋ）の事後確率が次の（２２６）式の複素ガウス分布であることが導かれる。

ここで、共分散行列~Σ（ｋ）及び平均＾~ｘ（ｋ）は次の（２２７）式〜（２３０）式及び（２３１）式〜（２３３）式で与えられる。但し、（２２７）式から（２２８）式への式変形の際に逆行列補題を用いた。

後に分かるように、Ｅステップでは、行列~Σ全体を計算する必要は実はなく、行列~ΣをＩ次の正方行列を要素とするＪ−１次の正方ブロック行列とみなしたときの対角ブロックのみを計算すれば良い。ｊ番目の対角ブロックをΣ_ｊ（ｋ）と書くことにすると、これは次の（２３４）式のように計算できる。明らかに、Ｊ＝２の場合には、（２３４）式は（１０１）式と等価である。

また、＾ｘ_ｊ（ｋ）を、＾~ｘ（ｋ）をＩ次元縦ベクトルを要素とするサイズ（Ｊ−１）×１のブロック行列とみなしたときの（ｊ，１）ブロックと定義すると、これは次の（２３５）式のように表される。

＾ｘ_ｊ（ｋ）及びΣ_ｊ（ｋ）は次の（２３６）式及び（２３７）式のようにも表されることに注意する。

ただし、上式においては、ｊ＝１，・・・，Ｊ−１である。また、Ｅ（・｜ｙ（ｋ））は、隠れ変数の事後確率ｐ（ｘ_１（ｋ），・・・，ｘ_Ｊ−１（ｋ）｜ｙ（ｋ），Θ）に関する期待値（（２３８）式）である。

これにならい、ｊ＝Ｊに対しても、＾ｘ_Ｊ（ｋ）及びΣ_Ｊ（ｋ）を次の（２３９）式〜（２４３）式及び（２４４）式〜（２５１）式のように定義する。

よって、以下の（２５２）式及び（２５３）式は、ｊ＝１，・・・，Ｊに対して成立する。

以下、Ｑ関数の詳細な導出について説明する。観測信号ベクトルｙ（ｋ）と隠れ変数ｘ_１（ｋ），・・・，ｘ_Ｊ−１（ｋ）との同時分布の対数は、次の（２５４）式〜（２５６）式のように書き下される。ただし、（２５４）式〜（２５６）式では定数項は無視してあり、したがって等号はその両辺が定数の差を除いて等しいことを表している。

ここで、次の各式（（２５７）式〜（２６２）式）が成立することに注意する。

したがって、Ｑ関数は次式（（２６３）式〜（２６６）式）で与えられる。（２６４）式は（２０９）式に他ならない。

ここで、（２６６）式における絶対値の二乗｜・｜^２は、要素ごとに取るものとする。また、等号「＝」の上に「ｃ」を付した記号は、両辺がΘに依らない定数の差を除いて等しいことを表す（先と異なる意味で記号が用いられていることに注意されたい）。

［第５の実施形態の変形例１］
次に、第５の実施形態の変形例１について説明する。本第５の実施形態の変形例１では、信号パラメータ分析部における処理の変形例について説明する。第５の実施形態では、信号パラメータ分析部は、不動点法に基づいて同時無相関化行列Ｐの推定値＾Ｐを更新していた。これに対し、第５の実施形態の変形例１に掛かる信号パラメータ分析部は、不動点法以外の任意の最適化法（例えば、勾配法、自然勾配法、ニュートン法など）に基づいて、同時無相関化行列Ｐの推定値＾Ｐを更新する。

例えば、勾配法の場合、本第５の実施形態の変形例１に係る信号パラメータ分析部は、同時無相関化行列Ｐの推定値＾Ｐを次の（２６７）式により更新する（これは、（２１５）式から直ちに導かれる）。

ただし、信号パラメータ分析部における同時無相関化行列Ｐの推定値＾Ｐの最初の更新の際には、（２６７）式の右辺において、＾Ｐを＾Ｐ^（０）で置き換える。μは学習率と呼ばれる正の数であり、例えばμ＝０．０５などと定めればよい。

また、自然勾配法の場合、本第５の実施形態の変形例１に係る信号パラメータ分析部は、同時無相関化行列Ｐの推定値を次の（２６８）式により更新する（これも（２１５）式から直ちに導かれる）。

ただし、信号パラメータ分析部における同時無相関化行列Ｐの推定値＾Ｐの最初の更新の際には、（２６８）式の右辺において、＾Ｐを＾Ｐ^（０）で置き換える。μは、学習率と呼ばれる正の数であり、例えば、μ＝０．０５などと定めればよい。自然勾配法によれば、空間の計量が適切に考慮されるため、一般に勾配法よりも安定な更新が可能であることが知られている。

［第５の実施形態の変形例２］
次に、第５の実施形態の変形例２について説明する。本第５の実施形態の変形例２では、信号パラメータ分析部における処理の変形例について説明する。第５の実施形態や第５の実施形態の変形例１では、信号パラメータ分析部は、Ｑ関数に基づいて同時無相関化行列Ｐの推定値＾Ｐを更新していた。これに対し、本第５の実施形態の変形例２に係る信号パラメータ分析部は、直接、対数尤度関数Ｌ_２（Θ）に基づいて、同時無相関化行列Ｐの推定値＾Ｐを更新する。

本第５の実施形態の変形例２に係る信号パラメータ分析部は、同時無相関化部からの同時無相関化された観測信号ベクトルｙ´（ｋ）と、空間共分散行列更新部からの同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊと、分散パラメータ更新部からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、（信号パラメータ分析部における同時無相関化行列Ｐの推定値＾Ｐの最初の更新の際には、更に、初期化部からの同時無相関化行列Ｐの初期推定値＾Ｐ^（０）と、）に基づいて、同時無相関化行列Ｐの推定値＾Ｐを次式（（２６９）式）により更新する（Ｅは単位行列を表すものとする）。

ただし、信号パラメータ分析部における＾Ｐの最初の更新の際には、（２６９）式の右辺において、＾Ｐを＾Ｐ^（０）で置き換える。なお、信号パラメータ分析部は、（２６９）式を複数回連続して適用することにより、同時無相関化行列Ｐの推定値＾Ｐを更新してもよい。

（２６９）式の対数尤度関数Ｌ_２（Θ）に基づく同時無相関化行列Ｐの推定値＾Ｐの更新則の導出について説明する。（２０１）式の対数尤度関数を同時無相関化行列Ｐの複素共役Ｐ^＊に関して偏微分すると次の（２７０）式を得る。

（２７０）式が零行列に等しいと置いて、両辺のｖｅｃを取ると、次式（（２７１）式）を得る。

（２７１）式の右辺の行列は、次の（２７２）式及び（２７３）式のように変形できる。

よって次式（（２７４）式及び（２７５）式）を得る。

（２７５）式に（２７６）式を代入することにより、（２６９）式の不動点法による同時無相関化行列Ｐの推定値＾Ｐの更新則が得られる。

このように、不動点法に基づくことで、勾配法や自然勾配法とは異なり学習率の調整が不要であるため、安定な最適化が可能である。

また、第５の実施形態の変形例１と同様に、本第５の実施形態の変形例２においても、不動点法に限らず任意の最適化法（勾配法、自然勾配法、ニュートン法など）により同時無相関化行列Ｐの推定値＾Ｐを更新してよい。以下、そのような例について説明する。

勾配法の場合、本第５の実施形態の変形例２に係る信号パラメータ分析部は、観測信号ベクトル作成部からの観測信号ベクトルｙ（ｋ）と、空間共分散行列更新部からの同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊと、分散パラメータ更新部からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、（信号パラメータ分析部における同時無相関化行列Ｐの推定値＾Ｐの最初の更新の際には、更に、初期化部からの同時無相関化行列Ｐの初期推定値＾Ｐ^（０）と、）に基づいて、同時無相関化行列Ｐの推定値＾Ｐを次の（２７７）式により更新する（これは、（２７０）式より直ちに導かれる）。

ただし、信号パラメータ分析部における同時無相関化行列Ｐの推定値＾Ｐの最初の更新の際には、（２７７）式の右辺において、＾Ｐを＾Ｐ^（０）で置き換える。なお、信号パラメータ分析部は、（２７７）式を複数回連続して適用することにより同時無相関化行列Ｐの推定値＾Ｐを更新してもよい。また、信号パラメータ分析部は、（２７７）式に（２７６）式を代入することにより得られる、同時無相関化された観測信号ベクトルｙ´（ｋ）を用いて表した式を、（２７７）式の代わりに用いても良い。μは、学習率と呼ばれる所定の正の数であり、例えば、μ＝０．０５などと定めればよい。この勾配法の場合には、（２７７）式の１回の適用につき逆行列計算が１回のみで行列乗算は不要（対角行列に対する逆行列計算はカウントしていない）であるため、不動点法の場合と比べてより少ない計算量で同時無相関化行列Ｐの推定値＾Ｐを更新できる。

また、自然勾配法の場合、本第５の実施形態の変形例２に係る信号パラメータ分析部は、観測信号ベクトル作成部からの観測信号ベクトルｙ（ｋ）と、空間共分散行列更新部からの同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊと、分散パラメータ更新部からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、（信号パラメータ分析部における同時無相関化行列Ｐの推定値＾Ｐの最初の更新の際には、更に、初期化部からの同時無相関化行列Ｐの初期推定値＾Ｐ^（０）と、）に基づいて、同時無相関化行列Ｐの推定値＾Ｐを次の（２７８）式により更新する（これも（２７０）式より直ちに導かれる）。

ただし、信号パラメータ分析部における同時無相関化行列Ｐの推定値＾Ｐの最初の更新の際には、（２７８）式の右辺において、＾Ｐを＾Ｐ^（０）で置き換える。なお、信号パラメータ分析部は、（２７８）式を複数回連続して適用することにより同時無相関化行列Ｐの推定値＾Ｐを更新してもよい。また、信号パラメータ分析部は、（２７８）式に（２７６）式を代入することにより得られる、同時無相関化された観測信号ベクトルｙ´（ｋ）を用いて表した式を、（２７８）式の代わりに用いても良い。μは、学習率と呼ばれる所定の正の数であり、例えば、μ＝０．０５などと定めればよい。また、自然勾配法によれば、空間の計量が適切に考慮されるため、一般に勾配法よりも安定な更新が可能であることが知られている。

［第５の実施形態の変形例３］
次に、第５の実施形態の変形例３について説明する。第５の実施形態の変形例３では、第３の実施形態にならい、周波数ごとに第５の実施形態を適用する。

第５の実施形態の変形例３においては、Ｊ個の音源信号が混在する状況において、それぞれ異なる位置でマイクロホンにより取得されたＩ個（Ｉは２以上の整数）の観測信号ｙ_ｉ（ｔ）が信号分析装置に入力されるものとする。観測信号ｙ_ｉ（ｔ）は、周波数分解され、周波数ごとに第５の実施形態に記載の方法で信号分離がなされ、最後に時間領域に変換される。

図５及び図６を用いて、第５の実施形態の変形例３に係る信号分析装置の構成及び信号分析処理の処理手順について説明する。

周波数領域変換部３１０（ステップＳ４０、ステップＳ４１）及び時間領域変換部３３０（ステップＳ４３）での処理は、第３の実施形態と同じである。以下、相違点である周波数領域信号分離部３２０（ステップＳ４２）での処理について説明する。

周波数領域信号分離部３２０は、周波数領域変換部３１０から（２７９）式に示す時間周波数領域の観測信号を受け取り、ｆごとに第５の実施形態の方法により周波数領域での信号分離を行うことにより、（２８０）式に示す時間周波数領域の源信号成分の推定値を計算し、時間領域変換部３３０に渡す（ステップＳ４２）。

すなわち、ｆ＝１に対しては、（２８１）式に示す観測信号を入力として、第５の実施形態の方法により（２８２）式の源信号成分の推定値を得る。

そして、ｆ＝２に対しては、（２８３）式に示す観測信号を入力として、第５の実施形形態の方法により（２８４）式の源信号成分の推定値を得る。

以下、同じ処理を各ｆに対して繰り返し、最後に、ｆ＝Ｆに対しては、（２８５）式に示す観測信号を入力として、第５の実施形態の方法により（２８６）式の源信号成分の推定値を得る。

こうして得られた時間周波数領域の源信号成分の推定値（（２８７）式参照）を時間領域変換部３３０に渡す。

［第５の実施形態の効果］
本発明の効果を確認するために、従来の方法及び第５の実施形態を用いた音源分離実験について説明する。

図１４は、従来の方法及び第５の実施形態を用いた音源分離実験の実験条件を示す図である。図１５は、従来の方法及び第５の実施形態を用いた音源分離実験結果（real time factor）を示す図である。図１６は、従来の方法及び第５の実施形態を用いた音源分離実験結果（音源分離性能）を示す図である。

８秒の英語音声に、図１４に示す実験室で測定された残響インパルス応答を畳み込んで残響音声を生成し、複数人の残響音声を加算することで観測信号を生成した。このとき、観測信号に対して、各方法を用いて信号分離を行った場合の性能を図１５及び図１６に示す。図１５は、real time factor（観測信号長を１としたときの、処理に要した時間）である。また、図１６は、音源分離性能である信号対歪み比（値が大きいほど音源分離性能が高い）である。

図１５及び図１６に示すように、本第５の実施形態に係る方法によれば、音源数が３個以上である場合にも、４２０倍以上の計算時間削減が実現でき、従来の方法を上回る音源分離性能を得ることができた。

このように、本第５の実施形態によれば、従来の方法と比べて計算時間を大きく削減できるのに加え、従来の方法よりも高い信号分離性能を実現できる。本第５の実施形態は、信号数が２個でない場合にも適用でき、信号数が３個以上の場合にも適用できる。

［第６の実施形態］
次に、第６の実施形態に係る信号分析装置の構成、処理の流れ及び効果について説明する。第６の実施形態では、第４の実施形態の変形例を実施するための具体的な手法について説明する。

本第６の実施形態では、第５の実施形態と同様に、一般にＪ個の源信号が混在し、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊが未知の状況で、分散パラメータｖ_ｊ（ｋ）と、空間共分散行列Ｒ_ｊをモデル化するパラメータである同時無相関化された空間共分散行列Λ_ｊ及び同時無相関化行列Ｐと、源信号成分ｘ_ｊ（ｋ）と、を同時に推定することにより、信号分離を実現する方法について説明する。

なお、第５の実施形態は、ＥＭアルゴリズムに基づいていたのに対し、本第６の実施形態は補助関数法に基づいている。本第６の実施形態では、それぞれ異なる位置で取得されたＩ個の観測信号ｙ_ｉ（ｋ）が信号分析装置に入力されるものとする。

［信号分析装置の構成及び処理］
図１７及び図１８を用いて、第６の実施形態に係る信号分析装置の構成及び信号分析処理の処理手順について説明する。図１７は、第６の実施形態に係る信号分析装置の構成の一例を示す図である。図１８は、第６の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。

図１７に示すように、信号分析装置６０１は、初期化部（図示しない）、観測信号ベクトル作成部６１０、同時無相関化部６３０、シングルチャネルウィーナーフィルタ部６４０、観測共分散行列更新部６５０、分散パラメータ更新部６６０、空間共分散行列更新部６７０、信号パラメータ分析部６２０、収束判定部（図示しない）及び同時無相関化逆変換部６８０を有する。

まず、信号分析装置の各部の概要について説明する。初期化部（図示しない）の処理（ステップＳ８１）、観測信号ベクトル作成部６１０の処理（ステップＳ８２，Ｓ８３）、同時無相関化部６３０の処理（ステップＳ８４）及びシングルチャネルウィーナーフィルタ部６４０の処理（ステップＳ８５）については、第５の実施形態と同様であるから、説明を省略する。

観測共分散行列更新部６５０は、分散パラメータ更新部６６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）（ただし例外として、観測共分散行列更新部６５０における最初の処理の際には、初期化部からの分散パラメータｖ_ｊ（ｋ）の初期推定値＾ｖ_ｊ ^（０）（ｋ））と、空間共分散行列更新部６７０からの同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊ（ただし例外として、観測共分散行列更新部６５０における最初の処理の際には、初期化部からの同時無相関化された空間共分散行列Λ_ｊの初期推定値＾Λ_ｊ ^（０））と、を受け取って、同時無相関化された観測共分散行列Ｕ´（ｋ）を次式（２８８）式により更新する（ステップＳ８６）。

ただし、観測共分散行列更新部６５０における最初の処理の際には、上式において、＾ｖ_ｊ（ｋ）を＾ｖ_ｊ ^（０）（ｋ）に、＾Λ_ｊを＾Λ_ｊ ^（０）に、それぞれ置き換える。

分散パラメータ更新部６６０は、シングルチャネルウィーナーフィルタ部６４０からの同時無相関化された源信号成分ｘ_ｊ´（ｋ）の推定値＾ｘ_ｊ´（ｋ）と、観測共分散行列更新部６５０からの同時無相関化された観測共分散行列Ｕ´（ｋ）と、空間共分散行列更新部６７０からの同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊ（ただし例外として、分散パラメータ更新部６６０における最初の処理の際には、初期化部からの同時無相関化された空間共分散行列Λ_ｊの初期推定値＾Λ_ｊ ^（０））と、を受け取って、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）を次式（（２８９）式）により更新する（ステップＳ８７）。

ただし、分散パラメータ更新部６６０における最初の処理の際には、上式において、＾Λ_ｊを＾Λ_ｊ ^（０）で置き換える。

空間共分散行列更新部６７０は、シングルチャネルウィーナーフィルタ部６４０からの同時無相関化された源信号成分ｘ_ｊ´（ｋ）の推定値＾ｘ_ｊ´（ｋ）と、観測共分散行列更新部６５０からの同時無相関化された観測共分散行列Ｕ´（ｋ）と、分散パラメータ更新部６６０からの分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、を受け取って、同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊを次の（２９０）式により更新する（ステップＳ８８）。

信号パラメータ分析部６２０は、第５の実施形態の変形例２と同様の処理を行って、同時無相関化行列Ｐの推定値＾Ｐを更新する（ステップＳ８９）。

図示しない収束判定部は、収束したかどうかの判定を行う（ステップＳ９０）。収束判定部が収束していないと判定した場合（ステップＳ９０：Ｎｏ）、同時無相関化部６３０での処理（ステップＳ８４）に戻って、処理が継続される。

一方、収束判定部が収束したと判定した場合（ステップＳ９０：Ｙｅｓ）、同時無相関化部６３０における処理（ステップＳ９１）、シングルチャネルウィーナーフィルタ部６４０における処理（ステップＳ９２）、及び、同時無相関化逆変換部６８０における処理（ステップＳ９３）が、第５の実施形態と同様にして行われることにより、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）が計算され、出力される。

［信号分析装置の処理の背景］
次に、第６の実施形態に係る信号分析装置の背景について説明する。第５の実施形態で述べたように、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊが未知の状況で、（１８０）式により、源信号成分ｘ_ｊ（ｋ）の推定値＾ｘ_ｊ（ｋ）を計算する（すなわち信号分離を実現する）ためには、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊを推定する必要がある。

補助関数法に基づいて、（１８１）式の最適化問題を解くことにより、分散パラメータｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊを最尤推定することができる。

この補助関数法では、補助変数Φ及び次式を満たす補助関数Ｌ_１ ⁻（Ψ，Φ）に基づいて、補助関数Ｌ_１ ⁻（Ψ，Φ）を補助変数Φについて最大化することによりΦを更新するステップと、補助関数Ｌ_１ ⁻（Ψ，Φ）が減少しないようにパラメータΨを更新するステップと、が交互に反復される（（２９１）式参照）。

本第６の実施形態では、補助関数として（２９２）式のＬ_１ ⁻（Ψ，Φ）を用いる。

ただし、右辺における「定数」は、パラメータΨ及び補助変数Φによらない定数である。（２９２）式が（２９１）式の補助関数の条件を満たすことは、参考文献３及び参考文献４と同様に示すことができる。なお、補助変数Φは、Ｕ（ｋ）及びＧ_ｊ（ｋ）からなり、Ｕ（ｋ）は正定値エルミート行列であり、Ｇ_ｊ（ｋ）は（２９３）式を満たすものとする。

補助関数Ｌ_１ ⁻（Ψ，Φ）を補助変数Φについて最大化することによりΦを更新するステップでは、分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）及び空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊに基づいて、Ｕ（ｋ）及びＧ_ｊ（ｋ）を（２９４）式及び（２９５）式により更新する。Ｕ（ｋ）は観測信号ベクトルｙ（ｋ）の共分散行列に他ならないから、観測共分散行列と呼ぶ。

補助関数Ｌ_１ ⁻（Ψ，Φ）が減少しないようにパラメータΨを更新するステップでは、上で更新した補助変数に基づいて、次の（２９６）式及び（２９７）式により分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）と、空間共分散行列Ｒ_ｊの推定値＾Ｒ_ｊと、を更新する。

ここで、Ａ＃Ｂは、（２９８）式で定義されるＡとＢとの幾何平均である。

また、行列Ｃ_ｊ及び行列Ｄ_ｊは、それぞれ（２９９）式及び（３００）式により計算される。

（２９４）式〜（３００）式の更新式の導出は、参考文献３及び参考文献４と同様にできる。更に、補助変数Ｇ_ｊ（ｋ）を代入消去することにより、補助変数Ｇ_ｊ（ｋ）が現れない更新則を導くことができる。この場合の更新則（１回の反復）を以下の（３０１）式〜（３０５）式に示す。明らかに、Ｊ＝２の場合には、（３０１）式〜（３０５）式による更新は、（１５６）式〜（１５９）式による更新と等価である。

しかしながら、上記の方法では、（３０２）式〜（３０４）式に含まれる、逆行列計算Ｕ（ｋ）^−１をサンプル点ごとに行う必要があった。そのため、上記の方法は膨大な計算量を要する。

これに対し、本第６の実施形態では、第５の実施形態と同様に、Ｊ個の源信号成分ｘ_１（ｋ），・・・，ｘ_Ｊ（ｋ）の同時無相関化に基づいて、逆行列計算Ｕ（ｋ）^−１を少ない計算量で計算できる対角行列に対する逆行列計算に帰着させることにより、少ない計算量での信号分離を実現する。具体的には、本第６の実施形態では、第５の実施形態と同様に、Ｊ個の空間共分散行列が（１９７）式〜（１９９）式の形に表されるという制約条件の下、同時無相関化行列Ｐ及び同時無相関化された空間共分散行列Λ_１，・・・，Λ_Ｊを最尤推定する。この最尤推定は、（２００）式の最適化問題として定式化される。本第６の実施形態では、補助関数法に基づいて、上述の最尤推定を実現する。この方法では、補助関数法の１回の反復を適用することにより分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）及び同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊを更新するステップと、同時無相関化行列Ｐの推定値＾Ｐを更新するステップと、が交互に反復される。

補助関数法に基づいて分散パラメータｖ_ｊ（ｋ）の推定値＾ｖ_ｊ（ｋ）及び同時無相関化された空間共分散行列Λ_ｊの推定値＾Λ_ｊを更新するステップは、（３０６）式の補助関数Ｌ_２ ⁻（Θ，Φ）に基づいている。

ただし、右辺の「定数」はパラメータΘ及び補助変数Φによらない定数である。

先と同様に補助変数Ｇ_ｊ（ｋ）を代入消去すると、更新則は下記の（３０７）式〜（３１１）式のように与えられる。ただし、（３１０）式の右辺における斜体のＤは、行列の非対角成分を０で置き換える演算子である。

これは、更に次の（３１２）式〜（３１５）式に書き直せる。

同時無相関化行列Ｐの推定値＾Ｐの更新則の導出については、第５の実施形態の変形例２と同様であるから説明を省略する。以上で、第６の実施形態における更新則が導かれた。

［第６の実施形態の変形例］
第６の実施形態では対数尤度に基づいて同時無相関化行列Ｐの推定値＾Ｐを更新するが、補助関数に基づいて同時無相関化行列Ｐの推定値＾Ｐを更新してもよい。詳細は、第５の実施形態及び第５の実施形態の変形例１と同様だから省略する。

［第６の実施形態の効果］
第６の実施形態によれば、同時無相関化に基づき、計算量を大幅に削減できる。また、ＥＭアルゴリズムは反復回数が少ない場合には必ずしも高い信号分離性能を実現できないが、第６の実施形態によれば、補助関数法に基づき、反復回数が少ない場合にも高い信号分離性能を実現できる。本第６の実施形態は、信号数が３個以上でも適用できる。

［第７の実施形態］
本第７の実施形態では、観測信号に含まれる源信号の総数ＬがＬ≧３の場合でも少ない計算量で信号分離を可能にする方法を示す。第５の実施形態や第６の実施形態でもそのような方法を示したが、本実施形態ではもう一つの実現方法を示す。

第５の実施形態や第６の実施形態は、「すべての時間周波数点において観測信号は、Ｌ個すべての源信号を含む」というモデル化（（１７９）式参照）の下、Ｌ個すべての源信号の同時無相関化（Ｌ個すべての源信号の空間共分散行列の同時対角化）に基づいて、Ｌ≧３の場合でも少ない計算量で信号分離を実現していた。この同時無相関化（同時対角化）は、「Ｌ個すべての源信号の空間共分散行列が１つの同時無相関化行列によって同時対角化できる」という制約条件（（１９７）式〜（１９９）式参照）に基づいていた。しかしながら、実際にはこの制約条件は現実に即していない場合もあり、そのような場合には、第５の実施形態や第６の実施形態ではモデル化精度が低下し、結果として信号分離性能が低下するおそれがあった。なお、第５の実施形態や第６の実施形態では、Ｌ個すべての源信号の同時無相関化（Ｌ個すべての源信号の空間共分散行列の同時対角化）に基づいているため、同時無相関化される源信号の最大個数（同時対角化される源信号の空間共分散行列の最大個数）をＪとすれば、Ｊ＝Ｌである。したがって、第５の実施形態や第６の実施形態では、これらをともにＪで表していた。

一方、本第７の実施形態では、源信号の時間周波数領域におけるスパース性に着目する。源信号は、時間周波数領域におけるスパース性、すなわち、「少数の時間周波数成分にエネルギーが集中する」という性質を持つことが多い（例えば、参考文献５「O. Yilmaz and S. Rickard, “Blind separation of speech mixtures via time-frequency masking”, IEEE Transactions on Signal Processing, vol. 52, no. 7, pp. 1830−1847, Jul. 2004.」参照）。この場合、各時間周波数点において観測信号は、Ｌ個すべての源信号を含むのではなく、Ｌ個よりも少ない源信号しか含まないとみなせることが多い。

そこで、本第７の実施形態では、「各時間周波数点において観測信号は高々２個の源信号しか含まない」というモデル化の下、各時間周波数点において観測信号に含まれる高々２個の源信号の同時無相関化（各時間周波数点において観測信号に含まれる高々２個の源信号の空間共分散行列の同時対角化）に基づいて、Ｌ≧３の場合でも少ない計算量で信号分離を実現する。一般化固有値問題に基づいて２個以下の源信号は常に一つの同時無相関化行列によって同時無相関化できる（２個以下の源信号の空間共分散行列は常に一つの同時無相関化行列によって同時対角化できる）から、本第７の実施形態では、第５の実施形態や第６の実施形態とは異なり、源信号の空間共分散行列に特段の制約条件を課す必要がない。

したがって、本第７の実施形態によれば、「Ｌ個すべての源信号の空間共分散行列が１つの同時無相関化行列によって同時対角化できる」という制約条件が現実に即していない場合でも高精度なモデル化が可能であり、結果として、そのような場合でも高い信号分離性能を実現可能である。なお、本第７の実施形態では、同時無相関化される源信号の最大個数Ｊ（同時対角化される源信号の空間共分散行列の最大個数Ｊ）はＪ＝２である。

本第７の実施形態では、Ｌ個（Ｌは２以上の自然数）の源信号が混在する状況において、それぞれ異なる位置で取得されたＩ個（Ｉは２以上の自然数）の観測信号ｙ_ｉ（ｔ）が信号分析装置に入力されるものとする。

［信号分析装置の構成及び処理］
図５及び図６を用いて、本第７の実施形態に係る信号分析装置の構成及び信号分析処理の処理手順について説明する。図５は、本第７の実施形態に係る信号分析装置の構成の一例を示す図である。図６は、本第７の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。

図５に示すように、信号分析装置３０１は、周波数領域変換部３１０、周波数領域信号分離部３２０、時間領域変換部３３０からなる。

まず、信号分析装置３０１の各部の概要について説明する。周波数領域変換部３１０は、入力された観測信号ｙ_ｉ（ｔ）を受け取り（ステップＳ４０）、短時間フーリエ変換などにより時間周波数領域の観測信号ｙ_ｉ（ｋ，ｆ）を計算する（ステップＳ４１）。ここで、ｋは時間フレームのインデックスであり、ｆは周波数ビンのインデックスである。ｙ_ｉ（ｋ，ｆ）をｙ_ｉ ^ｆ（ｋ）とも表す。

周波数領域信号分離部３２０は、ｆごとに信号分離を行って（ステップＳ４２）、時間周波数領域の源信号成分ｘ_ｊ（ｋ，ｆ）の推定値＾ｘ_ｊ（ｋ，ｆ）を得る。＾ｘ_ｊ（ｋ，ｆ）を＾ｘ_ｊ ^ｆ（ｋ）とも表す。

時間領域変換部３３０は、時間周波数領域の源信号成分ｘ_ｊ（ｋ，ｆ）の推定値＾ｘ_ｊ（ｋ，ｆ）に基づいて、逆短時間フーリエ変換などにより、時間領域の源信号成分ｘ_ｊ（ｔ）の推定値＾ｘ_ｊ（ｔ）を計算し（ステップＳ４３）、出力する。

次に、周波数領域信号分離部３２０についてより詳細に説明する。図１９及び図２０を用いて、周波数領域信号分離部３２０の構成及び周波数領域信号分離処理の処理手順について説明する。図１９は、周波数領域信号分離部３２０の構成の一例を示す図である。図２０は、周波数領域信号分離処理の処理手順の一例を示すフローチャートである。

図１９に示すように、周波数領域信号分離部３２０は、初期化部（図示しない）、観測信号ベクトル作成部７１０、信号パラメータ分析部７２０、同時無相関化部７３０、信号存在事後確率更新部７４０、シングルチャネルウィーナーフィルタ部７５０、事後共分散行列更新部７６０、混合重み更新部７７０、分散パラメータ更新部７８０、空間共分散行列更新部７９０、収束判定部（図示しない）及び源信号成分推定部８００を有する。これら周波数領域信号分離部３２０の各部での処理は、周波数ビン毎に行う。

周波数領域信号分離部３２０は、ｆ←０により周波数ビンのインデックスｆを初期化する（ステップＳ１０１）。周波数領域信号分離部３２０は、ｆ←ｆ＋１により周波数ビンのインデックスｆを１増やす（ステップＳ１０２）。

図示しない初期化部は、混合重みα_Ｂ（ｆ）の初期推定値＾α_Ｂ ^（０）（ｆ）と、分散パラメータｖ_ｊ（ｋ，ｆ）の初期推定値＾ｖ_ｊ ^（０）（ｋ，ｆ）と、空間共分散行列Ｒ_ｊ（ｆ）の初期推定値＾Ｒ_ｊ ^（０）（ｆ）と、を設定する（ステップＳ１０３）。これらは、例えば乱数に基づいて設定すればよい。

観測信号ベクトル作成部７１０は、周波数領域変換部からの時間周波数領域の観測信号ｙ_ｉ（ｋ，ｆ）を受け取り（ステップＳ１０４）、観測信号ベクトルｙ（ｋ，ｆ）を（３１６）式により作成する（ステップＳ１０５）。

ここで、本第７の実施形態における観測信号ベクトルｙ（ｋ，ｆ）のモデル化について説明する。本実施形態では、各時間周波数点（ｋ，ｆ）において観測信号ベクトルｙ（ｋ，ｆ）は高々２個の源信号成分しか含まないと仮定する。各時間周波数点（ｋ，ｆ）において観測信号ベクトルｙ（ｋ，ｆ）に含まれる源信号成分は、そのような源信号成分のインデックスの全体である集合Ｂ（ｋ，ｆ）により表すことができる（ただし、集合Ｂ（ｋ，ｆ）の元の個数＃Ｂ（ｋ，ｆ）は＃Ｂ（ｋ，ｆ）≦２を満たすものとする）。このとき、観測信号ベクトルｙ（ｋ，ｆ）は、Ｂ（ｋ，ｆ）のすべての元ｊに対する源信号成分ｘ_ｊ（ｋ，ｆ）の和として、（３１７）式によりモデル化される。

例えば、ある時間周波数点（ｋ，ｆ）でＢ（ｋ，ｆ）＝｛１，３｝である場合、（ｋ，ｆ）では、（３１７）式は、（３１８）式のようになる。

Ｂ（ｋ，ｆ）が取りうる値（集合）の全体を、集合｛１，・・・，Ｌ｝の空でない部分集合からなる集合系（すなわち、冪集合２^{｛１，・・・，Ｌ｝}の部分集合であって空集合を元として含まないもの）Ｈで表す。例えば、Ｂ（ｋ，ｆ）が＃Ｂ≦２なる任意の空でない集合Ｂ⊂｛１，・・・，Ｌ｝になり得る場合、Ｈは、（３１９）式及び（３２０）式で与えられる。

（３１９）式及び（３２０）式の集合系ＨをＨ_ｍａｘで表す。本第７の実施形態では、集合系Ｈは、あらかじめ与えられており、Ｈ⊂Ｈ_ｍａｘであるとする。また、＃Ｂ＝２なる集合Ｂ∈Ｈが少なくとも一つはあるものとする。

すなわち、集合系Ｈの元（集合）Ｂが存在して、＃Ｂ＝２であるとする。

信号パラメータ分析部７２０は、空間共分散行列更新部７９０からの空間共分散行列Ｒ_ｊ（ｆ）の推定値＾Ｒ_ｊ（ｆ）（ただし、例外として、信号パラメータ分析部７２０での最初の処理の際には、初期化部からの空間共分散行列Ｒ_ｊ（ｆ）の初期推定値＾Ｒ_ｊ ^（０）（ｆ））を受け取って、Ｌ個の源信号成分のうちの任意の２個ｘ_ｊ（ｋ，ｆ），ｘ_ｌ（ｋ，ｆ）（ｊ，ｌは１以上Ｌ以下の相異なる任意の自然数）に対して、同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）と対角行列Λ_{ｊ，｛ｊ，ｌ｝}（ｆ）を計算する（ステップＳ１０６）。ここで、同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）は２個の源信号成分ｘ_ｊ（ｋ，ｆ），ｘ_ｌ（ｋ，ｆ）を同時無相関化する行列（すなわち、２個の源信号成分ｘ_ｊ（ｋ，ｆ），ｘ_ｌ（ｋ，ｆ）の空間共分散行列を同時対角化する行列）であり、添え字である集合｛ｊ，ｌ｝は、同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）によって同時無相関化される源信号成分ｘ_ｊ（ｋ，ｆ），ｘ_ｌ（ｋ，ｆ）のインデックスを表している。また、対角行列Λ_{ｊ，｛ｊ，ｌ｝}（ｆ）は、源信号成分ｘ_ｊ（ｋ，ｆ）の空間共分散行列を同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）によって対角化した際の対角行列である。

具体的には、信号パラメータ分析部７２０は以下の処理を行う。まず、２個の源信号成分ｘ_ｊ（ｋ，ｆ），ｘ_ｌ（ｋ，ｆ）（ｊ，ｌは、１≦ｊ＜ｌ≦Ｌなる任意の２つの自然数）の空間共分散行列の推定値に対する一般化固有値問題＾Ｒ_ｌ（ｆ）ｐ＝λ＾Ｒ_ｊ（ｆ）ｐを解くことにより、（ただし、例外として、信号パラメータ分析部７２０での最初の処理の際には、これらの源信号成分の空間共分散行列の初期推定値に対する一般化固有値問題＾Ｒ_ｌ ^（０）（ｆ）ｐ＝λ＾Ｒ_ｊ ^（０）（ｆ）ｐを解くことにより、）一般化固有値λ_{１，｛ｊ，ｌ｝}（ｆ），・・・，λ_{Ｉ，｛ｊ，ｌ｝}（ｆ）と、一般化固有値λ_{１，｛ｊ，ｌ｝}（ｆ），・・・，λ_{Ｉ，｛ｊ，ｌ｝}（ｆ）にそれぞれ対応する一般化固有ベクトルｐ_{１，｛ｊ，ｌ｝}（ｆ），・・・，ｐ_{Ｉ，｛ｊ，ｌ｝}（ｆ）と、を求める。ただし、一般化固有ベクトルｐ_{１，｛ｊ，ｌ｝}（ｆ），・・・，ｐ_{Ｉ，｛ｊ，ｌ｝}（ｆ）は、（３２１）式を満たすように取るものとする。

なお、信号パラメータ分析部７２０での最初の処理の際には、（３２１）式において、＾Ｒ_ｊ（ｆ）を＾Ｒ_ｊ ^（０）（ｆ）で置き換える。信号パラメータ分析部７２０は、次に、一般化固有ベクトルｐ_{１，｛ｊ，ｌ｝}（ｆ），・・・，ｐ_{Ｉ，｛ｊ，ｌ｝}（ｆ）を列ベクトルとする行列を同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）とし、Ｉ次単位行列を対角行列Λ_{ｊ，｛ｊ，ｌ｝}（ｆ）とし、一般化固有値λ_{１，｛ｊ，ｌ｝}（ｆ），・・・，λ_{Ｉ，｛ｊ，ｌ｝}（ｆ）を対角要素とする対角行列を対角行列Λ_{ｌ，｛ｊ，ｌ｝}（ｆ）とする。

同時無相関化部７３０は、観測信号ベクトル作成部７１０からの観測信号ベクトルｙ（ｋ，ｆ）と、信号パラメータ分析部７２０からの同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）と、を受け取って、源信号成分ｘ_ｊ（ｋ，ｆ），ｘ_ｌ（ｋ，ｆ）（ｊ，ｌは、１≦ｊ＜ｌ≦Ｌなる任意の２つの自然数）について同時無相関化された観測信号ベクトルｙ´_{｛ｊ，ｌ｝}（ｋ，ｆ）を次の（３２２）式により計算する（ステップＳ１０７）。

ここで、添え字である集合｛ｊ，ｌ｝は、同時無相関化される源信号成分ｘ_ｊ（ｋ，ｆ），ｘ_ｌ（ｋ，ｆ）のインデックスを表す。

信号存在事後確率更新部７４０は、信号パラメータ分析部７２０からの同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）及び対角行列Λ_{ｊ，｛ｊ，ｌ｝}（ｆ）と、同時無相関化部７３０からの同時無相関化された観測信号ベクトルｙ´_{｛ｊ，ｌ｝}（ｋ，ｆ）と、混合重み更新部７７０からの混合重みα_Ｂ（ｆ）の推定値＾α_Ｂ（ｆ）（ただし、例外として、信号存在事後確率更新部７４０における最初の処理の際には、初期化部からの混合重み＾α_Ｂ（ｆ）の初期推定値＾α_Ｂ ^（０）（ｆ））と、分散パラメータ更新部７８０からの分散パラメータｖ_ｊ（ｋ，ｆ）の推定値＾ｖ_ｊ（ｋ，ｆ）（ただし、例外として、信号存在事後確率更新部７４０における最初の処理の際には、初期化部からの分散パラメータｖ_ｊ（ｋ，ｆ）の初期推定値＾ｖ_ｊ ^（０）（ｋ，ｆ））と、を受け取って、集合系Ｈの各元（集合）Ｂに対して、信号存在事後確率γ_Ｂ（ｋ，ｆ）を（３２３）式により更新する（ステップＳ１０８）。

なお、集合系Ｈの定義によっては、（３２３）式において、ただ一つの元からなる集合（一元集合）｛ｊ｝に対する同時無相関化行列Ｐ_｛ｊ｝（ｆ）、対角行列Λ_{ｊ，｛ｊ｝}（ｆ）及び同時無相関化された観測信号ベクトルｙ´_｛ｊ｝（ｋ，ｆ）が現れることがあるが、これらは二元集合に対するＰ_{｛ｊ，ｌ｝}（ｆ）、Λ_{ｊ，｛ｊ，ｌ｝}（ｆ）及びｙ´_{｛ｊ，ｌ｝}（ｋ，ｆ）を用いて、（３２４）式、（３２５）式及び（３２６）式により定義される。

ここで、ρ：｛１，．．．，Ｌ｝→｛１，．．．，Ｌ｝は、「すべてのｊ∈｛１，．．．，Ｌ｝に対してρ（ｊ）≠ｊ」なる所定の関数である。このように定義する理由は、次の通りである。

すなわち、１個の源信号成分ｘ_ｊ（ｋ，ｆ）の同時無相関化とは、単なるｘ_ｊ（ｋ，ｆ）の無相関化に他ならず、それは２個の源信号成分ｘ_ｊ（ｋ，ｆ），ｘ_ρ（ｊ）（ｋ，ｆ）の同時無相関化により当然実現されるからである。

シングルチャネルウィーナーフィルタ部７５０は、同時無相関化部７３０からの同時無相関化された観測信号ベクトルｙ´_{｛ｊ，ｌ｝}（ｋ，ｆ）と、信号パラメータ分析部７２０からの対角行列Λ_{ｊ，｛ｊ，ｌ｝}（ｆ）と、分散パラメータ更新部７８０からの分散パラメータｖ_ｊ（ｋ，ｆ）の推定値＾ｖ_ｊ（ｋ，ｆ）（ただし、例外として、シングルチャネルウィーナーフィルタ部７５０における最初の処理の際には、初期化部からの分散パラメータｖ_ｊ（ｋ，ｆ）の初期推定値＾ｖ_ｊ ^（０）（ｋ，ｆ））と、を受け取って、同時無相関化された事後平均ｍ´_{ｊ，｛ｊ，ｌ｝}（ｋ，ｆ）（ｊ，ｌは１以上Ｌ以下の相異なる任意の自然数）を（３２７）式により更新する（ステップＳ１０９）。

事後共分散行列更新部７６０は、信号パラメータ分析部７２０からの対角行列Λ_{ｊ，｛ｊ，ｌ｝}（ｆ）と、分散パラメータ更新部７８０からの分散パラメータｖ_ｊ（ｋ，ｆ）の推定値＾ｖ_ｊ（ｋ，ｆ）（ただし、例外として、事後共分散行列更新部７６０における最初の処理の際には、初期化部からの分散パラメータｖ_ｊ（ｋ，ｆ）の初期推定値＾ｖ_ｊ ^（０）（ｋ，ｆ））と、を受け取って、同時無相関化された事後共分散行列Σ´_{｛ｊ，ｌ｝}（ｋ，ｆ）（ｊ，ｌは１≦ｊ＜ｌ≦Ｌなる任意の２つの自然数）を次の（３２８）式により更新する（ステップＳ１１０）。

混合重み更新部７７０は、信号存在事後確率更新部７４０からの信号存在事後確率γ_Ｂ（ｋ，ｆ）を受け取って、（３２９）式により混合重みα_Ｂ（ｆ）の推定値＾α_Ｂ（ｆ）（Ｂは集合系Ｈの任意の元（集合））を更新する（ステップＳ１１１）。

分散パラメータ更新部７８０は、信号パラメータ分析部７２０からの対角行列Λ_{ｊ，｛ｊ，ｌ｝}（ｆ）と、信号存在事後確率更新部７４０からの信号存在事後確率γ_Ｂ（ｋ，ｆ）と、シングルチャネルウィーナーフィルタ部７５０からの同時無相関化された事後平均ｍ´_{ｊ，｛ｊ，ｌ｝}（ｋ，ｆ）と、事後共分散行列更新部７６０からの同時無相関化された事後共分散行列Σ´_{｛ｊ，ｌ｝}（ｋ，ｆ）とを受け取って、分散パラメータｖ_ｊ（ｋ，ｆ）の推定値＾ｖ_ｊ（ｋ，ｆ）を（３３０）式により更新する（ステップＳ１１２）。

なお、集合系Ｈの定義によっては、（３３０）式において、一元集合｛ｊ｝に対する同時無相関化された事後平均ｍ´_{ｊ，｛ｊ｝}（ｋ，ｆ）及び同時無相関化された事後共分散行列Σ´_｛ｊ｝（ｋ，ｆ）が現れることがあるが、これらは（３３１）式及び（３３２）式により定義される。

空間共分散行列更新部７９０は、信号パラメータ分析部７２０からの同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）と、信号存在事後確率更新部７４０からの信号存在事後確率γ_Ｂ（ｋ，ｆ）と、シングルチャネルウィーナーフィルタ部７５０からの同時無相関化された事後平均ｍ´_{ｊ，｛ｊ，ｌ｝}（ｋ，ｆ）と、事後共分散行列更新部７６０からの同時無相関化された事後共分散行列Σ´_{｛ｊ，ｌ｝}（ｋ，ｆ）と、分散パラメータ更新部７８０からの分散パラメータｖ_ｊ（ｋ，ｆ）の推定値＾ｖ_ｊ（ｋ，ｆ）と、を受け取って、空間共分散行列Ｒ_ｊ（ｆ）の推定値＾Ｒ_ｊ（ｆ）を（３３３）式により更新する（ステップＳ１１３）。

図示しない収束判定部は、収束したかどうかの判定を行う（ステップＳ１１４）。収束判定部が収束していないと判定した場合（ステップＳ１１４：Ｎｏ）、信号パラメータ分析部７２０での処理に戻って、処理が継続される。

一方、収束判定部が収束したと判定した場合（ステップＳ１１４：Ｙｅｓ）、源信号成分推定部８００での処理が行われる。源信号成分推定部８００は、信号パラメータ分析部７２０からの同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）と、信号存在事後確率更新部７４０からの信号存在事後確率γ_Ｂ（ｋ，ｆ）と、シングルチャネルウィーナーフィルタ部７５０からの同時無相関化された事後平均ｍ´_{ｊ，｛ｊ，ｌ｝}（ｋ，ｆ）と、を受け取って、源信号成分ｘ_ｊ（ｋ，ｆ）の推定値＾ｘ_ｊ（ｋ，ｆ）を（３３４）式により推定し（ステップＳ１１５）、出力する。

周波数領域信号分離部３２０は、ｆ＝Ｆであるかどうか判定する（ステップＳ１１６）。ｆ＝Ｆでないと判定された場合（ステップＳ１１６：Ｎｏ）、周波数領域信号分離部３２０でのｆ←ｆ＋１の処理に戻って、処理が継続される。

一方、ｆ＝Ｆであると判定された場合（ステップＳ１１６：Ｙｅｓ）、周波数領域信号分離部３２０での処理は終了する。

［信号分析装置の処理の背景］
次に、本第７の実施形態に係る信号分析装置の処理の背景について説明する。本第７の実施形態では、観測信号ベクトルｙ（ｋ，ｆ）は、（３１７）式でモデル化される。

本第７の実施形態では、源信号成分ｘ_ｊ（ｋ，ｆ）が平均０、共分散行列Ｓ_ｊ（ｋ，ｆ）＝ｖ_ｊ（ｋ，ｆ）Ｒ_ｊ（ｆ）の複素ガウス分布に従うと仮定する。このとき、Ｂ（ｋ，ｆ）が値Ｂ∈Ｈを取るという条件下での観測信号ベクトルｙ（ｋ，ｆ）の確率分布は、次の（３３５）式で与えられる。

ここで、Φ_１はすべてのモデルパラメータをまとめて表したものである。また、集合Ｂ（ｋ，ｆ）の確率分布は次の（３３６）式で与えられる。

ただし、α_Ｂ（ｆ）は、次の（３３７）式の制約条件を満たす。

（３３８）式及び（３３６）式より、観測信号ベクトルｙ（ｋ，ｆ）の確率分布は、次の（３３８）式で与えられる。

以下、α_Ｂ（ｆ）を混合重みと呼ぶ。パラメータΦ_１は、混合重みα_Ｂ（ｆ）、分散パラメータｖ_ｊ（ｋ，ｆ）及び空間共分散行列Ｒ_ｊ（ｆ）からなる。パラメータΦ_１は、次の（３３９）式の尤度関数の最大化により推定できる。

参考文献６「R. Ikeshita, M. Togami, Y. Kawaguchi, Y. Fujita, and K. Nagamatsu, “Local Gaussian model with source-set constraints in audio source separation”, IEEE International Workshop on Machine Learning for Signal Processing, Sept. 2017.」には、ＥＭアルゴリズムに基づいて上式の尤度関数を最大化することにより、パラメータΦ_１を推定する方法が開示されている。この従来技術では、ＥステップとＭステップを交互に反復することにより、パラメータΦ_１を推定する。

Ｅステップでは、信号存在事後確率γ_Ｂ（ｋ，ｆ）、事後平均ｍ_ｊ，Ｂ（ｋ，ｆ）及び事後共分散行列Σ_Ｂ（ｋ，ｆ）を次の（３４０）式〜（３４４）式により計算する。

ここで、信号存在事後確率γ_Ｂ（ｋ，ｆ）は、観測信号ベクトルｙ（ｋ，ｆ）が与えられたという条件下での集合Ｂ（ｋ，ｆ）の条件付き確率分布である事後分布として、次の（３４５）式で定義される。

ここで、＾Φ_１はパラメータΦ_１の現在の推定値を表す。また、事後平均ｍ_ｊ，Ｂ（ｋ，ｆ）及び事後共分散行列Σ_Ｂ（ｋ，ｆ）とは、観測信号ベクトルｙ（ｋ，ｆ）と集合Ｂ（ｋ，ｆ）が与えられたという条件下での源信号成分ｘ_ｊ（ｋ，ｆ）の条件付き確率分布である事後分布（複素ガウス分布であることが示される）の平均及び共分散行列として定義される。この事後分布は次の（３４６）式のように表される。

Ｍステップでは、混合重みα_Ｂ（ｆ）の推定値＾α_Ｂ（ｆ）、分散パラメータｖ_ｊ（ｋ，ｆ）の推定値＾ｖ_ｊ（ｋ，ｆ）及び空間共分散行列Ｒ_ｊ（ｆ）の推定値＾Ｒ_ｊ（ｆ）を、次の（３４７）式〜（３４９）式の更新則により更新する。

しかしながら、この従来技術では、（３４０）式の信号存在事後確率γ_Ｂ（ｋ，ｆ）の更新、（３４２）式の事後平均ｍ_ｊ，Ｂ（ｋ，ｆ）の更新、及び（３４４）式の事後共分散行列Σ_{｛ｊ，ｌ｝}（ｋ，ｆ）の更新において、行列式計算、逆行列計算及び行列乗算を時間周波数点ごとに行う必要があった。具体的には、（３４０）式における（３５０）式の形の行列式計算と、（３４０）式及び（３４２）式及び（３４４）式における（３５１）式の形の逆行列計算と、（３４４）式における（３５２）式の形の行列乗算と、を時間周波数点（ｋ，ｆ）ごとに行う必要があった。

行列式計算、逆行列計算及び行列乗算はいずれもＯ（Ｉ^３）のオーダーの計算量を必要とする行列演算であり、かつ時間周波数点の数は通常、非常に大きいため、この従来技術では膨大な計算量を要していた。

いま仮に、（３５０）式〜（３５２）式における共分散行列＾ｖ_ｊ（ｋ，ｆ）＾Ｒ_ｊ（ｆ）及び＾ｖ_ｌ（ｋ，ｆ）＾Ｒ_ｌ（ｆ）がともに対角行列であるとすると、（３５０）式〜（３５２）式における行列式計算、逆行列計算及び行列乗算は単なる対角要素に対するスカラー演算に帰着するため、少ない計算量で計算できる。しかしながら、通常、源信号成分（ベクトル）の異なる要素は互いに相関を持つため、共分散行列＾ｖ_ｊ（ｋ，ｆ）＾Ｒ_ｊ（ｆ）及び＾ｖ_ｌ（ｋ，ｆ）＾Ｒ_ｌ（ｆ）は対角行列ではない。

そこで、本第７の実施形態では、一般化固有値問題に基づいて２つの空間共分散行列＾Ｒ_ｊ（ｆ），＾Ｒ_ｌ（ｆ）を同時対角化することを考える。λ_{１，｛ｊ，ｌ｝}（ｆ），・・・，λ_{Ｉ，｛ｊ，ｌ｝}（ｆ）を空間共分散行列＾Ｒ_ｌ（ｆ），＾Ｒ_ｊ（ｆ）に対する一般化固有値問題（（３５３）式）の一般化固有値とし、ｐ_{１，｛ｊ，ｌ｝}（ｆ），・・・，ｐ_{Ｉ，｛ｊ，ｌ｝}（ｆ）を一般化固有値λ_{１，｛ｊ，ｌ｝}（ｆ），・・・，λ_{Ｉ，｛ｊ，ｌ｝}（ｆ）に対応する一般化固有ベクトルとする。

ただし、一般化固有ベクトルｐ_{１，｛ｊ，ｌ｝}（ｆ），・・・，ｐ_{Ｉ，｛ｊ，ｌ｝}（ｆ）は、次の（３５４）式を満たすように取る。

ｐ_{｛ｊ，ｌ｝}（ｆ）を一般化固有ベクトルｐ_{１，｛ｊ，ｌ｝}（ｆ），・・・，ｐ_{Ｉ，｛ｊ，ｌ｝}（ｆ）を列ベクトルとする行列、Λ_{｛ｊ，ｌ｝}（ｆ）を一般化固有値λ_{１，｛ｊ，ｌ｝}（ｆ），・・・，λ_{Ｉ，｛ｊ，ｌ｝}（ｆ）を対角要素とする対角行列とすると、次の（３５５）式が成り立つ。

（３５５）式は、行列Ｐ_{｛ｊ，ｌ｝}（ｆ）が２つの空間共分散行列＾Ｒ_ｊ（ｆ），＾Ｒ_ｌ（ｆ）を同時対角化する（２つの源信号成分ｘ_ｊ（ｋ，ｆ），ｘ_ｌ（ｋ，ｆ）を同時無相関化する）ことを示している。そこで、行列Ｐ_{｛ｊ，ｌ｝}（ｆ）を同時無相関化行列と呼ぶ。

（３５５）式を空間共分散行列＾Ｒ_ｊ（ｆ），＾Ｒ_ｌ（ｆ）について解いて、（３４０）式、（３４２）式及び（３４４）式に代入し、（３５６）式〜（３５８）式と定義する（ｊ，ｌは１以上Ｌ以下の互いに異なる任意の整数）と、（３６２）、（３６３）式及び（３６４）式（後述）を得る。

（３５６）式及び（３５７）式を事後平均ｍ_{ｊ，｛ｊ，ｌ｝}（ｋ，ｆ）及び事後共分散行列Σ_{｛ｊ，ｌ｝}（ｋ，ｆ）について解き、（３４８）式及び（３４９）式に代入すると（３６６）式及び（３６８）式を得る。

これで、本第７の実施形態に係る信号分析装置における更新則が導出された。

このように、本第７の実施形態に係る信号分析装置では、信号パラメータ分析部７２０は、２個の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列＾Ｒ_ｌ（ｆ），＾Ｒ_ｊ（ｆ）に対するｐ_{ｉ，｛ｊ，ｌ｝}（ｆ）^Ｈ＾Ｒ_ｊ（ｆ）ｐ_{ｍ，｛ｊ，ｌ｝}（ｆ）＝０（ｉ，ｍは１以上Ｉ以下の互いに異なる任意の整数）を満たす、それぞれ異なるＩ個（Ｉは２以上の整数）の一般化固有ベクトルｐ_{１，｛ｊ，ｌ｝}（ｆ），・・・，ｐ_{Ｉ，｛ｊ，ｌ｝}（ｆ）を列ベクトルとする行列である同時無相関化行列Ｐ_{｛ｊ，ｌ｝}（ｆ）またはそのエルミート転置Ｐ_{｛ｊ，ｌ｝}（ｆ）^Ｈを、Ｉ個の位置で取得される観測信号による観測信号ベクトルについて２個の源信号に対応する成分を無相関化するためのパラメータとして得る。

［第７の実施形態の実施例］
本実施例では、集合系Ｈの具体例について説明する。本実施例では、集合系Ｈを（３５９）式とする。

これは、各時間周波数点において観測信号ベクトルは、１番目の源信号成分をつねに含み、２番目〜Ｊ番目の源信号成分のうちのちょうど１つを含むことを意味する。例えば、１番目の源信号成分は背景雑音のようなスパースでない源信号成分をモデル化し、２番目〜Ｊ番目の源信号成分は音声のようなスパースな源信号成分をモデル化する。

信号パラメータ分析部７２０、同時無相関化部７３０、信号存在事後確率更新部７４０、シングルチャネルウィーナーフィルタ部７５０、事後共分散行列更新部７６０、混合重み更新部７７０、分散パラメータ更新部７８０、空間共分散行列更新部７９０及び源信号成分推定部８００における処理の更新式を（３６０）式〜（３７０）式に示す。本第７の実施形態では、Ｈの元である集合｛１，２｝，｛１，３｝，・・・，｛１，Ｊ｝は、それぞれインデックス２，３，・・・，Ｊと１対１に対応しているから、簡単のため以下の記法を用いる。すなわち、ｙ´_{｛１，ｊ｝}（ｋ，ｆ）をｙ´_ｊ（ｋ，ｆ）と書き、Ｐ_{｛１，ｊ｝}（ｆ）をＰ_ｊ（ｆ）と書き、γ_{｛１，ｊ｝}（ｋ，ｆ）をγ_ｊ（ｋ，ｆ）と書き、＾α_{｛１，ｊ｝}（ｆ）を＾α_ｊ（ｆ）と書き、Λ_{ｌ，｛１，ｊ｝}（ｆ）をΛ_ｌｊ（ｆ）と書き、ｍ´_{ｌ，｛１，ｊ｝}（ｋ，ｆ）をｍ´_ｌｊ（ｋ，ｆ）と書き、Σ´_{｛１，ｊ｝}（ｋ，ｆ）をΣ´_ｊ（ｋ，ｆ）と書く。ただし、ｊ＝２，・・・，Ｊとし、ｌ＝１，ｊとする（下記の更新式でも同様）。なお、本実施例では、信号パラメータ分析部７２０は、Ｌ個の源信号成分のうちの任意の２個ｘｊ（ｋ，ｆ），ｘｌ（ｋ，ｆ）（ｊ，ｌは１以上Ｌ以下の相異なる任意の自然数）に対して、同時無相関化行列Ｐ｛ｊ，ｌ｝（ｆ）と対角行列Λｊ，｛ｊ，ｌ｝（ｆ）を計算する必要はなく、（３６０）式にあるように、ｘ１（ｋ，ｆ），ｘｊ（ｋ，ｆ）（ｊは２以上Ｌ以下の任意の自然数）に対して、同時無相関化行列Ｐ｛１，ｊ｝（ｆ）（すなわちＰｊ（ｆ））と対角行列Λｌ，｛１，ｊ｝（ｆ）（すなわちΛｌｊ（ｆ））を計算すれば十分である。

［第８の実施形態］
本第８の実施形態では、第３の実施形態に係る信号分析装置の他の構成例について説明する。図２１は、第８の実施形態に係る信号分析装置の構成の一例を示す図である。図２２は、第８の実施形態に係る信号分析処理の処理手順の一例を示すフローチャートである。図２１に示すように、第８の実施の形態に係る信号分析装置８０１は、周波数領域変換部３１０、周波数領域信号分離部３２０−１及び時間領域変換部３３０を有する。

信号分析装置８０１では、周波数領域信号分離部３２０−１において、スペクトル基底ｔ_ｎ（ｆ）と、アクティベーション基底ｗ_ｎ（ｋ）と、分配関数ｚ_ｎｊと、同時無相関化された空間共分散行列Λ_ｊ（ｆ）と、同時無相関化行列Ｐ（ｆ）との推定値とを計算する。

スペクトル基底ｔ_ｎ（ｆ）と、アクティベーション基底ｗ_ｎ（ｋ）と、分配関数ｚ_ｎｊとは、分散パラメータｖ_ｊ（ｋ，ｆ）をモデル化するパラメータである。分散パラメータｖ_ｊ（ｋ，ｆ）は、これらのパラメータを用いて（３７１）式によりモデル化される。

ここで、ｎは、スペクトル基底およびアクティベーション基底のインデックスである。Ｎは、スペクトル基底およびアクティベーション基底の個数である。

周波数領域信号分離部３２０−１は、初期化部（図示しない）、観測信号ベクトル作成部８１０、同時無相関化８３０、シングルチャネルウィーナーフィルタ部８４０、観測共分散行列更新部８５０、分散パラメータ更新部８６０、空間共分散行列更新部８７０、信号パラメータ分析部８２０、収束判定部（図示しない）及び同時無相関化逆変換部８８０を有する。

［信号分析装置の構成及び処理］
図２１及び図２２を用いて、第８の実施形態に係る信号分析装置の構成及び信号分析処理の処理手順について説明する。

図示しない初期化部は、分散パラメータｖ_ｊ（ｋ，ｆ）と、スペクトル基底ｔ_ｎ（ｆ）と、アクティベーション基底ｗ_ｎ（ｋ）と、分配関数ｚ_ｎｊと、同時無相関化された空間共分散行列Λ_ｊ（ｆ）と、同時無相関化行列Ｐ（ｆ）と、の初期推定値を計算する初期化を行う（ステップＳ１２１）。初期推定値は、上付きの（０）を付して、^ｖ_ｊ ^（０）（ｋ，ｆ）などと表す。

周波数領域変換部３１０は、観測信号を取得し（ステップＳ１２２）、観測信号を周波数領域に変換する（ステップＳ１２３）ことにより、時間周波数領域の観測信号ｙ_ｉ（ｋ，ｆ）を計算する。

観測信号ベクトル作成部８１０は、観測信号ベクトルｙ（ｋ，ｆ）を第７の実施形態と同様にして作成する（ステップＳ１２４）。

同時無相関化部８３０は、同時無相関化された観測信号ベクトルｙ´（ｋ，ｆ）を（３７２）式により計算し、更新する（ステップＳ１２５）。

ここで、＾Ｐ（ｆ）は信号パラメータ分析部８２０からの同時無相関化行列の推定値である。

観測共分散行列更新部８５０は、同時無相関化された観測共分散行列Ｕ´（ｋ，ｆ）を（３７３）式により更新する（ステップＳ１２６）。

ここで、＾ｖ_ｊ（ｋ，ｆ）は、分散パラメータ更新部８６０からの分散パラメータの推定値である。＾Λ_ｊ（ｆ）は、空間共分散行列更新部８７０からの同時無相関化された空間共分散行列の推定値である。

分散パラメータ更新部８６０は、次の手順により分散パラメータの推定値＾ｖ_ｊ（ｋ，ｆ）を更新する（ステップＳ１２７）。

まず、分散パラメータ更新部８６０は、ｎ番目の基底に対応する同時無相関化された空間共分散行列Σ_ｎ（ｆ）を（３７４）式により更新する。（３７４）式において、＾ｚ_ｎｊは、分配関数の推定値である。

次に、分散パラメータ更新部はスペクトル基底の推定値＾ｔ_ｎ（ｆ）を（３７５）式により更新する。（３７５）式において、＾ｗ_ｎ（ｋ）はアクティベーション基底の推定値である。

次に、分散パラメータ更新部８６０は、アクティベーション基底の推定値＾ｗ_ｎ（ｋ）を（３７６）式により更新する。

次に、分散パラメータ更新部８６０は、分配関数の推定値＾ｚ_ｎｊを（３７７）式により更新する。

最後に、分散パラメータ更新部８６０は、分散パラメータの推定値＾ｖ_ｊ（ｋ，ｆ）を（３７８）式により更新する。

空間共分散行列更新部８７０は、（３７９）式により、同時無相関化された空間共分散行列の推定値＾Λ_ｊ（ｆ）を更新する（ステップＳ１２８）。

ここで、＃は、行列の幾何平均であり、Ａ_ｊ（ｆ）およびＢ_ｊ（ｆ）は、（３８０）式および（３８１）式で計算される行列である。

信号パラメータ分析部８２０は、同時無相関化行列の推定値＾Ｐ（ｆ）を（３８２）式により更新する（ステップＳ１２９）。

図示しない収束判定部は、収束したかどうかの判定を行う（ステップＳ１３０）。収束判定部が収束していないと判定した場合（ステップＳ１３０：Ｎｏ）、同時無相関化部８３０での処理（ステップＳ１２５）に戻って、処理が継続される。

一方、収束判定部が収束したと判定した場合（ステップＳ１３０：Ｙｅｓ）、シングルチャネルウィーナーフィルタ部８４０は、同時無相関化された源信号成分の推定値＾ｘ´_ｊ（ｋ，ｆ）を（３８３）式により計算する（ステップＳ１３１）。

シングルチャネルウィーナーフィルタ部８４０は、同時無相関化された源信号成分の推定値＾ｘ´_ｊ（ｋ，ｆ）を更新する（ステップＳ１３２）。

同時無相関化逆変換部８８０は、源信号成分の推定値＾ｘ_ｊ（ｋ，ｆ）を（３８４）式により計算する（ステップＳ１３３）。

時間領域変換部３３０は、第７の実施形態と同様の処理を行って、時間領域の源信号成分の推定値＾ｘ_ｊ（ｔ）を計算し、出力する（ステップＳ１３４）。

［第８の実施形態の効果］
第３の実施形態では、音の方向性のみを基に信号処理を行っていた。これに対して、本第８の実施形態では、音の構造（パーツ）も考慮して信号分析を行うため、音源分離の精度を向上させることができる。

そして、第８の実施形態では、パーミテーション問題が原理的に起こらないという効果を奏する。本第８の実施形態では、全ての周波数において、どのようなパーツがあるかを考慮した処理を行うため、周波数ごとに独立に処理しないことから、パーミテーション問題が原理的に起こらない。

［システム構成等］
また、図示した各装置の各構成要素は機能概念的なものであり、必ずしも物理的に図示のように構成されていることを要しない。すなわち、各装置の分散・統合の具体的形態は図示のものに限られず、その全部又は一部を、各種の負荷や使用状況等に応じて、任意の単位で機能的又は物理的に分散・統合して構成することができる。さらに、各装置にて行われる各処理機能は、その全部又は任意の一部が、ＣＰＵ及び当該ＣＰＵにて解析実行されるプログラムにて実現され、あるいは、ワイヤードロジックによるハードウェアとして実現され得る。

また、本実施形態において説明した各処理のうち、自動的に行われるものとして説明した処理の全部又は一部を手動的に行うこともでき、あるいは、手動的に行われるものとして説明した処理の全部又は一部を公知の方法で自動的に行うこともできる。この他、上記文書中や図面中で示した処理手順、制御手順、具体的名称、各種のデータやパラメータを含む情報については、特記する場合を除いて任意に変更することができる。すなわち、上記学習方法及び音声認識方法において説明した処理は、記載の順にしたがって時系列に実行されるのみならず、処理を実行する装置の処理能力あるいは必要に応じて並列的にあるいは個別に実行されてもよい。

［プログラム］
図２３は、プログラムが実行されることにより、信号分析装置１，２０１，３０１，４０１，５０１，６０１が実現されるコンピュータの一例を示す図である。コンピュータ１０００は、例えば、メモリ１０１０、ＣＰＵ１０２０を有する。また、コンピュータ１０００は、ハードディスクドライブインタフェース１０３０、ディスクドライブインタフェース１０４０、シリアルポートインタフェース１０５０、ビデオアダプタ１０６０、ネットワークインタフェース１０７０を有する。これらの各部は、バス１０８０によって接続される。

メモリ１０１０は、ＲＯＭ１０１１及びＲＡＭ１０１２を含む。ＲＯＭ１０１１は、例えば、ＢＩＯＳ（Basic Input Output System）等のブートプログラムを記憶する。ハードディスクドライブインタフェース１０３０は、ハードディスクドライブ１０９０に接続される。ディスクドライブインタフェース１０４０は、ディスクドライブ１１００に接続される。例えば磁気ディスクや光ディスク等の着脱可能な記憶媒体が、ディスクドライブ１１００に挿入される。シリアルポートインタフェース１０５０は、例えばマウス１１１０、キーボード１１２０に接続される。ビデオアダプタ１０６０は、例えばディスプレイ１１３０に接続される。

ハードディスクドライブ１０９０は、例えば、ＯＳ（Operating System）１０９１、アプリケーションプログラム１０９２、プログラムモジュール１０９３、プログラムデータ１０９４を記憶する。すなわち、信号分析装置１，２０１，３０１，４０１，５０１，６０１の各処理を規定するプログラムは、コンピュータ１０００により実行可能なコードが記述されたプログラムモジュール１０９３として実装される。プログラムモジュール１０９３は、例えばハードディスクドライブ１０９０に記憶される。例えば、信号分析装置１，２０１，３０１，４０１，５０１，６０１における機能構成と同様の処理を実行するためのプログラムモジュール１０９３が、ハードディスクドライブ１０９０に記憶される。なお、ハードディスクドライブ１０９０は、ＳＳＤ（Solid State Drive）により代替されてもよい。

また、上述した実施形態の処理で用いられる設定データは、プログラムデータ１０９４として、例えばメモリ１０１０やハードディスクドライブ１０９０に記憶される。そして、ＣＰＵ１０２０が、メモリ１０１０やハードディスクドライブ１０９０に記憶されたプログラムモジュール１０９３やプログラムデータ１０９４を必要に応じてＲＡＭ１０１２に読み出して実行する。

なお、プログラムモジュール１０９３やプログラムデータ１０９４は、ハードディスクドライブ１０９０に記憶される場合に限らず、例えば着脱可能な記憶媒体に記憶され、ディスクドライブ１１００等を介してＣＰＵ１０２０によって読み出されてもよい。あるいは、プログラムモジュール１０９３及びプログラムデータ１０９４は、ネットワーク（ＬＡＮ（Local Area Network）、ＷＡＮ（Wide Area Network）等）を介して接続された他のコンピュータに記憶されてもよい。そして、プログラムモジュール１０９３及びプログラムデータ１０９４は、他のコンピュータから、ネットワークインタフェース１０７０を介してＣＰＵ１０２０によって読み出されてもよい。

以上、本発明者によってなされた発明を適用した実施形態について説明したが、本実施形態による本発明の開示の一部をなす記述及び図面により本発明は限定されることはない。すなわち、本実施形態に基づいて当業者等によりなされる他の実施形態、実施例及び運用技術等はすべて本発明の範疇に含まれる。

１，２０１，３０１，４０１，５０１，６０１，８０１ 信号分析装置
１０，２１０，４１０，５１０，６１０，７１０，８１０ 観測信号ベクトル作成部
２０，２２０，４２０，５２０，６２０，７２０，８２０ 信号パラメータ分析部
３０，２３０，４３０，５３０，６３０，７３０，８３０ 同時無相関化部
４０，２４０，５４０，６４０，７５０，８４０ シングルチャネルウィーナーフィルタ部
５０，２８０，５８０，６８０，８８０ 同時無相関化逆変換部
２５０，５５０，７６０ 事後共分散行列更新部
２６０，４４０，５６０，７８０，８６０ 分散パラメータ更新部
２７０，４５０，５７０，７９０，８７０ 空間共分散行列更新部
３１０ 周波数領域変換部
３２０ 周波数領域信号分離部
３３０ 時間領域変換部
６５０，８５０ 観測共分散行列更新部
７４０ 信号存在事後確率更新部
７７０ 混合重み更新部
８００ 源信号成分推定部

Claims (20)

1. Ｊ個（Ｊは２以上の整数）の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列になるような行列である同時無相関化行列Ｐまたは／及び、そのエルミート転置Ｐ^Ｈを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するためのパラメータとして得る信号パラメータ分析部
   を有することを特徴とする信号分析装置。
2. 請求項１に記載の信号分析装置であって、
   前記信号パラメータ分析部は、前記空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列に近づくように求めた行列である同時無相関化行列Ｐまたは／及びそのエルミート転置Ｐ^Ｈを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するためのパラメータとして得ることを特徴とする信号分析装置。
3. 請求項１に記載の信号分析装置であって、
   前記信号パラメータ分析部は、不動点法、勾配法、自然勾配法、または、ニュートン法に基づいて、同時無相関化行列Ｐの推定値を更新することを特徴とする信号分析装置。
4. 請求項１に記載の信号分析装置であって、
   前記信号パラメータ分析部が得た前記同時無相関化行列Ｐのエルミート転置Ｐ^Ｈを、入力されたＩ個の観測信号による観測信号ベクトルに掛けることにより、Ｊ個の源信号に対応する成分が無相関化された同時無相関化済観測信号ベクトルを得る同時無相関化部
   をさらに有することを特徴とする信号分析装置。
5. 請求項４に記載の信号分析装置であって、
   前記同時無相関化部が得た前記同時無相関化済観測信号ベクトルの各要素に、
   前記同時無相関化行列Ｐに基づく行列Ｐ^ＨＲ_ｊＰ（ｊは１以上Ｊ以下の整数）の対角成分に基づくシングルチャネルウィーナーフィルタを適用して同時無相関化済源信号成分推定値によるベクトルを得るウィーナーフィルタ部
   をさらに有することを特徴とする信号分析装置。
6. 請求項５に記載の信号分析装置であって、
   前記ウィーナーフィルタ部が得た前記同時無相関化済源信号成分推定値によるベクトルに、前記同時無相関化行列Ｐのエルミート転置Ｐ^Ｈの逆行列を掛けることにより、源信号成分の推定値を得る同時無相関化逆変換部
   をさらに有することを特徴とする信号分析装置。
7. ２個の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_１，Ｒ_２に対するｐ_ｉ ^ＨＲ_２ｐ_ｍ＝０（ｉ，ｍは１以上Ｉ以下の互いに異なる任意の整数）を満たす、それぞれ異なるＩ個（Ｉは２以上の整数）の一般化固有ベクトルｐ_１，・・・，ｐ_Ｉを列ベクトルとする行列である同時無相関化行列Ｐまたはそのエルミート転置Ｐ^Ｈを、前記Ｉ個の位置で取得される観測信号による観測信号ベクトルについて２個の源信号に対応する成分を無相関化するためのパラメータとして得る信号パラメータ分析部
   を有することを特徴とする信号分析装置。
8. 請求項７に記載の信号分析装置であって、
   前記信号パラメータ分析部が得た前記同時無相関化行列Ｐのエルミート転置Ｐ^Ｈを、入力されたＩ個の観測信号による観測信号ベクトルに掛けることにより、２個の源信号に対応する成分が無相関化された同時無相関化済観測信号ベクトルを得る同時無相関化部
   をさらに有することを特徴とする信号分析装置。
9. 請求項８に記載の信号分析装置であって、
   前記同時無相関化部が得た前記同時無相関化済観測信号ベクトルの各要素に、
   前記同時無相関化行列Ｐに基づく行列Ｐ^ＨＲ_ｊＰ（ｊは１または２）の対角成分に基づくシングルチャネルウィーナーフィルタを適用して同時無相関化済源信号成分推定値によるベクトルを得るウィーナーフィルタ部
   をさらに有することを特徴とする信号分析装置。
10. 請求項９に記載の信号分析装置であって、
    前記ウィーナーフィルタ部が得た前記同時無相関化済源信号成分推定値によるベクトルに、前記同時無相関化行列Ｐのエルミート転置Ｐ^Ｈの逆行列を掛けることにより、源信号成分の推定値を得る同時無相関化逆変換部
    をさらに有することを特徴とする信号分析装置。
11. Ｊ個（Ｊは２以上の整数）の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列になるような行列Ｐを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するための同時無相関化行列として、前記同時無相関化行列のエルミート転置を、前記観測信号ベクトルに掛けることにより、Ｊ個の源信号に対応する成分が無相関化された同時無相関化済観測信号ベクトルを得る同時無相関化部
    を有することを特徴とする信号分析装置。
12. 請求項１１に記載の信号分析装置であって、
    前記同時無相関化部が得た前記同時無相関化済観測信号ベクトルの各要素に、
    前記同時無相関化行列Ｐに基づく行列Ｐ^ＨＲ_ｊＰ（ｊは１以上Ｊ以下の整数）の対角成分に基づくシングルチャネルウィーナーフィルタを適用して同時無相関化済源信号成分推定値によるベクトルを得るウィーナーフィルタ部
    をさらに有することを特徴とする信号分析装置。
13. 請求項１２に記載の信号分析装置であって、
    前記ウィーナーフィルタ部が得た前記同時無相関化済源信号成分推定値によるベクトルに、前記同時無相関化行列Ｐのエルミート転置Ｐ^Ｈの逆行列を掛けることにより、源信号成分の推定値を得る同時無相関化逆変換部
    をさらに有することを特徴とする信号分析装置。
14. Ｊ個（Ｊは２以上の整数）の混在する源信号から得た同時無相関化された共分散行列に基づいて、源信号の空間的特性をモデル化するパラメータである空間共分散行列を求める空間共分散行列計算部、または／及び前記同時無相関化された共分散行列に基づいて、源信号のサンプル点ごとの分散をモデル化するパラメータである分散パラメータを求める分散パラメータ計算部
    を有することを特徴とする信号分析装置。
15. 請求項１４に記載の信号分析装置であって、
    補助関数法を用いて、前記空間共分散行列または／及び分散パラメータを求めることを特徴とする信号分析装置。
16. 信号分析装置が実行する信号分析方法であって、
    Ｊ個（Ｊは２以上の整数）の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列になるような行列である同時無相関化行列Ｐまたは／及び、そのエルミート転置Ｐ^Ｈを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するためのパラメータとして得る工程
    を含んだことを特徴とする信号分析方法。
17. 信号分析装置が実行する信号分析方法であって、
    ２個の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_１，Ｒ_２に対するｐ_ｉ ^ＨＲ_２ｐ_ｍ＝０（ｉ，ｍは１以上Ｉ以下の互いに異なる任意の整数）を満たす、それぞれ異なるＩ個（Ｉは２以上の整数）の一般化固有ベクトルｐ_１，・・・，ｐ_Ｉを列ベクトルとする行列である同時無相関化行列Ｐまたはそのエルミート転置Ｐ^Ｈを、前記Ｉ個の位置で取得される観測信号による観測信号ベクトルについて２個の源信号に対応する成分を無相関化するためのパラメータとして得る工程
    を含んだことを特徴とする信号分析方法。
18. 信号分析装置が実行する信号分析方法であって、
    Ｊ個（Ｊは２以上の整数）の混在する源信号それぞれの空間的特性をモデル化する空間共分散行列Ｒ_ｊ（ｊは１以上Ｊ以下の整数）について、Ｐ^ＨＲ_ｊＰがすべて対角行列になるような行列Ｐを、それぞれ異なるＩ個（Ｉは２以上の整数）の位置で取得される観測信号による観測信号ベクトルについてＪ個の源信号に対応する成分を無相関化するための同時無相関化行列として、前記同時無相関化行列のエルミート転置を、前記観測信号ベクトルに掛けることにより、Ｊ個の源信号に対応する成分が無相関化された同時無相関化済観測信号ベクトルを得る工程
    を含んだことを特徴とする信号分析方法。
19. 信号分析装置が実行する信号分析方法であって、
    Ｊ個（Ｊは２以上の整数）の混在する源信号から得た同時無相関化された共分散行列に基づいて、源信号の空間的特性をモデル化するパラメータである空間共分散行列を求める工程
    を含んだことを特徴とする信号分析方法。
20. コンピュータを、請求項１〜１５のいずれか一つに記載の信号分析装置として機能させるための信号分析プログラム。

Images