**1. Introduction**

Loss of life and property remains an unavoidable consequence of major earthquakes. Studies of the consequences of major earthquakes have attempted to analyze the damage and make recommendations for reducing loss in the event of future earthquakes throughout history [1–3]. Liquefaction-induced lateral displacement is one of the most prevalent and damaging of these effects. It can cause enormous blocks of soil to move by a few millimeters to 10 m or more, inflicting substantial damage to lifeline networks, buried utilities, and a variety of other subsurface and civil engineering projects. Liquefaction-induced lateral displacement is most common on gentle slopes built on loose sand with a groundwater table close to the surface of the ground; however, open faces such as stream channels can also be susceptible [4].

**Citation:** Ahmad, M.; Amjad, M.; Al-Mansob, R.A.; Kami ´nski, P.; Olczak, P.; Khan, B.J.; Alguno, A.C. Prediction of Liquefaction-Induced Lateral Displacements Using Gaussian Process Regression. *Appl. Sci.* **2022**, *12*, 1977. https://doi.org/ 10.3390/app12041977

Academic Editors: Andrea Chiozzi, Elena Benvenuti and Željana Nikoli´c

Received: 31 December 2021 Accepted: 8 February 2022 Published: 14 February 2022

**Publisher's Note:** MDPI stays neutral with regard to jurisdictional claims in published maps and institutional affiliations.

**Copyright:** © 2022 by the authors. Licensee MDPI, Basel, Switzerland. This article is an open access article distributed under the terms and conditions of the Creative Commons Attribution (CC BY) license (https:// creativecommons.org/licenses/by/ 4.0/).

Various approaches have been presented to estimate the magnitude of lateral displacement to date, and from the technical perspective, they can be classified as: (1) numerical analysis based on finite element or finite difference approaches (e.g., Finn et al. [5], Liao et al. [6] and Arulanandan et al. [7], (2) simplified analytical methods, e.g., Newmark [8], Towhata et al. [9], and Kokusho and Fujita [10], (3) empirical methods based on either laboratory testing set or analytical methods of lateral spreading case history records (e.g., Hamada et al. [11] and Youd et al. [12]) and (4) machine learning approaches (e.g., Wang and Rahman [13]). These different approaches are reviewed herein, with particular emphasis on empirical models and soft computing techniques.

### *1.1. Finite Element Analysis*

To simulate different aspects of liquefaction and lateral spreading, including seismic loads, rapid loss of shear strength, redistribution of pore water pressure, and soil softening, Liao et al. [6] reported that very complex finite element and finite difference approaches are required. Very intricate numerical techniques, large computer skills, and extensive resources are necessary to create a realistic three-dimensional simulation inside the real-time domains. A number of well-known finite element method (FEM) and finite difference (FD) software programmes are used for liquefaction-induced lateral displacement assessments and earthquake soil dynamic analysis. Other finite element models provided by Hamada et al. [14] and Orense and Towhata [15] to determine the lateral ground deformations generated by earthquakes. Gu et al. [16,17] estimated liquefaction deformation using a planar strain model. It successfully anticipated the pattern of displacements at a wildlife site in California [17], but overestimated the magnitude of displacements by around 30%.

### *1.2. Simplified Analytical Models*

### 1.2.1. Sliding Block Model

Newmark [8] proposed a model based on strategy of sliding block on a frictional sloping surface that predicted seismically induced ground deformations by integrating accelerations above the sliding block's yield acceleration to obtain its velocities. The angle of inclination and the factor of safety over sliding are associated to yield acceleration. When the driving force (seismic acceleration) equal to or greater than resisting force (yield acceleration), block will begin to slide. The total cumulative resulting deformation is then determined by integrating the sliding block velocity. Yegian et al. [18] used Newmark's concept to introduce their model for predicting the permanent ground displacement expressed as

$$D = N\_{eq} T^2 a\_p f\left(\frac{a\_y}{a\_p}\right) \tag{1}$$

where *D* is the lateral ground deformation, *Neq* denote cycles number equivalent to uniform base motion, *T* denotes time interval (s), *a<sup>y</sup>* denotes yield acceleration (*g*), *a<sup>p</sup>* denotes peak acceleration (*g*), and *f* denotes dimensionless function that depends on base motion. Baziar et al. [19] also used Newmark's concept, assuming an equivalent sinusoidal base acceleration record, to propose their model for predicting the permanent ground displacement expressed as:

$$
\log D = 1.46 \log I\_a - 6.642a\_y + 1.546\tag{2}
$$

*D* denotes lateral ground deformation (cm), *I<sup>a</sup>* presents arias intensity (m/s), and *a<sup>y</sup>* presents yield acceleration (*g*).

### 1.2.2. Minimum Potential Energy Model

This model was proposed by Towhata et al. [9] depending on the results of shaking table testing. The final position of soil layers was found by the principle of minimal potential energy, using the Lagrangian equations of motion, and assuming the variation of lateral ground deformation with depth as a sine function and with neglecting inertial effects during dynamic loading. Tokida et al. [20] used the same principle to establish equations for predicting the maximum lateral displacement at the center of a slide as:

$$D = 1.73 \times 10^{-5} L^{1.99} H^{0.298} T^{-0.275} \theta^{0.963} \text{ for } (10 \text{ } m \le L \le 100 \text{ } m) \tag{3}$$

$$D = 1.29 \times 10^{-5} L^{1.99} H^{0.28} T^{-0.243} \theta^{0.995} \text{ for } (100 \, m < L \le 1000 \, m) \tag{4}$$

where *D* is horizontal displacement (m), *L* = length of slide (m), *H* represents average thickness of liquefied layer (m), *T* represents average thickness of liquefied surface layer (m), and *θ* is the slope of ground surface express in percentage.

### 1.2.3. Shear Strength Loss and Strain Re-Hardening Model

Bardet et al. [21] reported that in 1997, Byrne proposed a method to find the final position of a liquefying slope using the finite element software tool Fast Lagrangian Analysis of Continua. In liquefaction region, it is assumed that the liquefied material is initially free of shear, and is subjected to isotropic pressure. After such immediate melting of liquefied soil, the shear stress (*τ*) was supposed to rise with shear strain unless reached a certain residual shear strength (*τST*). Although the liquefied soil regains shear strength, the shear modulus was supposed to take a constant value *GLIQ*. The final position of the slope is determined using the dynamic equation of motion.

### 1.2.4. Viscous Models

Hadush et al. [22] reported that Aydan [23] considered the liquefied subsoil to act as a visco-elastic object and used an upgraded Lagrangian numerical approach to find the deformation velocities for the liquefied soil sub-layers. They also proposed a numerical method based on cubic interpolated pseudoparticles for liquefaction-induced lateral displacement analysis in the context of fluid dynamics. Liao et al. [6] reported that Hamada et al. [24] recommended to use viscous models to estimate the liquefaction-induced lateral displacement. Kokusho and Fujita [10] studied the role of water film in lateral flow failure during earthquakes, on the basis of field survey results collected from Niigata (1964) earthquake. It was reported that the water films produced under the fine soil sub layers did not actually have shear resistance, and a significant factor for the large lateral flow displacement.

### *1.3. Empirical Models and Soft Computing Techniques*

Hamada et al. [11] provided a preliminary relation of measuring horizontal ground displacement in meter's relying on 60 case histories, majority of them are obtained in Niigata and Noshiro, Japan. It can be seen from Table 1 that the equations are very common and easy to apply that contains only two parameters of site geometry and not considered seismic and geotechnical parameters but it has been suggested for limited dataset making it insufficiently broad to be extended to additional lateral displacement sites.

Youd and Perkins [25] suggested "liquefaction severity index" (*LSI*) to estimate maximum horizontal ground displacement generated by an earthquake. The *LSI* (inches) was calculated using distance to seismic energy source (*R*) (km); and moment magnitude (*Mw*) with maximum range of horizontal ground displacement as 2.5 m. This model assumes that the value of *LSI* depends on only seismic parameters (*R*, *Mw*). At the time, the proposed equation drew the attention of engineers. Although this method may have been useful for assessing lateral spreads inside the western United States, but lacks applicability and hence didn't receive widespread use.

Bardet et al. [21] used multiple linear regression (MLR) for developing relation to estimate lateral ground deformation for free face and sloping ground situations, respectively, utilizing data gathered by Bartlett and Youd [26,27], including three kinds of input variables:


3. Geotechnical parameters (in percent)—average mean particle size within *T*<sup>15</sup> (*D*5015, mm) and averaged fines contents in *T*<sup>15</sup> (*F*15)

The MLR approach was used to create the Youd et al. [12] model, presented in 1992, was built on upgraded results of Bartlett and Youd [26,27] for estimating lateral ground displacement, *D<sup>H</sup>* (m). As indicated in Table 1 the models include free face and sloping ground conditions equations. This model gained attraction amongst geotechnical engineers due to its utilization of a huge dataset from various earthquakes, as well as geometry of site, geotechnical data, and seismic characteristics. Although, it does have some limits in terms of applications. For example, the free face equation was used when 5 ≤ *W* ≤ 20%, Jafarian and Nasri [28] gathered the latest dataset of liquefaction-induced lateral ground deformation based on uncertainties of different boreholes, which outperformed Hamada et al. [11] Kanibir [29], Al Bawwab [30], Javadi et al. [31], Youd et al. [12], and Baziar and Azizkandi [32] models.

**Table 1.** Empirical and machine learning approaches for liquefaction-induced lateral displacement.


Note: *N*160*<sup>s</sup>* : (*N*1)<sup>60</sup> value corresponds to Js, Js is the lowest factor of safety below water table using simplified approach; *θ*: larger slope of either ground surface or the base of liquefied soil (%); *H*: thickness of liquefied zone (m); *R*\*: modified source distance factor that is a function of earthquake magnitude.

Soft computing is made up of a variety of techniques that function together, such as: artificial neural network, genetic algorithm, neuro-computing etc. Wang and Rahman [13] reported that new area of machine learning has arisen for handling decisions, modeling, and control issues. Baziar and Ghorbani [33] and Wang and Rahman [13] both used artificial neural networks (ANN) to estimate horizontal ground displacement. Javadi et al. [31] computed lateral displacement for free face and sloping ground using genetic programming (GP) using upgraded case data from Youd et al. [12]. In comparison to the MLR approach, the proposed GP approach has some advantages. Table 1 shows the proposed equations for free

face and sloping ground. Adaptive neuro-fuzzy inference system (ANFIS) based approach was suggested by Jadanian [34] using 426 case histories data and shows an improvement to the Youd et al. [12], Kanibir [29], Bardet et al. [21], and Rezania et al. [4] models.

Soft computing methodologies are more accurate than analytical formulas, according to all of these studies. The findings revealed that the ML models mentioned above are capable of obtaining the experimental observations with acceptable accuracy. However, this field continues to be further explored.

The Gaussian process regression (GPR) approach has been successfully applied in many domains, but its application in geotechnical engineering is limited based on literature surveys. Considering the improved performance of GPR, it is, however, used for the first time in this study to predict the liquefaction-induced lateral spread displacement for free face condition. To demonstrate the efficacy of the proposed GPR-based model, the results are compared with various well-known models for calculating the *DH*.

### **2. Gaussian Process Regression**

Gaussian process regression (GPR) is one of the appropriate and newly-proposed methods that have been employed for various machine learning examples. GPR is a stochastic, non-parametric technique for addressing complicated and non-linear challenges. GPR assumes that the target variable m is determined as follows:

$$m = f(n(k)) + \varepsilon \tag{5}$$

where *f* represents unidentified functional dependency, *n* represents the number input parameters, and ε represents Gaussian noise with variance *σ<sup>a</sup>* 2 . It's a method of indicating precedence straight over function space. The mean and covariance of a Gaussian distribution are matrices and vectors, respectively. The GPR model can determine the prediction distributions, which is similar to ensuring input knowledge [35]. The GPR approach is based on the idea that surrounding data informs neighbours.

GPR makes use of a number of kernel functions. A restriction of GPR regression is the selection of a suitable kernel function. Pearson VII kernel function (*PUK*) is utilized for GPR proposed model in this work.

$$PUK = \left( 1 / \left[ 1 + \left( 2\sqrt{\left\| \mathbf{x}\_i - \mathbf{x}\_j \right\|}^2 \sqrt{2^{\left(1/\omega\right)} - 1} / \sigma \right)^2} \right]^{\omega} \right) \tag{6}$$

where *ω* and *σ* are the Person's width, and peak tailing factor, respectively.

### **3. Case-History Database**

The case-study dataset used for this work was collected using three sources (Chu et al. [36], Youd et al. [12], and Cetin et al. [37]) which contains a total of 247 records of lateral displacement related to free face ground conditions.

The input parameters chosen by Youd et al. [12] have been largely acknowledged amongst researchers as a full and acceptable set for controlling lateral displacement. As a result, several other scholars have chosen the same characteristics as important indicators (e.g., Javadi et al. [31]; Jafarian and Nasri [28]; Baziar and Saeedi Azizkandi [32]). In addition, with inclusion of a ground's intensity measure, peak ground acceleration (PGA, *amax*) is employed in the present study to increase data set, making it more competent and effective in accounting for earthquake causes. By considering the causative fault types of all earthquakes, Sadigh et al. [38] employed attenuation equation to predict the PGA.

In this research, the following seven key parameters have been used to evaluate lateral displacement: earthquake magnitude (*M*), peak ground acceleration (*amax*, *g*), horizontal distance to seismic energy source (*R*, km), average particle size in *T*<sup>15</sup> (*D*5015, mm), average fines material (particles < 0.075 mm) in *T*<sup>15</sup> (*F*15, %), accumulative thickness of saturated layers with adjusted SPT number (*N*1)<sup>60</sup> < 15 (*T*15, m), free-face ratio (*W*, %), while the output is liquefaction-induced lateral displacement (*DH*, m).

In this work, training datasets are based on 80% of the data available (198 sets of data in free face characteristics). The testing dataset has been used to evaluate the proposed models' prediction abilities. The 49 historical records data are used as testing datasets in this study. The training and testing datasets were partitioned depending on statistical features of the datasets, such as mean and standard deviation. The model efficiency is enhanced by the statistical consistency of the training and testing datasets, which makes it easier to evaluate them. Table 2 shows the evaluation metrics of input and output variables in training and testing datasets for free face. Summary of liquefaction-induced lateral ground deformation database is presented in Appendix A, Table A1.


**Table 2.** Statistical parameters for free-face condition.

### **4. Correlation Analysis**

Correlation coefficients (*ρ*) have been used to test the significance of the relation between different factors (see Table 3). The equation for *ρ* is as:

$$\rho(u,v) = \frac{cov(u,v)}{\sigma\_{\mathcal{U}}\sigma\_{v}}\tag{7}$$

where *cov* indicates covariance, *σ<sup>u</sup>* represents standard deviation of *u*, and *σ<sup>v</sup>* defines standard deviation of *v*. |*ρ*| > 0.8 signifies a strong relation among *u* and *v*, values from 0.3–0.8 represents a moderate relationship, while |*ρ*| < 0.30 signifies a weak relation [39]. According to Song et al. [40], a relation is considered as "strong" if |*ρ*| > 0.8. *M*, *amax*, *R*, *D*5015, *F*15, *T*<sup>15</sup> and *W* have moderate to weak relations, as seen in Table 3. As a result, no variables from the lateral displacement estimation model were eliminated. Table 3 reveals that the correlation coefficient has a maximum absolute value of 0.761 and there is no "strong" link between different pairs of components.

**Table 3.** Correlation between parameters.


### **5. Construction and Evaluation of Prediction Model**

Figure 1 illustrates the prediction model's creation process. In this case, 80% and 20% of the dataset were chosen as training and test sets, respectively, based on statistical integrity. Second, the predictive model was constructed using the trial-and-error approach based on training set utilizing the optimum hyperparameters configurations. Iterative method was utilized to find optimum values for the hyperparameters after setting them to random values (within a reasonable range). The values of the key kernel parameters, omega (*ω*) and sigma (*σ*) are 0.4 while noise is 0. 3in the GPR model after multiple trials. Finally, the testing data was used to evaluate the proposed GPR model's performance using four common evaluation metrics: coefficient of determination (*R* 2 ), coefficient of correlation (*r*), mean absolute error (*MAE*), root mean square error (*RMSE*), ratio of the root mean square error (*RSR*) to the standard deviation of measured values, and Nash–Sutcliffe coefficient (*NSE*). The *R* <sup>2</sup> and *NSE* values that are higher, and *RSR* values that are lower, imply that proposed model's prediction accuracy is better. Waikato Environment for Knowledge Analysis software was used throughout the whole calculation process The Pearson VII function-based kernel [41] was employed in this study for the GPR model.

**Figure 1.** The flowchart for GPR based model to predict liquefaction induced lateral displacement.

The generated model's performance was assessed using *R* 2 , *r*, *MAE*, *RMSE*, *RSR*, and *NSE*.

$$\mathcal{R}^2 = 1 - \frac{\sum\_{i=1}^n (\mathbf{x}\_i - y\_i)^2}{\sum\_{i=1}^n (\mathbf{x}\_i - \overline{\mathbf{x}})^2} \tag{8}$$

$$r = \frac{\sum\_{i=1}^{n} (\mathbf{x}\_i - \overline{\mathbf{x}})(y\_i - \overline{y})}{\sqrt{\sum\_{i=1}^{n} (\mathbf{x}\_i - \overline{\mathbf{x}})^2} \sqrt{\sum\_{i=1}^{n} (y\_i - \overline{y})^2}} \tag{9}$$

$$MAE = \frac{1}{n} \sum\_{i=1}^{n} (x\_i - y\_i) \tag{10}$$

$$RMSE = \sqrt{\frac{1}{n} \sum\_{i=1}^{n} (x\_i - y\_i)^2} \tag{11}$$

$$RSR = \frac{\sqrt{\sum\_{i=1}^{n} (\mathbf{x}\_i - \mathbf{y}\_i)^2}}{\sqrt{\sum\_{i=1}^{n} (\mathbf{x}\_i - \overline{\mathbf{x}})^2}} \tag{12}$$

$$NSE = 1 - \frac{\sum\_{i=1}^{n} (x\_i - y\_i)^2}{\sum\_{i=1}^{n} (x\_i - \overline{x})^2} \tag{13}$$

where *n* denotes the set of data points, *x<sup>i</sup>* and *y<sup>i</sup>* denotes the actual and estimated output of data's *i*th sample, respectively; *x* and *y* represents the mean actual and estimated output of the dataset, respectively. The *r* value varies from −1 to 1. A perfect distribution between actual and estimated values is represented by value of *r* equal to 1, whereas a value of 0 shows no relation [42]. For *MAE* = 0, the model's value is perfectly aligned with the real value, and the model is deemed "ideal." The *MAE* value is between 0 and +∞. The mean squared difference between outputs and targets is termed as *RMSE*, and its value ranges from 0 to +∞. The *NSE* scale ranges from -∞ to 1, with 1 representing the ideal match. A strong relation is indicated by an *NSE* score of more than 0.65 [43,44]. The RSR ranges from a perfect 0 to a significant positive number. A smaller *RSR* indicates low RMSE, indicates that the model is more predictive. The *RSR* and *NSE* categorization ranges are shown in Table 4 as very good, good, adequate, and inadequate [44].

**Table 4.** Statistical indicators for model performance evaluation.


### **6. Result and Discussion**

### *6.1. Performance of GPR Model*

The GPR model's efficiency were assessed using coefficient of correlation (*r*), mean absolute error (*MAE*), root mean square error (*RMSE*), ratio of root mean square error (*RSR*) and Nash-Sutcliffe coefficient (*NSE*). The trend line for GPR in training and testing phases has been drawn by comparing the observed regression in Figure 2 scatter plot, and the GPR findings have the maximum inclination to the line *y* = *x* (i.e., *R* <sup>2</sup> = 0.9402 in training and *R* <sup>2</sup> = 0.894 in testing phases). Table 5 shows clearly that for the training model, *r* = 0.9697, *MAE* = 0.3403, *RMSE* = 0.5597, *RSR* = 0.248 and *NSE* = 0.938. Whereas for the testing model *r* = 0.9455, *MAE* = 0.5443, *RMSE* = 0.8438, *RSR* = 0.387 and *NSE* = 0.851. The trend line for GPR in training and testing phases has been drawn by comparing the observed regression in Figure 2 scatter plot, and the GPR findings have the maximum inclination to the line *y* = *x* (i.e., *R* <sup>2</sup> = 0.9402 in training and *R* <sup>2</sup> = 0.894 in testing phases).

**Figure 2.** Scatter plot presenting the measured *D<sup>H</sup>* values versus the predicted *D<sup>H</sup>* (**a**) training and (**b**) testing model.


**Table 5.** Performance statistics of GPR model in comparison with extra available models in literature.

Note: - represents that this performance statistic is not included in the reference.

The performance of the developed GPR model was compared to the evolutionary polynomial regression (EPR), MLR and ANN models in literature based on the *R* 2 , *MAE*, and *RMSE* criteria and the results are summarized in Table 5. Furthermore, in terms of the *MAE* and *RMSE* statistical measures in training, the lowest value was found for GPR (*MAE* = 0.3403 m, *RMSE* = 0.5597 m) compared to EPR (*MAE* = 0.537 m, *RMSE* = 1.003 m), ANN (*MAE* = 0.702 m, *RMSE* = 1.074 m), and MLR (*MAE* = 0.81 m, *RMSE* = 1.24 m). Whereas the prediction results in the testing, the MAE and RMSE values was found less for GPR (*MAE* = 0.544 m, *RMSE* = 0.8438 m) compared to EPR, ANN, and MLR except the MAEs values of EPR (=0.291 m) and MLR (=0.43 m) models. The superiority may be owing to the fact that the GPR model excellently captures the nonlinear relationships between lateral displacement and its influencing factors. It can therefore be concluded that, based on statistical indices, the GPR model had the best results. However, due to the use of different number of datasets, a comparison between these results is unwarranted. A project that uses different datasets is needed to gives generalized model to geotechnical earthquake engineering.

### *6.2. Sensitivity Analysis*

The sensitivity results of the GPR model were examined using Yang and Zang's approach for determining the impact of input variables on *DH*. This strategy, which has been used in a number of research [44–47], is as follows:

$$r\_{\rm ij} = \frac{\sum\_{k=1}^{n} (\boldsymbol{\pi}\_{ik} \times \boldsymbol{\pi}\_{ok})}{\sqrt{\sum\_{k=1}^{n} \boldsymbol{\pi}\_{ik} ^ 2 \sum\_{k=1}^{n} \boldsymbol{\pi}\_{ok} ^ 2}} \tag{14}$$

*xik* and *xok* are the actual and estimated variables, respectively, and *n* represents number of datasets (i.e., 198 data samples). For each input variable, the *rij* value varies from zero to one, with the greatest *rij* values indicating the most efficient output factor (i.e., *DH*). Figure 3 displays the *rij* values for all input factors. The accumulative thickness of saturated layers with adjusted SPT number, *T*<sup>15</sup> (*rij* = 0.843) has the greatest effect on the *DH*. Furthermore, Table 3 shows that the accumulative thickness of saturated layers with adjusted SPT number, *T*<sup>15</sup> has the highest *ρ* of 0.518 in all other parameters validating the sensitivity analysis results.

**Figure 3.** Sensitivity analysis of input variables.

### **7. Conclusions**

The GPR model was used to estimate liquefaction-induced lateral displacement in this work. The predictive model was built using seven input parameters and one output parameter. Performance measures such as R<sup>2</sup> , *r*, *RMSE*, *MAE*, *RSR*, *NSE*, and visual inspection such as scatter plots were used to assess the effectiveness of the developed model. This study's findings can be summarized as follows:


The GPR approach has the advantage of becoming easily modified as new data becomes available, reducing need for expertise and time to modify an existing design aid or equation and/or suggest a new equation.

**Author Contributions:** Conceptualization, M.A. (Mahmood Ahmad); methodology, M.A. (Mahmood Ahmad), M.A. (Maaz Amjad); software, M.A. (Maaz Amjad) and M.A. (Mahmood Ahmad); validation, R.A.A.-M., P.K., P.O., M.A. (Maaz Amjad) and B.J.K.; formal analysis, M.A. (Mahmood Ahmad) and M.A. (Maaz Amjad); investigation, P.K., P.O., A.C.A. and B.J.K.; resources, P.K.; data curation, M.A. (Maaz Amjad) and M.A. (Mahmood Ahmad); writing—original draft preparation, M.A. (Maaz Amjad); writing—review and editing, M.A. (Maaz Amjad) and M.A. (Mahmood Ahmad); supervision, M.A. (Mahmood Ahmad), R.A.A.-M. and P.K.; project administration, M.A. (Mahmood Ahmad), A.C.A.; funding acquisition, R.A.A.-M., P.K., P.O. and A.C.A. All authors have read and agreed to the published version of the manuscript.

**Funding:** This research received no external funding.

**Institutional Review Board Statement:** Not applicable.

**Informed Consent Statement:** Not applicable.

**Data Availability Statement:** The data used to support the findings of this study are included within the article.

**Conflicts of Interest:** The authors declare no conflict of interest.

### **Appendix A**


**Table A1.** Summary of Liquefaction-Induced Lateral Ground Deformation Database.

**Table A1.** *Cont.*


**Table A1.** *Cont.*


**Table A1.** *Cont.*



**Table A1.** *Cont.*

### **References**

