On Calculating the Inverse of Separation Matrix
in Frequency-Domain Blind Source Separation
Hiroshi Sawada, Shoko Araki, Ryo Mukai, and Shoji Makino
NTT Communication Science Laboratories, NTT Corporation,
2-4 Hikaridai, Seika-cho, Soraku-gun, Kyoto 619-0237, Japan{sawada, shoko, ryo, maki}@cslab.kecl.ntt.co.jp

Abstract. For blind source separation (BSS) of convolutive mixtures,the frequency-domain approach is eﬃcient and practical, because theconvolutive mixtures are modeled with instantaneous mixtures at eachfrequency bin and simple instantaneous independent component analysis(ICA) can be employed to separate the mixtures. However, the permutation and scaling ambiguities of ICA solutions need to be aligned toobtain proper time-domain separated signals. This paper discusses theidea that calculating the inverses of separation matrices obtained by ICA
is very important as regards aligning these ambiguities. This paper also
shows the relationship between the ICA-based method and the timefrequency masking method for BSS, which becomes clear by calculating
the inverses.

1

Introduction

With acoustical applications of blind source separation (BSS), such as solving
a cocktail party problem, signals are generally mixed in a convolutive mannerwith reverberations. Let s1 , . . . , sN be N source signals and x1 , . . . , xM be M
sensor observations. Then, the convolutive mixture model is formulated asxj (t) =

N 

k=1

hjk (l) sk (t − l), j = 1, . . . , M,(1)lwhere t represents time and hjk (l) represents the impulse response from source
k to sensor j. If we consider sounds mixed in a practical room situation, impulse
responses hjk (l) tend to have hundreds or thousands of taps even with an 8 kHz
sampling rate. This makes the convolutive BSS problem very diﬃcult comparedwith the BSS of simple instantaneous mixtures.
A practical approach for such convolutive mixtures is frequency-domain BSS[1-10], where we apply a short-time Fourier transform (STFT) to the sensor
observations xj (t). Then, the convolutive model (1) can be approximated as an
instantaneous mixture model at each frequency:
xj (f, t) =N

hjk (f )sk (f, t), j = 1, . . . , M,

k=1

J. Rosca et al. (Eds.): ICA 2006, LNCS 3889, pp. 691–699, 2006.
c Springer-Verlag Berlin Heidelberg 2006


(2)

692

H. Sawada et al.

where f represents frequency, t is now down-sampled with the distance of theframe shift, hjk (f ) is the frequency response from source k to sensor j, andsk (f, t) is a frequency-domain time-series signal of sk (t) obtained with an STFT.
The vector notation of (2) isx(f, t) =N

hk (f )sk (f, t),

(3)

k=1

where x = [x1 , . . . , xM ]T and hk = [h1k , . . . , hM k ]T .
Once we assume instantaneous mixtures at each frequency, and also if thenumber of sources N does not exceed the number of sensors M , we can apply standard instantaneous independent component analysis (ICA) [11] to the
mixtures x(f, t) to obtain separated frequency components:y(f, t) = W(f ) x(f, t),(4)where y = [y1 , . . . , yN ]T is the vector of separated frequency components andW = [w1 , . . . , wN ]H is an N × M separation matrix [1-6]. However, the ICA
solution has the permutation and scaling ambiguities. We need to align theseambiguities to obtain proper time-domain separated signals.
Various studies have tried to solve these permutation and scaling problems
because they constitute a critical issue. Some of these studies have attempted
to solve the problems by using information obtained from (4), i.e the separation
matrix W and/or the separated signals y, or by imposing some constraints on
W. By contrast, we believe that the inverses of separation matrices W provide
useful information for solving these problems. The main topic of this paper is to
discuss the procedures for solving these problems by calculating the inverses.
There is also a frequency-domain BSS method that is based on time-frequency(T-F) masking [7-10]. It does not employ a standard ICA to separate the mixtures, and can be applied even if the number of sources N exceeds the numberof sensors M . The method relies on the sparseness of source signals. It classiﬁes
the mixtures x(f, t) based on spatial information extracted from them. As thesecond topic of this paper, we show a link between the ICA-based method andthe T-F masking method. The link becomes clear once we have the decomposition (6) of mixtures by calculating the inverse of W. Based on the link, we seethat some of the techniques used in solving the permutation problem can also
be used for classifying the mixtures in the T-F masking method, and vice versa.

2

Calculating the Inverses of Separation Matrices

Figure 1 shows the ﬂow of ICA-based frequency-domain BSS that we consider
in this paper. The inverse of separation matrix W is represented as[a1 , · · · , aN ] = W−1 , ai = [a1i , . . . , aM i ]T ,(5)

On Calculating the Inverse of Separation Matrix in Frequency-Domain BSS

693

which we call basis vectors obtained by ICA, because the mixture x(f, t) isrepresented by their linear combination by multiplying W−1 and (4):x(f, t) =Nai (f )yi (f, t).

(6)

i=1

The basis vectors provide the key information with which to solve the permutation and scaling problems as shown in the following sections. If W is not square,
we use the Moore-Penrose pseudoinverse instead of the inverse. It is not diﬃcult
to make W invertible by using an appropriate ICA procedure, such as whitening
followed by unitary transformation (e.g. FastICA [11]).

3

Solving the Permutation Problem

Various methods have been proposed for solving the permutation problem:
1. making the separation matrices W(f ) smooth along frequencies f [1, 2],
2. maximizing the correlation of separated signal envelopes |yi | [3],3. analyzing the directivity patterns calculated from W(f ) [4],4. manipulating basis vectors ai (f ) [5, 6].
The third and fourth methods utilize the same information because W and ai
are related to each other through the inversion (5). However, the fourth method
is easier to apply when there are more than two sources [5, 6], because basis
vectors ai directly represent estimated mixing system information (6). Thissection describes how to utilize this information for solving the permutation
problem.
3.1

Assumption and Basic Idea

If ICA works well, we obtain separated signals yi (f, t) that should be close tosource signals sk (f, t) up to the permutation and scaling ambiguities. If we compare (3) and (6), we see that the basis vectors ai (f ), which are obtained by ICA
and the subsequent inversion of W, should be close to the vectors hk (f ), again,
up to the permutation and scaling ambiguities. The use of diﬀerent subscripts,
i and k, indicates the permutation ambiguity.

STFT

ICA

ISTFT

Inverse

SourcePermutation& Scaling

SensorFig. 1. Flow of ICA-based frequency-domain BSS

Fig. 2.
model

Direct-path

694

H. Sawada et al.

The method presented here assumes a direct-path model (Fig. 2) for the vector
hk = [h1k , . . . , hM k ]T , even though in reality signals are mixed in a multi-pathmodel (1). This simpliﬁed model is expressed in the frequency domain:hjk (f ) ≈ λjk · e−j 2πf τjk ,(7)

where τjk and λjk ≥ 0 are the time delay and attenuation from source k tosensor j, respectively. Since we cannot distinguish the phase (or amplitude) of
sk (f, t) and hjk (f ), these two parameters can be considered to be relative (thisfact causes the scaling ambiguity). Thus, without loss of generality, we normalize
them and align the scaling ambiguity byτjk = (djk − dJk )/c,M

2
j=1 λjk = 1,

(8)
(9)

where djk is the distance from source k to sensor j (Fig. 2), and c is the propagation velocity of the signal. Normalization (8) makes τJk = 0, i.e. the relative
time delay is zero at a reference sensor J.
By following the normalizations (8) and (9), the scaling ambiguity of basis
vectors ai is aligned by the operationai ←ai −j arg(aJi )e||ai ||(10)

which makes arg(aJi ) = 0 and ||ai || = 1. Now, the task as regards the permutation problem is to determine a permutation Πf that relates the subscript i andk with i = Πf (k), and to estimate parameters τjk , λjk that make the model (7)match the aji (f ) element of a basis vector. This can be formulated so as to ﬁnd
Πf , τjk and λjk that minimize the cost function:
J =

N 
M


|aji (f ) − λjk · e−j 2πf τjk |2 , i = Πf (k),(11)

f ∈F k=1 j=1

where F is the set of frequencies that we have to consider.
3.2

Clustering Frequency-Normalized Basis Vectors

If we consider the frequency range where spatial aliasing does not occur:
F = {f : −π < 2πf τjk < π, ∀ j, k}(12)

we can introduce the frequency normalization technique [6] to minimize the costfunction (11). Let dmax be the maximum distance between the reference sensor
J and any other sensor. Then the relative time delay is bounded by
max |τjk | ≤ dmax /c
jk

(13)

On Calculating the Inverse of Separation Matrix in Frequency-Domain BSS

and therefore the frequency range F can be expressed withc
F = {f : 0 < f <
}.
2dmax

695

(14)

The frequency normalization technique [6] removes frequency dependence fromthe elements of scale-normalized basis vectors (10):arg[aji (f )]āji (f ) ← |aji (f )| exp j.
(15)
4f c−1 dmax
The rationale of dividing the argument by 4f c−1 dmax is discussed in [6]. With
this operation, the cost function (11) is converted to
J¯ =

N 
M


|āji (f ) − h̄jk |2 , i = Πf (k)(16)

f ∈F k=1 j=1

where

π c · τjk·]
2 dmaxis a frequency-normalized model. In a vector notationh̄jk = λjk · exp[−j

J¯ =

N


||āi (f ) − h̄k ||2 , i = Πf (k),

(17)

(18)

f ∈F k=1

where āi = [ā1i , . . . , āM i ]T and h̄k = [h̄1k , . . . , h̄M k ]T . Because the model h̄k
do not depend on frequency, J¯ can be minimized eﬃciently by a clustering
algorithm that iterates the following two updates until convergence:Πf ← argminΠ

N||āΠ(k) (f ) − h̄k ||2 , for each f ∈ F ,

(19)

k=1

h̄k ←



āΠf (k) (f ), h̄k ← h̄k /||h̄k || , for each k = 1, . . . , N .

(20)

f ∈F

The ﬁrst update (19) optimizes the permutation Πf for each frequency f withthe current model h̄k . The second update (20) calculates the most probablemodel h̄k with the current permutations. This set of updates is very similar tothat of the k-means algorithm [12].
After the algorithm has converged, the permutation ambiguity in each frequency
bin is aligned byak (f ) ← aΠf (k) (f ), yk (f, t) ← yΠf (k) (f, t), k = 1, . . . , N.

(21)

In addition to aligning the permutations, the method estimates the model
parameters by
2 dmaxarg(h̄jk ), λjk = |h̄jk | .
τjk = − ·(22)πcFrom these parameters and sensor array geometry, we can perform source
localization, such as direction-of-arrival (DOA) estimation.

696

4

H. Sawada et al.

Solving the Scaling Problem

The ultimate goal as regards the scaling problem is to recover each source sk (t),i.e. multichannel blind deconvolution. However, this is very diﬃcult with coloredsource signals, such as speech. A feasible goal [13, 14] is simply to recover theobservation of each source k at a reference sensor J(23)
l hJk (l) sk (t − l) .
If we consider the frequency-domain counterpart of the above discussion, there
is no practical way to recover the amplitude and phase of sk (f, t) blindly, butthere is a feasible way to recover those ofhJk (f )sk (f, t)(24)

instead [3, 13]. We use this criterion for the scaling problem.
Calculating the inverse (5) and obtaining the linear combination form (6) ofx(f, t) provides an instant solution to the scaling problem. If ICA works well andthe permutation ambiguity is solved, we obtain separated signals yk (f, t) thatshould be close to source signals sk (f, t), now only up to the scaling ambiguity. If
we compare (3) and (6), we see that ak (f )yk (f, t) should be close to hk (f )sk (f, t)
and therefore aJk (f )yk (f, t) should be close to hJk (f )sk (f, t). Thus the scalingalignment can be performed simply byyk (f, t) ← aJk (f )yk (f, t) .

(25)

In other words, there is no scaling ambiguity to be considered in (6) if we do notdiscriminate between ai (f ) and yi (f, t).

5

A Link to the Time-Frequency Masking Method

This section reveals a link between the ICA-based method and the timefrequency (T-F) masking method. The link becomes clear by the linear combination form (6) of x(f, t) obtained by the inverse (5) of the ICA separation
matrix W(f ).
Let us explain the T-F masking method, in which we assume the sparsenessof source signals, i.e. at most only one source is active for each time-frequency
slot (f, t). Based on this assumption, the mixture model (3) can be simpliﬁed asx(f, t) = hk (f )sk (f, t), k ∈ {1, . . . , N }.

(26)

where k depends on each time-frequency slot (f, t). Then, the method classiﬁes the observation vectors x(f, t), ∀ f, t into N clusters C1 , . . . , CN so that thek-th cluster contains observation vectors in which the k-th source is the only active source. After the classiﬁcation, time domain separated signals yk (t) are obtained by applying an inverse STFT (ISTFT) to the following classiﬁed frequency
components
xJ (f, t) x(f, t) ∈ Ck ,yk (f, t) =(27)
0otherwise.

On Calculating the Inverse of Separation Matrix in Frequency-Domain BSS

697

In the classiﬁcation, the spatial information expressed in x(f, t) is extracted
and used. Typically, the phase diﬀerence normalized with frequency and/or theamplitude diﬀerence between two sensors: x2 (f, t) 
arg[x2 (f, t)/x1 (f, t)](28)
and/or 
2πf
x1 (f, t) are calculated for the classiﬁcation [7-9]. However, these papers presented only
cases with two sensors. Recently, we proposed a new technique for using all theinformation of more than two sensors [10]. The technique used there is similar
to that presented in Sec. 3. Thus, we consider that various techniques for classifying observation vectors x(f, t) in the T-F masking method can be used toclassify basis vectors ai (f, t) for solving the permutation problem in the ICAbased method, and vice versa.
Let us discuss this relationship in the following. If the sparseness assumption
is satisﬁed, the linear combination form (6) obtained by ICA is reduced tox(f, t) = ai (f )yi (f, t), i ∈ {1, . . . , N }.

(29)

where i depends on each time-frequency slot (f, t). If we compare (26) and (29),
we see that hk (f )sk (f, t) should be close to ai (f )yi (f, t) for each time-frequency
slot (f, t). Thus, the spatial information expressed in observation vectors x(f, t)
with the sparseness assumption (26) is the same as that of basis vectors ai (f, t) up
to the scaling ambiguity. Therefore, we can use the same techniques for extracting
spatial information from observation vectors x(f, t) and basis vectors ai (f, t).
The normalization formulas (10) and (15) and the clustering procedure of(19) and (20) can be used not only for the ICA-based method but also the T-F
masking method. Of course, we need to replace ai (f ) with x(f, t) and modify(19) for a standard clustering algorithm such as the k-means algorithm [12].
Figure 3 shows the ﬂows of both methods in accordance with this idea.

STFT

ICA

Basis vectors

Permutation

ISTFT

Normalize& clusterObservation vectors
STFT

T-F masking

ISTFT

Normalize& clusterFig. 3. Flows of ICA-based method (above) and time-frequency (T-F) masking method(below)

698

6

H. Sawada et al.

Experimental Results

We have performed experiments with the conditions shown in Fig. 4. We used
16 combinations of three speeches to evaluate the separation performance. The
system did not have to know the sensor geometry for solving the permutation
problem, but just the maximum distance dmax = 4 cm between the reference
sensor and any other sensor. The computational time was less than 3 seconds
for 3-second speech mixtures, meaning that real-time processing was possible.

Reverberation time:
130 ms

a

0.8m
1m
1.35m
heightc

Room size:
4.45 m × 3.55 m × 2.50 m

b 1.35mheight0.8m

SIR / SDR (dB)1mheightSampling rate:
8000 Hz
STFT frame:
size = 128 ms
shift = 32 ms

15

Sources:
3-second speeches

SIR
SDR

10

5

4cm
4cm

3.2cm

Fig. 4. Experimental conditions

0

ICA

T−F masking
Methods

Fig. 5. Separation performance

Figure 5 shows the average signal-to-interference ratio (SIR) and signal-todistortion ratio (SDR), whose detailed deﬁnitions can be found in [10]. Basically,
the SIR indicates how well the mixtures are separated into the sources, andthe SDR indicates how close each separated signal is to the observation of thecorresponding source at the reference sensor. Since the number of sensors was
suﬃcient for the number of sources in this case, ICA-based method worked better
than T-F masking as shown in Fig. 5. We have also already obtained results withanother setup where the number of sensors was insuﬃcient (N = 4, M = 3) andthe T-F masking method still worked [10].

7Conclusion

In the ICA-based frequency-domain BSS, the permutation and scaling ambiguityof the ICA solution should be aligned. Once we have the form (6) by calculating
the inverses of ICA separation matrices W, the scaling ambiguity does not haveto be considered. To align the permutation ambiguity, we can exploit the mixing
system information represented in basis vectors ai , as Sec. 3 presents an eﬃcientmethod. The form (6) clariﬁes the relationship between the ICA-based method
and the T-F masking method. The same technique as that presented in Sec. 3
can be used in the T-F masking method for clustering observations x(f, t) andextracting the mixing system information.

On Calculating the Inverse of Separation Matrix in Frequency-Domain BSS

699

References
1. Smaragdis, P.: Blind separation of convolved mixtures in the frequency domain.
Neurocomputing 22 (1998) 21–34
2. Parra, L., Spence, C.: Convolutive blind separation of non-stationary sources.
IEEE Trans. Speech Audio Processing 8 (2000) 320–327
3. Murata, N., Ikeda, S., Ziehe, A.: An approach to blind source separation based on
temporal structure of speech signals. Neurocomputing 41 (2001) 1–24
4. Saruwatari, H., Kurita, S., Takeda, K., Itakura, F., Nishikawa, T., Shikano, K.:
Blind source separation combining independent component analysis and beamforming. EURASIP Journal on Applied Signal Processing 2003 (2003) 1135–1146
5. Mukai, R., Sawada, H., Araki, S., Makino, S.: Frequency domain blind source
separation for many speech signals. In: Proc. ICA 2004 (LNCS 3195). SpringerVerlag (2004) 461–469
6. Sawada, H., Araki, S., Mukai, R., Makino, S.: Blind extraction of a dominantsource signal from mixtures of many sources. In: Proc. ICASSP 2005. Volume III.
(2005) 61–64
7. Aoki, M., Okamoto, M., Aoki, S., Matsui, H., Sakurai, T., Kaneda, Y.: Sound
source segregation based on estimating incident angle of each frequency component of input signals acquired by multiple microphones. Acoustical Science andTechnology 22 (2001) 149–157
8. Rickard, S., Balan, R., Rosca, J.: Real-time time-frequency based blind source
separation. In: Proc. ICA2001. (2001) 651–656
9. Yilmaz, O., Rickard, S.: Blind separation of speech mixtures via time-frequency
masking. IEEE Trans. Signal Processing 52 (2004) 1830–1847
10. Araki, S., Sawada, H., Mukai, R., Makino, S.: A novel blind source separation
method with observation vector clustering. In: Proc. 2005 International Workshop
on Acoustic Echo and Noise Control (IWAENC 2005). (2005)
11. Hyvärinen, A., Karhunen, J., Oja, E.: Independent Component Analysis. John
Wiley & Sons (2001)
12. Duda, R.O., Hart, P.E., Stork, D.G.: Pattern Classiﬁcation. 2nd edn. Wiley Interscience (2000)
13. Matsuoka, K., Nakashima, S.: Minimal distortion principle for blind source separation. In: Proc. ICA 2001. (2001) 722–727
14. Takatani, T., Nishikawa, T., Saruwatari, H., Shikano, K.: Blind separation of binaural sound mixtures using SIMO-model-based independent component analysis.
In: Proc. ICASSP 2004. Volume IV. (2004) 113–116