Statistical characterization of scattering delay in synthetic aperture radar imaging

Distinguishing between the instantaneous and delayed scatterers in synthetic aperture radar (SAR) images is important for target identification and characterization. To perform this task, one can use the autocorrelation analysis of coordinate-delay images. However, due to the range-delay ambiguity the difference in the correlation properties between the instantaneous and delayed targets may be small. Moreover, the reliability of discrimination is affected by speckle, which is ubiquitous in SAR images, and requires statistical treatment. Previously, we have developed a maximum likelihood based approach for discriminating between the instantaneous and delayed targets in SAR images. To test it, we employed simple statistical models. They allowed us to simulate ensembles of images that depend on various parameters, including aperture width and target contrast. In the current paper, we enhance our previously developed methodology by establishing confidence levels for the discrimination between the instantaneous and delayed scatterers. Our procedure takes into account the difference in thresholds for different target contrasts without making any assumptions about the statistics of those contrasts.

1. Introduction. Signal processing algorithms for synthetic aperture radar (SAR) imaging require a model for signal propagation and a model for scattering about the target. For example, standard SAR assumes dispersionless propagation with the speed of light and a point scatterer with constant instantaneous reflectivity. These models are deterministic. In addition, stochastic treatment may be justified for certain imaging scenarios. For example, physical characteristics of a turbulent medium (e.g., density, velocity, etc.) are typically considered random fields [22,24]. Accordingly, if radar signals propagate through such a medium, the resulting SAR image is described in statistical terms [20,11,12].
Stochastic approach can also be used to describe the scattering of SAR signals about the target. A detailed stochastic treatment of instantaneous extended targets can be found in [23]. The goal of the current study is to address the scattering delay and its detection in SAR. For scatterers with delayed response, a stochastic model for SAR imaging has been built in our work [14]. It assumes the deterministic propagation with constant speed as in standard SAR, while scattering about both instantaneous and delayed targets is described in stochastic framework. Hereafter, we extend the results of [14] by introducing confidence levels for the detection of targets with delayed response.
A delayed component in scattering may carry valuable information about the properties of the scatterer, such as its internal structure and characteristic size. The main difficulty in detecting the scattering delay is to separate it from the propagation delay, which is at the core of SAR reconstruction. The authors of [10] propose to interpret the scattering delay as a third dimension added to target reflectivity and SAR image (on top of two spatial coordinates). Then, a point scatterer in space ∼ δ(z − z 0 ) used in standard SAR is replaced with a point scatterer in space-time ∼ δ(z − z 0 )δ(t − t 0 ). Subsequent analysis in [10] focuses on the resulting coordinate-delay point spread function (PSF). The aforementioned difficulty in separating scattering delay from the propagation delay manifests itself via slow decay of PSF along certain directions in the space of its arguments, known as ambiguity directions. This effect is called the range-delay ambiguity.
However, a completely deterministic treatment like that of [10] does not take into account the stochastic effects in scattering, and hence cannot be applied directly to distributed SAR targets [23]. A key manifestation of stochasticity in scattering is speckle, which may be thought of as strong and rapid variations of the amplitude and phase of a SAR image while the target parameters of interest remain smooth. Speckle is common in images of most natural and man-made targets when illumination is coherent [15], which is the case for SAR. In the presence of speckle, weak and slow variations along the ambiguity directions in the SAR image can be undetectable, as demonstrated in [14]. This means that the discrimination between the instantaneous and delayed targets becomes unreliable, unlike in the deterministic case considered in [10].
A standard approach to problems of this kind is two-fold. To increase the reliability of classification one can increase the sample, i.e., the amount of data supplied to the discrimination functional. In [14], we have demonstrated the advantages of a bigger sample size. To quantify the reliability of classification outcomes, one needs to employ the confidence levels, which is the primary focus of the current work. Specifically, we refine the discriminating functional for coordinate-delay SAR images, introduce the confidence levels for it, and demonstrate the performance of the discrimination procedure for various system and target parameters. Confidence levels are crucial for applications that rely on the analysis of SAR images because in practice, it is often impossible to obtain additional images of the same target under similar conditions. Section 2 presents the coordinate-delay SAR imaging procedure, builds the corresponding imaging operator, and analyzes its properties in terms of the point spread function. In Section 3, we introduce models for instantaneous and delayed scatterers and analyze the autocorrelation properties of the resulting coordinate-delay SAR images. Section 4 presents two models of radar targets to be used in the discrimination problems. A binary classification procedure and its extension that uses confidence levels are introduced in Sections 5 and 6.1, respectively. In Section 6.2, we analyze the cumulative distribution functions for the values of the discrimination functional and relate them to the quality of discrimination by the original and extended classifier. In the same section, we also introduce the confidence levels and the corresponding threshold values for the discrimination functional in the case of known target contrast. The generalization of confidence levels to all target contrasts is presented in Section 6.3. To assess the efficiency of the proposed approach to discrimination between the instantaneous and delayed targets, we use the Monte-Carlo simulation procedure described in [14]. It lets us build ensembles of sampled SAR images and analyze the statistics of the discriminating functional. Section 7 presents the results of simulation in terms of achievable discrimination quality given a certain confidence level. Section 8 discusses possible future work. Additional bibliography can be found in [10,14].
2. Coordinate-delay SAR image and point spread function.
2.1. Assumptions, data collection geometry, and antenna waveform. Let (x 1 , x 2 , x 3 ) be Cartesian coordinates with vertical direction x 3 and origin in the target area on the ground (flat horizontal surface). Let the platform trajectory (i.e., synthetic aperture) be an arc of a circle of radius L in the horizontal plane at altitude H above the ground, with the center of the circle at (0, 0, H). Denote R = √ H 2 + L 2 and let θ be the incidence angle measured from the vertical. Then H = R cos θ and L = R sin θ. Denote the angular extent of synthetic aperture by ϕ T and parametrize the trajectory by the angle ϕ, |ϕ| ϕ T /2 1: (1) In practice, synthetic aperture is a finite set of locations x n = x (ϕ n ) on the trajectory that correspond to the chosen discrete arguments {ϕ n } in (1). Construction of a SAR image involves summation over n, but we assume that the spacing between the angles {ϕ n } is fine enough so as to enable a sufficiently accurate replacement of this summation with integration over the continuous variable ϕ, see formula (9).
The interrogating pulse will be taken as a narrow-band linear chirp (frequency modulated signal) with the carrier frequency ω 0 , bandwidth B, duration τ , and rate α > 0 related by (2) 2ατ = B ω 0 and Bτ 1.
Using the notation χ τ for the characteristic function: we can write the expression for the chirp waveform as follows: We denote the incident and scattered fields by u i and u s , respectively. For the signal emitted by the antenna at a location x , we have u i (t, x ) = P (t), whereas the scattered field recorded by the antenna at the same point is u s x (t) ≡ u s (t, x ); this notation assumes the so-called start-stop approximation, i.e., we ignore the antenna motion during transmit and receive (see [13,Chapter 6] for more detail). Additionally, we assume that the signals emitted from different x don't interfere. We will use z = (z 1 , z 2 , 0) for target coordinates; note that with the antenna trajectory as in (1), the coordinates z 2 and z 1 are usually called range and azimuth (or crossrange), respectively. For any x = x (ϕ) and for any z in the target area, we introduce a shortcut notation R z = |z − x | for the propagation distance. The propagation attenuation factor will be denoted by K(R z ), so that the incident and scattered fields obey the following relations: where c is the speed of light. Expression (4) relies on the first Born approximation [8,13], also referred to as the single-scattering approximation. For a typical signal wavelength λ = 2πc/ω 0 , we will assume where max |z | is the characteristic size of the target. While the fields vary on the scale of λ, the attenuation coefficient varies on a much longer scale of R, e.g., K(R z ) ∼ R −2 z . Hence, given (5), one can simplify (4) as follows:

2.2.
Coordinate-delay SAR image. The goal of SAR imaging is to build an approximate reconstruction of the reflectivity function of the target. In the current formulation, the coordinate-delay reflectivity function ν(t z , z ) defines a relation between the incident and scattered fields; this relation is local in space but distributed in delay time t z : The lower limit of the integral in (7) accounts for the causality of scattering. Combining (7) with (6) and dropping a constant factor K(R) 2 , we obtain the following linear model for the field scattered by a distributed non-instantaneous target: In (8) and below, the integrals without limits will assume integration over the entire real axis. The form (8) implies that ν does not depend on x , the property called angular coherence, which typically requires that the synthetic aperture is far away from the target and the angle subtended by it and centered at the target is small, see Section 2.1. The matched filter is obtained by taking the complex conjugate of the signal waveform (3) and denoted by P (. . .). For an image coordinate y = (y 1 , y 2 , 0), we introduce R ϕ y = |y − x (ϕ)| and build the coordinate-delay image I(t y , y ) by applying the matched filter to the received signal: Replacing x and R z in (8) with x (ϕ) and R ϕ z = |z − x (ϕ)|, respectively, and substituting the result into (9), we obtain the expression for the imaging operator that relates the SAR image and the reflectivity: with the kernel given by (11) W (t y , y ; t z , z ) = Function W in (10) is the point spread function (PSF) of the imaging operator in the following sense: W (t y , y ; t d , z d ) coincides with the image I(t y , y ) due to a space-time point scatterer where δ(. . .) is the Dirac delta function and a = const. In particular, we can see from (11) that max |W (t y , y ; t z , z )| = 1 is attained when y = z and t y = t z .
In order to calculate W for the collection geometry and signal waveform defined in Section 2.1, we take the Taylor expansions of R ϕ y and R ϕ y up to the second term in ϕ, such that according to (1), we have For future convenience, we define while T 0 is obtained from (13) by setting ϕ = 0: We also introduce the following function: Then, the result of integration in (11) under assumptions (2) and (13) can be expressed as follows (see details in [14]): , where the sinc function is defined as sinc ξ = sin ξ/ξ.

2.3.
Properties of the PSF and the range-delay ambiguity. We can see that the list of arguments of W in (16) can be contracted to (17) W (t y , y ; t z , z ) ≡ W (t y − t z , y − z ).
The range-delay ambiguity can be observed in the behavior of W of (17) along the direction defined by In the space of arguments of W , formula (18) defines a family of straight lines that we call ambiguity lines. Consider first the case where ϕ T is so small that the quadratic in ϕ term in (13) can be dropped. In this case, the second argument of Φ in (16) is zero, and (15) reduces to We can see that W lin of (19) is constant along the ambiguity lines (18), (14). In other words, a delayed scatterer (12) is indistinguishable from an instantaneous scatterer at a point z d on the same ambiguity line, i.e., a scatterer with (20) where e 2 is a unit vector in the downrange direction. Note that by setting t y = t z = 0 in (14) and (19) we reduce the latter to the standard SAR formula for the imaging kernel, leading to the following expressions for the range and azimuthal (i.e., crossrange) resolution: Obviously, in the coordinate-delay settings (10) and (19), the value of ∆ R does not completely characterize imaging in range because of the range-delay ambiguity.
Returning to the expression for W given by (13)- (16), we notice that the term Φ(0, ·) yields unambiguous resolution in the range coordinate. Using the explicit form for the marginal function and C and S are the Fresnel integrals [9], we can derive the following asymptotic relations: The width of the main lobe in range due to (23) can be evaluated from |v 2 | √ 360 ∼ 20, which corresponds to the distance of At the same time, for the scatterer (12) and with z d given by (20), we have (25) I(0, z d ) It is possible to interpret formulae (22) and (25) as follows: since Φ(0, 0) = 1 and |Φ(0, v 2 )| < 1 for any v 2 = 0, we can always discriminate between the pair of delayed and instantaneous scatterers given by (12) and (20), respectively. Namely, any coordinate-delay point (t d , z d ) can be tested for containing a delayed scatterer (12) by checking whether this point corresponds to a maximum of |I| on the ambiguity line (18) passing through it. Similarly to the resolution of standard SAR described by (19) and (21), we assume that the location of this maximum can be determined accurate to the width of the main lobe of (25). This yields (26) ϕ 2 T ω 0 t d 40 as a detectability condition for a delayed return due to the point scatterer (12). Additionally, we introduce the parameter κ to characterize the ratio between two range scales, ∆ R and ∆ U , see (21) and (24): such that relation (26) can be rewritten as On the way to more realistic setups, we are going to include a homogeneous background and consider a certain range of delay times rather than a fixed delay t d as in (12). These changes introduce a new effect, called speckle, into the consideration. In the following sections, we will formulate the metrics of detectability of a delayed return in the presence of speckle.
3. Statistical coordinate-delay models of distributed radar scatterers. The deterministic description of SAR imaging presented in Section 2 is appropriate for the space-time point scatterers given by (12). However, once the support of ν becomes nonsingular in any argument, a statistical description appears to be a proper way of describing interference of multiple scatterers within one resolution cell (see details in [14,23,15]).
In particular, a homogeneous instantaneous reflectivity, or the background, is modelled by (28) background: where µ b (z ) is a two-dimensional circular Gaussian white random field: . In (29), . . . denotes statistical averaging and σ 2 b is a positive deterministic constant describing the statistically averaged reflectivity [23]. When a random reflectivity function (28)-(29) is substituted into the imaging operator (10), the resulting image I b (t y , y ) is a stationary circular Gaussian random field, as confirmed in numerous experiments [26,23,15]. For the lower moments of such field, we have . The suggested physical model behind this behavior is that each resolution element is populated by a large number of uncorrelated point scatterers such that each image pixel is a result of interference of individual returns with homogeneously distributed phases [23,15].
In addition to the spatially homogeneous model (28)-(29), we introduce two models of inhomogeneous scatterers. First, we want to model a scatterer that is small in size and exhibiting a certain range of delays. Such scatterer may be representative of an opening into some cavity, e.g., a manhole or window in a wall. Multi-path reflections or structural dispersion inside this enclosed space will result in a range of response delays; if there are many such paths (or many electromagnetic cavity eigenmodes), then the responses with different delay times can be considered essentially uncorrelated, similarly to the scatterers at different locations in the model (28)-(29). Hence, in order to describe a delayed scatterer (henceforth called tscatterer) located at z d (cf. (12)), we modify the model (28)-(29) as follows: where µ t (t) will describe the random behavior of the reflectivity in time. The averaged reflectivity will be characterized by the product of a positive constant σ 2 t and a non-negative dimensionless function of dimensionless time F t (Bt/2), yielding µ t (t) as a non-stationary one-dimensional circular Gaussian white random process with the following properties: . We require certain properties of F t (ζ), in particular, that F t (ζ) = 0 for ζ < 0 from the causality considerations (cf. (7)) and the integrability to satisfy a sufficient condition for the existence of the process (32) [1]. Note that the moments of the image I t (t y , y ) due to such reflectivity function will obey the same relations as I b in (30).
Similarly to (31)-(32), we define an instantaneous inhomogeneous scatterer with the support on a straight line drawn in the range direction: where µ s (s) is another one-dimensional inhomogeneous circular Gaussian white random process: . The support of the scatterer in (33)-(34) is related to that in (31)-(32) through the ambiguity relation illustrated by (20). This means that images due to these two scatterers, I t (t y , y ) and I s (t y , y ), may resemble each other, and the problem of detection of a delayed scatterer may be formulated as a problem of discrimination between these two cases. In this context, it makes sense to assume that F s (ζ) = F t (ζ).
Introduce the dimensionless coordinates η d , ζ d , and ψ d with the origin at certain (t d , z d ) as follows: This way, the coordinate ψ d is aligned with the ambiguity line (18). We take t d = 0 for simplicity and present the second order statistics of images due to the scatterer models (28), (31), and (33), for η d = η d = 0 and ζ d = ζ d in the following concise form: . In (36), the following notations are used.
• For α ∈ {b, s, t} denoting the scatterer types given by (28), (33), and (31), respectively, formula (36) is obtained directly by substituting the expressions for ν α (t z , z ) into (10) and (16). We can choose H α to be dimensionless with sup |H α | ∼ 1 as follows: such that the normalizing coefficients are given by • We introduce the noise term I n to represent the receiver noise and processing errors. We do so by allowing α = n in (36) and formally setting so that all noise terms are uncorrelated. At the same time, the value of σ 2 n K n is a separate problem parameter as specified in (43) below.
In Fig. 1, we plot expectations for image intensities, |I t | 2 and |I s | 2 , in the plane y 1 − z d1 = 0 (cf. (18)) for functions with different values of ζ max . This is done by setting ψ d = ψ d in (36) and (37). As ζ max increases, the parallelogram-shaped level lines of |H t | and |H s | stretch in vertical and horizontal directions, respectively, which is in agreement with the analysis made in [14]. The size of the parallelograms in the direction along the ambiguity lines is determined by the width of the main lobe of Φ(0, κζ), see (23). It can be argued that the shapes in one column of Fig. 1 differ substantially from the respective shapes in the other column if this width is smaller than the support of F t and F s , i.e., (40) κζ max 20, cf. (26 ). Hence, the "difference" between the plots of |I t | 2 and |I s | 2 can be increased by increasing either ζ max or κ (or both). It appears quite feasible to apply traditional image processing techniques, such as edge detection and segmentation [5,2,28,18], to the shapes in Fig. 1 in order to determine the type of the scatterer and its parameters, such as ζ max . However, the intensities of actual images look dramatically different from their statistical averages because of the speckle, and while in practice there is typically only a single image acquisition of the scene of interest, computation of statistical averages from the empirical data is ruled out. For images with speckle, such as simulated in Fig. 2, the mere detection (let alone classification) of the target in certain cases looks problematic. The goal of the next sections is to quantify our ability to distinguish between the t-scatterer and s-scatterer as defined in (31)-(34) in the presence of speckle.

4.
Two models for coordinate-delay SAR images. Using the scatterer types described in Section 3, we are going to build models of radar targets to be used in discrimination problems. For simplicity, we assume that there are two possible configurations of scatterers in the target: where ν b , ν t , and ν s are defined in (28), (31), and (33), respectively, and z d in (41a) and (41b) is the same. The names "s-model" and "t-model" are intended to match the terms "s-scatterer" and "t-scatterer" introduced in Section 3, see (31) and (33). The coordinate-delay SAR images resulting from substitution of (41) into (10) are Figure 1. Plots of |I t | 2 and |I s | 2 for different values of ζ max and κ, see (39) and (27). The dashed lines passing through the origin indicate the ambiguity direction, see (14), (18). For the middle row of plots, the condition κζ max 20 (see (40)) is satisfied, and the difference in the orientation of the parallelogram-shaped level lines is more apparent than for the top and bottom rows. then given by either (42a) I s-model (t y , y ; z d ) = I b (t y , y ) + I n (t y , y ) + I s (t y , y ; z d ) or (42b) I t-model (t y , y ; z d ) = I b (t y , y ) + I n (t y , y ) + I t (t y , y ; z d ), where the terms I n are described by (36) with (38).
For delta-correlated scatterers as in Section 3, the correlation of the image described by (10) is determined by the properties of the imaging kernel W . In particular, the correlation of the image rapidly decreases across the ambiguity lines (18) owing to the sinc term in (16). This allows us to simplify the presentation of the correlation function of images given by (42) by specifying a discrete set of ambiguity lines with large enough spacing between them such that the image values at different lines can be considered uncorrelated.
We assume that the terms I b , I n , I s , and I t in (42a) or (42b) are independent; hence, the moments of the total image are sums of the moments of the corresponding components. In turn, due to the Gaussianity, the moments of each component can be calculated using formula (36). Consequently, formulae (36)-(38) provide a complete description of the statistics of the total image for the arguments corresponding to one and the same ambiguity line. Fig. 2 shows examples of simulated coordinate-delay SAR images due to the targets (41a) (right panels) and (41b) (left panels) with F s and F t given by (39). The relative scatterer intensities, or contrasts, are defined as follows: Note that for Fig. 2, we have chosen p n = 0, i.e., no noise component, while the target contrasts q s = q t = q take three different values. It can be seen that the visible shape features distinguishing the two types of scatterers in Fig. 1 appear much less prominent even for a high contrast of q = 0.8 in the top row of Fig. 2, and practically disappear for the lower contrasts. We will see that the value of q is very important for the effectiveness of the discrimination algorithms described in Sections 5 and 6.1. The effect of the value of p n has not been as prominent, and we always set it to p n = 0.1.

5.
Detection of the delayed response. In this work, we reduce the problem of detection of the delayed response to discrimination between the scenario (41a) involving only instantaneous scatterers and (41b) that includes a component ν t with a scattering delay. Assume that we have some observation data u s x (ϕ) (t), see (8). Using (9), we can build a coordinate-delay image I(t y , y ). We also assume that we can identify candidates for z d as locations of a sharp increase of the image intensity along the range direction at t y = 0, see Figs. 1,2. 1 After that, the neighborhood of each candidate location z d goes through the discrimination procedure described below. This procedure attributes the apparent inhomogeneity at z d to one of the two classes in (41).
We define the sampling procedure for the coordinate-delay SAR image as follows. Let {ζ m } be a discrete set of values of ζ d , see (35), for some z d . These values will be bounded by ζ max used in (39). In addition, we will introduce another parameter, ζ min = 3π, to cut off the transitional effects due to the behavior of F s (ζ) and F t (ζ) given by (39) in the vicinity of ζ = 0. Thus, this set will be defined according to (44) ζ m = πm, where m ∈ N, ζ min ζ m ζ max .
Each value of ζ m defines an ambiguity line passing through a neighborhood of (0, z d ), and, according to the discussion in Section 4, the spacing of π between the adjacent Figure 2. Simulated coordinate-delay SAR images with different contrasts. We use κ = 2.5 and ζ max = 5π, which corresponds to the middle row in Fig. 1. To build each image, we sample ζ d and ψ d with a step of π each, and for every ζ d (i.e., each ambiguity line) generate the multivariate circular Gaussian image components according to (36) with the weights σ 2 α K α calculated via (43). The total images are then computed according to the appropriate expression in (42). Simulation of multivariate normal random variables is performed using the MATLAB c function mvnpdf.
values of ζ m allows us to consider the values on different ambiguity lines uncorrelated.
For each m in (44), we choose N m values of ψ mj , 1 j N m ; these values will play the role of ψ d and ψ d for a given ζ d = ζ m in (36). Then, I mj will denote the coordinate-delay SAR image sampled in a neighborhood of (0, z d ). Using to index the triples of terms on the right hand sides of (42a) and (42b), respectively, we can represent the second order statistics for the corresponding sums with the help of (36) in the following form: (46) Remember that the statistical averages in (46) are unavailable in a practical setting. Instead, we will use the actual data I mj I mj and for each line of (46) build an objective function for optimization with the unknown scatterer intensities as the design variables. For each of the two scenarios in (42), i.e., for S and T in (45), our discrimination algorithm will seek a set of values for unknowns σ 2 α that maximizes the probability density of the image with the statistics described by (46). Then, we will choose the model that yields the larger of the two maxima. Essentially, this is a maximum likelihood (ML) based procedure [19,23].
The probability density of the sampled image {I mj } for either of the two models is calculated as follows. For each m we create a real-valued vector r m of dimension 2N m : (47) r m = (Re I m1 , Im I m1 , Re I m2 , Im I m2 , . . . Re I mNm , Im I mNm ) T .
Then, due to the circular Gaussianity and independence of all I α , formula (46) can be recast as where each individual 2 × 2 block on the right hand side is given by Then, we extend formulae (50) by including the data from multiple ambiguity lines given by a set of ζ m . In a simplified treatment suggested in Section 4, we consider the data for different ζ m uncorrelated. Then, for the full dataset vector R that combines all r m -vectors (47): Inverse Problems and Imaging Volume 14, No. 3 (2020), 511-533 we have The actual dataset vector Q representing a sampled image has the same structure as R of (51a): where each vector q m corresponds to image values taken at a certain ambiguity line. We will then consider (cf. (50), (51b)) as functions of the unknown scatterer intensities {σ 2 α } that enter M (m) via (49) for each of the models in (41). The functions p s-model (Q) and p t-model (Q) are called the likelihood functions [19]. The discrimination procedure solves two optimization problems formulated as follows: (54)p s = max subject to σ 2 b , σ 2 n , σ 2 s , σ 2 t 0. The resultingp s andp t yield the maximum likelihood (ML) values for the corresponding scatterer models. It is common to consider the logarithm of the likelihood rather than the likelihood itself. Accordingly, we introduce (55) l = logp t − logp s , and the classification based on the comparison of the two maxima [14] is performed as follows: (56) if l > 0 then the target is classified as a delayed scatterer (41b), (42b) else the target is classified as an instantaneous scatterer (41a), (42a).
6. Statistical characterization of observations. if l > l + then the target is classified as a delayed scatterer (41b), (42b) elseif l < l − then the target is classified as an instantaneous scatterer (41a), (42a) else the classification result is uncertain.
As compared to algorithm (56), we have introduced two classification thresholds, l − and l + , to be defined in Section 6.2, instead of a single threshold l = 0. Accordingly, we have three classification outcomes: s-model, t-model, and uncertain, instead of the two outcomes in algorithm (56). The second fundamental reason for possible misclassification is that formulae (53) yield a nonzero probability density for any model and any data Q, so a certain fraction of errors is inevitable regardless of the algorithm. The quality of the classification is characterized by the frequency of errors. Suppose that we have obtained a representative ensemble of sampled images of a target described by the s-model and another such ensemble for the t-model. Executing either of the algorithms (56) or (57) on each image in these ensembles, we can evaluate the performance of the classification by means of the confusion matrices as in Table 1. The rows named "input: s" and "input: t" denote the models (42a) and (42b), respectively, whereas the columns correspond to the outcomes of the particular classification algorithm. The ideal confusion matrix in Table 1(a) will have r s = r t = 0, whereas for Table 1(b) this will be r s = r t = r s = r t = 0.
The frequency of classification errors depends on several factors. System parameters, e.g., bandwidth, aperture width, etc., form one group. Another group contains parameters of the target, such as its contrast. The roles of these groups of parameters have been investigated in [14]. Ultimately, the classification quality depends on the discrimination algorithm. The choice of the algorithm and its settings may depend on the specific application. For example, a wide gap between l − and l + in (57) should decrease r s and r t in Table 1(b) at the cost of a large fraction of uncertain outcomes, i.e., large values of r s and r t . In Section 6.3, we discuss a procedure whereby the classification errors can be kept below a specified level. 6.2. Confidence levels for classification with a given target contrast. A standard approach to controlling the estimation errors for noisy measurements includes confidence intervals or levels [19]. In parameter estimation problems, a confidence interval is built around the measured value of a certain parameter to indicate a possible range for the true value of this parameter. The boundaries of such an interval are determined from an ensemble of measurements of the parameter of interest or a probability distribution function representing it, such that only a small percentage of outliers, say 5%, falls beyond this interval. Similarly, for a classification problem, an individual measurement can be assigned a numerical characteristic that will express the certainty that this observation falls into (or beyond) a specific category [19]. For the procedure described in Section 5, the value of l defined by (55) can play the role of such parameter. Table 1. Confusion matrices: (a) for classification using algorithm (56); (b) for classification using algorithm (57). The rows correspond to the models in (42), whereas the columns indicate the outcomes of a particular classification algorithm. The entries are relative frequencies of the events calculated for two ensembles with the same contrast, i.e., q s = q t = q.
(a) output: s output: t input: Yet in the case of SAR imaging, building an ensemble of observations to study the statistical properties of the discrimination procedure is not realistic, as explained in Section 1. In [14], we introduced a Monte-Carlo procedure that simulates ensembles of sampled coordinate-delay SAR images of instantaneous and delayed targets, see (42). 2 We used those ensembles to evaluate the efficiency of algorithm (56) for different target contrasts. In the current work, we extend the approach of [14] to define the confidence levels for target classification. In the simplest setting, the simulated ensembles of sampled SAR images represent two scenarios in (42) with equal target contrasts (43): In addition to contrast, each scenario has a set of associated parameters, such as κ, the values of {ζ m } used for sampling, etc. The output of simulation is an ensemble of datasets Q of type (52) that we use in lieu of the actual measurements. While the discrimination procedure does not "know" which of the two target models in (42) and what contrast were used to generate a given dataset Q, we can associate the outcomes of the procedure, and in particular, the set of values ofp s andp t , with the type of underlying model and the values of its parameters. We will describe the statistics of these outcomes with the help of a cumulative distribution function (cdf), 3 which for a real-valued random variable ξ and a given argument x yields the probability that ξ < x: In the context of discrimination between the two types of scatterers, the random variable will be l defined in (55), and we will use the following notations: for the ensemble generated from the s-model with q s = q, and for the ensemble generated from the t-model with q t = q. While the target contrast is explicitly specified as the second argument of cdf in (60), other system and target parameters affecting the probability in (59) will be considered fixed until Section 7. Note that the subscript at a cdf in (60) corresponds to the rows in the confusion matrices in Table 1, whereas the choice of the model in the optimization problem is denoted by the lower index inp s andp t , see (53) and (54). Figure 3(a) plots cdfs of logp s and logp t for a pair of ensembles of sampled images that differ only by the type of the actual inhomogeneous scatterer (note that tiny values of likelihood are not unusual for multidimensional data as in (53)). For the same data, Figure 3(b) plots cdf s (·, q) and cdf t (·, q), see (60), for q s = q t = q = 0.5. This plot clearly shows the separation between these two ensembles, such that most of the values of l are negative for the ensemble generated from the s-model and positive for the ensemble generated from the t-model. This implies that the discrimination results by algorithm (56) are correct in most cases (remember that l is calculated from the observations by a procedure that has no access to the underlying value of contrast or model type). We can establish the following relation between the curves in Fig. 3(b) and the values in Table 1: and, similarly, For example, the value of cdf t (0; q) yields the fraction of targets in the ensemble built from the t-model that algorithm (56) incorrectly classifies as s-targets.
To introduce confidence levels, we choose a small value, say p = 0.05, i.e., 5%, as a threshold for admissible classification errors. In other words, our goal is to make sure that (62) r s p and r t p, see Table 1(b). Then, we define two values, l − and l + , implicitly as solutions to the following equations: (63) cdf t (l − ; q) = p and cdf s (l + ; q) = 1 − p.
The cdfs in (63) are nondecreasing in their first argument, but may be discontinuous. Although this does not present a major obstacle to subsequent considerations, we will assume for simplicity that all cdfs are continuous and monotonic; in this case, solutions l − and l + always exist and unique for p < 1.
We will consider first the case where l − < l + , which is equivalent to as shown in Fig. 3(c) for p = 0.05. For an ensemble of datasets Q generated from the t-model, the frequency of the cases l(Q) < l − will be equal to p. If this value of l − is used as the lower threshold in algorithm (57), with the above dataset as the input, we will also have (64a) Considering an ensemble generated from the s-model, we obtain, in a similar way, the following: (64b) r s = P l(Q) l + | s-model = p.
From relations (64), we see that using the interval (l − , l + ) defined by (63) in algorithm (57), we can keep the rate of classification errors, in particular, r s and r t in Table 1(b), at the predefined level as stated in (62). The rate of uncertain outcomes from algorithm (57) can be expressed as follows: Relations (63) and (65) are illustrated in Figs. 3(c) and 4. The case where from (63) we obtain l − l + (e.g., for the ensembles in Fig. 3(c), this will happen for p = 0.2) can be interpreted as follows: the separation between the two ensembles of values of l is so good that the error level of p can be guaranteed with no need for a confidence interval. In this case, we can use algorithm (56) with any l * ∈ [l + , l − ] as a single threshold. Alternatively, we can find l * as a solution to cdf s (l * , q) + cdf t (l * , q) = 1.
Finally, we should note that taken alone, the definitions of thresholds in (63) can be seen as a way of excluding either extremely large positive or extremely large negative values of l. However, when the thresholds defined in (63) are used in algorithm (57), it is a neighborhood of l = 0 that gets thrown away. This highlights Separation between the graphs of cdf s (· ; q) and cdf t (· ; q) for different values of κ and q. The thick colored vertical bars indicate the percentage of uncertain classifications for the ensembles generated from the s-model (the left set of bars in each plot) and t-model (the right set of bars), see also (65) and Fig. 3(c).
the difference between the problems of parameter evaluation and classification. For the latter, once the rate of classification errors has been fixed at p, see (64), the quality of classification is determined by the percentage of uncertain outcomes, i.e., the values of r s and r t in Table 1(b). From Fig. 4, we can see that as either κ or q increase, the curves of cdf s (· ; q) and cdf t (· ; q) become better separated, the intervals between l − and l + shrink, and the above percentages decrease, as expected.
6.3. Generalization to all target contrasts. The confidence intervals introduced in Section 6.2 depend on the target contrasts q s and q t defined by (43). The latter values should be considered unavailable to the image processing algorithm. Hence, the definitions of l − and l + in (63) should be modified in order to make them independent of target contrasts.
One way of achieving this goal is to use prior information about the target contrasts. For example, we can assume that the probability distribution of the target contrast is known. This means that we can consider q to be a random variable with known probability, and instead of ensembles with a given value of q used in Section 6.2 generate a pair of ensembles, one for the s-model and one for the tmodel, with the given statistics of target contrasts (in general, this statistics can be different for s-target and t-target models). Then, cdf s and cdf t built from these ensembles should replace the cdfs of (60) in definitions (63).
An alternative approach that uses no prior information about the contrast is to take the minimal l − and, correspondingly, maximal l + , over the entire range of target contrasts: With l − and l + redefined as in (66), the procedure (57) should perform with the classification error rates r s and r t not exceeding p for ensembles generated from any probability distributions of contrasts q.
Definitions (66) are used in the numerical experiments described in Section 7. In order to actually compute the maximum and minimum in (66), we have taken the values of q in the range from 0.0 to 0.9 with the step of 0.1. 7. Simulation results. Discrimination between the instantaneous and delayed targets hinges upon our ability to resolve the range-delay ambiguity (see equations (20) and (25)) in the presence of clutter and noise. The quality of discrimination depends on the system and target parameters. In Section 6, the only variable parameter of the model was the contrast q of (58). In this section, we explore the dependence of the discrimination quality on the parameters κ and ζ max that determine the threshold for having the range-delay ambiguity resolved, see (40).
Figs. 5 and 6 show the dependence of the off-diagonal entries of the confusion matrices in Table 1 on ζ max and κ, respectively, for two different values of the target contrast q. The lower half in each color panel represents the second row in Table 1(a) or 1(b), with the colors denoting the individual entries. The upper half corresponds to the first rows in Table 1(a) or 1(b); for clarity of presentation, this part is flipped vertically with respect to the bottom half. The dashed vertical lines are drawn at κζ max = b Φ , where b Φ ≈ 23 [14] is the first local minimum of |Φ(0, ·)|, see (22) and (40).
As expected, the discrimination quality improves with the increase of ζ max , see Fig. 5. A less expected effect that can be observed in Fig. 6 is the saturation of the fraction of reliable classifications for κ 0.4; this may require further attention. Introduction of the confidence level successfully keeps the number of classification errors below p. This, on the other hand, makes a number of correct classifications deemed uncertain.
8. Discussion. The goal of analyzing the scattering delay is to enhance the amount of information supplied by a radar imaging system as compared to standard SAR. At the same time, the proposed methodology uses tools from image classification and pattern recognition. Hence, future developments of this work may come from solving two completely different classes of problems.
In the field of radar imaging, one possibility for the next step is to consider a wider class of functions F s and F t as compared to the indicator functions (39) used in this work. Another option is to explore the stability of the discrimination method to the incorporation of highly coherent components in the received signal: this problem was addressed by means of time-frequency analysis in [4,7]. Additional steps that can improve the applicability and performance of the discrimination procedure are suggested in [14]. Further, the choice of the contrast parameter q for setting the confidence levels, see (58) and (66), may not always be optimal from the standpoint of applications. For example, in a different setting we may be interested in detecting the cases where the dimensionless delay ζ max of (39) exceeds a certain threshold value. Such a problem will require significant modification to the classification algorithm (57) and definitions of the confidence levels (66).
The most noticeable developments in the area of image classification and pattern recognition are currently related to the advances in the artificial intelligence (AI) [3]. The concepts of deep learning and multi-layer convolution neural network (CNN) have received a wide recognition because of their demonstrated efficiency in image classification tasks [17,21]. Yet introducing elements of AI into the analysis of coordinate-delay SAR images may be complicated for several reasons. First, these images are expensive to build, and we cannot expect to be able to obtain the Figure 5. Dependence of the discrimination quality on ζ max , see (39), for two different target contrasts. The left column of plots corresponds to algorithm (56) and Table 1(a), and the right column corresponds to algorithm (57) and Table 1(b). The dashed vertical lines are drawn at κζ max = b Φ ≈ 23, cf. (40). The lower and upper parts of the colored panels represent ensembles generated from the s-model and t-model, respectively. The percentage of correct classifications is shown in green (two different shades are used to distinguish between the ensemble types), incorrect in red, and uncertain in yellow (only the right column of plots).
training sets as massive as those with optical images. Second, the appearance and properties of "signal" and "noise" in coordinate-delay SAR images, see Fig. 2, are very different from those in photography. Hence, besides the convolution and activation operations (i.e., nonlinearity) that are the building blocks of an image classification CNN, we may want to use the transformations that take into account the correlation properties of the images given by (36)-(37). Reports about successful application of deep learning to the problems of automated target recognition in standard SAR images are encouraging [27,6], but at the same time the scarcity of the real data and the difficulties in augmenting it with modelled data are recognized as a major problem [6,25].
As a combination of these two directions, we can apply the modern classification techniques to the entire output of the optimization problems (54). This means that in addition to the minimum values used in the classifier (56), we will take into account the arguments of the minima, i.e., the minimizing scatterer intensities: (p;σ 2 α , α ∈ S) s-model × (p;σ 2 α , α ∈ T ) t-model , see (45) and (54). The resulting parameter space is 8-dimensional, which is hard to process without assistance from some classification algorithm. In our initial trials involving a linear classifier (see, e.g., Figure 6. Dependence of discrimination quality on κ, see (27). The notations are the same as in Fig. 5. [16, Chapter 4]), we did not observe any significant improvements as compared to the method (56) that uses only two out of the eight parameters. This topic may require more attention in the future.