Kronecker PCA based robust SAR STAP
Abstract
This paper proposes a spatio-temporal decomposition for the detection of moving targets in multiantenna SAR. As a high resolution radar imaging modality, SAR detects and localizes non-moving targets accurately, giving it an advantage over lower resolution GMTI radars. Moving target detection is more challenging due to target smearing and masking by clutter. Space-time adaptive processing (STAP) is often used to remove the stationary clutter and enhance the moving targets. In this work, it is shown that the performance of STAP can be improved by modeling the clutter covariance as a space vs. time Kronecker product with low rank factors. Based on this model, a low-rank Kronecker product covariance estimation algorithm is proposed, and a novel separable clutter cancelation filter based on the Kronecker covariance estimate is introduced. The proposed method provides orders of magnitude reduction in the required number of training samples, as well as improved robustness to corruption of the training data. Theoretical properties of the proposed estimation algorithm are established showing significant reductions in training complexity under the spherically invariant random vector model (SIRV). Finally, an extension of this approach incorporating multipass data (change detection) is presented. Simulation results and experiments using the Gotcha SAR GMTI challenge dataset are presented that confirm the advantages of our approach relative to existing techniques.
I Introduction
The detection (and tracking) of moving objects is an important task for scene understanding, as motion often indicates human related activity [29]. Radar sensors are uniquely suited for this task, as object motion can be discriminated via the Doppler effect. In this work, we propose a spatio-temporal decomposition method of detecting ground based moving objects in airborne Synthetic Aperture Radar (SAR) imagery, also known as SAR GMTI (SAR Ground Moving Target Indication).
Radar moving target detection modalities include MTI radars [29, 11], which use a low carrier frequency and high pulse repetition frequency to directly detect Doppler shifts. This approach has significant disadvantages, however, including low spatial resolution, small imaging field of view, and the inability to detect stationary or slowly moving targets. The latter deficiency means that objects that move, stop, and then move are often lost by a tracker.
SAR, on the other hand, typically has extremely high spatial resolution and can be used to image very large areas, e.g. multiple square miles in the Gotcha data collection [34]. As a result, stationary and slowly moving objects are easily detected and located [11, 29]. Doppler, however, causes smearing and azimuth displacement of moving objects [25], making them difficult to detect when surrounded by stationary clutter. Increasing the number of pulses (integration time) simply increases the amount of smearing instead of improving detectability [25]. Several methods have thus been developed for detecting and potentially refocusing moving targets in clutter. Our goal is to remove the disadvantages of MTI and SAR by combining their strengths (the ability to detect Doppler shifts and high spatial resolution) using space time adaptive processing (STAP) with a novel Kronecker product spatio-temporal covariance model, as explained below.
SAR systems can either be single channel (standard single antenna system) or multichannel. Standard approaches for the single channel scenario include autofocusing [14] and velocity filters. Autofocusing works only in low clutter, however, since it may focus the clutter instead of the moving target [14, 29]. Velocity filterbank approaches used in track-before-detect processing [25] involve searching over a large velocity/acceleration space, which often makes computational complexity excessively high. Attempts to reduce the computational complexity have been proposed, e.g. via compressive sensing based dictionary approaches [26] and Bayesian inference [29], but remain computationally intensive.
Multichannel SAR have the potential for greatly improved moving target detection performance [11, 29, 18]. Standard multiple channel configurations include spatially separated arrays of antennas, flying multiple passes (change detection), using multiple polarizations, or combinations thereof [29]. Disadvantages to these approaches include the higher data rate created by collecting multiple channels and the fact that multiple passes involve long delays, registration issues, and having to fly the same orbit more than once [29].
I-a Previous Multichannel Approaches
Several techniques exist for using multiple radar channels (antennas) to separate the moving targets from the stationary background. SAR GMTI systems have an antenna configuration such that each antenna transmits and receives from approximately the same location but at slightly different times [34, 11, 29]. Along track interferometry (ATI) and displaced phase center array (DPCA) are two classical approaches [29] for detecting moving targets in SAR GMTI data, both of which are applicable only to the two channel scenario. Both ATI and DPCA first form two SAR images, each image formed using the signal from one of the antennas. To detect the moving targets, ATI thresholds the phase difference between the images and DPCA thresholds the magnitude of the difference. A Bayesian approach using a parametric cross channel covariance generalizing ATI/DPCA to channels was developed in [29]. Space-time Adaptive Processing (STAP) learns a spatio-temporal covariance from clutter training data, and uses these correlations to filter out the stationary clutter while preserving the moving target returns [11, 18].
A second configuration uses phase coherent processing of the signals output by an antenna array for which each antenna receives spatial reflections of the same transmission at the same time. This contrasts with the above configuration where each antenna receives signals from different transmissions at different times. In this second approach the array is designed such that returns from different angles create different phase differences across the antennas [18, 31, 27, 23, 9]. In this case, the covariance-based STAP approach, described above, can be applied to cancel the clutter [31, 18, 23].
In this paper, we focus on the first (SAR GMTI) configuration and propose a covariance-based STAP algorithm with a customized Kronecker product covariance structure. The SAR GMTI receiver consists of an array of phase centers (antennas) processing pulses in a coherent processing interval. Define the array such that is the radar return from the th pulse of the th channel in the th range bin. Let . The radar data is complex valued and is assumed to have zero mean. Define
(1) |
The training samples, denoted as the set , used to estimate the SAR covariance are collected from representative range bins. The standard sample covariance matrix (SCM) is given by
(2) |
If is small, may be rank deficient or ill-conditioned [29, 18, 21, 22], and it can be shown that using the SCM directly for STAP requires a number of training samples that is at least twice the dimension of [33]. In this data rich case, STAP performs well [29, 11, 18]. However, with antennas and time samples (pulses), the dimension of the covariance is often very large, making it difficult to obtain a sufficient number of target-free training samples. This so-called “small large ” problem leads to severe instability and overfitting errors, compromising STAP tracking performance.
By introducing structure and/or sparsity into the covariance matrix, the number of parameters and the number of samples required to estimate them can be reduced. It has been noted [8, 18, 11] that the spatiotemporal clutter covariance is low rank in general, indicating that the clutter lives in a spatiotemporal subspace of dimension . This reduces the number of parameters describing the covariance matrix from to . Hence, a common approach to STAP clutter cancelation [18, 31] is to estimate a low rank clutter subspace from and use it to estimate and remove the clutter component in the data [2, 18]. We call these methods Low Rank STAP (LR-STAP). Efficient algorithms, including some involving subspace tracking, have been proposed [3, 36]. Other methods adding structural constraints such as persymmetry [18, 9], and robustification to outliers either via exploitation of the SIRV model [17] or adaptive weighting of the training data [15] have been proposed. Fast approaches based on techniques such as Krylov subspace methods [19, 24, 30, 35] and adaptive filtering [12, 13] exist. All of these techniques remain sensitive to outlier or moving target corruption of the training data, and generally still require large training sample sizes [29]. In addition, to the best of our knowledge, none of these techniques explicitly incorporate the known spatio-temporal structure of the data into the covariance estimator. The contribution of this paper is to apply covariance estimation techniques designed to exploit spatio-temporal structure in order to significantly reduce the number of training samples required as well as to provide a degree of robustness to corrupted training data.
We exploit the explicit space-time arrangement of the covariance by modeling the clutter covariance matrix as the Kronecker product of two smaller matrices
(3) |
where is rank 1 and is low rank. In this setting, the matrix is the “temporal (pulse) covariance” and is the “spatial (antenna) covariance,” both determined up to a multiplicative constant.
Kronecker product covariances arise in a variety of applications, including MIMO radar [40], geostatistics [38], recommendation systems [1], multi-task learning [4], and genomics [42]. A rich set of algorithms and associated performance guarantees exist for estimation of covariances in Kronecker product form, including iterative maximum likelihood [39, 38], noniterative L2 based approaches [39], sparsity promoting methods [38, 43], and robust ML SIRV based methods [22]. Many of these methods have been shown to achieve significant reductions in the number of training samples required for estimation, in line with the reduction in the number of parameters in the Kronecker covariance model [39, 38, 43].
In this paper, an iterative L2 based algorithm is proposed to directly estimate the low rank Kronecker factors from the observed sample covariance. Convergence and symmetric positive semidefiniteness of the estimator is established. Theoretical results indicate significantly fewer training samples are required, and it is shown that the proposed approach improves robustness to corrupted training data. Critically, robustness allows significant numbers of moving targets to remain in the training set. We then introduce the Kron STAP filter, which projects away both the spatial and temporal clutter subspaces. This projects away a higher dimensional subspace than does LR-STAP, thereby achieving improved noise and clutter cancelation. We note that this algorithm differs significantly from the set of methods known as Kron PCA, which involves modeling the covariance as a sum of Kronecker products [37, 21, 20].
To summarize, the main contributions of this paper are: 1) the exploitation of the inherent Kronecker product spatio-temporal structure of the clutter covariance; 2) the introduction of the low rank Kronecker product based Kron STAP filter; 3) an algorithm for estimating the spatial and temporal clutter subspaces that is highly robust to outliers due to the additional Kronecker product structure; and 4) theoretical results demonstrating improved signal-to-interference-plus-noise-ratio; and 5) an extension to multipass STAP .
The remainder of the paper is organized as follows. Section II, presents the multichannel SIRV radar model. Our low rank Kronecker product covariance estimation algorithm and our proposed STAP filter are presented in Section IIIwith an extension to the case of moving target detection with multiple passes . Section IV gives theoretical performance guarantees and Section V gives simulation results and applies our algorithms to the Gotcha dataset.
In this work, we denote vectors as lower case bold letters, matrices as upper case bold letters, the complex conjugate as , the matrix Hermitian as , and the Hadamard (elementwise) product as .
Ii SIRV Data Model
Let be an array of radar returns from an observed range bin across channels and pulses. We model as a spherically invariant random vector (SIRV) with the following decomposition [41, 31, 18, 16]:
(4) |
where is Gaussian sensor noise with and we define . The signal of interest is the sum of the spatio-temporal returns from all moving objects, modeled as non-random, in the range bin. The return from the stationary clutter is given by where is a random positive scalar having arbitrary distribution, known as the texture, and is a multivariate complex Gaussian distributed random vector, known as the speckle. We define and note that is determined by the arrangement of stationary scatters on the ground. The means of these components of are zero. The resulting clutter plus noise () covariance is given by
(5) |
The ideal (no calibration errors) random speckle is of the form [29, 11]
(6) |
where . The representation (6) follows because the antenna configuration in SAR GMTI is such that each antenna receives signals emitted at different times at the same points in space [29, 34]. The representation (6) gives a clutter covariance of
(7) |
where
(8) |
depends linearly on the spatial covariance function of the clutter reflectivity, which in turn depends on the spatial characteristics of the clutter in the region of interest [11]. While in SAR GMTI is not exactly low rank, it is approximately low rank in the sense that significant energy concentration in a few principal components is observed over small regions [5].
Due to the long integration time and high cross range resolution associated with SAR, the returns from the general class of moving targets are more complicated. However, if we restrict to targets having constant doppler shift (proportional to the target radial velocity) within a range bin, the return has the form
(9) |
where is the target’s amplitude, , the depend on doppler shift and the platform speed and antenna separation [29], and depends on the target, , and its cross range path. The unit norm vector is known as the steering vector. For sufficiently large , will be small and the target will lie outside of the SAR clutter spatial subspace. Furthermore, as observed in [14], for long integration times the return of a moving target is significantly different from that of uniform stationary clutter, implying that moving targets generally lie outside the temporal clutter subspace [14] as well.
In practice, the signals from each antenna have gain and phase calibration errors that vary slowly across angle and range [29]. It was shown in [29] that in SAR GMTI these calibration errors can be accurately modeled as constant over small regions. Let the calibration error on antenna be and , giving an observed return and a clutter covariance of
(10) |
implying that the in (3) has rank one.
Ii-a Space Time Adaptive Processing
Let the vector be a spatio-temporal “steering vector” [18], that is, a matched filter for a specific target location/motion profile. For a measured array output vector define the STAP filter output , where is a vector of spatio-temporal filter coefficients. By (4) and (9) we have
(11) |
The goal of STAP is to design the filter such that the clutter is canceled ( is small) and the target signal is preserved ( is large). For a given target with spatio-temporal steering vector , we say that the filter is an optimal clutter cancellation filter if it maximizes the SINR (signal to interference plus noise ratio), defined as the ratio of the power of the filtered signal to the power of the filtered clutter and noise [18]
(12) |
where is the clutter plus noise covariance in (5).
It can be shown [11, 18] that, if the clutter covariance is known, under the SIRV model the optimal filter for targets at steering vector is given by the filter
(13) |
where
(14) |
Since the true covariance is unknown, we consider filters of the form
(15) |
and use the measurements to learn an estimate of the best .
Traditionally, the sample covariance (2) has been used to learn the clutter covariance [18]. For the inverse of the sample covariance matrix can be used to reliably estimate the optimal filter (14). Generally, in STAP and a regularized inverse of the sample covariance is often used as an approximation to (14). For this a dimensionality reduction method called clutter subspace processing can be used, giving an alternative filter that approximates (14) by projecting onto a low dimensional subspace. This approach is effective when the clutter subspace is of low rank . Most STAP techniques were developed for classical GMTI radars, for which the covariance is low rank by Brennan’s rule [8]. This approach is also valid for SAR GMTI since the clutter covariance is also low rank [29, 11].
In clutter subspace processing a clutter subspace is estimated using the span of the top principal components of the clutter sample covariance [11, 18]. The corresponding clutter cancelation filter is given by the matrix that projects onto the space orthogonal to the estimated clutter subspace:
(16) |
Since the sample covariance requires a relatively large number of training samples, obtaining sufficient numbers of target free training samples is a practical problem [29, 18]. In addition, if low amplitude moving targets are accidentally included in training, the sample covariance will be corrupted and partially cancel moving targets as well, which is especially problematic in online STAP implementations [29, 3]. The STAP approach discussed below mitigates these problems as it directly takes advantage of the inherent space vs. time Kronecker structure of the clutter covariance .
Iii Kronecker STAP
Iii-a Kronecker Subspace Estimation
In this section we develop a subspace estimation algorithm that accounts for spatio-temporal covariance structure and has low computational complexity.
Following the approach of [39, 21, 37, 22], we fit the low rank Kronecker product model (10) to the sample covariance matrix subject to , where the goal is to estimate . The estimation of the parameters and in (10) is performed by minimizing the following objective function
(17) |
For a matrix define to be its block submatrices, i.e. . Also, let where is the permutation operator such that for any matrix .
The invertible Pitsianis-VanLoan rearrangement operator maps matrices to matrices and, as defined in [37, 39] sets the th row of equal to , i.e.
(18) | ||||
The unconstrained (i.e. ) objective in (17) is shown in [39, 37, 21] to be equivalent to a rearranged rank-one approximation problem, with a global minimizer given by
(19) |
where is the first singular component of .
When the low rank constraints are introduced, a closed-form solution of (17) is no longer available. An alternating minimization algorithm is derived in Appendix A and is summarized by Algorithm 1. In Algorithm 1, denotes the matrix obtained by truncating the Hermitian matrix to its first principal components, i.e.
(20) |
where is the eigendecomposition of , and the (real and positive) eigenvalues are indexed in order of decreasing magnitude. The objective (17) is not convex, but since it is an alternating minimization algorithm, Algorithm 1 gives monotonic convergence of the objective 17 to a local minimum [7]. In practice, we typically initialize LR-Kron with either where are from the unconstrained estimate (19). Monotonic convergence then guarantees that LR-Kron improves on this simple closed form estimator.
Iii-B Robustness Benefits
Besides reducing the number of parameters, Kronecker STAP enjoys several other benefits arising from associated properties of the estimation objective (17).
The clutter covariance model (10) is low rank, motivating the PCA singular value thresholding approach of classical STAP. This approach, however, is problematic in the Kronecker case because of the way low rank Kronecker factors combine. Specifically, the Kronecker product has the SVD [28]
(21) |
where and are the SVDs of and respectively. The singular values are . As a result, a simple thresholding of singular values is not equivalent to separate thresholding of the singular values of and and hence won’ t necessarily adhere to the space vs. time structure.
For example, suppose that the set of training data is corrupted by inclusion of a sparse set of moving targets. By the model (9), the th moving target gives a return (in the appropriate range bin) of the form
(22) |
where are unit norm vectors.
This results in a sample data covariance created from a set of observations with , corrupted by the addition of a set of rank one terms
(23) |
Let and . Let be the eigenvalues of , , and let be the maximum eigenvalue of . Assume that moving targets are indeed in a subspace orthogonal to the clutter subspace. If , performing rank PCA on will result in principal components of the moving target term being included in the “clutter” covariance estimate.
If the targets are approximately orthogonal to each other (i.e. not coordinated), then . Since the smallest eigenvalue of is often small, this is the primary reason that classical LR-STAP is susceptible to moving targets in the training data [29, 18].
On the other hand, Kron-STAP is significantly more robust to such corruption. Specifically, consider the rearranged corrupted sample covariance:
(24) |
This also takes the form of a desired sample covariance plus a set of rank one terms. For simplicity, we ignore the rank constraints in the LR-Kron estimator, in which case we have (19)
(25) |
where is the first singular component of . Let be the largest singular value of . The largest singular value will correspond to the moving target term only if the largest singular value of is greater than . If the moving targets are uncoordinated, this holds if for some , . Since models the entire clutter covariance, it is on the order of the total clutter energy, i.e. . In this sense Kron-STAP is much more robust to moving targets in training than is LR-STAP.
Iii-C Kronecker STAP Filters
Once the low rank Kronecker clutter covariance has been estimated using Algorithm 1, it remains to identify a filter , analogous to (16), that uses the estimated Kronecker covariance model. If we restrict ourselves to subspace projection filters and make the common assumption that the target component in (4) is orthogonal to the true clutter subspace, then the optimal approach in terms of SINR is to project away the clutter subspace, along with any other subspaces in which targets are not present. If only target orthogonality to the joint spatio-temporal clutter subspace is assumed, then the optimal STAP filter is the projection matrix:
(26) |
where are orthogonal bases for the rank and subspaces of the low rank estimates of and , respectively, obtained by applying Algorithm 1. This is the Kronecker product equivalent of the standard STAP projector (16).
Additional information is available, however. Specifically, by (9), no moving target should lie in the same spatial subspace as the clutter. We thus propose the spatial-only filter (which we call spatial-only Kron STAP)
(27) |
to cancel as much of the clutter “subspace leakage” as possible while minimizing target cancelation. This leakage is due to noise and covariance estimation errors. Furthermore, as noted in Section II, if the dimension of the clutter temporal subspace is sufficiently small relative to the dimension of the entire temporal space, moving targets will have temporal factors () whose projection onto the clutter temporal subspace are small. Under these assumptions, it is thus near-optimal to project away both the temporal and spatial clutter subspaces. We thus propose a Kronecker STAP filter of the following form:
(28) |
We denote by Kron-STAP the method using LR-Kron to estimate the covariance and (28) to filter the data. Our clutter model has spatial factor rank (10), implying that the defined in (28) projects the array signal onto a dimensional subspace. This is significantly smaller than the dimensional subspace onto which (26) and unstructured STAP project the data. As a result, much more of the clutter that “leaks” outside the primary subspace can be canceled, thus allowing lower amplitude moving targets to be detected.
Iii-D Multipass STAP
In surveillance applications, it is often of interest to determine what, if anything, has changed in a scene between a reference time and a later time , e.g. disappearance/appearance of parked vehicles, or the appearance of vehicle footprints [29, 2, 6, 32]. When SAR is used for such change detection applications, the radar platform will generally fly past the scene and form a “reference” image at time , and then at time fly a path as close as possible to the original and form a new “mission” image. These images are then compared and changes detected. However, moving targets will almost always be detected as changes, along with the changes in the stationary scene background [29]. When changes of background are of primary interest, moving targets may in fact mask changes in the stationary scene due to displacement and smearing. Hence, it is advantageous to identify moving targets in both scenes prior to or parallel to background change detection. In addition, it may be of interest to detect moving targets in the imagery for their own sake [29]. We thus exploit the additional scene information arising from having two images to better estimate the clutter subspace, and follow STAP with subsequent noncoherent change detection.
Our Kronecker STAP based change detection approach concatenates the spatial channels of both registered phase histories (), forming a “ channel phase history”
(29) |
Since two images are involved with potentially different calibration errors, the clutter subspace is of rank 2. Thus, a rank 2 spatial clutter subspace and a low rank temporal subspace are estimated using LRKron and projected away via the KronSTAP filter. This two pass procedure is easily extended to handle multiple () passes of the radar sensor.
Iv SINR Analysis
For a STAP filter matrix and steering vector , the data filter vector is (15) [18]. With a target return of the form , the filter output is given by (11), and the SINR by (12).
Define to be the optimal SINR, achieved at (14).
Suppose that the clutter has covariance of the form (10). Assume that the target steering vector lies outside both the temporal and spatial clutter subspaces as per above and [18]. Suppose that LR-STAP is set to use principal components. Suppose further that Kron STAP uses 1 spatial principal component and temporal components, so that the total number of principal components of LR-STAP and Kron STAP are equivalent.
Under these assumptions, if approaches zero the SINR achieved using LR-STAP, Kron STAP or spatial Kron STAP with infinite training samples achieves [18] .
We analyze the asymptotic convergence rates under the finite sample regime. Define the SINR Loss as the loss of performance when using the estimate (corresponding to ) as the filter instead of :
(30) |
Let , be the eigenvalues of . Under the Kronecker model, we have
(31) |
since only has one nonzero singular value.
Theorem IV.1 (Lr-Stap Sinr [18]).
For large , the expected SINR Loss of LR-STAP is
Under the Kronecker model we have
(34) |
We now turn to Kron STAP. Note that the Kron STAP filter can be decomposed into a spatial stage (filtering by ) and a temporal stage (filtering by ):
(35) |
where and (28). Under the idealized model in this section, either the spatial or the temporal stage is sufficient to project away the clutter subspace. We assume the naive estimator
(36) |
for the spatial subspace (). This is equivalent to approximating the sample spatial covariance as rank 1. The analysis of [18] thus applies with and , except some of the samples are correlated. Using the Kronecker structure of the covariance it is trivial to show (for the SIRV distribution) that the worst case occurs when all the clutter temporal correlations are all , in which case reduces to an iid sample SCM with Gaussian noise variance and we can directly obtain the following via Theorem IV.1
Theorem IV.2 (Kron STAP SINR).
For large and using the estimator (36), the expected SINR Loss of Kron STAP using the estimator (36) for the spatial subspace satisfies
(37) |
where .
In the small regime this becomes
(38) |
Since by (7) , the gains of using Kron STAP can be quite significant.
Finally we consider the case where errors occurred in estimating the spatial covariance, either due to subspace estimation error or to having a rank greater than one, e.g., due to small calibration errors. Specifically, suppose the estimated (rank one) spatial subspace is , giving a Kron STAP spatial filter . Suppose further that spatial filtering of the data is followed by the temporal filter based on the temporal subspace estimated from the training data. Define the SINR loss from using an estimate of as
(39) |
where is the maximum achievable SINR given that the spatial filter is fixed at . Then it is shown in Appendix B that the expected SINR Loss of the temporal Kron STAP stage is given by the following theorem.
Theorem IV.3 (Kron STAP (temporal stage) SINR).
Suppose that a value for the spatial subspace estimate (with ) and hence is fixed. Let the steering vector for a constant Doppler target be per (9), and suppose that is fixed and is arbitrary. Then for large
(40) | ||||
In the small regime this becomes
(41) |
Note that in the regime relevant when , , where is the first singular vector of . This gives and if is indeed rank one. Hence, can be interpreted as quantifying the adverse effect of mismatch between and its estimate. To avoid cancelation of the moving targets, it is necessary that , and since in the ideal large sample regime all the clutter is removed by the temporal stage, can be smaller than . Hence this slower SINR convergence rate in on a smaller amount of cancelation than the spatial stage (since should be small) is still faster than that of LR-STAP in general.
V Numerical Results
V-a Dataset
For evaluation of the proposed Kron STAP methods, we use measured data from the 2006 Gotcha SAR GMTI sensor collection [34]. This dataset consists of SAR passes through a circular path around a small scene containing various moving and stationary civilian vehicles. The example images shown in the figures are formed using the backprojection algorithm with Blackman-Harris windowing as in [29]. For our experiments, we use 31 seconds of data, divided into 1 second (2171 pulse) coherent integration intervals.
As there is no ground truth for all targets in the Gotcha imagery, target detection performance cannot be objectively quantified by ROC curves. We rely on non ROC measures of performance for the measured data, and use synthetically generated data to show ROC performance gains. In several experiments we do make reference to several higher amplitude example targets in the Gotcha dataset. These were selected by comparing and analyzing the results of the best detection methods available.
V-B Simulations
We generated synthetic clutter plus additive noise samples having a low rank Kronecker product covariance. The covariance we use to generate the synthetic clutter via the SIRV model was learned from a set of example range bins extracted from the Gotcha dataset, letting the SIRV scale parameter in (5) follow a chi-square distribution. We use , , , and , and generate both training samples and a set of testing samples. The rank of the left Kronecker factor , , is 1 as dictated by the spatially invariant antenna calibration assumption and we chose based on a scree plot, i.e., was the location of the knee of the spectrum of . Spatio-temporal Kron-STAP, Spatial-only Kron-STAP, and LR-STAP were then used to learn clutter cancelation filters from the training clutter data. The learned filters were then applied to testing clutter data, the mean squared value (MS Residual) of the resulting residual (i.e. ) was computed, and the result is shown in Figure 1 as a function of . The results illustrate the much slower convergence rate of unstructured LR-STAP. as compared to the proposed Kron STAP, which converges after sample. The mean squared residual does not go to zero with increasing training sample size because of the additive noise floor.
To explore the effect of model mismatch due to spatially variant antenna calibration errors (), we simulated data with a clutter spatial covariance having rank 2 with non-zero eigenvalues equal to 1 and . The STAP algorithms remain the same with , and synthetic range bins containing both clutter and a moving target are used in testing the effect of this model mismatch on the STAP algorithms. The STAP filter response, maximized over all possible steering vectors, is used as the detection statistic. The AUC of the associated ROC curves is plotted on the left in Figure 2 as a function of the number of training samples. Note again the poor performance and slow convergence of LR-STAP, and that spatio-temporal Kron-STAP converges very quickly to the optimal spatial Kron-STAP performance, and more slowly converges to a superior performance as the temporal filter estimate converges.
Finally, we repeat the AUC vs. sample complexity experiment of the previous paragraph with 5% of the training data having synthetic moving targets with random Doppler shifts. The results are shown in Figure 3. As predicted by the theory in Subsection III-B, the Kronecker methods remain largely unaffected by the presence of corrupting targets in the training data, whereas significant losses are sustained by LR-STAP. This confirms the superior robustness of the proposed Kronecker structured covariance in our Kron STAP method.
V-C Gotcha Experimental Data
In this subsection, STAP is applied to the Gotcha dataset. For each range bin we construct steering vectors corresponding to 150 cross range pixels. In single antenna SAR imagery, each cross range pixel is a Doppler frequency bin that corresponds to the cross range location for a stationary target visible at that SAR Doppler frequency, possibly complemented by a moving target that appears in the same bin. Let be the matrix of steering vectors for all 150 Doppler (cross range) bins in each range bin. Then the SAR images at each antenna are given by and the STAP output for a spatial steering vector and temporal steering (separable as noted in (9)) is the scalar
(42) |
Due to their high dimensionality, plots for all values of and cannot be shown. Hence, for interpretability we produce images where for each range bin the th pixel is set as . More sophisticated detection techniques could invoke priors on , but we leave this for future work.
Shown in Figure 7 are results for several examplar SAR frames, showing for each example the original SAR (single antenna) image, the results of spatio-temporal Kronecker STAP, the results of Kronecker STAP with spatial filter only, the amount of enhancement (smoothed dB difference between STAP image and original) at each pixel of the spatial only Kronecker STAP, standard unstructured STAP with (similar rank to Kronecker covariance estimate), and standard unstructured STAP with . Note the significantly improved contrast of Kronecker STAP relative to the unstructured methods between moving targets (high amplitude moving targets marked in red in the figure) and the background. Additionally, note that both spatial and temporal filtering achieve significant gains. Due to the lower dimensionality, LR-STAP achieves its best performance for the image with fewer pulses, but still remains inferior to the Kronecker methods.
To analyze convergence behavior, a Monte Carlo simulation was conducted where random subsets of the (bright object free) available training set were used to learn the covariance and the corresponding STAP filters. The filters were then used on each of the 31 1-second SAR imaging intervals and the MSE between the results and the STAP results learned using the entire training set were computed (Figure 4). Note the rapid convergence of the Kronecker methods relative to the SCM based method, as expected.
Figure 4 (right) shows the normalized ratio of the RMS magnitude of the 10 brightest filter outputs for each ground truthed target to the RMS value of the background, computed for each of the STAP methods as a function of the number of training samples. This measure is large when the contrast of the target to the background is high. The Kronecker methods clearly outperform LR-STAP.
V-D Multipass Kron STAP
Representative two pass Kronecker STAP results are shown in Figure 5, comparing to two pass LR-STAP and to standard (gain calibrated) incoherent change detection. For the STAP methods, noncoherent change detection is performed following filtering by reforming each image (via maximum steering vectors as in the previous section) and subtracting the resulting pixel magnitudes. It can be seen that additional clutter cancelation capabilities can be gained by using Kronecker STAP on multiple passes.
As in the single pass case, Figure 6 shows relative RMSE convergence results and the normalized RMS ratio between targets and background. Again, Kron STAP outperforms the other methods, and both STAP methods outperform standard incoherent change detection.
Vi Conclusion
In this paper, we proposed a new method for clutter rejection in high resolution multiple antenna synthetic aperture radar systems with the objective of detecting moving targets. Stationary clutter signals in multichannel single-pass radar were shown to have Kronecker product structure where the spatial factor is rank one and the temporal factor is low rank. Exploitation of this structure was achieved using the Low Rank KronPCA covariance estimation algorithm, and a new clutter cancelation filter exploiting the space-time separability of the covariance was proposed. The resulting clutter covariance estimates were applied to STAP clutter cancelation, exhibiting significant detection performance gains relative to existing low rank covariance estimation techniques. As compared to standard unstructured low rank STAP methods, the proposed Kronecker STAP method reduces the number of required training samples and enhances the robustness to corrupted training data. These performance gains were analytically characterized using a SIRV based analysis and experimentally confirmed using simulations and the Gotcha SAR GMTI dataset.
Appendix A Derivation of Algorithm 1
We have the following objective function:
(43) |
To derive the alternating minimization algorithm, fix (symmetric) and minimize (43) over low rank :
(44) |
where is the th element of and denotes the complex conjugate of . This last minimization problem (A) can be solved by the SVD via the Eckart-Young theorem [10]. First define
(45) |
and let be the eigendecomposition of . The eigenvalues are real and positive because is positive semidefinite (psd) Hermitian if is psd Hermitian [39]. Hence by Eckardt-Young the minimizer of the objective (A) is
(46) |
Similarly, minimizing (43) over with fixed positive semidefinite Hermitian gives
(47) |
where now describes the eigendecomposition of
(48) |
Iterating between computing and completes the alternating minimization algorithm.
By induction, initializing with either a psd Hermitian or and iterating until convergence will result in an estimate of the covariance that is psd Hermitian since the set of positive semidefinite Hermitian matrices is closed.
Appendix B Proof of Theorem iv.3
After the spatial stage of Kron STAP projects away (27) the estimated spatial subspace (where ) the remaining clutter has a covariance given by
(49) |
By (9), the steering vector for a (constant Doppler) moving target is of the form . Hence, the filtered output is
(50) | ||||
Let and define . Then
(51) |
where and
(52) | ||||
which are proportional to and respectively. The scalar is small if is accurately estimated, hence improving the SINR but not affecting the SINR loss. Thus, the temporal stage of Kron STAP is equivalent to single channel LR-STAP with clutter covariance and noise variance .
Given a fixed , Algorithm 1 dictates (48), (47) that
(53) | |||
which is thus the low rank approximation of the sample covariance of
(54) |
Since , is an SIRV (Gaussian random vector scaled by ) with
(55) |
Furthermore, which is Gaussian with covariance . Thus, in both training and filtering the temporal stage of Kron STAP is exactly equivalent to single channel LR STAP. Thus to prove Theorem IV.3 it is straightforward to apply the analysis in the proof of Theorem IV.1 [18], with the noise variance in training effectively being times the noise variance in testing.
References
- [1] G. I. Allen and R. Tibshirani, “Transposable regularized covariance models with an application to missing data imputation,” The Annals of Applied Statistics, vol. 4, no. 2, pp. 764–790, 2010.
- [2] Y. Bazi, L. Bruzzone, and F. Melgani, “An unsupervised approach based on the generalized Gaussian model to automatic change detection in multitemporal SAR images,” Geoscience and Remote Sensing, IEEE Transactions on, vol. 43, no. 4, pp. 874–887, 2005.
- [3] H. Belkacemi and S. Marcos, “Fast iterative subspace algorithms for airborne STAP radar,” EURASIP Journal on Advances in Signal Processing, vol. 2006, 2006.
- [4] E. Bonilla, K. M. Chai, and C. Williams, “Multi-task Gaussian process prediction,” in NIPS, 2007.
- [5] L. Borcea, T. Callaghan, and G. Papanicolaou, “Synthetic aperture radar imaging and motion estimation via robust principal component analysis,” SIAM Journal on Imaging Sciences, vol. 6, no. 3, pp. 1445–1476, 2013.
- [6] F. Bovolo and L. Bruzzone, “A detail-preserving scale-driven approach to change detection in multitemporal sar images,” Geoscience and Remote Sensing, IEEE Transactions on, vol. 43, no. 12, pp. 2963–2972, 2005.
- [7] S. Boyd and L. Vandenberghe, Convex optimization. Cambridge university press, 2009.
- [8] L. Brennan and F. Staudaher, “Subclutter visibility demonstration,” Adaptive Sensors Incorporated, Tech. Rep. RL-TR-92-21, 1992.
- [9] E. Conte and A. De Maio, “Exploiting persymmetry for CFAR detection in compound-Gaussian clutter,” in Radar Conference, 2003. Proceedings of the 2003 IEEE. IEEE, 2003, pp. 110–115.
- [10] C. Eckart and G. Young, “The approximation of one matrix by another of lower rank,” Psychometrika, vol. 1, no. 3, pp. 211–218, 1936.
- [11] J. H. Ender, “Space-time processing for multichannel synthetic aperture radar,” Electronics & Communication Engineering Journal, vol. 11, no. 1, pp. 29–38, 1999.
- [12] R. Fa and R. C. De Lamare, “Reduced-rank stap algorithms using joint iterative optimization of filters,” Aerospace and Electronic Systems, IEEE Transactions on, vol. 47, no. 3, pp. 1668–1684, 2011.
- [13] R. Fa, R. C. de Lamare, and L. Wang, “Reduced-rank STAP schemes for airborne radar based on switched joint interpolation, decimation and filtering algorithm,” Signal Processing, IEEE Transactions on, vol. 58, no. 8, pp. 4182–4194, 2010.
- [14] J. R. Fienup, “Detecting moving targets in SAR imagery by focusing,” Aerospace and Electronic Systems, IEEE Transactions on, vol. 37, no. 3, pp. 794–809, 2001.
- [15] K. Gerlach and M. Picciolo, “Robust, reduced rank, loaded reiterative median cascaded canceller,” Aerospace and Electronic Systems, IEEE Transactions on, vol. 47, no. 1, pp. 15–25, January 2011.
- [16] G. Ginolhac, P. Forster, F. Pascal, and J.-P. Ovarlez, “Performance of two low-rank STAP filters in a heterogeneous noise,” Signal Processing, IEEE Transactions on, vol. 61, no. 1, pp. 57–61, Jan 2013.
- [17] G. Ginolhac, P. Forster, J. P. Ovarlez, and F. Pascal, “Spatio-temporal adaptive detector in non-homogeneous and low-rank clutter,” in IEEE International Conference on Acoustics, Speech, and Signal Processing, 2009, pp. 2045–2048.
- [18] G. Ginolhac, P. Forster, F. Pascal, and J. P. Ovarlez, “Exploiting persymmetry for low-rank Space Time Adaptive Processing,” Signal Processing, vol. 97, pp. 242–251, 2014.
- [19] J. Goldstein, I. S. Reed, and L. Scharf, “A multistage representation of the Wiener filter based on orthogonal projections,” Information Theory, IEEE Transactions on, vol. 44, no. 7, pp. 2943–2959, Nov 1998.
- [20] K. Greenewald and A. Hero, “Robust kronecker product pca for spatio-temporal covariance estimation,” Signal Processing, IEEE Transactions on, vol. PP, no. 99, pp. 1–1, 2015.
- [21] K. Greenewald, T. Tsiligkaridis, and A. Hero, “Kronecker sum decompositions of space-time data,” in Proceedings of IEEE CAMSAP, 2013.
- [22] K. Greenewald and A. Hero, “Regularized block toeplitz covariance matrix estimation via kronecker product expansions,” in Proceedings of IEEE SSP, 2014.
- [23] A. Haimovich, “The eigencanceler: adaptive radar by eigenanalysis methods,” Aerospace and Electronic Systems, IEEE Transactions on, vol. 32, no. 2, pp. 532–542, April 1996.
- [24] M. Honig and J. Goldstein, “Adaptive reduced-rank interference suppression based on the multistage Wiener filter,” Communications, IEEE Transactions on, vol. 50, no. 6, pp. 986–994, June 2002.
- [25] J. K. Jao, “Theory of synthetic aperture radar imaging of a moving target,” Geoscience and Remote Sensing, IEEE Transactions on, vol. 39, no. 9, pp. 1984–1992, 2001.
- [26] A. S. Khwaja and J. Ma, “Applications of compressed sensing for sar moving-target velocity estimation and image compression,” Instrumentation and Measurement, IEEE Transactions on, vol. 60, no. 8, pp. 2848–2860, 2011.
- [27] I. Kirsteins and D. Tufts, “Adaptive detection using low rank approximation to a data matrix,” Aerospace and Electronic Systems, IEEE Transactions on, vol. 30, no. 1, pp. 55–67, Jan 1994.
- [28] C. V. Loan and N. Pitsianis, “Approximation with kronecker products,” in Linear Algebra for Large Scale and Real Time Applications. Kluwer Publications, 1993, pp. 293–314.
- [29] G. Newstadt, E. Zelnio, and A. Hero, “Moving target inference with Bayesian models in SAR imagery,” Aerospace and Electronic Systems, IEEE Transactions on, vol. 50, no. 3, pp. 2004–2018, July 2014.
- [30] D. A. Pados, S. Batalama, G. Karystinos, and J. Matyjas, “Short-data-record adaptive detection,” in Radar Conference, 2007 IEEE. IEEE, 2007, pp. 357–361.
- [31] M. Rangaswamy, F. C. Lin, and K. R. Gerlach, “Robust adaptive signal processing methods for heterogeneous radar clutter scenarios,” Signal Processing, vol. 84, no. 9, pp. 1653 – 1665, 2004, special Section on New Trends and Findings in Antenna Array Processing for Radar. [Online]. Available: http://www.sciencedirect.com/science/article/pii/S016516840400091X
- [32] K. I. Ranney and M. Soumekh, “Signal subspace change detection in averaged multilook SAR imagery,” Geoscience and Remote Sensing, IEEE Transactions on, vol. 44, no. 1, pp. 201–213, 2006.
- [33] I. Reed, J. Mallett, and L. Brennan, “Rapid convergence rate in adaptive arrays,” Aerospace and Electronic Systems, IEEE Transactions on, vol. AES-10, no. 6, pp. 853–863, Nov 1974.
- [34] S. M. Scarborough, C. H. Casteel Jr, L. Gorham, M. J. Minardi, U. K. Majumder, M. G. Judge, E. Zelnio, M. Bryant, H. Nichols, and D. Page, “A challenge problem for SAR-based GMTI in urban environments,” in SPIE, E. G. Zelnio and F. D. Garber, Eds., vol. 7337, no. 1. [Online]. Available: http://link.aip.org/link/?PSI/7337/73370G/1, 2009, p. 73370G.
- [35] L. Scharf, E. Chong, M. Zoltowski, J. Goldstein, and I. S. Reed, “Subspace expansion and the equivalence of conjugate direction and multistage Wiener filters,” Signal Processing, IEEE Transactions on, vol. 56, no. 10, pp. 5013–5019, Oct 2008.
- [36] M. Shen, D. Zhu, and Z. Zhu, “Reduced-rank space-time adaptive processing using a modified projection approximation subspace tracking deflation approach,” Radar, Sonar Navigation, IET, vol. 3, no. 1, pp. 93–100, February 2009.
- [37] T. Tsiligkaridis and A. Hero, “Covariance estimation in high dimensions via kronecker product expansions,” IEEE Trans. on Sig. Proc., vol. 61, no. 21, pp. 5347–5360, 2013.
- [38] T. Tsiligkaridis, A. Hero, and S. Zhou, “On convergence of kronecker graphical lasso algorithms,” IEEE Trans. Signal Proc., vol. 61, no. 7, pp. 1743–1755, 2013.
- [39] K. Werner, M. Jansson, and P. Stoica, “On estimation of covariance matrices with kronecker product structure,” IEEE Trans. on Sig. Proc., vol. 56, no. 2, pp. 478–491, 2008.
- [40] K. Werner and M. Jansson, “Estimation of kronecker structured channel covariances using training data,” in Proceedings of EUSIPCO, 2007.
- [41] K. Yao, “A representation theorem and its applications to spherically-invariant random processes,” Information Theory, IEEE Transactions on, vol. 19, no. 5, pp. 600–608, 1973.
- [42] J. Yin and H. Li, “Model selection and estimation in the matrix normal graphical model,” Journal of Multivariate Analysis, vol. 107, 2012.
- [43] S. Zhou, “Gemini: Graph estimation with matrix variate normal instances,” The Annals of Statistics, vol. 42, no. 2, pp. 532–562, 2014.