1

Self-Supervised Multisensor Change Detection

Sudipan Saha, Member, IEEE, Patrick Ebel, and Xiao Xiang Zhu, Fellow, IEEE

Abstract—Most change detection methods assume that pre- To exploit the potential of deep learning while not using any
change and post-change images are acquired by the same sensor. training label or additional unlabeled images, transfer learning
However, in many real-life scenarios, e.g., natural disaster, it based CD methods are popular that reuse a pre-trained network
is more practical to use the latest available images before and
after the occurrence of incidence, which may be acquired using for bi-temporal feature extraction and comparison [1].
different sensors. In particular, we are interested in the combi- A striking feature of satellite data is its variability, in terms
nation of the images acquired by optical and Synthetic Aperture of different sensors. Images captured using passive optical sen-
Radar (SAR) sensors. SAR images appear vastly different from sor are quite similar to the natural images studied in the com-
arXiv:2103.05102v3 [[Link]] 23 Jan 2022

the optical images even when capturing the same scene. Adding to puter vision. However, images captured by the active sensors,
this, change detection methods are often constrained to use only
target image-pair, no labeled data, and no additional unlabeled e.g., Synthetic Aperture Radar (SAR) are remarkably different
data. Such constraints limit the scope of traditional supervised from the optical images [9], [10], [11]. While optical sensors
machine learning and unsupervised generative approaches for use wavelengths near to visible light (approx. 1 micron), SAR
multi-sensor change detection. Recent rapid development of self- uses a wavelength of 1 cm to 1 m. Moreover, optical sensors
supervised learning methods has shown that some of them can rely upon the natural illumination (e.g., sun) to create the
even work with only few images. Motivated by this, in this work
we propose a method for multi-sensor change detection using only brightness observed by the sensor, while the SAR sensors
the unlabeled target bi-temporal images that are used for training carry their own illumination source, in the form of radio
a network in self-supervised fashion by using deep clustering and waves transmitted by an antenna. Moreover, satellite images
contrastive learning. The proposed method is evaluated on four are captured with different number of spectral bands (one to
multi-modal bi-temporal scenes showing change and the benefits few hundreds), different spatial resolutions (few cm/pixel to
of our self-supervised approach are demonstrated.
Kms/pixel), different polarization. While this vast variation
Index Terms—Change Detection, Deep Learning, Self- provides an opportunity for detailed Earth observation, it is not
supervised learning, Multisensor analysis. trivial to use same set of methods for images from different
sensors. Due to this reason, most existing CD methods assume
I. I NTRODUCTION that the pre-change and post-change images are acquired using
the same sensor. The temporal frequency at which same sensor
Our Earth is rapidly changing, both due to natural and can image same place depends on the revisit period of the
man-made causes. Satellite image based change detection satellite on which the sensor is mounted. However, better the
(CD) is generally used to monitor the temporal evolution spatial resolution, more close the satellite is to the Earth, more
of the dynamic Earth [1], [2], [3], [4], [5], [6], [7]. CD time it takes to revisit same place. This is a hindrance in the
ingests bi-temporal images as input and segregates all pixels use of same-sensor CD in time-bound applications, e.g., fast
as changed/unchanged. CD is a crucial step for several ap- response for disaster management and precision agriculture.
plications, including disaster management, urban monitoring, Using different sensors may allow us to obtain temporal
forestry, glacier monitoring, and precision agriculture. Con- sequences with better temporal frequency without sacrificing
sidering the variation of applications, rarity of occurrences of spatial resolution. However, it is not trivial to process multi-
some change-inducing incidents (e.g., natural disasters), and sensor bi-temporal images as they are affected by the spectral
large geographic variation, it is imprudent to assume that large- characteristics of the sensors. Moreover, different sensors
scale training datasets corresponding to all such tasks can be capture different type of information, making their comparison
ever collected. Thus, there is a significant inclination in the often challenging [12]. The difficulty of this problem is further
CD literature towards methods that can process the target bi- accentuated by the fact that we are interested to detect change
temporal region-of-interest without using any training label without using any labeled training data or any abundant pool
or any additional pool of unlabeled image. Motivated by its of unlabeled data.
excellent performance in computer vision, researchers have The emergence of deep learning has seen many such prob-
applied deep learning to satellite image change detection [8]. lems solved that were thought to be very challenging in the
past [13], [14]. Self-supervised learning has shown remarkable
Sudipan Saha is with Data Science in Earth Observation, Techni-
cal University of Munich, Taufkirchen/Ottobrunn, Germany. E-mail: sudi- success recently, even when only few images are available
[Link]@[Link] [15]. Intrigued by this, in this paper we explore the challenging
Patrick Ebel is with Data Science in Earth Observation, Technical Univer- problem of change detection between optical and SAR images,
sity of Munich, Taufkirchen/Ottobrunn, Germany.
Xiao Xiang Zhu is with Remote Sensing Technology Institute, German the disparity between which is evident in Figure 1. We exploit
Aerospace Center (DLR), Weßling, Germany and also with Data Science in recent developments in the self-supervised learning and deep
Earth Observation, Technical University of Munich, Taufkirchen/Ottobrunn, clustering to propose a method for challenging SAR-optical
Germany.
Published in IEEE Transactions on Geoscience and Remote Sensing. Please CD where one of the bi-temporal images is acquired by an
visit the journal’s page for the official version. optical sensor, while the other is acquired by a SAR sensor.
 2

The proposed method requires only the bi-temporal target

scene (where change is to be detected), no training label, and
no additional unlabeled data. The target bi-temporal scene is
typically large, few hundred pixels by few hundred pixels.
Smaller bi-temporal patches (e.g., 64×64) are extracted from it
to train a two-branch network, similar to the Siamese network
[16]. Each branch of the network has a projection module and
a predictor. Projection modules learn features unique to optical
and SAR data without sharing weights, while predictors share
the weight. The output of the predictors are used to estimate (a) (b)
deep clustering loss for both images separately. Moreover,
Fig. 1. Visual contrast for Las Vegas between (a) optical image (prechange)
considering prior probability of changed pixels are much and (b) SAR image (postchange) . Optical and SAR images emphasize
less than the unchanged ones, a temporal consistency loss is different properties of the target area, thus performing CD on them is
proposed that ensures that pixels in the same location at two challenging.
different time tend to get same label. To ensure that this does
not lead the network to learn a trivial solution, a contrastive differencing, i.e., change vector analysis (CVA) [17]. A num-
loss is used. By combination of these losses, the proposed ber of superpixel and spatial neighborhood based variants of
method learns useful semantic feature from the multi-sensor CVA have been proposed, e.g., Parcel Change Vector Analysis
(SAR-optical) bi-temporal target scene and after training the (PCVA) [18] and Robust Change Vector Analsis (RCVA)
network predictions can be compared for change detection. [19]. Most deep learning based unsupervised change detection
The contributions of this paper are as follows: methods use transfer learning. [1] proposed deep change
1) We propose a self-supervised learning method for vector analysis (DCVA), a CD framework that combines ideas
change detection in a bi-temporal scene where one from CVA with feature extraction based on pre-trained neural
image is captured by the optical sensor and the other by networks. In nutshell, a deep model that has been trained
the SAR sensor. The proposed method, only exploiting for some other task is reused to obtain pixelwise bi-temporal
the available target unlabeled scene, effectively absorbs deep features from the target scene. Bi-temporal deep features
several concepts from the recent self-supervised learn- are then compared to obtain deep change hypervectors for
ing literature, e.g., deep clustering, augmented view, each pixel in the scene that are analyzed based on magnitude
Siamese network, and contrastive learning. By effec- (`2 norm) to identify the changed pixels. While [20] shows
tively exploiting these concepts and modifying them that sensor-specific pre-trained network is more suitable for
appropriately for the target multi-sensor bi-temporal transfer learning, [5] advocates models trained on ImageNet
data, proposed method is able to train a network that [21] for transfer learning in CD. There is another class of
is further used for bi-temporal comparison and change unsupervised CD methods that pre-classifies some pixels with
detection. high confidence as changed/unchanged using some traditional
2) We show the versatility of the self-supervised learning approach and further uses those confident samples for training
on spatio-temporal satellite data that is very different a CD model [22].
from typical computer vision images. Even though some It is not trivial to process multi-sensor bi-temporal images
form of aerial images (e.g., drone images) are often as they are affected by differences in the spatial resolution
studied in the computer vision, we stress that our satellite and differences in the spectral characteristics of the sensors.
data (both optical and SAR) are significantly different Due to this, there are very few works that can work in
from the typical aerial images. the setting where pre-change and post-change images have
3) We experimentally show the efficacy of the pro- different spatial resolution [23], [24] or bands with different
posed method on four different bi-temporal multi-sensor spectral characteristics [25]. Moreover, those works deal with
scenes. only minor variation in spatial or spectral characteristics. [23]
The rest of this paper is organized as follows. Related proposed a cycle-consistent generative adversarial network
works are briefly discussed in Section II. Section III outlines based method to learn transcoding between multi-sensor multi-
the proposed method. Datasets and experimental results are temporal domain. However, their work assumes that a large
detailed in Section IV. Finally, we conclude the paper in (unlabeled) area corresponding to both sensors are available
Section V. as training data. [26] used a symmetric convolutional cou-
pling network (SCCN) and [27] used denoising autoencoder
II. R ELATED W ORK (DAE) for CD in multisensor images. Though those works
In this Section we briefly discuss existing works on unsu- considered optical-SAR images, they applied their methods
pervised CD (with focus on the multi-sensor CD) and self- on scenes with limited spatial complexity. While our work is
supervised learning. strongly motivated from the existing works on multi-sensor
CD [23], [24], it takes them a step further by considering the
A. Change Detection challenging scenario of optical-SAR CD in complex urban
Prior to the emergence of deep learning, most unsupervised scenes and furthermore by integrating recent developments in
change detection methods used the concept of pixelwise image self-supervised learning.
 3

B. Self-supervised learning We can extract a set of bi-temporal patches of size R0 × C 0

Considering the difficulty of collecting labeled data and (R0 < R and C 0 < C) from the images X1 , Z2 . In practice,
abundance of unlabeled data, machine learning researchers one training iteration involves only a batch of B patches from
have focused on developing unsupervised and self-supervised X1 , denoted as X = {x11 , ..., xB1 } and corresponding patches
deep learning methods in the recent past. Gidaris et. al. [28] from Z2 , denoted as Z = {z21 , ..., z2B }. xb1 and z2b are processed
used image rotation as a pre-text task to learn unsupervised separately with deep clustering loss, as detailed in Section
semantic feature. Several other pre-text tasks have been ex- III-C. Furthermore, considering that xb1 and z2b represent same
plored in the literature, e.g., relative patch prediction [29] and location at two different times and prior probability of change
image inpainting [30]. Deep clustering, i.e., joint learning of is less, a temporal consistency loss (see Section III-D) is
the parameters of deep network and the cluster assignment of formulated using each such pair. Furthermore, Z is shuffled
the resulting features, has also been shown to be effective to form negative samples Z 0 and a contrastive loss is used
for unsupervised representation learning [31]. Remarkably, between pairs from X and Z 0 , as outlined in Section III-E.
[15] have shown that above-mentioned unsupervised methods Proposed method is outlined in Figure 2.
learn useful semantic features even with a single-image input.
Contrastive methods function by bringing the representation A. Bi-temporal patches are multiple views of the same location
of different views of the same image (‘positive pairs’) closer
while spreading representations of different images (‘negative We recall from Section II-B that many self-supervised
pairs’) apart [32], [33], [34]. Boostrap your own latent [35] learning approaches build upon the concept of bringing closer
and its variant SiamSiam [16] eliminate the requirement of the representation of the multiple views of the same image.
negative pair by using multiple views of the same image. In Different views of the same image are generally obtained
more details, SiamSiam [16] ingests as input two randomly by different augmentation techniques, e.g., random crops.
augmented views of an image and processes it through a We argue that multi-sensor bi-temporal patches xb1 and z2b
Siamese architecture. Each Siamese branch consists of an can be similarly thought to be multiple views of the same
encoder and a prediction head. The encoders share weight location. They represent augmentation of the same place,
between two views. where the augmentation transformation is naturally caused by
Proposed method is strongly inspired from the above self- multisensor differences and other factors including weather
supervised methods. Like deep clustering [31], the proposed condition. Considering that the prior probability of change
method uses the concept of simultaneous representation learn- is less [38], most of the times such pair of patches xb1 and
ing and cluster/label assignment. The bi-temporal images can z2b represent the same information, but from the eyes of two
be considered to be views of same scene, like SiamSiam [16]. different viewers (sensors).
Like the contrastive methods, the proposed method uses the
idea of bringing closer the representation of positive pairs and B. Siamese representation
spreading apart the negative pairs. Like [15], the proposed
method works on single scene (a pair of images capturing Since bi-temporal patches can be seen as multiple views
same location at two different times). of the same location, we argue that semantic information can
Multi-temporal satellite image processing researchers have be captured from them by using a Siamese-like architecture.
also proposed self-supervised representation learning methods, Similar to [16], both branches of the two-branch network
e.g., deep clustering for multi-temporal segmentation [36] and have projection modules fopt and fsar for the optical and
learning by rearranging randomly shuffled time-series images SAR branch, respectively. Additionally, both branches have
[37]. Proposed method is related to them, using the concept prediction modules hopt and hsar for the optical and SAR
of deep clustering as in [36]. branch, respectively. However, unlike [16], the projection
modules fopt and fsar do not share weight. This is because
SAR and optical images are significantly different processed
III. P ROPOSED METHOD by two different projection modules using different sets of
Let X1 , Z2 be two images of size R×C taken over the same weights. However, the prediction modules hopt and hsar share
geographical region at time t1 and t2 , respectively. Without weights and henceforth simply denoted as h.
loss of generality we assume that the pre-change image X1 The projection and the prediction networks consists of L1
is acquired by optical sensor (RGB) and the post-change and L2 (generally L2 = 1) convolutional layers, respectively,
image Z2 is acquired by the SAR sensor. Since SAR image where L = L1 + L2 . The two projections compute a projected
is grayscale, same channel is replicated thrice to make it 3- representation from the optical and SAR images and project
channel like the optical input. We aim to detect changes from them to a common domain. In the ideal scenario, where the
the images X1 , Z2 in an unsupervised manner, i.e., without us- projectors have perfectly learned to project optical and SAR
ing any training labels and any additional unlabeled data pool. images into a common domain and the bi-temporal images do
Our goal is to divide the set of all pixels Ω into two subsets not show any change, the output generated for an input pair is
Ωc and ωnc corresponding to changed and unchanged pixels, expected to be identical. However, practically even in absence
respectively. Like most existing unsupervised CD methods [1], of any change, there are differences caused by multi-sensor
we assume that prior probability of occurrence of change is acquisition and other factors that are not trivial for projection
less compared to no change [38]. modules to mitigate.
 4

Shuffle 0
{xb1 }B
b=1 {z2b }B
b=1 {z2b }B
b=1 X1 X2

Deep trainable model Deep model

0
{y1b }B
b=1 {y2b }B
b=1 {y2b }B
b=1 Bitemporal deep feature extraction
Iterative refinement

0 Estimate G
{cb1 }B
b=1 {cb2 }B
b=1 {cb2 }B
b=1

Threshold G
Compute self-supervised losses
L1 , L1,2 , and L01,2
ωnc , Ωc

Fig. 2. Proposed unsupervised multi-sensor (optical-SAR) CD framework. The left hand side denotes the self-supervised training process while the right hand
side shows the CD process using already trained model.

Optical input Conv1opt Conv2opt Conv3opt Conv4opt Optical output

Conv5common
SAR input Conv1sar Conv2sar Conv3sar Conv4sar SAR output

Fig. 3. The network simplified architecture with L1 = 4 and L2 = 1. Optical and SAR inputs are processed separately and subsequently fed to a common
prediction layer.

All but last convolution layers are followed by ReLU to 1-dimensional label cb1,n by argmax classification. This
b
activation function. They are further followed by batch normal- is achieved by selecting the kernel/feature in y1,n that has
b
ization layer. We do not use any pooling layer, hence the size maximum value. If the k-th feature of y1,n is represented by
b
of the input is preserved in the output. While filters of spatial y1,n (k), then label cb1,n is obtained as following:
size 3 × 3 are used for all convolution layers for projection,
the prediction module uses 1 × 1 filter. The kernel number cb1,n = arg max y1,n
b
(k) (3)
of the final layer is K and can be thought of as K different k∈K
clusters/classes. Each pixel can be assigned to one of these K
clusters (detailed in Section III-C). The network architecture The rationale behind finding the highest activation of an input
is shown in Figure 3. pixel is that the pixels that obtain the highest activation in
the same feature are likely to have similar semantics, thus
C. Deep clustering belonging to the same group. While there are several possible
Deep clustering process involves the joint learning of the ways to define the pseudo-label, our approach more closely
parameters of deep network and the cluster assignment of follows the ones based on argmax classification of the final
the resulting features [31]. Deep clustering helps the network layer [39], [40]. Once the pixels are assigned to the K clusters,
to learn discriminative features that can identify different parameters of the deep network can be updated by using a loss
b
classes/clusters in the images. Considering processing of the between the feature y1,n and the cluster cb1,n . We use cross-
two images as independent process, deep clustering can be entropy loss as:
performed for each of them. The output obtained by the
network for a paired input patches xb1 and z2b is: `b1,n = crossentropy(y1,n
b
, cb1,n ) (4)
y1b = h(fopt (xb1 )) (1)
In practice, the loss term L1 is computed by taking mean
y2b = h(fsar (z2b )) (2) of `b1,n over all pixels in xb1 and all patches in the batch
y1b has same spatial dimension R0 × C 0 as xb1 and has kernel (b = 1, ..., B). L1 is used to adjust the weights of h and fopt .
number (or, feature dimension) K. The deep clustering process Similarly, L2 is computed from z2b (b = 1, ..., B) and used to
is performed over the pixels, i.e., each pixel is assigned to a modulate the weights of h and fsar .
cluster. Without loss of generality, we henceforth explain the While the deep clustering helps to learn representation for
b
deep clustering process in reference to a generic pixel y1,n each sensor separately, they do not ensure that the indepen-
b b
from y1 . The dimension of y1,n is K that can be converted dently learned features are aligned with each other.
 5

Algorithm 1 Self-supervised training for multisensor CD unsupervised setting in which our work is based on, we simply
1: Initialize W1 , ..., WL shuffle the batch of patches Z to Z 0 . Recall that X and Z
2: for i ← 1 to I do have location-wise paired patches. This implies that X and Z 0
3: Sample B patches from X1 , denoted as X = have unpaired patches. Thus they should be more dissimilar
{x11 , ..., xB
1} in comparison to the paired patches in Section III-D. We
0
4: Obtain corresponding B patches from Z2 , denoted as encourage features computed for xb1 and z2b to be dissimilar.
Z = {z21 , ..., z2B } This is achieved by computing (negative) absolute error loss
b0
5: Obtain Z 0 as random shuffling of Z for each input pixel xb1,n and z2,n :
6: for j ← 1 to J do 0 0
`b12,n = −||y1,n
b b
− y2,n ||1 (6)
7: for b ∈ B do
0 0
8: y1b = h(fopt (xb1 )) `b12,n has negative value. Ideally `b12,n should be encouraged to
9: be more and more negative. However in practice we note that,
10: y2b = h(fsar (z2b )) simply shuffling Z to Z 0 does not always ensure that X and Z 0
11: have semantically different patches. Even after shuffling they
0 0
12: y2b = h(fsar (z2b )) may have the semantically paired patches, however penalized
13: in this step for producing similar features. Thus to control
0
14: end for its impact, we penalize the network with `b12,n only when it
15: Calculate deep clustering losses L1 , L2 b
approaches 0, i.e., y1,n b0
and y2,n become too similar. This
16: Calculate temporal consistency loss L1,2 is achieved by computing the loss term L01,2 as mean of
17: Calculate contrastive loss L01,2 0
exponentials of `b12,n over all considered pixels for all patches
18: if i ≤ I1 then in the batch.
19: Use loss (L1 + L2 )/2 to modulate W1 , ..., WL
20: else F. Overall loss and network refinement
21: For each 3 consecutive iterations j, use L1 ,
He initialization process [41] is used to initialize all the
L1,2 , and L01,2 , respectively, to modulate W1 , ..., WL
trainable weights of the network W1 , ..., WL , corresponding
22: end if
to L layers. For updating of weights, we exploit stochastic
23: end for
gradient descent (SGD) mechanism with momentum [42].
24: end for
The training process is executed in two different steps of I1
and I2 epochs (summing to I). For each batch of data, J
D. Temporal consistency iterations are performed. For the first I1 epochs only sum
of deep clustering losses (L1 + L2 ) is used to modulate the
Recalling from Section III-B, multi-sensor bi-temporal network weights. For subsequent I2 epochs, in one training
patches xb1 and z2b are multiple views of the same location iteration L1 is used as loss function, in the following iteration
in absence of any change. In other words, in co-registered bi- L1,2 is used and in the following iteration L01,2 is used.
temporal images, pixels in the same spatial location generally The combination of three loss functions yield a balanced
tend to belong to the same object as changes have a low training process taking into account coherent cluster formation,
prior probability than the unchanged class. Thus, the features temporal feature consistency, and feature dissimilarity for
computed for the bi-temporal paired patches xb1 and z2b should unpaired patches. Alternatively, sum of L1 , L1,2 , and L01,2 can
be similar in most cases. For each input pixel xb1,n and z2,n
b
, also be used as aggregated loss function. The self-supervised
we compute absolute error (AE) loss as: mechanism for network training is shown in Algorithm 1.
`b12,n = ||y1,n
b b
− y2,n ||1 (5)
G. Change detection
A loss term L1,2 is computed by taking mean of `b12,n
over all Once the network is trained, it can be used to detect change
considered pixels for all patches in the batch. The proposed between X1 and Z2 . Since the network is fully convolutional,
temporal consistency only ensures that the pixels at same it enables us to obtain pixelwise feature vector of dimension
location however at two different time tend to have same label. K from X1 and Z2 . Similar to [1], the pixelwise change
This may lead a to a degenerate solution where all pixels information is captured by taking the magnitude (`2 norm) of
simply have same prediction for both times. Moreover, some difference of the feature vectors computed from pre-change
bi-temporal pairs xb1 and z2b may be indeed changed, however, and post-change pixels. Changed pixels (Ωc ) generate higher
penalized for producing dissimilar output in this step. difference magnitude in comparison to the unchanged ones
ωnc and they can be distinguished by using any suitable
E. Contrastive learning threshold determination scheme [43].
While Section III-D encourages the features computed for
IV. E XPERIMENTAL VALIDATION
paired patch xb1 and z2b to be similar, in this Section we
encourage the network to produce dissimilar feature for dif- A. Datasets
ferent input by employing concepts inspired from contrastive We use four paired optical (pre-change) - SAR (post-
learning. While we do not have negative samples under the change) images to validate the proposed method. Optical
 6

(a) (b) (c) (a) (b) (c)

(d) (e) (f)

(d) (e) (f) Fig. 6. CD results for the Abu Dhabi. CD maps: (a) Reference, (b) Proposed,
(c) FCC between reference and proposed (the correctly detected region are in
black, false alarms are in green, missed alarms are in pink), (d) PCVA, (e)
DCVA, (f) SCCN.

(g) (h) (i)

Fig. 4. CD results for Las Vegas. CD maps: (a) Reference, (b) Proposed,
(c) FCC between reference and proposed (the correctly detected region are
in black, false alarms are in green, missed alarms are in pink), (d) CVA, (e) (a) (b) (c)
RCVA, (f) PCVA, (g) DCVA, (h) Encoder-decoder, and (i) SCCN.

(d) (e) (f)

Fig. 7. Qualitative CD results for the Montpellier. Change detection maps:
(a) Reference, (b) Proposed, (c) FCC between reference and proposed (the
correctly detected region are in black, false alarms are in green, missed alarms
are in pink), (d) PCVA, (e) DCVA, (f) SCCN.
(a) (b) (c)

States (824 × 716 pixels) (Figure 4), Chongqing in China

(730 × 544 pixels) (Figure 5), Abu Dhabi (799 × 785 pixels)
(Figure 6), and Montpellier in France (426 × 451 pixels)
(Figure 7). Thus this provides us an opportunity to validate
the proposed method on geographically distributed complex
(d) (e) (f) urban scenes with large variation.
Fig. 5. CD results for the Chongqing. CD maps: (a) Reference, (b) Proposed,
(c) FCC between reference and proposed (the correctly detected region are in
black, false alarms are in green, missed alarms are in pink), (d) PCVA, (e) B. Compared methods
DCVA, (f) SCCN.
To verify the effectiveness of the proposed method, we
compare it to related unsupervised change detection methods:
images are acquired by the Sentinel-2 sensor and are taken 1) Change vector analyis (CVA) [17], [45], a classical
from the Onera Satellite Change Detection (OSCD) dataset difference-based unsupervised model for change detec-
[44]. They show 10 m/pixel spatial resolution. OSCD dataset tion.
is originally a single-sensor dataset consisting only Sentinel- 2) Robust change vector analysis (RCVA) [19] that mod-
2 images. Recalling the importance of multisensor CD (see ifies CVA by taking into account pixel neighborhood
Section I), we extend this dataset by collecting the post-change effects.
SAR Sentinel-1 images for the nearest available date as the 3) Parcel change vector analysis (PCVA) [18] that incopo-
post-change image in original OSCD dataset. Both Sentinel- rates notion of object (superpixels) in CVA.
2 and Sentinel-1 sensors are part of the European Space 4) Deep change vector analysis (DCVA) [1] that detects
Agency’s Copernicus program. change by comparing bi-temporal deep features ex-
The four scenes are collected over Las Vegas in United tracted using a pre-trained network. We used second
 7

convolution layer of pre-trained VGGNet [46] for feature a false color composition between the reference map and
extraction. the obtained result is shown in Figure 4(c). The proposed
5) Image-to-image transfer model based on an encoder- method can detect most of the changed objects with fewer
decoder network architecture that projects pre-change false alarms in comparison to the compared methods. In many
optical images into post-change SAR image [47]. CD cases, proposed method partly detects the changed object, thus
map can be obtained by difference of the simulated missing some objects only partially (shown in pink in Figure
pre-change SAR image (obtained as projection of pre- 4(c)). CVA (Figure 4(d)) performs poorly and incorrectly
change optical image) and the original post-change SAR detects most urban area as changed. Result obtained by RCVA
image. (Figure 4(e)) is similar to CVA. While PCVA (Figure 4(f)),
6) Denoising autoencoder (DAE) based joint feature extrac- DCVA (Figure 4(g)), encoder-decoder (Figure 4(h)), DAE, and
tion [27]. SCCN (Figure 4(i)) improve the result over CVA, proposed
7) SCCN [26] that first identifies some unchanged pixels method still outperforms them by large margin. Quantitative
and uses them to learn a coupled network. evaluation (Table II) clearly shows the superiority of the
While the methods 1-3 are not deep learning based, the proposed method over state-of-the-art unsupervised methods.
following ones are deep learning based. The methods 1-4 do This can be attributed to superior capability of the proposed
not have any explicit adaptation for multi-sensor input, while method to ingest multi-sensor multi-temporal images.
the method 5, 6, and 7 have. Further studies are conducted by varying different parame-
ters on Las Vegas image pair.
C. Experimental settings Training epochs I is varied with different values as tab-
ulated in Table III while setting K = 4. We observe clear
The proposed method and compared methods are fed with
improvement in performance from I = 1 to 2. Recalling from
pre-processed images and post-processed similarly. For the
Section III-F that for first I1 = 1 iterations only deep cluster-
proposed method, we use I = 5 (I1 = 1, I2 = 4), J = 50
ing loss is used, this shows that bi-temporal deep clustering
K = 4, L1 = 4, and L2 = 1. We show the architecture
itself is not sufficient to learn the correspondence between two
of the network in Table I. A relatively simple architecture is
images and the other losses (L1,2 and L01,2 ) are required. I = 2
used considering that the number of patches available to us is
onwards, we observe an increment in performance initially
very few compared to the images in typical computer vision
followed by performance getting saturated/dropping. Despite
datasets. Moreover, our target image has coarse resolution
variation in performance, proposed method outperforms all
(10 m/pixel) compared to natural images in computer vision.
compared methods for I = 3, 5, 10.
Spatial complexity in such coarse images can be handled by
Kernel number of last layer (K) is varied from 2 to 16 in
simpler architecture compared to those in computer vision.
multiplicative steps of 2 while fixing the I = 5. The variation
64 × 64 patches are used to train the model and patches
in performance is shown in Table IV. While performance
are extracted from the bitemporal scene with a stride of 32.
improves from K = 2 to K = 4, a gradual fall in performance
The actual number of training patches for a scene depends on
is observed henceforth. Increasing value of K is equivalent
the size of the particular scene. E.g., for the Las Vegas scene
to allowing the scene to be partitioned more classes. Since
(824 × 716 pixels), number of patches extracted is 504. For
the spatial area of the scene is fixed and not too large (only
optimization, Stochastic Gradient Descent method is used with
few hundred pixels by few hundred pixels), large number of
learning rate set to 0.001.
classes potentially leads the model to learn irrelevant classes,
impacting change detection performance.
TABLE I
S TRUCTURE OF THE NETWORK FOR PROCESSING ONE OF THE TWO INPUTS Thresholding is done using Otsu’s method [43], as it is
popular in unsupervised CD methods [19], [48]. However any
Layer Kernel number Kernel size Stride other suitable method can be used, as shown in Table V.
convolution 64 (3,3) 1
convolution 64 (3,3) 1 Results obtained by ISODATA method [49], [50] and adaptive
convolution 64 (3,3) 1 method [1] are similar to Otsu’s method [43].
convolution 64 (3,3) 1 Loss plot visualization in Figure 8 shows the interplay be-
convolution K (1,1) 1 tween different components of loss. L1 consistently decreases
(Figure 8(a)) except it rises for a while after epoch 1 when
We show result in terms of sensitivity (accuracy in percent-
L1,2 and L01,2 are introduced to the training process. L1,2 and
age computed over reference changed pixels) and specificity
L01,2 balances each other as shown in Figure 8(b).
(computed over reference unchanged pixels). In more details,
Projection layers fopt and fsar need to modeled indepen-
given true positive (TP), true negative (TN), false positive
dently by not sharing weights between them to capture the
(FP) and false negative (FN), sensitivity is TP/(TP+FN) and
different semantic properties of optical and SAR patches, as
specificity is TN/(TN+FP).
hypothesized in Section III-B. Here we test this hypothesis by
instead sharing the weights between fopt and fsar . For I = 5
D. Results and K = 4 the proposed method fails to detect most of the
Las Vegas: The reference CD map (ground truth) for Las changes. This shows that it is crucial to model the optical and
Vegas is shown in Figure 4(a). Figure 4(b) shows the result SAR patches differently.
obtained by the proposed method. For better visualization, Computation time requirement is not high. We tested our
 8

TABLE II
0.7
2.0
L1,2 C OMPARISON OF DIFFERENT METHODS ON L AS V EGAS
0.6 L01,2

0.5 1.5
Method Sensitivity Specificity
Loss

Loss
0.4

0.3
1.0 Proposed 50.28 88.06
0.2
0.5 CVA 9.64 77.13
0.1

0.0
RCVA 8.65 78.97
0 5000 10000 15000
Iteration
20000 25000 30000 0 5000 10000 15000
Iteration
20000 25000 30000
PCVA 23.60 67.92
DCVA 20.67 75.40
(a) (b) Encoder-decoder 46.30 68.58
DAE 39.07 83.72
Fig. 8. Evolution of the loss over training iterations for Las Vegas: (a) deep
SCCN 24.58 75.30
clustering loss L1 , (b) temporal consistency loss L1,2 and contrastive loss
L01,2

TABLE III
VARIATION OF RESULT FOR L AS V EGAS AS I IS VARIED
code on a machine equipped with a Quadro T2000 GPU,
which is a low end GPU. For processing Las Vegas dataset I Sensitivity Specificity
(training process over 5 epochs), it takes approx. 460 seconds. 0 (Just initialized) 32.21 85.30
1 61.70 53.40
Las Vegas scene is 824 × 716 pixels with 10 m/pixel reso- 2 52.20 79.82
lution and thus processing it is equivalent to processing an 3 46.30 86.52
approximate area of 8*7 = 56 sq. km in terms of geography. 5 50.28 88.06
10 41.79 88.80
Same sensor bi-temporal input can be ingested by the pro-
posed method, though designed for multi-sensor CD. For Las
Vegas prechange optical - postchange optical input, proposed TABLE IV
method can obtain a sensitivity of 64.74% and specificity of VARIATION OF RESULT FOR L AS V EGAS AS K IS VARIED
97.89%. However, we note that some characteristics of the
K Sensitivity Specificity
proposed method (e.g., temporal consistency loss) are designed 2 41.47 87.48
to reduce the representation gap of multi-sensor input, which 4 50.28 88.06
is less relevant in single-sensor input. So the proposed method 8 42.75 84.70
16 60.01 77.86
may not be the most suitable choice for single-sensor scenarios
as there are numerous existing CD techniques particularly
designed for same-sensor scenario [1]. TABLE V
Chongqing and Abu Dhabi: Reference CD map (ground VARIATION OF RESULT FOR L AS V EGAS AS THRESHOLD DETERMINATION
SCHEME IS VARIED
truth) for Chongqing is shown in Figure 5(a). Figure 5(b) and
5(c) show the result obtained by the proposed method and Thresholding Sensitivity Specificity
false color composition between the reference map and the Otsu 50.28 88.06
obtained result, respectively. Proposed method outperforms all ISODATA 50.48 87.95
Adaptive 50.19 83.65
compared methods, as can be observed in quantitative result
in Table VI). Similar result is obtained for Abu Dhabi (Figure
6 and Table VII). TABLE VI
Montpellier: Reference CD map (ground truth) for Mont- C OMPARISON OF DIFFERENT METHODS ON C HONGQING
pellier is shown in Figure 7(a). The proposed method (Figure Method Sensitivity Specificity
7(b)) outperforms most of the state-of-the-art methods includ- Proposed 36.17 83.26
ing PCVA (Figure 7(d)) and DCVA (Figure 7(e)), as shown CVA 40.70 48.82
in Table VIII. However, SCCN (Figure 7(f)) outperforms the RCVA 41.96 44.28
PCVA 35.15 56.76
proposed method. The performance of the proposed method DCVA 32.67 79.18
is relatively poor for Montpellier, which can be possibly Encoder-decoder 19.59 82.23
explained by: 1) smaller size of Montpellier scene, which DAE 17.58 82.60
implies less data to learn proposed self-supervised network SCCN 30.67 85.73
and 2) uniform (showing mostly urban areas) geospatial char-
acteristics of Montpellier scene in comparison to Las Vegas TABLE VII
and Chongqing that show complex distribution consisting of C OMPARISON OF DIFFERENT METHODS ON A BU D HABI
both urban and non-urban areas.
Method Sensitivity Specificity
Proposed 48.92 84.38
V. C ONCLUSIONS CVA 4.18 73.09
RCVA 5.75 73.47
This paper proposed a self-supervised learning based PCVA 13.30 65.24
method for CD in multi-sensor bi-temporal images where one DCVA 20.35 76.74
of the image is acquired by optical sensor and the other one Encoder-decoder 36.52 74.22
DAE 46.29 72.29
is captured by SAR sensor. The proposed method effectively
SCCN 9.79 83.88
utilizes several concepts from the self-supervised learning,
 9

TABLE VIII [6] S. T. Seydi and M. Hasanlou, “A new land-cover match-based change
C OMPARISON OF DIFFERENT METHODS ON M ONTPELLIER detection for hyperspectral imagery,” European Journal of Remote
Sensing, vol. 50, no. 1, pp. 517–533, 2017.
Method Sensitivity Specificity [7] M. Puhm, J. Deutscher, M. Hirschmugl, A. Wimmer, U. Schmitt, and
Proposed 43.05 81.88 M. Schardt, “A near real-time method for forest change detection based
CVA 9.09 74.60 on a structural time series model and the kalman filter,” Remote Sensing,
RCVA 8.76 72.31 vol. 12, no. 19, p. 3135, 2020.
PCVA 25.74 58.59 [8] J. E. Ball, D. T. Anderson, and C. S. Chan, “Comprehensive survey
DCVA 32.37 76.49 of deep learning in remote sensing: theories, tools, and challenges for
Encoder-decoder 46.10 75.52 the community,” Journal of Applied Remote Sensing, vol. 11, no. 4, p.
DAE 24.46 73.97 042609, 2017.
SCCN 50.57 97.18 [9] M. Hirschmugl, J. Deutscher, C. Sobe, A. Bouvet, S. Mermoz, and
M. Schardt, “Use of sar and optical time series for tropical forest
disturbance mapping,” Remote Sensing, vol. 12, no. 4, p. 727, 2020.
[10] N. Zhou, X. Li, Z. Shen, T. Wu, and J. Luo, “Geo-parcel-based change
e.g., deep clustering, Siamese network, multiple view, and detection using optical and sar images in cloudy and rainy areas,” IEEE
contrastive learning and operates under severe constraints, Journal of Selected Topics in Applied Earth Observations and Remote
i.e., nothing except the target scene is used and no labeled Sensing, 2020.
data or additional unlabeled image is used. Despite strong [11] U. I. Ahmed, B. Rabus, and M. F. Beg, “Sar and optical image fusion
for urban infrastructure detection and monitoring,” in Remote Sensing
difference in the input modalities and operating under stringent Technologies and Applications in Urban Environments V, vol. 11535.
constraints, it can identify a large fraction of the changed International Society for Optics and Photonics, 2020, p. 115350M.
pixels. Comparisons with the existing methods working under [12] M. Schmitt and X. X. Zhu, “Data fusion and remote sensing: An ever-
growing relationship,” IEEE Geoscience and Remote Sensing Magazine,
unsupervised scenario show that the proposed method brings vol. 4, no. 4, pp. 6–23, 2016.
significant improvement, especially when the target scene is [13] X. X. Zhu, D. Tuia, L. Mou, G.-S. Xia, L. Zhang, F. Xu, and
large. Potential improvement of the proposed method may F. Fraundorfer, “Deep learning in remote sensing: A comprehensive
review and list of resources,” IEEE Geoscience and Remote Sensing
be achieved by prior learning of clusters on the unrelated Magazine, vol. 5, no. 4, pp. 8–36, 2017.
domains/sensors and transferring them to target sensors on the [14] G. Camps-Valls, D. Tuia, X. X. Zhu, and M. Reichstein, Deep learning
fly [51]. Additionally, our future work will focus on extending for the Earth Sciences: A comprehensive approach to remote sensing,
the method to other application domains, e.g., comparison of climate science and geosciences. John Wiley & Sons, 2021.
[15] Y. M. Asano, C. Rupprecht, and A. Vedaldi, “A critical analysis of self-
biomedical images. supervision, or what we can learn from a single image,” arXiv preprint
arXiv:1904.13132, 2019.
ACKNOWLEDGEMENT [16] X. Chen and K. He, “Exploring simple siamese representation learning,”
arXiv preprint arXiv:2011.10566, 2020.
The work is jointly supported by the European Research [17] W. A. Malila, “Change vector analysis: an approach for detecting forest
Council (ERC) under the European Union’s Horizon 2020 re- changes with Landsat,” in LARS symposia, 1980, p. 385.
search and innovation programme (grant agreement No. [ERC- [18] F. Bovolo, “A multilevel parcel-based approach to change detection in
very high resolution multitemporal images,” IEEE GRSL, vol. 6, no. 1,
2016-StG-714087], Acronym: So2Sat), by the Helmholtz As- pp. 33–37, 2009.
sociation through the Framework of Helmholtz AI [grant [19] F. Thonfeld, H. Feilhauer, M. Braun, and G. Menz, “Robust change
number: ZT-I-PF-5-01] - Local Unit “Munich Unit @Aeronau- vector analysis (RCVA) for multi-sensor very high resolution optical
satellite data,” International Journal of Applied Earth Observation and
tics, Space and Transport (MASTr)” and Helmholtz Excellent Geoinformation, vol. 50, pp. 131–140, 2016.
Professorship “Data Science in Earth Observation - Big Data [20] S. Saha, F. Bovolo, and L. Bruzzone, “Building change detection in vhr
Fusion for Urban Research”(W2-W3-100) and by the German sar images via unsupervised deep transcoding,” IEEE Transactions on
Geoscience and Remote Sensing, 2020.
Federal Ministry of Education and Research (BMBF) in the
[21] J. Deng, W. Dong, R. Socher, L.-J. Li, K. Li, and L. Fei-Fei, “Imagenet:
framework of the international future AI lab ”AI4EO – Artifi- A large-scale hierarchical image database,” in 2009 IEEE conference on
cial Intelligence for Earth Observation: Reasoning, Uncertain- computer vision and pattern recognition. Ieee, 2009, pp. 248–255.
ties, Ethics and Beyond” (Grant number: 01DD20001). [22] H. M. Keshk and X.-C. Yin, “Change detection in sar images based
on deep learning,” International Journal of Aeronautical and Space
Sciences, pp. 1–11, 2019.
R EFERENCES [23] S. Saha, F. Bovolo, and L. Bruzzone, “Unsupervised multiple-change
[1] S. Saha, F. Bovolo, and L. Bruzzone, “Unsupervised deep change detection in vhr multisensor images via deep-learning based adaptation,”
vector analysis for multiple-change detection in vhr images,” IEEE in IGARSS 2019-2019 IEEE International Geoscience and Remote
Transactions on Geoscience and Remote Sensing, vol. 57, no. 6, pp. Sensing Symposium. IEEE, 2019, pp. 5033–5036.
3677–3693, 2019. [24] P. Zhang, M. Gong, L. Su, J. Liu, and Z. Li, “Change detection based
[2] A. Appice, N. Di Mauro, F. Lomuscio, and D. Malerba, “Empowering on deep feature representation and mapping transformation for multi-
change vector analysis with autoencoding in bi-temporal hyperspectral spatial-resolution remote sensing images,” ISPRS Journal of Photogram-
images,” in CEUR Workshop Proceedings, vol. 2466, 2019. metry and Remote Sensing, vol. 116, pp. 24–41, 2016.
[3] H. Chen, C. Wu, B. Du, L. Zhang, and L. Wang, “Change detection [25] M. Volpi, G. Camps-Valls, and D. Tuia, “Spectral alignment of multi-
in multisource vhr images via deep siamese convolutional multiple- temporal cross-sensor images with automated kernel canonical correla-
layers recurrent neural network,” IEEE Transactions on Geoscience and tion analysis,” ISPRS J Photogramm Remote Sens., vol. 107, pp. 50–63,
Remote Sensing, vol. 58, no. 4, pp. 2848–2864, 2019. 2015.
[4] F. Rahman, B. Vasu, J. Van Cor, J. Kerekes, and A. Savakis, “Siamese [26] J. Liu, M. Gong, K. Qin, and P. Zhang, “A deep convolutional coupling
network with multi-level features for patch-based change detection in network for change detection based on heterogeneous optical and radar
satellite imagery,” in 2018 IEEE Global Conference on Signal and images,” IEEE transactions on neural networks and learning systems,
Information Processing (GlobalSIP). IEEE, 2018, pp. 958–962. vol. 29, no. 3, pp. 545–559, 2016.
[5] A. Pomente, M. Picchiani, and F. Del Frate, “Sentinel-2 change detection [27] T. Zhan, M. Gong, X. Jiang, and S. Li, “Log-based transformation
based on deep features,” in IGARSS 2018-2018 IEEE International feature learning for change detection in heterogeneous images,” IEEE
Geoscience and Remote Sensing Symposium. IEEE, 2018, pp. 6859– Geoscience and Remote Sensing Letters, vol. 15, no. 9, pp. 1352–1356,
6862. 2018.
 10

[28] S. Gidaris, P. Singh, and N. Komodakis, “Unsupervised repre-

sentation learning by predicting image rotations,” arXiv preprint
arXiv:1803.07728, 2018.
[29] C. Doersch, A. Gupta, and A. A. Efros, “Unsupervised visual repre-
sentation learning by context prediction,” in Proceedings of the IEEE
international conference on computer vision, 2015, pp. 1422–1430.
[30] D. Pathak, P. Krahenbuhl, J. Donahue, T. Darrell, and A. A. Efros,
“Context encoders: Feature learning by inpainting,” in Proceedings of
the IEEE conference on computer vision and pattern recognition, 2016,
pp. 2536–2544.
[31] M. Caron, P. Bojanowski, A. Joulin, and M. Douze, “Deep clustering
for unsupervised learning of visual features,” in Proceedings of the
European Conference on Computer Vision (ECCV), 2018, pp. 132–149.
[32] T. Chen, S. Kornblith, M. Norouzi, and G. Hinton, “A simple frame-
work for contrastive learning of visual representations,” arXiv preprint
arXiv:2002.05709, 2020.
[33] X. Chen, H. Fan, R. Girshick, and K. He, “Improved baselines with mo-
mentum contrastive learning,” arXiv preprint arXiv:2003.04297, 2020.
[34] Y. Tian, C. Sun, B. Poole, D. Krishnan, C. Schmid, and P. Isola,
“What makes for good views for contrastive learning,” arXiv preprint
arXiv:2005.10243, 2020.
[35] J.-B. Grill, F. Strub, F. Altché, C. Tallec, P. H. Richemond,
E. Buchatskaya, C. Doersch, B. A. Pires, Z. D. Guo, M. G. Azar et al.,
“Bootstrap your own latent: A new approach to self-supervised learning,”
arXiv preprint arXiv:2006.07733, 2020.
[36] S. Saha, L. Mou, C. Qiu, X. X. Zhu, F. Bovolo, and L. Bruzzone,
“Unsupervised deep joint segmentation of multitemporal high-resolution
images,” IEEE Transactions on Geoscience and Remote Sensing, 2020.
[37] S. Saha, F. Bovolo, and L. Bruzzone, “Change detection in image time-
series using unsupervised lstm,” IEEE Geoscience and Remote Sensing
Letters, 2020.
[38] F. Bovolo and L. Bruzzone, “The time variable in data fusion: a change
detection perspective,” IEEE Geoscience and Remote Sensing Magazine,
vol. 3, no. 3, pp. 8–26, 2015.
[39] A. Kanezaki, “Unsupervised image segmentation by backpropagation,”
in 2018 IEEE international conference on acoustics, speech and signal
processing (ICASSP). IEEE, 2018, pp. 1543–1547.
[40] S. Saha, S. Sudhakaran, B. Banerjee, and S. Pendurkar, “Semantic
guided deep unsupervised image segmentation,” in International Confer-
ence on Image Analysis and Processing. Springer, 2019, pp. 499–510.
[41] K. He, X. Zhang, S. Ren, and J. Sun, “Delving deep into rectifiers:
Surpassing human-level performance on imagenet classification,” in
Proceedings of the IEEE international conference on computer vision,
2015, pp. 1026–1034.
[42] S. Ruder, “An overview of gradient descent optimization algorithms,”
arXiv preprint arXiv:1609.04747, 2016.
[43] N. Otsu, “A threshold selection method from gray-level histograms,”
IEEE transactions on systems, man, and cybernetics, vol. 9, no. 1, pp.
62–66, 1979.
[44] R. C. Daudt, B. Le Saux, A. Boulch, and Y. Gousseau, “Urban change
detection for multispectral earth observation using convolutional neural
networks,” in IGARSS 2018-2018 IEEE International Geoscience and
Remote Sensing Symposium. IEEE, 2018, pp. 2115–2118.
[45] L. Bruzzone and D. F. Prieto, “Automatic analysis of the difference
image for unsupervised change detection,” IEEE Transactions on Geo-
science and Remote sensing, vol. 38, no. 3, pp. 1171–1182, 2000.
[46] K. Simonyan and A. Zisserman, “Very deep convolutional networks for
large-scale image recognition,” arXiv preprint arXiv:1409.1556, 2014.
[47] Y. Xu, S. Xiang, C. Huo, and C. Pan, “Change detection based on auto-
encoder model for VHR images,” in MIPPR 2013: Pattern Recognition
and Computer Vision, vol. 8919. International Society for Optics and
Photonics, 2013, p. 891902.
[48] S. Saha, Y. T. Solano-Correa, F. Bovolo, and L. Bruzzone, “Unsuper-
vised deep transfer learning-based change detection for hr multispectral
images,” IEEE Geoscience and Remote Sensing Letters, 2020.
[49] T. Ridler, S. Calvard et al., “Picture thresholding using an iterative
selection method,” IEEE trans syst Man Cybern, vol. 8, no. 8, pp. 630–
632, 1978.
[50] M. Sezgin and B. Sankur, “Survey over image thresholding techniques
and quantitative performance evaluation,” Journal of Electronic imaging,
vol. 13, no. 1, pp. 146–165, 2004.
[51] W. Menapace, S. Lathuilière, and E. Ricci, “Learning to cluster under
domain shift,” arXiv preprint arXiv:2008.04646, 2020.