Introduction

We present a seismic hazard assessment of the TMVB that takes into account the seismic activity observed in the region from 1858 to 2014. Although this region in central Mexico is not as seismically active as some other parts of the country, such as the Pacific coast or the Gulf of California, there is evidence of major historical crustal earthquakes (M>7) occurring on the TMVB in the last 500 years (^{García-Acosta and Suárez; 1996}; ^{Suárez and Caballero-Jiménez; 2012}). As in the case of many other continental faults, the crustal faults in the TMVB responsible for these large earthquakes apparently lie dormant for several thousand years.

Although relatively infrequent, the occurrence of these shallow, crustal events in the vicinity of densely populated areas poses an important seismic hazard to the more heavily populated area of Mexico. The larger cities in the country, such as Mexico City, Guadalajara, Morelia and Jalapa lie on the TMVB. According to the latest census of the National Institute of Statistics and Geography in 2010 (http://www.inegi.org.mx/est/contenidos/proyectos/ccpv/), more than 43 million people (approximately 40% of the population in Mexico) live in cities and towns located on this geological province. The observations of the past 500 years indicate that crustal earthquakes of large magnitude occur throughout the TMVB, regardless of the dearth of small magnitude seismicity.

Traditionally, the seismic hazard studies conducted for cities on the TMVB are based only on the instrumental seismicity catalog. However, the background seismicity is very low and, at first glance, this suggests that the seismic hazard is moderate. The Federal Power Commission (CFE) publishes a manual for the design of civil works against earthquakes, which has become a national standard. CFE classifies the TMVB as a region of moderate hazard (https://www.scribd.com/doc/52197523/CFESismo-08). Nevertheless, the 1912 Acambay and the 1920 Jalapa earthquakes, together with historical reports of macroseismic data, suggest that although these crustal earthquakes have long recurrence periods, they may be of relatively large magnitudes and take place in close vicinity to major cities and towns. The purpose of this paper is to integrate the large historical earthquakes of the XIX^th century in order to estimate the seismic hazard in the TMVB. Our results indicate that the hazard estimate including the historical data is much larger than previously estimated from instrumental data.

Tectonic overview of the TMVB

The TMVB is a volcanic arc that spans Central Mexico from the Pacific coast to the Gulf of Mexico. It has a variable width that oscillates between 80 and 230 km (Figure 1). ^{Pardo and Suárez (1993}; ¹⁹⁹⁵⁾ showed that the oblique orientation of the TMVB relative to the subduction zone is due to the geometry and to the lateral changes in dip of the subducted Cocos plate beneath central Mexico (^{Suárez et al., 1990}). ^{Pardo and Suárez (1993}; ¹⁹⁹⁵⁾ speculated that to the north of the region where the seismic activity within the subducted slab ceases, the subducted Cocos plate bends sharply downwards. Recent results of the tomographic inversion of teleseismic waves observed by the dense network of temporary seismic stations of the MesoAmerican Subduction Experiment (MASE), ^{Pérez-Campos et al. (2008)} suggests that in central Mexico, the subducted Cocos plate dips at an angle of 75º into the mantle, reaching a depth of 120 km beneath the TMVB.

Figure 1 Extent of the TMVB reported by ^{Ferrari et al. (2012)}. 

The TMVB is affected mostly by extensional tectonics (e.g., ^{Ferrari et al., 2012}; ^{Suter et al., 1995}). ^{Suter et al. (1992}; ²⁰⁰¹⁾ attribute the extensional regime in the Morelia-Acambay fault system, in the central part of the TMVB, to isostatically compensated surface loads in response to the high elevation of the volcanic belt. ^{Ferrari et al. (2012)}, based on the heat flow measurements of ^{Ziagos et al. (1985)}, suggest that the origin of this extension is related to the hot and low-density mantle (~950 to 1000° C) that underlies the Moho beneath the TMVB. This buoyant thermal effect presumably promotes the uplift and subsequent extensional deformation of the volcanic arc.

This extensional regime is characterized by seismically active, east west oriented normal faults throughout the volcanic belt. Most of these faults exhibit pronounced scarps that are clearly mapped on the surface for distances as long as 50 km (Figure 2). Many of these faults are geologically active and cut alluvium and Quaternary scoria volcanoes, indicating recent deformation. ^{Suter et al. (1992)} estimate an average rate slip rate of 2 mm/yr. on these extensional faults. The Quaternary vertical slip rates of these faults have a mean of 0.07 mm/ yr and the bulk extension rate of the TMVB is estimated to be 0.2 ± 0.05 mm/yr (^{Suter et al. 2001}).

Figure 2 Distribution of the faults with Quaternary activity mapped on the TMVB. PVG represents the Puerto Vallarta graben; MG Mecatán Graben; SPC the San Pedro-Ceboruco Graben; ACH the Amatlán de Cañas Halfgraben; PBCM the Plan de Barrancas-Cinco Minas Graben; SM the San Marco fault; CR the Colima Rift; CG the Citlala Graben. The major cities are abbreviated as: GD: Guadalajara; LN: León; MR: Morelia; QR: Querétaro; MC: Mexico City; PB: Puebla; and JL: Jalapa. 

Some of these geologically active faults like the Acambay-Tixmadejé, the Venta de Bravo faults, and the faults bounding the Chapala graben, are located just a few kilometers away from densely populated cities, such as Mexico City, Guadalajara, Puebla or Morelia. Although many of these extensional crustal faults are clearly mapped on the surface, in some cases they are blind and show no surface expression. In the region of the 1920 earthquake, for example, there are no clearly identified faults mapped on the surface.

Historical and instrumental seismicity in the TMVB

During the past 100 years, two large crustal earthquakes in the TMVB have been recorded instrumentally. On November 19, 1912, an earthquake (M_w 6.9) took place near the city of Acambay rupturing one of the faults bounding the Acambay graben, about 80 km from Mexico City (Figure 3). This earthquake caused large destruction in the town of Acambay and in neighboring towns and villages (^{Urbina and Camacho, 1913}) and caused moderate damage in Mexico City (^{Suter, 2015}). A few years later, on January 4, 1920, an earthquake M_w 6.4 took place on the eastern part of the TMVB (Figure 3). According to ^{Suárez (1992)}, this event is one of the more deadly earthquake ever recorded in Mexico, just after the large 1985 Michoacán earthquake (M_w 8.1). More than 1,500 fatalities were estimated. Most of them drowned or buried under the debris flows caused by the landslides on the steep cliffs bounding the course of the Pescados river (^{Comisiones del Instituto Geológico Mexicano, 1922}). The most recent moderatesized earthquake in the TMVB took place on February 22, 1979. The Maravatio earthquake (m_b 5.3) occurred near the western end of the Acambay fault and caused only slight damage (^{Astiz, 1980}).

Figure 3 Spatial and temporal distribution of the more important historical and instrumental seismicity within the TMVB since the XVI^th century. 

The recurrence period of some segments of these crustal faults were estimated on the basis of paleoseismological studies in the vicinity of the Acambay graben. ^{Langridge et al. (2000)} conducted a study on the fault trace of the 1912 Acambay earthquake. Their results show the presence of three previous ruptures with an approximate recurrence time of 3,600 years. Similar paleo seismological studies conducted on the neighboring faults, on the southern boundary of the Acambay graben, indicate rupture intervals ranging from 600 to 10,000 years, depending on whether each fault breaks as a single unit or in separate segments (^{Langridge et al., 2013}; ^{Ortuño et al., 2015}). Thus, based on the paleo seismological data, a significant earthquake may be expected in this region of central Mexico every 300 to 600 years (^{Zúñiga et al., 2012}).

Prior to the instrumental period, several large earthquakes have been identified as crustal events in the TMVB, based on the damage and felt reports. The earliest reference is the December 27, 1568, earthquake (Figure 3). Fortunately, for an earthquake that occurred so early in the written historical record of Mexico, there are detailed descriptions of the damage and of the substantial and extended ground deformation caused by this earthquake (^{Suárez et al., 1994}). Practically all of the newly built Franciscan churches in this area, located to the southwest of the city of Guadalajara, were destroyed. Based on the area of large seismic intensity, ^{Suárez et al., (1994)} suggest a magnitude of M_w 7.0. Later, ^{Suter (2015)} suggested a magnitude M_w 7.2 for this earthquake interpreting the ground rupture reported by historical sources as evidence of the fault trace. There are contradictions in the historical reports regarding the year of this earthquake. ^{Suter (2015)} suggests that this earthquake actually took place in 1567 and not in 1568.

Based on the inversion scheme of seismic intensity data proposed by ^{Bakun and Wentworth (1997)}, ^{Suárez and Caballero Jiménez (2012}) and ^{Suárez et al. (2016)} determined the approximate epicentral location and the magnitude of three large earthquakes in the TMVB that occurred during the XIX^th century. The Santa Juliana earthquake on June 19, 1858, (M 7.6), the San Cristobal event of February 11, 1875, (M 7.0) and the Pinal de Amoles earthquake of November 26, 1887, (M 6.1) have numerous and detailed historical reports (^{García Acosta and Suárez, 1996}). The number and quality of the macroseismic data available was sufficient to allow for a quantitative determination of the source parameters using the inversion scheme of the observed Modified Mercalli Intensities (MMI) as a function of distance (Figure 3). ^{Singh et al. (1996)} interpreted the Santa Juliana earthquake as an in-slab event. ^{Suárez and Caballero-Jiménez (2012)} and ^{Suárez et al. (2016)} provide evidence that this event is a crustal event.

Probabilistic seismic hazard assessment

The seismic hazard analysis presented here was developed following the methodology proposed by ^{Cornell (1968)} and ^{Esteva (1968)}. To this end, we compiled a robust seismicity catalog identifying the location of the seismic sources of interest that allows us to compute the mean recurrence periods of earthquakes in the region, as a function of magnitude. The seismic hazard estimations were obtained using the computer program CRISIS 2015, an updated version of CRISIS 2007 (^{Ordaz et al., 2007}; https://ecapra.org/crisis-2007).

Seismicity catalog and its level of completeness

The seismicity catalog used was based on an exhaustive search of the seismic databases. Thus the data were culled from the catalogs of the International Seismological Centre (ISC), the U.S. Geological Survey (USGS), the Global Centroid Moment Tensor Catalog (CMT), and of the Servicio Sismológico Nacional (SSN) were consulted. Duplications and aftershocks were avoided.

A total of 567 historical and instrumental earthquakes were initially selected. A maximum depth cutoff of 33 km was defined to guarantee that the earthquakes occurred in the crust (Figure 4). Due to the fact that only moderate and large earthquakes contribute to large ground accelerations, only M ≥ 4 earthquakes were taken into account in this study. Furthermore, due to the dearth of seismic stations in the region, events M<4 are poorly located and there is a great uncertainty in their focal depth. In fact, many of the small earthquakes presumed to be in the TMVB crust are probably deeper, in-slab events within the subducted Cocos plate (green circles on Figure 4).

Figure 4 Epicentral locations of earthquakes observed on the TMVB sfromince 1858 to 2012. These events represent the selected seismic catalog that constitutes the basis of the seismic hazard analysis presented here. 

A selected dataset of 43 events satisfied the selection criteria and constitute the basis of this analysis. The construction and operation of modern seismic networks, such as the national network of the SSN and the Worldwide Standard Seismic Network (WWSSN) substantially improved the detection capability of small earthquakes. Since 1964 the number of located small earthquakes has increased dramatically. ^{Zuñiga et al. (2016}) propose that the seismicity catalog of Mexico is complete for events approximately m_b > 3.6 since 1964. This observation is based on a completeness analysis of the Mexican catalog. The temporal changes and threshold levels observed in the Mexican seismicity catalog are similar to the variations observed worldwide (e.g., ^{Zúñiga and Wyss 1995}; ^{Zúñiga et al. 2005}; ^{Michael 2014}). In a conservative manner, we consider here that the seismicity catalog of the TMVB is complete for M ≥ 4 since 1964.

On the other hand, the population on the TMVB has been relatively dense since the XIX^th century. As the historical catalog demonstrates, it is very likely that the damage and effects of moderate sized earthquakes would have been reported in the historical accounts. Thus we qualitatively suggest that the catalog is complete for earthquakes M ≥ 6 since the second part of the XIX^th century.

Seismicity models

The more frequently used statistical model to represent the rate of occurrence of earthquakes is the Gutenberg-Richter relation. This relation establishes that the number of earthquakes is a function of the magnitude. This frequency magnitude relation is expressed as:

1

This relation between magnitude M and the logarithm of the number of earthquakes N is parameterized by the intercept of a and the slope b of the linear fit (^{Ishimoto and Ida, 1939}; ^{Gutenberg and Richter, 1944}). A more conservative expression of earthquake occurrence is obtained by associating a maximum magnitude, M_max, to the distribution (^{Aki, 1965}). The modified seismicity model is expressed as:

where, λ_c = N/A is the exceedance rate of the minimum magnitude earthquake (M_c) in the catalog. β =1/(M_p-M_c) = b ln(10). N is the number of earthquakes with magnitudes M ≥ M_c, the parameter A is the duration of the seismicity catalog and M_p is the average magnitude of all the events (^{Aki, 1965}).

Three seismicity models were constructed based on the compiled seismicity catalog. The first one, Model A, is a traditional parametric model following the modified Gutenberg-Richter relation and using only the instrumental catalog (green curve on Figure 5). The exceedance rate for a cut-off magnitude λ_c and the β value were computed using the maximum likelihood method proposed by ^{Aki (1965)}. The resulting seismicity parameters for this model are λ_c = 3.333 and β = 0.750. These values are similar to those reported by ^{Ordaz et al., (2012)}, which are the current values used to quantify the seismic hazard of the TMVB.

Figure 5 The exceedance rates of earthquakes expected in the TMVB are shown as a function of magnitude. Model A (red curve) is obtained from the fit of instrumental data according to a conventional Gutenberg-Richter distribution and the maximum likelihood method. Model B (blue curve) results from merging the historical and instrumental data. This semi parametric approach appears to be a conservative description of the exceedance rates of earthquakes observed in the TMVB since the 16th century. Finally, seismicity Model C (green curve) is obrtained by fusing the exceedance rates reproduced by the fit of instrumental data (4 ≤ M < 6) with the exceedance rates of earthquakes M > 6 observed from the seismic catalogs. 

Two other non-parametric models are proposed that integrate the historical and the instrumental seismic data (Figure 5). The seismicity Model B (blue line on Figure 5) represents an effort to describe the seismic activity of the TMVB in a semi-parametric manner, using both the instrumental and the historical earthquakes of the catalog. In this approach, two different Gutenberg-Richter distributions are estimated: one of the curves corresponds to the instrumental data and the other to the historical seismicity. Thus considering the assumption that the seismicity catalog is complete for magnitudes M≥4 since 1964 and for M≥6 since 1858, we estimate one Gutenberg-Richter linear function for earthquakes 4≤M<6 earthquakes, and a separate one for earthquakes M≥6 events. As we have not observed any earthquake M≥6 on the TVMB since 1964, we propose to fuse both fits into one curve by adding the exceedance rate of instrumental earthquakes to historical events with M>6. The resulting seismic parameters are: λ_c = 0.033 and β = 1.282

Finally, a third model (Model C and red curve on Figure 5) results from fusing the exceedance rates obtained from the fit of instrumental data (4≤M<6) with the exceedance rates of M≥6 earthquakes observed from the catalog. This model represents the upper bound of peak ground accelerations PGA that would be expected in this geological province. The mean recurrence times of earthquakes of a certain magnitude that may occur anywhere in the TMVB, according to the three seismicity models, are shown on Figure 6. Notice that for earthquakes M_w>5.5, the recurrence time between Model A, using only instrumental data, is substantially shorter than for the two other models, which include historical data.

Figure 6 Mean return periods of crustal earthquakes expected on the TMVB as a function of magnitude. These estimations are computed using the three different seismicity models discussed. For earthquakes M>5.5, the results from the traditional model that uses only instrumental data (green bars) begin to differ significantly with the results predicted by the semi parametric models (bars in blue and red), which also consider historical seismicity. 

Figure 7 Seismic hazard map of the TMVB for a return period of 500 years. Peak ground accelerations are estimated thinking the TMVB as a whole and assuming the semi-parametric seismicity model B as the preferred model describing seismicity rates in the volcanic arc. Such estimations might represent upper boundaries of PGA, however, they seem to fit the few observations available to date. 

Seismic hazard estimations of the TMVVB

We propose that model B is a conservative estimation of the seismicity in the TMVB. Using this model, we estimate the seismic hazard of the TMVB and a uniform seismic hazard (USH) map of the central portion of Mexico. In order to estimate the USH of the TMVB, three types of earthquake sources are considered: 1) shallow crustal seismicity occurring on the TMVB; 2) in-slab earthquakes occurring within the subducted Cocos plate and; 3) subduction zone earthquakes along the Middle American Trench (MAT).

Attenuation models

Due to the lack of strong motion records of crustal earthquakes in the TMVB, it was necessary to use Ground Motion Prediction Equations (GMPE) from other regions of the world with similar geological characteristics. We use the ground motion prediction equations reported by ^{Abrahamson et al., (2014)}, which describe the attenuation of the spectral response values in the horizontal and vertical components, for crustal earthquakes in tectonically active regions.

In the case of in-slab earthquakes, the parameters reported by ^{Ordaz et al., (2012)} and the GMPE’s proposed by ^{García et al., (2005)} were utilized. On the other hand, the seismicity rate in the Middle American Trench was determined using a traditional Gutenberg- Richter distribution. This was obtained compiling information of 288 thrust earthquakes reported in the Global Centroid Moment Tensor (CMT) database (Dziewonsky et al., 1981). Earthquakes from this catalog were selected from January 1^st, 1976 to December 31, 2014, within the geographic coordinates 15ºN ≥ latitude ≤ 24ºN latitudes and 106ºW ≤ longitude ≥ 92º W. All the events have magnitudes Mw ≥ 4.7 and focal depths of less than 40 km. In this case, the GMPE’s estimated by ^{Arroyo and Ordaz, (2010)}. for subduction zone earthquakes were utilized.

In computing ground motion expected for the return periods of 100, 500, 2500 and 10,000 years, we assumed a Poisson model to describe the probabilities of an earthquake occurrence. In other words, the seismicity rates on each seismic source are time-independent. The resulting seismicity parameters for the various seismic sources are summarized on Table 1.

Figure 8 Uniform hazard map of central Mexico for a return period of 500 years. This model is built considering subduction thrust-faulting earthquakes occurring in the Middle American Trench, normal-faulting events within the subducted Cocos plate and shallow crustal earthquakes in the TMVB. The seismicity parameters used to create this model are reported in Table 1. 

Table 1 Seismicity parameters used in this study to assess the seismic hazard of the TMVB. 

Discussion of results and conclusions

In this study, we propose that a non-traditional, semi-parametric seismicity model, which incorporates instrumental and historical data, represents better the seismicity of the TMVB and is better suited to assess the seismic hazard. For example, according to the seismicity Model B, the average recurrence time of a M≥7 earthquake on the TMVB is approximately 150 years. In contrast, the mean recurrence time estimated only on the instrumental catalog (Model A) is 12,000 years. The historical earthquake record of the past 160 years indicates that two earthquakes M≥7 have taken place on the TMVB. The mean recurrence time predicted by Model B of earthquakes M≥6.5 is of 66 years. The seismicity catalog shows that there are four events in this magnitude range in the last 160 years. Again, this is in good agreement with the predictions made by Model B and clearly shows that using only instrumental data, recurrence times of moderate and large earthquakes in the TMVB are considerably underestimated.

If we also consider the data available from paleo seismological studies, a similar picture emerges (Table 2). Studies performed in segments of crustal faults located on the TMVB indicate that the mean return periods of large earthquakes (M_w>7.0) are in the order of 3,600 to 12,500 years (Table 2). The return period predicted by instrumental data is 2,600 years. However, it should be stressed that the return periods estimated from paleo seismological studies refer only to short segments of the fault system in the TMVB. Clearly, considering the whole system of crustal faults observed in the TMVB, the return period predicted by Model A would be an underestimation.

Table 2 Paleoseismological studies carried out in segments of crustal faults located on the TMVB. 

Due to the epistemic uncertainties stemming from the limited data available, it was necessary to generate scenarios representing upper and lower bounds of the seismic hazard in the region. These models estimate PGA’s on hard soil of the TMVB for the more densely populated locations located on this geological province based on the uniform seismic hazard assessment presented above (Tables 3 and 4).

Table 3 Peak ground accelerations estimated for basin fill of the TMVB according to three different models describing its seismicity rate. 

Table 4 Uniform seismic hazard estimations expected in Mexico City, Guadalajara, Morelia and Jalapa, according to three seismicity models of the TMVB, for different return periods. 

It is worth to compare the results obtained here for the more important cities on the TMVB with the values of PGA proposed by ^{Ordaz (2004)} for a return period of 500 years (Table 5). The PGA values estimated by ^{Ordaz (2004)} are almost identical to those estimated from Model A in this study. This is not surprising, considering that both use a very similar instrumental catalog. However, when comparing with the results of Model B, the PGA’s estimated using only the instrumental seismicity are 18 to 56% smaller than those predicted in the model using the historical catalog. The largest differences are in the cities of Morelia and Guadalajara, due to the proximity of large historical earthquakes (Table 5).

Table 5 Comparison of the peak ground accelerations expected in Mexico City, Guadalajara, Morelia and Jalapa, for a return period of 500 years, according to ^{Ordaz (2004)} and this article. 

In conclusion, our results suggest strongly that the estimations of PGA’s may be severely underestimated using the limited instrumental catalog. Important differences exist in the uniform seismic hazard estimations with the incorporation of historical seismicity into the seismicity catalog. The seismic hazard determined with the use of the historical data, albeit only available for a short period, agree better with the seismic record and with the paleo seismological observations available. Clearly, this is a first approach to the important issue of estimating seismic hazard in this highly populated region of Mexico. It is evident that as more seismic, geological and geodetic data become available, hazard estimates should be improved. In particular, future work should consider including the local site response of the larger cities on the TMVB, particularly in Mexico City, in order for these results to be incorporated in the building codes of central Mexico.