Edge detection of remote sensing image based on Grünwald-Letnikov fractional difference and Otsu threshold

1.   Introduction

Remote sensing images include all kinds of geomorphic information, such as roads, farmland and buildings. The edge detection of these different geomorphic features is of great significance for topographic survey, urban planning, precision agriculture and military analysis. However, with the development of remote sensing technology, the resolution of remote sensing images is improving, the presentation of geomorphic information is becoming more and more abundant, and the difficulty of identifying and extracting edge information is also increasing.

In image processing, edges are the collections of pixels between adjacent areas of an image. They reflect the discontinuity of image features and contain a lot of structural information. Edge detection is an important part of image processing as well as a basic means of image segmentation ^[1], object detection ^[2] and information extraction ^[3]. Currently, the methods used for remote sensing image edge detection can be roughly divided into three categories. Fist, there is edge detection based on gradient: C. Xu et al. proposed multispectral image edge detection via Clifford gradient ^[4], S. Amstutz et al. proposed edge detection using topological gradients ^[5], J. T. Tang et al. proposed edge detection based on curvature of the gravity gradient tensor ^[6], and V. B. S. Prasath et al. proposed multiscale gradient maps augmented fisher information ^[7]. Second, there is edge detection based on transformation: H. H. Zhao et al. proposed optimal Gabor filter–based edge detection of high spatial resolution remotely sensed images ^[8], W. C. Lin et al. proposed edge detection in medical images with quasi high–pass filter ^[9], G. B. Chen et al. proposed a road identification algorithm for remote sensing images based on wavelet transform and a recursive operator ^[10], and A. Isar and C. Nafornita et al. proposed robust edge detection base on Hyperanalytic wavelet ^[11]. Third, there is edge detection based on machine learning: M. Han et al. proposed an extreme learning machine based on cellular automata of edge detection for remote sensing images ^[12], L. Huang et al. proposed edge detection in UAV remote sensing images using a method integrating Zernike moments with clustering algorithms ^[13], Z. Qu et al. proposed image fusion using deep neural networks for image edge detection ^[14], and X. G. Zheng proposed a robot image edge detection method based on a Gaussian positive–negative radial basis neural network ^[15]. The first category edge detectors use mostly differentiation to detect the discontinuity of gray value. Commonly used first order differential detectors mainly include the Sobel detector ^[16] and the Canny detector ^[17], and second order differential detectors include the Laplacian detector ^[18]and the Laplacian of Gaussian (LoG) detector ^[19]. These algorithms are easy to calculate, but they have some disadvantages such as discontinuous edge extraction and poor anti–noise performance. In view of the problems of these detectors, many researchers have made corresponding improvements. ^[20] introduced a quantum improved Sobel edge detection algorithm with non–maximum suppression and double threshold techniques for the novel enhanced quantum representation method. Based on Canny and LoG detectors, ^[21] developed a new efficient edge detector (CLoG) used for the detection of edge areas in digital images. The detection effect of these improved algorithms is indeed improved, but the detection results tend to mistakenly detect fake edges in the presence of excessive noise or artifacts ^[22].

In recent years, there has been the application of fractional differentials in image enhancement ^[23], image data enhancement ^[24], image denoising ^[25], image recovery ^[26], image encryption ^[27]and data hiding ^[28], extracting image edge features ^[29]. More and more widely, the application of fractional differentials in image edge detection has gradually become an important research topic. Because a fractional difference can maintain the low–frequency contour features of the smooth region of the image and enhance the high–frequency edge and texture details of the region with large or insignificant grayscale changes ^[25], it has achieved good results in edge detection of medical images ^[30], grayscale images ^[31], indexed images ^[32], true color images ^[33] and pseudo–color ^[34]. However, its application in the edge detection of remote sensing images is almost absent and needs further research. In order to extract rich detail features, we try to apply the fractional differential method to edge detection of remote sensing images.

Table 1 shows the comparison of some fractional difference mask operators in the last three years. There are three classical definitions of fractional derivatives: Riemann–Liouville (R–L)^[35], Caputo–Fabrizio(C–F)^[36,37] and Grünwald–Letnikov (G–L)^[38], which are the theoretical basis for the construction of masks. Considering that existing mask operators of fractional difference basically only pay attention to changes in two directions ($X$–direction and $Y$–direction), we try to construct mask operators in four directions for the edge detection research of remote sensing images. In addition, to better describe the variation, we also add a regulating parameter to the mask operator. Therefore, the main contributions of this paper are as follows:

1) A fractional difference mask operator with two parameters is constructed in four directions($X$-direction, $Y$-direction, diagonal direction and inverse diagonal direction) based on G-L fractional derivative.

2) Based on fractional difference mask and Otsu threshold, edge detection of remote sensing images from the NWPU VHR-10 dataset ^[40] is carried out.

3) The performance of the proposed method is compared with other detectors from the quantitative aspects of binary graph similarity (BGS) and edge pixel ratio (EPR).

4) The influence of two parameters and threshold on edge detection results is discussed.

The rest of the paper is organized as follows. Section 2 describes the related mathematical theory in detail (Grünwald–Letnikov fractional difference, fractional convolution mask and Otsu threshold). Section 3 introduces the operation procedure of the edge detection algorithm, sources of remote sensing images and performance measures for evaluating edge detection operation. The simulation results and discussions are illustrated in Section 4. Finally, conclusion and future scope is summarized in Section 5.

2.   Theoretical basis

2.1.   G–L fractional difference

The G–L fractional derivative is defined as

where $\gamma$ is the order of fractional differentiation, $\gamma \in R$ (and may be a fraction), $h$ is the number of steps, and $\Gamma(k) = (k-1)!$. So, the fractional partial derivative of $f(x, y)$ is defined as

When $h = 1$, the G–L fractional difference can be expressed as ^[41]

Further, the difference of fractional partial derivative can be written as

When $\gamma = 1$, the differences of the first partial derivative are

When $\gamma = 2$, the differences of the second partial derivative are

2.2.   Fractional convolution mask

The integer–order mask operators in horizontal direction and inverse diagonal direction are

and

where $\rho$ is called the regulating parameter and is a natural number.

The differential operator $\chi = \rho (f(x+1, y)-f(x-1, y))$ is constructed with $IOMH$. We consider using the two smallest difference values to approximate the gradient value and construct the fractional differential mask operators. Because $f(x+1, y)-f(x-1, y) = 2\frac{\partial f(x+1, y)}{\partial x}$,

Then, the fractional convolution mask in the $X$–direction is

Rotate the $FCMX$ symmetrically to obtain the fractional convolution mask in the $Y$–direction, that is,

The differential operator $\chi = \rho (f(x+1, y-1)-f(x-1, y+1))$ is constructed with $IOMI$. For symmetry, transform $\chi$ to

Because $f(x+1, y+1)-f(x+1, y-1) = 2\frac{\partial f(x+1, y+1)}{\partial y}$, $f(x+1, y+1)-f(x-1, y+1) = 2\frac{\partial f(x+1, y+1)}{\partial x}$, $f(x+1, y-1)-f(x-1, y-1) = 2\frac{\partial f(x+1, y-1)}{\partial x}$, $f(x-1, y+1)-f(x-1, y-1) = 2\frac{\partial f(x-1, y+1)}{\partial y}$,

Then,

and

The fractional convolution mask in the inverse diagonal direction is

Rotate the $FCMI$ symmetrically to obtain the fractional convolution mask in the diagonal direction, that is,

Add $FCMX$, $FCMY$, $FCMI$ and $FCMD$ to obtain the final fractional convolution mask:

2.3.   Otsu threshold

The basic idea of the maximum interclass variance method ^[42,43,44] is to divide the image into background and foreground and calculate the best threshold to distinguish the two pixel regions, so as to maximize the distinction between the two types of pixels. Assume that the gray–scale image $G(x, y)$ has $\kappa$ gray-scale levels, and $N_m (m = 1, \cdots, \kappa)$ represents the pixel number in each gray level. The pixel number of $G(x, y)$ is $N = \sum\limits_{m = 1}^\kappa N_m$, so the occurrence probability of gray-scale level $m$ is $P_m = N_m/N$.

If the threshold $\varpi$ is used as the limit to divide the image into two regions, the foreground pixel region is $A$, and the background pixel region is $B$, then the probabilities of foreground $A$ and background $B$ are, respectively, $P_A = \sum\limits_{m = 1}^\varpi P_m$, $P_B = \sum\limits_{m = \varpi+1}^\kappa P_m$.

Let the total gray mean of $G(x, y)$ be $s_0 = \sum\limits_{m = 1}^\kappa mP_m$, and then the gray means of $A$ and $B$ are, respectively, $s_A = \sum\limits_{m = 1}^\varpi mP_m$, $s_B = \sum\limits_{m = \varpi+1}^\kappa mP_m$. The interclass variance of $A$ and $B$ is

When $\varepsilon(\varpi)$ is the maximum value, $\varpi^*$ is the optimal threshold, that is, Otsu threshold.

3.   Preparation for experiment

3.1.   Edge detection algorithm based on FCM and Otsu

Suppose the remote sensing image is $RSI(R, G, B)$, and the process of edge detection using the FCM and Otsu threshold is as follows:

Step 1. $RSI(R, G, B)$ is preprocessed, and gray–scale image $G(x, y)$ is obtained.

Step 2. Image $E(x, y)$ is obtained by convolution of image $G(x, y)$ with $FCM$, using the appropriate values of $\rho$ and $\gamma$.

Step 3. By calculating the Otsu threshold $\varpi^*$ of image $G(x, y)$, image $E(x, y)$ is binarized to obtain the edge detection image $F(x, y)$.

In this process, the selection of parameters $\rho$ and $\gamma$ determines the effect of edge detection, so they need to be adjusted for the best effect in the experiment.

3.2.   Data source

The remote sensing images used in this study are all from the NWPU VHR-10 dataset. This very–high–resolution (VHR) remote sensing image dataset was constructed by ^[40] from Northwestern Polytechnical University (NWPU). This dataset contains a total of $800$ VHR remote sensing images, which were cropped from Google Earth and the Vaihingen data set. The images cover 10 geographical scenes, including airplane, ship, storage tank, baseball diamond, tennis court, basketball court, ground track field, harbor, bridge and vehicle. Several of them were selected for edge detection research in this paper.

3.3.   Performance indicators

In order to evaluate the performance ^[45] of the edge detection operator in this paper, we choose the following indicators:

1) Suppose $B_1(x, y)$ and $B_2(x, y)$ are $M\times N$ binary graphs, $C = \{(i, j)\in(M\times N)\mid B_1(i, j) = B_2(i, j)\}$, and $num(C)$ is the number of elements in set $C$. The BGS is calculated as follows:

2) Suppose $B(x, y)$ is an $M\times N$ binary graph, and $S$ is the number of edge pixels. The EPR is calculated as follows:

4.   Results and analysis

4.1.   Comparison of edge detectors

We used the proposed algorithm to detect the edge of an image in the NWPU VHR-10 dataset, and we compared the results with some edge detectors: Roberts, Sobel, LoG, Laplacian, Canny, ^[33] and ^[34]. The original image $I_1$, Otsu image and the edge images extracted by the seven detectors are shown in Figure 1. $I_1$ has a pixel size of $1073\times 705$ and Otsu threshold of $0.5725$. In the FCM image, $\rho = 4.5, \gamma = 0.7$.

Compared with the other seven detectors, the edge image extracted by the FCM detector has richer texture layers and more obvious details. In order to compare these eight detectors quantitatively, we calculated the BGS value, EPR value and required time, respectively, as shown in Table 2.

As can be seen from Table 2, the FCM detector extracts more edge information, and the edge image extracted by the FCM detector is more similar to the Otsu image, but the time spent is moderate.

4.2.   Influence of threshold

In order to investigate the influence of thresholds, for the original image $I_2$, we made edge detection images extracted by the FCM detector with different thresholds, as shown in Figure 2. $I_2$ has a pixel size of $1259 \times 820$ and Otsu threshold of $0.4941$. In the FCM image, $\rho = 3, \gamma = 0.5$. The BGS values of edge images extracted by the FCM detector with different thresholds relative to the Otsu image and the EPR values of edge images extracted by FCM detector with different thresholds are shown in Table 3.

The difference presented by the images is not obvious, but it can be seen from Table 3 that BGS increases, and EPR decreases with the change of thresholds. The method proposed in this paper uses Otsu to determine the threshold, avoids the subjectivity of artificial assignment, and achieves a balance between BGS and EPR.

4.3.   Influence of parameters $\rho$ and $\gamma$

In order to investigate the influence of the two parameters $\rho$ and $\gamma$ in the FCM detector, edge images corresponding to $\rho = 3$ and different $\gamma$ values are presented in Figure 3, and edge images corresponding to different $\rho$ values and $\gamma = 0.6$ are presented in Figure 4. $I_3$ has a pixel size of $974 \times 762$ and Otsu threshold of $0.4196$. Figures 3 and 4 show that when $\rho$ is constant and $\gamma$ is in the interval $[0, 1]$ or when $\gamma$ is constant and $\rho$ is larger, more background details and texture information are extracted.

The BGS and EPR values obtained by different values of the parameters $\rho$ and $\gamma$ are listed in Tables 4 and 5, respectively. As can be seen from Tables 4 and 5, if $\gamma$ remains unchanged, both BGS and EPR values increase with the increase of $\rho$. If $\rho$ is unchanged, the BGS and EPR values increase as $\gamma$ increases from $-0.5$ to $0.5$, and they decrease as $\gamma$ increases from $0.5$ to $1.7$. These features are consistent with the rules shown in Figures 3 and 4, so the values of BGS and EPR can reflect the extraction of edge information to a certain extent. The larger BGS or EPR is, the richer the extracted edge information.

When $\rho = 3$, the line graphs of BGS and EPR with respect to $\gamma$ on the interval $[-1.5, 2.5]$ are shown in Figure 5. Both BGS and EPR are symmetric with respect to $\gamma = 0.5$. When $\gamma = 0.5$, the line graphs of BGS and EPR with respect to $\rho$ are shown in Figure 6. As $\rho$ increases, both BGS and EPR increase.

Figure 7 presents the surface diagrams of BGS and EPR with respect to two parameters and more intuitively shows how BGS and EPR change as the parameters change. Obviously, the shape of the two surface diagrams is very similar, that is, the changes of BGS and EPR with respect to the two parameters are almost the same.

5.   Conclusions

Based on the Grünwald–Letnikov fractional difference formula, a method for constructing a fractional convolution mask is proposed. Then, the method is combined with the Otsu threshold for edge detection of remote sensing images. We performed both visual and quantitative comparison (BGS and EPR) with existing edge detectors and found that the proposed method can better extract edge information of remote sensing images, reflect more texture details and does not need smoothing processing. However, this study also has some shortcomings: First, the size of the mask is only $3\times 3$, and other sizes of mask structures are not considered; second, the detection effect of the image with noise is not good. Haze and clouds seriously reduce the visibility of remote sensing images and are the main source of noise in remote sensing images, so our next work is to remove haze ^[46] and clouds ^[47] from remote sensing images in combination with other mathematical theories, and then better carry out edge detection or target recognition.

Acknowledgments

We are thankful to anonymous reviewers and my colleagues (Ju Wu and Fang liu) for their help. This work is funded by the National Natural Science Foundation of China (61872304); National Natural Science Foundation of China (11502121); The Project of State Administration of Science, Technology and Industry for National Defence (20zg6108); The Application Basic Research Plan Project of Sichuan Province (No. 2021JY0108); The Science and Technology Research Plan Project of Sichuan Province (No. NJFH20-003).

Conflict of interest

All authors contributed to the study conception and design. Material preparation, data collection and analysis were performed by Bin Wu and Hua Kong. The first draft of the manuscript was written by Chao Chen, and all authors commented on previous versions of the manuscript. All authors read and approved the final manuscript. The authors declare there is no conflict of interest.