JP5406866B2 - 音源分離装置、その方法及びプログラム - Google Patents

音源分離装置、その方法及びプログラム Download PDF

Info

Publication number
JP5406866B2
JP5406866B2 JP2011036559A JP2011036559A JP5406866B2 JP 5406866 B2 JP5406866 B2 JP 5406866B2 JP 2011036559 A JP2011036559 A JP 2011036559A JP 2011036559 A JP2011036559 A JP 2011036559A JP 5406866 B2 JP5406866 B2 JP 5406866B2
Authority
JP
Japan
Prior art keywords
sound source
signal
occupancy
estimating
separation
Prior art date
Legal status (The legal status is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the status listed.)
Active
Application number
JP2011036559A
Other languages
English (en)
Other versions
JP2012173584A (ja
Inventor
章子 荒木
智広 中谷
Current Assignee (The listed assignees may be inaccurate. Google has not performed a legal analysis and makes no representation or warranty as to the accuracy of the list.)
Nippon Telegraph and Telephone Corp
Original Assignee
Nippon Telegraph and Telephone Corp
Priority date (The priority date is an assumption and is not a legal conclusion. Google has not performed a legal analysis and makes no representation as to the accuracy of the date listed.)
Filing date
Publication date
Application filed by Nippon Telegraph and Telephone Corp filed Critical Nippon Telegraph and Telephone Corp
Priority to JP2011036559A priority Critical patent/JP5406866B2/ja
Publication of JP2012173584A publication Critical patent/JP2012173584A/ja
Application granted granted Critical
Publication of JP5406866B2 publication Critical patent/JP5406866B2/ja
Active legal-status Critical Current
Anticipated expiration legal-status Critical

Links

Images

Landscapes

  • Obtaining Desirable Characteristics In Audible-Bandwidth Transducers (AREA)
  • Circuit For Audible Band Transducer (AREA)

Description

本発明は信号処理の技術分野に関し、特に1個以上の音源信号がノイズとともに混在している観測信号から、それぞれの音源に対応する分離信号を推定する音源分離技術に関する。
マルチチャネルウィーナーフィルタを用いた音源分離装置が従来技術として知られている(図1参照)。この従来技術の詳細は、例えば、非特許文献1に記載されている。以下、概略を説明する。K個の音源10(k)(但し、kは音源のインデックス番号であり、k=1,…,K)それぞれから同時に発せられる音源信号s(k)(t)を2個の収音手段(第一収音手段2及び第二収音手段2)で、ある収録時間内(例えば、5秒間)に観測する。この観測状況を状況Xとする。この状況Xの下、第一収音手段2で収音された観測信号をx(t)とし、第二収音手段2で収音された観測信号をx(t)とし、観測された観測信号をX(t)=[x(t),x(t)]とする。収音手段とは例えばマイクロホンのことであり、「」は行列の転置を表し、tを離散時刻とし、t=1,…,Tとする。周波数領域変換部11は、観測信号X(t)を周波数領域に変換することで、周波数毎の時系列信号である観測信号ベクトルX(n,f)=[x(n,f),x(n,f)]に変換する。以降、断りのない場合、観測信号とは、周波数毎の時系列信号である観測信号ベクトルX(n,f)=[x(n,f),x(n,f)]を意味し、時間領域の観測信号の場合、それを明記する。周波数領域への変換は、例えば短時間フーリエ変換を用いれば良い。なお、フレーム数をNとし、nはフレームに対応する時刻を表し、n=1,2,…,Nとする。また、fはサンプリング周波数fをF等分した離散点(周波数ビン)を表す。例えば、f∈{0,(1/F)f,…,((F−1)/F)f}である。
ここで、観測信号ベクトルは、以下の式(1)で表されると仮定する。
Figure 0005406866
但し、c(k)(n,f)=[c(1,k)(n,f),c(2,k)(n,f)]はk番目の音源信号s(k)(n,f)のソースイメージを表し、各成分c(1,k)(n,f),c(2,k)(n,f)は、それぞれ音源信号s(k)(n,f)の第一収音手段2、第二収音手段2における観測値を表す。つまり、ソースイメージc(k)(n,f)は、各収音手段における観測値x(n,f),x(n,f)中の音源信号s(k)(n,f)に基づく信号である。また、h(k)=[h(1,k),h(2,k)であり、各成分h(1,k),h(2,k)は、それぞれ音源10(k)から第一収音手段2、第二収音手段2までの周波数応答を示し、s(k)(n,f)は音源信号s(k)(t)を周波数領域で表現した信号であり、以降、断りのない場合、音源信号とは、周波数毎の時系列信号であるs(k)(n,f)を意味し、時間領域の音源信号の場合、それを明記する。
音源信号の観測時間内においては、音源10(k)、第一収音手段2及び第二収音手段2は固定されており、また、K個の音源10(k)は全て、異なる位置に配置されているとする。すなわち、周波数応答h(k)は時刻nに依存せずに、kの値によって異なる値をとるものと仮定する。
従来技術では、ソースイメージc(k)(n,f)の相関行列
Figure 0005406866
が音源信号s(k)(n,f)の時間周波数毎の分散v(k)(n,f)と、時刻nによらない周波数毎のM行M列の空間相関行列B(k)(f)を用いて、
(k)(n,f)=v(k)(n,f)B(k)(f) (2)
(k)(f)=h(k)(f)(h(k)(f))
とモデル化できると仮定する。但し、「」はエルミート転置を表す。そして、空間相関行列B(k)(f)のクラスタリングを行なうことで、分散v(k)(n,f)と空間相関行列B(k)(f)を推定する(詳細は非特許文献1参照)。さらに、推定した分散v(k)(n,f)と空間相関行列B(k)(f)に基づくマルチチャネルウィーナーフィルタW(k)(n,f)を用いて、音源分離を達成する。すなわち、空間相関行列B(k)が、音源10(k)の位置毎に固有の値を持つことを利用して音源分離を行なう。具体的には、図1において、まずモデルパラメタ初期化部121において、分散v(k)(n,f)と空間相関行列B(k)(f)に適当な初期値を与える。
次に、モデルパラメタ推定部122において、分散v(k)(n,f)と空間相関行列B(k)(f)を周波数毎に推定する(詳細は非特許文献1参照)。この処理は周波数毎に独立に行なわれるため、各パラメタのインデックス(k)と、実際の音源10(k)とが一致しない場合がある。そこで、パーミュテーション解決部123において、各パラメタのインデックスkと実際の音源10(k)との対応を全ての周波数において一致させる。
音源分離部13では、得られたパラメタを用いて、ウィーナーフィルタ
Figure 0005406866
を計算する。但し、
Figure 0005406866
とし、R(k)(n,f)は式(2)により与えられる。そして音源分離部13はさらにソースイメージc(k)(n,f)の推定値である分離信号
c^(k)(n,f)=W(k)(n,f)X(n,f)
を計算し出力する。最後に時間領域変換部14にて周波数領域の分離信号c^(k)(n,f)を時間領域の分離信号c^(k)(t)に変換する。
N. Q. K. Duong, E. Vincent, and R. Gribonval, "Under-determined reverberant audio source separation using a full-rank spatial covariance model", IEEE Transactions on Audio, Speech and Language Processing, 2010, vol. 18, no. 7, pp. 1830-1840.
従来技術は、よく響く(=残響の多い)部屋において観測された信号も高精度・高品質で分離できることが報告されているが、モデルパラメタ推定のための繰り返し計算の収束が遅い。また、従来法は、音源数Kが既知である必要があり、音源数保持部124にて音源数Kの情報を与える必要がある。
本発明は、収束性の改善や音源数Kが未知の場合などを考慮した、高精度・高品質で音源分離可能な音源分離技術を提供することを目的とする。
上記の課題を解決するために、本発明の第一の態様に係る音源分離技術によれば、1個以上の音源信号がノイズとともに混在している観測信号から、それぞれの音源に対応する分離信号を推定する。複数の収音手段で観測された観測信号を周波数領域の信号に変換し、観測信号と、各音源が各時間周波数において観測信号にどの程度寄与するかを表す音源占有度と、を用いて、空間パラメタを推定し、観測信号と音源占有度とを用いて、音源信号の分散と空間相関行列とを含むモデルパラメタを推定し、音源信号の分散と空間相関行列に基づくウィーナーフィルタを生成し、観測信号と、空間パラメタまたはモデルパラメタと、を用いて、音源占有度を推定し、ウィーナーフィルタを用いて、観測信号からソースイメージを推定し分離信号を生成し、分離信号を時間領域の信号に変換する。
上記の課題を解決するために、本発明の第二の態様に係る音源分離技術によれば、1個以上の音源信号がノイズとともに混在している観測信号から、それぞれの音源に対応する分離信号を推定する。複数の収音手段で観測された観測信号を周波数領域の信号に変換し、観測信号と音源占有度とを用いて、音源信号の分散と空間相関行列とを含むモデルパラメタを推定し、音源信号の分散と空間相関行列に基づくウィーナーフィルタを生成し、観測信号と、モデルパラメタと、を用いて、音源占有度を推定し、音源占有度を用いて、有効音源を推定し、ウィーナーフィルタを用いて、観測信号からソースイメージを推定し有効音源に対応する分離信号を生成し、分離信号を時間領域の信号に変換する。
本発明は、音源信号s(k)(n,f)に関してスパース性の仮定を導入し、各音源が各時間周波数(n,f)において観測信号X(n,f)にどの程度寄与するかを表す音源占有度M(k)(n,f)を導入する(音源占有度については参考文献1及び2参照)。音源占有度M(k)(n,f)に基づく音源分離手法と従来のマルチチャネルウィーナーフィルタW(k)(n,f)による手法とのハイブリッド構成とすることで、収束性の改善や音源数が未知の場合などを考慮した、高精度・高品質な音源分離を行なうことができるという効果を奏する。
[参考文献1]H. Sawada, S. Araki, and S. Makino, “A two-stage frequency-domain blind source separation method for underdetermined convolutive mixtures”, in Proc. WASPAA2007, 2007, pp.139-142
[参考文献2]S. Araki, T. Nakatani, and H. Sawada, “Blind sparse source separation for unknown number of sources using Gaussian mixture model fitting with Dirichlet prior”, in Proc. ICASSP'09, 2009, pp.33-36
音源分離装置1の機能ブロック図。 音源分離装置2の構成例を示す機能ブロック図。 音源分離装置2の処理フローを示す図。 音源分離装置1及び2の収束性を示す図。 音源分離装置2の有効音源推定結果を示す図。 音源数が既知の音源分離装置1及び音源数が未知の音源分離装置2の音源分離性能を示す図。 音源分離装置3の構成例を示す機能ブロック図。 音源分離装置3の処理フローを示す図。 音源分離装置4の構成例を示す機能ブロック図。 音源分離装置4の処理フローを示す図。
<本実施形態のポイント>
本実施形態では、空間相関行列B(k)(f)のクラスタリングを行なう際に、各音源信号s(k)(n,f)に関してスパース性の仮定を導入する。すなわち、音源信号s(k)(n,f)は稀にしか大きな値をとらず、s(k)(n,f)とs(k’)(n,f)(但しk≠k’とする)は各時間周波数(n,f)で互いに重ならない、ということを仮定する。これは、互いに異なる音声信号などで確認される性質である。これを仮定すると前記式(1)は、
X(n,f)=c(z(n,f))(n,f)+e(z(n,f))(n,f)=h(z(n,f))(f)s(z(n,f))(n,f)+e(z(n,f))(n,f) (3)
となる。ここで、z(n,f)は時間周波数(n,f)において占有的(支配的)な音源のインデックス番号を表し、e(z(n,f))(n,f)は、この占有的な音源以外の成分(占有的な音源のソースイメージc(z(n,f))(n,f)から見ると雑音成分のソースイメージ)である。
さらに、観測信号X(n,f)が与えられた条件下で、z(n,f)=k番目の音源信号が占有的であるという事後確率
Figure 0005406866
を導入する。このM(k)(n,f)を、音源占有度と呼ぶ。音源占有度の導入により、有効音源や有効音源数の推定が可能になったり、音源占有度M(k)(n,f)の概念を利用する収束の早い音源分離手法と組み合わせることが可能となる。このため、音源数Kが未知であっても動作する、収束の速い音源分離手段を提供することが可能となる。
以下、本発明の実施形態について、説明する。
<第一実施形態に係る音源分離装置2>
図2及び図3を用いて第一実施形態に係る音源分離装置2を説明する。音源分離装置2は、周波数領域変換部21と、空間パラメタ推定部221と、音源占有度推定部25と、モデルパラメタ推定部222と、パーミュテーション解決部223と、音源分離部23と、時間領域変換部24と、有効音源推定部26と、を備える。
K個の音源10(k)それぞれから同時に発せられる音源信号s(k)(t)をM個の第m収音手段2で(但し、Kを1以上の整数とし、k=1,2,…,Kとし、Mを2以上の整数とし、m=1,2,…,Mとする)、ある収録時間内(例えば、5秒間)に観測する。この観測状況を状況Xとする。この状況Xの下、第m収音手段2で収音された観測信号をx(t)とし、観測された観測信号をX(t)=[x(t),x(t),…,x(t)]とする。収音手段とは例えばマイクロホンのことである。音源分離装置2は、第m収音手段2を介して、K個の音源信号s(k)(t)がノイズとともに混在している観測信号X(t)を取得し、この観測信号から、それぞれの音源10(k)に対応する分離信号c^(k)(t)を推定する。なお、本実施形態では、説明を簡単にするために収音手段を2個(M=2)とし、第一収音手段2及び第二収音手段2としているが、2個に限定されるものではない。
<周波数領域変換部21>
周波数領域変換部21は、M個の収音手段2で観測されたT個の時間領域の観測信号X(t)=[x(t),x(t),…,x(t)]を受け取り(但し、t=1,2,…,T)、これを周波数毎の時系列信号である周波数領域の観測信号X(n,f)=[x(n,f),x(n,f),…,x(n,f)]に変換し(s21)、空間パラメタ推定部221、モデルパラメタ推定部222及び音源占有度推定部25に出力する。以降、断りのない場合、観測信号とは、周波数毎の時系列信号である観測信号ベクトルX(n,f)を意味し、時間領域の観測信号の場合、それを明記する。周波数領域への変換は、例えば短時間フーリエ変換を用いれば良い。なお、本実施形態では、Nをフレームの個数とし、nをフレームに対応する時刻のインデックスとしてn=1,2,…,Nとし、Fを周波数ビンの個数とし、fを周波数ビンに対応する周波数のインデックスとしてf=1,2,…,Fとする。
<空間パラメタ推定部221と音源占有度推定部25>
空間パラメタ推定部221は、観測信号X(n,f)を受け取り、観測信号X(n,f)と、各音源10(k)が各時間周波数(n,k)において観測信号X(n,f)にどの程度寄与するかを表す音源占有度M(k)(n,f)と、を用いて、空間パラメタφを推定し(s221)、音源占有度推定部25に出力する。
音源占有度推定部25は、観測信号X(n,f)と空間パラメタφとを受け取り、これらの値を用いて、音源占有度M(k)(n,f)を推定し(s251)、空間パラメタ推定部221へ出力する。例えば、以下のようにして空間パラメタφ、音源占有度M(k)(n,f)を推定する。
観測信号ベクトルX(n,f)を空間パラメタφについてモデル化し、音源占有度M(k)(n,f)を用いて観測信号ベクトルX(n,f)を分類する。なお、音源占有度M(k)(n,f)を用いて、収束が早く、かつ、音源数Kが未知の場合にも動作する手法を用いることができることがポイントである。ここでの空間パラメタφとしては、例えば、ノルムで正規化された空間ベクトルh ̄(k)(f)=h(k)(f)/‖h(k)(f)‖(例えば、参考文献1参照)やマイク間位相差A(n,f)=arg(x(n,f)/x(n,f))(例えば、参考文献2参照)を用いることができる。
本実施形態では、ノルムで正規化された観測信号X ̄(n,f)=X(n,f)/‖X(n,f)‖を、ノルムで正規化された空間ベクトルh ̄(k)(f)=h(k)(f)/‖h(k)(f)‖を用いて
Figure 0005406866
とモデル化する。但し、σ(k)(f)はX ̄(n,f)の分散を表し、本実施形態での空間パラメタは、φ={h ̄(k)(f),σ(k)(f),α(k)=p(z(n,f)=k)である。なお、α(k)は、ある音源のインデックスkが占有的である音源のインデックスz(n,f)である確率を表す。また、式(4)において、音源数Kが既知の場合はK’=Kとし、音源数Kが未知の場合は想定される音源数より十分大きな値をK’として用いる(例えばK’=10)。よって、以下において、k=1,2,…,K’とする。
空間パラメタφの推定は、対数尤度関数
Figure 0005406866
を最大とする空間パラメタφを推定することで行なう。ここではEMアルゴリズムを用いて、空間パラメタφの推定を行なう。補助関数Q(φ)は、
Figure 0005406866
により与えられ、空間パラメタφは、以下のE−stepとM−stepを繰り返して計算することにより、高速に推定することができる。なお、φ’は現在の更新までで得られている空間パラメタである。
(空間パラメタφ及び音源占有度M(k)(n,f)の初期化)
まず、空間パラメタφ及び音源占有度M(k)(n,f)を初期化する(s2211)。例えば、全てのk、n、fに対して、以下のように初期値を与え、初期化する。
Figure 0005406866
また、i=0とする。次に、更新回数iを1ずつ増やしながら(S2212)、収束するまで以下の計算を行なう。なお、二回目以降の処理においては、一回前の処理において生成された空間パラメタφと音源占有度推定部25から受け取った最新の音源占有度M(k)(n,f)を用いて以下の処理を行なう。
(M−step(空間パラメタφの更新)(S2214))
空間パラメタ推定部221において、M−step(空間パラメタφの更新)を行なう。以下の行列
Figure 0005406866
を求め、さらに、この行列Rの最大固有値に対応する固有ベクトルを求め、その固有ベクトルを新たなh ̄(k)(f)とすることによって、h ̄(k)(f)を更新する。
また、X ̄(n,f)の分散(σ(k)(f))
Figure 0005406866
として求め、σ(k)(f)を更新する。
また、ある音源のインデックスkが占有的である音源のインデックスである確率α(k)=p(z(n,f)=k)を、ある音源の音源占有度M(k)(n,f)の全ての時間周波数(n,f)に対する平均値として求める。つまり、α(k)
Figure 0005406866
として求め、α(k)を更新する。
空間パラメタ推定部221は、推定し、更新した空間パラメタφを音源占有度推定部25に出力する。
(E−step(音源占有度M(k)(n,f)の更新)(S251))
音源占有度推定部25は、受け取った最新の空間パラメタφと、観測信号X(n,f)を用いて、
Figure 0005406866
を求め、音源占有度M(k)(n,f)を更新する。
M−stepとE−stepを収束するまで(例えば、更新回数iが20に達するまで、もしくは、Q関数の値の変化量が十分小さくなるまで)繰り返す(s2213)。
収束すると、空間パラメタ推定部221は、音源占有度推定部25から受け取った最新の音源占有度M(k)(n,f)をモデルパラメタ推定部222に出力する。
<モデルパラメタ推定部222と音源占有度推定部25>
モデルパラメタ推定部222は、観測信号X(n,f)と音源占有度M(k)(n,f)とを受け取り、これらの値を用いて、音源信号s(k)(n,f)の分散v(k)(n,f)と空間相関行列B(k)(f)とを含むモデルパラメタΘ={θ(k)を推定し、分散v(k)(n,f)と空間相関行列B(k)(f)に基づくウィーナーフィルタW(k)(n,f)を生成する(s222)。なお、{a(k)は全てのkに関するaの集合を意味し、Θ={θ(k)={θ(1),θ(2),…,θ(K’)}である。
音源占有度推定部25は、観測信号X(n,f)とモデルパラメタΘとを受け取り、これらの値を用いて、音源占有度を推定し(s252)、モデルパラメタ推定部222に出力する。例えば、以下のようにしてモデルパラメタΘ、音源占有度M(k)(n,f)を推定する。
まず、ソースイメージc(k)(n,f)の相関行列R(k)(n,f)=c(k)(n,f)(c(k)(n,f))が音源信号s(k)(n,f)の分散v(k)(n,f)と、時刻nに依存しない空間相関行列B(k)(f)を用いて、
(k)(n,f)=v(k)(n,f)B(k)(f) (12)
とモデル化できると仮定する。そして、この空間相関行列B(k)(f)のクラスタリングを行なうことで、モデルパラメタΘを推定する。このために、本実施形態では、以下のようにソースイメージc(k)(n,f)と観測信号X(n,f)をモデル化する。
ソースイメージc(k)(n,f)を、平均0、分散v(k)(n,f)B(k)(f)の複素正規分布Nを用いて、
p(c(k)(n,f);B(k)(f),v(k)(n,f))=Nc(c(k)(n,f);0,v(k)(n,f)B(k)(f)) (13)
とモデル化する。
また、式(3)の観測信号X(n,f)が、ほぼz(n,f)番目のソースイメージc(z(n,f)))(n,f)のみから成り、それ以外の音源はほぼゼロとして観測されることを表現するために、観測信号X(n,f)を
Figure 0005406866
とモデル化する。なお、k’はz(n,f)を除く音源のインデックス番号とし、k’=1,2,…,(z(n,f)−1),(z(n,f)+1),…,K’であり、δはディラックのデルタ関数を表し、U(z(n,f))(n,f)は占有的な音源のソースイメージc(z(n,f))(n,f)に対する雑音成分e(z(n,f))(n,f)の相関行列であり、以下の式により表される。
Figure 0005406866
モデルパラメタ推定部222にて推定すべきモデルパラメタは、Θ={θ(k)=({v(k)(n,f)},{B(k)(f)})である。モデルパラメタΘの推定は、対数尤度関数
Figure 0005406866
を最大とするモデルパラメタΘを推定することで行なう。なお、Dは、全ての時間周波数(n,f)及び全ての音源のインデックスkに対する、観測信号X(n,f)の集合と、占有的な音源のインデックスz(n,f)の集合と、ソースイメージc(k)(n,f)の集合と、からなる集合を表し、∫dcは、隠れ変数c(k)(n,f)についての周辺化を意味する。なお、式(16)のL(n,f)は、
Figure 0005406866
である。本実施形態ではEMアルゴリズムを用いてモデルパラメタΘの推定を行なう。補助関数は、
Figure 0005406866
により与えられる。なお、Θ’は現在の更新までに得られているモデルパラメタである。
(モデルパラメタの初期化(s2221))
モデルパラメタ推定部222は、空間パラメタ推定部221から音源占有度M(k)(n,f)を受け取ると、まず、分散v(k)(n,f)と空間相関行列B(k)(f)を初期化する(s2221)。例えば、空間パラメタ推定部221で受け取った音源占有度M(k)(n,f)と観測信号X(n,f)の要素x(n,f)(第一収音手段2の観測値)を用いて、
Figure 0005406866
として初期化する。
また、更新回数i=0とする。以下、iを増やしながら(S2222)、E−stepとM−stepを収束するまで繰り返し、モデルパラメタΘは、この繰り返しにより更新されながら推定される。
(M−step(モデルパラメタΘの更新)(s2224))
モデルパラメタ推定部222において、
Figure 0005406866
を計算する。なお、Tr(A)は、行列Aの対角成分の和を返す処理(トレース)を意味する。E−stepとM−stepとの繰り返し処理において、初めて式(23)を計算する場合は、式(23)におけるM(k)(n,f)は、空間パラメタ推定部221の出力値である音源占有度M(k)(n,f)を用い、2回目以降は音源占有度推定部25の最新の出力値である音源占有度M(k)(n,f)を用いる。ここで、
Figure 0005406866
であり、式(28)におけるW(k)(n,f)が、音源分離のためのマルチチャネルウィーナーフィルタW(k)(n,f)であり、ソースイメージc(k)(n,f)の推定値である分離信号c^(k)(n,f)は、これを用いて式(26)で計算される。
さらに、モデルパラメタ推定部222は、
Figure 0005406866
を求める。この値は、後述する式(32)において用いる。
(E−step(音源占有度M(k)(n,f)の推定(s252))
式(19)のQ関数の中のp({c(k)(n,f)},z(n,f)|X(n,f))の項は、
Figure 0005406866
と表すことができ、p(z(n,f)=k|X(n,f))が音源占有度M(k)(n,f)に対応する。
音源占有度推定部25は、音源占有度M(k)(n,f)を以下の式により推定する。
Figure 0005406866
なお、c^(k)(n,f)とr(k)(n,f)は式(26)および式(27)で与えられる。また||A||=A−1Aとする。
以上を収束するまで(例えば、更新回数iが20に達するまで、もしくは、Q関数の値の変化量が十分小さくなるまで)繰り返す(s2223)。
収束すると、モデルパラメタ推定部222は、音源占有度推定部25から受け取った最新の音源占有度M(k)(n,f)と式(28)により生成した最新のウィーナーフィルタW(k)(n,f)を周波数毎に紐付けてパーミュテーション解決部223に出力する。
<パーミュテーション解決部223>
パーミュテーション解決部223は、音源占有度M(k)(n,f)とこれに紐付けられたウィーナーフィルタW(k)(n,f)を受け取り、音源毎にウィーナーフィルタW(k)(n,f)をまとめる(s223)。
空間パラメタφ及びモデルパラメタΘの推定は周波数毎に行なわれるため、各パラメタのインデックス番号kと、そのクラスタに対応する実際の音源10(k)とが一致しない場合がある。例えば、ある周波数fではk=1が音源10(1)に、k=2が音源10(2)に対応するが、別の周波数f’ではk=1が音源10(2)に、k=2が音源10(1)に対応する、というように、周波数毎に対応関係がばらばらになってしまうことが一般的である。これをパーミュテーションの問題という。そこで、パーミュテーション解決部223において、全ての周波数f=1,2,…Fで各パラメタのインデックス(k)と実際の音源10(k)とが完全に一対一対応するように整える。これは、例えば次のように行なわれる。
まず、各周波数fおよび各インデックスkにおいて得られた各音源占有度M(k)n,fを、
γ(k)(f)=[M(k)(1,f),・・・,M(k)(N,f)]
というベクトルとする。同じ音源であれば、音源占有度M(k)(n,f)は、全ての周波数で同期する性質があることを利用し、異なる周波数間でのベクトルγ(k)(f)とγ(k)(f’)の相関が全ての周波数で最大となるように、インデックスの番号を入れ替える。すなわち、ベクトルaとbの相関係数をρ(a,b)とした場合に、
Figure 0005406866
を最大とするkの配列Π(k)(f)を求める。ここで配列Π(k)(f)は、1,・・・,K’の整数が適切な順序で並んだ物であり、γ ̄(k)は、全ての周波数におけるインデックスkに対応するγ(k)(f)の平均値である。上記Jの最大化は、例えば以下の繰り返し演算により行なうことができる。
Figure 0005406866
これにより全ての周波数で各パラメタのインデックスkと音源10(k)との対応関係を揃えることができる。
パーミュテーション解決部223は、音源毎にまとめた音源占有度M(k)(n,f)を有効音源推定部26に出力する。さらに、音源毎にまとめた音源占有度M(k)(n,f)に紐付けられたウィーナーフィルタW(k)(n,f)を音源分離部23に出力する。
<有効音源推定部26>
有効音源推定部26は、音源毎にまとめられた音源占有度M(k)(n,f)を受け取り、この音源占有度M(k)(n,f)を用いて、有効音源を推定し(s26)、有効音源を音源分離部23に出力し、有効音源数Kを音源分離装置2の出力値として出力する。音源のインデックスk=1,2,…,K’のうち、有効音源に対応するインデックスの集合を{k}とし、その有効音源数Kとする。例えば、以下の方法で有効音源を推定する。本実施形態で利用する音源占有度M(k)(n,f)は、各時間周波数(n,f)におけるk番目の音源10(k)の占有度を表しているため、パーミュテーション問題を解決した後の音源占有度M(k)(n,f)の平均値を、各インデックスkについて求めれば、占有度の高い音源を求めることが可能となる。よって、
Figure 0005406866
を計算し、p(z(n,f)=k)の値が予め設定した閾値thより大きいとき、インデックスkを占有度の高い音源10(k)のインデックスとして判定し、そのkの集合{k}を出力する。また有効音源数Kも出力する。
<音源分離部23>
音源分離部23は、ウィーナーフィルタW(k)(n,f)と有効音源の集合{k}を受け取り、有効音源の集合{k}に対応するマルチチャネルウィーナーフィルタW(k)(n,f)を生成する。さらに、音源分離部23は、観測信号X(n,f)を受け取り、マルチチャネルウィーナーフィルタW(k)(n,f)を用いて、式(26)により観測信号X(n,f)からソースイメージc(k)(n,f)を推定した、有効音源に対応する分離信号c^(k)(n,f)を生成し(s23)、時間領域変換部24に出力する。
<時間領域変換部24>
時間領域変換部24は、周波数領域の分離信号c^(n,f)を受け取り、これを時間領域の分離信号c^(t)に変換し(s24)、この値を音源分離装置2の出力値として出力する。なお、時間領域への変換は、周波数領域変換部21で用いた変換方法に対応するものであればよい。
<効果>
このような構成とすることで、各パラメタの収束の早く、音源数が未知の場合にも動作する、高精度・高品質な音源分離を行なうことができる。
<シミュレーション結果>
第一実施形態の効果を調べるため、従来技術(非特許文献1参照)及び第一実施形態の音源分離装置で音源分離を行なった。実験にて、マイクロホン数は2、音源数は2または3とした。サンプリング周波数は8kHz、マイクロホンの間隔は4cmである。
図4は、音源数Kが既知の場合(2または3)に、4通りの音声の組合せについて、信号全体の歪みの尺度SDR(Signal to distortion ratio)を評価し、その平均を求めたものである。図4において、HB1は第一実施形態においてモデルパラメタΘの更新を1回のみにした場合の性能を、HB50は従来技術においてモデルパラメタの更新を50回にした場合の性能を示す。なお、図4及び図6において、実験時の残響時間は250msまたは400msとしている。第一実施形態は、学習回数1回のみにもかかわらず、従来技術で、学習を50回行なった場合よりも高い性能を示すことが分かる。これより、第一実施形態は、少ないモデルパラメタΘの更新で高い性能を示すことから、その収束性の早さが示された。
図5は、音源数未知の場合に、K’=8として実施例を用いた場合の重み係数p(z(n,f)=k)(式(37)参照)をプロットしたものである。この結果より、有効音源推定部26において有効音源及びその数の推定が可能であることが分かる。
図6は、音源数Kが未知の場合に、4通りの音声組合せについてSDRを評価し、その平均を求めた物である。図6において、HB(K given),HB(K unknown)はそれぞれ、正しい音源数Kを従来技術のシステムに与えた場合と、音源数未知の条件下でK’=8として第一実施形態を適用した場合の性能を示している。第一実施形態は、音源数未知の場合でも、音源数既知の従来技術と同程度の分離性能を示すことがわかる。
<他の変形例>
本実施形態においては、ソースイメージ及び観測信号のモデルとしてそれぞれ(13)、(14)を用いたが、それぞれ他の適切なモデルを用いることも可能である。
本実施形態においては、有効音源の推定をモデルパラメタΘの推定の後に行なったが、これを空間パラメタφの推定後に行ない、モデルパラメタΘの推定は、推定された音源数Kの音源に対してのみ行なってもよい。この場合、空間パラメタφの推定後に、音源占有度M(k)(n,f)とそれに紐付けられる空間パラメタφをパーミュテーション解決部223の入力とし、パーミュテーション解決部223は、音源毎に音源占有度M(k)(n,f)(と空間パラメタφ)をまとめ、有効音源推定部26に出力する。有効音源推定部26は、上述の方法により、有効音源を推定し、有効音源に対応するインデックスの集合を{k}と有効音源数Kと、それに対応する音源占有度M(k)(n,f)と空間パラメタφをモデルパラメタ推定部222に出力する。
本実施形態においては、有効音源推定部26において、有効音源数Kを推定しているが、予め有効音源数Kが利用者等により与えられている場合には、有効音源推定部26を設けなくともよい。その場合には、各部は、図示しない音源数保持部から記憶されている有効音源数Kを取得する。このような構成の場合には、各パラメタの収束の早く、かつ、高精度・高品質な音源分離を行なうことができる。
なお、空間パラメタφ、モデルパラメタΘ及び音源占有度M(k)(n,f)の初期値は上述した値以外の値であってもよい。例えば、各パラメタが取りうる値をランダムに設定してもよい。
<第二実施形態>
図7及び図8を用いて第二実施形態に係る音源分離装置3を説明する。第一実施形態と異なる部分についてのみ説明する。音源分離装置3は、周波数領域変換部21と、音源占有度推定部35と、モデルパラメタ推定部322と、パーミュテーション解決部223と、音源分離部23と、時間領域変換部24と、有効音源推定部26と、を備える。空間パラメタ推定部221を備えていない点、及び、音源占有度推定部35とモデルパラメタ推定部322における処理内容が第一実施形態とは異なる。
<モデルパラメタ推定部322と音源占有度推定部35>
モデルパラメタ推定部322は、観測信号X(n,f)を受け取り、観測信号X(k)(n,f)と音源占有度M(k)(n,f)を用いて、音源信号の分散v(k)(n,f)と空間相関行列B(k)(f)とを含むモデルパラメタΘを推定し、音源信号の分散v(k)(n,f)と空間相関行列B(k)(f)に基づくウィーナーフィルタW(k)(n,f)を生成する(s322)。
音源占有度推定部35は、観測信号X(n,f)とモデルパラメタΘとを受け取り、これらの値を用いて、音源占有度M(k)(n,f)を推定し(s352)、モデルパラメタ推定部222に出力する。
第一実施形態とは異なり、空間パラメタ推定部がないため、音源占有度推定部35は、観測信号X(n,f)と空間パラメタφとを用いて、音源占有度M(k)(n,f)を推定する必要がなく、音源占有度M(k)(n,f)の初期値をモデルパラメタ推定部322において与える(s3221)。例えば、全てのk、n、fに対して、M(k)(n,f)=1として初期値を与える。
他の処理については第一実施形態と同様である。
<効果>
第二実施形態の音源分離装置3は、収束の高速化の効果はなくなるが、音源数が未知の場合にも動作し、高精度・高品質な音源分離を行なうことができる。
<第三実施形態>
図9及び図10を用いて第三実施形態に係る音源分離装置4を説明する。音源分離装置2と異なる部分についてのみ説明する。音源分離装置4は、周波数領域変換部21と、空間パラメタ推定部421と、音源占有度推定部25と、モデルパラメタ推定部422と、パーミュテーション解決部223と、音源分離部23と、時間領域変換部24と、有効音源推定部26と、を備える。
第一実施形態では、空間パラメタφの推定を十分収束するまで行なってから、モデルパラメタΘの推定を行なったが、本実施形態では、空間パラメタφとモデルパラメタΘの推定を、それぞれ1回ずつ更新しながら、全体としての最適化を行なう点が異なる。
<空間パラメタ推定部421>
空間パラメタ推定部421は、観測信号X(n,f)を受け取り、観測信号X(n,f)と、各音源信号s(k)(n,f)に対応するソースイメージc(k)(n,f)が各時間周波数(n,k)において観測信号X(n,f)にどの程度寄与するかを表す音源占有度M(k)(n,f)と、を用いて、空間パラメタφを推定し(s421)、音源占有度推定部45に出力する。
第一実施形態のs2211において説明した空間パラメタφ及び音源占有度M(k)(n,f)の初期化とs2221において説明したモデルパラメタΘの初期化を繰り返し処理に先立ち行なう。
以下、iを増やしながら(S2222)、空間パラメタ推定部421と音源占有度推定部25とモデルパラメタ推定部422における処理を繰り返し、空間パラメタφ及びモデルパラメタΘは、この繰り返しにより更新されながら推定される。
空間パラメタ推定部421において、空間パラメタφの推定(s2214)は第一実施形態と同様の処理により行なわれ、空間パラメタφは音源占有度推定部25に出力される。
音源占有度推定部25は、第一実施形態と同様の処理により、音源占有度M(k)(n,f)を推定し(s251)、空間パラメタ推定部421へ出力する。空間パラメタ推定部421は、受け取った音源占有度M(k)(n,f)をモデルパラメタ推定部422に出力する。
<モデルパラメタ推定部422>
モデルパラメタ推定部422は、第一実施形態と同様の処理により、モデルパラメタΘを推定し、これに基づくウィーナーフィルタW(k)(n,f)を生成する(s422)。なお、本実施形態では、空間パラメタ推定部421と音源占有度推定部25とモデルパラメタ推定部422における処理をひとまとめとして処理として繰り返すため、モデルパラメタ推定部422においてiを増やす必要はない。また、モデルパラメタ推定部422は、モデルパラメタΘを受け取るたびに収束しているか否かを判定し(s4223)、収束すると、モデルパラメタ推定部222は、音源占有度推定部25から受け取った最新の音源占有度M(k)(n,f)と式(28)により生成した最新のウィーナーフィルタW(k)(n,f)を周波数毎に紐付けてパーミュテーション解決部223に出力する。収束していない場合には、第一実施形態と同様の処理により、モデルパラメタΘを推定し、音源占有度推定部25に出力する。
音源占有度推定部25は、第一実施形態と同様の処理により音源占有度M(k)(n,f)を推定し(s252)、モデルパラメタ推定部422に出力する。さらに、モデルパラメタ推定部422は受け取った音源占有度M(k)(n,f)を空間パラメタ推定部421に出力する。
<効果>
このような構成とすることで第一実施形態と同様の効果を得ることができる。
<プログラム及び記録媒体>
上述した音源分離装置は、コンピュータにより機能させることもできる。この場合はコンピュータに、目的とする装置(各種実施形態で図に示した機能構成をもつ装置)として機能させるためのプログラム、またはその処理手順(各実施形態で示したもの)の各過程をコンピュータに実行させるためのプログラムを、CD−ROM、磁気ディスク、半導体記憶装置などの記録媒体から、あるいは通信回線を介してそのコンピュータ内にダウンロードし、そのプログラムを実行させればよい。
<その他の変形例>
本発明は上記の実施形態及び変形例に限定されるものではない。例えば、上述の各種の処理は、記載に従って時系列に実行されるのみならず、処理を実行する装置の処理能力あるいは必要に応じて並列的にあるいは個別に実行されてもよい。その他、本発明の趣旨を逸脱しない範囲で適宜変更が可能である。また、各パラメタは、図示しない記憶部等に格納され、各部はこの記憶部から各パラメタを取得する構成としてもよい。
2,3,4 音源分離装置
収音手段
10 音源
21 周波数領域変換部
23 音源分離部
24 時間領域変換部
25,35,45 音源占有度推定部
26 有効音源推定部
221,421 空間パラメタ推定部
222,322,422 モデルパラメタ推定部
223 パーミュテーション解決部

Claims (7)

  1. 1個以上の音源信号がノイズとともに混在している観測信号から、それぞれの音源に対応する分離信号を推定する音源分離装置であって、
    複数の収音手段で観測された前記観測信号を周波数領域の信号に変換する周波数領域変換手段と、
    前記観測信号と、各音源が各時間周波数において観測信号にどの程度寄与するかを表す音源占有度と、を用いて、空間パラメタを推定する空間パラメタ推定手段と、
    前記観測信号と前記音源占有度とを用いて、前記音源信号の分散と空間相関行列とを含むモデルパラメタを推定し、前記音源信号の分散と空間相関行列に基づくウィーナーフィルタを生成するモデルパラメタ推定手段と、
    前記観測信号と、前記空間パラメタまたは前記モデルパラメタと、を用いて、前記音源占有度を推定する音源占有度推定手段と、
    前記ウィーナーフィルタを用いて、前記観測信号からソースイメージを推定し分離信号を生成する音源分離手段と、
    前記分離信号を時間領域の信号に変換する時間領域変換手段と、を備える、
    音源分離装置。
  2. 請求項1記載の音源分離装置であって、
    前記音源占有度を用いて、有効音源を推定する有効音源推定手段と、をさらに備え、
    前記音源分離手段は、前記有効音源に対応する分離信号を生成する、
    音源分離装置。
  3. 1個以上の音源信号がノイズとともに混在している観測信号から、それぞれの音源に対応する分離信号を推定する音源分離装置であって、
    複数の収音手段で観測された前記観測信号を周波数領域の信号に変換する周波数領域変換手段と、
    前記観測信号と、各音源が各時間周波数において観測信号にどの程度寄与するかを表す音源占有度を用いて、前記音源信号の分散と空間相関行列とを含むモデルパラメタを推定し、前記音源信号の分散と空間相関行列に基づくウィーナーフィルタを生成するモデルパラメタ推定手段と、
    前記観測信号と前記モデルパラメタとを用いて、前記音源占有度を推定する音源占有度推定手段と、
    前記音源占有度を用いて、有効音源を推定する有効音源推定手段と、
    前記ウィーナーフィルタを用いて、前記観測信号からソースイメージを推定し前記有効音源に対応する分離信号を生成する音源分離手段と、
    前記分離信号を時間領域の信号に変換する時間領域変換手段と、を備える、
    音源分離装置。
  4. 1個以上の音源信号がノイズとともに混在している観測信号から、それぞれの音源に対応する分離信号を推定する音源分離方法であって、
    複数の収音ステップで観測された前記観測信号を周波数領域の信号に変換する周波数領域変換ステップと、
    前記観測信号と、各音源が各時間周波数において観測信号にどの程度寄与するかを表す音源占有度と、を用いて、空間パラメタを推定する空間パラメタ推定ステップと、
    前記観測信号と、前記空間パラメタと、を用いて、前記音源占有度を推定する第一音源占有度推定ステップと、
    前記観測信号と前記音源占有度とを用いて、前記音源信号の分散と空間相関行列とを含むモデルパラメタを推定し、前記音源信号の分散と空間相関行列に基づくウィーナーフィルタを生成するモデルパラメタ推定ステップと、
    前記観測信号と、前記モデルパラメタと、を用いて、前記音源占有度を推定する第二音源占有度推定ステップと、
    前記ウィーナーフィルタを用いて、前記観測信号からソースイメージを推定し分離信号を生成する音源分離ステップと、
    前記分離信号を時間領域の信号に変換する時間領域変換ステップと、を備える、
    音源分離方法。
  5. 請求項4記載の音源分離方法であって、
    前記音源占有度を用いて、有効音源を推定する有効音源推定ステップと、をさらに備え、
    前記音源分離ステップは、前記有効音源に対応する分離信号を生成する、
    音源分離方法。
  6. 1個以上の音源信号がノイズとともに混在している観測信号から、それぞれの音源に対応する分離信号を推定する音源分離方法であって、
    複数の収音ステップで観測された前記観測信号を周波数領域の信号に変換する周波数領域変換ステップと、
    前記観測信号と、各音源が各時間周波数において観測信号にどの程度寄与するかを表す音源占有度を用いて、前記音源信号の分散と空間相関行列とを含むモデルパラメタを推定し、前記音源信号の分散と空間相関行列に基づくウィーナーフィルタを生成するモデルパラメタ推定ステップと、
    前記観測信号と前記モデルパラメタとを用いて、前記音源占有度を推定する音源占有度推定ステップと、
    前記音源占有度を用いて、有効音源を推定する有効音源推定ステップと、
    前記ウィーナーフィルタを用いて、前記観測信号からソースイメージを推定し前記有効音源に対応する分離信号を生成する音源分離ステップと、
    前記分離信号を時間領域の信号に変換する時間領域変換ステップと、を備える、
    音源分離方法。
  7. 請求項1から請求項3の何れかに記載の音源分離装置として、コンピュータを機能させるためのプログラム。
JP2011036559A 2011-02-23 2011-02-23 音源分離装置、その方法及びプログラム Active JP5406866B2 (ja)

Priority Applications (1)

Application Number Priority Date Filing Date Title
JP2011036559A JP5406866B2 (ja) 2011-02-23 2011-02-23 音源分離装置、その方法及びプログラム

Applications Claiming Priority (1)

Application Number Priority Date Filing Date Title
JP2011036559A JP5406866B2 (ja) 2011-02-23 2011-02-23 音源分離装置、その方法及びプログラム

Publications (2)

Publication Number Publication Date
JP2012173584A JP2012173584A (ja) 2012-09-10
JP5406866B2 true JP5406866B2 (ja) 2014-02-05

Family

ID=46976515

Family Applications (1)

Application Number Title Priority Date Filing Date
JP2011036559A Active JP5406866B2 (ja) 2011-02-23 2011-02-23 音源分離装置、その方法及びプログラム

Country Status (1)

Country Link
JP (1) JP5406866B2 (ja)

Families Citing this family (5)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP6463904B2 (ja) * 2014-05-26 2019-02-06 キヤノン株式会社 信号処理装置及び音源分離方法及びプログラム
CN105989851B (zh) * 2015-02-15 2021-05-07 杜比实验室特许公司 音频源分离
JP7112269B2 (ja) * 2018-07-09 2022-08-03 日本放送協会 方向別収音装置及びプログラム
CN110111808B (zh) * 2019-04-30 2021-06-15 华为技术有限公司 音频信号处理方法及相关产品
CN113362848B (zh) * 2021-06-08 2022-10-04 北京小米移动软件有限公司 音频信号处理方法、装置及存储介质

Family Cites Families (2)

* Cited by examiner, † Cited by third party
Publication number Priority date Publication date Assignee Title
JP4138290B2 (ja) * 2000-10-25 2008-08-27 松下電器産業株式会社 ズームマイクロホン装置
JP4787777B2 (ja) * 2007-03-13 2011-10-05 日本電信電話株式会社 信号分離装置、信号分離方法、信号分離プログラム、記録媒体

Also Published As

Publication number Publication date
JP2012173584A (ja) 2012-09-10

Similar Documents

Publication Publication Date Title
JP7191793B2 (ja) 信号処理装置、信号処理方法、及びプログラム
US20210089967A1 (en) Data training in multi-sensor setups
US9357298B2 (en) Sound signal processing apparatus, sound signal processing method, and program
WO2016100460A1 (en) Systems and methods for source localization and separation
US10192568B2 (en) Audio source separation with linear combination and orthogonality characteristics for spatial parameters
WO2016152511A1 (ja) 音源分離装置および方法、並びにプログラム
JP5337072B2 (ja) モデル推定装置、音源分離装置、それらの方法及びプログラム
US10373628B2 (en) Signal processing system, signal processing method, and computer program product
JP5406866B2 (ja) 音源分離装置、その方法及びプログラム
CN106031196B (zh) 信号处理装置、方法以及程序
JP2022529912A (ja) 深層フィルタを決定するための方法および装置
Nesta et al. Robust Automatic Speech Recognition through On-line Semi Blind Signal Extraction
KR101243897B1 (ko) 신호의 시간 지연 및 감쇄 추정에 기반한 반향 환경에서의 암묵 음원 분리 방법
Sheeja et al. CNN-QTLBO: an optimal blind source separation and blind dereverberation scheme using lightweight CNN-QTLBO and PCDP-LDA for speech mixtures
GB2510650A (en) Sound source separation based on a Binary Activation model
Duong et al. Spatial covariance models for under-determined reverberant audio source separation
JP5387442B2 (ja) 信号処理装置
JP4946330B2 (ja) 信号分離装置及び方法
Bavkar et al. PCA based single channel speech enhancement method for highly noisy environment
Mirzaei et al. Under-determined reverberant audio source separation using Bayesian Non-negative Matrix Factorization
Miyazaki et al. Theoretical analysis of parametric blind spatial subtraction array and its application to speech recognition performance prediction
Li et al. Low complex accurate multi-source RTF estimation
Ukai et al. Multistage SIMO-model-based blind source separation combining frequency-domain ICA and time-domain ICA
JP4714892B2 (ja) 耐高残響ブラインド信号分離装置及び方法
JP7126659B2 (ja) 信号処理装置、信号処理方法及び信号処理プログラム

Legal Events

Date Code Title Description
A621 Written request for application examination

Free format text: JAPANESE INTERMEDIATE CODE: A621

Effective date: 20130213

A977 Report on retrieval

Free format text: JAPANESE INTERMEDIATE CODE: A971007

Effective date: 20130925

TRDD Decision of grant or rejection written
A01 Written decision to grant a patent or to grant a registration (utility model)

Free format text: JAPANESE INTERMEDIATE CODE: A01

Effective date: 20131022

A61 First payment of annual fees (during grant procedure)

Free format text: JAPANESE INTERMEDIATE CODE: A61

Effective date: 20131101

R150 Certificate of patent or registration of utility model

Ref document number: 5406866

Country of ref document: JP

Free format text: JAPANESE INTERMEDIATE CODE: R150

S531 Written request for registration of change of domicile

Free format text: JAPANESE INTERMEDIATE CODE: R313531

R350 Written notification of registration of transfer

Free format text: JAPANESE INTERMEDIATE CODE: R350