IM-based seismic reliability assessment of the pre-code masonry building stock in metropolitan area of Lisbon

IM-based seismic reliability assessment of the pre-code masonry building stock in metropolitan area of Lisbon Vasco Bernardo a,*,1 , Alfredo Campos Costa b,2 , Paulo B. Lourenço a,3 a ISISE, Department of Civil Engineering, University of Minho, 4800058 Guimar˜ aes, Portugal b National Laboratory for Civil Engineering, Av. do Brasil 101, 1700075 Lisbon, Portugal ARTICLE INFO Keywords: Pre-code masonry buildings Seismic probabilistic approach SAC/FEMA formulation Reliability analysis Code calibration ABSTRACT Earthquakes have a long history of causing catastrophic damage to communities, resulting in structural collapses, loss of life, and economic turmoil. To enable informed decision-making and reduce the impact of these events in earthquake-prone regions, seismic risk studies provide relevant information to support stakeholders in the implementation of effective risk-based policies. The present work addresses the seismic reliability of the pre-code masonry building stock in the Metropolitan Area of Lisbon, which is the region of Portugal that faces the highest seismic risk due to the coexistence of moderate to high seismic hazard and highest demographic-economic exposure. The adopted general framework combines several hazard studies developed for the region under investigation and a synthetic database of masonry buildings representative of the pre-code building stock in Lisbon. Through analytical–numerical probabilistic approaches, new second-order hazard solutions with structural dependency are derived for the mean annual frequency of limit-state exceedance, which can be integrated into national application documents for Eurocode 8. In light of these results, the reliability assessment of the building stock is conducted in several Local Administrative Units by means of an improved SAC/FEMA formulation. The study represents the first comprehensive investigation of its kind in this region, providing essential information to define appropriate target safety level for code calibration and support future risk studies. 1. Introduction Over the past few decades, a notable emphasis has been placed on the seismic behavior of buildings, motivated by the concern to preserve built cultural heritage and safeguard human lives. This concern is particularly prominent in the case of unreinforced masonry (URM) structures, which are highly vulnerable when subjected to seismic ground motions (e.g., [1–4]). In this regard, the performance-based earthquake engineering (PBEE) framework has a paramount importance in the earthquake engineering community to define more appropriate mitigation policies and manage the post-earthquake functionality of structures. This is well demonstrated through the approach embraced by the Pacific Earthquake Engineering Research (PEER) Center on the estimation of the mean annual frequency (MAF) of exceeding a given limit-state (LS) defined in terms of system-level decision variables (e.g., economic losses, fatalities, repair cost, and downtime of facilities) [5]. The PEER’s PBEE methodology addresses four main stages: (i) hazard analysis –probabilistic seismic hazard analysis (PSHA) [6] to estimate the MAF for different intensity measures (IMs); (ii) structural analysis – probabilistic-based seismic response of the structure (e.g., numerical analysis –static or dynamic) for a certain IM to predict the Engineering Demand Parameter (EDP) [7–9]; (iii) damage analysis –estimation of physical damage at the component or system-level (structural or nonstructural) as a function of the nonlinear structural performance (expressed by the EDP) and definition of fragility curves [10] to described the relationship between the damage measures DM (e.g., damage levels related to repair measures or occupancy safety) and EDP; (iv) loss analysis –losses estimation by converting the damage information from (iii) into decision variables (DV) [5]. The simplification of such approaches have been provided by SAC/FEMA methodology [11], which introduced a comprehensive framework based on closed-form solutions for calculating the MAF of LS exceedance. * Corresponding author. E-mail address: [email protected] (V. Bernardo). 1 0000-0002-5071-8155. 2 0000-0001-6869-9186. 3 0000-0001-8459-0199. Contents lists available at ScienceDirect Structural Safety journal homepage: www.elsevier.com/locate/strusafe https://doi.org/10.1016/j.strusafe.2024.102514 Received 19 January 2024; Received in revised form 27 April 2024; Accepted 12 July 2024 Structural Safety 111 (2024) 102514 Available online 14 July 2024 0167-4730/© 2024 The Author(s). Published by Elsevier Ltd. This is an open access article under the CC BY license ( http://creativecommons.org/licenses/by/4.0/ ). The SAC/FEMA format offers a straightforward expression to convolve seismic hazard with structural response. It enables the MAF estimation of exceeding a LS associated with structural response while considering both aleatory and epistemic uncertainties (e.g., [12–14]). In literature, some studies can be found employing and expanding upon the initial SAC/FEMA methodology (e.g., Kazanti et al. [15], Vamvatsikos et al. [16], Bernardo et al. [17], Jalayer [9], FajFar and Dolsek [18], Franchim and Pinto [19]). It is noteworthy to emphasize the contribution of Vamvatsikos and Dolsek [20] on the improvement of the original SAC/FEMA expression to include time-dependent capacity over the lifetime of the structure. Furthermore, the work of Vamvatsikos [21] stands out, as it significantly improved the accuracy of the original closed-form SAC/FEMA solution. This goal was achieved by replacing the original first-order approximation for seismic hazard with a secondorder fit in log–log space to capture the hazard curvature, which was identified as the main accuracy error in the original SAC/FEMA formulation (Aslani and Miranda [22]; Bradley and Dhakal [23]). Rom˜ ao et al. [24] also provided alternative closed-form solutions to improve the accuracy of the original SAC/FEMA MAF formulation by employing nonlinear models in log–log domain to represent the seismic hazard distribution over the whole range of IMs. Taking advantage of the second-order improvement of SAC/FEMA probabilistic format, reliability analysis was conducted here on the precode masonry building stock in the Metropolitan Area of Lisbon (MAL). The vast majority of these buildings are characterized by regular geometry and composed of rubble stone masonry and/or brick masonry supporting the timber floors. Ancient buildings built after the 1755 Lisbon earthquake (Mw =8.5) display the first anti-seismic features in Europe, namely with the inclusion of the internal wooden cage –Pombaline cage [25]. This construction technique was adulterated over the years with the adoption of low-quality materials and marked by the abandonment of the Pombaline cage. Later, the timber floors were gradually replaced by concrete slabs. These typologies were only designed to support gravity loads since the first design code for buildings safety against earthquakes (RSCCS) emerged in 1958 [26]. Further details about the building stock in the MAL region and its structural characterization can be found in Bernardo et al. [27]. The MAL is considered the region of Portugal with highest seismic risk (e.g., Sousa [28], Silva et al. [29]) given the coexistence of a moderate to high seismic hazard, high population density and high building stock exposure. In particular, the seismic risk of MAL was evaluated by Campos Costa et al. [30] and Sousa et al. [31]. The former indicates that approximately 36 % of the total economic risk corresponds to URM buildings, while the latter highlighted that total losses in MAL region are mainly affected by seismic events at short distances. Despite the valuable insights provided by these studies, they relied on expert judgment to assess the seismic vulnerability of the building inventory. Recently, Bernardo et al. [32] evaluated the seismic vulnerability of the pre-code masonry building stock by means of numerical analysis on several representative archetypes. Making use of this synthetic database, the present work comprises the main following steps: (i) estimation of the fundamental period for 18.000 representative masonry buildings in the MAL region; (ii) computation of probabilistic seismic hazard analysis (PSHA) for the region under investigation combining different seismic source zones, including the European seismic hazard reference model [33] for consideration in the upcoming version of Eurocode 8; (iii) first-mode response estimation of the entire building database located in the various Local Administrative Units (LAU-parishes) of the MAL region; (iv) derivation of new second-order hazard coefficients for the MAF estimation based on the structural response for different building typologies and return periods; (v) computation of the reliability index in the various parishes considering aleatory and epistemic uncertainties through an improved SAC/FEMA approach. Based on such outcomes, this study provides new second-order hazard coefficients for the MAF estimation of pre-code buildings located in different regions of the MAL. These results can be included in the National Application Documents (NAD) of Eurocode 8. It is important to empathize that the current version of the NAD of the Eurocode 8 express the hazard function in terms of PGA through a first-order power law, with unique coefficients for the entire MAL region and without any structural response dependency. Furthermore, the results of reliability analysis carried out may serves as valuable support for regional mitigation strategies and define new target safety level for code calibration. Regarding this last, the current standards, namely the Eurocode 8, defines the structural safety levels based only on iso-hazard. Therefore, this study also plays an important role in this regard, as it allows to define reliability indices for the LAU (civil parishes) in the region, accounting for the building stock exposure. 2. Formulation for the mean annual frequency estimation The mean annual frequency (MAF) –λ–of a performance measure exceeding a specified threshold, can be described by the Total Probability Theorem [34] through: λ(DV) =∫∫∫G(DV|DM)|dG(DM|EDP)||dG(EDP|IM)||dλ(IM)| (1) This equation includes the conditional probabilities of the four stages introduced in Section 1, where G(X|Y)–generic notation –is the complementary cumulative distribution function for X given Y, and dλ(IM)is the absolute value of the derivate of the λof exceeding a certain IM. Eq. (1) can be decoupled by assuming Markovian independence of the four stages, resulting in the MAF of reaching or exceeding a damage state DM, or a specified demand level EDP defined by: λ(DM) =∫∫G(DM|EDP)|dG(EDP|IM)||dλ(IM)| (2) λ(EDP) = ∫G(EDP|IM)|dλ(IM)| (3) Taking IM in Eq. (1) as the input variable of a seismic fragility function regarding DM, the conditional probability G(DM|IM)results in Eq. (4). Replacing Eq. (4) into Eq. (3), the MAF of reaching or exceeding DM for a specified IM is expressed by Eq. (5). G(DM|IM) = ∫G(DM|EDP)dG(EDP|IM)(4) λ(DM) = ∫G(DM|IM)|dλ(IM)| (5) Performance measures as EDP, DV or DM can express the limit state (LS), with corresponding conditional probabilities defined by a general fragility function FR,LS, yielding the conditional probability of exceeding a LS for a given ground motion intensity level IM =x. Note that, a lognormal cumulative distribution function is commonly used to derive fragility functions [35], resulting in: FR,LS =P[LS|IM =x] = Φ[1 βLS ln(x IMLS)] (6) where, Φis the standard normal cumulative distribution function, IMLS the median value of the lognormal distribution of IM, and βLS is the logarithmic standard deviation. Hence, the MAF of exceeding a specified LS for a certain IM is given by Eq. (7), which represents the convolution between the fragility and seismic hazard for a given IM =x.FR,LS can be determined through nonlinear structural analyses and performing probabilistic seismic hazard analysis [9,36]. λLS =∫FR,LS(x)|dλIM(x)| (7) V. Bernardo et al. Structural Safety 111 (2024) 102514 2 Given the wide use of Eq. (7), some authors (e.g., Cornell et al. [37] and Jalayer [9]) have derived analytical closed-form solutions for the LS probability using a relationship between the ground motion IM and seismic demand (D) through a first-order power-law approximation: D=a(IM)b(8) where the regression parameters (a,b) can be estimated by linear regression in log-space from the results obtained by probabilistic seismic demand analysis (e.g., [7,9,20]). Considering FR,LS of a specific LS in terms of structural capacity (C), FR,LS =FR,C, and assuming both Dand Clognormally distributed, with βD|IM corresponding to the dispersion of median seismic demand over the range of IM, and βCthe dispersion value of median seismic capacity, results in: FR,LS =FR,C=P[D≥C|IM =x] = Φ⎡ ⎢ ⎣ ln(D)− ln(C)  β2 D|IM +β2 C √⎤ ⎥ ⎦(9) Regarding the seismic hazard term λIM in Eq. (7), a first-order power-law model was also proposed by Cornell [38] to provide the MAF of exceeding a particular IM: λIM =P[IM ≥x] = k0(IM)−k1=k0exp(− k1lnIM)(10) where k0and k1are positive coefficients representing the intercept and slope, respectively, of the linear relationship of seismic hazard in the log–log space. Finally, substituting Eq. (10) and Eq. (6) into Eq. (7), the result of the integral is the well-known expression for probabilistic seismic demand analysis [39]: λLS =P[D≥LS] = λIM(IMLS)exp(1 2 k2 1 b2β2 D|IM)(11) where IMLS represents the ground motion intensities associated with the LS. If the LS is defined in terms of capacity, substituting Eq. (10) and Eq. (9) into Equation Eq. (7) yields: λLS =P[D≥C] = λIM(IMC)exp[1 2 k2 1 b2(β2 D|IM +β2 C)](12) with IMCthe corresponding ground motion intensities to the structural capacity. Eq. (12) is the basis of the demand and capacity factor design method, employing a probabilistic performance-based approach and used in SAC/FEMA. Some authors (e.g., [22,23]) pointed out accuracy issues in the original SAC/FEMA formulation. These concerns are mainly related to the seismic hazard using a first-order approximation, which may lead to inaccuracies, especially when dealing with significant curvature in hazard functions. To overcome this limitation, Vamvatsikos [21] proposed a new shape for the seismic hazard model, alternative to Eq. (10) and employed in this study, through a second-order power law approximation: λIM =P[IM ≥x] = k0exp(− k2ln2IM −k1lnIM)(13) wherein k0,k1>0 and k2≥0 represent, respectively, the intercept, slope, and curvature of the seismic hazard. As a result, the closed-form solution of Eq. (11) can be written as: λLS =P[D≥LS] =  p √k1−p 0[λIM(IMLS)]pexp(1 2pk2 1 b2β2 D|IM)(14) p=1 1+2k2β2 D|IM (15) where pis a correction factor to account for the contribution of the hazard function curvature. Note that, for p=1, i.e., non-curvature of the seismic hazard model (k2=0), Eq. (14) simplifies to Eq. (11). Further details regarding performance evaluation with second-order hazard approximation and comparison with former SAC/FEMA solutions are discussed in Vamvatsikos [40]. 3. Characterization of the Metropolitan Area of Lisbon (MAL) The Metropolitan Area of Lisbon (MAL), Fig. 1, comprises 18 municipalities and 211 parishes, covering a total area of 2,957.5 km 2 and holding a population density of approximately 950 people/km 2 . The MAL region is characterized by the highest population density and economic concentration in Portugal, contributing with around 34 % to the GDP (2020) of the country. Although MAL has not experienced strong earthquakes in recent years, its history is marked by seismic occurrences that resulted in significant damage to masonry structures [41], such as the 1755 Lisbon earthquake (Mw =8.5), the 1909 Benavente earthquake (Mw =6.3), and the 1969 Algarve earthquake (Mw = 7.8) [42]. 3.1. Masonry building sock and representative archetypes The MAL region comprises 434,600 residential buildings and 1,423,654 dwellings, with masonry buildings accounting for approximately 35 % of the building stock. The subsequent analysis focuses on residential pre-code URM buildings up to five stories high, which constitute most of the building stock in the MAL region. For this purpose, representative archetypes were derived from an extensive collection of original detailed drawings (see Bernardo et al. [27]), which allow to characterize the geometric parameters summarized in Table 1 (plan dimensions Lx and Ly; ground and upper floor stories height H 0 and H n ; openings ratio OR subscript: front (ORF) and back (ORB) facade; interior walls density IWD; walls thickness Th subscript: facades (1), lateral side (2), interior (3), partition (4); average walls thickness reduction on the façade AWTR). Based on this information, archetypes with multiple indexes were derived, see Table 2. For convenience, only buildings with five stories height are shown, however, buildings with one to four stories follow the same geometry layout. Further details can be found in Bernardo [43]. 3.2. Soil classes and geographic distribution in MAL The MAL region is composed of different soil layers with different dynamic properties, which can produce significant modifications on the frequency energy contents of the surface seismic ground motion. A prior study conducted by Campos Costa et al. [30] in the MAL region classified several soil columns units from a geotechnical survey into three main types of soils (hard, intermedium and soft). This classification is based on the soil classes defined in EN 1998-1-1 (EC8) according to the mean shear wave velocity (Vs,30) of the soil profile units. A similar classification was adopted in the present study (see Table 3). Considering the spatial resolution in terms of parishes, the most representative soil layer in each parish (see Fig. 1) was assigned by weighting the total construction area of the pre-code masonry building stock located in the same ground. The contribution of others typological classes was excluded for this purpose. It was concluded that around 28 %, 42 % and 30 % of the masonry building stock is located in a hard, intermediate and soft soil, respectively. The amplification of the ground motion at the surface was simplified by considering the Portuguese national annex of the current EC8. This requires a non-dimensional parameter soil factor S to indirectly account for the local site conditions, and for each soil class, computed from: ag≤1.0:S=Smax (16) V. Bernardo et al. Structural Safety 111 (2024) 102514 3 1.0<ag<4.0:S=Smax −(Smax −1) 3(ag−1)(17) ag≥4.0:S=1.0 (18) where ag(m/s 2 ) is the spectral ground acceleration and the values of Smax are provided in Table 3. 4. Probabilistic seismic hazard analysis The definition of the hazard curves for the subsequent analysis were derived from the probabilistic seismic hazard analysis based on the classic Cornell approach [44]. The annual frequency λof exceeding a hazard level hfor a certain region considering the contribution of N seismic sources zones is given by: λ(H>h) = ∑ N k ν k∫∫P(H>h|m,r)kfM(m)kfR(r)kdrdm (19) where P(H>h|m,r)kis the conditional probability of exceeding hin a source zone k, characterized by attenuation law(s) that described the amplitude of seismic ground motion at a specified location for an event of magnitude Moccurring at a distance Rfrom the source; fMand fRare, respectively, the probability density functions for M(following a Gutenberg-Richter model) and R, both considered independent variables; ν kis the mean annual rate of occurrence of an earthquake with M>mmin, assuming a Poissonian process. Regarding the seismic source zones (seismogenic zones), i.e., regions with similar seismological/tectonic/geological attributes, different models are proposed for mainland Portugal based on the historical and instrumental earthquakes (main events only). A detailed discussion of different models is presented in Carvalho et al. [45]. The seismogenic regions considered in the present study follow three models, which were equally weighted in the PSHA given the lack of other information (Fig. 2): (i) Model A, adopted in the development of the Portuguese national annex of EC8 [46]; (ii) Model B proposed by the recent SERA project [33]; and (iii) Model C in the scope of ERSTA project [47]. Due to the seismotectonic characteristics of Portugal mainland, the seismic activity is distinguished by scenarios with offshore epicenters (interplate) or events that mainly occur inland (intraplate). The latter are characterized by moderate magnitudes, short source-site distance and higher frequency content, having the potential to result in more significant losses in the MAL region [28]. In particular, the Lower Tagus Valley Fault (LTVF), see Fig. 2, with a length of approximately 110 km (SW-NE), is the most relevant seismic event for the region under interest. Fig. 1. Geographic location of (left) and distribution of pre-code building stock (center) and soil classes (right). Table 1 Summary of the statistical properties for the geometric parameters collected [27]. Moments L x [m] L y [m] IWD[-] H 0 [m] H n [m] OR F [-] OR B [-] Th 1 [m] Th 2 [m] Th 3 [m] Th 4 [m] AWTR [m] Mean µ12.6 12.1 0.054 3.23 3.01 0.23 0.21 0.47 0.34 0.21 0.14 0.11 Std. deviation σ 5.00 4.1 0.01 0.42 0.24 0.08 0.08 0.14 0.11 0.05 0.02 0.06 Table 2 Representative archetypes and size of the URM building stock in the MAL region (plan dimension indicated by L x £L y ). A1 A2 A3 B1 B2 B3 C1 C2 C3 Table 3 Soil classes and amplification factors. Soil Class Bedrock Hard Intermediate Soft Shear wave velocity Vs,30 [m/s] >800 360–800 180–360 <180 Max. soil coefficient Smax [-] 1.00 1.35 1.60 2.00 V. Bernardo et al. Structural Safety 111 (2024) 102514 4 For this reason, only the seismic region integrating the LTVF were considered next. This region was described by the truncated GutenbergRichter law with the median parameters summarized in Fig. 2 [45]:mmax is the maximum magnitude; aand bare, respectively, the intercept and slope of the logarithmic relationship between the number of earthquakes and their magnitudes. For the seismic ground motion attenuation laws, the spectral acceleration at bedrock according to Ambraseys et al. [48], Atkinson &Boore [49], Akkar &Bommer [50] and Chiou & Youngs [51] was adopted. These are suitable for stable continental crust in Portugal mainland. The attenuation laws were also equally weighted in a logic tree. The computation of the PSHA was performed in the OpenQuake engine platform [52,53].Fig. 3 shows the mean hazard map, the response spectra for 10 % probability of exceedance in 50 years (return Fig. 2. Seismic sources and parameters adopted in the Gutenberg-Richter law: m max,aand b. From left to right: Model A, B and C, respectively. Fig. 3. Seismic hazard for the MAL region –bedrock (top) and surface (bottom). From left to right: mean hazard map, response spectra (RP =475 years) and hazard curves. V. Bernardo et al. Structural Safety 111 (2024) 102514 5 period RP =475 years) and the hazard curves in terms of peak ground acceleration (PGA). The multiple gray lines plotted define the median results of the PSHA in each parish evidencing the spatial variation in seismicity for the MAL region. For instance, statistics moments for spatial variability of the maximum acceleration response spectra (RP = 475 years; T ≈0.137 s) are summarized: bedrock with μ =0.32, σ = 0.03, skewness = − 1.70; surface with μ =0.37, σ =0.07, skewness = −0.37. Higher density of spectral curves is close to maximum acceleration values (negative skewness), namely for seismic hazard at bedrock. Comparing the results of Fig. 3 with the current version of the national annex of EC8, the mean values of PGA herein computed are between 0.11 g and 0.21 g (RP =475 years; bedrock). On the other hand, EC8 prescribes a global PGA value of 0.17 g for the seismic zonation that includes the MAL region. Therefore, EC8 seems to underestimate the seismic hazard in the surrounding areas of LTVF and overestimate the seashore NW regions. This brief observation highlights the importance to adopt a higher spatial resolution for the definition of new seismic zones in next code generations. Summing up, Fig. 4 shows the boxplot of PSHA results in terms of PGA for different RP with regards to the 16th, 50th and 84th quantile. The results were estimated in all parishes, where the spatial variability is also depicted by the box-whisker diagram with median values on the MAL region represented by fisheye marks. 5. Intensity-based structural response The estimation of the structural response is expressed in terms of IM as a function of the first-mode spectral acceleration for 5 % damping Sa(T1), which is widely used in literature for structures with short periods and mainly governed by the first-mode vibration and lower contribution of higher-modes. For this purpose, a large synthetic database generated by Bernardo et al [54] and comprising 18.000 representative masonry buildings was used. The classes of structures were represented with multiple index buildings (see Table 2) based on Monte Carlo simulations in which pre-defined statistics on material and geometrical characteristics were considered. The computation of the fundamental period T1is based on the Jacobi inverse algorithm implemented in the research version of TreMuri [55]. Fig. 5 depicts the empirical cumulative distribution function (ECDF) and the variability of T1for the entire database grouped by the number of stories. The estimated value for the coefficient of variation (cov) of T1 is approximately 0.30, independently of the number of stories, as shown in box-whisker plot. In what concerns the structural response in terms of Sa(T1),Fig. 5 also shows the spatial variability in the MAL region and for different levels of PGA, corresponding to the range of RP between 10 to 10.000 years. The soil amplification effects in the estimated values of Sa(T1)are also depicted. Note that, the first-mode response of each structure was treated per se, so the estimation of Sa(T1) values considers the entire database of buildings located in the different parishes and subjected to different PGA levels depending on the local hazard. 6. Estimation of the mean annual frequency 6.1. Second-order hazard fitting The mean annual frequency (MAF) estimation considers the results of the previous section and employs Eq. (14). First, a second-order power-law fit in log–log space is adopted to capture the hazard curve in Sa(T1) domain and extract the coefficients k0,k1and k2. Note that, these polynomial coefficients are conditioned to the pivot point (PP) required to estimate the MAF. For instance, if it is of interest the MAF of limit sate exceedance for the 475-years RP (rate per year λ0.00211), the fit should be centered in the corresponding PP and its vicinity. Indeed, the strict selection of the PP can reproduce a broad hazard domain, however, it is naturally dependent on the hazard shape function. Fig. 6 exemplifies the 2nd order hazard fitting for a given building (T1=0.24 s) in the database located in two distinct regions defined by PGA475(g) = [0.11, 0.21]. Three different PP were selected corresponding to the values of Safor RP (years) =[95, 475, 2000]. As can be noticed in Fig. 6, the fitting in Safor RP =475 years almost captures the entire domain in the upper hazard curve. The same observation is not valid for the lower hazard curve, since RP <100 years cannot be described by the fitting at the chosen PP. On the other hand, by defining the PP in Safor RP =95 and 2000 yrs., allows to capture the associated lower and higher hazard points, respectively. Hence, a local fitting is essential in the region of interest to estimate the polynomial regression coefficients with accuracy. Fig. 6 also shows the fitted coefficient for the mean curve that accounts the dispersion in response and spatial variability in the MAL region for RP =475 years. Box-diagrams of Fig. 7 and Fig. 8 show the distribution of the estimated regression coefficients from Eq. (13) (k 0 , k 1 and k 2 ) aggregated by number of stories, corresponding to the median PSHA. Therefore, the analytical hazard curves fitted corresponds to the IM-based median responses in Sa(T1) for RP =475 years. Disaggregation by parishes for the case of three-stories high buildings is also depicted in the same figures. Table 4 and Table 5 summarize the mean regression coefficients for the MAL region by considering different seismic hazard quantiles (16th, 50th and 84th). The results provided correspond to the fitting on Safor the MAF of 95-475and 975-years RP, which are usually required in the standards for the verification of the damage limitation, life safety and near collapse, respectively. To support a discussion of the above results, the corresponding timebased vulnerability curves (IM-response) are plotted (semi-log scale) in Fig. 9 for the median hazard. By comparing the results for Sa(T1= 0.01s), an approximation of the PGA indicated by the dashed line in the plot, it can be noticed that buildings with 2 and 3 stories tend to be more vulnerable for the return periods selected, i.e., the rate per year expected is higher in these structures for the same value of demand acceleration. It was found that these results clearly show a first-mode frequencydependence intercepting the spectral acceleration curves in the range of maximum acceleration, approximately between 0.1 s and 0.3 s. One the Fig. 4. Summary of seismic hazard results (PGA) for the MAL region (bedrock). V. Bernardo et al. Structural Safety 111 (2024) 102514 6 other hand, buildings with 1 and 5 stories are less vulnerable for the opposite reason. It is worth mentioning that the range of maximum Sa defined in the current Portuguese national annex of EC8 is located between periods of 0.1 s and 0.25 s, in line with the values obtained here. The comparison of the values of Table 4 and Table 5 with the ones proposed in EC8 is not straightforward, since the standard describes the hazard function through a first-order power law approximation (zero curvature, k2=0) and as a function of the PGA. A value of k1=2.5 is prescribed in the standard for the entire Portugal mainland (onshore seismic action). Naturally, the k1values proposed in Table 4 and Table 5 are higher since they depend on the values of Sa(T1). 7. Mean annual frequency of limit-state exceedance The limit state of our interest to be verified corresponds to the structural capacity Sacomputed from the nonlinear static pushover analysis, performed on the buildings database (further details in Bernardo et al [54]). The output is summarized in the lognormal probability plot depicted in Fig. 10 for different number of stories, showing that the data set follows a lognormal distribution hypothesis with the corresponding statistic moments. The median capacity Saclearly decreases with the number of stories. The capacity dispersion values result in coefficient of variation (cov) between 0.20 and 0.28. These results can be Fig. 5. Structural capacity/response. Left: ECDF (top) and box-whisker plot of T1 (bottom); center and right: S a(T1)in the MAL region for different RP and cutsection for 475-RP –bedrock (top) and surface (bottom). Fig. 6. Hazard 2nd order power-law fit –First row: pivot points in S a(T1=0.24s) for RP (years) =[100, 475, 2000]; Second row: batch fitting (RP =475 years) in several parishes and for buildings with 1, 3 and 5 stories, respectively. V. Bernardo et al. Structural Safety 111 (2024) 102514 7 used in Eq. (12) and Eq. (14) to compute the mean annual frequency of limit state exceedance λLS. However, the λLS estimation in the present study treated the capacity independently for each building of the synthetic database. Additionally, epistemic uncertainty with cov equal to 0.25 was adopted in the Saof each structure. This value is recommended by FEMA P-58 [56] for medium quality and completeness of the analytical models that describe the structural capacity. Based on the previous assumptions and employing the hazard results of Section 6.1, the reliability index βwas computed using the inverse standard normal cumulative density function β= − Φ−1(λLS).Fig. 11 depicts the empirical cumulative density function (ECDF) of βvalues for different buildings located in different regions according to the MAL hazard distribution. The results were computed for the 16th 50th and 84th hazard quantiles for 475-years RP –grey curves. The mean ECDF shows the βdistribution for the mean hazard in the MAL region. The vertical dashed-line of β=2.5 was defined as the requirements for the life safety verification (RP =475 yrs.) of existing structures [57]. As shown, the probability of the building stock having a β<2.5 (at bedrock) is around 0.17, 0.38, 0.50, 0.64 and 0.78, for buildings with 1 to 5 stories, respectively. These probabilities increase from 0.20 to 0.90 Fig. 7. Box-diagram of regression coefficients (RP =475 years; bedrock) –first row: building stock aggregated; second row: disaggregation by AML parishes (threestories buildings). Fig. 8. Box-diagram of regression coefficients (RP =475 years; surface) –first row: building stock aggregated; second row: disaggregation by AML parishes (threestories buildings). V. Bernardo et al. Structural Safety 111 (2024) 102514 8 by considering the soil specific conditions at surface. To provide information regarding the overall βdistribution in the MAL region and prioritize the parishes for seismic mitigation strategies, the weighted averages βwwere obtained based on the number-of-stories βcoefficients and buildings floors areas in each region according to Census data. Fig. 12 shows the βwmap for the median hazard and highlights the parishes with βwvalues lower than 2.5 in red, considering the hazard results at bedrock and surface. The results reflect the hazard distribution in the MAL and the soil amplification effect at surface in regions with moderate to high hazard, e.g., downtown Lisbon and the surrounding areas, west and south coastal areas, Tagus Valley (NE) and Setubal district (SE) regions. The estimated number of existing pre-code masonry buildings with βw<2.5 represents around 80 % of the total building stock. The disaggregation of βwin the priority regions (βw<2.5) was performed to evaluate the most influencing buildings typologies, i.e., number-of-floors cluster, by normalizing with respect to the βwand total buildings floors areas in each region. The results are presented in the ECDF depicted in Fig. 12. Disaggregation based on archetypes was not considered due to the similar structural responses across different geometry layouts, reflecting lower discrepancies in the ratio of masonry walls area to floor area in plan, as discussed in Bernardo et al.[32]. Looking at the median results at bedrock in Fig. 12, can be concluded that buildings with 4-, 5-, 3-stories high, respectively, have a greater influence on the βwvalues. This observation is related to the lower β coefficients (lower seismic capacity) associated with these buildings’ typologies, and their high density in these specific regions. On the other hand, the median results at surface are slightly different since the soil effect seems to amplify the response of 1-, 2-, 3-stories high buildings and attenuate the remaining typologies. Structures with 3 to 5 stories represent the most vulnerable typologies in that case. Naturally, this overall disaggregation needs to be stratified independently by parish to effectively target mitigation measures at the most vulnerable typologies of buildings in those regions. 8. Final comments and conclusions The present paper focused on the seismic reliability assessment of the pre-code masonry building stock in the Metropolitan Area of Lisbon (MAL) region, employing structural numerical models of representative buildings and the local seismic hazard. The main purpose of the study is i) to provide new second-order hazard solutions with structural dependency for the mean annual frequency of limit-state exceedance, with the aim to replace the first-order hazard solutions proposed by the current EC8; ii) provided valuable information for the development of national application documents for Eurocode 8 supported through a reliability-based code calibration, to ensure an appropriate target safety level. The results can also be useful to identify the most vulnerable regions in the MAL and support further risk studies. This study represents the first comprehensive investigation of its kind in the region based on reliability analysis and adopting a probabilistic framework. The first part of the study conducts a probabilistic seismic hazard analysis in the region under investigation by combining three different seismogenic models. Considering the range of fundamental periods in the building database, hazard functions at bedrock and at surface were then derived in terms of first-mode spectral acceleration. The Table 4 Summary of the mean regression coefficients (k 0 , k 1 and k 2 ) in the MAL region (bedrock) by number of stories and different seismic hazard quantiles (16th, 50th and 84th). No. of stories 1 2 3 4 5 16th 50th 84th 16th 50th 84th 16th 50th 84th 16th 50th 84th 16th 50th 84th k 0 (10 ¡5 )15.0 3.1 7.9 5.1 4.6 22.8 4.2 3.9 29.7 5.0 3.8 51.9 8.7 3.6 10.2 k 1 1.74 3.76 3.52 3.09 3.81 3.11 3.21 3.80 2.81 2.98 3.75 2.02 2.30 3.56 2.90 k 2 0.02 0.42 0.36 0.27 0.46 0.31 0.30 0.45 0.26 0.25 0.43 0.06 0.13 0.37 0.25 k 0 (10 ¡5 )0.2 1.3 11.1 0.8 2.8 19.8 1.1 5.3 27.0 1.0 5.5 26.4 0.7 4.0 22.4 k 1 4.46 4.43 3.37 4.19 4.27 3.21 3.99 3.53 2.84 3.98 3.30 2.75 4.00 3.15 2.63 k 2 0.47 0.57 0.34 0.46 0.58 0.34 0.43 0.35 0.28 0.42 0.32 0.26 0.41 0.28 0.23 k 0 (10 ¡5 )0.2 1.2 10.5 0.7 3.0 18.6 1.1 4.9 26.3 0.9 4.8 25.0 0.5 3.9 21.0 k 1 4.62 4.39 3.38 4.35 4.13 3.25 4.10 3.66 2.88 4.09 3.54 2.82 4.18 3.39 2.71 k 2 0.49 0.55 0.35 0.49 0.53 0.36 0.45 0.41 0.29 0.44 0.37 0.29 0.44 0.33 0.26 Table 5 Summary of the mean regression coefficients (k 0 , k 1 and k 2 ) in the MAL region (surface) by number of stories and different seismic hazard quantiles (16th, 50th and 84th). No. of stories 1 2 3 4 5 16th 50th 84th 16th 50th 84th 16th 50th 84th 16th 50th 84th 16th 50th 84th k 0 (10 ¡5 )7.2 5.3 23.0 4.3 13.8 40.4 4.1 7.6 50.1 4.2 6.3 75.5 5.9 5.4 13.3 k 1 2.83 3.62 3.18 3.49 3.39 3.01 3.44 3.61 2.72 3.36 3.63 1.98 2.91 3.47 2.84 k 2 0.22 0.42 0.32 0.37 0.39 0.31 0.36 0.44 0.26 0.34 0.43 0.06 0.25 0.38 0.26 k 0 (10 ¡5 )0.6 2.8 20.2 1.7 5.9 35.8 2.3 9.5 46.0 2.0 9.0 44.1 1.3 7.2 36.6 k 1 4.33 4.16 3.27 4.08 4.07 3.10 3.89 3.56 2.75 3.89 3.46 2.67 3.90 3.32 2.56 k 2 0.47 0.53 0.34 0.47 0.56 0.34 0.44 0.41 0.28 0.43 0.38 0.26 0.41 0.33 0.23 k 0 (10 ¡5 )0.5 2.6 19.3 1.6 5.1 34.1 2.2 9.8 45.1 1.9 10.3 42.6 1.2 7.5 34.9 k 1 4.49 4.22 3.27 4.23 4.06 3.13 3.98 3.36 2.78 3.98 3.16 2.73 4.07 3.04 2.63 k 2 0.50 0.57 0.35 0.50 0.58 0.36 0.45 0.36 0.29 0.45 0.31 0.29 0.45 0.28 0.26 V. Bernardo et al. Structural Safety 111 (2024) 102514 9