A Space-Time Adaptive Processing Method Based on Sparse Bayesian Learning for Maneuvering Airborne Radar

Space-time adaptive processing (STAP) is an effective technology in clutter suppression and moving target detection for airborne radar. Because airborne radar moves at a constant acceleration, and there is a lack of independent and identically distributed (IID) training samples caused by the heterogeneous environment, using the conventional STAP methods directly cannot ensure a good performance. To eliminate these effects and improve the performance of clutter suppression, a STAP method based on a sparse Bayesian learning (SBL) framework for uniform acceleration radar is proposed here. This paper introduces the signal model of the uniform acceleration radar. To promote the sparsity, a generalized double Pareto (GDP) prior is introduced into our method, and the estimation of hyper parameters via expectation maximization (EM) is given. The effectiveness of the proposed method is demonstrated by simulations.


Introduction
Space-time adaptive processing (STAP) [1,2] is a two-dimensional adaptive filtering technique, which is widely used in the airborne early warning (AEW) radar. It plays an important role in suppressing strong clutter so that weak, slow-moving targets can be detected. Whether the STAP can achieve a good clutter suppression performance depends on the estimation accuracy of the clutter covariance matrix (CCM) of cells under test (CUT), which is usually estimated by the independent and identically distributed (IID) training samples.
Unfortunately, it is known from the Reed-Mallet-Brennan (RMB) [3] rule that the accurate estimation of the CCM requires at least twice the system degrees of freedom (DOF), and it is quite difficult to obtain in the nonhomogeneous environment. Therefore, the performance of the traditional STAP methods degrades severely.
To solve this problem, many algorithms have been proposed to ameliorate the performance of clutter suppression.
From the start of the STAP technique application, reduced-dimension and reducedrank methods [4][5][6][7][8], including the extended factored approach (EFA) [9] and multistage wiener filter (MWF) [10], have been developed and widely used because of the relatively fewer requirements of IID training samples and faster convergence rate compared with the full STAP. The EFA method is easy to implement in practice because the DOF of radar systems can be reduced using linear transformations. However, this method still needs many training samples, especially for the systems with a large DOF. The MWF method has the advantage of a small calculation when the suitable clutter rank is selected. It is worth noting that regardless of whether the selected rank is higher or lower than the actual rank, it will greatly reduce the clutter suppression performance of the MWF algorithm. Direct data domain (DDD) algorithms [11,12] can cleverly avoid clutter inhomogeneity with just a training sample. The benefit of the DDD algorithm is achieved at the cost of aperture loss. Furthermore, KA algorithms [13,14] have been proposed, which take full advantage of parameterized models and a priori knowledge. The performance of KA algorithms relies heavily on the accuracy of a priori knowledge, and this information is hard to obtain in practice.
Over the past twenty years, the sparse recovery STAP (SR-STAP) methods [15][16][17][18] have been proposed to estimate the CCM by using the sparsity of space-time snapshot clutter data in the angle-Doppler domain, which can achieve a good performance. These kinds of methods refine the angle-Doppler two-dimensional plane into a finite grid point, and then use the information of the grid points to obtain an accurate representation of a clutter space-time profile of CUT data. The CCM is reconstructed with the obtained clutter space-time profile and the space-time dictionary matrix. The most important advantage of the SR-STAP methods is that it greatly reduces the need for training samples. However, in practical application, these SR-STAP algorithms are particularly sensitive to regularization parameters, so it is necessary to choose proper regularization parameters. Otherwise, the performance of sparse recovery will be affected, resulting in the degradation of clutter suppression.
The sparse Bayesian learning (SBL) algorithms proposed by Tipping [19] have been introduced to sparse recovery in the single measurement (SMV) case, which has drawn much attention. While retaining the advantage of sparse recovery, there is no need to choose regularization parameters. Even in the case of strong atomic correlation in the dictionary, a good performance can still be achieved. In order to utilize the information of the adjacent samples, Wipf extended it to multiple measurements (MMV) [20]. Most recently, SBL was introduced in the STAP with the MMV named M-SBL-STAP algorithm [21], which shows the potential of the SBL in the STAP.
In this paper, for the case of accelerated motion of the radar platform, the direct application of the traditional STAP method cannot achieve the equal performance of clutter suppression. Thus, we will improve the performance of the STAP method with the sparse Bayesian framework. Unlike most M-SBL-STAP methods that use Gaussian distribution priors and Laplace distribution priors, we used the generalized double Pareto (GDP) [22] prior instead in order to improve the performance of the M-SBL-STAP methods. Compared with the priors mentioned above, the GDP prior has a stronger sparsity.
The main contributions of this paper are listed as follows: 1.
Based on the sparse Bayesian learning using the generalized double Pareto prior, a hierarchical Bayesian model built for the STAP is proposed here.

2.
Improving the performance of clutter suppression for the radar platform that moves at a constant acceleration.

3.
Compared with other STAP methods, the simulation results and performance analysis are given to demonstrate the effectiveness of the proposed method.
The work is organized as follows: Section 2 introduces the data model and Section 3 gives the details of the key theory and implements the steps of the proposed algorithm. Section 4 provides the simulation results and performance analyses of the proposed algorithm. Conclusions are drawn in Section 5.
Notations: Boldface uppercase and boldface lowercase letters are used to represent matrices and vectors, respectively; ⊗ and denote the Kronecker product and the Hadamard product, respectively; (·) T , (·) * and (·) H denote the transposition, conjugate, and conjugatetransposition, respectively; (·) −1 denotes matrix inversion and |·| denotes the absolute value operation; · F denotes the Frobenius norm and · 2,0 stands for a mixed norm. E(·) is the expectation operation.

Signal Model and Problem Formulation
This paper considers the airborne radar system equipped with a uniform linear array (ULA) of N elements, and this radar system transmits M pulses at a constant pulse repetition frequency (PRF) f r during a coherent processing interval (CPI). The height of the radar platform is H. The azimuth and elevation angles are denoted as θ and ϕ, respectively. As shown in Figure 1, since this paper considers the accelerated motion of the platform and assumes that the direction of the speed and the acceleration are consistent, the airborne radar platform that moves along the X axis and its velocity and acceleration are v and a, respectively. v is often expressed as spectively. As shown in Figure 1, since this paper considers the acceler platform and assumes that the direction of the speed and the accelerat the airborne radar platform that moves along the X axis and its velo tion are v and a , respectively. v is often expressed as presence of the acceleration, it is worth noting that the temporal steerin ent from that of a conventional uniform motion radar, but the spatial the same. Here, we give a brief explanation. First, we assume that a sc located on the Neglecting the effect of range ambiguity, a general model for the received data of the lth range gate x l ∈ MN × 1 can be expressed as: where c and n denote the clutter component and white Gaussian noise component, respectively. The Nc denotes the number of the clutter patch, which is evenly distributed in each clutter range ring. For the n− th clutter patch, a n is the complex amplitude and s n is the space-time steering vector. The azimuth and elevation angles corresponding to the n− th clutter patch are θ n and ϕ n , respectively. s n is obtained by the Kronecker product of the temporal steering vector s t n and the spatial steering vector s s n . Due to the presence of the acceleration, it is worth noting that the temporal steering vector is different from that of a conventional uniform motion radar, but the spatial steering vector is the same. Here, we give a brief explanation. First, we assume that a scattering point P is located on the l− th range gate and its direction vector can be expressed as: Then, the delay of receiving the echo from the scattering point P caused by the m− th pulse is given as: where c is the speed of light. We can clearly infer from Equation (3) that the Doppler frequency has changed because of acceleration. So, the s t n and s s n can be expressed as follows: where f t = 2v cos θ n cos ϕ n /(λ f r ) f t = a cos θ n cos ϕ n /(λ f 2 r ) f s = d cos θ n cos ϕ n /λ (5) Hence, assuming that the clutter component and the noise component are mutually uncorrelated, the clutter-plus-noise covariance matrix (CNCM) can be written as: where p n denotes the clutter patch power and σ 2 n is the noise power.
In order to obtain the best performance for target detection under the criterion of the maximum SINR, the optimal STAP weight vector of the lth range gate can be given by: where s t is the space-time steering vector of the target.
Combining the sparse recovery techniques with the STAP to make full use of the sparse characteristics of the clutter data, the requirements of the training samples can be significantly reduced. In the application of sparse recovery, first we discretize the angle-angular two-dimensional plane uniformly into K = N s M d grid points, and the N s = η s N(η s > 1) and M d = η d M(η d > 1) represent the number of points in the spatial and temporal domains, respectively. Then, each grid point is associated with a spacetime steering vector b k , (k = 1, . . . , K) and all of these vectors form the dictionary matrix B = [b 1 , . . . , b K ] ∈ C MN×K . Finally, the signal model in the STAP can be represented as where X = [x 1 , . . . , x L ] ∈ C MN×L is the clutter data matrix consisting of L IID training samples. A = [a 1 , . . . , a L ] ∈ C K×L denotes the sparse coefficient matrix. Each column in A is independent and shares the same covariance matrix. N = [n 1 , . . . , n L ] ∈ C MN×L denotes the zero-mean Gaussian noise matrix. According to the sparse recovery algorithms, it can be inferred that the problem of representing X is equivalent to estimating A by solving the following objective function: where ε is the noise error allowance.

Proposed Method
In this section, a STAP method based on the sparse Bayesian learning (SBL) framework using GDP prior for the uniform acceleration radar is presented here.
To begin, the likelihood of measurement X has the following form: Each column in A follows a zero-mean Gaussian distribution: where Γ = diag([γ 1 , γ 2 , . . . , γ K ] T ) and γ k is the variance for k− th row of A. 0 denotes the zero vector.
In [23], the scholars introduce the following hierarchical hyper-priors over Γ, where each hyper-parameter γ k in Γ follows the independent Gamma distribution with the new parameter ξ k : Due to the introduction of the parameter, we need to model the new parameter ξ = [ξ 1 , . . . , ξ K ] as independent Gamma distribution: where h is a small constant and the value can be adjusted. It has been demonstrated in [23] that the marginal probability distribution of a l follows the GDP distribution and an appropriate choice of h can promote the sparsity of the coefficient matrix A.
The prior of A can be expressed as: Based on the Bayesian theorem, the posterior density of A can be obtained by combining the likelihood and prior as follows: The posterior p(A|X, Γ , σ 2 ) is a multivariate Gaussian distribution and is modulated by the hyper-parameters Γ, ξ and σ 2 . By using the classic expectation maximization (EM) method, we can find the hyper-parameters mentioned above.
In the E-step, according to the Equation (17), covariance and mean are given as: In the M-step, we update the hyper-parameters by maximizing the expectation of J(Γ, ξ), which is defined as follows: The expansion of Equation (20) is marked as Equation (21). In Formula (21), µ k and Σ kk denotes the kth row of µ and the kth element on the diagonal of Σ, respectively.
The derivative of J(Γ, ξ) with respect to γ k is as follows: Setting Equation (22) to zero to obtain the update formula of γ k : Similarly, the derivative of J(ξ, γ) with respect to ξ k is given as: Setting Equation (24) to zero to obtain the update formula of ξ k : In addition to the update of the hyper-parameters, we also need to update the noise σ 2 at each iteration, which follows the procedure in [21]: where µ l denotes the lth column of µ. It is worth noting that µ l is different from µ k in Equation (21). The M-GDP-STAP algorithm is shown in Table 1. Table 1. Pseudocode for M-GDP-STAP algorithm.
Step 1: Input the clutter data X and dictionary matrix B.
Step 5: Repeat step 3 and step 4 until a stopping criterion is satisfied.
Step 6: Estimate the CNCM R using the formula where µ l,k denotes the lk− th element of µ and β is the load factor.

Performance Assessment
In this section, a simulation of the results is provided to demonstrate the effectiveness of the proposed method. Ward's clutter model is adopted here, and the clutter-to noise ratio is set to be 40 dB. The 400th range gate is chosen to be the CUT. We set the numbers of spatial and temporal grid points to N s = 4N and M d = 4M, respectively. The parameter h is set to be 0.1. The main parameters of a side-looking airborne radar system are listed in Table 2. In order to evaluate the clutter suppression performance of the proposed method, we used the metric of improvement factor (IF), which is a widely used indicator to assess the performance of our proposed method and other reference methods, defined as follows: The proposed method is compared with the multiple orthogonal matching pursuit (M-OMP-STAP), the multiple focal underdetermined system solver (M-FOCUSS-STAP), and the iterative adaptive approach (IAA-STAP). In addition, we also added the traditional STAP algorithm based on the angle Doppler compensation (ADC-STAP) for comparison. Figures 2 and 3 plot the IF versus the normalized Doppler frequency curves under ideal and non-ideal conditions, respectively. For the non-ideal condition, we added amplitude error (standard deviation 0.03) and the phase random error (standard deviation 2) to the simulation. All simulation results of the IF are averaged in over 100 independent Monte Carlo trials.
under ideal and non-ideal conditions, respectively. For the non-ideal condit ed amplitude error (standard deviation 0.03) and the phase random error (s viation 2) to the simulation. All simulation results of the IF are averaged in dependent Monte Carlo trials. As shown in Figures 2 and 3, it can be clearly seen that the propo achieves better clutter suppression performance than other methods. In Fi the proposed method and the M-SBL achieve a near-optimal performance. the proposed method is slightly better than the M-SBL because of the existen tude and phase errors. In the simulations, the running time of the M-GDP than that of the M-SBL because the latter has a slow convergence.
To better illustrate the performance, Figures 4a-f and 5a-f plot the Capo different STAP methods under the ideal condition and non-ideal condition, r As shown in Figures 2 and 3, it can be clearly seen that the proposed method achieves better clutter suppression performance than other methods. In Figure 2, both the proposed method and the M-SBL achieve a near-optimal performance. In Figure 3, the proposed method is slightly better than the M-SBL because of the existence of amplitude and phase errors. In the simulations, the running time of the M-GDP is much less than that of the M-SBL because the latter has a slow convergence.
To better illustrate the performance, Figures 4a-f and 5a-f plot the Capon spectra of different STAP methods under the ideal condition and non-ideal condition, respectively. the proposed method and the M-SBL achieve a near-optimal performance. In Figure 3, the proposed method is slightly better than the M-SBL because of the existence of amplitude and phase errors. In the simulations, the running time of the M-GDP is much less than that of the M-SBL because the latter has a slow convergence.
To better illustrate the performance, Figures 4a-f and 5a-f plot the Capon spectra of different STAP methods under the ideal condition and non-ideal condition, respectively.  From Figure 4, in the ideal case, we note that that the clutter spectrum obtained by the M-SBL and the M-GDP are the closest to the ideal spectrum both in location and power, which indicates that the exact solution and power for clutter can both be obtained by these two algorithms. The clutter power cannot be recovered directly by the M-OMP and the M-FOCUSS from Figure 4b,c, respectively. From Figure 4d, the M-IAA overestimates the noise level, which leads to weak targets that cannot be detected. the M-SBL and the M-GDP are the closest to the ideal spectrum both in location and power, which indicates that the exact solution and power for clutter can both be obtained by these two algorithms. The clutter power cannot be recovered directly by the M-OMP and the M-FOCUSS from Figure 4b,c, respectively. From Figure 4d, the M-IAA overestimates the noise level, which leads to weak targets that cannot be detected. From Figure 5, in the non-ideal case, we note that that the clutter spectrum obtained by the M-SBL and the M-GDP are the closest to the ideal spectrum both in the location and power, which indicates that the exact solution and power for clutter can both be obtained by these two algorithms. The clutter power cannot be recovered directly by the M-OMP and the M-FOCUSS from Figure 5b,c, respectively. From Figure 5d, the M-IAA From Figure 4, in the ideal case, we note that that the clutter spectrum obtained by the M-SBL and the M-GDP are the closest to the ideal spectrum both in location and power, which indicates that the exact solution and power for clutter can both be obtained by these two algorithms. The clutter power cannot be recovered directly by the M-OMP and the M-FOCUSS from Figure 4b,c, respectively. From Figure 4d, the M-IAA overestimates the noise level, which leads to weak targets that cannot be detected.
From Figure 5, in the non-ideal case, we note that that the clutter spectrum obtained by the M-SBL and the M-GDP are the closest to the ideal spectrum both in the location and power, which indicates that the exact solution and power for clutter can both be obtained by these two algorithms. The clutter power cannot be recovered directly by the M-OMP and the M-FOCUSS from Figure 5b,c, respectively. From Figure 5d, the M-IAA overestimates the noise level, which leads to weak targets that cannot be detected.
We compared the average running time of different algorithms and the results are listed as follows: A total of 100 independent Monte Carlo trials were taken to obtain the average running time. Although the M-SBL and the M-GDP have shown a great performance in IF curves and clutter spectrums, the running time of the M-GDP is much less than that of the M-SBL. Although the running time of some algorithms is far less than that of the M-GDP, their performance cannot be accepted, such as the M-OMP, the M-FOCUSS and the M-IAA. The average running time is shown in Table 3.

Conclusions
In this paper, in order to improve the clutter suppression for the radar platform that moves at a constant acceleration, we study the model of sparse Bayesian learning and introduce the generalized double Pareto distribution prior into the SBL-based STAP framework. Compared with the classical sparse-based STAP method, the GDP prior can promote the sparsity of the sparse coefficient matrix, which will enhance the sparse recovery. At the end of this paper, the simulation results and the performance analysis are given to demonstrate the effectiveness of the proposed method.
We outline a few future research directions. For airborne radar with a higher DOF, the computational complexities of the SBL-based algorithms will be greatly improved; thus, the corresponding fast algorithm will be one of the directions the research will take. In addition, the derivation of algorithms suitable for two-dimensional array radars is also a direction worthy of research.

Conflicts of Interest:
The authors declare no conflict of interest.