3D Buildings Change Detection from Aerial and Satellite Stereo Imagery Using Kullback–Leibler Divergence Algorithm

Monitoring city sprawls is considered an essential subject in urban planning. The most important object in the urban areas is the building; therefore, finding an automatic method for detecting the changes in the buildings is considered a priority task for researchers to consider the changes in a district, especially for assessing the damages during disasters and updating geo-database. However, using 2D images to detect changes is ineffective because of the various imaging environments and the parameters of the sensors. Furthermore, during the change detection process, it is difficult to distinguish between the building and other objects due to the similarity of spectral properties. Therefore, it is necessary to use stereo images for DSM generation and then find the changes. This paper proposes a Kullback–Leibler divergence (KLD) algorithm to detect urban area changes based on stereo imagery. Two DSMs have been obtained through the photogrammetric process using two different sensors. The first data set is based on the stereo aerial imagery captured in 2012 and the second stereo is from the worldview-2 sensor captured in 2017. Before applying the KLD algorithm, the aerial’s DSM, which has an original resolution of 0.3m, was resampled to 1 m to make it similar to the satellite’s DSM. Three study areas have been selected for the algorithms test, located in Erbil-Iraq. The assessment shows that the KLD detected changes better than other methods after removing the small fragments through the post-processing step. For the evaluation, the confusion matrix has been determined for each study area. The analysis demonstrates that the overall accuracy for the three study areas where 89.3%, 91.1% and 88.9 %, respectively.


INTRODUCTION
Monitoring Change detection of buildings is considered to be very significant for the decision-makers. It provides necessary information that helps to assess the changes in various applications such as disaster assessment, urban planning, cadastral map updating, and others. Although manually detecting the changes is considered very accurate and reliable but is considered to be very time-consuming and expensive. Therefore, the trend is currently focusing on the automatic and semiautomatic methods since it has proven the result rapidly. Generally, aerial or satellite imagery change detection can be classified into 2D and 3D change detection [1]. The common change detection techniques were dominated by 2D image-based algorithms based on pixel intensity analysis of the imagery using radiometric values [1][2][3][4][5][6]. In addition to applying changes detection to the urban area, it has been widely used to monitor the changes in the vegetation and observer the cultivation [7,8]. On the other hand, various studies have been conducted using DSMs for change detection based on using 3D features of the objects [9][10][11][12][13]. Also, the combination of the 2D and 3D features has been approved to give good results in change detection to detect the changes by monitoring the newly constructed objects Huang [14].
However, from the above studies, it can be noticed many limitations will lead to minimising the detection accuracies, such as threshold values. It has been noticed that selecting the optimal threshold value is considered a challenge in the image detection method, and it has directly related to the amount of detected change in the pixels. Furthermore, important information will be omitted during the change detection process, due to the spatial relations between the neighbourhood of the pixel have not been taken into consideration [15]. However, an alternative method to solve the above limitation is to use object-based and information-based techniques. For instance, the information-based algorithm, the Kullback-Leibler Divergence method, is based on detecting the similarity between two probability density functions(PDF), which can be used in change detection [16].
Based on the literature, none of the studies is focused on using DSM in change detection using two different sensors (i.e. detecting the changes between the aerial and satellite data). Therefore, this study addresses the changes from different sensors using Kullback Leibler Divergence.
The structure of the paper is as follows, the introduction in this section. The applied methodology and used datasets are described in section 2. In section 3, the experiments have been illustrated, and the result evaluation has also been presented. In section 4, the discussion of the applied method has been presented. Finally, the conclusion and future improvement have been presented in section 5.

STUDY AREA AND DATASETS
An area of 5x5 km square with significant urban development in Erbil city is selected and acquired for the study. Within the study area, the methodology is evaluated in three different test areas: Area 1 (Ankawa), Area 2 (Dream City) and Area 3 ((sarbasty (32 park)). There are two groups of data in the data set: the first one is stereo satellite imagery from the world view-02 (WV02) sensor, captured in 2017 with a 0.5m GSD resolution as shown in Figure 1, while the second set is aerial imagery of Erbil city acquired in 2012 with a 0.10m GSD resolution by using an Ultracam Xp sensor. The Aerial images were provided by Salahaddin University's Geomatics (Surveying) Engineering Department.

Measuring Ground Control Points (GCPs)
For the purpose of processing the satellite and aerial images, Also for evaluating the accuracy of generated DSMs, a set of GCPs was measured by using a Leica-1200 GPS instrument. Through the use of the GPS post-processing kinematics (PPK) technique, the GCPs were observed. During the observation, a measurement was taken consciously at the base station for seven hours, and the required observation for each rover station was approximately three minutes.
The GCPs were acquired based on the projection Universal Transverse Mercator (UTM) zone 38N and datum WGS 84 ellipsoid. Afterwards, the GPS data was imported to the Leica geo-office software for post-processing purposes, and they were processed according to the nearest Continuously Operating Reference Station (CORS), which was based in Iraq Survey Erbil (ISER) as a source of GPS corrections. The result relative to the base station was thirty GCPs within millimetres accuracy.

The Reference Data
It is necessary to collect accurate ground truth, or reference data in order to assess the accuracy of a method, and also to determine a method's suitability and efficiency (Tian, 2013). In this research, the reference data was prepared for each study area by manual compilation of the buildings using summit evolution software. The stereo Aerial and satellite images were added and processed to extract boundaries of all changed and unchanged buildings, then the digitized reference data are imported to the GIS software (ESRI, Arc Map 10) for area computation of each changed and unchanged building. The computed reference data (shown in Figure 2 where the changed buildings are highlighted in Green) will be used for evaluating the accuracy of change detection algorithms.
a-test area one

2.4.Applying Change Detection
In this research, the KL-Divergence change detection algorithm has been used to identify the changes. The flow chart of the proposed building change detection methodology is presented in Figure 3. Initially, the stereo aerial and satellite images are used to generate the DSMs. Afterwards, the DSMs are used to apply the KL-Divergence algorithm for detecting 3D building changes. Figure 3 General workflow of information-based KL-Divergence method that has been implemented in this research.

Data Processing and DSM Generation
In this study, ERDAS Imagine software is used to process the stereo aerial and satellite images. The DSMs are generated by the image matching process using the least square matching algorithm. The least squares matching (LSM) algorithm minimizes gray value differences between the reference window and search window through an adjustment process, in which geometric (location, size, and shape of the search window) and radiometric (pixel gray values) corrections are determined for each matching window based on the affine transformation parmaters [17].
Generally, DSM co-registration is necessary before applying any change detection method to avoid false change detection caused by image displacement [1]. However, In the framework of this research, no image registration techniques were required, because the data sets have been derived using the same GCPs during the processing of the stereo images either from aerial or satellite sensors, for the purpose of producing DSMs and ortho-imagery.

Applying Kullback Leibler Divergence
KL divergence is a fundamental equation of information theory that measures a statistical similarity between the shapes of two probability distribution functions (PDFs), which has been applied widely in theory and practice fields [18]. The higher the similarity between the two PDFs, the lower the divergence. Assuming P(x) and Q(x) are two probability density functions. The KL-Divergence can assess the divergence between these two probability density functions as follow : In case if the two probability density functions are close to each other, the Kullback-Leibler divergence is small. In contrast, it is larger if there is a great deviation between them [19]. The KLD is always nonnegative and equals zero only if the two distributions are identical. In this research, for applying the KL-Divergence algorithm on generated DSMs, a specific code written in the MATLAB programming platform has been used. First, the Probability Density Functions (PDF)s of the pixel values are generated; then, the KL-Divergence of the PDFs are calculated as a measure of change, Afterwards, The obtained result after applying KLD was converted into the binary change map. However, some virtual changes that do not belong to the buildings have appeared as object differences also, which theys are removed through the following steps.
Step 1) Elimination of the Vegetation Since the main focus of this research is only on detecting the changes in the buildings, therefore the elimination of the vegetation from the 3D change detection scheme is essential. To separate the changes in the buildings from vegetation, the NDVI masks have been computed from the near-infrared and visible light reflected from stereo satellite images as it is given in equation (2). NDVI = (NIR -RED) / (NIR + RED) (2) Where: NIR and VIS stand for the spectral reflectance measurements acquired in the near-infrared and RED bands , respectively. Then NDVI mask was subtracted from the binary change map. Later, the vegetation was removed from the change map. But Spurious change in altitude can also be caused by other land covers or DSM computation errors displayed around the buildings that affected Building change extraction results. Therefore, more adaptive post-processing steps are required to keep only the real changes in the buildings.
Step 2) Applying Mathematical Morphology Spurious detected changes that are displayed as objects around the buildings need to be removed by applying some mathematical morphology(MM) operations, which are known as an opening followed by a closing operation. MM is a non-linear process that is commonly used in image analysis which is based on modifying the geometrical shapes within the image, rather than pixel values [20]. The basic morphological operations that are introduced are erosion; dilation; opening and closing. Erosion and dilation are the simplest operations in mathematical morphology. In practice, dilations and erosions are usually employed in pairs [21]. The result of the dilation and erosion combination is called opening or closing depending on which is applied first. In this research, MATLAB software has been used for the purpose of applying (opening and closing) operations.

EXPERIMENT AND RESULT EVALUATION
In this section, the proposed methodology in section 2 for 3D building change detection is performed. It is started with generation of DSM models from both aerial and satellite imagery followed by applying KLdivergence algorithm for building change detection and ended by evaluate/assess the produced results. Hence, results and their assessment have been reported in this paper as two following main sub-section.

DSM generation Evaluation
Photogrammetric approach has been used to generate DSM from both aerial and satellite imagery, Figure (4) presents the generated DSM result for three selected study areas. Basically, generation procedure is done by the image matching process using least square (LS) matching algorithm. Furthermore, a particular strategy has been conducted in setting the main parameters of the software such as: identifying the window size and correlation coefficient for image matching threshold value determination. It is important to mention that the topographic category type has been set as a (flat area) due to the location of the study area which is an urban area without rolling or hills. (note: I haven't included the details of setting those parameters as it was mentioned in the thesis because of the no. of pages limitation) As a result, two DSMs with 0.3m and 1m grid cell size from aerial photograph and satellite image respectively have been generated for all selected study areas (Figure 4 (a, c and e) presents the aerial-DSM for test area one, two and three respectively, while (b, d and f) presents the satellite-DSMs for the corresponding areas). basically, both generated DSM as a multiple raster datasets should be stored with the same cell resolution. Hence, the produced DSM from aerial images were resampled to 1 m grid cell size in order to have the same cell resolution with the one that generated from satellite image. As this study is focusing on 3D Building changes (i.e height change detection), Nearest neighbor technique is the method that has been conducted in resampling procedure because this method will keep the elevation values from the input raster dataset. That means the closest center in the output raster product is determined by the center from the input raster data, which in turn means no changing in the elevation values will produce.  The detailed point accuracy information is shown in Table 1.

KL-divergence Evaluation
KL-divergence algorithm applied on generated DSMs via conducting matlab code. Initially, pixel value's PDFs are generated followed by measuring of change through calculating the KL-divergence of those PDFs. Figure (5  Those obtained results have converted into the binary change map in order to apply other necessary procedures such as subtracting the NDVI mask to eliminate the vegetation, and also apply mathematical morphology operation to remove pseudo changes that could appear around building. As a result, the noises coming from DSM failures were considerably reduced as can be seen from figure (6).
a-test area one b-test area two c-test area three Figure 6 The final change map after applying KLD and morphological operations: (a) test area one, (b) test area two and (c) test area three.
In term of accuracy assessment, the effectiveness of the proposed method for the detection of the position and size of the changed buildings and the overall change situation has been evaluated by using confusion matrix for the object based assessment as presented in the Table 2. Basically, the obtained change map is assessed with respect to the reference (ground truth) data. In this assessment method, the changed buildings are treated as single objects without consideration of their shape and size. Only the effectiveness of the detection of distinct changed areas (buildings) in the change map is considered in the evaluation (Tian et al., 2014). Concerning building change detection accuracy, in some cases, the number of correctly detected building are more important. In the object-based evaluation, the confusion matrix is based on objects instead of pixels. Therefore, four parameters are measured to assess change detection result.
1. True Positivethe number of changed buildings correctly detected as changed.
2. False Positivethe number of unchanged buildings incorrectly detected as changed.
3. False Negativethe number of changed buildings incorrectly detected as unchanged 4. True Negativethe number of unchanged buildings correctly detected as unchanged.
Several certain quality measures can be derived from the above parameters to assess the performance of a change detection algorithm. Basically, the following quality measures are the most commons:

1) Branching Factor (BF)
The branching factor is the ratio of 'False Positive' and 'True Positive'. It represents the rate at which the incorrect objects are classified by the algorithm as changes.

2) Miss Factor (MF)
The miss factor is the ratio of 'False Negative' and 'True Positive'. It represents the rate at which the change objects in the reference data are missed by the algorithm in classification.

3) Completeness
Completeness is a quantity measure that defines the proportion of the total number of correctly classified changed objects to the total number of changed objects in the reference data.
Correctness is the fraction of total number of objects correctly identified by the algorithm as changed objects to the total number of changed objects in the reference data. Unlike completeness, it is a significant quality measure shows the correctness of the change detection.
The ideal values of completeness and correctness are 100%, which means that all changes are correctly identified, and it is achieved when FN = 0 in case of completeness or FP = 0 in case of correctness.

4) Quality Percentage (QP)
The quality percentage is the combination of both completeness and correctness to provide quality percentage for an accuracy measure. The ideal value of QP should be 100%.

5) Overall Accuracy (OA) and Kappa Accuracy (KA):
Overall accuracy is computed by dividing the total correct (i.e., the sum of the major diagonal) by the total number of objects (N) in the error matrix (Tian, 2013).
Overall accuracy alone may create ambiguity, because it takes into account the similar no change objects in both the reference data and the change map while the number of these pixels is presumably very high in most cases, even if one of the classes suffers adversely in classification. Therefore it is suggested to account for other quality assessment parameters for accuracy assessment such as (Varade, 2011).
However, overall accuracy that only takes into consideration the diagonal elements for the ratio in the confusion matrix from table 2, Kappa accuracy considers the non-diagonal elements. Therefore, for the confusion matrix, it is a measure of overall agreement.
represents the hypothetical probability of the agreement between the extracted result and reference data; N is the total number of observations in confusion matrix (Tian, 2013).
As a result, figure (7) presents the counted change and unchanged buildings in the reference data. Then table (3) and table (4) respectively presents the generated confusion matrix for the tree tested areas and the quality measures that used for analysis purpose. Figure 7 The red lines are the reference data ,overlain on aerial images to count the number of buildings for accuracy assessment

DISCUSSION
The proposed methodology performed to evaluate 3D building change detection from satellite and aerial photography confirm its effectiveness. Results of both generated DSMs and KL-divergence have detailed discussed in this section. Interestingly, Table 1, showed that the Aerial-DSM has better accuracy in term of RMSE, especially for the (test area one which is equal to 0.2275m and test area three =0.1826m), than the satellite-DSM, for the same areas RMSE were (0.9399m and 0.4762m) respectively. Furthermore, it is clear that the Aerial-DSM have a better accuracy, this is because that the aerial images have a higher GSD resolution than satellite images which is 0.10m and 0.50m respectively which in turn means that more detailed information about shape of the building with their edges are identified from aerial DSMs than the satellite-DSMs. On the other hand, table (4) records best result for test area 3 as the overall accuracy was 88.88 %. It is obvious that KLD algorithms was able to detect the changed buildings and constructed geometrical building shapes close to the original shapes. However, in the object based assessment method the changed buildings are treated as single objects without consideration of their shape and size. Moreover, an overall accuracy of (89.291%, 91.145% and 88.889 %) are achieved for test areas 1, 2 and 3 respectively. This in turn mean that KLD algorithm both can sufficiently distinguish between changed and unchanged buildings. However, when it has been compared with reference data it showed that are not sufficiently robust to be implemented alone. As in order to increase the accuracy and the rate of truly building change detection, further post-processing steps will need like removing unwanted error that may cause due to the DSM computation error and/or due to the different nature of the DSM sources.

CONCLUSION
This paper proposes an efficient algorithm for change detection in an urban area using different optical stereo datasets. The datasets are composed of multi-temporal DSMs extracted from aerial and satellite imagery. For detecting the changes, the KLD algorithm has been used. The KLD algorithm successfully detected all the height changes in the area, including the changes in the vegetation. To retain the buildings only, the NDVI index mask has been sued to remove the vegetation. Furthermore, the experiment results show that the mathematical morphology can be applied to improve the result by removing the false objects around the buildings that were produced due to vegetation mask submission and height subtraction. By comparing the result obtained from KLD with other methods such as (Pixle-Based image differencing), it can be noticed that the suggested method can achieve better results [22]. However, the algorithm mainly depends on the quality of the used DSMs and the building, in case the building's features are more accurate the errors will be less. Therefore, it is expected that the accuracy of change detection can be increased, in case a better stereo-matching algorithm is used for DSM generation such as Semi global matching. It is also possible to extend the suggested algorithm to detect the changes caused due to disasters such as earthquake effects and forest changes since it depends on introducing the height in addition instead of radiometric information.