Landslide Hazard Assessment Method Considering the Deformation Factor: A Case Study of Zhouqu, Gansu Province, Northwest China

Abstract

Landslides are geological disasters that can cause great damage to natural and social environments. Landslide hazard assessments are crucial for disaster prevention and mitigation. Conventional regional landslide hazard assessment results are static and do not take into account the dynamic changes in landslides; thus, areas with landslides that have been treated and stabilized are often still identified as high-risk areas. Therefore, a new hazard assessment method is proposed in this paper that combines the deformation rate results obtained by interferometric synthetic aperture radar (InSAR) with the results of conventional hazard assessments to obtain the hazard assessment level while considering the deformation factor of the study area, with Zhouqu, Gansu Province, selected as the case study. First, to obtain the latest landslide inventory map of Zhouqu, the hazard assessment results of the study area were obtained based on a neural network and statistical analysis, and an innovative combination of the deformation rate results of the steepest slope direction from the ascending and descending data were obtained by InSAR technology. Finally, the hazard assessment level considering the deformation factor of Zhouqu was obtained. The method proposed in this paper allows for a near-term hazard assessment of the study area, which in turn enables dynamic regional landslide hazard assessments and improves the efficiency of authorities when conducting high-risk-area identification and management.

1. Introduction

A landslide is a movement of mass of soil (earth or debris) or rock down a slope under the influence of gravity and influence of natural (snowmelt, abnormally high rainfall, earthquakes) or anthropogenic stresses (earthworks, vibrations, deforestation, exploitation of materials or aquifers [1,2,3,4,5]. Landslides are ubiquitous phenomena, and they are present all over the world [6], and some areas are more frequently affected, such as mountain slopes subject to increasing anthropisation [7,8,9], high seismicity areas [10], or regions with intense or long periods of low-to medium-intensity rainfalls [11,12,13,14]. Landslides are geological hazards that cause serious damage to the natural and social environments [15]. Landslides cause approximately 1000 deaths and $4 billion worth of economic losses annually [16]. Landslides and related processes were documented, causing more than 61,000 deaths worldwide between 1900 and 2009 [17]. Therefore, to avoid additional human and economic losses caused by landslides, it is important to carry out landslide hazard assessment work, detect high-risk areas, and take appropriate measures to mitigate and manage landslides before they occur.

A landslide hazard assessment is the process of ranking the land surface into zones based on the current or potential landslide sensitivity, hazard, and risk; such assessments are a key step in landslide investigations and landslide risk management [1,18]. Landslide hazard assessment is usually based on the spatial and temporal probability of landslide occurrences and is performed following three main steps: (i) the creation of a landslide inventory map, (ii) a landslide susceptibility analysis, and (iii) a landslide hazard analysis [19,20,21]. Landslide inventory maps play a crucial role in landslide hazard assessments, and the quality and completeness of their investigation affect the reliability of hazard assessments [3,22,23,24]. Different methods and techniques (i.e., aerial photograph interpretation, field survey, high-resolution imagery, and LiDAR interpretation [25]) can be used to build a landslide database. Landslide susceptibility represents the potential occurrence of phenomena based on predisposing factors (topographical, geological, geotechnical, soil, and land-cover conditions) for different slope failure surfaces [26,27]. Landslide hazard analysis relies on the two previous steps and considers the triggering factors (rainfall or earthquakes) [18]. Conventional regional hazard assessment results have contributed greatly to governmental land use management and planning; however, most of these hazard assessment results are static. Most previous regional landslide hazard studies were based on landslide inventory maps, and the correlations between landslides and genetic elements have been statistically analyzed based on predisposing factors in the study area as a regional hazard assessment and hazard classification method [28,29,30,31]. Other scholars have used specific software to calculate safety factors on representative slope sections with stability problems and then to complete a quantitative landslide hazard assessment, but this is only for specific areas and can only be applied to medium-to-low depth landslides and is not universally applicable [32]. There are some scholars have developed a physical-based model based on the assessment of landslides induced by climatic events, which can be applied to different types of landslides and can perform landslide hazard assessments at the meso or broader scale, but the results remain constant over a long period of time and still do not allow for the assessment of landslide hazard over a specific time period [21]. Some scholars have estimated the number of landslides in the study area based on geoclimatic models created from characteristic annual precipitation values and relief indices, which in turn determines the landslide hazard risk. However, the assessment results also have a high degree of uncertainty [33]. In the face of growing economic conditions and abundant human activities, and the loss of risk memory and over-confidence in different protection systems [20,34,35,36,37], the wide distribution of medium-high-risk areas prevents the authorities from effectively utilizing and developing resources, and high-risk areas containing landslides that have been controlled are thus still included in the medium-high-risk level areas when conducting hazard assessments. Therefore, conventional hazard assessment methods cannot meet the current needs of landslide risk management and prevention worldwide.

The development of InSAR technology has made it possible to assess the hazards associated with landslides in their recent stages. Most of the previous hazard assessment methods based on InSAR technology were combined with the one-dimensional deformation rates of radar line-of-sight (LOS) directions acquired by InSAR [38]; these data can also be used to perform landslide hazard assessments in the study area over time. However, due to the different deformation directions of landslides and their positive and negative deformation magnitudes, it is impossible to accurately compare the deformation conditions of different slopes in the same area. Moreover, SAR data from different perspectives are sensitive to different slope aspects, and considering only a single data point can also lead to missing deformation information [39,40]. Therefore, in this paper, we propose a landslide hazard assessment method combining the deformation rate results of different orbit SAR data images in the face of global and take Zhouqu as a case study. First, a latest landslide inventory map of Zhouqu is detected and obtained; a susceptibility assessment is carried out based on the neural network method, followed by a hazard assessment; the deformation rate level (ascending and descending data) of the steepest slope direction is obtained based on InSAR technology and effectively combined; and finally, the hazard assessment levels of Zhouqu are obtained while considering the deformation factor based on the distribution of the deformation rate level in the study area and the hazard assessment level. The landslide hazard assessment method provides landslide hazard assessment results for a specific time period, which can be used in any area worldwide that suffers from landslide hazards, and provides a reference for authorities to control and manage high-risk areas, provides a feasible solution for obtaining accurate ground deformation rate levels during processing, and uses a matrix approach to effectively combine the results generated in this paper. The whole process can provide a theoretical basis and reference for the readers to conduct relevant work.

2. Study Area

The study area, Zhouqu County, is located in the middle reaches of the Bailong River in the southern region of Gansu Province, northwestern China, with a length of approximately 99.4 km from east to the west and a width of approximately 88.8 km from north to south for a total area of approximately 3010 km² and a latitude and longitude coordinate range of 103°51′39″E–104°45′31″E, 33°13′06″N–34°01′00″N, as shown in Figure 1. Zhouqu is located at the intersection area of the eastern edge of the Tibetan Plateau, the western flank of the Western Qinling Mountains, and the Minshan Mountain Range; this mountainous area has experienced extensive tectonics and erosion. The area is characterized by overlapping mountains, steep valleys, gullies, narrow valley channels, steep slopes, and rapid streams. The overall topography is high in the northwest and low in the southeast, with complex topography and large differences in elevation (Figure 1c); the area contains typical temperate monsoon zone alpine valley landforms [41]. Figure 1c shows the distribution of rivers in Zhouqu and historical landslides in Gansu Province. The rivers in Zhouqu belong to the Jialing River system in the Yangtze River basin, and the three main rivers are the Bailong River, Gongba River, and Boyu River. The historical landslide database was obtained from the Geographic Data Sharing Infrastructure, Resource and Environment Science and Data Center (http://www.resdc.cn accessed on 5 October 2022), and landslides have mainly been distributed in the northern and central areas of Zhouqu. Historically, frequent earthquakes around the area and human activities gradually intensifying with economic and social development have revived numerous ancient landslides in Zhouqu and caused extensive losses of life and property, with many new landslides [42].

The stratigraphy in northern and southern Zhouqu belongs to different subdivisions of the Songpan-Ganzi stratigraphy and Qinling stratigraphy, respectively, and the overall orientation is consistent and trends NW-NNW. The region is influenced by geological tectonic activity that makes the Quaternary strata and magmatic rocks exhibit mostly unconformable contacts. The study area is located within two different geotectonic units. It is bounded by the Yangbuliangzi-Dainian line, and the southern part of the study area belongs to the northeastern part of the Songpan-Ganzi fold system, which is less active and less developed in terms of folding and fracturing than the northern region. The northern part belongs to the east-west fold belt of the Qinling Mountains; this area is strongly active and has developed strike faults. During the long-term geological and tectonic development of the region, roughly parallel extrusion zones formed along the northwest spreading direction (Figure 2a).

Zhouqu is located in a transition zone from the northern subtropical zone to the northern temperate zone and is influenced by atmospheric circulation and topography, with obvious vertical climate zoning and two distinctive dry and wet seasons. The climate in the southwestern part of the county is warm and humid, while the northeastern area is cool and dry, and the temperatures in the valleys are significantly higher than those in the mountains. The average annual rainfall in the territory is 434.0 mm, with a maximum daily rainfall of 96.3 mm. The precipitation distribution varies greatly across the county, with more precipitation falling in the southwest area than in the northeast area and more precipitation occurring in the mountains than in the river valleys (Figure 2b) [43].

3. Research Methods

The technical route of this study is shown in Figure 3 and consists of four steps. First, a landslide inventory of the study area was drawn. On the basis of the landslide inventory, a neural network model was used to assess the susceptibility and consider the landslide time probability and event probability factors to obtain the hazard assessment level. Based on InSAR technology, the deformation rate level of the steepest slope direction was obtained, and, finally, the hazard assessment level was combined with the deformation rate level to obtain the landslide hazard assessment level considering the deformation factor of the study area. Among them, the landslide database is the basis of our entire research, so we needed to map the latest landslides as much as possible before conducting the work, which has a great impact on the accuracy of the assessment results. Secondly, if single orbit SAR data are used, the deformation of some regions will be ignored. Therefore, SAR data fusion of different orbits and different wavelengths can better determine the ground deformation and improve the accuracy of the assessment results.

3.1. Drawing the Landslide Inventory Map

First, SAR image data covering Zhouqu County were collected, and the Sentinel-1A data provided by ESA were used in this study; these data included the ascending Sentinel-1A data from 11 January 2020 to 24 January 2022 for view 118 and the descending Sentinel-1A data from 11 January 2020 to 5 February 2022 for view 63. The coverage and specific parameter information of these data are shown in Figure 1b and

Table 1. Basic parameters of Sentinel-1radar satellite image used in the study.

Radar Satellite	Sentinel-1A
Orbital direction	Ascending/descending
Band	C
Radar wavelength (cm)	5.6
Spatial resolution (m)	5 × 20
Revisit period (d)	12
Polarization mode	VV
Incidence angle (°)	36.9/39.2
Acquisition dates	11 January 2020 to 24 January 2022/ 11 January 2020 to 5 February 2022
Number of images	118/63

. Then, the different orbit SAR data were processed based on time-series InSAR technology to obtain the potential deformation zones in the study area, visual interpretation was performed based on the optical images, and, finally, field verification and calibration steps were performed to determine the landslide locations and boundaries. The landslides acquired by InSAR were supplemented with the historical landslide database to obtain the latest landslide inventory map of Zhouqu.

The time-series InSAR method used for landslide detection in this study was the Stacking-InSAR method; this method is effective in reducing atmospheric disturbances [44,45] and is essentially a weighted average of the phase maps obtained by D-InSAR, superimposing the phases of landslide deformation over the entire study time scale. The atmosphere is spatially correlated and temporally uncorrelated, so the atmospheric phase is diluted in the Stacking-InSAR results, while the deformation information is retained, thereby improving the efficiency of landslide detection. Stacking InSAR has a strong perception of areas that have been in deformation for a long time, so it is often used for landslide identification. This technology reduces atmospheric interference rather than removing it, which can effectively save processing time. The principle is easy to understand, and the operation is simple.

3.2. Hazard Assessment Methods

First, correlations were established between the distribution of landslides in the area and geographic environmental factors, such as lithology and faults. Then, the landslide susceptibility assessment results were obtained based on the neural network method. Finally, the landslide event and temporal probability factors were considered to obtain the hazard assessment level of the study area [46].

3.2.1. Susceptibility Assessment

The methods used for this susceptibility assessment were certainty factors and neural networks. The certainty factors were proposed by Shortliffe E H (1976) and further improved by Heckerman (1986) to analyze and study the correlation of factors affecting the occurrence of an event [47,48]. The range of values was [−1, 1], with positive values indicating high correlations and high geohazard probabilities and negative values indicating low correlations and lower likelihoods of geohazards occurring. The utilized calculation formula is expressed as follows:

$$C_{\mathrm{F}}=\{\begin{array}{c}\frac{P{P}_{\mathrm{a}}-P{P}_{\mathrm{s}}}{P{P}_{\mathrm{a}}(1-P{P}_{\mathrm{s}})},P{P}_{\mathrm{a}}\ge P{P}_{\mathrm{s}}\\ \frac{P{P}_{\mathrm{a}}-P{P}_{\mathrm{s}}}{P{P}_{\mathrm{s}}(1-P{P}_{\mathrm{a}})},P{P}_{\mathrm{a}}<P{P}_{\mathrm{s}}\end{array}$$

(1)

where PP_a is the conditional probability of an event (geohazard) occurring in the factor classification data a. In practical studies, this probability is usually expressed as the ratio of the number (or area) of geohazards in the factor classification a to the area ratio of the data classification a. PP_s is the ratio of the total number (or area) of hazards in the whole study area to the total area of the study area.

The correlation between each factor category and landslide hazard can be obtained by a deterministic model, and the same number of landslide and nonlandslide data points were selected for the neural network hazard assessment. The selected neural network method was the multilayer perceptron model (MLP) [49]. The MLP is a multilayer feedforward neural network, a predecessor of deep neural networks that can be solved quickly by setting a small number of hidden layers for simple problems and accurately by setting multiple hidden layers for complex problems to transform the network into a deep neural network that is flexible, versatile, and widely used.

3.2.2. Hazard Assessment

The temporal probability of a landslide occurrence can usually express the frequency of recurrent landslide events, while the frequency indicates the number of events of a specific induced event type (e.g., rainfall, earthquake) at a specific time [19,50]. In this study, the amount of rainfall is the main consideration, and the landslide hazard index H equation is shown below:

$$H=P(s)*P(t)*P(e)$$

(2)

where P(e) is the event probability, expressed as the probability of occurrence of geohazard events; P(s) is the spatial probability, referring to the geohazard susceptibility; and P(t) is the temporal probability, mainly considering the frequency of occurrence of triggering factors such as rainfall.

3.3. Classification of the Deformation Rate Level

SBAS-InSAR processing was performed on the basis of the Stacking-InSAR results to effectively remove atmospheric interference, reduce the spatial and temporal decoherence, and improve the coherence [51]. At the same time, this technique can obtain the annual deformation rate of landslides and extract landslide time series curves, which is convenient for studying the historical deformation state of landslides. The deformation information obtained from the SBAS-InSAR characterizes a one-dimensional projection of the actual surface motion in three dimensions along the satellite LOS direction. To overcome the differences in the deformation rate results caused by the surface geometry and satellite view and to involve InSAR results in the hazard assessment, the deformation rates of slopes with different orientations had to be compared. Therefore, the obtained deformation rates in the LOS direction were all projected to the same direction—the steepest slope direction—to enrich the deformation velocity information [52,53]. The equation is shown as follows:

$$\{\begin{array}{c}{V}_{SLOPE}=\frac{V_{LOS}}{C}\\ C=\cos \beta \end{array}$$

(3)

where V_LOS is the velocity measured by the satellite along the LOS direction, C is the scale factor of the actual three-dimensional surface displacement and the displacement measured by InSAR, and β is the angle between the radar satellite LOS direction and the direction of the steepest slope. However, the method has some limitations, and here, we referred to the solutions adopted by Herrera et al. (2013) and Bianchini et al. (2013) [39,40]: (1) When β approaches 90°, C approaches 0, and V_SLOPE tends to infinity; therefore, V_SLOPE cannot be set higher than 3.33 times V_LOS. The threshold value of C is set to 0.3; when −0.3 < C < 0, C = −0.3; when 0 < C < 0.3, C = 0.3. (2) The target point with a positive value of V_SLOPE is excluded. A positive value of the target point indicates that the slope body moves upward along the direction of the steepest slope, but this is not consistent with reality.

Based on the above equations, the V_SLOPE deformation results of the two kinds of data in the ascending and descending orbits were obtained. Due to the different sensitivities of different satellite views to landslides of different orientations, the deformation information is missing from the individual data points, which in turn led to some differences in the V_SLOPE deformation results obtained from different data points. Therefore, to better characterize the deformation conditions in the study area within the studied time frame, satellite data from different angles obtained within similar study periods were combined. In this study, the V_SLOPE deformation rate level combination matrix for the ascending and descending orbit data was introduced, as shown in

Table 2. V_SLOPE (acquisition by ascending and descending data) deformation rate level combination matrix.

	V4	V3	V2	V1
V4	4	4	4	4
V3	4	3	3	3
V2	4	3	2	2
V1	4	3	2	1

; this matrix reflects all the deformation information contained in the two different datasets. The classification scheme of the deformation rate refers to the classification rules used in Zhou et al. (2022) [38]: V1 (0–2 mm/mouth), V2 (2–4 mm/mouth), V3 (4–6 mm/mouth), and V4 (>6 mm/mouth).

3.4. Hazard Assessment Methods Considering the Deformation Factor

Conventional hazard assessments have certain shortcomings. If a landslide has been managed and its deformation is thus effectively controlled, its hazard assessment level should be reduced; at the same time, if human activities are abnormally active and are sufficient to destroy the original equilibrium of the slope, it will increase the probability of slope instability, and the hazard assessment level should thus be increased. Therefore, the slope deformation process is always undergoing dynamic changes. If we want to reduce the landslide hazard, pinpoint high-risk areas, and provide authorities with a reasonable control plan, the current deformation states of landslides must be considered. The matrix method used by Zhou et al. (2022) [38] was referenced in this study to combine the deformation rate level of the combined V_SLOPE (acquisition by both ascending and descending data) with the hazard assessment level, and the assessment matrix is shown in

Table 3. V_SLOPE deformation rate level combined with the landslide hazard assessment level matrix.

	V4	V3	V2	V1
H4	4	4	3	2
H3	4	3	2	2
H2	4	3	2	1
H1	4	3	1	1

.

4. Dynamic Landslide Hazard Assessment in Zhouqu

4.1. Landslide Inventory Map

Stacking-InSAR-based technology was used to process the ascending and descending Sentinel-1A data for the two-year period from 2020 to 2022. The visual interpretation results of potential landslides based on optical images and the field verification are shown in Figure 4 and Figure 5. Figure 4 shows the landslide detection results obtained based on the ascending data. Fifty-two landslides were detected, and their distribution is shown in Figure 4a; these landslides were mainly concentrated in the northern and southeastern areas of Zhouqu. Figure 4b shows the optical images, boundary extent, and Stacking-InSAR results of four typical landslides with very distinct deformation streaks. After on-site verification, the front edge of landslide L16 (Suoertou) was found to be inhabited by a wide range of residents, and the road and wall showed a large numbers of cracks (Figure 4c). Landslide L24 showed partial collapse at the front edge and obvious cracks at the back edge, with obvious signs of deformation (Figure 4d). Figure 5 shows the landslide detection results based on the descending data. A total of 23 landslides were detected, of which 6 landslides were detected from the ascending data. The distribution of these landslides is shown in Figure 5a; they are mainly concentrated in northeastern Zhouqu. A typical landslide optical image, boundary extent, and InSAR results are shown in Figure 5b. The site photos are shown in Figure 5c,d. The locations where the photos were taken are shown in Figure 5b. The photos of landslides L14 (Xieliupo) and L16 (Mentouping) detected from the descending data exhibit obvious signs of deformation and landslide boundaries, and some local landslides are distributed on the surfaces of these landslides. Large tension cracks exist on the back edge of landslide L14, and multiple gully distributions also exist on the front edge of landslide L16. The typical landslides described above are all historically documented [42].

The landslides detected by InSAR technology were supplemented with the historical landslide database to obtain the latest landslide inventory map of Zhouqu. This map includes a total of 106 landslides with a total area of approximately 35 km²; the distribution of these landslides is shown in Figure 6. The yellow borders are landslides detected from the ascending SAR data (A), the red borders are landslides detected from the descending SAR data (D), and the purple borders are landslides in the historical landslide database not detected by InSAR technology (H). The distribution of the landslide area within each village in Zhouqu is shown in Figure 7, among which Pingding, Chengguan, and Dongshan landslides accounted for the highest percentages and Baleng, Bacang, Boyu, and Changang landslides in four villages were less distributed; other landslides accounted for a relatively high average percentage.

4.2. Landslide Hazard Assessment

4.2.1. Susceptibility Assessment

The selection of landslide susceptibility assessment factors is diversified, and appropriate factors can be selected according to historical landslide survey data or relevant literature in the study area. We can also select as many factors related to landslide occurrence as possible according to literature references, and judge and retain more appropriate factors for this study in the process of landslide susceptibility assessment. In this paper, we have referred to several historical related studies in the study area [43,54], and selected 9 factors, including elevation, slope, aspect, distance from river, distance from fault, distance from road, land cover [55], lithology, and rainfall to assess the susceptibility of Zhouqu, among which the elevation data are obtained from the SRTM 30-m DEM, the slope and aspect data are obtained from the DEM, and the river and road data are obtained from the National Geographic Information Resources Catalog Service System (www.webmap.cn accessed on 5 October 2022) (

Table 4. Summary of data in this study.

Data	Source
Historical landslides	http://www.resdc.cn accessed on 5 October 2022
Sentinel-1A	https://search.asf.alaska.edu/#/ accessed on 5 March 2022
Elevation	SRTM 30 m (http://gdex.cr.usgs.gov/gdex/ accessed on 5 October 2022)
Slope	Generated by DEM
Aspect	Generated by DEM
River	www.webmap.cn accessed on 5 October 2022
Fault	Ref [43]
Road	www.webmap.cn accessed on 5 October 2022
Lithology	Ref. [43]
Land cover	Ref. [55]
Precipitation	Ref. [43]

). These factors are all highly relevant to landslide generation. Geomorphological factors include elevation, slope, aspect, and land cover, geological factors include lithology and faults, environmental factors include rivers and rainfall, and human activities include roads [56]. Certainty factors were used to calculate the degree of correlation (CF value) between each classification level and landslides in Zhouqu for the 9 factors. Based on the landslide dimensions in Zhouqu, all raster data in the calculation process were resampled to a 30-m resolution. The classification of each factor, the area of each class of landslide, and the CF value results are shown in Figure 2 and Figure 8 and

Table 5. CF values obtained for different classes of the seven factors influencing landslides in Zhouqu.

Factors	Classification	Category Area (km2)	Landslides Area (km2)	CF
Elevation	<1500 m	73.4643	3.5352	0.76988
	1500–2000 m	409.806	15.723	0.70842
	2000–2500 m	835.794	14.265	0.33002
	2500–3000 m	989.172	1.1637	−0.89876
	3000–3500 m	527.346	0	−1.00000
	>3500 m	180.735	0	−1.00000
Slope	<15°	245.678	7.2675	0.61836
	15–30°	914.25	15.6267	0.33101
	30–45°	1438.03	10.4841	−0.36871
	45–60°	385.91	1.2636	−0.71762
	60–75°	31.842	0.045	−0.87835
	>75°	0.6075	0	−1.00000
Aspect	Plane	0.1107	0.0018	0.29617
	North	401.039	1.9134	−0.58792
	Northeast	435.179	4.2876	−0.14467
	East	419.045	8.5392	0.44074
	Southeast	343.301	7.1505	0.45310
	South	376.195	4.6467	0.06979
	Southwest	374.297	3.7287	−0.13508
	West	354.901	2.7738	−0.32288
	Northwest	312.249	1.6452	−0.54470
Distance from river	<200 m	360.988	5.2317	0.20892
	200–400 m	336.267	6.3063	0.39131
	400–600 m	319.784	5.535	0.33951
	600–800 m	304.063	4.3821	0.20441
	800–1000 m	286.675	3.3129	0.00495
	>1000 m	1408.54	9.9189	−0.39039
Distance from fault	<500 m	459.375	10.3554	0.49556
	500–1000 m	369.251	4.8051	0.11765
	1000–1500 m	307.184	2.6424	−0.25417
	1500–2000 m	258.879	3.9915	0.25711
	2000–2500 m	217.983	2.6874	0.06800
	>2500 m	1403.65	10.2051	−0.37047
Distance from road	<200 m	99.5868	2.736	0.58819
	200–400 m	90.7038	2.9367	0.65232
	400–600 m	86.9958	2.8377	0.65498
	600–800 m	84.2985	2.6154	0.63667
	800–1000 m	82.3068	2.7306	0.66097
	>1000 m	2572.43	20.8305	−0.29826
Lithology	Semi-hard limestone	639.29	6.6528	−0.09606
	Sandstone	1213.58	3.3966	−0.75874
	Powdery clay	35.0595	0.0918	−0.77434
	Slate	245.839	11.3859	0.76045
	Granite	10.6371	0.045	−0.63481
	Hard limestone	142.073	0.8721	−0.46909
	clastic rock	648.323	12.2382	0.65383
	Conglomerate	81.5085	0	−1
Land cover	Cropland	0.0603	0	−1
	Forest	288.133	16.0047	0.80220
	Shrub	2047.54	3.528	−0.85163
	Grassland	16.6779	0.1017	−0.47262
	Water	654.945	15.0147	0.50418
	Snow/Ice	2.6793	0.0009	−0.97112
	Barren	0.0117	0	−1
	Impervious	1.8612	0	−1
	Wetland	4.4055	0.0369	−0.27394
Precipitation	<500 mm	457.53	12.474	0.58493
	500–600 mm	673.176	13.3902	0.42677
	600–700 mm	732.11	7.3233	−0.13147
	700–800 mm	578.009	1.2105	−0.8196
	800–900 mm	575.486	0.2844	−0.9575

, respectively.

The covariance among the landslide factors can affect the performance of the evaluation model, and it is thus necessary to check whether covariance problems exist among the selected factors before performing neural network modeling. The tolerance and variance inflation factor (VIF) are two commonly used covariance detectors. When the VIF is greater than or equal to 5 or the tolerance is less than or equal to 0.2, the indicator is considered to exhibit covariance. As shown in

Table 6. Covariance test results of the landslide factors.

Factors	Tolerances	VIF
Elevation	0.416	2.406
Slope	0.863	1.158
Aspect	0.848	1.179
Distance from river	0.731	1.368
Distance from fault	0.816	1.226
Distance from road	0.766	1.305
Lithology	0.700	1.428
Land cover	0.581	1.722
Precipitation	0.483	2.071

, the minimum tolerance and maximum VIF are 0.416 and 2.406, respectively; thus, there is no covariance between the factors, and both can participate in the model calculation.

In this paper, the CF value of each factor evaluation level was used as the input data of the model, and the susceptibility index (landslide: 1, nonlandslide: 0) was used as the output terminal of the model to perform susceptibility modeling. The landslide transcription data were randomly divided into two parts; 70% of the landslide sample units were used for model training, and the remaining 30% were used to verify the model performance. To avoid an imbalance in the model training process, the same amounts of data were randomly selected for the neural network analysis between the nonlandslide samples and landslide samples for susceptibility modeling. The results of the conventional hazard assessment are shown in Figure 9, and the natural interruption method was used to classify the susceptibility level. Among them, the high-susceptibility area accounts for 14.83% of the area of Zhouqu, and these regions are concentrated mainly in the northern and middle-eastern parts of Zhouqu. The medium-susceptibility area accounted for 10.37%, the low-susceptibility area accounted for 17.25%, and the very-low-susceptibility area accounted for 57.55%. The ROC curve (Figure 10) was used for the accuracy assessment; the area under the curve (AUC) was 0.913 > 0.5, and the accuracy of the training results was thus superior.

4.2.2. Hazard Assessment

In Zhouqu, the amount of rainfall directly affects the frequency of geological disasters. According to the collected rainfall data, May–September is the rainfall-prone period in Zhouqu and is also the high season for landslide disasters. Therefore, the large-scale and wide-impact landslide geological hazards in the past 5 years were collected (

Table 7. Occurrence times and corresponding rainfall totals of relatively large-scale landslides that occurred in Zhouqu in the past 5 years.

Landslide Name	Date of Landslide	Rainfall of the Month (mm)
Jiangdingya	12 July 2018	145.8
Hena	July 2018	145.8
Yahuokou	19 July 2019	153.92
Beishangou	17 August 2020	296.42

), among which the smallest monthly rainfall found to induce landslides was 145.8 mm and the largest was 296.42 mm. By counting the monthly rainfall totals in the study area from May to September in 2018–2022, the number of occurrences of rainfall below 145.8 mm was found to be 5, and the probability of occurrence was 0.2; therefore, the temporal probability P(t) value in the landslide hazard calculation formula was set to 0.2. The spatial probability P(s) is the susceptibility value, and the event probability P(s) is the ratio of the landslide area to the total landslide area in each susceptibility classification. The risk assessment results of Zhouqu are shown in Figure 11; as the figure shows, the high-risk area accounted for 6.60% of the area of Zhouqu, mainly concentrated in the Zhouqu County area and some northern and central-eastern areas. The medium-susceptibility area accounted for 5.08%, the low-susceptibility area accounted for 3.81%, and the very-low-susceptibility area accounted for 84.50%.

4.3. Landslide Hazard Assessment Results Considering the Deformation Factor

The deformation rate results of the SBAS-InSAR processing method based on Stacking-InSAR are shown in Figure 12. The red points (negative values) indicate that the target is moving away from the satellite along the LOS direction; the blue points (positive values) indicate that the target is moving toward the satellite along the LOS direction; the green points indicate relative stability. Figure 12a shows the SBAS-InSAR results obtained based on the ascending SAR data; the annual deformation rate reaches maximum values of −148 mm/yr (negative) and 68 mm/yr (positive). Figure 12b shows the SBAS-InSAR results obtained based on the descending SAR data, with the annual deformation rates reaching up to −102 mm/yr (negative) and 107 mm/yr (positive).

The deformation results obtained in the LOS direction of the ascending and descending SAR data were projected to the deformation results obtained in the steepest slope direction according to the method described in Section 3.3; the results are shown in Figure 13, and both data volumes were reduced by approximately 10%. Different incidence angles have different sensitivities to different slope orientations, so the two V_SLOPE results show partial differences, but the main deformation zones are roughly the same, and the deformation classes of the deformation zones are also basically the same. Figure 13c shows the ratio of the number of data points corresponding to each deformation rate in the V_SLOPE results of the ascending and descending SAR data, and the ratios of the different deformation rate levels to the total number are approximately the same, thus laterally confirming the accuracy of the InSAR and V_SLOPE results.

Based on the use of the ordinary kriging method to interpolate the V_SLOPE results obtained from ascending and descending data, the deformation rate levels in Zhouqu (Figure 14) were obtained according to the combination matrix of deformation grades shown in Table 2. The deformation rate levels were divided into four levels, high, medium, low, and very low, among which the high-deformation-rate level accounted for 0.17%, the medium-deformation-rate level accounted for 2.09%, the low-deformation-rate level accounted for 33.88%, and the very-low-deformation-rate level accounted for 63.86%. The deformation rate classification results characterized the deformation conditions during the studied timeframe in Zhouqu.

In Table 3, the results were combined with the matrix to obtain a hazard assessment level map of landslides in Zhouqu while considering the deformation factor (Figure 15): high-hazard areas spanned 0.27%, medium-hazard areas covered 2.98%, low-hazard areas covered 11.04%, and very-low-hazard areas covered 85.71%. Among them, the medium- and high-risk areas were mainly located around Zhouqu County and in the central-eastern region. Compared to the conventional hazard assessment results, the high-hazard areas indicated by the hazard assessment results considering the deformation factor were reduced, and the current areas with relatively large deformation and high hazard degrees are retained. One of these areas, containing the Jiangdingya landslide located near Nanyu that occurred in July 2018 and caused great economic losses, experienced significant deformation before the landslide occurred [42] and was indicated as a high-hazard area in the conventional hazard assessment results. However, this landslide was well-managed after its occurrence, as the deformation has been stable within the current studied time frame, and the hazard assessment level considering the deformation factor was assessed as a low hazard assessment level, which is more suitable for the actual situation. Moreover, according to the hazard assessment level deformation factor, the high-risk areas can be precisely located, and the corresponding landslides can be targeted directly according to those high-hazard areas. As shown in Figure 16 (the area in the red border is shown in Figure 15), the results obtained from the ascending SAR data detected based on SBAS-InSAR indicated a maximum deformation rate of −148 mm/yr. A few randomly selected feature points are plotted in the corresponding time-series curves. The results show that both landslides (the Suoertou and Zhongpaicun landslides) were in a continuous deformation phase during the studied time frame, and their accumulated deformation reached a maximum of 350 mm. By providing recent hazard assessment results to the relevant management organizations, appropriate mitigation measures could be taken for these high-risk areas. It is also possible to directly manage landslides in high-risk areas and control their deformation rates to reduce their danger levels. This assessment method saves resources, improves overall efficiency, and is more conducive to protecting people’s lives and property.

Landslide evolution is a dynamic change process, and the landslide danger may differ among different regions. The time-series InSAR technique can effectively monitor the deformation rates of surface slopes in different states and can track the regional slope deformation state, thus providing key basic data and the possibility to dynamically evaluate landslide hazards.

5. Discussion

The current study faced some limitations. InSAR technology can obtain ground-deformation information in the studied time range, thus providing key basic data for the dynamic tracking of landslide hazards. In this study, the deformation information obtained by InSAR was involved in the hazard assessment calculation; we ensured the integrity of the deformation information as much as possible, but this information was inevitably affected by complex atmospheric disturbances, the topographic relief, and the vegetation in the mountainous area, causing some errors to arise in the deformation rate level distribution of the study area and affecting the accuracy of the hazard assessment results. This requires enhancing the capability of InSAR technology processing, ensuring less interference from atmospheric and other errors to reduce the on-site verification workload, and improving the overall efficiency.

The hazard assessment method introduced in this paper can also realize dynamic landslide hazard assessments. Based on InSAR technology, the ground deformation rates at different times of the same year can be obtained and combined with the landslide hazard assessment results to obtain landslide hazard assessment results characterizing different periods; then, corresponding preventive and control measures can be taken according to the identified high-risk areas in different periods. However, due to the different vegetation coverage conditions in different periods, the coherence gaps in the interferograms obtained by InSAR processing are large, causing large errors when comparing the processing results. Therefore, using SAR satellite data with different wavelengths and incidence angles for processing can increase the accuracy of the results.

Landslide risk evaluations are the basis of risk management, and reliable risk evaluation results can provide more scientific and reasonable technical support for the formulation of disaster management departments. The method proposed in this paper only performs regional hazard assessment without considering the vulnerability of people and buildings; there may not be any people living or any infrastructure distributed in the identified high-risk areas, making the hazard assessment results much less useful. Therefore, landslide risk assessments should be carried out on the basis of the assessment methods introduced in this paper while also taking into account the vulnerability of the building distribution and human activities in the study area to serve the disaster management authorities better.

6. Conclusions

This paper proposes a new method for landslide hazard assessment by combining the ground deformation rate levels obtained from different orbit SAR data with the conventional landslide hazard assessment. The method can obtain the distribution of landslide hazard levels in the study area over a period of time and revise the results obtained with conventional landslide hazard assessment methods. The novel method takes into account the deformation state presented after landslide management and maintenance, as well as the revival of ancient landslides and the generation of new landslides caused by human activities or other factors, thus making the landslide hazard assessment results more accurate and valuable. This paper takes Zhouqu County as a case study, and the main results are summarized as follows:

(1)

First, a landslide detection was carried out in the whole territory of Zhouqu by Stacking-InSAR technology. Through an on-site verification, we confirmed the 52 landslides detected from the ascending SAR data and 23 landslides detected from the descending SAR data, among which 6 landslides were jointly detected. Complementing these detected landslides with the historical landslide database, we formed a landslide inventory map including 106 landslides.

(2)

The landslide susceptibility assessment of Zhouqu was carried out by using neural network modeling on the basis of a landslide inventory map. The landslide hazard assessment results of Zhouqu were obtained after considering the event probability and temporal probability of the triggering factor (rainfall); in the results, the high-risk area accounted for 6.60% of the area of Zhouqu.

(3)

Based on Stacking-InSAR, the deformation rate results along the LOS direction for the whole Zhouqu area were obtained by using the SBAS-InSAR technique. To compare the deformation conditions of different slope directions and enrich the deformation information, the deformation rate in the LOS direction was projected to the steepest slope direction, and the results obtained from the ascending and descending data were combined with the matrix-based method to form the deformation rate level results of Zhouqu. The landslide hazard assessment levels considering the deformation factor of Zhouqu were obtained by combining those results with the landslide hazard assessment levels. Among the results, the medium-high-hazard areas accounted for 3.25% of the study area and had locations directly corresponding to the high-hazard Suoertou and Zhongpaicun landslides.

This landslide hazard assessment method is applicable to any country or region in the world that needs landslide prevention and control, which is more convenient for the authorities to target in high-risk areas and landslides, improve governance efficiency, and ensure the safety of people’s property. When using this method for landslide hazard assessment, it is necessary to select SAR satellite data with appropriate wavelength and resolution according to the historical landslide size and vegetation coverage in the study area, which depends on whether the area is covered by SAR satellites with different bands, and whether the covered SAR satellites contain different orbits, which will affect the accuracy of the ground deformation rate level. At the same time, when conducting conventional landslide hazard assessment, it is also necessary to select appropriate assessment factors according to the geological and geomorphic conditions of the study area. At present, the dynamic evaluation of landslide hazard can be realized through InSAR technology. What needs to be done in the future is to obtain the landslide hazard evaluation results of different regions in real time, dynamically and automatically, and realize the automatic alarm function under the premise of launching SAR satellites of different bands and various revisit periods. This is undoubtedly an innovation which can further reduce casualties and realize a better life for mankind.