O EURASIP Journal on

Advances in Signal Processing

a SpringerOpen Journal

RESEARCH Open Access

Generalized generating function with tucker decomposition and alternating least squares for underdetermined blind identification

Fanglin Gu1*, Hang Zhang1, Wenwu Wang2 and Desheng Zhu1

Abstract

Generating function (GF) has been used in blind identification for real-valued signals. In this paper, the definition of GF is first generalized for complex-valued random variables in order to exploit the statistical information carried on complex signals in a more effective way. Then an algebraic structure is proposed to identify the mixing matrix from underdetermined mixtures using the generalized generating function (GGF). Two methods, namely GGF-ALS and GGF-TALS, are developed for this purpose. In the GGF-ALS method, the mixing matrix is estimated by the decomposition of the tensor constructed from the Hessian matrices of the GGF of the observations, using an alternating least squares (ALS) algorithm. The GGF-TALS method is an improved version of the GGF-ALS algorithm based on Tucker decomposition. More specifically, the original tensor, as formed in GGF-ALS, is first converted to a lower-rank core tensor using the Tucker decomposition, where the factors are obtained by the left singular-value decomposition of the original tensor's mode-3 matrix. Then the mixing matrix is estimated by decomposing the core tensor with the ALS algorithm. Simulation results show that (a) the proposed GGF-ALS and GGF-TALS approaches have almost the same performance in terms of the relative errors, whereas the GGF-TALS has much lower computational complexity, and (b) the proposed GGF algorithms have superior performance to the latest GF-based baseline approaches.

Keywords: Blind identification; Generalized generating function; Tensor decomposition; Tucker decomposition; Underdetermined mixtures

1. Introduction

Blind identification (BI) of linear mixtures has recently attracted intensive research interest in many fields of signal processing including blind source separation (BSS). This work is devoted to BI of underdetermined mixtures with complex sources. Underdetermined mixtures are commonly encountered in many practical applications, such as in the radio communication context, where the reception of more sources than sensors becomes increasingly possible with the growth in reception bandwidth. In these applications, one often has to also deal with complex sources. One reason is that the communication signals are usually complex-valued such as the quadrature-amplitude modulation (QAM) signal,

* Correspondence: gu.fanglin@gmail.com

1College of Communication Engineering, PLA University of Science &

Technology, Nanjing 210007, People's Republic of China

Fulllist of author information is available at the end of the article

and minimum-shift keying (MSK) signal. Another reason is that frequency domain methods are often used for blind separation or identification from convolutive mixtures due to its computational efficiency [1,2], while the objective functions used in the frequency domain are usually defined on complex-valued variables.

A large number of methods for BI of underdetermined mixtures start from the assumption that the sources are sparse by nature (i.e., in its own domain such as the time domain) or could be made sparse in another domain (e.g., a transform domain). A predefined transform such as short-time Fourier transform (STFT) or a learned transform using, e.g., simultaneous codeword optimization (SimCO), is usually applied to sparsify the data [3,4] if the signal by nature is not sparse. Due to the sparsity of the sources, the scatter plot typically shows high signal values in the directions of the mixing vectors, which can be localized by using some clustering techniques [5,6]. It should

ringer

© 2013 Gu et al.; licensee Springer. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.

be noted that although some signals such as speech signals have some degree of sparsity in one domain or another, many other signals such as the majority of communication signals do not possess such a property. Hence, it is necessary to develop BI methods for the underdetermined mixtures that do not impose any sparsity constraint on the sources.

To this aim, many methods for BI of underdetermined mixtures turn to the use of various decomposition methods based on different data structures such as correlation [7,8] and higher-order cumulant [9-14] matrices. The main idea of these algorithms is to construct a tensor based on the cumulants of the observations and then to estimate the mixing matrix by the decomposition of such a tensor. This is notably the case for second-order blind identification of underdetermined mixtures (SOBIUM) [7], fourth-order blind identification of underdetermined mixtures (FOBIUM) [9], fourth-order-only blind identification (FOOBI) [10], FOOBI-2 [10], and blind identification of mixtures of sources using Redundancies in the daTa Hexacovariance matrix (BIRTH) [11,12] algorithms, which use second-order statistics tensors and fourth- and sixth-order cumulant tensors, respectively. A family of the methods named blind identification of over-complete mixtures of sources (BIOME) is proposed in [13], based on the even-order cumulants of the observations. However, all the methods proposed in [7-14] exploit only the statistical information contained in the data measured by second-order or higher-order statistics.

In order to exploit statistical information more effectively, a family of BI approaches was proposed in [15-18] by exploiting the statistical information with the characteristic function (CAF) or generating function (GF). In these works, the authors showed that the mixing matrix can be estimated up to trivial scaling and permutation indeterminacies by decomposing the tensor composed of partial derivatives of the GF. It is worth mentioning that the algorithms in [15-17] have been only applied to BI problems involving real-valued sources. In [18], the CAF approach was extended to the case of mixtures of complex-valued sources, which often occurs in digital communications. However, extra effort is required to obtain the correct real and imaginary combination of the mixing matrix since the real and imaginary parts of the mixing matrix are treated separately, leading to an increased computational cost due to the increased dimension of the matrix that needs to be processed. In this paper, we propose the Generalized Generating Function (GGF) to exploit the statistical information carried on the complex random variable. We show that the proposed GGF can exploit the statistical information carried on complex random variables in a more effective way than the GF presented in [18] due to the algebraic structure adopted

by GGF (as detailed in Algebraic structure based on generalized generating function). Furthermore, a simple method for the mixing matrix estimation is derived based on tensor decomposition where the tensor is composed of the Hessian matrices of the GGF of the observations.

The remainder of this paper is organized as follows. In Problem formulation the BI problem is formulated and relevant assumptions are presented. In Algebraic structure based on generalized generating function we firstly generalize the definition of the GF for complex-valued random variables and then derive the corresponding core equation for BI. In Blind identification based on tensor decomposition, the GGF-ALS and GGF-TALS approaches are developed for the estimation of the mixing matrix. In the GGF-ALS algorithm, the mixing matrix is estimated by directly decomposing the tensor, constructed from the Hessian matrices of the second GGF of the observations, using the alternating least squares (ALS) algorithm. In the GGF-TALS algorithm, the Tucker decomposition is firstly applied to convert the original tensor to a lower-order core tensor, then the mixing matrix is obtained by decomposing the core tensor with the ALS algorithm. Furthermore, the factors of Tucker decomposition are obtained by the left singular vectors of the original tensor's mode-3 matrix. Computer simulations are used to illustrate the performance of the proposed GGF approaches in Simulations and analysis. Finally, the paper is concluded in Conclusions.

2. Problem formulation

Considering the following linear mixture model

z (t) = As (t) + w(t) (1)

where the stochastic vector z(t) e CQ represents the observation signals, s(t) e CP contains the unobserved source signals, and w(t) e C Q denotes additive noise. From now on, the noise w(t) is simply ignored for convenience, except when running computer experiments. The unknown mixing matrix A e C Q x P characterizes the way that the sources are acquired by the sensors. BI aims to estimate the mixing matrix from the observations based on the assumption that the source signals are statistically independent. The mixing matrix obtained may in turn be used to estimate the original source signals from the observations. In addition, we make the following assumptions:

(i). The mixing matrix A is of full (row) rank.

(ii).The number p of sources is known.

(iii). The number of sensors is smaller than the number of sources, i.e., Q < P.

3. Algebraic structure based on generalized generating function

3.1. Core equation based on generalized generating function

For a real stochastic vector x e R Q, the GF 0x(u) obtained by dropping the term of the square root of (-1) in the exponent of a CAF is defined as

&(u) = E[exp(uTx)],u e RQ,

fz(u) = E[exp(uHz)], u e CQ

Note that the statistical information exploited by GF/ GGF is related to the number of processing points. Theoretically, a complete statistical description of the probability density function requires the evaluation of the GF/GGF at all (infinitely many) possible processing points. However, this often becomes computationally infeasible. In practice, such statistical information is obtained approximately by the evaluation of GF/GGF at a finite number of processing points. Hence, in comparison with the GF presented in [18], the GGF defined in (4) can exploit the statistical information carried on the complex variables more effectively when the number of the processing points stays the same, thanks to the incorporation of the imaginary part of the exponent to the function. Furthermore, as compared with the use of the GF in (3), using the GGF in (4) offers a simpler way for the estimation of the mixing matrix due to the exploitation of an elegant algebraic structure.

Now, replacing z by its model and neglecting the noise contribution yield

fz(u) = E[exp(uHAs)] = f(AHu), u e CQ

Defining 0z(u) = log fz(u), which is often referred to as the 'second' GGF, and using the source independence

property, the second GGF of the observations can be rewritten as

pz(u) = ELa

where u e R Q is an arbitrary vector referred to as a processing point [15], and E[ ] denotes an expectation operator. Nevertheless, both the observation vector z and the mixing matrix A discussed in this paper belong to the complex field. Hence, a definition of GF for complex variables is required. One such definition has been presented in [18] as

& (R(u), 3(u)) = E[ exp(R(uHz))]

= E[exp(RT(u)R(z) + 3T(u)3(z))],u e CQ,

where R() and 1() denote taking the real and imaginary parts from their arguments (i.e., complex-valued vectors) to form a real-valued vector of the same dimension. It is actually defined by assimilating C to R2. Thus the GF of a complex variable in (3) is defined as a function of the real and imaginary parts. In this paper, we generalize the definition of GF for real stochastic vector in (2) to the following complex form

Consequently, by calculating the derivative of the conjugate gradient of 0z(u) with respect to u (more details can be found in Appendix 1), we can obtain the following core equation for the Hessian matrix tyz(u),

Vz (u) = AVS(AH u)AH

Vz (u) = Vut [Vu. pz (u)] = d^r

d(pz (u)

where (■) denotes the conjugate operator. It is necessary to point out that Vs(AHu) is a diagonal matrix (more details can be found in Appendix 2).

3.2. Estimating ^z(u)

In this subsection, we discuss how to consistently estimate the Hessian matrix tyz(u). Under the ergodicity assumption, the mean value of a random variable can be estimated by a time average. Hence, we can estimate the GGF of the observation vector as

f z (u) = „Y. T=1 exp(uHz(i))

The conjugate gradient qz(u) = d fz(u)/d u* which is a Q x 1 vector can be estimated by

4 (u)= L 1exK uHz(t^ z(t)

Similarly, the gradient Zz(u) = d qz(u)/5 u which is a Q x Q matrix can be estimated by

Zz(u) =3 ^ == 7E 1 exp (uH z(i))z(i)zH(i)

Based on the above analysis, the Hessian matrix tyz(u) of the second GGF 0z(u) can be obtained as

Vz (u) =

Zz (u) C (u& H(u) f z (u) f 2(u)

4. Blind identification based on tensor decomposition

In this section, the GGF-ALS and GGF-TALS algorithms are developed for the estimation of the mixing matrix. In the GGF-ALS algorithm, the mixing matrix is estimated by decomposing the tensor which is formed of the Hessian matrices of the GGF of the observations. An improved version of the GGF-ALS algorithm, i.e., the GGF-TALS algorithm, is also developed, where the

Tucker decomposition is firstly employed to convert the original tensor as used in the GGF-ALS algorithm to a lower-rank core tensor, and the mixing matrix is then estimated by decomposing the core tensor with the ALS algorithm.

4.1. The GGF-ALS algorithm

Evaluating the Hessian matrices %z(u) of GGF at a series of processing points ui, u2,..., uK and using Equation (7), one can obtain the following joint diagonalization (JD) problem

Vz (ui) = (AHui) AH Vz K) = A^s(AHuK )AH

mijk = iatla*iidkl-, which we write as

M = YH=1 al°a/°dl'

minllM-V~~ ap °a*°dpl

A,D II ¿—ip=ir p

line search (ELS) [2i] and extrapolating search direction (ESD) [22], are proposed to accelerate the rate of convergence of the ALS. Hence, the ALS is chosen here to compute the CAND.

To a large extent, the practical importance of tensor decomposition stems from its uniqueness properties. It is clear that the tensor decomposition can only be unique up to a permutation of the rank-i terms and scaling of the factors of the rank-i terms. Therefore, we consider the tensor decomposition (ii) as essentially unique if any other matrix pair A and D that satisfies (ii) is related to A and D via

A = A PAi; D = DPA2

in which %s(AHuk) is diagonal, k =i,...,K. The problem we need to address is to estimate the mixing matrix A based on the set {%z(ui), •••, %z(uK)}. For the determined/overde-termined case, it is obvious that JD methods, such as the AC-DC method [i9], can be used to estimate the mixing matrix. However, this method does not work when Q < P i.e., in the underdetermined case.

As shown in [7], the JD problem (9) can be seen as a particular case of the parallel factor (PARAFAC) decomposition, also known as canonical decomposition (CAND), of the third-order tensor M e cQxQxK built by stacking the K matrices %x(uk) along the third mode. Specifically, the tensor M e CQxQxK is built by stacking %z(ui), M,..., %z(uk) as follows: (M)^ = (%(uk))j, i =,-,Q, j = i,...,Q, k = i,...,K. Define a matrix D e C KxP by (D)kl = (%s(AHuk)ll, l = i,...,P, k = i,...,K. Then we have

with Ai, A2 e C PxP being diagonal matrices, satisfying AiAiA2 = I, and P e RPxP being a permutation matrix.

The k-rank. The Kruskal rank or k-rank of a matrix A, denoted by ka, is the maximal number A such that any set of A columns of A is linearly independent.

Theorem 1. The tensor decomposition of (11) is essentially unique if [23]

2ka + KD > 2(P + i)

where ° denotes the tensor outer product, and al and dl represent the lth column of A and D, respectively. In this way, the mixing matrix A can be estimated by solving the following problem. Given the third-order tensor M e CQxQxK, we can compute its CAND with P components of the rank-one tensors that best approximates M ,

where llllF is the Frobenius norm.

Several algorithms exist for the computation of tensor decomposition. The standard way for computing the tensor decomposition is by using an 'ALS' algorithm [20]. Several improved versions, such as the enhanced

We call a property generic when it holds with probability one. Generically, the mixing matrix is of full rank and of full k-rank when the parameters it involves are drawn from continuous probability densities. Hence, in practice, ka = min(Q, P) and kd = min(K, P).

In summary, we come to the following conclusion: when Q > P, P > 2, then the generic essential uniqueness is guaranteed for K > 2; when Q < P and if K > P, then the generic essential uniqueness is guaranteed for P < 2Q - 2, if K < P, then the generic essential uniqueness is guaranteed for P < Q - i+K/2.

4.2. The GGF-TALS algorithm

In order to ensure the GGF-ALS algorithm has robust performance, a large value is often chosen for K in the tensor M e CQxQxK. Nevertheless, this will lead to a heavy computational load. To reduce the computational complexity, the Tucker decomposition [23-25] is firstly applied to represent the tensor M as a lower-rank core tensor, whose size is much smaller than the tensor M . Then the ALS algorithm is used to perform the CAND of the core tensor. In this way, the computational complexity can be reduced dramatically.

Matricization. Matricization, also known as unfolding or flatting, is the process of turning an N-way tensor into a matrix. The mode-n matricization of a tensor F e c/ix/2"x/N is denoted by Fw which arranges the mode-n fibers to be the columns of the resulting

matrix, that is, the tensor element (¿i, to the matrix element (¿n, j) where

j = 1 + -1)Jk

-, iN) is mapped

k=1 k^n

Jk = ![ I»

(FxH U)

L =MV^

Tucker decomposition is a form of higher-order principal component analysis (PCA). It decomposes a tensor into a core tensor multiplied by a matrix along each mode as shown in Figure 1. Thus, in the three-way case where F e C/1X/2X/3, we have

F = Gx1A(1) x2A(2) x3A(3) + U

= TJ1 T/2 T/3 g.

z^/, =1 Z^/, =1 Z^/3 =1«

(1) (2) (3) , a 'oa °a:

/1/2/3 /1 /2 /3

where GeCJ 1X/2X/3 is the core tensor; A(1) e C1 X/ 1, A(2) e Cl2X/2, A(3) e C/3X/3 are the factor matrices (which are usually column unitary for real-valued data) and can

be thought of as the principal components in each mode; and U e C/1X/2X/3 represents errors or noise.

With the help of Tucker decomposition, the tensor M, which is composed of the Hessian matrices tyz(uk) of the second GGF of the observations, can be compressed. Since the mixing matrix A is of full rank and the Hessian matrices tys(AHuk) of the second GGF of the sources are diagonal, the mode-n (n = 1, 2, 3) matrices of tensor M are M(1) e C Q x QK, M(2) e C Q x QK, and M(3) e CK x QQ respectively. Meanwhile, the n-rank of

the M(1) and M

The «-rank. Let F be an Nth-order tensor of size I1 x I2 — x IN. Then the n-rank of F , denoted by rankn(F), is the column rank of F(n). In other words, if we let rn = rankn(F) for n = 1,...,N, then we can say that F is a rank r1, • • •, rN tensor.

The n-mode product. The n-mode (matrix) product of a tensor F e C/1X/2"'X/N with a matrix U/Xln is denoted by Fxn U and is of size I1 x ••• xIn_ 1 x/xIn + 1 x • IN. Elementwise, we have

rank1(M) = rank2(M) = rank( M(1)) = rank ( M(2^ = Q

On the other hand, if we assume rank3(M) = L and L < K, then M is a rank-(Q,Q,L) tensor. Thus, the Tucker decomposition of tensor M is the so-called Tucker1 decomposition [23]

M=T x11 x21 x3 G

where T e cqxqxl is the core tensor, I e R Q x Q is an identity matrix, and G e C K x L is a column-unitary matrix. This is equivalent to a standard two-dimensional PCA since

M(3) = G x T(3).

where T(3) e C L x QQ is the mode-3 matrix of core tensor T. It is obvious that (18) corresponds to the PCA of M(3). Therefore, G e C K x L consists of the L-leading left singular vectors of M(3). Since K > L and G is a column-unitary matrix, it is straightforward to derive

T(3) = GH x M(3)

Therefore, the core tensor can be obtained by

T = MX11 x21 X3 GH (20)

Because the first and second factors of the Tucker decomposition in (17) are identity matrices, the core tensor T e cqxqxl is also a symmetric tensor [23], as

Figure 1 The Tucker decomposition of a three-way tensor f. For data compression, the size of the core tensor g is smaller than that of f (/1 < J1, /2 < J2, and /3 < J3), and A(1), A(2), and A(3) are thin matrices.

for the tensor M e CQxQxK and the CAND of the core tensor is as follows

T = 1 7°dI;

where al and dl are the lth column of A and D, respectively. The CAND process can also be implemented by the ALS algorithm. Nevertheless, our goal is to estimate the mixing matrix A, whereas the CAND of core tensor T e cQxQxL only conduces A. Hence, it is necessary to derive the mixing matrix A based on A. Since the first and second factors of Tucker decomposition (i7) are identity matrices, it is straightforward to derive

A = IA= A

4.3. Computational analysis

In this subsection, we aim at giving an insight into the numerical complexity of the proposed algorithm. For the GGF-ALS algorithm, the ALS algorithm is directly used to decompose the tensor M e CQxQxK ; therefore, the computational complexity is O(3PKQ2 + QKP2 + Q2P2) per iteration. For the GGF-TALS algorithm, the computational complexity is dominated by the Tucker decomposition of the original tensor M e CQxQxK (realized by the singular-value decomposition (SVD) of the mode-3 matrix M(3)) and the decomposition of the core tensor T e cQxQxL using the ALS algorithm. Therefore, the computational complexities for these two operations per iteration are O(Q6) and O(3PLQ2 + QLP2 + Q2P2), respectively. Table i shows the computational complexity for the GGF-ALS, and GGF-TALS methods.

It is worth mentioning that a large value is usually chosen for the number of processing points K in order to accumulate sufficient statistical information from data, while the rank of the core tensor order L is often chosen to be much smaller than K. For example, K = i00, L = 8 are chosen typically in our simulations as shown in the next section. The SVD is required to be calculated only once, whereas the ALS algorithm needs multiple iterations to achieve convergence, for example, i,000 iterations as in our experiments. Moreover, the number of sources P and the number of sensors Q are usually small, for example, P = 4, Q = 3 in our simulations. For these reasons, it can be readily derived

Table 1 Computational complexity of GGF-ALS, GGF-TALS methods

Method Operation Complexity

GGF-ALS ALS (one iteration) o(3pkq2 + qkp2 + q2p2)

GGF-TALS SVD o(q6)

ALS (one iteration) o(3plq2 + qlp2 + q2p2)

from Table i that the computational cost of the GGF-TALS algorithm is, in practice, much lower than that of the GGF-ALS algorithm (note that the complexity of SVD becomes negligible in this case).

5. Simulations and analysis

In this section, simulations are provided to illustrate the performance of the proposed GGF approaches for underdetermined mixtures of complex sources. The performance of the tested algorithms is evaluated and compared in terms of the relative error performance index (PI) versus the sample size and the signal-to-noise ratio (SNR) of the observations. Here the relative error PI is defined as [7] PI = E{||A-A1|/||A||}, in which the norm is the Frobenius norm and A represents the optimally ordered and scaled estimate of the mixing matrix A.

The experiments refer to the scenario that P = 4 narrowband source signals are received by a uniform circular array (UCA) with Q = 3 identical sensors of radius Ra Considering a free space propagation model, the entries of the mixing matrix A are given by

aqp = exp(2 jaq cos(Op) cos(+ cos(Op) sin(

where aq = (RJA)cos(2n(q - i)/Q), flq = (Ra/A)sin(2n(q - i)/Q), and j = -\/-i. We have Ra/A = 0.55. The direction of arrival (DOA) of the different sources are given by 0i = 3n/i0, O2 = 3n/i0, O3 = 2n/5, O4 = 0 and & = 7n/i0, 02 = 9n/i0, = 3n/5, 04 = 4n/5. The sources are unit-variance 4-QAM with a uniform distribution, shaped by a raised cosine pulse shaping filter with a roll-off p = 0.3. All sources have the same symbol duration T = 4Te, where Te is the sample period. The observations are contaminated by additive zero-mean complex Gaussian noise.

First, we compare the performance of the GGF-ALS algorithm with that of the GGF-TALS algorithm and investigate the influence of the rank of the core tensor L on the performance of the GGF-TALS algorithm. To this end, the following two simulation experiments are conducted. In the first simulation, we evaluate the performance of the GGF-TALS algorithm with different core tensor rank L for a fixed number of processing points K, and compare it with GGF-ALS with the same number of processing points. In this simulation, K is chosen to be i00 for the GGF-ALS algorithm, both the real and imaginary parts of processing points are randomly drawn from [-i; i], the SNR of the observations ranges from 0 to 25 dB, and the number of samples is 4,000. The core tensor rank for the GGF-TALS algorithm is chosen as L = i0,8,5, respectively. The threshold value described in (i2) to stop the ALS algorithm is i0-5, and i00 Monte Carlo experiments are run.

The average performance of the GGF-TALS algorithm versus the core tensor rank for a fixed number of processing points is shown in Figure 2. We can see that the

Figure 2 Performance of GGF-TALS versus core tensor order for a fixed number of processing points. The results for GGF-ALS are also shown for comparison.

performance curves of the GGF-TALS and that of the GGF-ALS almost coincide when the rank of the core tensor is larger than 8, whereas the performance of the GGF-TALS algorithm slightly deteriorates when the rank of the core tensor is less than 8. This indicates that the GGF-TALS algorithm maintains the identification accuracy offered by GGF-ALS despite the fact that the core tensors used by GGF-TALS have a much lower rank than that of the original tensor used by GGF-ALS.

In the second experiment, we investigate the performance of the GGF-TALS algorithm with different

number of processing points K when the rank of the core tensor L is fixed. The simulation conditions are the same as those in the previous simulation except the number of processing points for the GGF-ALS algorithm which are 20, 40, and 100 respectively, and the rank of the core tensor for the GGF-TALS which is 8. The average performance of the GGF-TALS versus the number of processing points for the fixed number of the core tensor rank is shown in Figure 3.

We can see that the average performance of GGF-TALS is consistent with that of GGF-ALS when both algorithms use the same number of processing points. For instance, when the number of processing points is 100, the performance of GGF-TALS is close to the performance of GGF-ALS, so is for K = 40. Therefore, the GGF-TALS algorithm consistently offers similar performance to GGF-ALS when the number of processing points is varied.

Second, we compare the performance of the GGF-ALS and GGF-TALS algorithms with that of the latest GF-based BI algorithm. Here, the LEMACFAC-2 in [18] is chosen as the baseline algorithm. The number of processing points and the value of processing points for the GGF-ALS, GGF-TALS, and LEMACFAC-2 are the same. Specifically, the number of processing points is 100, both the real and imaginary parts of processing points are randomly drawn from [_1; 1], and the rank of the core tensor is 8. The threshold value described in (12) to stop the LM algorithm [26] exploited in LEMACFAC-2 method is also 10_5.

Figure 4 shows the PI of the tested algorithms as a function of the SNR when 4,000 samples are used. It can be seen from this figure that the performance of the GGF-ALS is consistent with that of the GGF-TALS, and both perform better than LEMACFAC-2. Figure 5 shows

—e— GGF-ALS(K=20) —B— GGF-ALS(K=40) 0 GGF-ALS(K=100)

- 0 - GGF-TALS(K=20,L=8) -B - GGF-TALS(K=40,L=8)

- » - GGF-TALS(K=100,L=8)

SNR{dB)

Figure 3 The performance of GGF-TALS versus the number of processing points. The results for GGF-ALS are also shown for comparison.

- O - LEMACFAC-2 in [18]

The proposed GGF-TALS(K=1 00,L=8) —0— The proposed GGF-ALS(K=100)

ftftOOOOOOOtfOOOOO»

SNR(dB)

Figure 4 The performance of the tested algorithms versus the SNR of the observations.

■ e - LEMACFAC-2 in [26]

The proposed GGF-TALS(K=100,L=8) — The proposed GGRALS(K=100)

-O-----

-Ba__~ ^ ~ 1

-----"<>

3000 4000 5000 6000 7000 800 samples

Figure 5 The performance of the tested algorithm versus the size of the data set.

the PI of the tested algorithms as a function of the number of data samples N when the SNR is equal to 20 dB. Again, the LEMACFAC-2 underperforms the GGF-ALS and the GGF-TALS algorithms, and the GGF-ALS and the GGF-TALS have almost the same performance. This confirms that the GGF algorithms perform better than the method based on the GF defined in (3) in exploiting the statistical information carried on the complex variables when using the same number of processing points.

Appendix 1

In this appendix, we show the computational details of the core equation in Equation 7. First, the differentiation of (6) with respect to u gives

9^z(u) y^P d<Psp (a^u)

-/p=1 3u*

3ft, ( a» 3( aHu

v^P r s \ P ! \ P ! V^P

^ 3u* =^

3( a]>

Second, the differentiation of (23) with respect to uT gives

dft (u) _ y^P tfll) ^ d(^ T

9uT9u* ^ p=1.

Z—/p=1 Z^ p=1

m ajfu)^aJfuV p duT p dVs( a£u) 9( a«u

9( a£u) a(aHu

d<PsA a>

9( a^u] d( aHu

H apap

In a more compact form, we can obtain the core equation

Vz (u) = u)AH

6. Conclusions

We have presented two algorithms, namely, GGF-ALS and GGF-TALS, for blind identification from underdetermined mixtures of complex sources using the second GGF of the observations. In the GGF-ALS algorithm, the mixing matrix is estimated by directly using the ALS algorithm to decompose the tensor constructed from the Hessian matrices of the second GGF of the observations. The GGF-TALS algorithm is an improved version of the GGF-ALS algorithm, where the Tucker decomposition is first used to convert the original tensor into a lower-rank core tensor, and the mixing matrix is then estimated by applying the ALS algorithm to the core tensor. Simulation results have shown that (a) the proposed GGF-ALS and GGF-TALS approaches have almost the same performances in terms of the relative errors, whereas the GGF-TALS has a much lower computational complexity, and (b) the proposed GGF algorithms have superior performance to the latest GF-based BI approaches, since the GGF algorithms can exploit the statistical information carried on complex variables in a more effective way.

Appendix 2

Let <s(u) denote the characteristic function of the source signals s(t). Due to the statistical independence of the elements of s(t) = [si(t), •••, sP(t)]T, we have

<s(u) = 0si (~i) x 0s2 (~2) x • x <s„ (~p) (26)

where (uj) = exp(|wjjj , j = i, •••, P. Therefore, for the second characteristic function ps(u), we have

fts(u) = <Ps, (~1) + Vs2 (~2) + + <PSp (uP)

Consequently, the Hessian matrix %s(u) can be easily obtained as

Vs( u) = VaT

= diag

9 dûF

d(PSl (~1) d(PS2 (U2) • d<ps„ (UP) ÔK^Sj ' d~2du2 ' ' dUpdup

where diag(-) represents an operator forming a diagonal matrix by assigning its arguments to the entries of the main diagonal.

Competing interests

The authors declare that they have no competing interests. Acknowledgments

This work is supported in part by NaturalScience Foundation of China under grant 61001106 and NationalProgram on Key Basic Research Project of China under grant 2009CB320400, and Engineering and PhysicalScience Research Councilof the UK under grant EP/H050000/1. The authors would like to thank the anonymous referees for their constructive comments for improving this paper and Dr. Mark Barnard for proofreading the manuscript.

Author details

1College of Communication Engineering, PLA University of Science & Technology, Nanjing 210007, People's Republic of China. 2Department of Electronic Engineering, University of Surrey, Guildford GU2 7XH, UK.

Received: 2 November 2012 Accepted: 21 June 2013 Published: 1 July 2013

References

1. B Chen, AP Petropulu, Frequency domain blind MIMO system identification based on second- and higher order statistics. IEEE Trans. Signal Process. 49(8), 1677-1688 (2001)

2. K Rahbar, JP Reilly, A frequency domain method for blind source separation of convolutive audio mixtures. IEEE. Trans. Speech Audio Process

13(5), 832-844 (2005)

3. A Aissa-EI-Bey, N Linh-Trung, K Abed-Meraim, A Belouchrani, Y Grenier, Underdetermined blind separation of nondisjoint sources in the time-frequency domain. IEEE Trans. Signal Process. 55(3), 897-907 (2007)

4. W Dai, T Xu, W Wang, Simultaneous codeword optimization (SimCO) for dictionary update and learning. IEEE Trans. Signal Process.

60(12), 6340-6353 (2012)

5. P Bofill, M Zibulevsky, Underdetermined blind source separation using sparse representations. Signal Process. 81, 2353-2362 (2001)

6. C Fevotte, SJ Godsill, A Bayesian approach for blind separation of sparse sources. IEEE Trans Audio Speech Lang Process 14(6), 2174-2188 (2006)

7. LD Lathauwer, J Castaing, Blind identification of underdetermined mixtures by simultaneous matrix diagonalization. IEEE Trans. Signal Process.

56(3), 1096-1105 (2008)

8. P Tichavsky, Z Koldovsky, Weight adjusted tensor method for blind separation of underdetermined mixtures of nonstationary sources. IEEE Trans. Signal Process. 59(3), 1037-1047 (2011)

9. A Ferreol, L Albera, P Chevalier, Fourth-order blind identification of underdetermined mixtures of sources (FOBIUM). IEEE Trans. Signal Process. 53(5), 1640-1653 (2005)

10. LD Lathauwer, J Castaing, JF Cardoso, Fourth-order cumulant-based blind identification of underdetermined mixtures. IEEE Trans. Signal Process. 55(6), 2965-2973 (2007)

11. L Albera, A Ferreol, P Comon, P Chevalier, Sixth order blind identification of under-determined mixtures (B/RTH) of sources. Proceedings of the 4th International Symposium on Independent Component Analysis and Blind Signal Separation (Nara, 2003), pp. 909-914

12. ALF de Almeida, X Luciani, P Comon, Blind identification of underdetermined mixtures based on the hexacovariance and higher-order cyclostationarity. IEEE Workshop on Statistical Signal Processing (Cardiff, UK, 2009), pp. 669-672

13. L Albera, A Ferreol, P Comon, P Chevalier, Blind identification of over-complete mixtures of sources (BIOME). Lin. Algebra Appl. 391,1-30 (2004)

14. A Karfoul, L Albera, G Birot, Blind underdetermined mixture identification by joint canonical decomposition of HO cumulants. IEEE Trans. Signal Process. 58(2), 638-649 (2010)

15. A Yeredor, Blind source separation via the second characteristic function. Signal. Process. 80(5), 897-902 (2000)

16. E Eidinger, A Yeredor, Blind MIMO identification using the second characteristic function. IEEE Trans. Signal Process. 53(11), 4067-4079 (2005)

17. P Comon, M Rajih, Blind identification of underdetermined mixtures based on the characteristic function. Signal. Process. 86(9), 2271-2281 (2006)

18. X Luciani, ALF de Almeida, P Comon, Blind identification of underdetermined mixtures based on the characteristic function: the complex case. IEEE Trans. Signal Process. 59(2), 540-553 (2011)

19. A Yeredor, Non-orthogonal joint diagonalization in the least-squares sense with application in blind source separation. IEEE Trans. Signal Process. 50(7), 1545-1553 (2002)

20. LD Lathauwer, BD Moor, J Vandewalle, On the best rank-1 and rank-(r,, R2,..., rn) approximation of higher-order tensors. SIAM J. Matrix Anal. Appl. 21, 1324-1342 (2000)

21. D Nion, LD Lathauwer, An enhanced line search scheme for complex-valued tensor decompositions: application in DS-CDMA. Signal Process.

88(3), 749-755 (2008)

22. Y Chen, D Han, L Qi, New ALS methods with extrapolating search directions and optimal step size for complex-valued tensor decompositions. IEEE Trans. Signal Process. 59(12), 5888-5898 (2011)

23. TG Kolder, BW Bader, Tensor decompositions and applications. SIAM Rev. 51 (3), 455-500 (2009)

24. G Bergqvist, EG Larsson, The higher-order singular value decomposition: theory and an application. IEEE Signal Process. Mag. 27(3), 151-154 (2010)

25. G Zhou, A Cichocki, Fast and unique tucker decompositions via multiway blind source separation. Bull. Pol. Acad. Sci. 60(3), 389-405 (2012)

26. P Comon, X Luciani, ALF de Almeida, Tensor decomposition, alternating least squares and other tales. J. Chemometr. 23, 393-405 (2009)

doi:10.1186/1687-6180-2013-124

Cite this article as: Gu et al.: Generalized generating funetion with tueker deeomposition and alternating least squares for underdetermined blind identifieation. EURASIP Journal on Advances In Signal Processing 2013 2013:124.

Submit your manuscript to a SpringerOpen journal and benefit from:

7 Convenient online submission 7 Rigorous peer review 7 Immediate publication on acceptance 7 Open access: articles freely available online 7 High visibility within the field 7 Retaining the copyright to your article

Submit your next manuscript at 7 springeropen.com