Introduction

Intraslab seismicity defines the geometry of the subducting plates and sheds light on its thermal structure. In central Mexico, this seismicity ceases well before reaching the active volcanic front (^{Pardo and Suárez, 1995}). For example, a seismicity section from the middle America trench passing through the Nevado de Toluca volcano (N21°E), reveals that the intraslab earthquake activity ends at a depth of about 70 km, ~ 240 km from the trench and ~ 80 km from the volcano (Figure 1). A similar section from the trench through Popocatépetl shows that these events reach a depth of about 55 km and cease to occur ~ 240 km from the trench about 80 km from the volcano (Figure 1). Until now no intraslab earthquake has been reliably located whose epicenter lies near these active volcanoes. Here we report an intraslab earthquake which occurred on 19 May 2015 about 38 km west of the Nevado de Toluca at a depth of 98 km. The event provides first evidence from an earthquake that the Cocos plate is at a depth of ~ 100 km beneath the Nevado de Toluca volcano. The recordings also provide evidence of high attenuation of shear waves as they propagate through the mantle-wedge asthenosphere between the subducting Cocos plate and the continental crust.

Towns and cities in the Mexican altiplano are exposed seismic hazard from interplate earthquakes on the subduction thrust, crustal earthquakes in the Mexican Volcanic Belt (MVB), and intraslab earthquake in the subducted Cocos plate. The 19 September 2017 (M_w7.1) Morelos-Puebla earthquake is the most recent example of the intraslab event which destroyed many villages and towns in the epicentral region. It was the second most destructive earthquake in the history of Mexico City, next only to the 1985 (M _w8.0) interplate earthquake. The earthquake was located to the SSE of Mexico City at an epicentral distance of 114 km from Ciudad Universitaria (18.41 °N, 98.71 °W; depth H = 57 km). It is the closest, reliably located, intraslab earthquake (large or small) to Mexico City. One of many questions raised in the wake of the destruction left by the earthquake is whether a similar earthquake can occur at closer distance to Mexico City albeit at somewhat greater depth (^{Singh et al., 2018}). In this context, the detection and analysis of small intraslab earthquakes below the Central Mexican Volcanic Belt, such as the 19 May 2015 event, attains further importance.

Location and focal mechanism

Although the earthquake was small, it was reasonably well recorded by an accelerographic station in the epicentral region (TEJU, Figure 1 b), and many broadband stations of the Servicio Sismológico Nacional (SSN, Mexican National Seismological Service; ^{Pérez-Campos et al., 2018}; SSN, 2018) and the newly-completed Valley of Mexico network (^{Quintanar et al., 2018}). Figure 2 shows EW velocity seismograms at the closest 4 stations with clear signal. The accelerogram at TEJU was high-pass filtered at 1 Hz and then integrated to obtain the velocity trace shown in the figure. From the location of these stations (Figure 1 b) and the (S-P) times on the seismograms, it is clear that the event was relatively deep. We examined recordings at 31 stations within an epicentral distance of 360 km and found useful phase data at 19 stations. S-phase was read on the transverse component. The earthquake was located using 19 P-phase and 6 S-phase arrival times, and the crustal model used by the SSN. The depth, H, was varied between 0 and 200 km. A sharp minimum in the rms residual (0.31 s) occurs at H = 98 km. The final solution is: 19.092 °N, 100.131 °W; H = 98.2 km; origin time 03:06:53.9; azimuthal gap = 1580; and estimated errors in latitude, longitude, and H of 3.9 km, 3.3 km, and 3.3 km, respectively. The epicenter of the 19 May 2015 earthquake is 38 km to the west of Nevado de Toluca.

Figure 1 (a) Tectonic and location map (modified from ^{Singh et al., 2018}). Focal mechanisms and epicenters of earthquakes from Global CMT catalog with some exceptions. Event 1 is the intraslab Morelos-Puebla earthquake of 19 September 2017 (Mw7.1). Color of the focal mechanism is keyed to the depth of the event. Black stars: intraslab earthquakes studied here. Black dashed lines: depth contours of the Benioff zone as compiled by ^{Ferrari et al. (2012)}. Thick green dashed line defines the observed limit of moderate- and large-magnitude intraslab seismicity. Aftershock areas of large interplate earthquakes are shown in magenta color. Shaded area: Trans Mexican Volcanic Belt (^{Ferrari et al., 2012}). (b) Enlarged view of the rectangular area marked in (a), showing location of four stations whose data is analyzed in this study and the epicenters of 2015 (Mw3.3) and 2016 (M4.1) earthquakes. NT. V.: Nevado de Toluca volcano. P. V.: Popocatépetl volcano. (c) Section along line AA' shown in (a) and (b). The seismicity within ± 50 km of AA' and the stations are projected on the section. USL: ultra slow-velocity layer at the top of subducting oceanic crust shown in orange (^{Song et al., 2009}; ^{Pérez-Campos et al., 2008}). Brick color: mantle-wedge asthenosphere. 

Figure 2 EW velocity seismograms at the four closest stations with clear signal. Visual examination of the seismograms suggests a relative high-frequency deficiency at TEJU and MAVM as compared to ARIG and TLIG. 

The first motions could be read at 9 stations, which yields a normal-faulting focal mechanism (Figure 3): fault plane 1 with strike 316°, dip 46°, rake -90°; fault plane 2 with strike 116°, dip 46°, rake -90°. We note that the solution is not well constrained.

Figure 3 Normal-faulting focal mechanism which fits the first-motion data. 

The SSN catalog lists an earthquake on 10 November 2016 (M4.1), ~ 45 km south of the 2015 epicenter. We relocated this event following the same procedure as for the 2015 event. The event was located using 15 P-phase and 7 S-phase arrival times. The minimum rms residual (0.39 s) occurs at H = 68 km. The final solution is: 18.687 °N, 100.120 °W; H = 67.6 km; origin time 08:23:50.7; azimuthal gap = 660; and estimated errors in latitude, longitude, and H of 1.8 km, 2.5 km, and 4.1 km, respectively.

The earthquakes of 2016 and 2015 are 250 km and 290 km from the trench, respectively (Figure 1). The intraslab seismicity limit, defined by the Global CMT locations (shown in dashed green line in Figure 1) ends about 220 km from the trench. We note, however, that the teleseismic locations of Mexican earthquakes are systematically shifted by ~ 30 km to the NE (^{Singh and Lermo, 1985}; ^{Hjôrleifsdóttir et al., 2016}). If we take into account this mislocation, the intraslab seismicity limit defined by previous moderate and large earthquakes may be ~ 190 km from the trench. The events of 2015 and 2016 demonstrate that low-magnitude intraslab seismicity does persist beyond 190 km.

The focal depth of the 2016 earthquake agrees with the 60 km-depth contour of the Benioff zone. This is also true for the 19 September 2017 (M_w7.1; H = 57 km) earthquake (Figure 1). Generally, the hypocenters of well-located small/moderate intraslab earthquakes in central Mexico are less than 15 km below the plate interface (e.g., ^{Pacheco and Singh, 2010}). If we assume that the closest distance of the plate interface from the 2016 and 2015 events is 15 km then the depth of the slab increases by ~ 40 km within ~ 40 km, a dramatic change in dip over a short distance (Figure 1 c). Generally, the dip of T axis in subducted slabs coincides with the dip of the slab (see, e.g., ^{Pardo and Suárez, 1995} and Pacheco and Singh, 2009 for the Mexican subduction zone). Thus, the horizontal orientation of T-axis of the 2015 earthquake is not consistent with a steeply-dipping slab. We recall, however, that the focal mechanism is not well constrained. Clearly, locating more intraslab events in the region is critical to mapping the detailed geometry of the subducting Cocos plate as it changes its dip and plunges in the asthenosphere.

Source Spectrum and Q

The source acceleration spectrum of the earthquake, f2M˙O(f) was estimated from the analysis of the S wave at the four stations. The method is described elsewhere (^{Singh et al., 1999}); only a brief description is given here. The geometrical mean of Fourier acceleration spectral amplitude, A(f,R), of the two horizontal components of ground motion at a station may be written as:

where , M˙O(f), the source displacement spectrum or the moment rate spectrum, tends to M ₀ as f - 0, R = hypocentral distance, ρ = density in the focal region (taken as 3.2 g/cm³), β = shear-wave velocity in the focal region (taken as 4.68 km/s), and Q(f) = quality factor, which includes both anelastic absorption and scattering. Following ^{García et al. (2004)}, the geometrical spreading term, G(R), in Equation (1) was taken as 1 /R. Appropriate Q(f) for the earthquake is the one derived by ^{García et al. (2004)} for intraslab earthquakes in Mexico; i.e., Q(f)=251f ^0.58. However, the acceleration spectra, A(f), at the four stations, shown in Figure 4, suggest that Q for paths to ARIG and TLIG is very different from that for paths to TEJU and MAVM (Figure 1). To further appreciate this difference, we first used Q(f) = 251f^0.58 to correct spectra at all four stations and computed .The median source acceleration spëéti170m and ± one standard deviation curves are shown in Figure 5 a. We note that the standard deviation is extremely large at high frequencies (f>8 Hz), strongly suggesting that the same Q is not adequate for all paths.

Figure 4 Median horizontal Fourier spectra of S-wave at stations TEJU, MAVM, ARIG, and TLIG. Comparison of shapes and amplitudes of the spectra at higher frequencies suggest much greater attenuation of high-frequency S waves along the paths from source to TEJU and MAVM. 

Figure 5 Median source acceleration spectrum, f2M˙O(f), and ± one standard deviation curves. (a) Source spectrum using recordings from all four stations and Q(f) = 251f^0.58. (b) Source spectrum using recordings from ARIG and TLIG and Q(f) = 251f^0.58. (c) Source spectrum using recordings from TEJU and MAVM and Q(f) = 251f^0.58. (d) Source spectrum using recordings of TEJU and MAVM and Q = 350. In (a) to (d) the theoretical curve (continuous line) corresponds to Brune ω2-source model with M0 = 1.14x10¹⁴ Nm (Mw3.3) and fc = 8.4 Hz. 

Based on the similarity of the spectra in Figure 4, we then formed two groups. One group consisted of ARIG and TLIG and the other group included TEJU and MAVM. We computed the source acceleration spectrum for ARIG and TLIG group using Q(f) = 251f ^0.58 (Figure 5b). We interpreted the median spectrum by Brune a>²-source model (Brune, 1970) and obtained an estimate of seismic moment (M ₀ ) and corner frequency (f_c). The spectrum is well fit with M₀ = 1.14x10¹⁴ Nm (M_w3.3) and f _c = 8.4 Hz (Figure 5b). The stress drop (Δσ) computed using the Brune model (Brune, 1970) is 5.5 MPa. As expected, the theoretical spectrum for a Brune source model with M ₀ = 1.14x10¹⁴ Nm (M_w3.3) and f_c = 8.4 Hz does not fit the median spectrum of the TEJU and MAVM group at high frequencies (f > 5 Hz) if Q(f) = 251f ^0.58 is used for attenuation correction in equation 1 (Figure 5c).

If we assume that a>²-source model is appropriate for the earthquake, and M ₀ and f _c obtained for the ARIG and TLIG group correctly describes the source, then Q = 350 (1.0 < f < 20 Hz) for paths to TEJU and MAVM roughly fit the observed median spectrum (Figure 5 d). We note that this Q is 0.25 and 0.37 of Q = 251f^0.58 at f= 20 and 10 Hz, respectively. The two Qs are equal at f = 1.8 Hz. ^{Shapiro et al. (2000)} reported that shear waves traversing Popocatépetl volcano were highly attenuated. However, no such attenuation was discernible in the case of Nevado de Toluca. Thus, the relatively low Q (high attenuation) along the path to TEJU and MAVM is probably a consequence of shear-wave propagation through the mantle-wedge asthenosphere between the subducting Cocos plate and the continental crust. In any case, the paths to TEJU and MAVM for this earthquake do not cross the volcanic edifice (Figure 1).

Low Q in the mantle wedge in Mexico has been previously reported by ^{Singh et al. (2006)}. Q of the mantle wedge below Japan has been mapped in great detail (see e.g., ^{Takanami et al., 2000}; ^{Tsumura et al., 2000}), thanks to high density of seismic stations and their favorable location with respect to intraslab seismicity. Low Q of the mantle wedge implies large attenuation of seismic waves reaching the back arc. An impressive demonstration of this comes from the extensive recordings of Japanese inslab earthquakes. A rather low Q for S waves (-60) in the mantle wedge is required to explain recorded seismic waves on the western seaboard of Japan, across the volcanic front, as compared to eastern Japan where the high-frequency seismic waves propagate through the subducted slab and the continental lithosphere with much less attenuation.

Conclusions

Intraslab earthquakes are linked to metamorphic dehydration reactions in the minerals of subducting oceanic crust and upper mantle (e.g., ^{Kirby et al., 1996}; ^{Hacker et al., 2003}). Thus, a lack of intraslab seismicity in the "warm" subducting Cocos plate (age near the trench is ~5.6 my) well before the volcanic front in central Mexico may be attributed to release of fluid mostly occurring before the slab reaches this distance. The existence of active volcanoes suggests some fluid release which causes partial melting of the mantle-wedge asthenosphere above and results in the volcanic activity. The slab, however, may be too hot (T > 600 °C; ^{Manea et al. , 2004}) for brittle rupture and, hence, the lack of seismicity. In this context, the detection of the 2015 intraslab is important as it provides the first evidence from an earthquake that the plate is at a depth of ~ 100 km in central Mexico. It was a normal-faulting earthquake even though the focal mechanism is not well constrained. The seismograms of the 2015 earthquake also provide evidence of high attenuation of shear waves propagating, partially, in the mantle-wedge asthenosphere.

Perhaps the intraslab seismicity beneath active volcanoes in central Mexico at M ~3 level is abundant but has not been documented in the past because of poor instrumentation. New seismic stations in the region should provide high-quality recordings of such events and lead to firm support to the preliminary conclusions reached in this study that are based on the analysis of just one event. Mapping of small intraslab earthquakes below the MVB is important to (1) accurately delineate the geometry of the Benioff zone, (2) understand the likelihood of moderate and large intraslab earthquakes beneath the Mexican Volcanic Belt, (3) estimate ground motions from postulated intraslab earthquakes, and (4) assess impact of intraslab earthquakes on the seismic hazard to cities and town located in the region.