Ground-motion scenarios consistent with PSH deaggregation for Tehran , capital city of Iran

This study presents the results of probabilistic seismic hazard (PSH) deaggregation for 5 %-damped 0.2 and 1.0 s spectral accelerations, corresponding to mean return periods (MRPs) of 50, and 475 yr for Tehran city. The aim of this paper is to quantify the dominant events that have the most contribution on ground-motion exceedance from the above mentioned hazard levels. The scenario earthquakes are characterized by bins of magnitude ( M), source-to-site distance ( R), and epsilon ( ε). The results reveal that for Tehran city, the hazard is mainly controlled by local seismicity. Generally, as the spectral acceleration period increase, the contribution of larger and more distant scenario earthquakes to the overall seismic hazard increase.


Introduction
One of the advantages of probabilistic seismic hazard assessment (PSHA) is that it allows computation of the mean annual rate of exceedance at a particular site based on aggregated hazard of all possible range of magnitudes from potential seismic sources, occurring at many different source-tosite distances.This integrating nature of PSHA leads to lose the concept of a single event threat for the site (McGuire, 1995).The so-called dominant events (Bazzurro and Cornell, 1999) or design earthquakes (McGuire, 1995) are required to evaluate many characteristics of the ground motion such as dynamic time history analysis, duration, nonstationary of motion and critical pulses (Cramer and Petersen, 1996).Design earthquakes can be chosen based on deterministic assumptions, such as the magnitude of the worst historical event reported, and its best guessed location derived from known active structures.On the other hand, usually for the common case of a site, high-frequency hazard is dominated by frequent local earthquakes and low-frequency hazard is dominated by large magnitude earthquakes over a wide range of distances, therefore the selection of a unique earthquake scenario does not guarantee an equal seismic protection for the site under analysis (Sousa and Costa, 2009).In order to derive a seismic scenario consistent with the results of PSHA for a site and determine the relative contribution of events to the overall seismic hazard, the concept of deaggregation was introduced.The deaggregation process separates the contributions to the mean annual rate of exceedance (MRE) of a specific ground-motion value at a site due to scenarios of given magnitude M, distance R, and, often, ε, the groundmotion error term.The ε value is defined as the number of standard deviations by which the (logarithmic) ground motion generated by a given M-R pair deviates from the median value estimated by a prediction equation (Barani et al., 2009).The seismic hazard deaggregation process has been performed by Bernreuter (1992) for the first time, and following this study, deaggregation methods and applications were extensively discussed and improved by Chapman (1995), McGuire (1995), Cramer and Petersen (1996), Bazzurro and Cornell (1999), Harmsen et al. (1999), Sokolov (2000), Harmsen (2001), Harmsen and Frankel (2001), Montilla et al. (2002), Halchuk and Adams (2004), etc.In order to perform deaggregation of seismic hazard for Tehran city, probabilistic seismic hazard assessment is evaluated for the study area encompassed by 49.5-54 • E longitudes and 34-37 • N latitudes.Selecting this area allows to get involved in assessments and calculations all of the elements that could affect the seismicity of Tehran.The megacity of  Tehran with highly dense population (about 8.5 million people) and political and economic centralization has not experienced a damaging earthquake since 1830, although it is located in an active seismic zone and surrounded by active faults with evidence of historical destructive earthquakes.The major aim of this study is the deaggregation of the hazard results in terms of magnitude, distance and epsilon, to investigate earthquake occurrences that have the most contribution to the resulting ground-motion hazard in Tehran.In this regard, 20 area seismic sources are delineated in the study region and the quantified hazard values in terms of the peak ground acceleration (PGA) over bedrock, are mapped for 63 % and 10 % probabilities of exceedance in 50 years on grid intervals of 0.1 • , using the attenuation relation of Boore and Atkinson (2008) NGA.Then, based on obtained results of PSHA, deaggregation of the hazard is performed in terms of magnitude-distance pairs, and epsilon for Tehran city.

Tectonic setting
The active deformation of the Iranian plateau, as a tectonically active part of Alpine-Himalayan tectonic belt, is controlled by the Arabian-Eurasian convergence.Shortening is accommodated by distributed faulting and it makes regions with different seismotectonic features.The study area covers the Alborz in the north and northern part of Central Iran in the south (Fig. 1).The Alborz is a narrow mountainous area of approximately 100 km wide, curving along the South Caspian Sea.The mean elevation of the Alborz decreases from +3000 m above the mean sea level in the inner belt to −28 m at the Caspian Sea shoreline in the north, and to approximately +1500 m north of Central Iran in the north of Tehran (JICA and CEST, 2000).Despite the relatively high elevation of Alborz Mountains, the Bougure anomaly has very small negative values and the crustal thickness is less than 35 km along the Caspian coast; therefore, it must be considered as mountains without root whose uplift is due to thrusting of allochthonous masses over each other in a compressional tectonic regime (Dehghani and Makris, 1983).Fault plane solutions of earthquakes in the study Table 1.Beginning year of complete reporting of earthquakes in the Alborz-Azerbaijan and Central-East Iran seismotectonic provinces (Mirzaei et al., 1997a).Table 2. M s − m b relationships for Alborz-Azerbaijan and Central-East Iran seismotectonic provinces (Mirzaei et al., 1997a).
region indicate both thrust and strike-slip faulting, though the thrust faulting is prominent (Fig. 1).Because the Alborz seems to be rootless, it is probable that many active faults are present in the basement rocks and do not necessarily have surface expressions (JICA and CEST, 2000).Based on recent global positioning system (GPS) studies, the rate of northsouth shortening throughout the Alborz is 5±2 mm yr −1 , and the left-lateral shear of the overall belt occurs at a rate of 4 ± 2 mm yr −1 (Vernant et al., 2004a).Central Iran as a part of the Central-East Iran seismotectonic province is located in the south of the Alborz and forms a range-and-basin system with a lower altitude than the Alborz, with active mountain-bordering thrust and some long strike-slip faults (JICA and CEST, 2000).Central Iran has much in common with Alborz in its structure and history; thus the sedimentary record is so similar in the two regions that many rock-unit names of Paleozoic, Mesozoic and Tertiary formations can be applied in both the Alborz and a large part of Central Iran (Stocklin, 1974).In addition, both regions were folded and faulted with similar intensity by the Alpine diastrophism in Tertiary time (Stocklin, 1974;Mirzaei et al., 1998).The north part of Central Iran in the vicinity of the Alborz approximately follows the trend of the Alborz and shows strong tectonic activity with surface deformation during both historic and present time.In general, active faults in Central Iran are parallel to regional structures and are related to the overall northeasterly-directed compressional stress (JICA and CEST, 2000).Recent GPS data taken from a station south of Pishva fault in northern Central Iran, in the study area, indicates northward and eastward velocities with a rate of 11.6 and 2.8 mm yr −1 , respectively (Vernant et al., 2004b).

Probabilistic seismic hazard assessment
According to traditional methodology of PSHA (Cornell, 1968), at first, potential seismic sources should be delineated and then seismicity parameters in each source will be determined.For regions with insufficient earthquake data, seismicity studies in potential seismic sources using the traditional method brings a noticeable uncertainty into results.To overcome such shortcomings, the traditional method has been modified by Shi et al. (1992); therefore, inhomogeneity of seismicity in time and space could properly be regarded.
In this study, the modified methodology of PSHA is applied.Three major steps of modified method are (Shi et al., 1992;Shabani and Mirzaei, 2007) 1. delineation of seismotectonic provinces and evaluation of seismicity parameters (b value, annual mean occurrence rate (λ), and maximum possible magnitude, M max ) in each seismotectonic province; 2. determination of potential seismic sources; estimation of M max in the potential seismic sources and evaluation of spatial distribution function for each magnitude interval in each potential seismic source; and 3. dividing the interest region into a series of grid points and assessment of seismic hazard for every grid point, using characteristics of seismic activity in seismotectonic provinces and potential seismic sources through an earthquake ground-motion attenuation relationship.
Seismotectonic province is defined as an area that under the present-day geodynamic regimes has a comparable tectonic setting and unified seismicity pattern (Ye et al., 1995;Mirzaei et al., 1998).According to Mirzaei et al. (1998), territory of Iran can be divided into five major seismotectonic provinces based on all available geophysical, geological, tectonic and earthquake data.The study area is mainly located in the Alborz-Azerbaijan seismotectonic province in the north, and a small part of the Central-East Iran seismotectonic province in south (Fig. 1).
A uniform catalog of earthquakes containing historical and instrumental events covering the period from 400 BC to 2010 is used.The earthquake database is mainly compiled from ISC and USGS/NEIC for the modern instrumental time period , and the catalog of earthquakes provided by Ambraseys and Melville (1982) is the basic source of parameters for the historical (before 1900) and early instrumental  time periods.According to Mirzaei et al. (1997a), the years in which there was first complete reporting of earthquakes in the Alborz and Central-East Iran seismotectonic provinces are given in Table 1.The catalog Table 3. Estimation of seismicity parameters from recorded earthquakes in the Alborz-Azerbaijan and Central-East Iran seismotectonic provinces (Mirzaei et al., 1997b   of earthquakes has been made uniform using the relationships between M s and m b defined by Mirzaei et al. (1997a) for Alborz-Azerbaijan and Central-East Iran seismotectonic provinces (Table 2).Since we encountered an incomplete earthquake catalog in the study region, the procedures introduced by Kijko and Sellevoll (1992), which permit incorporation of magnitude uncertainty to estimate seismicity parameters from incomplete data files, are applied to the uniform catalog of earthquakes for estimating the seismicity parameters.Results of the estimation of seismicity parameters in the Alborz-Azerbaijan and Central-East Iran seismotectonic provinces are presented in Table 3.
Figure 2 shows epicenter distributions and potential seismic source zones in Tehran region.A total of 20 area sources and their corresponding M max are chosen after Mirzaei et al. (1999).Seismic sources have been delimited mainly based on the fault extent, seismogenic crust (a part of the earth crust in which large earthquakes usually originate), and mechanism of earthquake faulting or a type of active faults.The estimation of maximum magnitude in potential seismic source zones is usually according to the features of seismic activity and tectonic analogy.It implies that structures with analogous tectonic setting are capable of generating the same size earthquakes, and absence of earthquake record for a structure is not evidence for absence of earthquake occurrence on it (Mirzaei et al., 1999).
Generally, in seismic hazard analysis, annual mean occurrence rate of earthquakes in potential seismic sources is calculated by definite integration on magnitude-frequency relation.This method cannot reflect spatial inhomogeneity

Zone
M max 4.0 < M s ≤ 4.5 4.5 < M s ≤ 5.0 5.0 < M s ≤ 5.5 5.5 < M s ≤ 6.0 6.0 < M s ≤ 6.5 6.5 < M s ≤ 7.0 7.0 < M s ≤ 7.5 of seismicity and reveal a realistic activity rate of small and large-magnitude earthquakes in the potential seismic sources, especially for insufficient earthquake data.To fulfill this problem the concept of spatial distribution function for specified magnitude intervals has been introduced by Shi et al. (1992).In this method the annual mean occurrence rate of earthquakes in a seismotectonic province should be allocated to each magnitude interval in the corresponding potential seismic sources in order to properly reflect the inhomogeneity of seismicity in time and space, and to avoid the underestimation of a potential hazard of large magnitude earthquakes.Different kinds of seismological, tectonic, and geophysical data can be used to indicate the possible future earthquake activities in the interest region, providing basis for evaluation of spatial distribution function.In the study region, four controlling factors (Shabani and Mirzaei, 2007) are considered for evaluation of the spatial distribution function: 1. the reliability of delineation of potential seismic sources; 2. tectonic setting of potential seismic sources; 3. structural elements in potential seismic sources; 4. characteristics of seismic activity in potential seismic sources.
For calculation of spatial distribution function based on the controlling factors, the method of equal weight summation given by Yan (1993) is used.In this method for the selected factor (k) and each magnitude interval (m j ) in the l-th potential source, a distribution coefficient (W lm j k ) is assigned.Then in each seismotectonic province the distribution coefficient is normalized to obtain the factor load (Q lm j k = W lm j k W lm j k ).In next step, loads of controlling factors in each potential source are used to define total load (R lm j = Q lm j k ).Eventually the total load is normalized in each province to obtain the spatial distribution function for the j -th magnitude interval in the l-th potential seismic sources (f lm j = R lm j R lm j ) (for more details of controlling factors and spatial distribution function calculation steps see Shabani and Mirzaei, 2007).
Table 4 shows M max in potential seismic sources and spatial distribution function for different magnitude intervals."Background earthquake" for each seismotectonic environment is considered.In the concept of background earthquake (background seismicity), small and moderate-sized earthquakes may occur in the defined area randomly.We consider threshold magnitudes of M s = 5.5 and 6.0 for background seismicity of Central-East Iran and Alborz-Azerbaijan seismotectonic provinces, respectively.
For the l-th potential seismic source in the seismotectonic province, the annual mean occurrence rate of the j -th magnitude interval is as follows (Gao, 1988;Shi et al., 1992;  Mirzaei et al., 1997b): where λ l,m j and f l,m j are the annual occurrence rate and spatial distribution function of the j -th magnitude interval in the l-th potential seismic source, respectively, and λ is the annual mean occurrence rate of the seismotectonic province; β = b ln 10 and b is the b value in the frequency-magnitude relationship of the seismotectonic province; m j is the central value of the j -th magnitude interval; sh is the hyperbolic sine function; M is the width of magnitude interval; and M min and M max are the minimum magnitude that can affect the engineering site (usually M s = 4.0) and the maximum expected magnitude in the seismotectonic province, respectively.
In the regions with high seismic activity, any seismic hazard assessment needs a set of attenuation equations which could appropriately represent the variable parameters of the strong motion in that region.For the Iranian plateau, there is a limited number of reliable strong ground-motion data to perform a specific attenuation relationship, especially in the near-source areas of large earthquakes.For example the Iranian strong-motion data have a total of only 4 recordings for magnitudes greater than 7.0 and distances less than 30 km (Shoja-Taheri et al., 2010).Due to this reason, applying worldwide attenuation models for similar tectonic environment is appropriate.We have adopted for this study the attenuation relationship for rock (V s > 750 m s −1 ) developed by Boore and Atkinson NGA model (2008).The Next Generation Attenuation (NGA) project has resulted in the publication of five new ground-motion models to predict PGA, PGV, and response spectral ordinates for periods 0.0 up to 10.0 s (Abrahamson and Silva, 2008;Boore and Atkinson, 2008;Campbell and Bozorgnia, 2008;Chiou and Youngs, 2008;Idriss, 2008).The NGA database is one of the largest databases of uniformly processed strong ground-motion recorded in shallow crustal earthquakes in active tectonic regions, including data recorded in the 1978 Tabas and 1990 Manjil earthquake that occurred, respectively, in central-east Iran and Alborz-Azerbaijan seismotectonic provinces of Iran.Thus, these new equations might also be applied as ground-motion prediction in other parts of the world with the similar tectonic environment, such as the region encompassing the Iranian plateau.Shoja-Taheri et al. ( 2010) carried out a comprehensive study on comparison of the three NGA equations (Boore and Atkinson, 2008;Campbell and Bozorgnia, 2008;Chiou and Youngs, 2008) as representative of all the NGA models, with the strong-motion data recorded in Iran, and concluded that the three NGA models are generally applicable to the presently available Iranian dataset.Shoja-Taheri et al. ( 2010) also mentioned that the NGA model of Boore and Atkinson (2008) is somewhat more applicable for Iran than the other two models, since some of the details on the near source effects in it (e.g., effects of hanging wall and footwall, distance of the upper part of the fault rupture to the surface, and so on) are excluded.
Seismic hazard assessment has been performed for a grid of points with 0.1 • interval in longitude and latitude in the study area.It is noteworthy that, we used β values calculated for seismotectonic provinces for each seismic source delimited within them.PGAs over bedrock, for mean return periods (MRPs) of 50 and 475 yr, have been estimated and seismic zoning maps of the Tehran region was determined using EZ-FRISK (version 7.43) computer program (Fig. 3).The results show a good agreement with active faults and seismicity of the region.Sites with the highest level of PGA are mainly adaptable to potential seismic sources number 1, 2, 7, 8, and 10.Within these seismic sources several large earthquakes are documented during historical and instrumental time periods (Fig. 2).The maximum PGAs are estimated to be 0.3 g and 0.12 g for 475 and 50 yr return periods, respectively.The Tehran city is placed in relatively high risk area, as well.Comparison of the seismic macrozonation hazard map of Iran, Standard 2800 (BHRC, 2005), with results of this work shows good compatibility.Generally, as it is obvious from Fig. 3, the PGA value in the Alborz area is higher than the Central Iran area.

Seismic hazard deaggregation
PSH deaggregation is used to identify the distribution of earthquake scenarios that contribute to exceedance of a given spectral acceleration (S a ) level, in terms of magnitude (M), distance (R), and the deviation parameter in the attenuation equation (ε).In other words, if a ground motion of amplitude A occurs at the site of interest, certain magnitudes, distances, and ε values are more likely to have caused that amplitude than others, and the deaggregation shows these relative contributions.Deaggregation results could change with the spectral ordinate and return period, thus more than one single event may dominate the hazard especially if multiple sources affect the hazard at the site.These results can provide useful information on better defining the design scenario and selecting corresponding time histories for seismic design.In most cases, as the probability decreases, the hazard sources closer to the site dominate.Larger, more distant earthquakes contribute more significantly to hazard for longer periods than shorter periods (Halchuk and Adams, 2004).Spectral acceleration of 5 % damping, S a (T ), where T is the period in seconds, has been calculated in order to perform deaggregation of the seismic hazard for Tehran city using the EZ-FRISK (version 7.43) code.Bins of width 0.4 in magnitude, 10 km in distance, and 0.2 in ε are selected.EZ-FRISK code reports the deaggregation results in terms of probability density function (PDF), which is obtained by dividing the probability mass function (PMF) contribution of each bin by the bin's size, thus the PDF representation is independent of the bin's amplitude (Bazzurro and Cornell, 1999).The hazard can be simultaneously deaggregated in different types of bins.EZ-FRISK presents the result of seismic hazard deaggregation in terms of 1-D M, R and ε bins and 2-D M-R bins.A detailed review of different deaggregation representations (PDF and PMF) and procedures (1-D, 2-D and 3-D deaggre-gation, geographical deaggregation in terms of latitude and longitude) can be found in Bazzurro and Cornell (1999).
Usually, both mean and modal values of M, R and ε are reported for identifying the dominant event, because each of them has its advantages and disadvantages.The mean values (M, R and ε) are defined unambiguously and do not depend on the bin size, but they may correspond to a scenario that is not realistic, if there are two or more sources with significant contribution to the hazard.The modal values represent the most likely event that may induce the specified (or larger) acceleration level at the site, but it should be mentioned that these values depend on the dimension of the bins.
Figures 4a through d show the M-R deaggregation plots of the 5 %-damped linear elastic S a (T ) values at periods of 0.2 and 1.0 s which correspond to MRPs of 50 and 475 yr.The results of computed controlling earthquakes have been summarized in Table 5.The 50 yr S a (0.2 s) hazard for Tehran (Fig. 4a) is controlled by almost nearby earthquakes of moderate magnitude.Larger and more distant events dominate the S a (1.0 s) value for the same MRP.In particular, scenario earthquakes ranging from a magnitude M =6.2 at R = 35 km for S a (0.2 s), to a magnitude M = 6.5 at R = 48 km for S a (1.0 s) has been determined to represent the relatively likely ground motion having a 63 % chance of being exceeded in 50 yr.For less probable, extreme ground motions with a 10 % chance of being exceeded in 50 yr, the scenario earthquakes ranging from a magnitude M = 6.4 at R = 25 km for S a (0.2 s), to a magnitude M = 6.7 at R = 31 km for S a (1.0 s) are represented as dominant events (Figs.4c  through d).It is worthwhile to mention again that a dominant event or a design earthquake is defined as a scenario which has the most contribution to the mean annual rate of exceedance of hazard from a specified level.A "dominant event" is not necessarily "the worst scenario" that can affect the site, and it is obvious that because of the nature of the PSHA approach, no single event will ever be able to fully describe the seismic threat at the site.
Figure 5 shows the variation of mean and modal values of M, R, and ε with spectral periods of 0.0, 0.2, 0.5, 1.0, and 2.0 s.The acceleration values deaggregated correspond to MRPs of 50 yr (left-hand panels) and 475 yr (right-hand panels).As a comparison, modal values of both 2-D and 1-D PDFs are displayed.We denote the mode of the 2-D distribution as M * 2-D − R * 2-D , while that of the 1-D one is indicated by M * , R * , and ε * .For an MRP of 475 yr, mean magnitude and distance tend to increase with period.This trend was expected due to the fact that increasing of the mean and modal magnitude and distance with increasing response period is a general behavior for uni-modal deaggregation plots (Harmsen et al., 1999).Unlike M and R, the modal values seem the same constant value for some spectral periods.Although it should be noticed that the modal values are the central values of the M, R, and ε bins used in the calculations; for example, R * = 25 km means that R * can change in range of 20 to 30 km.An analogous behavior can be observed analyzing variation of mean and modal magnitude and distance values for an MRP of 50 yr.Generally, M * 2-D and R * 2-D are lower than or equal to M * and R * , and show more stability with increasing period.For an MRP of 475 yr, variation of R is slightly with period and diagrams of mean and modal distances seem close together.On the contrary, the existence of large difference between mean and modal values of R for an MRP of 50 yr, especially for periods of greater than 0.5 s reveals that more than one seismic source zones contribute significantly to the site hazard for those spectral periods (Montilla et al., 2002).Analyzing ε behavior reveals that both mean and modal values of ε tend to increase as the MRP increases.For an MRP of 50 yr, the dominant ground motions are within 0.5 σ -0.7 σ of the median, while for an MRP of 475 yr, they are within 1.1 σ -1.3 σ .

Conclusions
PSHA was evaluated for an area encompassed by 49.5-54 • E longitudes and 34-37 • N latitudes in order to perform seismic hazard deaggregation for Tehran city.Generally, the trend of the deaggregation results with spectral period (T ) shows that the mean rate of exceeding long-period spectral accelerations is controlled by earthquakes that are larger in size and farther from the site than those dominating the PGA and short-period spectral acceleration hazard.It was also observed that the larger the MRP, the greater the contribution from closer and higher magnitude events and, as the MRP increases, both the mean and modal values of ε tend to increase.Particularly, since Tehran is located in a highseismicity region, nearby seismicity (20-30 km) has often the major contribution to the hazard.For an MRP of 50 yr, controlling earthquakes have the magnitude range between 6 and 6.4.Events with magnitude ranging from 6.4 to 6.8 have the most contribution to the hazard for an MRP of 475 yr and periods longer than 0.2 s, while the smaller magnitude ranging from 6 to 6.4 are dominant for shorter periods and PGA.Regarding deaggregation results it could be concluded that events with most contribution to corresponding hazard levels are affected by the potential seismic sources number 6 and 7, although geographic deaggregation is necessary for determining more precisely.

Figure 1 .
Figure 1.Main active faults, main stress directions and representative focal mechanisms in the study area.Thick line indicates the border of seismotectonic provinces.

Fig. 1 .
Fig. 1.Main active faults, main stress directions and representative focal mechanisms in the study area.Thick line indicates the border of seismotectonic provinces.

Figure 2 .
Figure 2. Map of seismicity and seismic source zones of Tehran region (modified after Mirzaei et al., 1999).

Figure 3 .Figure 3 .Fig. 3 .
Figure 3. Probabilistic seismic hazard map of Tehran region for peak horizontal acceleration on firm-rock site conditions for 63% (a) and 10% (b) probability of exceedance in 50 years.

Fig. 5 .
Fig. 5. Variation of mean and modal values of M, R, and ε with spectral period for Tehran."STDV" in the label on the y-axis of the charts at the bottom of the figure stands for standard deviation.

Table 4 .
Spatial distribution functions for different magnitude intervals and M max in the potential seismic sources and background seismicity sources of the study region.B.S.1 indicates background seismicity of Alborz province and B.S.2 indicates background seismicity of Central Iran province.

Table 5 .
M-R deaggregation results for 0.2 and 1.0 s spectral acceleration values corresponding to MRPs of 50 and 475 yr.