IMPA versus Cloud Analysis and IDA: Different Methods to Evaluate Structural Seismic Fragility

: Well-known methods for seismic performance assessment, such as incremental dynamic analysis (IDA), multi-stripes analysis (MSA) and the cloud method, involve nonlinear response time-history analyses to characterize the relationship between the chosen damage measure versus intensity measure. Over the past two decades, many authors have proposed simpliﬁed procedures or nonlinear static approaches to develop fragility. In these procedures, the capacity of the system is evaluated by nonlinear static procedures (i.e., the capacity spectrum method (CSM), the N2 method, modal pushover analysis (MPA)) and the demand is derived by response spectra. In addition to the familiar ones, incremental modal pushover analysis (IMPA) is a novel nonlinear static procedure proposed in recent years, and it is used in this research to present an IM -based fragility estimation. The accuracy and effectiveness of different methods to assess vulnerability are investigated by comparing fragility curves derived by MPA-based cloud analysis, IMPA and cloud analysis against IDA. The comparison gives valuable insights on the inﬂuence of scaling on different sets of records; however, a more extended validation is needed to conﬁrm the obtained results and draw more general conclusions. Results arise from two relatively small bins of record motions differing by ranges of Joyner-Boore distance and scattered in a range of magnitude are presented.


Introduction
Performance-based earthquake engineering (PBEE) procedures allow the prediction and evaluation of the probabilistic seismic performance of bridges and buildings in terms of system-level decision variables, such as loss of use, repair cost and casualties. In the United States, the first generation of PBEE assessment and design procedures for buildings (SEAOC Vision 2000, FEMA 273, ATC-40 [1][2][3][4]) took significant steps toward achieving performance-based earthquake engineering. Since then, the Pacific Earthquake Engineering Research Center (PEER) has been working on developing a more robust methodology that involves four stages: hazard analysis, structural analysis, damage analysis, and loss analysis [5]. In the third stage, damage analysis, fragility functions describe the conditional probability of component, element or system to be damaged for a given intensity measure. The first attempt to determine fragility curves can be dated back to 1975, when the Seismic Design Decision Analysis (SDDA) procedure was proposed in the US [6]. Further developments [7,8] were initially applied in the field of the nuclear industry to define a probabilistic relationship between an intensity measure representing seismic input and a damage measure representing the seismic failure of a component of a nuclear power plant. From then on, several methods to estimate fragility (expert-based, experimental, analytical, hybrid, empirical) have been developed by researchers worldwide, relying on different assumptions and restrictions to overcome prevalent intrinsic uncertainties. However, due to the extremely high subjectiveness, lack of data and other drawbacks typical of expertbased, empirical and experimental methods, the common practice has aroused its interest in analytical and hybrid methods during the last two decades.
Among analytical approaches to derivate fragility curves (probabilistic seismic demand model, elastic spectral analysis, nonlinear static analysis, linear or nonlinear time history analysis [9][10][11][12][13][14][15][16][17]), incremental dynamic analysis (IDA) is a parametric analysis method developed in 1998 and deeply discussed in 2002 [18,19]. IDA became a worldwide method used by engineers and researchers, and it is still widespread. However, the introduction of uncertainties due to an excessively coarse description of seismic input with varying intensity is inevitable. Eventually, many authors have pointed out that a simple amplitude scaling of ground motion records is one of the main shortcomings in IDA, together with its high computational demand [20][21][22][23]. In IDA, a certain number of inputs amplitude scaled to define IM = im, then NL-THA is performed and DCR LS so determined are used to define the distribution of DCR LS |IM = im. This process is repeated by varying the scale factor to define the seismic response in a whole range of seismic intensities. According to previous research [24][25][26], the median response of a structure subjected to scaled ground motion records is comparable to that of a structure subjected to unscaled earthquake ground motion records. However, because a single IM is a highly simplified description of groundmotion severity, the value of DCR LS from different ground-motion with IM = im may be different, and thus the relationship DCR LS~f (IM) is probabilistic. A common limitation in current databases is the lack of strong ground motion records covering high-intensity intervals at specific periods of the structure [27,28]. Thus, an excessive scaling to fit within high-intensity intervals may occur, biasing the structural response [29], generating a false correlation between IMs and EDPs, and increasing uncertainties in the structural response.
In contrast to IDA and MSA (multiple-stripe analysis), the cloud method [26,30,31] involves nonlinear analysis of the structure subjected to a sample of different distance/intensity combined unscaled as-recorded ground motion, which may reduce the number of analyses, uncertainties in seismic input with intensity, computational effort in defining a seismic fragility curve [32], and is based on a regression in the logarithmic space of structural response versus seismic intensity.
To respond to the need for simplified, faster and/or approximate methods, from the late 1990 s, many studies have been published regarding the use of pushover analysis procedures to assess seismic vulnerability, mostly on bridges. In this case, the capacity of the system is evaluated by using nonlinear static procedures (NSPs) (i.e., the capacity spectrum method [33][34][35][36], the N2 method [37][38][39], and modal pushover analysis [20]), while the demand is estimated by response spectra. To assess the reliability of these analytical procedures, different authors have compared developed fragility curves to those obtained by nonlinear time history analysis [34,36,38].
Among NSPs, a novel procedure called incremental modal pushover analysis (IMPA) has been proposed in recent years by Bergami and his co-workers [40]. IMPA requires the execution of modal pushover analysis (MPA) and the evaluation of structural performance within a range of different seismic intensity levels to develop a multimodal capacity curve in terms of base shear versus top displacement. This approach is suitable for performing a displacement-based design procedure and structural analysis of existing structures, yet authors have not suggested the analytical estimation of fragility.
This paper aims to evaluate the reliability of structural fragility derived by the methods mentioned, advancing an IM-based derivation of structural fragility, strikingly similar to IDA, based on IMPA. It is known that IDA has a small sensitivity to record-to-record variability compared to other methodologies. Nevertheless, the results indicate that, amidst its slightly higher sensitivity, IMPA has the advantage of requiring considerably smaller computational effort to perform the structural analysis.
The authors argue that scaling response spectra at a higher range of intensities might introduce less uncertainties than a simple amplitude scaling of ground motions. Further steps of this research will address how the uncertainties in the seismic input affect the reliability of IMPA versus IDA seismic fragility for strong ground motions.
In the following paragraph, analyses are carried out on a real RC frame belonging to a school building located in Norcia (Italy). Each nonlinear dynamic and static procedure is briefly presented, including a step-by-step computational procedure of IMPA. Finally, the reliability of pushover-based estimation of seismic vulnerability is assessed by comparing these curves to those obtained by IDA.

Choice of Engineering Demand Parameter and Intensity Measure
The first-mode spectral acceleration S a (T 1 ,ξ = 5%) is commonly used as an intensity measure (IM) parameter [20,21]. Shome et al. [24] stated that the nonlinear response of an MDOF structure dominated by the first mode of vibration depends on the "intensity" of the records at the first period of vibration, while magnitude and distance play a minor role in it. For these reasons, the 5% damped spectral acceleration at the structure's first-mode period S a (T 1 ,ξ = 5%) or simply S a is adopted as the IM in this work, since the structure selected as the case study is dominated by the first mode of vibration (structure's first-mode period of vibration T 1 = 0.62 s and the modal mass participation at first-mode is 82%, see Section 3.1 for more details).
In the literature, various engineering demand parameters (EDPs) have been proposed [24], somehow representative of the structure's local or global damaged state. In this study, the critical demand to the capacity ratio for the desired limit state (LS), denoted as DCR LS [23,41], is assumed to be the EDP. It represents the demand-to-capacity ratio which brings the system closer to the onset of limit state (herein, the life-safety limit state). The weakest-link formulation is adopted to evaluate the DCR LS (Equation (1)), which means that if the demand-to-capacity ratio D jl /C jl is equal to or higher than unity in just one element, then the structure attains the expected limit state for the lth mechanism.
where N mech and N e are the numbers of the considered potential mechanism of failure and the number of the elements taking part in the lth mechanism, respectively. D jl and C jl (LS) are the demand and the limit state capacity, respectively, evaluated for the jth element of the lth mechanism.
In particular, in this work, only a ductile failure mechanism in columns and beams is considered as a potential failure mechanism (N mech = 1). Therefore, in this deformationbased critical DCR LS , the demand D is expressed in terms of maximum chord rotation in the jth component. Instead, the capacity C in terms of chord rotation is evaluated according to guidelines reported in the Commentary [42] of NTC 2018. Namely, the chord rotation for life-safety limit state is defined as 3 4 of that corresponding to near-collapse limit state θ u , evaluated according to Equation C8.7.2.5 of the Commentary (Equation (2)).
where θ y is the yield chord rotation, φ u and φ y are the ultimate and yield curvature, respectively, L pl is the plastic hinge length, and L v is the shear length.

Record Selection
The PEER Next Generation Attenuation (NGA)-West 2 Project [43] has been used to define the main database of 210 ground motions. In addition, two ground motions representing the 2016 Norcia earthquake, with epicentral distances of 4.6 km and 26.9 km, respectively, have been extracted from the Italian Accelerometric Archive [44] and included in the aforementioned database.
A set of 36 as-recorded ground motions listed in Table 1 has been defined with an average shear wave velocity to a depth of 30 m (V s,30 ) falling between 213 m/s and 724 m/s, therefore corresponding to the types of mass B' or C' (according to Eurocode 8 [45]) which are mixed into the set. This latter comprehends highly scattered values of magnitude M w , ranging from 5 to 7.5, and Joyner-Boore distance included between 0 km and 50 km. The set includes about 50% of near-fault (19 records), ranging from an epicentral distance of 0-10 km, and 50% of far-field records (17 records), or records with an epicentral distance greater than 10 km. The selection exhibits a prevalence of three fault mechanisms: normal, reverse, and strike-slip. Since the frame model used in this study is 2D (see Section 3), two orthogonal directions of the same seismic event are avoided. Record selection also comprehends a wide range of IM and distributed values of DCR LS , with at least one-third of the values greater than 1 [23]. The original set of 36 records was split into two subsets depending on the epicentral distance as shown in Figures 1 and 2, and these were studied separately to investigate the different effects of near-fault versus far-field ground motions. It is well known that the proximity to the fault renders the same ground motions (NF) different from ordinary (FF) ground motions [46]. The near-fault records selection avoids including pulse-like ground motions. Impulsive signals have been identified using the open-source algorithm proposed by Shahi and Baker [47,48]. This can identify pulses at arbitrary orientations using continuous wavelet transforms of two horizontal orthogonal components of a ground motion to identify the orientation that may contain a pulse.

Cloud-Based Analysis
The cloud-based Analysis (CA) is particularly suitable to assess structural fragility both for the simplicity of its formulation and for the low required computational effort. Conversely, it is extremely sensitive to the record selections and based on a few simplifying assumptions, such as fixed standard error of the regression [23,31,32,39,41,49].
CA adopts a linear regression model in the logarithm scale to fit the pairs of demand to capacity ratio (DCR LS ) and IM, where DCR LS are calculated through nonlinear analysis. The regression-based probability model describes the DCR LS for a given IM level and can be evaluated by Equations (3) and (4): where E[lnDCR LS |IM] is the expected value for the natural logarithm of DCR LS given IM, and η DCRLS|IM and σ lnDCRLS|IM are the median and logarithmic standard deviation for DCR LS given IM, respectively. The constants lna and b are the linear least square regression coefficients. Finally, the structural fragility obtained based on the CA is (Equation (5)): where Φ is the standard Gaussian cumulative distribution function. In this work, two different methodologies are adopted to find the relationships of IM versus DCR LS for the structure under investigation. Namely, in one case, time-history analyses are employed to evaluate the demand D jl (demand of the jth element of lth mechanism) at each time step. Within this paper, this approach is named dynamic cloud analysis, or shortly D-CA.
In the other case, the demand D jl is computed by adopting the modal pushover analysis (MPA). This approach is referred to as MPA-CA. The modal pushover analysis [50,51] is a nonlinear static procedure based on static analysis of the structure subjected to lateral forces distributed over the building height according to nth modal shape. Chopra and Goel [51] showed that this procedure is accurate enough for practical application. The MPA procedure used in this work adopts the capacity spectrum method (CSM), a nonlinear static analysis procedure to assess the seismic vulnerability of buildings originally proposed by Freeman [52]. The procedure permits finding a correlation between earthquake ground motions and building performance [53] (ATC, 1982) comparing a response spectrum (representing structure demand) and a pushover curve (representing building capacity) by an iterative procedure. This latter aims to the definition of the performance point (PP), which represents the state of maximum inelastic displacement of a building for a given seismic event. To plot the two curves in the same chart, RS and pushover curves are transformed into an acceleration displacement response spectrum (ADRS). The whole procedure employed can be summarized in the following steps (see Figure 3):

Incremental Dynamic Analysis (IDA)
In IDA, a nonlinear structural model is subjected to a set of scaled ground motion records (accelerogram aλ), each scaled to multiple levels of a monotonic scalable intensity measure such as Sa, PGA, PGV (herein IM = Sa(T1,ξ = 5%))."As-recorded" unscaled time histories are scaled by using a non-negative scale factor (λ) to obtain a scaled accelerogram Comme letter to
Evaluate the capacity via pushover curves in terms of base shear V b,n versus top displacement u r,n for the MDOF structure subjected to lateral forces profile proportional to nth mode shape; 3.
Determine maximum demand in terms of top displacement via the capacity spectrum method: • Convert the pushover curve of the nth mode shape to a capacity curve in the ADRS format by (Equations (6)- (9)): where M tot is the total mass of the structure, φ n is the nth natural vibration mode, φ n,r is the amplitude of the nth natural vibration mode at the roof of the structure, and Γ n and α n are the modal participation factor and modal mass of the nth mode, respectively; • Convert 5% damped response spectrum from the standard pseudo-acceleration S a versus the period of vibration T format to the ADRS format by (Equation (10)): where S D (T) is the displacement spectrum; • Plot demand and capacity diagrams together in the ADRS space. Determine the bilinear capacity curve. Iteratively determine the displacement demand for the nth mode shape. In this step, the dynamic analyses of a sequence of equivalent linear systems with successively updated values of equivalent viscous damping are involved; • Reduce the elastic spectrum by the equivalent viscous damping (Equation (11)); η = 10 5 + υ eq (11) • Determinate the performance point or the maximum expected demand in terms of top displacement;

4.
Convert displacement demand found in step 3 to global top displacement and individual component of local deformation (i.e., interstory drift) for the nth mode shape; 5.
Evaluate maximum demand to capacity ratio values according to Equation (1). If two or more modes of vibration are considered, combine the local modal responses according to the square-root-of-sum-of-squares (SRSS); 6.
Estimate parameters of the linear regression model in the logarithm scale to fit the pairs of demand to capacity ratio (DCR LS ) and IM; 7.
Draw structural fragility curve according to Equation (5).

Incremental Dynamic Analysis (IDA)
In IDA, a nonlinear structural model is subjected to a set of scaled ground motion records (accelerogram a λ ), each scaled to multiple levels of a monotonic scalable intensity measure such as S a , PGA, PGV (herein IM = S a (T 1 ,ξ = 5%))."As-recorded" unscaled time histories are scaled by using a non-negative scale factor (λ) to obtain a scaled accelerogram a λ , in which amplitudes are scaled without changing the frequency content of signals. The output of the analysis is represented by a collection of IDA curves, which are a plot of the recorded DCR LS (DM) against S a (T 1 ,ξ = 5%) (IM), all parameterized on the same IMs and DM [19].
Among all the analytical methods to develop the fragility based on IDA, the following is one of the simplest proposed [54]: In an EDP-based interpretation of the fragility (Equation (12)), the conditional probability of exceeding a limit state given an IM, herein the spectral acceleration, is equal to the probability of the demand to capacity ratio of exceeding 1 for a given Sa.
However, it is possible to express the fragility also as the complementary cumulative distribution function or "IM-based fragility" (Equation (13)). Incremental dynamic analysis is well suited to be represented by IM-based derivation of fragility. This interpretation expresses the seismic fragility as the probability of spectral acceleration values-denoted as S a DCR=1 and defined by intercepting all the IDA curves with the DCR LS = 1-to be smaller than a given value. DCR LS = 1 represents the threshold of a limit state LS and the intersection provides the empirical distribution of the random variable (IM), to which a model such as the lognormal appearing in Equation (12) can be fitted.
In this Equation (13), φ denotes the standard normal (Gaussian) cumulative distribution (CDF) of two-parameters (median or log of mean η and standard deviation β) estimated by the second-moment method or "METHOD A" (Equation (14)) proposed by Porter [54].

Incremental Modal Pushover Analysis (IMPA)
IMPA is a novel nonlinear static procedure proposed by Bergami and others first for buildings [40,55,56], and later adjusted also for bridges [57][58][59]. This procedure takes advantage of the simplicity of static analysis, but at the same time it grants the definition of the seismic demand for a certain range of intensity levels by scaling down response spectra. Conceptually, the procedure to find the maximum expected demand for the jth element of lth mechanism D jl is the same as that describe for MPA-CA in sub-Section 2.3.1, yet in IMPA procedure response spectra are scaled to multiple levels of a chosen monotonic scalable intensity measure as S a , PGA, PGV (herein IM = S a (T 1 ,ξ = 5%)) ( Figure 4). For each intensity level, the performance point (P.P.) can be determined and the demand measure D jl combined if two or more modes of vibration are considered, to define a multimodal DCR LS . The output of the analysis can be represented as a collection of "IMPA" curves, which are a plot of the recorded DCR LS against IM, wholly similar to IDA curves. The intersection of IMPA curves and the chosen threshold of the limit state (DCR LS = 1) generates the empirical distribution of the random variable (IM) for the probabilistic model of the fragility aforementioned in Section 2.3.1 (Equations (12) and (13)).

Frame Description
The transverse frame modelled and analysed in this study comes from an actual school building in Norcia (Italy, 42.7941 • North latitude, 13.0963 • East longitude). The building, which originally consisted of three aligned blocks, was later joined into a single complex during the various refurbishment works. The building consists of a reinforced concrete (RC) frame structure with a footprint of 12.80 × 59.80 m and a maximum height (from the foundations), corresponding to the roof beams, of about 16 m. The building consists of a one-floor basement, a ground floor, three storeys and an attic above ground. The inter-story height is 3.50 m for the basement floor and ground floor, 3.30 m for the other three floors, and 2.5 m for the attic ( Figure 5). Built in 1962, the school has survived a variety of seismic events before the 6.5 magnitude central Italy earthquake that caused severe damage to structural and especially non-structural parts in 2016. In accordance with the construction methods of the time of construction, the structure was designed using 2D models schematizing the reinforced concrete frames in the transverse direction of the building. Although the legislation of the time did not explicitly require it, the designer also took into account the seismic action by applying an acceleration of 0.07 g.  The frame under study belongs to one of the two lateral blocks ( Figure 5) and is a two-bay (5.65 m and 5.56 m span) regular cross frame. It is considered to be fixed at the base, and basement and soil-structure interaction has not been considered.
The geometry of the frame, the column and beam geometry, and the reinforcement details are shown in Figures 5 and 6. The permanent structural load (G1) and permanent non-structural load (G2) are calculated as G1 + G2 = 5.1 kN/m 2 (from ground floor to 2nd floor), G1 + G2 = 4.1 kN/m 2 (3rd floor), G1 + G2 = 4.22 kN/m 2 (roof beams), the live load is taken as Q1 = 3 kN/m 2 (from ground floor to 2nd floor), Q1 = 1 kN/m 2 (3rd floor), Q2 = 1.8 kN/m 2 (roof beams) and taken as concentrated gravity on the columns at the edge of each floor. Each floor was assigned a seismic mass equal to 1/6 of the total mass of one of the three original blocks of the building. Fundamental and second periods of the frame have been evaluated as T 1 = 0.62 s and T 2 = 0.21 s, respectively.

FE Model Description
The nonlinear FE model of the analysed frame was developed in the OpenSEES platform [60]. To account for the nonlinearity, "Beam With Hinges Element", already available in the OpenSEES library, was used to model columns and beams. This element adopts a lumped plasticity formulation with plastic hinges at the end of the element connected by an elastic link. This means that all nonlinearities are concentrated at the ends of the elements and can be only propagated along the length L p of the plastic hinge, unlike the distributed plasticity formulation where they may spread along the whole element. Therefore, the length of the plastic hinge L p plays an important role in avoiding the concentration of strain at the element ends. In this work, in good agreement with the L p evaluated by the equation proposed by Priestley and Park [61], it is assumed to be equal to the cross-section height. The two-point Gauss integration was used on the element interior, while two-point Gauss-Radau integration was applied over lengths of 4L p at the element ends. A total of six integration points were used [62,63]. To account for non-linearity, a four-point moment-curvature relationship was assigned to the element ends. The four-point law accounting for crack, yield, ultimate failure and collapse state (80% of ultimate failure) was evaluated through the software Response2000 [64], which can simulate nonlinear sectional behaviour by assuming a suitable law for the material. Due to the dependence of the sectional response on the applied axial load, it is assumed to be zero for beams, while it is estimated for the column considering their area of influence. Shear failure was not considered in the model. The Newton line search method was considered as the solution algorithm for the time-history analyses, which increased the effectiveness of the Newton-Raphson algorithm by introducing line search to solve the nonlinear residual equation. The tolerance and maximum number of iterations used were OpenSees default values [62]. Newmark integrator has been used and the convergence test was the normal displacement increment. The tolerance of the test is equal to 10 −7 and the number of maximum iterations is 50. Rayleigh damping is adopted to account for energy dissipation.
In IDA, the scale factor λ was chosen to scale the spectral acceleration at the fundamental period S a (T 1 ,ξ = 5%), which was scaled to IM = a λ ∈[0.1 g, 0.8 g] with ∆ aλ = 0.1 g. Similarly, in IMPA, the response spectra were scaled multiple times to obtain scaled spectral acceleration at the fundamental period S a (T 1 ,ξ = 5%) equal to IM= S aλ ∈[0.1 g, 0.8 g] with ∆S aλ = 0.1 g. The mass damping coefficient and the stiffness damping coefficient of the Rayleigh damping are evaluated by considering the first and the second natural frequency of the case study. The percentage of critical damping is equal to 5%.

Nonlinear Static Analysis Results
The pushover analysis has been performed in displacement control to reach a target displacement of 350 mm. The number of steps to reach the target displacement is 350 calculation steps with an increment of 1 mm each step. Figure 7 shows the capacity curves obtained applying two load distributions proportional to the first and second modal shapes, respectively. The capacity curve for the first load distribution reaches a maximum base shear of 492 kN corresponding to a top displacement of about 220 mm. The sequence of the plastic hinges activation with the relative calculation step is shown in Figure 8. A plastic hinge is activated when the reinforcement of the section reaches and exceeds the yield point: the first plastic hinges activated in the columns rather than in the beams, particularly in the upper stories which, according to Eurocode 8, have weak column-strong beam connections.

Nonlinear Dynamic Analysis Results
A total number of six subsets of near-fault records and six subsets of far-field records were examined. The subsets were randomly generated from the main set of 210 records to comply with the general rules for record selection explained in Section 2.2. The mean and standard deviation of the normal distributions for magnitude, epicentral distance and spectral acceleration of each record selection were given in Table 2. As expected, the results show that IDA and IMPA are less dependent on record selection, with the mean values of fragility curves ranging from 0.463 g to 0.525 g and from 0.387 to 0.432, respectively. In contrast, D-CA and MPA-CA show greater dependence on record selection, with mean values between 0.479 g and 0.724 g and between 0.321 g and 0.621 g, respectively ( Figure 9). It can be pointed out that IMPA provides the most conservative results for all the datasets studied, as shown in Table 3. Moreover, this methodology seems to be the more accurate with respect to IDA in estimating vulnerability for the far-field record selections.  Finally, Figure 8 shows the interstory drift at different steps of the analysis, with an interval of 50 steps. It shows a large concentration of interstory drift in the second and third stories, while the fourth one is moving almost rigidly.

Nonlinear Dynamic Analysis Results
A total number of six subsets of near-fault records and six subsets of far-field records were examined. The subsets were randomly generated from the main set of 210 records to comply with the general rules for record selection explained in Section 2.2. The mean and standard deviation of the normal distributions for magnitude, epicentral distance and spectral acceleration of each record selection were given in Table 2. As expected, the results show that IDA and IMPA are less dependent on record selection, with the mean values of fragility curves ranging from 0.463 g to 0.525 g and from 0.387 to 0.432, respectively. In contrast, D-CA and MPA-CA show greater dependence on record selection, with mean values between 0.479 g and 0.724 g and between 0.321 g and 0.621 g, respectively ( Figure 9). It can be pointed out that IMPA provides the most conservative results for all the datasets studied, as shown in Table 3. Moreover, this methodology seems to be the more accurate with respect to IDA in estimating vulnerability for the far-field record selections.  The two subsets of records studied in-depth (Sel. 1 NF and Sel. 7 FF (Table 1)) comprise two different records of the real seismic event to which the case study was exposed in 2016, the 6.5 magnitude Central Italy earthquake, which, as mentioned, caused severe damage to structural and especially non-structural component of the school. Figures 10 and 11 shows IDA curves and IMPA curves: the curves of IDA referring to the Norcia 2016 earthquake exceed the threshold of DCR LS for the intensities of S a equal to 0.5 g and 0.57 g for NF and FF selection, respectively. The results of D-CA show that, as expected, the DCR LS threshold for the SLV limit state is exceeded for both the near-fault record and far-field record, with DCR LS values of 2.57 and 1.01, respectively ( Figure 12). Similarly, even higher DCR LS are obtained from the cloud analysis based on nonlinear static analysis ( Figure 13).   To check the consistency of modal and multimodal IMPA and pushover-based cloud, the developed fragility curves (Figure 14) are compared with those of IDA. The accuracy of the prediction of the different fragility models with respect to IDA is quantified by normalized root-mean-square deviation (RMSD). It is evaluated according to the following Equation (15), assuming that the values derived by Equation (12) (IDA) are the reference ones: where n is the number of points, and y i . andŷ i are the predicted and reference probability of exceeding the considered limit state (LS), respectively.  The following Table 4 shows the comparison in terms of the percentage variation of the median, 16% and 84% fractiles of the fragility curves with respect to IDA and the absolute values of the standard deviation of each method. Regarding the selection of NF records, the cloud method appears solid in estimating the 50%, 16% and 84% fractiles compared to IDA with the smallest normalized root-mean-square deviation of 8%. However, when it comes to far-field records, IMPA appears to be the most accurate methodology for estimating fragility. The inclusion of two or more vibration modes in the assessment of the maximum multimodal DCR LS does not seem to be essential, as it leads to very conservative results.

Conclusions
This paper compares fragility curves obtained by various known static and dynamic nonlinear procedures. Incremental modal pushover analysis (IMPA) is proposed as an alternative to IDA, which is currently considered the most reliable method, to determinate IMPA curves and thus seismic fragility. Similarly, both MPA and NL-THA are used to determine capacity in the well-known cloud method. For this study, it was necessary to perform a relatively small number of nonlinear time histories using two different data sets. These differ in the range of Joyner-Boore distance (R jb ) and are scattered in a range of magnitude.
The following conclusions can be drawn from this study that is limited to a simplified 2D frame model and a small set of records: • Comparison of fragility curves shows that, in the case of methodologies distinguished by scaling (in terms of accelerograms (IDA) or response spectra (IMPA)), near-fault records and far-field record selections have led to nearly equivalent results. In contrast, the results in terms of fragility when using records without scaling, i.e., in D-CA and MPA-CA, show clear differences in the whole range of intensities; • Fragility curves that only consider the contribution of the first mode in determining DCR LS have led to more accurate results in relation to IDA, so the inclusion of higher mode contributions does not seem to be essential for low to medium buildings (up to nine stories [20]); • A total number of 12 subsets have been extracted from the main 210 set of records and exanimated, but results are not fully reported in this paper. The results have shown that D-CA leads to a smaller vulnerability than MPA-CA and IMPA in all selection and for the whole range of intensities. Methodologies based on the pushover analysis, on the contrary, have led to more conservative results, especially for 16% and 50% fractiles; • IDA shows less sensitivity to record-to-record variability. It should be noted, however, that IMPA, despite its slightly greater sensitivity, has the advantage of a large reduction in the computational effort required to perform the structural analysis. In IMPA, the total time required relates mostly to the post-processing phase, which is no different for small 2D frames or more complex 3D buildings.
A more comprehensive validation is needed to confirm the obtained results and draw more general conclusions.

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