1-Q-40

1-Q-40
絶対値の観測のみを用いた 2 つの複素信号の相関係数推定∗
○宮部滋樹 (筑波大), 小野順貴 (NII/総研大), 牧野昭二 (筑波大)
1
はじめに
本稿では，以上のような絶対値 Yi , i = 1, 2 のみが
ウィーナーフィルタやスペクトル減算といった雑音
除去，線形アレー信号処理のポストフィルタやエコー
抑圧，近年では非負値行列分解を用いた音源分離な
ど，本来複素数値で観測される時間周波数領域の信
観測される条件で相関係数 |ρ| を求める問題について
議論する。ここで，言うまでもなく分散 σi2 , i = 1, 2
の推定は標本平均で求めることができる。
[
]
2
σi2 ← E yi (n)
(8)
号の振幅のみを用いる非線形な信号処理が求められ
ることがある。振幅の処理が用いられる理由は，線形
ここで yi (n), n = 1, · · · , N は確率変数 Yi の N 個の
処理よりも強力な雑音抑圧が求められる場合のほか，
標本を表し，以下では確率変数を大文字で，標本を標
位相の信頼性が低い場合など様々である。このよう
本インデックス n を引数とする小文字で表すものと
な非線形な信号処理の多くは統計的な分析に基づい
する。しかし相互相関については絶対値観測の相互
ており，相互相関などの統計量の推定が重要となる。
相関に一致しない。
このような統計量は，振幅の統計量が求められる場
E [Y1 Y2 ] ̸= σ1 σ2 |ρ|
合だけでなく，元の複素数信号の統計量が求められる
(9)
場合もある。しかし，信頼できる位相の観測ができな
そのため以下のような絶対値の積の平均は相関係数
い場合などでは，振幅の観測のみを用いて複素数の
|ρ| の良い推定を与えない。
相互相関を推定するのは難しい。本稿では，位相情報
を隠れ変数と見なした確率モデルによって，振幅の観
|ρ| ←
測のみから複素相関を推定する手法について述べる。
2
N
1 ∑ y1 (n) y2 (n)
N n=1
σ1 σ2
(10)
従って，複素変数の相関係数を絶対値のみから推定す
問題設定
る問題は自明ではない。
いま，相関のある 2 つの零平均複素確率変数 Xi ∈
C, i = 1, 2 があるものとする。また，２つの確率変
3
数はそれぞれ以下のような統計量を持つものとする。
3.1
E [Xi ] = 0
[ ]
E Xi2 = 0
[
]
2
E |Xi | = σi2
E
[X1 X2∗ ]
= σ1 σ2 ρ
(1)
(2)
(3)
(4)
∗
ここで {·} は複素共役を，E [·] は期待値を表す。ま
た式 (2) は，Xi の偏角 Θi の周辺確率密度分布が一様
分布であることを示す。
Θi = ∠Xi
1
pΘi (θi ) =
, − π ≤ θi < π
2π
複素正規分布を仮定した相関推定
確率モデル
本稿では，確率変数 xi , i = 1, 2 が式 (1)–(4) を満
たす単純な仮定として，以下のような零平均 2 変量
複素正規分布に従うと仮定した確率モデルについて
議論する。
(
)
σ22 |x1 |2 +σ12 |x2 |2 −2σ1 σ2 Re[ρ∗ x1 x∗
2]
exp −
σ12 σ22 (1−|ρ|2 )
(
)
pX1 ,X2 (x1 , x2 ; ρ) =
2
π 2 σ12 σ22 1 − |ρ|
(11)
ここで
(5)
(6)
しかし複素変数 Xi , i = 1, 2 は観測されず，代わりに
xR
i = Re [xi ]
(12)
xIi = Im [xi ]
(13)
として，極座標系への変数変換
その絶対値 Yi のみが観測されるものとする。
Yi = |Xi |
(7)
I
R
I
dxR
1 dx1 dx2 dx2 = y1 y2 dy1 dy2 dθ1 dθ2
(14)
これは偏角 Θi が失われて観測されることに相当する。
∗
Estimation of correlation coeﬃcient between to complex signals only with observation of absolute values. by MIYABE, Shigeki (University of Tsukuba), ONO, Nobutaka (NII/Sokendai), MAKINO, Shoji
(University of Tsukuba)
日本音響学会講演論文集
- 735 -
2014年9月
により，観測される絶対値 Yi と観測されない偏角 Θi
また，完全データ {Y1 , Y2 , Θ} の同時確率密度関数が
の同時確率密度が以下のように与えられる。
以下のように与えられる。
pY1 ,Y2 ,Θ1 ,Θ2 (y1 , y2 , θ1 , θ2 ; ρ)
(
)
σ22 y12 +σ12 y22 −2σ1 σ2 |ρ|y1 y2 cos(θ1 −θ2 −∠ρ)
y1 y2 exp −
σ12 σ22 (1−|ρ|2 )
(
)
=
2
π 2 σ12 σ22 1 − |ρ|
pY1 ,Y2 ,Θ (y1 , y2 , θ; ρ)
= pΘ|Y1 ,Y2 (θ|y1 , y2 ; ρ) pY1 ,Y2 (y1 , y2 ; |ρ|)
(
)
σ 2 y 2 +σ 2 y 2 −2σ σ2 y1 y2 cos(θ−∠ρ)
2y1 y2 exp − 2 1 1 2σ2 σ2 11−|ρ|
2
)
1 2(
(
)
=
2
2
2
πσ1 σ2 1 − |ρ|
(15)
これを偏角 θ1 , θ2 で周辺化することにより，観測され
(19)
る絶対値 Y1 , Y2 の尤度が以下のような 2 変量レイリー
以上のように求められた事後密度と同時密度を用
分布 [1] として求められる。
いると，尤度を繰り返し計算で最適化する EM アル
pY1 ,Y2 (y1 , y2 ; |ρ|)
∫ π ∫ π
pY1 ,Y2 ,Θ1 ,Θ2 (y1 , y2 , θ1 , θ2 ; ρ) dθ1 dθ2
=
−π −π


4y1 y2
2
|ρ|
y
y
(
) I0 
( 1 2 )
=
2
2
σ12 σ22 1 − |ρ|
σ1 σ2 1 − |ρ|


2 2
2 2
y
+
y
σ
σ
1 2 
)
· exp − 2 1(
(16)
2
2
2
σ1 σ2 1 − |ρ|
ここで Iν (·) は次数 ν ∈ R の第 1 種変形ベッセル関数
である。相互相関パラメタ ρ は尤度の中では絶対値
を取った相関係数 |ρ| の形で現れ、確率が共分散の偏
角に依存しないことがわかる。従ってこの確率モデル
を用いた最尤法により，相関係数 |ρ| を推定すること
ができる。しかしこの最尤推定は解析的に解くこと
ができない。
3.2
EM アルゴリズムによる最尤推定
ゴリズムのための Q 関数は以下のように求めること
ができる。
Q (|ρ| ; |¯
ρ|)
∫ π
=
pΘ|Y1 ,Y2 (θ|y1 , y2 ; ρ¯) log pY1 ,Y2 ,Θ (y1 , y2 , θ; ρ)
−π
σ22 y12
=−
+ σ12 y22 − 2σ1 σ2 y1 y2 λ
)
(
2
σ12 σ22 1 − |ρ|
(
)
2
− log 1 − |ρ| + const.
(20)
ここで λ は以下のように求められる |¯
ρ| を固定した条
件での cos(θ − ∠ρ) の期待値積分を表す。
∫ π
λ=
cos (θ − ∠ρ) pΘ|Y1 ,Y2 (θ|y1 , y2 ) dθ
−π
(
)
¯ 1 y2
I1 σ σ2|ρ|y
¯ 2)
1 2 (1−|ρ|
)
= (
(21)
2|ρ|y
¯ 1 y2
I0 σ σ 1−|ρ|
2
¯ )
1 2(
観測できない偏角 θ1 (n), θ2 (n) と観測 y1 (n), y2 (n)
また，|¯
ρ| を固定した条件で Q 関数を最大化するパラ
を併せて完全データとして扱うことにより，式 (16)
メタ推定 |ρ| は，以下の偏微分を用いた極大点の推定
の最尤推定を繰り返し最適化により求める EM アル
により得られる。
∂Q (|ρ| ; |¯
ρ|)
=0
∂ |ρ|
ゴリズムを導出する。
まず，以下のように定義される θ を隠れ変数として
扱うこととする。
これを式 (8) の σi2 , i = 1, 2 の推定を用いて解くと，
Θ = Θ 1 − Θ2
(17)
最適な |ρ| の推定は以下のようになる。
|ρ| =
式 (6) の，θi の周辺分布が一様であるという仮定を用
いると，隠れ変数 θ の事後確率密度 p(θ|y1 , y2 ) が以
下のように求められる。
(23)
以上より，EM アルゴリズムによる相関係数 |ρ| の
を与え，E ステップの更新
)
(
1 (n)y2 (n)
I1 2|ρ|y
σ1 σ2 (1−|ρ|2 )
)
λ (n) ← (
2|ρ|y1 (n)y2 (n)
I0 σ σ 1−|ρ|2
)
1 2(
= pΘ|Θ2 ,Y1 ,Y2 (θ + θ2 |θ2 , y1 , y2 ; ρ)
pY1 ,Y2 ,Θ1 ,Θ2 (y1 , y2 , θ + θ2 , θ2 ; ρ)
pΘ2 (θ2 ) pY1 ,Y2 (y1 , y2 ; |ρ|)
(
)
1 y2 cos(θ−∠ρ)
exp 2|ρ|y
σ1 σ2 (1−|ρ|2 )
(
)
=
1 y2
2πI0 σ σ2|ρ|y
2
1 2 (1−|ρ| )
y1 y2 λ
σ1 σ2
推定は以下のようになる。最初に |ρ| に適当な初期値
pΘ1 |Y1 ,Y2 (θ|y1 , y2 ; ρ)
=
日本音響学会講演論文集
(22)
(18)
(24)
と M ステップの更新
|ρ| =
- 736 -
E [y1 (n) y2 (n) λ (n)]
σ1 σ2
(25)
2014年9月
1
Estimated correlation coefficient
Estimated correlation coefficient
1
0.8
0.6
0.4
0.2
0
0
Estimation with angle
Absolute correlation
Proposed
0.2
0.4
0.6
0.8
True correlation coefficient
0.8
0.6
0.4
0.2
Estimation with angle
Absolute correlation
Proposed
0
0
1
0.2
0.4
0.6
0.8
True correlation coefficient
1
Fig. 1 Estimation results of correlation coeﬃcients
between two pseudo random variables of complex
Fig. 2 Estimation results of correlation coeﬃcients
of the non-Gaussian data generated from pseudo
normal distribution.
random variables following complex generalized normal distributions.
を繰り返すことにより，|ρ| の推定は局所最適解に収
束する。ここで式 (24) の E ステップの更新において，
て生成した。ここではまず 2 つの複素一般化正規分
第 1 種変形ベッセル関数は引数が大きい場合に浮動
布に従う疑似乱数 [2] を生成し，次にそれらの線形混
小数点例外 Inf を生じてしまうため、このような場合
合によって相関を持たせたデータを生成した。相関係
には以下の漸近近似で置き換える必要がある。
(
)
2|ρ|y1 (n)y2 (n)
I1 σ σ 1−|ρ|2
)
1 2(
2 |ρ| y1 (n) y2 (n)
(
(
) ≫0
) ≈ 1 for
2
1 (n)y2 (n)
σ1 σ2 1 − |ρ|
I0 2|ρ|y
σ1 σ2 (1−|ρ|2 )
(26)
数の大きさは混合比によって変化させた。複素一般化
正規分布の密度関数は以下で与えられる。
(
)c
ci exp − |xσi |
( )i
pXi (xi ) =
2πσi2 Γ c2i
(27)
ここで形状母数 ci > 0 は 2 つの変数で異なったも
数値シミュレーション
4
4.1
のを用い，c1 = 0.5, c2 = 0.8 とした。また σi2 は，
2 変量正規分布の相関係数推定
2 変量複素正規分布の疑似乱数に対する相関係数推
定の数値シミュレーションを行う。ここで各次元の分
散の推定は自明であるため、各次元の分散が 1 の複
素ガウス変数の相互相関の絶対値を 0 から 1 まで変
えて、近似最尤推定と EM アルゴリズムによる推定
を比較した。各トライアルの標本数は 100、EM アル
ゴリズムは式 (10) の絶対値の相関を初期値としての
E[|Xi |2 ] = 1 となるように調整した。これらの形状
母数は，以下のように与えられる複素カートシス [3]
がおよそ 7.4, 2.3 の場合に相当し，ガウス性の面では
音声の短時間フーリエ分析と類似したものとなって
いる。
[
]
4
2
E |x|
|E [x]| ( [ 2 ])2
E |x|
−2
Kurt[x] = ( [
])2 −
2
E |x|
イタレーション 20 回で学習した。
=
シミュレーション結果を図 1 に示す。EM アルゴリ
Γ (2/c) Γ (6/c)
2
Γ (4/c)
−2
(28)
ズムの結果は高い精度であることがわかる。近似最尤
実験結果を Fig. 2 に示す。正規乱数に対する結果
推定の結果は、低い相関のもとでは推定精度が低い
と同様に真の相関係数の周辺に推定地が分布してい
ことがわかる。また EM アルゴリズムの結果は、や
るため，非ガウス性によるモデルミスマッチの影響は
や分散が大きくなるものの，真の相関係数の周辺に
大きくないことがわかる。
推定地が分布している。
4.2
4.3
正規分布に従わない確率変数の相関係数推定
音響データの時間周波数領域の相関係数推定
次に、マイクロホンアレーで観測した音声混合の
次に，正規分布を仮定したモデルのミスマッチを評
短時間フーリエ分析の振幅スペクトルに対する相関
価するため，非ガウス的なデータを疑似乱数を用い
係数の推定を行う。観測信号のサンプリング周波数
日本音響学会講演論文集
- 737 -
2014年9月
Estimated correlation coefficient
1
[2] M. Novey, T. Adali and A. Roy, “A
Complex Generalized Gaussian Distribution—
0.8
Characterization, Generation, and Estimation,”
IEEE Trans. Signal Processing, 58(3), 2010.
0.6
0.4
[3] S. Javidi, D. P. Mandic, C. C. Took, and A. Cichocki,“ Kurtosis-based blind source extraction
of complex non-circular signals with application
0.2
in EEG artifact removal in real-time, ”Frontiers
in Neuroscience, vol.5, no.105, 25 pages, 2011.
0
0
Amplitude correlation
Proposed
0.2
0.4
0.6
0.8
Correlation coefficient estimated with angle
doi:10.3389/fnins.2011.00105
1
Fig. 3 Estimation results of correlation coeﬃcients
between two channels of speech mixture.
は 48 kHz、窓幅は 4096 サンプル、フレームシフトは
2048 サンプル、フレーム数は 80、EM アルゴリズム
は近似最尤推定を初期値としてのイタレーション 20
回で学習した。
シミュレーション結果を図 3 に示す。厳密に定常な
分布に従う疑似乱数の場合と違い、EM アルゴリズム
の推定精度は低下していることがわかる。これは非
定常性などの，さまざまな要因が関係していると思
われるが，更なる原因究明が必要である。
5
おわりに
本稿では，複素変数の絶対値のみの観測から，EM
アルゴリズムを用いて元の複素変数の相関係数を求
める手法を提案した。元の複素変数は 2 変量複素正
規分布に従うと仮定し，偏角の差を隠れ変数として
扱うことにより，EM アルゴリズムに基づく最尤推定
を定式化した。数値シミュレーションの結果，提案手
法は観測データが正規分布に従わない場合にも相関
係数をある程度の精度で推定できることを確認した。
しかし，実際の音声の観測では相関推定の精度が低
下することがわかった。原因としては観測の非定常性
などが考えられるが，更なる検討が必要である。
参考文献
[1] F. J. Lopez-Martinez, D. Morales-Jimenez, E.
Martos- Naya, J. F. Paris, “On the bivariate
Nakagami-m cumulative distribution function:
closed-form expression and applications,” IEEE
Trans. Commun., vol. 61, no. 4, pp. 1404.1414,
2013.
日本音響学会講演論文集
- 738 -
2014年9月

Download Report