Regional seismic landslide susceptibility assessment considering the rock mass strength heterogeneity

Abstract Rock mass strength parameters are critical when applying the permanent-displacement analysis to landslide susceptibility assessment. Previous studies use limited measured data to assign a constant value for each lithology. Since the rock mass strength is spatially heterogeneous, assigning a constant value to each lithology cannot effectively reflect the spatial heterogeneity of rock mass strength, reducing the reliability of seismic landslide susceptibility assessment. Here we present an analysis of a well-documented seismic landslide inventory and the lithological environment, and develop an empirical model to quantify the rock mass strength heterogeneity. The rock mass strength heterogeneity was subsequently taken as an input of the permanent-displacement analysis to assessing regional seismic landslide susceptibility. Results on two case studies show that considering the rock mass strength heterogeneity would significantly improve seismic landslide susceptibility assessment. The findings of this study are valuable for earthquake emergency rescue and post-earthquake land plans.


Introduction
A strong earthquake may instantly trigger a large number of landslides that can be catastrophic, resulting in life loss, fatalities, or economic damage.For example, the 2008 Wenchuan earthquake-induced 197,481 seismic landslides in an area of $110,000 km 2 (Xu et al. 2014).Apart from seismic landslides, strong earthquakes also generate numerous discontinuity sets prone to slope instability.Heavy precipitation after the quake causes the slope stability to deteriorate, presenting serious threats to lives and properties (Guzzetti et al. 2006;Robinson et al. 2017).For instance, after the Wenchuan earthquake, the landslide frequency increased substantially compared with that before the earthquake (Huang and Li 2014).Thus, timely and reliably assessing the seismic landslide susceptibility is essential to prevent and mitigate the effects of seismic landslides (Tian et al. 2019;Segoni et al. 2020;Liu et al. 2021).
Statistical analysis and permanent-displacement analysis are two commonly used methods for seismic landslide susceptibility assessment (Pourghasemi et al. 2018;Ibrahim et al. 2021;Zou et al. 2021).The statistical analysis establishes a mathematical model by exploring the relationship between landslides and geo-environmental factors from seismic landslide inventory.The established mathematical model is then applied to assess the seismic landslide susceptibility for the whole area (Huang and Zhao 2018;Reichenbach et al. 2018;Merghadi et al. 2020).Its reliability relies on the quality and completeness of the landslide inventory.A high-quality and complete landslide inventory, however, is challenging or even impossible to produce, which limits the statistical analysis applications (Du et al. 2020;Chen et al. 2021).In contrast, the physically-based permanent-displacement analysis does not require the landslide inventory, and thus can be used to assess seismic landslide susceptibility in various scenarios (Wieczorek et al. 1985;Jibson 2011;Liu et al. 2018;Kumar et al. 2021).It takes a slope as a sliding-block with a certain critical acceleration.The critical acceleration represents the minimum acceleration that causes the sliding-block to slide.Once an earthquake's ground motion acceleration exceeds the block's minimum acceleration, the block initiates sliding deformation (also called permanent seismic displacement).Numerous studies have validated that the displacement is an important indicator to describe the seismic landslide susceptibility (Salvatici et al. 2018;Martino et al. 2020).To date, this method has been widely applied to landslide susceptibility assessment for many earthquakes, such as the 1994 Northridge earthquake (Jibson 1993;Jibson et al. 2000), 1999 Chichi earthquake (Tsai et al. 2019(Tsai et al. ), 2008 Wenchuan earthquake (Godt et al. 2008;Wang et al. 2013;Ma and Xu 2019), 2014 Ludian earthquake (Chen et al. 2019;Huang et al. 2020), and 2017 Jiuzhaigou earthquake (Chen et al. 2020).
When applying the permanent-displacement analysis to seismic landslide susceptibility assessment, critical acceleration of the slope is a prerequisite, which can be calculated based on the rock mass shear strength of the slope.Rock mass shear strength can be described by different types of strength models (e.g.Mohr-Coulomb, Hoek-Brown, Barton).Among them, the Mohr-Coulomb strength parameters have received much attention in permanent-displacement analysis (Jibson 1993;Gallen et al. 2015;Shinoda and Miyata 2017), which can be obtained from the laboratory test on selected samples (i.e.direct shear or triaxial tests) (Labuz and Zang 2012).However, the sample in laboratory is too small to reflect the widely developed scale-related structural features (e.g.fractures, joints, and discontinuities) (Hoek and Brown 1980;1997).Compared to the laboratory test, in-situ test can obtain reliable rock mass strength parameters, while conducting such a test at a regional scale is time-consuming, labour-intensive, or even impractical.Besides, the limited measurement data cannot represent the real rock mass strength parameters in a region.Therefore, obtaining reliable rock mass strength parameters remains an immense challenge, particularly at a regional scale.
Previous studies usually assigned a constant value for each lithology in a region based on the outcropping lithologies (Jibson et al. 2000;Godt et al. 2008;Wang et al. 2017;Salinas-Jasso et al. 2019;Kumar et al. 2021).This assignment ignores the variability of rock mass strength in space, resulting in uncertainty in seismic landslide susceptibility.To tackle this issue, researchers tried to obtain reliable rock mass strength parameters from different perspectives (Wang et al. 2008;Wang and Rathje 2015;Lee and Park 2016;Nguyen and Kim 2020).For example, the Monte Carlo method was applied to simulate the variability of rock mass strength parameters by assuming that the rock mass strength fits a specific distribution (e.g.normal, uniform) (Murphy and Mankelow 2004;Lee and Park 2016).Rodr ıguez-Peces et al.
(2020) designed a logic tree for descripting the uncertainty in rock mass strength parameter for the displacement analysis method.Although these methods could improve the reliability of the rock mass strength parameters, they require a large number of measurements, which limits the application of these methods.Given the availability of seismic landslide inventories, rock mass strength parameters can be inverted using ground motion parameters and seismic landslides (Gallen et al. 2015;Bunn et al. 2020;Huang et al. 2020;Townsend et al. 2020).For example, Shinoda et al. (2019) used different combinations of rock mass strength parameters to calculate permanent seismic displacements and determined the most appropriate rock mass strength parameters (i.e.effective cohesion and friction angle) by comparing the calculated displacement with the observed landslides.Jeandet et al. (2019) used Coulomb mechanics and relief constraints to explain the landslide size distribution, finding that reasonable low ranges of Mohr-Coulomb shear strength properties can effectively model seismic landslide size distribution.However, these methods also assumed a constant strength value for a lithology and assigned an average or residual strength parameter for each lithology based on the geological map.Although this assumption makes the operation easy, the variability within a geological unit is typically ignored due to practice constraints (Dreyfus et al. 2013;Rodr ıguez-Peces et al. 2020).
Rock mass properties depend not only on lithology but also on rock mass structure and weathering characteristics.Studies shown that, compared with lithology differences, rock mass structure (e.g.fractures, joints, and discontinuities) and weathering characteristics have a more substantial influence on rock mass properties (Schmidt and Montgomery 1995;Alberti et al. 2022;Medwedeff et al. 2022).Influenced by the lithological environment, the rock mass structure and weathering characteristics exhibit a strong degree of spatial heterogeneity (Baroni et al. 2010;Bicocchi et al. 2019;Fusco et al. 2021).For example, Medwedeff et al. (2022) evaluated landscape-scale patterns in near-surface mechanical characteristics and their relationship to environmental factors that affect rock strength by collecting in-situ measurements.They found that shallow bedrock strength is more dependent on the degree of chemical and physical weathering, rather than the mineral and textural differences between the metamorphic lithologies.Obviously, assigning a constant value for each lithology cannot reflect the spatial heterogeneity of rock mass strength in complex environment.Therefore, to obtain a reliable rock mass strength parameter, its spatial heterogeneity should be given reasonable consideration.
Here, we proposed a novel method to estimate the rock mass strength heterogeneity and applied it as the input of the permanent displacement analysis method for regional seismic landslide susceptibility assessment.We considered lithological environmental factors known to affect rock mass strength to estimate the rock mass strength heterogeneity.An empirical model was developed based on a well-documented seismic landslide inventory of the 2008 Wenchuan earthquake.The estimated spatial heterogeneity of rock mass strength and the permanent-displacement analysis were then jointly used for regional seismic landslide susceptibility assessment.Finally, two earthquakes with similar geo-environmental conditions were selected for model verification.

Study areas
This study selected the Wenchuan and Ludian earthquakes (Figure 1) as study cases.The Wenchuan Mw 7.9 earthquake occurred on May 12, 2008, with the epicentre located southwest of Yingxiu, Wenchuan (103.4E, 31.0N), at a focal depth of 14 km.The seismic rupture of this earthquake was located in the Longmenshan Fault on the eastern edge of the Qinghai-Tibet Plateau.By September 25, 2008, the Wenchuan earthquake resulted in 69,227 people dead, 17,923 people missing, and 374,643 people injured (Wang 2008).This earthquake also resulted in numerous mountain fractures, seismic landslides, and other secondary disasters, forming many new potential geological hazard sites, which adds to the difficulty of post-disaster recovery (Huang et al. 2012).This article selected a rectangular area (240 km Â 95 km) around the epicentre as the study area to evaluate seismic landslide susceptibility (Figure 1(a)).Following the Wenchuan earthquake, a Mw 6.1 earthquake occurred on August 3, 2014, with the epicentre located near Ludian (103.3E, 27.1 N), at a focal depth of 12 km.The seismic rupture of this earthquake was the northwest-oriented Baogushan-Xiaohe Fault, which is a secondary strike-slip fault of the Zhaotong-Ludian Fault (Figure 1(b)).By August 8, 2014, the Ludian earthquake resulted in 617 people dead, 112 people missing, and 3,143 people injured (Hong et al. 2016).Given that the Ludian and Wenchuan earthquakes have similar characteristics of geo-environment, the former was selected to verify the transferability of the empirical model of spatial heterogeneity developed based on the latter earthquake.

Data preparation
Data sources used in this article include landslide inventories, a digital elevation model (DEM), a geological map, ground motion parameters, fault, and river databases (Table 1).These data were converted into raster data with a spatial resolution of 90 m and 30 m for Wenchuan and Ludian, respectively, and re-projected to a unified spatial coordinate system (UTM-Zone 48, WGS84 datum).

Seismic landslide inventory
Seismic landslide inventory records the number, location, and size of landslides triggered by earthquakes (Harp et al. 2011).A complete landslide inventory is critical for seismic landslide susceptibility assessment and is also an important reference to verify the reliability of assessment results (Chen et al. 2020;2021).After the Wenchuan earthquake, different landslide inventories were compiled (Dai et al. 2011;Gorum et al. 2011;Xu et al. 2014).This study selected the landslide inventory released by Xu et al. (2014), compiled based on pre-and post-earthquake high-resolution remote sensing imagery with visual interpretation, due to its high-quality and completeness.The inventory includes 197,481 seismic landslides, affecting an area of $110,000 km 2 .The distribution of landslides triggered by the Wenchuan earthquake shows that seismic landslides were closely related to geo-environmental factors (Yin et al. 2009;Qi et al. 2010;Xu et al. 2014).For example, influenced by river distribution, seismic landslides mainly developed along both sides of rivers, and most seismic landslides distributed near faults, especially the seismogenic fault.
Given that only the landslide source area represents the occurring area, accumulation areas should be excluded from the inventory.This study defined landslide source areas as grid cells with elevations above the median elevation of each landslide (Jibson et al. 2000).In the Wenchuan earthquake, landslide source areas comprised nearly 1061.88 km 2 , accounting for 5.09% of the total area (Figure 2(a)).For the Ludian earthquake, this article selected the landslide inventory compiled by Xu et al. (2014) for the same consideration and extracted landslide source areas using the same method for Wenchuan (see Figure 2(b)).
Seismic landslides fall into deep-seated landslides and shallow landslides according to the landslide thickness.Since the developed method is mainly applicable to shallow landslides, deep-seated landslides need to be removed from the landslide inventory.However, precise landslide thickness is extremely difficult to obtain, making it challenging to distinguish shallow landslides from deep-seated landslides.This article, therefore, calculated the landslide thickness t using an area-volume empirical model: where V ls and A ls represent the landslide volume and the landslide area, respectively.The model parameters mainly obtained from an ordinary regression on data sets from the Longmen Shan landslide site measurements (Parker et al. 2011), which has been widely used to quantify the volume of the Wenchuan seismic landslide (Li et al. 2014;Fan et al. 2018).A landslide is taken as a deep-seated landslide if its thickness is above 10 m (Roering et al. 2003;Dou et al. 2015) and will be removed.For Wenchuan, the calculated landslide thickness ranges from 0.4 m to 46 m, with 93% shallow landslides.Finally, for Wenchuan and Ludian areas, after removing deep- seated landslides from the landslide inventories, the average thickness values of shallow landslides are 3.54 m and 1.80 m, respectively.This information will be used as the depth information in the subsequent permanent-displacement analysis.

Topography data
DEM is one of the most used topography data in landslide susceptibility assessment.This study used a Shuttle Radar Topography Mission DEM (Van Zyl 2001), which can be downloaded from the official website of United States Geological Survey (USGS).Topography-related factors, including slope angle and topography curvature, were extracted from the DEM (Figure 3).

Geological map and rock mass shear strength parameters
Shear strength parameters (c 0 and u 0 ) are assigned based on the mapped geologic units.The distribution of geological units is mainly derived from the 1:250,000 geological map published by China Geological Survey.This study aims to explore the influence of spatial heterogeneity of rock mass on seismic landslide susceptibility assessment, not the produce of accurate rock mass shear strength parameters.Thus, we assigned the rock mass shear strength parameters according to Ma and Xu (2019).
Similarly, for the Ludian earthquake area, the lithology was grouped into three lithological groups and assigned representative shear strength for each lithological group (Chen et al. 2019).Table 2 shows rock mass strength parameters assigned in different lithological groups in Wenchuan and Ludian, and Figure 4 shows the lithological group distributions in Wenchuan and Ludian.

Ground motion parameters
Ground motion parameters are important input parameters in permanent seismic displacement calculation, including peak ground acceleration (PGA), peak ground velocity, and Arias intensity.The ground motion parameters released by ShakeMap (Wald et al. 2005) are considered to be the best estimates of regional ground motions, which are widely used in seismic landslide susceptibility assessment.ShakeMap can provide several ground motion parameters, including PGA, peak ground velocity, and spectral acceleration at different periods.The PGA was selected as the ground motion parameter, since it is an important input parameter for the empirical displacement model selected in this article.Figure 5 shows the PGA maps of the Wenchuan and Ludian earthquakes.

Fault and river data
This article also considered the influence of fault and river distribution on rock mass properties.Specifically, distance to fault and river density were used to represent the fault and river distribution.Distance to fault was calculated based on fault distribution using Euclidean distance tool in QGIS and were then reclassified into six categories using the natural breakpoint classification algorithm, as shown in Figure 6(a,b).
For river density, this study selected river distribution with the 5th order, and the river density was calculated using density calculation tool in QGIS and was then reclassified into six categories, as shown in Figure 6(c,d).

Methodology
3.2.1.Principle of permanent-displacement analysis 3.2.1.1.Infinite-slope model.The infinite-slope model takes a potential landslide as a sliding rigid block resting on an inclined plane, which is a commonly used model in permanent-displacement analysis and can effectively simulate thinner landslides in a stiffer material, such as shallow landslides (Jibson 2011).Since shallow landslides are the main type of earthquake-induced landslides (Keefer 1984;2002), the infinite-slope model is suitable for modelling seismic landslides and can be used for seismic landslide susceptibility assessment (Wieczorek et al. 1985).Figure 7 shows an infiniteslope model under static conditions

Critical acceleration.
In infinite-slope model, the resistance force (Fr) of the sliding-block is determined by the rock mass properties, which are usually expressed by the Mohr-Coulomb parameter (i.e. the effective cohesion (c 0 ) and the internal friction angle (u 0 )).
Where, c and c w denote the unit weight of the sliding-block and water, respectively.m and t are the saturation and depth of the sliding-block, respectively.S is the contact area between the sliding-block and the sliding surface.h represents the slope angle of the sliding surface.
Under the earthquake condition, the effective sliding force (Fs 0 ) of the slidingblock is mainly composed of a gravity and the seismically driven load along the downward direction, which can be calculated as follow: where a is the ground motion acceleration and g is the acceleration of gravity, usually taken 9.80 m/s 2 .Suppose that the sliding-block is in the limit equilibrium state (Fr ¼ Fs 0 ), the ground motion acceleration at this time is the critical acceleration of the sliding-block (a c ), which is the minimum ground motion acceleration required to make the sliding-block slide.By combining Eq. (3) and Eq. ( 4), the critical acceleration can be calculated as follows (Newmark 1965): The above equation can also be further simplified as a function of the static safety factor of the sliding-block and the inclination angle of sliding surface.
where FS is the static safety factor of the sliding-block, which can be expressed as: The critical acceleration of the sliding-block depends on 1) the physical-mechanical parameters of the material of sliding-block (i.e. the effective cohesion, the effective friction, and the unit weight of the rock mass), 2) the inclination angle of the sliding surface, 3) the depth and saturation of the sliding-block.Here, the physical-mechanical parameters of the rock mass can be obtained using in-situ measurements or assigned empirical values based on geological unit.The inclination angle of the sliding surface is usually replaced by the topography slope.The depth (t) and the saturation (m) of the sliding-block are two important parameters that affect the sliding force of the sliding-block, which are challenging to obtain at the regional scale.Stateof-the-art often assigned constant values to t and m: In this article, the average landslide thickness in landslide inventory is used as the depth (t) parameter of the sliding-block for Wenchuan and Ludian earthquakes, and the saturation parameter of the sliding-block (m) is set to 0 since there was no heavy precipitation before Wenchuan and Ludian earthquakes.
3.2.1.3.Permanent seismic displacement calculation.Permanent seismic displacements can be calculated in two ways.The first way requires a complete time record of ground motion acceleration, and the displacement is calculated by quadratic integration of the time record of ground motion acceleration above the critical acceleration.Since the time record of ground motion acceleration contains the earthquake characteristics (e.g.magnitude, source mechanism), the calculated displacements are reliable.However, a complete time record of ground motion acceleration is limited, especially at a regional scale, which makes it relatively difficult to apply.The second way uses an empirical displacement model to calculate the displacement.The empirical displacement model can be established based on a large number of time records of ground motion acceleration from global/regional.For the latter one, the time record of ground motion acceleration was replaced with ground motion parameters, such as PGA.Compared to the first method, the second method has been widely used in practice, and thus we compared four benchmark empirical displacement models and selected the best one to calculate the permanent seismic displacement.3.2.2.Spatial heterogeneity analysis of rock mass strength 3.2.2.1.Lithological environment analysis.Lithological environment refers to the external geo-environment in which the outcropping lithology is located, which is one of the main causes leading to the spatial heterogeneity of rock mass strength.Accordingly, the spatial heterogeneity of rock mass strength can be quantified based on the lithological environment difference.Considering the importance, availability, and generalization of the lithological environment, we selected three major factors to represent the lithological environment, including river, fault, and topography curvature.
The river represents the hydrological impact of geo-environment.The seasonal fluctuation of river level can affect rock mass properties.Besides, the erosion and scouring of the river are easy to form high-steep terrains that in turn accelerate the deformation and destruction of the slope.Distance to river or river density is commonly used to reflect the river distribution.This study selected the river density as a hydrological factor.The fault distribution reflects the regional tectonic environment.Influenced by tectonic activity, the rock mass in the near-fault often develops dense structural surfaces, leading to low rock mass strength, and such region is prone to form seismic landslides.In this article, distance to fault was selected as the tectonic factor.Topography curvature reveals the degree of concavity and convexity of the terrain.In abrupt topographic change areas, such as concave and convex slopes, the stress redistribution and stress cracks are easily formed, which reduces the rock mass strength on the slope.Therefore, the topography curvature was also selected as a lithological environmental factor.

Multivariate analysis.
Earthquakes and their associated seismic landslides provide an opportunity to quantify the rock mass strength at a regional scale (Jeandet et al. 2019;Shinoda et al. 2019).Under earthquake condition, the earthquake provides a regional shear-strength test where the earthquake generates a measurable forcing (ground motion acceleration) and an observable response (seismic landslide) which depend on the rock mass strength (Gallen et al. 2015).Hence, a well-documented seismic landslide inventory, including location, volume, source area, and types of all seismic landslides, is critical in quantifying the rock mass properties.
Lithologic environmental factors contribute differently to the rock mass strength, affecting the spatial distribution of seismic landslides.Therefore, the spatial distribution of seismic landslides synthetically reflects the differences of lithologic environments.This article adopted logistic regression analysis (LRA) (Mousavi et al. 2011;Xie et al. 2018;Cao et al. 2019) to estimate contributions of different lithologic environmental factors to the rock mass strength based on the seismic landslide inventory.
The dependent variable in LRA is binary (e.g.0/1), where Y ¼ 0 represents the nonlandslide occurrence and Y ¼ 1 represents the landslide occurrence.The landslide occurrence is influenced by various geo-environmental factors, which can be expressed as: where x i ði ¼ 1, :::nÞ are the geo-environmental factors, including the factors related to rock mass properties, topography factors, and earthquake-related factors; and b i ði ¼ 1, :::nÞ are the coefficients of x i ði ¼ 1, :::nÞ, representing the contribution of each geo-environmental factor to the seismic landslide occurrence.
Odds ratio is the ratio of landslide occurrence probability (pðY ¼ 1Þ) to non-landslide occurrence probability (pðY ¼ 0Þ ¼ 1 À pðY ¼ 1Þ).The natural logarithm of the odds ratio is a linear expression of the independent variable x i ði ¼ 1, :::nÞ : The coefficient b i ði ¼ 1, :::nÞ is derived using the maximum likelihood estimation.If b i ði ¼ 1, :::nÞ is positive, the geo-environmental factor positively contributes to landslides, while the opposite is a negative contribution.
To estimate the rock mass strength heterogeneity, a linear equation was used to establish the spatial correlation between the rock mass strength heterogeneity and the lithological environment, which can be express as follows: where f is the comprehensive influence of different lithological environmental factors on the rock mass strength.A high value of f represents a high influence of the lithological environment on rock mass strength, indicating a low rock mass strength.
3.2.2.3.Empirical model of the rock mass strength heterogeneity.We selected the landslide inventory compiled by Xu et al. (2014), which recorded most of landslides induced by the Wenchuan earthquake, to develop an empirical model of rock mass strength heterogeneity.Specifically, to develop an empirical model of rock mass strength heterogeneity, a data set that includes an equal number of landslide points and non-landslide points was required.Landslide points were randomly selected in the landslide source areas, and non-landslide points were randomly selected in areas far from the landslide source areas (at least 1 km away).According to the established data set, the regression coefficient of each geo-environmental factor was calculated, see Table 3.A higher coefficient indicates a higher contribution of the geo-environmental factor to landslides.Table 3 shows that all selected factors positively contribute to seismic landslides except the topography curvature and distance to fault.All geo-environmental factors are statistically correlate to landslides, with significance being 5% confidence level.Among them, the lithology group and PGA have the greatest contributions to the seismic landslide.Besides, since the hard lithology group tends to form high and steep topography, its contribution to seismic landslides is higher than that of the low lithology group.
According to the lithological environment and the weight coefficients calculated using LRA, an empirical model was developed to estimate the rock mass strength heterogeneity, which can be expressed as follows: where f cur , f fault and f river represent the topography curvature, distance to fault, and river density, respectively.It is worth noting that different lithologic environmental factors have different effect on rock mass in different regions, which needs to be reasonably determined according to local conditions.Since the empirical model is developed based on the landslide induced by the Wenchuan earthquake, it is only suitable for areas with similar geo-environments to the Wenchuan earthquake, such as Ludian earthquake.

Procedure of seismic landslide susceptibility assessment
Critical acceleration positively correlates with the rock mass strength and plays an essential role in permanent seismic displacement calculation.Therefore, the rock mass strength heterogeneity can also affect the critical acceleration distribution.To incorporate the rock mass strength heterogeneity in the permanent-displacement calculation, the rock mass strength heterogeneity was normalized to [0.5 1] using a linear equation (Eq.( 12)), and the critical acceleration can be modified with the normalized heterogeneity as follows (Eq.( 13)): where f m is the rock mass strength heterogeneity after normalization; a c is the critical acceleration without considering the rock mass strength heterogeneity, and a cM is the modified critical acceleration considering the rock mass strength heterogeneity.
Combining yield acceleration and ground motion parameters, permanent seismic displacements can be quickly obtained.Based on the distribution of permanent seismic displacements, the earthquake affected areas can be categorized into different seismic landslide susceptibility levels.

Experimental results
To evaluate the performance of the developed method, this article selected three evaluation metrics: percentage of ground failures captured (%GFC), percentage of the where Dn represents permanent seismic displacement in centimetre, and T represents critical displacement.
In general, a high %GFC means that calculated permanent seismic displacement is more reliable, which can capture more observed earthquake-induced landslides.However, a high %GFC often corresponds to a high %Ls Pred, in which permanent seismic displacement overestimates landslide susceptibility.To balance %GFC and %Ls Pred, this study introduced the %Efficiency to reflect the validity of permanent seismic displacement.A high %Efficiency corresponds to a high %GFC and low %Ls Pred: As the instabilities frequently occur with a displacement larger than 1 cm (Bray and Travasarou 2007), we set the critical displacement threshold to 1 cm.The area with a displacement larger than 1 cm will be predicted as a seismic landslide.8(a), fault distribution controls the spatial distribution of rock mass strength, compared to topographic curvature and river distribution.Influenced by the fault distribution, the rock mass in the near-fault area corresponds to a lower strength, which is conducive to the occurrence of landslides.Figure 8(b) shows the critical acceleration calculated without considering the rock mass heterogeneity in a previous study (Ma and Xu 2019).The critical acceleration ranges from 0.01 g to 1.18 g, with most regions having a critical acceleration greater than 0.2 g. Figure 8(c) shows the modified critical acceleration considering the rock mass strength heterogeneity in Wenchuan.Critical acceleration ranges from 0.01 g to 1.64 g (Figure 8(c)).Considering the influence of complex lithologic environmental factors on the rock mass strength, the area with high critical acceleration in Figure 8(b) is greatly reduced.For example, the critical acceleration in the near fault area has lower value compared to that far away from the fault.Besides, the difference of critical acceleration in Figure 8(c) is more evident than that presented in Figure 8(b), which is meaningful to subsequently divide regional landslide susceptibility.
A suitable empirical displacement model is essential to calculate the regional permanent seismic displacement.To this end, this article compared four empirical displacement models developed based on medium and strong earthquakes, see Table 4.
The calculated displacements were then compared with the observed seismic landslides of the Wenchuan earthquake. Figure 9 shows that the permanent seismic displacements calculated using different empirical displacement models with and without considering the rock mass strength heterogeneity.This article used a uniform colour scale for display and superimposed seismic landslides of the Wenchuan earthquake.As shown in Figure 9, different empirical displacement models can obtain different displacements, which corresponds to different seismic landslide assessment result.Specifically, without considering the rock mass strength heterogeneity, the displacement calculated with RS09 has the best prediction accuracy in capturing 71.49% of the observed landslides, while the displacement calculated with Jib07A the lowest prediction accuracy in capturing only 38.41% of the observed landslides.Consideration of the rock mass strength heterogeneity substantially improved the prediction accuracy and validity of displacements calculated by four empirical displacement models.Among them, the prediction accuracy of displacement calculated with Jib07A has improved significantly from 38.41% to 46.26%, and the validity of displacement has increased from 26.42% to 31.37%.Besides, the prediction accuracy of  Mw 5-7.9 ln ðDÞ ¼ 5:52 À 4: (Rathje and Saygili 2009) Mw 5-7.9 ln ðDÞ ¼ 4:89 À 4:85 ac displacement calculated with RS09 also increased from 71.49% to 76.62%, and the validity increased from 42.13% to 47.26%.In summary, among four models, the displacement calculated with Jib07A and Jib07B underestimates the regional seismic landslide susceptibility.By contrast, the RS08 and RS09 have higher prediction accuracy, and thus they are more applicable to the permanent seismic displacement calculation of the Wenchuan earthquake.Compared to RS08, RS09 has higher predictive accuracy as it considers the effect of earthquake magnitude.Therefore, we chose RS09 to calculate the permanent seismic displacement.

Landslide susceptibility assessment
We divided the calculated displacement into three susceptibility levels: low susceptibility (0 < D < 1cm), moderate susceptibility (1 cm D < 5cm), and high susceptibility (5 cm D), according to the USGS division criteria (Dreyfus et al. 2013).Figure 10 presents two landslide susceptibility maps that consider with/without the rock mass strength heterogeneity.We can see that compared to Figure 10(a), areas with moderate susceptibility level in Figure 10(b) are reduced significantly, especially in the northeast of the study area (e.g.subarea A).Although this area has a dense distribution of faults, it is still of high rock strength after integrating other lithologic environments, and therefore corresponds to a low susceptibility level.In contrast, areas with high susceptibility level in Figure 10 5 shows the seismic landslide susceptibility level and the observed seismic landslide developed in each susceptibility level of Wenchuan.After considering the rock mass strength heterogeneity, areas with high susceptibility level increased by 2.74%, and the observed landslides in high susceptible areas increased by 7.37%, compared with that ignoring the rock mass strength heterogeneity.Results suggest that when the rock mass strength heterogeneity is ignored in the permanent seismic displacement calculation, the calculated displacements may underestimate the seismic landslide susceptibility.Therefore, it is necessary to consider the rock mass strength heterogeneity in permanent-displacement analysis for seismic landslide assessment.

Case study: Ludian earthquake
To validate the transferability of the empirical model developed in the Wenchuan earthquake, this study selected Ludian which has a similar geo-environment to Wenchuan, as another case study.

Permanent seismic displacement calculation
Figure 11 shows the permanent seismic displacements calculated using different empirical displacement models with and without considering the spatial heterogeneity of the rock mass strength.Without considering the spatial heterogeneity of the rock mass strength, the displacement calculated with the RS08 model has the best prediction accuracy and validity, with the %GFC and %Efficient being 67.74% and 55.40%, respectively; while the Jib07B has the worst prediction accuracy and validity, with the %GFC and %Efficient being 52.93% and 44.16%, respectively.Unlike the Wenchuan earthquake, when the earthquake magnitude is considered in the empirical displacement model, the prediction accuracy and validity of the calculated displacement are significantly reduced (58.74% vs. 52.93%,67.74% vs. 56.02%).This may be due to insufficient consideration of ground motion acceleration for moderate earthquake magnitudes when constructing the empirical displacement equation.With considering the spatial heterogeneity of rock mass strength, the prediction accuracy and the validity of the displacements with different empirical regress models were significantly improved, especially for the RS08 model, where the prediction accuracy increased from 67.74% to 82.66%, and the validity increased from 55.40% to 67.30%.Results indicate that considering the spatial heterogeneity of rock strength in the permanent seismic displacement calculation, the calculated displacement can effectively predict seismic landslides.To this end, with considering the spatial heterogeneity of rock mass strength, the RS08 was selected to calculate the permanent seismic displacement in Ludian.

Landslide susceptibility assessment
Similar to the first case study, the calculated permanent seismic displacements in Ludian were divided into low, moderate, and high susceptibility levels, as shown in  Figure 12.We can find that considering the spatial heterogeneity of rock mass strength, the area with moderate to high susceptibility levels increases to some extent (e.g.subarea B) compared to that without considering spatial heterogeneity of the rock mass strength.Table 6 shows the seismic landslide susceptibility level and the observed seismic landslide developed each susceptibility level of the Ludian earthquake.After considering the spatial heterogeneity of the rock mass strength, areas fallen in high susceptibility level increased by 1.83%, and the observed landslides in high susceptible areas increased by 14.34%, compared with that without considering the spatial heterogeneity of rock mass strength.The result also indicates that when the rock mass strength heterogeneity is ignored in the permanent seismic displacement calculation, the calculated displacements may underestimate the seismic landslide susceptibility.Besides, as the result in Ludian earthquake, the established empirical model is also applicable to the seismic landslide susceptibility assessment with similar scenarios to Wenchuan, such as Ludian.

Discussions
Prior studies demonstrated that assigning rock mass strength is critical in permanentdisplacement analysis, and the estimation of spatial heterogeneity of rock mass strength is critical in making a significantly reliable regional seismic landslide susceptibility assessment (Dreyfus et al. 2013).Although the importance of spatial heterogeneity of rock mass strength has been recognized (Rathje and Saygili 2009;Gallen et al. 2015;Qi and Li 2018;Rodr ıguez-Peces et al. 2020), few studies have considered the influence of spatial heterogeneity of rock mass strength in permanent-displacement analysis due to the difficulty of obtaining detailed rock strength measurements with spatial location.Many scholars have quantified rock strength heterogeneity from different perspectives.For example, Gischig et al. (2016) quantifies the spatial heterogeneity of post-earthquake rock strength based on pre-earthquake damage of slopes and the amplification effect of ground shaking; Cudjoe (2020) quantifies the heterogeneity of rock strength based on engineering rock classification methods using various quantifiers (e.g.joints, fissures).These methods are mainly for specific scenario areas or single slopes, which are not universal for regional scale.On the contrary, this article considers the spatial heterogeneity of rock mass strength, which is caused by lithological environment, and quantifies the spatial heterogeneity based on the lithological environment difference.To this end, we first examined the influence of lithological environment difference on the rock mass strength heterogeneity.Secondly, we proposed a novel method to invert the rock mass strength heterogeneity from a well-documented landslide inventory.Results on two cases show that there is significant heterogeneity in rock mass strength even in the same lithology.Besides, the permanent-displacement model that considers the rock mass strength heterogeneity can more reliably assess the seismic landslide susceptibility, and the susceptibility levels are in line with the observed landslides.
One of the common methods for calculating regional permanent seismic displacements is using empirical displacement models.Various empirical displacement models have been established based on datasets of seismic events in different regions.For example, Delgado et al. (2020) compared and analyzed the applicability of several empirical displacement models in Spain.Results indicate that the selection of an appropriate empirical displacement model depends on the earthquake magnitude, or the critical acceleration used to calculate the displacement.In the same manner, this study compared four empirical displacement models in the permanent-displacement analysis method and selected the empirical displacement model with the best performance by compared with the observed landslides.This study again verifies that different empirical displacement models significantly affect the regional seismic landslide susceptibility assessment.On the other hand, we also found that for moderate to high earthquake magnitude, reasonable consideration of magnitude in the empirical displacement model can improve the reliability of seismic landslide susceptibility assessment, while reduce its reliability for moderate to low earthquake magnitude.
Although this study improves the reliability of the permanent displacement-based approach for assessing seismic landslide susceptibility by considering the spatial heterogeneity of rock strength, some limitations of this study must be acknowledged.
The empirical model to estimate the spatial heterogeneity of rock mass strength was developed under specific geo-environmental condition and earthquake event.When applying the developed method to other regions, the similarity of the geo-environmental conditions and earthquake events with this study must carefully check.Second, source mechanism of earthquake also partially controls the seismic landslide distribution.This article did not consider these two factors due to the data availability.Third, we only established an empirical model of the rock mass strength heterogeneity on a single earthquake event that limits its application to a certain area.To broaden the scope of application, we may consider regional/global earthquake events to further study the rock mass heterogeneity.

Conclusions
This study proposes a permanent displacement model that takes into account the rock mass strength heterogeneity for regional seismic landslide susceptibility assessment.First, we developed a novel method to quantify the rock mass strength heterogeneity that results from different lithology environmental factors, and an empirical model was designed to estimate the rock mass strength heterogeneity based on the landslide inventory of Wenchuan earthquake.Second, we consider the rock mass strength heterogeneity in the calculation of critical acceleration of the permanent-displacement model.Finally, combining the critical acceleration and ground motion parameters, the permanent seismic displacement and seismic landslide susceptibility were derived.Results showed that the spatial heterogeneity significantly affects regional rock mass strength that in turn influences the seismic landslide susceptibility.The permanent-displacement model that considers the rock mass strength heterogeneity is more suitable for assessing seismic landslide susceptibility than that without considering the rock mass strength heterogeneity.

Figure 4 .
Figure 4. Maps of lithological groups of rock mass strength in (a) Wenchuan and (b) Ludian earthquakes.

Figure 6 .
Figure 6.Fault and river data used in Wenchuan and Ludian earthquakes.(a) Distance to fault and (c) river density in Wenchuan; (b) Distance to fault and (d) river density in Ludian.

Figure 7 .
Figure 7. Infinite-slope model under static conditions, where t, m, h, c 0 , u 0 c w and c represent depth, saturation, weight, slope angle, effective cohesion, effective friction angle, unit weight of the water, and unit weight of the sliding-block, respectively.

4. 1 .
Case study: Wenchuan earthquake 4.1.1.Permanent seismic displacement calculation Figure 8(a) shows the estimated rock mass strength heterogeneity using the empirical model developed in Wenchuan.As shown in Figure

Figure 8 .
Figure 8. Spatial heterogeneity of the rock mass strength and the critical acceleration in Wenchuan.(a) the spatial heterogeneity of rock mass strength in Wenchuan estimated using the empirical model, (b) the critical acceleration in Wenchuan calculated without considering the spatial heterogeneity of rock mass strength, and (c) the critical acceleration modified after considering the spatial heterogeneity of rock mass strength.

Table 4 .
List of empirical displacement models used in this study for permanent seismic displace- (b) are increased significantly compared to Figure 10(a), see subarea B. Table

Figure 9 .
Figure 9. Permanent seismic displacements in Wenchuan calculated using different empirical displacement models with and without considering the spatial heterogeneity of rock mass strength.

Figure 10 .
Figure10.Seismic landslide susceptibility assessment of Wenchuan earthquake using the permanent-displacement analysis method without considering the spatial heterogeneity of rock mass strength (a) and with considering the spatial heterogeneity of rock mass strength (b).

Figure 11 .
Figure 11.Permanent seismic displacements in Ludian calculated using different empirical displacement models with and without considering the spatial heterogeneity of rock mass strength.

Figure 12 .
Figure12.Seismic landslide susceptibility assessment of Ludian earthquake using the permanentdisplacement analysis method without considering the spatial heterogeneity of rock mass strength (a) and with considering the spatial heterogeneity of rock mass strength (b).

Table 1 .
Data sources used in Wenchuan and Ludian.

Table 2 .
Rock mass strength parameters assigned in different lithological groups in Wenchuan and Ludian.

Table 5 .
Seismic landslide susceptibility levels and the observed seismic landslide developed in each susceptibility level of the Wenchuan earthquake.

Table 6 .
Seismic landslide susceptibility levels and the observed seismic landslide developed in each susceptibility level of the Ludian earthquake.Without considering the spatial heterogeneity of rock mass strength