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

October 16-19, 2005, New Paltz, NY

ON REAL AND COMPLEX VALUED 1 -NORM MINIMIZATION FOR OVERCOMPLETE
BLIND SOURCE SEPARATION
Stefan Winter∗, Hiroshi Sawada, Shoji Makino
NTT Communication Science Laboratories,
NTT Corporation,
2-4 Hikaridai, Seika-cho, Soraku-gun,Kyoto 619-0237, Japan{wifan,sawada,maki}@cslab.kecl.ntt.co.jp
2. OVERCOMPLETE BSS

ABSTRACT

We will consider a convolutive mixing model with N speech sources si (t) (i = 1 . . . N ) and M (M < N ) sensors that yield
linearly mixed signals xj (t) (j = 1 . . . M ). The mixing can be
described by

A maximum a-posteriori approach for overcomplete blind source
separation based on Laplacian priors usually involves 1 -norm minimization. It requires different approaches for real and complex
numbers as they appear for example in the frequency domain. In
this paper we compare a combinatorial approach for real numbers with a second order cone programming approach for complex
numbers.
Although the combinatorial solution with a proven minimum
number of zeros is not theoretically justiﬁed for complex numbers,
its performance quality is comparable to the performance of thesecond order cone programming (SOCP) solution. However, it
has the advantage that it is faster for complex overcomplete BSS
problems with low input/output dimensions.

xj (t) =N ∞hji (l)si (t − l)(1)i=1 l=1

where hji (t) denotes the impulse response from source i to sensorj.
Instead of solving the problem in the time-domain, we choose
a narrowband approach in the frequency domain by applying ashort-time discrete Fourier transform (STDFT). Following [1] we
can approximate the mixing process in the frequency domain as
X(f, τ ) = H(f )S(f, τ ) ,1. INTRODUCTION(2)

where X ∈ KM , H ∈ KM ×N , S ∈ KN , K = C and τ denotes
the time frame.
This reduces the problem from convolutive to instantaneous
mixtures in each frequency bin f . For simplicity we will omit
the frequency and time dependence. Switching to the frequency
domain has the additional advantage that the sparseness of thesources is increased [2]. This is very important to ensure good separation performance since the separation is built on the assumptionof sparse source signals.
In complete BSS the mixing matrix H is square and (assuming
full rank) invertible. Therefore the BSS problem can be solved
by either inverting an estimate of the mixing matrix or directly
estimating its inverse and solving (2) for S.
However, this approach does not work for overcomplete BSS
where the mixing matrix is not invertible. An optimal solution ina statistical sense can be obtained by maximizing the a-posteriori(MAP) P (S|X, H) of the sources S given the mixtures X and themixing matrix H.

The high quality separation of speech sources is an important prerequisite for further processing such as speech recognition in environments with several active speakers. If the underlying mixing
process is unknown it is called blind source separation (BSS). We
consider the special case of overcomplete BSS where the numberof sources exceeds the number of sensors and take into accountconvolutive mixtures for reverberant environments.
If the sources are assumed to be statistically independent and
Laplacian distributed (which is typical for sparse sources) then amaximum a-posteriori (MAP) based separation approach leads toa two-stage process. In the ﬁrst step the mixing matrix must be
estimated. It is then used in the second step to separate the sources
by constrained 1 -norm minimization. Since we work in the frequency domain to account for convolutive mixtures we have todeal with complex numbers.
In this paper we investigate the difference between real andcomplex valued 1 -norm minimization and its implication for theovercomplete BSS of convolutive mixtures.
(3)
S = [ S1 . . . SN ]T = arg max P (S|X, H).
S
In Sec. 2 we ﬁrst explain the context of overcomplete BSS in
which we consider 1 -norm minimization. In Sec. 3.1 and Sec. 3.2
If we assume source signals whose spectral components have stawe provide a detailed description of real and complex valued 1 tistically independent phases and amplitudes with uniform and Laplanorm minimization before we consider their differences. The consequen- cian distributions, respectively, the cost function results in
ces of these differences for practical applications are presented inSec. 4min|Si |, i = 1, . . . , N
s.t.
HS = X(4)S∗
The author is on leave from the Department of Multimedia Communi-

for each time instance τ . |Si | denotes the amplitude of Si .

cation and Signal Processing, University Erlangen-Nuremberg.

0-7803-9154-3/05/$20.00 ©2005 IEEEi86

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

x(t)

S
T
D
F
T

S(f, τ )X(f, τ )

BSR
BMMR

I S
n
v T
e D
r
s F
e T

s(t)

The correct placement of the zeros can be determined by combinatorially testing all possibilities and accepting the one with thesmallest 1 -norm. As a simple example let us consider
1 0.6 −0.6
1
H=
,
X=
.
(8)
0 0.8
0.8
0.5

H

According to the dimensions of the problem at least one element of
the solution S must be zero. The 1 -norm of the possible solutions
are
⎡ 
 ⎤
−1 

1
1 0.6


⎣
⎦ = 1.25
(9)
0.5
0 0.8




0
1
⎡
⎤


1
1 −0.6 −1⎣ 0⎦ = 2
(10)




0.5
0
0.8
1
⎡
⎤


0
 
 
−1 
⎣
⎦ = 1.6(11)
1
0.6 −0.6




0.5
0.8
0.8
1

Figure 1: Overall separation system
This can be interpreted as a two-stage process as utilized forexample in [1, 2, 3]. In the ﬁrst step, which we call blind mixing
model recovery (BMMR), we estimate the mixing matrix H. Several methods have already been proposed for estimating the mixing
matrix in the BMMR step. One example is hierarchical clustering
as described in detail in [1]. Once the mixing matrix is known, thesecond step, which we call blind source recovery (BSR), is basically a constrained 1 -norm minimization to estimate the sources.
We assume that we already know the mixing matrix and concentrate in the following section on BSR by 1 -norm minimization.
Before we use the inverse STDFT to obtain time-domain signals, we apply the minimal distortion principle [4] to eliminate
scaling ambiguity. The overall separation system is depicted inFig. 1.

where |·|p denotes the p -norm. The chosen solution would be theone corresponding to (9).
This combinatorial method is based on the shortest-path algorithm [3] and the 0 -norm that basically counts the number of
non-zero elements. The combinatorial method stands in contrast
to the approach in [8] where conditions are given for which the0 -norm can be calculated by an p -norm with 0 < p ≤ 1.

3. CONSTRAINED 1 -NORM MINIMIZATION
Equation (4) cannot be solved analytically for real or complex
numbers. After an appropriate transformation numerical solutions
can be found in the area of mathematical programming, which
deals in general with the optimization of constrained functions [5].

3.2. 1 -norm minimization of complex valued problems
If complex numbers are involved then (6) can no longer be applied
because such numbers possess a continuous phase in contrast toa discrete phase of real numbers. Thus we cannot use algorithms
that solve linear programming problems for complex valued problems. However, 1 -norm minimization problems (4) with complex
numbers (K = C) can be transformed to second order cone programming (SOCP) problems in the following way.
Equation (4) is equivalent to3.1. 1 -norm minimization of real valued problems
If we had to consider only real valued problems (K = R) we could
apply linear programming (LP) [6] which solves problems of the
form

S
 = X,
 Si ≥ 0, i = 1, ..., N

min 
cT S,
s.t. H
(5)








 ∈ RM ×N and X
 ∈ RM . For K = R (4) can ∈ RN , H
where c, S
be transformed into (5) by separating positive and negative values
by
 + 




1
H
 ← S− ,  ← ← X. (6)
S
c←
, H
, X
1
−H
S

min t ∈ R

i1 , . . . , iM ∈ {1, . . . , N }s.t.

X = HS, and |S|1 ≤ t(12)N
By decomposing t =i=1 ti , ti ∈ R, the second constraint|S|1 ≤ t can be expressed byHere 1 stands for an unity matrix with appropriate dimensions.
S+ and S− are derived from S by setting all negative values orpositive values, respectively, to zero.
Although powerful algorithms for linear programming exist,
they are still time consuming. Depending on the dimensions of
the problem we can obtain a faster combinatorial algorithm if we
use a certain property of the solution. It can be shown that the
N -dimensional vector S, that solves (4), contains at least N − M
zeros if the columns of H are normalized [3, 7]. The normalization
can be assumed for BSS due to the scaling ambiguity.
The lower limit for the number of zeros can be considereda constraint imposed by the MAP estimation and can easily be
fulﬁlled by setting N − M elements of the solution to zero. Then
we only have to determine the remaining M elements. Assumingthat we know where to place the zeros the remaining elements arefound by multiplying the inverse of the quadratic matrix built by
the remaining mixing vectors hi with the constraining vector X[hi1 . . . hiM ]−1 X,

October 16-19, 2005, New Paltz, NY|S|1 =

N  (Si ) T (Si )  ≤ 1t = 1[t1 · · · tN ] = t
2

(13)

i=1

where (·) and (·) denote the real and imaginary parts, respectively. Thus we can rewrite (4) as (Si )  ≤ ti , ∀is.t.
X = HS, (14)
min 1t(Si ) 2t
By deﬁning⎡⎤
t1⎢ (S1 ) ⎥⎢ (S ) ⎥
1
⎥
⎢
⎥
⎢
..
⎥ ∈ R3N ,

S=⎢
⎥
⎢
.
⎥
⎢
⎥
⎢ tN
⎣ (S ) ⎦N
(SN )

(7)

87

⎤
1
⎢ 0 ⎥
⎢ 0 ⎥
⎥
⎢
⎢ . ⎥
3N
⎥
.

c=⎢
⎢ . ⎥∈R
⎥
⎢
⎢ 1 ⎥
⎣ 0 ⎦
⎡

0

(15)

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

 =
X



(X)
(X)



October 16-19, 2005, New Paltz, NY

1.27

∈ R2M

(16)

1.265
1.26

0
⎢ (h1 )⎢ −(h )
1
⎢
⎢
..
 =⎢
H
⎢
.
⎢
0
⎢
⎣ (h )
N

−(hN )⎤T0(h1 ) ⎥(h1 ) ⎥⎥⎥..
⎥ ∈ R2M ×3N
⎥
.
⎥
0
⎥
(h ) ⎦

1 -norm

⎡

1.255
1.25

1.245
1.24

(17)

1.235
0

0.05

0.1

N

(hN )

(α)

0.15

0.2

0.25

Figure 2: Smooth slope

we can write





 s.t. X
 =H
 S,  (Si )  ≤ ti
min 
cT S


(S
)
i

S
2

∀i

This means that if we have a combinatorial solution we can limit
our search for the minimum 1 -norm solution to the nullspace
N (H), i.e.

(18)




 
min Scomb + S

The second constraint in (18) can be interpreted as a second ordercone for each i.
Equation (18) describes an SOCP problem [9], which can besolved numerically for example with SeDuMi [10].

(24)

1

with

3.3. Analysis of real and complex valued 1 -norm minimization ∈ N (H)S

 = (1 − H− H)z,S⇔z ∈ CN(25)

where H− is an arbitrary generalized inverse of H. For N = 3
and M = 2 we can express the combinatorial solution and thenullspace without loss of generality by

In contrast to the real valued 1 -norm minimization problem wherea minimum number of zeros can be guaranteed theoretically in the
optimal solution, the number of zeros cannot be predicted withcomplex valued problems as the following simple example shows.
Let


√4
1 0.6
1
17
,
X
=
.
(19)
H=
0.8+j0.6
√
0.5
0 0.8


Scomb =[ h1h2 ]−1 X
0


(26)

17



Then the 1 -norm of the solution obtained by SOCP is given by 0.227 + 0.040i |Ssocp |1 =  0.511 − 0.091i  = 1.23.
(20)
 0.481 + 0.015i 

N (H) = α1

It does not contain any zeros as we could expect with real numbers
yet solves (4). In comparison, the 1 -norm of the optimal combinatorial solution is given by
⎡
⎤


0

−1 
 ⎥
⎢
√4
0.6
|Scomb |1 = ⎣1
⎦ = 1.24 (21)
17


√
0.5
0.8 0.8+j0.6


17
1

This observation reveals a very important difference from realvalued problems and prevents the theoretical justiﬁcation of a procedure similar to the combinatorial approach in Sec. 3.1. To better
explain this difference between real and complex numbers we takea look at a general solution based on a combinatorial solution andthe nullspace N (H) of H.
Even though the combinatorial solution Scomb does not necessarily minimize the 1 -norm, it fulﬁlls together with the SOCP
solution Ssocp
X = HScomb = HSsocp .

(22)

 = Ssocp − Scomb , (22) becomes
By deﬁning the difference S.
HScomb = HScomb + 
HS

(23)

=0

88

[ h1

h2 ]−1 h3
1


,

α∈C

(27)

With (26) and (27), the function to be minimized (24) can be written as



 
= |f11 (H, X) + αf12 (H, X)|(28)
Scomb + αS
1+|f21 (H, X) + αf22 (H, X)|+|f31 (H, X) + αf32 (H, X)|.
Here fij is a summand that only depends on H and X, which areconstant for any given problem. If only real values are involved
then (28) describes a piecewise linear function depending on α
whose slope can only change a limited number of times in a discrete manner.
However, once complex numbers are involved, their imaginary
part results in an inherent 2 -norm, which leads to smooth slopes
as they appear with second order or higher polynomials. This behavior becomes obvious in (18). There the 1 -norm is changed
from the sum of absolute values of real numbers to the sum ofthe 2 -norms of the real and imaginary part. The introduction ofthe 2 -norm explains the different behavior of complex valued 1 norm minimization compared with its real counterpart. An example is shown in Fig. 2, where the dependence of the 1 -norm
on α is shown (here only the dependence on the real part of α isshown). The combinatorial solution that minimizes the 1 -norm
is given there for α = 0. However, this is not the solution of (4),which is rather obtained for α = 0.

Table 2: Detailed performance comparison
SignalOriginal
Estimated
CS
SOCP
CS
SOCP
combination
1
0.938 0.913
0.824 0.895
0.958 0.948
0.895 0.876
20.898 0.955
0.842 0.891
3
0.919 0.912
0.847 0.876
4
Average
0.928 0.932
0.852 0.884

SOCP

5

10

15

Number of sources

20

4 mixtures

2

10

0

10

Combinatorial
SOCP
5

10

15

Number of sources

20

0

10

Combinatorial
SOCP
5

10

15

Number of sources

20

5 mixtures

2

10

0

10

Combinatorial
SOCP
5

10

15

Number of sources

20

5. CONCLUSION
Although the combinatorial solution with at least N − M zeros
is not theoretically justiﬁed for complex numbers, its performance
quality is comparable to that of the SOCP solution. It has the advantage that it is faster for overcomplete BSS problems with low
input/output dimensions.

Even though the combinatorial solution (CS) with a minimum number of zeros in Sec. 3.1 cannot be justiﬁed theoretically for complex numbers, in practice its performance is comparable to that of
the SOCP solution. We tested both approaches with the estimatedas well as the original mixing matrix and combinations of threeout of four different speech signals. The experimental conditions
can be found in Table 1. As a performance measure we used thecorrelation coefﬁcientτCombinatorial
0

10

3 mixtures2

10

Figure 3: Comparison of computational complexity

4. PERFORMANCE COMPARISON

xTsi (t)yi (t + τ )|xsi (t)|2 |yi (t)|2

2 mixtures2

10

Computation time (sec) Computation time (sec)Table 1: Experimental conditions
Source directions
50◦ , 90◦ , 120◦
Sensor distance
40 mm
7.4 seconds
Source signal length
Reverberation time TR 130 ms
Sampling frequency fs 8 kHz
von Hann
Window type
Filter length1024 points
256 points
Shifting interval

cor = max

October 16-19, 2005, New Paltz, NY

Computation time (sec) Computation time (sec)2005 IEEE Workshop on Applications of Signal Processing to Audio and Acoustics

6. REFERENCES[1] S. Winter, H. Sawada, S. Araki, and S. Makino, “Overcomplete BSS for convolutive mixtures based on hierarchicalclustering,” in Proc. ICA 2004, Sept. 2004, pp. 652–660.
[2] P. Boﬁll and M. Zibulevsky, “Blind separation of more
sources than mixtures using sparsity of their short-time
fourier transform,” in Proc. ICA 2000, June 2000, pp. 87–
92.
[3] F. Theis, “Mathematics in independent component analysis,” Ph.D. dissertation, University of Regensburg, 2002.
[Online]. Available: http://homepages.uni-regensburg.de/˜thf11669/phdthesis.html[4] K. Matsuoka and S. Nakashima, “Minimal distortion principle for blind source separation,” in Proc. ICA 2001, Dec.
2001, pp. 722–727.
[5] [Online]. Available: http://www.mathprog.org/[6] M. S. Lewicki and T. J. Sejnowski, “Learning overcomplete
representations,” Neural Computation, vol. 12, no. 2, pp.
337–365, 2000. [Online]. Available: citeseer.nj.nec.com/
lewicki98learning.html[7] I. Takigawa, M. Kudo, and J. Toyama, “Performance analysisof minimum 1 -norm solutions for underdetermined source
separation,” IEEE Trans. Signal Processing, vol. 52, no. 3,
pp. 582–591, Mar 2004.
[8] D. M. Malioutov, M. Cetin, and A. S. Willsky, “Optimalsparse representations in general overcomplete bases,” inProc. ICASSP 2004, 2004, pp. 793–796.
[9] A. Pruessner, M. Bussieck, S. Dirkse, and A. Meeraus,“Conic programming in GAMS,” INFORMS Annual
Meeting, Atlanta, Oct 2003. [Online]. Available: http://www.gams.com/presentations/present˙conic.pdf[10] J. Sturm, “Using SeDuMi 1.02, a MATLAB toolbox for
optimization over symmetric cones,” Optimization Methods
and Software, vol. 11-12, pp. 625–653, 1999, special
issue on Interior Point Methods. [Online]. Available:
http://fewcal.kub.nl/sturm/software/sedumi.html(29)

between the output signal yi (t) and the corresponding ﬁltered version xsi of the original speech signal at the microphones. τ determines the shift of the signal vector. The results are shown in Table
2. A subjective evaluation of the separated sources supports thisresult.
As expected, the SOCP solution yields better results than the
combinatorial solution. However, better estimation of the mixing
matrix reduces the advantage of the SOCP solution.
While the difference in performance quality might be negligible in practical applications, the computational complexity reveals great differences. The combinatorial solution has a low initial computational complexity but grows exponentially with the input dimension N . On the other hand, the SOCP solution has a highcomputational complexity even for low input dimensions N , but it
grows only according to√
N log(1/) .
(30)
O denotes the precision of the numerical algorithm [10]. Figure
3 illustrates this fact and shows on a logarithmic scale the timerequired by the two approaches to separate the sources in one frequency bin for different numbers of sources and sensors.
One reason for the big difference in the initial computational
complexity can be found in the reusability of previous results. For
overcomplete BSS in the frequency domain the minimum 1 -norm
solution must be calculated several times with the same mixing
matrix. SOCP always starts from scratch and therefore cannot
proﬁt from the reuse of earlier results. In contrast, the combinatorial solution is built on the inverses of selected mixing vectors.
Once they are calculated they can be reused as long as the mixing
matrix does not change.

89