Remote Sensing Image Change Detection Based on NSCT-HMT Model and Its Application

Traditional image change detection based on a non-subsampled contourlet transform always ignores the neighborhood information’s relationship to the non-subsampled contourlet coefficients, and the detection results are susceptible to noise interference. To address these disadvantages, we propose a denoising method based on the non-subsampled contourlet transform domain that uses the Hidden Markov Tree model (NSCT-HMT) for change detection of remote sensing images. First, the ENVI software is used to calibrate the original remote sensing images. After that, the mean-ratio operation is adopted to obtain the difference image that will be denoised by the NSCT-HMT model. Then, using the Fuzzy Local Information C-means (FLICM) algorithm, the difference image is divided into the change area and unchanged area. The proposed algorithm is applied to a real remote sensing data set. The application results show that the proposed algorithm can effectively suppress clutter noise, and retain more detailed information from the original images. The proposed algorithm has higher detection accuracy than the Markov Random Field-Fuzzy C-means (MRF-FCM), the non-subsampled contourlet transform-Fuzzy C-means clustering (NSCT-FCM), the pointwise approach and graph theory (PA-GT), and the Principal Component Analysis-Nonlocal Means (PCA-NLM) denosing algorithm. Finally, the five algorithms are used to detect the southern boundary of the Gurbantunggut Desert in Xinjiang Uygur Autonomous Region of China, and the results show that the proposed algorithm has the best effect on real remote sensing image change detection.


Introduction
Image change detection is the process of identifying changes in land cover through analyzing remote sensing images acquired in the same geographical location at different times [1]. It is widely used in the fields of video surveillance [2], medical diagnosis [3], land use [4], and natural disaster detection [5]. Changes in land use and oasis coverage have a direct impact on the human environment and ecological processes. Information regarding changes in land use and coverage are important for natural resource management and scientific decision-making [6,7], and are a critical component of research in the areas of the environment, forestry, hydrology, agriculture, geography, and ecology. Using remote sensing images for land use and oasis cover change detection can quickly and accurately obtain change information and development trends. However, the rapid updating of land cover data

Generation of the Difference Image
The ways to generate a difference map include the difference method and ratio method, and both can effectively construct difference images from optical remote sensing images, but the difference method is sensitive to image quality and the spectral characteristics of the objective conditions, resulting in a high rate of missed detection and false positives. The ratio method includes the logarithmic ratio method and mean ratio method. Among them, the advantage of the logarithmic ratio method is the ability to convert multiplicative speckle noise to additive noise, and the background information of the difference image obtained by the logarithmic transformation is relatively flat. The disadvantages are that the logarithmic ratio algorithm compresses the variation range of the difference image, and has the characteristics of enhancing the low intensity pixels to weaken the high intensity pixels, and cannot reflect the real change trends to the maximum extent [23]. The mean ratio method, on the other hand, does not have such problems. The mean value rule can effectively enhance the contour of the change area and changes in a small area, and can also prevent the loss of change information. In order to take full account of the relevant context, the difference map constructed by the mean value ratio depends only on relative changes in image intensity. It can truly reflect the changes in the region and retain more details. Therefore, in this paper, we use the mean ratio method proposed by Inglada et al. to generate the difference image [24].
Assume F and G are images of the same location taken at different times and the size of these images are r pixels × c pixels, and they have been calibrated by ENVI.
The mean-ratio operation formula is as follows: where µ 1 (i, j) and µ 2 (i, j) respectively represent the local mean value of F and G at the point (i, j) in the window size of ω × ω, ω = 1, 3, 5, 7, ...

A Denoising Method Based on the NSCT-HMT Model
After obtaining the mean difference map, this paper uses the non-sampling contourlet HMT (NSCT-HMT) model proposed by Wang et al. to denoise the difference image [25]. The original NSCT-HMT model is mainly used to process the optical images. In this paper, we apply it to remote sensing image denoising. The overview of the basic rationale of the adopted methodologies: 1.
The mean difference map was dealt with the NSCT decomposition; 2.
After the NSCT decomposition, the coefficients of NSCT are fitted by Gauss mixture model; 3.
Using the HMT model to estimate the parameters; 4.
The Monte-Carlo method is used to generate random white noise images and balance the variance of the image in the NSCT domain; 5.
The EM algorithm is used to estimate the NSCT coefficient of the denoised image; 6.
Using the NSCT synthesis to get denoised image.

NSCT Transform and the Hidden Markov Model
Da et al. proposed a kind of non-sampling contourlet transform (NSCT) that was based on the contourlet transform [26]. The NSCT does not have a sampling process in the image decomposition and reconstruction process, and it is composed of two-dimensional non-sampled pyramid filter banks and non-down-sampling filter banks. The result of this treatment is that the NSCT transform has translation invariance in addition to being multi-resolution, multi-scale, and anisotropic. Although the NSCT coefficients do not obey a Gauss distribution, they follow a zero mean mixed Gauss distribution. The NSCT coefficients can be modeled using a mixture of two Gauss functions. They have the characteristics of persistence and aggregation, and have a strong dependence on the neighborhood The hidden Markov model (HMT) is a double stochastic process consisting of a hidden Markov chain with a certain number of states and a set of random functions. It can be divided into a hidden layer and an observation layer. The hidden layer cannot be directly observed, and includes a certain number of states. The observation layer (actual observation) can have some observations corresponding to the state. The observation vector is generated by a sequence of states with probability density distribution. It is used to express various states by a probability density distribution.

Establishment of NSCT-HMT Model
In theory, a signal needs to be accurately represented by a series of Gaussian distributions. However, in the actual image processing application, it is possible to use the 2~4 Gaussian mixture to approximate the image [27,28]. The coefficients of NSCT can be modeled using a mixture of two Gauss functions. Formulas (2) and (3) are used to build the coefficients of the Gauss mixture model [25]: where N i represents the coefficient of a certain subband in each layer. s i indicates in which state the coefficients are selected. fs i (m) is the probability distribution function of the large and small states, and ∑ n m=1 f s i (m) = 1; f N i |s i (N i |s i = m) represents the probability distribution density function corresponding to the Gauss model when the coefficients are in a specific state. µ im and σ 2 im represent the mean and variance of the Gauss distribution, respectively.
The probability distribution function of N 1 in each direction sub-band node can be expressed as P N 1 (m), where m is the number of hidden states. State transition probabilities can be expressed as The above parameters are represented by a formal θ N−HMT of four tuples [25]: where P represents the number of coefficients and M represents the number of states. The NSCT-HMT model can be represented by the four tuple of parameters as shown in Formula (4).

Image Denoising Based on NSCT-HMT Model
In order to calculate of the model parameters, we use the HMT-EM algorithm proposed by Fan et al. to estimate the parameters [29], and the local optimal solution is obtained by iteration. A noisy image can be understood as the superposition of original image information and noise. If we assume that the NSCT coefficient of the difference image is y, then y is composed of the NSCT coefficient x of the original image and the NSCT coefficient n of the noise. So the image denoising problem can be transformed into the estimation of the NSCT coefficient x of the original image in the case of the NSCT coefficient y of the noisy image.
Image denoising based on NSCT-HMT model is as follows: 1. By Equation (4), the parameters of the model can be expressed as: where j, k, i respectively refer to the scale, direction, and coefficient. m indicates the hidden state of the NSCT coefficient. 3.
The NSCT-HMT parameter model of the denoised image can be expressed as: The above parameters are used to estimate the NSCT coefficients of the denoised image. When the state S j,k,i is constant, the NSCT coefficients obey the Gauss distribution. It is assumed that the noise is Gauss white noise with a mean value of 0, and it is independent of the NSCT coefficient.
Then the NSCT coefficient of the Bayes denoised image [26] can be represented as: 4. The conditional probability obtained by the EM algorithm [26] can be expressed as:

5.
The NSCT coefficient of the denoised image can be estimated as:

FLICM Clustering
In this paper, we use the FLICMC clustering algorithm proposed by Krindis et al. [31]. This algorithm modifies the objective function of the traditional FCM algorithm, and introduces the fuzzy factor G ki, : where d ij is a spatial euclidean distance between the i pixel and the neighborhood j pixel. µ kj indicates the degree of ownership of the j pixel to class k. x j represents the neighborhood pixels near the center i pixel in the local window. v k is the clustering center of class k. The objective function of the FLICM algorithm is shown below: where x i is the center pixel of local window. v k and µ kj respectively represent the class k clustering centers and fuzzy membership matrix. They can be expressed as:

Implementation Steps
The implementation steps are as follows: Step 1: The original remote sensing images are calibrated by ENVI software; Step 2: The difference image is obtained by using Equation (1) on the calibrated images; Step 3: The NSCT transform is used to transform the difference image, and the transformed coefficients are modeled using a hidden Markov tree. Equation (6) is used to obtain the coefficients of the original image. Finally, Equation (10) is used to estimate the NSCT coefficients of the denoised image.
Step 4: The Inverse NSCT transform is used to obtain the difference image after denoising.
Step 5: After the above treatment, Equation (12) is used to cluster the difference map, and the final change detection results are obtained.
The algorithm flow is shown below as Figure 1.

Implementation Steps
The implementation steps are as follows: Step 1: The original remote sensing images are calibrated by ENVI software; Step 2: The difference image is obtained by using Equation (1) on the calibrated images; Step 3: The NSCT transform is used to transform the difference image, and the transformed coefficients are modeled using a hidden Markov tree. Equation (6) is used to obtain the coefficients of the original image. Finally, Equation (10) is used to estimate the NSCT coefficients of the denoised image.
Step 4: The Inverse NSCT transform is used to obtain the difference image after denoising.
Step 5: After the above treatment, Equation (12) is used to cluster the difference map, and the final change detection results are obtained.
The algorithm flow is shown below as Figure 1.

The Selection of the Window Size for SAR Image
To illustrate the effect of window size on change detection results, we examine the relationship between window size ω × ω and percentage correct classification (PCC).
As can be seen in Figure 2, when the window size is 3 × 3, PCC can get the optimal solution. The reason is that if the selection is too large, the noise removal is relatively clean, but it will lead to the

The Selection of the Window Size for SAR Image
To illustrate the effect of window size on change detection results, we examine the relationship between window size ω × ω and percentage correct classification (PCC).
As can be seen in Figure 2, when the window size is 3 × 3, PCC can get the optimal solution. The reason is that if the selection is too large, the noise removal is relatively clean, but it will lead to the loss of edge and detail information. If selected too small, more details will be retained, but the denoising effect is poor, so for the SAR images in this paper, ω = 3. loss of edge and detail information. If selected too small, more details will be retained, but the denoising effect is poor, so for the SAR images in this paper, ω = 3.

SAR Image Data Description
In order to verify the effectiveness and to illustrate the practicality of the proposed method, we use two sets of real SAR images in [32]. The reason why we choose these two sets of SAR images is that these two groups of images are well known and often used for comparison.

Real Remote Sensing Image of the Bern Data Set
The first sets of real SAR image data of Bern, Switzerland were obtained by ERS-2 in April 1999 and May 1999, and are shown below in Figure 3a,b. The images are 301 × 301 pixels in size and have a gray level of 256. Changes are mainly caused by a flood, and the reference map is shown in Figure 3c.

SAR Image Data Description
In order to verify the effectiveness and to illustrate the practicality of the proposed method, we use two sets of real SAR images in [32]. The reason why we choose these two sets of SAR images is that these two groups of images are well known and often used for comparison.

Real Remote Sensing Image of the Bern Data Set
The first sets of real SAR image data of Bern, Switzerland were obtained by ERS-2 in April 1999 and May 1999, and are shown below in Figure 3a,b. The images are 301 × 301 pixels in size and have a gray level of 256. Changes are mainly caused by a flood, and the reference map is shown in Figure 3c. loss of edge and detail information. If selected too small, more details will be retained, but the denoising effect is poor, so for the SAR images in this paper, ω = 3.

SAR Image Data Description
In order to verify the effectiveness and to illustrate the practicality of the proposed method, we use two sets of real SAR images in [32]. The reason why we choose these two sets of SAR images is that these two groups of images are well known and often used for comparison.

Real Remote Sensing Image of the Bern Data Set
The first sets of real SAR image data of Bern, Switzerland were obtained by ERS-2 in April 1999 and May 1999, and are shown below in Figure 3a,b. The images are 301 × 301 pixels in size and have a gray level of 256. Changes are mainly caused by a flood, and the reference map is shown in Figure 3c.

SAR Image Data Experimental Results and Analysis
To evaluate the effectiveness of the proposed algorithm, we use the above three real remote sensing image data sets for change detection. The algorithm proposed in this paper is compared with the MRF-FCM [11], NSCT-FCM [12], PA-GT [13] and PCA-NLM [15] algorithms. Results of the qualitative analysis are shown in Figures 5 and 6.
As can be seen from Figure 5f, there is no isolated clutter in the resulting graph, which means that the algorithm proposed in this paper can effectively suppress non-changing clutter. It can be seen from Figure 6f that the resulting image effectively suppresses noise effectively and preserves the details of the change. In

SAR Image Data Experimental Results and Analysis
To evaluate the effectiveness of the proposed algorithm, we use the above three real remote sensing image data sets for change detection. The algorithm proposed in this paper is compared with the MRF-FCM [11], NSCT-FCM [12], PA-GT [13] and PCA-NLM [15] algorithms. Results of the qualitative analysis are shown in Figures 5 and 6.
As can be seen from Figure 5f, there is no isolated clutter in the resulting graph, which means that the algorithm proposed in this paper can effectively suppress non-changing clutter. It can be seen from Figure 6f that the resulting image effectively suppresses noise effectively and preserves the details of the change. In

SAR Image Data Experimental Results and Analysis
To evaluate the effectiveness of the proposed algorithm, we use the above three real remote sensing image data sets for change detection. The algorithm proposed in this paper is compared with the MRF-FCM [11], NSCT-FCM [12], PA-GT [13] and PCA-NLM [15] algorithms. Results of the qualitative analysis are shown in Figures 5 and 6.
As can be seen from Figure 5f, there is no isolated clutter in the resulting graph, which means that the algorithm proposed in this paper can effectively suppress non-changing clutter. It can be seen from Figure 6f   The results are shown in Table 1. As can be seen in Table 1, in comparison with MRF-FCM, NSCT-FCM, PA-GT and PCA-NLM, the proposed algorithm has the lowest number of false positives, the lowest total number of errors, and the highest PCC and the Kappa index. As for the running time, although the proposed algorithm takes a bit longer than the NSCT-FCM, PA-GT and PCA-NLM algorithms, the detection accuracy is much higher than theirs.

The General Applicability of Verification Glgorithm for SAR Image
The five algorithms are used to process the experimental data of 30 groups. The average results are shown in Table 2.  The results are shown in Table 1. As can be seen in Table 1, in comparison with MRF-FCM, NSCT-FCM, PA-GT and PCA-NLM, the proposed algorithm has the lowest number of false positives, the lowest total number of errors, and the highest PCC and the Kappa index. As for the running time, although the proposed algorithm takes a bit longer than the NSCT-FCM, PA-GT and PCA-NLM algorithms, the detection accuracy is much higher than theirs.

The General Applicability of Verification Glgorithm for SAR Image
The five algorithms are used to process the experimental data of 30 groups. The average results are shown in Table 2. For the running time, this paper takes the experimental data of 30 groups of 300 × 300 pixels images. The data for each group were run 20 times and then averaged. The results of each parameter are shown in Table 2. The following conclusions can be drawn from Table 2: 1.
Compared with the previous algorithms, the proposed algorithm can suppress clutter noise better while preserving detailed information from the original image.

2.
The proposed algorithm has a faster processing speed than MRF-FCM, but it is a little slower than NSCT-FCM, PA-GT and PCA-NLM. The complexity of the proposed algorithm needs to be improved.

3.
The proposed algorithm is suitable to detect SAR images.

Multi-Spectral Image Change Detection
Due to the fact that SAR image processing is intrinsically complex, with the presence of the speckle noise, it is difficult to analysis SAR images, and the multi-spectral images are easy affected by the presence of cloud cover or different sunlight conditions. The application to SAR and multi-spectral optical data should be different, but there are many algorithms that can be used to detect both the synthetic aperture radar images and the optical images. Just as the algorithms proposed in [13][14][15][16][17][18][19][20]. We are lucky enough to find an algorithm which can be used for both SAR and multi-spectral remote sensing image change detection.

The Selection of the Window Size for Multi-Spectral Remote Sensing Image
When the image spatial resolution is 30 m, the relationship between window size ω × ω and percentage correct classification (PCC) is as follows: as can be seen in Figure 7, when the window size is 1 × 1, PCC can get the optimal solution, so we choose ω = 1 when we are dealing with multi-spectral remote sensing images. For the running time, this paper takes the experimental data of 30 groups of 300 × 300 pixels images. The data for each group were run 20 times and then averaged. The results of each parameter are shown in Table 2. The following conclusions can be drawn from Table 2: 1. Compared with the previous algorithms, the proposed algorithm can suppress clutter noise better while preserving detailed information from the original image. 2. The proposed algorithm has a faster processing speed than MRF-FCM, but it is a little slower than NSCT-FCM, PA-GT and PCA-NLM. The complexity of the proposed algorithm needs to be improved. 3. The proposed algorithm is suitable to detect SAR images.

Multi-Spectral Image Change Detection
Due to the fact that SAR image processing is intrinsically complex, with the presence of the speckle noise, it is difficult to analysis SAR images, and the multi-spectral images are easy affected by the presence of cloud cover or different sunlight conditions. The application to SAR and multispectral optical data should be different, but there are many algorithms that can be used to detect both the synthetic aperture radar images and the optical images. Just as the algorithms proposed in [13][14][15][16][17][18][19][20]. We are lucky enough to find an algorithm which can be used for both SAR and multi-spectral remote sensing image change detection.

The selection of the Window Size for Multi-spectral Remote Sensing Image
When the image spatial resolution is 30 m, the relationship between window size ω × ω and percentage correct classification (PCC) is as follows: as can be seen in Figure 7, when the window size is 1 × 1, PCC can get the optimal solution, so we choose ω = 1 when we are dealing with multi-spectral remote sensing images. When the window size is 1 × 1, Equation (1) will be changed into: where X1 (i, j) and X2 (i, j) respectively represent the value of the pixel at the point ( i , j ).

Image Data Description
The set of real multi-spectral remote sensing image data are acquired by the Landsat Enhances Thematic Mapper Plus (ETM+) sensor of the Landsat-7 satellite, in an area of Mexico in April 2000 When the window size is 1 × 1, Equation (1) will be changed into: where X 1 (i, j) and X 2 (i, j) respectively represent the value of the pixel at the point (i,j).

Image Data Description
The

Multi-Spectral Image Data Experimental Results and Analysis
The evaluation indexes are false negatives (FN), false positives (FP), overall error (OE), percentage correct classification (PCC), the Kappa index and the running time (T). The following conclusions can be drawn from Figure 9 and in Table 3: 1. The PA-GT and PCA-NLM algorithms can achieve good results when dealing with SAR images. But, in dealing with real multi-spectral remote sensing images, the results are poor. 2. The proposed algorithm can be used for multi-spectral remote sensing image change detection.

Multi-Spectral Image Data Experimental Results and Analysis
The evaluation indexes are false negatives (FN), false positives (FP), overall error (OE), percentage correct classification (PCC), the Kappa index and the running time (T). The following conclusions can be drawn from Figure 9 and in Table 3: 1.
The PA-GT and PCA-NLM algorithms can achieve good results when dealing with SAR images. But, in dealing with real multi-spectral remote sensing images, the results are poor.

2.
The proposed algorithm can be used for multi-spectral remote sensing image change detection.

Multi-Spectral Image Data Experimental Results and Analysis
The evaluation indexes are false negatives (FN), false positives (FP), overall error (OE), percentage correct classification (PCC), the Kappa index and the running time (T). The following conclusions can be drawn from Figure 9 and in Table 3: 1. The PA-GT and PCA-NLM algorithms can achieve good results when dealing with SAR images. But, in dealing with real multi-spectral remote sensing images, the results are poor. 2. The proposed algorithm can be used for multi-spectral remote sensing image change detection.   The five algorithms are used to process the experimental data of 30 groups. The images are 500 × 500 pixels in size. The average results are shown in Table 4. The following conclusions can be drawn from Table 4: 1.
Compared with the previous four algorithms, the proposed algorithm has a good balance in the number of false negatives and false positives.

2.
The proposed algorithm is more suitable to detect multi-spectral images.

Acquisition of Real Remote Sensing Image
The Gurbantunggut Desert is located in the hinterland of the Junggar basin in northern Xinjiang, China. It is China's largest fixed and semi-fixed desert, and is located at 44 • 11 ~46 • 20 N and 84 • 31 ~90 • 00 E. Its acreage is about 4.88 × 10 4 km 2 . Precipitation occurs mainly in the spring and annually does not reach 150 mm. Its average annual evaporation potential is higher at more than 2000 mm, and the annual average temperature is 6~10 • C.
We choose a set of real multi-spectral remote sensing image data which are acquired by the Landsat Enhances Thematic Mapper Plus (ETM+) sensor of the Landsat-7 satellite, in an area of Gurbantunggut Desert in July 1999 and July 2007. The instrument's pixel resolution is 30 m. Figure 10a Landsat Enhances Thematic Mapper Plus (ETM+) sensor of the Landsat-7 satellite, in an area of Gurbantunggut Desert in July 1999 and July 2007. The instrument's pixel resolution is 30 m. Figure 10a,b shows channel 3 of the 1999 and 2007 images, respectively. They are registered using ENVI. The images are 570 × 570 pixels, and have a gray scale of 256. The image data is mainly composed of desert, vegetation, and street construction. Changes mainly result from the government and citizens combating desertification and planting vegetation and crops.

Analysis the Change of Oasis Cover
The above five algorithms are used to detect the change of Gurbantunggut data set. The results of the change of oasis cover in the study area are shown in Figure 11 and Table 5. The white pixels in Figure 11 represent the oasis cover change.

Analysis the Change of Oasis Cover
The above five algorithms are used to detect the change of Gurbantunggut data set. The results of the change of oasis cover in the study area are shown in Figure 11 and Table 5. The white pixels in Figure 11 represent the oasis cover change.  In the blue and purple areas, parts of the places have indeed changed, but the MRF-FCM and PCA-NLM algorithms did not detect the changed parts. In the red area, there is no change in this area, but NSCT-FCM and PA-GT algorithms have detected changes.

Conclusions
A remote sensing image change detection algorithm based on the NSCT-HMT denoising model is proposed in this paper. The proposed change detection method is tested for both optical and synthetic aperture radar satellite images. The mean-ratio operation is adopted to obtain the difference image, which will be denosing by the NSCT-HMT model. Then, using the FLICM algorithm, the difference image is divided into the change area and unchanged area. The experimental results show that compared with MRF-FCM, NSCT-FCM, PA-GT and PCA-NLM algorithms, the algorithm  In the blue and purple areas, parts of the places have indeed changed, but the MRF-FCM and PCA-NLM algorithms did not detect the changed parts. In the red area, there is no change in this area, but NSCT-FCM and PA-GT algorithms have detected changes.

Conclusions
A remote sensing image change detection algorithm based on the NSCT-HMT denoising model is proposed in this paper. The proposed change detection method is tested for both optical and synthetic aperture radar satellite images. The mean-ratio operation is adopted to obtain the difference image, which will be denosing by the NSCT-HMT model. Then, using the FLICM algorithm, the difference image is divided into the change area and unchanged area. The experimental results show that compared with MRF-FCM, NSCT-FCM, PA-GT and PCA-NLM algorithms, the algorithm proposed in this paper takes full account of the characteristics of the spatial neighborhood. It not only enhances the changes of the regional profile and the small area, but also can effectively suppress the clutter noise, and the final change detection accuracy is improved. Finally, we apply the algorithms to detect the oasis coverage of the southern Gurbantunggut Desert, and the results show that the proposed algorithm has a good detection ability.
In this paper, the multi-spectral optical data sets are acquired by the Landsat-7, so we did not list the relationship between the value of window size and the image spatial resolution in the paper, but we did some experiments, and we found that when the image spatial resolution is lower than 10 m, the value of window size must choose 1 × 1, or using other method such as the different method or log-ratio method, etc., to get the multi-spectral optical different image. There are many places in this paper that need to be improved. For example, how to reduce the complexity of the algorithm and how to further improve the accuracy of the algorithm. In our future investigations, additional work will be conducted on them.