(57)【特許請求の範囲】
【請求項１】
・Tは行列又はベクトル・の転置を表すとして、入力変数xτ及び出力変数yτの学習データに基づいて、ガウス混合モデルで表現された結合分布p(ATxτ,yτ;Θ)のガウス混合モデルパラメータΘ及び圧縮行列Aを求めるガウス混合モデルパラメータ計算装置であって、
入力されたガウス混合モデルパラメータΘに基づいて、EMアルゴリズムにおける負担率を計算する負担率計算部と、
入力されたガウス混合モデルパラメータΘ及び上記負担率に基づいてガウス混合モデルパラメータΘを更新するガウス混合モデルパラメータ更新部と、
入力されたガウス混合モデルパラメータΘ、上記負担率及び入力された圧縮行列Aに基づいて、上記入力された圧縮行列Aを勾配法を用いた最適化により更新する圧縮行列更新部と、
上記更新されたガウス混合モデルパラメータΘ及び上記更新された圧縮行列Aを上記入力されたガウス混合モデルパラメータΘ及び上記入力された圧縮行列Aとして、上記負担率計算部、上記ガウス混合モデルパラメータ更新部及び上記圧縮行列更新部の処理を繰り返す制御を行う制御部と、
を含むガウス混合モデルパラメータ計算装置。
【請求項２】
請求項１のガウス混合モデルパラメータ計算装置と、
上記繰り返し制御により最終的に生成された圧縮行列Aを用いて、情報推定の対象となる入力変数xτを圧縮したATxτを計算する次元圧縮部と、
上記繰り返し制御により最終的に生成されたガウス混合モデルパラメータΘを用いて、上記計算されたATxτが得られた下でのyτの期待値である^yτを計算する推定部と、
を含む情報推定装置。
【請求項３】
請求項１のガウス混合モデルパラメータ計算装置であって、
上記入力変数xτは目的音が含まれる音信号の音響特徴量であり、上記出力変数yτは事前信号雑音比である、
ガウス混合モデルパラメータ計算装置。
【請求項４】
強調の対象となる目的音が含まれる音信号を周波数領域信号に変換する周波数領域変換部と、
上記周波数領域信号に基づいて、強調の対象となる目的音が含まれる音信号の音響特徴量xτを抽出する音響特徴量抽出部と、
請求項３のガウス混合モデルパラメータ計算装置により最終的に生成された圧縮行列Aを用いて、上記抽出された音響特徴量xτを圧縮したATxτを計算する次元圧縮部と、
請求項３のガウス混合モデルパラメータ計算装置により最終的に生成されたガウス混合モデルパラメータΘを用いて、上記計算されたATxτが得られた下での事前信号雑音比yτの期待値である^yτを計算する推定部と、
上記計算された^yτを用いてウィナーフィルタを求めるウィナーフィルタ計算部と、
上記求まったウィナーフィルタを上記周波数領域信号に適用してフィルタリング後周波数領域信号を得るフィルタリング部と、
上記フィルタリング後周波数領域信号を時間領域信号に変換する時間領域変換部と、
を含む音強調装置。
【請求項５】
・Tは行列又はベクトル・の転置を表すとして、入力変数xτ及び出力変数yτの学習データに基づいて、ガウス混合モデルで表現された結合分布p(ATxτ,yτ;Θ)のガウス混合モデルパラメータΘ及び圧縮行列Aを求めるガウス混合モデルパラメータ計算方法であって、
負担率計算部が、入力されたガウス混合モデルパラメータΘに基づいて、EMアルゴリズムにおける負担率を計算する負担率計算ステップと、
ガウス混合モデルパラメータ更新部が、入力されたガウス混合モデルパラメータΘ及び上記負担率に基づいてガウス混合モデルパラメータΘを更新するガウス混合モデルパラメータ更新ステップと、
圧縮行列更新部が、入力されたガウス混合モデルパラメータΘ、上記負担率及び入力された圧縮行列Aに基づいて、上記入力された圧縮行列Aを勾配法を用いた最適化により更新する圧縮行列更新ステップと、
制御部が、上記更新されたガウス混合モデルパラメータΘ及び上記更新された圧縮行列Aを上記入力されたガウス混合モデルパラメータΘ及び上記入力された圧縮行列Aとして、上記負担率計算部、上記ガウス混合モデルパラメータ更新部及び上記圧縮行列更新部の処理を繰り返す制御を行う制御ステップと、
を含むガウス混合モデルパラメータ計算方法。
【請求項６】
請求項５のガウス混合モデルパラメータ計算方法の各ステップと、
次元圧縮部が、上記繰り返し制御により最終的に生成された圧縮行列Aを用いて、情報推定の対象となる入力変数xτを圧縮したATxτを計算する次元圧縮ステップと、
推定部が、上記繰り返し制御により最終的に生成されたガウス混合モデルパラメータΘを用いて、上記計算されたATxτが得られた下でのyτの期待値である^yτを計算する推定ステップと、
を含む情報推定方法。
【請求項７】
周波数領域変換部が、強調の対象となる目的音が含まれる音信号を周波数領域信号に変換する周波数領域変換ステップと、
音響特徴量抽出部が、上記周波数領域信号に基づいて、強調の対象となる目的音が含まれる音信号の音響特徴量xτを抽出する音響特徴量抽出ステップと、
次元圧縮部が、請求項５のガウス混合モデルパラメータ計算方法により最終的に生成された圧縮行列Aを用いて、上記抽出された音響特徴量xτを圧縮したATxτを計算する次元圧縮ステップと、
推定部が、請求項５のガウス混合モデルパラメータ計算方法により最終的に生成されたガウス混合モデルパラメータΘを用いて、上記計算されたATxτが得られた下での事前信号雑音比yτの期待値である^yτを計算する推定ステップと、
ウィナーフィルタ計算部が、上記計算された^yτを用いてウィナーフィルタを求めるウィナーフィルタ計算ステップと、
フィルタリング部が、上記求まったウィナーフィルタを上記周波数領域信号に適用してフィルタリング後周波数領域信号を得るフィルタリングステップと、
時間領域変換部が、上記フィルタリング後周波数領域信号を時間領域信号に変換する時間領域変換ステップと、
を含む音強調方法。
【請求項８】
請求項１又は３のガウス混合モデルパラメータ計算装置、請求項２の情報推定装置及び請求項４の音強調装置の何れかの装置の各部としてコンピュータを機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
この発明は、ガウス混合モデルパラメータを計算するための技術、情報を推定するための技術又は特定の音を強調するための技術に関する。
【背景技術】
【０００２】
ある入力xτ∈RD,τ={1,2,...,T}から、それに対応する出力yτ∈R×Ωを推定する問題を考える。D,Ωは、所定の正の整数とする。音源強調を例に挙げて説明すると、この問題は、スポーツフィールドである特定の競技音（e.g. キック音やホイッスル音）を収録する場合において、M≧1本のマイクロホンで観測した信号からD次元の音響特徴量xを抽出し、目的音を抽出するフィルタ設計するパラメータyを推定する問題である。
【０００３】
音源強調の代表的な手法であるウィナーフィルタリングを行う場合、観測から推定するべきパラメータyτ=(y1,τ,y2,τ,...,yΩ,τ)Tは事前信号雑音比（以下、事前SNRとも表記する。）である。ここで、・Tは行列又はベクトル・の転置を表す。以下ではウィナーフィルタリングと事前SNR の用い方について簡単に説明する。目的音（e.g. 収音したい特定の競技音）の複素スペクトルをSω,τ、雑音（e.g. 歓声）の複素スペクトルをNω,τ、時間と周波数のインデックスをそれぞれω∈{1,2,...,Ω}、τ∈{1,2,...,T}としたとき、観測信号は
Fω,τ=Sω,τ+Nω,τ
と記述することができる。目的音を強調する時刻τのウィナーフィルタGω,τは、
【０００４】
【数１】
0006370751.tif 000002
【０００５】
と求めることができる。また、yω,τは近似計算された事前信号雑音比である。このGω,τと観測信号Fω,τとを掛け合わせることにより、目的音を強調した信号Hω,τを得ることができる。
【０００６】
Hω,τ=Gω,τFω,τ
つまり、適切なウィナーフィルタを設計するためには、観測信号Fω,τから事前SNRを精度よく推定する必要がある。
【０００７】
推定の問題を解く一般的な手法として、最小二乗平均誤差推定（以下、MMSEとも表記する。）法がある。入力xτ（i.e. 音響特徴量）を得たもとでの出力yτ（i.e.事前SNR）のMMSE推定量^yτは、以下のように計算できる。
【０００８】
【数２】
0006370751.tif 000003
【０００９】
ここで、Θは、後述するxτとyτの結合分布のパラメータである。つまり、MMSE推定量^yτは、入力xτを得たもとでの出力yτの期待値である。出力yτの推定精度を向上させるためには、学習データから条件付き分布を最適化する必要がある。ベイズの定理より、条件付き分布の最適化は、xτとyτの結合分布p(xτ,yτ;Θ)のパラメータΘを最適化すればよい。
【００１０】
この結合分布には、例えばガウス混合モデル（Gaussian Mixture Model: GMM）を用いることができる。結合変数zτ=(xτT,yτT)Tの分布p(xτ,yτ;Θ)をGMMで表現する場合、以下のような形になる。
【００１１】
【数３】
0006370751.tif 000004
【００１２】
ここで、Tは転置、N(μ,Σ)は平均μ、共分散行列Σの正規分布、KはGMMの混合数、μkzはk番目の正規分布の平均ベクトル、Σkzはk番目の正規分布の共分散行列である。また、Θ={wk,μkz,Σkz}k=1Kである。さらに、μkz,Σkzは、以下のように書き下すことができる。
【００１３】
【数４】
0006370751.tif 000005
【００１４】
μkx,μkyはそれぞれk番目の正規分布のxとyの平均ベクトル、Σkxx,Σkyyはそれぞれk 番目の正規分布のxとyの共分散行列、Σkxy,Σkyxはそれぞれk 番目の正規分布のxとyの各次元の相関を表す相関行列である。
【００１５】
上述のGMMに基づく、入力xτからyτのMMSE推定量^yτは以下のように求めることができる。
【００１６】
【数５】
0006370751.tif 000006
【００１７】
結合分布にGMMを用いたMMSE推定法は、その計算の容易さから、声質変換（例えば、非特許文献１参照。）など、音響処理の分野で幅広く用いられている。
【先行技術文献】
【非特許文献】
【００１８】
【非特許文献１】T. Toda, A.W. Black, K. Tokuda. “Voice conversion based on maximum likelihood estimation of spectral parameter trajectory,” IEEE Transactions on Audio, Speech and Language Processing, Vol. 15, No. 8, pp. 2222-2235, Nov. 2007.
【非特許文献２】小西貞則, “多変量解析入門, 9 章主成分分析” pp.225-244, 岩波書店, 2010.
【発明の概要】
【発明が解決しようとする課題】
【００１９】
ところが、入力xτの次元Dが、学習データに対して大きすぎる場合（例えばT=3000個程度の学習データで、D=1024次元の入力を扱う場合）、結合分布p(xτ,yτ;Θ)を正確に学習することができない。このデータの次元に対して学習データ数が足りず、学習が困難になる現象は「次元の呪い」と呼ばれている。こういった状況では結合分布を学習する前に、入力xτの次元を削減しなくてはならない。
【００２０】
入力xτの次元を小さくする方法として、２つのアプローチがある。１つは、開発者の経験や直感に基づき、入力xτに用いる音響特徴量を人手で設計/精査することである。例えばシュート音やバッティング音を収音したいとき、突発性の音源であるという特徴が事前に既知であるので、ある帯域における時間的なパワー差分を使えば、事前SNRを推定できそうだろう。つまり、「音響特徴量xτとしてある周波数帯域のパワー時間差分だけを音響特徴量として用いる」と経験的に決定する。しかし、シュート音のように、「突発性」という、その音を特徴づける性質が明確であればよいが、人間の声や動物の鳴き声など、音を特徴づける性質が未知であることもある。こういった場合、強調/抑圧したい音源ごとに特徴量を人手で設計する必要があり労力を要する。また目的音の時間周波数構造が複雑な場合は、ルールや特徴量自体を設計することが困難である。
【００２１】
そこで２つ目のアプローチとして、大量の音響特徴量候補xτを圧縮行列A∈RD×Q（Q<D）を用いて次元圧縮する方法がある。つまり、人手で特徴量を精査するのをあきらめ、機械学習で有効な特徴量を探索する。そのため、xτとしては考え付く限りの大量の次元Dの音響特徴量を用い、ATxτとなる変換を行う。そして、結合分布p(ATxτ,yτ;Θ)を学習する。この圧縮行列Aの設計法の代表的なものとして、主成分分析がある（例えば、非特許文献２参照）。主成分分析では、入力xτの分散が最大となるような圧縮行列Aを設計する。
【００２２】
入力に対して圧縮処理を行うということは、入力xτの中から必要な情報を推定し、不要と判断された情報を捨てることに相当する。ところが、主成分分析の圧縮行列の設計には出力変数yτの情報を用いない。つまり、主成分分析では、出力yτを推定するために重要な情報を捨ててしまう可能性がある。そのため、主成分分析による次元圧縮は、GMMを用いたMMSE推定の枠組みから見たとき最適ではない。最適な圧縮行列Aを設計するためには、結合分布p(xτ,yτ;Θ)の最適化と同時に、圧縮行列Aを最適化する必要がある。
【００２３】
この発明の目的は、ガウス混合モデルパラメータを求めるときに圧縮行列Aの最適化を行うガウス混合モデルパラメータ計算装置、情報推定装置、音強調装置、これらの方法及びプログラムを提供することである。
【課題を解決するための手段】
【００２４】
この発明の一態様によるガウス混合モデルパラメータ計算装置は、・Tは行列又はベクトル・の転置を表すとして、入力変数xτ及び出力変数yτの学習データに基づいて、ガウス混合モデルで表現された結合分布p(ATxτ,yτ;Θ)のガウス混合モデルパラメータΘ及び圧縮行列Aを求めるガウス混合モデルパラメータ計算装置であって、入力されたガウス混合モデルパラメータΘに基づいて、EMアルゴリズムにおける負担率を計算する負担率計算部と、入力されたガウス混合モデルパラメータΘ及び上記負担率に基づいてガウス混合モデルパラメータΘを更新するガウス混合モデルパラメータ更新部と、入力されたガウス混合モデルパラメータΘ、上記負担率及び入力された圧縮行列Aに基づいて、上記入力された圧縮行列Aを勾配法を用いた最適化により更新する圧縮行列更新部と、上記更新されたガウス混合モデルパラメータΘ及び上記更新された圧縮行列Aを上記入力されたガウス混合モデルパラメータΘ及び上記入力された圧縮行列Aとして、上記負担率計算部、上記ガウス混合モデルパラメータ更新部及び上記圧縮行列更新部の処理を繰り返す制御を行う制御部と、を備えている。
【００２５】
この発明の一態様による情報推定装置は、上記ガウス混合モデルパラメータ計算装置と、上記繰り返し制御により最終的に生成された圧縮行列Aを用いて、情報推定の対象となる入力変数xτを圧縮したATxτを計算する次元圧縮部と、上記繰り返し制御により最終的に生成されたガウス混合モデルパラメータΘを用いて、上記計算されたATxτが得られた下でのyτの期待値である^yτを計算する推定部と、を備えている。
【００２６】
この発明の一態様による音強調装置は、強調の対象となる目的音が含まれる音信号を周波数領域信号に変換する周波数領域変換部と、上記周波数領域信号に基づいて、強調の対象となる目的音が含まれる音信号の音響特徴量xτを抽出する音響特徴量抽出部と、請求項３のガウス混合モデルパラメータ計算装置により最終的に生成された圧縮行列Aを用いて、上記抽出された音響特徴量xτを圧縮したATxτを計算する次元圧縮部と、請求項３のガウス混合モデルパラメータ計算装置により最終的に生成されたガウス混合モデルパラメータΘを用いて、上記計算されたATxτが得られた下での事前信号雑音比yτの期待値である^yτを計算する推定部と、上記計算された^yτを用いてウィナーフィルタを求めるウィナーフィルタ計算部と、上記求まったウィナーフィルタを上記周波数領域信号に適用してフィルタリング後周波数領域信号を得るフィルタリング部と、上記フィルタリング後周波数領域信号を時間領域信号に変換する時間領域変換部と、を備えている。
【発明の効果】
【００２７】
ガウス混合モデルパラメータを求めるときに圧縮行列Aの最適化が行われる。
【図面の簡単な説明】
【００２８】
【図１】ガウス混合モデルパラメータ計算装置の例を説明するためのブロック図。
【図２】情報推定装置の例を説明するためのブロック図。
【図３】音強調装置の例を説明するためのブロック図。
【図４】画像の判別を行うためのガウス混合モデルパラメータを計算するガウス混合モデルパラメータ計算装置の例を説明するためのブロック図。
【図５】ガウス混合モデルパラメータ計算方法の例を説明するための流れ図。
【図６】情報推定方法の例を説明するための流れ図。
【図７】音強調方法の例を説明するための流れ図。
【発明を実施するための形態】
【００２９】
以下、図面を参照して、この発明の一実施形態であるガウス混合モデルパラメータ計算装置、情報推定装置、音強調装置及びこれらの方法について説明する。
【００３０】
まず、ガウス混合モデルパラメータ計算装置及び方法について説明し、その後にガウス混合モデルパラメータ計算装置及び方法を用いた情報推定装置及び方法、音強調装置及び方法について説明する。
【００３１】
［ガウス混合モデルパラメータ計算装置及び方法］
ガウス混合モデルパラメータ計算装置は、図１に示すように、初期化部１１、負担率計算部１２、ガウス混合モデルパラメータ更新部１３、圧縮行列更新部１４及び制御部１５を例えば備えている。
【００３２】
ガウス混合モデルパラメータ計算装置が、図５に例示する各ステップの処理を行うことによりガウス混合モデルパラメータ計算方法が実現される。
【００３３】
ガウス混合モデルパラメータ計算装置は、通常のGMMの学習で用いられる期待値最大化（expectation-maximization: EM）アルゴリズムに、勾配法を用いた圧縮行列の最適化を入れ込んだ一般化EM（Generalized-EM: GEM）アルゴリズムを提供するものである。
【００３４】
ガウス混合モデルパラメータ計算装置１には、T個の入力変数xτの学習データx1,...,Tと出力変数yτの学習データy1,...,Tが入力される。入力の次元数はD、出力の次元数はΩとする。D,Ωは正の整数である。また、ガウス混合モデルパラメータ計算装置１には、圧縮後の次元数QとGMMの混合数K、また後述する勾配法のステップサイズεが入力される。Q、K、εの値は任意であり、解きたい問題によって設定を変更すべきだが、例えばQ=Ω、K=8、ε=10-2などに設定できる。
【００３５】
＜初期化部１１＞
初期化部１１は、GMMパラメータΘ={wk,μkν,Σkν}k=1Kと圧縮行列A∈RD×Qを初期化する（ステップＳ１１）。以下、GMMパラメータのことを、ガウス混合モデルパラメータとも表記する。
【００３６】
初期化部１１は、GMMパラメータΘはランダム初期化してもよいし、K-means アルゴリズムなどで初期化することもできる。また、初期化部１１は、圧縮行列Aはランダム初期化してもよいし、主成分分析などで初期化することもできる。
初期化されたGMMパラメータΘは、負担率計算部１２、ガウス混合モデルパラメータ更新部１３及び必要に応じて制御部１５に出力される。
【００３７】
初期化された圧縮行列Aは、負担率計算部１２、ガウス混合モデルパラメータ更新部１３、圧縮行列更新部１４及び必要に応じて制御部１５に出力される。
【００３８】
なお、GMMパラメータΘと圧縮行列Aは人手で初期化してもよい。この場合は、初期化部１１の処理に代えて、人手で初期化されたGMMパラメータΘと圧縮行列Aがガウス混合モデルパラメータ計算装置１に入力される。
【００３９】
＜負担率計算部１２＞
負担率計算部１２は、初期化されたGMMパラメータΘ、及び、初期化された圧縮行列Aを入力として、通常のEMアルゴリズムと同様に負担率γk,τを計算する（ステップＳ１２）。計算された負担率γk,τは、ガウス混合モデルパラメータ更新部１３及び圧縮行列更新部１４に出力される。
【００４０】
【数６】
0006370751.tif 000007
ここで、
【００４１】
【数７】
0006370751.tif 000008
であり、
【００４２】
【数８】
0006370751.tif 000009
である。また、
【００４３】
【数９】
0006370751.tif 000010
【００４４】
である。ここで、ξはATxτを一文字で表すために表記上導入した変数であり、μkξはそれぞれk番目の正規分布のATxτとyの平均ベクトル、Σkξξ,Σkyyはそれぞれk番目の正規分布のATxτとyの共分散行列、、Σkyξ,Σkξ yのはそれぞれk 番目の正規分布のATxτとyの各次元の相関を表す相関行列、wkはk番目の正規分布の混合比（重みパラメータ）である。
【００４５】
このようにして、負担率計算部１２は、入力されたガウス混合モデルパラメータΘに基づいて、EMアルゴリズムにおける負担率を計算する。
【００４６】
＜ガウス混合モデルパラメータ更新部１３＞
ガウス混合モデルパラメータ更新部１３は、初期化されたGMMパラメータΘ、初期化された圧縮行列A、及び、負担率γk,τを入力として、以下の式でGMMパラメータΘ={wk,μkν,Σkν}k=1Kを更新する（ステップＳ１３）。更新されたGMMパラメータΘは、圧縮行列更新部１４及び制御部１５に出力される。
【００４７】
【数１０】
0006370751.tif 000011
【００４８】
このようにして、ガウス混合モデルパラメータ更新部１３は、入力されたガウス混合モデルパラメータΘ及び上記負担率に基づいてガウス混合モデルパラメータΘを更新する。
【００４９】
＜圧縮行列更新部１４＞
圧縮行列更新部１４は、初期化された圧縮行列A、負担率γk,τ、及び、更新されたガウス混合モデルパラメータΘを入力として、圧縮行列A∈RD×Qを勾配法で更新する（ステップＳ１４）。勾配法の具体的な手段には、最急降下法や確率的最急降下法などを用いることができる。例えば最急降下法の場合は、圧縮行列Aの各行ベクトルa1,...,Dを例えば以下の式により更新する（ステップＳ１４１）。
【００５０】
【数１１】
0006370751.tif 000012
【００５１】
ここで、εはステップサイズである。勾配ベクトル∇adを算出する。まず、圧縮行列演算を以下のように書き下す。
【００５２】
【数１２】
0006370751.tif 000013
【００５３】
また、結合ガウス分布の性質を用いて、各ガウス分布の指数部分の二次形式をyτの項、xτの項、クロス項の3つに展開することで、勾配ベクトル∇adは以下のように計算できる。
【００５４】
【数１３】
0006370751.tif 000014
【００５５】
ただし、ΛkyξとΛkξξは以下の式で求まる分割行列である。なお、分割行列のサイズはΣkyξ,Σkξξに等しい。
【００５６】
【数１４】
0006370751.tif 000015
【００５７】
また、圧縮後の入力変数ATxτの分散を大きくしたいときは、勾配ベクトル∇adを以下のように計算してもよい。
【００５８】
【数１５】
0006370751.tif 000016
【００５９】
圧縮行列更新部１４は、以下に例示するステップＳ１４２からステップＳ１４５の処理において、圧縮行列Aの各基底が直交するように、最急降下法の各ステップで直行化する。
【００６０】
圧縮行列更新部１４は、圧縮行列Aを、ATAの最大固有値で割る（ステップＳ１４２）。
【００６１】
圧縮行列更新部１４は、A←(3/2)A-(1/2)AATAとする（ステップＳ１４３）。すなわち、(3/2)A-(1/2)AATAを計算して、新たなAとする。
【００６２】
圧縮行列更新部１４は、Aの各列を、ノルムが1 になるように正規化する（ステップＳ１４４）。
【００６３】
圧縮行列更新部１４は、Aが収束するまで、ステップＳ１４１からステップＳ１４４の処理を繰り返す（ステップＳ１４５）。
【００６４】
収束判定は、例えば、事前に指定した繰り返し回数に達するかどうかにより判定することができる。すなわち、圧縮行列更新部１４は、ステップＳ１４１からステップＳ１４４の処理を行った回数をカウントして、カウントされた回数が事前に指定した繰り返し回数に達した場合には、Aが収束したと判断し、ステップＳ１４５の処理を終了する。事前に指定した繰り返し回数は、例えば100に設定することができる。
【００６５】
また、収束判定は、ステップＳ１４２の処理開始前の圧縮行列AとステップＳ１４５の処理終了後の圧縮行列のフロベニウスノルムが閾値λ1以下であるかを判定することにより行ってもよい。すなわち、圧縮行列更新部１４は、ステップＳ１４２の処理開始前の圧縮行列AとステップＳ１４５の処理終了後の圧縮行列のフロベニウスノルムが閾値λ1以下である場合には、Aが収束したと判断し、ステップＳ１４５の処理を終了する。閾値λ1は、例えば10-4などに設定できる。
【００６６】
このようにして、圧縮行列更新部１４は、入力されたガウス混合モデルパラメータΘ、負担率及び入力された圧縮行列Aに基づいて、入力された圧縮行列Aを勾配法を用いた最適化により更新する。
【００６７】
＜制御部１５＞
制御部１５は、GMMパラメータΘと圧縮行列Aの更新が収束するまで、ステップＳ１２からステップＳ１４の処理を繰り返す（ステップＳ１５）。
【００６８】
収束判定は、例えば、事前に指定した繰り返し回数に達するかどうかにより判定することができる。すなわち、制御部１５は、ステップＳ１２からステップＳ１４の処理を行った回数をカウントして、カウントされた回数が事前に指定した繰り返し回数に達した場合には、GMMパラメータΘと圧縮行列Aの更新が収束したと判断し、ステップＳ１４５の処理を終了する。事前に指定した繰り返し回数は、例えば1000に設定することができる。
【００６９】
また、収束判定は、以下の式で表される対数尤度Lの上昇が閾値λ2以下であるかを判定することにより行ってもよい。すなわち、制御部１５は、以下の式で例えば表される対数尤度Lの上昇が閾値λ2以下である場合には、GMMパラメータΘと圧縮行列Aの更新が収束したと判断し、制御部１５の処理を終了する。閾値λ2は、例えば10-5などに設定できる。対数尤度Lの上昇とは、更新後の対数尤度Lと更新前の対数尤度Lとの差のことである。
【００７０】
【数１６】
0006370751.tif 000017
【００７１】
このようにして、制御部１５は、更新されたガウス混合モデルパラメータΘ及び更新された圧縮行列Aを、負担率計算部１２、ガウス混合モデルパラメータ更新部１３及び圧縮行列更新部１４へ新たな入力として、負担率計算部１２、ガウス混合モデルパラメータ更新部１３及び圧縮行列更新部１４の処理を繰り返す制御を行う。
【００７２】
このように、ガウス混合モデルパラメータを学習する際に圧縮行列Aの最適化を行うことにより、入力xτが高次元な場合にも、最適な結合分布をGMMで学習できるようになる。また、結合分布の尤度を最大化する圧縮行列Aが求まるため、このようにして求まったガウス混合モデルパラメータΘ及び圧縮行列Aを用いてMMSE 推定を行った場合には、その精度が向上する。
【００７３】
また、ガウス混合モデルパラメータを学習する際に圧縮行列Aの最適化を行うことにより、これまでは人手で経験的に行ってきた入力変数の設計/選択を自動化することができる。これにより、このようにして求まったガウス混合モデルパラメータΘ及び圧縮行列Aを用いたMMSE 推定を様々な分野で行うことができる。
【００７４】
［情報推定装置及び方法］
情報推定装置は、図２に示すように、次元圧縮部５１及び推定部５２を例えば備えている。
【００７５】
情報推定装置が、図６に例示する各ステップの処理を行うことにより情報推定方法が実現される。
【００７６】
情報推定装置には、情報推定の対象となる入力変数xτが入力される。
【００７７】
＜次元圧縮部５１＞
次元圧縮部５１は、収束後の圧縮行列Aを入力として、情報推定の対象となる入力変数xτを、ガウス混合モデルパラメータ計算装置１及び方法で推定された圧縮行列Aで圧縮する（ステップＳ５１）。圧縮された入力変数xτであるATxτは、推定部５２に出力される。
【００７８】
すなわち、次元圧縮部５１は、ガウス混合モデルパラメータ計算装置１及び方法における繰り返し制御により最終的に生成された圧縮行列Aを用いて、情報推定の対象となる入力変数xτを圧縮したATxτを計算する。
【００７９】
＜推定部５２＞
推定部５２は、入力変数xτを圧縮したATxτを入力として、MMSE推定量を例えば以下のように求める（ステップＳ５２）。
【００８０】
【数１７】
0006370751.tif 000018
【００８１】
すなわち、推定部５２は、ガウス混合モデルパラメータ計算装置１及び方法における繰り返し制御により最終的に生成されたガウス混合モデルパラメータΘを用いて、計算されたATxτが得られた下でのyτの期待値である^yτを計算する。
【００８２】
このように、ガウス混合モデルパラメータを学習する際に圧縮行列Aの最適化を行うガウス混合モデルパラメータ計算装置１及び方法で求まったガウス混合モデルパラメータΘ及び圧縮行列Aを用いてMMSE 推定を行うことにより、入力xτが高次元な場合にもMMSE 推定を行うことができる。また、MMSE 推定精度を向上させることができる。また、MMSE 推定を様々な分野で行うことができる。
【００８３】
［音強調装置及び方法］
音強調装置は、図３に示すように、マイクロホン２、周波数領域変換部３、音響特徴量抽出部４、情報推定装置５、ウィナーフィルタ計算部６、フィルタリング部７、時間領域変換部８、周波数領域変換部９１、重畳部９２、音響特徴量抽出部９３、事前信号雑音比計算部９４及びガウス混合モデルパラメータ計算装置１を例えば備えている。
【００８４】
音強調装置が、図７に例示する各ステップの処理を行うことにより音強調方法が実現される。
【００８５】
音強調装置及び方法の処理は、学習処理（ステップＳ９１からステップＡ６の処理）と強調処理（ステップＳ２からステップＳ８の処理）の２個の処理に分けることができる。まず、学習処理について説明する。学習処理では、GEMアルゴリズムを用いたガウス混合モデルパラメータの学習が行われる。学習処理は、強調処理に先立って行われてもよいし、強調処理と並行して行われてもよい。
【００８６】
学習処理は、強調処理に先立って行われる場合には、音強調装置は、周波数領域変換部９１、重畳部９２、音響特徴量抽出部９３、事前信号雑音比計算部９４及びガウス混合モデルパラメータ計算装置１を備えていなくてもよい。
【００８７】
目的音の学習データと雑音の学習データの時間波形を用意する。ここで、時間波形のサンプリングレートや量子化ビット数は任意であるが、例えばサンプリングレートは48kHz、量子化ビット数は16bit などに設定できる。
【００８８】
＜周波数領域変換部９１＞
周波数領域変換部９１は、目的音の学習データ及び雑音の学習データを短時間フーリエ変換（STFT）などを用いて周波数領域に変換し、それぞれ収音したい特定の音（目的音）の複素スペクトルSω,τと雑音の複素スペクトルNω,τを得る（ステップＳ９１）。得られたSω,τとNω,τは、重畳部９２及び事前信号雑音比計算部９４に出力される。ここで、フーリエ変換長は1024 点（サンプリング周波数48kHz で約22ms）、シフト長は512 点（サンプリング周波数48kHz で約11ms）などに設定できる。
【００８９】
＜重畳部９２＞
重畳部９２は、目的音が含まれる音信号である観測信号を模擬的に設計するために、以下の式に基づきSω,τとNω,τを重畳し、Fω,τを得る（ステップＳ９２）。得られたFω,τは、音響特徴量抽出部９３に出力される。
Fω,τ=Sω,τ+Nω,τ
【００９０】
＜音響特徴量抽出部９３＞
音響特徴量抽出部９３は、入力されたFω,τから音響特徴量xτを抽出する（ステップＳ９３）。この抽出された音響特徴量xτは、ガウス混合モデルパラメータ計算装置１に出力され、GEMアルゴリズムを実行するガウス混合モデルパラメータ計算装置１における入力変数xτになる。
【００９１】
音響特徴量xτに用いる特徴量は任意であるが、例えば48次元のメル周波数ケプストラム係数（MFCC）およびその一階差分と二階差分や、48 次元のメルフィルタバンク出力（MFBO）およびその一階差分と二階差分などを用いることができる。また、観測に用いたマイクロホンの数MがM>1であれば、ビームフォーミングを行って、方向別にMFCC やMFBOを求めることもできる。またその他にも、スペクトルフラックスやスペクトルセントロイドなど、D = 512程度の様々な音響特徴量を用いることができる。
【００９２】
＜事前信号雑音比計算部９４＞
事前信号雑音比計算部９４は、以下の式に基づいて、学習データから得られたSω,τとNω,τを用いて事前SNRを計算する（ステップＳ９４）。計算された事前SNRは、ガウス混合モデルパラメータ計算装置１に出力される。この計算結果がGEMアルゴリズムを実行するガウス混合モデルパラメータ計算装置１における出力変数yτとなる。
【００９３】
【数１８】
0006370751.tif 000019
【００９４】
ここで、フーリエ変換長が大きいと、事前SNRの次元Ωも大きくなるため、上の式の計算結果をメルフィルタバンクで圧縮してもよい。この場合、圧縮された計算結果がガウス混合モデルパラメータ計算装置１における出力変数yτとなる。メルフィルタバンクの個数はたとえば32 程度に設定できる。
【００９５】
＜ガウス混合モデルパラメータ計算装置１＞
ガウス混合モデルパラメータ計算装置１は、音響特徴量抽出部９３が計算した音響特徴量を入力変数xτとし、事前信号雑音比計算部９４が計算した事前信号雑音比を出力変数yτとして、［ガウス混合モデルパラメータ計算装置及び方法］の欄で説明した手法により、ガウス混合モデルパラメータΘ及び圧縮行列Aを計算する（ステップＳ１）。計算されたガウス混合モデルパラメータΘ及び圧縮行列Aは、情報推定装置５に出力される。
【００９６】
なお、ガウス混合モデルパラメータ計算装置１に、圧縮後の次元数QとGMMの混合数K、勾配法のステップサイズεが入力されるが、最終的に音強調をする場合には、これらのパラメータは、目的音のスペクトルの複雑さに応じてチューニングすべきである。スポーツフィールドでキック音を強調する場合、キック音のスペクトル形状は単純なため、例えば、Q=12、K=12、ε=10-2などに設定できる。
【００９７】
以上が、音強調装置及び方法の学習処理である。以下、音強調装置及び方法の強調処理について説明する。
【００９８】
＜マイクロホン２＞
M個のマイクロホン１を用いて、強調の対象となる目的音が含まれる、所望の音が収音される（ステップＳ２）。ここで、Mは1以上の任意の自然数である。収音された音信号は、周波数領域変換部３に出力される。
【００９９】
＜周波数領域変換部３＞
周波数領域変換部３は、音信号を周波数領域へ変換し、周波数領域信号Fω，τを得る（ステップＳ３）。すなわち、周波数領域変換部３は、強調の対象となる目的音が含まれる音信号を周波数領域信号に変換する。得られた周波数領域信号Fω，τは、音響特徴量抽出部４に出力される。
【０１００】
＜音響特徴量抽出部４＞
音響特徴量抽出部４は、入力された周波数領域信号Fω，τから音響特徴量xτを計算する（ステップＳ４）。すなわち、音響特徴量抽出部４は、周波数領域信号Fω，τに基づいて、強調の対象となる目的音が含まれる音信号の音響特徴量xτを抽出する。計算された音響特徴量xτは、情報推定装置５に出力される。
【０１０１】
＜情報推定装置５＞
情報推定装置５は、ガウス混合モデルパラメータ計算装置１及び方法により計算されたガウス混合モデルパラメータΘ及び圧縮行列Aに基づいて、［情報推定装置及び方法］の欄で説明した手法により、事前信号雑音比yτの期待値である^yτを計算する（ステップＳ５）。
【０１０２】
すなわち、情報推定装置５の次元圧縮部５１は、ガウス混合モデルパラメータ計算装置１及び方法により最終的に生成された圧縮行列Aを用いて、抽出された音響特徴量xτを圧縮したATxτを計算する（ステップＳ５１）。
【０１０３】
また、情報推定装置５の推定部５２は、ガウス混合モデルパラメータ計算装置１及び方法により最終的に生成されたガウス混合モデルパラメータΘを用いて、計算されたATxτが得られた下での事前信号雑音比yτの期待値である^yτを計算する（ステップＳ５２）。
【０１０４】
＜ウィナーフィルタ計算部６＞
ウィナーフィルタ計算部６は、計算された^yτを用いてウィナーフィルタを求める（ステップＳ６）。求まったウィナーフィルタGω,τは、フィルタリング部７に出力される。
【０１０５】
具体的には、ウィナーフィルタ計算部６は、以下の式に基づき、推定した事前SNR^yτからウィナーフィルタGω,τを計算する。
【０１０６】
【数１９】
0006370751.tif 000020
【０１０７】
ただし、学習処理において、事前信号雑音比計算部９４が事前SNR をメルフィルタバンク圧縮していた場合、上記式の計算の実行前に、事前SNR^yτを線形周波数領域にスプライン補間する。
【０１０８】
＜フィルタリング部７＞
フィルタリング部７は、求まったウィナーフィルタを周波数領域信号Fω,τに適用してフィルタリング後周波数領域信号Hω,τを得る（ステップＳ７）。具体的には、フィルタリング部７は、以下の式 に基づき、フィルタリング後周波数領域信号Hω,τを得る。
【０１０９】
Hω,τ=Gω,τFω,τ
得られたフィルタリング後周波数領域信号Hω,τは、時間領域変換部８に出力される。
【０１１０】
＜時間領域変換部８＞
時間領域変換部８は、入力されたフィルタリング後周波数領域信号Hω,τを時間領域信号に変換する（ステップＳ８）。
【０１１１】
具体的には、時間領域変換部８は、フィルタリング後周波数領域信号Hω,τを逆フーリエ変換などの手法を用いて時間領域に変換する。
【０１１２】
時間領域信号は、目的音が強調された信号として音強調装置から出力される。
【０１１３】
［変形例］
ガウス混合モデルパラメータ装置、情報推定装置及びこれらの方法の用途は、音源強調に限らない。ガウス混合モデルパラメータ装置、情報推定装置及びこれらの方法は、音響強調以外の用途に用いることができる。
【０１１４】
つまり、入力変数xτは音響特徴量に限らず、画像等のなんらかのデータから特徴量を抽出した結果であってもよいし、位置座標などの生データであってもよい。同様に出力変数yτも事前SNR に限らず、画像の属するクラス等を表す変数でもよいし、元信号の振幅スペクトル|Sω,τ|であってもよい。
【０１１５】
例えば、画像データに対する判別を行うMMSE推定を行うためのガウス混合モデルパラメータΘ及び圧縮行列Aを求める場合には、例えば以下のようにすればよい。判別対象の画像とそれ以外の画像の各データが学習データとして、図４の特徴量抽出部９５に入力される。特徴量抽出部９５は、コーナー検出やSHIFT（Scale-Invariant Feature Transform）特徴などを用いて、判別対象の画像とそれ以外の画像の各データから特徴量を抽出して入力変数xτとする。また、出力変数yτは、入力変数xτに対応する画像が判別対象の画像ならばy=1それ以外はy=0となるような2値変数とする。ガウス混合モデルパラメータ計算装置１は、これらの入力変数xτ及び出力変数yτに基づいて、［ガウス混合モデルパラメータ計算装置及び方法］の欄で説明した処理を行う。
【０１１６】
その他、この発明の趣旨を逸脱しない範囲で適宜変更が可能であることはいうまでもない。
【０１１７】
[プログラム及び記録媒体]
ガウス混合モデルパラメータ計算装置、情報推定装置、音強調装置及びこれらの方法において説明した処理は、記載の順にしたがって時系列に実行されるのみならず、処理を実行する装置の処理能力あるいは必要に応じて並列的にあるいは個別に実行されてもよい。
【０１１８】
また、ガウス混合モデルパラメータ計算装置、情報推定装置及び音強調装置における各処理をコンピュータによって実現する場合、ガウス混合モデルパラメータ計算装置、情報推定装置及び音強調装置が有すべき機能の処理内容はプログラムによって記述される。そして、このプログラムをコンピュータで実行することにより、その各処理がコンピュータ上で実現される。
【０１１９】
この処理内容を記述したプログラムは、コンピュータで読み取り可能な記録媒体に記録しておくことができる。コンピュータで読み取り可能な記録媒体としては、例えば、磁気記録装置、光ディスク、光磁気記録媒体、半導体メモリ等どのようなものでもよい。
【０１２０】
また、各処理手段は、コンピュータ上で所定のプログラムを実行させることにより構成することにしてもよいし、これらの処理内容の少なくとも一部をハードウェア的に実現することとしてもよい。
【産業上の利用可能性】
【０１２１】
この技術は、例えばサッカー等のスポーツ中継などで、ボールのキック音や選手の叫び声、ホイッスルなど、ある特定の競技音だけをクリアに抽出するために例えば用いることができる。もちろん、スポーツ中継以外の分野において、特定の音を強調するためにも用いることができる。
【符号の説明】
【０１２２】
１      ガウス混合モデルパラメータ計算装置
１１    初期化部
１２    負担率計算部
１３    ガウス混合モデルパラメータ更新部
１４    圧縮行列更新部
１５    制御部
２      マイクロホン
３      周波数領域変換部
４      音響特徴量抽出部
５      情報推定装置
５１    次元圧縮部
５２    推定部
６      ウィナーフィルタ計算部
７      フィルタリング部
８      時間領域変換部
９１    周波数領域変換部
９２    重畳部
９３    音響特徴量抽出部
９４    事前信号雑音比計算部
９５    特徴量抽出部
【図１】
0006370751.tif 000021
【図２】
0006370751.tif 000022
【図３】
0006370751.tif 000023
【図４】
0006370751.tif 000024
【図５】
0006370751.tif 000025
【図６】
0006370751.tif 000026
【図７】
0006370751.tif 000027