Analysis of a Pareto Mixture Distribution for Maritime Surveillance Radar

Abstract

The Pareto distribution has been shown to be an excellent model for X-band high-resolution maritime surveillance radar clutter returns. Given the success of mixture distributions in radar, it is thus of interest to consider the effect of Pareto mixture models. This paper introduces a formulation of a Pareto intensity mixture distribution and investigates coherent multilook radar detector performance using this new clutter model. Clutter parameter estimates are derived from data sets produced by the Defence Science and Technology Organisation's Ingara maritime surveillance radar.

1. Introduction

In recent years at the Australian Defence Science and Technology Organisation (DSTO), clutter mixture models have been analysed, since they can provide tighter fits to spiky clutter in the distribution’s upper tail region [1, 2]. The main application has been to improve the fit of the K-Distribution to real data, resulting in the KK-Distribution introduced in [1]. This mixture model involves two separate K-Distributions combined in the amplitude domain. The justification for this is that one component of the mixture models the Bragg/whitecap scatterers, while the other accounts for sea spikes. It was also shown in [1] that the KK-Distribution modelled real polarimetric high-resolution sea clutter returns better than a single K-Distribution, for the case of horizontally polarised clutter, which tends to be spikier than vertically polarised returns.

Recently, the Pareto distribution has been proposed as an alternative clutter model and has also been validated for high-resolution radar clutter returns [3, 4]. This intensity model has been found to fit both low and high grazing angle clutter returns, for X-band radar operating in a maritime surveillance environment. It has been found to match the performance of the K-Distribution and to be almost as accurate as the KK-model. There are examples where it matches the performance of the KK-Distribution exactly [4]. The advantage of using a Pareto clutter model is that it is simpler mathematically, resulting in more basic and efficient detection schemes [5]. Thus it is of interest to examine whether, in comparison to the success of the KK-Distribution, a Pareto mixture model yields any improvement in clutter modelling and detection performance.

The paper is structured by introducing the Pareto distribution in Section 2 and illustrating why it is becoming of interest to the radar signal processing community. Section 3 proposes a Pareto mixture model and describes the result of fitting it to DSTO’s high-resolution radar clutter sets. Section 4 derives the Neyman-Pearson optimal detector, based upon this clutter model. In addition, a generalised likelihood ratio test (GLRT) detector is produced. Finally, Section 5 gives some examples of detector performance curves, with parameters estimated from real clutter returns.

2. The Pareto Distribution

As maritime surveillance radar resolution has improved over the years, the backscattering from the sea surface, known as sea clutter, has been found to deviate significantly from the traditional Gaussian model [6]. As a consequence of this, many models of clutter have been proposed and analysed, including the Lognormal [7], Weibull [8] and the K-Distribution [9, 10]. The KK-Distribution [11, 12] was proposed to see whether more accurate clutter fits could be achieved, especially for the very spiky horizontally polarised clutter returns. This model proved to be very successful but has a number of serious shortcomings. The first is that the model is a 5-parameter family, although [1] reduces this to 4, based upon DSTO X-band sea clutter data. Previous models studied in radar signal processing have been two parameter distributions. The complexity of the K- and KK-Distribution clutter models also results in difficult mathematical expressions for the coherent multilook detectors. A good example of this is the detector introduced in [13], which shows the KK-Distribution coherent multilook detector is analytically difficult to work with.

The Pareto distribution [14, 15] is named after the Italian economist Vilfredo Pareto (15 July 1848–19, August 1923) [16, 17] and is a power law probability distribution that has been found to be an excellent model of long-tailed phenomena [18]. Its original application was in the modelling of income over a population [19]. It has been used in the modelling of actuarial data; an example is in excess of loss quotations in insurance [20]. Its usefulness as a model for long-tailed distributions has resulted in applications to a number of diverse areas including physics, hydrology, and seismology [18]. It has also found suitable applications in engineering, such as internet teletraffic modelling [21]. From a radar perspective, its simplicity has resulted in the construction of much simpler detection schemes [5, 22, 23].

In the intensity domain, the Pareto distribution has a power law density given by 𝑓_𝛼p,𝛽p(𝑡)=𝛼^p𝛽^𝛼p_p𝑡^𝛼p+1,(1) for 𝑡≥𝛽_p. The parameter 𝛼_p>0 is called the distribution’s shape parameter, which governs the tail behaviour. Parameter 𝛽_p>0 is the scale parameter, which determines where the distribution’s support begins. Its mean is given by 𝔼(𝑋)=𝛼_p𝛽_p/(𝛼_p−1), while its variance is var(𝑋)=𝛽_p²𝛼_p/(𝛼_p−1)²(𝛼_p−2). Hence for finite first and second moments, it is required that 𝛼_p>2.

Figure 1 shows an example of correlated Pareto intensity clutter, generated using the compound Gaussian model defined in [5]. The figure shows the clutter in both the time and frequency domain. In this case, the Pareto shape parameter is 𝛼_p=25, while its scale parameter is 𝛽_p=1. The compound Gaussian model has a clutter covariance matrix Σ that is Toeplitz, with factor 𝜅=0.01, resulting in weak clutter returns. This means the (𝑖,𝑗)th entry of Σ is given by 𝜅^{|𝑖−𝑗|}, for 𝑖 and 𝑗 in the set {1,2,…,𝑁} with 𝑁 the number of rows/columns of Σ.

(a)

(b)

Figure 1 

An example of correlated Pareto clutter returns, where the correlation is small. (a) shows the clutter in the time domain, while the (b) is the result of a fast Fourier transform (FFT) applied to the same clutter.

Figure 2 shows the same Pareto distribution simulated with stronger clutter returns (Toeplitz factor 𝜅=0.99). In this case we see the correlation having a significant effect on the clutter returns.

(a)

(b)

Figure 2 

Correlated Pareto clutter, where the correlation is strong.

These figures show the difficulty of performing target detection in a spiky clutter environment. Both Figures 1 and 2 illustrate that many small targets would easily be missed because of the spikiness of the clutter. In addition, false detections are likely for the same reason. Hence it is important to try to find detectors which allow increased probability of detection, even if there are only minor gains in performance. Small gains are significant for detection of small targets in spiky clutter environments.

3. A Pareto Mixture Model

The extension of the Pareto model (1) to a mixture model can be based upon the definition of the KK-Distribution in [1]. The Pareto mixture distribution, in the intensity domain, is defined by 𝑓_{(𝛼pp1,𝛼pp2,𝛽pp,𝑘)}(𝑡)=(1−𝑘)𝛼^pp1𝛽_pp^𝛼pp1𝑡^𝛼pp1+1+𝑘𝛼^pp2𝛽_pp^𝛼pp2𝑡^𝛼pp2+1,(2) where 𝛽_pp is the distribution’s scale parameter (so that 𝑡≥𝛽_pp) and 𝛼_pp1 and 𝛼_pp2 are two shape parameters. These three parameters are all assumed to be positive real numbers. The mixing coefficient is 𝑘∈[0,1]. The standard Pareto distribution density (1) is recovered by setting 𝑘 to zero or one in (2). In contrast to the KK-Distribution, the Pareto mixture was chosen to have a single-scale parameter, and two shape parameters. The idea behind this was that since the scale parameter determines the support, it would be useful to have only one such parameter and allow for greater distributional fitting through two shape parameters. In MATLAB, an algorithm was created to optimise the fit of the Pareto mixture to real data, in a similar vein to MATLAB’s automatic Pareto fit function gpfit, using maximum likelihood estimation. This algorithm fitted the density (2) to individual data sets based upon the four distributional parameters (𝛼_pp1,𝛼_pp2,𝛽_pp, and 𝑘).

Figure 3 shows an example of the fit of a Pareto mixture distribution to real data, using empirical distribution functions. The data used for this purpose is the DSTO Ingara maritime sea clutter trials data, which is fully polarised high-resolution X-band clutter returns. Details of the data and radar can be found in [1, 24]. This data was collected during a trial in 2004, in the Southern Ocean, and has been well documented. A description of the radar and clutter analysis can be found in [1], while the fitting of the Pareto distribution to this data is described in [4]. Key points about the trial are that the radar operated in a circular spotlight mode, scanning the same patch of sea surface at different azimuth angles. The Ingara radar is a fully polarimetric, operating in X-band [24]. No targets were present in the trial, so that the returns are pure clutter.

Figure 3 

Clutter fit example, showing the empirical distribution functions 𝐹(𝑥) for real Ingara data compared to the fitted Pareto, Pareto mixture (PP), K- and KK-Distributions, as a function of intensity 𝑥. Here the latter two coincide almost exactly, while the Pareto mixture provides the best fit.

This figure shows the Ingara data corresponding to run 34690, at azimuth angle of 255^∘, with horizontal polarisation, as considered in [4]. The fitted Pareto parameters are 𝛼_p=3.191 and 𝛽_p=0.0047, while the Pareto mixture parameters are 𝛼_pp1=3.6637,𝛼_pp2=7.0837,𝛽_pp=0.0083, and 𝑘=0.5402. For comparison, a K- and KK-Distribution were also fitted, using the parameter estimation guidelines in [1]. The fitted K-Distribution has scale parameter 𝑐=74.45 and shape parameter 𝜈=2.863, while the KK-Distribution has scale parameters 𝑐₁=71.45 and 𝑐₂=344.39, the same shape parameter and mixing coefficient 𝑘=0.01. In this case the K- and KK-Distributions are extremely close together, while the Pareto mixture (denoted PP in Figure 3) fits the data more closely. Under magnification, the KK-Distribution is closer to the data than the K-Distribution, especially in the upper tail region. The Pareto mixture clearly fits the data better than a single Pareto distribution. This was observed in all the Ingara data examined, implying that the Pareto mixture model gave a better fit to the data.

4. Detector Performance

Neyman-Pearson detectors [15] can be constructed using the Pareto mixture distribution as an intensity clutter model. The formulation of the basic detection problem is outlined briefly, following the set up in [5]. The radar return is denoted 𝐳, which is a vector of length 𝑁, corresponding to the number of looks. The coherent multilook detection problem is stated as the statistical hypothesis test 𝐻₀∶𝐳=𝐜 against 𝐻₁∶𝐳=𝑅𝐩+𝐜, where all complex vectors are 𝑁×1. The vector 𝐩 is the Doppler steering vector, whose components are given by 𝐩(𝑗)=𝑒^{𝑗2𝜋𝑖𝑓𝐷}, for 𝑗∈{1,2,…,𝑁}, where 𝑓_𝐷 is the normalised target Doppler frequency (−0.5≤𝑓_𝑑≤0.5) and is assumed to be completely known. The complex random variable 𝑅 represents the target, and |𝑅| is its amplitude. The clutter vector 𝐜 will be modelled as a spherically invariant random process (SIRP). This means we assume that the clutter takes the form 𝐜=𝑆𝓖, where 𝓖 is a zero mean complex Gaussian process with covariance matrix Σ and 𝑆 is a univariate nonnegative random variable with density 𝑓_𝑆. Useful references on the SIRP approach include [25, 26]. As in [5], it will be assumed that Σ is semi-positive definite, implying that a Cholesky factorisation can be applied to its inverse. Consequently, there exists a matrix 𝐴 such that Σ⁻¹=𝐴^𝐻𝐴, where 𝐴^𝐻 is the Hermitian transpose. This enables a whitening approach to be applied to the original hypothesis test, since SIRPs are unaffected by a linear transformation [25]. As a result of this, the original hypothesis test can be re-expressed as 𝐻₀∶𝐫=𝐧 versus 𝐻₁∶𝐫=𝑅𝐮+𝐧, where 𝐫=𝐴𝐳,𝐧=𝐴𝐜 and 𝐮=𝐴𝐩. The advantage of this is that the clutter process 𝐧=𝑆𝐴𝓖 is, when conditioned on 𝑆, a complex Gaussian process because 𝐴𝓖 is multidimensional complex Gaussian process, with zero mean but covariance the 𝑁×𝑁 identity matrix. This enables the corresponding densities to be written down more simply. In particular, since 𝑓_𝐧∣𝑆(𝐱)=(1/𝜋^𝑁𝑠^2𝑁)𝑒^{−𝑠−2‖𝐱‖2}, it can be shown that 𝑓_𝐧(𝐱)=1𝜋^𝑁ℎ_𝑁‖𝐱‖²,(3) where ℎ_𝑁(𝑥)=^∞0𝑠^−2𝑁exp−𝑥𝑠²𝑓^𝑆(𝑠)𝑑𝑠(4) is called the characteristic function [25].

For the case of a fixed target model (so that 𝑅 is constant), the density under 𝐻₁ can be shown to be 𝑓_𝐻1(𝐫)=1𝜋^𝑁ℎ_𝑁‖𝐫−𝑅𝐮‖².(5) To extend this to the case where 𝑅 is unknown the GLRT uses the approximation: ‖𝐫−𝑅𝐮‖²≈‖𝐫‖²−||𝐮^𝐻𝐫||²‖𝐮‖²(6) in the density under 𝐻₁.

In order to examine the performance of detectors based upon a Pareto mixture distribution assumption, it is necessary to specify the SIRP that generates the desired marginal intensity distributions. To generate a Pareto mixture SIRP, one uses the characteristic function: ℎ_𝑁(𝑝)=(1−𝑘)𝜁𝑝,𝛼_pp1,𝛽_pp,𝑁+𝑘𝜁𝑝,𝛼_pp2,𝛽_pp,𝑁,(7) where 𝜁 is defined by 𝜁(𝑝,𝛼,𝛽,𝑁)=𝛽^𝛼Γ(𝑁+𝛼)Γ(𝛼)[𝑝+𝛽]^𝑁+𝛼.(8) It can be shown that the SIRP with characteristic function (7) has Pareto mixture marginal intensity distributions. The Neyman-Pearson optimal detector [15] is given by the ratio of densities under the respective hypotheses: 𝐿(𝐫)=𝑓^𝐻1(𝐫)𝑓_𝐻0(𝐫)^𝐻1≷_𝐻0𝜏,(9) where 𝐿(𝐫) is the likelihood function and 𝜏 is the detection threshold. The notation (9) means we reject the null hypothesis 𝐻₀ (return is pure clutter) if the likelihood exceeds the threshold (set by the probability of false alarm). When we reject 𝐻₀ we accept the alternative hypothesis 𝐻₁, implying that there is strong evidence of a target present in the radar return 𝐫. The densities (3) and (5), together with (7), can be applied to (9) to construct the optimal detector explicitly. The GLRT is obtained explicitly by applying (6) to (5) prior to it being applied to (9).

5. Detector Performance Considerations

Using the DSTO Ingara data sets, the performance of the optimal detector, and GLRT, for Pareto mixture clutter models was assessed. Due to the fact that the Ingara data is pure clutter, a synthetic Gaussian target model has been used throughout, as in [5]. Pareto and Pareto mixture clutter models have had their parameters estimated from Ingara data sets. The detector performance curves have then been generated using Monte Carlo simulation, using 10⁶ runs in each case.

The comparison of detectors based upon the Pareto and Pareto mixture clutter models produced some very interesting results. In most cases it was observed that the mixture distribution assumption did not improve detector performance significantly. This is illustrated in Figure 4, which is a detection performance, plotting the probability of detection (Pd) as a function of the target signal to clutter ratio (SCR). In this example, the mixture distribution has 𝛼_pp1=20,𝛼_pp2=100,𝛽_pp=0.05, and 𝑘=0.8. The single Pareto distribution has 𝛼_p=2.347 and 𝛽_p=0.001277. Additionally, 𝑁=10 looks and the normalised Doppler frequency is 𝑓_𝐷=0.5. The false alarm probability has been set to 10⁻⁶. The simulated clutter has been generated with a Toeplitz covariance matrix with 𝜅=0.4. This example corresponds to a typical horizontally polarised Ingara data set. Figure 4 shows the optimum detector result as well as the GLRT detector. As can be observed, although the mixture model provided a better clutter fit, there is only a minor gain in detection performance.

Figure 4 

Comparison of detector performance curves for Pareto and Pareto mixture models. The probability of detection (Pd) is plotted as a function of the signal to clutter ratio (SCR). The optimal detectors are denoted Optimal (PP) for the mixture model and Optimal (P) for the standard Pareto counterpart. The GLRT for the Pareto mixture is denoted GLRT (PP), while the same for a single Pareto model is GLRT (P). The figure shows that the mixture model provides only a very small detection gain.

Figure 5 is a second example of detector performance. Here, the same clutter parameters and normalised Doppler frequency have been used as for Figure 4. However, the number of looks is 𝑁=5 and the Toeplitz factor has been increased to 0.7. Additionally, the false alarm probability has been reduced to 10⁻⁷. Comparing the pairs of decision rules, we again observe that the mixture model does not provide a significant improvement in detection performance.

Figure 5 

A second example of detector performance, with same legend as for Figure 1. In this case, the clutter strength has been increased (Toeplitz factor of 0.7) and the number of looks decreased. As in the previous figure, the detection gain provided by a mixture model is very small.

Varying the number of looks, normalised Doppler frequency and Toeplitz covariance factor did not significantly alter the results. Hence it appears that although the KK-mixture distribution provides significant detection performance relative to that provided by a detector based upon a single K-Distributed clutter assumption [11], the same is not repeated when the Pareto assumption is applied. This implies a single Pareto distribution is sufficient from a detection performance perspective.

6. Conclusions

Although a Pareto mixture model can provide a better fit to the upper tail region in real radar clutter returns, this has been observed to not significantly improve detection performance, which is the primary application of radar. Increasing the number of degrees of freedom in a model adds computational cost, slowing the radar’s ability to perform its operation in real time. Hence it is computationally better to use a nonmixture model when designing radar detection schemes, under the basic Pareto clutter model assumption.