Change Detection Techniques with Synthetic Aperture Radar Images: Experiments with Random Forests and Sentinel-1 Observations

: This work aims to clarify the potential of incoherent and coherent change detection (CD) approaches for detecting and monitoring ground surface changes using sequences of synthetic aperture radar (SAR) images. Nowadays, the growing availability of remotely sensed data collected by the twin Sentinel-1A/B sensors of the European (EU) Copernicus constellation allows fast mapping of damage after a disastrous event using radar data. In this research, we address the role of SAR (amplitude) backscattered signal variations for CD analyses when a natural (e.g., a ﬁre, a ﬂash ﬂood, etc.) or a human-induced (disastrous) event occurs. Then, we consider the additional pieces of information that can be recovered by comparing interferometric coherence maps related to couples of SAR images collected between a principal disastrous event date. This work is mainly concerned with investigating the capability of different coherent/incoherent change detection indices (CDIs) and their mutual interactions for the rapid mapping of “changed” areas. In this context, artiﬁcial intelligence (AI) algorithms have been demonstrated to be beneﬁcial for handling the different information coming from coherent/incoherent CDIs in a unique corpus. Speciﬁcally, we used CDIs that synthetically describe ground surface changes associated with a disaster event (i.e., the pre-, cross-, and post-disaster phases), based on the generation of sigma nought and InSAR coherence maps. Then, we trained a random forest (RF) to produce CD maps and study the impact on the ﬁnal binary decision (changed/unchanged) of the different layers representing the available synthetic CDIs. The proposed strategy was effective for quickly assessing damage using SAR data and can be applied in several contexts. Experiments were conducted to monitor wildﬁre’s effects in the 2021 summer season in Italy, considering two case studies in Sardinia and Sicily. Another experiment was also carried out on the coastal city of Houston, Texas, the US, which was affected by a large ﬂood in 2017; thus, demonstrating the validity of the proposed integrated method for fast mapping of ﬂooded zones using SAR data.

The rest of the paper is organized as follows. Section 2 summarizes the fundamentals of coherent and incoherent CD methods. Section 3 introduces some synthetic change detection indices (CDIs) used to feed an RF. The overall CD strategy is then detailed. Available SAR images and ancillary data are shown in Section 4. Experimental results are presented in Section 5. Discussion, conclusions, and future perspectives are finally presented in Sections 6 and 7.

Fundamentals of Coherent and Incoherent Change Detection Methods
This section introduces the fundamentals of SAR change detection, distinguishing between incoherent and coherent approaches.
Let us consider two co-registered complex-valued SAR images X 1 and X 2 acquired at times t 1 and t 2 , respectively, and let X k = [X 1 (P), X 2 (P)] T be the pixel pair data vector extracted from SAR data values corresponding with a generic k-th pixel P of the imaged scene. The CD problem may generally be modeled as a binary classification [1], where each pixel is mapped into the set of possible labels Ω = {ω u , ω c } related to the unchanged and changed classes. Such a problem can effectively be formulated in a hypothesis-testing framework. In this context, the goal is to determine whether the pixel pair is a realization of a null hypothesis H 0 (unchanged scene) or an alternative hypothesis H 1 (changed scene). Under these two hypotheses, we can write that [83]: (1) and: H 1 X 1 (P) = X 1k = C 0k + n 1k X 2 (P) = X 2k = C 1k + n 2k (2) where the term C 0k is the complex reflectivity at the given pixel P that is common to both SAR images under the H 0 hypothesis, and C 1k is the complex reflectivity of the second image under the assumption that a change happened. Note that j = √ −1 and the phase term ϕ k accounts for a possible displacement of the second flight track for the first image related to a height topography change. Finally, the terms n 1k and n 2k account for noise decorrelation sources affecting the two available SAR images. For homogeneous natural distributed target environments, the observed SAR scenes can be modeled as consisting of many discrete scatterers with a uniform random phase and identically-distributed random amplitudes [57]. The observed complex radar backscattered signal can then be modeled as a zero-mean complex circular Gaussian process with the following probability density function (pdf) [11]: where Q is the covariance matrix of X k that, under the two hypotheses, is given by: where σ 2 c = E |C ok | 2 = E |C 1k | 2 , σ 2 n1 = E |n 1k | 2 and σ 2 n2 = E |n 2k | 2 . Note that E[·] stands for the statistical expectation extraction operation that, under the hypothesis that the signals are ergodic, can effectively be implemented with the spatial averaging procedure Remote Sens. 2022, 14, 3323 5 of 32 · N performed over a set of neighbor SAR pixels (i.e., with a multi-looking process [84]). The term σ c can also be expressed as a function of the complex correlation coefficient γ between the two interfering SAR images [85], as follows: The amplitude of the complex correlation coefficient is commonly referred to as coherence, and it assumes values between 0 and 1. Coherence is usually factorized into the product of several different contributions [84,[86][87][88]: γ =γ geo ·γ Dopp ·γ mis−reg ·γ vol ·γ temp ·γ thermal (7) Interested readers can refer to [86,88] to have a more detailed description of the different decorrelation sources.
Assuming perfect knowledge of the unknown parameters that are present in the covariance matrices of Equations (4) and (5), a formal way to solve a CD problem is to compute the log-likelihood ratio: considering a population made by N SAR pixels located neighboring a selected SAR pixel. After some mathematical manipulations and considering that the complex radar backscattered signal has the pdf of Equation (3), the decision variable Z can then be rewritten as: which is referred to in the literature as the clairvoyant detector. Note that the operator Tr(·) is the trace of a matrix. This detector gives the optimal detection performance given the perfect knowledge of the unknown parameters in the covariance matrices Q 0 and Q 1 . The expressions for the probability of detection and false alarm for the log-likelihood clairvoyant detector have also been derived in [85] as a function of a decision threshold. Nonetheless, the values of the covariance matrices are generally unknown, and thus, in practice, sub-optimal tests are implemented where the unknown parameters are measured directly through the data with their maximum likelihood estimates (MLEs) to form a generalized likelihood test (e.g., see [71]). Due to the complexity of these optimal tests and the need to estimate the mentioned model parameters, several approximate unsupervised/supervised CD approaches have been proposed in the literature (among others, see [14,45,[89][90][91][92][93][94]). Overall, they can be distinguished in the two general categories of the methods based exclusively on the analysis of the amplitude of SAR backscattered returns (incoherent approaches) and those that also consider the phase of the radar backscattered signals (coherent strategies).
The following sub-sections provide the readers with additional details on these two categories of CD methods.

Incoherent CD Approaches
Despite their intrinsic differences, the incoherent CD methods can be seen as made by three main steps: (1) pre-processing of the SAR images; (2) comparison of the pair of SAR images; and (3) automatic extraction of the changed areas by proper thresholding algorithms.
The first step consists of the radiometric calibration of the available SAR images [57,95] and the subsequent application of adaptive filters that reduce the impact of speckle noise in the SAR images. Among others, the most used de-noised filters are Lee [32], Kuan [48], Remote Sens. 2022, 14, 3323 6 of 32 and Frost [46]. Jointly, a spatial multi-look operation [96] is typically applied to enhance the SNR and reduce the images' spatial resolution as a side effect. The consequences of local topography can also be corrected [97], and the SAR images could then potentially be georeferenced for subsequent analyses.
The second step consists of comparing the two pre-processed SAR images. To this aim, we observe that signature variations in SAR images are due to texture and speckle, where texture is a measure of the intrinsic spatial variability of the backscattered coefficient (in the absence of noise), and the speckle is a multiplicative term. Mathematically, the received power can be written as [44,57]: where I is the backscattered power, and T and s are the texture and speckle noise random terms, respectively. When dealing with SAR images, the ratio operator r = X 2 N X 1 N is preferred to the image difference. Indeed, the ratio image pdf is given by [36]: where r is the true (unknown) value of the images' ratio. Accordingly, the first-order statistics depend only on the relative change and not on the single average values, separately. Moreover, the mean amplitude values ratio is more robust to calibration errors, typically represented by multiplicative terms that moderately vary from one image to another and then they are almost cancelled out in the image ratio. Furthermore, the ratio is reported on a log scale to represent significantly high and shallow intensity values. Accordingly, the detection image is obtained as: The third step consists of analyzing the log-ratio image and identifying a proper threshold that allows discriminating between the two changed and unchanged classes. This problem can be solved in the context of Bayesian decisions. A general algorithm for the threshold decision was proposed by Kittler and Illingworth (KI), which derives the best threshold from being used through a parametric estimation based on the assumption that the two classes are Gaussian distributed [65]. Some extensions of the KI method that assume a more reliable distribution of the two classes (e.g., a generalized Gaussian distribution) have also been proposed [35].

Coherent CD Approaches
Coherent change detection methods rely on using the interferometric cross-correlation factor (see Equation (7)). The statistics of the sample coherence have extensively been addressed in the literature [88,98]. In particular, the pdf of the sample coherenceγ, given the true coherence valueγ , is given by [84,98]: (13) where N is the effective number of looks (ENL) and 2 F 1 (·) is the Gauss hypergeometric function. It is worth mentioning that the coherence is a biased estimator: this means that the measured sample coherenceγ obtained by computing the expectation values using spatial averaging operations has an average value different from the actual coherence valuê γ . In particular, the first moment of the spatial coherence is given by [98]: Remote Sens. 2022, 14, 3323 7 of 32 and, consequently, the coherence bias is E[γ] −γ > 0. An approximate expression for the coherence bias, which could be considered to "calibrate" the coherence estimates for CD purposes, can straightforwardly be obtained in the complex domain by observing that the interferometric phase can be encoded in this domain as a point in a unitary circle. The mean resultant lengthχ of the circular data representing the interferometric phases into the N-points averaging multi-look window is obtained as: As demonstrated in ref. [99,100], the computed mean resultant lengthχ is an alternative measure of the coherence withχ ≤γ, and the two coherence factors, both biased, are related to one another as:χ The mean value of the mean resultant length of the circular data vector representing the phase of the scatterers within the multi-looking spatial window, under the suitable assumption that the latter is distributed with a Von-Mises probability density function [100] of concentration parameter κ, can be calculated as follows [101]: and O(·) on the right-hand side of Equation (17) is the big-O Landau symbol. Hence, taking into account these mathematical relations, an approximate estimate of the bias for the mean resultant length is: Equation (18) can straightforwardly be interpreted and shows that, as expected, the coherence bias reduces as the number of looks increases and the bias is more significant in low-coherence regions. Accordingly, to estimate the unbiased value of the coherence, the correction factor in Equation (18) might eventually be applied.
Using interferometric coherence values, changed and unchanged regions can be discriminated against one another by thresholding them. Specifically, using Equation (16), for a given detection threshold T, the probabilities of correct detection of a changed pixel and the probability of false alarm were derived in [85]: (20) To improve the detection performance, the sample coherenceγ can also be further averaged over a window and/or other alternative adaptations of the detector can be used, such as the ordered statistic coherence detector and the censored mean level detector [102]. It is worth noting that the coherence detector is sensitive not only to abrupt changes due to a disaster event (such as a fire, a flood, the effects of extreme weather conditions or man-made disturbances), which are the goals of change detection techniques, but also to other temporal decorrelation effects that are not directly related to the main event under investigation. Thus, ambiguities in change detection remain where the temporal decorrelation caused by the background environment is dominant to those components strictly related to the natural/human-induced events. For instance, this happens over vegetated areas, likely affected by wind or seasonal changes. When decorrelation caused by a specific event is coupled with temporal decorrelation from wind or rain, isolating the two different decorrelation sources becomes challenging. Tracking coherence differences between couples of SAR images has been exploited to study various phenomena (e.g., [103][104][105][106]). Some investigations address the problem of developing proper temporal decorrelation models that can be adopted to discriminate the background changes due to temporal decorrelation from those ascribed to the primary disastrous event. To overcome some limitations in low coherence regions, the work [107] proposed an alternative maximum likelihood coherence change detector that lowers the false alarm probability in low coherent areas considering the clutter-to-noise power ratio of the SAR images. Specifically, the new estimator is expressed as follows: where σ 2 n1 and σ 2 n2 are the mean noise power values of the two interfering SAR images. Over recent years, new methods have also emerged that consider both amplitude and phase information for the unsupervised extraction of surface changes from a sequence of co-registered SAR images, by extending the analyses from two single SAR images to a set of multi-temporal SAR images, e.g., see [39,92,108]. Specifically, the abrupt changes in the radar returns in a time series of SAR data can be detected using a maximum-likelihood approach by assuming that the transition times (TTs) are single or multiple. TTs can vary from one region to another, especially for long time-series of data. In this context, as a further extension of GRLT to multi-pass SAR images, the work [92] proposed a method to detect changes in sequences of dual-pol Sentinel-1 multi-pass SAR data by considering both point targets and distributed targets on the terrain. The automatic identification of the transition times is also a problematic task to achieve. Nowadays, the level of maturity of such methodologies is high; however, their use in the context of rapid mapping just after and/or between a disastrous event is still challenging and poses the base for subsequent developments, where the help of artificial intelligence (AI) frameworks can be beneficial for the identification and quantification of damaged regions.

Introduction to AI-Aided Change Detection Methods
Another emerging class of CD methods can overcome the limitations of traditional CD methods while working with a short-term sequence of multi-pass SAR data. It is represented by the methods based on AI. Machine Learning (ML) methods are increasingly applied to study and analyze data because of their ability to find relationships between target and input variables, such as those arising in satellite remote sensing of geophysical parameters. Over the last 20 years, the launch of a new constellation of RS satellites enabled the observation of the Earth system with an improved time-space granularity. Thus, unprecedented large amounts of data have been produced [109,110].
In this context, the RS community started the development of new advanced methodologies through the adoption of ML-based frameworks. Two factors drove the adoption of those frameworks: 1.
The possibility of combining information from multiple RS domains for the solution of a specific task and overcoming the limitations of traditional methodologies that, on the contrary, only allow the use of information from a particular domain [111][112][113].

2.
The development of methodologies capable of ingesting and analyzing a large amount of data automated and extrapolating useful new information [110].
Different unsupervised and supervised [31] ML methodologies have been explored in the RS CD task context.
Finally, in recent years, the combinations of unsupervised and supervised ML methodologies, namely called semi-supervised, have been proposed [126,[128][129][130][131]. The scope of these methods is to integrate the merit of both unsupervised and supervised methods to extrapolate information from RS data and construct new helpful knowledge. Semi-supervised methodologies are particularly useful in CD scenarios where it is very challenging to collect information for labeling the class to be predicted.
A complete explanation of new advanced ML-based RS methodologies and future trends is outside the scope of this work. Interested readers are referred to [80,113,132] for a detailed digression on such technologies.

Proposed Multi-Temporal SAR Change Detection Strategy
Rapid damage assessment mapping and analysis after a disastrous event exploiting SAR data stacks can efficiently be carried out using CD methods, profiting from the availability of short sequences of SAR images, encompassing the selected primary disastrous event under investigation, and regularly collected before, between, and after the event itself. SAR technology's current and future trend is to operate in a synergic way with data collected by constellations of SAR sensors that are jointly in orbit and capable of collecting data in different portions of the electromagnetic spectrum and with short revisiting times (weekly basis or less). A significant step forward in this direction has been represented by the sending into orbit, starting from 3 April 2014, of the twin sensors Sentinel-1A/B of the ESA Copernicus Constellation [133], which have long been mapping Earth's surface regularly every six days (in the most favorable case), permitting extensive analyses on the state of Earth' surface and the environment. Collected SAR data are freely available and have been fostering the development of new studies and processing methods for rapidly estimating changed areas. This section proposes and describes a joint multi-pass coherent/incoherent CD strategy that allows analyzing a short-term sequence of calibrated SLC images and fast computing damaged areas using incoherent/coherent change detection indices (CDIs). These indices are then cooperatively used to automatically extract changed areas using a random forest (RF) classifier.
First, in Section 3.1, we describe the set of synthetic incoherent/coherent indices used in our work. Subsequently, in Section 3.2, we shortly introduce the fundamentals of RF, which employ the values of the selected CDIs as inputs for binary classification (unchanged/changed). Finally, Section 3.3 details the steps of the proposed multi-pass CD strategy.

Coherent Change Detection Indices
When dealing with SAR multi-pass CD in the occurrence of a disastrous event that occurred at a known time, namely t ev , the easiest way to obtain a reliable change map is to use one single pair of SAR images collected at times t pre and t post in such a way to encompass the principal disastrous event, i.e., t pre < t ev < t post , with one single image before the event and one single image after the event. Incoherent CD methods estimate and subsequently threshold the log-ratio image map (see Equation (12)).
However, the information on the changes that occurred, ascribed to the disastrous event (e.g., a fire, an inundation, etc.), is also enclosed in all SAR data pairs that encompass the main event, i.e., more than one pre/post-SAR image can jointly be used to circumvent better the effects of noise and other disturbances in SAR images, which could lead to inaccurate decisions while using one single SAR data pair. This is especially true when coherent CD methods are exploited. Indeed, the coherence (see Equation (6)) measures the similarity between two SAR images. Thus, necessarily tracking of coherence changes requires the availability of at least three SAR images: two images are collected before the event at times t 1 and t 2 and used to estimate the (background) coherence of the scene γ pre just before the event occurred, and one SAR image is collected just after the occurrence of the event at the time t 3 , to compute the coherence change between the event, namely γ ev , see Figure 1. In this way, common coherence variations due to other causes unrelated to the main natural or man-made disastrous event can be isolated. In particular, the coherence ratio (CR) or the normalized coherence difference ratio (ND) of the pre-and post-coherence images can be calculated. The coherence ratio has long been explored for building damage assessment [134][135][136]; it is defined as follows: Considering the three-pass configuration of Figure 1, the pre-and co-event coherence maps are and , respectively. Inherited from studies of the building damage assessment, it is also the normalized coherence ratio, which is defined as follows: It is worth remarking that both the used change detection indices (CDIs) of Equations (22) and (23) depend exclusively on three SAR images; however, one image is common for the coherence estimation. If the values of the selected CDIs are associated with the common reference image, a time-series of these indices could also be generated and employed for unsupervised CD analyses, automatically tracking, for instance, the correct transition times where significant changes occurred on the imaged scene.

Temporal Decorrelation Models: Implications for Change Detection
In the previous Section 3.1., we presented two simple, coherent change indices based on the ratio and the normalized difference of a pair of coherence maps formed from a triplet of SAR images collected before and between a disastrous event. In this sub-section, we want to discuss the benefits of these two CDIs, considering the interferometric coherence properties and especially the contributions due to volumetric and temporal changes. Equation (7) shows that coherence depends on several factors. With constellations of SAR sensors characterized by a narrow orbital tube, such as the Sentinel-1A/B Copernicus constellation, the geometric decorrelation is somewhat limited and, in any case, it can be estimated and compensated for [86,88] when the topography of the area is known with suf- Considering the three-pass configuration of Figure 1, the pre-and co-event coherence maps areγ pre =γ 1,2 andγ ev =γ 2,3 , respectively. Inherited from studies of the building damage assessment, it is also the normalized coherence ratio, which is defined as follows: It is worth remarking that both the used change detection indices (CDIs) of Equations (22) and (23) depend exclusively on three SAR images; however, one image is common for the coherence estimation. If the values of the selected CDIs are associated with the common reference image, a time-series of these indices could also be generated and employed for unsupervised CD analyses, automatically tracking, for instance, the correct transition times where significant changes occurred on the imaged scene.

Temporal Decorrelation Models: Implications for Change Detection
In the previous Section 3.1., we presented two simple, coherent change indices based on the ratio and the normalized difference of a pair of coherence maps formed from a triplet of SAR images collected before and between a disastrous event. In this sub-section, we want to discuss the benefits of these two CDIs, considering the interferometric coherence properties and especially the contributions due to volumetric and temporal changes. Equation (7) shows that coherence depends on several factors. With constellations of SAR sensors characterized by a narrow orbital tube, such as the Sentinel-1A/B Copernicus constellation, the geometric decorrelation is somewhat limited and, in any case, it can be estimated and compensated for [86,88] when the topography of the area is known with sufficient details. For our analyses, the most important contributions of the coherence are the temporal and volumetric components. To address these signals, we exploit in this work the temporal coherence model proposed in [137] that incorporates the effects of temporally correlated changes and those that randomly occur (i.e., temporally uncorrelated) in the case of disastrous events, which are due to differences in the arrangement and dielectric properties of the scatterers on the ground.
The work [137] shows that temporal decorrelation for zero-spatial baseline and longtemporal baseline can be expressed as follows: where ∆T is the interferometric temporal baseline, τ g and τ v are the unique characteristic time constants for the ground and volume contributions, γ v t_rand and γ g t_rand are the temporally uncorrelated (random) complex cross-correlation changes for volume and ground layers, µ is the ground-to-volume scattering ratio and ω represents the polarimetric scattering mechanism. Using Equation (24), we can express the ratio coherence of Equation (21) as follows: geo,pre ·γ Dopp,pre ·γ mis−reg,pre ·γ vol,pre ·γ temp,pre ·γ thermal,prê γ geo,ev ·γ Dopp,ev ·γ mis−reg,ev ·γ vol,ev ·γ temp,ev ·γ thermal,ev For a narrow sensors' orbital tube (e.g., for Sentinel-1A/B SAR sensors) the variations of the geometric decorrelation effects into the two coherence maps are reduced by the coherence ratio, which is only sensitive to perpendicular baseline double difference (i.e., the difference between the perpendicular baseline of the considered InSAR data pair), and the thermal noise components almost cancel out. Assuming limited mis-registration errors and almost aligned Doppler spectra between the three involved SAR images, the most prominent contributions that persist in the coherence ratio operator are the volumetric and the temporal decorrelation effects. Considering Equation (25), the most relevant contributions are thus: where β represents the residual terms accounting for the other decorrelation artefacts that do not perfectly cancel out by the coherence ratio. When the ground contributions are prevalent, and we assume the availability of a set of SAR data acquired regularly, we can generate the pre-and co-event interferometric data pairs with the same temporal baseline, i.e., ∆T pre = ∆T ev . In this case, Equation (26) particularizes as: Similarly, if the prevalent contribution to the coherence is the volumetric component, we have: Accordingly, the coherence ratio results more sensitive than the single coherence measurement to time-uncorrelated (random) coherent components that are beneficial for coherent change detection purposes.
Let us now focus on the normalized coherence ratio estimator. Equation (24) shows that the coherence model's ground and volumetric components are two additive terms. Accordingly, the coherence difference in a somewhat homogenous region can better discriminate the two separate ground and volumetric contributions separately. For small values of the topographic phase ϕ k ∼ = 0, i.e., for flat terrains, and considering that the prevalent contributions to the coherence are the volumetric and temporal decorrelation terms, we respectively have: Equations (27)- (30) clarify that, depending on the physical phenomena under investigation, the volumetric or the ground decorrelation terms can become predominant and be emphasized by the coherence ratio and the normalized coherence difference estimators. For example, if a fire happens in a forest, we can assume that most of the changes are associated with the canopy volume consumed by the fire. In contrast, the ground is less disturbed, especially for low to mid-fire severities.
The analyses proposed in this sub-section demonstrate that the coherence ratio is preferable to the single co-event coherence map because it is more sensitive to the random time-uncorrelated changes. Moreover, in somewhat homogenous regions, when the primary event under investigation does not drastically change the structural properties of the imaged scene, the normalized coherence difference can be superior because it allows isolating the predominant time-uncorrelated coherence changes. In-depth analyses are still required to calculate the pdfs of these two operators to statistically determine the probability of false alarm and the detection probability of such two synthetic coherence detection indices. This is a matter for future research.

Extension to the Multi-Pass Case
Let us now extend our analysis to the case when M pre SLC SAR images are available before the primary event, collected at times t 1 , t 2 , . . . , t M pre , and additional SLC SAR images are available after the primary event, collected at times Information on occurred changes related to the main event that happened at the time t ev can be extracted by singularly analyzing all SAR image pairs between the primary event, see Figure 2. In particular, given the SAR data As imaged by the satellite sensors, the effects on the ground related to the primary event can be detected and extracted from every single log-ratio image {LR l } Θ 1 . Accordingly, we could generally have Θ different looks for the same scene. The time average of these different looks allows one to obtain a temporal multi-look log-ratio image which is less affected by noise. Mathematically, multi-look log-ratio image map is obtained as follows: where the summation is performed pixel-by-pixel. Note that random changes and systematic subtle multiple changes can happen in the whole t 1 , t M pre +M post time interval; accordingly, the temporal multi-look log-ratio image will tend to filter out the small random amplitude fluctuations that are present in the single maps and emphasize the major systematic contributions due to the primary event. On the other hand, if the selected time window is too large and the scene rapidly comes back to original pre-event conditions (and/or in a new state after t ev that has a reduced memory of the disastrous event), the temporal multi-look LR image will be sensitive to changes different from those we want to extract. For this reason, the proposed methodology is adequate for short-term time-series of data, with a maximum time span (temporal baseline) of one/two months (or less) between the main event, and it looks appropriate for fast damage assessment. where the summation is performed pixel-by-pixel. Note that random changes and systematic subtle multiple changes can happen in the whole time interval; accordingly, the temporal multi-look log-ratio image will tend to filter out the small random amplitude fluctuations that are present in the single maps and emphasize the major systematic contributions due to the primary event. On the other hand, if the selected time window is too large and the scene rapidly comes back to original pre-event conditions (and/or in a new state after ev t that has a reduced memory of the disastrous event), the temporal multi-look LR image will be sensitive to changes different from those we want to extract. For this reason, the proposed methodology is adequate for short-term timeseries of data, with a maximum time span (temporal baseline) of one/two months (or less) between the main event, and it looks appropriate for fast damage assessment. If dual-pol SAR data are available, the incoherent temporal multi-looking of Equation (31) can be computed for the two polarizations, separately. In addition, with dual-pol VV/VH Sentinel-1 SAR data, the sum and difference of incoherent temporal multi-looking log-ratio images can be computed and used as two separate, dependent synthetic incoherent CDIs; see the experimental results presented in Section 5.
The joint exploitation of coherent and incoherent information in a unique framework is attractive. In a multi-pass scenario with SAR images, some interferometric data pairs with the same temporal baseline before the event and co-event can be selected. The relevant coherence ratio and normalized coherence difference ratio images can be computed. Using interferometric SAR data pairs with the same temporal baseline partially compensates for the effects of temporal decorrelation not related to the primary event.
Indeed, coherence decorrelation is at the first order the same in pre-and in-between-event coherence maps. Nonetheless, as the temporal baseline increases, the single coherence maps will be more decorrelated, and the performance of the coherence detectors diminishes. For this reason, the coherence analyses should be retained to shorter time windows. This statement is supported by the experimental results shown in Section 5, demonstrating that the relative weight/importance of longer baseline coherence in the binary decision between changed/unchanged regions drastically reduces as the interferometric SAR temporal baseline increases. The developed method relies on the computation and subsequent utilization of the following two sets of coherence ratio and normalized coherence differences: If dual-pol SAR data are available, the incoherent temporal multi-looking of Equation (31) can be computed for the two polarizations, separately. In addition, with dual-pol VV/VH Sentinel-1 SAR data, the sum and difference of incoherent temporal multi-looking log-ratio images can be computed and used as two separate, dependent synthetic incoherent CDIs; see the experimental results presented in Section 5.
The joint exploitation of coherent and incoherent information in a unique framework is attractive. In a multi-pass scenario with M pre + M post SAR images, some interferometric data pairs with the same temporal baseline before the event and co-event can be selected. The relevant coherence ratio and normalized coherence difference ratio images can be computed. Using interferometric SAR data pairs with the same temporal baseline partially compensates for the effects of temporal decorrelation not related to the primary event. Indeed, coherence decorrelation is at the first order the same in pre-and in-between-event coherence maps.
Nonetheless, as the temporal baseline increases, the single coherence maps will be more decorrelated, and the performance of the coherence detectors diminishes. For this reason, the coherence analyses should be retained to shorter time windows. This statement is supported by the experimental results shown in Section 5, demonstrating that the relative weight/importance of longer baseline coherence in the binary decision between changed/unchanged regions drastically reduces as the interferometric SAR temporal baseline increases. The developed method relies on the computation and subsequent utilization of the following two sets of coherence ratio and normalized coherence differences:  (33) Note that the families of used CDIs are composed of Θ elements, representing the number of pairs of InSAR data pairs between the primary event that can be formed from the available SAR scenes.

Random Decision Forest: Basic Rationale and Application to Change Detection Analyses
This section introduces Breiman's random forest (RF), an ensemble learning method [120] widely used for classification and regression tasks, and shows how the change detection indices presented above can jointly derive surface change maps.
In the context of ML methods, Breiman's random forest (RF) is an ensemble learning method [120] widely used for classification and regression tasks. These ensemble algorithms often address non-linear and dimensional problems. Specifically, differently to Maximum Likelihood Classification (MLC) classifiers that assume a normal data distribution, the RF method is based on the Classification and Regression Tree (CART) algorithm [120,138] that does not make any assumptions regarding the frequency distribution of data. To this aim, the RF method is very suitable for inferring information of input data from remotely sensed acquisitions, which have multi-modal distribution and very rarely the normal one. In Breiman's RF, the joint use of randomness in CART and bagging techniques is the key that allows this model to perform better than many other classifiers [139].
Specifically, given a training dataset S(X, y) representing the set of inputs and y ∈ R m the corresponding outputs where m is the size of training data and d the number of features, the RF learning is pursued by applying the CART algorithm on several L training subsets S k (X k , y k ), k = 1, . . . , L, produced by random sampling S(X, y), either with or without replacement.
In this framework, each decision tree is grown by the CART algorithm that works by recursively splitting the subsets of S k (X k , y k ), randomly sampling input features d, and choosing from among those as split nodes the ones that better optimize a specific impurity function (e.g., Gini or Entropy). The recursive subdivision process stops once the maximum depth is reached, or no further subdivisions can be performed. Within the decision tree structure, subsets which are not split are called leaf nodes. A class label designates each leaf node. There may be two or more leaf nodes with the same class label. Figure 3 shows a decision tree structure for a binary classification problem.
(32) (33) Note that the families of used CDIs are composed of elements, representing the number of pairs of InSAR data pairs between the primary event that can be formed from the available SAR scenes.

Random Decision Forest: Basic Rationale and Application to Change Detection Analyses
This section introduces Breiman's random forest (RF), an ensemble learning method [120] widely used for classification and regression tasks, and shows how the change detection indices presented above can jointly derive surface change maps.
In the context of ML methods, Breiman's random forest (RF) is an ensemble learning method [120] widely used for classification and regression tasks. These ensemble algorithms often address non-linear and dimensional problems. Specifically, differently to Maximum Likelihood Classification (MLC) classifiers that assume a normal data distribution, the RF method is based on the Classification and Regression Tree (CART) algorithm [120,138] that does not make any assumptions regarding the frequency distribution of data. To this aim, the RF method is very suitable for inferring information of input data from remotely sensed acquisitions, which have multi-modal distribution and very rarely the normal one. In Breiman's RF, the joint use of randomness in CART and bagging techniques is the key that allows this model to perform better than many other classifiers [139].
Specifically, given a training dataset S( X , y) representing the set of inputs and y ÎR m the corresponding outputs where m is the size of training data and d the number of features, the RF learning is pursued by applying the CART algorithm on several L training subsets S k ( X k , y k ), k = 1,..., L, produced by random sampling S( X , y) , either with or without replacement. In this framework, each decision tree is grown by the CART algorithm that works by recursively splitting the subsets of S k ( X k , y k ), randomly sampling input features d , and choosing from among those as split nodes the ones that better optimize a specific impurity function (e.g., Gini or Entropy). The recursive subdivision process stops once the maximum depth is reached, or no further subdivisions can be performed. Within the decision tree structure, subsets which are not split are called leaf nodes. A class label designates each leaf node. There may be two or more leaf nodes with the same class label. Figure 3 shows a decision tree structure for a binary classification problem. This process allows the construction of L decision tree classifiers c k (S k ) by maximizing their strength and jointly minimizing their correlation [120,140], forming a bagged predictor able to map X → y . Finally, RF prediction is computed as majority votes overall tree predictions (see Figure 4). This process allows the construction of L decision tree classifiers c k (S k ) by maximizing their strength and jointly minimizing their correlation [120,140], forming a bagged predictor able to map X ® y . Finally, RF prediction is computed as majority votes overall tree predictions (see Figure 4). RF models further produce additional information about the analyzed data: (i) the feature importance [141] used to perform the classification or regression task and (ii) the internal structure of the data in terms of similarity of different data samples to one another. The RF model determines the former during its training (tree by tree as the RF structure is formed), estimating each feature's error rate variability when permuted while all others remain unchanged. The latter is determined in terms of proximity matrix P(i, j) Î R m´m where the element (i, j) represents the fraction of trees in which the samples i , j fall in the same leaf node. The higher the proximity measure, the more similar the pair of samples.
As with many other ML methods also, RF uses different hyperparameters. Nowadays, hyperparameter optimization is one of the most important tasks to perform. In the recent literature [36][37][38][39][40], hyperparameter optimizations have extensively been addressed, demonstrating that the correct tuning of ML model hyperparameters is mandatory to obtain high generalization performance. Specifically, RF models depend mainly on three hyperparameters: (i) the number of decision trees ( n tree ) that compose the model, (ii) the number of features ( m try ) that will be considered in the CART splitting procedure, and (iii) the tree depth k n . Only the first hyperparameter has a robust analytical demonstration [142][143][144][145][146] of its impact on forest generalization performance. For the other two, there is no theory to guide the best choice; they are strictly dependent on the problem to solve, and the related dataset of measurements [144] and this aim has to be tuned appropriately. RF models further produce additional information about the analyzed data: (i) the feature importance [141] used to perform the classification or regression task and (ii) the internal structure of the data in terms of similarity of different data samples to one another. The RF model determines the former during its training (tree by tree as the RF structure is formed), estimating each feature's error rate variability when permuted while all others remain unchanged. The latter is determined in terms of proximity matrix P(i, j) ∈ R m×m where the element (i, j) represents the fraction of trees in which the samples i,j fall in the same leaf node. The higher the proximity measure, the more similar the pair of samples.

Proposed RF-Aided CD Strategy
As with many other ML methods also, RF uses different hyperparameters. Nowadays, hyperparameter optimization is one of the most important tasks to perform. In the recent literature [36][37][38][39][40], hyperparameter optimizations have extensively been addressed, demonstrating that the correct tuning of ML model hyperparameters is mandatory to obtain high generalization performance. Specifically, RF models depend mainly on three hyperparameters: (i) the number of decision trees (n tree ) that compose the model, (ii) the number of features (m try ) that will be considered in the CART splitting procedure, and (iii) the tree depth k n . Only the first hyperparameter has a robust analytical demonstration [142][143][144][145][146] of its impact on forest generalization performance. For the other two, there is no theory to guide the best choice; they are strictly dependent on the problem to solve, and the related dataset of measurements [144] and this aim has to be tuned appropriately.

Proposed RF-Aided CD Strategy
In this section, we describe the proposed RF CD framework. As illustrated in Figure 5, the proposed framework consists of three modules. In this section, we describe the proposed RF CD framework. As illustrated in Figure  5, the proposed framework consists of three modules.  Starting from a sequence of calibrated, co-registered, and geocoded SAR acquisitions, the first module consists of pre-processing data and calculating the incoherent and coherent change detection indices.  In this work, we treated change detection as a pixel-based binary classification task that uses 1 and 0 to indicate changed and unchanged pixels. Therefore, we used an RF model combining CDIs and a reference change mask in the second module to perform supervised learning.  The final module applies a spatial average with a moving window to the RF predicted binary change mask. Eventually, the binary change mask is retrieved.

Material
We processed three SAR datasets consisting of 15, 13, and 12 SAR images collected in the Single-Look-Complex (SLC) format and Terrain Observation with Progressive Scans SAR (TOPSAR), Interferometric Wide (IW) mode by the Sentinel-1A/B sensors over the areas of interest (AOIs) of the central-western sector of Sardinia Island, the centralwestern area of Sicilia Island regions in Italy, as well as the Southern-Eastern sector of Texas coastal area, the US. Table 1 shows the list of used SAR data. The AOIs include: i. The Montiferru region in Sardinia (see Figure 6A). The territory is mainly characterized by a mountain chain and some valleys located in its inner parts, particularly in the municipalities of Santu Lussurgiu, Cuglieri, and Scano Montiferro. The terrain elevation of the investigated area ranges from the sea level to the highest point of Monte Urtigu, about 1050 m a.s.l., located in the municipality of Santu Lussurgiu. The site is historically characterized by the Mediterranean climate, presenting dry summers, cold and wet winters, and intermediate conditions in spring and autumn. In recent decades, due to global warming, the region has also been facing alterations drastically in precipitation regimes, with the most considerable precipitation runoff decrease [147]. Starting from a sequence of calibrated, co-registered, and geocoded SAR acquisitions, the first module consists of pre-processing data and calculating the incoherent and coherent change detection indices.
In this work, we treated change detection as a pixel-based binary classification task that uses 1 and 0 to indicate changed and unchanged pixels. Therefore, we used an RF model combining CDIs and a reference change mask in the second module to perform supervised learning. The final module applies a spatial average with a moving window to the RF predicted binary change mask. Eventually, the binary change mask is retrieved.

Material
We processed three SAR datasets consisting of 15, 13, and 12 SAR images collected in the Single-Look-Complex (SLC) format and Terrain Observation with Progressive Scans SAR (TOPSAR), Interferometric Wide (IW) mode by the Sentinel-1A/B sensors over the areas of interest (AOIs) of the central-western sector of Sardinia Island, the central-western area of Sicilia Island regions in Italy, as well as the Southern-Eastern sector of Texas coastal area, the US. Table 1 shows the list of used SAR data. The AOIs include: i. The Montiferru region in Sardinia (see Figure 6A). The territory is mainly characterized by a mountain chain and some valleys located in its inner parts, particularly in the municipalities of Santu Lussurgiu, Cuglieri, and Scano Montiferro. The terrain elevation of the investigated area ranges from the sea level to the highest point of Monte Urtigu, about 1050 m a.s.l., located in the municipality of Santu Lussurgiu. The site is historically characterized by the Mediterranean climate, presenting dry summers, cold and wet winters, and intermediate conditions in spring and autumn. In recent decades, due to global warming, the region has also been facing alterations drastically in precipitation regimes, with the most considerable precipitation runoff decrease [147]. ii. The Sicilian Apennines, specifically the area of the "Madonie" (see Figure 6B). Within this area is situated the Parco delle Madonie, which is the second nature reserve in Sicily. Its 35,000 hectares are home to towering mountains (at 1979 m, the highest peak is Pizzo Carbonara), large expanses of woodland, and a flourishing variety of flora and fauna. In terms of flora, there are over 2600 different species of plants, many of which are endemic to the area. Specifically, at an altitude of 1500 m, the land is entirely covered by the Madonie Forest. Below, on the hillsides, the area is mainly characterized by crops, including the cultivation of wheat, olives, and fruits. The area incorporates several historic towns and villages such as Polizzi Generosa, Petralia Soprana and Sottana, Gangi and Castelbuono. iii. The Houston metropolitan area (see Figure 6C) is the fifth-most populous urban area in the USA. The region contains the city of Houston (the most significant economic and cultural center of the South). Its port (the second largest port in the United States and the 16th largest globally) leads the US international trade. The metropolitan area is in the Mexico Gulf Coastal Plains. Much of the urbanized area was built on forested land, marshes, and prairie.

Acquisition Dates SARDINIA SICILY TEXAS
Concerning the Italian sites, the characterizing phenomena are wildfires. Wildfires have affected the Italian areas between July and August of 2021. In the Sardinia area, on 23 July 2021, at around 5.30 pm CET, a wildfire started in the Montiferru region, specifically on the road that connects Bonarcado and Santu Lussurgiu. In the first phase of its propagation, the wildfire affected an area mainly covered by herbaceous vegetation and wooded pastures. Then, on 24 July, the fire quickly got out of control due to the region's complex topography, extreme weather conditions (temperature close to 40 • C), low relative humidity, and strong winds from S-SE. The wildfire spread towards N-NW, favored by the increasing slopes and higher amounts of fuel load and started burning the forest areas of Santu Lussurgiu. The wildfire spreading caused huge damage, devastating an area of about 20,000 ha. In Sicily, the wildfire started simultaneously between 6 and 12 August of 2021. In Concerning the Italian sites, the characterizing phenomena are wildfires. Wildfires have affected the Italian areas between July and August of 2021. In the Sardinia area, on 23 July 2021, at around 5.30 pm CET, a wildfire started in the Montiferru region, specifically on the road that connects Bonarcado and Santu Lussurgiu. In the first phase of its propagation, the wildfire affected an area mainly covered by herbaceous vegetation and wooded pastures. Then, on July 24, the fire quickly got out of control due to the region's complex topography, extreme weather conditions (temperature close to 40 °C), low relative humidity, and strong winds from S-SE. The wildfire spread towards N-NW, favored by the increasing slopes and higher amounts of fuel load and started burning the forest Fire and flooding perimeters were derived using the Copernicus Emergency Management Service (EMS), which provides geospatial information within hours or days after the catastrophic phenomenon with its Rapid Mapping component. Specifically, the Rapid Damage Assessment (RDA) module of the European Forest Fire Information System (EFFIS) has been used for wildfire events. By analyzing near real-time MODIS daily images at 250 m spatial resolution, the RDA provides the daily update of the perimeters of the burn areas in Europe. The process involves delineating the extent of wildfire events based on the semi-automatic classification of MODIS satellite imagery using ancillary spatial datasets such as CORINE Land Cover, the active-fire detection products from MODIS and VIIRS. The Global Flood Awareness (GloFAS) system has been used for flood events. GloFAS provides a floods map using flood-mapping algorithms developed by the Global Flood Mapping (GFM) consortium (https://wiki.c-scale.eu/GFM/PUM/ConsortiumMembers, accessed on 1 June 2022). Specifically, three individual flood-mapping algorithms operate in parallel, analyzing historical time-series of Sentinel-1 SAR intensity data [148][149][150][151]. Topography-derived indices are used to refine the classification of water bodies. Once each algorithm generates its own "observed flood extent" map, an ensemble-based approach is used to combine the maps into a single "consensus map" in which a pixel is accepted as "flooded" when at least most of the classification algorithms that compose the ensemble classify it accordingly.

Experimental Results
Temporal multi-looked (amplitude) maps and Coherence Changes Indexes (CCDIs) extracted from the available SAR images were used to detect changes due to the wildfire and flooding events characterizing the selected AOIs. The SAR images were preliminarily radiometrically calibrated [152] to extract from the digital data the maps of the radar backscatter, i.e., the sigma naught maps [153]. We assumed the 24 July 2021 date as that related to the primary fire event. We selected a time series of acquisitions with a temporal baseline of ±6, ±12, ±18 days, before (6, 12, and 18 July), during (24 July), and after (30 July, 5 and 11 August) the wildfire event, respectively. Each SAR image of the time series was independently post-processed by applying the de-speckling noise-filtering algorithm [154] and co-registered using Enhanced Spectral Diversity (ESD) with respect to the 24 July acquisition. The maps of sigma nought were finally converted from linear to decibel, and six sigma naught differences ∆σ 0 (dB) with respect to the 24 July maps were computed and geocoded. SAR data pre-processing was accomplished using the Sentinel Application Platform tool [155].
For the sake of representability, Figure 7 shows the maps of ∆σ 0 (dB) with a temporal baseline of ±6 days for the VH and VV polarization.
Once the sigma nought difference time series was computed, the incoherent temporal multi-looking step described in Section 3.3. was singularly applied to the VH and VV SAR sigma nought time series.
As shown in Figure 8, the temporal multi-looking operation allows us to highlight better the variations of amplitude SAR backscatter ∆σ 0 at VH and VV polarizations connected to changes in ground properties related to the wildfire event.
Finally, as described in Section 2, the sum and difference of two incoherent temporal multi-look maps obtained at the two available VV and VH polarization channels were used as additional CDIs. The relevant maps are shown in Figure 9A,B.
Interferometric SAR coherence change ND post−pre and ρ post−pre CCIs were also determined and used to detect changes related to the wildfire disaster event. Considering triplets of SAR images using temporal baselines of ±6,±12,±18,±30,±36, and ±42 days with respect to the 24 July 2021 primary event date acquisition, pre-and co-fire event interferometric pairs CCDI time-series were obtained by applying Equations (22) and (23). The use of these temporal baselines was chosen to highlight as much as possible the changes between the pre-and co-disaster coherence maps. Figure 10 shows the maps of coherence ratio and normalized difference coherence obtained using the interferometric coherence related to pairs of SAR acquisitions with a temporal baseline of ±6 days (24-30 July and 18-24 July) for VH and VV polarization. For the sake of representability, Figure 7 shows the maps of Ds 0 (dB) with a temporal baseline of ±6 days for the VH and VV polarization. Once the sigma nought difference time series was computed, the incoherent temporal multi-looking step described in Section 3.3. was singularly applied to the VH and VV SAR sigma nought time series.
As shown in Figure 8, the temporal multi-looking operation allows us to highlight better the variations of amplitude SAR backscatter 0 s D at VH and VV polarizations connected to changes in ground properties related to the wildfire event.   Once the sigma nought difference time series was computed, the incoherent temporal multi-looking step described in Section 3.3. was singularly applied to the VH and VV SAR sigma nought time series.
As shown in Figure 8, the temporal multi-looking operation allows us to highlight better the variations of amplitude SAR backscatter 0 s D at VH and VV polarizations connected to changes in ground properties related to the wildfire event.  It is worth highlighting that both for the sigma nought and coherence maps, changes depend on several natural phenomena that can act simultaneously (e.g., soil moisture and humidity changes due to weather conditions), and their use in vegetated areas, such as the investigated one, is more challenging due to the presence of volume scattering decorrelation, especially at C-and X-bands. As shown in Figures 8-10, the fire event signal is easily visible and in agreement with the fire perimeter provided by the EFFIS platform. The same processing methodology was used for the Sicily and Texas AOIs. The acquisition of the 8 August 2021 and 24 August 2017 was used as main event dates, respectively. Finally, as described in Section 2, the sum and difference of two incoherent temporal multi-look maps obtained at the two available VV and VH polarization channels were used as additional CDIs. The relevant maps are shown in Figure 9A,B. Interferometric SAR coherence change ND post-pre and r post-pre CCIs were also determined and used to detect changes related to the wildfire disaster event. Considering triplets of SAR images using temporal baselines of ±6, ±12, ±18, ±30, ±36, and ±42 days with respect to the 24 July 2021 primary event date acquisition, pre-and co-fire event interferometric pairs CCDI time-series were obtained by applying Equations (22) and (23). The use of these temporal baselines was chosen to highlight as much as possible the changes between the pre-and co-disaster coherence maps. Figure 10 shows the maps of coherence ratio and normalized difference coherence obtained using the interferometric coherence related to pairs of SAR acquisitions with a temporal baseline of ±6 days (24-30 July and 18-24 July) for VH and VV polarization. It is worth highlighting that both for the sigma nought and coherence maps, changes depend on several natural phenomena that can act simultaneously (e.g., soil moisture and humidity changes due to weather conditions), and their use in vegetated areas, such as the investigated one, is more challenging due to the presence of volume scattering decorrelation, especially at C-and X-bands. As shown in Figures 8-10, the fire event signal is easily visible and in agreement with the fire perimeter provided by the EFFIS platform. The same processing methodology was used for the Sicily and Texas AOIs. The acquisition of the 8 August 2021 and 24 August 2017 was used as main event dates, respectively.

Discussion on Random Forest Training and Research Outcomes
For the three selected case-study areas in Italy and the US, we applied the RF-aided strategy described in Sections 3.4 and 3.5. In this section, we discuss the results of the experiments that were carried out using an RF model trained using as inputs the set of synthetic incoherent (sigma naught) and coherent (coherence-based) change detection indices at VV and VH polarizations.
Supervised training was carried out using the fire and flood perimeters from EFFIS and Glo-FAS by the determination of a synthetic binary change mask, in which all pixels within the perimeters are identified as changed and labeled as one and all pixels outside the perimeters are identified as unchanged and labeled as 0; see Figure 11. In addition, to avoid the class imbalance problem, an additional preprocessing step was carried out. Specifically, the dataset of samples related to the Sardinia and Sicily AOIs presents a class imbalance of about 3:1 for the positive (1) class; see Figure 12. As we can see from Figure 12A, the dataset presents a class imbalance for the positive class (1). Imbalanced datasets pose a challenge for predictive modeling such as RF because they introduce a bias in the training that results in the poor predictive performance of the minority class with respect to the majority one [156][157][158][159]. In this way, the random-oversampling and -undersampling techniques [157,158] were used to make the training dataset imbalance less severe; see Figure 12B.
After the class balancing, the dataset was split into two subsets. For the former, 60% of the total dataset was used for training and cross validation, for the latter the remaining 40% was used for the testing. Using Optuna's [143] hyperparameter optimization framework, an RF structure composed of 209 decision trees was designed. Stratified K-Fold In addition, to avoid the class imbalance problem, an additional preprocessing step was carried out. Specifically, the dataset of samples related to the Sardinia and Sicily AOIs presents a class imbalance of about 3:1 for the positive (1) class; see Figure 12. In addition, to avoid the class imbalance problem, an additional preprocessing step was carried out. Specifically, the dataset of samples related to the Sardinia and Sicily AOIs presents a class imbalance of about 3:1 for the positive (1) class; see Figure 12. As we can see from Figure 12A, the dataset presents a class imbalance for the positive class (1). Imbalanced datasets pose a challenge for predictive modeling such as RF because they introduce a bias in the training that results in the poor predictive performance of the minority class with respect to the majority one [156][157][158][159]. In this way, the random-oversampling and -undersampling techniques [157,158] were used to make the training dataset imbalance less severe; see Figure 12B.
After the class balancing, the dataset was split into two subsets. For the former, 60% of the total dataset was used for training and cross validation, for the latter the remaining 40% was used for the testing. Using Optuna's [143] hyperparameter optimization frame- As we can see from Figure 12A, the dataset presents a class imbalance for the positive class (1). Imbalanced datasets pose a challenge for predictive modeling such as RF because they introduce a bias in the training that results in the poor predictive performance of the minority class with respect to the majority one [156][157][158][159]. In this way, the random-oversampling and -undersampling techniques [157,158] were used to make the training dataset imbalance less severe; see Figure 12B.
After the class balancing, the dataset was split into two subsets. For the former, 60% of the total dataset was used for training and cross validation, for the latter the remaining 40% was used for the testing. Using Optuna's [143] hyperparameter optimization framework, an RF structure composed of 209 decision trees was designed. Stratified K-Fold [160,161] used 10-fold to cross-validate the RF model in the training step. To evaluate the effectiveness of the proposed method, we decided to use only the data related to the Sardinia and Sicily AOIs to train, validate, and test the RF model, and subsequently assess it considering the new totally unseen AOI related to the Galveston Bay, Texas. It is worth noting the overestimation performance of random forest while the dependent training and test dataset were used, as recently claimed in the work of [162,163]. However, in our case this effect is limited because the training and the validation pixels were taken from two independent datasets. Figure 13 shows the RF predicted change masks. By a spatial comparison with the binary change mask used as truth, the RF model can identify the area related to the wildfire and flooding events very well. However, it is also clear that misclassification of regions is generated by the proposed method. This is due to the presence of speckle noise consistently in each AOI analyzed. We used the confusion matrix (CM) test to quantitatively examine the RF prediction performance. considering the new totally unseen AOI related to the Galveston Bay, Texas. It is worth noting the overestimation performance of random forest while the dependent training and test dataset were used, as recently claimed in the work of [162,163]. However, in our case this effect is limited because the training and the validation pixels were taken from two independent datasets. Figure 13 shows the RF predicted change masks. By a spatial comparison with the binary change mask used as truth, the RF model can identify the area related to the wildfire and flooding events very well. However, it is also clear that misclassification of regions is generated by the proposed method. This is due to the presence of speckle noise consistently in each AOI analyzed. We used the confusion matrix (CM) test to quantitatively examine the RF prediction performance. Confusion matrices of the tested AOIs are shown in Figure 14. From their analysis, it is straightforward to note the goodness of the proposed method. We want to remark that the probability values reported in the confusion matrices were calculated starting from available polygons of burned/changed areas; thus, errors due to incorrect knowledge of the true areas affected by changes can affect the reported estimates. Confusion matrices of the tested AOIs are shown in Figure 14. From their analysis, it is straightforward to note the goodness of the proposed method. We want to remark that the probability values reported in the confusion matrices were calculated starting from available polygons of burned/changed areas; thus, errors due to incorrect knowledge of the true areas affected by changes can affect the reported estimates.
ages of the areas.
Confusion matrices of the tested AOIs are shown in Figure 14. From their analysis, it is straightforward to note the goodness of the proposed method. We want to remark that the probability values reported in the confusion matrices were calculated starting from available polygons of burned/changed areas; thus, errors due to incorrect knowledge of the true areas affected by changes can affect the reported estimates.  In addition, starting from the information given by the confusion matrices, (i) precision, (ii) recall, and (iv) F-1 scores were determined. We report these scores in Table 2. Analyzing the results shown in Table 2, from the precision and recall scores of the classification of unchanged/changed pixels, it is evident that the RF model expresses very high change detection capacity in all the three analyzed areas. F1 scores are also very high. In particular, for the Texas case the F1 scores indicate that the proposed RF model is able to identify changes with high precision and accuracy in new unseen areas afflicted by a different phenomenon.
Finally, we determined how much each selected feature contributes to the classification of changed, unchanged pixels in the analyzed phenomenon using the RF feature importance statistic. Figure 15 shows that among the coherent and incoherent set of features, the incoherent temporal multi-look VH map significantly contributes to the final estimate. The sum follows this, and the difference in mutual sigma naught ratio maps between different polarization is also relevant. Features related to CCDIs also make a high contribution. The ranking of importance reported in Figure 15 further confirms that, for what attains the SAR backscattered returns, the cross-polarization channel is more sensitive than the co-polarization one to fire damage. This result is in good agreement with the findings of the work [78] where the authors demonstrated with theoretical observations and direct evidence with Sentinel-1 SAR data that the VH channel is more apt to detect fire disturbance scars and this outcome is ascribable to the associated reduction in the volumetric scattering contribution after a fire event. Furthermore, the results show that the sum of SAR backscattered signals (sigma nought) on the co-pol VV and cross-pol VH channels has a great importance in the final estimate of changed areas.
It is remarkable that the results of Figure 15 show that the coherent indices contribute to the final estimate of changed areas but with less importance than information coming from incoherent indices. This finding was expected, knowing that coherence is more sensitive to subtle spatial-varying changes in ground properties over time after a disastrous event.
the co-polarization one to fire damage. This result is in good agreement with the findings of the work [78] where the authors demonstrated with theoretical observations and direct evidence with Sentinel-1 SAR data that the VH channel is more apt to detect fire disturbance scars and this outcome is ascribable to the associated reduction in the volumetric scattering contribution after a fire event. Furthermore, the results show that the sum of SAR backscattered signals (sigma nought) on the co-pol VV and cross-pol VH channels has a great importance in the final estimate of changed areas. It is remarkable that the results of Figure 15 show that the coherent indices contribute to the final estimate of changed areas but with less importance than information coming from incoherent indices. This finding was expected, knowing that coherence is more sensitive to subtle spatial-varying changes in ground properties over time after a disastrous event. Summarizing, the importance rank of the different coherent layers for a fire-disturbed region, which is described in Figure 15, demonstrates that: (i) Short temporal InSAR baselines are preferred to long-baselines because coherence rapidly varies after a primary event and tends to achieve a new (random) state, not linked to the primary event under investigation, just a few days after the event itself. This finding is in accordance with the fact that temporal decorrelation is sensitive not only to random changes (linked to the event) but also to composite ground and volumetric changes that determine a systematic decay of the coherence over time (see the model in Section 3.2). (ii) The normalized coherence difference has generally an enhanced importance than the coherence ratio. This finding was also expected, and it agrees with theory (see Section 3.2). Indeed, the normalized coherence difference has the beneficial effect of including in a unique estimator the advantage of the coherence difference and coherence ratio to discriminate and better isolate the random coherent components. (iii) The co-pol and cross-pol channels have almost the same importance, with a slight marked preference versus the co-pol VV polarization.

Conclusions and Future Perspective
This paper addressed a review of the theory of existing coherent and incoherent CD methodologies used to study and monitor ground surface changes that arise when a disaster event (i.e., wildfires or floods) occurs, using sequences of synthetic aperture radar (SAR) images. Then, we focused on new trends of evolution of the research branch on Change Detection with SAR data, by addressing specifically the emerging ML-driven CD approaches. We investigated the potential of different synthetic coherent/incoherent CDIs and their mutual interactions for the rapid mapping of "changed areas" relying on the joint exploitation of SAR sigma naught and interferometric coherence maps.
A classifier based on RF was trained combining different information coming from coherent/incoherent CDIs in a unique corpus, over different AOIs characterized by wildfire and flood phenomena.
Specifically, the study of the Montiferru and Madonie wildfires showed that the combined use of different CDIs able to synthetically describe ground surface changes associated with a disaster event (i.e., the pre-, cross-, and post-disaster phases) and RF classifier is a powerful way to quickly identify areas related to wildfire events and quickly assess fire damage.
The proposed methodology was also tested considering a flood event that affected the south of Texas in August of 2017, in the Galveston Bay, using S-1 SAR data.
The goodness of RF CD performance was also verified using a quantitative analysis using confusion matrix statistics that show a very high rate in identify changed and unchanged areas. The achieved F-1 scores of about 0.9 highlight the potential of RF methods in remote sensing CD tasks.
As a result, our findings demonstrate that S-1 C-band SAR data can provide suitable information on fire and flood events over the most severely affected areas (the highest level of fire or flood damage), supporting how such systems can be a complementary source of data to optical ones, in the case of cloud covers or plumes due to extreme weather and fire events.
As a further development, we intend to extend our research on change detection based on sets of heterogeneous optical and radar images. The inclusion of optical features and the fusion with the SAR one in the RF forest classifier is expected to improve generalization performance (e.g., in flooded areas). Fine tuning and extensive exploitation of the proposed method on a large scale requires the selection and processing of several independent SAR datasets and a heterogenous family of disasters (e.g., floods, droughts, extreme events, forest and vegetation disturbance, urban changes, earthquakes, volcano eruptions, manmade changes, etc.). This perspective is worthwhile and is an issue to be investigated in the future.

Data Availability Statement:
The data presented in this study are available on request from the corresponding author.

Conflicts of Interest:
The authors declare no conflict of interest.