(57)【特許請求の範囲】
【請求項１】
音データを検出する音検出部と、
検出された音データに基づいて、方向及び時間ごとのＭＵＳＩＣスペクトルを演算するＭＵＳＩＣスペクトル演算部と、
各方向の音源の有無による状態、状態の遷移を記述する状態遷移モデル、及び音源有りの状態及び音源無しの状態に観測されるＭＵＳＩＣスペクトルを記述する観測モデルを定め、ＭＵＳＩＣスペクトルの時系列データに基づいて、前記観測モデル及び前記状態遷移モデルのモデルパラメータの事後分布を推定する、モデルパラメータ推定部と、
推定されたモデルパラメータの事後分布に基づいて、方向及び時間ごとの音源存在事後確率のパーティクルをサンプルすることによって、音源定位を行う音源定位部と、を備えた音源定位システム。
【請求項２】
前記観測モデルとしてガウス混合モデルを使用する請求項１に記載の音源定位システム。
【請求項３】
音データを検出するステップと、
検出された音データに基づいて、方向及び時間ごとのＭＵＳＩＣスペクトルを演算するステップと、
各方向の音源の有無による状態、状態の遷移を記述する状態遷移モデル、及び音源有りの状態及び音源無しの状態に観測されるＭＵＳＩＣスペクトルを記述する観測モデルを定め、ＭＵＳＩＣスペクトルの時系列データに基づいて、前記観測モデル及び前記状態遷移モデルのモデルパラメータの事後分布を推定するステップと、
推定されたモデルパラメータの事後分布に基づいて、方向及び時間ごとの音源存在事後確率のパーティクルをサンプルすることによって、音源定位を行うステップと、を含む音源定位方法。
【請求項４】
前記観測モデルとしてガウス混合モデルを使用する請求項３に記載の音源定位方法。
【請求項５】
音源定位を行うステップは、Ｐ個のパーティクルをサンプルするステップと、各パーティクルについて重みを計算するステップと、各パーティクルの重みを正規化するステップと、各パーティクルの重みを使用してリサンプリングを行うステップと、を含む請求項３または４に記載の音源定位方法。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、音源定位システム及び音源定位方法に関する。
【背景技術】
【０００２】
音源定位は、マイクロフォンアレイを使用した同時発話混合音声の分離、遠隔ロボットのオペレータへの音源方向提示、移動ロボットによる音源検出と位置推定などに使用される重要な要素技術である。
【０００３】
マイクロフォンアレイを使用した音源定位法には、ビームフォーミングに基づく方法とMUltiple SIgnal Classification（ＭＵＳＩＣ）に基づく方法がある。これらのうち、ＭＵＳＩＣに基づくＭＵＳＩＣ法は、雑音に対してロバストであり、音源数がマイクロフォンの数未満という条件下では、比較的安定して複数音源の定位を行うことができる(たとえば、特許文献１)。
【０００４】
通常のＭＵＳＩＣ法では、ＭＵＳＩＣスペクトルと呼ばれる音源到来評価関数に対して閾値を設定して、音源がある方向に存在するかどうか判定する。適切な閾値を定めるには、環境中の音源数や残響時間などを考慮する必要がある。したがって、動的に音環境が変化する状況において音源定位を行うには、従来技術においては、閾値を手動で設定する必要があった。換言すれば、動的に音環境が変化する状況において、ＭＵＳＩＣスペクトルの閾値を自動的に定めることのできる音源定位システム及び音源定位方法は開発されていなかった。
【先行技術文献】
【特許文献】
【０００５】
【特許文献１】特許４０９５３４８号
【発明の概要】
【発明が解決しようとする課題】
【０００６】
したがって、動的に音環境が変化する状況において、ＭＵＳＩＣスペクトルの閾値を自動的に定めることのできる音源定位システム及び音源定位方法に対するニーズがある。
【課題を解決するための手段】
【０００７】
本発明の第１の態様による音源定位システムは、音データを検出する音検出部と、検出された音データに基づいて、方向及び時間ごとのＭＵＳＩＣスペクトルを演算するＭＵＳＩＣスペクトル演算部と、を備えている。本システムは、各方向の音源の有無による状態、状態の遷移を記述する状態遷移モデル、及び音源有りの状態及び音源無しの状態に観測されるＭＵＳＩＣスペクトルを記述する観測モデルを定め、ＭＵＳＩＣスペクトルの時系列データに基づいて、前記観測モデル及び前記状態遷移モデルのモデルパラメータの事後分布を推定する、モデルパラメータ推定部と、推定されたモデルパラメータの事後分布に基づいて、方向及び時間ごとの音源存在事後確率のパーティクルをサンプルすることによって、音源定位を行う音源定位部と、をさらに備えている。
【０００８】
本態様の音源定位システムによれば、観測モデル及び状態遷移モデルのモデルパラメータの事後分布を推定し、推定されたモデルパラメータの事後分布に基づいて音源定位を行うので、動的に音環境が変化する状況においても、手動で閾値を設定する必要なく、ロバストに音源定位を行うことができる。
【０００９】
本発明の第１の態様の第１の実施形態による音源定位システムは、前記観測モデルとしてガウス混合モデルを使用する。
【００１０】
本実施形態によれば、ガウス分布を使用することで解析的な計算が可能となる。
【００１１】
本発明の第２の態様による音源定位システムは、音データを検出するステップと、検出された音データに基づいて、方向及び時間ごとのＭＵＳＩＣスペクトルを演算するステップと、を含む。本方法は、各方向の音源の有無による状態、状態の遷移を記述する状態遷移モデル、及び音源有りの状態及び音源無しの状態に観測されるＭＵＳＩＣスペクトルを記述する観測モデルを定め、ＭＵＳＩＣスペクトルの時系列データに基づいて、前記観測モデル及び前記状態遷移モデルのモデルパラメータの事後分布を推定するステップと、  推定されたモデルパラメータの事後分布に基づいて、方向及び時間ごとの音源存在事後確率のパーティクルをサンプルすることによって、音源定位を行うステップと、をさらに含む。
【００１２】
本態様の音源定位方法によれば、観測モデル及び状態遷移モデルのモデルパラメータの事後分布を推定し、推定されたモデルパラメータの事後分布に基づいて音源定位を行うので、動的に音環境が変化する状況においても、手動で閾値を設定する必要なく、ロバストに音源定位を行うことができる。
【００１３】
本発明の第２の態様の第１の実施形態による音源定位方法は、前記観測モデルとしてガウス混合モデルを使用する。
【００１４】
本実施形態によれば、ガウス分布を使用することで解析的な計算が可能となる。
【００１５】
本発明の第２の態様の第２の実施形態による音源定位方法において、音源定位を行うステップは、Ｐ個のパーティクルをサンプルするステップと、各パーティクルについて重みを計算するステップと、各パーティクルの重みを正規化するステップと、各パーティクルの重みを使用してリサンプリングを行うステップと、を含む。
【００１６】
本実施形態によれば、推定されたモデルパラメータの分布に基づいてパーティクルをサンプルすることによって、方向及び時間ごとの音源事後確率のパーティクルを簡単な手順で求めることができる。
【図面の簡単な説明】
【００１７】
【図１】本発明の一実施形態による音源定位システムの構成を示す図である。
【図２】Ｍ個のマイクロフォンからなるマイクロフォンアレイの構成を示す図である。
【図３】対数スケールのＭＵＳＩＣスペクトルの分布を示す図である。
【図４】ＶＢ−ＨＭＭの確率変数間の条件付き独立性を示すグラフィカルモデルを示す図である。
【図５】モデルパラメータ推定部による、モデルパラメータの分布を推定する手順を示す流れ図である。
【図６】音源定位部が、各方向ビンにおける音源存在事後確率を表す、P 個のパーティクルを求める手順を示す流れ図である。
【図７】オンラインの音源定位実験に使用した音源の配置を示す図である。
【図８】従来の音源定位システムによるオンライン音源定位の結果を示す図である。
【図９】本実施形態の音源定位システムによるオンライン音源定位の結果を示す図である。
【発明を実施するための形態】
【００１８】
図１は、本発明の一実施形態による音源定位システム１００の構成を示す図である。音源定位システム１００は、音検出部１０１と、ＭＵＳＩＣスペクトル演算部１０３と、モデルパラメータ推定部１０５と、音源定位部１０７とを含む。
【００１９】
音検出部１０１は、Ｍ個のマイクロフォンからなるマイクロフォンアレイであってもよい。
【００２０】
図２は、Ｍ個のマイクロフォン１０１１からなるマイクロフォンアレイ１０１の構成を示す図である。図２において、Ｍ＝８である。一例として、８個のマイクロフォン１０１１は水平面上に配置されており、音源定位システム１００は、音源が、水平面上のどの方向に存在するかを定める。一例として、方向の解像度は５°であり、３６０°／５°＝７２のいずれかの方向を定める。
【００２１】
たとえば、音検出部であるマイクロフォンアレイによって、Ｍチャネルの音響信号が得られる。また、周波数ビンごとに水平面上のＤ方向（Ｄ＝７２）からの伝達関数が与えられているとする。音源定位システム１００は、Ｎ個の音源方向を定める。ただし、同時に定位可能な最大音源数Ｎmaxは、マイクロフォンの数未満である。
Ｎ≦Ｎmax＜Ｍ
【００２２】
ＭＵＳＩＣスペクトル演算部１０３におけるＭＵＳＩＣ(Multiple Signal Classification)スペクトルの算出方法について説明する。より詳細は、文献[R. O. Schmidt, “Multiple Emitter Location and Signal Parameter Estimation,” IEEE Trans. on Antennas and Propagation, vol. 34, no. 3, pp. 276-280, 1986; P. Dan`es and J. Bonnal, “Information-Theoretic Detection of Broadband Sources in a Coherent Beamspace MUSIC Scheme,” in Proc. of IROS-2010, 2011, pp. 1976-1981.]に記載されている。ＭＵＳＩＣ法は、時間周波数領域において適用される。具体的には、サンプリング周波数１６０００[Ｈｚ]で、窓長１２[ｐｔ]、シフト幅１６０[ｐｔ]の短時間フーリエ変換を行っている。
【数１】
0005629249.tif 000002
をＭチャネル音響信号の時間フレームτ、周波数ビンωにおける入力Ｍチャネル音響信号の複素振幅ベクトルとする。各周波数ビンω、ΔＴ［ｓｅｃ］間隔の時刻ｔに対して、
(１)入力信号の自己相関行列
【数２】
0005629249.tif 000003
の計算、
【数３】
0005629249.tif 000004
の固有値分解、
(３)固有ベクトルと伝達関数を用いたＭＵＳＩＣスペクトルの計算を行う。
【００２３】
上記の(１)乃至(３)について以下に説明する。
（１）入力信号の自己相関行列の計算
入力信号の自己相関行列は、時間ΔＴで観測したサンプル値の相関として計算する。
【数４】
0005629249.tif 000005
ただし、
【数５】
0005629249.tif 000006
はエルミート転置、
【数６】
0005629249.tif 000007
は時刻ｔにおける時間フレームを表す。入力ベクトル
【数７】
0005629249.tif 000008
のＭ個の要素は各チャネルに対応する。
（２）固有値分解
【数８】
0005629249.tif 000009
を次のように固有値分解する。
【数９】
0005629249.tif 000010
ここで、
【数１０】
0005629249.tif 000011
は固有ベクトル、
【数１１】
0005629249.tif 000012
は固有値からなる対角行列である。
【数１２】
0005629249.tif 000013
は、
【数１３】
0005629249.tif 000014
と
【数１４】
0005629249.tif 000015
のＭ個の固有ベクトルで表せ、
【数１５】
0005629249.tif 000016
となる。固有値
【数１６】
0005629249.tif 000017
は降順に並べられているものとする。
【００２４】
入力信号にＮ個の音源が含まれる場合、固有値
【数１７】
0005629249.tif 000018
から
【数１８】
0005629249.tif 000019
までが、音源のエネルギに対応する大きな値を有する。それに対し、残りの固有値
【数１９】
0005629249.tif 000020
から
【数２０】
0005629249.tif 000021
まではマイクロフォンの観測ノイズなどに対応する小さな値を有する。ここで、重要な点は、
【数２１】
0005629249.tif 000022
から
【数２２】
0005629249.tif 000023
のノイズに対応する固有ベクトルは、音源方向に対応する伝達関数ベクトルと直交するという点である[R. O. Schmidt, “Multiple Emitter Location and Signal Parameter Estimation,” IEEE Trans. on Antennas and Propagation, vol. 34, no. 3, pp. 276-280, 1986.]。
(３)固有ベクトルと伝達関数を用いたＭＵＳＩＣスペクトルの計算
ＭＵＳＩＣスペクトルは以下の式にしたがって計算する。
【数２３】
0005629249.tif 000024
ただし、
【数２４】
0005629249.tif 000025
は方向ｄ、周波数ビンωに対応するＭ次元の伝達関数ベクトルである。これらの伝達関数は、マイクロフォンアレイを使用して事前に測定したものである。ここで、観測されうる最大の音源数はＮmaxである。そのため、
【数２５】
0005629249.tif 000026
から
【数２６】
0005629249.tif 000027
までの固有ベクトルは、音源方向ｄに対応する伝達関数
【数２７】
0005629249.tif 000028
と直交する。したがって、式(３)の分母は音源方向のｄに対して０となる。すなわち、式(3)のＭＵＳＩＣスペクトル
【数２８】
0005629249.tif 000029
は発散する。ただし、実際には、壁からの反射音などの影響でＭＵＳＩＣスペクトルは発散せず鋭いピークとして観測されることが多い。
【００２５】
つぎに、周波数ビンごとのＭＵＳＩＣスペクトルを以下の式にしたがって合算する。
【数２９】
0005629249.tif 000030
ここで、
【数３０】
0005629249.tif 000031
は、周波数ビンωにおける最大固有値である。本実施形態では、音声信号を対象とするため、
【数３１】
0005629249.tif 000032
とした。
【００２６】
つぎに、モデルパラメータ推定部１０５の機能について説明する。モデルパラメータ推定部１０５は、変分ベイズ隠れマルコフモデル(ＶＢ−ＨＭＭ)を使用する。
【００２７】
状態ベクトルとしてＤ次元の２値ベクトルを使用する。各次元のベクトル値は、その方向に音源が存在するか否かを示す。
【００２８】
また、ＭＵＳＩＣスペクトルがガウス分布にしたがう観測値であるとみなし、音源有りの場合及び音源無しの場合に対応するガウス分布からなるガウス混合分布を観測モデルとする。ガウス分布を使用する理由は、複数の周波数ビンの値を加算して対数とったＭＵＳＩＣスペクトルが近似的にガウス分布とみなせることと、ガウス分布を使用することで解析的な計算が可能となることである。
【００２９】
図３は、対数スケールのＭＵＳＩＣスペクトルの分布を示す図である。図３の横軸は、対数スケールのＭＵＳＩＣスペクトルを示す。対数スケールのＭＵＳＩＣスペクトルは、以下の式によって求める。
【数３２】
0005629249.tif 000033
図３の縦軸は、観測数を表す。図３に点線で示す音源が存在しない（オフ）のときのガウス分布は、狭くかつ値の小さなＭＵＳＩＣスペクトルの領域に形成され、図３に実線で示す音源が存在する（オン）ときのガウス分布は、広くかつ値の大きなＭＵＳＩＣスペクトルの領域を覆っている。
【００３０】
モデルパラメータ推定部１０５で使用される観測モデルは以下の式で表せる。
【数３３】
0005629249.tif 000034
ただし、
【数３４】
0005629249.tif 000035
はｘ＝ｙのとき
【数３５】
0005629249.tif 000036
そうでなければ、
【数３６】
0005629249.tif 000037
を表す。また、
【数３７】
0005629249.tif 000038
は、平均μ、精度λの正規分布の確率密度関数を表す。パラメータμとλには、正規・ガンマ分布を使用する。
【数３８】
0005629249.tif 000039
ただし，N(・|m,L-1)は平均 m, 精度L (分散 1/L)の正規分布であり、以下の式で表せる。
【数３９】
0005629249.tif 000040
G(・|a,b)は形状a, 尺度b のガンマ分布であり、以下の式で表せる。
【数４０】
0005629249.tif 000041
正規分布のβ、 ガンマ分布のa は事前分布の影響の大きさ(事前に観測したとみなす観測数)を表す。本実施形態では、事前情報よりも学習時に観測されるデータを重視するため、
【数４１】
0005629249.tif 000042
に設定する。
【００３１】
m0 は平均パラメータμの事前情報で与えられた平均値を示す。本実施形態では、およそ25程度が適している。あるいは、ＶＢ−ＨＭＭの学習に用いる観測値のサンプル平均でも問題ない。
【００３２】
b0 は精度パラメータλの事前情報で与えられた「ばらけ具合」を示す。実験的に500に設定している。あるいは、ＶＢ−ＨＭＭの学習に用いる観測値のサンプル分散程度でも問題ない。
【００３３】
図４は、ＶＢ−ＨＭＭの確率変数間の条件付き独立性を示すグラフィカルモデルを示す図である。ＶＢ−ＨＭＭは、状態遷移確率のパラメータθkや、観測確率のパラメータμ、λが数値ではなく、確率変数として扱われる点が通常のＨＭＭと異なる。モデルパラメータ推定部１０５は、これらのパラメータの確率分布を学習する。
【００３４】
モデルパラメータ推定部１０５で使用される状態遷移モデルは基本的に、各方向ビンd について、 前状態で音源がない場合
【数４２】
0005629249.tif 000043
と音源がある場合
【数４３】
0005629249.tif 000044
から、次状態で音源が出現する、継続する、消滅するといった遷移を考える。本実施形態においては、移動する音源についても考慮するために、以下の表1に示すように前状態の組み合わせから成る四つの場合を考える。すなわち、前時刻の同方向ビン
【数４４】
0005629249.tif 000045
に音源が存在するかどうかと, 前時刻の隣接方向ビン
【数４５】
0005629249.tif 000046
のいずれかに音源が存在するかによって分類する。たとえば、θ1 は前時刻に当該方向d及び隣接ビンd±1 に音源が存在しない状態から音源が出現する確率、θ2 は、前時刻に方向d に音源が存在しないが、隣接ビンd±1 には音源が存在したため、その音源が方向d に移動してきて
【数４６】
0005629249.tif 000047
となる確率を表す。
【表１】
0005629249.tif 000048
【００３５】
状態遷移確率は、以下の式によって表せる。
【数４７】
0005629249.tif 000049
ここで、
【数４８】
0005629249.tif 000050
は表1 に従って、 方向ビンd の周りの前状態の値
【数４９】
0005629249.tif 000051
によって条件k に合致するときに
【数５０】
0005629249.tif 000052
その他の場合は0 を返す条件識別関数である。 初期状態としては、 音源は存在しない、すなわちすべてのd に対して
【数５１】
0005629249.tif 000053
とする。
【００３６】
状態遷移パラメータである
【数５２】
0005629249.tif 000054
には、式(8)の共役事前分布としてベータ分布を用いる。
【数５３】
0005629249.tif 000055
ただし
【数５４】
0005629249.tif 000056
はパラメータｃ、ｄを持つβ分布の確率密度関数である。
【００３７】
モデルパラメータ推定部１０５におけるＶＢ−ＨＭＭの学習は、事後分布
【数５５】
0005629249.tif 000057
を以下のように因数分解可能な分布に近似して推定する。
【数５６】
0005629249.tif 000058
(・)1:T は、 時刻1 からT までの確率変数の集合を表す。一般的なＶＢ−ＨＭＭ の推論については、文献[M. J. Beal, “Variational Algorithms for Approximate Bayesian Inference,”Ph.D. dissertation, Gatsby Computational Neuroscience Unit, University Colledge London, 2003.] に記載されている。
【数５７】
0005629249.tif 000059
はそれぞれのk に対し、式(11)に示すパラメータ
【数５８】
0005629249.tif 000060
を持つベータ分布となり、
【数５９】
0005629249.tif 000061
は、 式(12), (13) で表されるパラメータ
【数６０】
0005629249.tif 000062
を持つ正規ガウス分布となる。
【数６１】
0005629249.tif 000063
ただし、 変数st,d,j は, st,d = 0 のとき、st,d,0 = 1 、また、st,d =1 のとき, st,d,1 =1 となる変数である。 式(12)、(13)に用いられる正規分布の十分統計量は
【数６２】
0005629249.tif 000064
と定義する。 また、
【数６３】
0005629249.tif 000065
は式(10) の分布による期待値演算子である。各時刻の状態変数と状態遷移の期待値
【数６４】
0005629249.tif 000066
は次のように計算する。
【数６５】
0005629249.tif 000067
ただし、α(st,d,j)とβ(st,d,j)はそれぞれ前向き、後ろ向き再帰式により計算される。
【数６６】
0005629249.tif 000068
ここで、
【数６７】
0005629249.tif 000069
である。遷移、観測確率の幾何平均は以下の式で表せる。
【数６８】
0005629249.tif 000070
ここで、
【数６９】
0005629249.tif 000071
は、以下の式で定義されるディガンマ関数である。
【数７０】
0005629249.tif 000072
【００３８】
式(14)及び(15)は、ともに添え字ｊ、ｋを動かしたとき総和が1になるように正規化されている。
【数７１】
0005629249.tif 000073
は、状態遷移の条件ｋに関する前向き確率である。
【００３９】
図５は、モデルパラメータ推定部１０５による、モデルパラメータの分布を推定する手順を示す流れ図である。
【００４０】
図５のステップＳ１０１０において、モデルパラメータ推定部１０５は、初期値を定める。初期値は、式(14)及び式(15)の値に対して、たとえば、以下手順によって定めることができる。
【００４１】
式(14) の左辺 <st,d,j> は、時刻t、方向ビンd において音源が存在しない場合は st,d,0= 1かつst,d,1=0 、音源が存在する場合は st,d,0=0 かつ st,d,1=1 となる2値変数の期待値である。観測値xt,dが適当な閾値 (たとえば、m0の値)を超える場合は、たとえば、
<st,d,1> = 0.8, <st,d,0> = 1-0.8=0.2 
のように設定する。0.8の代わりに1にしてもほぼ同様に動作する。
【００４２】
式(15) の左辺 <st,d,jfk(st,d)> も、xt,dが閾値を超えるか否かに従って算出する。この値は、st,d,1=0,1 の2通りと、k=1〜4のいずれかでfk(st,d)=1となる4通りの組み合わせで、計8通りの組み合わせが存在する。xt,dの閾値処理と、前時刻のxt-1,dの閾値処理、前時刻・近傍方向のxt-1,d±1の閾値処理結果から、表1を参照して、該当する k の値を求める。たとえば、前時刻xt-1,dは閾値以下だが、隣接ビンのxt-1,d+1では閾値を超えた場合は k = 2 にあたる。xt,dが閾値を超えた場合、<st,d,1f2(st,d)> = 0.8、その他の7通りの組み合わせでは、<st,d,jfk(st,d)>=(1-0.8)/7 のように設定する。
【００４３】
図５のステップＳ１０２０において、モデルパラメータ推定部１０５は、式(18)及び式(19)を使用して、遷移、観測確率の幾何平均を求める。
【００４４】
図５のステップＳ１０３０において、モデルパラメータ推定部１０５は、ステップＳ１０２０で求めた遷移、観測確率の幾何平均、式(16)及び式(17)を使用して、α(st,d,j)とβ(st,d,j)を計算する。
【００４５】
図５のステップＳ１０４０において、モデルパラメータ推定部１０５は、ステップＳ１０３０で求めたα(st,d,j)とβ(st,d,j)、式(14)及び式(15)を使用して、各時刻の状態変数と状態遷移の期待値を求める。
【００４６】
図５のステップＳ１０５０において、モデルパラメータ推定部１０５は、ステップＳ１０４０で求めた各時刻の状態変数と状態遷移の期待値、式(11)から式(13)を使用して、モデルパラメータの事後分布を推定する。
【００４７】
図５のステップＳ１０６０において、モデルパラメータ推定部１０５は、収束判定を行う。具体的に、式(12), (13) で算出されるパラメータβ, m, a, b の値の変動がなくなることによって収束判定を行う。収束していないと判断されれば、ステップＳ１０２０に戻る。収束していると判断されれば、処理を終了する。
【００４８】
つぎに、音源定位部１０７の機能について説明する。音源定位部１０７は、モデルパラメータ推定部１０５が推定した、モデルパラメータの事後分布に基づいて、パーティクルフィルタを使用して、複数音源の存在事後確率計算を行う。パーティクルフィルタの推定対象は, ＭＵＳＩＣスペクトルの時系列データが与えられたときの、各方向ビンにおける音源存在事後確率である。この分布を、P 個のパーティクルを用いて以下のように近似計算する。
【数７２】
0005629249.tif 000074
ただし,
【数７３】
0005629249.tif 000075
はパーティクルp の重み、
【数７４】
0005629249.tif 000076
は状態ベクトルの値である。
【００４９】
図６は、音源定位部１０７が、各方向ビンにおける音源存在事後確率を表す、Ｐ個のパーティクルを求める手順を示す流れ図である。
【００５０】
図６のステップＳ２０１０において、音源定位部１０７は、Ｐ個のパーティクルをサンプルする。
【００５１】
Ｐは以下のように定める。Ｐが大きいほど式(20)の近似は正確になるが、Ｐの大きさに比例した計算時間がかかる。そこで、Ｐを決める一般的な手続きとしては、実用的な近似を得るほど十分大きなＰを与え、このようなＰの処理時間が大きすぎる場合は、Ｐの大きさを適宜縮小する。本実施形態においては、近似結果が収束し、理速度も十分高速な値として、Ｐ＝５００とした。
【００５２】
Ｐ個のパーティクルのサンプルは、以下の式で表せる分布を使用して行う。
【数７５】
0005629249.tif 000077
ただし、xt,dが極大値をとるｄのとき、C(xt,d)=1でその他の場合は、C(xt,d)=0となる。上記の分布の重みには、以下の式で表せるマハラノビス距離を使用する。
【数７６】
0005629249.tif 000078
【００５３】
時刻t において，式(22)で計算される分布q は、計Ｄ個のビンに対してオン
【数７７】
0005629249.tif 000079
であるか、またはオフ
【数７８】
0005629249.tif 000080
であるかの確率を与える。サンプル方法としては、各dについて
a)     C(xt,d)=0 の場合、j=0 とする。つまり
【数７９】
0005629249.tif 000081
b)     C(xt,d)=1の場合、j=0,1 それぞれの場合の分布qの確率を参照する。たとえば、
【数８０】
0005629249.tif 000082
の場合、0〜1の区間から一様乱数を生成し、その値が0.8以下の場合は、
【数８１】
0005629249.tif 000083
、0.8を超えた場合は
【数８２】
0005629249.tif 000084
のように値を決定する。
【００５４】
図６のステップＳ２０２０において、音源定位部１０７は、以下の式にしたがって、各パーティクルについて重みwp を算出する。
【数８３】
0005629249.tif 000085
式(24)及び式(25)の状態遷移、観測確率は、モデルパラメータ推定部１０５で使用された式(6)及び式(8)の事後分布で積分消去することで計算することができる。分布の共役性を用いると, この積分計算は、次のように解析的に求まる。
【数８４】
0005629249.tif 000086
ただし、
【数８５】
0005629249.tif 000087
は、平均m, 精度l , 自由度n のStudent t-分布である。さらに、最大の音源数をNmax に抑えるため、状態ベクトル
【数８６】
0005629249.tif 000088
に存在する音源数がNmax を超える場合には観測確率は
【数８７】
0005629249.tif 000089
とする。
【００５５】
図６のステップＳ２０３０において、音源定位部１０７は、各パーティクルの重みwp が
【数８８】
0005629249.tif 000090
となるよう正規化する。
【００５６】
図６のステップＳ２０４０において、処理を終了するかどうか判断する。たとえば、スイッチの状態によって処理を終了するかどうか判断してもよい。処理を終了しないと判断した場合には、ステップＳ２０５０に進む。そうでなければ処理を終了する。
【００５７】
図６のステップＳ２０５０において、音源定位部１０７は、リサンプリングを行う。リサンプリングは、パーティクルp が持つ値
【数８９】
0005629249.tif 000091
を、そのパーティクルの重み
【数９０】
0005629249.tif 000092
に比例した確率で複製することによって行う。方法の一例としては、以下の処理を p’=1 〜 P まで繰り返す。
a)     0〜1 区間から一様乱数を生成 rp’とする。
b)     p=1〜P まで
ｉ．    rp’← rp’-wp
ｉｉ．  rp’< 0 になったら、
【数９１】
0005629249.tif 000093
として，p のループを抜ける。
ｉｉｉ．wp’← 1/P (リサンプル後の重みは全てのパーティクルで等しい)
c)     a) に戻る。
【００５８】
つぎに、評価実験について説明する。評価実験は、本実施形態による音源定位システム及び固定閾値を使用する従来の音源定位システムの性能を比較した。モデルパラメータ推定部１０５による、オフラインでのＶＢ−ＨＭＭでの学習は、一人の話者がマイクロフォンの周囲を発話しながら動く音響信号で行った。
【００５９】
図７は、オンラインの音源定位実験に使用した音源の配置を示す図である。二人の話者３０１及び３０３が、マイクロフォンアレイ１０１の周囲を移動しながら発話し、固定されたスピーカ２０１から音楽が再生されている。 オフライン及びオンラインで用いられた信号の長さはともに20 [sec] である。
【００６０】
パラメータの設定は以下のとおりである。
Nmax = 3, α0 = [1,1], β0 =1, a0 = 1, b0 = 500 
また、パーティクル数はＰ＝５００ とした。実験で使用した室内の残響時間はRT20 = 840 [msec] であった。
【００６１】
図８は、従来の音源定位システムによるオンライン音源定位の結果を示す図である。図の横軸は秒単位の時間を示し、図の縦軸は度単位の方向を示す。従来の音源定位システムの閾値はPthres = 23,25,27 に設定されている。図８(ａ)、(ｂ)、(ｃ)は、それぞれ、閾値が23,25,27の場合の音源定位の結果を示す。図８(ａ)、(ｂ)、(ｃ)においては、閾値を超えるビンを、音源が存在するとして黒く表示している。図８(ａ)、(ｂ)、(ｃ)において、固定されたスピーカ及び移動する話者が黒で表示されている。しかし、図８(ａ)、(ｂ)の実線の囲みで示すように、閾値を低く設定した場合は音源の誤検出が頻発する。
【００６２】
図９は、本実施形態の音源定位システムによるオンライン音源定位の結果を示す図である。図の横軸は秒単位の時間を示し、図の縦軸は度単位の方向を示す。本実施形態の音源定位システムの初期値はm0 = 23,25,27 に設定されている。図９(ａ)、(ｂ)、(ｃ)は、それぞれ、初期値が23,25,27の場合の音源定位の結果を示す。図９(ａ)、(ｂ)、(ｃ)においては、事後分布の音源存在確率が0.95 以上のビンを、音源が存在するとして黒く表示している。図９(ａ)、(ｂ)、(ｃ)において、固定されたスピーカ及び移動する話者が黒で表示されている。図８(ａ)、(ｂ)の実線の囲みの位置に対応する、図９(ａ)、(ｂ)の点線の囲みの位置において音源の誤検出は少ない。このことは、本実施形態の音源定位システムにおいて、学習の初期値にかかわらず音源の誤検出が少ないことを意味する。また、本実施形態の音源定位システムにおいて音源存在確率の閾値を0.95から1.00 まで動かして結果を検証し、これらの閾値に対してもロバストに同様の結果を示すことを確認した。この結果から、本実施形態の音源定位システムのモデルパラメータ推定部１０５によるオフライン学習及び音源定位部１０７によるオンライン定位の枠組みが、自動的に音源定位に適したパラメータに収束することが確認できた。さらに、今回の実験条件から、本実施形態の音源定位方法は、１音源のみで学習した場合であっても、複数音源に対して安定した音源定位結果をもたらすことが実証された。
【符号の説明】
【００６３】
１０1…音検出部、１０３…ＭＵＳＩＣスペクトル演算部、１０５…モデルパラメータ推定部、１０７…音源定位部
【図１】
0005629249.tif 000094
【図２】
0005629249.tif 000095
【図３】
0005629249.tif 000096
【図４】
0005629249.tif 000097
【図５】
0005629249.tif 000098
【図６】
0005629249.tif 000099
【図７】
0005629249.tif 000100
【図８】
0005629249.tif 000101
【図９】
0005629249.tif 000102