JP2007304445A - 周波数成分の修復・抽出方法、周波数成分の修復・抽出装置、周波数成分の修復・抽出プログラムならびに周波数成分の修復・抽出プログラムを記録した記録媒体 - Google Patents

Abstract

【課題】入力音響信号の消えた周波数の成分を修復しながら、元の音響信号の周波数の成分を分離抽出する。
【解決手段】入力音響信号に対し、周波数解析を行い、スペクトログラムｘ（ｉ，ｊ）を算出し（ｉは時間、ｊは周波数）、ｘ（ｉ，ｊ）の各位置について、ベクトル量ｘ^ｋ（ｋ＝１〜Ｋ）の各要素を計算し、観測値行列集合Ｘ^ｋを生成し、周波数成分状態行列Θ^Ｆの各位置について、要素（ｉ，ｊ）を除くその前後の要素と観測値行列を特徴ベクトルｆ（ｉ，ｊ）として出力し、各位置において、特徴ベクトルｆ（ｉ，ｊ）を与えた時に任意の周波数成分状態が出現する事後確率を計算し、ｘ（ｉ，ｊ）の各位置の事後確率を用いて事後同時確率を増大させるΘ^Ｆを推定し、このΘ^Ｆより入力音響信号の周波数成分を出力する。
【選択図】図１

Description

この発明は、複数の音響信号や雑音が混在した入力信号のスペクトログラムから各音響信号の成分を抽出すると同時に失われた成分を修復する周波数成分の修復・抽出方法、装置、プログラム、および記録媒体に関するものである。

観測された混合信号を分離して、混ざり合う前のオリジナルの信号を推定するブラインド音源分離技術がある。この技術は、複数の信号源から発せられる未知の原信号から混合された混合信号を分離し、復元する技術である。ブラインド音源分離技術の中の有効な方法に、統計的に独立な成分を分解抽出する独立成分分析法（ＩＣＡ：Ｉｎｄｅｐｅｎｄｅｎｔ Ｃｏｍｐｏｎｅｎｅｔ Ａｎａｌｙｓｉｓ）がある。この手法によれば、例えば、複数の音信号をＡ個のマイクロホン（Ａは整数とする）で受信し、各音源から到来する音信号同士が統計的に独立であることを利用して前記の受信信号を処理することで、マイクロホンと同数のＡ個もしくはＡ個以下の音源からの信号を分離することができる。

また、独立成分分析法に類似の方法として特許文献１に示す周波数成分抽出装置がある。この装置は、周波数成分抽出に対して、合理的な制約条件を導入し、その制約条件のものとでの周波数成分の状態の最適化を図ることによって、制約に反した周波数成分すなわち雑音の影響を軽減した周波数成分抽出を行うものである。

また、画像の修復技術として、非特許文献１に示すマルコフ確率場とギブス・サンプラーを使う方法があるが、条件付確率を手作業で設計しなければならず、そのまま音響信号に適用しても、精度の高い修復・抽出を行うことが難しいという問題があった。
特許第３２６２２０４号 Ｓ．Ｇｅｍａｎ ａｎｄ Ｄ．Ｇｅｍａｎ．Ｓｔｏｃｈａｓｔｉｃ ｒｅｌａｘａｔｉｏｎ，Ｇｉｂｂｓ ｄｉｓｔｒｉｂｕｔｉｏｎｓ ａｎｄ ｔｈｅ Ｂａｙｅｓｉａｎ ｒｅｓｔｏｒａｔｉｏｎ ｏｆ ｉｍａｇｅｓ．ＩＥＥＥ Ｔｒａｎｓａｃｔｉｏｎｓ ｏｎ Ｐａｔｔｅｒｎ Ａｎａｌｙｓｉｓ ａｎｄ Ｍａｃｈｉｎｅ Ｉｎｔｅｌｌｉｇｅｎｃｅ，ＰＡＭＩ−６（６）：７２１−７４１，１９８４

上記の独立成分分析法では、少なくとも、信号源の個数と等しい数以上の入力信号が必要であり、通信断、雑音などにより、消えてしまった信号を修復することは出来ないという問題があった。
また、上記の周波数成分抽出装置は、着目時間位置ｉの近傍として、前の時間位置ｉ−１のみしか利用しておらず、着目時間位置ｉの後の近傍時間ｉ＋１以後については、観測信号を利用できず、それだけ正しく周波数成分を修復・抽出することができなかった。

入力音響信号に対して周波数解析を行い、サウンドスペクトログラム（以下、スペクトログラムという）を計算して、時間×周波数の実数値行列を求め、実数値行列で表現されたスペクトログラムから、スペクトログラムの各位置（ｉ、ｊ）について、Ｋ個の特徴量ｘ^ｋ（ｉ、ｊ）（ｋ＝１、．．．、Ｋ）の各要素を計算し、観測値行列集合［Ｘ^１、Ｘ^２、．．．、Ｘ^Ｋ］を生成する。ここで（ｉ、ｊ）はスペクトログラムの時間ｉと周波数ｊの行列の要素、Ｋは特徴数である。周波数成分の存在の確からしさを示す離散値を表す行列（以下、周波数成分状態行列Θ^Ｆという）の上記各位置（ｉ、ｊ）について、その要素（ｉ、ｊ）を除くその前後の要素と上記観測値行列を特徴ベクトルｆ（ｉ、ｊ）として出力し、スペクトログラム上の各位置（ｉ、ｊ）において、上記特徴ベクトルｆ（ｉ、ｊ）を与えたときに、任意の周波数成分状態が出現する確率を計算し、上記スペクトログラムの各位置の上記事後確率を用いて、事後同時確率を増大させる周波数成分状態行列を推定し、上記周波数成分状態行列の推定値より入力音響信号の周波数成分を出力する。

着目時間位置ｉの前のみならず後の近傍時間位置ｉ＋１以後の観測信号をも利用しているため、種々の音響信号に対し、従来より正確に周波数成分を修復・抽出することが出来る。

実施例１
この発明のハードウェア構成例を図１に示し、この実施例で説明する信号等の具体的形態図を図２に示し、この発明の主要な処理の流れを図１３に示す。特徴抽出手段４は、スペクトログラム生成手段６と観測値行列生成手段１０とで構成されている。観測値行列生成手段１０は、第ｋ特徴量生成手段１２ｋと（ｋ＝１、．．．、Ｋ）と第ｋ観測値行列記憶手段１４ｋとで構成されている。

複数の音響信号や雑音が混在した入力音声信号がサンプリングされ、ディジタル信号化された状態で、入力端子２に入力される。この入力音声信号の具体的波形例を図２Ａに示す。入力された音声信号は、時間ｉについて、スペクトログラム生成手段６で公知の周波数解析方法を用いて、信号波形のパワースペクトルｘ（ｉ、ｊ）（ｉ＝１、．．．、Ｉ、ｊ＝１、．．．、Ｊ）を導出する。ここで、ｉは時間を表し、ｊは周波数の数を表し、公知の周波数解析方法とは、離散フーリエ変換（ＤＦＴ：Ｄｉｓｃｒｅｔｅ Ｆｏｕｒｉｅｒ Ｔｒａｎｓｆｏｒｍ）や高速フーリエ変換（ＦＦＴ：Ｆａｓｔ Ｆｏｕｒｉｅｒ Ｔｒａｎｓｆｏｒｍ）等が考えられる。このようにして、入力信号のサウンドスペクトログラム（以下、スペクトログラムという）が生成される（Ｓ３００）。

スペクトログラムの具体的形態例を図２Ｂに示し、その時間、周波数の位置（ｉ、ｊ）の関係の概念図を図３に示す。図２Ｂにおいて、黒の濃度が高い部分がパワー値が大きいことを表し、特に、濃度の強い黒でノイズ部分を表している。時間ｉ、周波数ｊの部分を位置（ｉ、ｊ）と表す。非特許文献１に示す技術では、画像の各ピクセルが１つの属性（白色、黒色、またはグレー色）で表現されていたのに対し、この発明では、各位置に複数の属性（パワー値、ピークの有無等）を特徴として与えることが出来る。このスペクトログラムの各位置（ｉ、ｊ）のパワースペクトルｘ（ｉ、ｊ）はそれぞれ第ｋ特徴量生成手段１２ｋ（ｋ＝１、．．．、Ｋ）に入力される。ただし、Ｋは特徴数を表す。第ｋ特徴量生成手段１２ｋでスペクトログラムｘ（ｉ、ｊ）から特徴量ｘ^ｋ（ｉ、ｊ）（ｋ＝１、．．．、Ｋ）を生成し、第ｋ観測値行列記憶手段１４ｋに特徴量ｘ^ｋ（ｉ、ｊ）を記憶する。

この処理をスペクトログラムの全ての時間について行い、第ｋ観測値行列記憶手段１４ｋに記憶する。その際に、各特徴量ｘ^ｋ（ｉ、ｊ）がｉ行ｊ列の要素となる第ｋ観測値行列Ｘ^ｋとして記憶し、生成する（ステップＳ３０２）。これらＫ個の観測値行列記憶手段１４ｋ（ｋ＝１、．．．、Ｋ）により、観測値行列集合［Ｘ^１、．．．、Ｘ^ｋ、．．．、Ｘ^Ｋ］を記憶する観測値行列集合記憶手段１６が構成される。

また、音声信号の周波数特性によるパワー損失を補填するためにパワースペクトルの周波数方向の正規化を行い、正規化により得られた値を要素とする入力信号の特徴を抽出する構成も考えられる。この場合は、図１に破線で示すように正規化手段８を設け、スペクトログラム生成手段６で計算されたパワースペクトルｘ（ｉ、ｊ）が正規化手段８により、正規化される。

正規化の具体的手法を説明すると、パワースペクトルｘ（ｉ、ｊ）の値が減衰した周波数帯域におけるパワー値を補正する。例えば、各位置（ｉ、ｊ）の値ｘ（ｉ、ｊ）を、その周波数軸方向の±Ｎ_Ｆの位置のパワースペクトルの最大値と最小値を使って０〜１の値を持つように、正規化する方法などが考えられる。具体的に説明すると、音声信号が周波数領域に変換された場合の、波形例を図４Ａに示す。
図４Ａにおいて、横軸は周波数、横軸はパワー値を表す。周波数の範囲ｊ±Ｎ_Ｆ
内で、最大レベルを「１」とし、最小レベルを「０」として、周波数ｊの値を０〜１の値で表す。このようにすれば、図４Ａに示す波形は図４Ｂに示すようになる。ピークになっている周波数をｊとすると、図４Ｂは図４Ａと同様、横軸が周波数方向、縦軸がレベル（パワー値）を表す。なお、この正規化パワースペクトルもディジタル値で表す。正規化方法はこの方法に限らず、様々な方法が考えられる。

そして、正規化手段８で正規化された場合の、第ｋ特徴量生成手段１２ｋの具体的構成例を図５に説明する。第ｋ特徴量生成手段１２ｋは例えば、量子化部１９ｊで位置（ｉ、ｊ）のパワー値をＤ個のビンで量子化される。即ち正規化されたパワーの取りうる最大値「１」と最小値「０」との範囲がＤ個の範囲（ビン）に分割され、量子化部１９ｊでは、これに入力されたパワースペクトルｘ（ｉ、ｊ）が、属するビンｄ（ｄ＝１、．．．、Ｄ）を「１」とし、ほかのＤ−１個のビンｄ’（ｄ≠ｄ’、ｄ’＝１、．．．、Ｄ）は「０」とする。例えば、図６に示すように、周波数ｊ１の場合は、ビンｄ１の範囲にあるので、ビンｄ１については「１」、それ以外のビンについては「０」とする。

量子化部１９１〜１９Ｊの各第ｄビン（ｄ＝１、．．．Ｄ）の値ｘ^ｄ（ｉ、１）、．．．、ｘ^ｄ（ｉ、ｊ）、．．．、ｘ^ｄ（ｉ、Ｊ）が各時間ｉごとに第ｄ観測値行列記憶手段１４ｄにそれぞれ入力され、観測値行列Ｘ^ｄ［ｘ^ｄ（ｉ、１）、．．．、ｘ^ｄ（ｉ、ｊ）、．．．、ｘ^ｄ（ｉ、Ｊ）］として出力される。

スペクトログラムの各位置の観測信号値の特徴を表す観測値行列としては例えば、図４Ｂに示す波形において、ピークとなっている周波数成分が存在する位置（ｉ、ｊ）のｘ（ｉ、ｊ）の特徴として、ｘ^ｋ（ｉ、ｊ）を「１」とし、ピークになっていないｘ（ｉ、ｊ）の特徴ｘ^ｋ（ｉ、ｊ）は「０」として、観測値行列Ｘ^ｋを２値で表す方法や、観測値行列Ｘ^ｋをスペクトログラムの各位置（ｉ、ｊ）におけるｘ（ｉ、ｊ）のパワー値を連続値として表す方法などが考えられる。これら特徴の数Ｋは複数が好ましい。

観測値行列生成手段１０よりの観測値行列集合Ｘ＝［Ｘ^１、．．．、Ｘ^Ｋ］は周波数成分修復・抽出手段１８に入力される。
周波数成分修復・抽出手段１８の具体的構成例を図７に示す。周波数成分修復・抽出手段１８は特徴ベクトル生成手段２２、特徴ベクトル記憶手段２４、事後確率計算手段３１、同時確率増大化手段３２、繰り返し処理判定手段３３、周波数成分出力手段３４とで構成されている。

観測値行列生成手段１０よりの観測値行列集合Ｘは特徴ベクトル生成手段２２に入力されると、スペクトログラムにおける全ての位置（ｉ、ｊ）について、位置（ｉ、ｊ）の観測値とその近傍の位置（ｉ’’、ｊ’’）の観測値および状態からなる観測特徴ベクトルｆｘを求める。ここで、図８Ａに示すようにスペクトログラムの周波数成分の数がＪ、単位時間数がＩとし、位置（ｉ、ｊ）を着目位置と呼び、その近傍の位置（ｉ’’、ｊ’’）とは着目位置（ｉ、ｊ）の周波数方向前後の各ｍ個の位置、時間方向前後の各ｎ個の位置で定義される矩形状の範囲に含まれる位置を表す。また、観測特徴ベクトルは各ｘ_１、．．．、ｘ^ｋのそれぞれの（ｉ、ｊ）の要素から構成し、［ｘ^１（ｉ、ｊ）、．．．、ｘ^ｋ（ｉ、ｊ）、．．．、ｘ^Ｋ（ｉ、ｊ）］となる。

図８Ａにスペクトログラムの観測値行列を示し、その位置（ｉ、ｊ）を中心にして、太線で囲っている矩形範囲（領域）が前記近傍を含む部分である。図８Ｂにその矩形範囲の拡大の詳細図を示す。横辺の位置数が２ｎ＋１、縦辺の位置数が２ｍ＋１である。この矩形内の位置の総数Ｎ_１は以下の式（１）で表すことが出来る。
Ｎ_１＝（２ｍ＋１）・（２ｎ＋１）・・・（１）
着目位置（ｉ、ｊ）に対する矩形範囲に含まれる各位置の特徴ベクトルを含む多次元観測特徴ベクトルｆｘを生成する。この多次元観測特徴ベクトルｆｘの次元数は位置の総数Ｎ_１と観測特徴ベクトルの次元数（特徴数）Ｋを乗算したＫ・Ｎ_１個の次元となる。
スペクトログラムの各位置（ｉ、ｊ）における周波数成分の有無を表す行列を周波数状態行列Θ^Ｆとするこれを図８Ｃに示す。

着目位置（ｉ、ｊ）に対する前記矩形範囲において、その位置（ｉ、ｊ）を除いた位置（ｉ’’、ｊ’’）に周波数成分の有無を表すベクトルを周波数成分特徴ベクトルｆｈとする。この周波数成分特徴ベクトルｆｈは図８Ｄに示すように、×印で示す着目位置（ｉ、ｊ）を除いた近傍の位置（ｉ’’、ｊ’’）の状態を要素とするＮ_１−１個の次元数Ｎ_２を持つ。なお、スペクトログラムの全位置（ｉ、ｊ）における周波数成分の有無を表す、つまり周波数成分状態の存在の確からしさを表す行列を、周波数状態行列Θ^Ｆといい、周波数状態行列Θ^Ｆを生成する（ステップＳ３０４）。この周波数状態行列Θ^Ｆは修復したい信号の周波数成分についての状態を表す行列である。周波数状態行列Θ^Ｆは例えば、各時間、各周波数における原信号の有無を表す２値、つまり原信号があれば「１」を示し、原信号がなければ「０」を示すことや、周波数成分の確率（０％〜１００％）を離散値で表したりすることが考えられる。また、着目位置（ｉ、ｊ）を中心にする前記矩形内の全近傍位置（ｉ’’、ｊ’’）に対する周波数成分状態θ^Ｆ（ｉ、ｊ）の全要素を表す周波数特徴ベクトルｆｈはＮ_２個である。

これにより、着目位置（ｉ、ｊ）における特徴ベクトルｆ（ｉ、ｊ）として以下の式（３）に示すものを生成する（ステップＳ３０６）。
ｆ（ｉ、ｊ）＝（ｆｘ、ｆｈ）・・・（３）
またｆ（ｉ、ｊ）の要素数（次元数）Ｎ_ａは、
Ｎ_ａ＝Ｎ_１×ｄ＋Ｎ_２・・・（４）
また、この実施例では更に、隠れ状態がある場合に、この発明を適用した例である。隠れ状態とは、各位置（ｉ、ｊ）にノイズがあるか否かを表すノイズ状態である。位置（ｉ、ｊ）に、ノイズがあれば、「１」とし、ノイズがなければ「０」として表示する。あるいは、ノイズのレベルを離散値で表したりすることが考えられる。

着目位置（ｉ、ｊ）を中心とする前記矩形範囲内の全ての近傍位置（ｉ’’、ｊ’’）における隠れ状態θ^Ｆ（ｉ、ｊ）の各要素よりなるベクトルを隠れ特徴ベクトルｆｎとする。この隠れ特徴ベクトルｆｎの要素の個数は特徴ベクトルｆｈと同様、Ｎ_２個である。なお、実施例１の以下の説明については一般性を持たせるため、隠れ状態を考慮した場合を説明する。

これにより、着目位置（ｉ、ｊ）における特徴ベクトルｆ（ｉ、ｊ）として、以下の式（５）に示すものを生成する。
ｆ（ｉ、ｊ）＝（ｆｘ、ｆｈ、ｆｎ）・・・（５）
よって隠れ状態を考慮した場合のｆ（ｉ、ｊ）の要素数Ｎ_ｂは、以下の式で表すことが出来る。
Ｎ_ｂ＝Ｎ_１×ｄ＋２×Ｎ_２・・・（６）
上述の処理で求められた特徴ベクトルｆ（ｉ、ｊ）は特徴ベクトル記憶手段２４に一旦記憶される。

周波数成分修復・抽出手段１８は、モデルパラメータ蓄積手段２８に予め格納された事後確率のモデルＰ（Θ^Ｆ、Θ^Ｎ│Ｘ）のためのパラメータと、特徴ベクトル記憶手段２４よりの特徴ベクトルｆ（ｉ、ｊ）を用いて、位置（ｉ、ｊ）における周波数成分状態θ^Ｆ（ｉ、ｊ）の推定値であるθ^〜Ｆ（ｉ、ｊ）と隠れ状態θ^Ｎ（ｉ、ｊ）の推定値であるθ^〜Ｎ（ｉ、ｊ）を計算する。（ステップＳ３０８）なお、事後確率のモデルＰ（Θ^Ｆ、Θ^Ｎ│Ｘ）のためのパラメータの詳細は後述する。ここで、Ｐ（状態Ａ│データＢ）はデータＢが与えられた時の、状態Ａになる確率を表す。

事後確率計算手段３１と同時確率増大化手段３２と繰り返し処理判定手段３３により以下の一連の処理を行う。例えばギブス・サンプラーなどの手法を用いて、同時確率Ｐ（Θ^Ｆ、Θ^Ｎ│Ｘ）を増大させる。推定値θ^〜Ｆ（ｉ、ｊ）、θ^〜Ｎ（ｉ、ｊ）を繰り返し求め、最終的に推定値をθ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）とする。なお、同時確率とは、周波数成分状態θ^Ｆ（ｉ、ｊ）、隠れ状態θ^Ｎ（ｉ、ｊ）の各要素の全てが同時に起こる確率を示し、隠れ状態θ^Ｎ（ｉ、ｊ）を用いない場合は、周波数状態θ^Ｆ（ｉ、ｊ）の各要素が同時に出現する確率である。

前記事後確率のモデルはマルコフ確率場によって表現される。マルコフ確率場は、格子などのグラフにおいて、ある格子点（ノード）における状態はその近傍のみによって決まるという考え方に基づいて作られる確率モデルである。なお、マルコフ確率場および上述のギブス・サンプラーの詳細は例えば、非特許文献１などに記載されている。この実施例の場合は、以下の式（７）のように定式化される。
Ｐ（θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）│θ^Ｆ（ｉ’、ｊ’）、θ^Ｎ（ｉ’、ｊ’）、Ｘ；ｉ≠ｉ’、ｊ≠ｊ’）＝
Ｐ（θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）│θ^Ｆ（ｉ’’、ｊ’’）、θ^Ｎ（ｉ’’、ｊ’’）、Ｘ；（ｉ’’、ｊ’’）∈Ｇ）１≦ｉ≦Ｉ、１≦ｊ≦Ｊ・・・（７）
式（７）の右辺の条件「θ^Ｆ（ｉ’’、ｊ’’）、θ^Ｎ（ｉ’’、ｊ’’）、Ｘ；（ｉ’’、ｊ’’）∈Ｇ」として、特徴抽出手段４で抽出された特徴ベクトルｆ（ｉ、ｊ）を用いる。つまり、位置（ｉ、ｊ）における周波数成分状態および隠れ状態の取る確率が以下の式（８）で表される。
Ｐ（θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）│ｆ（ｉ、ｊ））・・・式（８）

また、ギブス・サンプラーは、θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）を同時確率Ｐ（Θ^Ｆ、Θ^Ｎ│Ｘ）に従って、生成する方法である。全ての位置（ｉ、ｊ）に対して、同時に生成するのではなく、各位置（ｉ、ｊ）において式（８）の条件付確率に従って、θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）を繰り返し生成することで、同時確率を増大させる周波数成分状態行列Θ^Ｆと隠れ状態行列Θ^Ｎを求める方法である。

具体的処理の流れを、図９に示し、同時確率増大化手段３２と事後確率計算手段３１と、の具体的構成例とその他、関係のある手段を図１０に示す。
図１０に示すように、同時確率増大化手段３２は温度パラメータ生成手段４２、候補状態生成手段４７、遷移確率生成手段４８、一様乱数生成手段５０、比較手段５２、更新手段５４と、で構成されている。また事後確率計算手段３１は遷移前確率計算手段４４と遷移後確率計算手段４６とで構成されている。
まず、周波数成分状態行列生成手段４０で、例えば、Θ^Ｆ、Θ^Ｎの各要素にランダムな値を与えるなどして、初期推定値θ^Ｆ（ｉ、ｊ）（０）、θ^Ｎ（ｉ、ｊ）（０）を生成する。（ステップＳ２００）。

繰り返し処理判定手段３３で、繰り返し処理判定パラメータの１つである繰り返し回数ｔを「１」に初期化する（ステップＳ２０２）。温度パラメータ生成手段４２で温度降下関数Ｔ（ｔ）を以下の式（９）で計算する（ステップＳ２０４）。
Ｔ（ｔ）＝Ｃ／ｌｏｇ（１＋ｔ）・・・（９）
ただし、Ｔは温度パラメータであり、Ｃを定数とする。温度降下関数Ｔ（ｔ）は遷移確率生成手段４８に入力される。

次に、繰り返し処理判定手段３３において、繰り返し処理判定パラメータの時間ｉの初期値を「ｎ＋１」と設定し（ステップＳ２０６）、続いて、繰り返し処理判定パラメータの周波数ｊの初期値を「ｍ＋１」と設定する（ステップＳ２０８）。特徴ベクトル生成手段２２で特徴ベクトルｆ（ｉ、ｊ）を生成し（ステップＳ２１０）、遷移前確率計算手段４４と遷移後確率計算手段４６に入力される。特徴ベクトルｆ（ｉ、ｊ）に対する現在の周波数状態θ^Ｆ（ｉ、ｊ）（ｔ−１）および隠れ状態θ^Ｎ（ｉ、ｊ）（ｔ−１）の遷移前の出現確率（式（１０））をモデルパラメータ蓄積手段２８中のモデルパラメータを使用して、遷移前確率計算手段４４で求める（ステップＳ２１２、ステップＳ３０８）。
Ｑ＝Ｐ（θ^Ｆ（ｉ、ｊ）（ｔ−１）、θ^Ｎ（ｉ、ｊ）（ｔ−１）│ｆ（ｉ、ｊ））
・・・（１０）
ここで、［θ^Ｆ（ｉ、ｊ）（ｔ−１）、θ^Ｎ（ｉ、ｊ）（ｔ−１）］はそれぞれ周波数成分状態と隠れ状態の現在の状態を表す。この遷移前の出現確率Ｑは遷移確率生成手段４８に入力される。

そして、次の状態θ^Ｆ（ｉ、ｊ）（ｔ）、θ^Ｎ（ｉ、ｊ）（ｔ）の候補である新しい状態［φ^Ｆ（ｉ、ｊ）、φ^Ｎ（ｉ、ｊ）］を候補状態生成部４７でランダムに生成する（ステップＳ２１４）。
そして、特徴ベクトルｆ（ｉ、ｊ）に対するこれら新候補の周波数状態行列および隠れ状態行列の遷移後の出現確率Ｒを、モデルパラメータを用いて遷移後確率計算手段４６で、以下の式（１１）で求める（ステップＳ３０８）。
Ｒ＝Ｐ（φ^Ｆ（ｉ、ｊ）、φ^Ｎ（ｉ、ｊ）│ｆ（ｉ、ｊ））・・・（１１）
遷移後確率Ｒは遷移確率生成手段４８に入力される。

上述の遷移前の確率Ｑ、遷移後の確率Ｒ、温度降下関数Ｔ（ｔ）、とにより遷移確率αを遷移確率生成手段４８で、以下の式（１２）で求める（ステップＳ２１８）。
α＝ｍｉｎ［１、（Ｒ／Ｑ）^{１／Ｔ（ｔ）}］・・・（１２）
ここで、ｍｉｎ（ａ、ｂ）は、ａとｂのうち、小さい方の値を示す。遷移確率αは比較手段５２に入力される。

一方、一様乱数生成手段５０で一様乱数ｒを生成する（ステップＳ２２０）。ここで、ｒは０≦ｒ＜１を満たす実数である。遷移確率αと一様乱数ｒとを比較手段５２で比較して（ステップＳ２２２）、一様乱数ｒより遷移確率αの方が大きい場合は、比較手段５２が更新指示信号を生成して、更新指示信号が更新手段５４に入力され、以下の式（１３）（１４）に示すように、更新手段５４で次の状態へ遷移させる（ステップＳ２２４、Ｓ３１０）。
θ^Ｆ（ｉ、ｊ）（ｔ）＝φ^Ｆ（ｉ、ｊ）・・・（１３）
θ^Ｎ（ｉ、ｊ）（ｔ）＝φ^Ｎ（ｉ、ｊ）・・・（１４）

また、一様乱数ｒより遷移確率αの方が大きくない場合は、更新指示信号が生成されず、以下の（１５）（１６）に示すように、更新手段５４で次の状態へ遷移しない（Ｓ２２６）。
θ^Ｆ（ｉ、ｊ）（ｔ）＝θ^Ｆ（ｉ、ｊ）（ｔ−１）・・・（１５）
θ^Ｎ（ｉ、ｊ）（ｔ）＝θ^Ｎ（ｉ、ｊ）（ｔ−１）・・・（１６）
ステップＳ２２４、Ｓ２２６の後、ｊ＝Ｊ―ｍ＋１でなければ（ステップＳ２２８）、ｊを「１」インクリメントして、ステップＳ２１０に戻り、ステップＳ２２８でｊ＝Ｊ―ｍ＋１であり、ｉ＝Ｉ−ｎ＋１でなければ、ｉを「１」インクリメントして、ステップＳ２０８に戻り（ステップＳ２３４）、このようにして、ステップＳ２１０、Ｓ２１２、Ｓ２１４、Ｓ２１６、Ｓ２１８、Ｓ２２０、Ｓ２２２、Ｓ２２４、Ｓ２２６の処理を全ての時間（ｉ＝ｎ＋１、．．．、Ｉ−ｎ＋１）、全ての周波数（ｊ＝ｍ＋１、．．．Ｊ−ｍ＋１）で行うように、繰り返し処理判定手段３３が制御する（ステップＳ３１２）。

次に、状態の変化がある程度、収束した場合は、すなわち、ｉとｊの１巡中にステップＳ２２４で状態が更新される数が予め決められたしきい値以下になった場合、あるいはこの更新される状態数がほぼ一定値になった場合、繰り返し処理判定手段３３が、更新停止信号を生成し、更新停止信号が更新手段５４に入力され、更新処理が停止される（ステップＳ２３６）。また、繰り返し回数ｔが予め決められた回数Ｍａｘｔを超えても、繰り返し処理判定手段３３が更新手段５４の更新処理を停止させる。なお、しきい値は例えば、スペクトログラムの時間方向の位置の個数Ｉと周波数方向の位置の個数Ｊを用いて、およそＩ・Ｊ・１０００の値であることが好ましい。

このようにして更新が停止されると、同時確率が増大された周波数状態行列Θ^Ｆ’が周波数成分出力手段３４から出力される（ステップＳ３１４）。
このようにして、消滅した周波数成分が修復され、またノイズに埋もれた信号の周波数成分が取り出され、かつノイズ成分を除した周波数状態行列Θ^Ｆが得られる。
また、必要に応じて、この同時確率が増大され、推定された周波数状態行列Θ^Ｆ’（以下、周波数状態行列Θ^Ｆ’という）入力音響信号のスペクトログラムを用いて、抽出された周波数成分のみから構成された音響信号を合成し、出力する構成も考えられる。その場合、図７中の破線で示す周波数記憶手段３６と、音響信号合成手段３８が周波数成分修復手段１８に付加される。

周波数成分出力手段３４よりの周波数状態行列Θ^Ｆ’が周波数記憶手段３６に一旦記憶される。周波数状態行列Θ^Ｆ’は音響信号合成手段３８で例えば、図１中のスペクトログラム生成手段６で得られた入力音響信号のスペクトログラムと周波数状態行列Θ^Ｆ’の各成分同士の積をとることにより、存在する周波数成分の振幅を計算し、各時間フレーム（各時間ｉ）において存在する周波数成分の正弦波を生成し、和を取る方法などで実現する。出力される周波数成分信号の具体的形態図を図２Ｃに、ノイズ信号の具体的形態図を図２Ｄに示す。図２Ｃにおいては、太線で囲まれた部分がノイズ信号が取り除かれ、破線部分は修復された信号成分を示す。なお、縦軸が周波数方向、横軸は時間方向を示す。

実施例２
上述の実施例１では、予め、人間の手作業での試行錯誤によりマルコフ確率場における式（１１）に示す条件付確率を推定してモデルパラメータをモデルパラメータ蓄積手段２８に格納した。
実施例２では、学習モデル生成手段を設けて、式（１１）の条件付確率を計算するためのパラメータを学習モデル生成手段が予め学習して生成する。

図１１は実施例１のハードウェア構成に、確率モデル生成手段５９を付加させた実施例２のハードウェア構成を示す。確率モデル生成手段５９は学習音響信号入力部６０と学習周波数成分状態入力部６２と学習特徴ベクトル生成手段６３と確率モデル学習手段６４とモデルパラメータ蓄積手段２８とで構成されている。なお、その他の部分は実施例１と同様であり、それらは図１１では省略している。

まず、学習音響信号入力部６０から、学習音響信号が入力されると、特徴抽出手段４に入力されて、上述と同様の処理で、学習音響信号から学習スペクトログラムが生成され、その学習スペクトログラムから学習観測値行列が生成され、この学習観測値行列が学習特徴ベクトル生成手段６３に入力される。学習特徴ベクトル生成手段６３で、学習スペクトログラムの各位置（ｉ、ｊ）の特徴ベクトルｆ（ｉ、ｊ）を学習特徴ベクトルｆ（ｉ、ｊ）として生成する。この学習特徴ベクトルｆ（ｉ、ｊ）は確率モデル学習手段６４に入力される。

また、学習音響信号の学習用周波数状態行列Θ^Ｆが学習周波数成分状態入力部６２に入力され、これより各位置（ｉ、ｊ）に対応する学習周波数成分状態θ^Ｆ（ｉ、ｊ）が生成出力され、学習周波数成分状態θ^Ｆ（ｉ、ｊ）は確率モデル学習手段６４に入力される。確率モデル学習手段６４では、各学習特徴ベクトルｆ（ｉ、ｊ）とこれに対応する学習周波数成分状態θ^Ｆ（ｉ、ｊ）の組から、確率モデルＰ（θ^Ｆ（ｉ、ｊ）│ｆ（ｉ、ｊ））を学習する。

また、隠れ状態を考慮した場合は、原学習音響信号とノイズとを独立に生成し、原学習音響信号の一部をノイズで置き換えて、学習音響信号とし、これを学習音響信号入力部６０に入力し、原学習音響信号の周波数成分状態行列Θ^Ｆとノイズの隠れ状態行列Θ^Ｎとを学習周波数成分状態入力部６２に入力する。学習周波数成分状態入力部６２から、各位置（ｉ、ｊ）に対応する学習周波数成分状態θ^Ｆ（ｉ、ｊ）と隠れ状態θ^Ｎ（ｉ、ｊ）が生成出力され、確率モデル学習手段６４で、各学習特徴ベクトルｆ（ｉ、ｊ）とこれに対応する学習周波数成分状態θ^Ｆ（ｉ、ｊ）の組から、確率モデルＰ（θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）│ｆ（ｉ、ｊ））を学習する。

学習された確率モデルＰのパラメータをモデルパラメータ蓄積手段２８に蓄積し、周波数成分修復・抽出手段１８では、このモデルパラメータ蓄積手段２８中のモデルパラメータを使用して、特徴ベクトルｆ（ｉ、ｊ）と周波数成分状態θ^Ｆ（ｉ、ｊ）から事後確率を計算すればよい。

学習特徴ベクトルｆ（ｉ、ｊ）と学習周波数成分状態θ^Ｆ（ｉ、ｊ）の組み合わせを用いて、確率モデルＰを学習する手法として、パルツェンウィンドウ法がある。パルツェンウィンドウ法は確率密度関数を推定する方法の一つであり、特徴ベクトルが多値をとる場合に用いられる手法である。パルツェンウィンドウ法において、ベクトルＹが与えられたときの確率密度関数Ｐ（Ｙ）は、例えば、以下の式（１７）で表される。
Ｐ（Ｙ）＝Ｃ／ＭΣ^Ｍ _ｉ＝１ｅｘｐ［−（Ｙ−Ｙ_ｉ）・（Ｙ−Ｙ_ｉ）／２σ^２］・・・（１７）
と表現される。ここで、Ｃ、σは任意の定数となり、Ｍは学習データの個数、Ｙ_ｉは学習特徴ベクトルであり、Ａ・ＢはベクトルＡ、Ｂの内積を表す。

隠れ状態θ^Ｎ（ｉ、ｊ）を考慮しない場合は、学習特徴ベクトルｆ（ｉ、ｊ）と学習周波数成分状態θ^Ｆ（ｉ、ｊ）の全ての組み合わせについて、Ｐ（Ｙ）を計算し、全ての組み合わせにおけるＰ（Ｙ）の和が「１」となるように正規化することで、確率モデルＰ（θ^Ｆ（ｉ、ｊ）│ｆ（ｉ、ｊ））を学習することが可能である。

また、隠れ状態θ^Ｎ（ｉ、ｊ）を考慮する場合は、学習特徴ベクトルｆ（ｉ、ｊ）と学習周波数成分状態θ^Ｆ（ｉ、ｊ）と隠れ状態θ^Ｎ（ｉ、ｊの全ての組み合わせについて、Ｐ（Ｙ）を計算し、全ての組み合わせにおけるＰ（Ｙ）の和が「１」となるように正規化することで、確率モデルＰ（θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）│ｆ（ｉ、ｊ））を学習することが可能である。
なお、パルツェンウィンドウ法の詳細は「Ｄｏｎａｌｄ Ｆ．Ｓｐｅｃｈｔ．Ｇｅｎｅｒａｔｉｏｎ ｏｆ ｐｏｌｙｎｏｍｉａｌ ｄｉｓｃｒｉｍｉｎａｎｔ ｆｕｎｃｔｉｏｎｓ ｆｏｒ ｐａｔｔｅｒｎ ｒｅｃｏｇｎｉｔｉｏｎ ＩＥＥＥ Ｔｒａｎｓａｃｔｉｏｎｓ ｏｎ Ｅｌｅｃｔｒｏｎｉｃ Ｃｏｍｐｕｔｅｒｓ ＥＣ−１６（３）：３０８−３１９，１９６７」等に記載されている。

また確率モデルＰを学習する手法として、最大エントロピー法がある。最大エントロピー法では、学習データ中にいくつかの特性に着目して、モデル化を行うが、このために素性関数（ｆｅａｔｕｒｅ ｆｕｎｃｔｉｏｎ）と呼ばれる２値の関数を用い、特徴ベクトルｆ（ｉ、ｊ）が「１」もしくは「０」をとる場合に用いられる手法である。

最大エントロピー法によると、あるデータｄを与えた時に、そのデータｄが番号ｃで表現されるカテゴリに属する確率Ｐ（ｃ│ｄ）を以下の式（１８）で表現する。
Ｐ（ｃ│ｄ）＝１／Ｚ（ｄ）ｅｘｐ［Σ^Ｎｆｅ _ｉ＝１λ_ｉＦ_ｉ（ｄ、ｃ）］・・・（１８）
Ｆ_ｉ（ｄ、ｃ）はデータｄとカテゴリｃの組み合わせから計算されるＮｆｅ個の関数の一つであり、特徴関数と呼ばれ、λ_ｉは学習によって決まるパラメータである。Ｚ（ｄ）は正規化のための関数で以下の式（１９）で定義される。
Ｚ（ｄ）＝Σ^Ｎｃａｔ _ｃ＝１ｅｘｐ［Σ^Ｎｆｅ _ｉ＝１λ_ｉＦ_ｉ（ｄ、ｃ）］・・・（１９）
ここで、Ｎｃａｔはカテゴリ数とする。

以下、一般性を持たせるために、隠れ状態がある場合について説明する。周波数成分状態θ^Ｆ（ｉ、ｊ）、隠れ状態θ^Ｎ（ｉ、ｊ）はともに離散値をとるので、値の全ての組み合わせに対して番号を振り、これをｃとする。例えば、θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）がともに２値の場合はカテゴリ数Ｎｃａｔは４となる。

次にＦ_ｉ（ｄ、ｃ）を以下のように定義する。まず、データｄをスペクトログラム上の位置（ｉ、ｊ）とみなす。ｆ（ｉ、ｊ）の要素をｆ_ｅ（ｉ、ｊ）（ｅ＝１、２、．．．、Ｎｂ）とする。任意のｄ、ｃに対し、サイズＮｂのブロックＮｃａｔ個からなるＮｃａｔ×Ｎｂ個のＦ_１（ｄ、ｃ）を次のように定義する。つまり、Ｆ_１（ｄ、ｃ）におけるｃ番目のブロックのＮｂ個の要素それぞれをｆ_ｅ（ｉ、ｊ）で表し、それ以外の値を「０」とする。なお、最大エントロピー法の詳細は「Ａ ｍａｘｉｍｕｍ ｅｎｔｒｏｐｙ ａｐｐｒｏａｃｈ ｔｏ ａｄａｐｔｉｖｅ ｓｔａｔｉｓｔｉｃａｌ ｌａｎｇｕａｇｅ ｍｏｄｅｌｌｉｎｇ．Ｃｏｍｐｕｔｅｒ Ｓｐｅｅｃｈ ａｎｄ Ｌａｎｇｕａｇｅ」などに記載されている。

また、サポートベクターマシン（ＳＶＭ：Ｓｕｐｐｏｒｔ Ｖｅｃｔｏｒ ｍａｃｈｉｎｅ）等の識別関数の演算結果を例えばシグモイド関数（１／１＋ｅｘｐ（−ｘ））に入力して得られた値、つまり０〜１の実数を事後確率とみなして、事後確率のモデルを学習する手法も利用可能である。
サポートベクターマシンなどのパターン識別手法においては、識別関数ｇ（Ｘ）は以下の式（２０）のように表現される。
ｇ（Ｘ）＝α_０＋Σ^Ｍ _ｉ＝１α_ｉｘ_ｉ・・・（２０）
ここで、ＭはベクトルＸの次元、ｘ_ｉはベクトルＸの要素、α_１は係数、α_０は定数である。識別問題では、ｇ（Ｘ）の符号で、Ｘがあるカテゴリに含まれる否かを判断する。ここでは、ｇ（Ｘ）は一般に全ての実数をとるので、ｇ（Ｘ）をシグモイド関数に入力して得られた出力
Ｐ（Ｘ）＝１／１＋ｅｘｐ（−Ｃｇ（Ｘ））を擬似的に確率であるとみなす。なお、Ｃは定数である。
ここで、Ｘを特徴ベクトルｆ（ｉ、ｊ）、Ｘ_ｐを学習サンプルの特徴ベクトルとし、θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）のすべての組み合わせについて、Ｐ（Ｘ）を学習し、全ての組み合わせにおけるＰ（Ｘ）の和が「１」となるように、正規化することで、確率モデルＰ（θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）│ｆ（ｉ、ｊ））を得ることが可能となる。

実施例３
実施例３は先の実施例とは特徴ベクトルｆ（ｉ、ｊ）の生成が異なり、これに伴って、処理計算が異なる。この異なる部分について主として説明する。
スペクトログラム上の位置（ｉ、ｊ）において観測特徴ベクトルｆｘ（ｉ、ｊ）を特徴ベクトルｆ（ｉ、ｊ）とする。即ち、実施例３における特徴ベクトルｆ（ｉ、ｊ）は位置（ｉ、ｊ）とその近傍位置（ｉ’’、ｊ’’）の各観測値行列Ｘ^k（ｋ＝１、…、Ｋ）の各成分よりなるベクトルｆｘのみであり、周波数成分状態行列Θ^Fの位置（ｉ、ｊ）を中心とする近傍の位置の（ｉ’’、ｊ’’）の各成分のベクトルｆｈは特徴ベクトルｆ（ｉ、ｊ）に含まれない。

事後確率計算では特徴ベクトルｆ（ｉ、ｊ）（先の実施例と区別して書くとｆｘ(ｉ、ｊ）に対し、位置（ｉ、ｊ）における成分ベクトルがθ^F（ｉ、ｊ）である事後確率Ｐ２（θ^F（ｉ、ｊ）｜ｆ（ｉ、ｊ））、つまりＰ２（θ^F（ｉ、ｊ）｜ｆｘ（ｉ、ｊ））を求める。
ここで、（θ^F（ｉ，ｊ）の取りうるすべての値について、確率Ｐ２（θ^F（ｉ，ｊ）｜ｆ_ｘ）を計算し、最も高い確率を与える（θ^F（ｉ，ｊ）を周波数成分状態行列の（ｉ，ｊ）成分に設定する。

この場合、隠れ状態を導入する場合は特徴ベクトルｆ（ｉ、ｊ）として観測特徴ベクトルｆｘ（ｉ、ｊ）とするが、事後確率はＰ２（θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）｜ｆｘ（ｉ、ｊ））を計算すればよい。確率Ｐ２（θ^Ｆ（ｉ、ｊ）｜ｆｘ（ｉ、ｊ））やＰ２（θ^Ｆ（ｉ、ｊ）、θ^Ｎ（ｉ、ｊ）｜ｆｘ（ｉ、ｊ））は、例えば、上述の最大エントロピー法、パルツェン窓法、識別関数を利用する方法などを用いて、学習して事前に用意しておく。これらの手法の適用は、確率モデル学習手段６４におけるｆ（ｉ、ｊ）の代わりに、ｆｘ（ｉ、ｊ）を用いる方法で実現できる。
また、この方法により周波数成分状態行列の（ｉ、ｊ）成分を設定した場合には、同時確率増大化過程で用いられる温度降下関数Ｔ（ｔ）の定数Ｃは０が望ましい。

実験結果
以下にこの発明の実験結果を示す。
実験データとして、６種類の楽器音（ピアノ、バイオリン、フルート、トランペット、アルト（声）、マリンバ）および１種類の人工音を使って、それぞれ７、５秒間の音楽のフレーズを作成し、それぞれのスペクトログラムをとる。このスペクトログラムの周波数方向においてピークとなる位置を「１」、それ以外を「０」としたマトリクスを正解データ（周波数成分状態）とする。そして、この正解データそれぞれに１５０ｍｓｅｃ〜３００ｍｓｅｃの白色雑音を５箇所加えた音響データをテストデータとして使う。

サンプリング周波数は全てにおいて、１１、０２５Ｈｚであり、スペクトログラム作成には、窓幅１０２４のハニング窓をかけた高速フーリエ変換を用いて、時間軸上のフレームのオーバーラップは無しとし、単位はｄｂとした。
スペクトログラムから観測値行列Ｘを作成する方法は、各位置（ｉ、ｊ）においてスペクトログラムのパワー値を３個のビンに量子化し、パワー値に該当するビンのみ「１」とし、その他のビンには「０」を入れる方法に加え、位置（ｉ、ｊ）が周波数特性波形におけるピークになっていれば「１」、なっていなければ「０」を入れる方法を使った。

事後確率の学習方法として、最大エントロピー法を用いた。また近傍の位置の範囲は、周波数方向に±３位置（ｍ＝３）とし、時間方向に±３位置（ｎ＝３）とした。
上記の実験データで、周波数成分の修復および抽出精度を計測し、従来の技術の結果と、この発明の結果と、を比較した表を図１２に示す。図１２Ａ、Ｂともに周波数成分の抽出度、修復度、抽出度と修復度の調和平均を示したものである。括弧内の数値は、この発明の手法を適用せず、単純なスペクトログラムのピーク抽出による結果である。
図１２Ａは、７、５秒の音楽フレーズ全ての領域について、位置（ｉ、ｊ）単位で正解データと比べた結果であり、図１２Ｂは、図１２Ａと同様の評価をノイズをかけた領域（時間帯）のみについて行ったものである。単純ピーク抽出での精度を表している括弧内の数値と比較すると、この発明の方が大きな精度を得ることが出来ることがわかる。

以上の各実施形態の他、本発明である周波数修復・抽出方法は上述の実施形態に限定されるものではなく、本発明の趣旨を逸脱しない範囲で適宜変更が可能である。また、周波数修復・抽出方法において説明した処理は、記載の順に従って時系列に実行されるのみならず、処理を実行する装置の処理能力あるいは必要に応じて並列的にあるいは個別に実行されるとしてもよい。
また、この発明の周波数修復・抽出方法における処理をコンピュータによって実現する場合、周波数修復・抽出方法が有すべき機能の処理内容はプログラムによって記述される。そして、このプログラムをコンピュータで実行することにより、周波数修復・抽出方法における処理機能がコンピュータ上で実現される。

この処理内容を記述したプログラムは、コンピュータで読み取り可能な記録媒体に記録しておくことができる。コンピュータで読み取り可能な記録媒体としては、例えば、磁気記録装置、光ディスク、光磁気記録媒体、半導体メモリ等どのようなものでもよい。具体的には、例えば、磁気記録装置として、ハードディスク装置、フレキシブルディスク、磁気テープ等を、光ディスクとして、ＤＶＤ（Ｄｉｇｉｔａｌ Ｖｅｒｓａｔｉｌｅ Ｄｉｓｃ）、ＤＶＤ−ＲＡＭ（Ｒａｎｄｏｍ Ａｃｃｅｓｓ Ｍｅｍｏｒｙ）、ＣＤ−ＲＯＭ（Ｃｏｍｐａｃｔ Ｄｉｓｃ Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）、ＣＤ−Ｒ（Ｒｅｃｏｒｄａｂｌｅ）／ＲＷ（ＲｅＷｒｉｔａｂｌｅ）等を、光磁気記録媒体として、ＭＯ（Ｍａｇｎｅｔｏ−Ｏｐｔｉｃａｌ ｄｉｓｃ）等を、半導体メモリとしてＥＥＰ−ＲＯＭ（Ｅｌｅｃｔｒｏｎｉｃａｌｌｙ Ｅｒａｓａｂｌｅ ａｎｄ Ｐｒｏｇｒａｍｍａｂｌｅ−Ｒｅａｄ Ｏｎｌｙ Ｍｅｍｏｒｙ）等を用いることができる。

また、このプログラムの流通は、例えば、そのプログラムを記録したＤＶＤ、ＣＤ−ＲＯＭ等の可搬型記録媒体を販売、譲渡、貸与等することによって行う。さらに、このプログラムをサーバコンピュータの記憶装置に格納しておき、ネットワークを介して、サーバコンピュータから他のコンピュータにそのプログラムを転送することにより、このプログラムを流通させる構成としてもよい。

このようなプログラムを実行するコンピュータは、例えば、まず、可搬型記録媒体に記録されたプログラムもしくはサーバコンピュータから転送されたプログラムを、一旦、自己の記憶装置に格納する。そして、処理の実行時、このコンピュータは、自己の記録媒体に格納されたプログラムを読み取り、読み取ったプログラムに従った処理を実行する。また、このプログラムの別の実行形態として、コンピュータが可搬型記録媒体から直接プログラムを読み取り、そのプログラムに従った処理を実行することとしてもよく、さらに、このコンピュータにサーバコンピュータからプログラムが転送されるたびに、逐次、受け取ったプログラムに従った処理を実行することとしてもよい。また、サーバコンピュータから、このコンピュータへのプログラムの転送は行わず、その実行指示と結果取得のみによって処理機能を実現する、いわゆるＡＳＰ（Ａｐｐｌｉｃａｔｉｏｎ Ｓｅｒｖｉｃｅ Ｐｒｏｖｉｄｅｒ）型のサービスによって、上述の処理を実行する構成としてもよい。なお、本形態におけるプログラムには、電子計算機による処理の用に供する情報であってプログラムに準ずるもの（コンピュータに対する直接の指令ではないがコンピュータの処理を規定する性質を有するデータ等）を含むものとする。

また、この形態では、コンピュータ上で所定のプログラムを実行させることにより、周波数修復・抽出方法を構成することとしたが、これらの処理内容の少なくとも一部をハードウェア的に実現することとしてもよい。

この発明装置の実施例中の特徴抽出手段４のハードウェアの構成例を示すブロック図。 入力された音響信号、この信号のスペクトログラム、出力される周波数成分信号、ノイズ信号の具体的形態例を示す図。 スペクトログラムの各位置（ｉ、ｊ）の位置関係の概念図。 周波数領域に変換された入力信号の波形図であり図４Ａは正規化前の波形図であり、図４Ｂは正規化後の波形図である。 第ｋ特徴量生成手段１２ｋの具体的構成例を示すブロック図。 正規化された音響信号において、Ｄ個のビンに分割されることを示す図。 周波数成分修復・抽出手段１８の具体的構成例を示すブロック図。 図８Ａはスペクトログラムとその位置（ｉ、ｊ）とその近傍領域を示す図であり、図８Ｂはベクトルｆ（ｉ、ｊ）が生成される領域を示す拡大図であり、図８Ｃは状態位置の配置例とその位置（ｉ、ｊ）の近傍領域を示す図であり、図８Ｄは状態行列の特徴ベクトルｆｈおよびｆｎが生成される領域を示す拡大図である。 同時確率増大化手段３２の処理の流れの例を示すフローチャート。 同時確率増大化手段３２の具体的構成例を示すブロック図。 確率モデル生成手段５９の具体的構成例を示すブロック図。 従来の技術とこの発明の周波数成分の修復・抽出を行った実験結果を示す図であり、図１１Ａはスペクトログラムの全ての領域について行った場合の値であり、図１１Ｂはノイズをかけた領域について行った場合の値を示すものである。 この発明の全実施例に共通する主たる処理の流れの例を示すフローチャート。

Claims (15)

1. 入力音響信号の周波数成分を抽出するとともに失われた成分を修復する周波数成分修復・抽出方法であって、
   入力音響信号に対し、周波数解析を行い、サウンドスペクトログラム（以下、スペクトログラムという）ｘ（ｉ、ｊ）（ｉは時間、ｊは周波数を示す）を算出するスペクトログラム生成過程と、
   上記スペクトログラム生成過程で生成されたスペクトログラムｘ（ｉ、ｊ）の各位置（ｉ、ｊ）の値から、Ｋ（Ｋ＞１）個の特徴量ｘ^ｋ（ｉ、ｊ）（ｋ＝１、．．．、Ｋ）を算出し、それぞれの特徴に関して、ｘ^ｋ（ｉ、ｊ）をｉ行ｊ列の要素とする観測値行列Ｘ^ｋ（ｋ＝１、．．．、Ｋ）を生成する観測値行列生成過程と、
   周波数成分状態の存在の確からしさを表す周波数成分状態行列Θ^Ｆを生成する周波数成分状態行列生成過程と、
   上記観測値行列生成手段で生成されたそれぞれの観測値行列Ｘ^ｋ内の位置（ｉ、ｊ）を中心とする近傍の値（位置（ｉ、ｊ）を含む）から、位置（ｉ、ｊ）に関する特徴ベクトルｆ（ｉ、ｊ）を生成する特徴ベクトル生成過程と、
   上記特徴ベクトル生成過程で生成された特徴ベクトルｆ（ｉ、ｊ）から、スペクトログラム上の位置（ｉ、ｊ）の任意の周波数成分状態が出現する確率を計算する事後確率計算過程と、
   位置（ｉ、ｊ）に対して上記特徴ベクトル生成過程と、上記事後確率計算過程と、による処理を行った後、上記事後確率を用いて、事後同時確率が大きくなるように、上記周波数成分状態行列を変更する同時確率増大化過程と、
   全ての位置（ｉ、ｊ）に対する上記同時確率増大化過程の処理後、所定の条件により、上記特徴ベクトル生成過程と、上記事後確率計算過程と、上記同時確率増大化過程と、を更に繰り返しを行うかどうかを判定する繰り返し処理判定過程と、
   上記繰り返し処理判定過程が繰り返しを行わないと判定すれば、上記周波数成分状態行列の値に基づき入力音響信号の周波数成分を出力する周波数成分出力過程と、を有することを特徴とする周波数成分修復・抽出方法。
2. 請求項１に記載の周波数修復・抽出方法において、
   上記特徴ベクトル生成過程は、上記観測値行列生成手段で生成されたそれぞれの観測値行列Ｘ^ｋ内の位置（ｉ、ｊ）を中心とする近傍の値（位置（ｉ、ｊ）を含む）と、上記周波数成分状態行列生成過程で生成された周波数成分状態行列Θ^Ｆ内の上記位置（ｉ、ｊ）を中心とする近傍の位置（位置（ｉ、ｊ）を除く）とを含めたものから、位置（ｉ、ｊ）に関する特徴ベクトルｆ（ｉ、ｊ）を生成する過程であることを特徴とする周波数成分修復・抽出方法。
3. 請求項１に記載の周波数修復・抽出方法において、
   上記特徴ベクトル生成過程は、上記スペクトログラム上の位置（ｉ、ｊ）およびその前後の近傍の各位置（ｉ’’、ｊ’’）における観測値行列の各要素値ｘ（ｉ、ｊ）により上記特徴ベクトルｆ（ｉ、ｊ）として生成する過程であり、
   上記事後確率計算過程は、上記与えられた特徴ベクトルｆ（ｉ、ｊ）から周波数成分状態行列Θ^Ｆの上記位置（ｉ、ｊ）の成分θ^Ｆ（ｉ、ｊ）が出現する事後確率モデルに上記特徴ベクトルｆ（ｉ、ｊ）を適用して上記確率を計算する過程であることを特徴とする周波数修復・抽出方法。
4. 請求項２に記載の周波数修復・抽出方法において、
   上記特徴ベクトル生成過程は、上記スペクトログラム上の位置（ｉ、ｊ）およびその前後の近傍の各位置（ｉ’’、ｊ’’）における観測値行列の各要素値ｘ（ｉ、ｊ）、と、上記位置（ｉ、ｊ）の前後の近傍の各位置（ｉ’’、ｊ’’）における周波数成分状態θ（ｉ’’、ｊ’’）^Ｆとにより上記特徴ベクトルｆ（ｉ、ｊ）として生成する過程であり、
   上記事後確率計算過程は、上記与えられた特徴ベクトルｆ（ｉ、ｊ）から周波数成分状態行列Θ^Ｆの上記位置（ｉ、ｊ）の成分θ^Ｆ（ｉ、ｊ）が出現する事後確率モデルに上記特徴ベクトルｆ（ｉ、ｊ）を適用して上記確率を計算する過程であることを特徴とする周波数修復・抽出方法。
5. 請求項２または４何れかに記載の周波数修復・抽出方法において、
   学習音響信号からそのスペクトログラムの各位置の特徴ベクトルｆ（ｉ、ｊ）を学習特徴ベクトルｆ（ｉ、ｊ）として生成する学習特徴ベクトル生成過程と、
   上記各学習特徴ベクトルｆ（ｉ、ｊ）とその各位置（ｉ、ｊ）と対応する学習周波数成分状態θ（ｉ、ｊ）^Ｆの組を用いて、事後確率のモデルを学習し、その学習されたモデルのパラメータをモデルパラメータ蓄積手段に蓄積する確率モデル学習過程とを更に有し、
   上記事後確率計算過程は、上記蓄積されたパラメータを使って、上記特徴ベクトルｆ（ｉ、ｊ）と周波数成分状態θ（ｉ、ｊ）^Ｆから事後確率を計算する過程であることを、特徴とする周波数修復・抽出方法。
6. 請求項５に記載の周波数修復・抽出方法において、
   上記確率モデル学習過程はパルツェン窓（Ｐａｒｚｅｎ Ｗｉｎｄｏｗ）法を用いて、事後確率のモデルを学習する過程であることを特徴とする周波数修復・抽出方法。
7. 請求項５に記載の周波数修復・抽出方法において、
   上記観測値行列生成過程は、上記観測値行列の各要素の値を全て「０」もしくは「１」（以下、２値という）として各観測値行列を生成する過程であり、
   上記確率モデル学習過程は、最大エントロピー法を用いて、事後確率のモデルを学習する過程であることを特徴とする周波数修復・抽出方法。
8. 請求項５に記載の周波数修復・抽出方法において、
   上記確率モデル学習過程は、上記学習特徴ベクトルにつき識別関数を演算し、その演算結果をシグモイド関数に入力して得られた値を事後確率とみなして、事後確率のモデルを学習する過程であることを特徴とする周波数修復・抽出方法。
9. 請求項２又は４に記載の周波数修復・抽出方法において、
   上記特徴ベクトル生成過程は、上記スペクトログラム上の各位置（ｉ、ｊ）において、１つ以上の隠れ状態θ（ｉ、ｊ）^Ｎで構成される行列（以下、隠れ状態行列Θ^Ｎという）の上記位置（ｉ、ｊ）の要素を除く成分をも含めて上記特徴ベクトルｆ（ｉ、ｊ）を生成する過程であり、
   上記事後確率計算過程は、上記与えられた特徴ベクトルｆ（ｉ、ｊ）から周波数成分状態θ（ｉ、ｊ）^Ｆおよび隠れ状態θ（ｉ、ｊ）^Ｎが出現する事後確率モデルに、上記特徴ベクトルｆ（ｉ、ｊ）を適用して確率を計算する過程であり、
   上記同時確率増大化過程において、上記初期状態生成過程はスペクトログラム上の各位置について、隠れ状態行列の初期値も発生させる過程であり、
   上記要素推定過程は周波数成分状態および隠れ状態を更新し、事後確率を増大させることで、周波数成分状態行列Θ^Ｆおよび隠れ状態行列Θ^Ｎの各要素の値を推定する過程であることを特徴とする周波数修復・抽出方法。
10. 請求項５〜８のいずれかに記載の周波数修復・抽出方法において、
    上記特徴ベクトル生成過程は、上記スペクトログラム上の各位置（ｉ、ｊ）において、１つ以上の隠れ状態θ（ｉ、ｊ）^Ｎで構成される行列（以下、隠れ状態行列Θ^Ｎという）の上記位置（ｉ、ｊ）の要素を除く成分をも含めて上記特徴ベクトルを生成する過程であり、
    上記事後確率計算過程は、上記与えられた特徴ベクトルｆ（ｉ、ｊ）から周波数成分状態θ（ｉ、ｊ）^Ｆおよび隠れ状態θ（ｉ、ｊ）^Ｎが出現する事後確率モデルに、上記特徴ベクトルを適用して確率を計算する過程であり、
    上記学習特徴ベクトル生成過程は、上記学習音響信号より、隠れ状態行列から上記位置（ｉ、ｊ）の要素を除く成分をも含む学習特徴ベクトルを生成する過程であり、
    上記確率モデル学習過程は、上記学習特徴ベクトルｆ（ｉ、ｊ）とその各位置と対応する学習周波数成分状態θ（ｉ、ｊ）^Ｆと隠れ状態θ（ｉ、ｊ）^Ｎとの組を用いて事後確率のモデルを学習し、その学習されたモデルのパラメータをモデルパラメータ蓄積手段に蓄積する過程であり、
    上記同時確率増大化過程において、上記初期状態生成過程はスペクトログラム上の各位置について、隠れ状態行列の初期値も発生させる過程であり、
    上記要素推定過程は周波数成分状態および隠れ状態を更新し、事後確率を増大させることで、周波数成分状態行列Θ^Ｆおよび隠れ状態行列Θ^Ｎの各要素の値を推定する過程であることを特徴とする周波数修復・抽出方法。
11. 請求項１〜１０のいずれかに記載の周波数修復・抽出方法において、
    正規化手段が、上記スペクトログラムにおける各位置（ｉ、ｊ）のパワースペクトルを局所的に正規化し、正規化された値を上記観測値行列生成過程で処理する実数値行列とする過程を、有することを特徴とする周波数修復・抽出方法。
12. 請求項１〜１１のいずれかに記載の周波数修復・抽出方法において、
    音響信号合成手段が、上記周波数成分出力過程により出力された周波数成分で構成される音響信号を合成する過程を、有することを特徴とする周波数修復・抽出方法。
13. 入力音響信号の周波数成分を抽出するとともに失われた成分を修復する周波数成分修復・抽出装置であって、
    入力音響信号に対し、周波数解析を行い、サウンドスペクトログラム（以下、スペクトログラムという）ｘ（ｉ、ｊ）（ｉは時間、ｊは周波数を示す）を算出するスペクトログラム生成部と、
    上記スペクトログラム生成過程で生成されたスペクトログラムｘ（ｉ、ｊ）の各位置（ｉ、ｊ）の値から、Ｋ（Ｋ＞１）個の特徴量ｘ^ｋ（ｉ、ｊ）（ｋ＝１、．．．、Ｋ）を算出し、それぞれの特徴に関して、ｘ^ｋ（ｉ、ｊ）をｉ行ｊ列の要素とする観測値行列Ｘ^ｋ（ｋ＝１、．．．、Ｋ）を生成する観測値行列生成過程と、
    周波数成分状態の存在の確からしさを表す周波数成分状態行列Θ^Ｆを生成する周波数成分状態行列生成部と、
    上記観測値行列生成手段で生成されたそれぞれの観測値行列Ｘ^ｋ内の位置（ｉ、ｊ）を中心とする近傍の値（位置（ｉ、ｊ）を含む）から、位置（ｉ、ｊ）に関する特徴ベクトルｆ（ｉ、ｊ）を生成する特徴ベクトル生成部と、
    上記特徴ベクトル生成過程で生成された特徴ベクトルｆ（ｉ、ｊ）から、スペクトログラム上の位置（ｉ、ｊ）の任意の周波数成分が出現する確率を計算する事後確率計算部と、
    位置（ｉ、ｊ）に対して上記特徴ベクトル生成過程と、上記事後確率計算過程と、による処理を行った後、上記事後確率を用いて、事後同時確率が大きくなるように、上記周波数成分状態行列を変更する同時確率増大化部と、
    全ての位置（ｉ、ｊ）に対する上記同時確率増大化過程の処理後、所定の条件により、上記特徴ベクトル生成過程と、上記事後確率計算過程と、上記同時確率増大化過程と、を更に繰り返しを行うかどうかを判定する繰り返し処理判定部と、
    上記繰り返し処理判定過程が繰り返しを行わないと判定すれば、上記周波数成分状態行列の値に基づき入力音響信号の周波数成分を出力する周波数成分出力部と、を有することを特徴とする周波数成分修復・抽出装置。
14. 請求項１〜１２の何れかに記載した周波数修復・抽出方法の各過程をコンピュータに実行させるための周波数修復・抽出プログラム。
15. 請求項１４記載の周波数修復・抽出プログラムを記録したコンピュータ読み取り可能な記録媒体。

Images