InSAR time series and LSTM model to support early warning detection tools of ground instabilities: mining site case studies

Early alarm systems can activate vital precautions for saving lives and the economy threatened by natural hazards and human activities. Interferometric synthetic aperture radar (InSAR) products generate valuable ground motion data with high spatial and temporal resolutions. Integrating the InSAR products and forecasting models make possible to set up early alarm systems to monitor vulnerable areas. This study proposes a technical support to early warning detection tools of ground instabilities using machine learning and InSAR time series that is capable of forecasting regions affected by potential collapses. A long short-term memory (LSTM) model is tailored to predict ground movements in three forecast ranges (i.e., SAR observations): 3, 4, and 5 multistep. A contribution of the proposed strategy is utilizing adjacent time series to decrease the possibility of falsely detecting safe regions as significant movements. The proposed tool offers ground motion-based outcomes that can be interpreted and utilized by experts to activate early alarms to reduce the consequences of possible failures in vulnerable infrastructures, such as mining areas. Three case studies in Spain, Brazil, and Australia, where fatal incidents happened, are analyzed by the proposed early alert detector to illustrate the impact of chosen temporal and spatial ranges. Since most early alarm systems are site dependent, we propose a general tool to be interpreted by experts for activating reliable alarms. The results show that the proposed tool can identify potential regions before collapse in all case studies. In addition, the tool can suggest an optimum selection of InSAR temporal (i.e., number of images) and spatial (i.e., adjacent measurement points) combinations based on the available SAR images and the characteristics of the study area.


Introduction
Natural and man-made (anthropogenic) hazards threaten lives and the economy through various factors influencing a region's exposure and vulnerability.Natural hazards (e.g., geophysical, hydrological, climatological, meteorological, and biological) caused 10,492 fatalities in 2021, according to the Emergency Event Database (EM-DAT) (CRED 2022), and man-made hazards (i.e., caused by humans or close to human settlements) caused 7.93 billion US dollars global insured losses in 2020, reported by Statista (Statista Research Department 2022).The International Federation of Red Cross and Red Crescent Societies stated that "disasters happen when a community is not appropriately resourced or organized to withstand the impact" and "disasters, therefore, can and should be prevented."To do so, early warning detection (EWD) and monitoring tools can support the experts' decisions to enhance preparedness, response, and resilience (Holmes et al. 2012).
Disaster risk reduction plans are able to save lives and reduce economic damages by preparing early warning systems (Intrieri et al. 2021).According to the checklist launched by the Sendai Framework for Disaster Risk Reduction 2015-2030, four elements are involved in the framework to reduce disaster risk (Nations 2015): (i) disaster risk knowledge; (ii) detection, monitoring, analysis, and forecasting of hazards and possible consequences; (iii) warning dissemination and communication; and (iv) preparedness and response capabilities.In this work, we focus on the second element (i.e., strengthening disaster risk governance to manage disaster risk) to provide an EWD tool using remote sensing (RS) data to detect and monitor potential hazards by machine learning (ML)-based forecasting models.
RS offers a rich geodatabase of satellite images to monitor the state of stable and moving regions.Multiple interferometry SAR (InSAR) techniques have been proposed to monitor and detect ground deformations on a millimeter scale.These techniques have been employed for monitoring and post-event analysis of natural and anthropogenic hazards, such as subsidence (Herrera et al. 2009), seismic (Fentahun et al. 2021), landslide (Calò et al. 2014), flood (Chaabani et al. 2020), mining (Palamà et al. 2022), road infrastructure (Orellana et al. 2020), and other activities (Cigna and Tapete 2021).
As one of the most important economic and industrial infrastructures, mining areas are prone to ground instabilities.Despite the methodologies used to monitor mining areas, various remotely sensed data have been employed to monitor the movements in mining areas infrastructures, such as optical (Sengupta et al. 2018), SAR (Palamà et al. 2022), UAV (Rauhala et al. 2017), LiDAR (Caudal et al. 2017), and hyperspectral (He and Barton 2021) data.Regarding the performance of InSAR techniques, a wide range of studies has been performed over mining infrastructures to monitor small-to large-scale spatial movements (Perski et al. 2009;Plattner et al. 2010;Carlà et al. 2018;Kim et al. 2021).For instance, Silva Rotta et al. (2020) comprehensively investigated possible causes of the Brumadinho dam collapse using optical and InSAR data from Sentinel-1.In addition to the impact of changes in soil moisture and land use/cover, InSAR data was provided to evaluate the deformation over the tailing dam.Carlà et al. (2018) also utilized InSAR data integrated with ground-based radar data to identify the cause of a slope failure in an open-pit mine.A slope-accelerating creep was observed due to the short revisit time of Sentinel-1 satellites, which matched the source area of the failure.It should be noted that the period of monitoring and temporal sampling of time series using InSAR techniques varies based on the availability of the SAR dataset.
ML and deep learning (DL) approaches have been increasingly proposed to improve monitoring capability and reduce disaster risk by forecasting near future (e.g., 6 days to monthly intervals) instabilities (i.e., movements and deformations) using InSAR time series products.The idea of time series forecasting is based on training a model using univariate or multivariate time series (i.e., previously observed values of single or multiple parameters) to forecast future values (Hill et al. 2021).These methods mainly were focused on the time series products of InSAR techniques over different deformation case studies (e.g., landslide, subsidence, and slope failure) to warn of abnormal movements.For instance, Ma et al. (2020) proposed a deep convolutional neural network (DCNN) to detect deformations on reclaimed lands at the Hong Kong International Airport priorly.In addition to the accurately predicted results, the DCNN model could predict linear trends of the settlements on reclaimed lands and the buildings' seasonal pattern.A neural network model has also been presented to detect anomalies using a large InSAR dataset over the entire of Italy (Milillo et al. 2022).More than 170 million InSAR time series were used; the proposed approach could improve the anomalous point identification speed compared to an analytical model proposed by Raspini et al. (2018).Considering the performance of recurrent neural networks for time series prediction, the following works focused on long short-term memory (LSTM) to detect anomalies in InSAR time series.Hill et al. (2021) compared multiple structures of LSTM and conventional methods in forecasting InSAR ground motion data of the West Yorkshire coal mine subsidence.Several limitations of the models and the effect of signals nature on the forecasting performance were also highlighted.Kathirvel et al. (2021) compared the model proposed by Hill et al. (2021) to a multiscale attention mechanism LSTM architecture, which could outperform for predicting seasonal and non-seasonal deformation signals over the volcanic activity of Fogo Island, Cape Verde.In the case of multivariate LSTM forecasting, Radman et al. (2021) used environmental parameters (e.g., rainfall, groundwater, and area variations) along with InSAR ground deformations to predict the future behavior of land subsidence over an area by the Lake Urmia, Iran.Recently, LSTM has been utilized to detect and classify sinkholerelated anomalies in a supervised bidirectional architecture (Kulshrestha et al. 2022).The detected anomalies were classified based on a sudden step and sudden velocity changes in deformation time series.Bidirectional LSTM has also been used by Lattari et al. (2022) to obtain change points affected by volcanic activities over the Fernandina volcano area.Raspini et al. (2018) was again considered as the baseline, and a time-gated LSTM was designed to use sampling intervals as additional information during learning.The proposed architecture identified three classes of trend changes, including step, velocity, and step and velocity (Lattari et al. 2022).
Based on our best knowledge of the current state of the art, there is a gap in providing a general tool to exploit the capability of InSAR time series in detecting spatio-temporal anomalies.In addition, no research has considered the adjacent time series of InSAR measurement points.This study presents a framework of InSAR time series forecasting using the LSTM model and considering the neighbor time series.The outcome of the proposed framework prepares spatio-temporal-and ground instability-based input for early warning systems to be integrated by other parameters for activating alarms over high-risk regions with the minimum false alarm rate.The method is tested over three failure cases in mining areas in Spain, Brazil, and Australia.Eventually, the goals of this work are (i) presenting an approach to support early warning detection (EWD) tools using InSAR time series products to provide an interpretable platform for experts in natural and human-initiated risk assessment; (ii) preparing an LSTM architecture employed by tuned parameters to forecast short timesteps (e.g., 3, 4, and 5); (iii) proposing an early alert detector to detect four types of behaviors, including normal, outlier, noise, and early warning; and (iv) investigating on impact of InSAR temporal (i.e., minimum number of images) and spatial (i.e., adjacent InSAR measurement points) characteristics in detecting regions for potential collapse in near future.
The reminder of this paper is structured as follows: "Case studies and InSAR data" section provides information on three collapse cases and InSAR datasets utilized in this study.Afterwards, the methodology is explained in "Methodology" section, including LSTM model characteristics and the algorithm of EWD tool.The performance of proposed tool and LSTM are presented in "Results" section.In "Discussion" section, the contribution and limitations of this work are discussed, along with the several suggestions for future studies.Finally, the conclusion is in "Conclusion" section.

Case studies and InSAR data
Three cases of collapses over mining sites are considered in this study to evaluate the capability of the EWD tool.These cases are in Spain, Brazil, and Australia: Cobre Las Cruces, Córrego do Feijão Mine (Brumadinho), and Cadia Valley Operations (Cadia).In the following sections, the case studies and characteristics of the hazards are explained, along with the InSAR datasets.

Cobre Las Cruces, Spain
Cobre Las Cruces is the largest open-case copper mine in Europe, located in the municipalities of Gerena, Salteras, and Guillena in the Seville province of Spain (Fig. 1).It is an open pit mine using a hydrometallurgy process to obtain  Florencio)), a massive mining-induced landslide occurred (Fig. 1a), covering 675,000 m 2 , including the slopes of the open-pit mine, where the wave of toxic waste traveled 1300 m from the mine's waste management plant to the bottom of the mine (50% of the mine's waste management plant was collapsed, Fig. 1c), along with the northern slope of the mine (Rodríguez 2019).It was also reported that this collapse contained an estimated 15 million m 3 (Rodríguez 2019).Fortunately, the collapse occurred during the shift change in the mining work, and there was no personal misfortune, although it took several months to remove material deposited at the bottom of the open pit mine, rebuild the slopes, and resume the extraction and treatment of copper (Jesús Florencio).Based on local geologists' point of view (Jesús Florencio; José 2019), the movements on the Niebla-Posadas aquifer, where the piezometric level was forcibly reduced, caused a new geotechnical instability creating a gigantic discontinuity surface that was not considered at all in the technical calculations carried out for the project.

Brumadinho, Brazil
The Córrego do Feijão Mine is one of the complexes of an iron ore mining company, located in Brumadinho, Minas Gerais State, Brazil (Fig. 2).The inactive Dam I of Córrego do Feijão Mine collapsed on January 25, 2019, known as the Brumadinho dam disaster (Fig. 2a).The dam's height was 86 m and had a crest length of 720 m, a waste disposal area of 249,500 m 2 , and a volume disposed of 11.7 million m 3 (Porsani et al. 2019).In addition to the dams, the complex included dams, an administrative center, a dining hall, a maintenance office, a cargo terminal, and a miniature railway network for the transport of iron ore (Fig. 2a).The failure caused 270 people to die after the dam released a mudflow through the complex infrastructures.The dam was classified as having a low risk of severe potential damage and went through Fig. 2 The Córrego do Feijão Mine in Brazil.The overview of the mine site in a is as follows: (1) the mine, (2) the overburden dump, (3) water dam VI, (4) Dam I, (5) the stockpile area, (6) canteen, and (7) Vale's administrative center (Porsani et al. 2019;Robertson et al. 2019).The first observed deformation was recorded by a video camera after 18 s (b) and 6 min and 25 s (c) (Robertson et al. 2019).The mud wave affected the area after the collapse (Silva Rotta et al. 2020) (d) multiple field monitoring, including two local stability declarations, biweekly field inspections, piezometers and water-level indicators, ground-based radar, video monitoring system, siren alert system, and downstream population registration (clarifications regarding Dam I of the Córrego do Feijão Mine 2019).However, about 12 million m 3 of a tailing-mud mixture was released after the collapse, destroying 8.5 km up to the Paraopeba River, extending for more than 300 km along the bed of the Paraopeba River toward the São Francisco River with waves up to 30 m high (Fig. 2d) (Porsani et al. 2019).It was reported by an expert panel that "a slope failure within the dam starting from the crest and extending to an area just above the First Raising (the Starter Dam).The dam crest dropped, and the area above the toe region bulged outwards before the dam's surface broke apart," based on recorded videos (Robertson et al. 2019).It was also investigated that vertical subsidence, up to 30 cm, was moving within the past 12 months before the dam collapse.The subsidence occurred due to seepage erosion, saturating the tailings dam, and sediments removal from the fill (Silva Rotta et al. 2020).

Cadia mine, Australia
Cadia Valley Operations (Cadia), located in central western New South Wales, Australia (Fig. 3), owns Australia's largest underground gold mine and one of the world's largest gold and copper deposits.A section of the northern dam wall, which links two tailing storage areas, collapsed into the southern tailings dam on March 9, 2018 (Fig. 3a) (Petley 2018).Two 2.7 magnitude earthquakes in shallow depth of 10 km happened a day before the failure (Geoscience Australia 2021).Although there is no evidence that seismic events caused the collapse (Thomas et al. 2019), seismic activities can still be considered as the reason related to the collapse due to the shallow depth of earthquakes and active fault observed in the seismic catalogue (Geoscience Australia 2021).Additionally, construction activities took place before the day of failure due to satellite imageries, along with several cracks observed a few hours before the collapse (Thomas et al. 2019).According to Thomas et al. (2019), a mixture of transient and progressive displacement was observed at

InSAR datasets
In this study, InSAR time series were generated using SAR Sentinel-1 A/B interferometric wide swath (IW) single look complex (SLC) images at full resolution.Table 1 provides the number of images and the time span of Sentinel-1 images for each case study.According to the revisit plan of the satellite, the minimum temporal sampling was 6 days for Cobre Las Cruces and 12 days for Brumadinho and Cadia.Due to technical issues, the observation intervals were not equal in the time series.
The InSAR processing was performed using the Persistent Scatterer Interferometry chain of the Geomatics (PSIG) Research Unit of the Centre Tecnològic de Telecomunicacions de Catalunya (CTTC), described in Devanthéry et al. (2014).The main steps of the PSIG processing implemented in this study are (1) generation of consecutive interferograms; (2) selection of points based on the dispersion of amplitude; (3) estimation of the residual topographic error and subsequent removal from original single-look interferograms; (4) phase unwrapping of the consecutive interferograms, generating a set of N unwrapped phase images, which are temporally ordered in correspondence with the dates of the SAR images processed, referred as time series of deformation (TSD); and (5) geocoding of the results.
The output of the InSAR processing is a deformation map composed of a set of selected points, called persistent scatterers (PS), with information on the estimated lineof-sight (LOS) deformation velocity and the accumulated deformation at each Sentinel-1 image acquisition time.Thus, in addition to the dispersion of amplitude criterion, parameters such as the temporal consistency among points and coherence values of interferograms were considered to ensure that the analysis involves only well-processed points.Besides, distance between selected points and the number of edges connected to a fixed point were made to satisfy a given threshold.

Methodology
A supportive EWD tool provides spatio-temporal analyses to detect potential future movements.This study proposes a supporting mechanism for EWD systems to supply spatial and temporal states of unstable regions as an input to early alarm systems to reduce the false alarm rate of high-risk alarm activation in vulnerable areas.Thus, the proposed tool attempts to introduce regions characterized by various ranges of anomalies inside their InSAR time series, where potential failures can occur.The tool's outcomes can be employed as informative inputs to EWD systems to be integrated by other technical and decision-making factors to activate alarms on high-risk targets finally.
Our tool is ensembled based on two main elements: the forecasting model and the early alert (EA) detector.We selected the LSTM model regarding its ability and reliability in forecasting the InSAR time series (Hill et al. 2021;Lattari et al. 2022).The idea of employing a forecasting model is to find possible anomalies in the last time steps (e.g., 3, 4, and 5 observations) of the InSAR time series, which could indicate potential future movements.Thus, all sets of time series are divided into three parts, including training and test splits of the model and the anomaly period.Afterward, all InSAR time series plus corresponding predicted values and adjacent time series are the input of the EA detector.Each time series is evaluated based on multiple criteria considering the predicted values and adjacent time series to calculate an indicator, presenting the potentiality of a point for a significant collapse.The following sections comprehensively explain the forecasting model and EA detector.

LSTM model and magnitude of anomaly
The common configuration of a classic artificial neural network (ANN) is composed by input layer, hidden layers, and output layer.The input layer enters the system with the initial data and feeds the following layers (hidden layer).Artificial neurons in the hidden layer, which lies between the input and output layers, receive a set of weighted inputs and generate an output using an activation function.The final layer of neurons that generates results from ANN modeling is known as the output layer.The input-implicit-output layer is completely connected to conventional neural networks.The steps in the sequences are unrelated to one another.Therefore, time series prediction cannot be done with a conventional neural network (Chen and Chou 2012).The recurrent neural network (RNN), which permits network feedback (Sak et al. 2014), stores the previous information and uses it to inform the calculation of the current output.The problem of "gradient explosion" may be effectively solved using a special type of RNN, called LSTM (Hochreiter and Computation 1997).
Long-term memory and data discovery are guaranteed by the LSTM's forgetting gate, input gate, and output gate.For the input and removal of data transmission, the three gate functions provide a dependable nonlinear control mechanism.
The LSTM model forecasts three forward time steps, including 3, 4, and 5, which are input values for the EA detector.Figure 4 shows how a time series is subdivided into train, test, and anomaly parts.The last part refers to the last values of a time series, where a potential ground motion can be triggered.The LSTM is separately trained for each case study and forecasting range, and time series are parallelly entered to the models to get the best combination of parameters for each multistep prediction (i.e., 3, 4, and 5) for each study area.
Hyperparameters of the LSTM models are tuned for every multistep forecasting and case study.Table 2 contains

Early alert detection
As mentioned above, we select three multistep anomaly periods to detect potential significant movements, including 5, 4, and 3 SAR observations in every time series.Considering the next step observations in each time series and corresponding predicted values and adjacent time series, each InSAR observation can be classified as one of the following four types of alerts: normal, outlier, noise, and potential  anomaly.Figure 5 contains the visual example of the four types of InSAR time series behavior within the anomaly period.Figure 5a shows a normal point, in which the predicted value meets an interval of the real observation (the threshold is explained in this section).A normal point demonstrates that the observation follows the dominant trend of the time series, referring to non-anomaly behavior.An outlier alert is shown in Fig. 5b, indicating that the predicted value locates out of the interval of the real observation, however, inside the next observation interval.This time series may be characterized by an anomaly in the following steps, but Fig. 5b observation is out of this probable trend.Figure 5c also shows another example of observation, which does not display our demanded deformations, called noise.This alert identifies a point which is not in the interval of corresponding observation neither the mean of next observation of adjacent time series interval.Additionally, a noise locates in opposite side of the two observations.Finally, type 4 alert warns an early alarm for the possible significant movement.In this alert, similar to the noise alert, the predicted value is not inside of both observations' intervals; however, it is between two observations, which are along the same direction of the deformation/dynamic trend.
It should be noted that each predicted value is evaluated by the corresponding observations of its and adjacent time series.The adjacent time series are chosen from neighborhood time series in radiuses of 20, 30, and 45 m, where the impact of this parameter is also discussed in the "Results" section.To classify the observations, four conditions are employed based on a threshold (Eq.1), computed by the residual standard deviation (rstd) of the time series out of the anomaly period.An order 3 polynomial function is fitted to calculate the residual values in each time series.
where x and x est are the observation value and the estimated value.
Algorithm 1 elaborates the proposed procedure to detect potential time series warning an unstable movement.The entire set of time series is imported one by one through the first loop with the length of number of time series.X vector contains a full-length time series excluded the anomaly period, A p , where the rest of observations is in the X p vector.So, the length of the X is "n steps," where n is the length of a time series with the entire observations.The forecasting model (LSTM) predicts corresponding values of A p using X.This model has been so far tuned as explained in the previous section.Thus, P vector contains predicted values considering the chosen anomaly period (i.e., forward multisteps 3, 4, and 5).Then, observations of adjacent time series, X n , are found based on the given neighborhood distance: 20, 30, and 45 m.The rstd is also computed using Eq. 1.The second loop begins after the above steps, containing the detection of alert types for every observation in the A p .The first Additionally, a noise is not inside the interval of corresponding observation neither the mean of next observation of adjacent time series.Finally, an observation characterized by an early warning alert satisfies the noise detector condition; however, the observations follow a positive or negative slope.Moreover, each observation is labeled based on its step and alert type assigned, in such a way that 0: normal, − 1: outlier and noise, and 2: early warning.Then after, the labels multiply by the numerator of each step to emphasize the impact of recent anomalies or possible noises.Through all these steps, the difference of real observation and predicted values is computed.At the final stage, the summation of alerts, EA, and differences, D, of each time series are calculated as the indicators of EWD.
It is worth noting that the minimum and maximum ranges of EA values vary for the selected forecast ranges.Indeed, the lowest values of EA are calculated − 6, − 10, and − 15 for 3, 4, and 5 forecast ranges, respectively, where the highest values are 12, 20, and 30.In this study, we consider 6 (i.e., two consecutives of early alerts in the first and second time step) as the indicator of potential anomaly.This value can be changed to higher ones based on the decision of experts for finding more significant or recent anomalies.Therefore, histograms of the EA and D are generated to designate the population of time series in different intervals of EA and D. Since the scope of this study is to detect time series with significant anomalies, we mostly focus on time series groups with high EA values.

Results
Before analyzing the EWD results, the performance of the LSTM model is reported for the case studies.Then, histograms of magnitudes of differences between predictions and real observations and anomaly estimations are provided to enhance the interpretation of the EWD results, as well as maps of these results to demonstrate how the tool could accurately detect regions for potential collapse.Additionally, the evaluation of selected time steps and radiuses on the effectiveness of the proposed methodology is provided.It should be noted that several results are supplied in the Appendix to better convey the article's scope in the main body.3. The most accurate results were obtained for the Cadia dataset in all forecasting ranges.Additionally, the best parameters were selected through a hyperparameter tunning process using values in Table 2.Moreover, Fig. 6 shows a time series sample with the prediction values and neighbor time series in 20 m.In this time series, an anomaly was triggered, causing a significant difference between the predicted values (red time series) and the real observations (blue time series).Additionally, five neighbor time series were found in a 20-m radius, which has been demonstrated in this figure.Since an abrupt change characterized the real-time series before the anomaly period, the proposed method estimated an EA amount of 20-and 150-mm magnitude of differences, D.

Analysis of anomalies for Cobre las Cruces
Histograms of estimated EA and magnitude of differences were provided to indicate the distribution of the estimated values in various intervals, including noise or outliers, normal, and significant anomalies.The provided histograms help to find the number of points with high values of D and EA, indicating potential time series characterized by significant anomalies.Figure 7 includes the results of EA estimation histograms for all anomaly periods and radiuses of the Cobre las Cruces dataset, containing nine different combinations (see Appendix for several examples of histograms).The bin widths of EA histograms were selected 3, 4, and 7 to properly categorize EA values for detecting safe 374 Page 10 of 23 and unsafe (i.e., potential) time series.The x-axis refers to the bin widths' intervals, and the y-axis shows the percent of time series detected in each interval.As shown in the subfigures of Fig. 7, the range of three anomaly periods is not equal because a more extended forecast range analyzes more observations to detect early warnings.Figure 7 also shows that a high percentage of time series was clustered as safe targets (i.e., EA < 6 ), which was expected.
A more considerable EA value refers to continuous and recent anomalies regarding the potential time series.To  this end, this tool can trigger an alarm for those areas with high EA values for further analysis considering the fact that a recent and notable failure can happen.The following results of differences' magnitude and maps are presented to address the level of anomalies and spatial position of potential time series.
Figure 8 shows the estimated D, the sum of the differences among real observations and predicted values, using the LSTM model.These results are the outputs of histograms for each forecast range in the Cobre las Cruces dataset (see examples in Appendix).The bin width was selected as 50 mm to properly categorize time series characterized by considerable differences.Overflow and underflow bins were also chosen to maintain a consistent number of bins among the radiuses results.It was also concluded from Fig. 8 that all plots approximately meet a normal distribution shape.The maximum number of time series had safe amounts of D (i.e., bars in the middle of graphs).Time series with maximum and minimum D can provide valuable information about the degree of anomalies.For instance, the last bins at the left and right of the bar charts indicate the most significant difference among the InSAR measurements and predictions.Since the dataset was a descending pass, two or three last bins on the right of the graph are anomalies that move far from the satellite.
The results, as mentioned earlier, are also shown in Fig. 9, a map of Cobre las Cruces to demonstrate the spatial distribution of outcomes and how the tool was precise in detecting potential regions where failures happened.Figure 9 includes EA and D estimations of the Cobre las Cruces dataset for a radius of 20 m and a forecast range of 5.In "Impact of forecast range and radius in EWD" section, further explanation of the results of the proposed tool is provided considering the chosen radiuses and forecasting ranges in performance of the tool in EWD.
Figure 9 illustrates the proposed tool's capability to forecast the potential deformations over the Cobre las Cruces.A red polygon in Fig. 9 has highlighted the collapsed region.It is evident in Fig. 9a that a high concentration of potential targets was located inside the drawn area.Other potential points are also in other places far from the collapsed zone.Additionally, Fig. 9b shows the magnitude of D over the case study.Multiple potential targets (i.e., more than − 50 mm) were detected over the collapsed zone, illustrated as reddish colors.There is also a concentration of high values of D in the west of the pit, which may provide valuable information regarding the possible causes of the collapse.

Analysis of anomalies for Brumadinho
Figure 10 provides the result of EA values' histograms of the Brumadinho dataset, with similar characteristics to Fig. 8. Since the point density of this case study was not as large as Cobre's, the graphs' shape is different.An equal number of time series were detected as safe and unsafe points in most combinations, possibly due to the lack of PS detected in this area.As stated before, the time series with large EA values are crucial regions for further analysis.Over the Brumadinho region, a greater percentage of time series were characterized by high EA values.Additionally, the three selected radiuses' forecast range results indicate similar behavior.
Regarding the D results, we followed an identical methodology, as explained in "Analysis of Anomalies for Cobre las Cruces" section, with different overflow and underflow bins.
Figure 11 also meets a normal distribution shape with a bin width of 50 mm.Although a limited number of time series had the highest and lowest differences values, these are crucial outcomes that need more interpretation.In Fig. 12, their location is highlighted as red.In all three forecast ranges, the minimum and maximum D values are almost equal for Brumadinho and Cobre las Cruces areas.Since these values vary between 150 and 200 mm, significant displacements were triggered earlier than collapses over these areas.

Analysis of anomalies for Cadia mine
Like Cobre las Cruces and Brumadinho datasets, the histograms of 9 combinations of forecast ranges and radiuses were generated for the Cadia dataset to facilitate the interpretation of the estimated EA and D values.As shown in Fig. 13, the percentages of potential time series detected with forecast ranges of 4 and 5 are notably higher than 3.In detail, around 20% of points were characterized by significant anomalies (i.e., EA ≥ 6 ) in all spatial distances of 3 forecast ranges; however, more than 45% of time series were categorized in other forecast ranges.
On the other hand, the graphs of D values for the Cadia dataset follow a normal distribution shape, where −50 ≤ D ≤ 50 locates among the center bars in Fig. 14.Moreover, the minimum and maximum D values of the Cadia dataset are smaller than Cobre las Cruces and Brumadinho, indicating that a lower level of instant deformation occurred over the region.It can also be concluded from Fig. 15 that the propagation of potential points is less than Cobre las Cruces dataset.In the eastern parts of Fig. 15b and west of the purple ellipse (i.e., where the dam collapsed), several potential targets were detected, which were addressed as possible triggers in a previous study (Thomas et al. 2019).For instance, in the east of the dam, construction activities were reported 2 days before the failure, increasing the dam's height.The visual evidence presented the presence of a bulldozer, and continuous construction activities started several days before the hazard (Thomas et al. 2019).

Impact of forecast range and radius in EWD
Table 4 contains the percentage of time series with the possibility of failure, i.e., (EA ≥ 6 and D > 50mm) and (EA ≥ 6 and D < −50mm) , inside the surfaces of dams in Brumadinho and Cadia and open pit of Cobre las Cruces, where the collapses happened over the case studies for all combinations of 3 forecast ranges and 3 radiuses.In all case studies, the number of detected time series as potential targets increased by the forecasting range.Additionally, this increase has been seen in the length of the radius.For instance, around 3-5% of potential time series were detected by increasing the forecast range and the radius.In addition, Table 4 incorporates the percentage of possible failures within the areas surrounding the critical regions of the case studies.Upon reviewing the Fig. 10 The percentage of EA estimations for all combinations of forecasting ranges (i.e., 3, 4, and 5 time steps) and radiuses of the Brumadinho dataset.The bin widths of forecasting ranges were selected as 3, 4, and 7, respectively table, it becomes apparent that the values decrease in the Cadia case about 8-12% and 3-6% in Brumadinho, whereas the differences observed in the Cobre las Cruces results are relatively minor (less than 1%).This disparity is primarily attributed to the fact that the selected critical area in the Cobre las Cruces case study is approximately 4-5 times wider than that of Cadia and Brumadinho.As for the smaller differences observed in Brumadinho when compared to Cadia, it can be attributed to the limited number of scatterers in the land cover of the Brumadinho region, resulting in a density that is not substantially larger than that of the dam's surface.
The inclusion of this table enhances the interpretation by providing insights into the percentage of time series with potential failure possibilities across the entire set of PS points, rather than solely focusing on the area of failure.Furthermore, it enables the generation of a comparable standardized indicator, allowing for meaningful comparisons between different case studies, depending on the type of critical areas (e.g., dam or open pit).Moreover, it offers a perception of the extent to which the provided EA and D values are representative.
The forecast range and radius can be selected due to various factors, including the point density, the average temporal interval of SAR observations, and the amount of noise and outliers.For instance, the number of time series affected by noise and outlier observations in Cadia mine is less than in other datasets.The noise and outlier observations decreased the level of EA values in a forecast range.Additionally, the propagation of PS points in the Cadia dataset was not circular due to the dam and vegetation influencing the number of neighbor points.Moreover, the 6-day InSAR data was only feasible over the Cobre las Cruces dataset, supplying shorter time intervals.It provides a shorter day among the forecasting ranges to perform a more reliable prediction.Indeed, the minimum days of forecasting range 5 are 1 month, which can be 2 months or more for the Cadia and Brumadinho case studies.Also, the point density of the Cobre las Cruces is almost five times larger than other datasets, preparing more information for the ML forecasting model, which positively impacts training the model.Apart from the fact that it may generate more undesirable observations (e.g., outlier and noise), a more reliable spatial analysis can be achieved due to the high density of PS points.More adjacent time series can provide ancillary information to reduce the false alarm rate for detecting early alert observations.Therefore, it can be proposed that the forecast range and radius should be selected considering the characteristics of the InSAR dataset and the case study, such as point density, the average temporal interval of SAR observations, and the amount of noise and outliers.

Discussion
We provided an initial-level tool to support a technical procedure for EWD and early alarm systems.This tool can employ assumptions over specific areas to show how the integration of ML and InSAR is powerful over mining areas using temporal (i.e., forecast range) and spatial (i.e., adjacent time series) for detecting regions with the possibility of a near-future collapse.Experts can develop the proposed tool (e.g., geologists and risk assessors) to implement in their specific areas, considering the type of infrastructure and vulnerability.
In this study, we presented the impact of chosen distances for spatial investigation of the EWD.At first glance, all radiuses in the three forecast ranges have almost equal percentages of detected time series.It may indicate that the proposed tool is robust to detect potential anomalies, although it may refer to the poor impact of selected spatial analysis.However, the chosen radiuses can demonstrate at what distance a possible failure may happen.
A deeper analysis of the impact of selected radiuses can be investigated in future studies because this work provides technical support for presenting a tool to emphasize the capability of InSAR and ML for proposing potential regions for probable future failures.Indeed, the outcomes of this tool must be further analyzed by geologists, on-site engineers (miners in our case studies), and sustainability experts to finally label an area as a dangerous target.
Besides, we have proposed these forecast ranges due to the short revisit time of Sentinel-1.However, a forecast range of 5 SAR observations for non-European areas reaches almost 2 months, which could be categorized as a  13 The percentage of EA estimations for all combinations of forecasting ranges (i.e., 3, 4, and 5 time steps) and radiuses of the Cadia dataset.The bin widths of forecasting ranges were selected as 3, 4, and 7, respectively long period for an early warning.Additionally, this study's spatial analysis (i.e., selection of radiuses) was performed based on the Sentinel-1 proper 14 × 4m pixel size, providing spatio-temporal information from the adjacent areas.
Regarding the outcomes of the proposed EWD tool, the assumptions (i.e., potential time series) can be more conservative regarding the false alarm level of the infrastructure.It means a mine site must be closed regarding the false alarm level.For instance, a high range of EA and D values from conservative parameters-forecast range of 3 and 20-m radiusmay trigger an alarm for a probable failure over the case study in case of reliable and available short-term SAR images.
Despite the contribution proposed by the EWD tool, several limitations were causing unreliable information and increasing the false alarm rate.First, it is important to note that the availability of SAR images in both ascending and descending modes plays a significant role in fully harnessing the potential of InSAR outcomes.In our study, we focused solely on the descending mode of Sentinel-1 due to data availability.However, future investigations should explore the capabilities of both ascending and descending pass datasets in enhancing the Earth's surface displacement analysis.Incorporating the ascending pass data could provide additional displacement information beyond LOS measurements.Therefore, future studies should consider the integration of both passes to further improve the effectiveness of the Earth surface displacement monitoring tool.Second, there were limitations in obtaining enough measurement points over areas covered by vegetation and specular targets (e.g., water), where no PS point was detected.Over the Brumadinho mining site, the PSIG processing chain was not able to get an adequate number of PS points; however, using distributed scatterer (DS) technique could improve the analysis of Brumadinho case study since there were limited numbers of PSs in this area regarding the land cover.Given the inherent susceptibility of mining sites to substantial movements, the installation of permanent corner reflectors can greatly enhance the quality of consistent monitoring.By strategically placing corner reflectors or other persistent reflectors of SAR signals during the construction of dam embankments (e.g., dikes), a systematic monitoring approach can be established to track displacements and detect sudden movements within critical infrastructures.Third, accurate time series forecasting by ML models mainly needs long time series characterized by a few rates of noise.Since Sentinel-1 serves from 2015 and the InSAR time series cannot be generated by constant 6-or 12-day observations, the forecasting input dataset encounters irregularity in the time series.Additionally, as stated in various studies of InSAR time series, filtering observations known as noise or outliers is not recommended.In addition to the complexity of defining and detecting them, they may include information regarding abrupt deformations.Furthermore, the coverage of measurement points detected by InSAR techniques affects the quality of the spatial analysis.For instance, in a distance of 20 m for a PS point in the middle of InSAR products from Sentinel-1, 12 adjacent PS points must exist with a pixel size of 14 × 4m .However, the low density of measurement points, caused by the technical issues and lack of PS, decreased the number of adjacent time series.Additionally, the spatial analysis of the proposed tool is based on the availability of adjacent PS, which may affect the evaluation of single points.In fact, there would be single points (or isolated points) over the area of interest that indicate a hazardous condition, but they are not spatially analyzed in our proposed tool.On the other hand, it is worth noting that the temporal evaluation is performing on all PS with or without adjacent PS in the chosen radiuses.Finally, temporal sampling of SAR observations, Sentinel-1 in this study, can only offer the minimum of 6 or 12 days (note that before this study, the minimum temporal sampling of Sentinel-1 was 6 days and will be shorter in the future), which may not be adequate input for EWD systems for daily monitoring and risk assessment.Also, shorter temporal sampling can improve the quality of InSAR time series over areas like the Brumadinho site since shorter temporal baselines usually enhance the coherency of InSAR datasets.
It is worth noting that considerable contributions have been investigated in this study, which is not the primary goal for proposing an EWD tool.Therefore, these findings can be discussed in the subsequent studies of this work.For instance, noise and outlier detection can be further addressed due to their high impact on EWD.Outlier values are also recommended to be refined to normal alerts by performing an extra threshold based on the minimum standard deviation in a time series.Additionally, time series characterized by vertical jumps-phase unwrapping errors-were observed in the extremely high range of D values, which were greater and smaller than the overflow and underflow bins.

Conclusion
In this study, we have proposed a technical supporting tool to EWD systems on ground motion data using InSAR time series and ML model, LSTM, over three collapsed mining sites.It was concluded that this tool can provide technical support for experts to interpret and activate early alarms over vulnerable regions since the tool is adjustable considering the available data and characteristics of the case study.The proposed tool included a forecasting model and a novel methodology, EA detector, to find potential time series for a possible significant movement in close future.The LSTM model was trained and tested to accurately forecast the anomaly period, including 3, 4, and 5 time steps.Additionally, we presented a new contribution to InSAR time series forecasting, spatial forecasting.Indeed, adjacent time series in 20-, 30-, and 45-m distances was added to the tool to improve the reliability of analysis by decreasing the impact of noise and outlier in InSAR measurements.This spatio-temporal tool could indicate regions with the possibility of collapse before the event in three case studies, where significant collapses happened.The results indicated that the various combinations of forecast ranges and radiuses showed similar results over the Cadia mine site, although there were other conclusions over the Cobre las Cruces and Brumadinho.This suggests that the forecasting range and distance for adjacent must be selected by experts considering the limitations in available SAR images (i.e., temporal point of view) and the number of detected PS points (i.e., spatial point of view) over the study area.Also, the forecasting model could be undergone of further developments to improve the accuracy of InSAR time series prediction over those areas with short time series and forecast ranges.

Appendix
Fig. 16 Examples of EA and D histograms for three datasets with various radiuses and anomaly periods (APs)

Fig. 1
Fig. 1 Cobre Las Cruces case study in Spain.a The pre-event image of the site and the red circle shows where the collapse happened.A closer view of the aftermath landslide is shown in b on January 25,

Fig. 3
Fig. 3 The map of the Cadia mine site in (a): (1) Cadia Hill open pit; (2) Cadia East underground mine; (3) South Waste Rock Dump; (4) Northern tailing storage facility; and (5) Southern tailing storage facility (Thomas et al. 2019).The collapse area has been highlighted the parameters and values used for the hyperparameter tunning.Root mean square error (RMSE) and mean absolute error (MAE) are computed to assess the performance of the model.After finding the best parameters, the predicted values are computed to be compared with InSAR values to calculate the magnitude of difference among the predicted values and real values, indicating the beginning of anomalies for the following time steps.The difference of predicted steps and corresponding observations is summed as the anomaly amount for each time series, called D. Indeed, D demonstrates what amount of instability was detected inside the most recent InSAR displacements, compared to previous observations.

Fig. 4
Fig.4A time series sample divided to train, test, and anomaly period.The anomaly period is compared with predicted by the LSTM model.The proposed EWD procedure is implemented on anomaly period

Fig. 5
Fig. 5 Visual examples of four types of alerts, proposed in this study: a normal, b outlier, c noise, and d early warning observations

Fig. 6
Fig. 6 Example of a Cadia time series and adjacent time series, along with the prediction values of the forecast range of 5 and 20-m radius.The zoom-in figure broadens the real and adjacent time series and predicted period

Fig. 8 Fig. 9 ◂
Fig.8The percentage of D estimations for all combinations of forecasting ranges (i.e., 3, 4, and 5 time steps) and radiuses of the Cobre las Cruces dataset.The bin widths of forecasting ranges were selected at 50 mm

Fig. 11
Fig.11The percentage of D estimations for all combinations of forecasting ranges (i.e., 3, 4, and 5 time steps) and radiuses of the Brumadinho dataset.The bin widths of forecasting ranges were selected at 50 mm

Fig. 12
Fig. 12 The spatial distribution of estimated EA (a) and D (b) values for the Brumadinho dataset with the forecast range of 4 and 45 m radius.The yellow polygon shows the collapsed area ◂

Fig. 14 ◂374
Fig.14The percentage of D estimations for all combinations of forecasting ranges and radiuses of the Cadia dataset.The bin widths of forecasting ranges were selected at 50 mm

Table 1
Characteristics of InSAR datasets This column indicates the point density over the surface of dams for Brumadinho and Cadia and open pit of Cobre las cruces *

Table 2
Hyperparameters and values used to tune the LSTM models Table 3 contains the RMSE and MAE results of the three datasets' best-tuned parameters of the LSTM model.It also includes the performance of the model in three forecasting ranges.The RMSE and MAE values of training and test of each forecasting range and dataset are reported in Table

Table 3
The best configuration of LSTM parameters achieved by hyperparameter tunning and evaluation accuracies for three forecasting ranges in the case studies

Table 4
The percentage of detected potential time series with the possibility of failure over the surface of dams and critical areas (the upper cells in each row of case studies) and over the entire datasets of case studies (the bottom cells in each row of case studies shaded in gray)