Prediction of liquefaction-induced lateral spreading based on Neural network

2.1. Earthquake motion database

A compilation of seismic events has been meticulously collected from various countries, including China, the United States, Canada, Iran, Turkey, and Japan, sourced from the Pacific Earthquake Engineering Research Center [17]. This collection encompasses 27 significant seismic events, among which notable examples are the 1995 Kobe earthquake in Japan, the 1999 Chi-Chi earthquake in Taiwan, and the 1989 Loma Prieta earthquake, amounting to a total of 1960 seismic records. The seismic motion database used for the analysis in this research is detailed in Table 1. Concurrently, corresponding seismic parameters for each event have been systematically gathered and organized. These encompass variables such as PGA, $M_w$, $V_s^{30}$, FM and $R_{rup}$.

2.2. Newmark sliding block method

During the soil liquefaction, once liquefaction is triggered and the sliding surface begins to form, the strength of saturated sandy soil rapidly diminishes. As the soil strength continues to decrease, the residual shear strength of the soil remains within the liquefied site. When the shear forces induced by the earthquake, combined with the gravitational forces acting on the soil above the sliding surface, exceed the residual shear strength, the soil above the sliding surface initiates movement. This movement accumulates displacement continuously until the earthquake subsides. The motion pattern of the soil above the sliding surface can be described using the Newmark sliding block method. Specifically, by performing two integrations of the acceleration exceeding the $k_y$, the resulting displacement under seismic conditions can be determined. The ultimate sliding displacement is the cumulative result of multiple increments of sliding and is illustrated in Fig. 2. This approach has been employed by several scholars to predict lateral spreading induced by liquefaction, including Yang [16], Baziar [18], Taboada [19], Olson, and Johnson [20].

Fig. 2Diagram of dynamic permanent displacement calculated by Newmark sliding block method

From Fig. 2, it is evident that the $k_y$ is a pivotal parameter in the computation of displacement using the Newmark sliding block method. Therefore, to obtain a substantial dataset of liquefaction-induced lateral spreading calculations suitable for machine learning algorithms, this section assumes various values for the $k_y$, 0.01 g, 0.03 g, 0.05 g, 0.75 g, 0.1 g, 0.15 g, 0.2 g, 0.25 g, and 0.3 g. The selection range for the $k_y$ is determined based on the analysis of corresponding models using 23 well-documented cases of liquefaction-induced lateral spreading, as compiled by Yang et al. [12]. It is important to note that the computation of $k_y$ necessitates the establishment of a limit equilibrium analysis model based on parameters such as soil distribution, soil strength, soil density, and groundwater levels [16]. Consequently, by employing the $k_y$, the information encompassing the parameters is inherently considered.

Subsequently, using the collected dataset of 1960 seismic records, lateral spreading under different yield accelerations is computed. Note, before calculating the liquefaction-induced lateral spreading, records with a $k_y$ exceeding the PGA are excluded, as the calculated lateral spreading would be zero. In the end, a total of 8914 valid liquefaction-induced lateral spreading calculations are obtained.

3.1. Parameter selection

Based on the findings from the studies mentioned in references [8] and [14], this research selects six seismic parameters, PGA, $k_y$, $M_w$, $V_s^{30}$, FM, and $R_{rup}$, as input features, with $\mathrm{l}\mathrm{n}\left(D_n\right)$ as the output feature. The model used for training the lateral spreading parameter sensitivity analysis is expressed in the following functional form:

where, $D_n$ represents the lateral spreading in mm. PGA stands for peak ground acceleration, measured in units of g, $k_y$ denotes the yield acceleration, also in units of g, $M_w$ signifies the moment magnitude, $V_s^{30}$ represents the average shear wave velocity in the top 30 meters of soil, measured in m/s, FM corresponds to the focal mechanism, which includes three forms, i.e., Normal, Reverse, and Strike-Slip faults, and $R_{rup}$ signifies the rupture distance, measured in km.

Within the dataset of 1960 seismic records used in this research, the ranges of the six seismic parameters mentioned above are depicted in Table 2. It is important to note that the $k_y$ is assumed, and in the FM category, 1 corresponds to Normal fault, 2 to Reverse fault, and 3 to Strike-Slip fault. The distribution between PGA and $M_w$ is illustrated in Fig. 3(a), while the distribution among $R_{rup}$, $M_w$, and FM is presented in Fig. 3(b).

Fig. 3Distribution of relevant parameters of seismic records [21]

a) Distribution of PGA and ${\mathrm{M}}_{\mathrm{w}}$

b) Distribution of ${\mathrm{R}}_{\mathrm{r}\mathrm{u}\mathrm{p}}$, ${\mathrm{M}}_{\mathrm{w}}$ and FM

3.2. Neural network model

A neural network is a type of machine learning algorithm designed to emulate the structure and functionality of the human neural system. It comprises a hierarchical structure of multiple neurons and is utilized for prediction and classification tasks by learning patterns and features from input data. Neural networks optimize model parameters through forward and backward propagation algorithms to minimize prediction errors. A neural network model consists of an input layer, hidden layers, and an output layer. The input layer primarily receives and passes input features to the hidden layers. The hidden layers play a key role in feature extraction and transformation. Finally, the output layer provides the responses of the predictive model.

In this research, there are six input neurons and one output neuron. Additionally, according to the Universal Approximation Theorem, a single hidden layer is sufficient to describe continuous nonlinear functions [14]. Therefore, a single hidden layer is employed. Considering the nonlinear relationship between input and output features, the sigmoid function is chosen as the activation function to produce predicted values that closely approximate input features. The operation stops when the output signal meets the desired criteria; otherwise, error backpropagation is executed. This involves feeding back error values to the hidden layer and continuously adjusting network weights, allowing the hidden layer to gain strong nonlinear mapping capabilities and ultimately produce the desired output. Note, to minimize scale differences between various features and enhance data uniformity, the input features are scaled to a range between -1 and 1 according to Eq. (2):

where, $a=$ –1 and $b=$ 1, and $x_{min}$ and $x_{max}$ takes the minimum and maximum values of the various features that are mapped.

Based on this framework, 80 % of the data is randomly allocated for training, while the remaining 20 % forms the testing set. Bayesian regularization is employed for training the data. After multiple iterations, the regression achieves the best results with 10 hidden neurons. Thus, the artificial neural network structure employed in this research is depicted in Fig. 4.

The neural network structure model depicted in Fig. 4 was utilized for training using the dataset from Section 1. Fig. 5 illustrates the performance concerning six parameters related to lateral spreading. Impressively, across the training set, testing set, and the overall assessment, the predicted results consistently exhibit a robust correlation coefficient ($R^2$) of 0.94 or higher. This high and consistent $R^2$ strongly signifies the predictive reliability of the model.

Fig. 4The neural network structure for lateral spreading parameter sensitivity analysis model

Fig. 5Performance of sensitivity analysis model for liquefaction-induced lateral spreading

3.3. Sensitivity analysis

Considering the challenges in acquiring the six parameters for engineering applications, this study aims to enhance the applicability of the lateral spreading prediction model by conducting a sensitivity analysis. This analysis will compare the influence of the six parameters on the prediction outcomes and identify the three most sensitive parameters to further refine and optimize the model, enhancing its suitability for practical use.

Therefore, all parameters, i.e., PGA, $k_y$, $M_w$, $V_s^{30}$, FM, and $R_{rup}$, were multiplied by scaling factor of 1.2 and 1.5. Using this modified data, a secondary prediction was conducted with the model, and the results were contrasted with the original predictions. Here, taking the PGA as an example is pivotal for illustration. Initially, the PGA undergoes scaling, being amplified by factors of 1.2 and 1.5. Subsequently, leveraging the sensitivity analysis model from Fig. 5, predictions for lateral spreading are made considering the 1.2 and 1.5 scaled PGA while keeping other data constant. Finally, based on the initial lateral spreading predictions, the $R^2$ and RMSE (Root Mean Square Error) are computed separately. This process reveals the $R^2$ and RMSE for each parameter under varying scaling factors. By comparing $R^2$ and RMSE, it becomes possible to assess the sensitivity of the six parameters. This analysis yields the sensitivity of each input feature, as presented in Table 3.

From the variation of $R^2$ and RMSE values in Table 3. It is evident that among the six input parameters, $M_w$, $k_y$, and PGA yield the most significant influence on liquefaction-induced lateral spreading. As the energy intensity of an input earthquake is often described by $M_w$, and PGA, additionally, $k_y$ emerge as a critical factor within the Newmark sliding block frame, representing the input and yield acceleration, respectively. Thus, the three parameters significantly influencing the analysis of lateral spreading [1]:

where, $N$ represents the number of observed samples, $y_{tar,i}$ denotes the $i$th target value, and $y_{pre,i}$ denotes the $i$th predicted value. RMSE measures the magnitude of deviation between the target values and the predicted values. It reflects the discrepancy between the predicted and actual values.

4.1. Prediction model of liquefaction-induced lateral spreading

Building upon the original data and the input feature sensitivity analysis results from Table 3, the prediction model is constructed using $M_w$, $k_y$, and PGA as input features, with lateral spreading calculated values as the output feature. The functional form of the predictive model is expressed as shown in Eq. (4):

where, it can be observed that this research employs a neural network structure with three input neurons and one output neuron.

Similarly, to account for the nonlinear relationship between seismic parameters and liquefaction-induced lateral spreading, the Sigmoid function is chosen as the activation function for the hidden layer to produce predicted values that approximate the input features. Regarding dataset splitting, the data is randomly divided into an 80 % training set and a 20 % testing set. Bayesian regularization is employed for training the data. After multiple attempts, it was found that the optimal regression performance is achieved with seven hidden neurons.

Hence, employing the artificial neural network structure illustrated in Fig. 6, the model was trained to predict lateral spreading. The training results for the lateral spreading prediction model are depicted in Fig. 7. As seen in Fig. 7, whether for the training set, testing set, or overall analysis, the predicted results demonstrate a $R^2$ of 0.93 or higher. This signifies a high level of reliability in the predictive liquefaction-induced lateral spreading of the model.

Fig. 6Neural network structure for liquefaction-induced lateral spreading prediction model

Fig. 7Performance of model for liquefied-induced lateral spreading

4.2. Comparison with empirical lateral spreading prediction models

Upon comparing Fig. 5 and Fig. 7, it becomes evident that the lateral spreading prediction model depicted in Fig. 7 maintains predictive accuracy while halving the number of input parameters, thereby considerably streamlining the model usability. Furthermore, to objectively assess the reliability of the proposed lateral spreading prediction model in this research, a comparative analysis with two established empirical models [4]-[5] that employ identical input variables ($M_w$, $k_y$, and PGA) was conducted. This comparative analysis serves to offer additional insights into the model proposed in this research.

It is worth noting that the Jibson 2007 model [4] is applicable for 5.3 ≤ $M_w$ ≤ 7.6. Therefore, to align with the constraints outlined in references [4]-[5] and the proposed model in this research, the specified ranges for $M_w$, $k_y$, and PGA are as follows, 5.3 ≤ $M_w$ ≤ 7.6, 0.01 ≤ $k_y$ ≤ 0.30 g, 0.015 ≤ PGA ≤ 1.779 g, as depicted in Fig. 8.

1. Jibson 2007 model [4], referred to as J07:

2. Rathje and Saygili 2009 model [5], abbreviated as RS09:

In Eqs. (5-6), $D_n$ represents the lateral spreading in cm.

Fig. 8Comparison of the lateral spreading predictions model with the existing empirical prediction models [21]

a) Variation with respect to $M_w$ with $k_y=$ 0.15g, PGA = 0.3 g

b) Variation with respect to $k_y$ with $M_w=$ 7.0, PGA = 0.3 g

c) Variation with respect to PGA with $M_w=$ 7.0, $k_y=$ 0.15 g

From Fig. 8(a), it can be observed that for constant values of other seismic parameters, as the $M_w$ increases, the lateral spreading also increases. This is because a larger $M_w$ implies higher energy carried by seismic waves during the corresponding seismic event, making liquefaction-prone areas more susceptible to liquefaction, resulting in larger lateral spreading.

In Fig. 8(b), it is evident that with constant values for other seismic parameters, an increase in the $k_y$ results in a decrease in the lateral spreading. This observation aligns with the displacement calculation theory of Newmark [1], where, under unchanged conditions, a higher $k_y$ leads to a smaller Newmark-calculated displacement. Moreover, a higher $k_y$ signifies a greater residual shear strength of liquefied soil. Soil displacement begins only when the combined shear forces resulting from seismic and gravitational effects surpass the residual shear strength. As a result, the lateral spreading decreases.

Fig. 8(c) illustrates that different PGA have varying effects on liquefaction-induced lateral spreading. However, there is an overall trend of increasing lateral spreading with higher PGA, like the pattern observed for moment magnitude. This alignment with theoretical expectations and the calculation principles of Newmark [1] reinforces the consistency of the model.

The alignment of the proposed predictive model with theoretical analysis and its good agreement with existing empirical models highlight its reliability. This suggests that the proposed liquefaction-induced lateral spreading prediction model can effectively capture the intricate nonlinear relationship between input parameters and the output parameter, making it a dependable prediction tool.

5.1. Example of lateral spreading of liquefaction

A comprehensive dataset of liquefaction-induced lateral spreading well-cases has been meticulously gathered and organized for validation based on published literature. This dataset encompasses detailed information about soil profiles, soil properties, and Standard Penetration Test (SPT) values of potentially liquefiable soils. The collected data has been compiled into an extensive lateral spreading database, as presented in Table 4. This table incorporates various parameters, including (N₁)₆₀, which signifies the normalized and standardized SPT blow count under standard atmospheric pressure, and (N₁)_60-cs, representing the equivalent clean-sand normalized SPT blow count, accounting for the purity of the sand. Additionally, $\mathrm{l}\mathrm{n}\left(D_{n,true}\right)$ denotes the natural logarithm of the observed lateral spreading in the field, measured in millimeters (mm).

5.2. Application of lateral spreading prediction model

The determination of the yield acceleration for lateral spreading involves various factors, including soil layer distribution, soil strength, soil unit weight, and groundwater level. To achieve this, a limit equilibrium analysis model is set up using Slide software, utilizing the Morgenstern-Price method for limit equilibrium analysis. A critical step in this process is establishing the residual shear strength of the liquefiable soil, as emphasized by previous research [16].

To calculate the residual shear strength of the liquefiable soil, three different formulas are utilized, i.e., Olson and Johnson [20], Idriss and Boulanger [37], and Kramer and Wang [38]. These formulas are applied based on the standard penetration test (SPT) values of the liquefiable soil. By adjusting the dynamic load factor to achieve a safety factor of exactly 1.0, corresponding to the horizontal seismic coefficient, the yield acceleration of the site is obtained. The yield accelerations determined using different residual shear strength prediction models are summarized in Table 5. This process ensures that the analysis captures the specific characteristics of the liquefaction potential of the soil in question.

The application of the model proposed in this research is compared with the Jibson 2007 model [4] and the Rathje and Saygili 2009 model [5]. Using the residual shear strength predictions from the Olson and Johnson [20], Idriss and Boulanger [37], and Kramer and Wang [38] models, the lateral spreading for the 22 liquefaction cases described in Table 5 are predicted. These predictions are compared with observed lateral spreading data to assess the effectiveness of the proposed method. The predictive results are illustrated in Fig. 9.

To provide an objective assessment of the accuracy of the lateral spreading prediction method presented in this research, the evaluation is conducted using $R^2$ and RMSE. The results are detailed in Table 6. This comprehensive analysis aims to validate the reliability and accuracy of the approach outlined in this research when compared to established models and actual observations.

Upon examining the evaluation metrics presented in Table 6, it is evident that the proposed method in this research consistently outperforms the Jibson 2007 model and the Rathje and Saygili 2009 model across most cases. Notably, when utilizing the residual shear strength predictions by Kramer and Wang [38], the performance of the proposed method is marginally lower than that of the Jibson 2007 model. However, in the remaining scenarios, the proposed method surpasses both the Jibson 2007 model and the Rathje and Saygili 2009 model. This comparative analysis underscores the effectiveness and accuracy of the approach developed in this research for predicting lateral spreading.

Fig. 9Prediction of liquefaction-induced lateral spreading based on different residual shear strength models of liquefied soil

a) Olson and Johnson model [20]

b) Idriss and Boulanger model [37]

c) Kramer and Wang model [38]