Abdeldjalil Aïssa-El-Bey

Karim Abed-Meraim

Yves Grenier

Blind Separation of Underdetermined Convolutive Mixtures using their Time-Frequency Representation

Keywords: blind source separation, underdetermined/overcomplete representation, vector clustering, subspace projection, speech signals, convolutive mixture, time-frequency distribution, sparse signal decomposition/representation

This paper considers the blind separation of nonstationary sources in the underdetermined convolutive mixture case. We introduce, two methods based on the sparsity assumption of the sources in the time-frequency (TF) domain. The first one assumes that the sources are disjoint in the TF domain; i.e. there is at most one source signal present at a given point in the TF domain. In the second method, we relax this assumption by allowing the sources to be TF-nondisjoint to a certain extent. In particular, the number of sources present (active) at a TF point should be strictly less than the number of sensors. In that case, the separation can be achieved thanks to subspace projection which allows us to identify the active sources and to estimate their corresponding time-frequency distribution (TFD) values. Another contribution of this paper is a new estimation procedure for the mixing channel in the underdetermined case. Finally, numerical performance evaluations and comparisons of the proposed methods are provided highlighting their effectiveness.

I. INTRODUCTION

T HE OBJECTIVE of blind source separation (BSS) is to extract the original source signals from their mixtures and possibly to estimate the unknown mixing channel using only the information of the observed signal with no, or very limited, knowledge about the source signals and the mixing channel. The BSS problem arises in many fields of study including speech processing, data communication, biomedical signal processing, etc [START_REF] Cichocki | Adaptive Blind Signal and Image Processing[END_REF]. Most approaches to blind source separation assume the sources are statistically independent and thus are often seek solutions of separation criteria using higher-order statistical information [START_REF] Cardoso | Blind signal separation: statistical principles[END_REF] or using only second order statistical information in the case where the sources have temporal coherency [START_REF] Belouchrani | A blind source separation technique using second-order statistics[END_REF], are nonstationary [START_REF] Belouchrani | Blind source separation based on time-frequency signal representations[END_REF], or eventually are cyclostationary [START_REF] Abed-Meraim | Blind source separation using second order cyclostationary statistics[END_REF].

Although the plethora of existing BSS algorithms, the underdetermined case (UBSS for underdetermined blind source separation) where the number of sources is greater than the number of sensors remains relatively poorly treated especially in the convolutive case, and its resolution is one of the challenging problems of blind source separation. In the instantaneous mixture case, some methods exploiting the sparseness of the sources in certain transform domain have been proposed for UBSS [START_REF] Bofill | Underdetermined blind source separation using sparse representations[END_REF]- [START_REF] Abrard | A time-frequency blind signal separation method applicable to underdetermined mixtures of dependent sources[END_REF]. Other methods consider Manuscript received July 1, 2006; revised January 31, 2007. A. Aïssa-El-Bey, K. Abed-Meraim and Y. Grenier are with the TSI Department, ENST-Paris, 46 rue Barrault 75634, Paris Cedex 13, France. Email: {elbey, abed, grenier}@tsi.enst.fr similarly underdetermined mixtures of delayed sources [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF], [START_REF] Rosca | Generalized sparse signal mixing model and application to noisy blind source separation[END_REF]. All these methods proceed 'roughly' as follows: The mixtures are first transformed to an appropriate representation domain; the transformed sources are then estimated using their sparseness, and finally one recovers their time waveforms by source synthesis (for more information, see the recent survey work [START_REF] O'grady | Survey of sparse and nonsparse methods in source separation[END_REF]).

The UBSS methods for nonstationary sources have been proposed, given that these sources are sparse in the timefrequency (TF) domain [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF], [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF]. The first method uses timefrequency distributions (TFDs), whereas the second one uses a linear TFD. The main assumption used in these methods is that the sources are TF-disjoint. In other words, there is at most one source present at any point in the TF domain. This assumption is rather restrictive, though the methods have also showed that they worked well under a quasi sparseness condition, i.e. sources are TF-almost-disjoint.

In this paper, we focus on the UBSS in convolutive mixtures case and target the relaxation of the TF-disjoint condition by allowing the sources to be nondisjoint in the TF domain; that is, multiple sources are possibly present at any point in the TF domain. This case has been considered in [START_REF] Linh-Trung | Underdetermined blind source separation of non-disjoint nonstationary sources in time-frequency domain[END_REF] for the separation of instantaneous mixtures, in [START_REF] Rosca | Generalized sparse signal mixing model and application to noisy blind source separation[END_REF] for the deconvolution of single-path channels with non-zero delays, in [START_REF] Araki | A novel blind source separation method with observation vector clustering[END_REF] where a priori information about the location of the considered sources as well as an approximation of the filter impulse response are considered, and in [START_REF] Araki | Blind separation of more speech than sensors with less distortion by combining sparseness and ICA[END_REF] where binary TF-masking (or directivity pattern based masking [START_REF] Araki | Underdetermined blind separation of convolutive mixtures of speech with directivity pattern based mask and ICA[END_REF], [START_REF] Araki | Underdetermined blind separation of convolutive mixtures of speech by combining time-frequency masks and ICA[END_REF]) and ICA technique are jointly used. In particular, we limit ourselves to the scenario where the number of sources present at any point is smaller than the number of sensors. Under this assumption, the separation of TF-nondisjoint sources is achieved thanks to subspace projection. Subspace projection allows us to identify at any point the active sources, and then to estimate their corresponding TFD values.

The main contribution of this paper consists in two new algorithms for UBSS in the TF domain; the first one uses vector clustering while the other uses subspace projection. Another side contribution of the paper is an estimation method for the mixing channel matrix.

The paper is organized as follows. Section II-A formulates the UBSS problem, introduces the underlying TF tools, and states some TF conditions necessary for the separation of nonstationary sources in the TF domain. In Section III-A, we propose a new method for the blind estimation of mixing channel. Section III-B deals with the TF-disjoint sources. It proposes a cluster-based TF-CUBSS (Time-frequency convolutive underdetermined blind source separation) algorithm.

Section III-C proposes the subspace-based TF-CUBSS algorithm for TF-nondisjoint sources. Some comments and remarks on the proposed methods are provided in Section IV. Finally, the performance of the above methods are numerically evaluated in Section V while Section VI is devoted for the concluding remarks.

II. PROBLEM FORMULATION

A. Data model

Let s 1 (t), . . . , s N (t) be the desired sources to be recovered from the convolutive mixtures x 1 (t), . . . , x M (t) given by:

x(t) = K k=0 H(k)s(t -k) + η(t) (1)
where

s(t) = [s 1 (t), . . . , s N (t)]
T is the source vector with the superscript T denoting the transpose operation,

x(t) = [x 1 (t), . . . , x M (t)]
T is the mixture vector, η(t) is the observation noise, and H(k

) def = [h 1 (k), . . . , h N (k)] are M × N matrices for k ∈ [0, K]
representing the impulse response coefficients of the channel that satisfies:

Assumption 1: The channel is such that each column vector of

H(z) def = K k=0 H(k)z -k def = [h 1 (z), . . . , h N (z)]
is irreducible, i.e. the entries of h i (z) denoted h ij (z), j = 1, . . . , M , have no-common zeros ∀i. Moreover, any M column vectors of H(z) form a polynomial matrix H(z) that it full rank over the unit-circle, i.e. rank( H(f )) = M ∀f . The sources are nonstationary, that is their frequency spectra vary in time. Often, nonstationarity gives rise to more difficulties in a problem, however, in this case it actually offers certain diversity that allows us to achieve the BSS without using higher-order approaches by directly exploiting the additional information of this TF diversity across the spectra [START_REF] Belouchrani | Blind source separation based on time-frequency signal representations[END_REF]. In that case, we often use the powerful tool of time-frequency signal analysis which basic concept is introduced next.

B. Time-frequency distributions

TF signal processing provides effective tools for analyzing nonstationary signals, whose frequency content varies in time. This concept is a natural extension of both the time domain and the frequency domain processing that involves representing signals in a two-dimensional space, the joint TF domain, hence providing a distribution of signal energy versus time and frequency simultaneously. For this reason, a TF representation is commonly referred to as a time-frequency distribution (TFD).

Well-known TFD 1 and most used in practice is the shorttime Fourier transform (STFT):

S x (t, f ) ∞ -∞ x(τ )w(τ -t) e -j2πf τ dτ, ( 2
)
where w(t) is a windowing function and x(t) a given nonstationary signal. Note that the STFT is a linear TFD and thus has the advantage of simplicity compared to other nonlinear (quadratic) TFDs, e.g. Wigner-Ville and Cohen's class distributions [START_REF] Boashash | Time Frequency Signal Analysis and Processing: Method and Applications[END_REF].

C. TF conditions on the sources

In order to deal with UBSS, one often seeks for a sparse representation of the sources [START_REF] Bofill | Underdetermined blind source separation using sparse representations[END_REF]. In other words, if the sources can be sparsely represented in some domain, then their separation can be carried out in that domain by exploiting their sparseness.

1) TF-disjoint sources: Recently, there have been several UBSS methods, notably those in [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF] and [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF], in which the TF domain has been chosen to be the underlaying sparse domain. These two papers have based their solutions on the assumption that the sources are disjoint in the TF domain. Mathematically, if Ω 1 and Ω 2 are the TF supports of two sources s 1 (t) and s 2 (t) then the sources are said TF-disjoint if Ω 1 ∩ Ω 2 = ∅. However, this is a rather strict assumption. A more practical assumption is that the sources are almostdisjoint in the TF domain [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF], allowing some small overlapping in the TF domain, for which the above two methods (in [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF] and [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF]) also worked.

2) TF-nondisjoint sources: In this paper, we want to relax the TF-disjoint condition by allowing the sources to be nondisjoint in the TF domain.

This is motivated by a drawback of the methods in [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF], [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF]. Although these methods worked under the TF-almostdisjoint condition, they did not explicitly treat the TF regions (points) where the sources were overlapping. A point at the overlapping of two sources was assigned 'by chance' to belong to only one of the sources. As a result, the source that picks up this point will have some information of the other source while the latter loses some information of its own. The loss of information can be recovered to some extent by the interpolation at the intersection point using TF synthesis. However, for the other source, there is an interference at this point, hence the separation performance may degrade if no treatment is provided. If the number of overlapping points increases (i.e. the TF-almost-disjoint condition is violated), the performance of the separation is expected to degrade unless the overlapping points are properly treated. This paper will give such a treatment using subspace projection. Therefore, we will allow the sources to be nondisjoint in the TF domain; that is, multiple sources are allowed to be present at any point in the TF domain. However, instead of being inevitably nondisjoint, we limit ourselves by making the following constraint:

Assumption 2: The number of active sources (i.e. sources that overlap) at any TF point is strictly less than the number of sensors. In other words, for the configuration of M sensors, there exists at most (M -1) sources at any point in the TF domain. For the special case when M = 2, Assumption 2 reduces to the disjoint condition. Note that in [START_REF] Araki | Blind separation of more speech than sensors with less distortion by combining sparseness and ICA[END_REF]- [START_REF] Araki | Underdetermined blind separation of convolutive mixtures of speech by combining time-frequency masks and ICA[END_REF], the case of M overlapping sources has been treated thanks to additional strong assumptions that we do not consider in our work. More specifically, the channels are assumed to be of single-path with a given direction of arrival, and the sources are such that one of them is present alone at certain time instant and can be removed by binary masking.

III. TF-CUBSS ALGORITHM

In order to solve the UBSS problem in the convolutive case, we propose to identify first the impulse response of the channels (see the algorithm's diagram in Figure 1). This problem in overdetermined case is very difficult and becomes almost impossible in the underdetermined case without side information on the considered sources. In this work and similarly to [START_REF] Huang | A blind channel identification-based two-stage approach to separation and dereverberation of speech signals in a reverberant environment[END_REF], we exploit the sparseness property of the audio sources by assuming that from time to time only one source is present. In other words, we consider the following assumption:

x 2 (t) x 1 (t)

Diagram of proposed TF-CUBSS algorithm combining channel identification and UBSS technique in TF domain.

Assumption 3: There exists, periodically, time intervals where only one source is present in the mixture. This occurs for all source signals of the considered mixtures (see Figure 2). To detect these time intervals, we propose to use informationcriteria based testing for the estimation of the number of sources present in the signal (see Section III-A for more details).

A. Channel estimation

Based on assumption 3, we propose here to apply SIMO (Single Input Multiple Output) based techniques to blindly estimate the channel impulse response. Regarding the problem at hand, we have to solve 3 different problems: first, we have to select time intervals where only one source signal is effectively present; then, for each selected time interval one should apply an appropriate blind SIMO identification technique to estimate the channel parameters; finally, the way we proceed, the same channel may be estimated several times and hence one has to group together (cluster) the channel estimates into N classes corresponding to the N source channels. 1) Source number estimation: Let define the spatiotemporal vector:

x d (t) = [x T (t), . . . , x T (t -d + 1)] T = N k=1 H k s k (t) + η d (t), (3)
where H k are block-Sylvester matrices of size dM × (d + K):

H k =    h k (0) • • • h k (K) 0 . . . . . . 0 h k (0) • • • h k (K)    s k (t) def = [s k (t), . . . , s k (t -K -d + 1)
] T and d is a chosen processing window size. Under the no-common zeros assumption (Assumption 1) and for large window sizes (see [START_REF] Wax | Detection of signals by information theoretic criteria[END_REF] for more details), matrices H k are full column rank.

Hence, in the noiseless case, the rank of the data covariance matrix

R def = E[x d (t)x H d (t)
] is equal to min(p(d + K), dM ) where p is the number of sources present in the considered time interval over which the covariance matrix is estimated. In particular, for p = 1, one has the minimum rank value equal to (d + K). Therefore, our approach consists in estimating the rank of the sample averaged covariance matrix R over several time slots (intervals) and select those corresponding to the smallest rank value r = d + K. In the case where p sources are active (present) in the considered time slot, the rank would be r = p(d + K) and hence p can be estimated by the closest integer value to r d+K . The estimation of the rank value is done here by Akaike's criterion [START_REF] Wax | Detection of signals by information theoretic criteria[END_REF] according to:

r = arg min k 2 6 6 6 6 4 -2 log 0 B B B @ M d Q i=k+1 λ 1/(M d-k) i 1 M d-k M d P i=k+1 λ i 1 C C C A (M d-k)T s + 2k(2M d -k) , (4)
where λ 1 ≥ . . . ≥ λ M d represent the eigenvalues of R and T s is the time slot size. This criterion represents the maximum likelihood estimate of the system parameters (given here by the signal eigenvectors and eigenvalues of the covariance matrix) penalized by the number of free adjusted parameters under the asymptotic Gaussian distribution of the latter (see [START_REF] Wax | Detection of signals by information theoretic criteria[END_REF] for more details). Note that it is not necessary at this stage, to know exactly the channel degree K as long as d > K (i.e. an over-estimation of the channel degree is sufficient) in which case the presence of one signal source is characterized by: d < r < 2d . Histogram representing the number of time intervals for each estimated number of sources for 4 audio sources and 3 sensors in convolutive mixture case.

2) Blind channel identification:

To perform the blind channel identification, we have used in this paper the Cross-Relation (CR) technique described in [START_REF] Xu | A least-squares approach to blind channel identification[END_REF], [START_REF] Aïssa-El-Bey | Blind system identification using cross-relation methods: Further results and developments[END_REF]. Consider a time interval where we have only the source s i present. In this case, we can consider a SIMO system of M outputs given by:

x(t) = K k=0 h i (k)s i (t -k) + η(t), ( 5
)
where

h i (k) = [h i1 (k) . . . h iM (k)] T , k = 0, • • • , K.
From ( 5), the noise-free outputs x j (k), 1 ≤ j ≤ M are given by:

x j (k) = h ij (k) * s i (k), 1 ≤ j ≤ M, (6)
where " * " denotes the convolution. Using commutativity of convolution, it follows:

h il (k) * x j (k) = h ij (k) * x l (k), 1 ≤ j = l ≤ M. ( 7
)
This is a linear equation satisfied by every pair of channels. It was shown that reciprocally, the previous M (M -1)/2 cross-relations characterize uniquely the channel parameters.

We have the following theorem [START_REF] Xu | A least-squares approach to blind channel identification[END_REF]: Theorem 1: Under the no-common zeros assumption (Assumption 1), the set of cross-relations (in the noise free case): [START_REF] Linh-Trung | Underdetermined blind source separation of non-disjoint nonstationary sources in time-frequency domain[END_REF] where h (z) = [h 1 (z) . . . h M (z)] T is a M × 1 polynomial vector of degree K, is satisfied if and only if h (z) = αh i (z) for a given scalar constant α. By collecting all possible pairs of M channels, one can easily establish a set of linear equations. In matrix form, this set of equations can be expressed as:

x l (k) * h j (k) -x j (k) * h l (k) = 0, 1 ≤ l < j ≤ M,
X M h i = 0, (9)
where

h i def = [h i1 (0) . . . h i1 (K), . . . , h iM (0) . . . h iM (K)]
T and X M is defined by:

X 2 = [X (2) , -X (1) ], X n =      X n-1 0 X (n) 0 -X (1) . . . . . . 0 X (n) -X (n-1)      , ( 10
)
with n = 3, . . . , M and:

X (n) =    x n (K) . . . x n (0) . . . . . . x n (T s -1) . . . x n (T s -K -1)    . ( 11
)
In the presence of noise, equation ( 9) can be naturally solved in the least-squares (LS) sense according to:

h CR = arg min h i =1 h H i X H M X M h i ( 12
)
which solution is given by the least unit-norm eigenvector of matrix X H M X M . It is shown in [START_REF] Xu | A least-squares approach to blind channel identification[END_REF] that the noise term in the quadratic form [START_REF] Rosca | Generalized sparse signal mixing model and application to noisy blind source separation[END_REF] has a mean value proportional to the identity matrix. Consequently, the channel estimates remains unbiased under white additional noise assumption.

Remark: We have presented here a basic version of the CR method. In [START_REF] Ahmad | Proportionate frequency domain adaptive algorithms for blind channel identification[END_REF] an improved version of the method (introduced in the adaptive scheme) is proposed exploiting the quasi-sparse nature of acoustic impulse responses. Other channel estimation techniques in the overcomplete case, e.g. [START_REF] Winter | Overcomplete BSS for convolutive mixtures based on hierarchical clustering[END_REF], can be used as well at this stage.

3) Clustering of channel vector estimates: The first step of our channel estimation method consists in detecting the time slots where only one single source signal is 'effectively' present. However, the same source signal s i may be present in several time intervals (see Figure 2 and Figure 3) leading to several estimates of the same channel vector h i . We end up, finally, with several estimates of each source channel that we need to group together into N classes. This is done by clustering the estimated vectors using k-means algorithm [START_REF] Frank | The data analysis handbook[END_REF]. The i th channel estimate is evaluated as the centroid of the i th class. [START_REF] Wax | Detection of signals by information theoretic criteria[END_REF] to detect the number of source, then application of the blind identification algorithm in [START_REF] Xu | A least-squares approach to blind channel identification[END_REF], [START_REF] Aïssa-El-Bey | Blind system identification using cross-relation methods: Further results and developments[END_REF] followed by vector clustering. 2) Mixture STFT computation by [START_REF] Araki | Blind separation of more speech than sensors with less distortion by combining sparseness and ICA[END_REF] and noise thresholding by ( 16) 3) Vector clustering by [START_REF] Araki | Underdetermined blind separation of convolutive mixtures of speech by combining time-frequency masks and ICA[END_REF] and [START_REF] Boashash | Time Frequency Signal Analysis and Processing: Method and Applications[END_REF]. 4) Source STFT estimation by [START_REF] Huang | A blind channel identification-based two-stage approach to separation and dereverberation of speech signals in a reverberant environment[END_REF]. 5) Source TF synthesis by [START_REF] Griffin | Signal estimation from modified shorttime fourier transform[END_REF].

B. UBSS algorithm with TF-disjoint assumption

As we have seen before, the STFT is often used for speech/audio signals because of its low computational cost. Therefore, in this section we propose a new cluster-based TF-CUBSS algorithm using the STFT for convolutive mixture case. Note that the STFT is a particular form of wavelet transforms which have been used in [START_REF] Zibulevsky | Independent Component Analysis: Principles and Practice[END_REF] for the UBSS of image signals.

After transformation into the TF domain using the STFT, the model in ( 1) becomes (in the noiseless case):

S x (t, f ) = H(f )S s (t, f ), (13)
where S x (t, f ) is the mixture STFT vector, S s (t, f ) is the source STFT vector and

H(f ) = [h 1 (f ) . . . h N (f )]
is the channel Fourier Transform matrix. Under the assumption that all sources are disjoint in the TF domain, (13) reduces to

S x (t, f ) = h i (f )S si (t, f ), ∀(t, f ) ∈ Ω i , ∀i ∈ N , ( 14
)
where N = {1, . . . , N } and Ω i is the TF support of the i th source.

Consequently, two TF points (t 1 , f 1 ) and (t 2 , f 2 ) belonging to the same region Ω i (i.e. corresponding to the source signal s i ) are 'associated' with the same channel h i . It is this observation that is used to derive the separation algorithm summarized in Table I and detailed next. First, we compute the STFT of the mixtures, S x (t, f ), by applying (2) for each of the mixture in x(t), as follows:

S x i (t, f ) = m=(L-1)/2-1 m=-(L-1)/2 w(t -m)x i (m)e -j2πf m , i = 1, . . . , M, (15a) S x (t, f ) = [S x1 (t, f ), . . . , S x M (t, f )] T . ( 15b
)
where w(t) is a chosen window (in our simulations we chose Hamming window) of length L.

Then, we apply a noise thresholding procedure which mitigates the noise effect and reduces the computational cost as only the selected TF points are further treated by our algorithm.

In particular, for each frequency f 0 , we apply the following criterion for all the time points t k belonging to the frequencyslice (t, f 0 )

If S x (t k , f 0 ) max t { S x (t, f 0 ) } > 1 , then keep (t k , f 0 ), ( 16
)
where 1 is a small threshold (typically, 1 = 0.01). Then, the set of all selected points, Ω, is expressed by

Ω = N i=1 Ω i ,
where Ω i is the TF support of the source s i (t). Note that, the effects of spreading the noise energy while localizing the source energy in the time-frequency domain amounts to increasing the robustness of the proposed method with respect to noise (see Part IV of [START_REF] Boashash | Time Frequency Signal Analysis and Processing: Method and Applications[END_REF]). Hence, by equation ( 16), we would keep only time-frequency points where the signal energy is non-negligible, the other time-frequency points are rejected, i.e. not further processed, since considered to represent noise contribution only. Also, due to the noise energy spreading, the contribution of the noise in the source time-frequency points is relatively, negligible at least for moderate and high SNRs.

On the other hand, note that the noise thresholding as well as TF masking induce non-linear distortion in the reconstructed signal. Now, how this distortion affects the source estimates is an open problem that still raises many questioning and research works including those which try to mitigate this distortion in the TF domain, e.g. [START_REF] Rosca | Statistical inference of missing speech data in the ICA domain[END_REF].

After noise thresholding, the clustering procedure can be done as follows: For each TF point, we obtain the spatial direction vectors by:

v(t, f ) = S x (t, f ) S x (t, f ) , (t, f ) ∈ Ω, (17)
and force them, without loss of generality, to have the first entry real and positive. Next, we cluster these vectors into N classes {C i | i ∈ N } by minimizing the criterion:

v(t, f ) ∈ C i ⇐⇒ i = arg min k v(t, f ) - h k (f )e -jθ k h k (f ) 2 (18
) where h k (f ) is the Fourier Transform of the k th channel vector estimate (given by [START_REF] Rosca | Generalized sparse signal mixing model and application to noisy blind source separation[END_REF] and the proposed clustering procedure) and θ k is the phase argument of h k1 (f ) (this is to force the first entry to be real positive). The collection of all points, whose vectors belong to the class C i , now forms the TF support Ω i of the source s i (t).

Therefore, we can estimate the STFT of each source s i (t) by:

S si (t, f ) = b h H i (f ) b h i (f ) 2 S x (t, f ), ∀ (t, f ) ∈ Ω i , 0, otherwise, (19)
since, from ( 14), we have

b h H i (f ) b h i (f ) 2 Sx(t, f ) = b h H i (f )h i (f ) b h i (f ) 2 Ss i (t, f ) ≈ Ss i (t, f ), ∀ (t, f ) ∈ Ω i .

C. UBSS algorithm with TF-nondisjoint assumption

We have seen the cluster-based TF-CUBSS methods, using the STFT, as summarized in Table I. This method relies on the assumption that the sources were TF-disjoint, which led to the TF-transformed structure in [START_REF] Araki | A novel blind source separation method with observation vector clustering[END_REF]. The latter is no longer valid, when the sources are nondisjoint in the TF domain.

Under the TF-nondisjoint condition, stated in Assumption 2, we propose in this section an alternative method using subspace projection.

Recall that the first two steps of the cluster-based quadratic TF-CUBSS algorithm do not rely on the assumption of TFdisjoint sources (see Table I). Therefore, we can reuse these steps to obtain the channel estimation and all the TF points of Ω. Under the TF-nondisjoint condition, consider a TF point (t, f ) ∈ Ω at which there are J sources s α1 (t), . . . , s α J (t) present, with J < M where α 1 , . . . , α J ∈ N denote the indices of the sources present at (t, f ). Our goal is to identify the sources that are present at (t, f ), i.e. α 1 , . . . , α J , and to estimate the STFT of each of these contributing sources.

We define the following:

s = [s α 1 (t), . . . , s α J (t)] T , (20a) Hα (f ) = [h α 1 (f ), . . . , h α J (f )]. (20b)
Then, ( 13) is reduced to the following

S x (t, f ) = Hα (f )S s(t, f ). ( 21
)
Let

Hβ (f ) = [h β 1 (f ), . . . , h β J (f )
] and Q β (f ) be the orthogonal projection matrix onto the noise subspace of Hβ (f ) expressed by:

Q β (f ) = I -Hβ (f ) HH β (f ) Hβ (f ) -1 HH β (f ). (22)
We have the following observation:

Q β (f )h i (f ) = 0, i ∈ {β 1 , . . . , β J } Q β (f )h i (f ) = 0, i ∈ N \{β 1 , . . . , β J } . ( 23
)
Consequently, as S x (t, f ) ∈ Range{ Hα (f )}, we have

Q β (f )S x (t, f ) = 0, if {β 1 , . . . , β J } = {α 1 , . . . , α J } Q β (f )S x (t, f ) = 0, otherwise .
(24) If H(f ) has already been estimated by the method presented in Section III-A, then this observation gives us the criterion to detect the indices α 1 , . . . , α J ; and hence, the contributing sources at the considered TF point. In practice, to take into account noise, one detects the column vectors of Hα (f ) by minimizing:

{α 1 , . . . , α J } = arg min β 1 ,...,β J { Q β (f )S x (t, f ) } . ( 25
)
Next, TFD values of the J sources at TF point (t, f ) are estimated by:

S s(t, f ) ≈ H# α (f )S x (t, f ), (26)
where the superscript ( # ) represents the Moore-Penrose's pseudo-inversion operator.

In the simulation, the optimization problem of ( 25) is solved using exhaustive search. This is computationally tractable for small array sizes but would be prohibitive if M is very large.

Table II provides a summary of the subspace projection based TF-CUBSS algorithm using STFT. [START_REF] Wax | Detection of signals by information theoretic criteria[END_REF] to detect the number of source, then application of the blind identification algorithm in [START_REF] Xu | A least-squares approach to blind channel identification[END_REF], [START_REF] Aïssa-El-Bey | Blind system identification using cross-relation methods: Further results and developments[END_REF] followed by vector clustering. 2) STFT computation and noise thresholding.

3) For all selected TF points, detect the active sources by [START_REF] Aïssa-El-Bey | Blind system identification using cross-relation methods: Further results and developments[END_REF] and ( 25). 4) Source STFT estimation by [START_REF] Zibulevsky | Independent Component Analysis: Principles and Practice[END_REF]. 5) Source TF synthesis by [START_REF] Griffin | Signal estimation from modified shorttime fourier transform[END_REF].

IV. DISCUSSION

We discuss here certain points relative to the proposed TF-CUBSS algorithms and their applications.

1) Number of sources: The number of sources N is assumed known in the clustering method that we have used. However, there exist clustering methods [START_REF] Frank | The data analysis handbook[END_REF] which perform the class estimation as well as the estimation of the number N . In our simulation, we have observed that most of the time the number of classes is overestimated, leading to poor source separation quality. Hence, robust estimation of the number of sources in the UBSS case remains a difficult open problem that deserves particular attention in future works.

2) Number of overlapping sources: In the subspace-based approach, it is also possible to consider a fixed (maximum) value of J that is used for all TF points. Indeed, if the number of overlapping sources is less than J , we would estimate close-to-zero source STFT values. For example, if we assume J = 2 sources are present at a given TF point while only one source is effectively contributing, then we estimate one closeto-zero source STFT value. This approach increases slightly the estimation error of the source signals (especially at low SNRs) but has the advantage of simplicity compared to using information theoretic-based criteria for estimating the value of J .

3) Separation quality versus number of sources: Although we are in the underdetermined case, the number of sources N should not exceed too much the number of sensors. Indeed, when N increases, the level of source interference increases, and hence, the source quasi-disjointness assumption is illsatisfied. Moreover, for a large number of sources, the likelihood of having two closely spaced sources, i.e. such that the spatial directions h i and h j are 'close' to linear dependency, increases. In that case, vector clustering performance degrades significantly. In brief, sparseness and spatial separation are the two limiting factors against increasing the number of sources.

4) Overdetermined case: Our algorithm can be further simplified in the overdetermined case where M ≥ N . In that context, the algorithm can be reduced to the channel estimation step, the STFT computation and noise thresholding then source STFT estimation using the channel matrix pseudo-inversion at each frequency:

S s (t, f ) = H # (f )S x (t, f ).

V. SIMULATION RESULTS

In the simulations, we have considered an array of M = 3 sensors, that receives signals from N = 4 independent audio sources (3 speech signals corresponding to 2 men and 1 woman plus a guitar signal). The filter coefficients are chosen randomly and the channel order is K = 6. The sample size is T = 8192 samples (corresponding approximately to 1 second recording of speech signals sampled at 8 KHz). The separation quality is measured by the normalized mean squares estimation errors (N M SE) of the sources evaluated over N r = 200 Monte-Carlo runs and defined as:

N M SE i def = 1 N r Nr r=1 min α α s i,r -s i 2 s i 2 (27) N M SE i = 1 N r N r r=1 1 - s i,r s H i s i,r s i 2 ( 28
)
N M SE = 1 N N i=1 N M SE i . ( 29
)
where s i def = [s i (0), . . . , s i (T -1)], s i,r (defined similarly) represents the r th estimate of source s i and α is a scalar factor that compensate for the scale indeterminacy of the BSS problem.

In Figure 4, the top four plots represent the TF representation of the original source signals, the middle three plots represent the TF representation of the M mixture signals and the bottom four plots represent the TF representation of the source estimates by the subspace-based algorithm (Table II) using STFT of length 1024. Figure 5 represents the same disposition of signals but in the time domain.

In Figure 6, we compare the separation performance obtained by the subspace-based algorithm with J = 2 and the cluster-based algorithm (Table I). It is observed that subspacebased algorithm provides much better separation results than those obtained by the cluster-based algorithm. This is mainly due to the high occurrence of overlapping sources in the TF domain for this type of signals so that the 'TF-disjointness' assumption used by the TF-CUBSS algorithm is poorly satisfied. This can be observed also from Figure 4, where we can see that the TFD supports of the 4 audio sources are clearly overlapping.

In Figure 7, we present the performance of channel identification obtained by using SIMO identification algorithm (in this case we choose only the time intervals where only one source is present using AIC criterion) with SIMO and MIMO identification algorithms2 (in this case we choose the time intervals where we are in the overdetermined case; i.e. where p = 1 or p = 2). It is observed that SIMO based identification provides better results than those obtained by SIMO and MIMO identification algorithms. Indeed, the advantage of considering overdetermined MIMO system identification resides in the fact that the occurrence of MIMO (i.e. number of time intervals where this situation occurs as shown in Figure 3) is much higher than that of SIMO case. However, as we observe it, this does not compensate for the higher estimation error of MIMO systems compared to SIMO systems. The plot in Figure 8 (respectively in Figure 9) presents the separation performance when using the exact matrix H compared to that obtained with the proposed estimate H using the cluster-based method (respectively the subspacebased method). The observed performance loss is due to the channel estimation error which is relatively high for low SNRs and becomes negligible for high SNRs.

a) S s 1 (t, f ) (b) S s 2 (t, f ) (c) S s 3 (t, f ) (d) S s 4 (t, f ) (e) Sx 1 (t, f ) (f) Sx 2 (t, f ) (g) Sx 3 (t, f ) (h) S ŝ1 (t, f ) (i) S ŝ2 (t, f ) (j) S ŝ3 (t, f ) (k) S ŝ4 (t, f )
In Figure 10, we compare the performance obtained with the subspace-based method for J = 2 and J = 3. In that experiment, we have used M = 4 sensors and N = 5 source signals. One can observe that, for high SNRs, the case of J = 3 leads to a better separation performance than for the case of J = 2. However, for low SNRs, a large value of J increases the estimation noise (as mentioned in Section IV) Figure 12 illustrates the algorithm's performance when we consider long impulse response channels. More specifically, the plots represent the separation performance for channels of length 50, 100 and 200 respectively. The channel taps are generated randomly using Gaussian law. We observe a slight performance degradation when the channel order increases but the separation quality remains quite good.

In Figure 13, we compare the separation performance of our algorithm, Deville's algorithm in [START_REF] Albouy | Alternative structures and power spectrum criteria for blind segmentation and separation of convolutive speech mixtures[END_REF] and Parra's algorithm in [START_REF] Parra | Convolutive blind separation of nonstationary sources[END_REF] in the overdetermined case of 2 sensors and 2 speech signals of one man and one woman (selected among the four previous sources). The algorithms in [START_REF] Albouy | Alternative structures and power spectrum criteria for blind segmentation and separation of convolutive speech mixtures[END_REF], [START_REF] Parra | Convolutive blind separation of nonstationary sources[END_REF] separate the sources only up to an unknown filter and hence we use in this experiment the interference to signal ratio (ISR) criterion defined in [START_REF] Parra | Convolutive blind separation of nonstationary sources[END_REF] instead of the NMSE. We observe a significant performance gain in favor of the proposed method especially at high SNR values. Moreover, our method has the following advantages : (i) it can treat the underdetermined case, (ii) it estimates the sources up to a constant not to an unknown filter like in [START_REF] Albouy | Alternative structures and power spectrum criteria for blind segmentation and separation of convolutive speech mixtures[END_REF], [START_REF] Parra | Convolutive blind separation of nonstationary sources[END_REF], (iii) the proposed frame selection procedure does not involve any thresholding (the choice of an appropriate threshold value is a difficult problem as it is strongly dependent on the context) or ad-hoc selection of frequency range like in [START_REF] Albouy | Alternative structures and power spectrum criteria for blind segmentation and separation of convolutive speech mixtures[END_REF].

VI. CONCLUSION

This paper introduces new methods for the UBSS of TFdisjoint and TF-nondisjoint nonstationary sources in the convolutive mixture case using their time-frequency representations. The first proposed method has the advantage of simplicity while the second uses a weaker assumption on the source 'sparseness', i.e. the sources are not necessarily TF-disjoint, and proposes an explicit treatment of the overlapping points using subspace projection, leading to significant performance improvements. Simulation results illustrate the effectiveness of our algorithms in different scenarios.

  Fig. 1.Diagram of proposed TF-CUBSS algorithm combining channel identification and UBSS technique in TF domain.

4 Fig. 2 .

 42 Fig. 2. Time representation of 4 audio sources: this representation illustrates the audio signal sparsity (i.e. there exists time intervals where only one source is present).

Figure 3

 3 Figure3illustrates the effectiveness of the proposed method where a recording of 6 seconds of M = 3 convolutive mixtures of N = 4 sources is considered. The sampling frequency is 8 KHz and the time slot size is T s = 200 samples. The sources consist of 3 speech signals corresponding to 2 men and 1 woman plus a guitar signal. The convolutive channel is of order K = 6 and its coefficients are generated randomly using Gaussian law. One can observe that the case p = 1 (one signal source) occurs approximatively 10% of the time in the considered context.

  Fig. 3.Histogram representing the number of time intervals for each estimated number of sources for 4 audio sources and 3 sensors in convolutive mixture case.

(

Fig. 4 .

 4 Fig. 4. Simulated example (viewed in TF domain) for the subspace-based TF-CUBSS algorithm in the case of 4 speech sources and 3 sensors. The top four plots represent the original source signals, the middle three plots represent the 3 mixtures, and the bottom four plots represent the source estimates.

Fig. 5 .

 5 Fig. 5. Simulated example (viewed in time domain) for the subspace-based TF-CUBSS algorithm in the case of 4 speech sources and 3 sensors.

Fig. 6 .

 6 Fig. 6. Comparison between subspace-based and cluster-based TF-CUBSS algorithms : normalized MSE (NMSE) versus SNR for 4 speech sources and 3 sensors.

Fig. 7 .

 7 Fig.7. NMSE versus SNR for 4 audio sources and 3 sensors in convolutive mixture case : comparison of the performance of identification algorithme using only SIMO system and the algorithm using SIMO and MIMO system.

Fig. 8 .

 8 Fig. 8. Comparison, for the cluster-based TF-CUBSS algorithm, when the mixing channel H is known or unknown: NMSE of the source estimates.

Fig. 9 .

 9 Fig. 9. Comparison, for the subspace-based TF-CUBSS algorithm, when the mixing channel H is known or unknown: NMSE of the source estimates.

Fig. 10 .

 10 Fig. 10. Comparison between subspace-based and cluster-based TF-CUBSS algorithms: NMSE of the source estimates for different ranks of the projection subspace, for the case of 5 sources and 4 sensors.

Figure 11

 11 Figure 11 illustrates the rapid degradation of the separation quality when we increase the number of sources from N = 4 to N = 7. This confirms the remarks made in Section IV.Figure12illustrates the algorithm's performance when we consider long impulse response channels. More specifically, the plots represent the separation performance for channels of length 50, 100 and 200 respectively. The channel taps are generated randomly using Gaussian law. We observe a slight performance degradation when the channel order increases but the separation quality remains quite good.In Figure13, we compare the separation performance of our algorithm, Deville's algorithm in[START_REF] Albouy | Alternative structures and power spectrum criteria for blind segmentation and separation of convolutive speech mixtures[END_REF] and Parra's algorithm in[START_REF] Parra | Convolutive blind separation of nonstationary sources[END_REF] in the overdetermined case of 2 sensors and 2 speech signals of one man and one woman (selected among the four previous sources). The algorithms in[START_REF] Albouy | Alternative structures and power spectrum criteria for blind segmentation and separation of convolutive speech mixtures[END_REF],[START_REF] Parra | Convolutive blind separation of nonstationary sources[END_REF] separate the sources only up to an unknown filter and hence we use in this experiment the interference to signal ratio (ISR) criterion defined in[START_REF] Parra | Convolutive blind separation of nonstationary sources[END_REF] instead of the NMSE. We observe a significant performance gain in favor of the proposed method especially

Fig. 11 .

 11 Fig. 11. Comparison between subspace-based and cluster-based TF-CUBSS algorithms: NMSE versus number of sources.

Fig. 12 .

 12 Fig. 12. NMSE versus SNR for 4 audio sources and 3 sensors: Comparison, for the subspace-based TF-CUBSS algorithm, for different filter size K.

Fig. 13 .

 13 Fig. 13. ISR versus SNR for 2 audio sources and 2 sensors: Comparison between the subspace-based TF-CUBSS algorithm, Parra's algorithm and Deville's algorithm.

In fact, the STFT does not represent an energy distribution of the signal in the TF plane. However, for simplicity, we still refer to it as a TFD.

For the identification of MIMO system, we have used the subspace method[START_REF] Abed-Meraim | A subspace algorithm for certain blind identification problems[END_REF] for the equalization step followed by SOBI algorithm[START_REF] Belouchrani | A blind source separation technique using second-order statistics[END_REF] for the separation step.


Abdeldjalil Aïssa-El-Bey

Nguyen Linh-Trung
email: linhtrung@ieee.org

Karim Abed-Meraim

Adel Belouchrani
email: adel.belouchrani@enp.edu.dz

Yves Grenier
email: grenier@tsi.enst.fr

Underdetermined Blind Separation of Nondisjoint Sources in the Time-Frequency Domain

Keywords: blind source separation, underdetermined/overcomplete representation, spatial time-frequency representation, vector clustering, subspace projection, speech signals, sparse signal decomposition/representation

This paper considers the blind separation of nonstationary sources in the underdetermined case, when there are more sources than sensors. A general framework for this problem is to work on sources that are sparse in some signal representation domain. Recently, two methods have been proposed with respect to the time-frequency (TF) domain. The first uses quadratic timefrequency distributions (TFDs) and a clustering approach, and the second uses a linear TFD. Both of these methods assume that the sources are disjoint in the TF domain; i.e. there is at most one source present at a point in the TF domain. In this paper, we relax this assumption by allowing the sources to be TF-nondisjoint to a certain extent. In particular, the number of sources present at a point is strictly less than the number of sensors. The separation can still be achieved thanks to subspace projection that allows us to identify the sources present and to estimate their corresponding TFD values. In particular, we propose two subspace-based algorithms for TF-nondisjoint sources, one uses quadratic TFDs and the other a linear TFD. Another contribution of this paper is a new estimation procedure for the mixing matrix. Finally, then numerical performance of the proposed methods are provided highlighting their performance gain compared to existing ones.

I. INTRODUCTION

S OURCE SEPARATION aims at recovering multiple sources from multiple observations (mixtures) received by a set of linear sensors. The problem is said to be 'blind' when the observations have been linearly mixed by the transfer medium, while having no a priori knowledge of the transfer medium or the sources. Blind source separation (BSS) has applications in several areas, such as communication, speech/audio processing, and biomedical engineering [START_REF] Nandi | Blind estimation using higher-order statistics[END_REF]. A fundamental and necessary assumption of BSS is that the sources are statistically independent and thus are often sought solutions using higher-order statistical information [START_REF] Cardoso | Blind signal separation: statistical principles[END_REF]. If some information about the sources is available at hand, such as temporal coherency [START_REF] Belouchrani | A blind source separation technique using second-order statistics[END_REF], source nonstationarity [START_REF] Belouchrani | Blind source separation based on time-frequency signal representations[END_REF], or source cyclostationarity [START_REF] Abed-Meraim | Blind source separation using second order cyclostationary statistics[END_REF] then one can remain in the second-order statistical scenario.

The BSS is said to be underdetermined if there are more sources than sensors. In that case, the mixing matrix is not invertible and, consequently, a solution for source estimation must also be found even if the mixing matrix has been estimated. A general framework for underdetermined blind source separation (UBSS) is to exploit the sparseness, if it exists, of the sources in a given signal representation domain [START_REF] Bofill | Underdetermined blind source separation using sparse representations[END_REF]. The mixtures are then transformed to this domain; one may then, estimate the transformed sources using their sparseness, and finally recover their time waveforms by source synthesis. For more information on BSS and UBSS methods, see for example a recent survey [START_REF] O'grady | Survey of sparse and nonsparse methods in source separation[END_REF].

Recently, several UBSS methods for nonstationary sources have been proposed, given that these sources are sparse in the time-frequency (TF) domain [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF]- [START_REF] Barkat | Algorithms for blind components separation and extraction from the time-frequency distribution of their mixture[END_REF]. The first method uses quadratic time-frequency distributions (TFDs), whereas the second one uses a linear TFD. The main assumption used in these methods is that the sources are TF-disjoint; in other words, there is at most one source present at any point in the TF domain. This assumption is rather restrictive, though the methods have also showed that they worked well under a quasi sparseness condition, i.e. sources are TF-almost-disjoint.

In this paper, we want to relax the TF-disjoint condition by allowing the sources to be nondisjoint in the TF domain; that is, multiple sources are possibly present at any point in the TF domain. This case has been considered in [START_REF] Linh-Trung | Underdetermined blind source separation of non-disjoint nonstationary sources in time-frequency domain[END_REF] (which corresponds to part of this work) and in [START_REF] Rickard | Desprit -histogram based blind source separation of more sources than sensors using subspace methods[END_REF] for the parametric mixing matrix case. In particular, we limit ourselves to the scenario where the number of sources present at any point is smaller than the number of sensors. Under this assumption, the separation of TF-nondisjoint sources is achieved thanks to subspace projection. Subspace projection allows us to identify at any point the sources present, and hence, to estimate the corresponding TFD values of these sources.

The main contribution of this paper is proposing two subspace-based algorithms for UBSS in the TF domain; one uses quadratic TFDs while the other uses linear TFD. In line with the cluster-based quadratic algorithm proposed in [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF], we also propose here a cluster-based algorithm but using a linear TFD, which is not a block-based technique like the quadratic one. Therefore, its low cost computation is useful for processing speech and audio sources. Another contribution of the paper is a method of estimation for the mixing matrix.

The paper is organized as follows. Section II-A formulates the UBSS problem, introduces the underlying TF tools, and states some TF conditions necessary for the separation of nonstationary sources in the TF domain. Section III deals with the TF-disjoint sources. It reviews the cluster-based quadratic TF-UBSS algorithm [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF], and from that, proposes a cluster-based linear TF-UBSS algorithm. Section IV proposes two subspace-based TF-UBSS algorithms for TF-nondisjoint sources, using quadratic and linear TFDs. In this section, we propose also a method for the blind estimation of mixing matrix. There is some discussion of the proposed methods in Section V. The performance of the above methods are numerically evaluated in Section VI.

II. PROBLEM FORMULATION

A. Data model

Let s 1 (t), . . . , s N (t) be the desired sources to be recovered from the instantaneous mixtures x 1 (t), . . . , x M (t) given by:

x(t) = As(t), (1)
where

s(t) = [s 1 (t), . . . , s N (t)]
T is the source vector with the superscript T denoting the transpose operation,

x(t) = [x 1 (t), . . . , x M (t)]
T is the mixture vector, and A = [a 1 , . . . , a N ] is the mixing matrix of size M ×N that satisfies:

Assumption 1: The column vectors of A are pair-wise linearly independent. That is, for any index pair i, j ∈ N , where N = {1, . . . , N }, and i = j, we have a i and a j linearly independent. This assumption is necessary because if otherwise, we have a 1 = αa 2 for example, then the input/output relation (1) can be reduced to

x(t) = [a 1 , a 3 , . . . , a N ] [s 1 (t) + αs 2 (t), s 3 (t), . . . , s N (t)] T ,
and hence the separation of s 1 (t) and s 2 (t) is inherently impossible.

It is known that BSS is only possible up to some scaling and permutation. We take advantage of these indeterminacies to further assume, without loss of generality, that the column vectors of A all have unit norm, i.e. a i = 1 for all i ∈ N .

The sources are nonstationary, that is their frequency spectra vary in time. Often, nonstationarity imposes more difficulties on a problem, however, in this case it actually offers a solution: one can solve the BSS problem without using higher-order approaches by directly exploiting the additional information of this TF diversity across the spectra; this solution was proposed in [START_REF] Belouchrani | Blind source separation based on time-frequency signal representations[END_REF]. We defer to a little later making TF assumptions on the sources, and for now we introduce the concept of TF signal processing.

B. Time-frequency distributions

TF signal processing provides effective tools for analyzing nonstationary signals, whose frequency content varies in time. This concept is a natural extension of both the time domain and the frequency domain processing that involve representing signals in a two-dimensional space the joint TF domain, hence providing a distribution of signal energy versus time and frequency simultaneously. For this reason, a TF representation is commonly referred to as a time-frequency distribution (TFD).

The general class of quadratic TFDs of an analytic signal z(t) is defined as [START_REF] Boashash | Time Frequency Signal Analysis and Processing: Method and Applications[END_REF]:

ρ zz (t, f ) ∞ -∞ e j2πν(u-t) Γ(ν, τ ) × z(u + τ 2 )z * (u - τ 2 ) e -j2πf τ dν du dτ, (2)
where Γ(ν, τ ) is a two-dimensional function in the so-called ambiguity domain and is called the Doppler-lag kernel, and the superscript ( * ) denotes the conjugate operator. We can design a TFD with certain desired properties by properly constraining Γ.

When Γ(ν, τ ) = 1 we have the following famous Wigner-Ville distribution (WVD):

ρ wvd zz (t, f ) ∞ -∞ z(t + τ 2 )z * (t - τ 2 ) e -j2πf τ dτ. ( 3
)
The WVD is the most widely studied TFD. It achieves maximum energy concentration in the TF plane around the instantaneous frequency for linear frequency-modulated (LFM) signals. However, it is in general non-positive and it introduces the so-called "cross-terms" when multiple frequency laws (e.g. two LFM components) exist in the signals, due to the quadratic multiplication of shifted versions of the signals.

Another well-known TFD and most used in practice is the short-time Fourier transform (STFT):

S z (t, f ) ∞ -∞ z(τ )h(τ -t) e -j2πf τ dτ, (4)
where h(t) is a window function. Note that the STFT is a linear TFD1 , and its quadratic version, called the spectrogram (SPEC), is defined as:

ρ spec zz (t, f ) |S z (t, f )| 2 . ( 5
)
Clearly, from the definition, there is no cross-terms effect present in STFT, hence in the SPEC. However, these distributions have very low TF resolution in comparison with the WVD. The low cost of implementation for the STFT, hence for the SPEC, in comparison with that for the WVD and, together with the advantage of being free of cross-terms, justifies the fact that the STFT is most used in practice, especially for speech or audio signals. But when it comes to FM-like signals, the WVD is preferred.

To combine the high resolution of the WVD while using the free cross-term property of the SPEC, the masked Wigner-Ville distribution (MWVD) is derived so that:

ρ mwvd zz (t, f ) ρ wvd zz (t, f ) • ρ spec zz (t, f ). (6)
There are many other useful TFDs in the literature, notably those that give high TF resolution while effectively minimizing the cross-terms, for example the B distribution [START_REF] Barkat | A high-resolution quadratic time-frequency distribution for multicomponent signal analysis[END_REF]. However, we only introduce here the TFDs above since they will be used in the later sections.

C. TF conditions on sources

Now, as we have introduced the concept of TF signal processing as a useful tool for analyzing nonstationary signals, some TF conditions need to be applied to the sources. Note that the TF method in [START_REF] Belouchrani | Blind source separation based on time-frequency signal representations[END_REF] does not work for UBSS because the mixing matrix is not invertible. In order to deal with UBSS, one often seeks for a sparse representation of the sources [START_REF] Bofill | Underdetermined blind source separation using sparse representations[END_REF]. In other words, if the sources can be sparsely represented in some domain, then the separation is to be carried out in that domain to exploit the sparseness.

1) TF-disjoint sources: Recently, there have been several UBSS methods, notably those in [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF] and [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF], in which the TF domain has been chosen to be the underlaying sparse domain. These two papers have based their solutions on the assumption that the sources are disjoint in the TF domain. Mathematically, if Ω 1 and Ω 2 are the TF supports of two sources s 1 (t) and s 2 (t) then Ω 1 ∩ Ω 2 = ∅. This condition can be illustrated in Figure 1. However, this is a rather strict assumption. A more practical assumption is that the sources are almost-disjoint in the TF domain [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF], allowing some small overlapping in the TF domain, for which the above two methods also worked. 2) TF-nondisjoint sources: In this paper, we want to relax the TF-disjoint condition by allowing the sources to be nondisjoint in the TF domain; as illustrated in Figure 2.

Ω 2 Ω 1 f frequency t Fig. 2. TF nondisjoint condition: Ω 1 ∩ Ω 2 = ∅
This is motivated by a drawback of the method in [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF]. Although this method worked well under the TF-almost-disjoint condition, it did not explicitly treat the TF regions where the sources were allowed to have some small overlapping. A point at the overlapping of two sources was assigned 'by chance' to belong to only one of the sources. As a result, the source that picks up this point will have some information of the other source while the latter loses some information of its own. The loss of information can be recovered to some extent by the interpolation at the intersection point using TF synthesis. However, for the other source, there is an interference at this point, hence the separation performance may degrade if no treatment is provided. If the number of overlapping points increases (i.e. the TF-almost-disjoint condition is violated), the performance of the separation is expected to degrade unless the overlapping points are treated. This paper will give such a treatment using subspace projection. Therefore, we will allow the sources to be nondisjoint in the TF domain; that is, multiple sources are allowed to be present at any point in the TF domain. However, instead of being inevitably nondisjoint, we limit ourselves by making the following constraint:

Assumption 2: The number of sources that contribute their energy at any TF point is strictly less than the number of sensors.

In other words, for the configuration of M sensors, there exist at most (M -1) sources at any point in the TF domain. For the special case when M = 2, Assumption 2 reduces to the disjoint condition.

We also make another assumption on the TF conditioning of the sources.

Assumption 3: For each source, there exists a region in the TF domain, where this source exists alone. Note that, this assumption is easily met and hence not restrictive for audio sources and FM-like signals. Also, it should be noted that this last assumption is, however, not a restriction on the use of subspace projection, because it will only be used later for the estimation of the mixing matrix. If otherwise, the mixing matrix can be obtained by another method, for example the one in [START_REF] Lathauwer | ICA techniques for more sources than sensors[END_REF], then Assumption 3 can be omitted.

III. CLUSTER-BASED TF-UBSS APPROACH FOR DISJOINT SOURCES

A. Quadratic TFD approach

In this section, we review a method proposed in [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF] based on the idea of clustering; hence, it is now referred to as the cluster-based quadratic TF-UBSS algorithm. For a signal vector z(t) = [z 1 (t), . . . , z N (t)] T , the Spatial Time Frequency Distribution (STFD) matrix is given by [START_REF] Belouchrani | Blind source separation based on time-frequency signal representations[END_REF]:

D zz (t, f )    ρ z 1 z 1 (t, f ) . . . ρ z 1 z N (t, f ) . . . . . . . . . ρ z N z 1 (t, f ) . . . ρ z N z N (t, f )    , (7)
where, for i, j ∈ N , ρ z i z j (t, f ) is the quadratic cross-TFD between z i (t) and z j (t) as obtained by ( 2), but with the first z being replaced by z i and the second by z j . By definition, the STFD takes into account the spatial diversity. By applying the STFD defined in [START_REF] O'grady | Survey of sparse and nonsparse methods in source separation[END_REF] on both sides of the BSS model in [START_REF] Nandi | Blind estimation using higher-order statistics[END_REF], we obtain the following TF-transformed structure:

D xx (t, f ) = AD ss (t, f )A H (8)
where D ss (t, f ) and D xx (t, f ) are, respectively, the source STFD matrix and mixture STFD matrix. [START_REF] Barkat | Algorithms for blind components separation and extraction from the time-frequency distribution of their mixture[END_REF]; noise thresholding by [START_REF] Linh-Trung | Underdetermined blind source separation of non-disjoint nonstationary sources in time-frequency domain[END_REF].

2) Noise thresholding and auto-source point selection by [START_REF] Linh-Trung | Underdetermined blind source separation of non-disjoint nonstationary sources in time-frequency domain[END_REF].

3) Vector clustering by [START_REF] Rickard | Desprit -histogram based blind source separation of more sources than sensors using subspace methods[END_REF] and k-means algorithm; source TFD estimation by [START_REF] Boashash | Time Frequency Signal Analysis and Processing: Method and Applications[END_REF]. 4) Source TF synthesis by [START_REF] Boudreaux-Bartels | Time-varying filtering and signal estimation using Wigner distributions[END_REF].

Let us call an auto-source TF point a point at which there is a true energy contribution/concentration of source or sources in the TF domain, and a cross-source point a point at which there is a 'false' energy contribution (due to the cross-term effect of quadratic TFDs). Note that, at other points with no energy contribution, the TFD value is ideally equal to zero. Under the assumption that all sources are disjoint in the TF domain, there is only one source present at any auto-source point. Therefore, the structure of

D xx (t, f ) is reduced to D xx (t a , f a ) = ρ s i s i (t a , f a ) a i a H i , ∀(t a , f a ) ∈ Ω i , (9)
where Ω i denotes, hereafter, the TF support of source s i (t).

The observation [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF] suggests that for all (t a , f a ) ∈ Ω i , the corresponding set of STFD matrices {D xx (t a , f a )} will have the same principal eigenvector a i . It is this observation that leads to the general separation method using quadratic TFDs in [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF]. Indeed, [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF] proposed several algorithms and pointed out that the choice of the TFD should be made carefully in order to have a 'clean' (cross-term free) TFD representation of the mixture, and chose the MWVD as a good candidate. This algorithm is summarized in Table I, and further detailed below for later use.

1) STFD mixture computation and noise thresholding: The STFD of the mixtures using the MWVD is computed by the following:

D wvd xx (t, f ) k,l = ρ wvd x k x l (t, f ) (10a) D stft xx (t, f ) k,l = S x k (t, f ), for k = l, 0, otherwise, (10b)
D mwvd xx (t, f ) = D wvd xx (t, f ) D stft xx (t, f ) 2 (10c)
In [START_REF] Barkat | Algorithms for blind components separation and extraction from the time-frequency distribution of their mixture[END_REF], k, l ∈ N , and denotes the Hadamard product.

2) Noise thresholding and auto-source point selection: A 'noise thresholding' procedure is used to keep only those points having sufficient energy, i.e. auto-source points. One way to do this is: for each time-slice (t p , f ) of the TFD representation, apply the following criterion for all the frequency points f q belonging to this time-slice:

If D mwvd xx (t p , f q ) max f { D mwvd xx (t p , f ) } > 1 , keep (t p , f q ), ( 11
)
where 1 is a small threshold (typically, 1 = 0.05). This 'hard thresholding' procedure has been preferred to the 'soft thresholding' using power-weighting of [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF] as it contributes also to reducing the computation complexity. The set of all the auto-source points is denoted by Ω. Since sources are TFdisjoint, we have Ω = N i=1 Ω i . This partition is found in the following way:

3) Vector clustering and source TFD estimation: For each point (t a , f a ) ∈ Ω, compute its corresponding spatial direction a(t a , f a )

a(t a , f a ) = diag D stft xx (t a , f a ) diag D stft xx (t a , f a ) , ( 12
)
and force it, without loss of generality, to have the first entry real and positive.

Having the set of spatial direction {a(t a , f a )|(t a , f a ) ∈ Ω} one can cluster them into N classes using any unsupervised clustering algorithm (see [START_REF] Frank | The data analysis handbook[END_REF] for different clustering methods). The clustering algorithm used in [START_REF] Linh-Trung | Separating more sources than sensors using time-frequency distributions[END_REF] is rather sensitive due to the threshold in use; a robust method should be investigated, and this deserves another contribution. If the number of sources has been well estimated, one can use the so-called k-means clustering algorithm [START_REF] Frank | The data analysis handbook[END_REF] to achieve a good clustering performance. The output of the clustering algorithm is a set of N classes {C i |i ∈ N }. Also, the collection of all the points that correspond to all the vectors in the class C i forms the TF support Ω i of the source s i (t).

Then, estimate the TFD of the source s i (t) (up to a scalar constant) as:

ρwvd si (t, f ) = trace D wvd xx (t, f ) , (t, f ) ∈ Ω i , 0, otherwise. (13)
4) Source TF synthesis: Having obtained the source TFD estimate ρwvd s i (t, f ), the estimation of the source s i (t) can be done through a TF synthesis algorithm. The method in [START_REF] Boudreaux-Bartels | Time-varying filtering and signal estimation using Wigner distributions[END_REF] is used for TF synthesis from a WVD estimate, based on the following inversion property of the WVD [START_REF] Boashash | Time Frequency Signal Analysis and Processing: Method and Applications[END_REF]:

x(t) = 1 x * (0) ∞ -∞ ρ wvd x ( t 2 , f ) e j2πf t df ,
which implies that the signal can be reconstructed to within a complex exponential constant e jα = x * (0)/|x(0)| given |x(0)| = 0. It can be observed that in this version of the quadratic TF-UBSS algorithm, the STFD matrices are not fully needed as only their diagonal entries are used in the algorithm. This should be taken into account to reduce the computational cost.

B. Linear TFD approach

As we have seen before, the STFT is often used for speech/audio signals because of its low computational cost. Therefore, in this section we briefly review the STFT method in [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF], and propose simultaneously a cluster-based linear TF-UBSS algorithm using the STFT to avoid some of the drawbacks in [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF].

First, under the transformation into the TF domain using the STFT, the model in (1) becomes:

S x (t, f ) = AS s (t, f ), (14)
where S x (t, f ) is the mixture STFT vector and S s (t, f ) is the source STFT vector. Under the assumption that all sources are disjoint in the TF domain, ( 14) is reduced to [START_REF] Frank | The data analysis handbook[END_REF]; noise thresholding by [START_REF] Griffin | Signal estimation from modified shorttime fourier transform[END_REF] 2) Vector clustering by [START_REF] Zibulevsky | Independent Component Analysis: Principles and Practice[END_REF] and [START_REF] Wax | Detection of signals by information theoretic criteria[END_REF].

S x (t a , f a ) = a i S s i (t a , f a ), ∀(t a , f a ) ∈ Ω i , ∀i ∈ N . ( 15
)
3) Source STFT estimation by (21). 4) Source TF synthesis by [START_REF] Griffin | Signal estimation from modified shorttime fourier transform[END_REF]. Now, in [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF], the structure of the mixing matrix is particular in that it has only 2 rows (i.e. the method uses only 2 sensors) and the first row of the mixing matrix contains all 1. Then, ( 15) is expanded to

S x 1 (t a , f a ) S x2 (t a , f a ) = 1 a 2,i S s i (t a , f a ), which results in a 2,i = S x2 (t a , f a ) S x1 (t a , f a ) . ( 16
)
Therefore, all the points for which the ratios on the right-hand side of ( 16) have the same value form the TF support Ω i of a single source, say s i (t). Then, the STFT estimate of s i (t) is computed by:

Ŝs i (t, f ) = S x1 (t, f ), ∀(t, f ) ∈ Ω i , 0, otherwise.
The source estimate ŝi (t) is then obtained by converting Ŝs i (t, f ) to the time domain using inverse STFT [START_REF] Griffin | Signal estimation from modified shorttime fourier transform[END_REF]. Note that, the extension of the UBSS method in [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF] to more than two sensors is a difficult task. Second, the division on the right-hand side of ( 16) is prone to error if the denominator is close to zero. To avoid the above mentioned problems, we propose here a modified version of the previous method referred to as the cluster-based linear TF-UBSS algorithm. In particular, from the observation [START_REF] Lathauwer | ICA techniques for more sources than sensors[END_REF], we can deduce the separation algorithm as shown next, and summarized in Table II.

1) Mixture STFT computation and noise thresholding: Compute the STFT of the mixtures, S x (t, f ), by applying (4) for each of the mixture in x(t), as follows:

S x i (t, f ) = ∞ -∞ x i (τ )h(τ -t)e -j2πf τ dτ, i = 1, . . . , M, (17a) S x (t, f ) = [S x 1 (t, f ), . . . , S x M (t, f )] T . ( 17b
)
Since the STFT is totally free of cross-terms, a point with a nonzero TFD value is ideally an auto-source point. Practically, we can select all auto-source points by only applying a noise thresholding procedure as that in the cluster-based quadratic TF-UBSS algorithm. In particular, for each time-slice (t p , f ) of the TFD representation, apply the following criterion for all the frequency points f k belonging to this time-slice

If S x (t p , f k ) max f { S x (t p , f ) } > 2 , then keep (t p , f k ), ( 18
)
where 2 is a small threshold (typically, 2 = 0.05). Then, the set of all selected points, Ω, is expressed by

Ω = N i=1 Ω i ,
where Ω i is the TF support of the source s i (t). Note that, the effects of spreading the noise energy while localizing the source energy in the time-frequency domain amounts to increasing the robustness of the proposed method with respect to noise. Hence, by equation ( 18) (or equation ( 11)), we would keep only time-frequency points where the signal energy is significant, the other time-frequency points are rejected, i.e. not further processed, since considered to represent noise contribution only. Also, due to the noise energy spreading, the contribution of the noise in the source time-frequency points is relatively, negligeable at least for moderate and high SNRs.

2) Vector clustering and source TFD estimation: The clustering procedure can be done in a similar manner as in the quadratic algorithm. First, we obtain the spatial direction vectors by:

v(t a , f a ) = S x (t a , f a ) S x (t a , f a ) , (t a , f a ) ∈ Ω, (19)
and force them, without loss of generality, to have the first entry real and positive. Next, we cluster these vectors into N classes {C i | i ∈ N }, using the k-means clustering algorithm. The collection of all points, whose vectors belong to the class C i , now forms the TF support Ω i of the source s i (t). Then, the column vector a i of A is estimated as the centroid of this set of vectors:

âi = 1 |C i | (t,f )∈Ω i v(t, f ), (20)
where |C i | is the number of vectors in this class. Therefore, we can estimate the STFT of each source s i (t) by:

Ŝsi (t, f ) = âH i S x (t, f ), ∀ (t, f ) ∈ Ω i , 0, otherwise. (21)
since, from (15), we have

âH i S x (t, f ) = âH i a i S si (t, f ) ≈ S si (t, f ), ∀ (t, f ) ∈ Ω i .
Note that the STFT is a particular form of wavelet transforms which have been used in [START_REF] Zibulevsky | Independent Component Analysis: Principles and Practice[END_REF] for the UBSS of image signals.

IV. SUBSPACE-BASED TF-UBSS APPROACH FOR NONDISJOINT SOURCES

We have seen the cluster-based TF-UBSS methods, using either quadratic TFDs such as the MWVD or linear TFDs such as the STFT, as summarized in Table I or Table II, respectively. These methods relied on the assumption that the sources were TF-disjoint, which has led to the enabling TFtransformed structures in [START_REF] Yilmaz | Blind separation of speech mixtures via time-frequency masking[END_REF] or [START_REF] Lathauwer | ICA techniques for more sources than sensors[END_REF]. When the sources are nondisjoint in the TF domain, then these equations are no longer true.

Under the TF-nondisjoint condition, stated in Assumption 2, we propose in this section two alternative methods, one for quadratic TFDs and the other for linear TFDs, for the UBSS problem using subspace projection.

A. Subspace-based quadratic TF-UBSS algorithm

Recall that the first two steps of the cluster-based quadratic TF-UBSS algorithm do not rely on the assumption of TFdisjoint sources (see Table I). Therefore, we can reuse these steps to obtain the set of auto-source points Ω. Now, under the TF-nondisjoint condition, consider an auto-source point (t b , f b ) ∈ Ω such that there are K sources, K < M , present at this point. Our goal is to identify the sources present at (t b , f b ) and to estimate the energy each of these sources contributes.

Denote α 1 , . . . , α K ∈ N the indices of the sources present at (t b , f b ), and define the following:

s = [s α1 (t), . . . , s α K (t)] T , (22a) Ã = [a α1 , . . . , a α K ]. (22b)
Then, under Assumption 2, ( 8) is reduced to

D wvd xx (t b , f b ) = ÃD ss (t b , f b ) ÃH , ( 23
)
Consequently, given that D ss is of full rank, we have

Range {D xx (t b , f b )} = Range{ Ã}. ( 24
)
Let P be the orthogonal projection matrix onto the noise subspace of D wvd xx (t b , f b ). Then, from (24), we obtain:

P = I -VV H , ( 25
)
and

Pa i = 0, ∀ i ∈ {α 1 , . . . , α K } , Pa i = 0, ∀ i ∈ N \ {α 1 , . . . , α K } (26)
In ( 25), V is the matrix formed by the K principal singular eigenvectors of D xx (t b , f b ).

Assuming that A has been estimated by some method, the observation in (26) enables us to identify the indices α 1 , . . . , α K ; and hence, the sources present at (t b , f b ). In practice, to take into account the estimation noise, one can detect these indices by detecting the K smallest values from the set { Pa i | i ∈ N }, as mathematically expressed by:

{α 1 , . . . , α K } = arg min K { Pa i | i ∈ N } , ( 27
)
where min K denotes the minimization to obtain the K smallest values. The TFD values of the K sources at (t b , f b ) are estimated as the diagonal elements of the following matrix:

Dss (t b , f b ) ≈ Ã# D xx (t b , f b ) Ã# H , ( 28
)
where the superscript ( # ) is the Moore-Penrose's pseudoinversion operator.

Here, we propose also an estimation method for A by using Assumption 3. This assumption states that, for each source s i (t), there exists a TF region R i where s i (t) exists alone. In other words, R i contains all the single-source auto-source points of s i (t). Therefore, we can reuse the observation (9) in the TF-disjoint case, but for some TF regions, as below: [START_REF] Barkat | Algorithms for blind components separation and extraction from the time-frequency distribution of their mixture[END_REF].

D xx (t, f ) = ρ s i s i (t, f )a i a H i , ∀(t, f ) ∈ R i , ∀i ∈ N . The union of these regions, R = N i=1 R i , is detected by the following: If λ max {D wvd xx (t, f )} trace{D wvd xx (t, f )} -1 < 3 , then (t, f ) ∈ R, ( 29
)
2) Noise thresholding and auto-source point selection by [START_REF] Linh-Trung | Underdetermined blind source separation of non-disjoint nonstationary sources in time-frequency domain[END_REF].

3) Single-source auto-source point selection by (29); mixing matrix estimation by ( 30) and (31) 4) For all auto-source points, perform subspace-based TFD estimation of sources by ( 25), ( 27) and (28) 5) Source TF synthesis by [START_REF] Boudreaux-Bartels | Time-varying filtering and signal estimation using Wigner distributions[END_REF].

where 3 is a small threshold value (typically, 3 ≤ 0.1) and λ max {D wvd xx (t, f )} denotes the maximum eigenvalue of D wvd xx (t, f ). Then, we can apply the same vector clustering procedure as in Section III-A.3 to estimate A. In particular, we first obtain all the spatial direction vectors:

a(t, f ) = diag D stft xx (t, f ) diag D stft xx (t, f ) , ∀(t, f ) ∈ R. ( 30
)
Next, we cluster these vectors into N classes {D i |i ∈ N } using the k-means clustering algorithm. The collection of all points, whose vectors belong to the class D i , now forms the TF region R i of the source s i (t). Finally, the column vectors A are estimated as the centroid vectors of these classes as:

âi = 1 |D i | (t,f )∈R i a(t, f ), ∀i ∈ N (31)
where D i is the number of points in R i .

Table III gives a summary of the subspace-based quadratic TF-UBSS algorithm.

B. Subspace-based linear TF-UBSS algorithm

Similarly, we propose here a subspace-based linear TF-UBSS algorithm for TF-nondisjoint sources using STFT. We also use the first step of the cluster-based linear TF-UBSS algorithm (see Table II) to obtain all the auto-source points Ω. Under the TF-nondisjoint condition, consider an autosource point (t b , f b ) ∈ Ω at which there are K sources s α 1 (t), . . . , s α K (t) present, with K < M . Then, ( 8) is reduced to the following

S x (t b , f b ) = ÃS s(t b , f b ), ∀(t b , f b ) ∈ Ω ( 32
)
where à and s are as previously defined in (22). Let Q represent the orthogonal projection matrix onto the noise subspace of Ã. Then, Q can be computed by:

Q = I -Ã ÃH Ã -1 ÃH . ( 33
)
We have the following observation:

Qa i = 0, i ∈ {α 1 , . . . , α K } Qa i = 0, i ∈ N \{α 1 , . . . , α K } . ( 34
)
If A has already been estimated by some method, then this observation gives us the criterion to detect the indices 20) and (37), and k-means algorithm. 4) For all auto-source points, perform subspace-based TFD estimation of sources by (33), ( 35) and (36). 5) Source TF synthesis by [START_REF] Griffin | Signal estimation from modified shorttime fourier transform[END_REF]. α 1 , . . . , α K ; and hence, the contributing sources at the autosource point (t b , f b ). In practice, to take into account noise, one detects the column vectors of à minimizing:

{α 1 , . . . , α K } = arg min β 1 ,...,β K QS x (t, f ) | Ãβ (35)
where

Ãβ = [a β1 , . . . , a β K ].
Next, TFD values of the K sources at TF point (t, f ) are estimated by:

Ŝs (t, f ) ≈ Ã# S x (t, f ). (36)
Here we propose a method for estimating the mixing matrix A. This is performed by clustering all the spatial direction vectors in [START_REF] Zibulevsky | Independent Component Analysis: Principles and Practice[END_REF] as for the preview TF-UBSS algorithm. Then within each class C i we eliminate the far-located vectors from the centroid (in the simulation we estimate vectors v(t, f ) such that:

v(t, f ) -âi > 0.8 max v(t,f )∈Ωi v(t, f ) -âi , (37)
leading to a size-reduced class Ci . Essentially this is to keep the vectors corresponding to the TF region R i , which are ideally equal to the spatial direction a i of the considered source signal. Finally, the i th column vector of A is estimated as the centroid of Ci .

Table IV provides a summary of the subspace projection based TF-UBSS algorithm using STFT.

V. DISCUSSION

We discuss here certain points relative to the proposed TF-UBSS algorithms and their applications.

1) Number of sources: The number of sources N is assumed known in the clustering method (k-means) that we have used. However, there exist clustering methods [START_REF] Frank | The data analysis handbook[END_REF] which perform the class estimation as well as the estimation of the number N . In our simulation, we have observed that most of the time the number of classes is overestimated, leading to poor source separation quality. Hence, robust estimation of the number of sources in the UBSS case remains a difficult open problem that deserves particular attention in future works.

2) Number of overlapping sources: In the subspace-based approach, we have to evaluate the number K of overlapping sources at a given TF point. This can be done by finding out the number of non-zero eigenvalues of D wvd xx (t, f ) using criteria such as Minimum Description Length (MDL) or Akaike Information Criterion (AIC) [START_REF] Wax | Detection of signals by information theoretic criteria[END_REF]. It is also possible to consider a fixed (maximum) value of K that is used for all auto-source TF points. Indeed, if the number of overlapping sources is less than K, we would estimate close-to-zero source STFT values. For example, if we assume K = 2 sources are present at a given TF point while only one source is effectively contributing, then we estimate one close-to-zero source STFT value. This approach increases slightly the estimation error of the source signals (especially at low SNRs) but has the advantage of simplicity compared to using information theoretic-based criterion. In our simulation, we did choose this solution with K = 2 or K = 3.

3) Quadratic versus linear TFDs: We have proposed two algorithms using quadratic and linear TFDs. The one using the quadratic TFD should be preferred when dealing with FM-like signals and for small or moderate sample sizes. For audio source separation often the case the sample size is large, and hence, to reduce the computational cost one should prefer the linear TFD based UBSS algorithm. Overall, the quadratic version performs slightly better than the linear one but costs much more in computations.

4) Separation quality versus number of sources: Although we are in the underdetermined case, the number of sources N should not exceed too much the number of sensors. Indeed, when N increases, the level of source interference increases, and hence, the source disjointness assumption is ill-satisfied. Moreover, for a large number of sources, the likelihood of having two sources closely spaced, i.e. such that the spatial directions a i and a j are 'close' to linear dependency, increases. In that case, vector clustering performance degrades significantly. In brief, sparseness and spatial separation are the two limiting factors against increasing the number of sources. Figure 8 illustrates the performance degradation of source separation versus the number of sources.

VI. SIMULATION RESULTS

A. Simulation results of subspace-based TF-UBSS algorithm using STFT

In the simulations, we use a uniform linear array of M = 3 sensors. It receives signals from N = 4 independent speech sources in the far field from directions θ 1 = 15, θ 2 = 30, θ 3 = 45 and θ 4 = 75 degrees respectively. The sample size is T = 8192 samples. In Figure 3, the top four plots represent the TF representation of the original sources signal, the middle three plots represent the TF representation of the M mixture signals and the bottom four plots represent the TF representation of the estimate of sources by the subspacebased algorithm using STFT (Table IV). Figure 4 represents the same disposition of signals but in the time domain. In Figure 5, we compare the separation performance obtained by the subspace-based algorithm with K = 2 and the clusterbased algorithm (Table II). It is observed that subspace-based algorithm provides much better separation results than those obtained by the cluster-based algorithm.

In the subspace-based method, one first needs to estimate the mixing matrix A. This is done by the cluster-based method presented previously. The plot in Figure 6 represents the normalized estimation error of A versus the SNR in dB. Clearly, the proposed estimation method of the mixing matrix provides satisfactory performance, while the plot in Figure 7 presents the separation performance when using the exact matrix A compared to that obtained with the proposed estimate Â. Figure 8 illustrates the rapid degradation of the separation quality when we increase the number of sources from N = 4 to N = 7. This confirms the remarks made in Section V.

(a) S s 1 (t, f ) (b) S s 2 (t, f ) (c) S s 3 (t, f ) (d) S s 4 (t, f ) (e) Sx 1 (t, f ) (f) Sx 2 (t, f ) (g) Sx 3 (t, f ) (h) S ŝ1 (t, f ) (i) S ŝ2 (t, f ) (j) S ŝ3 (t, f ) (k) S ŝ4 (t, f )
In Figure 9, we compare the performance obtained with the subspace-based method for K = 2 and K = 3. In that experiment, we have used M = 4 sensors and N = 5 source signals. One can observe that, for high SNRs, the case of K = 3 leads to a better separation performance than for the case of K = 2. However, for low SNRs, a large value of K increases the estimation noise (as mentioned in Section V) and hence degrades the separation quality.

B. Simulation results of subspace-based TF-UBSS algorithm using STFD

In this simulation, we use a uniform linear array of M = 3 sensors with half wavelength spacing. It receives signals from N = 4 independent LFM sources, each has 256 samples, in the presence of additive Gaussian noise where the SNR=20 dB.

We compare the cluster-based (Table I) and the proposed subspace-based (Table III) TF-UBSS algorithms. Fig- ures 10-(a,d,g,j) represent the TFDs (using WVD) of the four sources. Figures e,h,k) show the estimated source TFDs using the cluster-based algorithm, whereas Figures f,i,l) are those obtained by the subspace-based algorithm.

From Figures 10-(b,e) we can see that the overlapping points between source s 1 (t) and source s 2 (t) were picked up by source s 2 (t) with the cluster-based algorithm. On the other hand, using the subspace-based algorithm, the intersection points have been redistributed to the two sources (Figure 10-(c,f)).

In general, the overlapping points in the nondisjoint case have been explicitly treated. This provides a visual performance comparison.

In Figure 11, we compare the statistical separation performance between the subspace-based algorithm and the clusterbased algorithm using STFD, evaluated over 1000 Monte-Carlo runs.

One can also notice that the gain here is smaller than the one obtained previously for audio sources. This is due to the fact that the overlapping region of the considered signals is smaller. This result confirms the previous visual observation with respect to the performance gain in favor of our subspacebased method.

VII. CONCLUSIONS

This paper introduces new methods for the UBSS of TFnondisjoint nonstationary sources using time-frequency repre- sentations. The main advantages over the proposed separation algorithms are, first, a weaker assumption on the source 'sparseness', i.e. the sources are not necessarily TF-disjoint, and second, an explicit treatment of the overlapping points using subspace projection, leading to significant performance improvements. Simulation results illustrate the effectiveness of our algorithms in different scenarios compared to those existing in the literature.

Fig. 1 .

 1 Fig. 1. Source TF-disjoint condition: Ω 1 ∩ Ω 2 = ∅ (when Ω 1 ∩ Ω 2 ≈ ∅, sources are said to be TF-almost-disjoint).

Fig. 3 .

 3 Fig. 3. Simulated example (viewed in TF domain) for the subspace-based TF-UBSS algorithm with STFT in the case of 4 speech sources and 3 sensors. The top four plots represent the original source signals, the middle three plots represent the 3 mixtures, and the bottom four plots represent the source estimates.

Fig. 4 .

 4 Fig. 4. Simulated example (viewed in time domain) for the subspace-based TF-UBSS algorithm with STFT in the case of 4 speech sources and 3 sensors. The top four plots represent the original source signals, the middle three plots represent the 3 mixtures, and the bottom four plots represent the source estimates.

Fig. 5 .

 5 Fig. 5. Comparison between subspace-based and cluster-based TF-UBSS algorithms using STFT: normalized MSE (NMSE) versus SNR for 4 speech sources and 3 sensors.

Fig. 6 .

 6 Fig. 6. Mixing matrix estimation: normalized MSE versus SNR for 4 speech sources and 3 sensors.

Fig. 7 .

 7 Fig. 7. Comparison, for the subspace-based TF-UBSS algorithm using STFT, when the mixing matrix A is known or unknown: NMSE of the source estimates.

Fig. 8 .

 8 Fig. 8. Comparison between subspace-based and cluster-based TF-UBSS algorithms using STFT: NMSE versus number of sources.

Fig. 9 .

 9 Fig. 9. Comparison between subspace-based and cluster-based TF-UBSS algorithms using STFT: NMSE of the source estimates for different sizes of the projector, for the case of 5 sources and 4 sensors.

Fig. 10 .

 10 Fig.10.Simulated example (viewed in TF domain) for the subspacebased TF-UBSS algorithm with STFT in the case of 4 LFM sources and 3 sensors. From left to right, the figures respectively represent the original source TF signatures, the estimated source TF signatures using the clusterbased algorithm, and the estimated source TF signatures using the subspacebased algorithm.

Fig. 11 .

 11 Fig. 11. Comparison between subspace-based and cluster-based TF-UBSS algorithms using STFD: normalized MSE (NMSE) versus SNR for 4 LFM sources and 3 sensors.

TABLE I CLUSTER

 I -BASED QUADRATIC TF-UBSS ALGORITHM USING STFD 1) Mixture STFD computation by

In fact, the STFT does not represent an energy distribution of the signal in the TF plane. However, for simplicity, we still refer to it as a TFD.

Abdeldjalil Aïssa-El-Bey was born in Algiers, Algeria, in 1981. He received the State Engineering degree from École Nationale Polytechnique (ENP), Algiers, Algeria, in 2003, the M.S. degree in signal processing from Supélec and Paris XI University, Orsay, France, in 2004. Currently he is working towards the Ph.D. degree at the Signal and Image Processing Department of École Nationale Supérieure des Télécommunications (ENST) Paris, France. His research interests are blind source separation, blind system identification and equalization, statistical signal processing, wireless communications and adaptive filtering.


CONFIDENTIALITE DES DONNEES VOUS CONCERNANT

Dans le cadre de la recherche impliquant la personne humaine à laquelle [nom du promoteur] vous propose de participer, un traitement de vos données personnelles va être mis en oeuvre pour permettre d'analyser les résultats de la recherche au regard de l'objectif de cette dernière. A cette fin, les données médicales vous concernant et les données relatives à vos habitudes de vie, ainsi que, dans la mesure où ces données sont nécessaires à la recherche, vos origines ethniques ou des données relatives à votre vie sexuelle, seront transmises au Promoteur de la recherche ou aux personnes ou sociétés agissant pour son compte, en France ou à l'étranger. Ces données seront identifiées par un numéro de code et/ou vos initiales ou les trois premières lettres de votre nom Vous pouvez également accéder directement ou par l'intermédiaire d'un médecin de votre choix à l'ensemble de vos données médicales en application des dispositions de l'article L 1111-7 du Code de la Santé Publique. Ces droits s'exercent auprès du médecin qui vous suit dans le cadre de la recherche et qui connaît votre identité.

Nous vous informons que vous serez inscrit dans le fichier national des personnes qui se prêtent à des recherches impliquant la recherche humaine prévu à l'article L.1121-16 du code de la santé publique. Vous avez la possibilité de vérifier auprès du ministre chargé de la santé l'exactitude des données vous concernant présentes dans ce fichier et la destruction de ces données au terme du délai prévu par le Code de la santé publique.

  [à préciser selon les cas]. Ces données pourront également, dans des conditions assurant leur confidentialité, être transmises aux autorités de santé françaises ou étrangères et à d'autres entités de [nom du Promoteur]. Ces données pourront également être utilisées lors de recherches ultérieures pour des fins scientifiques ; en cas de retrait de consentement et sauf précision de votre part, les données collectées jusqu'à cette date pourront être utilisées. Conformément aux dispositions de loi relative à l'informatique aux fichiers et aux libertés, le promoteur a procédé à une déclaration auprès de la Commission Nationale de l'Informatique et des Libertés (CNIL). Vous disposez d'un droit d'accès et de rectification. Vous disposez également d'un droit d'opposition à la transmission des données couvertes par le secret professionnel susceptibles d'être utilisées dans le cadre de cette recherche et d'être traitées. Les données sont stockées [précisez lieu et organisme]. Le promoteur procédera à la collecte de ces données [précisez les modalités].