Semi-empirical relationships to assess the seismic performance of slopes from an updated version of the Italian seismic database

Seismic performance of slopes can be assessed through displacement-based procedures where earthquake-induced displacements are usually computed following Newmark-type calculations. These can be adopted to perform a parametric integration of earthquake records to evaluate permanent displacements for different slope characteristics and seismic input properties. Several semi-empirical relationships can be obtained for different purposes: obtaining site-specific displacement hazard curves following a fully-probabilistic approach, to assess the seismic risk associated with the slope; providing semi-empirical models within a deterministic framework, where the seismic-induced permanent displacement is compared with threshold values related to different levels of seismic performance; calibrating the seismic coefficient to be used in pseudo-static calculations, where a safety factor against limit conditions is computed. In this paper, semi-empirical relationships are obtained as a result of a parametric integration of an updated version of the Italian strong-motion database, that, in turn, is described and compared to older versions of the database and to well-known ground motion prediction equations. Permanent displacement is expressed as a function of either ground motion parameters, for a given yield seismic coefficient of the slope, or of both ground motion parameters and the seismic coefficient. The first are meant to be used as a tool to develop site-specific displacement hazard curves, while the last can be used to evaluate earthquake-induced slope displacements, as well as to calibrate the seismic coefficient to be used in a pseudo-static analysis. Influence of the vertical component of seismic motion on these semi-empirical relationships is also assessed.


Introduction
Stability of a slope subjected to seismic loading depends on slope and earthquake characteristics. Attainment of limit conditions can be mainly ascribed to seismic-induced inertial forces acting in the unstable mass and/or to a decrease of the shear strength induced by pore water pressure build-up (Di Filippo et al. 2019) and cyclic degradation of shear strength parameters (Bandini et al. 2015). Inertial forces induce permanent displacements accumulating until the end of the seismic event, while a decrease in the shear strength can result in slope collapse occurring in post-seismic conditions.
Focusing on the effects of seismic-induced inertial forces, seismic performance of slopes can be assessed in terms of the permanent displacement attained at the end of the seismic event, through dynamic analyses (Rathje and Cho 2019), displacement-based analyses (Ji et al. 2020) and either limit-equilibrium (Saade et al. 2016) or limit-analysis (Ausilio et al. 2000) force-based computations. These methods primarily differ in complexity and in the representation of the seismic action. While dynamic analyses can be quite complex and time-consuming, in limit-equilibrium analyses seismic-induced inertial forces are reproduced in a quite crude fashion through the pseudo-static approach. Hence, the well-known displacement-based Newmark's method (1965) can be deemed a good compromise both in terms of accuracy and complexity, particularly for shallow sliding surfaces, for which the assumption of rigid behaviour is typically adequate. In contrast, in the presence of deep sliding surfaces deformabililty of the unstable soil mass should be considered (Rathje and Antonakos 2011;Tsai and Chien 2016).
After parametrically integrating a set of acceleration time histories, several empirical relationships based on the Newmark's method are available in the literature. Using these models, the seismic performance of a slope can be estimated through the permanent displacements induced by earthquake loading, computed as a function of some ground motion parameters and the slope yield seismic coefficient k y Fotopoulou and Pitilakis 2015): this latter is the seismic coefficient for which the pseudo-static factor of safety (F S ) of a given slope is equal to unity. The models mainly differ in the adopted seismic database and in the selected functional form. These semi-empirical relationships are typically used to develop site-specific displacement hazard curves when a fully-probabilistic approach is followed (Rathje and Saygili 2011;Rathje et al. 2014;Du and Wang 2016;Macedo et al. 2018;Wang and Rathje 2018;Macedo and Candia 2020). In this approach, the aleatory variability of both ground motion properties and resulting displacements is explicitly taken into account. However, due to the inherent complexity of this procedure, these variabilities are typically neglected or not treated rigorously, thus following either a deterministic or a semi-probabilistic approach , as already done in several studies (e.g. Jibson 2007) when a screen analysis (Stewart et al. 2003) is meant to be performed. Much more frequently adopted in common practice is the pseudo-static approach, where the inertial forces acting in the slope are represented via constant static-equivalent forces proportional to the self-weight of the unstable mass through the seismic coefficient k. In the framework of the performance-based design, the seismic coefficient can be calibrated against the seismic performance desired for the slope, typically expressed in terms of threshold permanent displacements, d y . This procedure has been widely adopted in the literature for applications to slopes (Bray and Travasarou 2009;Rampello et al. 2010;Biondi et al. 2011;Bray et al. 2017), earth dams (Seed 1979;Hynes-Griffin and Franklin 1984), solid-waste landfills , hillside residential and commercial developments (Stewart et al. 2003), mountain reservoirs (Veylon et al. 2017) and retaining walls (Biondi et al. 2014;Gaudio et al. 2018aGaudio et al. , b, 2021. In this approach, the seismic performance of the slope is still evaluated applying the force-based pseudo-static approach, in which the equivalent seismic coefficient k is the critical seismic coefficient that corresponds to a threshold or limit displacement d y . In this paper, well-known semi-empirical relationships are updated by performing a parametric integration of a new version of the Italian seismic database, this covering the time frame 1972-2017 and thus including the seismic records of recent destructive seismic events (2009 L'Aquila, Maugeri et al. 2011;2012Emilia, Mucciarelli and Liberatore 2014, de Nardis et al. 2014 Central Italy seismic sequences, Mollaioli et al. 2018;Luzi et al. 2019).
The article is organised as follows. The main characteristics and the most significant ground motion parameters of the new Italian seismic database are firstly presented and compared to those obtained in a previous version of the database (SISMA 2008;Scasserra et al. 2008) and to the ones computed with widely-used ground motion prediction equations (GMPEs). Then, semi-empirical models for predicting permanent displacements are developed. One (scalar) and two (vector) ground motion parameter empirical models are obtained, adapting the ones recently proposed by Rathje and Cho (2019) to the Italian seismicity. Furthermore, upper-bound semi-empirical relationships are derived, in the functional form proposed by Ambraseys and Menu (1988), thus linking permanent displacements d to the ratio between the yield and the maximum seismic coefficient, k y /k max , where k max = PGA/g and PGA is the peak ground acceleration. Influence of the vertical component of seismic motion on these relationships is also assessed. The seismic coefficient k, used in pseudo-static calculations, is therefore computed for given threshold values of permanent displacements d y assuming different subsoil classes and acceleration levels. Different semi-empirical relationships are finally provided where permanent displacements d are expressed as a function of other ground motion parameters in addition to the ratio k y /k max , such as PGA, the mean period T m , the significant duration D 5-95 (Tropeano et al. 2017) and the Arias intensity I A (Jibson 1993(Jibson , 2007Chousianitis et al. 2014), identifying the most convenient models to evaluate earthquake induced displacements.
The results of this study can be mainly used for hazard mapping or to perform screen analyses (Stewart et al. 2003) for slopes located on the national territory based on an updated version of the Italian seismic database.  Luzi et al. 2016b). This database, named ITACA + ESM_2017 in the following, includes 947 records of 207 seismic events with moment magnitude M w ≥ 4, peak ground accelerations PGA ≥ 0.05 g, epicentral distance R ep < 100 km, and focal depth z ip ≤ 45 km, recorded by 297 stations located throughout the national territory.

Seismic database
Values of the local magnitude M L are available for 203 seismic events, while the moment magnitude M w is known for 141 events. Epicentral distance R ep and focal depth z ip are known for all records, while Joyner and Boore distance, R JB , (Joyner and Boore 1981) is known for 309 records. For the records with magnitude greater than 6, Joyner and Boore distance R JB was computed through the empirical relationship proposed by Malagnini and Montaldo (2004), as a function of the epicentral distance R ep : Since this relationship is valid for superficial magnitude M s ≥ 6, that is not available, it was assumed M s = M w according to Idriss (1985), thus using Eq. (1) for M w ≥ 6.
Information about the subsoil class identified by the classification provided by Eurocode 8, Part 1 (CEN 2003) is available for each station. The subsoil class is attributed on the basis of the average shear wave velocity in the upper 30 m, V s,30 , or, in the absence of such data, on surface geological categorisation . The records were divided into five groups: rock or rock-like subsoil, with V s,30 > 800 m/s (class A), dense granular and stiff cohesive subsoil, with V s,30 = 360-800 m/s (class B), medium dense granular and medium stiff cohesive subsoil, with V s,30 = 180-360 m/s (class C), loose granular and soft cohesive subsoil, with V s,30 < 180 m/s (class D) and subsoils with V s values of class C or D and thickness of 5 to 20 m underlain by stiffer material with V s > 800 m/s (class E). When only geological data were available, subsoil categories were identified by symbols A * , B * , C * , D * , E * . About 13% of the records (123) were attributed to subsoil class A, 49.5% (469) to subsoil class B, 31% (294) to subsoil class C, 1.5% (14) to subsoil D, and 5% (47) to subsoil E. Due to the scarce number of records assigned to subsoil classes D and E, these were considered as a single subsoil group together with subsoil class C, as already done in previous works (e.g. Rampello et al. 2010), and identified as soft soils in the following.
The database provides three acceleration time histories for each record: two horizontal components, related to the North-South and East-West directions, and one vertical component, relative to the Up-Down direction.
Characteristics of the database are described in Table 1 and in Figs. 1, 2 and 3: these were obtained by processing all the records as well as using the ESM-flatfile_2017 . Ground motion parameters describing both horizontal and vertical acceleration time histories are given in the Appendix 1. For each subsoil class, Table 1 shows the minimum and maximum magnitudes M, the focal depth z ip and (1) R JB = −3.5525 + 0.8845 ⋅ R ep   170 (18%), 126 (13.3%) and 65 (6.9%), respectively. Focal mechanism is not available for 5 (0.5%) registrations. Figure 1 shows the distribution of the number of events (a) and records (b), as a function of the local magnitude M L , and of the number of events as a function of the focal mechanism (c) for the ITACA + ESM_2017 and for the SISMA 2008 (Scasserra et al. 2008) databases, the latter including 247 records of 89 Italian seismic events recorded by 101 stations in the period 1972-2002. For the seismic databases mentioned above, the maximum number of events is obtained for local magnitude M L in the range of 3.5 to 4.5 (102 for ITACA + ESM_2017 and 45 for SISMA), while the maximum number of records is obtained for M L = 4.5-5.5 (444 for ITACA + ESM_2017 and 96 for SISMA).
As expected, the number of seismic events and records of ITACA + ESM_2017 database is higher than the one of SISMA database, except for lower local magnitudes M L = 2.5-3.5, since ITACA + ESM_2017 database includes events with M w > 4 only.
The relative frequency distribution of moment magnitude M w (a), peak ground acceleration PGA (b), mean period T m (c) and significant duration D 5-95 (d) of the ITACA + ESM_2017 database is shown in Fig. 2 for the three subsoil classes: rock-like (A), stiff (B) and soft soils (C, D, E). For all the three subsoil classes, modal moment magnitude is between 4.5 and 5.5 (Fig. 2a), typical of most aftershock records of the main Italian seismic sequences (Tropeano et al. 2017). For the horizontal acceleration time histories, the most frequent peak acceleration ( Fig. 2b) falls between 0.05 and 0.1 g, while the modal value of mean period T m (Fig. 2c) is between 0.2 and 0.4 s. Significant duration D 5-95 (Fig. 2d) shows the modal value between 3 and 6 s for the rock-like and the stiff subsoil classes, and between 6 and 12 s for soft subsoils.
In Fig. 3 the moment magnitude M w is plotted against the epicentral (Fig. 3a) and Joyner and Boore (b) distances, while peak ground velocity PGV (c), Arias intensity I A (d), mean period T m (e) and spectral acceleration Sa(1.5T s ), computed at the degraded period 1.5T s (f), are plotted against peak ground acceleration PGA. A fundamental period of the slope T s = 0.19 s was assumed to compare Sa(1.5T s ) to the values recently presented by Rathje and Cho (2019). As expected, increasing values of peak ground velocity PGV, Arias intensity I A and spectral acceleration Sa(1.5T s ) are obtained for increasing PGA, while no remarkable influence on the mean period T m is detected.
The most significant ground motion parameters were compared with attenuation curves obtained from widely used ground motion prediction equations (GMPEs). For the sake of brevity, only the figures related to the horizontal components of seismic motion recorded on stiff soils (subsoil class B) during strong earthquakes (magnitudes M = 5.5-6.9) are shown in the following, unless otherwise specified.
The curves obtained from the GMPEs proposed by Sabetta and Pugliese (1996) (Fig. 4a) and Ambraseys et al. (2005) (Fig. 4b) are plotted in Fig. 4 together with the peak horizontal acceleration derived from the seismic records. In each graph, the median curves and median curves ± 1/2 σ of the GMPEs are plotted, where σ is the standard deviation of the error. The curves proposed by Ambraseys et al. (2005) are given for three different focal mechanisms (F N for normal faulting earthquakes, F T for thrust faulting earthquakes and F O for odd faulting earthquakes). Points related to events whose Joyner and Boore distance is calculated through the empirical relationship proposed by Malagnini and Montaldo (2004) are represented in grey. It turned out that points belonging to the ITACA + ESM_2017 database can be reliably predicted by the attenuation laws. Specifically, PGA values lay between the curves proposed by Sabetta and Pugliese (1996) for magnitudes M = 5-7, except for 15 out of 144 (about 10%) seismic events exceeding the curve M = 7 for epicentral distances R ep ≥ 8 km and just 1 plotting below the μ−σ/2 curve with M = 5. Similarly, most of the points (105 out of 126, equal to the 83%) representing the ITACA + ESM_2017 database plot between the curves proposed by Ambraseys et al. (2005) for magnitudes M = 6-7, in the whole range of considered Joyner and Boore distances (R JB = 0.7-80 km) and particularly for high values of R JB (≥ 20 km). The reduced, but not negligible amount of data points plotting outside the range located by the attenuation laws, 10% and 17% with respect to Sabetta and Pugliese (1996) and Ambraseys et al. (2005), respectively, can be attributed to the different adopted databases: Sabetta and Pugliese (1996) referred to an old Italian database, Ambraseys et al. (2005) used a European database with M w ≥ 5, while ITACA + ESM_2017 includes the most recent Italian seismic events. Figure 5 shows the comparison between the attenuation laws proposed by Sabetta and Pugliese (1996) (Fig. 5a), Kayen and Mitchell (1997) (Fig. 5b) and the values of Arias intensity I A derived from the records. Note that a different definition of I A is adopted in the equation proposed by Sabetta and Pugliese (1996), as follows: For the GMPE proposed by Sabetta and Pugliese (1996) the maximum Arias intensity of the horizontal components was considered, while for the GMPE provided by Kayen and Mitchell (1997) the sum of the two values of the horizontal components was used. Most of the points representing the ITACA + ESM_2017 database lie between the area enclosed by the attenuation curves in this case too: 119 out of 144 (about 83%) values of Arias intensity plot between the curves when comparing with the GMPE by Sabetta and Pugliese (1996), especially for epicentral distances R ep > 4 km (Fig. 5a). A bit worse agreement is obtained with the relationship by Kayen and Mitchell (1997) (Fig. 5b), as 23 points out of 107 (≈ 21%) plot below the range detected by the attenuation laws.
The computed values of mean period T m and significant duration D 5-95 are compared in Fig. 6 with the GMPEs proposed by Rathje et al. (2004) and Kempton and Stewart (2006), for magnitudes M = 6-6.5. These relationships were developed for the rock-like subsoil class only. The comparison is also made with the GMPEs provided by Tropeano et Comparison between the peak ground horizontal accelerations recorded on stiff soils (ITACA + ESM_2017 database) and those computed using the GMPEs proposed by: a Sabetta and Pugliese (1996) and b Ambraseys et al. (2005) al. (2017) that reworked these relationships using the Italian database SISMA 2008. Figure 6a shows that 8 out of 31 points (about 26%) representing the mean period T m of the ITACA + ESM_2017 database fall just out of the range detected by the median ± σ curves provided by Tropeano et al. (2017). On the contrary, the values of significant duration D 5-95 obtained from the records agree with all the median curves proposed by the mentioned Authors, also capturing the increase of the significant duration with the Joyner and Boor distance.
To consider records significant for triggering sliding phenomena, a screening criterion of seismic records was introduced to predict permanent displacements with the (a) (b)

Fig. 5
Comparison between the Arias intensities recorded on stiff soils (ITACA + ESM_2017 database) and those computed using the GMPEs proposed by: a Sabetta and Pugliese (1996) and b Kayen and Mitchell (1997) (a) (b) Fig. 6 Comparison between a the mean periods and b the significant durations recorded on rock-like soils, (ITACA + ESM_2017 database) and those computed using the GMPEs proposed by: a Rathje et al. (2004) and Tropeano et al. (2017); b Kempton and Stewart (2006) and Tropeano et al. (2017) 1 3 rigid sliding-block model (Newmark 1965). According to Keefer and Wilson (1989), seismic events to be excluded from the database are those recorded at a source-to-site distance above the upper-bound curve shown in Fig. 7. The ITACA + ESM_2017 database includes 5 magnitude-distance pairs plotting slightly above this upper-bound curve, corresponding to 7 horizontal acceleration time histories (2 for subsoil class A, 1 for class B and 4 for group C, D, E). However, these time histories were not excluded in the end, as this would not involve significant changes in the results of the sliding-block analyses shown in the following.

Evaluation of seismic-induced permanent displacements
Permanent displacements induced by seismic actions were evaluated through a parametric integration of the acceleration time histories included in the database described in Sect. 2, i.e. by computing the seismic-induced permanent displacements through a double-integration of the equation of relative motion, as discussed below, for different values of some parameters characterising the seismic input (e.g. PGA) and the seismic resistance of the slope (e.g. the yield seismic coefficient k y ). To this end, the rigid-block model sliding on a horizontal plane (Newmark 1965) was adopted ( Fig. 8) both accounting for and ignoring the vertical component of ground motion. If the time history of the vertical seismic coefficient, k v (t), is considered, the slope resistance to earthquake loading, represented by the yield seismic coefficient, k y , is not constant during the analysis as it becomes a function of time ( Fig. 9). In this study, a simple though useful relationship to evaluate the yield seismic coefficient k y (t) in the presence of k v (t) is adopted following the formulation proposed by Sarma and Scorer (2009): where k y(k v =0) is the seismic yield coefficient computed in the absence of the vertical component. Equation (3) indicates that the yield coefficient k y , evaluated considering the vertical component of ground motion, can be either higher or lower than k y(k v =0) , depending on

Fig. 7
Screening of significant seismic events for slopes the sign assumed by k v (t): k v > 0 is directed opposite to gravity, then reducing the available seismic resistance.
Permanent displacements were evaluated via double integration of the equation of relative motion of the rigid block, that can be written as (Michalowski and You 2000) where k h (t) is the time history of the horizontal seismic coefficient, g = 9.81 m/s 2 is the gravitational acceleration, T f is the duration of the seismic input and C is a shape factor that accounts for the shape of the sliding surface. Referring to the simple scheme of infinite slope, the permanent displacement accumulated in the direction of the planar sliding surface can be evaluated using the following expression for the shape factor C (Crespellani et al. 1998): where φ′ is the angle of shearing resistance and β is the angle of slope inclination with respect to the horizontal. It is worth mentioning that Eq. (5) holds also to soils with effective cohesion c′ ≠ 0. Typical and stable configurations are related to values of φ′ and β in the range 22°-28° and 5°-20°, respectively, thus providing values of the shape factor C ranging between 1.03 and 1.12, close to unity. Hence, in calculation of permanent displacement a shape factor C = 1 was assumed without loss in generality. Integration was carried out in the time intervals where the horizontal seismic coefficient, k h (t), is greater than the yield seismic coefficient, k y (t) (Fig. 9), and in all time instants when the relative velocity of the sliding block is positive. A null relative velocity was obtained at the end of the integration procedure by adding a set of n zeros to the original time history of horizontal seismic coefficient, k h (t), n being the initial number of points of k h (t), this not affecting any ground motion parameter characterising the seismic input, such as its amplitude, frequency content or strong motion duration.
The parametric integration was carried out considering both signs of the acceleration time histories, assuming the maximum computed displacement as the slope displacement related to that seismic input. Influence of the vertical component of ground motion was observed to be negligible on the semi-empirical models developed in the following.

One and two-ground motion parameters models
Permanent displacements of slopes obtained through Newmark-type computations are first expressed as a function of ground motion parameters for a given value of the yield seismic coefficient. Following Rathje and Cho (2019), one or two ground motion parameters (GM) are used to derive the following expressions: where d is the displacement in units of cm, A 0 , A 1 and A 2 are regression coefficients, GM 1 and GM 2 are the ground motion parameters, σ is the standard deviation of the model and t is the reciprocal value of the normal standard distribution related to a generic level of confidence (t = 0 for median curves, here adopted). Displacements computed in these semiempirical relationships are obtained ignoring the vertical component of seismic motion; moreover, only displacements greater than 1 cm are taken into account to derive the regression coefficients, thus limiting the number of records used to 120, that is similar to the number of records considered by Rathje and Cho (2019) (105). Although the permanent displacements and the computed equations are sensitive to the assumed yield seismic coefficient k y , for the sake of comparison, the same value k y = 0.12 as that adopted by Rathje and Cho (2019) is here assumed for the yield seismic coefficient in the integration of the Italian seismic database.
Computed displacements are plotted in Fig. 10 as a function of five ground motion parameters (PGA, PGV, I A , T m and Sa(1.5T s )), together with the relevant best-fit curves and the values of the standard deviation σ. Table 2 lists the regression parameters, together with the standard deviation of the model σ and the coefficient of determination R 2 . Increasing permanent displacements d are obtained for increasing values of all the considered seismic parameters. Among the one-parameter models, the smallest standard deviation (σ = 0.535) is obtained using I A , that describes intensity, frequency content and duration of the ground motion, thus resulting the most efficient parameter. The standard deviation of PGV model (σ = 0.581) is slightly larger than the one computed using the Arias intensity. However, the PGV model may result in a simpler use in that PGV can be evaluated using one of Single Ground Motion Parameter semi-empirical relationships, where permanent displacement is expressed against: a peak ground acceleration, b peak ground velocity, c Arias intensity, d mean period, and e spectral acceleration at 1.5T s (120 records) the several Ground Motion Prediction Equations (GMPEs) available in the literature. Moreover, GMPEs developed for Arias intensity I A may be affected by a higher epistemic uncertainty than that of GMPEs built for PGA or PGV, this resulting in a higher standard deviation of the model (Campbell and Bozorgnia 2012;Douglas 2012). The standard deviation related to PGA, Sa(1.5T s ) and T m models is larger than 0.7: this is because PGA and Sa(1.5T s ) do not account for the frequency content of ground motion, while the mean period T m does not consider its intensity.
Considering the two ground motion parameter models, the lowest values of standard deviation are computed when considering the I A -T m and the I A -PGV models, with values of σ = 0.403 and 0.347, respectively. The standard deviation related to other pairs of parameters (PGA-PGV; PGA-T m ; and PGA-I A ) are about equal to 0.5, similar to the value obtained for the single parameter I A model (σ = 0.535). The smaller standard deviation computed using the couples I A -T m and I A -PGV results from the information on the frequency content provided by either T m or PGV. Similar values of σ computed for the two parameters (PGA, I A ) and one parameter (I A ) models suggest that information about intensity provided by PGA is already represented by I A so that data scattering is not reduced by adding PGA.
Based on the computed standard deviation, the best one-parameter model is the I A model with a value of σ similar to those computed for most of the two-parameter models. The best two-parameter models are the I A -T m and I A -PGV models, characterised by values of σ about 25% and 35% smaller than that of the I A model. However, the main issue in using these three models consists in the lack of hazard maps for expected, site-specific values of I A and T m . Therefore, the PGA-PGV model may be preferred to the previous ones due to the fact that several ground motion models are available for PGA and PGV (e.g. Rathje and Cho 2019).
Regression parameters and standard deviations computed by Rathje and Cho (2019) are also listed in Table 2: the intercepts A 0 are always larger than those obtained in this study, due to the larger computed permanent displacements, except for what obtained in the PGA, T m model. This result can be attributed to the different adopted seismic database that is characterised by much higher earthquake magnitudes (M w > 6.5). Values of parameter A 1 obtained in this study are 26 to 61% higher than those computed by Rathje and Cho (2019)  for one-ground motion parameter models. Similarly, higher values of A 2 are evaluated in the two-GM parameter relationships (+ 67% at most), except for the PGA, I A equation (− 70%). A fair agreement is obtained in terms of the standard deviation σ of the models, this being either higher or lower than those by Rathje and Cho (2019) of about the same amount (maximum difference of about ± 35%). The semi-empirical relationships showed in this paragraph relate permanent displacements to some ground motion parameters, for a given yield coefficient of the slope k y ; different semi-empirical models can be obtained specifying different values of k y , thus referring to slopes characterised by different seismic resistances. These curves can then be used as a tool to develop site-specific displacement hazard curves within a fully-probabilistic approach. By contrast, semi-empirical curves proposed in the following paragraphs (Sects. 4.2 and 4.3) are meant to be adopted either to calibrate the seismic coefficient k or to estimate the maximum expected permanent displacement d, in the framework of a deterministic approach.

Evaluation of the seismic coefficient
Empirical relationships to evaluate earthquake-induced displacements can be extended introducing the yield seismic coefficient k y as an added independent variable. Permanent displacements can then be computed for slope-specific values of k y and compared with threshold values related to the level of seismic performance desired for the slope at hand. However, the pseudo-static approach is still the most common procedure adopted to evaluate the seismic stability of a slope, thanks to its simplicity and to the experience accumulated by engineers in using this approach in the framework of the limit-equilibrium force-based computations. In a pseudo-static analysis, inertial actions are represented by static-equivalent forces proportional to the self-weight of the unstable mass through the seismic coefficient k, whose selection is crucial. The pseudo-static approach provides an evaluation of a safety factor against sliding that is assimilated to a failure mechanism. Nevertheless, as mentioned above, the inertial effects mainly produce a progressive development of permanent displacements in the slope rather than a failure mechanism, due to the transient and cyclic nature of seismic actions. Using the pseudo-static approach for the evaluation of slope stability under seismic conditions then requires relating the maximum expected displacements to the seismic coefficient k and the safety factor F S . Such an equivalence can be obtained using relationships of the kind mentioned above between earthquake-induced displacements and given ground motion parameters. Specifically, upperbound relationships relating earthquake-induced displacements to the ratio k y /k max can be used, following Rampello et al. (2010). For each subsoil class and acceleration level, an exponential relationship between the permanent displacement d and the ratio k y /k max can be written in the form: where A is a non-dimensional coefficient defining the slope of the curve in a semi-logarithmic scale and B is the displacement for k y = 0. Assuming a log-normal distribution for displacements, the median curve corresponding to the 50th percentile is obtained (coefficients A and B), as well as the upper-bound curve corresponding to the 94th percentile, (coefficients A and B 1 > B), which is a value in the range of those commonly adopted in the literature (e.g. Whitman and Liao 1985;Rampello et al. 2010;Biondi et al. 2014). The slope A computed for the median curve was also used for the upper-bound curve (A 1 = A), thus assuming a constant standard deviation irrespective of the ratio k y /k max . This hypothesis was validated by performing preliminary analyses in which negligible differences were obtained considering the dependence of σ on the ratio k y /k max .
In deriving these upper-bound relationships, incompleteness of data used for the parametric integration should be considered, in that the available acceleration time histories were recorded during the last 40-50 years, a very short period if compared to the return periods typically considered when checking the safety of slopes subjected to seismic events. Therefore, the horizontal acceleration time histories included in the database ITACA + ESM_2017 were scaled to reach desired values of PGA equal to 0.05, 0.15, 0.25 and 0.35 g, following the Italian Building Code (Ministero delle Infrastrutture 2018), limiting the scale factor F h in the range 0.5-2. Such a narrow interval for F h allowed to preserve site effects; moreover, the frequency content and the strong-motion duration of the records were not affected by the scaling operation. The same scaling factor (F h = F v ) was adopted for the time histories of vertical acceleration. In the presence of both horizontal components (NS and WE), the vertical accelerograms were scaled twice, according to the scale factors F h used in each horizontal direction. Therefore, four different sets of acceleration time histories were obtained for each subsoil group (rock-like, stiff and soft soils).
Calculations of permanent displacements were carried out using the four sets of scaled time histories varying the yield seismic coefficient k y(k v =0) from 10 to 80% of the maximum seismic coefficient k max (k y /k max = 0.1-0.8).

Influence of the vertical component of ground motion
Seismic-induced permanent displacements were computed both in the presence and in the absence of the vertical component of ground motion, fitting the results with the upperbound curves given in Eq. (8). The computed upper-bound curves for soft soils (subsoil classes C, D, E) are shown in Fig. 11: for low to intermediate levels of peak acceleration, the curves obtained considering the vertical component provide slightly greater displacements, while for PGA = 0.35 g, upper-bound curves overlap. The observed differences are less than 10% so that they can be ignored, leading to the conclusion that the vertical component of ground motion does not affect substantially the upper-bound curves and thus the Fig. 11 Comparison between the upper-bound curves (94th percentile) calculated in the presence and absence of the vertical component of ground motion for soft soils choice of the seismic coefficient adopted in a pseudo-static analysis. Similar results were obtained for the other subsoil classes considered in this study (A and B). The same conclusion was drawn in some studies adopting different seismic inputs (Gazetas et al. 2009;Sarma and Scorer 2009). Then, the time histories of horizontal acceleration only are used in the computation presented in the following.

Upper-bound empirical curves
Seismic-induced permanent displacements are plotted in Fig. 12 as a function of the ratio k y /k max between the yield and the maximum seismic coefficients for the subsoil class B (stiff soils) only. The median (50th percentile, dashed lines) and the upper-bound curves (94th percentile, solid lines) are plotted in the figure, for each acceleration level. As expected, seismic-induced displacements decrease as the ratio k y /k max increases and the peak acceleration PGA decreases, this occurring for all subsoil groups. The upper-bound curves computed for the three subsoil groups are compared in Fig. 13. The maximum permanent displacements are invariably computed for soft subsoils: this result is not surprising and is to be attributed, for a fixed value of the peak ground acceleration PGA, to the typical frequency content and significant duration characterising the acceleration time histories recorded on this type of soils. Indeed, propagation of seismic waves through them causes both the mean period T m and the significant duration D 5-95 increase (see Fig. 2), this leading to a strong increase of the permanent displacement accumulated till the end of the seismic event. In contrast, the lowest displacements are calculated for the stiff or rock-like subsoils, depending on the peak ground acceleration PGA and k y /k max ratio. Specifically, for low acceleration levels (PGA = 0.15 and 0.05 g) the lowest displacements are obtained for the subsoil class A, while for high accelerations (PGA = 0.35 g) they are attained for stiff soils (subsoil class B). This result might appear counterintuitive, as local (i.e. stratigraphic) effects are supposed to produce upper-bound relationships for subsoil class B plotting above the one obtained for subsoil class A. However, the largely different sample size for the two subsoil classes A and B, equal to 390 and 2490 respectively, Fig. 13 Comparison amongst upper-bound permanent displacements obtained for the different subsoil groups caused a higher coefficient of variation γ ln = σ ln /μ ln for the subsoil class A, and therefore a greater shift from the median to the upper-bound curve. Conversely, for PGA = 0.25 g, negligible differences are observed between subsoil classes A and B. It is worth mentioning that the adequateness of the sample size n has been checked for all subsoil classes and acceleration levels by verifying that the residuals of ln(d), ε ln , follow a normal distribution in all considered cases.
Values of coefficient A are in the range of 7.24-7.76 with a mean value equal to 7.45, while B 1 ranges between 0.11 and 1.54 m, with an average of 0.63 m and a coefficient of variation of about 76% (Table 3). These values are compared with those computed by Rampello et al. (2010) using the database SISMA (Scasserra et al. 2008): while the slope A is quite similar, with a maximum difference of 8%, coefficient B 1 computed with database ITACA + ESM_2017 always attains lower values, up to − 72%, except for the soft subsoils at high peak ground accelerations PGA (= 0.35 and 0.25 g), for which an almost double value is obtained. This outcome resulted in lower permanent displacements except for the soft subsoils at high PGA. The observed differences come from the different databases adopted in calculations, as many changes occurred over the years in the attribution of subsoil classes to the recording stations ; moreover, the database used by Rampello et al. (2010) is updated to 2002, while that adopted in this study is updated to 2017. Table 3 also lists the values of the standard deviation of the model σ and of the coefficient of determination R 2 : a quite narrow range of values is detected, with σ = 0.91-1.15 and R 2 = 0.66-0.77.
Computation of coefficients A and B 1 was repeated adopting horizontal acceleration time histories recorded by stations not affected by any local effects, thus obtaining a smaller set of accelerograms. To this end, the descriptive data of the recording site available in the database ITACA + ESM_2017, together with its topographical and morphological conditions, were selected. Specifically, free-field housing and plain, centre of Valley and not classified morphological conditions were chosen, thus reducing the number of considered stations from 297 to 115. Table 4 shows the coefficients A and B 1 obtained using the reduced database. Slight variations in the coefficient A were computed, with a maximum difference of about 7% obtained for rock-like subsoils and PGA = 0.35 g, while changes in the coefficient B 1 were at most equal to about 35% for soft subsoils and PGA = 0.05 g. Similar results were obtained for standard deviation of the model σ and the coefficient of determination R 2 , the former ranging between 0.92 and 1.28, the latter between 0.60 and 0.78, similarly to what obtained with the complete database. Upper-bound relationships from the complete database are then considered in the following, due to the small differences.

Seismic coefficient
Assuming a ductile behaviour along the sliding surface, thus assuming a constant shear strength, the seismic coefficient k can be estimated as a fraction of the maximum seismic coefficient k max by applying the position k = η·k max with η ≤ 1. Specifically, k can be related to the desired upper-bound threshold displacement d y inverting Eq. (8) as follows (Fig. 14): where coefficient B 1 corresponds to the upper-bound curves (94th percentile) presented in the previous paragraph for different acceleration levels, to introduce a fixed level of conservativism in the estimation of k (Stewart et al. 2003). Hence, if a pseudo-static analysis is performed and limit equilibrium is reached (F S = 1), it is k = k y , this corresponding to a maximum permanent displacement d = d y . This procedure constitutes the link between the force-based pseudo-static approach and the permanent displacements computed via the semi-empirical relationship as the one given in Eq. (9). Computed values of coefficient η are given in Table 5, for selected values of threshold displacements d y , equal to 15, 5 and 2 cm. Threshold values of 15 and 5 cm correspond to moderate to negligible levels of damage assuming a ductile slope behaviour during sliding, in the absence and presence of manufactures (Idriss, 1985), while d y = 2 cm was considered acceptable by Wilson and Keefer (1985) for stiff soils and rocks, characterised by a brittle behaviour. A minimum safe value of η = 0.10 was assumed when calculations returned η < 0.10. This occurred for high admissible permanent displacements (d y = 15 cm) and low acceleration levels (PGA = 0.05 g).
It is worth noting that, since the deformation pattern of a slope may deeply differ from the simple sliding along a well-defined surface, as assumed in the analyses, threshold value of d y should be solely considered as an index of seismic performance of the slope.
Different coefficients η can be computed assuming different threshold displacements. The discussed procedure can also be adopted assuming other levels of conservativism, using different percentile to define the upper-bound curves reported in Eq. (8).

Further semi-empirical relationships
In the procedure discussed in Sect. 4.2, seismic-induced permanent displacements were related to the ratio between the yield and the maximum seismic coefficients, k y /k max , thus describing the seismic input by its amplitude only. The remarkable scattering of the data points shown in Fig. 12 reveals that permanent displacements do not only depend on the peak acceleration acting on the potential sliding mass, but on other parameters characterising the seismic motion as well. This led to the development of a number of semi-empirical relationships, where the ground motion is described by additional parameters such as the mean period T m , that reflects the frequency content of the seismic input, the significant duration D 5-95 , that quantifies the duration of the strong-motion phase, and the Arias intensity I A , that represents the energy content of the seismic event.

Seismic input described via k max , T m and D 5-95
Empirical relationships are used herein where the permanent displacement d (variable 1) or the normalised displacement d (variable 2) (Yegian et al. 1991): is expressed as a function of peak acceleration PGA, mean period T m and significant duration D 5-95 , to account for the influence of frequency content and duration of ground motion (Madiai 2009;Biondi et al. 2011;Tropeano et al. 2017). Permanent displacement d strongly depends on site conditions, while the normalised displacement d has the advantage to be almost site-independent thanks to the inclusion of 3 ground motion parameters (PGA, T m and D 5-95 ) explicitly accounting for the stratigraphic effects. Following Biondi et al. (2011), three semi-empirical relationships can be obtained (a, b and c) and used for both variables d and d , thus providing six different equations: 1a, 1b and 1c for permanent displacement d and 2a, 2b and 2c for the normalised displacement d , namely: where σ is the standard deviation of the considered functional form and t is the abovementioned reciprocal value of the normal standard distribution (t = 1.555 for 94th percentile, here adopted). Again, as in Sect. 4.2, a constant standard deviation σ was assumed in computations. Relationship (a) (Eq. 11) is equivalent to the one presented in Eq. (8), while relationship (b) (Eq. 12) complies with the following conditions: d → ∞ for k y /k max = 0 and d = 0 for k y /k max = 1 (Ambraseys and Menu 1988). Relationship (c) (Eq. 13) was formulated following Saygili and Rathje (2008).
Coefficients A, B, C and D characterising the different regressions are given in the Appendix 2, together with the relevant values of the standard deviation σ and the coefficient of determination R 2 . Type 2 relationships provide a better best-estimate of analysis results with higher values of R 2 , thanks to the dimensionless expression adopted for The maximum increase of R 2 was computed for soft soils, with a maximum increase of about 22% for PGA = 0.05 g.
The upper-bound curves computed using Eqs. (11)-(13) are plotted in Fig. 15. The results related to regressions (a), (b) and (c) are superimposed in the figures, for the permanent displacement d (1) (Fig. 15a, c, e) and the non-dimensional displacement d (2) (Fig. 15b, d, f), for each subsoil class and two acceleration levels (PGA = 0.35 and 0.15 g). Regressions (a) to (c) provide very similar upper-bound displacements, for all subsoil groups and peak accelerations PGA, with the maximum difference of about 30% obtained for stiff soils at high accelerations (PGA = 0.35 g) and for ratio k y /k max = 0.80, the latter involving very low permanent displacements d or d . Influence of PGA is mostly mitigated using the non-dimensional displacement d , reducing sensibly the maximum difference between the permanent displacements computed for PGA = 0.35 and 0.15 g (Fig. 15) from about 380 to 50%, the latter value obtained for subsoil class B.
No significant improvement in the evaluation of permanent displacements is observed from a to c semi-empirical curves, so that relationship 2a is advised due to its simple form. However, relationships of type 2 need an estimate of ground motion parameters T m and D 5-95 : to this end, the GMPEs proposed in Sect. 2 as well as in the literature can be used, this requiring evaluation of magnitude M and of epicentral or Joyner and Boore distance, R ep or R JB , respectively (Tropeano et al. 2017).
Influence of subsoil groups on permanent displacements is shown in Figs. 16 and 17, for relationships of types 1 (d) and 2 ( d ). Similarly to what observed by Biondi et al. (2011) and Tropeano et al. (2017), permanent displacements d depend on subsoil group (Fig. 16), as also discussed in Sect. 4.2.2. The largest differences among subsoil groups in the evaluation of permanent displacement d are obtained comparing stiff (class B) and soft soils (classes C, D and E), the latter providing higher permanent displacements (about + 70%) for PGA = 0.35 g and k y /k max = 0.45. Conversely, Fig. 17 shows that the dependence of upper-bound non-dimensional displacements on subsoil groups almost disappears irrespective of the peak ground accelerations, except for PGA = 0.35 g, where the seismic inputs recorded on rock-like subsoil result in a lower non-dimensional displacement d . In this case, a maximum difference of about 90% is computed for the non-dimensional displacement d of rock-like (A) and stiff (B) subsoils, at k y /k max = 0.8, for which very low values of d are attained. Nonetheless, the curves obtained for the other values of PGA overlap, thus providing the same estimate of the seismic performance.

Seismic input described via I A and k max
The energy content of the seismic input can be introduced in the semi-empirical relationships through the Arias intensity I A as well, as proposed by Jibson (1993Jibson ( , 2007. The relationships adopted in this study are developed for the permanent displacement d only, thus providing relationships 1d, 1e and 1f: (e) (f) Fig. 15 Comparison between upper-bound (94th percentile) curves obtained using relationships 1 (a, c and e) and 2 (b, d and f) for given values of peak ground accelerations and for all subsoil groups 1 3 where I A is expressed in m/s and d in cm. As in Jibson (1993Jibson ( , 2007, the curves were obtained without distinguishing between subsoil groups. Values of I A = 0.2-10.0 m/s and k y = 0.02-0.40 were considered by Jibson (1993), while values of I A = 0.002-5.451 m/s and k y = 0.005-0.28 were considered in this study. Computed coefficients A i , B i and C i , for i = 1d, 1e and 1f, are listed in Tables 6, 7 and 8 for relationships given in Eqs. (14), (15) and (16), respectively. These values are also compared to those computed by Jibson (1993Jibson ( , 2007. The comparison between the coefficients calculated using relationship 1d (Table 6) shows a similar dependence of permanent displacements d on Arias intensity I A (coefficient A 1d ) but a doubled influence of the yield coefficient k y (coefficient B 1d ); similar results were obtained by Madiai (2009) using an older version of the Italian seismic database (Scasserra et al. 2008). A coefficient R 2 = 0.674 was computed for this relationship. Conversely, a better agreement was obtained using the other two relationships (Tables 7 and 8), with higher values of the coefficient of determination R 2 > 0.81.
Influence of the subsoil group on relationships 1 for all considered peak ground acceleration levels Fig. 17 Influence of the subsoil group on relationships 2 for all considered peak ground acceleration levels

3
The upper-bound (94th percentile) relationships 1d (Eq. 14, Table 6) are plotted in Fig. 18a, while Fig. 18b shows the comparison between relationships of type 1f (Eq. 16, Table 8). With reference to Fig. 18a, a fair agreement was obtained for low values of the yield coefficient k y (= 0.02, 0.05 and 0.10), while lower values of upper-bound displacements were computed for higher values of k y (= 0.20, 0.25 and 0.28). By contrast, a better agreement is provided by relationships 1f (Fig. 18b) for high ratios k y /k max (= 0.60, 0.70 and 0.80), while a poor agreement is observed for low ratios k y /k max . Again, the observed differences can be attributed to the different databases considered in the analyses. Specifically, Jibson (1993) used 11 strong-motion records from California (10) and Iran (1), while in Jibson (2007) 2270 strong-motion horizontal components, recorded during 30 earthquakes occurred worldwide, were adopted.

Summary and conclusions
Seismic performance of slopes can be successfully assessed within the framework of the performance-based design, where the seismic-induced permanent displacement d is typically assumed as an index of seismic performance. Permanent displacements induced by earthquake loading can be computed through the well-known Newmark's method, where the slope is assimilated to a rigid-block sliding on a horizontal plane whenever the seismic coefficient k(t) exceeds the yield (or critical) seismic coefficient k y of the slope at hand. In this paper, Newmark-type computations have been performed through a parametrical integration of an updated version of the Italian seismic database. This includes seismic events recorded in the Italian territory in the time frame 1972-2017, thus including the destructive 2009 L'Aquila, 2012 Emilia and 2016 Central Italy earthquakes. It is seen that ground motion parameters characterising the new Italian seismic database, such as the peak ground acceleration PGA, the Arias intensity I A , the mean period T m and the significant duration D 5-95, can be reasonably computed using, on a first approximation, wellknown ground motion prediction equations (GMPEs) already available in the literature. It is understood that the new database can be used to calibrate new GMPEs as well.
Both one-and two-ground motion parameter, semi-empirical models have been developed following Rathje and Cho (2019), where the permanent displacements are expressed as a function of some parameters of the seismic motion, for a given yield seismic coefficient k y . As expected, two-ground motion parameter relationships are characterised by a higher coefficient of determination R 2 when accounting for the frequency content and the duration of the seismic motion, in addition to its intensity: among these, the best model, both in terms of accuracy and ease of use, proved to be the (PGA, PGV) model. These semi-empirical models can be used to evaluate seismic-induced displacements and, in turn, the seismic performance of a slope characterised by a given k y , through the comparison with threshold displacements d y .
The semi-empirical relationships can be modified to explicitly introduce the dependency of slope displacements on the ratio between the yield and the maximum seismic coefficient, k y /k max . These relationships can be used to either estimate the permanent displacement d induced by earthquake loading or to evaluate the seismic coefficient k to be adopted when using a pseudo-static approach. To this end, upper-bound (94th percentile) relationships linking the permanent displacement d to the ratio k y /k max have been proposed and values of the seismic coefficient have been computed assuming threshold displacements d y = 2, 5 and 15 cm, three subsoil groups (rock, stiff and soft soils) and four acceleration levels (PGA = 0.05, 0.15, 0.25 and 0.35 g). It is worth mentioning that the procedure laid out in this study to compute the equivalent seismic coefficient can be readily adopted for different values of threshold displacements d y and for different percentiles. Influence of the vertical component of ground motion on these upper-bound curves has been found to be negligible.
Semi-empirical relationships including additional ground motion parameters provide a better agreement with the permanent displacements computed via the sliding-block analyses: specifically, the best agreement is obtained using the non-dimensional variable d = d∕ PGA ⋅ T m ⋅ D 5−95 and relationship 2a (R 2 max = 0.86). However, accounting for the frequency content and duration of the strong motion phase requires an estimate of mean period T m and significant duration D 5-95 using the above-mentioned GMPEs, thus introducing further dispersion in the relationship. The seismic input can also be described by the Arias intensity I A of ground motion, thus providing an additional simplified tool while considering the energy released during the earthquake. The discussed model forms are already well-known in the literature; however, due to the fact that the target of the analyses is evaluating the order of magnitude of permanent displacements instead of their exact value, it is believed that, in this context, increasing the complexity of the analyses by proposing new model forms would not have improved their prediction substantially. Indeed, the updated semi-empirical relationships are meant to be used in the framework of a screen analysis, where sites with expected poor seismic performance are distinguished from those with acceptable seismic performance.
Clearly, some limitations arise in the study. First, the proposed semi-empirical relationships are developed for the Italian seismicity only: nonetheless, the same discussed procedure can be adopted in other high-seismicity areas. Moreover, computation of permanent displacements is based on the simplifying assumptions underlying the adopted rigid sliding-block model, such as the one related to a constant shear strength during seismic loading, this possibly leading to underestimate slope displacements. Conversely, neglecting soil deformability results in assuming an infinite wavelength of seismic motion, this leading to overestimate permanent displacements d at the end of the earthquake. The latter assumption explains why the discussed results can be deemed more appropriate for slopes characterised by shallow translational mechanisms.

Appendix 2 : Coefficients for the semi-empirical relationships of Sect. 4.3.1
Coefficients A, B, C and D characterising the different regressions are given in Tables  11, 12, 13, 14, 15, 16 with a double subscript that identifies the functional form (a, b, c) and the dependent variable (1, 2). The coefficients calculated for each empirical relationship are listed together with the standard deviation σ and the coefficient of determination R 2 , for each subsoil group and acceleration level.