Identify Landslide Precursors from Time Series InSAR Results

Meng Liu · Wentao Yang,5 · Yuting Yang · Lanlan Guo · Peijun Shi

Abstract Landslides cause huge human and economic losses globally.Detecting landslide precursors is crucial for disaster prevention.The small baseline subset interferometric synthetic-aperture radar (SBAS-InSAR) has been a popular method for detecting landslide precursors.However, non-monotonic displacements in SBAS-InSAR results are pervasive, making it challenging to single out true landslide signals.By exploiting time series displacements derived by SBAS-InSAR, we proposed a method to identify moving landslides.The method calculates two indices (global/local change index) to rank monotonicity of the time series from the derived displacements.Using two thresholds of the proposed indices, more than 96% of background noises in displacement results can be removed.We also found that landslides on the east and west slopes are easier to detect than other slope aspects for the Sentinel-1 images.By repressing background noises, this method can serve as a convenient tool to detect landslide precursors in mountainous areas.

Keywords Monotonously changing displacements · Moving landslides · SBAS-InSAR · Time series of deformation

1 Introduction

As a main type of geological hazards, landslides caused many casualties and damages to mountain communities in the world (Petley 2012; Froude and Petley 2018).Existing work shows that some of the catastrophic landslides have deforming precursors before their failures (Intrieri et al.2018; Fan et al.2019; Ouyang et al.2019; Qi et al.2021).Therefore, identifying deforming landslides is crucially important to recognize imminent dangers of landslides to avoid tragic disasters.

Field reconnaissance is a traditional way to map landslide hazards (Strom and Korup 2006).However, it is often labor intensive and time consuming and dangerous in mountains.In contrast, remote sensing is an ef ficient way to identify moving landslides in large mountain regions (Xu et al.2020).Pixel of fset tracking (POT) and interferometric synthetic aperture radar (InSAR) are two most popular remote sensing techniques to identify landslide precursors (Ventura et al.2011; Casu and Manconi 2016; Lacroix et al.2018).Although POT (with either optical or SAR images) has been widely used to monitor landslide deformation, it is dif ficult to detect subtle deformations of a few centimeters (Li et al.2020).

In contrast, InSAR technology can monitor landslide deformations of a few centimeters by calculating phase differences in SAR images (Dai et al.2020; Zhang et al.2022).There are a few InSAR methods regarding dif ferent pairing strategies, such as D-InSAR (dif ferential interferometric synthetic aperture radar) (Wang et al.2013), PSI (persistent scatterer interferometry) (Guéguen et al.2009), and SBASInSAR (small baseline subset interferometric synthetic aperture radar) (Berardino et al.2002; Bayer et al.2017).By performing dif ferential interferometry, D-InSAR uses two SAR images before and after the movement of a landslide to calculate the deformation between them, which is fast to process but is always inf luenced by atmospheric conditions and temporal decorrelations (Lu et al.2019).To overcome these problems in D-InSAR, multi-temporal InSAR (MTInSAR), such as PSI and SBAS-InSAR were developed by inversing time series deformations (Bayer et al.2017; Lu et al.2019; Roy et al.2022).

Fig.1 Location of the study area in northeastern Tibetan Plateau and elevation.GD = Guanggaishan-Dieshan.

The MT-InSAR technology has been extensively used to derive time series deformations of disastrous landslides in a retrospective way (Intrieri et al.2018; Ouyang et al.2019)or monitor deformation of well-known moving landslides(Bian et al.2022).In these studies, cumulative displacements of landslides are carefully selected by considering large displacements (usually > 100 mm) and changing monotonously with time.However, it is possible that there are many non-monotonously changing points that may also have large cumulative displacements in MT-InSAR results.In addition, there may be some monotonously changing displacements with minor displacements, which are landslides submerged by background noises.Therefore, identifying moving landslides is challenging from MT-InSAR derived time series results.By exploiting spatial patterns of deformations, space cluster strategies have been popular in deriving landslides from MT-InSAR results (Bianchini et al.2012; Lu et al.2019; Solari et al.2019; Zhang, Zhu,et al.2021).However, the space cluster methods ignored time series information in MT-InSAR results, which could make it challenging to detect small moving landslides of a few pixels.Few works used time series displacements to extract landslides (Urgilez Vinueza et al.2022).Based on the assumption that cumulative displacements of landslides should change monotonously with time, this work aims to propose an algorithm to extract true moving landslides by fully exploiting a time series of SBAS-InSAR derived displacements.The work was carried out in Gansu Province in western China.First, we derived a time series of cumulative deformation using the SBAS-InSAR.Then, we applied the proposed method to the derived time series deformations to f ilter out non-monotonous deformation areas.Finally, we compared our method with two traditional ones to showcase our improvements.

2 Methodology

This section has two subsections.In the f irst subsection,we introduce the study area with a topographic map and describe its geological conditions.In the second subsection,we describe the method used for this work.In particular, we detail the indices we proposed by explaining their def initions with examples.

2.1 Study Area

The study area is located on the northeastern edge of the Tibetan Plateau, in the middle and upper reaches of the Bailong River—a secondary tributary of the Yangtze River(Fig.1).With frequent geological hazards, the study area has been identif ied as one of the most landslide sensitive areas in China (Zhang et al.2018).The study area is characterized by steep slopes and deep-incised valleys with limited space for local development.It has a semiarid climate with an average annual rainfall of 300–400 mm.The multi-year mean summer (June–August) rainfall is > 70 mm per month,accounting for 70–80% of the total rainfall (Zhang, Meng et al.2021).The elevation of the study area ranges from 1023 to 4057 m above sea level.The topographic features have a great impact on the local climate.This is manifested in that areas with relatively high altitudes receive more rainfall, while valley areas receive less rainfall.

The study area is located on the southern section of the Qaidam-West Qinling block.Fold and fault structures are extremely developed in the area and neotectonic movement is active.The local active faults mainly extend along the NWW-SEE and NEE-SWW directions (Yu et al.2015).The area is cross-sectioned by the Guanggaishan-Dieshan (GD)thrust, which consists of the North GD fault and the South GD fault.The 2008 Wenchuan Earthquake had an ef fect on the study area by causing large cracks in sediments and activating many landslides (Zhang et al.2018).The lithology of the study area is composed of thick layers of sandstone,limestone, slate, conglomerate, and sandy conglomerate (Bai et al.2012).The wide distribution of layered soft and hard rocks is one of the most notable characteristics of the geological environment in this area.

2.2 Methods

The method of the work can be divided into two main parts(Fig.2).The f irst part is to derive cumulative displacements using the traditional SBAS-InSAR method by using SAR images from both the ascending and descending tracks.The second part is to use the proposed method to extract monotonous deformation pixels.

2.2.1 Inversing Deformation in the Line of Sight (LOS)Direction by SBAS-InSAR

SBAS-InSAR technology is based on setting the spatial and temporal baseline thresholds as a condition, and then generates interferometric pairs, thus mitigating decoherence phenomenon (Lanari et al.2007; Chen et al.2021).We used a total of 46 scenes of Sentinel-1 single look complex(SLC) images under the descending track and 59 scenes of the ascending track data from April 2020 to April 2022 to derive surface deformation using SBAS-InSAR,1https:// search.asf.alaska.edu/and the orbit correction was performed through the precise orbit f ile2https:// scihub.coper nicus.eu/corresponding to the time.By setting a temporal baseline threshold of 60 days and a spatial baseline threshold of 20%, 261 and 219 interferometric pairs were generated for the ascending and descending track images, respectively.Among them, the maximum spatial thresholds in the descending track and ascending track image pairs are 228 m and 218 m, respectively.

Adaptive f iltering functions were used for interferogram processing, and unwrapped the phase by using minimum cost f low (MCF) algorithm (Werner et al.2003; Pepe and Lanari 2006).After these processes, we used the SRTMDEM3https:// dwtkns.com/data for terrain correction (Zhang, Zhu, et al.2021)and corrected atmospheric errors through atmospheric spatiotemporal f iltering.We set ground control points (GCPs)based on the selection of a relatively stable region.To estimate and remove the remnant constant phase, these points were used for ref inement and ref lattening.We inversed the f irst deformation rate to f latten the resulting interferogram by selecting the linear model.Based on the result of the f irst deformation rate, we removed atmosphere phase delay by using temporal high-pass f ilter and spatial low-pass f ilter to separate the phase components.Finally, the singular value decomposition (SVD) method was used to obtain line of sight (LOS) deformation results from unwrapped phase(Berardino et al.2002; Chen et al.2021).The SBAS-InSAR operation steps of this work were implemented in the ENVI/SARscape package.1:4 multi-look operation was used for the range and azimuth direction, and the f inal output image resolution is 20 × 20 m.

2.2.2 Def initions of Global Change Index (GCI) and Local Change Index (LCI) and Their Applications to Extract Landslides

Based on the SBAS-InSAR derived cumulative displacements, we proposed two indices to quantify monotonicity of the displacement curve for each pixel.The f irst proposed index is the global change index (GCI).It can be calculated in the following steps using the deformation time series.For a given deformation time series, we start from the second deformation value to the last.We compare the second deformation value to the f irst value.If the second is smaller than the f irst value, we count 1.Then, we compare the third deformation value with its preceding ones and count the number of deformations that are larger than it.We repeat this process to the end of the deformation time series.Finally, we sum up these counted numbers in all above mentioned repeating processes to get the GCI for this pixel.The mathematical def inition of the GCI is:

The GCI def ines the overall monotonicity of time series.To further explain the definition of GCI, we take the descending track in this work as an example.There are 46 Sentinel-1 SAR images, that is, 46 displacement values in the SBAS-InSAR time series.Wheni=x, we count the number of displacements from the f irst (S 1) to the last ()that is larger thanS x.We then sum up all counts from the 2nd to the 46th to get the GCI for the pixel.For an ideal constantly decreasing time series of 46 measures, we have a maximum GCI of 1035 (1 + 2 + 3 + … + 45), which means every measurement on the time series of the SBAS-InSAR result is less than all previous measurements.On the contrary, the 0 value of GCI indicates that the deformation time series increases monotonically.

To consider local f luctuations in the time series, we proposed the local change index (LCI), which only compares any neighboring displacements in the time series.The LCI is the number of relations where the later measurement is larger than the earlier one in the SBAS-InSAR time series result.The mathematical def inition of the LCI is:

where,nis the number of total used SAR images for a given track.The LCI def ines local f luctuations of the cumulative displacements in the time series.In the displacement time series, we compare thexth displacement () and its last neighboring one (.We count the number of the relations where>Take the 46 Sentinel-1 SAR images in the descending track as an example, the ideally consistent decrease of the displacement will result in a maximum LCI of 45.If LCI equals 0, it means the f irst displacement is smaller than the neighboring second one and the displacement increases consistently.

In theory, the GCI ranges from 0 to ((n−1) ×n) / 2 and the LCI ranges from 0 ton–1 fornSAR images.The 0 value of both GCI and LCI means monotonously increasing displacement, and their theoretical maximums indicate monotonously decreasing displacements.The GCI measures an overall monotonicity of the displacement, whereas the LCI measures the local f luctuation of the displacement.

Figure 3 are three examples of the time series of displacements from the ascending track.There are 59 SAR images in the ascending track.In theory, the range of GCI and LCI are 0–1711 and 0–58.The low GCI and LCI in Fig.3 a show an overall increasing trend of the displacement.The GCI of Pixel 1 is far less than one half of the theoretical maximum, indicating a persistent increase of the displacements.The low LCI of Pixel 1 means that there is more increase than decrease for any two neighboring displacements in the time series.On the contrary, the high GCI and LCI of Pixel 3 (Fig.3 c) show an overall decreasing trend of the displacements.The moderate GCI and LCI of Pixel 2 (Fig.3 b) show a complex change of the time series.The GCI is larger than one half of the theoretical maximum, which indicates an overall down-trend of the displacement.The LCI equals one half of the theoretical maximum, that is, half increase and half decrease of values among all 58 neighboring displacements, which is the most f luctuated time series.

Fig.3 Time series of displacements with low, moderate, and high global change index (GCI) and local change index (LCI)

Fig.4 Histogram of the global change index (GCI) of the study area for the ascending ( a) and descending tracks ( b).Large GCI to the right end of the histograms indicate the monotonous decrease of the displacements with time, whereas small GCI near 0 means the time series of the displacements increases monotonously with time.Two percentile values of the histograms (3% and 97%) are shown with vertical dashed lines.

2.2.3 Evaluation of the Method

By using the last cumulative displacement image, we calculated the mean (μ) and standard deviation (σ) of the study area.We set two thresholds (μ ± σ and μ ± 2σ) to eliminate background noises.All pixels in the last cumulative displacement image that are not within these ranges will be removed.Finally, we compared the result f iltered by these thresholds with the result f iltered by our proposed two indices (GCI and LCI).We also used a typical landslide to illustrate the dif ferences among them.

3 Results

We f irst use statistical tools (histograms and scatterplots)to show the proposed indices (GCI and LCI), which lead to thresholds used later.Second, we present the dif ference between ground deformation maps with and without using our proposed method.Third, we selected some well-known landslides to validate our results.

3.1 Global and Local Change Indices (GCI and LCI)in SBAS-InSAR Results

We calculated the GCI and LCI using the SBAS-InSAR derived time series of displacements for both tracks.Figure 4 shows the histograms of the GCI for the ascending and descending tracks in the study area.There are 580,412 and 581,204 valid pixels with GCI for the ascending and descending tracks, respectively.Because the number of used SAR images for the two tracks are dif ferent, the range of their derived GCI values is also dif ferent.The GCI values for the descending track ranges from 20 to 1035 and for the ascending track ranges from 28 to 1705.The medians, modes, means, and standard deviations of the GCI values for the ascending and descending tracks are 1127 and 474, 1440 and 499, 1059 and 481, and 384.11 and 216.94, respectively.For the ascending track, the mean(1059) is far from the mode (1440), indicating a negative skewness distribution.For the descending track, the mean(481) is very close to the mode (499), indicating that the frequency of GCI values is close to a normal distribution.In this work, we used two percentile thresholds (3% and 97%) to extract true moving landslides in the study area.Tails of the distribution are where the time series of cumulative displacements changes monotonously and are more likely to be moving slopes.

Fig.5 Histogram of the local change index (LCI) for the descending( a) and ascending tracks ( b).Inset plots ( c), ( d), ( e), ( f) are tails of the histograms.The theoretical minimums of both tracks are 0, and the theoretical maximums of the two tracks are 58 and 45, respectively.The actual minimum LCIs of the two tracks are 12 and 11, whereas the actual maximum LCIs are 53 and 45, respectively.

Fig.6 Scatterplots to show the relationship between the global change index (GCI) and displacement for the ascending ( a) and descending tracks ( b).Colors of the points show the local change index (LCI) of each pixel.Vertical dashed lines are 3% and 97% percentiles of the GCI.

We also plotted the histogram of the LCI for the ascending (Fig.5 a) and descending tracks (Fig.5 b).Because the total number of used SAR images for the ascending and descending tracks are 59 and 46, the theoretical maximum LCIs are 58 and 45, respectively.The theoretical minimum LCIs are 0 for both tracks.For both tracks, the actual minimum values (12 and 11 for the ascending and descending tracks) are larger than the theoretical minimum(0), indicating that there are no monotonously increasing displacements with time.There are 12 or more episodes of decreases in the time series of displacements for all increasing cumulative displacements.In contrast, the actual maximum LCI value of the descending track equals its theoretical maximum (45), indicating that there are perfect monotonously decreasing time series of cumulative displacements.For the ascending track, the dif ference between the actual (53) and the theoretical maximums(58) of the LCI is smaller than the dif ference between the actual (12) and theoretical minimums (0), indicating that there are less f luctuations for the decreasing than increasing cumulative displacements for the ascending track.

Figure 6 shows the relationship between the GCI, LCI and displacement of all pixels in two scatterplots for both tracks.The GCI and displacements are negatively correlated, that is,positive displacements have smaller GCIs and negative displacements have larger GCIs for both tracks.The horizontal 0 mm displacement line and the vertical 3% and 97% GCI lines partitioned the panel into six zones.Points in zone I have small GCIs, positive displacements and some of them have relatively small LCIs (red color), indicating that these pixels moved monotonously towards the sensor in the LOS direction.Points in zone II have large GCIs and negative displacements.Purple points in zone II also have large LCIs,indicating that these pixels moved monotonously away from the sensor in the LOS direction.There are more high-value LCI (purple colored) points than low value (red colored)points.Most of the high LCI points are found in region II of both tracks with displacements < -50 mm, meaning that most of the monotonously decreasing time series have larger displacements.

Fig.7 Cumulative deformations in the line of sight (LOS) direction of the ascending track ( a) of the study area from April 2020 to April 2022 and the cumulative deformations f iltered by the global change index (GCI) and local change index (LCI) ( b).Insets are histograms of the displacements to show the frequency.

3.2 Filtering of the SBAS-InSAR Results by Global Change Index (GCI) and Local Change Index(LCI)

Figure 7 a shows the cumulative deformation in the LOS direction derived from the SBAS-InSAR for the ascending track.The SBAS-InSAR failed to work in void locations without cumulative deformations due to low coherence and a lack of valid InSAR measures.We show the deformations of the ascending track in the range from−200 mm to 200 mm.Most of the SBAS-InSAR derived cumulative deformations (about 93%) are within the range of ±50 mm.There are more negative deformation pixels(about 61%) than positive ones (about 32%).Negative displacements indicate the increase of distance between the monitored point and the SAR sensor onboard the Sentinel-1 satellite, meaning that the point moved away from the satellite in the LOS direction, whereas positive values in both maps indicate decreases of the distance and the point moved toward the satellite.

Although large parts of the region have deformations, it is unlikely that all these deformations are related to landslides.Using the 3% and 97% percentile thresholds in the GCI and LCI histograms, we f iltered out pixels within the GCI and LCI range of 3%–97%.We showed the remaining pixels with their cumulative displacements in Fig.7 b, which accounts for about 2.9% of the originally SBAS-InSAR derived pixels.For the remaining pixels, the proportion of pixels within the range of− 50 mm–50 mm falls to 55%.These remaining pixels (Fig.7 b)are supposed to have monotonously changing deformations.

Fig.8 Cumulative deformations in the line of sight (LOS) direction of the descending track ( a) of the study area from April 2020 to April 2022 and the cumulative deformations f iltered by the global change index (GCI) and local change index (LCI) thresholds ( b).Insets are histograms of the displacements to show the frequency.

Figure 8 a shows the original SBAS-InSAR derived cumulative displacements of the descending track.Similarly,deformations of most pixels (> 88%) are within −50 mm–50 mm for the descending result.There is also a large proportion (> 10%) of the pixels with deformations larger than 50 mm.Figure 8 b shows the GCI and LCI f iltered results,which account for 3.8% of the original pixels after GCI and LCI thresholding.In the original result, most displacements are within +/−50 mm, whereas proportions of pixels outside that range increased from 11 to 53%.In Particular,proportions of pixels with displacements between 50 and 100 mm increased signif icantly from 8 to 32% compared to the original results.Filtered by the GCI and LCI thresholds,there are decreases of pixels for all eight displacement categories (Table 1).For example, > 96% pixels were removed in categories −50 mm–50 mm of both tracks.Decreases in most other displacement categories are less than 90%.Minimum decrease of pixels occurred in the category of −150 mm to −100 mm for both tracks.Red pixels in Fig.8 b are located on top of the mountain, where the elevation is about 3,500 m.We checked these locations on Google Earth’s high spatial resolution images and found that they are probably related to freeze-thaw ef fects.

With 3% and 97% percentile thresholds of the GCI and LCI, we f iltered out about 97.1% and 96.2% of original SBAS-InSAR derived displacements for the ascending and descending tracks, respectively.We made colored scatterplots to show the relationship between the GCI, LCI, and displacements for the remaining pixels (Fig.9).All these remaining points are located in zone I and zone II.In zone I, points with the smallest LCIs tend to have larger displacements.In zone II, points with the largest LCIs have large displacements with negative signs.These f indings mean that:(1) there are noises in our SBAS-InSAR derived displacements; and (2) if a landslide moves very slowly, the noise may take over the deformation signal in the LOS direction,resulting in low LCI in zone II and high LCI in zone I.

3.3 Validation

To validate, we selected the deformation results of some known landslides to verify the results of GCI and LCI.We also compared our results with a traditional method.

3.3.1 Validation with Ground Result in a Typical RegionWe selected a typical area with known landslides in the descending track for further analysis (Fig.10).Among these known landslides, the Beishan landslide on the westmost side has partly collapsed in January 2021.Two landslides on the east side (Huanian landslide, Shuidi landslide) havebeen deforming slowly.The black polygons in Fig.10 a are pixels with GCI < 3% or > 97% and the red polygons are pixels with LCI by the same percentile thresholds.Although their spatial patterns are similar, these two types of polygons do not always overlap.There are many overlaps of the two types of polygons for the three known moving landslides(Beishan, Huanian, and Shuidi), the west two of which were validated by f ield reconnaissance in August 2022 (Fig.10 b and c).Figure 10 b is a photo taken near the Huanian Village and Fig.10 c was taken for the Beishan landslide.In addition,we selected four points to represent four types of pixels in this subregion.P1 meets the LCI threshold criteria but not the GCI.P2 meets the GCI threshold criteria instead of the LCI.P3 meets neither criterion and P4 meets both criteria.We plotted their time series displacements in Fig.10 d.Time series of P1, P2, and P3 show only minor displacements at the end of the study period compared to P4.These three points do not have a monotonously decreasing trend and their interannual f luctuations of displacements are large compared to their overall displacements.Instead, P4 shows a consistent decreasing trend, indicating a persistent move away from the sensor in the LOS direction.Note the local increase between some neighboring measures takes its toll in the LCI.Figure 10 d indicates that the combined use of the GCI and LCI is better than any single measurement.

Table 1 Numbers and percentages of the original and fi ltered pixels for ascending ( asc.) and descending ( desc.) tracks ( Disp.is cumulative displacements; Ori. pix/perc. is the number/percent-

Fig.9 Scatterplots to show the relationship between the global change index (GCI) and displacement for remaining pixels f iltered by the GCI and local change index (LCI) thresholds of the ascending ( a) and descending tracks ( b).Colors of the points show the LCI value of each pixel.

Fig.11 The cumulative displacement of a landslide under the descending ( a) and ascending tracks ( b).( c) and ( d) are images of the landslide from Google Earth in August 2019 and November 2021.For each (descending and ascending) track, four points were selected to show deformation time series f iltered by global change index(GCI) and local change index (LCI) ( e and f).P1 and P5 meet both threshold criteria.P2 and P6 meet GCI threshold criteria but not LCI.P3 and P7 meet LCI threshold criteria but not GCI.P4 and P8 meet neither.

By visual comparison, we found that most of our detected landslides are consistent with an existing landslide database (Dai et al.2023).To further validate our results, we divided the landslides detected by GCI into three categories by examining the monitoring results of the ascending and descending tracks.The f irst category is landslides that can be monitored by both the ascending and descending tracks but the GCI values have the opposite trend (GCI > 886 under descending and GCI < 296 under ascending / GCI< 133 under descending and GCI > 1,596 under ascending).Figure 11 shows a landslide of this category.Figure 11 a and b show that the deformation monitoring area in the descending and ascending tracks is not completely consistent.Figure 11 c and d show the development of landslide back wall from high spatial resolution Google Earth images.For instance, in Fig.11 c we can see that the collapse is not obvious.Yet in Fig.11 d the partial collapse becomes severe.We plotted their time series displacements in Fig.11 e and f.In the deformation time series under the descending and ascending orbits, P1 and P5 meeting the GCI and LCI can better represent the monotonic deformation.P4 and P8,which conform to neither GCI nor LCI, are highly volatile.The second category is landslides that can be monitored by both the ascending and descending tracks and the GCI values have the same trend (GCI > 886 under descending and GCI > 1,596 under ascending / GCI < 133 under descending and GCI < 296 under ascending).Figure 12 shows another landslide detected by the descending and ascending tracks.There were some cracks on the slope in August 2019.Compared with August 2019, we can see that there were many partial collapses on the slopes in December 2020.The displacement information of landslide time series is shown in Fig.12 e and f.The third category is landslides that can be monitored by either the descending or ascending track.Typical landslides of this category are shown in Fig.10.

Fig.12 The cumulative displacement of a landslide under the descending ( a) and ascending tracks ( b).( c) and ( d) are Google Earth images of the landslide in August 2019 and December 2020.For each(descending and ascending) track, four points were selected to show deformation time-series f iltered by global change index (GCI) / local change index (LCI) ( e and f).P1 and P5 meet both threshold criteria.P2 and P6 meet GCI but not LCI threshold criteria.P3 and P7 meet LCI but not GCI threshold criteria.P4 and P8 meet neither.

3.3.2 Comparison of the New Method with the Traditional Method

Figure 13 are two violin plots of the last displacement images from the descending and ascending tracks.With means (＆#x1d707;desc.= −11.44 and＆#x1d707;asc.= 6.79) and standard deviations (＆#x1d70e;desc.= 23.53 and＆#x1d70e;asc.= 29.99) from both tracks,we set two thresholds, that is, μ ± σ and μ ± 2σ, to f ilter out the deforming pixels not within these ranges.With the traditional method, 78.5%, 77.1% and 94.4%, 92.9% pixels were removed with the f irst (μ ± σ) and second (μ ± 2σ)thresholds for the descending and ascending tracks respectively.We used straightforward thresholds to eliminate all smaller deformations in the result.In contrast, many large deformations (for example, < μ – 2σ or μ + 2σ) that meet those two traditional threshold crietria are removed by our method due to non-monotonicity in their time series.

Fig.13 Violin plots of the last displacement images for the descending (a) and ascending (b) tracks.Traditional thresholds of μ ± σ and μ ± 2σ are shown as horizontal dotted/dashed lines.Displacement results f iltered by our proposed method is shown as grey spindles within in the plots.

We compared the spatial pattern of identif ied landslides in Fig.14.There are lots of noises (especially the southern part of Fig.14 a) with the μ ± σ threshold, whereas a large part of the real deforming landslide was eliminated by the μ ± 2σ threshold (Fig.14 e).In comparison, our proposed thresholds can not only repress most of the background noise but also preserve real deforming landslide pixels (Fig.14 c and f).We did not show the time series of cumulative displacements for the landslide in Figs.14 d–14f, because deformation time series derived by our method are all monotonic pixels regardless of their last cumulative displacements (as demonstrated in Figs.10, 11, 12).

4 Discussion

InSAR has been widely used in deriving retrospective deformations of occurred landslides or already-known slowmoving landslides.Moving landslides typically have more monotonic deformation time series.However, how to locate monotonously changing locations in SBAS-InSAR results remains challenging.This work solved the problem by quantifying the monotonicity of the time series of cumulative displacements.

4.1 Ef fects of the Method to Detect Landslides in SBAS-InSAR Results

It is dif ficult to identify moving landslides from time series of displacements derived by SBAS-InSAR.Although there are previous attempts to reduce background noises (Bianchini et al.2012; Barra et al.2017; Lu et al.2019; Solari et al.2019; Zhang, Zhu, et al.2021), few works explored time series of the information to extract moving landslides.Although Urgilez Vinueza et al.( 2022) used time series displacements to derive acceleration/deceleration of landslides,the linear regression models in their work require dense SAR images and could easily neglect local outliers, which is also important to recognize landslide pixels.Based on the assumption that displacements of moving landslides should change monotonously, this work fully exploited time series of displacements to quantify changing trends in time series displacements.Using 3% and 97% percentile thresholds in GCI and LCI, we extracted moving landslides by f iltering out > 96% non-monotonously changing displacements,which are unrelated to landslides.Distinct to other previous works that exploit spatial patterns (Bianchini et al.2012;Lu et al.2019; Solari et al.2019; Zhang, Zhu et al.2021),this proposed method is pixel-based, which can fully exploit time series deformations in single pixels.An advantage of this method is that it not only can detect large landslides covering many pixels but also works well for very small landslides of a few pixels.Compared to the original cumulative displacements derived from the SBAS-InSAR, we can focus on slopes that are more likely to be moving landslides.

With the proposed method, more pixels in minor deformation categories (for example, > 96% pixels with − 50 mm–50 mm displacements) and a smaller number of pixels in larger deformation categories (for example, < -50 mm or> 50 mm displacements) were removed.It indicates that this method is more ef ficient to repress noises with minor deformations.Most derived displacements in the original SBASInSAR results are within this range (− 50 mm–50 mm) and extremely slow-moving landslides within this deformation range may be submerged.This method could extract these landslides from background noises.In addition, there are also pixels in larger deformation categories that are removed(for example, 91.71% removed in the ascending track for 150–200 mm), indicating that not all large deformations are moving landslides.Monotonously changing deformations have always been used in existing literature to show moving landslides (Intrieri et al.2018; Ouyang et al.2019).By using this method, we can easily locate optimal pixels to demonstrate time series displacements of moving landslides.

4.2 Implications for Optimally Monitored Slopes by InSAR

Theoretical maximum GCI/LCI indicates rigorous monotonous decrease of displacements, whereas theoretical minimum GCI/LCI indicates rigorous monotonous increase of displacements.In this work, we observed that the actual maximum values of GCI and LCI are identical to their theoretical maximum values.This finding suggests the presence of strict monotonicity in certain decreasing time series displacements.These monotonously decreasing pixels are landslides that move away from the sensor in the LOS direction.Take the ascending track as an example, these are probably moving slopes with east aspects that are parallel to the LOS direction.For these east facing slopes, the slope angle should be smaller than the complement of the incidence angle, otherwise there are shadows in SAR images.Therefore, slopes with east aspects are optimal positions for the ascending track and slopes with west aspects are optimal positions for the descending track (Bianchini et al.2012;Zhang et al.2022).

On the other hand, minimum GCIs/LCIs are much larger than the theoretical minimums, meaning that there is no monotonous movement of landslides towards SAR sensors in the LOS direction.For the ascending track, landslides that move towards the sensor probably are located on west facing slopes, which should not be too steep to have layover ef fects on SAR images.The move of landslides on these west facing gentle slopes will result in a continuous decrease of distance between the moving block and the SAR sensors, which will lead to a continuous decrease of cumulative displacements(small GCI and LCI).Therefore, displacements with low GCIs indicate that north slopes are the second optimal slopes for the ascending track and south slopes are the second optimal slopes for the descending track.These speculations are substantiated by plotting pixels that are retained by this method with aspects (Fig.15).In addition, with time series displacements, InSAR could detect landslides with the size of a few pixels.Considering SAR imaging mode, translational landslides and heads and toes of rotational landslides in the LOS direction may be the easiest type to be detected.

Fig.15 Pixels retained by the proposed method plotted with aspects.Distributions of pixels with small (with a threshold of < 3%) and large (> 97%) GCIs on dif ferent aspects for the ascending track are plotted in ( a) and ( b).The optimal slopes for the ascending track are south-facing and the second optimal slopes are west-facing.Distributions of pixels with small (with a threshold of < 3%) and large (>97%) GCIs on dif ferent aspects for the descending track are plotted in( c) and ( d).The optimal and second optimal slopes for the descending track are north and east aspects, respectively.

4.3 A Caveat to Use Percentile Thresholds

To remove background noise, this work used 3% and 97%percentile thresholds in histograms of the GCI and LCI.Although the selection of thresholds is empirical, it is based on the assumption that the majority of the study area are stable and pixels at the tails of the LCI and GCI histograms are moving landslides.It is evident that stricter thresholds (for example, < 1% or > 99%) will remove more pixels and are more likely to identify moving landslides.In this work, we found that both histograms of the GCI are asymmetrical and applying the same tail thresholds may not be optimal.For example, as the actual minimum GCI is 28, deformations of some pixels with GCI smaller than the 3% threshold could be non-monotonous.Although the GCI thresholds are moderate criteria, they can remove 96% noises.In contrast, histograms of the LCI are more symmetrical and could further repress non-landslide pixels.Although the proposed method could repress a majority of noises with percentile thresholds, it should be used in regional studies instead of for individual slopes.This is because percentile thresholds depend on samples of the study area.Future studies should identify f ixed thresholds for the GCI and LCI to pick up truly monotonously changing pixels.

4.4 Comparisons to Pixel Of fset Tracking Methods

InSAR is very sensitive to detect centimeter deformations in the LOS direction (Zhang et al.2022).In particular, east and west aspects of gentle gradients (less than complementary of the SAR incidence angle) are optimal slopes for InSAR of ascending and descending tracks,respectively.However, the result that there is an absence of minor displacements (< 30 mm) with the largest LCI points in Fig.9 means that SBAS-InSAR could still have dif ficulty in identifying extremely slow-moving landslides even for its optimal slope aspects.The f inding that the minimums of the GCI and LCI are much larger than 0 indicates that the performance of InSAR deteriorates on second optimal aspects of the west and east aspects for the ascending and descending tracks, respectively.Therefore,monitoring extremely slow-moving landslides in other aspects would be more challenging.

The movement of landslides is three-dimensional,whereas InSAR is only capable of monitoring one-dimensional motions of a certain (LOS) direction (Shi et al.2018).This work further substantiates that mountainous terrains pose signif icant challenges for InSAR to monitor moving landslides.For example, there could be many omission errors by using InSAR to detect moving landslides because of the spatial conf iguration of SAR sensors and topography.In contrast, pixel of fset tracking (POT) methods are sensitive to monitor horizontal two-dimensional motions, although the accuracy of POT depends on the spatial resolution of used images (Stumpf et al.2017; Lacroix et al.2018).Previous work shows that POT is capable of detecting surface change of 1/30–1/10 pixels (Leprince et al.2007; Provost et al.2022).For the 30 cm resolution WorldView imagery,POT can detect surface change of 10 mm–30 mm in the horizontal direction.Higher spatial resolution images can be easily taken by airborne vehicles, which means that more subtle surface changes can be detected for all slope conf igurations.The above points indicate that joint use of InSAR and POT may be more ef ficient to monitor moving landslides in mountain environments.

5 Conclusion

Cumulative displacements of landslides in MT-InSAR results should change monotonously with time.However,selecting pixels with monotonous displacements from MT-InSAR results can be dif ficult.This work proposed a method to quantify monotonicity of the time series of displacements in MT-InSAR results that can also be used to detect moving landslides.Dif ferent from the frequently used space cluster-based methods, the proposed method fully considers time series information for each pixel,which ensures that small landslides of a few pixels could also be detected.Using this method, we found that large deformations (> 50 mm or < − 50 mm) may not be moving landslides, whereas most small deformations (− 50 mm–50 mm) are not landslides.This method can remove more than 96% displacements in the original SBAS-InSAR results.We further found that east and west aspects are the best slopes to show ideal time series of displacements for moving landslides for the ascending and descending tracks, respectively.Moving landslides on these slopes tend to show monotonously changing deformations in the SBAS-InSAR results.Due to the special mechanism of SAR imaging, monitoring of moving landslides on other slope aspects would be challenging.

Acknowledgment This work is supported by the Second Tibetan Plateau Scientif ic Expedition and Research Program (STEP, Grant No.2019QZKK0906).

Open AccessThis article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source,provide a link to the Creative Commons licence, and indicate if changes were made.The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material.If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder.To view a copy of this licence, visit http:// creat iveco mmons.org/ licen ses/ by/4.0/.