2013 IEEE Workshop on Applications of Signal Processing to Audio and Acoustics

October 20-23, 2013, New Paltz, NY

EMPLOYING MOMENTS OF MULTIPLE HIGH ORDERS FOR HIGH-RESOLUTIONUNDERDETERMINED DOA ESTIMATION BASED ON MUSIC
Yuya Sugimoto1 , Shigeki Miyabe1 , Takeshi Yamada1 , Shoji Makino1 and Fred Juang2
1

2

University of Tsukuba, 1-1-1 Tennodai, Tsukuba, Ibaraki, 305-8577 Japan
Georgia Institute of Technology, 777 Atlantic Dr. NW. Atlanta, GA 30332-0250 USA
Email: sugimoto@mmlab.cs.tsukuba.ac.jp, {miyabe, makino}@tara.tsukuba.ac.jp,
takeshi@cs.tsukuba.ac.jp, juang@ece.gatech.edu
ABSTRACT

Several extensions of the MUltiple SIgnal Classification (MUSIC)
algorithm exploiting high order statistics were proposed to estimate
directions of arrival (DOAs) with high resolution in underdetermined conditions. However, these methods entail a trade-off between two performance goals, namely, robustness and resolution, inthe choice of orders because use of high-ordered statistics increases
not only the resolution but also the statistical bias. To overcome
this problem, this paper proposes a new extension of MUSIC using
a nonlinear high-dimensional map, which corresponds to the joint
analysis of moments of multiple orders and helps to realize the bothadvantages of robustness and high resolution of low-ordered andhigh-ordered statistics. Experimental results show that the proposed
method can estimate DOAs more accurately than the conventionalMUSIC extensions exploiting moments of a single high order.
Index Terms— MUSIC, Higher order statistics, Microphonearray, Underdetermined DOA estimation, Subspace analysis
1. INTRODUCTION
In microphone array signal processing, direction of arrival (DOA)
estimation is essential in various applications such as speaker tracking and preprocessing for signal enhancement. MUltiple SIgnal
Classification (MUSIC) [1] is a representative high-resolution DOA
estimation method based on the subspace analysis of the observed
signals. However, the effectiveness of MUSIC is conditioned on
the dimensionality of the covariance matrix because it requires theestimation of the noise subspace that is orthogonal to the transfer
function vectors of the observed signals. With N sound sources,
M (> N ) sensors are required and the DOA estimation performance degrades significantly as N approaches M .
For high-resolution DOA estimation in underdetermined conditions where the number of sources exceeds the number of microphones, several extensions of MUSIC were developed. In the early
90’s, extended MUSIC exploiting fourth-order cumulants was proposed [2], and recently this approach was generalized to cumulants
analysis of arbitrary even order, referred to as 2q-MUSIC [3]. Thesemethods improve the resolution of DOA estimation and can also estimate DOAs in underdetermined conditions by increasing the dimensionality of observation with higher-order cross cumulants. Asan alternative extension, we proposed the mapped MUSIC [4] using a high-dimensional nonlinear map of the observation, which isinspired by the nonlinear beamformer [5, 6, 7] for speech enhancement in a similar underdetermined condition. We showed that theThis work was supported by KAKENHI 23760325.

978-1-4799-0972-8/13/$31.00 ©2013IEEE

mapped MUSIC corresponds to the analysis of arbitrary even order cross moments. The mapped MUSIC achieves a performancein DOA estimation comparable to that of conventional extensionsof MUSIC based on cumulant analysis, but offers an advantage incomputational complexity due to its efficient algorithm in dimensionality expansion.
In this paper, we propose a method for joint analysis of moments of multiple orders to further improve the resolution of DOA
estimation in the MUSIC framework. Although use of higher order
statistics yields a higher resolution in the DOA estimate, the resultis less robust due to the greater bias in the higher-order statistics, as
examined in [4]. Thus the question of a strategy for trade-off between resolution and robustness arises. To explore the possibility ofrealizing both advantages of the high and low orders of statistics, weevaluate moments of multiple orders simultaneously. This goal was
found to be possible and can be accomplished by properly combining low-and-high-ordered moments in a nonlinearly expanded
space, as inspired by the mapped MUSIC method. Experimental results of speech DOA estimation demonstrate that the mapped MUSIC with the proposed map achieves a higher accuracy than other
MUSIC extensions.
2. PROBLEM STATEMENT
Throughout this paper, signals are expressed as complex amplitudes
in the time-frequency domain. Observed signals can be modeled asx(ω, t) = [x1 (ω, t), · · · , xM (ω, t)]T=N∑ai (ω)si (ω, t) + n(ω, t),(1)n(ω, t) = [n1 (ω, t), · · · , nM (ω, t)]T ,(2)ai (ω) = [a1,i (ω), · · · , aM,i (ω)] ,

(3)

i=1

T

where t = 1, · · · , L is a time frame index, ω is the angular frequency, M is the number of sensors, N is the number of soundsources, [·]T denotes transpose, si (ω, t) is the complex amplitudeof the i-th sound source, xj (ω, t) is the complex amplitude of thesignal observed with the j-th sensor, nj (ω, t) is the complex amplitude of the noise observed with the j-th sensor, and aj,i (ω) denotes
the transfer function from the i-th source to the j-th sensor.
In the signal model expressed by (1), each sensor observes a
mixture of source signals and noise. The problem in this paper is
to estimate the DOAs of the source signals s1 (ω, t), · · · , sN (ω, t)
by analyzing these observed signals and finding steering vec-

2013 IEEE Workshop on Applications of Signal Processing to Audio and Acoustics

tors b(ω; θ1 ), · · · , b(ω; θN ) similar to the transfer function vectors a1 (ω), · · · , aN (ω) from a set of steering vectors b(ω; θ)(||b(ω; θ)||2 = 1) of arbitrary directions θ.
1b(ω; θ) = √ [exp(−jωτ1 ), · · · , exp(−jωτM )] ,M

(4)where τi (i = 1, · · · , M ) are time delays from a reference point at
each sensor, given by a signal from direction θ.

(i = 1, · · · , N

3.1. DOA estimation algorithm of mapped MUSICThis section briefly reviews the DOA estimation with mapped MUSIC [4]. Mapped MUSIC maps the M -dimensional observed vector x(ω, t) onto M 0 -dimensional Euclidian space (M 0 ≥ M ) with0
a nonlinear map φ : CM → CM and conducts a similar analysis to MUSIC with mapped steering vectors φ(b(ω; θ)). To estimate DOAs properly with mapped MUSIC, the information about
the correlations between each microphone signal must be retained
after mapping. This requirement for the map φ is given by the following three conditions:
1. The magnitude relation of the norm is retained.
||x||2 ≥ ||y||2 → ||φ(x)||2 ≥ ||φ(y)||2
2. The origin remains intact.

(5)x = 0 → φ(x) = 0
3. The orthogonality between vectors is preserved.

(6)

where E[·] denotes the expectation of the argument, and [·]H denotes a complex conjugate transpose. The following equations describe the eigen decomposition of covariance matrix R(ω):H

R(ω) = V(ω)E(ω)V (ω),(9)

V(ω) = [v1 (ω), · · · , vM 0 (ω)], V (ω)V(ω) = IM 0 ,E(ω) = diag[e1 (ω), · · · , eM 0 (ω)],
e1 (ω) ≥ · · · ≥ eM 0 (ω),
H0

M = dim[φ(x(ω, t))],(12)

span[φ(ai (ω))|i = 1, · · · , N ] ' span[vj (ω)|j = 1, · · · , N 0 ],(13)
where span[·] denotes the space spanned by argument vectors, we
can define the noise subspace N (ω) in the mapped space;
(14)

(15)

|φH (b(ω; θ))v

k (ω)|

2.

(16)

Since the MUSIC score f (ω; θ) is defined for each narrowband,
it is necessary to integrate the information from all frequency bins
in order to obtain a broadband DOA estimation. Among various
average operations to integrate the narrowband MUSIC scores [8],we use the following geometric mean in this paper;
[
∏

]1

K

f (ω; θ)

,

(17)

ω

where K denotes the number of averaged frequency bins. EstimatedDOAs are given by finding peaks of f (θ).
3.2. Map for analysis of arbitrary even order momentsd

In [4], we proposed a map φd : CM → CM for an analysis ofcross moments of 2d-th order:
φd (x(ω, t)) ,

[ d
∏
{

xl⊗ ,

where vi (ω) (i = 1, · · · , M 0 ) are the eigenvectors associated witheach eigenvalue ei (ω) (i = 1, · · · , M 0 ) respectively, Ii denotes thei-dimensional identity matrix, diag[·] is a diagonal matrix consisting of elements within the argument vector, and dim[·] is the dimensionality of the argument. With N 0 denoting the number of theeigenvectors corresponding to the mapped directional source signals;

N (ω) , span[vk (ω)|k = N 0 + 1, · · · , M 0 ].

k=N 0 +1

(10)
(11)

1

f (ω; θ) , ∑M 0

(7)

(8)

k = N 0 + 1, · · · , M 0 ).

Since we assume b(ω; θi ) ∼= ai (ω) (i = 1, · · · , N ), we can findthe true sound source directions by searching the orthogonal projection of mapped steering vectors φ(b(ω; θ)) onto the mapped noise
subspace. Here, we define the following MUSIC score of mapped
MUSIC f (ω; θ), which is locally maximized at the direction θi :f (θ) ,We describe the DOA estimation algorithm of mapped MUSIC with the map φ satisfying (5)–(7). The covariance matrix ofφ(x(ω, t)) is expressed as
R(ω) = E[φ(x(ω, t))φH (x(ω, t))],Note that the value of N 0 depends on the property of the map φand the number of the source signal N in the original observation.
Although the orthogonality between the mapped transfer functions
φ(ai (ω)) and the mapped noise subspace N (ω) is retained by theconditions (5)–(7), the dimensionality N 0 of the mapped signal subspace is increased by the higher-dimensional map or the greaternumber of sources N . Thus we must choose N 0 appropriately tosatisfy the following condition:|φH (ai (ω))vk (ω)|2  1

3. MAPPED MUSICxH y = 0 → φH (x)φ(y) = 0

October 20-23, 2013, New Paltz, NYxl⊗
c1l (ω, t), · · · ,

l=1

x
x∗

d
∏

]T
xl⊗
cM d l (ω, t)

,

(18)

l=1

(if l is odd),
(if l is even)(19)

where [·]∗ denotes a complex conjugate. ckl is an index specifying the element index of vector x for calculating φd , and this is
the (k, l) element of the M d × d index matrix C that arranges drepeated-permutations of M in an arbitrary order as its row vectors:
1 1 ··· 1
 1 1 ··· 2 C , [ckl ]kl = 
,
..


.
M M ··· M(20)where [·]ij denotes a matrix consisting of the argument as its (i, j)element. We showed in our previous work [4] that the mapped MUSIC with the map φd is as accurate as the conventional MUSIC extensions based on cumulants analysis with less computational complexity.

2013 IEEE Workshop on Applications of Signal Processing to Audio and Acoustics

October 20-23, 2013, New Paltz, NY

4. PROPOSED METHOD

φd1 ,··· ,dm , φd1 ⊕ · · · ⊕ φdm ,{φd1 , · · · , φdm } = {φd | d = d1 , · · · , dm }di+1 > di (i = 1, · · · , m − 1),

(21)
(22)

where ⊕ denotes a direct sum. Mapped MUSIC with the map
φd1 ,··· ,dm gives the joint analysis of multiple moments because
the covariance matrix R⊕
d1 ,··· ,dm (ω) of φd1 ,··· ,dm contains the moments of multiple orders as
[
]
0
H
0
R⊕
d1 ,··· ,dm (ω) = E φd1 ,··· ,dm (x (ω, t))φd1 ,··· ,dm (x (ω, t))
Rd1 ,d1 (ω) Rd1 ,d2 (ω) · · · Rd1 ,dm (ω)
 Rd2 ,d1 (ω) Rd2 ,d2 (ω) · · · Rd2 ,dm (ω) 


=
,
..
..
..
..


.
.
.
.
Rdm ,d1 (ω) Rdm ,d2 (ω) · · · Rdm ,dm (ω)
(23)
[
]
0
Rdi ,dj (ω) , E φdi (x0 (ω, t))φH
dj (x (ω, t)) = [rpq ]pq(i, j = 1, · · · , m p = 1, · · · , M diq = 1, · · · , M ),∗ ()  dj
di
∏
∏
 xl⊗
 .
rpq = E 
xl⊗
cpl
cql
dj

l=1

Rotation direction

Window

20
95

150

20

340

Door

170

10

70

To utilize the advantages of high-resolution estimation with highordered moments and robust estimation with low-ordered moments,
we propose the new map for joint analysis of cross moments ofmultiple orders. The proposed map φd1 ,··· ,dm is given by a direct
sum of maps φd1 , · · · , φdm as

200

Room size
Reference point of angle 340 x 400 x 270
is displaced by 5∘
400

Figure 1: Recording environment.

map φd1 ,··· ,dm in each frequency bins as
x0 (ω, t) =wx(ω, t)
√ [
],
∑M
2dm
1
2dm
E M
i=1 |xi (ω, t)|(26)

where w is a parameter to adjust the L2dm -norm of the observed
signal. With large value of w, the proposed method analyzes thehigh-ordered moments more significantly, and vice versa.
5. EXPERIMENT

(24)
(25)

l=1

where x0 (ω, t) is an observed signal normalized properly (normalization rule is described later), and submatrices Rdi ,dj (ω)
gives (di + dj )-thcross moments. The dimensionality of∑order
m
di
R⊕
(ω)
is
M
, which is still larger than that of thed1 ,··· ,dm
i=1
highest-dimensional single cross moment matrix Rdm ,dm (ω). In
addition, by analyzing the low-ordered moments together with thehigh-ordered moments, the significance of the statistical bias in
the analysis of high-ordered moments is relaxed. With the increased dimensionality and reduced significance of statistical bias,
the mapped MUSIC with the proposed map φd1 ,··· ,dm achieves further improvement of DOA estimation performance.
Note that we must take care of the scaling of the observed signal for the sufficient joint analysis of the moments of multiple orders. As we find in (23)–(25), R⊕
d1 ,··· ,dm (ω) contains moments ofvarious orders, and it is obvious that the absolute values of the elements in the high-ordered moment matrices becomes significantly
larger than those of the low-ordered moment matrices when thenorm ||x(ω, t)|| of the observation is large. In contrast, the absolute values of the elements in the high-ordered moments becomesconsiderably small when the norm ||x(ω, t)|| is small. Since theeigen decomposition (9)–(12) takes the moments of the orders withlarge values into account more significantly, all the orders of moments must be weighted appropriately for the sufficient joint analysis of moments of multiple orders. Since such the proportion isdefined by the norm, we normalize the observed signal with theL2dm -norm based on the maximum degree dm within the proposed

This section describes an experiment to verify the effectiveness ofthe proposed method. We evaluated mapped MUSIC with the proposed map φ1,2,3 , MUSIC, mapped MUSIC (d = 2, 3) and 2qMUSIC (q = 2, 3) are implemented for comparison.
5.1. Experimental condition
We conducted an evaluation of the DOA estimation of each methodby using mixed signals that were created as convolutive mixturesof the Japanese speeches and the impulse responses, measured at
the environment shown in Fig.1. To emulate noisy observation, diffused pink noises [9] was superimposed to the observed signals. To
verify estimation performance for both overdetermined and underdetermined environments, we also evaluated the performance withtwo different number of sources. Table 1 shows the other experimental conditions.
Table 1: Experimental conditions
Sensor array configuration
Sound sources# of sound sources
3, 5# of sensors
4
Room size
Reverberation
Speech duration
Frame length
Frame shift length
Window function

Circular array with radius of 0.1 m
Speakers 1.5 m apart from array
Sampling frequency
16 kHzSNR
20 dB
3.4 m × 4 m × 2.7 m
T60 of 0.3 s
1s
512 samples
256 samples
Hamming window2013 IEEE Workshop on Applications of Signal Processing to Audio and Acoustics

20

15
RMSE [degree]

For the evaluation, we employed the root mean squared error(RMSE) between the estimated direction and the true sound direction. We evaluated 100 combinations of the positions of three or fivesources selected randomly from the directions {0◦ , 30◦ · · · , 330◦ }.
When the estimation score had fewer peaks than sources, we added
a penalty equal to the average error for all directions. This paper
employed N 0 = 15 in (13)–(16) as the number of source is three,and N 0 = 30 for five sources. We also employed normalization
parameter w = 2 in (26) for both of the number of sources. DOA
estimation with MUSIC when M ≤ N is performed by regarding
the one-dimensional subspace associate with the minimum eigenvalue as the noise subspace in every frequency bin.

October 20-23, 2013, New Paltz, NY

6. CONCLUSION
This paper proposed extended mapped MUSIC realizes to estimate
DOAs with higher resolution than other conventional MUSIC extensions. We reviewed mapped MUSIC with the map for moments
analysis of single arbitrary even order. We proposed new map as adirect sum of maps of multiple degrees, and showed that the proposed method corresponds to joint moments analysis of multipleorders taking advantage of both resolution with high-ordered moments and robustness with low-ordered moments. The experimentcomparing with the other MUSIC extensions based on single highorder statistics clarified the effectiveness of the joint analysis andexcellent DOA estimation performance by the proposed method.
7. REFERENCES[1] R. O. Schmidt, “Multiple emitter location and signal parameter estimation,” IEEE Trans. Antennas Propag., vol.AP-34,
no.3, pp.276–280, 1986.
[2] B. Porat and B. Friedlander, “Direction finding algorithms
based on higher order statistics,” IEEE Trans. Signal Process.,
vol.39, no.9, pp.2016–2024, 1991.
[3] P. Chevalier, A. Ferreol and L. Albera, “High resolution direction finding from higher order statistics: The 2q-MUSIC algo-

Mapped MUSIC
⣔ิ2
(d=2)
Mapped MUSIC
26.5
(d=3)
2q-MUSIC
27.4

(q=2)
5

2q-MUSIC
29.9
(q=3)

0

Proposed method⣔ิ6

Figure 2: Experimental results for 3 sources within 1 second.

5.2. Experimental results

40

30
RMSE [degree]

Figures 2–3 show the experimental results obtained under different number of sources. Throughout all the conditions, the proposed
method performs the best, and mapped MUSIC performs slightly
better than 2q-MUSIC. When we compare the results of mapped
MUSIC (d = 2) and that of (d = 3), the former shows the better
performance, and this tendency also appears as to 2q-MUSIC. As
examined in [3], these results suggest that the latter is biased greaterthan the former because of the shortage of snapshots, and the observation of one second is not enough to reduce the influence of statistical bias derived from the high-ordered statistics. While, the comparison between the proposed method and mapped MUSIC (d = 2)
denotes the effectiveness of joint analysis of moments of multipleorders. Although the high-ordered moments is inferior as for thesole use, we can confirm its usefulness to support the low-orderedmoments. Furthermore, the proposed method can be expected tokeep showing the best performance even if the time length becomeslarger because the information of the high-ordered becomes morerobust. From these discussion, the effectiveness of the proposed
method utilizing the advantage of low-ordered and high-orderedmoments simultaneously is verified.

10

MUSIC
MUSIC

20

10

MUSIC
⣔ิ1
Mapped MUSIC
22.6
(d=2)
Mapped MUSIC
26.5
(d=3)
2q-MUSIC
27.4(q=2)
2q-MUSIC
29.9
(q=3)

0

Proposed method⣔ิ6

Figure 3: Experimental results for 5 sources within 1 second.
rithm,” IEEE Trans. Signal Process., vol.54, no.8, pp.2986–
2997, 2006.
[4] Y. Sugimoto, et al., “Underdetermined DOA estimation by the
non-linear MUSIC exploiting higher-order moments,” Proc.
IWAENC, 2012.
[5] H. Saruwatari, et al., “Speech Enhancement Using Nonlinear Microphone Array Based on Complementary Beamforming,” IEICE Trans. Fundam., vol.E82-A, no.8, pp.1501–1510,
1999.
[6] H. Kamiyanagida, et al., “Direction of Arrival Estimation Using Nonlinear Microphone Array,” IEICE Trans. Fundam.,
vol.E84-A, no.4, pp.999-1010, 2001.
[7] S. Miyabe, et al., “Kernel-based nonlinear independent component analysis for underdetermined blind source separation,”
Proc. ICASSP, pp.1641–1644, 2009.
[8] M. Wax, T.-J. Shan and T.Kailath, “Spatio-temporal analysis
by eigenstructure methods,” IEEE Trans. Acoust. Speech Signal Process., vol.32, pp.817–827, 1984.
[9] I. A. McCowan and H. Bourlard, “Microphone array postfilter based on noise field coherence,” IEEE Trans. Speech Audio Process., vol.11, no.6, pp.709–716, 2003.