SAGE-Based Algorithm for Direction-of-Arrival Estimation and Array Calibration

Most existing array processing algorithms are very sensitive to model uncertainties caused by the mutual coupling and sensor location error. Tomitigate this problem, a novel method for direction-of-arrival (DOA) estimation and array calibration in the case of deterministic signals with unknown waveforms is presented in this paper. The analysis begins with a comprehensive perturbed array output model, and it is effective for various kinds of perturbations, such as mutual coupling and sensor location error. Based on this model, the Space Alternating Generalized Expectation-Maximization (SAGE) algorithm is applied to jointly estimate the DOA and array perturbation parameters, which simplifies the multidimensional search procedure required for finding maximum likelihood (ML) estimates. The proposed method inherits the characteristics of good convergence and high estimation precision of the SAGE algorithm. At the same time, it forms a unified framework for DOA and array perturbation parameters estimation in the presence of mutual coupling and sensor location error. The simulation results demonstrate the effectiveness of the algorithm.


Introduction
The problem of estimating direction-of-arrival (DOA) of multiple narrowband signals plays an important role in many fields, including radar, wireless communications, seismology, and sonar.The performance of most existing DOA estimation methods relies crucially on perfect knowledge of the array manifold.In practice, however, the array manifold is often affected by the mutual coupling, gain/phase uncertainty, sensor location error, and so forth.Without array manifold calibration, the performance of DOA estimation may deteriorate significantly.
To mitigate this problem, various array calibration methods have been proposed.Those methods formulate the array perturbations with unknown parameters and estimate those parameters together with the source directions, so as to realize array calibration and DOA estimation.Most of the existing methods only focus on certain type of array imperfection, such as mutual coupling [1][2][3][4], gain/phase uncertainty [5][6][7], and sensor location error [8,9].The DOA estimation problem in the presence of more than one type of array perturbation has also been studied in a small amount of literatures [10,11].However, the method in [10] is computationally much expensive as it requires multidimensional parameter searching, while the method in [11] needs to know the directions of the calibration sources.
Among all popular DOA estimation methods, the maximum likelihood (ML) approach provides the best asymptotic performance and remains stable in scenarios involving small numbers of snapshots, coherent signals, and low signal-tonoise ratios (SNR).The main drawback of the ML approach is the high computational complexity caused by optimization of the likelihood function.To reduce this difficulty, the Expectation-Maximization (EM) algorithm has been derived for both deterministic and stochastic signal models with known noise covariance structure [12,13].The Space Alternating Generalized EM (SAGE) algorithm is a variation of the widely used EM algorithm, which updates subsets of parameters sequentially in one iteration.Through the augmentation scheme specified by the EM or SAGE algorithm, the complicated multidimensional search involved in maximizing likelihood functions can be simplified to one-dimensional search.It was proved in [14] that SAGE converges faster than EM while retaining the advantages of numerical simplicity and stability due to its flexible augmentation scheme.Both the EM algorithm and the SAGE algorithm are earlier applied 2 International Journal of Antennas and Propagation to DOA estimation problem in accurately calibrated arrays.However, because the directional information can hardly be extracted directly when array perturbation exists, there is little literature using EM or SAGE algorithm to estimate DOA in perturbed arrays.In [15], we present a comprehensive perturbed output model, which makes the relationship between the array output and the perturbation parameters clearer.The sparse Bayesian method and EM algorithm are applied to estimate DOA and perturbation parameters in the case of stochastic signals.
In this paper, we focus on deterministic and unknown narrowband signals.And a single type of array imperfection is considered, with mutual coupling and sensor location error treated as typical examples.Based on the comprehensive perturbed array output model proposed in [15], we derive an SAGE-based algorithm, which applies the SAGE algorithm to jointly estimate DOA and perturbation parameters.Compared to the existing method, the proposed algorithm can achieve higher estimation precision; at the same time, it follows a unified framework to address the DOA estimation problem in the presence of array imperfections, with typical perturbations of mutual coupling and sensor location error taken into consideration.
This paper is outlined as follows.The comprehensive perturbed array output model is described in Section 2. The SAGE algorithm for DOA estimation and array calibration is developed in Section 3. Numerical results are presented and discussed in Section 4. Section 5 concludes this work.

Perturbed Array Output Formulation
Consider that  unknown and deterministic narrowband signals impinge onto an  element uniform linear array from directions of  = [ 1 , . . .,   ]  , where  < .We take the first array sensor as the reference.The perturbed array output is given as follows: where the array responding matrix A  () = [a  ( 1 ), . .., a  (  )] and a  () is the perturbed responding vector.The observe vector x() = [ 1 (), . . .,   ()]  is sampled at time instances  =  1 , . . .,   , where  represents the number of snapshots.The signal vector s() = [ 1 (), . . .,   ()]  , and   () is the waveform of the th signal.The noise vector k() = [V 1 (), . . ., V  ()]  is independent, identically complex, and normally distributed with zero mean and covariance matrix  2 I, where  is an unknown noise spectral parameter.The two types of typical array perturbations, including mutual coupling and sensor location error, are concerned in this paper, and the array geometry is assumed to be linear to simplify notations.The sensor location error is also assumed to exist along the array axes and thus does not destroy the linear geometry of the array.The perturbed array responding vector has diverse expressions in the case of different array imperfections.When mutual coupling is present, the vector is a  () = Ca() [2], with C representing the mutual coupling matrix.The mutual coupling matrix can be written more explicitly as C = toeplitz([1,  1 , . . .,   , 0  (−−1)×1 ]  ), where   ∈ C is the coupling coefficient of the two sensors displaced by  − 1 times the interelement spacing of the ULA, and   is very small when  >  and is, thus, neglected.The perturbed array responding vector in the presence of sensor location error is a  () = [ 2( 1 + d1 ) sin / , . . .,  2(  + d ) sin / ]  , with d1 , . . ., d being the location errors of the  sensors [8].In this paper, we take the first array sensor as the reference; thus  1 = 0 and d1 = 0.
Although the directional information of the incident signals is still reserved in the perturbed array outputs, it can hardly be extracted directly based on the array output formulation in (1), as the structure of A  () is not available beforehand due to the unknown array imperfections.In order to highlight the perturbation-free signal components and also make the relationship between the array output and the perturbation parameters clearer, we have established the following comprehensive formulation of the perturbed array output in [15]: Equation ( 2) adapts to the typical array perturbations, including mutual coupling and sensor location error, with c standing for a column vector consisting of the array perturbation parameters, c ∈ C ×1 for mutual coupling, and c ∈ R (−1)×1 for sensor location error.We use  to denote the dimension of c uniformly for notational convenience; that is, c = [ 1 , . . .,   ]  , and  =  − 1 in the case of sensor location error.The explicit formulation of A  () varies with respect to the array perturbation type, and the following equation holds for Q() and c: where Ψ is a function of c and independent of the signal directions, while Φ() is independent of c and relies on the signal directions; they diverge largely according to the array perturbation type.However, for each of the typical perturbations, it holds for is concluded by taking the differentiations of both sides of (3) with respect to   .In order to explain (2) and ( 3) more explicitly, we itemize the expressions of c, Ψ, Φ(), and Q().
Mutual Coupling.Equation ( 2) can be rewritten as follows in the case of mutual coupling: By combining (3) and ( 4), it can be concluded that , and G  = Ψ/  = C/  contain nonzero elements of 1 only on the ± diagonals.
Sensor Location Error.Equation ( 2) can be rewritten as follows in the case of sensor location error: The first equation in ( 5) is concluded via first-order Taylor expansion under the assumption of small sensor location errors, and Ψ/  have their ( + 1,  + 1)th element being the only nonzero value of 1.

SAGE-Based Algorithm for DOA Estimation and Array Calibration
It is well known that both EM algorithm and SAGE algorithm have been applied to the question of DOA estimation without array imperfections.And they can achieve the global convergence.The EM algorithm is a general numerical method for finding maximum likelihood estimates which is characterized by simple implementation and stable convergence.The SAGE algorithm is a generalized form of the EM algorithm, which allows a more flexible optimization scheme and converges faster than the EM algorithm.Instead of estimating all parameters at once in EM, SAGE breaks up the problem into several smaller ones and uses EM to update the parameter subset associated with each reduced problem.In this section, based on the existing comprehensive perturbed array output model, we propose an SAGE-based algorithm to jointly estimate DOA and array perturbation parameters.At the same time, for comparing the performance of the SAGE and EM algorithms in array self-calibration problem, the EM-based algorithm is also considered.
Based on the comprehensive model of the perturbed array output in (2), we construct the augmented data as follows for separating different signal components of array output: where . G , and (  ) diverge largely according to the array perturbation type.
Both EM and SAGE algorithms update each unknown parameter by iterative method.Each iteration consists of an E-step and an M-step.Given the estimate of the ( − 1)th iteration Ω (−1) , the th iteration consists of the following steps.

E-Step. Calculate
where the superscript () stands for the iteration index and E[•] denotes expectation operator.Equation ( 11) is equivalent to computing the following conditional expectations: For emphasizing the dependence on array perturbation parameters, we use A  (c (−1) ,  (−1) ) to denote the perturbed responding matrix in (15), which is different in above context.

M-
Step.The unknown parameters Ω are updated by maximizing F(Ω () ; Ω (−1) ), which can be realized by setting its differentiations of those parameters to zero.Consider 1) , )     a  (c (−1) , ) As can be seen from ( 15) to (18), the update processes of   and {  (  )}  =1 are mainly dependent on the corresponding augmented data {y  (  )}  =1 , and it does not matter with {y 1 (  ), . . ., y −1 (  ), y +1 (  ), . . ., y  (  )}  =1 .We know that the parameters c and  2 exist in all signal components model, so the update processes are rested on {y 1 (  ), . . ., y  (  )}  =1 .The iteration process can be divided into two parts: sequential update   and {s  (  )}  =1 for different  values and unified update c and  2 by using all of augmented data.The formulation of Q  () is very complicated, the iteration strategy of those parameters given above is not usable directly.In appendix, we carry out some deeper analysis to simplify the implementation of the iteration.

Simulation Results
The simulations in this section illustrate the performance of our proposed algorithm.The proposed algorithm gives a unified framework for DOA estimation in the presence of mutual coupling and sensor location error.In the first experiment, the scenario when mutual coupling exists is considered.For comparison, we also apply S-S method [2] and Y-L method [3] to the same batch of data, as well as the curves of the Cramer-Rao lower bound (CRLB) [16].In the second experiment, the scenario when sensor location error exists is considered.The maximum likelihood method [8] (denoted by ML method), IMTAM [8], and CRLB are also implemented for the performance comparison.
We consider a nominally uniform linear array of 8 sensors with interelement spacing equaling half a wavelength of the where  denotes the number of Monte Carlo experiments.
We set  = 500 in the simulation.θ() and  () are the estimated and true directions in the th simulation.Similarly, the RMSE of the mutual coupling coefficients is used for the array calibration precision evaluation, which is defined as where c keeps constant in each scenario, ĉ() is the estimated coupling coefficient vector in the th simulation,  is a tuning factor introduced to enhance the sense of the RMSE,  = ‖‖ 2 for mutual coupling, and  =  −1 for sensor location error.

Statistical Performance in the Presence of Mutual Coupling.
For each Monte Carlo experiment, we suppose that the coupling coefficient vector includes two coupling coefficients.Each coefficient is a random complex data.The range of real and imaginary parts of coefficients is 0.1∼0.6.Figure 1 shows the RMSE of DOA estimates of different methods against  input SNR.The results illustrate that both SAGE-based and EM-based algorithms can lead to significant improvement in estimation accuracy, especially for low levels of SNR, and they improve with a similar speed as that of the CRLB when the SNR increases.Figure 2 shows the RMSE of mutual coupling coefficients versus input SNR.The results illustrate that both S-S method and Y-L method cannot achieve high estimation precision, while the SAGE-based and EM-based algorithms can achieve the better performance.They are always close to the CRLB. Figure 3 shows the average number of iterations needed for convergence versus input SNR.Since each iteration of the SAGE-based and EM-based algorithms requires similar computations, their total computational costs are determined by the number of iterations.It can be observed that the number of iterations needed by the SAGE-based algorithm is always smaller than that by the EM-based algorithm over the whole range of SNR.

Statistical Performance in the Presence of Sensor Location Error.
In the second experiment, we suppose that the sensor locations are not accurately calibrated.They depart from the nominal positions by 0, 0.08, 0.12, 0.16, 0.04, −0.12, −0.08, and −0.16 times the half-wavelength; that is, c = [0.08,0.12, 0.16, 0.04, −0.12, −0.08, −0.16]  .The nominal sensor positions are used to initialize the algorithm; that is, c [0] = [0, 0, 0, 0, 0, 0, 0]  .The initial DOA estimate is given by  Figure 4 shows RMSE of DOA estimates of different methods against SNR.As can be shown in Figure 4, the SAGE-based and EM-based algorithms perform similarly over the entire SNR range.They are always close to CRLB.Both of them have significant improvement in estimation accuracy than ML method and IMTAM method.Figure 5 shows RMSE of sensor location error of different methods against SNR.We can easily observe that ML method fails to estimate sensor location error effectively (no array calibration result is available from IMTAM), while the SAGEbased and EM-based algorithms provide reasonable results.Figure 6 shows the average number of iterations needed for convergence versus each SNR.As can be shown in Figure 6, as the SNR increases, the numbers of iterations required by both algorithms are reduced quickly, and the SAGE-based algorithm gets a faster convergence rate over the entire SNR range.

Conclusion
In this paper, the DOA and array perturbation parameters estimation problem in the case of deterministic signals with unknown waveforms is studied.An SAGE-based method for DOA estimation and array calibration in the presence of mutual coupling and sensor location perturbation is proposed.First, a comprehensive array output model that is applicable to two typical array perturbations is introduced.Based on this model, we construct the augmented data needed by the SAGE algorithm and establish the augmented data unified likelihood for  snapshots.Then, the E-step and M-step of the SAGE algorithm for DOA and perturbation parameters estimation are derived.The implementation of the proposed method follows the standard SAGE iterations and thus has guaranteed convergence and high estimation precision.For comparison, the EM-based algorithm is also studied, which applies the EM algorithm to jointly estimate DOA and perturbation parameters.The simulation results show that the EM-based algorithm can obtain higher It can be concluded that ( ) = a(  )(2/) sin   , c = [ 1 , . . .,   ]  = [ d2 , . . ., d ] , and G , have their ( + 1,  + 1)th element being the only nonzero value of 1. k  () is a portion of noise extracted from k(), and k  () ∼ N(0,    2 I  ), 0 ≤   ≤ 1.

2 Figure 1 :
Figure 1: RMSE of DOA estimates of different methods against SNR.incident signals.Suppose that two equal-power independent signals impinge onto the array from directions of −9 ∘ and 25 ∘ .The snapshot number is fixed at 200.The algorithm is terminated if the increase in the data likelihood function is smaller than 10 −4 or the maximal number of iterations has arrived.The maximal number of iterations in the simulation is set to be 200.The average root-mean-square error (RMSE) of the DOA estimates is considered for statistical direction estimation precision evaluation, which is defined as

Figure 2 : 2 Figure 3 :
Figure 2: RMSE of coupling coefficients of different methods against SNR.

2 Figure 6 :
Figure 6: The average number of iterations needed for convergence against SNR.