説明

音源パラメータ推定装置と音源分離装置とそれらの方法とプログラム

【課題】音源モデルパラメータが予め与えられていなくとも音源パラメータと一緒に音源モデルパラメータも推定できる音源パラメータ推定装置を提供する。
【解決手段】音源モデルパラメータ更新部は、音源パワー特徴量と音源パワーパラメータと音源占有度と、音源モデル記憶部に記憶された音源パワーパラメータの事前確率密度関数と音源パワー特徴量のモデルとを入力として音源モデルパラメータを更新する。音源占有度更新部は、音源位置特徴量と音源パワー特徴量と各音源の更新された音源パワーパラメータと音源位置パラメータと、音源モデルパラメータと、音源モデル記憶部に記憶された音源パワーパラメータの事前確率密度関数と音源パワー特徴量のモデルとを入力として各音源の音源占有度を更新する。

【発明の詳細な説明】
【技術分野】
【0001】
この発明は、複数の音源が同時に生成した音響信号が混ざって複数のマイクロホンで収音された観測信号から、各音源の音源パラメータを推定する音源パラメータ推定装置と、その音源パラメータに基づいて各音源を分離する音源分離装置と、それらの方法とプログラムに関する。
【背景技術】
【0002】
図7に、非特許文献1に開示された従来の音源パラメータ推定装置900の機能構成を示す。音源パラメータ推定装置900は、音源モデル記憶部920と、特徴抽出部910と、音源パワーパラメータ更新部930と、音源位置パラメータ更新部940と、音源占有度更新部950と、を具備する。
【0003】
音源モデル記憶部920は、複数の音源信号それぞれの音源パワー時系列全体の状態を表す音源パワーパラメータの事前確率密度関数p(q(l))((l)はl番目の音源)と、その音源パワーパラメータq(l)が与えられた場合の各音源信号の各時間周波数点(n,k)における事後確率密度関数である音源パワー特徴量のモデルf=βq(l),n,k(S)とを記憶する。特徴抽出部910は、複数の音源信号を複数(m個)のマイクロホンで収音した時間領域の信号を時間周波数領域信号に変換した観測信号x(m)n,kを入力として、各時間周波数点における音源位置特徴量An,kと音源パワー特徴量Xn,kを抽出する。
【0004】
音源パワーパラメータ更新部930は、Ns個の音源ごとに音源パワー特徴量Xn,kと観測信号が得られた下での占有的な音源の事後確率密度関数である音源占有度^M(l)n,kと、音源パワー特徴量のモデルf=βq(l),n,k(S)と音源パワーパラメータの事前確率密度関数p(q(l))とを入力として、音源パワーパラメータ^q(l)を更新する。音源位置パラメータ更新部940は、音源位置特徴量An,kと音源占有度^M(l)n,kを入力として、各音源の音源位置パラメータ^φ(l)を更新する。
【0005】
音源占有度更新部950は、音源位置特徴量An,kと音源パワー特徴量Xn,kと各音源の更新された音源パワーパラメータ^q(l)と音源位置パラメータ^φ(l)と、音源パワー特徴量のモデルf=βq(l),n,k(S)と音源パワーパラメータの事前確率密度関数p(q(l))とを基に各音源の音源占有度^M(l)n,kを更新する。
【0006】
音源パラメータ推定装置900の音源パラメータ推定技術を用いた音源分離装置は、音源パラメータ推定装置900が出力する音源パワー特徴量と、更新した音源占有度と音源パワーパラメータと、各音源信号の各時間周波数点における音源パワー特徴量のモデルとを入力として複数の音源のそれぞれの音源分離信号を、最小自乗誤差推定により求める音源分離部を更に備える。
【先行技術文献】
【非特許文献】
【0007】
【非特許文献1】Tomohiro Nakatani, Shoko Araki, Takuya Yoshioka, Masakiyo Fujimoto, “Multichannel source separation based on source location cue with log-spectral shaping by hidden Markov source model,” Proc. of Interspeech-2010,pp. 2766-2769, Sep., 2010.
【非特許文献2】中谷智広、荒木章子、吉岡卓也、藤本雅清、“DOAクラスタリングと音声の対数スペクトルHMMに基づく音源分離”日本音響学会2010年秋季研究発表会講演論文集、pp.577-580, 9月, 2010年.
【発明の概要】
【発明が解決しようとする課題】
【0008】
従来の音源パラメータ推定装置は、音源モデル記憶部が記憶する音源パワーパラメータの事前確率密度関数と音源パワー特徴量のモデルのそれぞれの挙動を制御するパラメータである音源モデルパラメータが事前に与えられなければ音源パラメータの推定を行うことができなかった。
【0009】
この発明はこの課題に鑑みてなされたものであり、音源モデルパラメータが事前に与えられていない場合でも、音源パラメータの一部としてその他の音源パラメータと一緒に音源モデルパラメータをも推定できる音源パラメータ推定装置と音源分離装置と、それらの方法とプログラムを提供することを目的とする。
【課題を解決するための手段】
【0010】
この発明の音源パラメータ推定装置は、音源モデル記憶部と、特徴抽出部と、音源パワーパラメータ更新部と、音源位置パラメータ更新部と、音源モデルパラメータ更新部と、音源占有度更新部と、を具備する。音源モデル記憶部は、複数の音源信号それぞれの音源パワー時系列全体の状態を表す音源パワーパラメータの事前確率密度関数と、その音源パワーパラメータが与えられた場合の各音源信号の各時間周波数点における音源パワー特徴量のモデルとを記憶する。特徴抽出部は、複数の音源信号を複数のマイクロホンで収音した時間領域信号を時間周波数領域信号に変換した観測信号を入力として、各時間周波数点における音源位置特徴量と音源パワー特徴量を抽出する。音源パワーパラメータ更新部は、音源パワー特徴量と、観測信号が得られた下での占有的な音源の事後確率密度関数である音源占有度と、音源パワーパラメータの事前確率密度関数と各音源信号の音源パワー特徴量のモデルと、音源パワー特徴量のモデルと音源パワーパラメータの事前確率密度関数のそれぞれの挙動を制御するパラメータである音源モデルパラメータと、を入力として各音源の音源パワーパラメータを更新する。音源位置パラメータ更新部は、音源位置特徴量と音源占有度を入力として、各音源の音源位置パラメータを更新する。音源モデルパラメータ更新部は、音源パワー特徴量と音源パワーパラメータと音源占有度と、音源モデル記憶部に記憶された音源パワーパラメータの事前確率密度関数と音源パワー特徴量のモデルとを入力として音源モデルパラメータを更新する。音源占有度更新部は、音源位置特徴量と音源パワー特徴量と各音源の更新された音源パワーパラメータと音源位置パラメータと、音源モデルパラメータと、音源モデル記憶部に記憶された音源パワーパラメータの事前確率密度関数と音源パワー特徴量のモデルとを入力として各音源の音源占有度を更新する。
【発明の効果】
【0011】
この発明の音源パラメータ推定装置によれば、音源モデルパラメータが事前に与えられていない場合でも、音源パラメータの一部としてその他の音源パラメータと一緒に音源モデルパラメータをも推定することができるので、事前に統計的性質を与えることができない多様な音源に対しても最適な音源パラメータを与えることができる。
【0012】
その音源パラメータ推定装置を用いたこの発明の音源分離装置は、誤差の少ない音源分離信号を出力することが可能である。
【図面の簡単な説明】
【0013】
【図1】この発明の音源パラメータ推定装置100の機能構成例を示す図。
【図2】音源パラメータ推定装置100の動作フローを示す図。
【図3】この発明の音源分離装置200の機能構成例を示す図。
【図4】音源分離装置200の動作フローを示す図。
【図5】実験結果を示す図であり、(a)は混合前の音声から推定した結果、(b)は混合音から推定した結果を示す図である。
【図6】実験結果を示す図であり、分離後の信号のケプストラム歪みと混合音の数との関係を示す、(a)は混合条件1の混合音に対して処理した結果、(b)と(c)は残響のある別の環境で収録した混合音(混合条件2と3)に対して処理した結果を示す図である。
【図7】従来の音源パラメータ推定装置900の機能構成を示す図。
【発明を実施するための形態】
【0014】
以下、この発明の実施の形態を図面を参照して説明する。複数の図面中同一のものには同じ参照符号を付し、説明は繰り返さない。実施例の説明の前にこの発明の基本的な考えについて説明する。
【0015】
〔この発明の基本的な考え〕
この発明では、予め音源モデルパラメータが与えられていなくとも音源パラメータの推定が行えるようにするために、l番目の音源信号の音源パワー特徴量の時系列全体を{S(l)n,k}と表した時に、その同時確率密度関数を音源モデルパラメータψ(l)を含めてモデル化すると共に、複数の音が混在した信号の音源パワー特徴量Xn,kが与えられたときに、各音源の音源占有度、音源パワー特徴量のモデルと音源パワーパラメータの事前確率密度関数、および音源パワーパラメータに基づいて各音源の音源モデルパラメータを更新する音源モデルパラメータ更新部を有する点で新しい。
【0016】
最初に説明に用いる記号について説明する。観測信号には、Ns個の音源信号が重畳しており、その音源信号をNm本のマイクロホンで収音する。m番目のマイクロホンから収音した収音信号を短時間フーリエ変換等を用いて周波数領域の信号に変換した観測信号をX(m)n,kと表記する。nはn番目の時間つまりフレーム番号、kはk番目の周波数つまりビン番号であり、n番目の時間及びk番目の周波数に対応する時間周波数点を参照する場合に、時間周波数点(n,k)と表記する。なお、記号^の位置や添え字の表記とその位置は、式中の表記が正しい。
【0017】
この発明では、l番目の音源信号の音源パワー時系列全体{S(l)n,k}の同時確率密度関数を次式に示すようにモデル化する。
【0018】
【数1】

【0019】
ここで、q(l)はl番目の音源の音源パワー時系列全体の状態を表す音源パワーパラメータを表す。以下では全ての音源のq(l)をまとめてq=[q(1),…,q(Ns)]とも表記する。ψ(l)はl番目の音源の音源モデルパラメータ全体を現す。全ての音源のψ(l)をまとめてψ=[ψ(1),…,ψ(Ns)]とも表記する。
【0020】
また、βq(l),n,k,ψ(l)(S)は音源パワー特徴量のモデルであり、音源パワーパラメータq(l)と音源モデルパラメータψ(l)が与えられた下で各時間周波数点(n,k)の音源信号の音源パワーがS(l)n,kとなる確率密度関数である(式(3))。式(1)の総和演算は、q(l)が離散値ではなく連続値をとる場合にはq(l)に関する積分演算に置き換えて表現されるものとする。式(2)において、音源の状態が既知のもとでは、異なる時間周波数点における音源パワーS(l)n,kは相互に独立であるという仮定を導入している。
【0021】
また、この発明では式(4)に示すように、各時間周波数点(n,k)において最も大きなエネルギーを持つ音源信号(以下、占有的な音源信号と称する)の音源パワーS(l)n,kは、観測信号の音源パワーと一致すると仮定する。
【0022】
【数2】

【0023】
また、占有的ではない音源lに関しては、S(l)n,k≦Xn,kの関係を持つと仮定する。すると、各音源信号の状態が既知の条件の下で、観測信号の音源パワーXn,kの事後確率密度関数は次のように表現できることが知られている(参考文献:S. J. Rennie, J.R. Hershey, and P. A.01sen, “Hierarchical variational loopy belief propagation for multi-talker speech recognition,” Proc. ASRU-2009, pp. 176-181 2009.)。
【0024】
【数3】

【0025】
この発明では、更に上式は次のように分解可能であると仮定している。
【0026】
【数4】

【0027】
また、この発明では音源位置特徴量から音源位置パラメータφ^(l)を推定するため、音源位置特徴量のモデルp(An,k;φ)を導入する。音源位置特徴量のモデルp(An,k;φ)は、各音源信号のエネルギーは異なる時間周波数点にわたり疎に分布していると仮定し、その時間周波数点において占有的な音源の音源位置のみに依存して決まると仮定する。
【0028】
一般的に、全ての音源の音源位置パラメータφ(l)をまとめてφ=[ψ(1),…,ψ(Ns)]と表すと、音源位置特徴量のモデルp(An,k;φ)、つまり観測信号の音源位置特徴量の確率密度関数は、混合分布として式(8)に示すように展開することができる。
【0029】
【数5】

【0030】
式(8)において、Zn,kは時間周波数点(n,k)において占有的な音源の番号を表す確率変数であり、Zn,k=lは、l番目の音源が占有的な音源である場合を示す。また、p(Zn,k=l)は、l番目の音源が時間周波数点(n,k)において占有的な音源になる事前確率密度関数を表している。更に、以降の説明では次の表記を用いることにする。
【0031】
【数6】

【0032】
γφ(l),n,k(A)は、時間周波数点(n,k)において占有的な音源の番号がlの場合に、音源位置特徴量An,kが得られる確率密度関数を表す。これは、l番目の音源の音源位置パラメータφ(l)のみに依存するものとする。具体的なγφ(l),n,k(A)やφ(l)の定義については後述する。
【0033】
式(8)のもと、γφ(l),n,k(A)が定義されている場合、音源位置パラメータφ(l)と占有的な音源の番号に関する事前確率密度関数p(Zn,k=l)が与えられれば、音源位置特徴量のモデルp(An,k;φ)は一意に定めることができる。逆に、音源位置特徴量An,kが観測された場合に、最尤推定などの方法に従い、音源位置パラメータと占有的な音源の番号に関する事前確率密度関数p(Zn,k=l)やその事後確率密度関数を推定することができる。
【0034】
以上の定義に従うと、完全データの確率密度関数は式(10)に示すように導出される。
【0035】
【数7】

【0036】
式(10)において、qが音源パワーパラメータ、ψが音源モデルパラメータ、φが音源位置パラメータであり、これらのパラメータがパラメータ推定の対象である。この発明では、次の対数尤度関数を最大化する値として、音源パワーパラメータqと音源モデルパラメータψと音源位置パラメータφを推定する。
【0037】
【数8】

【0038】
式(12)で、確率変数Zn,kは隠れ変数として扱われる。隠れ変数を含む対数尤度関数の最大化には、例えば、期待値最大化アルゴリズムなどを用いることができる。期待値最大化アルゴリズムでは、音源パワーパラメータの推定値^qと音源位置パラメータの推定値^φと音源モデルパラメータの推定値^ψに基づき、観測信号が得られた下での占有的な音源の番号の事後確率密度関数^M(l)n,k=p(Zn,k|An,k,Xn,k,^q;^φ,^ψ)をも同時に推定する必要がある。この発明では、この関数の値を音源占有度と称し、この値も音源パラメータに含めて考える。
【0039】
以上説明した考えで、音源パワー特徴量のモデルβq(l),n,k,ψ(l)(S)と、音源位置特徴量のモデルp(An,k;φ)の両者を考慮しながら最適な音源パラメータを推定することで音源パラメータの推定誤差を減らすことができる。また、音源位置特徴量のモデルp(An,k;φ)(式(8))と、音源パワー特徴量のモデル(式(7))に、占有的な音源の番号を表す変数Zn,kを共有化することで、2つの特徴量を考慮しながら音源パラメータ推定の計算を簡単にすることができる。
【0040】
以上、述べたように、この発明の音源パラメータ推定方法によれば、l番目の音源信号の音源パワー時系列全体{S(l)n,k}の同時確率密度関数を音源モデルパラメータψ(l)を含めてモデル化すると共に、複数の音が混在した信号の音源パワー特徴量Xn,kが与えられたときに、各音源の音源占有度、音源パワー特徴量のモデルと音源パワーパラメータの事前確率密度関数、および音源パワーパラメータに基づいて各音源の音源モデルパラメータを更新する音源モデルパラメータ更新部を有することで、予め音源モデルパラメータが与えられていなくとも音源パラメータの推定が行える。
【実施例1】
【0041】
図1にこの発明の音源パラメータ推定装置100の機能構成例を示す。その動作フローを図2に示す。音源パラメータ推定装置100は、特徴抽出部910と、音源モデル記憶部90と、音源の数に対応した数の音源パワーパラメータ更新部60〜60Nsと、音源パワーパラメータ更新部60〜60Nsと同じ数の音源位置パラメータ更新部940〜940Nsと、音源占有度更新部70と、音源モデルパラメータ更新部80と、を具備する。その各部の機能は、例えばROM、RAM、CPU等で構成されるコンピュータに所定のプログラムが読み込まれて、CPUがそのプログラムを実行することで実現されるものである。
【0042】
音源パラメータ推定装置100は、従来技術で説明した音源パワーパラメータ推定装置900に対して音源モデルパラメータ更新部80を備える点で異なる。特徴抽出部910と音源位置パラメータ更新部940〜940Nsとは、参照符号から明らかなように音源パワーパラメータ推定装置900と同じものである。
【0043】
各機能部の動作を説明する。
【0044】
〔特徴抽出部〕
特徴抽出部910は、複数の音源信号を複数(m本)のマイクロホンで収音した時間領域信号を時間周波数領域信号に変換した観測信号x(m)n,kを入力として、各時間周波数点(n,k)における音源位置特徴量An,kと音源パワー特徴量Xn,kを抽出する。
【0045】
音源パワー特徴量Xn,kは、例えば、1本目のマイクロホンが収音した信号の対数パワースペクトルを音源パワー特徴量として抽出する場合には式(13)に示すように計算される。音源パワー特徴量Xn,kは、音源占有度更新部70と音源パワーパラメータ更新部60〜60Nsに入力される。
【0046】
【数9】

【0047】
音源位置特徴量An,kは、一般に各時間周波数点における異なるマイクロホン間での信号の位相差や強度比などに表れる。したがって、音源位置特徴量An,kは、信号の位相差や強度比を異なるマイクロホンペアごとにまとめて出来るベクトルであったり、そこから更に何らかの特徴抽出を行った結果の値として抽出される。例えば、2本のマイクロホンで収音した信号の位相差を音源位置特徴量An,kとして抽出する場合、式(14)に示すように計算される。
【0048】
【数10】

【0049】
〔音源モデル記憶部〕
音源モデル記憶部90は、複数の音源信号それぞれの音源パワー時系列全体の状態を表す音源パワーパラメータq(l)の事前確率密度関数p(q(l)(l))と、その音源パワーパラメータq(l)が与えられた場合の各音源信号の各時間周波数点における音源パワー特徴量のモデルβq(l),n,k,ψ(l)(S)とを記憶する。(S)は音源パワー特徴量Xn,kを表す変数である。
【0050】
〔音源占有度更新部〕
音源占有度更新部70は、音源位置特徴量An,kと音源パワー特徴量Xn,kと各音源の更新された音源パワーパラメータ^q (l)と音源位置パラメータ^φ (l)と、音源モデルパラメータ^ψ (l)と、音源モデル記憶部90に記憶された音源パワーパラメータの事前確率密度関数p(q(l)(l))と音源パワー特徴量のモデルβq(l),n,k,ψ(l)(S)とを入力として上記各音源の音源占有度を更新する。
【0051】
音源占有度更新部70は、Σl^M(l)=1となるように、音源占有度^M (l)n,kを、例えば乱数で初期化する(ステップS70)。若しくは、従来技術で説明した音源パラメータ推定装置900と同じ方法を用いても良い。その後、音源パワーパラメータ更新部60〜60Nsと音源占有度更新部70と音源位置パラメータ更新部940〜940Nsと音源モデルパラメータ更新部80の各処理を収束するまで繰り返す。
【0052】
〔音源パワーパラメータ更新部〕
音源パワーパラメータ更新部60〜60Nsは、音源パワー特徴量Xn,kと、各音源lごとに初期化された音源占有度^M (l)n,kと、音源モデル記憶部90に記憶された音源パワーパラメータ^q (l)の事前確率密度関数p(q(l)(l))と音源パワー特徴量のモデルβq(l),n,k,ψ(l)(S)と音源モデルパラメータ^ψ (l)とを入力として、音源パワーパラメータ^q (l)を式(15)に示すように更新(M-step)する(ステップS60)。
【0053】
【数11】

【0054】
〔音源位置パラメータ更新部〕
音源位置パラメータ更新部940〜940Nsは、各音源lごとに初期化された音源占有度^M (l)n,kと、音源位置特徴量An,kを入力として音源位置パラメータ^φ (l)を、式(17)に示すように更新(M-step)する(ステップS940)。
【0055】
【数12】

【0056】
〔音源モデルパラメータ更新部〕
音源モデルパラメータ更新部80は、音源パワー特徴量Xn,kと音源パワーパラメータ^q (l)と音源占有度^M (l)n,kと、音源モデル記憶部90に記憶された音源パワーパラメータの事前確率密度関数p(q(l)(l))と音源パワー特徴量のモデルβq(l),n,k,ψ(l)(S)とを入力として音源モデルパラメータ^ψ (l)を式(15′)に示すように更新する(ステップS81)。
【0057】
【数13】

【0058】
そして、音源占有度更新部70は、各音源lごとに更新された音源パワーパラメータ^q(l)と音源モデルパラメータ^ψ(l)と音源位置特徴量An,kと音源パワー特徴量Xn,kと、音源モデル記憶部90に記憶された音源パワーパラメータの事前確率密度関数p(q(l)(l))と音源パワー特徴量のモデルβq(l),n,k,ψ(l)(S)と、を入力として、音源占有度^M (l)n,kを式(18)に示すように更新(E-step)する(ステップS71)。
【0059】
【数14】

【0060】
ステップS60〜ステップS71の処理は、収束が得られるまで繰り返される(ステップS72のno)。より具体的な音源位置特徴量のモデル及び、音源パワー特徴量のモデルを用いた実施例2を次に説明する。
【実施例2】
【0061】
先ず、特徴抽出部910は、式(14)に基づきマイク間位相差を、音源位置特徴量An,kとして抽出する。また、各音源lに由来する観測信号のマイク間位相差は、各周波数ごとに異なる平均値μ(l)k、分散σ(l)kのガウス分布に従うと仮定する。すると式(9)は以下のように定義できる。
【0062】
【数15】

【0063】
但し、φ(l)k=[μ(l)k(l)k]は、音源位置パラメータφ(l)のうち周波数kのみに関する部分を取り出したものであり、φ(l)は全ての周波数kについてφ(l)を集めたφ(l)=[φ(l)1,…, φ(l)Nk]である。N(・)は、ガウス分布の確率密度関数を表す。
【0064】
一方、特徴抽出部910は、式(13)に基づき、どれか一つのマイクロホン信号の対数パワースペクトルを音源パワー特徴量Xn,kとして抽出するものとする。さらに、音源パワーパラメータq(l)は、q(l)={q(l)0,q(l)1,…}のように各時刻の状態を表す状態系列に分解され、一次のマルコフ過程に従い状態遷移が各時刻で起こると仮定する。
【0065】
但し、音源パワーパラメータq(l)0は隠れマルコフモデルの初期状態を表す。更に、式(3)で定義される各時間周波数点(n,k)におけるS(l)n,kの事後確率密度関数は、その時刻の状態q(l)のみに依存するガウス分布に従うと仮定する。これを数式で表すと次のようになる。
【0066】
【数16】

【0067】
ここで、π(l)i=p(q(l)0=i)は、隠れマルコフモデルの初期状態がiである事前確率、α(l)i,j=p(q(l)n=j|q(l)n-1=i)は、隠れマルコフモデルが状態iから状態jへ移る状態遷移確率、βi,n,k,ψ(l)(S)=p(S(l)n,k=S|q(l)n=i;ψ(l))=N(S(l)n,k(l)i,k(l)i,k)は、隠れマルコフモデルの状態iにおける出力の確率密度関数であり、μ(l)i,k及びσ(l)i,kはその平均と分散である。
【0068】
この定義に基づくと、音源モデルパラメータ^ψ(l)は、全てのi,j,k,lに対するπ(l), α(l)i,j(l)i,k(l)i,kで構成される。この実施例では、全て若しくは一部の音源モデルパラメータが音源パラメータの一部として期待値最大化アルゴリズムにより推定される。
【0069】
以上の仮定の下、図2で説明済みの期待値最大化アルゴリズムのM-step1は、各音源lごとに、音源パワーパラメータ更新部60〜60Nsが式(22)を満たす状態時系列^q (l)=[^q (l)0,…,^q (l)Ns]を、Viterbiアルゴリズムを用いて更新する。
【0070】
【数17】

【0071】
また、M-step2は、各音源lごとに、音源位置パラメータ更新部940〜940Nsが、全ての周波数kで、φ(l)k=[μ(l)k(l)k]を次のように更新する。
【0072】
【数18】

【0073】
また、M-step3は、音源lごとに、音源モデルパラメータ更新部80が音源モデルパラメータ^ψ (l)を更新する。まず、π(l)(l)を、i,jに関するπ(l)(l) i,jの集合とすると、π(l)(l)は以下のように更新される。
【0074】
【数19】

【0075】
上記の更新は、隠れマルコフモデルの学習のための既知の方法で容易に実現することが可能である。一方、μ(l)(l)を全てのi,kに関するμ(l)i,k(l) i,kの集合とすると、μ(l)(l)の更新は、以下のように実現される。
【0076】
【数20】

【0077】
上記の更新を実現する一つの方法として、準ニュートン法、共役勾配法などに代表される逐次最大化アルゴリズムを上げることができる。これには、一般的に知られる多くのアルゴリズムを適用することができる。一方、上記の関数中に含まれるlog(βq(l)n,n,k,ψ(l)(Xn,k))は解析的な扱いが容易であるが、log(ρq(l)n,n,k,ψ(l)(Xn,k))は積分演算が含まれているため解析的な扱いが複雑になり、逐次最大化アルゴリズムに必要とされる計算を比較的複雑にしてしまう問題がある。これを回避する一つの方法は、log(ρq(l)n,n,k,ψ(l)(Xn,k))を解析的な扱いが比較的容易な関数で近似することである。以下では、これについて少し詳しく説明する。
【0078】
まず、x=(Xn,k(l)i,k)(σ(l) i,k)-1/2とおき、f(x)= log(ρq(l)n,n,k,ψ(l)(Xn,k))と等価な以下の正規表現に書き換えておく。
【0079】
【数21】

【0080】
すると、上記の関数f(x)の一次導関数df(x)/dxは、以下のどちらか一方の関数で近似できる。
【0081】
【数22】

【0082】
これらの関数は、どちらも解析的な扱いが容易であり、特にμ(l)i,k(l) i,kの更新に用いることで、例えば、以下の2つの利点を得ることができる。
【0083】
その1は、μ(l)i,kの更新に関してf(x)の一次導関数をg1(x)で近似することで、高速な収束が得られる繰り返し推定法の一つであるニュートン法に関しても大域的な収束が保証される。その2は、g2(x)は、log(ρq(l)n,n,k,ψ(l)(Xn,k))の一次導関数と同じ数学的な形式をしているため、μ(l)i,k(l) i,kの更新においてf(x)の一次導関数を、g2(x)で近似することで、複雑な計算を行わずに逐次最大化アルゴリズムを適用することができる。
【0084】
したがって、例えばg1(x)を用いると、ニュートン法によるμ(l)i,kの1回の更新は以下のように実現することができる。
【0085】
【数23】

【0086】
ここでν(l)iは、q(l)n=iを満たす時間nの集合を表す。g1′(x)は、g1(x)の一次導関数を表す。
【0087】
一方、g2(x)を用いると、σ(l) i,kの更新は以下のように実現することができる。
【0088】
【数24】

【0089】
ただし、κ(l) i,kは以下の値をとるものとする。
【0090】
【数25】

【0091】
また、音源占有度更新部70が行うE-stepは、音源占有度を式(32)に示すように更新する。
【0092】
【数26】

【0093】
ここで、事前に学習した隠れマルコフモデルがない場合の音源モデルパラメータ^ψ (l)の初期化について他の方法を説明する。音源パワー特徴量の時系列を混合ガウスモデルでモデル化し、その結果得られた混合ガウス分布の分布パラメータψ′のうち混合比α′からπ(l)(l)を、各ガウス分布の平均μ′i,kと分散σ′ i,kからμ(l)i,k,σ(l) i,kを定める。このとき、混合ガウスモデルは以下の形をとる。
【0094】
【数27】

【0095】
音源パワー特徴量Xn,kが与えられた条件下で、混合ガウス分布の分布パラメータを定めるためには、例えば、期待値最大化アルゴリズムなどのように、一般的に知られている方法を適用することができる。その結果得られた分布パラメータを元に、音源モデルパラメータ^ψ (l)を以下のように定めることができる。
【0096】
【数28】

【0097】
これにより、観測された音源パワー特徴量の分布を近似的に表現する音源モデルパラメータの初期化を行うことができる。
【0098】
〔変形例1〕
次に、音源モデルパラメータ更新部80が更新する音源モデルパラメータψ(l)の一つであるμ(l)i,kの更新方法の変形例を説明する。例えば、μ(l)i,kの事前確率密度関数p(μ(l)i,k)=N(μ(l)i,k;~μ (l)i,k,~σ (l) i,k)が与えられていると仮定する。この事前確率密度関数の分布パラメータとしては、上記して説明した例で述べた方法などに基づき初期化された音源モデルパラメータの値を用いて、~μ(l)i,k=^μ(l)i,k,~σ(l) i,k=^σ(l) i,kと定めることが効果的であることが実験により確認されている。
【0099】
この事前確率密度関数を用いると^μ(l)i,kの更新式(29)は、事後確率最大化基準に基づき以下のように修正される。
【0100】
【数29】

【0101】
ここで、ρは事前確率密度関数の重みを調整するコントロールパラメータ(>0)であり、事前確率密度関数を信頼する程度に基づいて自由に定めることができる。
【0102】
〔音源分離装置〕
図3に、この発明の音源分離装置200の機能構成例を示す。その動作フローを図4に示す。音源分離装置200は、上記した音源パラメータ推定装置100と、音源分離部95と、を具備する。
【0103】
音源分離部95は、音源パラメータ推定装置100が出力する音源パワー特徴量Xn,kと、更新した音源占有度^M (l)n,kと音源パワーパラメータ^q(l)nと音源モデルパラメータ^ψ(l)と、各音源信号の各時間周波数点における音源パワー特徴量のモデルβq(l),n,k,ψ(l)(S)と、を入力として複数の音源のそれぞれの音源分離信号^S(l)n,kを最小自乗誤差推定により求める。
【0104】
音源分離の方法は次式によって行う。
【0105】
【数30】

【0106】
〔確認実験〕
この発明の音源分離性能を評価する目的で確認実験を行った。実験条件を説明する。観測信号を30組用意し、全ての観測信号において音源数はNs=2とした。各観測信号は、それぞれ2人の男性の発話、2人の女性の発話、若しくは1名の女性と1名の男性の発話の混合音で構成した。
【0107】
標本化周波数は16kHzとした。各観測信号に含まれる2つのマイクロホン信号は、各話者の発話に関するマイク間時間差がそれぞれ±1.5ミリ秒になるように、計算機上で信号を加算して合成した(混合条件1)。音源モデルパラメータの初期値は、上記した事前に学習した隠れマルコフモデルがない場合の初期化方法で初期化した値を用いた。そして、変形例1で説明した事前確率密度関数を利用した。各隠れマルコフモデルの状態数は4とした。
【0108】
実験結果を図5と図6に示す。図5の縦軸は、混合前の音声から推定した音源パワー特徴量に関する隠れマルコフモデルのそれぞれについての各状態の出力分布の平均のパワー(dB)、横軸は周波数(kHz)である。図5(a)は混合前の音声から推定した結果、図5(b)は混合音から推定した結果を示す。
【0109】
混合音中から推定されたパラメータ(図5(b))は、混合前の音声から推定されたパラメータと酷似しており、この発明のパラメータ推定精度の信頼性の高さを証明している。
【0110】
図6に、分離後の信号のケプストラム歪みと混合音の数との関係を示す。縦軸はケプストラム歪み(dB)、横軸は混合音の数(数が増えるほど観測信号が長くなる)である。太い実線(□)で示す特性は、この発明で音源モデルパラメータを混合前の各音声信号から学習した場合を示す。実線(○)で示す特性は、この発明で音源モデルパラメータを観測信号のみから推定した場合を示す。一点鎖線(△)で示す特性は、非特許文献1の方法を用いた場合を示す。
【0111】
図6(a)は、混合条件1の混合音に対して処理した結果、図6(b)と(c)は、残響のある別の環境で収録した混合音(混合条件2と3)に対して処理した結果を示している。全ての場合において、この発明の音源分離方法は、非特許文献1の方法よりもケプストラム歪みの小さな音源分離を実現している。また、音源モデルパラメータを事前学習していない場合でも、事前学習している場合に相当する性能が実現できていることが確認できる。
【0112】
上記装置における処理手段をコンピュータによって実現する場合、各装置が有すべき機能の処理内容はプログラムによって記述される。そして、このプログラムをコンピュータで実行することにより、各装置における処理手段がコンピュータ上で実現される。
【0113】
この処理内容を記述したプログラムは、コンピュータで読み取り可能な記録媒体に記録しておくことができる。コンピュータで読み取り可能な記録媒体としては、例えば、磁気記録装置、光ディスク、光磁気記録媒体、半導体メモリ等どのようなものでもよい。具体的には、例えば、磁気記録装置として、ハードディスク装置、フレキシブルディスク、磁気テープ等を、光ディスクとして、DVD(Digital Versatile Disc)、DVD-RAM(Random Access Memory)、CD-ROM(Compact Disc Read Only Memory)、CD-R(Recordable)/RW(ReWritable)等を、光磁気記録媒体として、MO(Magneto Optical disc)等を、半導体メモリとしてEEP-ROM(Electronically Erasable and Programmable-Read Only Memory)等を用いることができる。
【0114】
また、このプログラムの流通は、例えば、そのプログラムを記録したDVD、CD-ROM等の可搬型記録媒体を販売、譲渡、貸与等することによって行う。さらに、このプログラムをサーバコンピュータの記録装置に格納しておき、ネットワークを介して、サーバコンピュータから他のコンピュータにそのプログラムを転送することにより、このプログラムを流通させる構成としてもよい。
【0115】
また、各手段は、コンピュータ上で所定のプログラムを実行させることにより構成することにしてもよいし、これらの処理内容の少なくとも一部をハードウェア的に実現することとしてもよい。

【特許請求の範囲】
【請求項1】
複数の音源信号それぞれの音源パワー時系列全体の状態を表す音源パワーパラメータの事前確率密度関数と、その音源パワーパラメータが与えられた場合の各音源信号の各時間周波数点における音源パワー特徴量のモデルと、を記憶した音源モデル記憶部と、
上記複数の音源信号を複数のマイクロホンで収音した時間領域信号を時間周波数領域信号に変換した観測信号を入力として、各時間周波数点における音源位置特徴量と音源パワー特徴量を抽出する特徴抽出部と、
上記音源パワー特徴量と、上記観測信号が得られた下での占有的な音源の事後確率密度関数である音源占有度と、上記音源パワーパラメータの事前確率密度関数と上記各音源信号の音源パワー特徴量のモデルと、上記音源パワー特徴量のモデルと上記事前確率密度関数のそれぞれの挙動を制御するパラメータである音源モデルパラメータと、を入力として上記各音源の音源パワーパラメータを更新する音源パワーパラメータ更新部と、
上記音源位置特徴量と音源占有度を入力として、上記各音源の音源位置パラメータを更新する音源位置パラメータ更新部と、
上記音源パワー特徴量と上記音源パワーパラメータと上記音源占有度と、上記音源モデル記憶部に記憶された音源パワーパラメータの事前確率密度関数と音源パワー特徴量のモデルとを入力として上記音源モデルパラメータを更新する音源モデルパラメータ更新部と、
上記音源位置特徴量と音源パワー特徴量と各音源の更新された音源パワーパラメータと音源位置パラメータと、上記音源モデルパラメータと、上記音源モデル記憶部に記憶された音源パワーパラメータの事前確率密度関数と音源パワー特徴量のモデルとを入力として上記各音源の音源占有度を更新する音源占有度更新部と、
を具備する音源パラメータ推定装置。
【請求項2】
請求項1に記載した音源パラメータ推定装置において、
上記音源パワー特徴量の時系列は隠れマルコフモデルに従い、上記各音源信号のマルコフモデルの状態が既知の条件下における上記音源パワーパラメータの事後確率密度関数は共分散行列が対角行列で表されるガウス分布でモデル化されており、当該ガウス分布の平均と共分散行列の対角要素を上記音源モデルパラメータに含むことを特徴とする音源パラメータ推定装置。
【請求項3】
請求項1又は2に記載した音源パラメータ推定装置と、
その音源パラメータ推定装置が出力する音源パワー特徴量と、当該音源パラメータ推定装置が更新した音源占有度と音源パワー特徴量と音源パワーパラメータと音源モデルパラメータと、上記各音源信号の各時間周波数点における音源パワー特徴量のモデルと、を入力として複数の音源のそれぞれの音源分離信号を、最小自乗誤差推定により求める音源分離部と、
を具備する音源分離装置。
【請求項4】
複数の音源信号を複数のマイクロホンで収音した時間領域信号を時間周波数領域信号に変換した観測信号を入力として、各時間周波数点における音源位置特徴量と音源パワー特徴量を抽出する特徴抽出過程と、
上記音源パワー特徴量と、上記観測信号が得られた下での占有的な音源の事後確率密度関数である音源占有度と、音源モデル記憶部に記憶された音源パワーパラメータの事前確率密度関数と上記各音源信号の音源パワー特徴量のモデルと、上記音源パワー特徴量のモデルと上記事前確率密度関数のそれぞれの挙動を制御するパラメータである音源モデルパラメータと、を入力として上記各音源の音源パワーパラメータを更新する音源パワーパラメータ更新過程と、
上記音源位置特徴量と音源占有度を入力として、上記各音源の音源位置パラメータを更新する音源位置パラメータ更新過程と、
上記音源パワー特徴量と上記音源パワーパラメータと上記音源占有度と、上記音源モデル記憶部に記憶された音源パワーパラメータの事前確率密度関数と音源パワー特徴量のモデルとを入力として上記音源モデルパラメータを更新する音源モデルパラメータ更新過程と、
上記音源位置特徴量と音源パワー特徴量と各音源の更新された音源パワーパラメータと音源位置パラメータと、上記音源モデルパラメータと、上記音源モデル記憶部に記憶された音源パワーパラメータの事前確率密度関数と音源パワー特徴量のモデルとを入力として上記各音源の音源占有度を更新する音源占有度更新過程と、
を備える音源パラメータ推定方法。
【請求項5】
請求項4に記載した音源パラメータ推定方法において、
上記音源パワー特徴量の時系列は隠れマルコフモデルに従い、上記各音源信号のマルコフモデルの状態が既知の条件下における上記音源パワーパラメータの事後確率密度関数は共分散行列が対角行列で表されるガウス分布でモデル化されており、当該ガウス分布の平均と共分散行列の対角要素を上記音源モデルパラメータに含むことを特徴とする音源パラメータ推定方法。
【請求項6】
請求項4又は5に記載した音源パラメータ推定方法で抽出した音源パワー特徴量と、当該音源パラメータ推定方法で更新した音源占有度と音源パワー特徴量と音源パワーパラメータと音源モデルパラメータと、上記各音源信号の各時間周波数点における音源パワー特徴量のモデルと、を入力として複数の音源のそれぞれの音源分離信号を、最小自乗誤差推定により求める音源分離過程と、
を備える音源分離方法。
【請求項7】
請求項1乃至3の何れかに記載した音源パラメータ推定装置又は音源分離装置としてコンピュータを機能させるためのプログラム。

【図1】
image rotate

【図2】
image rotate

【図3】
image rotate

【図4】
image rotate

【図5】
image rotate

【図6】
image rotate

【図7】
image rotate


【公開番号】特開2012−173592(P2012−173592A)
【公開日】平成24年9月10日(2012.9.10)
【国際特許分類】
【出願番号】特願2011−36713(P2011−36713)
【出願日】平成23年2月23日(2011.2.23)
【出願人】(000004226)日本電信電話株式会社 (13,992)
【Fターム(参考)】