Remote sensing image change detection method

Abstract

The present invention provides a remote sensing image change detection method, characterized by comprising: acquiring old time-phase data and new time-phase data, wherein the new time-phase data is a new time-phase image, and the old time-phase data comprises vector data and an old time-phase image; for these different types of old time-phase data, using corresponding change detection strategies to perform change detection, so as to obtain two change detection results; establishing an error matrix to score a plurality of accuracy indicators; on the basis of score values and weights corresponding to the respective accuracy evaluation indicators, using an accuracy coupling evaluation model to compute coupling evaluation values of the change detection results; and using the change detection result having the higher coupling evaluation value as a final change detection result. The present invention combines image-to-image and vector-to-image change detection methods, and fully utilizes historical result data. The coupling evaluation model is used to comprehensively score the two change detection results, and to evaluate the change detection results from the aspects of numerical errors and error spatial distribution, which has the advantages of high detection accuracy and high efficiency and improves the intelligence level of change detection.

Description

A remote sensing image change detection method Technical Field

[0001] The present invention relates to the field of image processing technology, and in particular to a remote sensing image change detection method. Background Art

[0002] Change detection is the process of identifying changes in the state of a ground object or phenomenon through multiple observations. It has been widely used in fields such as land surveys, disaster assessment, urban research, and ecosystem monitoring. Depending on the data source used, change detection methods can be divided into two categories: image-to-image and vector-to-image.

[0003] Image-to-image methods use time-series remote sensing imagery as a data source, employ a specific similarity measure to measure differences between images, and perform change detection based on a certain threshold. Classic image-to-image methods include image difference, image ratio, change vector analysis, and post-classification comparison. Currently, a large number of classic methods have emerged, including image difference, change vector analysis, statistical methods, and post-classification comparison. Among them, the image difference method is based on image difference maps. While simple, its disadvantage is that it is susceptible to noise and other factors, resulting in low detection accuracy. The K-means algorithm is a common clustering method used for image change detection because it is simple and has good convergence, but it is sensitive to noise and prone to falling into local optimality. Some people have proposed the PCA-K-means method based on PCA for change detection, but this method will lose a lot of image detail information. The Fuzzy C-means (FCM) algorithm is a clustering algorithm based on fuzzy theory and is also widely used in image change detection. However, it has the problem that the initial value of the cluster center has a great influence on the clustering result, making it difficult to obtain the global optimal solution, and the detection accuracy still needs to be improved.

[0004] Vector-imagery methods utilize vector maps from the past and remote sensing images from the present for change detection. Existing vector-imagery methods can be categorized into two types: post-classification comparison and category-feature methods. The post-classification comparison method overlays the vector map onto the remote sensing image, extracting patches of unchanged category as training samples. This method then performs supervised classification on the remote sensing image and compares the vector map categories with the classification results to detect changes. Category-feature methods overlay the vector map onto the remote sensing image, extracting features of patches of each category within the image, analyzing the patterns of these patches by category, identifying patches that have changed, and integrating the change results across categories to detect changes. Existing research on vector-imagery methods is relatively limited, and a mature methodology has yet to be established. Furthermore, existing category-feature methods are poorly effective at detecting locally changed patches, failing to fully leverage existing vector data and meet the needs of intelligent change detection.

[0005] In addition, existing technologies often use the Kappa coefficient as an indicator of accuracy. A larger value indicates higher accuracy, effectively evaluating detection results to a certain extent and being widely used in accuracy assessment research. However, the Kappa coefficient is a coefficient based on probabilistic events. While it can provide a numerical estimate of accuracy, it lacks information about the spatial representation of error. This quantitative analysis, performed solely through the Kappa coefficient, lacks a coupled indicator evaluation system that comprehensively considers both numerical and spatial error in the accuracy of change detection.

[0006] Summary of the Invention

[0007] In response to the technical problems existing in the prior art, the present invention provides a remote sensing image change detection method. According to the type of old phase data in the image library, a corresponding change detection strategy is selected. In the target area, change detection is performed on the new and old phase data. The method has both image-image and vector-image change detection methods, making full use of historical results data. A coupling evaluation model is introduced to comprehensively evaluate the two change detection results, and the one with a higher coupling evaluation value is selected as the final detection result. This method has the advantages of high detection accuracy and high efficiency, and improves the automation and intelligence of change detection.

[0008] The present invention provides a remote sensing image change detection method, comprising the following steps:

[0009] Old phase data and new phase data are acquired, wherein the old phase data includes first type data and second type data, and the new phase data is a new phase image.

[0010] For the two types of old phase data, corresponding change detection strategies are respectively adopted to perform change detection on the old phase data and the later phase data to obtain a first change detection result and a second change detection result.

[0011] Multiple precision evaluation indicators are introduced, and an error matrix is ​​established to score multiple precision indicators. According to the score value and weight corresponding to each precision evaluation indicator, the precision coupling evaluation model is used to calculate the coupling evaluation value of the first change detection result and the second change detection result, and the change detection result with the highest coupling evaluation value is taken as the final change detection result.

[0012] Specifically, the first type of data is vector results, and the second type of data is old phase images.

[0013] Specifically, when the old phase data is a vector product, a first change detection strategy is adopted to obtain a first change detection result; when the old phase data is an old phase image, a second change detection strategy is adopted to obtain a second change detection result.

[0014] Specifically, the first change detection strategy includes the following steps:

[0015] The vector results are aligned with the new phase image to obtain the aligned new phase image, wherein the vector results include feature category information and old phase feature patches; under the constraints of the vector results, the aligned new phase image is segmented to obtain new phase feature patches; the feature information of the old phase feature patches and the new phase feature patches are extracted respectively to obtain a first feature map and a second feature map; the feature distances of corresponding patches in the first feature map and the second feature map are calculated by two change detection sub-algorithms to obtain a first initial change index set and a second initial change index set; the first initial change index set and the second initial change index set are weightedly combined to obtain a final change index set of the patches.

[0016] The Bayesian algorithm is used to obtain the change threshold of each object category. The patches are classified into two categories according to the final change index of the patches using the change threshold, and the changed patches are extracted to obtain the first change detection result.

[0017] Specifically, the change detection sub-algorithm includes a cosine distance algorithm and a Euclidean distance algorithm.

[0018] Specifically, the final change index calculation method is:

[0019] Among them, D(O1,O2) is the final change index, d r (O1, O2) is the rth initial change index calculated by the rth change detection sub-algorithm; w ir is the weight of the rth change detection sub-algorithm, ∑w r =1.

[0020] Specifically, the second change detection strategy includes the following steps:

[0021] The old phase image and the new phase image are preprocessed, and the grayscale value difference between the two is calculated to form a grayscale difference map; multiple texture features of the old phase image and the new phase image are extracted respectively, and the difference matrix of multiple texture feature statistics is calculated to form a texture difference map; the grayscale difference map and the texture difference map are fused to obtain the final difference map; the improved fuzzy C-means method is used to perform binary classification on the final difference map to obtain the second change detection result.

[0022] Specifically, the accuracy evaluation indicators include Kappa coefficient and accuracy rate.

[0023] Specifically, the method for calculating the coupling evaluation value includes the following steps:

[0024] The Kappa coefficient and the score corresponding to the accuracy are Z-score standardized to obtain the standardized score value Z k and Zc .

[0025] The entropy weight method is used to calculate the standardized score value Z s and Z h The weight N s and N h .

[0026] According to Z s and Z h and its weight N s and N h , build a coupling evaluation model for accuracy.

[0027] The coupling evaluation values ​​of the first change detection result and the second change detection result are respectively calculated according to the coupling evaluation model, and the change detection result with the higher coupling evaluation value is used as the final change detection result.

[0028] Specifically, the coupling evaluation model is: Q=Z s N s +N h (AZ h ), where A is the standardized index value Z h The maximum distance of the interval is rounded to an integer value.

[0029] The beneficial effects of the present invention are:

[0030] (1) The present invention selects the corresponding change detection strategy according to the type of old phase data in the image library, and has both vector-image and image-image change detection methods, making full use of historical results data. It also introduces a coupled evaluation model to comprehensively evaluate the accuracy of change detection results from two aspects: numerical error and error spatial distribution. The highest accuracy is selected as the final change detection result, thereby improving the accuracy of change detection.

[0031] (2) In the vector-image change detection strategy of the present invention, the initial change index of the pattern obtained by calculating each change detection sub-algorithm is combined by weighting to obtain the final change index of the pattern, which can avoid the error caused by the interference of the changing pattern by a single algorithm, thereby improving the accuracy of change detection.

[0032] (3) In the image-to-image process, the difference map is constructed by fusing the grayscale features and texture features of the two-phase images, which solves the problem of insufficient information when constructing the difference map based on a single type of feature. In the iterative process of the fuzzy C-means algorithm, the algorithm convergence speed is accelerated by modifying the membership of the pixel point closest to the class center, thereby improving the convergence speed of the fuzzy C-means algorithm and thus improving the accuracy and efficiency of change detection. BRIEF DESCRIPTION OF THE DRAWINGS

[0033] In order to more clearly illustrate the embodiments of the present invention or the technical solutions in the prior art, the following briefly introduces the drawings required for use in the embodiments. Obviously, the drawings described below are only some embodiments of the present invention. For ordinary technicians in this field, other drawings can be obtained based on these drawings without paying any creative work.

[0034] FIG1 is a schematic diagram of a flow chart of a remote sensing image change detection method according to an embodiment of the present invention;

[0035] FIG2 is a technical flow chart of a first change detection strategy according to an embodiment of the present invention;

[0036] FIG3 is a technical flow chart of a second change detection strategy according to an embodiment of the present invention;

[0037] FIG4 is a technical flow chart of binary classification of the final difference map G in the second change detection strategy according to an embodiment of the present invention. DETAILED DESCRIPTION

[0038] The following will clearly and completely describe the technical solutions in the embodiments of the present invention in conjunction with the accompanying drawings. Obviously, the described embodiments are only part of the embodiments of the present invention, not all of the embodiments. It should be noted that as long as no conflict arises, the various embodiments and features of each embodiment of the present invention can be combined with each other, and the resulting technical solutions are all within the scope of protection of the present invention.

[0039] Referring to FIG1 , the present invention provides a remote sensing image change detection method, comprising the following steps:

[0040] Old phase data and new phase data are acquired, wherein the old phase data includes first type data and second type data, and the new phase data is a new phase image.

[0041] In an embodiment of the present invention, the first type of data is vector data, and the second type of data is old phase images.

[0042] For the two types of old phase data, corresponding change detection strategies are respectively adopted to perform change detection on the old phase data and the later phase data to obtain a first change detection result and a second change detection result.

[0043] In an embodiment of the present invention, when the old phase data is a vector product, a first change detection strategy is adopted to obtain a first change detection result; when the old phase data is an old phase image, a second change detection strategy is adopted to obtain a second change detection result.

[0044] In an embodiment of the present invention, referring to FIG2 , the first change detection strategy includes the following steps:

[0045] Step 1: Register the vector result with the new phase image to obtain the registered new phase image. The vector result includes the feature category information and the old phase feature patches.

[0046] In an embodiment of the present invention, manual registration correction of the old phase vector result includes moving the overall position of the vector contour and adding, deleting and moving vector nodes, so that the registered vector result accurately matches the contour in the new phase image.

[0047] Step 2: Under the constraints of the vector results, the registered new phase image is segmented to obtain the new phase ground feature patches.

[0048] In the embodiment of the present invention, the vector result of the old phase is represented as V, and the ground feature category set in the vector result is Ω={Ω1,Ω2,…,Ω p}, where Ωp (1≤p≤q) is the pth land feature category; q is the total number of land feature categories.

[0049] The labeled watershed algorithm is a fast image segmentation method that can obtain continuous single-pixel edges. This method first generates a gradient image, then uses a gradient threshold to obtain a labeled region representing the interior of the feature. Flooding is then performed outward from the labeled region. When floodwaters from different labeled regions merge, a watershed is set to separate them. The watershed corresponds to an edge in the image. In an embodiment of the present invention, the edges of the vector results are first mapped onto the new phase image to form a priori watershed. The new phase image is then segmented using the traditional labeled watershed algorithm, subject to the constraints of the prior watershed.

[0050] The new phase image is represented as R, which is divided into n patches, that is, R = {R1, R2, ..., R n}, where R m (1≤m≤n) is the mth patch; the feature categories of all patches in the old period are set as T={T1,T2,…,T n}, where T m (1≤m≤n) is the patch R m The feature category in the old phase vector results, and T m ∈Ω.

[0051] In an embodiment of the present invention, vector-result-constrained image segmentation of the new phase preserves the boundaries of the original vector result while also segmenting newly created edges in the new phase image. For example, a vector result contains two image patches, R1 and R2. A local change occurs within patch R1. After image segmentation under vector-result constraints, patch R3 is added within R1. In the vector result, patches R1 and R3 have the same feature category.

[0052] In this embodiment of the present invention, under the constraints of the vector output V, the labeled watershed algorithm is used to segment the new phase image. The labeled watershed algorithm is a fast image segmentation method that can obtain continuous single-pixel edges. This method first generates a gradient image, then uses a certain gradient threshold to obtain a labeled region representing the interior of the feature. Flooding is then performed outward based on the labeled region. When floodwaters from different labeled regions merge, a watershed is set to separate them. The watershed corresponds to the edge in the image.

[0053] In an embodiment of the present invention, the edges of the vector results are first mapped onto the new phase image to form a priori watershed; then, under the constraints of the prior watershed, the new phase image is segmented using a traditional labeled watershed algorithm.

[0054] The new phase image is represented as R, which is divided into n patches, that is, R = {R1, R2, ..., R n}, where R m (1≤m≤n) is the mth patch, and the set of land feature categories of all patches in the old period is T={T1,T2,…,T n}, where T m (1≤m≤n) is the patch R m The feature category in the old phase vector results, and T m ∈Ω.

[0055] Step 3: Extract the feature information of the old-phase ground feature patches and the new-phase ground feature patches respectively to obtain the first feature map and the second feature map.

[0056] In an embodiment of the present invention, the pre-trained ResNet50 is used for feature extraction. ResNet50 includes a convolutional layer conv1 and four convolutional blocks conv2_x-conv5_x. The four convolutional blocks conv2_x-conv5_x include 3, 4, 6, and 3 residual blocks respectively. The structure of each residual block is shown in Figure 3, and the first feature map and the second feature map are obtained.

[0057] Step 4: Calculate the feature distances of corresponding spots in the first feature map and the second feature map using two change detection sub-algorithms to obtain a first set of initial change indices and a second set of initial change indices.

[0058] In the embodiment of the present invention, two change detection sub-algorithms, namely the cosine distance and the Euclidean distance algorithm, are used to respectively calculate the characteristic distances of the corresponding patches in the first feature map O1 and the second feature map O2 (i.e., the old phase feature patches and the new phase feature patches in the vector results at the same geographical location), and obtain a plurality of initial change indices. r}, O2={y1,y2,..,yr}, the specific steps are:

[0059] The cosine distance is used as the first change detection sub-algorithm. The feature distance of the corresponding patches between two feature maps is calculated according to the following formula to obtain the first initial change index set d1(O1, O2):

[0060] Using Euclidean distance as the second change detection sub-algorithm, the feature distance of corresponding patches between two feature maps is calculated according to the following formula to obtain the second initial change index set d2(O1, O2):

[0061] Among them, r=1,2,...,n, n is the number of spots in the feature map, x r is the rth spot in the first feature map O1, x r ∈O1;y r is the rth spot in the second feature map O2, y r ∈O2; d1(O1,O2) is the first feature in the corresponding patch, feature x r To the second feature feature map feature y r The cosine distance set, that is, the first initial change index set; d2(O1, O2) is the feature x in the first feature feature map in the corresponding pattern r To the second feature feature map feature y r The Euclidean distance set of , that is, the second initial change index set.

[0062] Step 5: Perform a weighted combination of the first initial change index set and the second initial change index set to obtain a final change index set of the patch.

[0063] The two change detection sub-algorithms are used to calculate the first initial change index set d1(O1, O2) and the second initial change index set d2(O1, O2) of the patches, respectively. By weighted combination, the final change index of the corresponding patches between the first feature map and the second feature map is obtained. The calculation formula is:

[0064] D(O1,O2) is the final change index, d i (O1, O2) is the rth initial change index set calculated by the rth change detection sub-algorithm; w r is the weight of the rth change detection sub-algorithm, ∑w r =1.

[0065] In the embodiment of the present invention, the weights of the two change detection sub-algorithms are both 0.5, and the final change index set of the image spots is obtained. The specific calculation formula is:

[0066] Step 6: Use the Bayesian algorithm to obtain the change threshold of each object category, use the change threshold to classify the patches into two categories according to the final change index corresponding to the patches, extract the changed patches, and obtain the first change detection result.

[0067] In the embodiment of the present invention, D(O1, O2)={z1, z2, ..., z r} is the set of the final change index of the spots, and r is the number of spots in the feature map. Compare the size of each final change index in D(O1,O2) with the obtained change threshold, and then divide the corresponding spots into change categories ω c With unchanged class ω u .

[0068] Assume that the conditional probability density function P(x|ω) of the changed class and the unchanged class c ), P(x|ω u ) obeys Gaussian distribution:

[0069] in For the change class ω c and unchanged class ω u The mean and standard deviation of The total density function P(x) of the spot change index can be expressed as: P(x)=P(ω c )×P(x|ω c )+P(ω u )×P(x|ω u ) (2)

[0070] Under the assumptions, the process of solving the change threshold can be transformed into estimating the change category ω c With unchanged class ω u The process of estimating the mean and standard deviation can be achieved by applying the expectation maximization (EM) algorithm. Each iteration consists of two steps: the expected value step and the expectation maximization step. The expected value step directly estimates the expected value of the probability density from the observed data based on the current value of the parameter to be estimated. The expectation maximization step updates the parameter estimate by maximizing the expected value. These two steps are performed alternately throughout the iterative process until convergence.

[0071] Apply EM algorithm to the two sub-Gaussian models for unchanged category ω u Parameter estimation, the calculation formulas for the prior probability, mean and standard deviation are:

[0072] Where t and t+1 are the current estimate and the estimate used in the next iteration respectively, and r is the number of patches. c The estimation method of the mean and standard deviation is the same as ω u .

[0073] When the difference between the prior probability, mean and standard deviation calculated between two adjacent iterations is less than the given threshold (ε=1×10 -8 ), the iteration terminates.

[0074] Obtain the mean of the unchanged class and the changed class according to the EM algorithm Standard deviation Then, the change threshold T is calculated according to the Bayesian minimum error probability theory: P(T|ω u )P(ω u )=P(T|ω c )×P(ω c ) (4)

[0075] By sorting out formula (4), we can get the solution formula of the change threshold T as follows:

[0076] Where μ u 、μ c , δ u , δ c 、P(ω u ) and P(ω c ) have been calculated using formula (3). The change threshold can be calculated using formula (5). Using the change threshold, the final change index is binary segmented to generate the results of change and unchanged semantics, and the patches of the change class are extracted as the first detection result.

[0077] In an embodiment of the present invention, referring to FIG3 , the second change detection strategy includes the following steps:

[0078] Step 1: Preprocess the old phase image and the new phase image, and calculate the grayscale value difference between the two to form a grayscale difference map.

[0079] In the embodiment of the present invention, pre-processing operations such as radiation correction and geometric correction are performed on the old phase image I1 and the new phase image I2 (both with a size of M×N) to be detected.

[0080] Calculate the grayscale difference image E0 of the preprocessed images I1 and I2: E0(i,j)=|I1(i,j)-I2(i,j)|

[0081] Where, 1≤i≤M, 1≤j≤N, I1(i,j) and I2(i,j) are the grayscale values ​​of the corresponding pixels in the two temporal images I1 and I2.

[0082] Step 2: Extract the four texture features of the old phase image and the new phase image, namely energy, contrast, correlation and entropy; and calculate the difference matrix F of the four texture feature statistics of images I1 and I2 to form a texture difference map F'.

[0083] Energy: reflects the uniformity and compactness of the texture in the image. The greater the difference in the values ​​in the gray-level co-occurrence matrix, the greater the energy value.

[0084] Where R(i,j) is the normalized gray-level co-occurrence matrix, and (i,j) represents the matrix pixel coordinates. Energy reflects the uniformity and compactness of texture in an image. The greater the difference in the gray-level co-occurrence matrix values, the greater the energy value.

[0085] Contrast: Describes the brightness contrast of pixel values ​​in an image, which is visually manifested as the clarity and vividness of the image color.

[0086] Correlation: describes the similarity between rows and columns in the matrix, reflecting the consistency of image texture.

[0087] Among them, μ x , μ y , Represents m x =∑ j R i,j and m y =∑ i R i,j The mean and standard deviation of .

[0088] Entropy: describes the average amount of information in an image and reflects the complexity of the texture in the image.

[0089] Taking the old phase image I1 as an example, a window of size m×n is slid across the image, and the four texture feature statistics within each window are calculated as the texture feature value of the central pixel in the window. Each movement distance is 1 pixel, and 0° scanning, 45° scanning, 90° scanning, and 135° scanning are performed in sequence. When the window is moved across the entire image, the texture data in the four scanning directions are obtained and superimposed, and the energy matrix F of image I1 is finally obtained. 11 , contrast matrix F 12 , the correlation matrix F 13 and entropy matrix F 14 , the size is M×N. Similarly, the energy matrix F of image I2 is obtained 21 , contrast matrix F 22 , the correlation matrix F 23 and entropy matrix F24 .

[0090] Taking the energy matrix as an example, calculate the energy difference matrix F1 of images I1 and I2: F1(i,j)=|F 11 (i,j)-F 21 (i,j)|

[0091] Similarly, we obtain the contrast difference matrix F2, the correlation difference matrix F3, and the entropy difference matrix F4. By superimposing the four characteristic statistic difference matrices, we obtain the texture feature difference matrix F. In order to make the values ​​of the grayscale difference map E0 and the texture difference map F' have the same value range, it is necessary to normalize the texture difference map D to obtain the normalized texture difference map F': The calculation method of the texture difference map F' is:

[0092] Among them, F max is the maximum value in the difference matrix F.

[0093] Step 3: Fuse the grayscale difference map E0 and the texture difference map F' to obtain the final difference map G, which is calculated as follows:

[0094] The final difference map G contains the grayscale value information and texture feature information of the two temporal remote sensing images I1 and I2;

[0095] Step 4: Use the improved fuzzy C-means method to perform binary classification on the final difference map G to obtain the second change detection result.

[0096] In an embodiment of the present invention, referring to FIG. 4 , step 4 specifically includes the following steps:

[0097] (1) Initialize the parameters of fuzzy C-means; fuzzy weighting index m = 2, maximum number of iterations T = 30, minimum error for stopping iteration ε = 0.0001, initial number of iterations t = 1, total number of clusters c = 2;

[0098] (2) The final difference map G is used as the input image of the improved fuzzy C-means method. The total number of pixels is n. c1 and c2 represent the changed class and unchanged class in the final difference map respectively. v1 and v2 represent the cluster centers of the two classes respectively. The cluster centers are initialized: v1 0 =X max , v2 0 =X min , where X max and X minare the maximum and minimum grayscale values ​​of the final difference map respectively; select C1 points with the smallest distance ratio, the distances of these pixels to the center of the changed class are significantly smaller than the distances to the center of the unchanged class, and are used to represent the number of points belonging to the changed class; select C2 points with the largest distance ratio, the distances of these pixels to the center of the changed class are significantly larger than the distances to the center of the unchanged class, and are used to represent the number of points belonging to the unchanged class.

[0099] (3) Calculate the membership matrix u of all pixels 1k (i=1,...,c):

[0100] Among them, p k represents the kth pixel, v i is the cluster center of the i-th class, v j is the cluster center of the jth class, the total number of clusters c = 2, and n is the total number of image pixels.

[0101] (4) Calculate the cluster centers V1 of all pixels in the final difference map to the two categories c1 and c2 respectively t-1 and V2 t-1 The Euclidean distances M1 and M2 and the distance ratio K of all pixels to the centers of the two categories are calculated as follows:

[0102] (5) Sort all the values ​​in the distance ratio K and select C1 pixels with the smallest distance as the change class, and modify the membership u of the change class pixels 1k =1,u 2k = 0; select the C2 pixel points with the largest distance ratio as the unchanged class, and modify the membership u of the unchanged class points 1k =0,u 2k =1.

[0103] (6) Calculate the change class center v1 t and unchanged class center v2 t , the calculation formula is:

[0104] And update the grayscale values ​​of the C1 and C2 pixels selected in step (5) as their respective cluster center values.

[0105] (7) If t>T (T represents the preset threshold of the number of iterations, and t represents the number of iterations that have been performed), or (where the minimum error for stopping iteration is ε=0.0001) and Then stop the iteration and get the final membership matrix u ik (i=1,2); otherwise t=t+1, increase the values ​​of C1 and C2, and jump to step (3).

[0106] (8) According to the final membership matrix u ik (i=1,2), judge the category attributes of all pixels: when u 1k >u 2k (k=1,2,…,n), then pixel k belongs to the changed class, otherwise pixel k belongs to the unchanged class;

[0107] (9) Obtain a second change detection result.

[0108] In an embodiment of the present invention, multiple precision evaluation indicators are introduced, and an error matrix is ​​established to score the multiple precision indicators. According to the score values ​​and weights corresponding to each precision evaluation indicator, the precision coupling evaluation model is used to calculate the coupling evaluation value of the first change detection result and the second change detection result, and the change detection result with the highest coupling evaluation value is used as the final change detection result.

[0109] In an embodiment of the present invention, an error matrix is ​​established to calculate the Kappa coefficient and the OA accuracy rate to obtain corresponding scoring values.

[0110] The calculation formula for the accuracy OA score is:

[0111] Among them, U accurate Indicates the number of pixels of the correct change class, which is the sum of the diagonal of the confusion matrix, U total is the total number of pixels.

[0112] The formula for calculating the score value of the Kappa coefficient is:

[0113] Among them, U represents the number of pixels in the change class; x tt Represents the value on the diagonal of the t-th category in the confusion matrix; s represents the total number of categories; x t∑ represents the sum corresponding to the column where the tth category is located; x ∑t Represents the sum corresponding to the row of category t.

[0114] Then, the Kappa coefficient and the corresponding score values ​​of the accuracy rate were imported into SPSS software and Z-score standardization was performed to obtain the standardized score value Z s and Z h and its weight N s and N h ; Then Z s and Z h Input it into the data_entropy.xls file, call the entropy weight method command in MATLAB software, and calculate the index value Z s and Z hThe weights are 0.6370 and 0.3630, and they are substituted into the simulation accuracy coupling evaluation model Q=Z s N s +N h (AZ h ), where A is the normalized index value Z h The maximum distance of the interval is rounded to an integer value. That is, we get: Q=0.6370Z s +0.3630(3-Z h )

[0115] The coupling evaluation values ​​Q of the first change detection result and the second change detection result are respectively calculated according to the coupling evaluation model, and the change detection result with the higher coupling evaluation value is taken as the final change detection result.

[0116] It should be noted that the terms used herein are only for describing specific embodiments and are not intended to limit the exemplary embodiments according to the present invention. As used herein, unless the context clearly indicates otherwise, the singular form is intended to include the plural form. In addition, it should be understood that when the terms "comprise" and / or "include" are used in this specification, they indicate the presence of features, steps, operations, devices, components and / or combinations thereof.

[0117] The above embodiments are only used to illustrate the present invention, and are not intended to limit the present invention. Ordinary technicians in the relevant technical field may make various changes and modifications without departing from the spirit and scope of the present invention. Therefore, all equivalent technical solutions also fall within the scope of the present invention. The scope of patent protection of the present invention should be defined by the claims.

Claims

1. A remote sensing image change detection method, characterized in that: The following steps are involved: Acquire old phase data and new phase data, wherein the old phase data includes first type data and second type data, and the new phase data is a new phase image; For the two types of old phase data, corresponding change detection strategies are respectively adopted to perform change detection on the old phase data and the later phase data to obtain a first change detection result and a second change detection result; Multiple precision evaluation indicators are introduced, and an error matrix is ​​established to score multiple precision indicators. According to the score values ​​and weights corresponding to each precision evaluation indicator, the precision coupling evaluation model is used to calculate the coupling evaluation value of the first change detection result and the second change detection result, and the change detection result with a high coupling evaluation value is taken as the final change detection result.

2. The method according to claim 1, characterized in that The first type of data is a vector product, and the second type of data is an old phase image.

3. The method according to claim 2, characterized in that When the old phase data is a vector product, a first change detection strategy is adopted to obtain a first change detection result; when the old phase data is an old phase image, a second change detection strategy is adopted to obtain a second change detection result.

4. The method according to claim 3, characterized in that The first change detection strategy comprises the following steps: Registering the vector result with the new phase image to obtain the registered new phase image, wherein the vector result includes the ground object category information and the ground object spots of the old phase; Under the constraints of vector results, the registered new phase image is segmented to obtain the new phase ground feature patches; Extracting feature information of old-phase ground feature patches and new-phase ground feature patches respectively to obtain a first feature map and a second feature map; Calculate the feature distances of corresponding spots in the first feature map and the second feature map by using two change detection sub-algorithms to obtain a first initial change index set and a second initial change index set; The first initial change index set and the second initial change index set are weightedly combined to obtain a final change index set of the image patch. The Bayesian algorithm is used to obtain the change threshold of each object category. The spots are classified into two categories according to the change threshold and the final change index of the spots is used to extract the change spots to obtain the first change detection result.

5. The method according to claim 4, characterized in that The change detection sub-algorithms include a cosine distance algorithm and a Euclidean distance algorithm.

6. The method according to claim 5, characterized in that The calculation method of the final change index is: Among them, D(O1,O2) is the final change index, d r (O1, O2) is the rth initial change index calculated by the rth change detection sub-algorithm; r is the weight of the rth change detection sub-algorithm, ∑w r =1.

7. The method according to claim 3, characterized in that The second change detection strategy comprises the following steps: Preprocess the old phase image and the new phase image, and calculate the gray value difference between the two to form a gray value difference map; Extracting multiple texture features of old phase images and new phase images respectively, calculating the difference matrix of multiple texture feature statistics, and forming a texture difference map; Fuse the grayscale difference map and the texture difference map to obtain the final difference map; The improved fuzzy C-means method is used to perform binary classification on the final difference map to obtain the second change detection result.

8. The method according to claim 1, characterized in that The accuracy evaluation indicators include Kappa coefficient and accuracy rate.

9. The method according to claim 8, characterized in that The method for calculating the coupling evaluation value comprises the following steps: The Kappa coefficient and the score corresponding to the accuracy are standardized by Z-score to obtain the standardized score value Z k and Z c ; The standardized score value Z is calculated using the entropy weight method. s and Z h The weight N s and N h ; According to Z s and Z h and its weight N s and N h , build a coupled evaluation model for accuracy; The coupling evaluation values ​​of the first change detection result and the second change detection result are respectively calculated according to the coupling evaluation model, and the change detection result with a high coupling evaluation value is taken as the final change detection result.

10. The method according to claim 8, characterized in that The coupling evaluation model is: Q = Z s N s +N h (A-Z h ), where A is the standardized index value Z h The maximum distance of the interval is rounded to an integer value.

Images