Large strike-slip earthquakes are marked generally by long surface ruptures and relatively narrow deformation zones (Xu et al., 2006; Choi et al., 2018; Antoine et al., 2022; Nurminen et al., 2022). The displacement on a strike-slip fault is always estimated by empirical relations between magnitude and coseismic displacement based on surface rupture measurements of large earthquake cases (Wells and Coppersmith, 1994; Cheng et al., 2020; Shaw, 2023). However, detailed mapping of surface ruptures of recent strike-slip earthquakes demonstrates that the coseismic displacement of surface ruptures occurs on the main fault plane and associated secondary structures (e.g., Antoine et al., 2022; Liu-Zeng et al., 2024). The displacement on the secondary structures is commonly accommodated by the distributed ruptures in a wide zone (Figure 1). In some cases, the role of distributed displacement is more than that on the main fault (Antoine et al., 2022; Liu-Zeng et al., 2024). This type of distributed displacement caused a wider surface rupture zone and affected the distribution of coseismic damages and earthquake-triggered geological hazards. It is urgent to assess the distributed displacement hazard of strike-slip faults for seismic safety of large linear engineering and seismic disaster prevention and mitigation.

Increasing studies in the past 20 years have been conducted on evaluating the displacement from active strike-slip faults (Lee and Trifunac, 1995; Peter, 2010; Inoue et al., 2020; Nurminen et al., 2020; Visage et al., 2023). A Probabilistic fault displacement hazard analysis (PFDHA) framework was first proposed based on permanent fault displacements and the classical probabilistic seismic hazard analysis was applied to assess the seismic risk of the Yucca Mountain nuclear waste disposal project in the United States (Youngs et al., 2003). When more attention focused on the distribution of surface displacement, six types of regression curves were used to fit the coseismic surface slip distribution of historical earthquakes and the spatial distribution pattern of surface displacement of large earthquakes along their seismogenic fault strike were explored (Wesnousky, 2006; Wesnousky, 2008). Since many related studies had been carried out in America, researchers tried to conduct similar method in other regions, using the actual seismic data in Japan with the PFDHA method to establish a surface-rupturing fault displacement prediction model which is applicable to Japan (Takao et al., 2013). But these studies do not consider the effect of distributed ruptures on the spatial distribution of coseismic displacement along a strike-slip fault. A unified method for calculating on-fault and off-fault displacement was presented based on fault sections with different geometry and structures. The off-fault displacement is similar to distributed displacement and may include some deformation that does not generate surface ruptures (Antoine et al., 2022). In addition, the distribution curve of coseismic displacements perpendicular to the fault strike is inconsistent with the specific coseismic offset data of large earthquakes (Petersen et al., 2011). Further work tried to mainly focus on distributed ruptures, the distance from the main fault was taken into account in a probability model to predict the displacement distribution of distributed rupture crossing the fault traces using surface rupture data from five strike-slip earthquakes in northern America. These five earthquakes span a range of magnitudes between M W 6.4 and 7.3, four of which are from M W 7.1 to 7.3 (Rodriguez Padilla and Oskin, 2023). However, the applicability of this model to strike-slip faults in other areas is unknown.

In western China, the increasing seismic potential of many strike-slip faults poses a high risk of the displacement of linear engineers across active faults. Especially the 2022 Menyuan earthquake ( M W 6.6) offset the railway tunnel (Li et al., 2023), stimulating the displacement hazard of distributed ruptures along strike-slip faults. The probability of permanent displacement was estimated across faults of the second line of West-East natural gas transmission pipeline by using the potential focal region parameters of China’s ground motion parameter zoning map (Zhao et al., 2008). Based on the surface-rupturing data of strike-slip faults in China, both parabolic and elliptic prediction equations related to surface displacement and surface rupture length were presented. Similarly, the permanent displacement risk curve of surface-rupturing zones of the Zemuhe active fault zone was obtained combined with the practical application of PFDHA method (Jin, 2019). However, these earlier studies only use the uniform models and do not consider the application of these models in China. Also, the distributed ruptures along a strike-slip fault are not taken into account.

In the Tibetan Plateau, several large earthquakes along large strike-slip faults have been well studied, and associated detailed surface ruptures have been mapped (Figure 2) (Li et al., 2012; Ren and Zhang, 2019; Yuan et al., 2021; Ren et al., 2022; Li et al., 2023; Liu-Zeng et al., 2024), providing an opportunity to construct the model of displacement hazard of distributed ruptures from strike-slip faults in China. In this study, we collected the recent surface rupture data of five strike-slip earthquakes in the Tibetan Plateau and built the probability model to predict the displacement hazard of distributed rupture crossing the fault trace when a certain magnitude occurs. This study not only proposes a workflow for the probabilistic displacement hazard of distributed ruptures from strike-slip faults but also supports the seismic hazard assessment of linear engineering crossing strike-slip faults in the Tibetan Plateau.

The Tibetan Plateau, controlled by the Cenozoic Indian-Eurasian collision, forms the most intense area of tectonic activity in China. Accompanied by the uplift of the plateau, the material began to the outward extrusion and produced many large strike-slip faults, such as the Kunlun, Altyn Tagh, Haiyuan, Ganze-Yushu-Xianshuihe faults (Molnar and Tapponnier, 1978; Xu et al., 2005; Ren et al., 2013). These faults undergo the late Quaternary slip rate of ∼5–10 mm/yr and have ruptured in serval large earthquakes in the past 20 years (Kirby et al., 2007; Ren et al., 2013; Yao et al., 2019).

We choose five strike-slip earthquakes in the northern Tibetan Plateau: the 1997 Mani ( M W 7.5), 2010 Yushu ( M W 6.9), 2014 Yutian ( M W 6.9), 2021 Maduo ( M W 7.4), and 2022 Menyuan ( M W 6.6) earthquakes (Figure 2) (Table 1). The original datasets are mainly obtained by scanning published papers or from relevant experts and researchers. When the specific methods collecting raw data may vary from person to person, field measurement and optical image correlation are usually chosen in current studies. In addition to the basic data of the five earthquakes in TableX, since they happened in different regions and seasons along with other possible factors, the actual geological and structural characteristics and the surface condition of every earthquake are significantly different, resulting in different surface fractures.

The 1997 Mani earthquake ruptured the left-lateral Mani fault, a possible westernmost segment of the Kunlun fault (Shan et al., 2006; Ren and Zhang, 2019). This earthquake occurred in November 1997 when the temperature was below freezing point and the epicentral area was covered by permafrost. The rupture zone of the Baixue Lake section is located in the lacustrine plain and passes through a series of large-scale ice water alluvial fans, and the Chaoyang-Shuangduan Lake section created a fault-plug pond formed by seismic steep ridges blocking the flow of melting ice and snow (Ren and Zhang, 2019).

The 2010 Yushu earthquake occurred on the Yushu segment of the Ganz-Yushu fault (Li et al., 2012; Sun et al., 2012). Located in the alpine region at a high altitude, Longbao Lake and the surrounding water were frozen during this earthquake, leading to a unique surface rupture trace in such frozen regions: ice cracks. The ice fissure zone is distributed in a planar shape and large scale, accompanied by sand liquefaction, showing a normal component dipping the Longbao Lake (Sun et al., 2012).

The 2014 Yutian earthquake mainly ruptured the left-lateral strike-slip Altyn Tagh fault (ATF) system with normal-slip (Li et al., 2016; Yuan et al., 2021). The earthquake vibrations near the Xiaoerkule Lake caused the instability of the diluvial fan, and the landslide and graben system formed after the slide of the salty lake’s shoreline due to gravity, which belongs to the associated shallow deformation. This earthquake was located in the inner flow area of the Tibetan Plateau, in the middle of the rupture zone a thick loose sedimentary cover was formed based on the lake, which would magnify the seismic rupture in the shallow surface and enlarge the surface rupture width (Yuan et al., 2021).

The 2021 Maduo earthquake mainly ruptured the left-lateral strike-slip fault with a partly normal slip, the southeastern branch of the Kunlun fault zone, the Jiangcuo fault (Ren et al., 2022). This earthquake ruptured across a long-range area with complex landform, including many valleys, lakes, wetlands, grasslands, mountains and some sand dunes, etc (Pan et al., 2021; Ren et al., 2022). In the rupture section across Mustan Bridge, north of Huanghe Township and the valley of Yellow River, as the rupture zone is mainly distributed along the valley, a large number of tensile fractures parallel to the river channel are formed due to slope instability, accompanied by sandblasting water and sand liquefaction (Liu-Zeng et al., 2024).

The 2022 Menyuan earthquake mainly ruptured the left-lateral strike-slip fault with a slight thrust component, the middle segment of Qilian-Haiyuan fault (Xue et al., 2022; Li et al., 2023). This earthquake happened in winter when the Liuhuanggou River was frozen, causing the ice river surface and flood plain to break in the nature of thrust with the phenomenon of tension cracks on the bank slope and slope slide under the effect of gravity from slope instability on the right bank of the river, which may amplify the actual left-lateral offsets (Li et al., 2023).

With the continuous development of western China, a large number of pipeline projects have been built or are being built in the Tibetan Plateau region. The damages of strike-slip faulting can be found in these five earthquakes selected in this study: Maduo earthquake caused the serious collapse of Mustan Bridge as part of the Gongyu Expressway, which is a major traffic route (Ren et al., 2022); The dislocation of the Lanzhou-Urumqi High-speed Railway’s Daliang Tunnel caused by Menyuan earthquake directly disrupted the traffic in this section for up to 18 months (Li et al., 2023); Yushu earthquake caused serious damage to the national and provincial trunk roads and transportation infrastructure in the Yushu area, with an affected area of 20,000 square kilometers (Li et al., 2012). Therefore, predicting the distribution of surface rupture caused by strike-slip earthquakes in the Tibetan Plateau region can help maintain the safety of transportation arteries in the region and provide reference for avoidance locations for transportation arteries currently being built or planned to be built in the region, in order to reduce potential damage in the future.

Rodriguez Padilla and Oskin (2023) proposed a mathematical approach to estimate fault displacement based on data from detailed surface-ruptures strike-slip earthquakes: P S > S 0 | x , M W = P r u p t u r e | x P S > S 0 | x , r u p t u r e , M W . (1)

Equation 1 shows the probability per unit area of finding a distributed rupture that accommodates a displacement that exceeds a displacement threshold at a given distance from the principal fault. The equation is produced by the joint probability of two parts: the former P r u p t u r e | x is the probability of rupture occurring per unit area at a distance x from the fault; The latter P S > S 0 | x , r u p t u r e , M W is the probability of finding a displacement that exceeds the threshold at a given distance from the fault for a given earthquake magnitude, given the presence of a rupture. In the following paragraphs, the process of establishing these two parts will be carried out separately before being combined into the final equation (Rodriguez Padilla and Oskin, 2023).

P r u p t u r e | x can be calculated by studying the spatial distribution of the rupture density, which can be given by the inverse power law (e.g., Padilla et al., 2022): V x = V 0 x + x f r x f r − γ . (2)

In which V 0 is the rupture density at the origin in number of ruptures per unit area ( 1   m 2 , x f r is a normalized factor which is related to the uncertainty of the location of the principal fault trace in meters, the exponent γ is the slope of the decay of rupture density with distance for values of x ≫ d in log-log space or scaling exponent. V x is the probability of a rupture occurrence per unit 1   m 2 .

With reference to the inverse power law used in Equation 2, we try to analyze the relationship between the displacement measurements and the distance to the principal fault trace: λ x = β x + x S x S − n . (3)

In which λ is the mean of the displacement at every distance bin, β is the average displacement at the origin, x is the location away from the principal fault in meters, x S is a normalized factor which we set as 1   m in this equation, and n is the slope of the relationship between mean displacement and distance in log-log space or scaling exponent.

To deeply understand the meanings of above parameters, we put them back in actual situation to help with this. The parameter V 0 describes the amount of the surface ruptures per unit area on the principal fault trace, so as the parameter β which describes the average displacement value on the principal fault trace. These two parameters are used to describe the surface-ruptures and the displacement on the principal fault. Unlike V 0 and β , the parameter γ and n are both the slope of the best-fitting curve in log-log space which can describe the state of decay.

By comparing several distribution functions, it is found that the population of displacement measurements can be better described by exponential distributions, and the distribution can be described as follows: f S | x = 1 λ e − S λ . (4)

Combining Equations 3, 4 yields: f S | x = 1 β x + x S x S n e − S β x + x S x S n . (5)

This equation is a probability density function (PDF) of displacement measurements with distance x from the principal fault trace. To solve the probability of finding a displacement that exceeds S 0 with a given earthquake magnitude: P S > S 0 | x , r u p t u r e , M W = ∫ S 0 S max 1 β x + x S x S n e − S β x + x S x S n d S = − e − S β x + x S x S n | S 0 S max (6)

Since the model focuses on the displacements on the distributed ruptures which do not tend to be of large value, the evaluation limits the threshold S 0 ≪ S max to appropriately predict the probability of distributed displacements above the threshold S 0 . Completing the integration of Equation 6 with this application can yield: P S > S 0 | x , r u p t u r e , M W = e − S 0 β x + x S x S n . (7)

Combining Equations 2, 7 yields the final model: P S > S 0 | x , M W = V 0 x + x f r x f r − γ e − S 0 β x + x S x S n . (8)

This model is presented to calculate the probability per unit area of finding a distributed rupture that accommodates a displacement that exceeds a displacement threshold at a given distance from the principal fault (Rodriguez Padilla and Oskin, 2023). This probability model can quantify the displacement hazard of surface ruptures and estimate the expected displacement distribution caused by strike-slip faults, to provide reference data for the design, evaluation, and maintenance of engineering structures and lifelines located near or across strike-slip faults.

The entire probability model involves the calculation of multiple unknown parameters, which will all influence the final results. Therefore, after clarifying the model-building process, a brief analysis of the characteristics of each parameter is needed, and then the expected results of the theoretical model should be summarized as a reference for further research.

The origin rupture density V 0 and the origin average displacement β reflect the rupture and displacement conditions along each principal fault trace of the earthquake, and further affect the establishment of their respective models. A larger V 0 and β will result in higher P S > S 0 at the same location. The influence of parameter V 0 on P S > S 0 is greater than that of parameter β when S 0 = 0.01 , 0.1 m, the probability curve with a larger V 0 gradually dissociates above other curves with the increase of distance x ; When the threshold S 0 is large, the influence of parameter β on P S > S 0 becomes prominent, which actually reflects the magnitude dependence of parameter β : the small magnitude earthquake event tend to have a small β , and its probability curve gradually dissociates below other curves with the increase of distance x .

For the magnitude dependence of parameter β , the equation was proposed between moment magnitude M W and mean displacement β to explore the empirical relationship between them (Wells and Coppersmith, 1994): log 10 β = b M W − a , (9) where parameters a and b are regression coefficients derived from the best fit of the data using the least square method, and all fitted data used are strike-slip earthquakes. When S 0 = 0.5 m, the P S > S 0 of the smallest magnitude earthquake is significantly smaller than that of other earthquakes. In other words, for strike-slip earthquakes with smaller magnitude, the probability of large rupture occurring at a certain distance from the fault will be significantly lower than that of earthquakes with larger magnitude, which also matches our prediction about the characteristics of the surface rupture displacement distribution for strike-slip earthquakes.

According to the definition, the values of the exponents γ and n obtained by fitting are affected by the rupture density and average displacement of each earthquake respectively, but there is a reasonable range of values. After obtaining the respective parameters and models of each earthquake, a theoretical model is established. The parameters of the theoretical model are obtained by combining all the seismic data used and fitting through Equations 2, 3, and the range of the exponents γ and n obtained by the joint fitting will be affected by all the seismic data used.

The reasonable outputs of theoretical model should also be summarized after the analysis of the characteristics of each parameter. The reasonable results obtained after the establishment of the theoretical model are consistent. The rupture density curve from the fit of Equation 2 shows a downward trend, with certain fluctuations at the end; The scatter-fitting curve of displacement measurements from Equation 3 also showed a downward trend, and the scatters concentrated around the best fitting curve with a little isolated from the entirety.

Then comes to the final model output by Equation 8, it can be seen from the curve that the theoretical P S > S 0 decreases with the increase of fault-perpendicular distance, and the slope keeps increasing. Set with the input of moment magnitude M W = 7 and threshold S 0 = 0.1 m, P S > S 0 has dropped to about 1/1,000 when the distance reaches about 100 m, and at 1,000 m it has dropped to less than 1/1,000, reaching a relatively small risk probability value. Such results indicate that there is a very low probability for the surface ruptures in the area far away from the principal rupture trace to accommodate a large value of displacement.

In this study, we introduce the approach of Rodriguez Padilla and Oskin (2023) to the surface rupture dataset of large strike-slip earthquakes in the Tibetan Plateau and then build the displacement probability distribution model. By analyzing the model building process and output results of the theoretical model, we can summarize the reasonable results of the theoretical model and use it as a reference to guide the subsequent model-building work in China.

Detailed analysis of the surface rupture dataset indicates that the coseismic surface rupture zone of a large strike-slip earthquake consists primarily of principal ruptures along the main fault and secondary ruptures caused by the secondary fault or fault branch (Figure 1). In addition, minor fractures occurred outside of the main fault zone. These fractures are related to non-tectonic factors such as landslides, and gravity instability under seismic ground motion (Liu-Zeng et al., 2024). These minor fractures are generally located in the local environment and are not directly produced by fault displacement. Therefore, we should remove these minor fractures, and only tectonic ruptures are involved in our analysis.

After preliminary sorting of the five original earthquake data (map documents and statistical tables of displacement measurements at the rupture), it is found that the original data will record the surface rupture traces with event particularity caused by the earthquake itself due to structural characteristics and some other causes (Li et al., 2012; Li et al., 2016; Ren and Zhang, 2019; Pan et al., 2021; Xue et al., 2022). Since the model is expected to be applied to a wide range of the Tibetan Plateau, to reduce the impact of non-tectonic factors (such as seismic motion) on the model results, we formulate criteria for the screening of the original data combined with the rupture characteristics of the five typical earthquakes on the Tibetan Plateau, hoping to get more consistent data to obtain a more uniform result and model.

For the convenience of data storage and search, we store all the data in a standardized format after data filtering. After comparing the data processing methods of standardized empirical databases that have been generated in recent years (e.g., Ancheta et al., 2013; Chiou et al., 2008; Nurminen et al., 2022), we chose to refer to the Fault Displacement Hazards Initiative (FDHI) database hosted and maintained by the Natural Hazards Risk and Resilience Research Center at the University of California, Los Angeles, to conduct standardized processing of surface rupture data form these five earthquakes in China (Sarmiento, 2021). The buildup of the model in this study mainly needs two parts of data: surface rupture maps and surface rupture displacement measurement, which are sorted into ESRI Shapefile (*.SHP format) and table form, respectively, after standardized processing.

Based on the original data collected from these five earthquakes, we first screened the surface rupture traces based on the following four criteria (Figure 3):

(1) Through image comparison, the rupture traces around the water body are generally interpreted as resulting from the instability of the water body under seismic ground motion. These ruptures strike along the boundary of the water body and do not have a similar strike to the main fault. Such fractures are removed.

After completing the standardized filtering of the original surface rupture traces maps, the displacement measurements located on the removed rupture traces are also deleted to ensure the uniformity of the data used. Then we unified the format of surface rupture displacement measurements according to that of FDHI (Sarmiento, 2021). Currently, the displacement measurements for each earthquake as recorded by various scholars are typically documented and reported separately in the original data, including horizontal slip and vertical component (e.g., Li et al., 2012; Li et al., 2016; Pan et al., 2021; Li et al., 2023). It is important to note that we focus on the horizontal slip of a strike-slip rupture and the influence of other factors causing vertical slip in a local area is not considered in this study. Therefore, all original displacement measurements should be classified before calculation and organized into a uniform representation of net horizontal slip.

After all original data have been preprocessed, we conduct the following work such as quantifying the fracture density and describing the variation in displacement measurements at the fracture point with the distance from the main fault. We integrated the FDHI database’s definitions of principal rupture and the location of principal faults in each seismic event, and consequently, a simplified principal rupture trace was delineated (Figure 4). A singular, continuous fault trace was defined along the densely distributed and consistent principal ruptures. This approach can enhance the convenience of data processing in the subsequent model development by utilizing a limited number of rupture trace lines.

After pre-processing of the original surface rupture data of these five earthquakes, we built our model based on these data. Firstly, the spatial distribution of surface rupture density is studied based on Equation 2. When all the rupture traces have been separated into segments at a 1 m interval, we calculated the smallest distance from each segment to the principal rupture trace, and fit the attenuation of surface rupture density for each earthquake based on these data (Figure 5).

When the surface rupture density distribution for each earthquake has been obtained, we fit each decay with an affine-invariant ensemble sampler for Markov chain Monte Carlo (Goodman and Weare, 2010; Foreman-Mackey et al., 2013), setting the prior ranges of the values of the three unknown parameters in Equation 2 to estimate the maximum-likelihood values by using the actual surface rupture data of each earthquake (Table 2).

We follow a similar method to estimate the attenuation of displacement distribution with the fault-perpendicular distance in Equation 3: With the smallest distance of each displacement measurement point to the principal rupture trace calculated, we calculate the attenuation of the surface rupture displacement measurements for each earthquake (Figure 6).

When the decay of displacement measurements with fault-perpendicular distance for each earthquake has been obtained, we set the prior ranges of the values of the two unknown parameters in Equation 3 to estimate the maximum-likelihood values by using the actual surface rupture data of each earthquake (Table 3). Then, we build the probability model for each of the five earthquakes after obtaining the parameters for the final model (Figure 7).

As the individual models are built, we fit the parameters of the general model based on all the surface rupture data from these five earthquakes according to the above method (Table 4). About the value of β , according to Equation 9 we give the input of M W . Using the parameters in Table 3, we build the general model for the Tibetan Plateau and set different values of M W , S 0 to observe the rationality of the results (Figure 8). In the establishment of the model, we imported the code into Jupyter Notebook, enabling end users to input the displacement threshold S 0 and earthquake moment magnitude M W to obtain outputs for P S > S 0 (Figure 9).

Based on a systematic analysis of surface ruptures from the five strike-slip earthquakes in the Tibetan Plateau in the past decades, we construct a general probabilistic model for displacement hazards associated with distributed ruptures from strike-slip faults. The following conclusions can be drawn:

(1) This study summarizes a workflow for building localized models of displacement hazard for active faults based on surface ruptures from large earthquake cases (Figure 12).