Jump Detection In Blurred Regression Surfaces

We consider the problem of detecting jump location curves of regression surfaces when they are spatially blurred and contaminated pointwise by random noise. This problem is common in various applications, including equi-temperature surface estimation in meteorology and oceanography and edge detection in image processing. In the literature, most existing jump-detection methods are developed under the assumption that there is no blurring involved, or that the blurring mechanism described by a point spread function (psf) is completely specified. In this article, we propose four possible jump detectors, without imposing restrictive assumptions on either the psf or the true regression surface. Their theoretical and numerical properties are studied and compared. We also propose a new quantitative metric for measuring the performance of a jump detector. A data-driven bandwidth selection procedure via the bootstrap is suggested as well. This article has supplementary material online.


INTRODUCTION
Surface estimation is an important problem in many scientific areas, including image processing, geology, meteorology, oceanography, and so forth. In many surface estimation problems, the true surfaces have jumps and discontinuities, which are often called edges of the related image intensity surface in the image processing literature (see Gonzalez and Woods 1992;Qiu 2005). Jumps in a surface often represent outlines of objects (e.g., the boundary locations of image objects), or places where the structure of a two-dimensional process changes abruptly (e.g., places where an equi-temperature surface in high sky or deep ocean changes dramatically). Therefore, it is important to detect the jump locations of regression surfaces accurately from the observed data. In practice, however, the observed data are usually not exactly the same as the true surfaces, but are degraded versions of their true surfaces. There are many different sources of degradations. For instance, in aerial reconnaissance, astronomy, and remote sensing, image intensity surfaces are often degraded by atmospheric turbulence, aberrations of the optical system, or relative motion between a camera and an object. Among different degradations, pointwise degradation (or random noise) and spatial degradation (or blur) are the most common in practice. As a demonstration, Figures 1(a)-1(c) shows a true image with a step edge at x = 0.5, a blurred version, and a blurred-noisy version, respectively. The cross-section of the blurred-noisy version at y = 0.5 is shown in Figure 1(d). From these plots, it can be seen that noise alters individual gray levels in a random pattern but the alteration by blur is spatially correlated. From Figures 1(b)-1(d), it can also be seen that the step edge of the original image is still visible in its blurred or blurrednoisy version, although the edge is tapered by some smoothing process that is often unknown in practice and such a smoothing process can vary over location as well (e.g., atmospheric turbulence often causes spatially variant blur in satellite images because both its direction and intensity are usually different at different locations). Other types of degradations involve chro-matic or temporal effects. See Bates and McDonnell (1986) for a detailed discussion about formation of various degradations. This article focuses on jump detection in regression surfaces when both spatial blur and pointwise noise are present in the observed data.
In this article, we propose four jump detectors in cases when the observed data contain both spatial blur and pointwise noise. The idea behind our proposed methods can be explained intuitively as follows. First, although spatial blur would alter the regression surface structure, especially around the true jump location curves (JLCs), the blurred regression surface would change most dramatically at the true JLCs along their normal directions (see Figures 1(b)-1(d)). Therefore, appropriate estimators of the surface gradients can be used for jump detection. Second, spatial blur is a smoothing process; it smooths the regression surface at the true jump locations. If we can recover the jump structure to a certain degree when removing noise, it may help detect jumps. Our proposed jump detectors make use of these properties of spatial blur in different ways. One major feature of these jump detectors is that they do not require restrictive assumptions on either the point spread function (psf) that describes the blurring mechanism or the true regression surface. Their theoretical and numerical properties are studied and compared. Further, in this article, we propose a new quantitative metric for measuring the performance of a jump detector. A data-driven bandwidth selection procedure via the bootstrap is suggested as well.
The remainder of the article is organized as follows: in the next section, our proposed jump detectors are described in detail. Some of their statistical properties are presented in Section 3. Their numerical performance is evaluated in Section 4. Several remarks conclude the article in Section 5. Some technical details are given in a supplementary file.

PROPOSED JUMP DETECTORS
In this section, we first describe our proposed jump detectors in detail in Section 2.1, and then propose a new metric for measuring the jump detection performance in Section 2.2. Based on the proposed performance measure, a bootstrap bandwidth selection procedure is suggested in that subsection as well.

Four Jump Detectors Based on Local Kernel Smoothing
In the image processing literature, a commonly used model for describing the relationship between a true regression surface f and its observed degraded version Z is as follows: where G{f }(x, y) = R 2 g (u, v; x, y)f (x − u, y − v) dudv denotes the convolution between a two-dimensional psf g (see Section 3 for a detailed mathematical description) and the true regression function f , ε(x, y) is the pointwise noise at (x, y), and is the design space of the surface. In model (1), it is assumed that the true regression function f is degraded spatially by g and pointwise by ε, that the spatial blur is linear in the sense that G{a 1 f 1 + a 2 f 2 } = a 1 G{f 1 } + a 2 G{f 2 }, where a 1 and a 2 are two arbitrary constants and f 1 and f 2 are any two images, and that the pointwise noise is additive. In most references, it is further assumed that the psf g, which describes the spatial blurring mechanism, is location invariant. That is, g (u, v; x, y) does not depend on (x, y). See Hall and Qiu (2007) for a related discussion.
For simplicity, let us assume that the design space is = [0, 1] × [0, 1], {(x i , y j ) = (i/n, j/n), i, j = 1, 2, . . . , n} are equally spaced design points in , and {(x i , y j , Z ij ), i, j = 1, 2, . . . , n} follow the model where {ε ij } are iid random errors with mean 0 and unknown variance σ 2 . For a given point ( where h n ∈ (0, 1/2) is a bandwidth parameter, we consider its circular neighborhood: In this neighborhood, let us consider the following local linear kernel (LLK) smoothing procedure: where K ij := K((x i − x)/h n , (y j − y)/h n ) and K is a circularly symmetric bivariate density kernel function defined on the unit disk centered at the origin. Let the solution to {a, b, c} of the problem (3) be denoted as { a(x, y), b(x, y), c(x, y)}. Then, ( b(x, y), c(x, y)) is an estimator of the gradient of the true regression surface f at (x, y). Now, we divide O n (x, y) into two halves, denoted as U n (x, y) and V n (x, y), along the direction perpendicular to ( b(x, y), c(x, y)). More specifically, we define The change in the values of f from V n (x, y) to U n (x, y) should be relatively large if (x, y) is on a JLC. Let be the Nadaraya-Watson local constant kernel (LCK) estimators of f (x, y), constructed from the observations in the two one-sided neighborhoods U n (x, y) and V n (x, y), respectively. Then, f LCK,+ (x, y) − f LCK,− (x, y) would be a good measure of the change in values of f along the direction from V n (x, y) to U n (x, y). Our first jump-detection criterion is defined by its standardized version and the corresponding jump detector is denoted as LCK. Theoretically, it can be checked that the jump-detection criterion LCK n (x, y) has the following two properties (proofs are given in the supplementary material): (i) If f is continuous around (x, y), then we have (ii) If f has a jump at (x, y), then we have Therefore, the criterion LCK n (x, y) contains useful information for jump detection. Our proposed jump detector LCK can be summarized as follows.
• At a given design point (x i , y j ), for i, j = 1, 2, . . . , n, solve the minimization problem (3) and obtain the solution vector { b(x i , y j ), c(x i , y j )}. • Divide O n (x i , y j ) into two halves U n (x i , y j ) and V n (x i , y j ) along the direction perpendicular to the estimated gradient vector { b(x i , y j ), c(x i , y j )}. Compute the Nadaraya-Watson LCK estimators of f (x, y) from the observations in U n (x i , y j ) and V n (x i , y j ), respectively, as described in (4). • The design point (x i , y j ) is flagged as a detected jump point if where Z 1−α n is the (1 − α n )th quantile of the standard normal distribution and α n is a significance level.
In practice, σ is often unknown, and it needs to be estimated from the observed data. To this end, it can be estimated by the conventional kernel smoothing approach in a chosen region in which the regression surface is relatively smooth. Or, σ 2 can be estimated by the residual mean squares of the jump-preserving surface estimation procedure suggested by Qiu (2004). Note that the jump detector LCK does not take into account the possible blurring in the observed data. As illustrated in Figure 1(b), blurring alters the image most significantly around a JLC (see its formal definition in Qiu 2005) and least significantly at places where the intensity surface is straight. So, in cases when blurring is present and (x, y) is on a JLC, if a design point in O n (x, y) is closer to the line that divides O n (x, y) into U n (x, y) and V n (x, y), then it is more likely that the corresponding observation Z ij has blurring involved. Thus, it should receive smaller weight in the related weighted average. To address this issue, we suggest using a univariate kernel function L for assigning such weights. By this idea, our second jump-detection criterion is defined by L ij := L(d ij /h n ), L is a univariate increasing density kernel function with support [0, 1], and d ij is the Euclidean distance from the point (x i , y j ) to the line that divides O n (x, y) into U n (x, y) and V n (x, y). The corresponding jump detector is labeled as LC2K, where LC2K denotes local constant smoothing with two kernel functions. The jump-detection criterion LC2K n (x, y) has similar asymptotic results to (5) and (6). Also, a numerical algorithm of the jump detector LC2K can be developed in a similar way to that of the jump detector LCK. In cases when the true regression surface f is quite steep but continuous around a given point (x, y), by the jump detectors LCK and LC2K, (x, y) can be falsely detected as a jump point because both their jump-detection criteria LCK n (x, y) and LC2K n (x, y) are good estimators of the directional derivative of f along the gradient direction at (x, y). To overcome this limitation, in the one-sided neighborhoods U n (x, y) and V n (x, y), instead of computing the Nadaraya-Watson LCK estimators (see (4)), we consider fitting local planes by solving the following LLK smoothing problems: and The solutions to a 0 of (8) and (9) are denoted as f LLK,+ (x, y) and f LLK,− (x, y), respectively, and they have the expressions for s 1 , s 2 = 0, 1, 2, and w ij (x, y) is defined to be the same as w ij (x, y) except that U n (x, y) needs to be replaced by V n (x, y) throughout. Then, we define our third jump-detection criterion by It can be checked that, in cases when f is continuous at (x, y), the value of |LLK n (x, y)| would be small, even when f is steep around (x, y), because the slope effect has been accommodated in fitting the local planes in (8) and (9). See, for instance, Qiu (2004) for a related discussion.
Obviously, the jump-detection criterion LLK n (x, y) has not taken into account the blurring effect, as discussed above about the jump detector LCK. By combining the ideas of the jumpdetection criteria LC2K n (x, y) and LLK n (x, y), we define the fourth jump-detection criterion by where f LL2K,+ (x, y) and f LL2K,− (x, y) are, respectively, the solutions to a 0 of the following local weighted least-square prob-lems: For the jump-detection criteria LLK n (x, y) and LL2K n (x, y), we can derive their asymptotic results similar to those in (5) and (6). Their corresponding jump-detection procedures are denoted as LLK and LL2K, respectively, and are briefly summarized as follows.
• At a given design point (x i , y j ), for i, j = 1, 2, . . . , n, solve the minimization problem (3) and obtain the solu- Compute f LLK,+ (x, y) and f LLK,− (x, y) (or, f LL2K,+ (x, y) and f LL2K,− (x, y)) from (8) and (9) (or (10) and (11)). • The design point (x i , y j ) is flagged as a detected jump point if the jump detector LLK is used and or, the jump detector LL2K is used and where σ is an appropriate estimator of σ .

Selection of the Bandwidth
The proposed jump detectors described in Section 2.1 have the bandwidth parameter h n involved. In one-dimensional cases, Gijbels and Goderniaux (2004) demonstrated that this parameter plays an important role in their one-dimensional jump-detection procedure. Based on our numerical experience, this is also true for the current two-dimensional procedure. In this subsection, we propose a procedure for choosing h n in two-dimensional cases.
The two-dimensional bandwidth selection problem is more complicated than its one-dimensional counterpart. To choose the bandwidth properly in two-dimensional cases, we first need to choose a metric for measuring the distance between the pointset S of the true jump points and its estimator S n by a jump detector with a given bandwidth h n . Qiu (2002) investigated this performance measure problem carefully, and suggested the following performance measure: where 0 ≤ ω ≤ 1 is a weight, denotes the entire design space, A \ B denotes the set of points in A but not in B, and |A| denotes the number of design points in the point set of A. Clearly, d Q is a weighted average of the false positive rate and the false negative rate of the related jump detector, and the weight ω reflects the relative importance of the two rates. In an application, if the two rates are equally important, then ω can be simply chosen as 0.5. It should be pointed out that d Q can be computed only when the pointset of the true jump points S is known. In mathematics, a popular metric for measuring the distance between two pointsets A and B is the following Hausdorff distance: where d E (s 1 , s 2 ) denotes the Euclidean distance between two points s 1 and s 2 . Qiu (2002) demonstrated that the above Hausdorff distance had two major limitations: (i) it was sensitive to individual points in the related pointsets, and (ii) it required extensive computations. As a comparison, the measure d Q is easier to compute and more robust to individual points. However, we find that d Q has its own limitations. For instance, assume that the true JLC is the line y = 0.5 and the detected jumps are those design points on the line y = 0.5 + 1/n 1 , where n 1 is a positive integer. Then, we would expect that the related jump detector performs better when n 1 gets larger. However, it can be checked that the value of d Q does not depend on n 1 in such a case, which is intuitively unreasonable. To avoid this drawback, in this article, we propose the following performance measure: where 0 ≤ ω ≤ 1 is a weight, and d E ((x , y ), S) denotes the Euclidean distance from the point (x , y ) to the pointset S, that is, the minimum of all pairwise distances between (x , y ) and (w, z) for (w, z) in S. Obviously, d KQ is a weighted average of two averages, where the first average in (12) is the average Euclidean distance from the individual detected jump points to the set of true jump points and the second average is the average Euclidean distance from the individual true jump points to the set of detected jump points. By using d KQ , the drawback of the measure d Q mentioned above is obviously overcome.
To demonstrate the strengths and limitations of the three performance measures discussed above, let us consider the following two examples. Without loss of generality, in both examples, ω is fixed at 0.5 and the sample size is fixed at n = 100. In the first example, it is assumed that there is a single JLC in the design space [0, 1] × [0, 1], and it is the line parallel to the x-axis at y = 0.5. We further assume that the detected jump points include all design points with y = 0.5 and a point at (0.5, 0.5 + r), where 0 ≤ r ≤ 0.5 is a constant. These detected jump points are shown in Figure 2 In this example, it is obvious that d H ( S n , S) does not reflect the jump detection performance well when r is large, and the other two performance measures are more reasonable to use. Next, in the second example, let us assume that the detected jump points are those design points on the line at y = 0.5 + t (shown in Figure 2(c)), where 0 ≤ t ≤ 0.5 is a constant, and the remaining setup is the same as those in the first example. In this example, it can be checked that the three performance measures are When t changes from 0 to 0.5, these performance measures are shown in Figure 2( However, we should point out that d KQ is not a mathematical metric and it is proposed for selecting bandwidth due to the reason that the mathematically well-defined Hausdorff distance is highly sensitive to individual points as demonstrated in Figure 2. Thus, the Hausdorff distance still needs to be adopted for characterizing asymptotic properties of jump detectors in Section 3, but d KQ is used in Section 4 to evaluate numerical performance in finite-sample scenarios. In practice, S can be replaced by its discrete version S * = {(x i , y j ) : d((x i , y j ), S) ≤ 1/(2n)} for the purpose of calculating d KQ . In simulations, the point set S is usually known; thus, h n can be chosen by minimizing d KQ ( S n , S; h n ). In practice, however, S is often unknown. In such cases, we propose a bootstrap procedure described below for choosing h n .
2.2.1 Bootstrap Procedure For Choosing h n .
1. Apply a jump detector to the original dataset {(x i , y j , Z ij ), i, j = 1, 2, . . . , n}. The set of detected jump points is denoted as S n . 2. Construct an estimator G{f }(x, y) of the blurred regression function G{f }(x, y) using a jump-preserving surface estimation method, and define the residuals by r ij = Z ij − G{f }(x i , y j ). 3. Draw with replacement a bootstrap sample of size n 2 , denoted as { r * ij , i, j = 1, 2, . . . , n}, from the residual set { r ij , i, j = 1, 2, . . . , n}, and then obtain a bootstrap dataset

The jump detector with the same bandwidth as that in
Step 1 is then applied to the bootstrap dataset, and a set of detected jump points from this bootstrap dataset can be obtained. 5. Steps 3 and 4 are repeated N times. The sets of the detected jump points from the N bootstrap datasets are denoted as S * 1n , S * 2n , . . . , S * Nn , respectively. Then, the bootstrap estimator of d KQ ( S n , S; h n ) is defined by The optimal bandwidth is then approximated by the minimizer of d BT KQ ( S n , S; h n ) with respect to h n . In Step 2 of the above procedure, G{f }(x, y) can be constructed by any jump-preserving surface estimator proposed in the literature (see Hillebrand and Müller 2007;Qiu 2007). In all numerical examples presented in this article, the local piecewise constant kernel estimator suggested in Qiu (2009) is used.

STATISTICAL PROPERTIES
We discuss some statistical properties of the proposed jump detectors. In the literature, the psf g (u, v; x, y) in model (2) is often assumed to be a density function (i.e., a nonnegative function with a unit integration on its support), since it is believed that the blurring process does not change the mass (Bates and McDonnell 1986) of the regression surface. This conventional assumption is also adopted here. Further, we assume g is symmetric about both x-axis and y-axis, which in practice is often satisfied by common types of blur such as out-of-focus blur, motion blur, or Gaussian blur. More specifically, we assume that, for any (x, y) ∈ , (i) g (u, v; x, y) ≥ 0, for all (u, v) ∈ R 2 , (ii) R 2 g (u, v; x, y) , y), J S ) ≤ h n }, and J C S,h n = \ J S,h n . From Theorem 3.1, each of the four proposed jump detectors provides a consistent estimator of S in the sense that it converges almost surely to S in the Hausdorff distance, after certain small regions around the singular points in S and around the border of the design space are excluded. It should be pointed out that the assumption used in the theorem that the blurring extent ρ n converges to 0 when n increases seems restrictive; but, it is actually quite flexible for the reason given below. Let r n = nρ n . Then, r n denotes the blurring extent in terms of the number of rows/columns of design points. The assumptions in Theorem 3.1 allows r n to increase to infinity when n increases. Namely, when the resolution of the observed regression surface increases, our proposed methods allow the number of rows/columns of design points involved in the blurring at a given point to increase at a certain rate, which should be flexible enough for most applications. The proof of Theorem 3.1 is given in the supplementary material.

NUMERICAL STUDIES
In this section, we present some numerical examples concerning the performance of the four proposed jump detectors. Throughout this section, if there is no further specification, the kernel function K used in (3) is chosen to be the truncated two-dimensional Gaussian density function 1/(2π − 3π exp(−0.5))[exp(−(x 2 + y 2 )/2 − exp(−0.5))]I x 2 +y 2 ≤1 , the two-dimensional kernel function used in constructing each Table 1. Numerical results for example 1 based on 100 replicated simulations. Each entry is (optimal bandwidth, bandwidth chosen by the proposed bootstrap procedure), and the averaged value of d KQ ( S n , S; h n ) when h n is chosen to be the optimal bandwidth  jump-detection criterion (see e.g., (4)) is chosen to be the Epanechnikov kernel function (2/π )(1 − x 2 − y 2 )I x 2 +y 2 ≤1 , and the one-dimensional kernel function L used in (7) and (10) is chosen to be (1/1.194958) exp(x 2 /2)I 0≤x≤1 , which is proportional to the reciprocal of the truncated one-dimensional Gaussian kernel function. We choose these kernel functions because the Epanechnikov kernel function is a standard choice in the statistical literature and the truncated Gaussian kernel functions are commonly used in the image processing literature.
In model (2), the psf is chosen to be g (u, v; x, y) = 3 πρ 2 n (x, y) and the additive random errors ε ij are generated from the distribution N (0, σ 2 ). The above psf is circularly symmetric at each (x, y) with the blurring extent controlled by ρ n (x, y). In this example, we simply choose ρ n (x, y) ≡ 0.02. In such cases, g (u, v; x, y) is location invariant. Next, we try to measure the performance of the related jump detectors quantitatively using the proposed performance measure d KQ ( S n , S; h n ), in which ω is fixed at 0.5 in this section. Without loss of generality, the bandwidth h n is assumed to take the value of k/n, where k < n is a positive integer. We consider two sample sizes 100 × 100 (i.e., because it provides a reasonable approximation to the minimizer of min h n E(d KQ ( S n , S; h n )). Average bandwidths computed over 100 replications chosen by our proposed bootstrap procedure are also presented in the table.
From Table 1, it can be seen that: (i) the jump detectors based on the LLK smoothing (i.e., LLK and LL2K) perform better than their counterparts based on the LCK smoothing (i.e., LCK and LC2K), (ii) generally speaking, bandwidth for each jump detector should be chosen larger for noisier data, (iii) performance of the jump detectors based on the LCK smoothing does not change much when noise level or sample size changes, (iv) the jump detectors based on the LLK smoothing perform better when the sample size gets larger, (v) the bootstrap procedure tends to select the bandwidth slightly larger than the optimal value, and (vi) the bootstrap procedure gives close-tooptimal bandwidths when the sample size increases. The result (i) is intuitively reasonable because the regression surface is steep in certain regions in this example (see Figure 3(a)), the jump detectors based on the LLK smoothing can accommodate such a case well, but the jump detectors based on the LCK smoothing would have many false jump detections due to the steep regression surface. The result (iii) is also caused by the false jump detections of the jump detectors based on the LCK smoothing.
In cases when n = 100 and σ = 0.3, the detected jumps by the four jump detectors with their optimal bandwidths presented in Table 1 are shown in Figure 4, where the results are from the simulation whose d KQ ( S n , S; h n ) value is the median of the 100 d KQ ( S n , S; h n ) values computed from the 100 replications. From the figure, it can be seen that (i) the detected jumps in both plots (a) and (b) contain some false detections in the region where the true regression surface is steep, (ii) both the jump detectors LCK and LC2K fail to detect some true jump points, and (iii) the detected jumps in plot (d) seem to be a little more variable than those in plot (c). All these results are consistent with those in Table 1 and with our theoretical justifications discussed before. This example demonstrates that the jump detectors LLK and LL2K can accommodate the slope effect of the true regression surface well, while the jump detectors LCK and LC2K cannot.
Next, we consider another example in which the true regression function is if 0.2 < x ≤ 0.7 (x − 0.7) 2 + 1.25, if 0.7 < x ≤ 1, and the psf is defined by (13) with ρ n (x, y) = 0.1(1 − y). Therefore, the true regression function f 2 has two JLCs, both are parallel to the y-axis at x = 0.2 and x = 0.7, and the jump size is larger (i.e., 2.0) for x ≤ 0.2 and smaller (i.e., 1.0) for x > 0.7. The psf is location-dependent in this example, with blurring extent increasing when y gets smaller. The blurred version of f 2 is shown in Figure 5(a). A noisy and blurred version of f 2 with n = 100 and σ = 0.3 is shown in Figure 5(b).
In the same format as Table 1, the optimal bandwidths, the corresponding d KQ ( S n , S; h n ) values, and the bandwidths chosen by the bootstrap procedure, based on 100 replicated simulations, are presented in Table 2. From the table, it can be seen that (i) the four jump detectors all perform better when the sample size increases, (ii) the jump detectors LLK and LL2K are more sensitive to the noise level in the sense that their d KQ ( S n , S; h n ) values increase more rapidly as the noise level increases, compared to the jump detectors LCK and LC2K, (iii) for all four jump detectors, a larger bandwidth is chosen when the noise level is larger, and (iv) the bandwidths chosen by the bootstrap procedure are all quite close to the optimal bandwidths. The result (ii) is consistent with the theoretical result that the LLK estimators would have a larger variability than the LCK estimators, if their bandwidths are similar (see Qiu 2005, chap. 2). The results (i) and (iii) are intuitively reasonable.
Simulation results shown in Figures 5(c)-5(f) have the same setup as that in Figure 4. From the plots, it can be seen that (i) all four jump detectors detect the larger jumps at x = 0.2 reasonably well, (ii) the jump detectors LCK and LC2K that are based on the LCK smoothing detectors perform better than the jump detectors LLK and LL2K that are based on the LLK smoothing when detecting the smaller jumps at x = 0.7, (iii) by comparing plots (c) and (d), it seems that the jump detector LC2K is slightly more robust to the blurring than the jump detector LCK, and (iv) by comparing plots (e) and (f), the jump detector LL2K seems more robust to the blurring than the jump detector LLK. This example demonstrates that (i) the jump detectors LC2K Figure 7. A test image with spatially varying blur involved.  Figure 7 by the jump detectors LCK, LC2K, LLK, and LL2K, respectively. and LL2K using two kernel functions in their construction seem more robust to the blurring, compared to the other two jump detectors, and (ii) the results by the jump detectors LCK and LC2K seem less variable, compared to the jump detectors LLK and LL2K.
Next, we consider a real test image shown in Figure 6(a), which is synthetic aperture radar (SAR) image of an area near the Thetford forest in England. This image has 250 × 250 pixels with gray levels in the range [0, 255]. Its blurred and noisy version by the psf in (13) with ρ n (x, y) ≡ 0.02 and by the iid noise from N (0, 5 2 ) is shown in Figure 6(b). The detected jumps by the four jump detectors LCK, LC2K, LLK, and LL2K are presented in Figures 6(c)-6(f), respectively. In all four jump detectors, α n is fixed at 0.001, and their bandwidths are chosen by the bootstrap procedure to be 3/250, 3/250, 11/250, and 9/250. From the plots, it can be seen that (i) all four jump detectors detect the major JLCs reasonably well, (ii) the JLCs detected by LCK and LC2K are thinner and more straight than those detected by LLK and LL2K; but, some false jumps are detected by them here and there, (iii) the jump detector LLK fails to detect certain jumps (e.g., the jumps located in the upperright portion of the image), and (iv) the jump detector LL2K detects more true jumps than the jump detector LLK, which is consistent with our analysis in Section 2. Given that the noise level in the SAR image is relatively high, we suggest using jump detectors based on the LCK smoothing (i.e., LCK or LC2K) in this example.
We conclude this section by considering an image shown in Figure 7. From the figure, it can be seen that objects within the focus (e.g., the pen) are quite sharp but objects in distance (e.g., the child's face) are blurred. So, this image provides us an example with spatially varying blur. By the jump-preserving image restoration procedure suggested by Qiu (2004), an estimate of σ is computed to be 4.558. The detected jump points by LCK, LC2K, LLK, and LL2K are shown in Figures 8(a)-8(d), respectively. From the plots, it can be seen that both LCK and LLK fail to detect some jumps in the blurred region, which is especially true for LLK. The detected jump points by LC2K and LL2K look reasonably good. Because the blur is quite severe in this example, we suggest using the jump detector LC2K or the jump detector LL2K.

CONCLUSIONS
We have presented four jump detectors for detecting jumps in blurred regression surfaces. All four jump detectors are based on local constant or LLK smoothing, two of which (i.e., LC2K and LL2K) take the possible blurring in the observed data into consideration while the other two do not. Also, a new quantitative metric for measuring the performance of a jump detector is proposed, which overcomes the major limitations of the Hausdorff distance and the distance metric proposed by Qiu (2002). A data-driven bandwidth selection procedure via bootstrap is suggested as well in the article.
Numerical examples presented in the previous section show that these jump detectors have their own strengths and limitations when handling various different situations. These results can be summarized as follows. (i) The jump detectors based on LCK smoothing (i.e., LCK and LC2K) are more robust to noise and their detected JLCs tend to be thinner. However, they cannot accommodate the slope effect of the true regression surface well and would detect many false jumps at places where the true regression surface is steep. (ii) The jump detectors based on LLK smoothing (i.e., LLK and LL2K) can accommodate the slope effect of the true regression surface well; but they are more sensitive to noise. (iii) The jump detectors using two kernel functions (i.e., LC2K and LL2K) are more robust to blurring; however, their detected jumps tend to have a larger variability. Based on these results, we provide the following guidelines for the use of these jump detectors in practice. (i) The two jump detectors based on the LCK smoothing (i.e., LCK and LC2K) are recommended in cases when the noise level in the observed image is relatively high. (ii) The two jump detectors based on the LLK smoothing (i.e., LLK and LL2K) should be used in cases when the observed image contains some parts at which the image intensity surface is steep. (iii) The jump detectors using two kernel functions in their construction (i.e., LC2K and LL2K) are recommended in cases when the observed image has a quite severe blur involved. Much future research is needed to find an appropriate way to combine the major strengths of the four jump detectors and avoid their major limitations at the same time.

SUPPLEMENTARY MATERIALS
Supplemental.pdf: This pdf file contains proofs for Theorem 1. ComputerCodesAndData.zip: This zip file contains Fortran source code of our proposed jump-detection methods and the two test images used in the article.