1 Introduction

Landslides are one of the most hazardous natural events, inflicting massive casualties, injuries, and property loss (Garcia-Delgado et al. 2022). Landslides have killed over 50,000 people in the twenty-first century (Mizutori 2018), and the global economic loss caused by landslides amounts to billions of dollars per year (Kjekstad and Highland 2009; Ha et al. 2020). Rainfall-triggered shallow landslides are unstable soil masses with small volumes (ranging from 101 to 105 m3) on the slope (Liu et al. 2022b). They are densely distributed across upstream catchments, contributing plenty of deposits to the channel (Bordoni et al. 2021), evolving into catastrophic debris flows and causing considerable damage to farmlands and civil engineering infrastructures, and, in some cases, causing massive casualties (Pradhan et al. 2019; Grima et al. 2020). The ability of susceptibility assessment methodologies to predict the occurrence of shallow landslides or identify regions highly susceptible to failure is beneficial for land use planning, infrastructure development, community resilience, and disaster risk reduction, among others (Dou et al. 2019; Camera et al. 2021).

Empirical, data-driven, and deterministic techniques have been commonly employed in previous studies to predict shallow landslide susceptibility (Ghosh et al. 2011; Daneshvar 2014; Feng et al. 2020; Bordoni et al. 2021). Among the various methods, expert empirical methods are developed by incorporating experience and judgment to determine the spatial limits of a hazard (Hess et al. 2017). These methods, however, have faced criticism due to their heavy dependence on expert experience and subjective judgment, resulting in the potential for controversial results (Nguyen et al. 2020a). The use of data-driven approaches has gained popularity in landslide susceptibility research owing to the advancements in data-mining algorithms and GIS tools. These approaches involve establishing relationships between past landslide occurrences and the controlling factors of landslides (Zêzere et al. 2017; Lin et al. 2021). In contrast to expert empirical approaches, data-driven methods provide a more objective and quantitative approach (Yong et al. 2022). Furthermore, deterministic approaches are commonly employed to investigate the mechanisms of landslide formation and estimate landslide susceptibility (Liu et al. 2016; Murgia et al. 2022). These physical models are mainly divided into two categories—pseudo-static analysis methods based on infinite slope theory (Pack et al. 2005; Keles and Nefeslioglu 2021) and finite element methods (Vandromme et al. 2015; Zhu et al. 2022; Qin et al. 2023). The former focus on the ratio of forces resisting slope destabilization to forces driving slip formation (Arnone et al. 2016) and are easy to apply and extend to regional scales with greater universality. The latter are a refined approach that can better simulate slope stability but are limited in universality due to high parameter requirements and calibration (Baum et al. 2008; Keles and Nefeslioglu 2021; Murgia et al. 2022; Yong et al. 2022). Despite the limitations, studies have continued to explore physical models since doing so opens up new opportunities and research trajectories (Bordoloi and Ng 2020; Murgia et al. 2022).

The impact of vegetation on slope stability of potential shallow landslide area carries substantial practical engineering implications (Liu et al. 2016; Murgia et al. 2022). Its research provides scientific support for assessing landslide risks, designing engineering structures, restoring ecosystems, and managing land resources (Tron et al. 2014; Feng et al. 2020). Consequently, researchers have conducted extensive research on how vegetation influences shallow landslide susceptibility (Hao et al. 2020; Hao et al. 2022; Liu et al. 2022b; Murgia et al. 2022). The effects of vegetation on slope stability have been the subject of numerous studies. These effects can be broadly categorized into two categories: positive effects of vegetation provided by root reinforcement (Bischetti et al. 2005; Arnone et al. 2016), root–water uptake (Zhu and Zhang 2015, 2019), and vegetation interception (Kozak et al. 2007; Bordoloi and Ng 2020), and negative effects of vegetation provided by root preferential flow (Shao et al. 2015; Wang et al. 2022), vegetation weight loading (Reichenbach et al. 2018; Zou et al. 2021), and wind loading (Wang et al. 2019; Zhuang et al. 2022). The major source of the root infiltration effect is the plant roots below the ground surface, which act as a conduit for precipitation infiltration and soil evaporation and hence promote soil infiltration (Wang et al. 2020; Liu et al. 2022a). According to Jotisankasa and Sirirattanachat (2017), increasing the quantity of root content in the soil leads to increased soil permeability. Qin et al. (2022) also confirmed this trend via experiments, revealing that elevating the root area density (RAD) from less than 0.05 mm2/cm2 to higher than 0.05 mm2/cm2 increases the saturated permeability (Ks) of root–soil mass by up to 4−15 times. They discovered that when a landslide is induced by rainfall, the hydraulic conductivity of the root–soil mass above the slip surface is 14−17 times greater than that of the soil below. This leads to accelerated rainfall infiltration and the formation of high pore-water pressure areas at the ends of root channels. Moreover, root–water uptake during transpiration can induce matric suction, contributing to improved slope stability and can be retained after a certain duration of rainfall, as demonstrated in studies by Ng et al. (2014), Zhu and Zhang (2015, 2019). According to the study conducted by Kim et al. (2013), even though the proportion of interception loss in the total observed rainfall is relatively small during an extreme rainfall event that causes landslides, it could still influence the onset of landslides through changes in the temporal distribution of rainfall. Furthermore, wind loads transferred through vegetation trunks can also influence the stability of slopes, indirectly applied to the soil mass (Wang et al. 2013; Zhuang et al. 2022). Consequently, vegetation has both positive and negative effects on slope stability, posing a significant challenge in modeling the stability of shallow slopes when considering these factors in an integrated manner (Bordoloi and Ng 2020; Liu et al. 2022b; Murgia et al. 2022).

Extensive research conducted over decades has focused on the effects of vegetation on slope failure, which can be quantified using numerical models (Wu and Sidle 1995; Schwarz et al. 2010; Liu et al. 2016; Murgia et al. 2022). The reinforcement effect of the root system, which enhances resistance to shear stress, can be quantified by additional soil cohesion (Wu and Sidle 1995). Chirico et al. (2013) and Schwarz et al. (2010) used the Wu model (WM) and the fiber bundle model (FBM) to quantify root reinforcement. Vegetation weight loads are often quantified using stand volume, as exemplified by Zou et al. (2021), who provided a model incorporating the reinforcement effect of the root system and vegetation weight loading through the finite element strength reduction method. Despite extensive discussions on the influence of vegetation on slope stability, there is still a research gap in identifying the positive and negative impacts of vegetation on a regional scale in transition zones with extremely rugged terrain.

This article introduces a novel physical-based methodology that focuses on the physical processes associated with vegetation to estimate shallow landslide susceptibility at a regional scale. This method encompasses the positive impacts of vegetation, such as rainfall interception and root system reinforcement, and the negative effects such as vegetation weight loads and preferential flow of root systems. Ultimately, the application of this method resulted into the development of a shallow landslide susceptibility map for a typical extreme topographic transition belt, serving as a tool for disaster prevention and mitigation in the region.

2 Study Area

The Dadu River Basin is located in the western part of Sichuan Province, China. Its latitudes range from 28° 13′ N to 33°17′ N, and its longitudes range from 100° 05′ E to 103° 46′ E, as shown in Fig. 1. The area is characterized by steep slopes, towering mountains, and valleys that form a V shape. The region, situated between the Songpan fold system and the Yangzi platform, is dominated by active neotectonic movement and complex geological structures (Zhang 2013). These conditions contribute to the region’s complex lithological change from the Sinian to the Quaternary (Zou et al. 2021). The region exhibits two distinct climates in the north and south, namely, a dry mountainous plateau climate characterized by sparse rainfall, and a humid and rainy monsoon climate (Jiang et al. 2022). Due to the impacts of active tectonics, enormous topographic relief, and complicated hydrological conditions, the region is also a typical meeting place for landslides (Xiong et al. 2020; Qin et al. 2023) (Fig. 1). According to incomplete records, there are over 1000 mountain disasters of different types and magnitudes from 2000 to 2020, resulting in direct economic damages of tens of millions of dollars and loss of lives (Ding et al. 2007; Ba et al. 2011). This region also exhibits significant vertical and horizontal differentiation of vegetation due to its diverse climatic and hydrological conditions as well as its steep terrain (Li et al. 2020; Zou et al. 2021). With elevations from 275 to 7445 m, the area has a diverse range of vegetation types, including broad-leaved forest, coniferous forest, shrub, meadow, grassland, and alpine vegetation. Consequently, it presents an ideal opportunity for exploring the effects of vegetation on shallow landslides at a regional scale.

Fig. 1
figure 1

The study area and spatial distribution of historical landslides in the Dadu River Basin. The historical landslides were identified in this study. Photographs by Hu Jiang, in September 2020 and April 2022

Full size image

3 Materials and Methods

This article outlines a physical-based approach to predicting the spatial likelihood of landslides at a regional scale. The methodology centers on the failure mechanism of vegetated slopes and explores the influence of plants on shallow landslides. Specifically, it considers the effects of vegetation canopy interception, vegetation weight loads, the improvement of vegetation root systems for drainage, and the reinforcement of vegetation root systems. The methodology is divided into four main components: (1) parameter quantification; (2) model construction; (3) model evaluation; and (4) landslide susceptibility mapping. Figure 2 presents the flowchart illustrating this process.

Fig. 2
figure 2

Flowchart of the proposed method

Full size image

3.1 Inventory of Landslides

Field investigation and remote sensing image interpretation were employed to gather landslide samples. These samples help determine the degree of enhancement of the model when considering the influence of vegetation on the occurrence of shallow landslides. Due to the scope of the research region, remote sensing image interpretation was used as a supplementary method to identify the spatial distribution of landslides that are inaccessible and located far from the main road network. Some interpretation markers (Aksoy and Ercanoglu 2012; Hungr et al. 2014), such as image features (shape, hue, and texture) and basic landslide features (landslide body, landslide wall, and landslide boundary), were fully taken into account to ensure data reliability. The image interpretation results, which form part of the shallow landslide inventory, require that the area of a single landslide is less than 10,000 m2 (Shao et al. 2022; Hwang et al. 2023). Finally, a total of 2334 landslides caused by different rainfall events were identified for model performance assessment (Fig. 1). The validation dataset is made of past landslide events, with the premise that these landslide spots and their surrounding regions are likely to be areas prone to shallow landslides owing to the spatial autocorrelation of the geographical environment (Dou et al. 2015; Tiranti et al. 2019).

3.2 Quantifying Landslide Susceptibility Considering the Positive and Negative Impacts of Vegetation on Landslides Occurrence

Many studies have shown that vegetation exerts a twofold role in shallow landslides occurring (Liu et al. 2022b): (1) On the positive side, vegetation improves soil shear strength through the root system reinforcement, which is conducive to maintaining slope stability (Kim et al. 2013); (2) On the negative side, vegetation weight loads increase total slope weight and sliding force, and the existence of root system, especially for the coarse roots, significantly changes soil porosity and provides preferential flow for water infiltration, which are detrimental to slope stability (Wang et al. 2020). This section presents a method for establishing a regional-scale landslide susceptibility assessment that depends on the physical formation mechanism of vegetated slope failure triggered by rainfall. This method integrates multiple models, namely the infinite slope model, vegetation canopy interception model, root system reinforcing model, vegetation weight load model, and subsurface flow model, to evaluate the impacts of plants on shallow landslides. Table 1 provides an overview of the sources and spatial resolution of data used in this study. To maintain the consistency of spatial resolution of different parameters and the suitability of the precise geomorphic characteristics described in the Digital Elevation Model (DEM) for evaluating shallow landslide susceptibility, the collected geodata were downscaled to 30 m spatial resolution with the raster processing tools of geographic information system (Jiang et al. 2022).

Table 1 Sources and spatial resolution of the datasets used in this study
Full size table

3.2.1 Quantifying the Parameters Related to the Positive Effects of Vegetation

The positive impacts of plants on slope stability are categorized (Liu et al. 2022b; Murgia et al. 2022) as (1) above-surface shading effects provided by the tree canopy, which reduces splash and infiltration of rainfall on the surface or soil mass; and (2) under-surface reinforcing effects provided by plant root systems, which affects soil shear strength. This study quantified the positive effects of plants on slope stability through the interception of rainfall by vegetation and the reinforcement of root systems, employing physical–mechanical and hydrological mechanisms.

3.2.1.1 Quantifying the Reinforcement Impacts of the Plant Root System Using Extra Root–Soil Mass Cohesion (C root)

The reinforcement effect provided by the plant root systems is a significant force contributing to the resistance against slope failure for vegetated slopes (Feng et al. 2020). The extra cohesion (Croot) offered by the plant root systems can be calculated using Eqs. 1−3 developed and refined by Wu et al. (1979).

$$ C_{{{\text{root}}}} \cong k^{\prime}T_{{{\text{root}}}} {\text{RAR}}, $$
(1)
$$ k^{\prime} = \cos \theta \tan \phi + \sin \theta , $$
(2)
$$ {\text{RAR}} = A_{{\text{r}}} /A_{{{\text{rs}}}} , $$
(3)

where Croot is the supplementary cohesion that represents reinforcing impacts of the plant root systems, and its unit is kPa; k′ can be considered as the correction coefficient and commonly specified as 1.2 (Wu et al. 1979); Troot denotes mean tensile strength per unit root cross-sectional area; θ and \(\phi\) are the shear distortion and internal friction angle, respectively; RAR denotes root area ratio; Ar denotes root cross-sectional area; Ars denotes rooted soil area. The parameters are shown in Table 2.

Table 2 Average tensile strength per average root cross-sectional area of different vegetation
Full size table
3.2.1.2 Quantifying the Effective Rainfall Reaching the Ground After Considering the Vegetation Canopy Interception and Evapotranspiration

The effective rainfall reaching the ground can be calculated (Terink et al. 2015) by subtracting evapotranspiration and vegetation canopy interception from total rainfall. The reference evapotranspiration (RETr) can be quantified with Eq. 4—the refined Hargreaves’s formula (Hargreaves and Samani 1985; Allen et al. 1998).

$$ {\text{RET}}_{{\text{r}}} = 0.0023 \times 0.408 \times {\text{Ra}}(T_{{{\text{avg}}}} + 17.8) \times {\text{TD}}^{0.5} , $$
(4)

where Ra represents extraterrestrial radiation calculated using those equations taking the day of the year and the latitude of the study region (Allen et al. 1998) as input parameters; Tavg represents average daily air temperature (Fig. 3); TD represents daily temperature range (Fig. 4).

Fig. 3
figure 3

Average daily air temperature (Tavg) in the Dadu River Basin

Full size image
Fig. 4
figure 4

Daily temperature range (TD) in the Dadu River Basin

Full size image

Canopy storage, or rainfall interception, can be computed using Eq. 5 developed by Kozak et al. (2007) and Terink et al. (2015).

$$ {\text{Scan}}_{\max } = 0.935 + 0.498{\text{LAI}} - 0.00575{\text{LAI}}^{2} , $$
(5)

where Scanmax denotes maximal canopy storage, which varies with plant type; LAI denotes the leaf-area index and can be calculated using Eq. 6 (Terink et al. 2015).

$$ {\text{LAI}} = {\text{LAI}}_{\max } \cdot \frac{{\log (1 - {\text{FPAR}})}}{{\log (1 - {\text{FPAR}}_{\max } )}}, $$
(6)

where LAI denotes the leaf-area index, and LAImax denotes the maximum leaf-area index determined by vegetation type (Sellers et al. 1996; Terink et al. 2015), as shown in Table 3. The fraction of photosynthetically active radiation (FPAR) can be computed using Eq. 7 (Sellers et al. 1996; Peng et al. 2012).

$$ {\text{FPAR}} = \min \left( {\frac{{({\text{SR}} - {\text{SR}}_{\min } )({\text{FPAR}}_{\max } - {\text{FPAR}}_{\min } )}}{{({\text{SR}}_{\max } - {\text{SR}}_{\min } )}} + {\text{FPAR}}_{\min } ,\;0.95} \right) $$
(7)

with

$$ {\text{SR}} = \frac{{1 + {\text{NDVI}}}}{{1 - {\text{NDVI}}}}, $$

where FPARmax and FPARmin are equal to 0.95 and 0.001, respectively. NDVI is the normalized difference vegetation index (Fig. 5).

Table 3 Maximum leaf-area index (LAImax) for different vegetation types
Full size table
Fig. 5
figure 5

Normalized difference vegetation index (NDVI) in the Dadu River Basin

Full size image

Then, the effective rainfall can be determined by integrating Eqs. 8 and 9.

$$ {\text{Scan}}_{t} = {\text{Scan}}_{t - 1} + P_{t} , $$
(8)
$$ {\text{Pe}}_{t} = \max (0,\;{\text{Scan}}_{t} - {\text{Scan}}_{{\max ,{\text{t}}}} ), $$
(9)

where Scant and Scant−1 denote the canopy storage at time t and t−1, respectively. Pt is the total rainfall at time t, and Pet is the effective precipitation after vegetation interception at time t, namely the throughfall.

After determining the effective rainfall, a portion of the canopy storage can be deducted due to atmospheric evaporation, and the canopy storage is refreshed according to Eqs. 10 and 11.

$$ {\text{Scan}}_{t} = {\text{Scan}}_{t} - {\text{Pe}}_{t} , $$
(10)
$$ {\text{Scan}}_{t} = {\text{Scan}}_{t} - \min (1.5{\text{RET}}_{{{\text{r}},{\text{t}}}} ,\;{\text{Scan}}_{t} ), $$
(11)

where RETr,t is the reference evapotranspiration. For a comprehensive procedure for updating the canopy storage, please consult the previous study conducted by Terink et al. (2015). In light of the activation effect of extreme rainfall events on landslide occurrence and the performance of identifying past disaster sites, this research chose the annual average of 3-h maximum rainfall (AAMR-3H) (Fig. 6) from 2000 to 2020 as the input rainfall of the model, from multiple pre-experiments with different hourly intervals.

Fig. 6
figure 6

Annual average of the 3-h maximum rainfall (AAMR-3H) in the Dadu River Basin

Full size image

3.2.2 Quantifying the Parameters Associated with the Negative Impacts of Plants

Scholars have highlighted that plants can also have negative impacts on the occurrence of shallow landslides. These negative impacts include vegetation weight loads (Zou et al. 2021) and the enhancement of plant root systems on drainage (Shao et al. 2015), which enhance the depiction of the vegetated slope failure mechanism in the proposed model. To qualify these negative impacts, vegetation weight loads (Wveg) and topographic wetness index (TWI) are calculated as follows.

3.2.2.1 Vegetation Weight Loads (W veg)

The weight of vegetation, which increases the downward-sliding force of vegetated slopes, is one of the triggering causes of shallow landslides. Vegetation weight loads can be computed using Eqs. 12) and 13 (Zou et al. 2021).

$$ W_{{{\text{veg}}}} = \frac{{V_{{{\text{wood}}}} \times \rho_{{{\text{wood}}}} \times g \times \rho_{{{\text{tree}}}} \times A_{{{\text{slope}}}} }}{{B_{{{\text{slope}}}} }}, $$
(12)
$$ V_{{{\text{wood}}}} = \left\{ {\begin{array}{*{20}c} {0.39\left( {D + 3.5} \right)^{2} \left( {0.48 + L} \right)/10,000\;(D > 10\;{\text{cm}}),} \\ {0.7854L\left( {0.2\left( {D + 0.45L} \right)} \right)^{2} /10,000\; \, (D \le 10\;{\text{cm}}),} \\ \end{array} } \right. $$
(13)

where Wveg denotes the vegetation weight loads (N/m); D and L denote tree diameter and tree height, respectively; Vwood and \(\rho_{{{\text{wood}}}}\) denote wood volume and wood density, respectively; g denotes gravitational acceleration; Aslope and Bslope denote slope area (m2) and slope width (m), respectively; and \(\rho_{{{\text{tree}}}}\) denotes tree density (1/m2). Although there are limitations, for example, the area is not fully covered by vegetation, the quantitative parameters for the entire Dadu River region were assigned based on the vegetation zoning map (Fig. 7). According to the dominant vegetation types in different vegetation zones, the relevant vegetation canopy morphology parameters were assigned.

Fig. 7
figure 7

Vegetation zonation in the Dadu River Basin

Full size image
3.2.2.2 Topographic Wetness Index (TWI)

Numerous studies suggest that the existence of a plant root system enhances soil permeability, leading to increased soil saturation and reduced soil shear strength (Li et al. 2021). The preferential flow facilitated by the plant root system is considered detrimental to slope stability from a hydrological perspective. Therefore, when quantifying parameters related to infiltration processes in this study, the effects of plant roots on the hydraulic properties are considered for the root–soil mass rather than solely the fill soil mass.

The topographic wetness index (TWI) (Pack et al. 2005) was incorporated into the proposed model to quantify the m, the saturated part of the soil layer, under the support of Darcy’s law and field survey findings that areas with higher soil moisture content tend to be in areas of confluence (Beven and Kirkby 1979). The TWI is also an essential parameter for typical infinite slope models, such as SHALSTAB and SINMAP, across a specific period in disaster-prone places. The TWI is calculated using the following three presumptions: (1) shallow lateral subsurface flow follows topographic gradients; (2) lateral discharge at each point is in equilibrium with a steady-state recharge P; and (3) the capacity for lateral flux at each point is Tsinθ, which is hydraulic conductivity multiplied by soil depth, h (Fig. 8). The following are the computational formulas:

$$ {\text{TWI}} = \min \left( {\frac{{PA_{{{\text{up}}}} }}{TB\sin \theta },\;1} \right), $$
(14)
$$ T = \left( {k_{{{\text{s}},{\text{soil}}}} + k_{{{\text{s}},{\text{root}}}} } \right)h, $$
(15)

where TWI is equal to m, the ratio of the water height above the failure plane to the potentially unstable soil depth; ks=ks,soil+ks,root denotes the saturated conductivity of root–soil mass, which can be assigned according to Table 4. P and T are effective precipitation and soil transmissivity, respectively; P/T may be treated as the hydrologic ratio combining climate and hydrogeological indicators. Aup and B denote the upslope drainage area and the width of the outflow boundary, respectively; sinθ denotes head gradient (Fig. 9), and Aup/(bsinθ) can be treated as the topographic ratio to capture the effects of topographic convergence on concentrating runoff and elevating pore pressures, showing that (1) the greater the upslope drainage area in respect to the cell width, the higher the m; and (2) the larger the gradient of slopes, the higher flow velocity of the subsurface flow and the smaller the corresponding m value.

Fig. 8
figure 8

Soil depth in the Dadu River Basin

Full size image
Table 4 Hydraulic properties of various vegetations’ root–soil combination
Full size table
Fig. 9
figure 9

Head gradient (sinθ) in the Dadu River Basin

Full size image

3.2.3 Infinite Slope Stability Model Accounting for the Hydrological and Physical–Mechanical Impacts of Vegetation (ISM-HPM-VEG)

A novel physical-based hybrid model for assessing shallow landslide susceptibility on a regional scale is proposed (Fig. 10). This model combines various components such as the infinite slope model (ISM), vegetation canopy interception model, root system reinforcing model, vegetation weight loading model, and subsurface flow model. The aim was to comprehensively assess the positive and negative impacts of plants on the occurrence of shallow landslides. The positive effects include vegetation canopy interception and reinforcement of root systems, whereas the negative impacts involve plant gravity loading and preferential flow through root systems. A safety factor (SF) can be used, based on the infinite slope model, to quantify and categorize the level of landslide susceptibility. The formula is as follows.

$$ {\text{SF}} = \frac{{\left( {C_{{{\text{soil}}}} + C_{{{\text{root}}}} } \right)b + \left( {\left( {W_{{{\text{soil}}}} + W_{{{\text{veg}}}} } \right)\cos \beta - \gamma_{{\text{w}}} mhb\cos^{2} \beta } \right)\tan \phi }}{{\left( {W_{{{\text{soil}}}} + W_{{{\text{veg}}}} } \right)\sin \beta }}, $$
(16)
$$ W_{{{\text{soil}}}} = hb\cos \beta \left[ {\left( {1 - m} \right)\gamma_{{\text{m}}} + m\gamma_{{{\text{sat}}}} } \right], $$
(17)
$$ \gamma_{{{\text{sat}}}} = \gamma_{{\text{w}}} \left( {\theta_{{\text{s}}} - \theta_{{\text{i}}} } \right) + \gamma_{{\text{m}}} , $$
(18)

where C = Csoil + Croot denotes the integral cohesion of root–soil mass; φ and b denote the internal friction angle and bottom length of a calculation unit, respectively; h denotes the depth of potentially unstable soil layer; Wsoil and Wveg denote unstable soil mass weight and vegetation weight loads, respectively; γ and γw denote soil unit weight and water unit weight, respectively; β denotes slope angle; γm and γsat are the soil unit weight above the water table and saturated soil unit weight, respectively. For an explanation of other parameters, please refer to Sect. 3.2.

Fig. 10
figure 10

A conceptual diagram of the infinite slope model

Full size image

In accordance with the findings of Zou et al. (2021), the geotechnical model in this study is comprised of three layers: the root–soil layer, the soil layer, and the bedrock layer. The depth search procedure, similar to the classical TRIGRS model (Baum et al. 2008), sequentially scans for the minimum FS value associated with the depth, starting from shallow depths and progressing to deeper ones until reaching the upper surface of the bedrock layer. In cases where the search depth falls within the root–soil layer, the reinforcement of soil strength by the root system is taken into account. For depths within the soil layer, the determination of soil strength parameters is based on the specific soil type, as stipulated in the table documenting the strength parameters of different soil types in Zou et al. (2021).

3.3 Indicators for Reliability Analysis

To evaluate the model’s effectiveness and reliability, some indicators, such as prediction accuracy (ACC), precision, recall, and F1 score (Oommen et al. 2018; Zhang et al. 2019), were calculated by Eqs. 19–23.

$$ {\text{Precision}} = {\text{TP}}/({\text{TP}} + {\text{FP}}), $$
(19)
$$ {\text{Recall}} = {\text{TP/}}({\text{TP}} + {\text{FN}}), $$
(20)
$$ {\text{Specificity}} = {\text{TN}}/({\text{TN}} + {\text{FP}}), $$
(21)
$$ F1\;{\text{score}} = 2\;*\;{\text{PPV}}\;*\;{\text{TPR}}/({\text{PPV}} + {\text{TPR}}), $$
(22)
$$ {\text{ACC}} = \left( {{\text{TP}} + {\text{TN}}} \right)/({\text{TP}} + {\text{FN}} + {\text{FP}} + {\text{TN}}), $$
(23)

where ACC denotes prediction accuracy. TP and TN represent test sites number with correct predictions, whereas FP and FN represent test sites number with incorrect predictions (Table 5).

Table 5 Confusion matrix
Full size table

Additionally, the area under the curve (AUC) was calculated using the receiver operating characteristic (ROC) curve (Pradhan et al. 2016) and was used to assess model performance, as illustrated in Fig. 11. The highest value of AUC is 1.0, and a higher AUC means higher accuracy. This study determined the mean ROC by using inventoried landslide sites and non-landslide points of the same size derived from 100 random events with a uniform distribution across non-landslide areas (Hess et al. 2017). The purpose of random selection is to achieve an even spatial distribution of the chosen non-landslide samples within the Dadu River Basin.

Fig. 11
figure 11

Area under the curve (ROC) of different models

Full size image

4 Results

The performance of the model was assessed using several indicators (Table 5). According to Table 6, the proposed model (ISM-HPM-VEG) has a better performance with an AUC of 0.850 and an ACC of 0.764, demonstrating that it can serve for predicting shallow landslide occurrence in the extreme topographic transition belt, by taking into account physical failure mechanisms of the vegetated slopes. Furthermore, the performance of the traditional pseudo-static analysis method (TPSA) was also evaluated, and its output results were used to compare with the proposed ISM-HPM-VEG model. The comparison results indicate that ISM-HPM-VEG outperforms TPSA by 9.93% for ACC and 5.59% for AUC and that it is effective and reliable.

Table 6 Performance of the infinite slope stability model accounting for the hydrological and physical–mechanical impacts of vegetation (ISM-HPM-VEG) and traditional pseudo-static analysis (TPSA) method
Full size table

The outputs of the aforementioned models were used to reclassify the susceptibility levels into five classes (very low, low, medium, high, and very high) taking 1.50, 1.25, 1.00, and 0.50 as the breaking points. Figure 12 depicts the relative proportion of different susceptibility levels of each model. The region with the very low susceptibility level has the highest percentage (55.60%) in the proposed model (ISM-HPM-VEG), with the remaining 17.12%, 14.44%, 12.50%, and 0.34% of areas falling into the low, medium, high, and very high susceptibility levels, respectively. Susceptibility maps were then drawn on a GIS platform (Figs. 13, 14). The spatial consistency of the two landslide susceptibility maps was analyzed using Spearman’s rank correlation coefficients (Torizin et al. 2022). The results implied that the spatial distribution characteristics of the susceptibility levels derived from the two models are highly consistent, with strong spatial similarity, and the correlation coefficient is 0.78 at the 0.05 confidence level. It indicated that the ISM-HPM-VEG is capable of assessing shallow landslide occurrence in the transition zone with rapidly changing terrain, and the outputs are trustworthy. Moreover, the distribution of landslide susceptibility levels revealed that the steep bank slopes on both sides of the river in the middle section of the Dadu River Basin, which has experienced long-term river erosion, are more likely to be classified as having high or very high susceptibility levels. The northern plateau zone, with flatter topography, tends to be classified as having low or very low susceptibility levels.

Fig. 12
figure 12

Proportion of susceptibility levels determined by the infinite slope stability model accounting for the hydrological and physical–mechanical impacts of vegetation (ISM-HPM-VEG) and traditional pseudo-static analysis (TPSA) method

Full size image
Fig. 13
figure 13

Susceptibility map using the infinite slope stability model accounting for the hydrological and physical–mechanical impacts of vegetation (ISM-HPM-VEG)

Full size image
Fig. 14
figure 14

Susceptibility map using the traditional pseudo-static analysis (TPSA) method

Full size image

5 Discussion

Based on the ACC and AUC metrics (Table 6), the ISM-HPM-VEG exhibits superior spatial recognition ability compared to the TPSA in the comprehensive diagnosis of shallow landslide occurrences. But the percentage of historical landslide sites occurring in areas with high or very high susceptibility levels by the ISM-HPM-VEG is lower than by the TPSA (Fig. 15), indicating a less discriminative capacity of ISM-HPM-VEG in identifying disaster-prone areas compared to the TPSA model. This phenomenon can be attributed to the fact that the ISM-HPM-VEG incorporates both the negative and positive effects of vegetation on landslides, leading to fewer regions being classified as having high or very high susceptibility levels. This occurs when the positive effect of vegetation significantly impacts slope stability. The ISM-HPM-VEG also achieves a higher precision score compared to the TPSA model. Furthermore, in extreme topographic transition zones where flat areas suitable for building construction are scarce, it is crucial to identify areas with very low and low landslide susceptibility that are suitable for human settlements. In that sense, the ISM-HPM-VEG can offer valuable and better support for land use planning and risk management.

Fig. 15
figure 15

Proportion of historical shallow landslides occurring in areas with different susceptibility levels by the infinite slope stability model accounting for the hydrological and physical–mechanical impacts of vegetation (ISM-HPM-VEG) and traditional pseudo-static analysis (TPSA) method

Full size image

By investigating the relationship between FS value and slope, the findings showed that slope has a critical impact on landslides occurring in the study area. The steeper the slope is, the greater the susceptibility level because of the larger shear stress caused by gravity (Pradhan et al. 2019; Araújo et al. 2022). Areas with slopes ranging from 30° to 40° are most prone to landslides (Fig. 16). This is largely consistent with the results of Qin et al. (2023) for the Dadu River Basin, where they found that the slope intervals with the highest likelihood of vegetated slope destabilization are (1) from 26.0° to 37.0° for arbor trees; (2) from 28.0° to 39.0° for shrubs; and (3) from 31.0° to 39.0° for meadows.

Fig. 16
figure 16

The number of grids with different susceptibility levels by slope group

Full size image

To further investigate the model’s sensitivity to vegetation-related parameters, this research set four representative values of m (that is, increased by 25%, 50%, 75%, and 100% of m) and Wveg (that is, increased by 25%, 50%, 75%, and 100% of vegetation weight loads), and the counterpart susceptibility maps were obtained, as shown in Fig. 17. As depicted in Fig. 18, the proportion of areas with high or very high susceptibility within the Dadu River Basin increased significantly when increasing the value of m, especially for the areas with very high susceptibility level, which have a 92.77% to 1122.62% increase. In contrast, the proportion of areas with low and very low susceptibility levels decreased significantly. The high sensitivity of m is attributed to the input parameter’s ability to represent the physical and hydrological mechanisms underlying landslide occurrence, such as rainfall, vegetation–rainfall interception, evaporation, infiltration, and slope destabilization mechanisms, indicating that there is a strong link between the aforementioned hydrological processes adopted in the derivation process of the safety factor and the triggering mechanism of shallow landslides.

Fig. 17
figure 17

Susceptibility maps based on different settings for the representative values of m

Full size image
Fig. 18
figure 18

Proportion of the different susceptibility levels after setting four representative values of m and Wveg. The area above the blue curve indicates the proportion of the region with a landslide susceptibility level above medium, which represents risk areas

Full size image

Unsurprisingly, the proportion of locations with high or very high susceptibility rises with increasing Wveg. This implies the proposed model’s ability to represent the negative impact of the vegetation weight loading on slope stability.

Through statistical analysis of the proportion of various landslide susceptibility levels under different ratios of Wveg to Wsoil (PDD) (Fig. 19), additional examination was conducted on the effect of vegetation gravity load on the stability of a vegetated slope. The results indicated that as the ratio of Wveg to Wsoil increased in landslide-prone regions (that is, regions with medium and above susceptibility levels), the PDD increased and then decreased. When the ratio of Wveg to Wsoil is between 0.001 and 0.01, the PDD reaches a maximum of 53.61%. In areas not prone to landslides, the trend is the opposite. This suggests that when the ratio of Wveg to Wsoil is comparatively small (that is, less than 0.01), an increase in vegetation gravity load is detrimental to slope stability. Generally, an increase in plant biomass will exacerbate the instability of a planted slope. However, when the ratio of Wveg to Wsoil is high (that is, greater than 0.01), an increase in vegetation biomass tends to be advantageous for slope stability. This may be because the positive effects of vegetation (for example, root reinforcement) limit or offset the negative effects of increasing vegetation gravity load on slope stability. In mountainous regions, the ratio of Wveg to Wsoil may play a crucial role in directing disaster prevention and reduction policies.

Fig. 19
figure 19

Proportion of areas with varying susceptibility levels for different ranges of parameter values. The area above the blue curve indicates the proportion of the region with a landslide susceptibility level above medium, which represents risk areas

Full size image

Figure 19 shows that an increase in the value of C (C = Csoil + Croot) is associated with a reduction in the overall percentage of areas that are susceptible to disasters. This suggests that strong reinforcement from roots greatly enhances slope stability, particularly when the C value exceeds 75kPa. Under this circumstance, the proportion of disaster-prone areas reaches a record low of 19.43%.

Although the ISM-HPM-VEG outperforms the TPSA in terms of predicting shallow landslide susceptibility in the study area, the ISM-HPM-VEG model’s input parameters are more complex due to the consideration of the positive and negative impacts of plants on slope failure, implying that the ISM-HPM-VEG is relatively suitable for those areas with counterpart soil surveys by which soil physical–mechanical parameters (including cohesion, internal friction angle, bulk weight) are recorded and easily accessible, as well as soil hydraulic parameters (such as saturation conductivity, saturated or stagnant volumetric water content, and so on). The next step will therefore explore methodologies that are more transferable to different places when incorporating vegetation-related factors to further reveal the regional-scale physical failure mechanisms of vegetated slopes.

Currently, at the regional scale, several physical models based on infinite slope stability have considered some of the vegetation impacts, but not comprehensively. For example, the GIS-TISSA model (Escobar-Wolf et al. 2021) provides an estimation of vegetation surcharge with varying heights and employs the root cohesion converted from the NDVI (Jr-Chuan et al. 2006) as an input parameter to capture the reinforcing impact of the root system on slope stability, but ignores the influence of vegetation interception; and the TRIGRS model (Kim et al. 2010) is characterized by its powerful capacity to resolve complex storm processes in complicated landscapes, by taking infiltration and hydraulic properties as input parameters that vary over the grid area, but is still lack of consideration of some effects of vegetation on landslide occurrences, such as the positive impact of vegetation canopy interception and root systems reinforcement. The comprehensive consideration of vegetation-related impacts is essential for assessing landslide susceptibility, particularly in extreme topographic transition belts with significant horizontal and vertical vegetation variabilities. In this regard, the proposed methodology is novel by incorporating the aforementioned positive and negative effects of vegetation for the regional-scale susceptibility assessment of shallow landslides, thereby improving the physical representation of vegetated slope failure mechanisms in models predicting the regional-scale spatial likelihood of shallow landslide occurrence.

6 Conclusion

A physical-based methodology is proposed to predict the spatial likelihood of shallow landslides in a transition belt characterized by rapidly changing terrain. The methodology integrates the infinite slope model with other models that consider the positive and negative impacts of plants on the occurrence of shallow landslides. The application results, across the entire Dadu River Basin, showed that (1) the proposed model is effective and capable of assessing landslide susceptibility in the region, after accounting for the positive and negative impacts of plants on shallow landslide occurrence; (2) the model’s performance is improved by 9.93% in ACC and 5.59% in AUC over the TPSA model; and (3) although the effects of plants on shallow landslide occurrence are positive and negative, the positive effects of vegetation are more important in the study area. Lastly, further exploration is needed to understand how the regional representation of specific vegetation types’ parameters affects the predictive capabilities of the model. The susceptibility study presented in this article may provide a scientific support to local governments for the prevention and mitigation of shallow landslides.