Change Detection of High Spatial Resolution Images Based on Region-Line Primitive Association Analysis and Evidence Fusion

: Change detection (CD) remains an important issue in remote sensing applications, especially for high spatial resolution (HSR) images, but it has yet to be fully resolved. This work proposes a novel object-based change detection (OBCD) method for HSR images that is based on region–line primitive association analysis and evidence fusion. In the proposed method, bitemporal images are separately segmented, and the segmentation results are overlapped to obtain the temporal region primitives (TRPs). The temporal line primitives (TLPs) are obtained by straight line detection on bitemporal images. In the initial CD stage, Dempster–Shafer evidence theory fuses the multiple items of evidence of the TRPs’ spectrum, edge, and gradient changes, and obtains the initial changed areas. In the reﬁning CD stage, the association between the TRPs and their contacting TLPs in the unchanged areas is established on the basis of the region–line primitive association framework, and the TRPs’ main line directions (MLDs) are calculated. Some changed TRPs omitted in the initial CD stage are recovered by their MLD changes, thereby reﬁning the initial CD results. Di ﬀ erent from common OBCD methods, the proposed method considers the change evidence of TRPs’ internal and boundary information simultaneously via information complementation between TRPs and TLPs. The proposed method can signiﬁcantly reduce missed alarms while maintaining a low level of false alarms in OBCD, thereby improving total accuracy. In our experiments, our method is superior to common CD methods, including change vector analysis (CVA), PCA-k-means, and iterative reweighted multivariate alteration detection (IRMAD), in terms of overall accuracy, missed alarms, and Kappa coe ﬃ cient. method for very high resolution images that fuses the results from di ﬀ erent methods, including CVA, IRMAD, and iterative slow feature analysis (ISFA) through D–S evidence theory. Bitemporal HSR images were stacked for image segmentation. CVA, IRMAD, and ISFA were applied to obtain three candidate change maps. Then, these change maps were fused on the basis of D–S theory to obtain the ﬁnal CD result. This method explores the uncertainty among di ﬀ erent change results and generates a change map that is more accurate than that produced by individual methods. CVA PCA–k-means ﬁxed


Introduction
Change detection (CD) identifies differences in an object or phenomenon by observing it at different times [1]. CD is an important technical step in many applications, such as damage assessment, The existing research verifies the feasibility and application value of OBCD. However, OBCD only considers objects' (segments') "internal" information, and may thus lose important CD evidence along object boundaries. Typically, changes in artificial objects, e.g., roads and buildings, are often related to the disappearance or appearance of structural straight lines that strongly denote the objects' shapes and configurations. Some lines are located on object boundaries and are highly informative in delineating objects when their tones are similar to the surrounding background. Thus, using such structural information, i.e., straight lines, is helpful in OBCD.
In previous studies [33,34], Wang et al. proposed the region-line primitive association framework (RLPAF), which uses structural lines in OBIA. In this technical framework, straight edge lines act as additional object primitives, i.e., line primitives (LPs), for succeeding feature extraction and analysis. The RLPAF conducts information complementation between segments, i.e., region primitives (RPs) and LPs, expands OBIA's feature set, and enhances OBIA performance. It has been validated in the extraction of several spatial objects, e.g., raft cultivation areas and photovoltaic panels, from HSR images [35,36].
In the present study, a novel OBCD method originating from the idea of the RLPAF is proposed. In this method, pieces of change evidence of temporal region primitives (TRPs) are fused on the basis of evidence theory for the initial CD. The change evidence of the temporal line primitives (TLPs) inside the TRPs and that along the TRPs' boundaries are then supplemented such that they offer auxiliary CD evidence for the "weakly" changed TRPs omitted in the initial CD. The experiments show that the proposed method greatly reduces missed alarms while maintaining low false alarm rates, and is thus obviously superior to other methods. This remainder of this paper is organized as follows. Section 2 introduces the technical route of the multitemporal region-line association modeling and feature extraction originated and extended from the RLPAF. It also introduces the implementation of the proposed method. Section 3 is the experimental part, which includes method analysis and comparison. Section 4 is the discussion section. Finally, Section 5 summarizes this study and introduces future work.

Methodology
As shown in Figure 1, the proposed method is composed of two main stages: the creation of TRPs and TLPs, and CD based on RPs and LPs. The CD stage is further composed of substages of initial detection based on evidence fusion and CD refinement based on the RLPAF.
The bitemporal images are segmented separately, and the segmentation results are overlapped to obtain TRPs. The initial detection fuses the change evidence of an object's spectrum, edge, and gradient features without consideration of TLPs; hence, this process may lose weakly changed regions. Thus, in the refining CD stage, the RLPAF-based relationships between the TRPs and the TLPs are established within the unchanged regions of the initial CD. Then, missed alarms are revised on the basis of the RLPAF rules, and the entire CD process is then completed.

Object Primitive Creation and Change Feature Extraction
As shown in Figure 2, the bitemporal images are first segmented separately using the hard boundary-constrained segmentation method [37]. Then, the bitemporal RPs are overlapped into a single layer. The bitemporal spectral, gradient, edge, and shape features are calculated to obtain the TRPs. In implementing OBCD, TRPs' bitemporal spectral and gradient features are quantized to form the frequency histograms-that is, the change feature vectors of TRPs. The histogram of edge pattern distribution describes the frequencies of different edge pixel distributions [38] and forms the TRPs' change feature vectors of edges. For multiband images, feature vectors of individual bands are stacked to obtain the whole feature vectors. Meanwhile, Burns' straight line detection method [39], which can extract low contrast lines, obtains straight lines from bitemporal images. Furthermore, line direction, length, and intensity are calculated to obtain the TLPs. The bitemporal images are segmented separately, and the segmentation results are overlapped to obtain TRPs. The initial detection fuses the change evidence of an object's spectrum, edge, and gradient features without consideration of TLPs; hence, this process may lose weakly changed regions. Thus, in the refining CD stage, the RLPAF-based relationships between the TRPs and the TLPs are established within the unchanged regions of the initial CD. Then, missed alarms are revised on the basis of the RLPAF rules, and the entire CD process is then completed.

Object Primitive Creation and Change Feature Extraction
As shown in Figure 2, the bitemporal images are first segmented separately using the hard boundary-constrained segmentation method [37]. Then, the bitemporal RPs are overlapped into a single layer. The bitemporal spectral, gradient, edge, and shape features are calculated to obtain the TRPs. In implementing OBCD, TRPs' bitemporal spectral and gradient features are quantized to form the frequency histograms-that is, the change feature vectors of TRPs. The histogram of edge pattern distribution describes the frequencies of different edge pixel distributions [38] and forms the TRPs' change feature vectors of edges. For multiband images, feature vectors of individual bands are stacked to obtain the whole feature vectors. Meanwhile, Burns' straight line detection method [39], which can extract low contrast lines, obtains straight lines from bitemporal images. Furthermore, line direction, length, and intensity are calculated to obtain the TLPs.

Initial Change Detection
For the evidence fusion-based OBCD, the structural similarities (SSIMs) [40] of the spectral, gradient, and edge features of TRPs at two time points are calculated first. These SSIM measures, namely, change evidence, are fused via D-S evidence theory to obtain the initial changes.

Initial Change Detection
For the evidence fusion-based OBCD, the structural similarities (SSIMs) [40] of the spectral, gradient, and edge features of TRPs at two time points are calculated first. These SSIM measures, namely, change evidence, are fused via D-S evidence theory to obtain the initial changes.

Feature Similarity Measure
SSIM comprehensively consider the mean, variance, and covariance of two variables and are defined as: where X and Y are two variables, and µ X , µ Y ; σ X , σ Y ; σ 2 X , σ 2 Y ; and σ XY are their mean, standard deviation, variance, and covariance, respectively. C1 and C2 are constants to prevent instability when the denominator is close to 0. A TRP's spectral SSIM is calculated using Equation (1), where X and Y are the bitemporal spectral feature vectors of the TRP. The TRP's gradient and edge SSIMs are obtained in the same way, in which X and Y are the bitemporal gradient and edge feature vectors, respectively.

Change Detection by Evidence Fusion
D-S evidence theory is a tool to model inaccurate and uncertain information [41]. We define a non-empty set U, i.e., the discriminative framework, that consists of a series of mutually exclusive and exhaustive elements. Proposition A in the question domain belongs to the power set 2 U . We define the basic probability assignment function (BPAF) m : 2 U → [0, 1] in 2 U and let: where m(A) represents the confidence of subset A in U with the current evidence. If A ⊂ U, then m(A) denotes a determined belief in A; if A = U, then m(A) denotes an uncertain assignment; if A ⊆ U and m(A) > 0, then A is called a focal element of m. D-S evidence theory combines different pieces of evidence with an orthogonal sum. We let m 1 , m 2 , . . . , m n be n BPAFs in 2 U ; their orthogonal sum is expressed as: and is further defined as: where: and k is the conflict degree of the evidence. Equation (4) is called Dempster's combination rule. BPAF is constructed on the basis of spectral, gradient, and edge SSIMs. The discrimination framework U for CD is defined as: where Y represents the changed classes, and N represents the unchanged classes. On the basis of this equation, BPAF is specifically defined as: where S i is the SSIM of a feature, and α i is the trust degree of evidence in the discrimination frame.

Change Detection Refinement Using RLPAF
The initial CD stage might lose some changed TRPs. These TRPs are rechecked on the basis of the relationships of their RLPAF with their contacted TLPs.
Given a region Q and its contacted straight line L, we define the direction operator set as: which denotes a subset extracted from Q that is located above, on, or below line L, as illustrated in Figure 3a. The topology operator set is as follows: Top(L, Q) = In(L, Q), Touch(L, Q), Out(Q, L), Pro j(L, Q) , where the first three operators denote a subset extracted from L, which is contained in, touched by, or outside Q. On the basis of these definitions, the region-line mutual conversion operators are defined, and they include region-to-region, region-to-line, and line-to-region conversions [33,34].  Some changed TRPs are marked as unchanged because of their low BPAF, which indicates that multiple instances of change evidence only from the TRPs might be insufficient. To this end, we consider TRPs' MLD changes for method refinement. Given that the line detection and directional sector division might be inaccurate, the refining stage relaxes the CD threshold for the compensation of possibly changed TRPs. This processing is more robust than a direct MLD comparison.
Algorithm 1 is the pseudocode of the two CD stages. The spectral, gradient, and edge BPAFs are calculated and fused via D-S theory to obtain the BN of each TRP (P). If P's BN is less than threshold T, we divide [−π/2, π/2] into 16 intervals and allocate the lines. Then, we select the average angles of the first and second sectors that contain the most lines as the first and second main line directions. Bitemporal MLDs are extracted for TRPs to compare the changes in line direction. Given that touched lines are involved in the MLD calculation, pieces of change evidence located inside objects and on object boundaries are comprehensively considered in OBCD.
Some changed TRPs are marked as unchanged because of their low BPAF, which indicates that multiple instances of change evidence only from the TRPs might be insufficient. To this end, we consider TRPs' MLD changes for method refinement. Given that the line detection and directional sector division might be inaccurate, the refining stage relaxes the CD threshold for the compensation of possibly changed TRPs. This processing is more robust than a direct MLD comparison.
Algorithm 1 is the pseudocode of the two CD stages. The spectral, gradient, and edge BPAFs are calculated and fused via D-S theory to obtain the B N of each TRP (P). If P's B N is less than threshold T, then P is included in the output results {P C }; otherwise, P's MLD 1 and MLD 2 are calculated. If MLD 1 is not equal to MLD 2 , then threshold T is relaxed to T 1 , and P is redetermined by T 1 for its belongings. These steps are looped until all the elements in {P} have been processed. Then, {P C } is output as the final CD result.

Algorithm 1. Two-stage change detection
Input: TRPs {P}, TLPs {L 1 }, and {L 2 }, Change threshold T, Scaling factor S Output: Changed TRPs {P C } For each P within {P}{ Calculate its spectral BPAF, gradient BPAF, and edge BPAF and fuse them to obtain B N If P's B N < T, put P to {P C } Else Obtain P's bitemporal MLD 1 and MLD 2 using its contacted lines extracted from

Experimental Results and Analysis
The proposed method was implemented and integrated into the RSFinder software system [34], which was developed by the authors on the basis of Microsoft Visual C++2010 and Geospatial Data Abstraction Library (GDAL) 1.8 library. Different experimental images captured by different sensors were selected for method validation and comparison to verify the characteristics and advantages of the proposed method. Table 1 Figure 4 shows the original bitemporal images of the three experimental areas.  For ALOS, the single-band panchromatic and pansharpened multispectral images using the Gram-Schmidt method [42] were used for method validation. For Gaofen-2, the pansharpened multispectral images were used. For image segmentation, the scale inputs for each experimental area were varied to minimize oversegmentation errors and avoid undersegmentation errors. For straight line detection, the gradient difference was fixed to 1, and the line length was 10.

Experimental Procedure
The parameters for the initial CD were determined by a few sample areas. Change threshold T was slightly different according to the image resolutions and ground field distributions, as shown in Table 1. Other parameters were fixed, and they included the following: (1) Trust degrees 1 , 2 , and 3 . These parameters affected the trust degree of each feature during evidence fusion. They were equal to 0.35, 0.85, and 0.65.
(2) Scaling factor S. This parameter affected the extraction of the weakly changed region during the refining CD stage. An excessively small value cannot detect weakly changed areas, whereas excessively large ones will bring obvious false alarms. The parameter was tuned and fixed to 1.5.
Three well-known CD methods, namely, CVA, PCA-k-means, and IRMAD, were selected for method comparison. CVA detects changes according to the magnitude and direction of multitemporal spectral vector changes. The CD threshold of CVA was manually tuned to 0.6 in our experiments, and optimal results were obtained. PCA-k-means has two parameters, namely, nonoverlapping blocks H (H was 3 in our experiments) and the dimensions of S (S was 3 in our experiments) of the eigenvector space. IRMAD is an iterative, multivariate CD method. Its maximum number of iterations was fixed to 50.
To conduct method evaluation, we used four common indicators: false alarm (FA ), missed alarm (MA), overall accuracy (OA), and Kappa coefficient. For ALOS, the single-band panchromatic and pansharpened multispectral images using the Gram-Schmidt method [42] were used for method validation. For Gaofen-2, the pansharpened multispectral images were used. For image segmentation, the scale inputs for each experimental area were varied to minimize oversegmentation errors and avoid undersegmentation errors. For straight line detection, the gradient difference was fixed to 1, and the line length was 10.
The parameters for the initial CD were determined by a few sample areas. Change threshold T was slightly different according to the image resolutions and ground field distributions, as shown in Table 1. Other parameters were fixed, and they included the following: (1) Trust degrees α 1 , α 2 , and α 3 . These parameters affected the trust degree of each feature during evidence fusion. They were equal to 0.35, 0.85, and 0.65.
(2) Scaling factor S. This parameter affected the extraction of the weakly changed region during the refining CD stage. An excessively small value cannot detect weakly changed areas, whereas excessively large ones will bring obvious false alarms. The parameter was tuned and fixed to 1.5.
Three well-known CD methods, namely, CVA, PCA-k-means, and IRMAD, were selected for method comparison. CVA detects changes according to the magnitude and direction of multitemporal spectral vector changes. The CD threshold of CVA was manually tuned to 0.6 in our experiments, and optimal results were obtained. PCA-k-means has two parameters, namely, nonoverlapping blocks H (H was 3 in our experiments) and the dimensions of S (S was 3 in our experiments) of the eigenvector space. IRMAD is an iterative, multivariate CD method. Its maximum number of iterations was fixed to 50.
To conduct method evaluation, we used four common indicators: false alarm (FA), missed alarm (MA), overall accuracy (OA), and Kappa coefficient. where: TP is the number of change image objects correctly detected, FP is the number of unchanged image objects incorrectly detected as changed ones, TN is the number of unchanged image objects correctly detected, and FN is the number of changed image objects incorrectly detected as unchanged ones.

Method Performance Analysis
Reference maps were obtained by using the visual interpretation of the changed TRPs in the experimental areas. Manmade object changes, which could be visually identified, were included in the reference maps. The increase and decrease in vegetation, water, and other natural objects, which were visually identifiable (normally with the changed areas larger than 100 pixels), were marked as changes. TRPs with changed areas larger than 40% were classified as changed segments. The same and unified standards were adopted in all the experiments. Two sets of statistical results were obtained by counting the segments or pixels, as shown in Tables 2-4. Overall, the results were similar. Therefore, the succeeding analyses are mainly based on the segment-based results.
3.2.1. Experimental Area 1 Figure 5 shows the results of experimental area 1. From this figure, we can see that the proposed method had the lowest missed alarms and false alarms, and that it obtained the best results in experimental area 1. This evaluation was further verified, as shown in Table 2.
As shown in Figure 5a,b, changes mainly occurred in the building and road sections in this area. CVA caused a large number of false alarms and missed alarms, as shown in Figure 5g. PCA-k-means and IRMAD had fewer false alarms, but more serious missed alarms compared with CVA, as shown in Figure 5h,i. The proposed method had the least false alarms and missed alarms, as presented in Figure 5l. The multifeature evidence fusion made the method robust to false alarms in the initial CD, and the region-line associated analysis detected some weak changes in the refining stage, resulting in the method's superior performance.
As shown in Table 2, CVA's OA and Kappa were the lowest. PCA-k-means and IRMAD had fair FAs, but their MAs were high. IRMAD was the best among the three methods used for comparison, with its OA, FA, and Kappa values being equal to 93.65%, 4.13%, and 0.47, respectively. These measures of the proposed method were 97.68%, 0.87%, and 0.77, and were thus better than those of IRMAD. CVA had the best MA of 28.77% among the three methods, whereas the proposed method's MA was 26.03%, which was better than that of CVA. In general, the proposed method had the highest OA and Kappa and lowest MA and FA among all the methods, indicating its superior performance. Figure 5k illustrates the results of directly relaxing the scale factor without considering MLD changes. Directly relaxing caused many false alarms while reducing missed alarms, as observed in Figure 5. As shown in Table 2   By contrast, the MLD-based refining significantly reduced MA and prevented a significant FA increase. As shown in Table 2, MA dropped from 32.88% to 26.03%, whereas FA only increased from 0.63% to 0.87%. Accordingly, OA and Kappa increased from 97.52% to 97.68% and from 0.74 to 0.77, respectively. All these results indicated the MLD constraint, effectively preventing possible false alarms and increasing the overall method accuracy relative to the case of direct threshold relaxation.  Figure 6 shows the results of experimental area 2, and Table 3 presents the detection precision in area 2. From Figure 6 and Table 3, we can see that the proposed method also obtained the best results in this experimental area.
In this area, CVA performed the worst with serious missed alarms and false alarms, as shown in Figure 6. PCA-k-means had few missed alarms, but serious false alarms. IRMAD had few false alarms, but serious missed alarms. Similar to the results in area 1, IRMAD was the best among the three methods used for comparison in this area, with the OA, FA, and Kappa values being equal to 82.64%, 10.99%, and 0.22, respectively. PCA-k-means had the best MA of 62.66% among the three methods.
The false alarms of the proposed method increased in this area, but it still performed the best among all the methods, with its OA reaching 95.32%, its MA and FA respectively reaching 12.53% and 3.30%, and its Kappa reaching 0.79, indicating its superior performance. Figure 6k illustrates the results of directly relaxing the scale factor without considering MLD changes. As shown in Figure 6 and Table 3, directly relaxing caused many false alarms while reducing missed alarms, OA dropped from 94.39% to 87.65%, and Kappa dropped from 0.71 to 0.58 in area 2.
By comparison, the MLD-based refining significantly reduced MA and prevented significant FA increase. In area 2, MA decreased by 21.41%, whereas FA only increased by 1.77%. Accordingly, OA and Kappa increased from 94.39% to 95.32% and from 0.71 to 0.79, respectively. All these results indicated the MLD constraint, effectively preventing possible false alarms and increasing the overall method accuracy.  Table 4 are the results and detection precision of experimental area 3, respectively. In this area, changes occurred in the buildings and some water sections. CVA and IRMAD had serious missed alarms and false alarms, as shown in Figure 7 and Table 4. PCA-k-means had few false alarm errors, but had obvious missed alarms. The proposed method still obtained the best results among all the methods.  Figure 7k illustrates the results of directly relaxing the scale factor without considering MLD changes. Similar to the results of areas 1 and 2, directly relaxing caused many false alarms while reducing missed alarms. By contrast, the MLD-based refining significantly reduced MA and prevented significant FA increase. Compared with the result of initial detection in this area, the refinement caused OA to decrease from 94.91% to 94.58%. Meanwhile, the refinement caused the OA  Figure 7k illustrates the results of directly relaxing the scale factor without considering MLD changes. Similar to the results of areas 1 and 2, directly relaxing caused many false alarms while reducing missed alarms. By contrast, the MLD-based refining significantly reduced MA and prevented significant FA increase. Compared with the result of initial detection in this area, the refinement caused OA to decrease from 94.91% to 94.58%. Meanwhile, the refinement caused the OA measured by pixels to increase from 96.08% to 96.16%. OA in this area varied slightly. Kappa significantly increased from 0.55 to 0.62 (by segments) and 0.67 (by pixels). All these results indicated the MLD constraint, effectively preventing possible false alarms and increasing the overall method accuracy relative to the case of direct threshold relaxation. Figures 8-10 show some enlarged areas, which illustrate the characteristics of the proposed method. Initial CD caused some missed alarms in the three experimental areas.

Further Analyses
Remote Sens. 2019, 11, x FOR PEER REVIEW 16 of 23 measured by pixels to increase from 96.08% to 96.16%. OA in this area varied slightly. Kappa significantly increased from 0.55 to 0.62 (by segments) and 0.67 (by pixels). All these results indicated the MLD constraint, effectively preventing possible false alarms and increasing the overall method accuracy relative to the case of direct threshold relaxation. Figures 8-10 show some enlarged areas, which illustrate the characteristics of the proposed method. Initial CD caused some missed alarms in the three experimental areas. As shown in Figure 8, the buildings and surroundings in the panchromatic image were difficult to distinguish, resulting in missed alarms in the initial CD. However, the MLD changed significantly from vertical to near horizontal with the emergence of buildings. Thus, the missed building was recovered in the refining stage, as shown in Figure 8g,h.

Further Analyses
Missed and false alarms existed in all the comparative methods, as shown in Figure 9. Given that the bitemporal gradient and edge features were similar in some changed areas, the proposed method lost some changes in the initial CD stage, as shown in Figure 9g. However, because many horizontal lines emerged along with the building construction, the MLD changed significantly in the second image. Thus, many missed alarms were compensated for in the refining stage, as shown in Figure 9h. Figure 10 shows a changed area from bare land to buildings. Tones were close between the bitemporal images, easily causing errors in the initial CD stage. However, lines were cluttered in the first image. In the second image, MLD was approximately horizontal along with the construction of buildings. Given the obvious MLD changes, the missed TRPs were recovered in the refining stage, as shown in Figure 10g,h.
Meanwhile, CVA had obvious false alarms in Figures 8 and 9, whereas obvious missed alarms were found in Figure 10 with the same CD settings. PCA-k-means had obvious false alarms in Figure  9 and missed alarms in Figures 8 and 10. IRMAD had obvious missed alarms in Figure 8 and obvious false and missed alarms in Figures 9 and 10. The proposed method showed the most stable performance among all the methods.

Discussion
In this section, method efficiency, parameters, and possible improvements are further discussed.

Method Efficiency
The main steps of the method include TRP and TLP creation, feature similarity calculation, CD by evidence fusion, and CD refinement using RLPAF. The time consumption of each step is shown in Table 5. The experiments were carried out on Windows 7 64-bit OS with a CPU (Intel Core i7-4790, 3.60 GHz), RAM (8 GB), and a GPU (NVIDIA GT 630, 2 GB). Table 5. Method efficiency (unit: seconds). As shown in Figure 8, the buildings and surroundings in the panchromatic image were difficult to distinguish, resulting in missed alarms in the initial CD. However, the MLD changed significantly from vertical to near horizontal with the emergence of buildings. Thus, the missed building was recovered in the refining stage, as shown in Figure 8g,h.
Missed and false alarms existed in all the comparative methods, as shown in Figure 9. Given that the bitemporal gradient and edge features were similar in some changed areas, the proposed method lost some changes in the initial CD stage, as shown in Figure 9g. However, because many horizontal lines emerged along with the building construction, the MLD changed significantly in the second image. Thus, many missed alarms were compensated for in the refining stage, as shown in Figure 9h. Figure 10 shows a changed area from bare land to buildings. Tones were close between the bitemporal images, easily causing errors in the initial CD stage. However, lines were cluttered in the first image. In the second image, MLD was approximately horizontal along with the construction of buildings. Given the obvious MLD changes, the missed TRPs were recovered in the refining stage, as shown in Figure 10g,h.
Meanwhile, CVA had obvious false alarms in Figures 8 and 9, whereas obvious missed alarms were found in Figure 10 with the same CD settings. PCA-k-means had obvious false alarms in Figure 9 and missed alarms in Figures 8 and 10. IRMAD had obvious missed alarms in Figure 8 and obvious false and missed alarms in Figures 9 and 10. The proposed method showed the most stable performance among all the methods.

Discussion
In this section, method efficiency, parameters, and possible improvements are further discussed.

Method Efficiency
The main steps of the method include TRP and TLP creation, feature similarity calculation, CD by evidence fusion, and CD refinement using RLPAF. The time consumption of each step is shown in Table 5. The experiments were carried out on Windows 7 64-bit OS with a CPU (Intel Core i7-4790, 3.60 GHz), RAM (8 GB), and a GPU (NVIDIA GT 630, 2 GB).
TRP and TLP creation took the most time among these steps, because the bitemporal images needed to be segmented separately and straight lines were detected in both times. In the CD refinement step using RLPAF, the modeling and calculation of MLDs took the second-longest amount of time. However, image complexity influenced the method efficiency. Compared with that in area 1 with a single-band image, the feature similarity calculation of the multiband images in areas 2 and 3 took a significantly longer time. Meanwhile, the time consumption of the CD refinement step in area 3 was significantly longer than those in areas 1 and 2, because the TLPs in area 3 were more densely distributed.

Parameter Influence
To analyze the influences of parameters on method accuracy, we modified the change threshold, scaling factor, and the trust degrees of the spectral, edge, and gradient features in every step while keeping the other parameters fixed. Figure 11 shows the influence when the trust degree changes. Kappa increased gradually and reached the maximum values when the spectral trust degree increase from 0.20 to 0.35, as shown in Figure 11(a1-a3). The results indicated that the spectral feature improved the method accuracy when its weight (trust degree) was within a certain range. However, it reduced the method accuracy when the weight was unbalancedly large. The gradient and edge features also showed the same performance, as illustrated in Figure 11b,c. This evidence indicated that none of these individual features could be fully trusted, and that synthetical feature utilization was thus necessary. TRP and TLP creation took the most time among these steps, because the bitemporal images needed to be segmented separately and straight lines were detected in both times. In the CD refinement step using RLPAF, the modeling and calculation of MLDs took the second-longest amount of time. However, image complexity influenced the method efficiency. Compared with that in area 1 with a single-band image, the feature similarity calculation of the multiband images in areas 2 and 3 took a significantly longer time. Meanwhile, the time consumption of the CD refinement step in area 3 was significantly longer than those in areas 1 and 2, because the TLPs in area 3 were more densely distributed.

Parameter Influence
To analyze the influences of parameters on method accuracy, we modified the change threshold, scaling factor, and the trust degrees of the spectral, edge, and gradient features in every step while keeping the other parameters fixed. Figure 11. Influences of different trust degrees of features on method accuracy. (a1-a3) Influences of spectral trust degree in three experimental areas; (b1-b3) gradient trust degree; (c1-c3) edge trust degree. Figure 11 shows the influence when the trust degree changes. Kappa increased gradually and reached the maximum values when the spectral trust degree increase from 0.20 to 0.35, as shown in Spectral features were generally more vulnerable to temporal changes than gradient and edge features. Kappa obtained the largest values when the trust degrees of spectral, gradient, and edge features reached 0.35, 0.85, and 0.65, respectively, as shown in Figure 11. These results indicated that the spectral features took a less important role than the gradient and edge features in the proposed method, and were helpful in improving the method's robustness. Figure 12 shows the method accuracy curves when modifying the change threshold and scaling factor. Kappa first increased and then decreased along with the increase in change threshold, as shown in Figure 12(a1-a3). Many missed alarms existed when the change threshold was small, and clear false alarms appeared with large thresholds. Similar Kappa increases and decreases were found when the scaling factors increased, as shown in Figure 12(b1-b3). In the stage of CD refinement using RLPAF, the change threshold was relaxed to detect the weakly changed TRPs. However, the unchanged TRPs were considered changed when the change threshold was excessively relaxed, resulting in additional false alarms and reduced method accuracy. In this study, the method inputs were tuned using small pieces of training areas and were validated as appropriate. Figure 11a1-a3. The results indicated that the spectral feature improved the method accuracy when its weight (trust degree) was within a certain range. However, it reduced the method accuracy when the weight was unbalancedly large. The gradient and edge features also showed the same performance, as illustrated in Figure 11b,c. This evidence indicated that none of these individual features could be fully trusted, and that synthetical feature utilization was thus necessary. Spectral features were generally more vulnerable to temporal changes than gradient and edge features. Kappa obtained the largest values when the trust degrees of spectral, gradient, and edge features reached 0.35, 0.85, and 0.65, respectively, as shown in Figure 11. These results indicated that the spectral features took a less important role than the gradient and edge features in the proposed method, and were helpful in improving the method's robustness.  Figure 12 shows the method accuracy curves when modifying the change threshold and scaling factor. Kappa first increased and then decreased along with the increase in change threshold, as shown in Figure 12a1-a3. Many missed alarms existed when the change threshold was small, and clear false alarms appeared with large thresholds. Similar Kappa increases and decreases were found when the scaling factors increased, as shown in Figure 12b1-b3. In the stage of CD refinement using RLPAF, the change threshold was relaxed to detect the weakly changed TRPs. However, the unchanged TRPs were considered changed when the change threshold was excessively relaxed, resulting in additional false alarms and reduced method accuracy. In this study, the method inputs were tuned using small pieces of training areas and were validated as appropriate.

Further Discussion on Method Performance
In CD engineering applications, missed alarms must be recovered by visually comparing the entire bitemporal image, whereas false alarms only need to be validated within the detected changed areas. Thus, the intention of the application is to prevent missed alarms more than false alarms. On the basis of RLPAF, straight line changes were introduced as new OBCD evidence, and weakly changed TRPs were extracted. Thus, missed alarms were significantly reduced when low levels of false alarms were maintained, thereby improving the method's practicability.
However, the MLD-based constraint still left a small number of false alarms in the natural areas, typically in vegetation, because trivial changes in vegetation appearance also caused line distribution changes. Further determining the specific types of changes may be useful. These issues could be solved by integrating image classification techniques, which could separate manmade and natural

Further Discussion on Method Performance
In CD engineering applications, missed alarms must be recovered by visually comparing the entire bitemporal image, whereas false alarms only need to be validated within the detected changed areas. Thus, the intention of the application is to prevent missed alarms more than false alarms. On the basis of RLPAF, straight line changes were introduced as new OBCD evidence, and weakly changed TRPs were extracted. Thus, missed alarms were significantly reduced when low levels of false alarms were maintained, thereby improving the method's practicability.
However, the MLD-based constraint still left a small number of false alarms in the natural areas, typically in vegetation, because trivial changes in vegetation appearance also caused line distribution changes. Further determining the specific types of changes may be useful. These issues could be solved by integrating image classification techniques, which could separate manmade and natural objects, and classify change types into the proposed method. To this end, deep learning-based techniques, which can perform automatic feature learning and feature expression, have shown better accuracy than traditional methods. Hence, deep learning-based techniques will be considered for a classification-involved OBCD in our future work.

Conclusions
In this study, a pioneering OBCD method for HSR images is proposed on the basis of the analysis of region-line primitive association and evidence fusion. The steps include obtaining TRPs by using multitemporal image segmentation, acquiring TLPs via multitemporal straight line detection, conducting multifeature initial OBCD using evidence theory, calculating MLD features through RLPAF, and CD refinement with MLD change analysis. The proposed method conducts OBCD on TRPs and TLPs, and is thus distinct from common OBCD methods. The advantages of the proposed method include the following aspects: (1) Multifeature fusion in the initial CD stage obtains fair method accuracy. Tone differences between bitemporal images generally exist after image radiation correction. Through multifeature fusion, the proposed method reduces such influence by introducing more robust change evidence than image tones. The results include few false alarms and improved method accuracy.
(2) RLPAF extends the OBIA feature set, such as the feature subsets of line and region-line association [33,34], and thus offers new and effective information for OBCD. In common OBCD methods, feature extraction and analysis are all based on region primitives-that is, segments. In such technical scenarios, change clues along the boundaries of objects are easily ignored and might cause missed detection errors. The internal and tangent TLPs related to the TRPs are synthetically involved in the proposed method. Therefore, change clues from the object's interior and boundary are considered to improve detection accuracy.
(3) In the refinement stage, CD is limited within the areas with distinctive MLD changes. Such region-line synthetical analysis significantly reduces the missed alarms and prevents a clear increase in false alarms. In real applications, missed alarms are more intended for prevention than false alarms. These characteristics improve its application values.
In all the experiments, the proposed method obtains the best accuracy measures, which verify its most stable performance among all the methods. However, a few false alarms exist typically in natural objects and vegetation areas. In the future, we will continually refine the method to improve its accuracy and practicability. Furthermore, the proposed method focuses on the change analysis of bitemporal images. Future studies will concentrate on the analysis of the changes in serial images to investigate the distribution patterns of changes in time and space, describe and explain the change rules, and predict future changes.