Deep Unsupervised Learning for 3D ALS Point Cloud Change Detectioniris.de-gelis@irisa.frIris de Gélis is also with the Chair of Data Science in Earth Observation as a Beyond Fellow in the International Future Lab AI4EO, Technical University of Munich (TUM), Germany. Her work is also partly funded by the CNES, Toulouse, France.

Inspired by Saha et al. 2022, we propose a self-supervised approach that does not require complementary data to train the feature extraction network. While in Saha et al. 2022, self-supervised learning idea is based on learning to extract similar features from very different SAR and optical acquisitions from a same scene, we thought the variation in 3D points distribution may also be an advantage. Let us note that even in unchanged parts, 3D PCs may have different distributions due to the various acquisition plans, sensors, weather conditions, etc. Although differences in distributions make the direct comparison of PCs impossible, this property can be an asset for training a network to predict similar attributes over an unchanged area regardless of distribution.

This is the idea of the self-supervised part. The network is trained using three different losses on an unlabeled training set from the same two campaigns of acquisition as the testing set. At each iteration, the back-propagation of the gradient is made using alternatively one of the three losses. Thereby, in each iteration, a batch of $\mathcal{B}$ tiles of the older PC, denoted as $\mathcal{X}_{1}=\left\{x_{1}^{1},\ldots,x_{1}^{\mathcal{B}}\right\}$, and the corresponding $\mathcal{B}$ tiles of the newer PC, $\mathcal{X}_{2}=\left\{x_{2}^{1},\ldots,x_{2}^{\mathcal{B}}\right\}$, are independently given to the network, resulting in features $y$:

where $y_{1}^{b}$ and $y_{2}^{b}$ have the dimension $N_{1}^{b}\times K$ and $N_{2}^{b}\times K$, respectively. We recall that $N_{1}^{b}$ and $N_{2}^{b}$ are the number of points in the corresponding tiles. $K$ refers to the dimension of the output which is in practice the number of desired clusters (see the deep clustering loss below).

As for losses, we alternatively use three different terms, as illustrated in Figure 2. The first one is based on the deep clustering principle to force the network to learn discriminative features. Deep clustering relies on a pseudo-label assignment which will be used to train the network (Caron et al. 2018). In this study, pseudo-labels are obtained for each point by taking the argument of the maxima as the output of the network as in Saha et al. 2022. For example, for each point $i$ of the tile $x_{1}^{b}\in\mathcal{X}_{1}$, the corresponding pseudo-label $c_{1,i}^{b}$ is defined as:

where $K$ is the number of clusters, which is a hyper-parameter to fix. It can be associated with the number of semantic classes to segment in a single PC (note that this does not concern the number of classes of change between two PCs). However, intuitively, if $K$ is small, learned features will not be discriminative enough as large sets of points will be classified in the same class. On the contrary, with excessively high value, features will be too precise, and no generalization will be possible. In this study, it has been set empirically.

Based on these pseudo-labels, two deep clustering losses $\mathcal{L}_{1}$ and $\mathcal{L}_{2}$ are defined as the cross-entropy between $(c_{1},y_{1})$ and between $(c_{2},y_{2})$ respectively. The average of these two terms is taken to modulate weights:

However, with such losses, we observed that the network was collapsing and predicting all the points in a single cluster. To prevent this obvious solution, a weighting of the cross-entropy losses was done by applying the following weights:

for each cluster $k\in\left\{1,\ldots,K\right\}$, $C_{k}$ being the number of points in the cluster $k$. $\alpha$ is fixed to $K{\sum}_{h=1}^{K}C_{h}$ as done in the public implementation of KP-FCNN in the Torch-Points3D framework (Chaton et al. 2020). Weights are recomputed at each epoch. Intuitively, the deep clustering loss enables the network to learn discriminative features to be able to segment each point into clusters.

In addition to these clustering losses, we add a temporal consistency loss whose rule is to push the network to make similar predictions for tiles from different times but at similar places. As a matter of fact, even in unchanged areas, the point distribution in 3D point clouds differ from each other. Thus, by assuming that permanent changes between two dates are rare in proportion to the unchanged parts in urban areas, the temporal consistency loss enforces the network to make similar predictions for each point of the newer PC compared to the corresponding nearest point in the older PC. Therefore, both predictions ($y_{1}^{b}$ and $y_{2}^{b}$) are ordered before computing the loss:

The temporal consistency loss, $\mathcal{L}_{1,2}$, is then given by taking the mean of $l_{12,i}^{b}$ over all considered points ($p_{2i}^{b}\in\mathcal{X}_{2}$) of all tiles of the batch. Notice that this strategy of nearest point correspondence has already been successfully employed in supervised 3D PC change detection (de Gélis et al. 2023).

The third loss is a contrastive loss to encourage the network to produce dissimilar features for different tiles. As proposed in Saha et al. 2019, the contrastive loss is computed in a similar way to $\mathcal{L}_{1,2}$ by having previously randomly shuffled the batch $\mathcal{X}_{2}$ into $\mathcal{X}_{2}^{\prime}$ to obtain different tiles between $\mathcal{X}_{1}$ and $\mathcal{X}_{2}^{\prime}$. The loss $l_{12,i}^{b^{\prime}}$ is defined as follows:

Similarly to Saha et al. 2019, $\mathcal{L}_{1,2}^{\prime}$ is given by taking the mean of the exponential of the term $l_{12,i}^{b^{\prime}}$ over all considered points of all tiles in the batch $\mathcal{X}_{2}^{\prime}$:

where $N_{\mathcal{X}_{2}^{\prime}}$ is the number of points in the batch $\mathcal{X}_{2}^{\prime}$. Here, the exponential is added to avoid over-penalizing the network when $l_{12,i}^{b^{\prime}}$ is too far from 0. Indeed, even by shuffling $\mathcal{X}_{2}$, some areas can keep the same semantic, for example there might always be some ground points.

To summarize, the deep clustering loss makes it possible to learn discriminative features, the temporal consistency loss forces the network to predict similar features for similar areas regardless of the point distribution, and the contrastive loss avoids a trivial solution where all predictions are similar for both times. The overall process is given in Algorithm 1 and illustrated in Figure 2. This method is referred to as Self-Supervised Learning (SSL) in the following part.

Once a model is trained to perform a segmentation task, it can be used on both input PCs to extract features at different levels of abstraction and complexity depending on the layer. These extracted features can be used in order to highlight changes applying the Deep Change Vector Analysis (DCVA) principle, initially developed for 2D pixel change retrieving through deep features comparison (Saha et al. 2019). As shown in the comparison module of Figure 1, the point-wise change identification is realized by taking the magnitude of the difference ($\delta^{l}$) between feature vectors computed for each point of the newer PC, $p_{2i}\in\mathcal{P}_{2}$, with the nearest point of the older PC, $p_{1j}\in\mathcal{P}_{1}$. In other words, the feature difference $\delta^{l}$ is computed between features $\mathcal{F}_{1}^{l}$ and $\mathcal{F}_{2}^{l}$ for each PC, $\mathcal{P}_{1}$ and $\mathcal{P}_{2}$ respectively, according to the following equation:

with $f_{1j}^{l}\in\mathcal{F}_{1}^{l}$, $f_{2i}^{l}\in\mathcal{F}_{2}^{l}$, $p_{1j}\in\mathcal{P}_{1}$ and $p_{2i}\in\mathcal{P}_{2}$. The magnitude of the difference, also called deep feature magnitude coefficient, is obtained by taking the L2-norm of $\delta^{l}$. A threshold is applied on the deep feature magnitude coefficient to distinguish between changed and unchanged points. As in Saha et al. 2019, the threshold is selected using the unsupervised Otsu algorithm (Otsu 1979). Let us remark that Otsu thresholding is commonly used for change detection, for both 2D (Saha et al. 2019; Du et al. 2019; Zheng et al. 2021) and 3D (dos Santos et al. 2021; Marmol and Borowiec 2023). The trained network extracts similar features for two similar areas, thus the deep magnitude coefficient is close to zero in the unchanged part. The choice of the layer from which features are taken is a hyper-parameter to be set. Combined with SSL, the method is called SSL-DCVA.

We performed the same experiments on a real dataset consisting of two different ALS campaigns throughout the Netherlands. We extracted some tiles for the training and the test of the method among the publicly available dataset Actueel Hoogtebestand Nederland (AHN). This dataset consists of a total of four surveys throughout the entire country (Sande et al. 2010), allowing multi-date change extraction (Cserép and Lindenbergh 2023). The two last surveys (AHN3 and AHN4) have been semi-automatically annotated to assign a semantic label to each point. Five different classes are distinguished: ground, buildings, water, civil engineering structures (e.g., bridges) and clutter. This mono-date segmentation allowed us to manually derive a change ground truth for the test set. Concerning the changes, we decided to include in change areas the new or demolished building (in this case the building footprint on the ground points is marked as changed), new vegetation and new clutter. Since our focus is on 3D object changes, our manual annotation was ignoring changes in land cover if not characterized by modifications of the geometry. As such, changes like grassland to bare soil are only discernible using red-green-blue (RGB) color information, and do not imply 3D geometric changes. Besides, let us notice that the RGB data provided in AHN have not been acquired at the same date as LiDAR PCs implying numerous disagreement between 3D geometry and RGB information of changed scenes. Therefore, we explicitly discard RGB information.

As far as training PCs are concerned, they are sub-parts of tiles 31HN1_22, 31HN1_23, 31HZ1_04, from the divided AHN dataset (https://geotiles.nl/). Concerning the test set, qualitative analysis was performed on a sub-part of tile 37EN1_08 (see Figure 4(a,c)). This sub-part is selected to show sufficient change instances to properly evaluate the method. AHN data does not contain any change annotation. Thus, deriving change labels from available semantic labels is not obvious due to point cloud characteristics (e.g., no point to point corresponding). Thereby, we selected a surface of 12,400 m² from the tile 37EN1_08 and conducted a manual annotation. Even if here the change detection is only binary, some different types of changes are present in this area including new buildings, demolished buildings, new vegetation and clutter.

The density of AHN3 is about 12 points/m² while AHN4 is about 22 points/m². Height and planimetric stochastic errors are 5 cm. In addition to point coordinates, AHN data also includes LiDAR intensity and the number of returns. However, since here we focus on raw 3D PCs only, we rely solely on the 3D point coordinates to feed the network.

To assess our proposed self-supervised strategy to extract relevant deep features, we will compare our results with those obtained by our 3D PC DCVA adaptation using a network pre-trained on an annex task such as semantic segmentation using labels from a publicly available dataset as done in the 2D case in (Saha et al. 2019). Indeed, while datasets annotated according to the change are not common when dealing with 3D PCs, public datasets with a mono-date semantic annotation in urban environment are widely spread (Hackel et al. 2017; Roynard et al. 2018; Varney et al. 2020; Kölle et al. 2021). By choosing a public dataset as close as possible to the unlabeled change detection dataset to perform supervised training of the network, one could expect that extracted features will be consistent in unchanged areas and different in the case of changes. In our study, we use the Hessigheim 3D (H3D) Aerial LiDAR Survey to train the network for semantic segmentation. The H3D dataset consists of four different PCs at various dates and comes with labels related to 11 semantic classes that have been manually annotated (Kölle et al. 2021). In practice, the training is performed on H3D PCs acquired in March 2016 on behalf of the national mapping agency of Baden-Württemberg, Germany. This survey has a mean point density of about 20 points/m². In the following study, we refer to this method by Supervised Semantic Segmentation Training (SSST). Combined with DCVA, the method is called SSST-DCVA.

Additionally, even if no specific unsupervised deep method exists so far in the literature for 3D change detection, we decided to compare with the supervised Siamese KPConv network (de Gélis et al. 2023) trained on the simulated Urb3DCD dataset (de Gélis et al. 2021) and directly applied to AHN-CD testing set without any retraining. Thereby, this strategy referred as Siamese KPConv transfer is another unsupervised baseline as no label from the target dataset (i.e., AHN-CD) is used during the training. We also provide results achieved with a pure supervised baseline, training Siamese KPConv on the AHN-CD train set (with automatic annotation of change). To further benchmark our method against the state-of-the-art, we provide a comparison with Cloud-to-Cloud (C2C) (Girardeau-Montaut et al. 2005) and multi-scale model-to-model cloud comparison (M3C2) (Lague et al. 2013) distance-based methods. These methods constitute unsupervised baselines for 3D point-based binary change detection (Shirowzhan et al. 2019). In particular, to obtain final binary change information, a thresholding based on Otsu algorithm (Otsu 1979) is applied on the Hausdorff point-to-point distance computed in C2C. The traditional M3C2 method uses local surface normal and orientation to compute the 3D distance between two PCs (Lague et al. 2013). This method relies on statistical tests on distances between local surface normal and orientation features of the two PCs to automatically extract significant changes. In addition, we provide a comparison with an unsupervised machine learning method, namely a $k$-means, trained on 3D hand-crafted features proposed by Tran et al. 2018.

For the same computational reason as 2D images are divided into patches, original PCs are also divided into tiles to run deep learning experiments. As done in the supervised deep framework for 3D PC change detection from de Gélis et al. 2023, vertical 3D cylinders are chosen to make sure that whatever changes occur between the two dates, at least the ground is visible. Also, considering that the chosen back-bone model is KP-FCNN, which is relying on Kernel Point Convolution, a first sub-sampling rate has to be chosen (Thomas et al. 2019). Let us remark that the choice of the radius of cylinders and of the first sub-sampling rate has to be made in consideration of the scale of the sought changes. The thinner the first sub-sampling rate is, the more points from the original PC will be given to the network to extract useful features. The larger the radius is, the more context will be taken into account. Thereby a compromise between a large radius size and a thin sub-sampling rate has to be done due to memory limitations. We set the parameters following previous studies (de Gélis et al. 2023) and use the same first sub-sampling rate value. As far as the radius is concerned, we had to change the value for the AHN-CD dataset to fit the memory constraints.

Concerning SSL, the training is realized using cylinders of 20 m in radius and a first sub-sampling rate of 0.5 m. Then, the best results are obtained using 6 clusters for the deep clustering loss, and after 15 epochs of training. One hundred cylinders are used for each epoch, with a batch size of 10.

For the SSST part, KP-FCNN is trained using cylinders of 10 m in radius with a first sub-sampling rate of 0.2 m. A total of 6,000 cylinders are used for the training at each epoch. The batch size is 10. The 11 H3D classes are fused into 7 classes better suited for the target datasets. It requires around 90 epochs to converge.

Networks are optimized using a Stochastic Gradient Descent with a momentum of 0.98. The learning rate is set at 0.01 and decreases exponentially. These settings were chosen on the basis of the initial publication of KP-FCNN back-bone (Thomas et al. 2019). Regardless of the first sub-sampling rate, the final results are given at the original resolution. Indeed, a reprojection step is performed with regard to the nearest neighbor in order to obtain a final result for each point of the original PC.

For the choice of the layer to take features for the DCVA, several configurations have been tested. Knowing that KP-FCNN has 9 layers, the best results, reported here, are obtained using the 7th layer for SSST and the 8th for the SSL strategy.

Following the thresholding step, a cleaning of isolated predictions is realized to spatially smooth the results. Note that this cleaning step is systematically applied to ensure fair comparison.

Experiments were conducted with a single Graphic Processing Unit (GPU) (Nvidia Tesla V100 SXM2 16 GB).

Quantitative results are given in Table 1. The mean of accuracy (mAcc), the mean of intersection over union (IoU) and the IoU for both changed and unchanged classes are given. Corresponding qualitative results on the manually annotated test set are presented in Figure 3. To complete the qualitative analysis of the results, a larger scene has been visually inspected to understand the behavior of our methods in multiple conditions (see Figure 4).

As can be seen in Table 1, SSL-DCVA outperforms other methods, including SSST-DCVA. It is worth noting that despite its simplicity (no training is required), C2C provides relevant results. However, the point-to-point distance seems limited in places where for example trees have been replaced by a building of approximately the same height (see regions of interest in Figure 3(i)) or where a new building replaced an old one as in the top of zoom 1 where buildings in AHN3 and AHN4 are very different (Figure 4(a,d,h)). Conversely to C2C, M3C2 provides inconsistent results here (see Figure 3(h)): the ground elevation has changed slightly between the two acquisitions, so almost all areas are marked as changed. Notice that even when removing ground points for metric computation, the M3C2 method is lagging behind other methods (still about 10% of mAcc behind the SSST-DCVA algorithm evaluated under the same conditions). Moreover, in this study, we aim at detecting changes in object semantics (new buildings, demolition, new vegetation, etc.), so a change in the ground height is not of interest to us and has not been marked as changed in the ground truth. To further explain the relative poor results of M3C2, we recall that this method was originally developed in a geoscience context to detect changes at different scales, including centimetric (Lague et al. 2013). Thereby, it might be inadequate for urban environments, as already noticed in de Gélis et al. 2021. A distance-based method may not distinguish between topographic and semantic changes as long as the geometry of objects has changed. Learning-based methods (see Figure 3(d-g)) seem to also retrieve small changed objects, which is not possible with distance-based methods without including too many changes.

By looking at the results of the Siamese KPConv change detection network with transfer onto AHN dataset from Urb3DCD simulated dataset, we can see that performances are quite similar to the SSST-DCVA but are overtaken by SSL-DCVA method. In particular, some differences with the ground truth are visible in the demolished area and at object boundaries. This is probably due to the difference of building types present in the selected area of AHN dataset. Indeed, Urb3DCD dataset contains buildings from a french city center different from train and test areas, as for example, AHN3 data (time 1) contains a glasshouse. The same problem occurs with SSST-DCVA methods, since H3D PCs have different resolution and quality than AHN PCs. This shows the advantage of training directly on a dataset with similar properties to the test set and using recent developments in self-supervised learning. However, when compared to the supervised Siamese KPConv network, unsupervised methods can still largely be improved. The main differences of SSL-DCVA with the ground truth are visible on small objects such as vehicles, road signs or vegetation (see Figure 3(c) and (d)). Furthermore, as can be seen in the buildings on the left side of Figure 3 and right side of Figure 4(c), some omissions remain on new buildings with a flat roof. When looking at the mono-date segmentation of the PC realized before the DCVA step, one can see that flat roofs are classified in the same class as ground so, when comparing features, no changes are highlighted. This raises the difficulty of late-fusion change identification. Indeed, errors in the feature extraction step are propagated in the comparison step. Finally, some false detections are visible on the ground, forming a large trapezium (see Figure 3(d)). This is due to changes in the orientation of the ground surface.

Our method encounters difficulties in unchanged vegetated areas (see the top of zoom 2 in Figure 4(i)) certainly because of the complexity of LiDAR data in such areas with a high variation of point distribution even without changes in the semantics of objects. This results in a mixture of points predicted as changed and unchanged. Furthermore, these vegetated areas may have grown, and the acquisition not realized in the same season implies some differences on the 3D representation of trees. Note that the same problems occur with the other learning-based methods (see the top of zoom 2 in Figure 4(j,k). Looking at Figure 4(l), we can observe that C2C method is not better in this zoom where the vegetation has been removed. Indeed, in AHN3 some points are acquired from the ground to the top of the tree canopy thanks to the LiDAR sensor, thereby the point-to-point distance is not an efficient indicator for changes. Finally, it seems that SSL-DCVA (as well as SSST-DCVA) is more prone to commission than omission changes, while C2C shows the opposite behaviour (see Figure 3(d,e,i) or Figure 4(i,j,l)). From the user point of view, we believe that it is better to obtain more commissions than omissions. Indeed, as changed parts are rare compared to unchanged ones in general, it is faster to check errors and correct changed predictions than unchanged ones.

The point-to-point nearest neighbor correspondence lacks precision in the presence of occlusion in the 3D PCs. Indeed, due to the geometry of acquisition, some occlusions may appear in PCs, these hidden parts may not be similar in the two compared PCs leading to difficulties when comparing points in the DCVA part.

Once again SSL-DCVA seems more interesting than SSST-DCVA when looking at training time. SSST-DCVA takes about 17 hours to train on H3D dataset, while SSL-DCVA only requires about 9 minutes to train (see Table 1). The DCVA part on the manually annotated test set takes about 40 seconds.

Following our experimental assessment, we have identified two difficulties faced by our method. The first one is related to the hypothesis of rare changes required for the temporal consistency loss (Equation 8). Although already used in the literature (Saha et al. 2022), this assumption should be verified considering the training set of the studied dataset, whereas the test set should contain enough changes so that the thresholding operation is valid. For example, this assumption prevents us from applying the presented method on the Urb3DCD dataset (de Gélis et al. 2021) because it contains a high number of changed objects in the training set. Then, the other issue with our method comes from the DCVA part which relies on a point-to-point comparison based on the nearest point. This point comparison is not optimal in occluded parts as well as in dense urban areas. Let us note that this issue has been already mentioned when describing C2C misclassifications (see region of interest in Figure 3(i)). Even if DCVA comparison relies on multiple deep features, avoiding the problem when the two points being compared have different latent embeddings, the problem remains when the latent embedding of the two points under comparison is similar, meaning the same class is predicted by the back-bone network (whether trained by transfer learning or self-supervision). Indeed, in this case, the deep magnitude coefficient computed from deep features will be similar.

Finally, based on the state-of-the-art in unsupervised 2D image change detection, we proposed to adapt the SSL strategy developed in Saha et al. 2022 to 3D PCs change detection task. Our study showed the possibilities offered by SSL to tackle this task. However, there is still room for improvement. Indeed, different SSL strategies have been already developed in the literature for 3D point clouds understanding (Sauder and Sievers 2019; Xie et al. 2020; Alliegro et al. 2021; Chen et al. 2021; Zhang et al. 2021b) or 2D image change detection (Leenstra et al. 2021; Cai et al. 2021; Chen and Bruzzone 2022; Saha et al. 2022). SSL is vast, and many different strategies can be elaborated to train a neural network to extract interesting features. Further SSL studies can be conducted to tackle 3D PCs change detection task based on the existing literature in 3D PCs understanding and/or 2D image change detection. Besides, further development of pre-text tasks directly designed for 3D PCs change detection task would be relevant. For example, a possible improvement would be to incorporate a change-related task in the SSL training of the network, so that the learned features are directly related to change. Which will enable us to get rid of the nearest point comparison, one of the possible difficulties of our method.