(57)【要約】【課題】音響信号により適した距離尺度に基づく、拡散性雑音の存在下での雑音抑圧時及び音声強調時等に用いるモデルパラメータの推定技術を提供する。【解決手段】観測信号を、非定常の拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、モデル化するときのモデルパラメータは拡散性雑音の共分散行列を特定するモデルパラメータと、ノイズフリー信号の共分散行列を特定するモデルパラメータとを含み、モデル推定装置は、観測信号から得られる共分散行列と、拡散性雑音の共分散行列とノイズフリー信号の共分散行列との和との板倉-齋藤距離が小さくなるように、モデルパラメータを更新するモデル推定部を含む。【選択図】図１
2015135437.tif 000002
【特許請求の範囲】
【請求項１】
観測信号を、非定常の拡散性雑音と前記拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、モデル化するときのモデルパラメータは拡散性雑音の共分散行列を特定するモデルパラメータと、ノイズフリー信号の共分散行列を特定するモデルパラメータとを含み、
観測信号から得られる共分散行列と、拡散性雑音の共分散行列とノイズフリー信号の共分散行列との和との板倉-齋藤距離が小さくなるように、前記モデルパラメータを更新するモデル推定部を含む、
モデル推定装置。
【請求項２】
観測信号を、非定常の拡散性雑音と前記拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、１つ以上のフレームからなるブロックの番号をt∈{1,…,T}とし、ブロック内のフレームの番号をu∈{1,…,U}とし、ブロックtごとの観測信号ytuの共分散行列をΦytとし、ブロックtごとのノイズフリー信号xtuの共分散行列をΦxtとし、ブロックtごとの拡散性雑音vtuの共分散行列をΦvtとし、
観測信号ytから得られる共分散行列Φytと、Φxt＋Φvtとの板倉-齋藤距離の、すべてのブロックtについての総和が小さくなるように、前記共分散行列Φxtと前記共分散行列Φvtとを更新するモデル推定部を含み、
前記ノイズフリー信号は音源信号である、
モデル推定装置。
【請求項３】
フレームごとの観測信号を非定常の拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、１つ以上のフレームからなる時間区間をブロックとして、
複数ブロック分の観測信号から得られる共分散行列と、複数ブロック分の拡散性雑音の共分散行列と複数ブロック分のノイズフリー信号の共分散行列との和と、の板倉-齋藤距離が小さくなるように、前記複数ブロック分の拡散性雑音の共分散行列を特定するパラメータと前記複数ブロック分のノイズフリー信号の共分散行列を特定するパラメータとを更新するモデル推定部を含み、
前記ノイズフリー信号は、音源から発せられる直接音と初期反射成分とからなる信号であるクリーン信号に残響が重畳された信号であり、
前記複数ブロック分のノイズフリー信号の共分散行列を特定するパラメータは、フレームごとの前記ノイズフリー信号に含まれる残響をK次の自己回帰過程によりモデル化したときの自己回帰過程の予測係数行列Gk(k=1,2,…,K)と、フレームごとの前記クリーン信号の共分散行列と、を含む、
モデル推定装置。
【請求項４】
請求項１または２または３記載のモデル推定装置であって、
前記拡散性雑音の共分散行列は、フレームごとに、行列のなす線型空間における部分空間Vに属するとしてモデル化されたものである、
モデル推定装置。
【請求項５】
請求項２記載のモデル推定装置であって、
観測信号をモデル化するときのモデルパラメータの集合をΘとし、前記モデル推定部は、モデルパラメータの集合Θが与えられたときの前記観測信号ytuの集合Yと前記音源信号xtuの集合Xとの同時分布の対数log(p(X,Y;Θ))の前記集合Xに関する期待値が最大となるように、前記共分散行列Φxtと前記共分散行列Φvtとを更新する、
モデル推定装置。
【請求項６】
請求項５記載のモデル推定装置であって、
前記モデル推定部は、
前記観測信号ytuと現在の前記共分散行列Φxtと前記共分散行列Φvtとを用いて、
【数１１０】
2015135437.tif 000112
により、観測信号ytuが与えられたときの音源信号xtuの事後確率の平均μx|ytuと共分散行列Φx|ytとを更新する事後確率更新部と、
行列・のエルミート転置を・Hとし、観測信号を取得する際に用いるマイクロホンの個数をMとし、行列・の対角成分の和をTr[・]とし、空間的な特性を表すパラメータを表す信号コヒーレンス行列をBとし、前記平均μx|ytuと前記共分散行列Φx|ytとを用いて、
【数１１１】
2015135437.tif 000113
を求め、
【数１１２】
2015135437.tif 000114
により信号パワーφxtを更新する信号パワー更新手段と、
前記^Φxtと前記信号パワーφxtとを用いて、
【数１１３】
2015135437.tif 000115
により、前記信号コヒーレンス行列Bを更新する信号コヒーレンス行列更新手段と、
前記信号パワーφxtと前記信号コヒーレンス行列Bとを用いて、
【数１１４】
2015135437.tif 000116
により前記共分散行列Φxtを更新する信号共分散行列更新手段と、
前記部分空間Vの正規直交基底を{Qd}1≦d≦Dとし、m∈{1,…,M}とし、n∈{1,…,M}とし、m番目のマイクロホンとn番目のマイクロホンとの距離をLmnとし、音速をcとし、前記観測信号ytuに対応する周波数をfとし、sinc関数をsinc(・)とし、(m,n)成分として
【数１１５】
2015135437.tif 000117
を持つ雑音コヒーレンス行列をΓとし、前記平均μx|ytuと前記共分散行列Φx|ytと前記観測信号ytuとを用いて、
【数１１６】
2015135437.tif 000118
を求め、前記部分空間Vが
【数１１７】
2015135437.tif 000119
なる性質を持つ場合、
【数１１８】
2015135437.tif 000120
により、前記共分散行列Φvtを更新し、前記部分空間Vとして定数雑音コヒーレンスモデルを用いる場合、
【数１１９】
2015135437.tif 000121
により前記共分散行列Φｖtを更新する雑音共分散行列更新手段と、を含む、
モデル推定装置。
【請求項７】
請求項３記載のモデル推定装置であって、
観測信号をモデル化するときのモデルパラメータの集合をΘとし、前記モデル推定部は、モデルパラメータの集合Θが与えられたときの前記観測信号の集合Yと前記ノイズフリー信号の集合Ｒとの同時分布の対数log(p(R,Y;Θ))の前記集合Rに関する期待値が最大となるように、前記モデルパラメータの各々を更新する、
モデル推定装置。
【請求項８】
請求項７記載のモデル推定装置であって、
前記ブロックの番号をt∈{1,…,T}（Tはブロックの総数）とし、ブロック内のフレームの番号をu∈{1,…,U}とし、観測信号ytuのブロックtごとの共分散行列をΦytとし、クリーン信号xtuのブロックtごとの共分散行列をΦxtとし、拡散性雑音vtuのブロックtごとの共分散行列をΦvtとし、行列・のエルミート転置を・Hとし、
【数１２０】
2015135437.tif 000122
とし、
k∈{1,…,K}とし、i∈{1,…,T}とし、単位行列をIとし、所定の遅延をΔとし、j∈{1,…,T}とし、行列FはT2個のM×M行列をブロックに持つブロックテプリッツ行列であり、前記Fの(i,j)番目のブロックは
【数１２１】
2015135437.tif 000123
であるとして、
前記複数ブロック分のノイズフリー信号の共分散行列はF-H~ΦxF-1である、
モデル推定装置。
【請求項９】
請求項８記載のモデル推定装置であって、
前記モデル推定部は、
行列・の転置を・Tとし、~yu:=[yTuT…y1uT]Tとし、t'∈{1,…,T}とし、現在の前記予測係数行列Gkを用いて前記Fを更新し、現在の前記共分散行列Φxtを用いて前記~Φxを更新し、前記Fと前記~Φxとを用いて、
【数１２２】
2015135437.tif 000124
により、~Φrを更新し、現在の前記共分散行列Φvtを用いて前記~Φvを更新し、前記~Φrと前記~Φvと前記観測信号ytuとを用いて、
【数１２３】
2015135437.tif 000125
を求め、前記~μur|yを用いて、
【数１２４】
2015135437.tif 000126
により、観測信号ytuが与えられたときのノイズフリー信号rtuの事後確率の平均μtur|yを更新し、前記~Φrと前記~Φvとを用いて、
【数１２５】
2015135437.tif 000127
を求め、前記~Φr|yを用いて、
【数１２６】
2015135437.tif 000128
により、観測信号ytuが与えられたときのブロックtおよびブロックt’におけるノイズフリー信号rtuおよびrt'uの事後共分散行列Φtt'r|yを更新する事後確率更新部と、
~G:=[G1T…GKT]Tとし、行列・の対角成分の和をTr[・]とし、空間的な特性を表すパラメータを表す信号コヒーレンス行列をBxとし、前記平均μtur|yと前記共分散行列Φtt'r|yとを用いて、
【数１２７】
2015135437.tif 000129
を求め、前記~Gと前記Ψtt'とを用いて、
【数１２８】
2015135437.tif 000130
を求め、前記^Φxtと現在の信号コヒーレンス行列Bxとを用いて、
【数１２９】
2015135437.tif 000131
により信号パワーφxtを更新する信号パワー更新手段と、
前記^Φxtと前記信号パワーφxtとを用いて、
【数１３０】
2015135437.tif 000132
により、前記信号コヒーレンス行列Bxを更新する信号コヒーレンス行列更新手段と、
前記信号パワーφxtと前記信号コヒーレンス行列Bxとを用いて、
【数１３１】
2015135437.tif 000133
により前記共分散行列Φxtを更新する信号共分散行列更新手段と、
前記共分散行列Φvtは、行列のなす線型空間における部分空間Vに属するとしてモデル化されたものであり、部分空間Vの正規直交基底を{Qd}1≦d≦Dとし、m∈{1,…,M}とし、n∈{1,…,M}とし、m番目のマイクロホンとn番目のマイクロホンとの距離をLmnとし、音速をcとし、前記観測信号ytuに対応する周波数をfとし、sinc関数をsinc(・)とし、(m,n)成分として
【数１３２】
2015135437.tif 000134
を持つ雑音コヒーレンス行列をBvとし、前記平均μtur|yと前記共分散行列Φtt'r|yとを用いて、
【数１３３】
2015135437.tif 000135
を求め、前記部分空間Vが
【数１３４】
2015135437.tif 000136
なる性質を持つ場合、
【数１３５】
2015135437.tif 000137
により、前記共分散行列Φvtを更新し、前記部分空間Vが定数雑音コヒーレンスモデルである場合、
【数１３６】
2015135437.tif 000138
により、前記共分散行列Φvtを更新する雑音共分散行列更新手段と、
前記Ψtt'と前記信号パワーφxtとを用いて、
【数１３７】
2015135437.tif 000139
を求め、前記~Gを用いて、前記予測係数行列Gkを更新する予測係数行列更新手段と、を含む、
モデル推定装置。
【請求項１０】
請求項２または３のモデル推定装置により推定されたモデルパラメータと、前記観測信号ytuとを用いて、前記観測信号ytuが与えられたときの前記ノイズフリー信号の事後確率の平均を求める平均計算部と、
前記事後確率の平均を時間領域の信号に変換し、時間領域における前記ノイズフリー信号の推定値を求める時間領域変換部と、を含む
雑音抑圧装置。
【請求項１１】
請求項３または８または９のいずれかのモデル推定装置により推定されたモデルパラメータと、前記観測信号ytuとを用いて、前記観測信号ytuが与えられたときの前記ノイズフリー信号の事後確率の平均μr|ytuを求める平均計算部と、
請求項３または８または９のいずれかのモデル推定装置により推定された前記予測係数行列Gkと前記平均計算部で求めた事後確率の平均μtur|yとを用いて、前記観測信号ytuが与えられたときの前記クリーン信号の事後確率の平均^xtuを計算する強調音声生成部と、
前記平均^xtuを時間領域の信号に変換し、時間領域における前記クリーン信号の推定値を求める時間領域変換部と、を含む
音声強調装置。
【請求項１２】
観測信号を、非定常の拡散性雑音と前記拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、モデル化するときのモデルパラメータは拡散性雑音の共分散行列を特定するモデルパラメータと、ノイズフリー信号の共分散行列を特定するモデルパラメータとを含み、
観測信号から得られる共分散行列と、拡散性雑音の共分散行列とノイズフリー信号の共分散行列との和との板倉-齋藤距離が小さくなるように、前記モデルパラメータを更新するモデル推定ステップを含む、
モデル推定方法。
【請求項１３】
観測信号を、非定常の拡散性雑音と前記拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、１つ以上のフレームからなるブロックの番号をt∈{1,…,T}とし、ブロック内のフレームの番号をu∈{1,…,U}とし、ブロックtごとの観測信号ytuの共分散行列をΦytとし、ブロックtごとのノイズフリー信号xtuの共分散行列をΦxtとし、ブロックtごとの拡散性雑音vtuの共分散行列をΦvtとし、
観測信号ytから得られる共分散行列Φytと、Φxt＋Φvtとの板倉-齋藤距離の、すべてのブロックtについての総和が小さくなるように、前記共分散行列Φxtと前記共分散行列Φvtとを更新するモデル推定ステップを含み、
前記ノイズフリー信号は音源信号である、
モデル推定方法。
【請求項１４】
フレームごとの観測信号を非定常の拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、１つ以上のフレームからなる時間区間をブロックとして、
複数ブロック分の観測信号から得られる共分散行列と、複数ブロック分の拡散性雑音の共分散行列と複数ブロック分のノイズフリー信号の共分散行列との和と、の板倉-齋藤距離が小さくなるように、前記複数ブロック分の拡散性雑音の共分散行列を特定するパラメータと前記複数ブロック分のノイズフリー信号の共分散行列を特定するパラメータとを更新するモデル推定ステップを含み、
前記ノイズフリー信号は、音源から発せられる直接音と初期反射成分とからなる信号であるクリーン信号に残響が重畳された信号であり、
前記複数ブロック分のノイズフリー信号の共分散行列を特定するパラメータは、フレームごとの前記ノイズフリー信号に含まれる残響をK次の自己回帰過程によりモデル化したときの自己回帰過程の予測係数行列Gk(k=1,2,…,K)と、フレームごとの前記クリーン信号の共分散行列と、を含む、
モデル推定方法。
【請求項１５】
請求項１３のモデル推定方法により推定されたモデルパラメータと、前記観測信号ytuとを用いて、前記観測信号ytuが与えられたときの前記ノイズフリー信号の事後確率の平均を求める平均計算ステップと、
前記事後確率の平均を時間領域の信号に変換し、時間領域における前記ノイズフリー信号の推定値を求める時間領域変換ステップと、を含む
雑音抑圧方法。
【請求項１６】
請求項１４のモデル推定方法により推定されたモデルパラメータと、前記観測信号ytuとを用いて、前記観測信号ytuが与えられたときの前記ノイズフリー信号の事後確率の平均μr|ytuを求める平均計算ステップと、
請求項１４のモデル推定方法により推定された前記予測係数行列Gkと前記平均計算ステップで求めた事後確率の平均μtur|yとを用いて、前記クリーン信号の事後確率の平均^xtuを計算する強調音声生成ステップと、
前記平均^xtuを時間領域の信号に変換し、時間領域における前記クリーン信号の推定値を求める時間領域変換ステップと、を含む
音声強調方法。
【請求項１７】
請求項１〜９の何れかのモデル推定装置、または、請求項１０の雑音抑圧装置、または、請求項１１の音声強調装置として、コンピュータを機能させるためのプログラム。
【発明の詳細な説明】
【技術分野】
【０００１】
本発明は、複数のマイクロホンで観測された、拡散性雑音を含む観測信号を用いて、拡散性雑音を含まない信号であるノイズフリー信号を推定する雑音抑圧技術、及び、雑音及び残響を含まないクリーン信号を推定する音声強調技術に関する。さらに、観測信号を、拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとして、雑音抑圧時または音声強調時等に用いるモデルパラメータの推定技術に関する。
【背景技術】
【０００２】
＜拡散性雑音下における雑音抑圧及び音声強調＞
拡散性雑音とは、多数の雑音源、または空間的な広がりをもつ雑音源により生じる、様々な方向から到来する雑音を指す。多数の雑音源による拡散性雑音の例としては、食堂での人々の話し声や食器の音があり、空間的な広がりをもつ雑音源による拡散性雑音の例としては、電車内での車体の振動による雑音がある。このような拡散性雑音は、ビームフォーミングのような、従来の雑音抑圧技術による抑圧が困難であり、雑音抑圧技術の適用範囲を大幅に制限する要因となっていた。また、拡散性雑音は、少数の点音源から生じる雑音と比べて、モデル化がより難しく、実環境における音声強調技術の応用範囲を大きく制限する要因となっていた。
【０００３】
観測信号を、拡散性雑音と拡散性雑音を含まない信号との和でモデル化する場合、モデル化したときのモデルパラメータとして、拡散性雑音の共分散行列Φvt及び音源信号の共分散行列Φxtを含み、共分散行列Φvt及びΦxtの推定精度を上げることで、雑音抑圧性能及び音声強調性能を向上させることができる。
【０００４】
信号と拡散性雑音の空間的な性質を利用して、共分散行列Φvt及びΦxtを推定する従来技術として非特許文献１が知られている。
【０００５】
非特許文献１では、次式で定義されるコスト関数を最小化することにより、共分散行列Φvt及びΦxtを推定する。
【０００６】
【数１】
2015135437.tif 000003
【０００７】
ここで、DEU(Φyt,Φxt+Φvt):=||Φyt-(Φxt+Φvt)||F2は、観測信号を用いて得られる観測信号の共分散行列Φytと共分散行列のモデルΦxt+Φvtとのユークリッド距離の平方であり、||・||Fはフロベニウスノルムである。なお、A:=Bは、「Aの意味するところをBと定義する」ことを意味し、式中「:=」を
【０００８】
【数２】
2015135437.tif 000004
【０００９】
とも記載する。
【先行技術文献】
【非特許文献】
【００１０】
【非特許文献１】N. Ito, E. Vincent, N. Ono, and S. Sagayama, “General algorithms for estimating spectrogram and transfer functions of target signal for blind suppression of diffuse noise,” in Proc. IEEE International Workshop on Machine Learning for Signal Processing (MLSP), 2013.
【発明の概要】
【発明が解決しようとする課題】
【００１１】
上記のように、非特許文献１では、簡単のため、式(a1)のように、ユークリッド距離の平方DEUという単純な尺度により、ΦytとΦxt+Φvtとの間の距離を測っている。しかしながら、より高い音声品質、ならびに、より高精度な音声認識の実現のためには、このような尺度は必ずしも最適ではない。なぜなら、音声品質の向上であれ、音声認識精度の向上であれ、人間の聴覚も、音声認識の特徴量(例えばMFCC: mel - frequency cepstral coefficients)も、ともに対数的であるため、ユークリッド距離の平方ではなく、対数パワー領域での誤差を小さくする必要があるからである。実際、例えば、振幅が大きい信号に多少の歪みが生じても、聴感上さほど気にならないのに対し、振幅が小さい信号に少しでも歪みが生じると、敏感に知覚される。
【００１２】
本発明は、このような背景に鑑みてなされたものであり、音響信号により適した距離尺度に基づく、拡散性雑音の存在下での雑音抑圧技術、音声強調技術、雑音抑圧時及び音声強調時等に用いるモデルパラメータの推定技術を提供することを目的とする。
【課題を解決するための手段】
【００１３】
上記の課題を解決するために、本発明の第一の態様によれば、観測信号を、非定常の拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、モデル化するときのモデルパラメータは拡散性雑音の共分散行列を特定するモデルパラメータと、ノイズフリー信号の共分散行列を特定するモデルパラメータとを含み、モデル推定装置は、観測信号から得られる共分散行列と、拡散性雑音の共分散行列とノイズフリー信号の共分散行列との和との板倉-齋藤距離が小さくなるように、モデルパラメータを更新するモデル推定部を含む。
【００１４】
上記の課題を解決するために、本発明の他の態様によれば、観測信号を、非定常の拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、１つ以上のフレームからなるブロックの番号をt∈{1,…,T}とし、ブロック内のフレームの番号をu∈{1,…,U}とし、ブロックtごとの観測信号ytuの共分散行列をΦytとし、ブロックtごとのノイズフリー信号xtuの共分散行列をΦxtとし、ブロックtごとの拡散性雑音vtuの共分散行列をΦvtとし、モデル推定装置は、観測信号ytから得られる共分散行列Φytと、Φxt＋Φvtとの板倉-齋藤距離の、すべてのブロックtについての総和が小さくなるように、共分散行列Φxtと共分散行列Φvtとを更新するモデル推定部を含む。ノイズフリー信号は音源信号である。
【００１５】
上記の課題を解決するために、本発明の他の態様によれば、フレームごとの観測信号を非定常の拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、１つ以上のフレームからなる時間区間をブロックとして、モデル推定装置は、複数ブロック分の観測信号から得られる共分散行列と、複数ブロック分の拡散性雑音の共分散行列と複数ブロック分のノイズフリー信号の共分散行列との和と、の板倉-齋藤距離が小さくなるように、複数ブロック分の拡散性雑音の共分散行列を特定するパラメータと複数ブロック分のノイズフリー信号の共分散行列を特定するパラメータとを更新するモデル推定部を含む。ノイズフリー信号は、音源から発せられる直接音と初期反射成分とからなる信号であるクリーン信号に残響が重畳された信号である。複数ブロック分のノイズフリー信号の共分散行列を特定するパラメータは、フレームごとのノイズフリー信号に含まれる残響をK次の自己回帰過程によりモデル化したときの自己回帰過程の予測係数行列Gk(k=1,2,…,K)と、フレームごとのクリーン信号の共分散行列と、を含む。
【００１６】
上記の課題を解決するために、本発明の他の態様によれば、観測信号を、非定常の拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、モデル化するときのモデルパラメータは拡散性雑音の共分散行列を特定するモデルパラメータと、ノイズフリー信号の共分散行列を特定するモデルパラメータとを含み、モデル推定方法は、観測信号から得られる共分散行列と、拡散性雑音の共分散行列とノイズフリー信号の共分散行列との和との板倉-齋藤距離が小さくなるように、モデルパラメータを更新するモデル推定ステップを含む。
【００１７】
上記の課題を解決するために、本発明の他の態様によれば、観測信号を、非定常の拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、１つ以上のフレームからなるブロックの番号をt∈{1,…,T}とし、ブロック内のフレームの番号をu∈{1,…,U}とし、ブロックtごとの観測信号ytuの共分散行列をΦytとし、ブロックtごとのノイズフリー信号xtuの共分散行列をΦxtとし、ブロックtごとの拡散性雑音vtuの共分散行列をΦvtとし、モデル推定方法は、観測信号ytから得られる共分散行列Φytと、Φxt＋Φvtとの板倉-齋藤距離の、すべてのブロックtについての総和が小さくなるように、共分散行列Φxtと共分散行列Φvtとを更新するモデル推定ステップを含む。ノイズフリー信号は音源信号である。
【００１８】
上記の課題を解決するために、本発明の他の態様によれば、フレームごとの観測信号を非定常の拡散性雑音と拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化するものとし、１つ以上のフレームからなる時間区間をブロックとして、モデル推定方法は、複数ブロック分の観測信号から得られる共分散行列と、複数ブロック分の拡散性雑音の共分散行列と複数ブロック分のノイズフリー信号の共分散行列との和と、の板倉-齋藤距離が小さくなるように、複数ブロック分の拡散性雑音の共分散行列を特定するパラメータと複数ブロック分のノイズフリー信号の共分散行列を特定するパラメータとを更新するモデル推定ステップを含む。ノイズフリー信号は、音源から発せられる直接音と初期反射成分とからなる信号であるクリーン信号に残響が重畳された信号である。複数ブロック分のノイズフリー信号の共分散行列を特定するパラメータは、フレームごとのノイズフリー信号に含まれる残響をK次の自己回帰過程によりモデル化したときの自己回帰過程の予測係数行列Gk(k=1,2,…,K)と、フレームごとのクリーン信号の共分散行列と、を含む。
【発明の効果】
【００１９】
本発明によれば、従来技術よりも推定精度よく雑音抑圧時及び音声強調時等に用いるモデルパラメータの推定することができるという効果を奏する。
【図面の簡単な説明】
【００２０】
【図１】第一実施形態に係るモデル推定装置の機能ブロック図。
【図２】第一実施形態に係るモデル推定装置の処理フローの例を示す図。
【図３】第二実施形態に係る雑音抑圧装置の機能ブロック図。
【図４】第二実施形態に係る雑音抑圧装置の処理フローの例を示す図。
【図５】第二実施形態に係る雑音抑圧装置の効果を説明するための図。
【図６】第三実施形態に係るモデル推定装置の機能ブロック図。
【図７】第三実施形態に係るモデル推定装置の処理フローの例を示す図。
【図８】第四実施形態に係る音声強調装置の機能ブロック図。
【図９】第四実施形態に係る音声強調装置の処理フローの例を示す図。
【図１０】第四実施形態に係る音声強調装置の効果を説明するための図。
【図１１】評価指標による客観評価の結果を示す図。
【図１２】第五実施形態に係る雑音抑圧装置の機能ブロック図。
【図１３】第五実施形態に係る雑音抑圧装置の処理フローの例を示す図である。
【発明を実施するための形態】
【００２１】
以下、本発明の実施形態について説明する。なお、以下の説明に用いる図面では、同じ機能を持つ構成部や同じ処理を行うステップには同一の符号を記し、重複説明を省略する。以下の説明において、テキスト中で使用する記号「^」、「~」等は、本来直後の文字の真上に記載されるべきものであるが、テキスト記法の制限により、当該文字の直前に記載する。式中においてはこれらの記号は本来の位置に記述している。
【００２２】
＜第一実施形態＞
本実施形態では、観測信号から拡散性雑音抑圧時に用いるモデルパラメータを推定する。なお、拡散性雑音抑圧とは，上述の拡散性雑音を抑圧する技術である。
【００２３】
＜拡散性雑音抑圧の問題設定＞
本実施形態では、特に断りのない限り、各信号を短時間フーリエ変換（short-time Fouriertransform:STFT）などの時間周波数領域で表現する。フレームの総数をIで表し、フレームの番号をi∈{1,…,I}で表す。また、本実施形態では、簡潔さのため、周波数ビンの番号の表示を省略するが、この省略が混乱をもたらすおそれは小さい。本実施形態では、周波数ビン毎に独立に処理を行うからである。
【００２４】
M個のマイクロホンによる観測信号をyi∈CM、音源信号をxi∈CM、拡散性雑音をvi∈CMにより表す。ただし、Cは複素数の全体の集合である。この記法を用いて、拡散性雑音抑圧の定義を言い換えると、「Y:={yi}Ii=1が与えられたときに、X:={xi}Ii=1を推定する問題」ということになる。ここで、記法{yi}Ii=1は、{yi|i=1,…,I}を表す。ここで、xi、vi、yiは、いずれも、M個のマイクロホン位置における各種信号を並べたM次元のベクトルである。  本実施形態では、観測信号を収録する音場は、静止した単一の点音源に起因する音源信号と、拡散性雑音とからなり、観測信号yiは、次式で表されるとモデル化する。
yi=xi+vi    (b1)
つまり、観測信号yiを、拡散性雑音viと拡散性雑音を含まない信号（ノイズフリー信号）である音源信号xiとの和でモデル化する。ここで、簡単のため、音源信号xiと拡散性雑音viとについて、以下の仮定を置く。
・時間的独立性:{xi}Ii=1は独立な系列である。すなわち、任意の相異なるi、jに対し（i∈{1,…,I}、j∈{1,…,I}、i≠j）、音源信号xiと音源信号xjとは独立である。また、{vi}Ii=1も独立な系列である。
・相互独立性:{xi}Ii=1と{vi}Ii=1とは、互いに独立である。すなわち、任意のi、jに対し（i∈{1,…,I}、j∈{1,…,I}）、音源信号xiと拡散性雑音vjとは独立である。
・ガウス性:音源信号xiと拡散性雑音viとは、平均が0の複素ガウス分布に従う。
ここで、平均μ∈CM、共分散行列Σ∈CM×Mの複素ガウス分布の確率密度関数は、次式で与えられる。
【００２５】
【数３】
2015135437.tif 000005
【００２６】
ただし、det・は行列・の行列式、・Hは行列・のエルミート転置、exp[・]は自然対数の底数の・乗である。
【００２７】
さらに、短時間定常性をモデル化するため、I個のフレームを、連続するフレームからなる、T個のブロックに区分し、各ブロックにおいて、音源信号xi及び拡散性雑音viの共分散行列は一定と仮定する。ブロックの番号をt=1,…,T、各ブロックt内のフレームの番号をu=1,…,Uとするとき、iと(t,u)との対応
【００２８】
【数４】
2015135437.tif 000006
【００２９】
すなわち、一般に、iをUで割ったときの商および剰余をそれぞれqおよびrとするとき、
【００３０】
【数５】
2015135437.tif 000007
【００３１】
である。このように定められたt、uを用いて、xi→xtu、vi→vtuのように番号を振りなおすとき、短時間定常性は次のように表すことができる。
・短時間定常性:音源信号xtuおよび拡散性雑音vtuの共分散行列は、ブロックの番号tには依存するが、各ブロック内でのフレームの番号uには依存せず、Φxt、Φvt∈CM×Mと置ける。
【００３２】
なお、フレームの総数Iは、Uの整数倍とは限らないため、端数は最後のブロックで調整するものとする。すなわち、最後のブロックは、I-U(T-1)個のフレームからなる。
【００３３】
上述の仮定のもと、観測信号ytuが与えられたときの、音源信号xtuの最大事後確率推定量を導くことができる。これは、観測信号ytuが与えられたときの音源信号xtuの事後確率p(xtu|ytu;Θ)を最大化するxtuと定義される。そこで、この事後確率を導くため、まず同時確率を導くと、
【００３４】
【数６】
2015135437.tif 000008
【００３５】
となる。式(b5)およびベイズの定理より、事後確率は、
【００３６】
【数７】
2015135437.tif 000009
【００３７】
ここで、
【００３８】
【数８】
2015135437.tif 000010
【００３９】
であり、
【００４０】
【数９】
2015135437.tif 000011
【００４１】
は、xtuに依存しない定数を除いて等しいことを表す。よって、
【００４２】
【数１０】
2015135437.tif 000012
【００４３】
であるから、xtuの最大事後確率推定量は、(b8)で与えられる。なお、最大事後確率推定量(b10)は、最小平均二乗誤差推定量であるマルチチャネルウィーナーフィルタ(参考文献１〜４参照)と一致している。
(参考文献１) N.Q.K. Duong, E. Vincent, and R. Gribonval, “Under-determined reverberant audio source separation using a full-rank spatial covariance model”, IEEE Trans. ASLP, Sep. 2010, vol. 18, no. 7, pp. 1830-1840.
(参考文献２)H.L. Van Trees, "Optimum Array Processing", John Wiley & Sons, NewYork, 2002.
(参考文献３)K.U. Simmer, J. Bitzer, and C. Marro, “Post-filtering techniques”, in
Microphone Arrays, M. Brandstein and D.Ward, Eds., pp. 39-60. Springer, Berlin Heidelberg, 2001.
(参考文献４)S. Doclo and M. Moonen, “GSVD-based optimal filtering for single and multimicrophone speech enhancement”, IEEE Trans. SP, Sep. 2002, vol. 50, no. 9,pp. 2230-2244. 
これは、ガウス分布の平均と最頻値が一致することに起因する。
【００４４】
マルチチャネルウィーナーフィルタ(b8)の設計は、共分散行列ΦxtおよびΦvtの推定に帰着する。この推定問題は、観測信号ytuを用いて推定できる観測信号ytuの共分散行列を
【００４５】
【数１１】
2015135437.tif 000013
【００４６】
と置くと、次のように定式化できる。
【００４７】
【数１２】
2015135437.tif 000014
【００４８】
ここで、D(・,・)は何らかの距離尺度である。また、Σuは、t=1,…,T-1に対しては、ΣUu=1を意味し、t=Tに対しては、ΣI-U(T-1)u=1を意味する。(b12)において、共分散行列ΦxtとΦvtとは共に未知である。したがって、(b12)は劣決定の逆問題であり、制約条件がない場合、(b12)の最適解は無数に存在する。実際、^Φxt、^Φvtが最適解ならば、Φxt+Φvt=^Φxt+^Φvtを満たす任意の共分散行列Φxt、Φvtは最適解である。したがって、この逆問題を解くためには、共分散行列ΦxtとΦvtとの性質を利用して、これらを適切にモデル化する必要がある。
【００４９】
＜共分散行列を用いた信号と雑音のモデル化＞
共分散行列ΦxtとΦvtとをモデル化するために、本実施形態では、信号と雑音の空間的な性質を利用する。
【００５０】
クリーン信号の信号源が、静止した単一の点音源であるという仮定の下では、音源信号の空間的な特性は、時間によらず一定と考えられる。この空間的な特性の時間的な不変性に着目して、次式のような共分散行列Φxtのモデルを考えることができる。
【００５１】
【数１３】
2015135437.tif 000015
【００５２】
ここで、φxt∈R(Rは実数の全体の集合)は、スペクトルに対応する、時間依存のパラメータであり、B∈CM×Mは、空間的な特性を表す、時間非依存のパラメータである。以下、φxt及びBをそれぞれ信号パワー及び信号コヒーレンス行列ともいう。
【００５３】
非特許文献１でも(b13)と類似の共分散行列のモデルが用いられているが、本実施形態における音源信号の共分散行列のモデルでは、(b13)において、信号コヒーレンス行列Bに何の制約も課さないのに対し、非特許文献１に記載の音源信号の共分散行列のモデルでは、(b13)において、信号コヒーレンス行列Bがランク１行列hhHであるという制約を課している点が異なる（ｈはステアリングベクトル）。ランク１行列の逆行列は定義できないため、非特許文献１の音源信号の共分散行列のモデルに対しては、音源信号をガウス分布でモデル化することができず、したがって、音響信号に適した尺度である板倉-齋藤距離を用いることができない。これに対し、本実施形態では、信号コヒーレンス行列Bに上記のような制約を課さないため、板倉-齋藤距離を適用することができる。
【００５４】
音源信号とは異なり、拡散性雑音は、不特定多数の音源により生じるため、たとえそれらの音源がすべて静止していたとしても、各フレームで異なる音源がアクティブとなることにより、雑音全体の空間的な特性は時間的に変化する。また、多数の話者が同時に発話する場合のように、拡散性雑音は、スペクトルも時間変化することが多い。このように、拡散性雑音は、一般にスペクトルも空間的な特性も時間的に変化するため（非定常であるため）、時間的な不変性に基づかずにモデル化する必要がある。これに対し、非特許文献１及び参考文献５では、等方性という空間的な不変性に着目して、拡散性雑音をモデル化することが提案された。
（参考文献５）N. Ito, "Robust Microphone Array Signal Processing against Diffuse Noise", Ph.D. thesis, the University of Tokyo, 2012.
拡散性雑音は、不特定多数の音源に起因するため、方向依存性は小さく、等方的とみなせる。そこで、この等方性のモデルとして、非特許文献１及び参考文献５では、「拡散性雑音の二点間のクロススペクトルは、その二点間の方向には依存せず、その距離のみで決まる」と仮定する。この仮定のもとで、「Φvtは、M×Mエルミート行列全体がなす実ベクトル空間Hの、ある低次元部分空間Vに属する」ことが示せる（非特許文献１及び参考文献５参照）。本実施形態においても非特許文献１及び参考文献５と同じ拡散性雑音のモデルを仮定する。すなわち、我々の拡散性雑音モデルは、次式で表される。
Φvt∈V    (b14)
Vは、行列がなすベクトル空間Hの部分空間であるから、行列部分空間と呼ぶ。換言すると、Vは、行列のなす線型空間における部分空間である。
【００５５】
行列部分空間Vの具体形は、アレイ配置や雑音場に関する付加的な仮定などによって、複数提案されている。まず、任意のアレイ配置に適用可能で、雑音に関する付加的な仮定も置かない、最も一般的な場合が、実数値雑音共分散（real-valued noise covariance）モデル
【００５６】
【数１４】
2015135437.tif 000016
【００５７】
である。また、結晶型（参考文献５及び６参照）と呼ばれる、特別なクラスのアレイ配置の場合には、ブラインド雑音無相関化（blind noise decorrelation:BND）モデル
【００５８】
【数１５】
2015135437.tif 000017
【００５９】
が適用できる。ここで、Uは、既知の定数ユニタリ行列である。
（参考文献６）N. Ito, H. Shimizu, N. Ono, and S. Sagayama, “Diffuse noise suppression using crystal-shaped microphone arrays”, 2011, vol. 19, no. 7, pp. 2101-2110.
また、波長に比べ十分距離の大きい2点間では、拡散性雑音は無相関である、と仮定すると、マイクロホン間の距離を大きくした極限で、空間的無相関雑音（spatially uncorrelated noise）モデル
【００６０】
【数１６】
2015135437.tif 000018
【００６１】
が成り立つ。最後に、マイクロホン間の相対的な位置関係が既知で、3次元実空間R3の任意の方向から、スペクトルの等しい雑音が、平面波として到来するという仮定の下では、定数雑音コヒーレンス（fixed noise coherence）モデル
【００６２】
【数１７】
2015135437.tif 000019
【００６３】
が成立する。ここで、Γはコヒーレンス行列と呼ばれ、その(m,n)成分は、
【００６４】
【数１８】
2015135437.tif 000020
【００６５】
で与えられる。ここで、sinc(・)はsinc関数である。また、Lmnはm番目とn番目のマイクロホンの距離、cは音速であり、fは観測信号ytuに対応する周波数であり、Lmnとしては実測値を用い、cとしては、例えば、近似値c=340m/sを用いればよい。
【００６６】
これらの雑音モデルの比較については、参考文献５を参照されたい。本実施形態の方法は、これらの全ての雑音モデルに適用可能である。
【００６７】
これらの拡散性雑音の共分散行列のモデルは、拡散性雑音の共分散行列が時間的に変化することを許しているため、非定常の拡散性雑音も扱うことができるという大きな利点がある。従来の多くの雑音抑圧技術では、定常雑音か、少数の点音源から生じる方向性雑音を仮定することが多かった。しかしながら、駅・空港、乗り物、オフィス、工場、カフェ、居酒屋、繁華街といった、さまざまな実環境において、非定常の拡散性雑音が存在する。定常雑音か、少数の点音源から生じる方向性雑音を仮定する従来技術では、このような雑音は扱えないため、このことが、雑音抑圧技術の実環境での応用範囲を大きく制限していた。これに対し、本実施形態における雑音の共分散行列のモデル化では、雑音の空間的な性質に着目することで、このような非定常の拡散性雑音もモデル化できる。
【００６８】
本実施形態における雑音の共分散行列のモデル化に基づくモデル推定の方法として、すでに非特許文献１において、ユークリッド距離を評価尺度とする手法が提案されていた。しかし、既述の通り、非特許文献１では音源信号のモデルに制約を設けているため板倉‐齋藤距離を最小化することはできない。つまり、非定常の拡散性雑音を扱えるようなモデルにおいて、板倉-齋藤距離を最適化する手法の導出が可能であるかどうかは、知られていなかった。
【００６９】
＜拡散性雑音抑圧の従来技術＞
従来技術の課題についてより詳細に説明する。
【００７０】
拡散性雑音抑圧は、音源信号の共分散行列Φxtと拡散性雑音の共分散行列Φvtとを推定する問題に帰着させることができる。すなわち、これらを推定することができれば、マルチチャネルウィーナーフィルタ(b8)により、観測信号ytuを用いて、音源信号xtuを推定することができる。この推定値の精度は、共分散行列ΦxtとΦvtとの推定精度に依存するため、共分散行列ΦxtとΦvtとを正確に推定することが非常に重要である。
【００７１】
非特許文献１では、この問題を解くために、次式で定義されるコスト関数を用いていた。
【００７２】
【数１９】
2015135437.tif 000021
【００７３】
簡単のため、フレーム数Iが、ブロックあたりのフレーム数Uで割り切れ、端数が生じない場合で説明した。さらに、残響はフレーム長に比べて十分短いなどの仮定の下、前述の雑音モデルの下で、(b20)の最小化により共分散行列ΦxtとΦvtとを推定するための反復アルゴリズムを導くことができる。
【００７４】
しかし、前述の通り、ユークリッド距離の平方ではなく、対数パワー領域での誤差を小さくする必要がある。
【００７５】
＜本実施形態のポイント＞
本実施形態では、非特許文献１のように、ユークリッド距離の平方の総和(b20)を最小化するのではなく、最尤法により音源信号の共分散行列Φxtと拡散性雑音の共分散行列Φvtとを推定する。後述のように、この最尤推定は、次式のコスト関数の最小化と等価である。
【００７６】
【数２０】
2015135437.tif 000022
【００７７】
このコスト関数は、観測信号ytから得られる共分散行列Φytと、Φxt＋Φvtとの板倉-齋藤距離の、すべてのブロックtについての総和を表し、本実施形態では、この値が小さくなるように共分散行列Φxtと共分散行列Φvtとを更新する。ここで、DISは、次式のマルチチャネル板倉-齋藤ダイバージェンスである（参考文献７〜９参照）。
【００７８】
【数２１】
2015135437.tif 000023
【００７９】
（参考文献７）H. Sawada, H. Kameoka, S. Araki, and N. Ueda, “Multichannel extensions of non-negative matrix factorization with complex-valued data”, IEEE Trans. ASLP, May 2013, vol. 21, no. 5, pp. 971-982.
（参考文献８）K. Yoshii, R. Tomioka, D. Mochihashi, and M. Goto, “Infinite positive semidefinite tensor factorization for source separation of mixture signals”, in Proc. International Conference on Machine Learning (ICML), Jun. 2013, pp. 576-584.
（参考文献９）B. Kulis, M. Sustik, and I. Dhillon, “Low-rank kernel learning with Bregman matrix divergences”, Journal of Machine Learning Research, Feb. 2009, vol. 10, pp. 341-376.
【００８０】
つまり、本実施形態では、観測信号yｔを拡散性雑音vｔと拡散性雑音を含まない信号であるノイズフリー信号xｔとの和でモデル化し、観測信号yｔから得られる共分散行列Φｙtと、拡散性雑音の共分散行列Φvtとノイズフリー信号の共分散行列Φxtとの和（Φxt＋Φvt）と、の板倉-齋藤距離が小さくなるように各モデルパラメータ（ここではΦvtとΦxt）を更新する。本実施形態では、拡散性雑音の共分散行列を特定するモデルパラメータはΦvtそのものであり、ノイズフリー信号の共分散行列を特定するモデルパラメータはΦxtそのものである。拡散性雑音の共分散行列Φvtは、観測信号yｔを拡散性雑音vｔと音源信号xｔとの和でモデル化したときの拡散性雑音vｔの統計的性質を表すモデルパラメータであり、音源信号の共分散行列Φxtは上記モデル化における音源信号xｔの統計的性質を表すモデルパラメータであるともいえる。
【００８１】
また、簡単のため、フレーム数Iが、ブロックあたりのフレーム数Uで割り切れ、端数が生じない場合で説明した。スカラーk>0に対して、DIS(kA,kB)=DIS(A,B)が成り立つことから、DISは、上記の音響信号の対数性を反映した、適切な尺度となっていることが分かる。なお、(b22)は、NTT(登録商標)の板倉・齋藤(参考文献１０参照)により提案され、音声処理における有効性が示されてきた、板倉-齋藤ダイバージェンスを、複数チャネルの場合へ拡張した尺度とみなせる。
(参考文献１０)F. Itakura and S. Saito, “Analysis synthesis telephony based on the maximum likelihood method”, in Rep. 6th International Congress on Acoustics, 1968, pp. C-17-C-20.
【００８２】
＜本実施形態における共分散行列推定の原理＞
前述のように、本実施形態では、最尤法によりΘ:={{Φxt}Tt=1,{Φvt}Tt=1}={{φxt}Tt=1,B,{Φvt}Tt=1}を推定する（(b13）参照)。目的関数である対数尤度は、次式により与えられる。
【００８３】
【数２２】
2015135437.tif 000024
【００８４】
以下では、簡単のため、フレーム数Iが、ブロックあたりのフレーム数Uで割り切れ、端数が生じない場合で説明するが、端数が生じる場合への拡張は容易である。この場合、上式は、以下のように単純化できる。
【００８５】
【数２３】
2015135437.tif 000025
【００８６】
したがって、J(Θ)の最大化は、(b21)の最小化と等価である。拘束条件(b14)に注意して、Θを推定する問題は、以下のように定式化される。
【００８７】
【数２４】
2015135437.tif 000026
【００８８】
J(Θ)は、log det(Φxt+Φvt)および(Φxt+Φvt)-1なる項を含むため、各パラメータに関するJ(Θ)の偏微分は複雑であり、(b27)の最適解を閉形式で求めることは困難である。一方、音源信号xtuを隠れ変数とみなし、expectation-maximization(EM)アルゴリズム(参考文献１１参照)を用いれば、効率的な最適化が可能である。
（参考文献１１）A.P. Dempster, N.M. Laird, and D.B. Rubin, “Maximum likelihood from incomplete data via the EM algorithm”, Journal of the Royal Statistical Society: Series B (Methodological), 1977, vol. 39, no. 1, pp. 1-38.
EMアルゴリズムは、次のEステップとMステップを交互に反復するものであり、(b27)の局所解への収束が保証される。
・Eステップ:現在のパラメータの推定値Θ'と(b10)とを用いて、音源信号xtuの事後確率p(xtu|ytu;Θ')を計算する。本実施形態における(b10)はガウス分布であり、平均μx|ytuおよび共分散行列Φx|ytにより完全に決定されるため、Θ'と(b8)(b9)とを用いて、推定値(μx|ytu)'および(Φx|yt)'を計算すれば十分である。ただし、(・)'は、現在のパラメータの推定値Θ'を用いて計算することを示す。
・Mステップ:Eステップで得た事後確率p(xtu|ytu;Θ')に関する、同時分布の対数logp(xtu,ytu;Θ)の期待値として定義されるQ関数
【００８９】
【数２５】
2015135437.tif 000027
【００９０】
が最大となるようΘを更新する。但し、<・>p(x_tu|y_tu;Θ')（下付添え字のp(x_tu|y_tu;Θ')は、p(xtu|ytu;Θ')を表す）は、事後確率p(xtu|ytu;Θ')を用いて計算した、音源信号xtuに関する期待値を表す。換言すると、モデルパラメータの集合Θが与えられたときの観測信号ytuの集合Yと音源信号xtuの集合Xとの同時分布の対数log(p(X,Y;Θ))の集合Xに関する期待値が最大となるように、共分散行列Φxtと共分散行列Φvtとを更新する。
【００９１】
このEMアルゴリズムにより、パラメータの推定値Θ'は、(b27)の局所最適解に収束することが保証されている。
【００９２】
Mステップにおけるパラメータの更新式を導くため、Q関数の具体形を求める。(b5)を(b28)に代入すれば、
【００９３】
【数２６】
2015135437.tif 000028
【００９４】
ここで、定数項を無視し、なおかつ、<・>p(x_tu|y_tu;Θ')を<・>と略記した。なお、Tr[・]は、行列・の行列の対角成分の和を表す。(b10)より、(b29)における<xtuxtuH>および<(ytu-xtu)(ytu-xtu)H>は、Eステップで計算した(μx|ytu)'および(Φx|yt)'を用いて、下記のように計算できる。
【００９５】
【数２７】
2015135437.tif 000029
【００９６】
【数２８】
2015135437.tif 000030
【００９７】
Q関数(b29)を、各パラメータに関して最大化することにより、更新式を導ける。
【００９８】
まず、Q関数(b29)の信号パワーφxtに関する偏微分を0と置くことにより((b13)参照)、
【００９９】
【数２９】
2015135437.tif 000031
【０１００】
これを信号パワーφxtについて解くと、
【０１０１】
【数３０】
2015135437.tif 000032
【０１０２】
ただし、
【０１０３】
【数３１】
2015135437.tif 000033
【０１０４】
と置いた。次に、信号コヒーレンス行列Bに関する偏微分より、次式を得る((b13)参照)。
【０１０５】
【数３２】
2015135437.tif 000034
【０１０６】
これをBについて解くと、次式を得る。
【０１０７】
【数３３】
2015135437.tif 000035
【０１０８】
拡散性雑音の共分散行列Φvtの更新式は、行列部分空間Vの具体形に依存する。以下では、まず、
【０１０９】
【数３４】
2015135437.tif 000036
【０１１０】
なる性質を持つクラスの行列部分空間Vに対しては、統一的な更新式を導出できることを示す。
【０１１１】
このクラスは、実数値雑音共分散モデルVreal、ブラインド雑音無相関化モデルVBND、空間的無相関雑音モデルVuncorを含む。次に、上記のクラスに含まれない定数雑音コヒーレンスモデルVcohに対する更新式を、別途導出する。
【０１１２】
行列部分空間Vが(b40)を満たすとき、(Φvt)-1∈Vに注意すると、Q関数(b29)の共分散行列Φvtに依存する項は、次のように変形できる。
【０１１３】
【数３５】
2015135437.tif 000037
【０１１４】
ただし、P[・]は、Hの標準的な内積<A,B>:=Tr[AB]に関する、行列部分空間Vへの正射影を表す。(b41)におけるP[・]内の項は、一般に、行列部分空間Vに属する成分と、これに直交する成分との両方を持つが、(Φvt)-1∈Vとの内積により後者は消えるため、(b41)のようになる。(b41)を最大化する共分散行列Φvt∈Vを求めるために、拘束条件Φvt∈Vを一旦忘れて、(b41)をΦvtで偏微分すると、次式を得る。
【０１１５】
【数３６】
2015135437.tif 000038
【０１１６】
ただし、
【０１１７】
【数３７】
2015135437.tif 000039
【０１１８】
である。P[・]の定義より明らかなように、(b42)は、確かにΦvt∈Vを満たす。
【０１１９】
上記では、偏微分を用いて更新式(b42)を導いたが、以下のように、より直観的に導くこともできる。(b41)を-Uで割れば、共分散行列Φvtに依存しない定数項を除いて、下記のマルチチャネル板倉-齋藤ダイバージェンス(参考文献７〜９参照)に一致する。
【０１２０】
【数３８】
2015135437.tif 000040
【０１２１】
したがって、(b41)の最大化は、(b45)の最小化と等価である。DIS(・,・)は非負であり、二つの引数が等しいときに限り0となる。P[^Φvt]は共分散行列Φvtの実行可能領域である行列部分空間Vに属するから、(b45)はΦvt=P[^Φvt]のときに最小となる。
【０１２２】
次に、行列部分空間Vとして定数雑音コヒーレンスモデルを用いる場合、すなわち、
【０１２３】
【数３９】
2015135437.tif 000041
【０１２４】
と表される場合の共分散行列Φvtの更新式を導出する。ここで、φvtは未知の雑音パワーであり、Γは既知と仮定する雑音コヒーレンス行列である。Q関数(b29)に(b46)を代入して、未知の雑音パワーφvtに関して偏微分すると、(b35)の導出と同様に、次式を得る。
【０１２５】
【数４０】
2015135437.tif 000042
【０１２６】
以上で導出したEMアルゴリズムの一回の反復を以下にまとめる。
[Eステップ]
【０１２７】
【数４１】
2015135437.tif 000043
【０１２８】
[Mステップ]
【０１２９】
【数４２】
2015135437.tif 000044
【０１３０】
なお、Mステップにおける、行列部分空間Vへの直交射影は、次式により計算できる。
【０１３１】
【数４３】
2015135437.tif 000045
【０１３２】
ここで、Qdは、行列部分空間Vの正規直交基底であり、DはV行列部分空間の次元である。なお、Qdの具体形については、非特許文献１、参考文献５を参照されたい。
【０１３３】
＜本実施形態に係るモデル推定装置１００＞
図１は第一実施形態に係るモデル推定装置１００の機能ブロック図、図２はその処理フローの例を示す図である。モデル推定装置１００は、周波数領域変換部１１０、事後確率更新部１２０、パラメータ更新部１３０、パラメータ保持部１４０を含む。また、事後確率更新部１２０、パラメータ更新部１３０、パラメータ保持部１４０をまとめて、モデル推定部１５０と称する。
【０１３４】
＜周波数領域変換部１１０＞
周波数領域変換部１１０は、M個のマイクロホンで取得した、時間領域の観測信号~yτ(m)(m=1,…,M)を受け取り、これらの値を用いて、短時間フーリエ変換などの時間周波数変換により、周波数領域の観測信号ytを計算し（ｓ１１０）、事後確率更新部１２０及びパラメータ更新部１３０に出力する。ここで、M>1であり、yの上の〜は時間領域における表現であることを表し、τは時間領域におけるサンプルの番号である。また、前述の通り、フレームは、ブロックの番号tとブロック内でのフレームの番号uとによって指定し、周波数ビンの番号は省略している。
【０１３５】
＜事後確率更新部１２０＞
事後確率更新部１２０は、周波数領域の観測信号ytuを受け取り、パラメータ保持部１４０から現在のパラメータ(要は、EMアルゴリズムの反復処理によって更新する中で、現在のパラメータであることを意味し、「直近に求めたパラメータ」「最新のパラメータ」と言い換えてもよい)の集合の推定値Θ'を取り出す。これらの値を用いて、観測信号ytuが与えられたときの音源信号xtuの事後確率p(xtu|ytu;Θ')を更新し、パラメータ更新部１３０に出力する。上述の通り、ガウス分布の場合、実際には、事後確率の平均μx|ytuおよび共分散行列Φx|ytを更新して出力すれば十分である。
【０１３６】
以下、事後確率更新部１２０における処理を、詳細に説明する。図１のように、事後確率更新部１２０は、事後平均更新手段１２１と、事後共分散行列更新手段１２２とを含む。なお、事後確率更新部１２０における初めての処理に先立ち、パラメータの集合Θの初期値を設定し（ｓ１）、パラメータ保持部１４０に用意しておく。この初期値は、例えば、非特許文献１により推定した、信号のステアリングベクトル^h、音源信号のパワースペクトログラム^φxt、雑音共分散行列^Φvtを用いて、以下のように計算すればよい。
【０１３７】
【数４４】
2015135437.tif 000046
【０１３８】
ここで、(b50)のεIM×Mの加算は、信号コヒーレンス行列Bの初期値が特異行列となってアルゴリズムが不安定になるのを防ぐ役割がある。εは正数で、例えば10-4とすればよく、IM×MはM×Mの単位行列を表す。
【０１３９】
事後平均更新手段１２１は、周波数領域の観測信号ytuと、現在のパラメータの集合の推定値Θ'に含まれる音源信号の共分散行列Φxtと拡散性雑音の共分散行列Φvtとを用いて、事後確率の平均μx|ytuを次式により更新して出力する（ｓ１２１）。
【０１４０】
【数４５】
2015135437.tif 000047
【０１４１】
事後共分散行列更新手段１２２は、現在のパラメータの集合の推定値Θ'に含まれる共分散行列ΦxtとΦvtとを用いて、事後確率の共分散行列Φx|ytを、次式により更新して出力する（ｓ１２２）。
【０１４２】
【数４６】
2015135437.tif 000048
【０１４３】
＜パラメータ更新部１３０＞
パラメータ更新部１３０は、観測信号ytuと、事後確率更新部１２０で計算した事後確率の平均μx|ytuと共分散行列Φx|ytとを受け取る。また、パラメータ保持部１４０から現在のパラメータの集合の推定値Θ'を取り出す。パラメータ更新部１３０は、これらの値を用いて、パラメータの集合Θを更新して、パラメータ保持部１４０に出力する。また、所定の更新回数を終了している場合には（ｓ２）、更新したパラメータの集合Θを最終的なパラメータの推定値として出力する。なお、更新の対象をパラメータの集合Θと呼び、その際に用いる、過去に求めた（更新した）パラメータの集合を推定値Θ'と呼ぶ。以下、パラメータ更新部１３０における処理を、詳細に説明する。図１に示すように、パラメータ更新部１３０は、信号パワー更新手段１３１、信号コヒーレンス行列更新手段１３２、信号共分散行列更新手段１３３、雑音共分散行列更新手段１３４を含む。
【０１４４】
信号パワー更新手段１３１では、事後確率の平均μx|ytu、共分散行列Φx|yt、信号コヒーレンス行列Bを用いて、信号パワーφxtを更新して出力する（ｓ１３１）。まず、
【０１４５】
【数４７】
2015135437.tif 000049
【０１４６】
により、^Φxtを計算する。次に、
【０１４７】
【数４８】
2015135437.tif 000050
【０１４８】
により、信号パワーφxtを更新する。
【０１４９】
信号コヒーレンス行列更新手段１３２は、信号パワー更新手段１３１で計算した^Φxtと信号パワーφxtとを用いて、信号コヒーレンス行列Bを次式により更新して出力する（ｓ１３２）。
【０１５０】
【数４９】
2015135437.tif 000051
【０１５１】
信号共分散行列更新手段１３３は、信号パワー更新手段１３１で計算したφxtと、信号コヒーレンス行列更新手段１３２で計算した信号コヒーレンス行列Bとを用いて、音源信号の共分散行列Φxtを次式により更新して出力する（ｓ１３３）。
【０１５２】
【数５０】
2015135437.tif 000052
【０１５３】
雑音共分散行列更新手段１３４は、事後確率の平均μx|ytuと共分散行列Φx|yt、および観測信号ytuを用いて、拡散性雑音の共分散行列Φvtを更新して出力する（ｓ１３４）。まず、
【０１５４】
【数５１】
2015135437.tif 000053
【０１５５】
により、^Φvtを計算する。次に、
【０１５６】
【数５２】
2015135437.tif 000054
【０１５７】
により、拡散性雑音の共分散行列Φvtを更新して出力する。ここで、P[^Φvt]は、QdをVの正規直交基底として、
【０１５８】
【数５３】
2015135437.tif 000055
【０１５９】
により計算され、コヒーレンス行列Γの(m,n)成分γmnは、
【０１６０】
【数５４】
2015135437.tif 000056
【０１６１】
により計算される。Lmnはm番目とn番目のマイクロホンの距離、cは音速、ｆは観測信号ytuの周波数である。＜拡散性雑音抑圧の問題設定＞で説明したとおり、本実施形態のモデル推定装置には、周波数ビンごとの観測信号が入力され、周波数ビンごとに処理が行われる。つまり、観測信号ytuはある周波数ビンの観測信号であり、ｆはその周波数ビンの周波数を表す。
【０１６２】
パラメータ更新部１３０は、上述の信号パワー更新手段１３１、信号コヒーレンス行列更新手段１３２、信号共分散行列更新手段１３３、雑音共分散行列更新手段１３４の各部の処理により得られたパラメータの集合Θ:={{Φxt}Tt=1,{Φvt}Tt=1}を出力する。
【０１６３】
＜パラメータ保持部１４０＞
パラメータ保持部１４０は、パラメータ更新部１３０での更新処理により得られたパラメータの集合Θを受け取り、保持し（ｓ１４０）、事後確率更新部１２０とパラメータ更新部１３０とにおける次回の処理の際に提供する。
【０１６４】
＜効果＞
このような構成により、従来技術よりも推定精度よく雑音抑圧時に用いるモデルパラメータの推定することができる。
【０１６５】
＜第二実施形態＞
第一実施形態と異なる部分を中心に説明する。
【０１６６】
本実施形態では、第一実施形態で推定したパラメータの集合Θ:={{Φxt}Tt=1,{Φvt}Tt=1}をを用いて、観測信号ytuに含まれる雑音を抑圧する雑音抑圧装置２００について説明する。
【０１６７】
図３は雑音抑圧装置２００の機能ブロック図、図４はその処理フローの例を示す図である。雑音抑圧装置２００は、周波数領域変換部１１０、モデル推定部１５０、平均計算部２６０、時間領域変換部２７０を含む。
【０１６８】
周波数領域変換部１１０、モデル推定部１５０の処理内容は第一実施形態で説明した通りである。
【０１６９】
モデル推定部１５０は、反復終了後のパラメータの集合Θ:={{Φxt}Tt=1,{Φvt}Tt=1}を出力し、平均計算部２６０での処理に供する。
【０１７０】
なお、図４において、事後確率更新部１２０、パラメータ更新部１３０及びパラメータ保持部１４０における処理をまとめて、モデル推定処理（ｓ１５０）という。
【０１７１】
＜平均計算部２６０＞
平均計算部２６０は、周波数領域変換部１１０から出力された観測信号ytuとモデル推定部から出力されたパラメータの集合Θ:={{Φxt}Tt=1,{Φvt}Tt=1}を用いて、式(b52)により、
観測信号ytuが与えられたときの音源信号xtuの事後確率の平均μx|ytuを得て出力する。
【０１７２】
＜時間領域変換部２７０＞
時間領域変換部２７０は、平均計算部２６０の出力である事後確率の平均μx|ytuを受け取り、この値に対し、逆短時間フーリエ変換（inverse STFT）などの時間周波数変換の逆変換を適用し、時間領域における音源信号の推定値である、時間領域の信号^xτ∈RMに変換し（ｓ２７０）、これを雑音抑圧装置２００の出力値として出力する。
【０１７３】
＜効果＞
以上のように、本実施形態では、非特許文献１におけるようなユークリッド距離の平方ではなく、より音響信号に適したマルチチャネル板倉-齋藤ダイバージェンスを用いることにより、対数的な音響信号に、より適した処理を行うことが可能である。さらに、本実施形態は、観測信号の確率的生成モデルに基づくため、確率的生成モデルに基づく他の音声強調手法（例:音源分離(参考文献１及び７参照)、残響除去（参考文献１２参照）との統合に向いている（つまり、第一実施形態で求めたモデルパラメータは拡散性雑音状況下の音源分離、残響除去等の技術に用いることができる）。
（参考文献１２）T. Nakatani, T. Yoshioka, K. Kinoshita, M. Miyoshi, and B.-H. Juang, “Speech dereverberation based on variance-normalized delayed linear prediction”, IEEE Trans. ASLP, Sep. 2010, vol. 18, no. 7, pp. 1717-1731.
このような統合により、様々な音響的事象を含む実環境において、広く適用可能な新しい音声強調の枠組みを構築することができると期待されるため、大きな利点である。
【０１７４】
＜実験結果＞
本実施形態の効果を確かめるため実験を行った。観測信号としては、REVERB challenge(参考文献１３参照)のデータベース中の、AMI_WSJ20-Array1-*_T10c030x.wav（*=1,…,8）を用いた。
(参考文献１３)K. Kinoshita, M. Delcroix, T. Yoshioka, T. Nakatani, E.A.P. Habets,R. Hab-Umbach, V. Leutnant, A. Sehr, W. Kellermann, R. Maas, S. Gannot, and B. Raj, "The REVERB challenge: A common evaluation framework for dereverberation and recognition of reverberant speech", in Proc. WASPAA, Oct. 2013.
【０１７５】
これは、雑音・残響を含む実環境において、一人の話者の音声を、8チャネルのマイクロホンアレイを用いて収録した信号である。データ長は2秒で、サンプリング周波数は16kHzであった。
【０１７６】
本実施形態と非特許文献１のそれぞれにより推定した共分散行列を用いて設計した、マルチチャネルウィーナーフィルタによる雑音抑圧性能を比較した。短時間フーリエ変換および逆短時間フーリエ変換のフレーム長、フレームシフト、窓は、それぞれ、1024点（64ms）、32点（2ms）、Hamming窓とした。また、ブロック長は32フレームとした。アルゴリズムの反復回数は20回とした。雑音抑圧性能の評価指標としては、ケプストラム距離（cepstrum distance:CD）を用いた（定義は参考文献１３参照）。その際、参照信号(雑音・残響をほぼ含まない信号)としては、ヘッドセットによる収録音AMI_WSJ20-Headset1_T10c030x.wavを用いた。ケプストラム距離が小さい程、参照信号に近く、雑音抑圧性能が高いことを示す。図５に、観測信号、参照信号、処理音のスペクトログラムとケプストラム歪み（括弧内の数字）を示す。スペクトログラムの横軸は時間（秒）、縦軸は周波数（kHz）を表す。(a)は観測信号、(b)は非特許文献１による処理音、(c)は本実施形態による処理音、(d)は参照信号である。比較のため、手法間の差異が顕著であった0-4kHzの範囲に限って図示する。(b)の非特許文献１に比べて、(c)の本実施形態による処理音は、より拡散性雑音を効果的に抑圧できたことが分かる。これは、本実施形態による処理音のケプストラム距離3.9dBが、非特許文献１のケプストラム距離4.6dBよりも小さいことからもわかる。
【０１７７】
＜第三実施形態＞
第一実施形態と異なる部分を中心に説明する。なお、必要に応じて、第一実施形態で定義した記号の一部を再定義する。
【０１７８】
本実施形態では、拡散性雑音下での音声強調技術、特に、残響の影響を受けた単一のノイズフリー音と、拡散性雑音とからなるとモデル化される観測信号を用いて、雑音と残響の影響が除去されたクリーン音成分を推定する際に用いるモデルパラメータ（クリーン信号の共分散行列、拡散性雑音の共分散行列、予測係数行列）を推定する。音源の位置は、未知であってもよく、また、それは静止しているとモデル化する。なお、音源がゆるやかに移動する場合には、音源が静止していると近似しうる複数のブロックに観測信号を分割し、ブロックバッチ処理を適用することにより、性能の低下を抑制することができる。
【０１７９】
＜拡散性雑音下での音声強調の問題設定＞
本実施形態では、特に断りのない限り、各信号を短時間フーリエ変換（short-time Fourier transform: STFT）などの時間周波数領域で表現する。フレームの総数をTで表し、フレームの番号をt∈{1,…,T}で表す。
【０１８０】
以下のように変数を定義する。
・yt∈CM:M個のマイクロホンによる観測信号
・rt∈CM:残響の影響を含む、ノイズフリー信号
・xt∈CM:残響の影響を含まない、クリーン信号、すなわち、音源から発せられる直接音と初期反射成分とからなる信号
・vt∈CM:拡散性雑音
・Gk∈CM×M:残響をモデル化する自己回帰過程の予測係数行列（k∈{1,…,K}はタップの番号、Kは予測次数）
【０１８１】
ここで、yt、rt、xt、vtは、M個のマイクロホン位置における各種信号を並べたM次元のベクトルである。
このとき、観測信号のモデルは、数学的には次式により表現できる。
【０１８２】
【数５５】
2015135437.tif 000057
【０１８３】
ここで、遅延Δは、直近のΔブロックを残響除去のための線形予測から除外することにより、残響に起因する時間相関のみを無相関化し、音声固有の相関を白色化するのを防ぐ効果がある（参考文献１２参照）。また、t≦0に対しては、rt=0と定める。第一実施形態も第三実施形態も、観測信号を、非定常の拡散性雑音と拡散性雑音を含まない信号（ノイズフリー信号）との和でモデル化する点は共通するが、第一実施形態では残響を考慮せずに、拡散性雑音と音源信号との和でモデル化するのに対し、本実施形態では残響を考慮して、拡散性雑音とクリーン信号に残響が重畳された信号rt（残響の影響を含む、(c2)参照）との和でモデル化する点が異なる。つまり、第一実施形態ではノイズフリー信号を音源信号と仮定しているのに対し、第三実施形態ではノイズフリー信号がクリーン信号に残響が重畳された信号であると仮定してモデル化している点が異なる。
【０１８４】
次に、統計的信号処理の枠組みで音声強調を行うために、上記のモデルを、確率的生成モデルに拡張する。本実施形態では、クリーン信号xt,ノイズフリー信号rt,拡散性雑音vt,観測信号ytを確率変数とみなし、予測係数行列Gkを確定的な未知パラメータとみなす。{xt}1≦t≦Tが与えられると(c2)により{rt}1≦t≦Tが定まり、この{rt}1≦t≦Tと{vt}1≦t≦Tが与えられると(c1)により{yt}1≦t≦Tが定まる。したがって、(c1)(c2)を生成モデル化することは、{xt}1≦t≦Tと{vt}1≦t≦Tの確率分布を定義することに帰着する。この確率分布として、本実施形態では、簡単のため、以下の仮定を満たすクラスを考える。
・時間的独立性:{xt}1≦t≦Tは独立な系列である。すなわち、任意の相異なるt,t'に対し（t∈{1,…,T}、t'∈{1,…,T}、t≠t'）、xtとxt'とは独立である。また、{vt}1≦t≦Tも独立な系列である。
・相互独立性:2つの系列{xt}1≦t≦T,{vt}1≦t≦Tは、互いに独立である。すなわち、任意のt,t'に対し（t∈{1,…,T}、t'∈{1,…,T}）、xtとvt'とは独立である。
・ガウス性:xtおよびvtは、平均がともに0、共分散行列がΦxtおよびΦvtの複素ガウス分布に従う。
【０１８５】
ここで、平均μ∈CM、共分散行列Σ∈CM×Mの複素ガウス分布の確率密度関数は、次式で与えられる。
【０１８６】
【数５６】
2015135437.tif 000058
【０１８７】
以上の仮定により、観測信号の生成モデルを立てることは、結局、クリーン信号xtおよび拡散性雑音vtの共分散行列ΦxtおよびΦvtのモデル化に帰着する。この詳細については、後述する。
【０１８８】
以上のモデル化では、各フレームで異なる共分散行列Φxt、Φvtを用いていたが、信号の短時間定常性をモデル化するために、モデルを拡張する。各フレームが、フレームシフトをU倍にすることにより得られるU個のサブフレームよりなると仮定する。共分散行列Φxt、Φvtは、フレームの番号tには依存するが、各フレーム内のサブフレームの番号u∈{1,…,U}には依存しないと仮定する。また、各信号は、xtu,rtu,vtu,ytu∈CMのように、tとuを用いて表記し、同一フレームtに属するサブフレームu∈{1,…,U}に対するこれらの信号は、互いに独立で同一の分布に従うもの（independent and identically distributed: i.i.d.)とする。(c1)(c2)は、次式のように書き改められる。
【０１８９】
【数５７】
2015135437.tif 000059
【０１９０】
上で定義した記法を用いて、本実施形態で扱う音声強調の問題を改めて述べると、観測信号Y:={ytu}1≦t≦T,1≦u≦Uが与えられたときに、クリーン信号X:={xtu}1≦t≦T,1≦u≦Uを推定する問題である。特に、そこでの主要な問題は、上記のモデルのパラメータ
Θ:={Θx,Θv,Θg}   (c6)
Θx:={Φxt}1≦t≦T    (c7)
Θv:={Φvt}1≦t≦T    (c8)
Θg:={Gk}1≦k≦K    (c9)
の推定であり、音声強調技術の性能は、その推定精度に大きく左右される。
【０１９１】
＜信号と拡散性雑音の共分散行列のモデル化＞
上述の通り、観測信号の生成モデルを立てることが、共分散行列ΦxtとΦvtのモデル化に帰着する。非特許文献１では、信号と拡散性雑音の空間的な性質が利用されている。
【０１９２】
クリーン音が、静止した単一の点音源から発せられるという仮定の下では、クリーン音の空間的な特性は、時間によらず一定である。ゆえに、クリーン信号のモデルは、第一実施形態と同じく(b13)を用いる。第一実施形態と同様に、本実施形態においても信号コヒーレンス行列Bは非特許文献１で用いているランク１行列であるという制約を課さない。
【０１９３】
拡散性雑音については、第一実施形態と同じく(b14)によりモデル化する。
【０１９４】
Vの具体形は、アレイ配置や雑音場に関する付加的な仮定などによって、複数提案されている。第一実施形態で説明したように、以下の四つのモデルが考えられる。
・実数値雑音共分散（real-valued noise covariance）モデル
【０１９５】
【数５８】
2015135437.tif 000060
【０１９６】
・ブラインド雑音無相関化（blind noise decorrelation: BND）モデル
【０１９７】
【数５９】
2015135437.tif 000061
【０１９８】
・空間的無相関雑音（spatially uncorrelated noise）モデル
【０１９９】
【数６０】
2015135437.tif 000062
【０２００】
・定数雑音コヒーレンス（fixed noise coherence）モデル
【０２０１】
【数６１】
2015135437.tif 000063
【０２０２】
ここで、Bvは雑音コヒーレンス行列であり、その(m,n)成分は、
【０２０３】
【数６２】
2015135437.tif 000064
【０２０４】
で与えられる。
【０２０５】
これらの雑音モデルの比較については、参考文献５を参照されたい。本実施形態の方法は、これらの全ての雑音モデルに適用可能である。
【０２０６】
＜本実施形態のポイント＞
上記の目的を実現するために、本実施形態では、最尤法によりモデルパラメータ（クリーン信号の共分散行列、拡散性雑音の共分散行列、予測係数行列）を推定する。この最尤法は、ユークリッド距離に基づく尺度と比べ、音声により適した尺度である、板倉-齋藤行列ダイバージェンス（参考文献７〜９）の最小化と等価である。すでに説明した通り、非特許文献１は残響が存在しない場合、すなわち本実施形態のモデルにおいてGk=0とした場合において、ユークリッド距離に基づいてモデルパラメータを推定するものであった。比較のため、非特許文献１と同じ条件であるGk=0の場合について述べると、本実施形態における最尤推定は、次式のコスト関数の最小化と等価である。
【０２０７】
【数６３】
2015135437.tif 000065
【０２０８】
ここで、DISは、次式で定義される。
【０２０９】
【数６４】
2015135437.tif 000066
【０２１０】
スカラーk>0に対して、DIS(kA,kB)=DIS(A,B)が成り立つことから、DISは、上記の音響信号の対数性を反映した、適切な尺度となっていることが分かる。
【０２１１】
Gk≠0の場合にも、同様のことが言える。なお、(c19)は、NTTの板倉・齋藤により提案され（参考文献１０参照）、音声処理における有効性が示されてきた、板倉-齋藤ダイバージェンスを、複数チャネルの場合へ拡張した尺度とみなせるため、本実施形態では、板倉-齋藤行列ダイバージェンスと呼ぶ。
【０２１２】
＜原理＞
＜目的関数＞
本実施形態では、パラメータの集合Θ:={{φxt}1≦t≦T,Bx,{φvt}1≦t≦T,{Gk}1≦k≦K}を推定するために、最尤法に基づき、次式の対数尤度を最大化する。
【０２１３】
【数６５】
2015135437.tif 000067
【０２１４】
すなわち、拘束条件(b14)に注意して、
【０２１５】
【数６６】
2015135437.tif 000068
【０２１６】
によりΘを推定する。
【０２１７】
前述の観測信号生成モデルに基づき、対数尤度(c20)の具体的な表式を導出する。{ytu}1≦t≦Tは、残響に起因する時間相関を持つため、独立な系列ではない。
【０２１８】
そこで、{ytu}1≦t≦Tを並べて拡張したベクトル
【０２１９】
【数６７】
2015135437.tif 000069
【０２２０】
の確率密度関数を導出することにする。まず、(c1)(c2)のモデルを、拡張されたベクトルを用いて書き直すと、次式のようになる。
~yu=~ru+~vu    (c23)
FH~ru=~xu    (c24)
ここで、~ru,~vu,~xuは、~yuと同様に定義する。つまり、
~ru:=[rTuT…r1uT]T∈CMT×1
~vu:=[vTuT…v1uT]T∈CMT×1
~xu:=[xTuT…x1uT]T∈CMT×1
である。また、F∈CMT×MTは、ブロックテプリッツ行列であって、T2個のM×M行列をブロックに持ち、(i,j)番目のブロックは次式で与えられる。
【０２２１】
【数６８】
2015135437.tif 000070
【０２２２】
ただし、Iは単位行列を表す。つまり、以下のように表される。
【０２２３】
【数６９】
2015135437.tif 000071
【０２２４】
(c24)より、FHは雑音除去オペレータ、F-H:=(FH)-1は残響付加オペレータと解釈できる。つまり、F-Hはクリーン信号に重畳される残響の特徴を表すオペレータ（パラメータ）ともいえる。ここで、FHは対角要素がすべて1の上三角行列であり、したがってdet(FH)=1であるから、FHは正則であることに注意する。次に、拡張されたベクトル~xu,~vuの確率密度関数は、前述の観測信号生成モデル、および~xu,~vuの定義より、次式で与えられる。
p(~xu;Θ)=NC(~xu;0,~Φx)    (c26)
p(~vu;Θ)=NC(~vu;0,~Φv)    (c27)
ここで、~Φx∈CMT×MTは、ブロック対角行列であって、T2個のM×M行列をブロックに持ち、i∈{1,…,T}番目の対角ブロックはΦxT-i+1に等しい。~Φvも同様に定義される。つまり、以下のように表される。
【０２２５】
【数７０】
2015135437.tif 000072
【０２２６】
よって、観測モデル(c23)(c24)、および~xu,~vuの確率密度関数(c26)(c27)より、対数尤度は、
【０２２７】
【数７１】
2015135437.tif 000073
【０２２８】
である。ここで、前述の観測信号生成モデルより、{ytu}1≦t≦Tは独立な系列ではないが、{~yu}1≦u≦Uは独立な系列であることに注意する。
【０２２９】
なお、(c28)は、次式のように変形できる。
【０２３０】
【数７２】
2015135437.tif 000074
【０２３１】
ここで、DISは前述の板倉-齋藤ダイバージェンスであり、
【０２３２】
【数７３】
2015135437.tif 000075
【０２３３】
である。したがって、J(Θ)=logp(Y;Θ)の最大化は、DIS(~Φy,F-H~ΦxF-1+~Φv)の最小化と等価である。よって、本実施形態では、~Φyと、F-H~ΦxF-1+~Φvとの板倉-齋藤距離が小さくなるように、共分散行列Φxtと共分散行列Φvtと予測係数行列Gkとを更新する。特に、Gk=0、すなわちF=Iの場合には、前述のように(c18)の最小化に帰着する。
【０２３４】
つまり、本実施形態においても、第一実施形態と同様に、観測信号yｔを拡散性雑音vtと拡散性雑音を含まない信号であるノイズフリー信号との和でモデル化する。ただし、本実施形態におけるノイズフリー信号は、音源から発せられる直接音と初期反射成分とからなるクリーン信号xtに残響が重畳された信号である。そして、本実施形態では、観測信号ytから得られる共分散行列~Φｙと、拡散性雑音の共分散行列~Φvとノイズフリー信号の共分散行列F-H~ΦxF-1との和（F-H~ΦxF-1+~Φv）と、の板倉-齋藤距離が小さくなるように各モデルパラメータ（ここではΦvtとΦxtとGk）を更新する。上述の通り、~Φv はΦvtにより特定される共分散行列なので、拡散性雑音の共分散行列を特定するモデルパラメータはΦvtである。また、~Φx はΦｘt により特定される共分散行列であり、F-HとF-1は予測係数行列Gkにより特定されるオペレータであるので、拡散性雑音を含まない信号(ノイズフリー信号)を特定するモデルパラメータはΦxtとGkである。なお、~Φvは観測信号を非定常の拡散性雑音とノイズフリー信号との和でモデル化したときの拡散性雑音の統計的性質を表しており、F-H~ΦxF-1はノイズフリー信号の統計的性質を表しているといえる。
【０２３５】
＜Q関数＞
R:={rtu}1≦t≦T,1≦u≦Uを隠れ変数とみなし、expectation-maximization(EM)アルゴリズムにより（参考文献１１）、パラメータの集合Θの更新式を導くことができる。パラメータの集合Θの現在の推定値をΘ'で表すと、Q関数は、
【０２３６】
【数７４】
2015135437.tif 000076
【０２３７】
で定義される。これを計算するために、まず、完全データの対数尤度関数を計算すると、
【０２３８】
【数７５】
2015135437.tif 000077
【０２３９】
次に、隠れ変数の事後確率を計算する。上式と、ベイズの定理により、
【０２４０】
【数７６】
2015135437.tif 000078
【０２４１】
となる。したがって、
【０２４２】
【数７７】
2015135437.tif 000079
【０２４３】
である。ただし、
【０２４４】
【数７８】
2015135437.tif 000080
【０２４５】
である。次式により、μr|ytuおよびΦr|ytt'を定義する。
【０２４６】
【数７９】
2015135437.tif 000081
【０２４７】
(c35)および(c39)より、Q関数は
【０２４８】
【数８０】
2015135437.tif 000082
【０２４９】
と書ける。ここで、
【０２５０】
【数８１】
2015135437.tif 000083
【０２５１】
と定義した。
【０２５２】
＜更新式の導出＞
式(c45)のQ関数を、各パラメータに関して最大化することにより、Mステップにおける更新式が得られる。信号パワーφxtに関する偏微分より、
【０２５３】
【数８２】
2015135437.tif 000084
【０２５４】
よって、
【０２５５】
【数８３】
2015135437.tif 000085
【０２５６】
次に、信号コヒーレンス行列Bxに関する偏微分より、
【０２５７】
【数８４】
2015135437.tif 000086
【０２５８】
よって、
【０２５９】
【数８５】
2015135437.tif 000087
【０２６０】
拡散性雑音の共分散行列Φvtの更新式については、雑音抑圧技術（第一実施形態参照）の場合と同様に、
【０２６１】
【数８６】
2015135437.tif 000088
【０２６２】
を得る。なお、第一実施形態のように、「V=VBND,Vreal,Vuncor」を「Vが
【０２６３】
【数８７】
2015135437.tif 000089
【０２６４】
なる性質を持つ場合」に置き換えてもよい。最後に、^Φxtが~Gを含むことに注意して、~G*に関する偏微分より（・*は・の複素共役を表す）、次式を得る。
~G=D-1N    (c58)
ただし、DおよびNは、次式により定義される。
【０２６５】
【数８８】
2015135437.tif 000090
【０２６６】
以上で導出したEMアルゴリズムの一回の反復を、以下にまとめる。
[Eステップ]
【０２６７】
【数８９】
2015135437.tif 000091
【０２６８】
[Mステップ]
【０２６９】
【数９０】
2015135437.tif 000092
【０２７０】
＜本実施形態に係るモデル推定装置３００＞
図６は第三実施形態に係るモデル推定装置３００の機能ブロック図、図７はその処理フローの例を示す図である。
【０２７１】
モデル推定装置３００は、モデル推定装置の機能構成例を示す。モデル推定装置３００は、周波数領域変換部１１０、事後確率更新部３２０、パラメータ更新部３３０及びパラメータ保持部１４０を含む。また、事後確率更新部３２０、パラメータ更新部３３０、パラメータ保持部１４０をまとめて、モデル推定部３５０と称する。図７を参照して、モデル推定装置３００の動作例を、手続きの順に従って説明する。周波数領域変換部１１０、及びパラメータ保持部１４０の処理内容は第一実施形態と同様である。
【０２７２】
＜事後確率更新部３２０＞
事後確率更新部３２０は、周波数領域変換部１１０で計算した、周波数領域の観測信号ytuを受け取り、パラメータ保持部１４０に保持されている、現在のパラメータの集合の推定値Θ'を取り出し、これらの値を用いて、観測信号ytuが与えられたときのノイズフリー信号rtuの事後確率の平均μr|ytuおよび共分散行列Φr|ytt'を更新し、更新した値をパラメータ更新部３３０に出力する。
【０２７３】
以下、事後確率更新部３２０における処理を、詳細に説明する。事後確率更新部３２０は、事後平均更新手段３２１と、事後共分散行列更新手段３２２とを含む。
【０２７４】
なお、事後確率更新部３２０における初めての処理に先立ち、パラメータの集合Θの初期値を設定し（ｓ１）、パラメータ保持部１４０に保持しておく。この初期値は、例えば、非特許文献１により推定した、信号のステアリングベクトル^h、信号のパワースペクトログラム^φxt、雑音共分散行列^Φvtを用いて、以下のように計算すればよい。
【０２７５】
【数９１】
2015135437.tif 000093
【０２７６】
ここで、εは正数（例えば10-4）であり、Bxの初期値が特異行列となってアルゴリズムが不安定になるのを防ぐ役割がある。
【０２７７】
事後平均更新手段３２１は、周波数領域の観測信号ytuと、現在のパラメータの集合の推定値Θ'とを用いて、事後確率の平均μr|ytuを更新して出力する（ｓ３２１）。まず、定義に従って、Fおよび~Φxを更新する（(c25)参照、~Φx∈CMT×MTは、ブロック対角行列であって、T2個のM×M行列をブロックに持ち、i∈{1,…,T}番目の対角ブロックはΦxT-i+1に等しい）。次に、
【０２７８】
【数９２】
2015135437.tif 000094
【０２７９】
により~Φrを更新する。次に、(c57)によりΦvtを求める。上述の通り
【０２８０】
【数９３】
2015135437.tif 000095
【０２８１】
なので、この処理は~Φvを更新することに相当する。つまり、(c57)によりΦvtを求めることで~Φvを更新する。次に、~Φrと~Φvを用いて
【０２８２】
【数９４】
2015135437.tif 000096
【０２８３】
により~μr|yuを更新して出力する。ここで、
【０２８４】
【数９５】
2015135437.tif 000097
【０２８５】
である。 (c43)の定義より、~μr|yuはμr|ytuを統合した行列であるので、~μr|yuを更新することはμr|ytuも同時に更新することを意味する。
【０２８６】
事後共分散行列更新手段３２２では、~Φrと~Φvとを用いて、
【０２８７】
【数９６】
2015135437.tif 000098
【０２８８】
により~Φr|yを更新する。(c44)の定義より、~Φr|yはΦr|ytt'を統合した行列であるので、~Φr|yを更新することはΦr|ytt'も同時に更新することを意味する（ｓ３２２）。
【０２８９】
＜パラメータ更新部３３０＞
パラメータ更新部３３０は、事後確率更新部３２０で計算した事後確率の平均μr|ytuと共分散行列Φr|ytt'、ならびに、観測信号ytuを受け取り、パラメータ保持部１４０から保持されている現在のパラメータの集合の推定値Θ'を取り出し、これらの値を用いて、パラメータの集合Θを更新し、パラメータの集合Θをパラメータ保持部１４０に格納する。また、所定の更新回数を終了している場合には（ｓ２）、更新したパラメータの集合Θを最終的なパラメータの推定値として出力する。以下、パラメータ更新部３３０における処理を、詳細に説明する。
【０２９０】
図６のように、パラメータ更新部３３０は、信号パワー更新手段３３１、信号コヒーレンス行列更新手段３３２、信号共分散行列更新手段３３３、雑音共分散行列更新手段３３４及び予測係数行列更新手段３３５を含む。
【０２９１】
信号パワー更新手段３３１では、事後確率更新部３２０から出力された事後確率の平均μr|ytuと共分散行列Φr|ytt'、および、パラメータ集合の現在の推定値Θ'を用いて、信号パワーφxtを計算して出力する（ｓ３３１）。まず、
【０２９２】
【数９７】
2015135437.tif 000099
【０２９３】
により、Ψtt'を計算する。次に、
【０２９４】
【数９８】
2015135437.tif 000100
【０２９５】
により^Φxtを更新して出力する。最後に、
【０２９６】
【数９９】
2015135437.tif 000101
【０２９７】
により、信号パワーφxtを更新して出力する。
【０２９８】
信号コヒーレンス行列更新手段３３２では、信号パワー更新手段３３１で計算した^Φxtと信号パワーφxtとを用いて、信号コヒーレンス行列Bxを、次式により更新して出力する（ｓ３３２）。
【０２９９】
【数１００】
2015135437.tif 000102
【０３００】
信号共分散行列更新手段３３３では、信号パワー更新手段３３１で計算した信号パワーφxtと、信号コヒーレンス行列更新手段３３２で計算した信号コヒーレンス行列Bxとを用いて、クリーン信号の共分散行列Φxtを次式により更新して出力する（ｓ３３３）。
【０３０１】
【数１０１】
2015135437.tif 000103
【０３０２】
雑音共分散行列更新手段３３４では、事後確率更新部３２０から出力された事後確率の平均μr|ytuと共分散行列Φr|ytt'、パラメータの集合の現在の推定値Θ'、および観測信号ytuを用いて、拡散性雑音の共分散行列Φvtを更新して出力する（ｓ３３４）。まず、
【０３０３】
【数１０２】
2015135437.tif 000104
【０３０４】
により、^Φvtを計算する。次に、
【０３０５】
【数１０３】
2015135437.tif 000105
【０３０６】
により、拡散性雑音の共分散行列Φvtを更新する。ここで、P[^Φvt]は、{Qd}1≦d≦DをVの正規直交基底として、
【０３０７】
【数１０４】
2015135437.tif 000106
【０３０８】
により計算され、信号コヒーレンス行列Bvの各成分bvmnは、
【０３０９】
【数１０５】
2015135437.tif 000107
【０３１０】
により計算される。
【０３１１】
予測係数行列更新手段３３５では、信号パワー更新手段３３１で計算したΨtt'および信号パワーφxtを用いて、予測係数行列Gkを更新して出力する（ｓ３３５）。まず、
【０３１２】
【数１０６】
2015135437.tif 000108
【０３１３】
によりDを更新する。次に、
【０３１４】
【数１０７】
2015135437.tif 000109
【０３１５】
によりNを更新する。そして、
~G←D-1N    (c80)
により~Gを更新する。(c52)の定義より、~GはGkを統合した行列なので、~Gを更新することはGkも同時に更新することを意味する。
【０３１６】
＜効果＞
このような構成により、従来技術よりも推定精度よく音声強調時に用いるモデルパラメータの推定することができる。
【０３１７】
＜第四実施形態＞
第三実施形態と異なる部分を中心に説明する。本実施形態では、第三実施形態に係るモデル推定装置によって推定されたモデルパラメータを用いて、雑音と残響の影響が除去されたクリーン信号を推定する音声強調装置４００について説明する。
【０３１８】
図８は第四実施形態に係る音声強調装置４００の機能ブロック図、図９はその処理フローの例を示す図である。音声強調装置４００は、周波数領域変換部１１０、モデル推定部３５０、平均計算部４５０、強調音声生成部４６０、時間領域変換部４７０を含む。
【０３１９】
なお、図９において、事後確率更新部３２０、パラメータ更新部３３０及びパラメータ保持部１４０における処理をまとめて、モデル推定処理（ｓ３５０）という。
【０３２０】
周波数領域変換部１１０及びモデル推定部３５０における処理は第三実施形態で説明した通りである。
【０３２１】
＜平均計算部４５０＞
平均計算部４５０は、周波数領域変換部１１０から出力された観測信号ytuとモデル推定部１５０から出力されたパラメータの集合Θを受け取り、式（c66）及び式(c43)により、観測信号ytuが与えられたときのノイズフリー信号rtuの事後確率の平均μr|ytuを得て出力する（ｓ４５０）。
【０３２２】
＜強調音声生成部４６０＞
強調音声生成部４６０は、平均計算部４５０から出力された事後確率の平均μr|ytuおよびモデル推定部１５０から出力されたパラメータの集合Θを受け取り、これらの値を用いて、次式により、強調音声^xtuを生成（計算）し（ｓ４６０）、時間領域変換部４７０に出力する。
【０３２３】
【数１０８】
2015135437.tif 000110
【０３２４】
ここで、FHは、
【０３２５】
【数１０９】
2015135437.tif 000111
【０３２６】
により、パラメータの集合Θに含まれる予測係数行列Gkを用いて構成される行列である。なお、強調音声^xtuは、観測信号ytuが与えられたときのクリーン信号の事後確率の平均に相当する。
【０３２７】
＜時間領域変換部４７０＞
最後に、時間領域変換部４７０は、強調音声^xtuに対し、逆短時間フーリエ変換（inverse STFT）などの時間周波数変換の逆変換を適用し、時間領域におけるクリーン信号の推定値である、時間領域の信号^xτ∈RMに変換し（ｓ４７０）、音声強調装置４００の出力値として出力する。
【０３２８】
＜効果＞
以上のように、本実施形態では、非特許文献１におけるようなユークリッド距離の平方ではなく、より音響信号に適した板倉-齋藤行列ダイバージェンスを用いることにより、対数的な音響信号に、より適した処理を行うことが可能である。さらに、本実施形態は、観測信号の確率的生成モデルに基づくため、確率的生成モデルに基づく他の音声強調手法（例:音源分離、参考文献１、７参照）との統合に向いている。このような統合により、様々な音響的事象を含む実環境において、広く適用可能な新しい音声強調の枠組みを構築することができると期待されるため、大きな利点である。
【０３２９】
＜実験結果＞
本実施形態の効果を確認するために、本実施形態の音声強調方法による、雑音および残響抑圧実験をおこなった。
【０３３０】
観測信号としては、参考文献１３記載のREVERB challengeのReal Dataデータベースに含まれる、8チャネルの実環境収録音AMIWSJ20-Array1-*T10c020c.wav（*=1,…,8はマイクロホンの番号）から、長さ2秒の区間を切り出したものを用いた。
【０３３１】
この観測信号は、残響時間RT60〜0.7sの残響、および雑音を含み、サンプリング周波数は16kHzであった。
【０３３２】
比較手法としては、第四実施形態の音声強調方法において、Eステップの1〜4を、トリビアルな更新則μr|ytu←ytuおよびΦr|ytt'←0で置き換えて、雑音モデルを無効にしたアルゴリズムを用いた。このアルゴリズムは、信号の空間相関をモデル化し、全てのチャネルを用いて信号分散を計算する点を除けば、参考文献１２の周波数領域NDLP法（variance-normalized delayed linear prediction)と等価である。
【０３３３】
評価指標としては、参考文献１３にならい、ケプストラム距離（cepstrum distance:CD）、対数尤度比（log-likelihood ratio: LLR)、周波数重み付き区分信号対雑音比frequency-weighted segmental signal-to-noise ratio: FWSegSNR）、および信号対残響変調エネルギー比（speech-to-reverberation modulation energy ratio: SRMR）を用いた。これらの指標の計算のための参照信号として、ヘッドセットによる録音AMI WSJ20-Headset1 T10c020c.wav を用いた。FWSegSNR とSRMRが大きいほど、また、CDとLLRが小さいほど、音声強調性能が高いとみなす。
【０３３４】
分析条件は、フレーム長を1024点、フレームシフトを256点（それぞれ64ms,16msに相当）、窓をハミング窓とした。予測次数はK=3、予測遅延はΔ=3、反復数は20とした。
【０３３５】
図１０にシミュレーション結果のスペクトログラムを示す。スペクトログラムの横軸は時間（秒）、縦軸は周波数（kHz）を表す。(a)は観測信号、(b)は比較手法による処理音、(c)は第四実施形態の音声強調法で、残響と雑音の両者を考慮したモデルを学習するが、マルチチャネルウィーナーフィルタによる雑音抑圧のみを行ったもの、(d)は第四実施形態の音声強調法による、雑音と残響の同時抑圧を行ったもの、(e)は参照信号である。差異が最も明白な0-4kHzの範囲を示した。(a)の観測信号に比べ、(b)の比較手法は、残響を抑圧できており、特に1sおよび1.7s付近で顕著であるが、雑音はほとんど抑圧できていない。(c)は、(a)に対して雑音を抑圧できており、特に信号の最初の非音声区間で顕著であるが、残響はほとんど抑圧できていない。(d)は、雑音と残響を効果的に抑圧できている。
【０３３６】
図１１に評価指標による客観評価の結果を示す。ラベル(a)-(d)は、図１０で定義した通りである。太字および斜体は、各指標について一番目および二番目に良い数値を表す。第四実施形態の音声強調方法により、CD,FWSegSNR,SRMRでは最も良い結果、LLRでも2番目に良い結果がえられた。これらの結果は、第四実施形態の音声強調方法の有効性を示している。
【０３３７】
＜第五実施形態＞
第四実施形態と異なる部分を中心に説明する。
【０３３８】
本実施形態では、第三実施形態に係るモデル推定装置によって推定されたモデルパラメータを用いて、雑音を除去したノイズフリー信号rtuを推定する雑音抑圧装置５００について説明する。
【０３３９】
図１２は第五実施形態に係る雑音抑圧装置５００の機能ブロック図、図１３はその処理フローの例を示す図である。雑音抑圧装置５００は、周波数領域変換部１１０、モデル推定部３５０、平均計算部４５０、時間領域変換部５７０を含む。
【０３４０】
周波数領域変換部１１０及びモデル推定部３５０及び平均計算部４５０における処理は第四実施形態で説明した通りである。
【０３４１】
＜時間領域変換部５７０＞
時間領域変換部５７０は、平均計算部４５０の出力である事後確率の平均μr|ytuを受け取り、この値に対し、逆短時間フーリエ変換（inverse STFT）などの時間周波数変換の逆変換を適用し、時間領域におけるノイズフリー信号の推定値である、時間領域の信号^rτ∈RMに変換し（ｓ５７０）、これを雑音抑圧装置５００の出力値として出力する。
【０３４２】
＜その他の変形例＞
本発明は上記の実施形態及び変形例に限定されるものではない。例えば、上述の各種の処理は、記載に従って時系列に実行されるのみならず、処理を実行する装置の処理能力あるいは必要に応じて並列的にあるいは個別に実行されてもよい。その他、本発明の趣旨を逸脱しない範囲で適宜変更が可能である。
【０３４３】
＜プログラム及び記録媒体＞
また、上記の実施形態及び変形例で説明した各装置における各種の処理機能をコンピュータによって実現してもよい。その場合、各装置が有すべき機能の処理内容はプログラムによって記述される。そして、このプログラムをコンピュータで実行することにより、上記各装置における各種の処理機能がコンピュータ上で実現される。
【０３４４】
この処理内容を記述したプログラムは、コンピュータで読み取り可能な記録媒体に記録しておくことができる。コンピュータで読み取り可能な記録媒体としては、例えば、磁気記録装置、光ディスク、光磁気記録媒体、半導体メモリ等どのようなものでもよい。
【０３４５】
また、このプログラムの流通は、例えば、そのプログラムを記録したＤＶＤ、ＣＤ−ＲＯＭ等の可搬型記録媒体を販売、譲渡、貸与等することによって行う。さらに、このプログラムをサーバコンピュータの記憶装置に格納しておき、ネットワークを介して、サーバコンピュータから他のコンピュータにそのプログラムを転送することにより、このプログラムを流通させてもよい。
【０３４６】
このようなプログラムを実行するコンピュータは、例えば、まず、可搬型記録媒体に記録されたプログラムもしくはサーバコンピュータから転送されたプログラムを、一旦、自己の記憶部に格納する。そして、処理の実行時、このコンピュータは、自己の記憶部に格納されたプログラムを読み取り、読み取ったプログラムに従った処理を実行する。また、このプログラムの別の実施形態として、コンピュータが可搬型記録媒体から直接プログラムを読み取り、そのプログラムに従った処理を実行することとしてもよい。さらに、このコンピュータにサーバコンピュータからプログラムが転送されるたびに、逐次、受け取ったプログラムに従った処理を実行することとしてもよい。また、サーバコンピュータから、このコンピュータへのプログラムの転送は行わず、その実行指示と結果取得のみによって処理機能を実現する、いわゆるＡＳＰ（Application Service Provider）型のサービスによって、上述の処理を実行する構成としてもよい。なお、プログラムには、電子計算機による処理の用に供する情報であってプログラムに準ずるもの（コンピュータに対する直接の指令ではないがコンピュータの処理を規定する性質を有するデータ等）を含むものとする。
【０３４７】
また、コンピュータ上で所定のプログラムを実行させることにより、各装置を構成することとしたが、これらの処理内容の少なくとも一部をハードウェア的に実現することとしてもよい。
【図１】
2015135437.tif 000140
【図２】
2015135437.tif 000141
【図３】
2015135437.tif 000142
【図４】
2015135437.tif 000143
【図６】
2015135437.tif 000144
【図７】
2015135437.tif 000145
【図８】
2015135437.tif 000146
【図９】
2015135437.tif 000147
【図１１】
2015135437.tif 000148
【図１２】
2015135437.tif 000149
【図１３】
2015135437.tif 000150
【図５】
2015135437.tif 000151
【図１０】
2015135437.tif 000152