The Sensitivity of Global Structural Parameters for Unreinforced Masonry Buildings Subjected to Simulated Ground Motions

Abstract

This research performs a parametric study based on Equivalent Single Degree of Freedom (ESDOF) models for simplified seismic analysis of unreinforced masonry (URM) structures. This is a necessary action due to the fact that it is not affordable to model and analyze populations of masonry buildings by using detailed continuum-based models during regional seismic damage and loss estimation studies. Hence, this study focuses on the sensitivity of major structural parameters of a selected idealized hysteretic model for URM buildings. The numerical models are subjected to region-specific simulated ground motion time histories generated using validated seismological parameters. The variations in dynamic analysis results are evaluated using statistical tools for major structural and seismological parameters. The results reveal that the strength factor is the most influential structural parameter, whereas magnitude and distance have a significant impact on the response of idealized URM models as seismological parameters. Furthermore, the specific seismic performance exhibiting limited ductility capacity and the narrow margin of safety between the initial state of inelastic behavior and the ultimate (collapse) state for URM buildings is verified by the statistical approaches employed in this study.

1. Introduction

One of the oldest types of construction is unreinforced masonry (URM). This construction type still constitutes a significant percent of the building stock in Asia and South America as well as in some parts of Europe and North America. By definition, URM structures have load-bearing walls made of different materials like brick, concrete block, stone and adobe, which have been combined by mortar in between. In general, this construction type is known to have poor and brittle behavior under seismic motion due to the low tensile strength capacity of the walls; although the seismic performances of masonry structures vary widely all around the world according to the local construction practices, material properties and workmanship. In addition, masonry buildings are rigid and heavy structures; hence, the load-bearing masonry walls attract higher levels of earthquake loads when compared to more flexible frame types of construction systems. This makes URM buildings highly vulnerable to seismic action, and there are numerous studies that have aimed to identify the seismic risk of this construction type in order to mitigate the consequences of this vulnerability [1,2,3,4,5,6,7]. However, there is a challenge in quantifying the seismic risk of URM type, especially in the case of regional damage and loss estimation studies that include populations of masonry buildings with different structural characteristics. A reliable seismic response can be obtained when continuum or finite mesh-based modeling is employed for masonry buildings. However, it is computationally not feasible to model and analyze numerous masonry buildings with complex approaches, even if they are represented in terms of predefined building sub-classes. Then, simplified but still reliable alternatives should be employed to achieve this task. The Equivalent Single Degree of Freedom (ESDOF) approach is the most commonly used tool for simplified seismic analysis of building structures. In using ESDOF models, the most critical issue is to attain major structural parameters that dominate the global behavior of the structure under concern. That is the only way to obtain reliable response statistics for a large population of URM buildings during regional seismic risk analysis.

A key step of regional seismic risk analysis is the definition of a representative set of ground motion records. Recently simulated ground motions are in use as alternatives to real records. A literature survey reveals the existence of alternative approaches to simulate region-specific ground motion time histories. Among these, deterministic approaches are more accurate but require detailed and well-defined input-model data, while stochastic ones are straightforward and more practical, providing reliable broadband ground motion records [8,9,10]. The simulated ground motion datasets provide a homogeneous distribution in terms of magnitude, source-to-site distances and soil information compared to real datasets. Thus, they are useful in the assessment of the seismic performance of structures under a wide set of regional input ground motions [11,12,13,14,15].

The primary objective of this study is to investigate the seismic performance of URM structures when subjected to region-specific simulated ground motion time histories generated using validated seismological parameters. For this purpose, in this study, a simulated ground motion dataset in Duzce (Turkiye) [11,12,15,16,17,18,19] generated by the stochastic finite-fault ground motion simulation approach of Motazedian and Atkinson [8] is employed. Next, through statistical methods, the aim is to analyze the effects of various seismological and structural characteristics on the performance levels of URM structures. To simplify the nonlinear structural analyses, multi-degree of freedom (MDOF) masonry structures are transferred into ESDOF masonry models. The limit states in terms of displacement of the idealized SDOF models are assessed through evaluation of the fundamental period, strength ratio, ductility, post-yielding ratio, post-capping ratio, residual strength, and degradation parameter. Next, sensitivity analyses are performed to observe changes in damage states (DSs) with respect to the earthquake magnitude, source-to-site distance and soil conditions for varying performance levels of the URM models. Finally, to evaluate relationships between DSs with seismological and structural parameters and to distinguish the most influential parameters on the estimation of DSs, alternative statistical tools, including discriminant analysis and multiple regression analysis, are employed. In this study, the authors make use of data-driven approaches, including discriminant analysis and multiple regression analysis, to evaluate relationships between DSs with seismological and structural parameters and to distinguish the most influential parameters on the estimation of DSs. Such approaches are very popular in civil and earthquake engineering research, and they have been used by many researchers [20,21,22,23,24,25,26]. Moreover, the recent advances in uncertainties regarding the sensitivity analysis of global structural parameters (especially for masonry structures) can be accessed in [27,28,29,30,31,32].

2. Simplified Modeling Approach for URM Buildings

Due to the complexity of modern structures in both elevation and plan, as well as the diversity of architectural considerations, the demand for structural modeling shifts toward complex and detailed solutions. Easy access to large-size and high-speed computing with current advances in computer technology also encourages the use of huge and detailed numerical models in structural earthquake engineering. However, detailed modeling and complex analysis always take significant time and effort. Moreover, it requires expertise to define numerous structural parameters in detailed models appropriately and also to interpret the huge amount of response data after the analysis. The situation becomes even more tangled if a large population of buildings should be modeled rather than an individual structure. Hence, a possible solution for this issue is to use simplified numerical models. There are many studies in the literature that have proposed simplified modeling approaches. Among these, the simplest and the most common approach is using an ESDOF model. Starting with the pioneering work of Biggs [33] and later with milestone studies like Saiidi and Sozen [34] and Qi and Moehle [35], many researchers employed this technique. It was also implemented in well-known and widely used international standards ATC-40 [36] and FEMA-273 [37] to promote nonlinear static procedures in seismic analysis of structures. This simplified model is basically obtained from the capacity curve of a complex structure after pushover analysis. The trade-off between simplicity and accuracy is based on two gross assumptions for the ESDOF approach. According to the first one, the real structure should be regular enough to ensure a dominant fundamental (first) mode that represents total response in a satisfactory manner. Second, the deformed shape of the structure is assumed to remain constant during vibrational motion. These assumptions reveal that the ESDOF idealization yields reliable seismic response values for regular and first-mode dominant structures. However, it should be kept in mind that except for midrise, large, and complex apartment buildings in urban regions, the majority of the masonry buildings, especially the rural ones, are low-rise, small-sized and nearly regular plan buildings, which may fit within the limitations of the approach. Hence, seismic response statistics obtained from the ESDOF approach are deemed to provide reasonable predictions of the maximum displacement response of actual MDOF masonry buildings.

The ESDOF modeling approach requires using a hysteretic model with global structural parameters to obtain an inelastic dynamic response under seismic excitation. There are many different hysteretic models, from the simplest to the most complex. The important issue is to select an optimum model that has sufficient model parameters to simulate the global seismic response of the structure in a satisfactory manner. To this extent, one of the most suitable models is the one that has been proposed by Ibarra et al. [38], called the Ibarra–Medina–Krawinkler model. This is a vertex-oriented hysteretic model which can simulate all types of degradation (i.e., stiffness, strength and pinching), which was enhanced by Lignos [39] with a modified envelope curve and renewed modes of deterioration. This hysteretic model was used in the OpenSees seismic analysis platform in terms of the moment versus chord rotation relationship. In order to use the hysteretic model in the ESDOF approach, these parameters have been converted to force and displacement, respectively, by using their basic structural definitions.

Figure 1 shows the envelope curve for the modified Ibarra–Medina–Krawinkler (MIMK) hysteretic model for symmetric response. There are three major strength parameters. These are effective yield strength (f_y), capping strength (f_c) and residual strength (f_r). Up to yield strength, the structure is assumed to behave in a linear elastic manner with elastic stiffness (k_e). The region between the yield and the capping strength with stiffness kp represents the post-yield inelastic behavior of the structure. The negative stiffness region after the capping strength (i.e., k_c) denotes an unstable response up to a highly reduced capacity (i.e., residual strength), which stands for the initiation of total structural collapse. There are also the corresponding displacement parameters in the model that can be listed as yield displacement (δ_y), capping displacement (δ_c), residual displacement (δ_r) and ultimate displacement (δ_u). For the post-yield and post-capping slopes, non-dimensional model parameters α_py and α_pc can be defined as follows:

$${\alpha }_{py}=k_p/k_e$$

(1)

$${\alpha }_{pc}=k_c/k_e$$

(2)

For the simulated cyclic behavior, the MIMK model has a cyclic degradation parameter (γ) that represents the amount of hysteretic energy dissipation through cycles with degrading behavior. The rate of strength degradation is taken constant as the default value in OpenSees (i.e., unity).

Considering the above discussions, seven parameters have been selected to define the global seismic response of the ESDOF model subjected to earthquake ground motions. These parameters are period (T), strength ratio (η), ductility (μ) and the aforementioned hysteretic model parameters α_py, α_pc, λ and γ. The first three structural parameters can be easily obtained through Equations (3)–(5), where m is the mass and g is the gravitational acceleration.

$$T=2\pi \sqrt{m/k_e}$$

(3)

$$\eta =f_y/mg$$

(4)

$$\mu ={\delta }_c/{\delta }_y$$

(5)

In the next phase, discrete parameters are attained for each of these structural and model parameters to classify URM building models according to their vibrational properties and structural conditions. The assigned values for each sub-class are given in

Table 1. Discrete values of structural and model parameters assigned for URM sub-classes.

Vibrational Properties	Very Rigid	Rigid	Average	Flexible	Very Flexible
Period (T)	0.05 (T1)	0.1 (T2)	0.2 (T3)	0.3 (T4)	0.4 (T5)
Seismic Conditions	Very poor	Poor	Moderate	Good	Very Good
Strength Ratio (ɳ)	0.1 (F1)	0.3 (F2)	0.5 (F3)	0.7 (F4)	0.9 (F5)
Ductility Ratio (µ)	2.0 (M1)	2.5 (M2)	3.0 (M3)	3.5 (M4)	4.0 (M5)
Post-yield Ratio (αs)		0 (Y2)	0.05 (Y3)	0.1 (Y4)
Post-capping Ratio (αr)		−0.4 (C2)	−0.3 (C3)	−0.2 (C4)
Residual Strength Ratio (λ)		0 (R2)	0.2 (R3)	0.4 (R4)
Model degradation parameter (γ)		200 (D2)	400 (D3)	800 (D4)

. The abbreviations of sub-classes for each parameter are also provided in parenthesis in the table. These abbreviations are later used to demonstrate the seismic response statistics of the parametric dynamic analyses.

In relation to vibrational properties, there are five sub-classes in terms of the fundamental period, i.e., very rigid to very flexible structures. In order to obtain these values, the empirical period formulations in different international standards and the studies based on ambient vibration tests of actual masonry buildings are considered. Starting from Uniform Building Code [40], the American Design guidelines employ an empirical period formulation in the form of a power function:

$$T=C_t{H}^x$$

(6)

where C_t and x are the constants that were determined for each construction type by using the available field data and based on expert opinion. In the final version of the ASCE standard series, i.e., ASCE/SEI 7-16 [41], the values for these constants have been provided as 0.049 and 0.75, respectively. Eurocode 8 [42] has a similar formulation for URM buildings.

$$T=0.05H^{0.75}$$

(7)

There are also previous studies that propose empirical period formulations for Turkish masonry buildings. In an early attempt, an international technical team performed field surveys on rural masonry buildings after major earthquakes in Turkiye in the 1970s and 1980s [43]. Their final report declared that Equation (8) could be used as a rule of thumb for rural masonry buildings in Turkiye, based on field measurements. In this equation, N is the number of stories.

$$T=0.05-0.06N$$

(8)

In another study that focused on masonry buildings in the Marmara region, the following empirical relationships were obtained for masonry buildings with timber floors (Equation (9)) and those with reinforced concrete floors (Equation (10)) [44].

$$T=0.039H$$

(9)

$$T=0.062H^{0.87}$$

(10)

Considering all these formulations and noting that most of the masonry building stock (about 90%) is composed of 1- or 2-story structures, whereas the remaining ones (especially in urban regions) have 3–6 stories, it would be representative to consider 0.05 s and 0.4 s for lower and upper bounds, respectively, where a period of 0.2 s could be a good candidate for the average value.

The second classification is for the strength ratio (η) parameter. As seen in Table 1, the parameter η takes five discrete values between 0.1 and 0.9. These values have been selected in accordance with some previous experimental work. Two reduced-scale URM buildings made of clay brick units were tested by Costley and Abrams [45]. The experimental capacity curves for the specimens in terms of base shear coefficient/moment versus drift are shown in Figure 2. During the tests, the base shear coefficient of the specimen with a smaller number of openings (left plot) was obtained as 1.2, whereas the other specimen had a η value of 0.7 (right plot). Benedetti et al. [46] tested 24 half-scale specimen buildings that represented non-engineered URM structures with stone and brick walls. He obtained the base shear coefficient, equivalent to parameter η, in the range of 0.1–0.3 from the tested specimens. Tomazevic et al. [47] conducted an experimental campaign by testing 1/5 scale multi-story URM building specimens constructed with different types of materials. According to the test results, parameter η (BSC in Figure 3) ranged from 0.5 to 1.9.

It should be mentioned that the above η values from different experimental programs were obtained in laboratory conditions, and they are generally higher than the strength ratio of actual URM buildings due to the differences in construction practice, material quality and workmanship. Verification of this statement can be achieved by referring to a previous study in which a real and full-scale URM building was tested to near collapse [48]. It was a typical two-story brick masonry building with some structural deficiencies (Figure 4a). These types of masonry buildings constitute a significant part of the Turkish building stock. The building was subjected to cyclic lateral loads up to near collapse. The corresponding cyclic base shear force versus top displacement plot is given in Figure 4b. It can be clearly seen that the parameter η is about 0.56 from the plot, lower than the expected value before the test. Hence, it can be concluded that all the aforementioned research has been taken into account during the selection of discrete values for the parameter η.

The next parameter is the ductility ratio (μ), which takes five discrete values between 2 and 4 in Table 1. The selection of these values for the parameter μ is again carried out thanks to previous experimental data. For instance, Zavala et al. [49] conducted a field static load test on a full-scale two-story brick masonry building. From the force–displacement relationship of the tested building, μ can be calculated as 2.8. Another experimental study was performed by Magenes and Penna [50], who tested three half-scale stone masonry building specimens on a shaking table. The parameter μ was approximately 5.0 from the cyclic force–displacement relationships obtained after the tests for the specimen with diagonal cracking failure mode. Lourenço et al. [51] conducted tests on half-scale concrete masonry building specimens. From the base shear coefficient versus drift curves of the URM specimens, the parameter μ can be obtained approximately as 5.5. Apart from the experimental data, field observations after major earthquakes have revealed that masonry structures have limited ductility up to collapse. When the elastic limit is exceeded by cracking or damage in masonry walls, masonry building rapidly reaches the ultimate state before it can exhibit considerable inelastic action. This is especially true if out-of-plane mechanisms dominate the global response. Considering this fact, the selected values are more conservative than the experimentally obtained ones.

For the remaining hysteretic model parameters, all the discrete values in Table 1 have been selected in accordance with the past studies from the literature in which the MIMK model had been employed. For each parameter, three discrete values are assigned for sub-classes of poor, moderate and good conditions. Overall, there are 2025 combinations of input parameter sets for the hysteretic model in terms of seismic conditions. When the sub-classes for vibrational properties are also included, it makes a total of 10,125 structural combinations for URM models to cover all range of inherent characteristics for this specific construction type.

3. Simulated Ground Motion Dataset

In this study, the simulated ground motion dataset of Karimzadeh et al. [11] is employed to evaluate a homogeneous variation of seismological parameters that are unavailable in real ground motion datasets. The dataset was generated for the Duzce region in Northwest Turkiye using the stochastic finite-fault ground motion simulation approach of Motazedian and Atkinson [8].

The stochastic finite-fault ground motion simulation technique based on a dynamic corner frequency method as implemented in EXSIM software models a series of sub-faults that comprise a rectangular fault domain [8]. Each sub-fault is modeled as a point source with a ω-2 spectrum [52]. The fundamental assumption is that the hypocenter from which rupture propagates radially is located on one of the sub-faults, and subsequently, the rupture travels to the centers of the remaining sub-faults. Finally, the contributions from all sub-faults with specific lengths and widths are summed in the time domain while accounting for their respective time delays. The dynamic corner frequency approach conserves the total radiated energy from the ruptured plane regardless of the sub-fault dimension.

To depict near-field and far-field motions, the dataset was simulated in two alternate zones.

The far-field region has 66 dummy stations, while the near-field region contains 280 dummy stations closer. A total of 344 fictitious seismic stations were investigated between 30.0 and 32.0 latitudes in the east and 40.0 and 42.0 longitudes in the north. It is noted that the simulated records are baseline corrected and bandpass filtered between 0.1 to 25 Hz by using the fourth-order Butterworth filtering type.

Simulations were performed for hypothetical events with moment magnitudes of Mw = 5.0, 5.5, 6.0, 6.5, 7.0, and the Duzce 1999 event (M = 7.1) using two local site conditions classified as soft and hard soils. In this study, for nonlinear time history analysis (NLTHA) of ESDOF models, 240 time history earthquake records are selected from the database by using a pseudo-random algorithm. Half of these records belong to soft soil type while the other half corresponds to hard soil condition. During the selection process, 20 records per scenario earthquake are selected. PGA is employed as the major ground motion intensity parameter because it correlates well with the seismic response of rigid masonry structures. The selected records cover a PGA range of 0.1 g to 1.0 g. The PGA values are evenly distributed throughout ten intervals to capture the entire seismic behavior.

4. Seismic Response Analysis for Masonry Buildings

This section focuses on the results obtained from the seismic response analysis of ESDOF models that represent masonry buildings subjected to site-specific synthetic ground motion records. As discussed in the previous sections, three structural parameters (T, η, μ) and four hysteretic model parameters (α_s, α_r, λ, γ) are employed to predict the global seismic response of URM structures. Each parameter has a specific number of discrete values with an average (typical) one and deviations from this average value for both favorable and unfavorable conditions, as seen in Table 1. When all the combinations of parametric discrete values are calculated, it adds up to 10,125 cases. Furthermore, there are three seismological parameters (M, RJB distance, and soil condition). Since 20 different earthquake stations (with different RJB distances) have been selected for the synthetic ground motions with six magnitude values (5.0, 5.5, 6.0, 6.5, 7.0, 7.1) and two soil conditions (hard site, soft site), the total number of synthetic ground motion records become 240. Hence, the number of the NLTHA for the generated ESDOF simulations (10,125) subjected to synthetic earthquake records (240) sum up to 2,430,000. This means a huge response database is obtained after a considerable number of dynamic analyses. This section examines the obtained seismic response data in three parts. First, the seismic response of simulations with the most unfavorable and favorable structural combinations, together with the typical (average) structural combination, are compared in terms of seismological parameters. In the second part, the effect of each structural and model parameter is assessed in an individual manner while the rest of the parameters are kept constant. Finally, damage state probabilities for simplified masonry models are obtained in terms of PGA and magnitude for both soil conditions separately.

In order to classify and discretize seismic response data, damage states (DS) are employed in this study. Damage states are specific structural performance states bounded by performance thresholds called limit states (LS). The exceedance or non-exceedance status of limit states determines the damage state that the structure experiences under the given loading intensity. In this study, three limit states (i.e., four damage states) have been defined in accordance with the employed hysteretic model in this study (Figure 5). Up to LS-1, the structure mostly behaves in the elastic range; therefore, this damage state (i.e., DS-1) can be regarded as “none/slight damage”. Exceedance of LS-1 indicates that the structure is in the inelastic range up to displacement at capping strength, i.e., LS-2. This damage state (DS-2) can be named “moderate damage”. After the exceedance of displacement at capping strength (LS-2), the capacity of the structure is significantly reduced, and damage evolves very rapidly up to displacement at residual strength (LS-3). This region of seismic behavior can be regarded as “extensive/heavy damage” (DS-3). Finally, the exceedance of displacement at residual capacity (LS-3) indicates the “collapse” of the structure (i.e., DS-4). After each seismic response analysis, the maximum displacement attained from the corresponding ESDOF model subjected to a specific synthetic ground motion record is compared with the limit state values to determine the damage state of the model with a certain parametric set of values. This means that for each dynamic analysis, there is a discrete performance level result in terms of DS. All the discussions related to seismic response analysis are carried out based on the final DS values.

All the seismic response analyses of ESDOF models were carried out using the earthquake analysis platform OpenSees. A simple inverted pendulum model was developed in which MIMK hysteretic model is employed for the inelastic cyclic behavior of the rotational hinge. Dynamic analyses are performed by using the Newmark numerical integration method.

In the first phase of seismic response evaluation after parametric analysis, the overall DS values of the ESDOF models with the most unfavorable, average (typical) and favorable combination of parametric values are compared under soft and hard soil conditions as shown in

Table 2. Overall DS value matrices for M-RJB pairs for (a) the most unfavorable combination in soft soil condition, (b) the most favorable combination in soft soil condition, (c) the most unfavorable condition in hard soil condition and (d) the typical (average valued) combination in soft soil condition.

(a)	Distance interval (km)					(b)	Distance interval (km)
M	0-5	5–10	10–15	15–20	20–25	M	0-5	5–10	10–15	15–20	20–25
7.1	DS-4	DS-4	DS-4	DS-4	N.A	7.1	DS-2	DS-2	DS-1	DS-1	DS-1
	(4.00)	(3.33)	(4.00)	(4.00)			(1.78)	(1.33)	(1.00)	(1.00)	(1.00)
7	DS-4	DS-4	DS-4	N.A	DS-3	7	DS-2	DS-1	DS-1	N.A	DS-1
	(4.00)	(4.00)	(3.33)		(3.00)		(1.33)	(1.00)	(1.00)		(1.00)
6.5	DS-4	DS-4	DS-4	N.A	DS-3	6.5	DS-1	DS-1	DS-1	N.A	DS-1
	(4.00)	(4.00)	(3.50)		(2.50)		(1.00)	(1.00)	(1.00)		(1.00)
6	DS-4	DS-4	DS-3	DS-2	DS-2	6	DS-2	DS-1	DS-1	DS-1	DS-1
	(3.20)	(3.25)	(2.71)	(1.50)	(1.50)		(1.20)	(1.00)	(1.00)	(1.00)	(1.00)
5.5	DS-3	DS-2	DS-1	DS-1	DS-1	5.5	DS-1	DS-1	DS-1	DS-1	DS-1
	(2.80)	(1.50)	(1.00)	(1.00)	(1.00)		(1.00)	(1.00)	(1.00)	(1.00)	(1.00)
5	DS-1	DS-1	DS-1	DS-1	DS-1	5	DS-1	DS-1	DS-1	DS-1	DS-1
	(1.00)	(1.00)	(1.00)	(1.00)	(1.00)		(1.00)	(1.00)	(1.00)	(1.00)	(1.00)
(c)	Distance interval (km)					(d)	Distance interval (km)
M	0-5	5–10	10–15	15–20	20–25	M	0-5	5–10	10–15	15–20	20–25
7.1	DS-4	DS-3	DS-4	DS-2	N.A	7.1	DS-3	DS-2	DS-2	DS-1	DS-1
	(3.85)	(2.66)	(4.00)	(2.00)			(2.85)	(2.00)	(2.00)	(1.00)	(1.00)
7	DS-4	DS-3	DS-4	N.A	DS-3	7	DS-3	DS-2	DS-1	N.A	DS-1
	(3.58)	(2.66)	(4.00)		(3.00)		(2.50)	(2.00)	(1.00)		(1.00)
6.5	DS-4	DS-4	DS-3	N.A	DS-2	6.5	DS-2	DS-2	DS-1	N.A	DS-1
	(3.50)	(4.00)	(2.50)		(1.50)		(1.90)	(1.50)	(1.00)		(1.00)
6	DS-4	DS-4	DS-3	DS-1	DS-1	6	DS-2	DS-1	DS-1	DS-1	DS-1
	(3.40)	(3.25)	(2.43)	(1.00)	(1.00)		(1.60)	(1.00)	(1.00)	(1.00)	(1.00)
5.5	DS-3	DS-2	DS-1	DS-1	DS-1	5.5	DS-1	DS-1	DS-1	DS-1	DS-1
	(2.80)	(1.50)	(1.00)	(1.00)	(1.00)		(1.00)	(1.00)	(1.00)	(1.00)	(1.00)
5	DS-1	DS-1	DS-1	DS-1	DS-1	5	DS-1	DS-1	DS-1	DS-1	DS-1
	(1.00)	(1.00)	(1.00)	(1.00)	(1.00)		(1.00)	(1.00)	(1.00)	(1.00)	(1.00)

. Overall DS values in each cell of the matrices in Table 2 are obtained by checking the exceedance/non-exceedance of limit states for that specific combination of magnitude and distance interval (the R_JB distance intervals are classified as 0–5, 5–10, 10–15, 15–20 and 20–25 in km). This means a discrete value is assigned to each damage state (i.e., 1 for DS-1, 2 for DS-2, 3 for DS-3 and 4 for DS-4), and then the average DS values for each M-R_JB pair is calculated, which is the number in parenthesis in each cell. For instance, a value of 1.0 means that all the related seismic response analyses have yielded displacement values not exceeding the yield displacement (i.e., all under elastic response), so the overall performance of the simulation is DS-1. However, when the value is 1.5, this means some of the analyses have yielded displacement values exceeding the yield displacement, so the overall performance of the simulation is regarded as DS-2. Accordingly, the matrices show the most unfavorable and favorable combinations of structural parameters as well as the typical (average valued) one for soft soil conditions (see Table 2a,b,d). In this manner, within this huge response dataset, the seismic response values at the parametric boundaries can be compared with each other and with the central values of parameters. The cells with the abbreviation “N.A.” mean there is no response data for that specific combination of M-R_JB pairs.

Referring to Table 2a,c, the most unfavorable combination of parameters is obtained as the structural simulation coded as T1-F1-M1-Y2-C2-R2-D2. This combination of parameters represents a very rigid masonry structure with poor strength and very limited ductility. This generic structural typology is also assumed to exhibit insufficient hysteretic energy dissipation with severe deterioration under earthquake loading. There are many examples of such non-engineered and deficient masonry buildings in the actual building stock, especially in rural regions. On the other hand, Table 2b gives the results for the most favorable combination of parameters, i.e., T5-F5-M5-Y4-C4-R4-D4. This specific combination of parameters stands for multi-story (flexible) masonry structures with high strength and ductility capacity, which can exhibit a considerable amount of hysteretic energy with minor deterioration during seismic action. This combination of parameters symbolizes well-designed, modern, urban-type masonry apartments. Finally, Table 2d yields overall DS values for the typical (average valued) structural parameter combination, abbreviated as T3-F3-M3-Y3-C3-R3-D3. This set of parameters represents common masonry construction with 2–3 stories in small cities with some deficiencies due to the used material, workmanship, or the construction technique, but not all. In all the matrices, variation with magnitude is given in column-wise order, whereas the variation with distance is provided in row-wise order. This means the uppermost cell has the most critical value in terms of seismic hazard (i.e., largest magnitude and closest distance). Each DS is given a different color code for better visualization of the results.

The first comparison can be carried out between unfavorable and favorable combinations of structural parameters under soft soil conditions (Table 2a,b). There seems to be a significant difference in seismic response, especially for M > 6 and R_JB < 15 km. In this range of seismological parameters, nearly all the ESDOF models with unfavorable combinations of parameters experience collapse (i.e., DS-4), whereas the ones with a favorable combination of parameters possess moderate damage (i.e., DS-2) at most. These response statistics are in accordance with the field observations after major earthquakes in Turkiye, such that rural and non-engineered masonry buildings are imposed to severe damage and collapse during even moderate earthquakes [53], whereas well-constructed and robust masonry buildings can stand still even after a major earthquake while the nearby reinforced concrete buildings have been razed to the ground [54]. The results also provide some insight into seismic risk evaluation of masonry buildings for earthquake loss mitigation strategies since it is revealed that making strengthening/demolition decisions for non-engineered masonry buildings close to active fault lines that have the potential to trigger large-magnitude earthquakes is a priority. It should also be added that for smaller magnitudes (M < 5.5) and far distances (R_JB > 10 km), the structural differences do not matter.

The second comparison from Table 2 is for the effect of local soil conditions for which the seismic response of ESDOF models with unfavorable combinations of parameters are compared for soft and hard sites (Table 2a,b). It is observed that under soft soil conditions, the DS values are slightly more critical, especially for cases of higher magnitude and closer distance. Just like in the previous case, as the magnitude becomes lower and the distance becomes longer, the effect of local soil conditions diminishes.

Finally, the response statistics obtained from the ESDOF models with typical (average) values of structural parameters (Table 2d) are compared with the two limiting cases in Table 2a,b for soft soil conditions. The results show that the central values of parameters yield seismic responses between the bounds as expected but much closer to the results obtained from the models with the most favorable combination without any case of collapse. Especially for M ≥ 6, the results are almost the same for all distance intervals. This may be attributed to the fact that masonry structures with an average strength (η = 0.5) and ductility (μ = 3) have the ability to withstand even major earthquakes without collapse, and they mainly satisfy the life safety performance criterion.

In the second phase of seismic response evaluation after parametric analysis, the sensitivity of each structural (T, η, μ) and model (α_py, α_pc, λ and γ) parameter with varying seismological (M, R_JB distance, site condition) parameters is assessed in an individual manner while the rest of the parameters are kept constant. In this paper, the sensitivity analyses are discussed in a detailed manner only for global structural parameters since the change in model parameters was observed not to be sensitive to DS for all M values and R_JB distance intervals. DS values have a slight tendency to decrease with model parameters by shifting from soft to hard site conditions. For more details regarding the sensitivity of model parameters, one can refer to Koc [55].

The first parameter to be investigated is T, with five discrete values, as mentioned in Table 1. The other structural and model parameters are kept constant at their average (central) values. The seismic response values in terms of DS that have been obtained from the dynamic analyses form a grid network with discrete T values (shown as blue dots in the figure), as seen in Figure 6. At each node, response data are clustered, and trend lines are drawn according to the number of data points at each related note of the grid. The plots at the top row show the variation with respect to M for soft and hard site conditions, whereas the ones at the bottom row are for the variation of R_JB distance with respect to both site conditions.

Considering the plots in Figure 7, it can be stated that DS values are sensitive to T for M values of 7.0–7.1 and R_JB distance interval of 0–5 km, whereas there is only a slight sensitivity for M = 6.5 and 5 < R_JB < 10 (in km). DS values tend to decrease as site conditions shift from soft to hard. For the rest of the magnitudes and distance intervals, T does not seem to be effective on DS values. It is also observed that DS values generally have a decreasing trend as the period increases. This is due to the fact that structures with shorter periods (i.e., more rigid systems) attract higher earthquake forces, and therefore they are more liable to seismic damage.

The second parameter to be investigated is η with five discrete values, as shown in Table 1. The other structural and model parameters are kept constant at their average values. The response data and the trend lines that show the variation of DS with η for different values of seismological parameters are shown in Figure 7. From the plots, it is observed that η has more influence on DS values when compared to T for all magnitudes except M = 5.0 and for all R_JB distances. As M increases and R_JB distance becomes closer, the sensitivity of DS on η also increases. The interval of 20–25 km seems to have higher sensitivity than the interval of 15–20 km unexpectedly. This is mainly due to the smaller number of records in the distance interval of 15–20 km in accordance with the selection of the earthquake stations in the case study region. Furthermore, just like in the case of parameter T, DS values tend to decrease slightly as site conditions shift from soft to hard.

As the final global structural parameter, the sensitivity of DS on μ is investigated in terms of seismological parameters as the remaining parameters are again kept constant at their average values (Figure 8). The trend lines in the plots reveal that DS values are slightly influenced by μ only for magnitudes M = 7.0 and 7.1 and for 0 < R_JB < 5 (in km). The decreasing trend on DS values with a shift from soft to hard site conditions that has already been obtained for parameters T and η is also valid for parameter μ.

In summary, among the seven global structural and model parameters that represent ESDOF models for masonry structures, η seems to have the highest influence on seismic response and damage, followed by parameters T and μ. Hysteretic model parameters were found to have no significant influence on damage states. In all cases, the sensitivity of DS increases with a larger magnitude and closer distance. In addition, site conditions also slightly influence seismic response in the sense that the damage states are more critical for soft sites rather than stiff sites.

In the third and the last phase of seismic response evaluation, damage state probabilities of the ESDOF models are calculated in terms of magnitude and PGA for both soft and stiff site conditions. All the response data from dynamic analyses are employed to determine the probabilities of exceeding each limit state as previously defined in Figure 5. In order to generate these plots, the probability density function (PDF) for each DS is obtained first from the relevant response data for PGA values ranging from 0.1 g to 1.0 g and six discrete M values (5.0, 5.5, 6.0, 6.5, 7.0, 7.1). Then, the limit state curves are derived by calculating the number of cases exceeding the related damage state, normalized by the total number of cases for each ground motion intensity level. Finally, this operation is repeated for all the discrete values of PGA and M to obtain the statistical values in Figure 9, which compares the limit state curves for both soil conditions. From the plots, it is observed that the site conditions do not have a significant influence on the limit states of ESDOF models except for a slight difference in the case of LS-1 for both soil conditions. Another important observation is that LS-2 and LS-3 are generally close to each other, so they can be treated as a single-limit state. This is an expected result due to the fact that masonry structures have limited ductility capacity and a narrow margin of safety between the initial state of inelastic behavior and the ultimate (collapse) state. Figure 9a reveals that there is a 50% probability that masonry structures generally start to experience inelastic action at PGA values of 0.4 g, whereas they are expected to be imposed to heavy damage or collapse by a probability of more than 50% as PGA values exceed 0.9 g. In terms of magnitude (Figure 9b), it can be stated that 25% of the ESDOF models that represent generic masonry structures are expected to receive damage (i.e., exceeding LS-1) when M = 6.0, whereas the percent of damaged building models increases to 60% when M = 7.0. In addition, the percentage of collapsed structures is nearly 25% when M = 7.0. The authors believe that this type of global seismic risk quantification is very valuable to estimate regional earthquake damage and loss and take precautions to mitigate the adverse effects of potential seismic actions.

5. Statistical Evaluation of Global Model Parameters

For statistical evaluation of the results of this study, two alternative methods are selected: discriminant analysis (DA) and multiple regression analysis (MRA). DA is employed to study the discriminating ability of the independent variables by classifying the data into various groups. In the next stage, MRA is used to predict DS, given the most influential independent parameters [56].

5.1. Discriminant Analysis (DA)

DA is a multivariate technique for classifying the relative weights of independent variables among different case groups. This analysis provides a linear discriminant relationship between dependent and independent variables as follows:

$$Z_i=W_1{X}_1+W_2{X}_2+\dots \dots +W_n{X}_n$$

(11)

where Z_i is the score of discriminant relationship i, W_n is the discriminant coefficient of independent variable n, and X_n is the nth independent variable.

Herein, the DA is employed to obtain the independent variables which discriminate the DSs effectively. The independent variables regarding seismological parameters include earthquake magnitude, local soil conditions at the recording site (Vs30) and RJB distance. The independent variables in terms of structural parameters are period, strength ratio, ductility ratio, post-yielding ratio, post-capping ratio, residual strength ratio and hysteretic model degradation. The discriminant analysis procedure is applied to the dynamic response analysis results. The statistical significance of DA can be evaluated in terms of the eigenvalue (λ_e) and Wilks’ lambda (Λ). According to Bruin [57], λ_e and Λ can be estimated through the canonical correlation (cc) as follows:

$${\mathsf{\Lambda }}_e={cc}^2/(1-{cc}^2)$$

(12)

$$\mathsf{\Lambda }=1-{cc}^2$$

(13)

For cases with cc approaching the value of 1.0, the highest discriminating relationship of DS with the independent variables is expected. Similarly, a Λ value close to 0 or a λ_e approaching infinity indicates a strong discriminating ability. The results of the DA are presented in

Table 3. Eigenvalues of the discriminant relationships.

Discriminant Relationships	Eigenvalue	% of Variance	Cumulative %	Canonical Correlation
1	1.299	96.7	96.7	0.752
2	0.036	2.7	99.4	0.188
3	0.008	0.6	100.0	0.087

,

Table 4. Coefficients of the independent variables for the best relationship.

Parameter	Relationship 1
Strength Ratio (η)	0.812
Magnitude (M)	−0.687
Distance (RJB)	0.513
Soil (S)	0.141
Period (T)	0.105
Ductility Ratio (µ)	0.027
Residual Strength Ratio (λ)	−0.015
Post-capping Ratio (αr)	0.003
Post-yielding Ratio (αs)	0.001

and

Table 5. Classification of the results between the groups of DS.

DS	Predicted Group Membership				Total
	DS-1	DS-2	DS-3	DS-4
DS-1	77.5	16.7	4.5	1.4	100.0
DS-2	4.6	56.3	26.4	12.7	100.0
DS-3	0.0	23.5	53.6	22.9	100.0
DS-4	0.0	7.8	24.0	68.2	100.0

.

In this study, three discriminant relationships are obtained, and the statistical results of these relationships are presented in Table 3. Accordingly, the first relationship should be selected as its eigenvalue is greater than 1.0, and the canonical correlation of the eigenvalue at the first relationship is the closest to 1.0. It is also observed that the first relationship has a 96.7% of discriminating ability as far as three continuous discriminant relationships are considered. Dividing the eigenvalue of the relationship to the sum of all eigenvalues determines the percentage of variance.

Table 4 presents the coefficients of independent variables for the best relationship. The independent variables include M, S, R_JB, T, η, µ, α_s, α_r, λ and γ, which stand for magnitude, soil condition, Joyner and Boore distance, period, strength ratio, ductility ratio, post-yielding ratio, post-capping ratio, residual strength ratio and degradation parameter, respectively. Results reveal that the strength ratio, magnitude of the earthquake and source-to-site distance are the most effective parameters for the URM structures in the order of the discriminating power.

The period is a crucial structural metric for identifying seismic behavior during earthquake excitation. However, it appears that the strength ratio of URM structures dominates more significantly. The effects of other structural characteristics, including ductility, residual strength ratio, post-yielding ratio, post-capping ratio and hysteretic model degradation parameter, are not significant for the considered URM structural models.

Next, classification between different DSs is performed. As observed in Table 5, 77.5% of DS-1 is accurately classified, whereas this value is 56.3%. 53.6% and 68.2%, respectively, for DS-2, DS-3 and DS-4. Overall, 71.7% of the original grouped cases are accurately classified. Findings also reveal that DS-1 and DS-4 are better classified compared to other damage states. This could be attributed to their less uncertain nature.

5.2. Multiple Regression Analysis

In this study, MRA, the most common form of simple linear regression, is employed to predict displacement values. The multiple regression relationship is similar to the discriminant relationship in terms of the dependent variable, which is identical for both analyses. MRA relies on a continuous dependent variable, while DA classifies items into discrete groups based on a discrete dependent variable. The equation of MRA is given as follows:

$$Y=a+c_1{K}_1+\dots \dots +c_n{K}_n$$

(14)

where Y corresponds to the value of the dependent variable, c_n is the coefficient of the n^th independent variable, K_n is the value of the n^th independent variable, and finally, a is a constant.

The regression results in terms of standardized coefficients of the independent variables are shown in

Table 6. Coefficients of the independent variables in MRA.

Parameter	Unstandardized Coefficients	Std. Error of the Unstandardized Coefficients
Constant (a)	0.316	0.005
Magnitude (Mw)	0.474	0.001
Soil Type (S)	−0.148	0.001
Distance (RJB)	−0.039	0
Period (T)	−0.653	0.003
Strength Ratio (η)	−1.632	0.002
Ductility Ratio (µ)	−0.044	0.001
Post-yielding Ratio (αs)	−0.146	0.011
Post-capping Ratio (αr)	−0.041	0.005
Residual Strength Ratio (λ)	0.072	0.003
Degradation Parameter (γ)	−1.78 × 10−5	0

. The multiple regression relationship of this study is developed using the unstandardized coefficients:

$$\begin{array}{l}{Y}^{\prime }\\ =0.316\\ +0.474M-0.148S-0.039R_{JB}-0.653T-1.632\mathrm{ɳ}-0.044\mathrm{\mu }-0.146{\alpha }_s-0.041{\alpha }_r\\ +0.072\lambda -0.0000178\gamma \end{array}$$

(15)

Results reveal that the strength ratio of the structure, the period of the structure and the magnitude of the earthquake are the most dominant parameters for this problem. The performance of the linear regression is evaluated through different model performance indicators, and the results are provided in

Table 7. Model summary of MRA.

R	R2	Adjusted R2	Std. Error of the Estimate
0.721	0.520	0.520	0.692

. To this end, the term R, as the multiple correlation coefficient, measures the strength of the relationship between the independent and dependent variables. The R value is obtained as 0.72. The term R² is the measurement of the variance in the dependent variable caused by the model’s independent variables. The adjusted R² is another version of R² which is modified for the number of independent variables in the model. This metric only increases if a new predictor makes the model more accurate. Both R² and adjusted R² values in this study are obtained as 0.52, which indicates a moderate correlation.

6. Conclusions

In this study, the seismic performance of unreinforced masonry (URM) structures is investigated by subjecting them to region-specific simulated ground motion time-series. The study utilizes a simulated ground motion dataset generated for Duzce, Turkiye, using a validated set of input parameters through the stochastic finite-fault ground motion simulation approach. Statistical methods are then employed to analyze the impact of various seismological and structural characteristics on the performance of URM structures. Masonry structures are idealized by using Single Degree of Freedom (SDOF) models. Limit states are defined in terms of displacement using selected structural properties as the explanatory variables. Sensitivity analysis is performed to evaluate relationships between damage states and seismological/structural parameters and to identify the most influential parameters in estimating damage states. The study employs data-driven methods, namely discriminant analysis and multiple regression analysis, in an effective manner to assess the seismic vulnerability and regional risk assessment of masonry structures.

The main findings of this study are the following:

• This study employs a significant number of structural and seismological parameters to assess the seismic performance measures of simplified URM models. It is a very challenging task to deal with the uncertainties arising from seismic demand due to ground motion variability and seismic capacity of URM structures with highly unpredictable seismic behavior. The variation in response statistics in such a parametric study is expected to be huge, and it may not be possible to obtain solid conclusions from the highly scattered results. However, this study differs from many others due to the fact that it utilizes a region-specific simulated ground motion set instead of real ground motion records extracted from different earthquakes all around the world. This consideration has two advantages: first, the ground motion variability is relatively reduced and second, the results are more consistent and interpretable for seismological parameters since region-specific simulations have been carried out. This is also reflected in the coefficient of determination (R2) statistics of the proposed regression formulation, which shows moderate correlation as a promising result for the seismic response statistics of URM structures;
• The statistical results reveal that deficient URM buildings close to the fault rupture (RJB < 15 km) are expected to collapse under a major earthquake with a magnitude M > 6.5. Superior URM buildings seem to experience slight to moderate damage, while typical URM buildings with average parameter values exhibit moderate to heavy damage under the same seismic conditions. The sensitivity of damage states is reduced with lower magnitude and longer distance;
• The influence of site conditions on damage levels of URM buildings does not seem to be very significant, but in favor of hard soil conditions compared to the soft ones;
• Among the structural parameters, the strength ratio seems to be the most dominant one regarding the seismic performance of URM models. Period and ductility seem to have relatively moderate effects on seismic performance;
• Damage state probability plots of URM models verify that masonry buildings have limited ductility capacity and a narrow margin of safety between the initial state of inelastic behavior and the ultimate (collapse) state. Hence, intermediate limit states are close to each other and difficult to identify. We conclude that using only two limit states (onset of damage and ultimate) for URM buildings during regional damage and loss estimation studies would be more appropriate.