Regional Recognition and Classification of Active Loess Landslides Using Two-Dimensional Deformation Derived from Sentinel-1 Interferometric Radar Data

Abstract

Identification and classification of landslides is a preliminary and crucial work for landslide risk assessment and hazard mitigation. The exploitation of surface deformation velocity derived from satellite synthetic aperture radar interferometry (InSAR) is a consolidated and suitable procedure for the recognition of active landslides over wide areas. However, the calculated displacement velocity from InSAR is one-dimensional motion along the satellite line of sight (LOS), representing a major hurdle for landslide type and failure mechanism classification. In this paper, different velocity datasets derived from both ascending and descending Sentinel-1 data are employed to analyze the surface ground movement of the Huangshui region (Northwestern China). With global warming, precipitation in the Huangshui region, geologically belonging to the loess basin in the eastern edge of Qing-Tibet Plateau, has been increasing, often triggering a large number of landslides, posing a potential threat to local citizens and natural and anthropic environments. After processing both SAR data geometries, the surface motion was decomposed to obtain the two-dimensional displacements (vertical and horizontal E–W). Thus, a classification criterion of the loess landslide types and failure mode is proposed, according to the analysis of deformation direction, velocities, texture, and topographic characteristics. With the support of high-resolution images acquired by remote sensing and unmanned aerial vehicle (UAV), 14 translational slides, seven rotational slides, and 10 loess flows were recognized in the study area. The derived results may provide solid support for stakeholders to comprehend the hazard of unstable slopes and to undertake specific precautions for moderate and slow slope movements.

1. Introduction

Chinese loess landslides, mainly located in the upper and middle sectors of the Yellow River in Northwest China, are recurrent and prevalent phenomena on the Loess Plateau [1]. A series of factors, involving geological, geomorphological, and hydro-geological properties of loess contributes to their erosion and thus to failure. According to previous statistics, 50% of large-scale landslides estimated have a direct relationship with rainfall [2,3]; in addition, more extreme rainfall events have given rise to increasing occurrence of landslides as feedback to environmental fluctuation, especially since the recent climate-change context [4,5]. Many researchers have reported that Qing-Tibet Plateau (QTP) has experienced warming and wetting trends, and the precipitation amount and the number of precipitation days have increased drastically during past decades [6,7]. For this reason, the upper reaches of the Yellow River-Huangshui Basin, located in the northeastern edge of QTP, became a disaster-prone region where numerous loess landslides occurred in recent years (Figure 1) [8].

Besides, under the implementation of the strategy for the infrastructure development of Western China, road and tunnel construction, mining activities, as well as cutting slopes involving the loess terrace were promoted, however very often triggering a large number of unstable slopes, due to incorrect land management. Moreover, other triggering factors including seismic activities and freeze–thawing effects also make the loess slopes unstable and more vulnerable to topple, fall, slides, and shallow flow, leading to numerous causalities, loss of economic resources, and even ecological and environmental issues [9,10,11]. Additionally, many studies have been conducted on hazard assessment and triggering mechanisms of landslides [12,13,14,15,16,17,18]. However, landslides inventory mapping (i.e., creation of a database with distribution, boundary, state of activity and type of phenomena) represents the primary and leading task for landslide hazard assessment and risk mitigation [19].

Over the past years, the Chinese government paid great attention to preventing soil loss and protecting the ecological environment; indeed, a full spectrum of loess-slope failures was completely inventoried via multi-source technologies. Traditional field survey, as a primary tool for defining the location and failure types, has been widely applied in landslide inventory mapping [20,21]. However, this kind of investigation is time consuming, costly, and mostly focused on individual landslides or limited areas. Recent studies have increasingly utilized high-resolution remote sensing (RS) images and unmanned aerial vehicle (UAV) photogrammetry to delineate the spectral, textural, morphometric, and topographic properties of loess landslides in a wide region [22,23]. Although the general condition of these hazards (i.e., main scarp, fissures, flank, toe, etc.) can be obtained from optical images, it represents a qualitative analysis of the state-of-the-art and cannot reflect and help mapping and evaluating instability [24]. Interferometric synthetic aperture radar (InSAR) and a series of multi-temporal InSAR (MT-InSAR) techniques enable the measurement of millimetric ground displacements and analyze the distribution of active slope, the definition of moving boundaries, and the computation of cumulative movements [25,26]. However, the definition of the slope movements, either considering the rate and the type, is still challenging for local authorities and administrators, when they are in charge of taking early-warning actions. Loess slopes present different failure modes and the displacements vary according to several factors, such as type of movement, material size, geotechnical characteristics, etc. Generally, shallow movements of loess deposits can be grouped into two main categories: slides and flow-like phenomena. Translational and rotational slides, as defined by Cruden and Varnes scheme [27], generally involve large slope sectors and a large amount of material, with usually moderate velocities. By contrast, loess flows or flowslides occur in a very short time and evolve, generating catastrophic consequences, even when the initial deformation is relatively small [28,29]. Therefore, it is necessary and paramount to further put forward classification of these active slopes derived from InSAR results, with the purpose of understanding landslide kinetics and taking corresponding precautions accordingly to different landslide forms.

Technically, differential synthetic aperture radar interferometry (DInSAR) sensors work regularly in side-looking geometry and the calculated displacement velocity is along the radar line of sight (LOS). For slightly small incidence angle (usually between 23° and 40°), the sensor is much more sensitive to vertical deformation than horizontal deformation [30]. Hence, the resulting datasets can estimate only a small component of the real motion of the landslide, i.e., the projection along the satellite LOS. This aspect may cause large errors of slope motion direction and misleading of the classification of slope movement types when using a one-dimensional component instead of real motion of the landslide. Fuhrmann et al. [31] assessed uncertainty of the horizontal and vertical motion and pointed out these errors are related to the magnitude of movement. Aiming at overcoming this limitation, the exploitation of multiple orbits, incidence angle, and different geodetic techniques has become increasingly popular. For example, the authors of [32,33,34,35,36] combined synthetic aperture radar (SAR) ascending and descending information to extract the vertical and horizontal components. Fialko [37] and Hu [38] calculated the three-dimensional displacement via integrating InSAR and Offset-tracking methods that can derive azimuth measurement. Vollrath et al. [39] proposed an effective fusion of InSAR and Global Position System (GPS) network to estimate 3-D movements. Even if each of the approaches has its strengths, there is no unique and well-developed strategy to recognize fully 3-D active landslides since different survey regions have various deformation features as well as different availability of SAR and external data [40].

In this paper, landslides were detected and classified according to the decomposition of interferometric vectors. Ascending and descending acquisitions of Sentinel-1 mission were combined to generate two-dimensional displacements and thus to detect and map active loess landslides. A further step of this work is represented by the identification of the type of loess landslides according to SAR data. The proposed procedure here presented includes four steps: (1) detection of unstable areas obtained by processing and interpreting Sentinel-1 data, based on the average deformation velocity of both orbit datasets; (2) 2D deformation distribution of the whole area, via decomposition from LOS displacement into vertical and east–west horizontal components; (3) classification scheme of loess slides according to different patterns of deformation recognized and ancillary information, where available; (4) landslide inventory map update by applying the proposed scheme at a regional scale. To further validate the proposed scheme, three case studies are presented, where DInSAR information was integrated with UAV and field surveys carried out in the selected areas. This proposed scheme could represent a significant boost for rapid mapping, i.e., in case of emergency following earthquakes or extreme rainfall events, for quickly identifying and classifying landslides. In this way, priorities of interventions may be established, as well as a step forward to the comprehension of loess slopes dynamic can be obtained.

2. Study Area

The study area is a tectonic margin zone, covering about 1.6 × 10⁴ km², between the Tibet Plateau and Chinese Loess Plateau, where thrust and strike-slip fault systems developed during the Mesozoic and Cenozoic (Figure 2), indicating the northeastward growth of the Tibet [41]. Simultaneously, under the influence of a series of tectonic activities and surface erosion, sharp gorges and narrow river catchments consisting of Yellow river, Huangshui river, and Datong river catchments formed in different ages (Figure 2) [42]. The stratigraphic setting of this region is made of Proterozoic gneisses and schists, Cambrian gray limestones and green basalts, and Mesozoic clastic sediments. The lithology is comprised of red mudstones, sandstones, siltstones, and conglomerates [43]. From the hydrological point of view, the Huangshui River, 374 km long, originates in the north part of Tibet at an elevation of 4000 m, and flows eastwards into the Yellow River in Gansu Province, becoming a first-order tributary. The study area is characterized also by very steep slopes and narrow valleys, with height difference up to 900 m [8]. The regional climate is characteristically plateau arid and continental semi-arid with 34 °C during summer and −20 °C during winter on average. The mean annual precipitation is 496.7 mm, concentrated between June and September. However, considering the climate fluctuation during recent years, the rainfall rate has increased almost twice from 2015 to 2019, according to the statistics of the meteorological station (Figure A1). This trend has induced an increase in the instability of the loess slopes of the Huangshui region, posing a serious threat to the local population and activities.

3. Dataset and Methodology

3.1. Dataset

In order to cover the whole study area, 35 ascending images and 60 descending images with vertically transmitted vertically received (VV) polarizations data were acquired, spanning the time interval from January 2018 to March 2019. Since the region of interest is located at the boundary of two descending bursts, extra concatenate sequential bursts from Frame 466,471 were processed (Figure 2). Precise Orbit Ephemerides (POD) from ESA (https://qc.sentinel1.eo.esa.int/) were used to remove orbit errors and improve co-registration accuracy. Besides, an external shuttle radar topography mission (SRTM) digital elevation model (DEM) of 30 m resolution, downloaded from US Geological Survey’s EROS (Earth Resources Observation and Science) Data Center, was used to remove the topographic contribution and eliminate the noise. In order to validate our proposed classification, optical remote sensing images were downloaded from Google Earth, and field UAV investigations were carried out on 30 November 2019, to generate surface texture characteristics and digital surface models (DSMs) of typical loess failures.

3.2. DInSAR Processing Chain

DInSAR is a powerful monitoring technology used to detect ground displacement from centimeters to millimeters by calculating differential phases between two consecutive periods [44,45]. However, limited by systematic noise, atmospheric phase delay, and geometrical decorrelation, phase values between the sensors and target cannot be reliable, thus, a series of multi-temporal approaches were developed including Permanent Scatterers (PS) [46], small baseline subset (SBAS) [47], stable point network (SPN) [48], Standford Method for Persistent Scatters (StaMPS) [49], Coherent Pixels Technique (CPT) [50], Permanent Scatter Pairs (PSP) InSAR, and the SqueeSAR [51]. Among these advanced algorithms, PS and SBAS are two main categories in terms of selecting analysis targets. The focus of this work is to scan wide areas, to spot unstable, localized zones which could correspond to slope instabilities or ground subsidence. Considering the intermittent behavior of loess landslides, the unfriendly radar scenario (soft rocky outcrops and boundaries with rapid evolution) and the small size of occurring instabilities, the SBAS approach turned out to be the most suitable monitoring solution. Four main steps were followed:

1 Generation of interferograms

After preprocessing of Sentinel-1 Terrain Observation with Progressive Scan images (updating POD orbit parameter, burst extracting, and coregistration, etc.), a multi-looking of 10 × 2 was applied to reduce the speckle noise and produce a 25 × 25 m ground resolution cell. A temporal baseline threshold of 90 days and a spatial baseline threshold of 120 m were set, to generate interferograms for ascending and descending datasets from paths 128 and 135 by using commercial differential software, (Gamma Remote Sensing, Bern, Switzerland, https://www.gamma-rs.ch/software).

2 Distributed Scatterers (DSs) Selection

Since most of the unstable slopes in the study area are located in a rural environment, selecting high quality and density scatters is a challenging but crucial step for landslide detection. Groups of pixels with homogeneous radar reflectivity were selected. Such DSs represent the average phase of several objects with different scattering properties [52].

3 Estimation of mean deformation map

The phase of DS points of each interferogram at time t_a, t_b consist of four main deformation compositions [47]:

$$\mathsf{\Delta }\phi (x,r)=\frac{4\pi }{\lambda }\left((d(t_b,x,r)-d(t_a,x,r))+\mathsf{\Delta }{\phi }_{top}(x,r)+\mathsf{\Delta }{\phi }_{atm}(x,r)+\mathsf{\Delta }{\phi }_{noise}(x,r)\right)$$

(1)

where λ is the radar wavelength, d(t_b,x,r), d(t_a,x,r), are the ground displacement of every distributed scatter pixel, Δφ_top is the topographic phase error obtained from SRTM, Δφ_atm is the atmospheric phase delay, Δφ_noise is the noise effect from systematic thermal noise, orbital and temporal decorrelation. In order to relieve topographical, atmosphere errors, and systematic noise in the time-series interferometric pairs, seeking a minimum norm least-squares using singular value decomposition (SVD) of DSs was utilized. Note that the atmospheric phase is mitigated by high-pass temporal and low-pass spatial filtering based on statistical patterns [53].

4 Decomposition of interferometric vectors

The combination of ascending and descending data is a standard procedure to obtain the horizontal E–W and vertical displacement components (Figure 3) [31,32,54]. The two-dimensional velocities can be obtained by implementing a linear conversion, as indicated by the following equations [32,37]:

$$\left(\begin{array}{l}{V}_{asc}\\ V_{des}\end{array}\right)=\left(\begin{array}{l}\cos {\theta }_{asc}\sin {\theta }_{asc}\\ \cos {\theta }_{des}\sin {\theta }_{des}\end{array}\right)\left(\begin{array}{l}{V}_{\mathrm{v}}\\ V_{\mathrm{h}}\end{array}\right)$$

(2)

$$V_{\mathrm{h}}=\frac{V_{asc}\cos {\theta }_{des}-V_{des}\cos {\theta }_{asc}}{\sin {\theta }_{asc}\cos {\theta }_{des}-\sin {\theta }_{des}\cos {\theta }_{asc}}$$

(3)

$$V_{\mathrm{v}}=\frac{V_{asc}\sin {\theta }_{des}-V_{des}\sin {\theta }_{asc}}{\cos {\theta }_{asc}\sin {\theta }_{des}-\cos {\theta }_{des}\sin {\theta }_{asc}}$$

(4)

where V_asc, V_dec is the velocity from ascending and descending geometries. V_v, V_h is the vertical and east–west values, θ_asc, θ_des are the incidence angles for the ascending and descending orbits (39.4° and 41.5°, respectively).

3.3. UAV Survey and Data Processing

To acquire high-resolution optical images and detailed topographic features, a field survey was performed on 30 November 2019 using a DJ-Phantom 4 Pro aircraft (Dajiang Baiwang technology company, Shenzhen, China). The flight missions were pre-programmed through the open source Altizure software (https://altizure.com). A total of 1209 photos were captured with 10 cm ground sampling distance (GSD) over the three selected cases. Ground control points (GCPs) were not used in our survey, since precise surface elevation data were not necessary for the aimed scale of analysis. A bundle adjustment process using image match algorithms was applied for raw data, and orthorectified images digital surface models (DSMs) were thus generated. All the processing procedure was performed by Pix4Dmapper photogrammetry software, (Pix4D company, Ecublens, Switzerland, https://www.pix4d.com/try-software).

3.4. Classification Scheme of Loess Landslide

A criterion to identify and classify landslide types occurring on loess slopes is here proposed, according to the analysis of deformation characteristics, such as deformation rates, main direction, and spatial distribution derived from DInSAR data (

Table 1. Classification scheme of loess landslides based on deformation characteristics. Note: H represents the horizontal direction; V represents the vertical direction; V_h represents the horizontal displacement velocities; V_v represents the vertical displacement velocities; High V_h represents high magnitudes velocities in the horizontal direction, High V_V, V_h illustrates high velocities in both vertical and horizontal directions.

Type	Major Movement Direction	Deformation Velocities	Deformation Location
Translational slides	H	High Vh	Whole landslide body
Rotational slides	V	High Vv, Vh,	Source area
Loess flow	H	High Vh	Parallel to the slope

). Three main types of landslides can be distinguished on loess formations in Tibet: translational slides, rotational slides, and flows. It is worth mentioning that this distinction follows the Cruden and Varnes [27] scheme and it is based on basic assumptions derived from literature. According to the above-proposed scheme, three typical loess landslides can be detected and classified. A scheme summarizing the general idea of the classification scheme is in Figure 4, whilst the whole technical flowchart is shown in Figure 5.

• Translational slides: the major displacement direction is conventionally horizontal, the characteristic deformation rate is low, with slow and intermittent movements; the horizontal deformation is typically dominant, much higher than the vertical one, along the whole landslide body. This sort of slide can be usually distributed in slope sectors characterized by gentle topographic gradient.
• Rotational slides: these events are characterized by a balance between horizontal and vertical displacements, however, their distribution in the landslide body is different. The vertical movement is located mostly in the upper part of slopes, whilst the rest of the landslide moves with a direction usually parallel to the slope. This failure is prone to develop in steep slopes with arc-shaped rupture surface under unsaturated loess layers. Field investigation shows that many cracks can be often seen in the crown and shear deformation area, representing the rapid, and even surge sliding velocities. Besides, a large number of fissures, typically present in the middle or toe part of the slope, demonstrate that the deformation direction changes from vertical to horizontal.
• Loess Flows: theoretically, horizontal velocities are predominant, even though a vertical component is present, depending on the steepness of the slope; generally, flows are characterized by a more rapid motion with a very short time interval. However, differently from transitional and rotational slides, loess rapid flows are difficult to detect by single deformation characteristics, since abrupt movement cannot be detected by satellite imagery, therefore here only slow flows or intermittent phases are considered.

4. Results

4.1. Mean Deformation Velocities along the LOS

Since ground-truth data are not available in the study area, the Minhe city buildings, in the center of the study area, were chosen as reference stable point. About 20 million distributed scatterers were extracted both from ascending (26,616,239) and descending (20,996,749) orbits, respectively. The maximum deformation velocities of the two orbits vary from −51 to 53 mm/year and −48 to 49 mm/year, respectively (Figure A2 and Figure A3). Note that the negative value indicates displacement away from the satellite view and positive values refer to displacements towards the sensors. Based on the standard deviation of deformation velocities, all the areas with displacement rate between 0–7 mm/year, 7–15 mm/year, >15 mm/year, were defined, respectively, as stable areas, relative strong deformation areas (RSDA), and strong deformation areas (SDA). Thus, 98% of the study area can be considered as stable, 1.8% can be classified as RSDA, and the remaining 0.2% is SDA (Figure A2) Comparing with the ascending data, the results obtained with the descending orbit (97.5% of the study area is stable, 2.2% is RSDA, and 0.3% can be classified as SDA) show a higher concentration of unstable slopes in the study area (Figure A3).

For evaluating the uncertainty of surface movements along LOS, the standard deviation of each distributed area was estimated under the assumption that a high standard deviation shows lower precision of the corresponding mean velocities and vice versa (Figure A4). The whole area has a general pattern of standard deviation below 5 mm/year as shown in regions A and C, indicating that the general deformation results are reliable. Only small regions (B from ascending and D from descending) have a high deviation value up to 8 mm/year, which may account for the low temporal coherence and abrupt deformation characteristic [54].

4.2. Horizontal and Vertical Deformation Velocities Distribution

According to decomposing Equations (3)–(5), the horizontal and vertical deformation velocity maps were generated. As can be seen in Figure 6, approximately 98.3% of the points show displacements rates between −7 and 7 mm/year, which indicated the overall stability of the study area. About 30 slopes, with total areas of 48 km² were identified as characterized by a clear horizontal motion. On the other hand, only 15 slopes, with total areas of 48 km², were recognized as characterized by a major vertical direction motion (Figure 7). Accordingly, the whole deformation scenario in the Huangshui basin is slightly more prone to shallow horizontal movement than to vertical motion, being an area susceptible to soil erosion induced by rainfall.

4.3. Classification Scheme Validation

In order to validate the proposed classification scheme, three different phenomena near MinHe, i.e., a translational slide, a rotational slide, and a flow, were selected for further UAV investigations, on 8 December 2019.

4.3.1. Translational Slide

A translational slide has been identified in the area of GaJia valley (Figure 6). From the UAV image (Figure 8), the landslide is 1100 m long, 760 m wide, and has 200 m height difference and 20° of average slope degree. The whole topography has low energy relief, separated by deep gullies in the central part, and some cultivated lands. In terms of spatial coverage and magnitude of deformation velocities, horizontal displacement is larger than vertical, with maximum velocity reaching 30 mm/year (Figure 8). According to the profile analysis of horizontal deformation, the upper section of the landslide (0–500 m) is larger than the lower section. Comparing with horizontal, the vertical deformation has two distinct increasing trends downward at 100–300 m and 400–600 m, where typical geomorphological features (e.g., scarps or step terrace) can be observed (Figure 9). According to the two-direction displacements, the real displacements assumed as a slow creep behavior along the slope, and according to the geomorphological features, the landslide seems to be shallow. This motion behavior is very frequent in translational slides, according to the previous classification of loess landslides [21,55]. In terms of hazards, these kinds of movements are usually triggered by abundant rainfalls, posing a potential threat to the local villagers living beneath the loess hill.

4.3.2. Rotational Slide

The rotational slide, selected for the validation analysis, has been identified near Jianzi village (Figure 7c No.2). Based on the visual interpretation of Figure 10a, the slope is divided into four sectors: Sector A is located on the top of the slope, with a length of 600 m and a width of the main crown of about 20 m (Figure 10a and Figure 11a,b). This sector is also characterized by the presence of loess debris that buried a road passing in the area (Figure 11b). Sector B is in the right flank of landslide, characterized by a higher slope angle. Sectors C and D, localized in the central part of the landslide, are characterized by the presence of numerous fissures and cracks. In terms of distribution and magnitude of deformation, vertical displacement velocities are predominant in this slope sector (Figure 10b), while horizontal deformation is higher in the two flanks and in the toe of the landslide (Figure 10c). From the profile analysis of B-B’, it is possible to observe that both vertical and horizontal deformational pattern have the same motion rate (Figure 10d), with maximum deformation rates located over the cracks identified in the field (2 m high) (Figure 11c). According to the synthesized real displacement, the rupture surface of this landslide can be shaped as deep-concave. The profile C-C’ (Figure 10e) highlights a main vertical component of the deformation. However, the motion trend of this part is discontinuous and can be divided into two independent parts. The first starts at 200 m and terminates at the 500 m area, and it is characterized by vertical deformation rate up to 16 mm/year. The lower part has a mean vertical motion displacement up to 12 mm/year. The area connecting Sections B and C is flat (the area estimated from distance 450–720 m), characterized by high horizontal motion along with the generation of large numbers of compressing fissures about 20 cm wide (Figure 11d). Additionally, a group of compressing fissures with NW–SE direction and extensive fissures along NE–SW direction can be observed near the flatform (Figure 10a), implying that this gentle landscape prevented upper loess (Section C) movement forward, and rebalanced sliding surface from depth to near-surface. In the cross-section C-C’, both two flanks present concave terrains (0–230m, 380–630 m) where horizontal deformation is the major motion, while the middle landform (from 230 to 380 m) has relative gentle slope characterized by vertical motions increasing from 10 to 17 mm/year (Figure 10f). Therefore, the integration of the two-dimensional deformation results, UAV images, profile analysis, and previous studies on the loess landslide failure mode [56], allowed to define this landslide as rotational slide with a strip-slip deformation slide rather than a long-term creep based on the previous conclusion.

4.3.3. Loess Flow

The flow-like landslide identified and selected for the validation analysis is located near the Chuang-Da expressway, in MinHe county (Figure 7c, No.3). The photo-interpretation and the field investigation showed that the sliding loess mass deposited in the lower part and blocked the gully in front of the slope (Figure 12a). Two-dimensional displacement InSAR measurements highlight that this slope has experienced both vertical and horizontal movement, particularly in the upper part of the body (Figure 12b,c). It can be also clearly seen from the profile analysis E-E’ that the deformation is increasing continuously proceeding downstream (Figure 12d). The vertical deformation increases suddenly from 2 to 12 mm/year, keeping such high velocity for about 300 m, as well as the horizontal movement. Besides, the cross-profile F-F’ illustrates that this western sector of the slope has both high horizontal and vertical displacement velocities, implying that deposited loess is susceptible to accumulate and flow along the channel (Figure 12e). It is worth noting that no clear horizontal displacement was detected in the toe, most likely due to very fast and large deformation, exceeding the monitoring capability of InSAR techniques. In order to analyze the slope behavior during August 2018, coherence and phase were compared in the time spans August 6–18 and August 18–30. In the first period, coherence shows a generally high value and the interferometric phase seems stable (Figure A5a,b). However, in the following 12 days, the lower part of the slope presents lower coherence and its corresponding interferometric phase seems noisier (Figure A5c,d). Combining with the remote sensing images from Google Earth (Figure A5e,f), the change in such value can be explained with a flow event that occurred during 18–30 August 2018.

5. Discussion

The occurrence of landslides in the Loess Plateau is increasing recently, triggered by more and more diffuse human activities and extreme weather conditions in Qing-Tibet Plateau. Therefore, the detection and the mapping of unstable slopes are significant tasks for stakeholders in landslide hazard assessment and risk mitigation. InSAR technology, with its large coverage and millimeter accuracy, is popularly applied in landslides monitoring and landslides inventory mapping. However, limited by a one-dimensional line of sight viewing geometry, the derived displacement results are not sufficient to capture fully slope movement, especially in the horizontal direction (N–S, and E–W), leading to difficult interpretations, estimation of displacement rates, and evolutions by end-users in charge of disaster prevention. The combination of descending and ascending acquisitions provides the opportunity to detect active phenomena and estimate the E–W and vertical components. Many authors have already used components of movement to characterize single slope [26,32,57]. In this work, a procedure based on the exploitation of both components was enhanced, and leveraging on the wide-area mapping capability of Sentinel-1, it was applied on a basin scale. The main idea is to exploit Sentinel-1 acquisition strategy and a simple post-processing approach to scan wide areas and to highlight unstable zones. Nevertheless, some slopes along the north–south direction are still difficult to detect owing to the inherent limitation of InSAR.

In addition, a criterion to classify loess landslides by two-dimensional deformation velocities, sliding direction, and spatial distribution is proposed. Such classification is based on the decomposition of movement vectors derived by DInSAR information, which allows to obtain vertical and horizontal components of the movement occurring on the ground surface. Taking into account the prevalence or the balance between the two forces, according to their rates and to their spatial distribution, as well as to the basic knowledge on landslide-type features, three main types of landslides occurring on loess slopes were semi-automatically identified. Translational slides are mostly characterized by horizontal movements, and at times parallel to the slope, thus horizontal and vertical displacements have motion rates very comparable, with a slight prevalence of the first; rotational slides present a major vertical displacement in the initiation area, while the rest of the landslide body moves parallel to the slope; flow-like landslides are characterized by generally higher velocities of displacement, very rarely with a major vertical component. In this way, a reliable landslide inventory map of the study area was obtained in a very rapid way.

In Figure 13 a landslide inventory map, generated through the application of the proposed scheme, is depicted. It is made of 14 translational slides, seven rotational landslides, and 10 loess flows, detected according to the decomposition of interferometric vectors and following the proposed procedure. No human activities (e.g., excavations, loading, deforestation) were detected from historical Google Earth image archive, thus, slope movements most likely are triggered by the precipitations, as also stated by previous studies [58,59]. In terms of deformation magnitude, translational slides range from values of 4 to 30 mm/year along the horizontal direction. The average displacement is about 17 mm, indicating a slow-moving trend. Rotational slides have displacements values ranging from 16 to 20 and 8 to 15 mm/year in the vertical and horizontal directions, respectively. Considering the formation mechanism, the location of vertical deformation is prevalent in the source area. Loess flows are characterized by horizontal displacement rate varying from 10 to 30 mm/year, while the average horizontal value is up to 23 mm/year, indicating a fast motion along the slope. It is worth noting that landslides occurring on loess slopes are generally characterized by continuous deformation trends, making our procedure, based on the average velocity information, very well-suited for the detection and the mapping of slow-moving phenomena; however, some loess slides or loess earth flows can be characterized by higher velocities, abrupt or intermittent motion, still representing a challenge for DInSAR approaches, especially during monitoring stages. Moreover, another main disadvantage influencing the disparity is the SAR geometric distortions, i.e., loess landslides located within “blind” areas caused by layover, shadow, or foreshortening of single ascending or descending observations [55]. This aspect would give rise to inconsistency in the deformation magnitude and direction from different observations, thus affecting the decomposition of interferometric vectors and classification results. Nonetheless, combining multi-sensor satellite SAR sensors with different incidence angle, flight direction, and resolution could further upgrade the reliability of the classification scheme [60].

Even though the achieved results are sensitive to side-look image geometry of SAR, the procedure presented here may demonstrate a very important tool for landslide management activities and for landslide inventory mapping, during early-warning phases, especially when dealing with large or inaccessible areas. Considering financial budgets, integration of DInSAR information with other sources may be an added value in order to have as more precise as possible information over entire regions, as proven in this work. However, DInSAR data have been repeatedly demonstrated to be a reliable, cost-effective, and time-saving source of information [18], since an in-situ investigation over a large area such as this would have required much more time and a sensibly higher budget. Additionally, Sentinel-1 SAR products are freely accessible through the Sentinels Scientific Data Hub (https://scihub.copernicus.eu/), thus providing the scientific community, as well as public and private companies, with consistent archives of openly available radar data, suitable for monitoring applications. It is worth mentioning that vertical movements detected in flat areas, attributable to subsidence phenomena (Figure 6b, Figure 7b No.5, Figure 7c No.4), were not considered in this work. However, it must be stated that this classification criterion is based exclusively on movement characteristics, rather than location of sliding/rupture surface, travel distance, volume of displaced mass, or failure mechanism [21]. This classification criteria may be an effective, low-cost, and rapid guide for the stakeholders to take specific actions, according to the detection of different phenomena.

6. Conclusions

The aim of this paper was to provide a rapid and reliable tool capable of generating landslide inventory maps, with basic classification of the type of loess landslides, by applying simple and transferable differential interferometric synthetic aperture radar (DInSAR) processing and post-processing procedures. Two-dimensional displacement velocities, both in horizontal and vertical directions, were extracted by computing and combining ascending and descending Sentinel-1 data, acquired over Huangshui region (Northwestern China). The classification scheme is based on basic knowledge of typical landslide mechanisms developing on loess formations. In particular, three types of landslides, very common and spread over loess slopes, were considered: translational slides, rotational slides, and flows. Translational slides are commonly characterized by prevalent horizontal and parallel to the slope movements, while rotational slides present a dominant vertical component in the source area. Flow-type landslides move parallel to the slope, with considerable velocities, generally higher than the abovementioned types. According to these basic rules, and after validating DInSAR-derived data via field unmanned aerial vehicles (UAV) surveys on three test-areas, a general landslide inventory map of the Huangshui region was generated. In this way, more than 30 landslides were detected and distinguished. The obtained results confirmed the reliability of SAR interferometry for large areas monitoring, and the analysis of landslide mechanisms through the use of post-processing procedure is a valid and consolidated way of estimating the prevalent motion of a failure. Starting from these assumptions, the classification scheme proposed here may represent a valid and transferable way of generating rapid landslide inventory maps containing various information, including those about the type of landslide. The gathering of all this information may provide substantial and precious support to land management activities, in order to take precise and rapid mitigation actions.