Small-Target Detection between SAR Images Based on Statistical Modeling of Log-Ratio Operator

The log-ratio (LR) operator is well suited for change detection in synthetic aperture radar (SAR) amplitude or intensity images. In applying the LR operator to change detection in multi-temporal SAR images, a crucial problem is how to develop precise models for the LR statistics. In this study, we first derive analytically the probability density function (PDF) of the LR operator. Subsequently, the PDF of the LR statistics is parameterized by three parameters, i.e., the number of looks, the coherence magnitude, and the true intensity ratio. Then, the maximum-likelihood (ML) estimates of parameters in the LR PDF are also derived. As an example, the proposed statistical model and corresponding ML estimation are used in an operational application, i.e., determining the constant false alarm rate (CFAR) detection thresholds for small target detection between SAR images. The effectiveness of the proposed model and corresponding ML estimation are verified by applying them to measured multi-temporal SAR images, and comparing the results to the well-known generalized Gaussian (GG) distribution; the usefulness of the proposed LR PDF for small target detection is also shown.


Introduction
During the last two decades, the problem of change detection [1] from at least two geo-registered synthetic aperture radar (SAR) images, collected at different times, has been recognized to be important for practical applications related to environmental monitoring [2][3][4], damage assessment [5,6], urban studies [7,8], and forest monitoring [9][10][11]. Many approaches [12,13] have been developed in the literature to deal with the problem of change detection. Based on these approaches, a brief and excellent review of state-of-the-art change detection techniques of SAR images has been given in [14]. Additionally, more recent methods, for instance, based on nonlocal means of denoising, mean-shift clustering, and multiple neural-network models, as well as a hierarchical method, have been reported in [7,[14][15][16], respectively.
Specifically, as one of the most important applications of change detection technology, the identification of small objects (i.e., objects occupying only a few pixels, such as vehicles) between two SAR amplitude or intensity images has also been discussed extensively [17][18][19], where small objects only appear in one of two SAR images, and are labelled as a changed class in a fixed scene.
Because of the multiplicative nature of speckle noise [1,20], it has been proven to be more effective to use the ratio operator than the difference operator [21,22] as a change-detection metric to compare two SAR temporal images [23][24][25]. In practice, to uncompress the range of variations of the image 1.
By transforming the image ratio into the log-scale domain, we derive the PDF of the LR statistics. Moreover, ML estimates [34] of the parameters in the LR PDF are derived analytically.

2.
The change detection of vehicles based on the measured data is assessed by combining the proposed parameter estimation and the CFAR technique. As a result, the usefulness of the proposed distribution model and the effectiveness of the corresponding parameter estimation are verified for small target detection.
The rest of this paper is organized as follows. In Section 2, the detailed derivations of the proposed model and ML parameter estimation of the LR statistics are given. Section 3 is devoted to the applications of the proposed model for small target detection; the experimental analysis of typical measured data is also given in this section. Section 4 concludes the paper.

LR Statistics and Distribution
Let us consider two equal-sized and co-registered SAR intensity images with M rows and N columns, I 1 = {I 1 (i, j), 1 ≤ i ≤ M, 1 ≤ j ≤ N} (reference image) and I 2 = {I 2 (i, j), 1 ≤ i ≤ M, 1 ≤ j ≤ N} (test image), acquired over the same geographical region at different times. The ratio of image intensity can be defined as the ratio of average intensity computed for a window of size w × w, as follows: where Ω i,j defines the neighbourhood around pixel (i, j) in a sliding window of size w × w. Furthermore, for a speckled, but otherwise homogeneous region, it is well-known that the statistics of the intensity ratio as shown in Equation (1) can be characterized by the density [31,32] where n is the number of looks, and Γ is the Gamma function. Equation (2) is still valid in the presence of texture, provided that the backscattering coefficient with one n-look pixel is constant [17,32]. The parameter ρ represents the coherence magnitude between the two images. The symbol τ denotes the true intensity ratio, which can be estimated using [31] where E(·) is the expectation. It should be noted that the expectations in Equation (3) are computed by averaging over each entire image. By transforming the intensity ratio into the log-scale domain, the PDF of the LR statistics X, X = log R, can be derived using Equation (2) to give (see Appendix A) As shown in Figure 1, the PDF defined by (4) has the desirable characteristic of being symmetrical about ln(τ). This can be proved mathematically as follows. where n is the number of looks, and G is the Gamma function. Equation (2) is still valid in the presence of texture, provided that the backscattering coefficient with one n -look pixel is constant [17,32]. The parameter r represents the coherence magnitude between the two images. The symbol τ denotes the true intensity ratio, which can be estimated using [31] ( ) ( ) where ( ) E  is the expectation. It should be noted that the expectations in Equation (3) As shown in Figure 1, the PDF defined by (4)   Assuming that two arbitrary variables, x 1 and x 2 , are symmetrical about ln(τ), then x 1 = ln(τ) + a and x 2 = ln(τ) − a, where a is a nonzero constant. Via (4), the PDF value of x 1 is τ + e ln (τ)+a e n ln (τ)+na τ + e ln (τ)+a 2 − 4τρ 2 e ln (τ)+a Likewise, p X (x 2 ) can be written as By multiplying the numerator and denominator of the right-hand side in Equation (5) by the same factor e −(2n+1)a , we obtain As seen from (7), the LR PDF shown in (4) is symmetrical around ln(τ). This symmetry is rather important for the simplification of change detection problems in practice. As is known, there might be positive and negative changed parts among the changed pixels [20,23,28]. Taking the small target detection as an example, we often do not know in advance what targets appear in which of the two images (sometimes the changed targets appear in both co-registered SAR images), thus it cannot be determined whether the target pixels should correspond to the positive or negative changed parts after LR operation. As a result, two thresholds are commonly needed to classify the changed and unchanged categories. However, the number of decision thresholds might decrease to one if it is assumed that the positive and negative changed parts are identically important, because of the symmetrical characteristic of the LR PDF.

ML Estimation
As shown in Equation (4), the PDF of random variable X is parameterized by three parameters, τ, n, and ρ, in which τ is easy to obtain by Equation (3). If τ is estimated by Equation (3) in advance, it can be regarded as known. Therefore, the central task of parameter estimation for the density shown in Equation (4) is focused on the estimates of the number of looks n, and the coherence magnitude ρ. Generally, the ML approach, the method of moments (MoM) [34], as well as the more recent method of log-cumulants (MoLC) [12,28,35], could be candidates for the parameter estimates of a known density. However, an analytical description of the moments of the LR statistics is not easy to acquire by Equation (4), if not impossible, thus MoM and MoLC are not realistic for estimating parameters n and ρ. The ML estimates of n and ρ are derived in the following.
For a sample set X 1 , X 2 , · · · , X m of m independent observations of the random variable X, all with the same distribution obeying Equation (4), the log-likelihood function can be derived as Hence, the ML estimatesn andρ correspond to parameters n and ρ, respectively, and can be written as l(n,ρ) = max n,ρ l(n, ρ) The previous maximization can be realized by applying the first-order derivatives of the log-likelihood function shown in Equation (8) with respect to n and ρ. Equating these derivatives to zero, we obtain the nonlinear equations as follows:  (10) where Ψ(·) is the so-called Psi or Digamma function [34], i.e., the logarithmic derivative of the Gamma function. The solution of the above equations corresponds to the ML estimation of the parameters n and ρ.

Applications to Small target Detection
In this section, we test our proposed distribution model and its parameter estimation approach by applying it to small target detection using measured multi-temporal SAR data. Our goals are to gauge the performance of the proposed distribution model of LR statistics and corresponding parameter estimation, and to determine the usefulness of the proposed distribution model for small target detection. For these purposes, a dataset of VHF-band CARABAS-II SAR images containing changed vehicles is used to test the techniques.

Data Description
The 24 public SAR images included in VHF-band CARABAS-II dataset [18]

Model Fitting
We test the proposed model and parameter estimation by applying it to the data described above. It is worth noting that the distribution of the LR image is expected to be close to Gaussian [33,36]. However, several studies, e.g., [28,33], have shown that a Gaussian approximation is not sufficiently accurate to characterize the LR statistics in practice. Instead, a more general and accurate parametric model, i.e., the GG model, was adopted to describe the statistical behavior of the LR images. Although the GG model is still empirical, it is an attractive candidate for modeling LR images, and has proven to be useful for change detection in SAR images [3,28,33]. Hence, to assess the modeling capability of the proposed model, we compare it with the GG model in this investigation. The GG distribution is characterized by the following PDF [37,38]

Model Fitting
We test the proposed model and parameter estimation by applying it to the data described above. It is worth noting that the distribution of the LR image is expected to be close to Gaussian [33,36]. However, several studies, e.g., [28,33], have shown that a Gaussian approximation is not sufficiently accurate to characterize the LR statistics in practice. Instead, a more general and accurate parametric model, i.e., the GG model, was adopted to describe the statistical behavior of the LR images. Although the GG model is still empirical, it is an attractive candidate for modeling LR images, and has proven to be useful for change detection in SAR images [3,28,33]. Hence, to assess the modeling capability of the proposed model, we compare it with the GG model in this investigation. The GG distribution is characterized by the following PDF [37,38] where γ = 1 σ Γ 3 c /Γ 1 c , in which σ is the variance of the distribution. The parameters µ and c are the mean and shape parameters of the distribution, respectively. The symbol |·| represents the absolute value. It should be emphasized that c = 2 yields the Gaussian density function, which indicates that the GG model encompasses the modeling ability of the Gaussian distribution. Besides, the GG PDF shown in Equation (11) is symmetrical around µ. The parameter estimation technique of the GG model can be found in [37].
Similar to conventional change detection applications, such as environmental monitoring [2][3][4] and damage assessment [5], speckle filtering with a finite-sized window on a SAR image pair is necessary before calculating the LR statistics for small target detection applications, because the presence of speckle can lead to poor performance of the detection algorithm [28,39]. However, an appropriate despeckling window size is not easy to determine, and has to be obtained empirically because no prior training pixels are available. As suggested in [18,19], a 5 × 5 window is feasible for the detection of small targets when the resolution of the SAR image is not larger than 1 m, because it represents a good compromise between speckle smoothing and detection of small targets. In addition, it is worth noting that speckle filtering is less important at the model fitting stage than at the small target detection stage, based on the fact that the theoretical models should accurately describe the LR statistics for any finite size of filtering window. For brevity, we report results using window sizes of 1 × 1 (i.e., non-filtered)  (1) is taken with a 5 × 5 window before calculating the LR statistics, as discussed in Section 3.3.
The parameter estimates of the proposed model and the GG model for two experiments with different window sizes are listed in Table 1. Based on these estimates, the fitting results are provided so that the performance difference in fitting the whole histogram and the histogram tails can be seen clearly (Figures 4 and 5). To clearly show the details of the fitting results, all figures are plotted both on linear and semi-log scales. It is clear from a visual inspection of these figures that almost perfect fitting results are obtained by the theoretical model presented in Equation (4) for the two SAR LR images. It is also clear, from a visual comparison between the estimated PDFs and image histograms, that the proposed model in Equation (4) agrees better with the given LR statistics than the GG model, despite the filtering window sizes, indicating the high precision of the fitting LR statistics using the proposed model. Furthermore, Figure 4b,d, along with Figure 5b,d, show less deviation from the tail of the histograms using the proposed model compared to using the GG model, implying better detection performance of the proposed model, because the histogram tail commonly contributes to false alarms.     Figure 2: (a,b) are the results for a 1 × 1 filtering window displayed on linear and semi-log scales, respectively; (c,d) are the results for a 5 × 5 filtering window displayed on linear and semi-log scales, respectively.  To quantitatively assess the fitting results and compare the different models, specific criteria as goodness-of-fit measures can be employed. In information theory, the Kullback-Leibler (KL) divergence is known to be an effective measure and can realize the global comparison of PDFs [40].
Given the theoretical PDF ( ) p w , and observed PDF ( ) q w , the KL divergence or relative entropy between two densities is defined as [40]: The approximated numerical calculation is given by where ( ) Q w is the value of probability of the observeddata, and ( ) ). In our case, the symmetrized KL divergence is adopted as When the observed density equals the theoretical density, KL D is zero. Otherwise, KL D is a positive value. The KL divergence reflects the overall similarity of the observed and theoretical densities. The smaller the value of KL divergence, the higher the similarity, indicating better fitting accuracy. Based on Equation (14), the KL values of the fitting results shown in Figures 4 and 5 are compared in Table 2. This table clearly confirms the conclusion that were drawn by visual inspection, i.e., the proposed model is in better agreement with the histograms of the LR statistics for different   Figure 3: (a,b) are the results for a 1 × 1 filtering window displayed on linear and semi-log scales, respectively; (c,d) are the results for a 5 × 5 filtering window displayed on linear and semi-log scales, respectively. Note:μ,σ andĉ are estimates of parameters µ, σ and c in Equation (11), respectively;τ,n andρ are estimates of parameters τ, n and ρ in Equation (4), respectively.
To quantitatively assess the fitting results and compare the different models, specific criteria as goodness-of-fit measures can be employed. In information theory, the Kullback-Leibler (KL) divergence is known to be an effective measure and can realize the global comparison of PDFs [40].
Given the theoretical PDF p(w), and observed PDF q(w), the KL divergence or relative entropy between two densities is defined as [40]: The approximated numerical calculation is given by where Q(w) is the value of probability of the observeddata, and P(w) denotes the theoretical probability value. Note that D(q||p) is not symmetrical (i.e., D(q||p) = D(p||q)). In our case, the symmetrized KL divergence is adopted as When the observed density equals the theoretical density, D KL is zero. Otherwise, D KL is a positive value. The KL divergence reflects the overall similarity of the observed and theoretical densities. The smaller the value of KL divergence, the higher the similarity, indicating better fitting accuracy.
Based on Equation (14), the KL values of the fitting results shown in Figures 4 and 5 are compared in Table 2. This table clearly confirms the conclusion that were drawn by visual inspection, i.e., the proposed model is in better agreement with the histograms of the LR statistics for different window sizes than with the GG model. This is to be expected, because the KL values are smaller for the proposed model than for the GG model. In summary, the proposed model shown in Equation (4) has better accuracy in describing the statistical behaviors of the LR operator and is superior to the known GG model. Moreover, the corresponding ML estimators of parameters are also effective. Additionally, unlike the empirical GG model, the proposed model is rigorously derived based on mathematics and SAR principles, making it more appropriate and controllable in practical applications.

Detection Application
For change detection problems, a popular decision technique is the Kittler-Illingworth (KI) unsupervised threshold selection criterion which has been developed under the Gaussian, GG, Nakagami, Log-Normal, and Weibull assumptions for modeling the statistical distributions of changed and unchanged classes in either the LR or ratio images [28,33]. As stated in [23,28], the KI method takes change detection in the ratio or LR images as a two-category classification problem, and performs well on the condition that the changed and unchanged parts have comparable histogram peaks. However, in small target detection, this condition generally cannot be satisfied because small targets only comprise a few changed pixels in a SAR image pair, leading to a monomodal LR image histogram, and the targets labelled as changed pixels fall into the tail part of the histogram. As a result, the KI method can result in several detection errors [23] and cannot be used for small target detection.
In contrast, the CFAR technique is well-known because of its adaptive ability of determining a threshold in the field of small target detection. As for change detection, because there might be positive and negative changed parts among the change pixels, two thresholds T 1 and T 2 (assuming T 1 > T 2 ) are often needed in practice. It is also reasonable to assume that the positive and negative changed parts contribute the same false alarm rate, when one does not know in advance what targets will appear in which of the two images.
Under a given value of false alarm probability, denoted by P f a , the corresponding CFAR thresholds for a specific density p X (x) is obtained from (15) Based on Equations (4) and (11), the threshold T 1 under different background assumptions (i.e., obeying the proposed model or being GG-distributed) can be accurately calculated by Equation (15) with the help of the numerical solution. Because of the symmetry of the two densities, T 1 and T 2 should be symmetrical about ln(τ) for the proposed model, and about µ for the GG distribution. Then, T 2 = 2 ln τ − T 1 in the assumption of the proposed model, and T 2 = 2µ − T 1 in the assumption of the GG distribution. Accordingly, for the test cell η in the LR image, the decision rule of detection can be adopted as target : Given the theoretical false alarm probability P f a = 10 −3 , using (16) for the image pairs shown in Figures 2 and 3, the corresponding CFAR detection results of the LR operator with the distribution assumptions of the proposed model and the GG model are given in Figures 6 and 7, respectively.
First, Figures 6b and 7b show that all targets are well detected by the LR CFAR detector based on the proposed model, and far fewer false alarms occur compared to CFAR detection using the GG model from a visual inspection. This demonstrates the effectiveness of the CFAR method based on the proposed model and allows us to complete the adaptive detection of small targets with a desirable false alarm probability, as expected.  show that all targets are well detected by the LR CFAR detector based on the proposed model, and far fewer false alarms occur compared to CFAR detection using the GG model from a visual inspection. This demonstrates the effectiveness of the CFAR method based on the proposed model and allows us to complete the adaptive detection of small targets with a desirable false alarm probability, as expected.
Second, if detections are found to exist within a radius of 10 m from a ground truth position, the target is declared to be found, and can be counted. All other detected changes that are not considered to be related to a target are regarded as false alarms. A whole performance comparison between the Second, if detections are found to exist within a radius of 10 m from a ground truth position, the target is declared to be found, and can be counted. All other detected changes that are not considered to be related to a target are regarded as false alarms. A whole performance comparison between the CFAR detections based on the proposed model and the GG model is provided in Table 3. In Table 3, the actual false alarm rate (FAR) is estimated by FAR = N c M × N − N t (17) where N c is the number of false alarm pixels, M × N is the size of the detected image, and N t is the number of the detected target pixels. From Table 3, on one hand, the actual FARs generated by both the GG model and proposed model do not rigorously match the theoretical P f a . This is mainly due to that both the GG model and proposed model exhibit some deviation from the measured histogram tails, as shown in Figures 4d  and 5d. On the other hand, it can be seen that both the CFAR methods based on the GG model and the proposed model are able to detect all targets easily, but the detection method using the GG model generates higher FAR than the proposed model in two experiments. In other words, in terms of actual FAR, the CFAR method based on the proposed model shows better target detection performance than the CFAR method based the GG model. This is in accord with the analysis in Section 3.2, since the proposed model exhibits higher modeling capacity of LR statistics than the empirical GG model as a whole, as well as the tail part of the LR histogram, which contributes to most of the false alarms.
In summary, all the results above suggest that the proposed LR PDF and its parameter estimation approach are promising for application in the detection of small targets between SAR images.

Conclusions
In this paper, we analytically derived the PDF of the LR operator, along with the ML estimates of parameters in the LR PDF. The experimental results using real data validate the effectiveness of the proposed model and ML estimation. Moreover, in combination with a commonly used CFAR decision technique, the proposed model and corresponding parameter estimation have excellent potential for applications in small target detection between SAR images. We expect to test more measured data in the future.