Root Sparse Bayesian Learning for Off-Grid DOA Estimation

The performance of the existing sparse Bayesian learning (SBL) methods for off-gird DOA estimation is dependent on the trade off between the accuracy and the computational workload. To speed up the off-grid SBL method while remain a reasonable accuracy, this letter describes a computationally efficient root SBL method for off-grid DOA estimation, where a coarse refinable grid, whose sampled locations are viewed as the adjustable parameters, is adopted. We utilize an expectation-maximization (EM) algorithm to iteratively refine this coarse grid, and illustrate that each updated grid point can be simply achieved by the root of a certain polynomial. Simulation results demonstrate that the computational complexity is significantly reduced and the modeling error can be almost eliminated.


I. INTRODUCTION
Direction-of-arrival (DOA) estimation of narrow-band sources is an important topic in array signal processing and has attracted tremendous interest in many fields, such as radar, sonar, mobile communications [1]. The conventional subspacetype DOA estimation methods, including MUSIC and ESPRIT, require a large number of snapshots to achieve their highresolution performance, and they may fail to work when signals are highly correlated or coherent due to multipath propagation. Recently, the emerging technique of sparse representation has given renewed interest to the DOA estimation problem [2]- [5]. The sparse representation methods exhibit many advantages compared to the conventional subspace-type ones, e.g., improved robustness to noise, limited number of snapshots and correlation of signals.
Sparse Bayesian learning (SBL) [6], [7] is one of the most popular approaches for the sparse signal recovery, which exploits the sparsity information with a sparse prior assumption for the signal of interest from a Bayesian perspective. Theoretical results show that SBL includes the l 1 -norm minimization method as a special case when a maximum a posteriori (MAP) optimal estimate is adopted with a Laplace signal prior [8], [9]. The excellent performance achieved by the SBL methods for DOA estimation relies crucially on the assumption that the true DOAs lie on (practically, close to) the sampling grid points. In practice, however, this assumption is usually unavailable. The gap between the true DOA and its nearest grid point is known as the off-grid gap.
A number of improved methods have been proposed to deal with the off-grid DOA estimation [10]- [12]. Yang et al. [11] applied a linear approximation to the true DOA and proposed an off-grid SBL method, where the off-grid gap is assumed to be uniformly distributed (noninformative). Using the sample covariance matrix, Zhang et al. [12] further gave an improved off-grid SBL method to reduce the effect of noise variance. The existing methods [11], [12] remain a major problem that their performance is dependent on the trade off between the accuracy and the computational workload. When a coarse grid is used, it may lead to a high modeling error; while, if a dense sampling grid is used, the massive involved computation could make this method unfavorable for real applications.
In this letter, we try to propose a computationally efficient root SBL method for off-grid DOA estimation, which adopts a coarse grid and considers the sampled locations in the coarse grid as the adjustable parameters. We utilize an expectationmaximization (EM) algorithm to iteratively refine this coarse grid, and illustrate that each updated grid point can be simply achieved by the root of a certain polynomial. In this way, the computational complexity is significantly reduced and the modeling error can be almost eliminated.

II. DATA MODEL
Consider K narrow-band far-field sources impinging on an M -element uniform linear array (ULA), where the distance between adjacent sensors is d. The K sources, s 1 (t), s 2 (t), . . . , s K (t), arrive at the array from distinct directions, θ 1 , θ 2 , . . . , θ K , with respect to the normal line of the array. The M × 1 array output vector y(t) is then given by where λ is the wavelength of the source, and n(t) = [n 1 (t), n 2 (t), . . . , n M (t)] T is an unknown noise vector.
With the definitions of Y [y(t 1 ), y(t 2 ), . . . , y(t T )], N [n(t 1 ), n(t 2 ), . . . , n(t T )] and S [s(t 1 ), s(t 2 ), . . . , s(t T )], we have In order to cast the problem of DOA estimation as a sparse representation problem, generally we let {θ i }K i=1 be a fixed sampling grid that uniformly covers the DOA range [−π/2, π/2], whereK denotes the grid number. If the grid is fine enough such that the true DOAs lie on (or, practically, close to) the grid, we can use the following model for Y: where Aθ [α(θ 1 ), . . . , α(θK)] andŜ is aK × T complex matrix whose the ith row corresponds to the signal impinging on the array from a possible source atθ i . However, the assumption that the true DOAs are located on the predefined spatial grid is not always valid in practical implementations.
To solve this problem, a linear approximation method is used in [11].
. . ,K}, is the nearest grid point to θ k , the steering vector α(θ k ) is approximated by the linearization: where b(θ n k ) = α ′ (θ n k ). By absorbing the approximation error into the noise, the observation model (3) can be rewritten as diag{·} represents the diagonal matrix operator, and ζ is a zero vector except that the n k -th element Obviously, the modeling error caused by off-grid gap can be alleviated by the new model (5), but it can not be fully eliminated. When a coarse grid is used, (5) may still lead to a high modeling error; while, if a dense sampling grid is used, the massive involved computation could make the linear approximation method unfavorable for real applications.

III. THE PROPOSED ROOT SBL METHOD
In this section, we try to propose a computationally efficient root SBL method for off-grid DOA estimation. To keep a low computational load, we adopt the original model (3) with a coarse grid (i.e.,K is small). To handle the modeling error in (3), we consider the sampled locations in the coarse grid as the adjustable parameters, and utilize an expectation-maximization (EM) algorithm to iteratively refine the grid. Note that the name is due to the fact that each updated grid point can be simply achieved by using the root of a polynomial.

A. Sparse Bayesian Formulation
Firstly, we address the sparse Bayesian model that is commonly used in SBL [6]. A typical SBL treatment ofŜ begins by assigning a non-stationary Gaussian prior distribution with a distinct inverse variance δ i for each row ofŜ. Letting δ = [δ 1 , δ 2 , . . . , δK] T and ∆ = diag(δ), we have whereŝ t denotes the tth column ofŜ. In order to obtain a twostage hierarchial prior that favors most rows ofŜ being zeros, the hyper-parameter δ i s are further modeled as independent Gamma distributions [6], [8], i.e., where ρ is a small positive constraint (e.g., ρ = 0.01 [8]). Under an assumption of circular symmetric complex Gaussian noises, we have where β σ −2 denotes the noise precision with σ 2 being the noise variance,θ [θ 1 ,θ 2 . . . ,θK] and y t denotes the tth column of Y. As β is usually unknown, we model it as a Gamma hyperprior p(β) = Γ(β; a, b), where we set a, b → 0 as in [6], [8] so as to obtain a broad hyperprior.

B. Bayesian Inference
As p(Ŝ, β, δ|Y;θ) cannot be explicitly calculated, an EM algorithm will be exploited to perform the Bayesian inference. The principle behind the EM algorithm is to repeatedly construct a lower-bound on the evidence function p(β, δ|Y;θ), or equivalently, ln p(Y, β, δ;θ) (E-step), and then optimize that lower-bound (M-step).
In the rest part, we will focus on the parameter update for θ. Ignoring terms in the logarithm independent thereof, we just have to maximize ln p(Y|Ŝ, β;θ) p(Ŝ|Y,β,δ;θ) To refine each sampled locationθ i , i = 1, 2, . . . ,K, or equivalently, its exponential form vθ i ( e −j2πd/λ sin(θi) ), we calculate the derivative of (17), with respect to vθ i , and then set it to zero: where a i , µ ti and γ ij denote the ith column, the ith element and the (i, j)th element of Aθ, µ t and Σ, respectively, y t−i y t − j =i µ tj a j , a ′ i da i /dvθ i , and (·) * stands for conjugate operation. Here, we use the fact that After some algebraic operations, (18) can be rewritten in a polynomial form: As the polynomial is of order M − 1, it has M − 1 roots in the complex plane. According to the definition of vθ i , the selected root for refining vθ i should be with a unit absolute value; however, due to the presence of noise, the roots may not be on the unit circle. In this case, the closest root to the unit circle is selected (which is denoted by z i ⋆ ), and the candidate point for the refined grid iŝ On the other hand, we note that when the original coarse grid uniformly covers the DOA range, the DOA estimates on the rough grid is near the true DOAs. Hence, we further screen the candidate from where it falls into.θ new ; otherwise, it is rejected and the corresponding grid point retains unchanged.   Actually, we do not need to refine everyθ i in each iteration, because anyθ i s corresponding to rows that have small entries can be safely discarded. In the practical implementation, we may set a threshold to select some proper active grid points as follows. Let f t be the Frobenius norm of the mean ofŝ t . Then, the indexes of the grid points that need to be activated can be selected by finding the first η maxima f t , where 1 ≤ η ≤ M . Note that if we use a small η, the selected grid points in the current iteration may miss some true DOAs. However, as the Frobenius norms f t s will vary with the iteration, the missing DOAs might be activated in the next iterations. Simulation results recommend to set η ≥ K for grid-refining, especially η = M if the number of sources K is not available.
Since η is much smaller than the grid numberK, the computation of (20) is negligible, and the most demanding steps for our method are in (14) and (15), whose computational complexity is of order O(MK 2 ) per iteration. It is the same as the OGSBI method [11]. However, it is worth noting that our method can work with a very coarse grid, which will bring a significant improvement for computational cost.

IV. SIMULATION RESULTS
In this section, we will present several simulation results to illustrate the performance of our proposed method. We will compare the proposed method to OGSBI [11] and the  Fig. 1 shows the RMSE of DOA estimation versus grid interval based on 200 Monte Carlo runs. It is seen that our method outperforms the state-of-the-art methods, especially when the grid interval is large. The reason is that the linear approximation adopted in OGSBI will lead to a high modeling error when a coarse grid is used [see (5)]; while our method can properly handle the modeling error through viewing the sampled locations in the coarse grid as the adjustable parameters. Fig. 2 shows the total CPU time versus grid interval based on 200 Monte Carlo runs. As can be seen from the figure, the computational times required by all the methods decrease as the grid gets coarser. Our method is much faster than the other methods, especially when the grid interval is large. Therefore, the results in Figs. 1 and 2 recommend to use a coarser grid with r = 4 • or r = 6 • , as it can give a fast DOA estimation but remain a reasonable accuracy.
The last simulation investigates the effect of active grid points on the DOA estimation performance. Consider the same scenario as in Simulation 1, except for SNR = 0 dB. Fig. 3 shows the RMSE of DOA estimation versus the number of active grid points based on 200 Monte Carlo runs. Fig. 4 shows the total CPU time versus the number of active grid points. It is shown that the choice of the number of active grid points does not affect the performance much, as long as it is larger than the number of signals K; however, it does bring a light computational cost, if a small value is chosen. Another observation is that if the grid is very coarse, knowing the exact value of K can improve the DOA estimation performance.

V. CONCLUSION
We have proposed a computationally efficient root SBL method that approaches the problem of DOA estimation with off-grid model error. Unlike the existing off-grid SBL method in [11] that applies a linear approximation to the true DOA, our new method considers the sampled locations in the coarse grid as the adjustable parameters, and utilizes an EM algorithm to iteratively refine the grid. We further illustrate that each updated grid point can be simply achieved by the root of the polynomial (19). Simulation results demonstrate that the computational complexity is significantly reduced and the modeling error can be almost eliminated.