Monitoring the Recent Activity of Landslides in the Mailuu-Suu Valley (Kyrgyzstan) Using Radar and Optical Remote Sensing Techniques

: Central Asian mountain regions are prone to multiple types of natural hazards, often causing damage due to the impact of mass movements. In spring 2017, Kyrgyzstan su ﬀ ered signiﬁcant losses from a massive landslide activation event, during which also two of the largest deep-seated mass movements of the former mining area of Mailuu-Suu—the Koytash and Tektonik landslides—were reactivated. This study consists of the use of optical and radar satellite data to highlight deformation zones and identify displacements prior to the collapse of Koytash and to the more superﬁcial deformation on Tektonik. Especially for the ﬁrst one, the comparison of Digital Elevation Models of 2011 and 2017 (respectively, satellite and unmanned aerial vehicle (UAV) imagery-based) highlights areas of depletion and accumulation, in the scarp and near the toe, respectively. The Di ﬀ erential Synthetic Aperture Radar Interferometry analysis identiﬁed slow displacements during the months preceding the reactivation in April 2017, indicating the long-term sliding activity of Koytash and Tektonik. This was conﬁrmed by the computation of deformation time series, showing a positive velocity anomaly on the upper part of both landslides. Furthermore, the analysis of the Normalized Di ﬀ erence Vegetation Index revealed land cover changes associated with the sliding process between June 2016 and October 2017. In addition, in situ data from a local meteorological station highlighted the important contribution of precipitation as a trigger of the collapse. The multidirectional approach used in this study demonstrated the e ﬃ ciency of applying multiple remote sensing techniques, combined with a meteorological analysis, to identify triggering factors and monitor the activity of landslides.


Introduction
Between 2005 and 2014, non-event related landslides were responsible for 4506 deaths in the 40 most mountainous countries of the world, excluding tens of thousands of deaths linked to landslides triggered by earthquakes [1]. Such natural disasters can often lead to catastrophic consequences, especially in mountainous areas where intense rainfall, high regional seismicity, weak geological materials, and steep slopes act as major triggers or predisposition factors. Therefore, it is essential to monitor landslides and to analyze related mechanisms in order to reduce and prevent the negative socioeconomic impacts linked to these natural hazards. The improvement of spatial and temporal resolution of satellite imagery, along with the development of novel techniques, makes it currently possible to cover large  [27]). The red pin shows the location of the town of Maily-Say, the black and blue lines indicate the faults and the waterbodies respectively.

Study Area
Kyrgyzstan is located in the middle of Central Asia and crossed by the geohazard-prone mountain range of the Tien Shan, with peaks reaching altitudes of 7000 m. According to the Ministry of Emergency Situations of the Kyrgyz Republic, more than 7% of the total territory, representing an area of approximately 15,000 km 2 , could potentially be affected by landslides. The most susceptible regions to slope failure are located along the northern and eastern rim of the Fergana Basin, which is also marked by the highest population density [6]. In this area, mass movements are concentrated in  [27]). The red pin shows the location of the town of Maily-Say, the black and blue lines indicate the faults and the waterbodies respectively.
In this study, we employ multiple distinct and complementary methods to understand the recent collapse of the Koytash landslide and the deformation affecting the neighboring Tektonik mass movement (see view shown in Figure 2). We detect terrain elevation changes related to both fast and slow displacements of the ground surface by using UAV-based imagery combined with radar and optical remote sensing techniques. In addition, we conducted a meteorological analysis to identify the triggering conditions, which led to slope instabilities.  [27]). The red pin shows the location of the town of Maily-Say, the black and blue lines indicate the faults and the waterbodies respectively.

Study Area
Kyrgyzstan is located in the middle of Central Asia and crossed by the geohazard-prone mountain range of the Tien Shan, with peaks reaching altitudes of 7000 m. According to the Ministry of Emergency Situations of the Kyrgyz Republic, more than 7% of the total territory, representing an

Study Area
Kyrgyzstan is located in the middle of Central Asia and crossed by the geohazard-prone mountain range of the Tien Shan, with peaks reaching altitudes of 7000 m. According to the Ministry of Emergency Situations of the Kyrgyz Republic, more than 7% of the total territory, representing an area of approximately 15,000 km 2 , could potentially be affected by landslides. The most susceptible regions to slope failure are located along the northern and eastern rim of the Fergana Basin, which is also marked by the highest population density [6]. In this area, mass movements are concentrated in a range of altitudes between 700 and 2000 m. These landslides, often characterized by deep and steep scarps, mobilize weakly consolidated sediments of Tertiary or Quaternary age, including loess deposits [7,19,28]. The Mailuu-Suu Valley is a former uranium mining area in the Jalal-Abad region, southern Kyrgyzstan at the north-eastern border of the Fergana Basin. This region is particularly prone to landslide hazards and, during the last 50 years, has experienced severe landslide disasters in the vicinity of numerous nuclear waste tailings. The high intensity of natural hazards in that area is due, on one hand, to its location between the mountain range and the basin marked by a relatively wet climate and, on the other hand, to the presence of soft geological materials near the surface [23] (see the geological map of Maily-Say in Figure A1). This type of sediment is prone to gravitational mass movements initiated by intense precipitations or strong seismic ground shaking [6,23,[29][30][31]. Due to the presence of radioactive wastes produced by past mining activities, this region is exposed to a high risk of environmental contamination and demands particular attention ( Figure 3). As a result, landslides represent a major threat to the local population of the small town of Maily-Say and of the regions downstream, not only through the direct impact, but also through indirect effects that lead to the creation of landslide dams that may block the river and cause the formation of a lake. Such lakes represent a further threat for downstream populations due to the potential failure of the landslide dam, causing the release of huge amounts of water resulting in floods. In the case of Maily-Say, the lake could also possibly reach the level of the nuclear waste tailings and thus cause their instability and the related mobilization of radioactive material. In 2003, Vandenhove et al. [30] identified 23 mine tailings and 13 nuclear landfills for a total volume of approximately 3 million m 3 (see tailing locations in the geological map in Figure A1). Due to its critical situation, Mailuu-Suu Valley was and still is the target area of several international risk assessment projects [4].
During the 1990s, three of the largest landslides of the Mailuu-Suu Valley (i.e., Koytash, Tektonik, and Isolith; the first two are shown in Figure 2), displaced more than 5 million m 3 of material [32]. These three landslides are located near the core of the central anticline and represent the greatest threat for the destabilization of some nuclear waste tailings [6,33]. According to Havenith et al. [6], the Mailuu-Suu landslides now extend over an area of approximately 10 km 2 . Havenith et al. [5] also showed that distant earthquakes (epicentral distance > 100 km), even of a magnitude less than 6, can induce slight displacements on landslide slopes in the Mailuu-Suu Valley.
The Koytash landslide (see location in Figure 3) has been extensively studied for many years as it has been active since the 50's [5,23,29,34,35]. Two elements indicate that this landslide is characterized by a roto-translational movement: First, the presence of a steeply dipping scarp at the crest of the landslide, showing slickensides, indicates sliding rather than a simple detachment as would be expected for a translational landslide. Second, the landslide body slides along a paleogenic limestone basement with a slope steeper than the ground surface. Although some translational sliding occurs above the limestone along the internal paleogenic clay layers, the base of the landslide is more than 50 m deep, which requires a rotational component. Furthermore, some anti-dip slope scarps in the lower part of the landslide also suggest some rotational movements. These observations are based on Electrical Resistivity Tomography (ERT) measurements completed by one of the authors (H.B. Havenith) together with the GEOPRIBOR team (led by I. Torgoev, Kyrgyzstan) [36]. Havenith et al. [5] evaluated its average annual speed at 60 cm/year and its volume at 5 million m 3 . In 2009, Schlögel [37] estimated the average speed of the movement at 1-2 mm/day and recognized Koytash as a major threat due to the possibility of a river dam. This roto-translational landslide was massively reactivated on the 22nd of April 2017 blocking the Mailuu-Suu River. Due to the damming of the valley, a lake formed with a depth of~15 m, almost reaching the base of the nearest tailing. A few weeks later, after the dredging of the dam to open a channel for the river, the dam partially failed and the lake progressively disappeared [38]. The rupture in spring 2017 also destroyed a few houses at the foothill of Koytash; these houses had already been evacuated due to the landslide risks [38]. Second, in addition to Koytash, we studied the large Tektonik landslide (perspective view in Figure 2 and map location in Figure 3). Tektonik is a complex (multi-rotational) and flow-like landslide with a length of more than one kilometer in an E-W direction. This landslide, initiated on 4 July 1992, displaced 2.5 million m 3 of material and caused significant environmental damages [6]. The main escarpment is located at an altitude of around 1350 m, while its lower part, divided into two lobes, joins the Mailuu-Suu River at an altitude of around 1000 m. A secondary scarp is located just below the main escarpment. Its activity is controlled by many factors, such as the loess cover in the upper part and the natural groundwater conditions that make it unstable [39]. The main risk for the surrounding population is related to its high river damming potential. This landslide first collapsed massively a few weeks after the earthquake of 15 May 1992, which struck the region with a magnitude of 6.2 [5,26,40,41]. The epicenter of this earthquake was located only 30 km SSE from the town of Maily-Say. Although the aforementioned earthquake certainly contributed to the development of the Tektonik landslide in combination with the rise of the groundwater table after spring precipitations, one must consider that loss of stability over longer terms was caused by the presence of underground mining galleries [5,42]. In 1994In , 2002In , and 2005, the massive movements of the Tektonik landslide repeatedly lead to the blocking of the Mailuu-Suu River, forming a dam and causing numerous upstream and later also downstream floods [5]. In 2017, massive failure occurred on Koytash landslide, while the Tektonik landslide had only slightly moved in its upper part. The materials and methods used to show the respective behaviors of the two mass movements are detailed in the following section.

Comparison of Multitemporal Digital Elevation Models
Multitemporal DEMs were used to study the evolution of the topography before and after the reactivation of the Koytash landslide (partly also for the Tektonik landslide). We collected all existing Second, in addition to Koytash, we studied the large Tektonik landslide (perspective view in Figure 2 and map location in Figure 3). Tektonik is a complex (multi-rotational) and flow-like landslide with a length of more than one kilometer in an E-W direction. This landslide, initiated on 4 July 1992, displaced 2.5 million m 3 of material and caused significant environmental damages [6]. The main escarpment is located at an altitude of around 1350 m, while its lower part, divided into two lobes, joins the Mailuu-Suu River at an altitude of around 1000 m. A secondary scarp is located just below the main escarpment. Its activity is controlled by many factors, such as the loess cover in the upper part and the natural groundwater conditions that make it unstable [39]. The main risk for the surrounding population is related to its high river damming potential. This landslide first collapsed massively a few weeks after the earthquake of 15 May 1992, which struck the region with a magnitude of 6.2 [5,26,40,41]. The epicenter of this earthquake was located only 30 km SSE from the town of Maily-Say. Although the aforementioned earthquake certainly contributed to the development of the Tektonik landslide in combination with the rise of the groundwater table after spring precipitations, one must consider that loss of stability over longer terms was caused by the presence of underground mining galleries [5,42]. In 1994In , 2002In , and 2005, the massive movements of the Tektonik landslide repeatedly lead to the blocking of the Mailuu-Suu River, forming a dam and causing numerous upstream and later also downstream floods [5]. In 2017, massive failure occurred on Koytash landslide, while the Tektonik landslide had only slightly moved in its upper part. The materials and methods used to show the respective behaviors of the two mass movements are detailed in the following section.

Comparison of Multitemporal Digital Elevation Models
Multitemporal DEMs were used to study the evolution of the topography before and after the reactivation of the Koytash landslide (partly also for the Tektonik landslide). We collected all existing DEMs for this region (i.e., Shuttle Radar Topography Mission (SRTM), SPOT, ASTER, and TanDEM-X DEMs), dating from 2000 to 2011, to characterize the topography prior to the collapse (Table 1). Due to the lack of current (post-failure) DEMs for this region, we had to build a new model based on unmanned aerial vehicle (UAV) images. These images were acquired on the 15th of August 2017, during a field campaign, shortly after the failure of Koytash. In total, 434 high-resolution images were processed using the Agisoft Photoscan ® software. After removing all possible artefacts (e.g., excess vegetation), this software triangulates a set of control points (i.e., corresponding features between a minimum of three images) to determine the depth of a point in space relative to the image's 2D plane. Altogether, the control points create a cloud of points that represent the skeleton on which the reconstruction of the DEM is based [43]. Before comparing the DEMs, they were georeferenced in the WGS84 UTM coordinate system, interpolated according to the same grid, and resampled to have an identical spatial extension corresponding to the area surveyed with the UAV. To georeference the DEMs, we collected, during a field campaign in August 2017, 24 Differential Global Positioning System (DGPS) in situ points (with an accuracy of~5 cm, see location of points in Figure 4). Out of the 24 measured DGPS points, we chose 4 stable positions (outside the failure area) that served as ground control points (GCPs) for the adjustment of the Z-axis. The DGPS altitudes measured for these 4 stable GCPs were compared with the corresponding altitudes for these same GCPs extracted from each DEMs. The relative mean differences in altitude were calculated between the DEMs and the GCPs (see Table A1). For each pixel, the mean differences obtained were added or subtracted from the initial altitudes. These steps were carried out with the ArcMap (10.5.1) version of ArcGIS ® (developed by ESRI). In order to detect mass displacements due to the collapse in April 2017, we subtracted the TanDEM-X elevation model (pre-landslide) from the UAV DEM (post-landslide). The TanDEM-X and UAV DEMs have vertical accuracies of 4 m and 0.62 m, respectively. This technique allowed us to highlight an uplift at the foot of Koytash (accumulation zone) and a subsidence at the crest of the landslide (depletion zone). from the initial altitudes. These steps were carried out with the ArcMap (10.5.1) version of ArcGIS ® (developed by ESRI). In order to detect mass displacements due to the collapse in April 2017, we subtracted the TanDEM-X elevation model (pre-landslide) from the UAV DEM (post-landslide). The TanDEM-X and UAV DEMs have vertical accuracies of 4 m and 0.62 m, respectively. This technique allowed us to highlight an uplift at the foot of Koytash (accumulation zone) and a subsidence at the crest of the landslide (depletion zone).

Radar Remote Sensing
The basic principles of InSAR rely on SAR systems sending radar signals to the Earth that are then reflected by the ground surface to the satellite. InSAR refers to several methods based on phase measurements [44,45]. The phase, measured for the return signal to the satellite, depends on the distance between the specific target and the satellite. The launch in 2014 of the first Earth Observation Satellite from the Copernicus program of the European Space Agency (ESA), Sentinel-1, brought new open access possibilities to monitor natural processes such as landslides. In this study, D-InSAR enables detecting displacements from the Copernicus Sentinel-1 sensors to the scattering target along the LOS of the satellite ( Figure 5) [46]. . Geometry of a SAR image and principle of D-InSAR during two successive passes (t1 and t2) of the satellite over the same area, in this case, respectively, before and after a landslide activation (or a reactivation as for Koytash) and generation of an interferogram; ∆R = difference of distance travelled by the phase φ; λ = wavelength of the sensor (modified from the Geoscience Australia [47]).
This method is based on the difference of phase recorded when a displacement occurs between two SAR acquisitions. The repeat cycle of the SAR satellite, in this case Sentinel-1, measures a phase value for each distinct pixel of the target area at two different times (e.g., every 6 days with the two operational Sentinel-1A and Sentinel-1B satellites). The phase difference, also called interferometric phase (∆φint), is expressed by where φ(t1) and φ(t2) are the SAR phase values at t1 and t2; ∆φdef is the phase difference related to the ground displacement in the LOS direction; ∆φorb is the orbital geometry variations of the satellite between two acquisitions due to the Earth curvature; ∆φtopo is the topographic component; ∆φatm is the atmospheric contribution (related to a potential modification of the signal propagation in the atmosphere); and ∆φnoise is the background noise resulting from changes in soil diffusion properties, changes in thermal properties of the atmosphere, and possible errors of co-registration [9]. As the atmospheric component acts at a regional scale, the impact of ∆φatm can be ignored when focusing on local displacements. In the case of D-InSAR, the topographic component is removed by using a preexisting DEM to leave only the phase difference associated to the displacement [48]. Although, a DEM sampled at SAR resolution would be necessary to completely remove the topography, lower resolution global DEMs, available in open access, are suitable when studying ground deformations (e.g., SRTM DEM [49][50][51][52][53][54][55]). The orbital contribution can be removed as we selected SAR images acquired along a same orbit, which we assume to be correctly georeferenced. The processing chain we used for the generation of interferograms is synthetized as a flowchart in the Appendix A, Figure A1. The open access Sentinel Application Platform (SNAP) software acquired via the Copernicus Open Access Hub provided by ESA, was used to process Sentinel-1 Single Look Complex (SLC) images. After being imported, the data was co-registered so that the Figure 5. Geometry of a SAR image and principle of D-InSAR during two successive passes (t 1 and t 2 ) of the satellite over the same area, in this case, respectively, before and after a landslide activation (or a reactivation as for Koytash) and generation of an interferogram; ∆R = difference of distance travelled by the phase ϕ; λ = wavelength of the sensor (modified from the Geoscience Australia [47]).
Measuring displacements with InSAR can be done using either descending or ascending geometry. When studying a roto-translational landslide such as Koytash, the movements are mainly perpendicular to the surface at the top and the foot of the landslide. As Sentinel-1 is a right-looking satellite, we chose the ascending geometry for which the LOS is almost parallel to the main rotational deformation and thus more or less perpendicular to the slope ( Figure 5). Thus, relatively to the landslide surface, a subsidence is a downward negative displacement, while an uplift is positive. However, velocities are projected according to the general downhill movement of the landslide.
This method is based on the difference of phase recorded when a displacement occurs between two SAR acquisitions. The repeat cycle of the SAR satellite, in this case Sentinel-1, measures a phase value for each distinct pixel of the target area at two different times (e.g., every 6 days with the two operational Sentinel-1A and Sentinel-1B satellites). The phase difference, also called interferometric phase (∆ϕ int ), is expressed by where ϕ(t 1 ) and ϕ(t 2 ) are the SAR phase values at t 1 and t 2 ; ∆ϕ def is the phase difference related to the ground displacement in the LOS direction; ∆ϕ orb is the orbital geometry variations of the satellite between two acquisitions due to the Earth curvature; ∆ϕ topo is the topographic component; ∆ϕ atm is the atmospheric contribution (related to a potential modification of the signal propagation in the atmosphere); and ∆ϕ noise is the background noise resulting from changes in soil diffusion properties, changes in thermal properties of the atmosphere, and possible errors of co-registration [9]. As the atmospheric component acts at a regional scale, the impact of ∆ϕ atm can be ignored when focusing on local displacements. In the case of D-InSAR, the topographic component is removed by using a pre-existing DEM to leave only the phase difference associated to the displacement [48]. Although, a DEM sampled at SAR resolution would be necessary to completely remove the topography, lower resolution global DEMs, available in open access, are suitable when studying ground deformations (e.g., SRTM DEM [49][50][51][52][53][54][55]). The orbital contribution can be removed as we selected SAR images acquired along a same orbit, which we assume to be correctly georeferenced.
The processing chain we used for the generation of interferograms is synthetized as a flowchart in the Appendix A, Figure A1. The open access Sentinel Application Platform (SNAP) software acquired via the Copernicus Open Access Hub provided by ESA, was used to process Sentinel-1 Single Look Complex (SLC) images. After being imported, the data was co-registered so that the master image, taken as a reference, and the slave images coincide. After this essential step, a preliminary interferogram was generated. Then, in order to preserve solely the phase component related to variations in the optical path, we subtracted the topographic contribution from the interferogram [56]. To do this, the DEM of the Shuttle Radar Topography Mission (SRTM, 30 m resolution), was subtracted from the preliminary product. During the last decades, many researchers used SRTM DEMs for their InSAR processing [49][50][51][52][53][54][55]. As the same topographic reference (DEM) is used for each pair of images, it does not affect the measurement of relative displacements of a few centimeters. Subsequently, the Goldstein Phase Filtering, Multilooking, and Range-Doppler terrain correction steps were applied to reduce noise, resize pixels in square format, and correct geometric distortions linked to the topography. The phase of the differential interferogram, initially expressed in radians, was then unwrapped to convert it into ground displacements. As a starting location for the unwrapping process, we chose a reference DGPS point (#23) situated in a stable and coherent area outside the boundaries of the landslides (see Figure 4a and Table A1). Finally, these displacements along the LOS (D LOS ) are calculated using the following equation, with λ = wavelength of the satellite (i.e., 5.55 cm for Sentinel-1 C-band radar). The dataset we used for this analysis consists of SAR satellite images captured between January 2016 and December 2017 over the Mailuu-Suu area ( Table 2). We used Sentinel-1A images taken in the Interferometric Wide swath mode (IW), whereas Sentinel-1B images were unavailable over the area. As the Koytash landslide is facing west, we used ascending geometry for the SAR images acquired in the Terrain Observation by Progressive Scans (TOPSAR) mode. In total, 50 SLC products were analyzed to estimate the displacement rate before, during, and after the collapse of Koytash. Table 2. Properties of the Sentinel-1A SAR images collected for the Mailuu-Suu area and used in this study.

Sentinel-1A Properties
Microwave band (wavelength) C-band ( [57]. FASTVEL was recently developed by TRE-Altamira and is available with restricted access on the ESA Geohazards Exploitation Platform (v2) [58]. This specific processing chain allows the creation of maps (in GeoTIFF format) of average ground displacement velocities during observation periods (in cm/year) based on the Persistent Scatterer Interferometry (PSI) technique [59]. The FASTVEL algorithm employs ENVironment SATellite (ENVISAT), Advanced SAR (ASAR), European Remote Sensing (ERS), as well as Sentinel-1 products to identify surface movements from a set of SAR images analyzing the stack of differential interferograms to extract terrain displacement information (e.g., extent and magnitude). Using the MTA mode, a map of the corrected topography (DEM error and reference DEM in meters) in GeoTIFF format (not used in this study) is produced. The main information of PSI products, in the LOS, in which each line of the database will represent one measurement point are delivered in a CSV file [60].

Optical Remote Sensing
We also studied the Koytash landslide using optical remote sensing by calculating the difference in vegetation indices (i.e., NDVI), allowing for the semiautomatic detection of (re)activated zones. We performed this complementary analysis because it enables to overcome certain limitations of D-InSAR, essentially linked to the significant loss of coherence for rapid movements [61,62]. To calculate the NDVIs, we used Pleaide-1A images coming from a single sensor which currently have the highest spatial resolution (Figure 6a). Therefore, this certifies that the bands used for image processing present the same wavelength and spectral response. In order to obtain the best acquisitions to perform the NDVI analysis, we pre-selected all Pleaide-1A images with less than 10% of cloud cover and chose the pre-and post-images closest to the date of the landslide failure. Accordingly, we compared an image from the 6th of June 2016 (pre-) with one taken on the 1st of October 2017 (post-landslide), both acquired in the descending mode and presenting a spatial resolution of 0.5 m (Figure 6b,c). We co-registered the raw images before comparing them to one another [23,25]. In addition, the topographic component and possible distortions were corrected (orthorectification). These steps are part of the geometric corrections that were necessary to perform a change detection analysis, obtain better results, and facilitate the interpretation of the Earth surface process [19,25,63].
NDVI is an indicator of vegetation health and biomass that is based on leaf absorption and reflection. During photosynthesis, live leaves mainly absorb incoming solar radiation within the red spectral band while reflecting near-infrared light. Within each spectral band, reflectances are ratios of the scattered over the incoming radiation. NDVI values are calculated using the following equation, where ρ NIR and ρ R are the spectral reflectance measurements, respectively, in the near-infrared band (NIR) and in the red band (R). Theoretically, NDVI values range between −1 and +1 [20,61]. Negative values correspond to areas with water, ice, or snow, which have higher ρ R than ρ NIR . Barren surfaces (rock or soil) have similar R and NIR reflectances resulting in NDVIs near 0. Positive values refer to areas with vegetation. Accordingly, low positive values (~0.2-0.5) represent shrubs or grassland while higher values approaching +1 (~0.6-0.9) represent luxuriant vegetation likely indicating slope stability [64,65]. This produced monochromatic images (grayscale). To improve the readability of the output image, we decided to modify the color scale (blue-green). Once the NDVIs were obtained for the entire master and slave images, we calculated the NDVI difference using the following equation, Computed in ArcMap, this method consists in subtracting, pixel-by-pixel, the NDVI values of the pre-landslide image from the post-landslide image. The result is a monochromatic image, this time with values ranging from −2 to +2. All the steps described hereinabove, except the latter, were achieved using the ENVI software (version 5.4.1) developed by Harris Geospatial. and chose the pre-and post-images closest to the date of the landslide failure. Accordingly, we compared an image from the 6th of June 2016 (pre-) with one taken on the 1st of October 2017 (postlandslide), both acquired in the descending mode and presenting a spatial resolution of 0.5 m (Figure  6b,c). We co-registered the raw images before comparing them to one another [23,25]. In addition, the topographic component and possible distortions were corrected (orthorectification). These steps are part of the geometric corrections that were necessary to perform a change detection analysis, obtain better results, and facilitate the interpretation of the Earth surface process [19,25,63].

Meteorological Data
In order to have a deeper understanding of the environmental conditions in which Koytash was massively reactivated during spring 2017, we performed a rudimentary meteorological analysis. The only in situ data available came from the Ak-Terek weather station, located about 30 km east of the Koytash landslide at an altitude of 1772 m. Despite the potential bias linked to the distance from our study site, we preferred to use in situ data (including rainfall, snow depth, and temperature), as previous studies highlighted the important underestimation of satellite rainfall data [66]. Since 1934, temperature, rain, and snow precipitation data have been recorded daily by the station. To identify and highlight potential triggers to the reactivation of the Koytash landslide, we examined the meteorological records between the 1st of January 2011 and the 31st of December 2017. Unfortunately, in situ rainfall data was missing between August 27th, 2014 and December 1st, 2015 as well as from January to June 2016. This lack of records might be explained either by intense snowfall precipitations or by an error due to a deficiency of the station. After collecting and compiling the significant weather-related data, we produced relevant graphs showing both the daily and monthly evolution of meteorological factors prior and after the collapse of Koytash in April 2017; they are shown at the end of the next section.

Differences between Pre-and Post-Failure Topography
In order to observe the effects of rapid movements on the topography of the Koytash landslide during its reactivation, we compared all previously mentioned DEMs. The smaller deformations prior to failure will be analyzed via D-InSAR processing in the next sub-section. Figure 7a compiles all DEMs, after the application of corrective factors (see Section 3.1), along the A-B profile of the main axis of the Koytash landslide. According to the DGPS in situ data, the TanDEM-X from 2011 represents the best topographic fit with a calculated relative mean difference of 4.55 m in altitude (see Table A1). Thus, this TanDEM-X was compared with the post-failure UAV DEM for this analysis. all DEMs, after the application of corrective factors (see Section 3.1), along the A-B profile of the main axis of the Koytash landslide. According to the DGPS in situ data, the TanDEM-X from 2011 represents the best topographic fit with a calculated relative mean difference of 4.55 m in altitude (see Table A1). Thus, this TanDEM-X was compared with the post-failure UAV DEM for this analysis. As shown in Figure 7b, the difference of DEMs revealed the displacement of sediments during the reactivation episode of April 2017. Several zones were identified: a depletion zone, located in the upper part of the Koytash landslide and indicating the removal of rocky material along the scarp, and an accumulation zone, situated near the toe of the landslide body where more than 30 m of rocks were uplifted. The two aforementioned zones were then mapped according to the degree of surface level change (i.e., minor and major, see Figure 7c). The hatched area on Figure 7c represents an artifact that could potentially be related to an excess of vegetation. The collapse along the scarp caused a downward displacement of the material resulting in an uplift and the accumulation of sediments at the foot of the landslide, which partially blocked the Mailuu-Suu River. This analysis also highlighted a minor depletion zone on the top and a small accumulation zone in the middle of the neighboring landslide, Tektonik. These results are in accordance with the roto-translational nature of the landslide movement as well as with the in situ observations described in Section 2.  As shown in Figure 7b, the difference of DEMs revealed the displacement of sediments during the reactivation episode of April 2017. Several zones were identified: a depletion zone, located in the upper part of the Koytash landslide and indicating the removal of rocky material along the scarp, and an accumulation zone, situated near the toe of the landslide body where more than 30 m of rocks were uplifted. The two aforementioned zones were then mapped according to the degree of surface level change (i.e., minor and major, see Figure 7c). The hatched area on Figure 7c represents an artifact that could potentially be related to an excess of vegetation. The collapse along the scarp caused a downward displacement of the material resulting in an uplift and the accumulation of sediments at the foot of the landslide, which partially blocked the Mailuu-Suu River. This analysis also highlighted a minor depletion zone on the top and a small accumulation zone in the middle of the neighboring landslide, Tektonik. These results are in accordance with the roto-translational nature of the landslide movement as well as with the in situ observations described in Section 2.

Landslide Displacements from D-InSAR Analysis
On the basis of the 50 available Sentinel-1A radar images, we created 46 differential interferograms. Three periods (pre-, intra-, and post-reactivation), characterized by an interval of 12 days, were considered: 2-14 August 2016, 18-30 March, 2017 and 9-21 August 2017. The three unwrapped differential interferograms are presented in Figure 8 to estimate displacements at pre-, intra-, and post-sliding periods. The phase unwrapping measures values initially expressed in radians [0, 2π] that were converted in D LOS (cm) measurements. In order to obtain consistent and comparable results between the three considered periods, the process of phase unwrapping was carried out according to the four known stable DGPS points (11,13,23, and 24; see Figure 4).

Landslide Displacements from D-InSAR Analysis
On the basis of the 50 available Sentinel-1A radar images, we created 46 differential interferograms. Three periods (pre-, intra-, and post-reactivation), characterized by an interval of 12 days, were considered: 2-14 August 2016, 18-30 March, 2017 and 9-21 August 2017. The three unwrapped differential interferograms are presented in Figure 8 to estimate displacements at pre-, intra-, and post-sliding periods. The phase unwrapping measures values initially expressed in radians [0, 2π] that were converted in DLOS (cm) measurements. In order to obtain consistent and comparable results between the three considered periods, the process of phase unwrapping was carried out according to the four known stable DGPS points (11,13,23, and 24; see Figure 4).  Figure 8a). From the perspective of the landslide surface, a subsidence is negative while an uplift is positive. The Koytash landslide presents, in its upper part, an oval-shape deformation field with a diameter (along the A-B profile) of approximately 400 m. This deformation, extending from west to east, is characterized by negative displacements. Moving away from the central point of this anomaly, for which the DLOS reached 1.1 cm of displacement away from the sensor over 12 days, the DLOS values tend towards 0 in the middle part of the landslide body. The displacements are negative inside the landslide boundaries and therefore represent a predominantly vertical movement of subsidence close to the scarp area.
The second unwrapped phase result (Figure 8b, t1) is based on SAR acquisitions taken one month before the recent reactivation. For Koytash, the conditions are similar to the differential interferogram of August 2016 with the presence of an ellipsoidal anomaly upstream of the landslide while the  Figure 8a). From the perspective of the landslide surface, a subsidence is negative while an uplift is positive. The Koytash landslide presents, in its upper part, an oval-shape deformation field with a diameter (along the A-B profile) of approximately 400 m. This deformation, extending from west to east, is characterized by negative displacements. Moving away from the central point of this anomaly, for which the D LOS reached 1.1 cm of displacement away from the sensor over 12 days, the D LOS values tend towards 0 in the middle part of the landslide body. The displacements are negative inside the landslide boundaries and therefore represent a predominantly vertical movement of subsidence close to the scarp area.
The second unwrapped phase result (Figure 8b, t 1 ) is based on SAR acquisitions taken one month before the recent reactivation. For Koytash, the conditions are similar to the differential interferogram of August 2016 with the presence of an ellipsoidal anomaly upstream of the landslide while the deformation field moved slightly upwards, following the shape of the main escarpment, with a slight increase in the displacement amplitude. At t 1 (March 2017), the maximum D LOS measured at this location reached −1.5 cm, while the foothills of the landslide are characterized by low but positive values, nearing 0. Previously negative, these values demonstrate the transition from a subsidence to an uplift at the toe of the landslide.
The third interferogram (Figure 8c, t 2 ), created from the SAR acquisitions of August 2017, differs from the previous ones. Indeed, more than two months after the reactivation of the Koytash landslide, the movement gradually became stationary. The deformation pattern visible on the previous interferograms is no longer observed as the collapse has caused a sharp decrease in its activity. The third interferogram (Figure 8c, t2), created from the SAR acquisitions of August 2017, differs from the previous ones. Indeed, more than two months after the reactivation of the Koytash landslide, the movement gradually became stationary. The deformation pattern visible on the previous interferograms is no longer observed as the collapse has caused a sharp decrease in its activity. The obtained DLOS values are between −0.17 cm in the depletion area near the scarp and 0.15 cm in the foothills.
VLOS maps showing the entire spectrum of measured velocities as well as velocities with values strictly greater than 1.5 cm/year are presented in Figure 9a  The DLOS values of the three differential interferograms (Figure 8a-c) as well as the VLOS values from FASTVEL (Figure 9a) were extracted at 39 points, every 35 m, along three A-B profiles separated by 40 m as presented in Figure 10. In Figure 10a, the A-B topography profiles show the altitude of the landslide slope before (TanDEM-X, 2011) and after the failure (UAV, 2017). These two DEMs help distinguish, from left to right along the A-B profile, the uplift (accumulation zone), the transition zone, and the scarp (depletion zone) of the landslide. The role of this changing topographic profile is to facilitate the planar interpretation of the corresponding VLOS and DLOS values in Figure 10b,c,  Figure 10. In Figure 10a, the A-B topography profiles show the altitude of the landslide slope before (TanDEM-X, 2011) and after the failure (UAV, 2017). These two DEMs help distinguish, from left to right along the A-B profile, the uplift (accumulation zone), the transition zone, and the scarp (depletion zone) of the landslide. The role of this changing topographic profile is to facilitate the planar interpretation of the corresponding V LOS and D LOS values in Figure 10b,c, respectively. As a reminder, the DEMs used to show the topography of the landslide in Figure 10a do not correspond to the SRTM DEM used to remove the topographic contribution from the interferometric phase (see Section 3.2). When interpreting these InSAR results, it is important to consider that we measured slow and nonlinear displacements changing over time. Although D LOS (cm) and V LOS (cm/year) are covering different periods of time and expressed using different units, they are both calculated based on Sentinel-1A images and thus have the same resolution. In Figure 10b At the toe of the landslide, the velocities recorded before the accumulation zone, are negative and do not exceed −0.5 cm/year, indicating that the area can be considered as stationary (with < velocities < . Correspondingly, the DLOS measured for this zone, at the three periods, are included between −0.2 and 0.2 cm. In the accumulation zone, the VLOS values increase linearly although the displacements are close to 0 cm/year for all three periods. Overall, as we move along the profile, from the accumulation zone to the depletion zone, we observe higher velocities and an increase in negative displacements at t0 and t1. Due to the roto-translational nature of Koytash, there are some slow translational movements in the transition zone at t0 and t1, combined to light downward displacements linked to the continuous creeping of the slope. In the depletion zone, where topographic differences are the highest, the VLOS values reach maxima, up to 4 cm/year, with a negative displacement peak of −1.5 cm in 12 days for t1. Additionally, the slopes of the DLOS points at each period of time (t0, t1, and t2) reveal the different acceleration intensities. Although all displacements are negative, we notice that the slope steepness deduced by the dots is greater for the t1 than for the t0 interferogram (see Figure 10c). Finally, after the collapse of Koytash (t2), we notice that the entire slope has more or less stabilized with displacement under 0.5 cm in both directions. At the toe of the landslide, the velocities recorded before the accumulation zone, are negative and do not exceed −0.5 cm/year, indicating that the area can be considered as stationary (with −λ 4π < velocities < λ 4π ). Correspondingly, the D LOS measured for this zone, at the three periods, are included between −0.2 and 0.2 cm. In the accumulation zone, the V LOS values increase linearly although the displacements are close to 0 cm/year for all three periods. Overall, as we move along the profile, from the accumulation zone to the depletion zone, we observe higher velocities and an increase in negative displacements at t 0 and t 1 . Due to the roto-translational nature of Koytash, there are some slow translational movements in the transition zone at t 0 and t 1 , combined to light downward displacements linked to the continuous creeping of the slope. In the depletion zone, where topographic differences are the highest, the V LOS values reach maxima, up to 4 cm/year, with a negative displacement peak of −1.5 cm in 12 days for t 1 . Additionally, the slopes of the D LOS points at each period of time (t 0 , t 1 , and t 2 ) reveal the different acceleration intensities. Although all displacements are negative, we notice that the slope steepness deduced by the dots is greater for the t 1 than for the t 0 interferogram (see Figure 10c). Finally, after the collapse of Koytash (t 2 ), we notice that the entire slope has more or less stabilized with displacement under 0.5 cm in both directions.

Mapping of Reactivated Zones and NDVI Analysis
The multitemporal series based on high-resolution optical satellite images is shown in Figure 11. This analysis helps observing the evolution of the landslides since the last inventory, based on the Quickbird image acquired in 2007, carried out by Schlögel et al. [23]. Figure 11 Figure 12 illustrates NDVI values pre-(a) and post-landslide reactivation (b) as well as the NDVI difference between the two periods (c). In Figure 12c, the bright zones represent areas which recorded a vegetation loss between 2016 and 2017, while the darker zones indicate a gain of vegetation or no visible change. The bright areas in Figure 12c may result from the reactivation process that occurred between the two acquisitions, causing the collapse of several existing landslides. The areas reactivated in 2017, mapped by comparing the 2016 and 2017 images are shown in red in Figure 12 to demonstrate the equivalence between these landslide limits and the results obtained semi-automatically by calculating the difference in NDVI. This analysis identified the collapse zones of both Koytash and Tektonik. The lighter areas outside the landslide limits are probably linked to high seasonal variability in the vegetation cover between the end of spring (in June) and the end of summer (in October) when cloud free images were available.
April 2017 (Figure 11e); and SPOT6 of 31 April 2017 (Figure 11f). Taking into consideration the seasonal effect, the deformations zones, visible on this sequence of images, seem to have remained generally stable until 2015, showing several zones of old-growth vegetation. In Figure 11e, we observe that on 9 April 2017, the Tektonik landslide had already started to slide, while Koytash collapsed on the 22nd of April 2017, blocking the Mailuu-Suu River. Although less spectacular than Koytash, the images of 9 April 2017 and 31 May 2017 reveal the failure of the upper part of Tektonik.  Figure 12c may result from the reactivation process that occurred between the two acquisitions, causing the collapse of several existing landslides. The areas reactivated in 2017, mapped by comparing the 2016 and 2017 images are shown in red in Figure 12 to demonstrate the equivalence between these landslide limits and the results obtained semiautomatically by calculating the difference in NDVI. This analysis identified the collapse zones of both Koytash and Tektonik. The lighter areas outside the landslide limits are probably linked to high seasonal variability in the vegetation cover between the end of spring (in June) and the end of summer (in October) when cloud free images were available.

Identification of the Meteorological Triggering Factors
A meteorological analysis was performed to identify the triggering factors contributing to the reactivation of Koytash. We calculated total cumulative precipitations from 2011 to 2016, as well as the total precipitations for 2017, to detect a potential variability in rainfall that might be at the origin of the collapse of the Koytash landslide. On average, since 2011, total annual precipitations reached 1460.5 mm, while only 877.4 mm where recorded in 2017. However, despite the fact that 2017 was drier and hotter than previous years, half of the annual rainfall recorded in 2017 fell between March

Identification of the Meteorological Triggering Factors
A meteorological analysis was performed to identify the triggering factors contributing to the reactivation of Koytash. We calculated total cumulative precipitations from 2011 to 2016, as well as the total precipitations for 2017, to detect a potential variability in rainfall that might be at the origin of the collapse of the Koytash landslide. On average, since 2011, total annual precipitations reached 1460.5 mm, while only 877.4 mm where recorded in 2017. However, despite the fact that 2017 was drier and hotter than previous years, half of the annual rainfall recorded in 2017 fell between March and April (i.e., 414.2 mm), just before the reactivation episode ( Figure 13). We also observed that the month of April 2017 recorded the fourth highest precipitation peak since January 2011. and April (i.e., 414.2 mm), just before the reactivation episode ( Figure 13). We also observed that the month of April 2017 recorded the fourth highest precipitation peak since January 2011. When focusing on the meteorological conditions in March and April 2017, we notice that the daily temperature gradually increases from March to April (Figure 14). Consequently, the snow cover decreased progressively during March and disappeared completely before mid-April, about ten days before the collapse. Additionally, precipitations increased starting from 27 March 2017. Three main peaks stand out, on April 3, 13, and 22, with daily precipitation greater than 40 mm. The action of water from the intense rain precipitations in April, especially the peak on 22 April 2017 (52 mm), as When focusing on the meteorological conditions in March and April 2017, we notice that the daily temperature gradually increases from March to April ( Figure 14). Consequently, the snow cover decreased progressively during March and disappeared completely before mid-April, about ten days before the collapse. Additionally, precipitations increased starting from 27 March 2017. Three main peaks stand out, on April 3, 13, and 22, with daily precipitation greater than 40 mm. The action of water from the intense rain precipitations in April, especially the peak on 22 April 2017 (52 mm), as well as the rapid snow melting, represent important triggering factors leading to the reactivation of the Koytash landslide. Figure 13. Comparison of precipitation and temperature for the year 2017 (dark gray) versus precipitation and temperature for a period from 2011 to 2016 (hatched) based on in situ data from the Ak-Terek station (located at a higher altitude than the Mailuu-Suu River slopes, at about 1800 m, where the temperature is about 4° C lower than in the Mailuu-Suu region).
When focusing on the meteorological conditions in March and April 2017, we notice that the daily temperature gradually increases from March to April ( Figure 14). Consequently, the snow cover decreased progressively during March and disappeared completely before mid-April, about ten days before the collapse. Additionally, precipitations increased starting from 27 March 2017. Three main peaks stand out, on April 3, 13, and 22, with daily precipitation greater than 40 mm. The action of water from the intense rain precipitations in April, especially the peak on 22 April 2017 (52 mm), as well as the rapid snow melting, represent important triggering factors leading to the reactivation of the Koytash landslide.

Discussion
In Central Asia, numerous studies conducted on mass movements have previously demonstrated the high landslide hazard level of Kyrgyzstan [6]. In 2015, as part of a large-scale geohazards analysis, Havenith et al. [27] performed an extensive landslide susceptibility mapping of the Tien Shan mountain range. Previous studies identified the Mailuu-Suu Valley as one of the regions particularly prone to slope instabilities [3]. Both studies also confirmed that the number of landslides is increasing in this region, with a predominance of the reactivation and enlargement of existing landslides. These mass movements are likely to grow or merge with neighboring landslides when reactivating.
In this study, we used different types of satellite data to investigate the recent reactivation of two of the largest landslides of the Mailuu-Suu Valley, Kyrgyzstan. The Koytash and Tektonik landslides represent a major threat for the local populations due to their recurrent past reactivations, the presence of numerous neighboring nuclear waste tailings, and the risk of damming the Mailuu-Suu River [34]. Although former remote sensing studies in the Mailuu-Suu Valley used satellite data to identify landslides through change detection methods [4,23], this study brings a new dimension by using modern techniques, such as InSAR analysis and UAV imagery, to measure displacement rates and detect deformation zones.
D-InSAR was used to detect slow displacement rates of Koytash and Tektonik one year before their reactivation. In 2015, Teshebaeva et al. [12] performed a D-InSAR analysis to detect slow-moving deep-seated landslides near the town of Uzgen, 90 km southeast of our study area. Similarly, they highlighted a strong correlation between deformations peaks and precipitations and revealed the continuous activity affecting the slope of the landslides in their area of interest. While they used an L-band Advanced Land Observing Satellite Phased Array type L-band SAR (ALOS PALSAR) data set, we processed the newly available C-band Sentinel-1 satellite images. Although, L-band radar sensors have a longer wavelength (i.e., λ = 23.8 cm), which partially penetrates surface vegetation and allows the measurement of significant D LOS over long periods, their temporal resolution is of 46 days [67]. Therefore, due to the low vegetation cover in our study area, we chose the open access C-band SAR images that have only 12 days of repeat cycle (over the Mailuu-Suu area), reducing potential temporal decorrelation in comparison to the former study. Our analysis showed that the upper part of Koytash was the most active with negative D LOS values showing an early stage of subsidence along the scarp. Furthermore, the differential interferograms suggest that the reactivations of Koytash and Tektonik were not sudden but resulted from a long-term deformation. This could partly be related to a geological phenomenon called creeping, characterized by a slow downward movement of the soil. It is interesting to note that, during the entire study period, the displacement magnitudes of Koytash were greater than those measured for Tektonik. This result justifies why only the upper part of Tektonik collapsed while Koytash was entirely reactivated. The calculation of V LOS also led to this conclusion. When looking at the lower part of Koytash, a transition to very small positive values was detected in March 2017. This may be linked to the fact that the displacements are oriented perpendicularly to the LOS of the satellite, making them difficultly detectable. A few months after the collapse, in August 2017, we noticed that the movement had slowed down, implying a gradual stabilization of the landslide. These results suggest that it took a maximum of four months for the landslide to stabilize after being active for at least nine months prior to the failure. Accordingly, it would be justified to perform a back-analysis of the last decades by using ERS data to produce supplementary interferograms and thus estimate when Koytash's movement was initiated.
One of the most constraining limitations encountered during our D-InSAR analysis was the 1D expression of D LOS inducing the systematic underestimation of intensities. Therefore, depending on the field properties, it is only possible to track positive and negative displacements along the LOS of the satellite while deformation in the orthogonal direction is not detectable. In ascending orbital geometry, the observed displacements and velocities are positive when approaching the satellite and negative when moving away from the satellite. Although the multitemporal analysis with the FASTVEL algorithm requires at least 25 images, it is only able to detect mean velocity patterns (or hotspot areas as defined in [59]) that are significantly changing over such a long temporal coverage. Moreover, the velocity rates and directions are evolving overtime in a complex and vegetated terrain meaning that only the relative spatial changes are comparable to detect hazardous areas. To fully apprehend the complex deformation behaviors and facilitate their interpretation, sophisticated multitemporal InSAR techniques with various sensors (multi-angle) [48], ideally combined with in situ measurements, enable to break down the displacements into 3D in order to measure the true direction of the ground movement. Fuhrmann & Garthwaite [68] demonstrated the advantages of combining LOS measurements from different viewing geometries. Resolving 3D surface displacements would provide an improved overview of the complex geomorphology and dynamic of landslides. It would also be interesting to exploit additional SAR processing chains such as the SBAS (Small BAseline Subset) [69,70] or PSI multi-time series [62,71] to go one step further in multitemporal InSAR analysis. Another frequently encountered error in D-InSAR analysis is related to phase unwrapping and often leads to the misinterpretation of results [72]. In the case of Koytash and Tektonik, several parameters, such as the topography of the region and significant variations in the deformation rates, made the phase unwrapping a complex procedure. The topographic component, including both the roughness and vegetation cover of the landslide, can create artifacts (e.g., layover or foreshortening) responsible for sudden variations in slopes. Thus, an incorrectly unwrapped interferometric phase induces inaccurate displacement rates due to possible phase jumps (2π) [11]. To overcome this limitation, in this study the unwrapping process was initiated from a reference point (DGPS) located in a stable and coherent area, outside the landslide boundaries. As demonstrated by Schlögel et al. [73] and Teshebaeva et al. [12], the use of L-band in combination to C-band SAR images could help obtain better results by reducing time decorrelation effects. L-band acquisitions may be more suitable for highly vegetated areas as well as larger displacements. Another limitation of D-InSAR is the loss of coherence when confronted to rapid displacements or changes in land cover (e.g., agriculture). Additionally, by increasing the frequency of SAR acquisitions, the temporal decorrelation decreases as well. In the case of Mailuu-Suu, the use of Sentinel-1A, with only 12 days of repeat cycle, enabled us to observe some precursors before the main reactivation. It is noteworthy that for certain well-covered areas, the combination of Sentinel-1A and Sentinel-1B allows shorter return periods of 6 days. However, the effectiveness of this technique is counterbalanced by its numerous limitations. The main sources of decorrelation are linked to meteorological conditions, spatial and temporal resolutions, magnitudes, velocities, and orientation of the movements, or even variations in soil properties.
Although radar analyses can only highlight slow displacements because the limit of detectability is a function of the wavelength of the sensor, abrupt topographic modifications as well as the evolution of the land cover can be underlined through change detection. This can be achieved using several techniques such as the difference in DEMs, NDVIs or with a multitemporal optical analysis. Variations in surface topography can be used to delimit the depletion and accumulation zones. Our results, based on TanDEM-X and UAV DEMs, revealed that the reactivation of Koytash was characterized by a combined rotational and translational movement. Indeed, soft rocks detached from the scarp (depletion zone) applied an intense, almost vertical, pressure on the material below which caused an uplift in the foot (accumulation zone). By contrast, due to the superficial deformation of the Tektonik landslide, we could assume that the displacement was more likely translational. The performance of this technique is based on the quality of the DEMs. Depending on their coordinate systems and their georeferencing, the comparison and the difference between the products will be more or less exploitable. In the case of the Mailuu-Suu Valley, DGPS points, collected in August 2017 during the field campaign, were used to adjust raw DEMs and to ensure the accuracy of the results. However, such corrections can sometimes lead to significant errors, which must be considered. Therefore, it is essential to test the overlapping compatibility of the input files. Moreover, even though DEMs created on the basis of drone acquisitions have highly superior resolution, they can only be collected on a local scale because they are extremely time-consuming. Nevertheless, high-resolution DEMs are widely used to estimate the volume of displaced geological material [18,74,75].
Studying changes in land cover using NDVI is frequently used to survey the evolution of landslides and has been extensively applied to landslide monitoring studies [19,20,23]. A NDVI analysis of the Mailuu-Suu region, based on five landslide inventories spanning the past 50 years (1962, 1984, 1996, 2002, and 2007), demonstrated the efficiency of using NDVI to detect new sliding activation [23]. In this study, we used vegetation indices to highlight deformation zones and verify if our findings coincided with the D-InSAR results. On the basis of the difference of NDVI map, we observed bright areas (i.e., soil denudation) and darker zones (i.e., growing vegetation). As expected, bright areas were identified within the landslides' boundaries, but also outside the collapsed zones. These might indicate the seasonal variability marked by a strong difference in the vegetation cover depending on the periods of acquisition of the two images (e.g., end of spring (June) and end of summer (October)). Normalization, pixel-by-pixel, of the images by the average NDVI value could possibly mitigate this seasonal effect. However, many parameters (e.g., humidity and anthropogenic changes in land use) must also be considered in order to distinguish the zones that are truly reactivated. Nonetheless, in the case of this study, NDVI was used as a complementary qualitative analysis to verify if the landslide boundaries coincided with important changes in vegetation.
The reactivation of Koytash and Tektonik was also observed in the multitemporal optical analysis showing the evolution of the landslides since 2007. Koytash's recent reactivation can be considered as the most important failure event since the first movements were detected in the 1960s [76]. This massive collapse led to the obstruction of the Mailuu-Suu River forming a dam and subsequently a lake. As the level of water rises, the lake becomes prone to contamination from an upstream uranium tailing and thus represents an important environmental risk for potential water pollution. In 2017, the level of the lake almost reached the toe of this nuclear waste tailing. The formation of this lake may also be a direct threat to the local populations due to related floods.
Considering the actual climatic context, this study highlighted the importance of understanding the major triggering role of meteorological factors affecting the (re)activation of landslides. To support the remote sensing analysis, and in order to better understand the influence of weather on natural hazards occurrence in Central Asia, we analyzed different potential triggers of the activation. The meteorological analysis performed in the frame of this study demonstrated the triggering role of rainfall combined with the rapid snow melt in the reactivation process as it had already been demonstrated for Kyrgyzstan [77] but also for other regions [78,79]. The significant concentration of rain as well as the sudden snowmelt contributed to the water saturation of the soil before Koytash's failure. These predominant weather settings may have conditioned the deformation just before the reactivation, by increasing the pore pressure and thus the mobility of the sediments. Both the Koytash and Tektonik landslides are located in the bed of the Mailuu-Suu valley (Figure 3a). Therefore, during rainfall and the melting of snow caps upstream, water flows through the valley following the Mailuu-Suu River. Furthermore, the landslides are located on steep slopes with overhanging mountains from where surface runoff also flows after the snow melts and intense rainfall. Due to the triggering role of precipitations, it would be justified to automate the daily calculation of meteorological statistics and develop an early warning system to alert the surrounding populations in case of an upcoming disaster. For example, in Africa, Monsieurs et al. [66] were able to determine a rainfall threshold over which landslide susceptibility increases considerably.
With a perpetual increasing number of landslides, studies at regional scale become more and more necessary. With the current technological advances in satellite data, we are now able to process larger amounts of data and span wider areas. However, local and detailed analyses remain crucial when considering landslides such as Koytash and Tektonik that are recurrently reactivated and represent an important threat to surrounding populations. Due to the great density of landslides and the related risks, the Mailuu-Suu region is and will remain a critical area requiring continuous monitoring.

Conclusions
The aim of this study was to investigate the recent behavior of two landslides-Koytash and Tektonik-in the Mailuu-Suu Valley, Kyrgyzstan. This study was divided into three main objectives: a comparison of multitemporal DEMs, a D-InSAR analysis, and an analysis by optical remote sensing. To better understand the geo-environmental conditions causing the reactivation of the Koytash and Tektonik landslides, we also performed a complementary meteorological analysis. We showed that intense precipitations as well as the rapid melting of an important snowcap contributed to the activation of these two and many other landslides in this region in 2017. However, by creating displacement and velocity maps, the InSAR analysis revealed signs of sliding activity long before the last rupture of Koytash and Tektonik. This indicates that these landslides result from long-term deformations partly related to creeping. Indeed, InSAR allows the measurement of small magnitude displacements that are often not visible in the field. Although this method is not able to detect rapid movements (e.g., the collapse of landslides), InSAR can be used to monitor large-scale areas and identify long-term precursory deformations. Another major limitation is that InSAR measures displacements along the LOS that cannot be interpreted in 3D. Nonetheless, due to the open access availability of SAR acquisitions with short repeat cycles (e.g., 6 days for Sentinel 1A & B), InSAR is a powerful tool to monitor areas prone to landslides and is thus playing a key role in risk management.
Unlike InSAR, NDVI is a change detection method capable of highlighting the rapid evolution of vegetation. The difference of NDVI, between 2016 and 2017, identified several areas of landslide related land cover change, including the failure of Koytash and Tektonik. Nevertheless, anthropogenic changes in land use or seasonal variations in vegetation cover, depending on the acquisition date, may lead to misinterpretations of the results. Therefore, as shown in this study, NDVI can be used as a complementary qualitative analysis for landslide monitoring, capable of confirming deformations zones identified through InSAR. Similarly to NDVI, the multitemporal optical analysis highlighted deformation zones. The use of high-resolution optical images allows the delimitation of landslide boundaries, enabling precise surface calculations, which is not always possible with other remote sensing techniques. By mapping the evolution of the landslides in the Mailuu-Suu Valley, this optical analysis revealed that 30% of total landslide surface was reactivated in 2017. Additionally to the previous methods, we created a post-landslide DEM using UAV imagery, only possible at a local scale. The comparison of the latter with a pre-landslide TanDEM-X DEM, allowed us to determine a maximal accumulation of 30.13 m and a maximal depletion of 41.4 m. This evolution of the topography demonstrated that the collapse of the upper part of Koytash pushed on the underlying layers, which transferred the pressure downslope, resulting in an uplift due to the push-up effect and sliding material aggradation. When studying recent reactivations, post-landslide DEMs are rarely available, and it is therefore necessary to create a new DEM. In conclusion, the combination of radar and optical data revealed the recent landslide activity while the meteorological analysis helped us to identify the triggering factors responsible for the collapse of numerous landslides in Kyrgyzstan in spring 2017. The inherent limitations of the remote sensing methods used in this work, despite their considerable advantages, show that it is useful and essential to combine each of the techniques in order to obtain complete and conclusive results. It is therefore crucial to pursue the study of mass movements, as well as the development of novel cutting-edge techniques, to better understand their triggering conditions and, eventually, anticipate their rupture.  Figure A1. Geological map projected on the TanDEM-X 12 m hillshade of the Mailuu-Suu target area, with the location of nuclear waste tailings (green) and landslides (red outlines), as mapped by Schlögel et al. [23].