Landslide Characterization Applying Sentinel-1 Images and InSAR Technique: The Muyubao Landslide in the Three Gorges Reservoir Area, China

: Landslides are a common natural hazard that causes casualties and unprecedented economic losses every year, especially in vulnerable developing countries. Considering the high cost of in-situ monitoring equipment and the sparse coverage of monitoring points, the Sentinel-1 images and Interferometric Synthetic Aperture Radar (InSAR) technique were used to conduct landslide monitoring and analysis. The Muyubao landslide in the Three Gorges Reservoir area in China was taken as a case study. A total of 37 images from March 2016 to September 2017 were collected, and the displacement time series were extracted using the Stanford Method for Persistent Scatterer (StaMPS) small baselines subset method. The comparison to global positioning system monitoring results indicated that the InSAR processing of the Muyubao landslide was accurate and reliable. Combined with the ﬁeld investigation, the deformation evolution and its response to triggering factors were analyzed. During this monitoring period, the creeping process of the Muyubao landslide showed obvious spatiotemporal deformation di ﬀ erences. The changes in the reservoir water level were the trigger of the Muyubao landslide, and its deformation mainly occurred during the ﬂuctuation period and high-water level period of the reservoir.


Introduction
Landslides are one of the most common types of natural hazards that cause serious economic losses, casualties, and damages to buildings, critical infrastructures and industrial settlements [1]. The impact of landslides is particularly high in rural mountainous areas, where land management is difficult to achieve, and risk mitigation is scarce due to the large distances, difficult logistics and harsh climatic conditions. As an example, on 2 May 2014, a large-scale landslide occurred in the Badakhshan Province in Northeastern Afghanistan due to continuous heavy rainfall. It resulted in nearly 2700 deaths and more than 300 houses were buried [2]. Again, a giant landslide occurred in the Himachal Pradesh region of Northern India on 13 August 2017 that killed 47 people and nearly a 300-m highway was completely destroyed.
Monitoring and early warning systems are effective methods to reduce the risk of landslides [3], but the application of them on each single, potentially unstable slope is often impossible in mountain areas, due to the high spatial frequency of affected slopes. In such cases, remote sensing may help to reduce the cost and time for the application of mitigation measures, because it can provide a preliminary low-cost assessment of the severity of the slope instability and allow for the prioritization of critical cases. The use of remote-sensing data is one of the hotspots of landslide research at present. Among all the physical manifestations of a mass movement, the surface deformation is the most intuitive and comprehensive to measure and use for hazard assessment. It is a critical indicator for developing landslide early warning systems. Therefore, the implementation of deformation monitoring is of great significance for landslide prevention and mitigation [4,5]. The majority of fatal landslide events occur in less developed countries, such as China, India, Nepal, Bangladesh, etc. [1]. Traditional deformation monitoring equipment, such as an inclinometer, global positioning system (GPS), etc., are well-suitable, but their high costs limit their applications in underdeveloped areas.
Synthetic Aperture Radar Interferometry (InSAR) is an effective surface deformation monitoring method with a wide range and high accuracy. Commercial radar imaging can be expensive and mainly used for landslide hazard prevention in developed areas, including landslide identification [6,7], monitoring and early warning [8][9][10] and risk assessment [11]. With the free distribution of Sentinel-1 satellite images, InSAR has become a low-cost landslide monitoring technique. Small Baselines Subset (SBAS) InSAR is a time series analysis method of InSAR and can achieve effective monitoring in rural and urban areas. It has been gradually applied in landslide-prone areas, which provides new technical means for landslide monitoring and early warning [12][13][14][15][16].
The Three Gorges Reservoir area (TGRA) has been impounded since 2003. The periodical scheduling of the reservoir water level from 145 m to 175 m has been carried out since 2009. Affected by heavy rainfall and reservoir water level scheduling, a large number of new landslides and reactivation of old landslides have occurred [17]. Until now, more than 500 landslides have undergone varying degrees of deformation. However, due to the high cost of dedicated landslide monitoring, the latter was carried out only in 254 landslides in the TGRA [18], which cannot meet the needs of disaster prevention and mitigation.
The Muyubao landslide located in the TGRA was taken as a case study in this paper. Based on the Sentinel-1 radar images, the SBAS InSAR analysis method was applied to extract the information of landslide movement. The characterization of its temporal and spatial deformation was conducted as well. Moreover, the evolution mechanism of the Muyubao landslide was analyzed with the combined consideration of its engineering geological conditions and triggering factors, including rainfall and fluctuation of reservoir water level. This study aims to explore the feasibility and effectiveness of Sentinel-1 images and InSAR technique in landslide monitoring and analysis in the TGRA as a case example from an underdeveloped landslide-prone region.

Geological Conditions
The Muyubao landslide occurred in Shazhenxi Town, Zigui County of Hubei Province on the bank of the Yangtze River and 56 km away from the Three Gorges Dam. The Muyubao landslide is delimited laterally by N-S-oriented trenches, while the upper portion is straight and smooth. It is of a chair-like shape, with smooth topography in the middle and lower parts and steeper terrain in the upper part (     The Muyubao landslide developed within a flysch formation dipping in the same direction and roughly parallel to the slope. The stratigraphy is made by the typical succession of quartz sandstone and sand-mudstone layers of the Jurassic Xiangxi formation (Figure 3). The sliding body is mainly composed of two parts. The surface layer is made up by loose deposits composed of alluvial sub-clay, gravel layers with clay, colluvial and eluvial sub-clay and collapsed block stones. The lower part of the sliding body is composed of highly disrupted layered Quartz sandstone, which is relatively intact in the western upper portion and weathered and fractured in the upper and lower portions.
Remote Sens. 2020, 12, x FOR PEER REVIEW 4 of 19 The Muyubao landslide developed within a flysch formation dipping in the same direction and roughly parallel to the slope. The stratigraphy is made by the typical succession of quartz sandstone and sand-mudstone layers of the Jurassic Xiangxi formation (Figure 3). The sliding body is mainly composed of two parts. The surface layer is made up by loose deposits composed of alluvial sub-clay, gravel layers with clay, colluvial and eluvial sub-clay and collapsed block stones. The lower part of the sliding body is composed of highly disrupted layered Quartz sandstone, which is relatively intact in the western upper portion and weathered and fractured in the upper and lower portions. The upper portion of the sliding body has a linear detachment surface in the profile map, with an inclination of about 25° and a thickness of 60~90 m. The front part is an uplift platform formed by shear sliding, with the thickness of 80~120 m, and the inclination angle of the rock layer is about 27° (Figure 3). The thickness values here refer to the central part of the landslide, as illustrated in the profile of Figure 3, which may decrease towards the boundaries. The landslide is still active, as it has experienced several reactivation events in the recent past [19,20]. This makes it a priority case for monitoring, since, in case the movement would develop towards sudden failure, it will endanger the lives and property of 140 households (500 people) in the landslide area and threaten the safety of the road and shipping of the Yangtze River.

Field Investigation of Landslide Deformation
After the reservoir impounding to 175 m in the TGRA in September 2009, the deformation of the Muyubao landslide began to accelerate. Macro-deformations appeared in the landslide mass. In 2011, the road surface across both the boundaries of the landslide was damaged, and the roadbed on the west boundary was broken in widths of 20 cm (Figure 4a). The cracks in the wall of the village buildings on the eastern side of the landslide continued to be stretched (Figure 4b). In the investigation of 2014, a local collapse was found on the west side of the landslide (Figure 4c  The upper portion of the sliding body has a linear detachment surface in the profile map, with an inclination of about 25 • and a thickness of 60~90 m. The front part is an uplift platform formed by shear sliding, with the thickness of 80~120 m, and the inclination angle of the rock layer is about 27 • (Figure 3). The thickness values here refer to the central part of the landslide, as illustrated in the profile of Figure 3, which may decrease towards the boundaries. The landslide is still active, as it has experienced several reactivation events in the recent past [19,20]. This makes it a priority case for monitoring, since, in case the movement would develop towards sudden failure, it will endanger the lives and property of 140 households (500 people) in the landslide area and threaten the safety of the road and shipping of the Yangtze River.

Field Investigation of Landslide Deformation
After the reservoir impounding to 175 m in the TGRA in September 2009, the deformation of the Muyubao landslide began to accelerate. Macro-deformations appeared in the landslide mass. In 2011, the road surface across both the boundaries of the landslide was damaged, and the roadbed on the west boundary was broken in widths of 20 cm (Figure 4a). The cracks in the wall of the village buildings on the eastern side of the landslide continued to be stretched (Figure 4b). In the investigation of 2014, a local collapse was found on the west side of the landslide (Figure 4c

Small Baselines Subset InSAR Analysis
With the development of the InSAR technique and the increasing of synthetic aperture radar (SAR) images, Ferretti et al. [21] proposed the method called Persistent Scatterer InSAR (PSInSAR) based on the differential InSAR technique. The main idea of the PSInSAR is to use multiple SAR images covering the same area to analyze the stability of the amplitude and phase and, then, to identify the pixels that are less-affected by the spatiotemporal decorrelation. To obtain accurate deformation information, the components of the phase need to be jointly analyzed and modeled to remove errors. In the PSInSAR method, a unique master image is selected from all images to generate interferograms. The strategy will result in the presence of long-space baseline interference pairs. It produces a low density of the target reflectors in case of nonurban areas. Hence, in rural areas with few artificial buildings, the permanent scattering pixels may be very sparse [22].

Small Baselines Subset InSAR Analysis
With the development of the InSAR technique and the increasing of synthetic aperture radar (SAR) images, Ferretti et al. [21] proposed the method called Persistent Scatterer InSAR (PSInSAR) based on the differential InSAR technique. The main idea of the PSInSAR is to use multiple SAR images covering the same area to analyze the stability of the amplitude and phase and, then, to identify the pixels that are less-affected by the spatiotemporal decorrelation. To obtain accurate deformation information, the components of the phase need to be jointly analyzed and modeled to remove errors. In the PSInSAR method, a unique master image is selected from all images to generate interferograms. The strategy will result in the presence of long-space baseline interference pairs. It produces a low Remote Sens. 2020, 12, 3385 6 of 20 density of the target reflectors in case of nonurban areas. Hence, in rural areas with few artificial buildings, the permanent scattering pixels may be very sparse [22].
The small baseline subset method was initially proposed to overcome the problem of decorrelation by making full use of the interferograms with both small temporal baselines and short perpendicular baselines [23,24]. In rural areas, plenty of pixels have no dominated scattering characteristics in the SAR image. The pixels showing a slow decorrelating filtered phase (SDFP) are widely distributed in natural terrain and maintained good coherence during a short time interval. The Stanford Method for Persistent Scatterer (StaMPS) was proposed in 2004 [25]. In the StaMPS SABS method, these kinds of pixels are identified and analyzed to monitor surface displacement. The SDFP pixels are selected through their phase characteristics. In order to reduce the calculation burden, an initial subset of pixels containing almost all SDFP pixels are firstly selected through amplitude analysis. The amplitude dispersion value is calculated as the indicator of phase stability. The wrapped phase is composed of spatially correlated phase and spatially uncorrelated look angle error. The spatially correlated phase, including ground deformation, elevation error and orbit error, is estimated by the bandpass filtering in the frequency domain. The spatially uncorrelated look angle error mainly consists of the spatially uncorrelated elevation error. It is estimated through its correlation with the perpendicular baseline. The residual, obtained by removing these two terms from the wrapped phase, gives an estimation of the decorrelation noise γx (Equation (1)), which indicates the stability of the pixel phase: where N is the number of interferograms, ψ x,i is the wrapped phase,φ u θ,x,i is the spatially uncorrelated look angle error and ψ x,i is the spatially correlated term. The final SDFP pixels are selected through the threshold analysis of γx. Then, these SDFP pixels can be decomposed with the three-dimensional phase unwrapping method [26]. The unwrapped phase on a given pixel can be expressed as: where ψ x,i is the unwrapped phase, and φ D, x,i and φ N,x,i are the phase components due, respectively, to ground deformation, topographic error, atmospheric disturbance, inaccurate orbit information and noise. In the StaMPS SBAS analysis, a theoretical framework for three-dimensional phase unwrapping was proposed. These different phase components were estimated through iterative filtering with consideration of their characteristics in spatial-temporal domains. The deformation phase can be obtained by removing the other components from the unwrapped phase. Afterwards, a time series of deformation can be obtained from the phase using singular value decomposition.

Time Series InSAR Processing of Muyubao Landslide
Sentinel-1 is a two-satellite constellation with the prime objective of land. It is the first satellite developed by the European Commission and the European Space Agency (ESA) for the Copernicus Programme. This satellite was launched in April 2014. After half a year of trial operation, it began to be gradually used after October 2014. One goal of the mission is to provide C-Band SAR data continuity following the retirement of ERS-2 and the end of the Envisat mission. In this study, 37 images acquired by Sentinel-1a, in the period from 10 March 2016 to 13 September 2017 ( Figure 5), were collected from the ESA [27]. We collected the precise orbit parameters of the Sentinel-1 in the ESA as well [28]. The angle of the orbit with the north is 12.16 • , and the average incidence angle of three sub-swaths are 34  Interferograms are generated applying the two-track differential method [29]. The SRTM (Shuttle Radar Topography Mission) DEM (digital elevation model) covering the study area is utilized to generate single-look differential interferograms. Temporal and perpendicular thresholds are, respectively, set to 90 days and 1000 m for generating interferograms. Considering that these thresholds will generate many interferograms and result in a heavy computational burden, each image was mostly combined to three subsequent images for interferometric processing. Data processing was carried out in COMET-LiCSAR [30]. As a result, 97 interferograms were produced as the input for the SBAS analysis ( Figure 6). The amplitude difference dispersion was applied to selected coherent targets, and it was set as 0.6. Subsequently, the time series displacement was obtained by applying the StaMPS SBAS method. Moreover, in order to reduce the effects of residual atmospheric artefacts, the obtained time series displacement within the Muyubao landslide has been referenced to a stable neighbor point (31°1′52.62″N, 110°29′57.48″E) [15,31]. The direction of displacement acquired by the SBAS method is the line-of-sight (LOS) of radar satellite. Before applying the measurement to the landslide analysis, the LOS displacements should be projected along real slope components to get actual vectors. Hilley et al. [32,33] proposed a widely used projection method; the LOS displacement can be projected into the direction of the steepest descent slope. In this case study, a field investigation showed that the slip direction was 16°. Hence, the LOS displacement was projected onto the downward direction of along the sliding direction. Interferograms are generated applying the two-track differential method [29]. The SRTM (Shuttle Radar Topography Mission) DEM (digital elevation model) covering the study area is utilized to generate single-look differential interferograms. Temporal and perpendicular thresholds are, respectively, set to 90 days and 1000 m for generating interferograms. Considering that these thresholds will generate many interferograms and result in a heavy computational burden, each image was mostly combined to three subsequent images for interferometric processing. Data processing was carried out in COMET-LiCSAR [30]. As a result, 97 interferograms were produced as the input for the SBAS analysis ( Figure 6). The amplitude difference dispersion was applied to selected coherent targets, and it was set as 0.6. Subsequently, the time series displacement was obtained by applying the StaMPS SBAS method. Moreover, in order to reduce the effects of residual atmospheric artefacts, the obtained time series displacement within the Muyubao landslide has been referenced to a stable neighbor point (31 • 1 52.62"N, 110 • 29 57.48"E) [15,31]. Interferograms are generated applying the two-track differential method [29]. The SRTM (Shuttle Radar Topography Mission) DEM (digital elevation model) covering the study area is utilized to generate single-look differential interferograms. Temporal and perpendicular thresholds are, respectively, set to 90 days and 1000 m for generating interferograms. Considering that these thresholds will generate many interferograms and result in a heavy computational burden, each image was mostly combined to three subsequent images for interferometric processing. Data processing was carried out in COMET-LiCSAR [30]. As a result, 97 interferograms were produced as the input for the SBAS analysis ( Figure 6). The amplitude difference dispersion was applied to selected coherent targets, and it was set as 0.6. Subsequently, the time series displacement was obtained by applying the StaMPS SBAS method. Moreover, in order to reduce the effects of residual atmospheric artefacts, the obtained time series displacement within the Muyubao landslide has been referenced to a stable neighbor point (31°1′52.62″N, 110°29′57.48″E) [15,31]. The direction of displacement acquired by the SBAS method is the line-of-sight (LOS) of radar satellite. Before applying the measurement to the landslide analysis, the LOS displacements should be projected along real slope components to get actual vectors. Hilley et al. [32,33] proposed a widely used projection method; the LOS displacement can be projected into the direction of the steepest descent slope. In this case study, a field investigation showed that the slip direction was 16°. Hence, the LOS displacement was projected onto the downward direction of along the sliding direction. The direction of displacement acquired by the SBAS method is the line-of-sight (LOS) of radar satellite. Before applying the measurement to the landslide analysis, the LOS displacements should be projected along real slope components to get actual vectors. Hilley et al. [32,33] proposed a widely used projection method; the LOS displacement can be projected into the direction of the steepest descent slope. In this case study, a field investigation showed that the slip direction was 16 • . Hence, the LOS displacement was projected onto the downward direction of along the sliding direction.

InSAR Results of Muyubao Landslide
As stated in Section 3.2, the displacement time series of the Muyubao landslide was extracted by applying the time series InSAR technique and projected onto the direction of the downward sliding direction. Figure 7 shows the average deformation speed of the Muyubao landslide. The strong deformation mainly occurred at the eastern upper portion of the landslide, with deformation rates up to 100 mm/year. In the direction perpendicular to the sliding direction, the deformation of the eastern side was slightly stronger than the western side of the Muyubao landslide.

InSAR Results of Muyubao Landslide
As stated in Section 3.2, the displacement time series of the Muyubao landslide was extracted by applying the time series InSAR technique and projected onto the direction of the downward sliding direction. Figure 7 shows the average deformation speed of the Muyubao landslide. The strong deformation mainly occurred at the eastern upper portion of the landslide, with deformation rates up to 100 mm/year. In the direction perpendicular to the sliding direction, the deformation of the eastern side was slightly stronger than the western side of the Muyubao landslide. Some points in different parts of the Muyubao landslide were selected to further analyze the spatial deformation characteristics. The locations and cumulative displacement time series curves are shown in Figures 7 and 8, respectively. It can be seen that the cumulative displacements of PSI-04 and PSI-05 show the active part of the landslide (Figure 7). The cumulative displacement of these two points from March 2016 to September 2017 exceeded 250 mm (Figure 8d,e). The deformations of the monitoring points in the middle part were smaller, with the cumulative displacement about 150 mm (Figure 8f,h). The deformation rates were the smallest (about 50 mm/y) at the western front part of the landslide (Figures 7 and 8i-k). Compared to the linear monitoring curves at the upper portion of the landslide, the monitoring curves near the lower portion showed slight step-like characteristics (PSI-01, PSI-02, PSI-06, PSI-09 and PSI-10), which were directly affected by periodic rainfall and reservoir scheduling [34,35]. In summary, the deformation of the Muyubao landslide seems to show obvious differences in space. The deformation on the eastern upper portion of the landslide was the strongest, followed by the middle part, and the western part showed the smallest deformation. Some points in different parts of the Muyubao landslide were selected to further analyze the spatial deformation characteristics. The locations and cumulative displacement time series curves are shown in Figures 7 and 8, respectively. It can be seen that the cumulative displacements of PSI-04 and PSI-05 show the active part of the landslide (Figure 7). The cumulative displacement of these two points from March 2016 to September 2017 exceeded 250 mm (Figure 8d,e). The deformations of the monitoring points in the middle part were smaller, with the cumulative displacement about 150 mm (Figure 8f,h). The deformation rates were the smallest (about 50 mm/y) at the western front part of the landslide (Figures 7 and 8i-k). Compared to the linear monitoring curves at the upper portion of the landslide, the monitoring curves near the lower portion showed slight step-like characteristics (PSI-01, PSI-02, PSI-06, PSI-09 and PSI-10), which were directly affected by periodic rainfall and reservoir scheduling [34,35]. In summary, the deformation of the Muyubao landslide seems to show obvious Remote Sens. 2020, 12, 3385 9 of 20 differences in space. The deformation on the eastern upper portion of the landslide was the strongest, followed by the middle part, and the western part showed the smallest deformation.

The Deformation Analysis of Muyubao Landslide
The Three Gorges Reservoir oscillates regularly between 145 m and 175 m, with annual cycles (Figure 9). It can be seen from Figure 9 that, from January to February 2017, the TGRA remained at a normal water storage level of 175 m. Under the uplift pressure of the reservoir water, the Muyubao landslide was in a state of continuous deformation, the southeast edge of which showed a larger displacement rate (Figures 10f and 11f). From March to April 2017, with the slow decline of the reservoir level, the landslide gradually stabilized as a whole, and the displacement rates of the monitoring points were less than 10 mm/month (Figures 10g and 11g). From May to June 2017, the TGRA was experiencing its rainy season, while the reservoir level began to decline rapidly. The displacement rates of all monitoring points of the landslide increased gradually, some of which reached 30 mm/month (Figure 11h). Due to the hysteresis characteristics of rainfall infiltration, the displacement rates of the upper potion were greater than the middle part, which grew gradually with the increase of rainfall infiltration capacity and infiltration depth (Figure 10h). The landslide showed slow deformation when the reservoir water kept at the level of 145 m from July to August 2017. The infiltration lag effect of the previous rainfall, combined with the large continuous rainfall in the reservoir area, triggered the strong deformation period of the year 2017. Some of the monitoring displacement rates reached 40 mm/month (Figure 11c,i).

The Deformation Analysis of Muyubao Landslide
The Three Gorges Reservoir oscillates regularly between 145 m and 175 m, with annual cycles (Figure 9). It can be seen from Figure 9 that, from January to February 2017, the TGRA remained at a normal water storage level of 175 m. Under the uplift pressure of the reservoir water, the Muyubao landslide was in a state of continuous deformation, the southeast edge of which showed a larger displacement rate (Figures 10f and 11f). From March to April 2017, with the slow decline of the reservoir level, the landslide gradually stabilized as a whole, and the displacement rates of the monitoring points were less than 10 mm/month (Figures 10g and 11g). From May to June 2017, the TGRA was experiencing its rainy season, while the reservoir level began to decline rapidly. The displacement rates of all monitoring points of the landslide increased gradually, some of which reached 30 mm/month (Figure 11h). Due to the hysteresis characteristics of rainfall infiltration, the displacement rates of the upper potion were greater than the middle part, which grew gradually with the increase of rainfall infiltration capacity and infiltration depth (Figure 10h). The landslide showed slow deformation when the reservoir water kept at the level of 145 m from July to August 2017. The infiltration lag effect of the previous rainfall, combined with the large continuous rainfall in the reservoir area, triggered the strong deformation period of the year 2017. Some of the monitoring displacement rates reached 40 mm/month (Figure 11c,i).
The tensile cracks at the upper portion of the landslide provides a passage for the infiltration of continuous heavy rainfall. Therefore, the monitoring points of the eastern upper portion showed the largest cumulative displacements and displacement rates, which gradually decreased towards the front part (Figure 10i).
Remote Sens. 2020, 12, x FOR PEER REVIEW 10 of 19 The tensile cracks at the upper portion of the landslide provides a passage for the infiltration of continuous heavy rainfall. Therefore, the monitoring points of the eastern upper portion showed the largest cumulative displacements and displacement rates, which gradually decreased towards the front part (Figure 10i).

The Reliability Analysis of InSAR Monitoring in the Muyubao Landslide
Interferometry techniques have a huge potential to monitor landslides, as demonstrated by worldwide examples [12,15,32]. The monitoring accuracy of InSAR is affected by various factors, such as the atmosphere and terrain. Therefore, it is necessary to evaluate the reliability of the time series displacement before it is applied for landslide analysis. In the InSAR monitoring of the Muyubao landslide, the reliability analysis of the results can be carried out from two aspects.
(a) Comparison of InSAR monitoring results at different locations. The monitoring point PSI-12 outside the Muyubao landslide was selected for comparative analysis. According to the InSAR monitoring results (Figure 12), the monitoring location of PSI-12 did not deform during the monitoring period, which is consistent with the results of the field investigation. While the monitoring points (PSI-01~PSI-11) within the landslide were deforming during the same period ( Figure 8). By comparing the InSAR results of the monitoring points within and outside the landslide, it demonstrates that the InSAR monitoring in this study can identify the deformation and nondeformation zones.

The Reliability Analysis of InSAR Monitoring in the Muyubao Landslide
Interferometry techniques have a huge potential to monitor landslides, as demonstrated by worldwide examples [12,15,32]. The monitoring accuracy of InSAR is affected by various factors, such as the atmosphere and terrain. Therefore, it is necessary to evaluate the reliability of the time series displacement before it is applied for landslide analysis. In the InSAR monitoring of the Muyubao landslide, the reliability analysis of the results can be carried out from two aspects.
(a) Comparison of InSAR monitoring results at different locations. The monitoring point PSI-12 outside the Muyubao landslide was selected for comparative analysis. According to the InSAR monitoring results (Figure 12), the monitoring location of PSI-12 did not deform during the monitoring period, which is consistent with the results of the field investigation. While the monitoring points (PSI-01~PSI-11) within the landslide were deforming during the same period (Figure 8). By comparing the InSAR results of the monitoring points within and outside the landslide, it demonstrates that the InSAR monitoring in this study can identify the deformation and nondeformation zones. Remote Sens. 2020, 12, x FOR PEER REVIEW 13 of 19 Figure 12. The InSAR monitoring results of PSI-12.
(b) Comparison of monitoring results with different techniques. GPS is a reliable monitoring technique of ground deformation, which has been widely used in various fields [36,37]. In order to verify the accuracy of InSAR monitoring, the comparison was carried out between the monitoring results of GPS and PSI points at similar locations of the Muyubao landslide (the locations are shown in Figure 7). Considering that the automatic GPS monitoring equipment was installed in June 2016, the monitoring results from June 2016 to September 2017 were used for comparison. The results of GPS-01-PSI-08 and GPS-02-PSI-07 are shown in Figures 13 and 14, respectively. The cumulative displacement of InSAR and GPS at similar locations are highly consistent in magnitude. It indicates that the InSAR monitoring of this study is reliable and can be utilized for landslide analysis.   (b) Comparison of monitoring results with different techniques. GPS is a reliable monitoring technique of ground deformation, which has been widely used in various fields [36,37]. In order to verify the accuracy of InSAR monitoring, the comparison was carried out between the monitoring results of GPS and PSI points at similar locations of the Muyubao landslide (the locations are shown in Figure 7). Considering that the automatic GPS monitoring equipment was installed in June 2016, the monitoring results from June 2016 to September 2017 were used for comparison. The results of GPS-01-PSI-08 and GPS-02-PSI-07 are shown in Figures 13 and 14, respectively. The cumulative displacement of InSAR and GPS at similar locations are highly consistent in magnitude. It indicates that the InSAR monitoring of this study is reliable and can be utilized for landslide analysis. (b) Comparison of monitoring results with different techniques. GPS is a reliable monitoring technique of ground deformation, which has been widely used in various fields [36,37]. In order to verify the accuracy of InSAR monitoring, the comparison was carried out between the monitoring results of GPS and PSI points at similar locations of the Muyubao landslide (the locations are shown in Figure 7). Considering that the automatic GPS monitoring equipment was installed in June 2016, the monitoring results from June 2016 to September 2017 were used for comparison. The results of GPS-01-PSI-08 and GPS-02-PSI-07 are shown in Figures 13 and 14, respectively. The cumulative displacement of InSAR and GPS at similar locations are highly consistent in magnitude. It indicates that the InSAR monitoring of this study is reliable and can be utilized for landslide analysis.

The Formation Mechanism of the Muyubao Landslide
Based on the analysis of geological conditions and monitoring data (Figure 8), the deformation of the Muyubao landslide showed a large spatial difference. Large deformation mainly occurred in the middle and upper portions, while the deformation of the uplift at the front part was smaller (Figures 7 and 10). The profile with the largest deformation is shown in Figure 15. It can be seen that the deformation of the entire section is relatively small in the period of July-September 2016. It became larger in the period of January-March 2017. At the front edge of the profile near the Yangtze River, the Muyubao landslide deformed slightly in both periods. A field investigation was conducted for further analysis of the landslide deformation mechanism.
By analyzing the results of the field investigation, it can be seen that the deformation and failure mode is controlled by the lithology condition and the structure of the slope. Under the long-term gravitational process, the slope was creeping down along the weak mudstone layer. The front stress accumulated gradually, which resulted in the bending and uplifting in the middle and front stratum of the slope. Once the potential slip surface of the bending part connects, the accumulated energy may suddenly release, which may cause collapse and high-speed sliding. In conclusion, the evolution process can be divided into four stages: slipping stage, bending and uplifting stage, local shear stage and completely connecting stage of the shearing surface (Figure 16a) [19].
The bending and uplifting in the front stratum limited the deformation of the front part, which occurred when the slope was creeping down along the lines of weakness. Due to the topography and geomorphology conditions, the uplift in the front increases the sliding resistance effect of the landslide. When rainfall occurs, the gullies distributed on the surface of the landslide provides favorable channels for the rapid discharge of surface water, thus reducing the threat of water seepage to the landslide stability. Meanwhile, in the front part, the sliding surface is of a gentle slope, and the rock mass is bent and dips inside, which forms a sliding resistance for the landslide as well. In short, under natural conditions, the geological characteristics are beneficial to the overall stability of the landslide. However, due to the combined influence of factors in the TGRA, the Muyubao landslide has been undergoing a continuous creep deformation. gravitational process, the slope was creeping down along the weak mudstone layer. The front stress accumulated gradually, which resulted in the bending and uplifting in the middle and front stratum of the slope. Once the potential slip surface of the bending part connects, the accumulated energy may suddenly release, which may cause collapse and high-speed sliding. In conclusion, the evolution process can be divided into four stages: slipping stage, bending and uplifting stage, local shear stage and completely connecting stage of the shearing surface (Figure 16a) [19].  The bending and uplifting in the front stratum limited the deformation of the front part, which occurred when the slope was creeping down along the lines of weakness. Due to the topography and geomorphology conditions, the uplift in the front increases the sliding resistance effect of the landslide. When rainfall occurs, the gullies distributed on the surface of the landslide provides

The Relationship between Landslide and Influencing Factors
For the landslide in the TGRA, most of their movements were influenced by the fluctuation of the reservoir water level [38]. Seen from the InSAR results of the Muyubao landslide ( Figure 17, March 2016-September 2017), the deformation mainly occurred during the fluctuation period and high-water level period of the reservoir (>170 m). From March to June 2016, when the reservoir water level dropped from 167 m to 145 m, the landslide deformed by about 70 mm. The relationship between the stability of the Muyubao landslide and the reservoir water level during this monitoring period was similar to that of many other reservoir landslides [15,39,40]. Such as the Lorenzo-1 and Rules Viaduct landslides of the Rules Reservoir in Southern Spain [15], the InSAR data shows that the three acceleration periods of both landslides are related to drawdown periods of the water level. Many similar cases are also studied in the TGRA [17,41]. When the reservoir level drops rapidly, the pore-water pressure within the landslide begins to dissipate, and the dissipation speed lags greatly behind the reservoir drawdown speed. This process will induce hydrodynamic pressures. High hydrodynamic pressure is no more balanced by water lateral confining in the basin, which increases the sliding force of the landslide, thus reducing its stability. From November 2016 to February 2017, when the reservoir was at a water level of higher than 170 m, the Muyubao landslide deformed by 25 mm. Therefore, it was still deforming during the period of high water level. The high water level increased the height of the groundwater level inside the reservoir landslide, which changed the saturation state of the soil and reduced the shear strength. At the same time, the uplift force formed by the reservoir water in the submerged sliding mass reduced the sliding resistance force of the landslide, thereby reducing the stability of the landslide. We also analyzed the relationship between the rainfall and acceleration of movement within the Muyubao landslide. It showed a deceleration in movement related to the rainfall peaks. Due to the large thickness and scale of the Muyubao landslide, it is difficult for rainfall to seep into the slip surface to promote landslide movement. The decorrelation between reservoir landslides and rainfall are also shown in the Lorenzo-1 and Rules Viaduct landslides and other landslides in the TGRA [42]. Based on the long-term GPS monitoring data of the Muyubao landslide, Deng et al. [19] and Wan et al. [20] analyzed the relationship between the landslide deformation and its inducing factors (reservoir water level and rainfall). These papers showed how water level variation and rainfall have influences on landslide displacement patterns, which is consistent with the conclusion of this paper. Moreover, it indicates that the application of InSAR data is reliable for analyzing the deformation mechanism of landslides. Viaduct landslides and other landslides in the TGRA [42]. Based on the long-term GPS monitoring data of the Muyubao landslide, Deng et al. [19] and Wan et al. [20] analyzed the relationship between the landslide deformation and its inducing factors (reservoir water level and rainfall). These papers showed how water level variation and rainfall have influences on landslide displacement patterns, which is consistent with the conclusion of this paper. Moreover, it indicates that the application of InSAR data is reliable for analyzing the deformation mechanism of landslides.

The Future Development of InSAR in Landslide Application
The InSAR technique is suitable for observing slow-moving landslides at a large scale. There are always spatial differences in the deformation of large-scale landslides [43,44]. The landslide deformation characteristics (direction, intensity, etc.) in different positions and periods are not the same. The displacement obtained by the InSAR technique was a one-dimensional LOS direction (East-West). Nowadays, the LOS displacement is often projected into the steep slope direction to

The Future Development of InSAR in Landslide Application
The InSAR technique is suitable for observing slow-moving landslides at a large scale. There are always spatial differences in the deformation of large-scale landslides [43,44]. The landslide deformation characteristics (direction, intensity, etc.) in different positions and periods are not the same. The displacement obtained by the InSAR technique was a one-dimensional LOS direction (East-West). Nowadays, the LOS displacement is often projected into the steep slope direction to make it consistent with the main sliding direction, which is convenient for deformation analyses. However, this method cannot meet the demands of accurate landslide monitoring. Therefore, it is necessary to develop an extraction method for landslide three-dimensional deformation based on the ascending and descending orbit images and the landslide deformation evolution.
In the current early warning system, the monitoring methods of surface deformation are mostly the global navigation satellite system (GPS, Beidou, etc.). The remote-sensing techniques, such as InSAR, will be widely applied in landslide monitoring and early warning in the future. The monitoring results of these techniques are not point data but a sort of areal deformation metrics that can provide much more deformation information. Therefore, we need to develop data mining algorithms to obtain more information from the massive motoring data of these novel techniques. It will significantly improve the effectiveness of the landside early warning system.

Conclusions
In this study, Sentinel-1 images and the StaMPS SBAS method were utilized to extract the time series displacement of the Muyubao landslide. During this monitoring period, the changes in the reservoir water level were the main triggering factor of the Muyubao landslide, and its deformation mainly occurred during the fluctuation period and high-water level period of the reservoir. This landslide also showed a significant spatial difference, which was the strongest on the eastern upper portion of the landslide, followed by the middle part, and the western front portion showed the slightest deformation.
The Muyubao landslide experienced a translational mode during its evolution process. The natural geological characteristics (formation mechanism, slope structure, topography and geomorphology) are beneficial to the overall stability of the landslide under natural conditions. However, when the external environment changed, the combined effects of reservoir water and rainfall contributed to the deformation of the Muyubao landslide.
The monitoring and characterization of the Muyubao landslide was carried out by applying Sentinel-1 images and the time series InSAR technique. Through field investigation and verification with GPS, we found that this method can effectively monitor slow-moving landslides and presents the advantage of high-density coverage of monitoring points. The InSAR monitoring results can provide us more information for comprehensive understanding of the target landslides. In short, the application of Sentinel-1 images and InSAR techniques in carrying out landslide monitoring is an effective and economical method.
Author Contributions: C.Z. processed the dataset, performed the analysis and wrote the original draft of this paper. Y.C. contributed to processing the dataset and performed part of the analysis. K.Y. and Y.W. supervised the interpretation of the results from the geomorphological point of view. F.C. supervised the analysis. X.S., F.C. and B.A. edited the manuscript. All authors contributed to paper writing and revision. All authors have read and agreed to the published version of the manuscript.