<figure><img src="image_1.png" /></figure>

ARTICLE

https://doi.org/10.1038/s41467-019-10605-2

OPEN

# The mechanism of tidal triggering of earthquakes at mid-ocean ridges

Christopher H. Scholz ID <fnref n="1" />,<fnref n="3" />, Yen Joe Tan ID <fnref n="1" />,<fnref n="3" /> & Fabien Albino ID ²

The strong tidal triggering of mid-ocean ridge earthquakes has remained unexplained because the earthquakes occur preferentially during low tide, when normal faulting earthquakes should be inhibited. Using Axial Volcano on the Juan de Fuca ridge as an example, we show that the axial magma chamber inflates/deflates in response to tidal stresses, producing Coulomb stresses on the faults that are opposite in sign to those produced by the tides. When the magma chamber's bulk modulus is sufficiently low, the phase of tidal triggering is inverted. We find that the stress dependence of seismicity rate conforms to triggering theory over the entire tidal stress range. There is no triggering stress threshold and stress shadowing is just a continuous function of stress decrease. We find the viscous friction parameter A to be an order of magnitude smaller than laboratory measurements. The high tidal sensitivity at Axial Volcano results from the shallow earthquake depths.

<footnote n="1">¹ Lamont-Doherty Earth Observatory, Columbia University, Palisades, NY 10964, USA. ² School of Earth Sciences, University of Bristol, Bristol B58 1RJ, UK.</footnote>

<footnote n="3">³These authors contributed equally: Christopher H. Scholz, Yen Joe Tan. Correspondence and requests for materials should be addressed to C.H.S. (email: scholz@ldeo.columbia.edu) or to Y.J.T. (email: yjt@ldeo.columbia.edu)</footnote>
Tidal triggering of earthquakes has been found to be elusive despite a long search¹⁻³. On the continents, the signal is so weak that a significant statistical correlation between earthquakes and tides can be detected only with very large catalogues⁴. In the oceans, where loading from ocean tides result in tidal stresses an order of magnitude larger than the solid earth tides, the tidal triggering signal has been found to be stronger. Cochran et al.⁵ used a large catalog of oceanic earthquakes to show that shallow thrust earthquakes may be found to correlate with maximum tidally generated Coulomb stresses when the tides are large enough. The earthquakes preferentially occur at low water when the normal stresses on low-angle thrust faults are reduced such that they become unclamped. Much stronger tidal triggering has been observed with ocean bottom seismometer networks in magmatic areas at mid-ocean ridges⁶⁻¹⁰. These are the most promising places to test theories of earthquake triggering. In these cases, however, even the most basic mechanism of the triggering is not understood because seismicity peaks at low tides, when it should be inhibited on normal faults. The most well-studied of these cases is at Axial Volcano on the Juan de Fuca ridge. We shall study this case, and later see if the results obtained there can also be applied to the other mid-ocean ridges.

At Axial Volcano, the rate of normal faulting earthquakes is maximum at low tide. We show that they are driven by the tidally induced inflation of the magma chamber. The seismicity rate is modulated by tidal stresses at all tidal phases in agreement with triggering models based on nucleation theory. These results show that the viscous friction parameter A must be much smaller than indicated in laboratory experiments.

## Results

**Tidal triggering of earthquake at Axial Volcano.** Axial Volcano, which is at the intersection of a mid-ocean ridge with a hotspot (Fig. 1), erupts on a decadal time scale. Each eruption is followed by a caldera collapse accompanied by thrusting on outwardly dipping ring faults, followed by a re-inflation period, at the latter stages of which the ring faults become reactivated in normal faulting¹¹⁻¹³. The best observations of tidal triggering were for the normal faulting earthquakes in the months prior to the 2015 eruption⁷.

At Axial Volcano, the ocean tides are very large (3 m) so that ocean loading dominates the solid earth tides and the vertical tidal stress dominates and is in phase with the ocean tides (Supplementary Fig. 1), so we need only consider the vertical

<figure><img src="image_2.png" /><figcaption>**Fig. 1** Cross-section of seismicity at Axial Volcano. East-west cross-section of seismicity in 3 months preceding the 2015 eruption at Axial Volcano. Red curve is the roof of the axial magma chamber¹⁵. The bathymetry is from a compilation of latest cruises, the most recent in 2010⁶². Inset, location map for Axial Volcano</figcaption></figure>

component in our analysis. Tension is taken as positive for tidal stresses, so the maximum tidal stress corresponds to the minimum water depth. To avoid ambiguity, in this paper we will refer to high and low tides in the conventional way as high and low water, recalling that low water produces tension and high water produces compression.

Figure 1 shows a cross-section view of the seismicity for the three months prior to the 2015 eruption, which illuminates the ring faults. This data set contains ~60,000 earthquakes with a magnitude of completeness ($M_c$) = 0.1¹⁴. Figure 2 shows a histogram of the seismicity plotted as a function of tidal phase, in which 0° is the maximum low tide. The correlation is obvious and requires no statistical treatment. It was first proposed that this was also a case of fault unclamping⁶,⁸,¹⁰, but when it was established that these earthquakes were dominated by normal faulting¹² this viewpoint became untenable. Both the seismicity trends in Fig. 1 and the focal mechanisms¹² indicate a mean fault dip of 67°. The reduction of vertical stress brought about by low tide will produce a Coulomb-stress change on such steeply dipping normal faults that inhibits their slip. It is, rather, the high tides that will produce a Coulomb stress on the faults that encourages slip. This seeming paradox is resolved by including the effect of the axial magma chamber on the distribution of stress.

**The response of the magma chamber.** The red curve in Fig. 1 delineates the roof of the axial magma chamber obtained from seismic imaging¹⁵. Inflation of the magma chamber drives the normal faulting on the ring faults. This is demonstrated in Fig. 3, where we show east-west cross-sections of a 3D model containing a magma chamber with dimensions defined by seismic imaging¹⁵,¹⁶. Figure 3a shows the Coulomb failure stress change, $\Delta CFS = \Delta\tau + \mu\Delta\sigma$, on 67° dipping faults that results from a magma chamber overpressure of 1 MPa ($\Delta\tau$ is the change in shear stress resolved on the fault in the slip direction, $\Delta\sigma$ is the change in normal stress on the fault plane, and $\mu$ is the friction coefficient). Positive (red) $\Delta CFS$ values encourage normal fault slip, negative ones (blue) inhibit it. The primary features in Fig. 3a are the zones of positive $\Delta CFS$ that correspond to the seismicity shown in Fig. 1. See Methods for details about the model.

Because the magma chamber is a soft inclusion, its presence will profoundly affect the stress field in its vicinity resulting from any external load. We simulate the response to tides by

<figure><img src="image_3.png" /><figcaption>**Fig. 2** Seismicity versus tidal phase. Histogram of earthquakes plotted vs. the phase of the vertical component of the tidal stress, in which 0° is the peak stress (tension is positive), which corresponds to the lowest ocean height</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>**Fig. 3** Coulomb stress on normal faults above the magma chamber. Distribution of Coulomb stress changes on 67° dipping normal faults near the axial magma chamber. Positive values favor fault slip, negative values inhibit it. **a** For an overpressure of 1 MPa within the magma chamber with friction coefficient of the faults $\mu = 0.8$. **b** For a decrease in vertical stress equivalent to a reduction in water level of 1 m. In **b** an effective friction $\mu' = 0.4$ is used. The bulk modulus of the rock is assumed to be $K_r = 55$ GPa, and the bulk modulus of the magma chamber is assumed to be $K_m = 1$ GPa. The heavy line in Fig. 3b indicates the fault upon which $\Delta$CFS was measured</figcaption></figure>

calculating the distribution of $\Delta$CFS on 67° dipping faults resulting from a reduction in vertical stress corresponding to a 1 m drop in the ocean tide. This is shown in Fig. 3b. The pattern is very similar to that of Fig. 3a, demonstrating how a low tide can stimulate activity on these faults. This pattern arises because the reduction of vertical stress causes the magma chamber, owing to its higher compressibility, to inflate relative to the surrounding rock, which produces a stress field congruent with that of Fig. 3a. This is superimposed on a uniform $\Delta$CFS from the tidal stress, which is negative in the case of a low tide. Likewise, high tides cause the magma chamber to deflate, which also produces Coulomb stresses opposite in sign to the tidal ones. Which component is larger determines whether normal faulting earthquakes are stimulated by the low tide or the high tide.

The relative expansion of the magma chamber depends inversely with $K_m/K_r$, the bulk modulus of the magma chamber relative to that of the surrounding rock, so this is the critical parameter that determines the behavior of the system. In the calculation of Fig. 3a we used $\mu = 0.8$. The lack of a phase shift between the tidal stress and the seismicity (Fig. 2) indicates that tidal loading must be under undrained conditions (see Supplementary Fig. 2 for a check on this assumption). Therefore, for calculations such as shown in Fig. 3b we use an effective friction $\mu' = (1 - B)\mu$, where B is Skempton's coefficient. The proper value to use for B is not well established: plausible values are between 0.5 and $1^{17}$, so we explored values of $\mu'$ from 0.4 to 0. In Fig. 3b we used $\mu' = 0.4$, $K_m = 1$ GPa and $K_r = 55$ GPa.

The ratio $K_m/K_r$ and $\mu'$ control the $\Delta$CFS on overlying normal faults that results from an applied vertical tidal stress. We use a metric, $\chi$, which normalizes the $\Delta$CFS by the vertical tidal stress and which, therefore, is time invariant. We define $\chi$ by the $\Delta$CFS on a 67° dipping fault, indicated by the bold line in Fig. 3b, averaged from the corner of the magma chamber to the surface, normalized by the applied vertical tidal stress. This parameter is shown in Fig. 4 as a function of $K_m$ for several values of $\mu'$ at a fixed $K_r = 55$ GPa. Positive $\chi$ values indicate that normal faulting

<figure><img src="image_5.png" /><figcaption>**Fig. 4** Properties of the magma chamber deformation system. The vertical axis $\chi$ is the average $\Delta$CFS on a 67° dipping normal fault from the tip of the magma chamber to the surface (bold line in Fig. 3b) resulting from an applied vertical tidal stress, normalized by that stress. The red area defines the conditions in which low tides encourage seismicity on normal faults and high tides discourage it, and the blue area vice versa. The plus sign indicates the conditions shown in Fig. 3b</figcaption></figure>

earthquakes will be favored by low tides, negative values by high tides and vice versa for thrust earthquakes. All conditions in Fig. 4 within the red region therefore favor normal faulting earthquakes on the low tide and inhibit them on the high tide, and within the blue region, vice versa. Normal faulting may be generated by low tides for values of $K_m < 8$ GPa, depending on the value of $\mu'$. The point indicated by the plus in Fig. 4 is the case illustrated in Fig. 3b. The bulk modulus of gas-free mid-ocean ridge (MOR) magma is 12 GPa$^{18}$, but at the pressure of the magma chamber
(~40 MPa) this value can be reduced by one to two orders of magnitude by the presence of volatiles¹⁹. Thus, at this pressure, a magma of $K_m = 1$ GPa would contain 2650 ppm CO₂ by weight¹⁸. This is greater than the highest values typically seen for CO₂ content of MOR magma²⁰, but this difference could easily be accounted for by the inclusion of exsolved H₂O. So, the values of $K_m$ that we find would promote normal faulting at low tide are realistic. Our illustrative example (Fig. 3b) indicates $\chi = 0.176$, a figure that will enter into the modeling calculations of the triggered seismicity in the next section.

**Comparison of the earthquake triggering with theory.** During the three month period prior to the eruption, the seismic moment release on the eastern ring fault indicated a slip magnitude approximately the same as observed geodetically¹², indicating that the faults are seismically coupled. In this case, there are two models that relate change in seismicity rate to a rapid change in driving stress. These are based on earthquake nucleation models²¹, one derived from the rate-state friction law²² and the other from subcritical crack growth due to stress corrosion²³,²⁴. The rate-state friction version is

$$ \frac{R}{r} = \exp\left[\frac{\Delta\text{CFS}}{A\sigma}\right], \qquad (1) $$

and the stress corrosion version is

$$ \frac{R}{r} = \left[ 1 + \left( \frac{\Delta \text{CFS}}{\Delta \tau} \right)^n \right], \qquad (2) $$

where R is the instantaneous seismicity rate, r is the background rate, here taken as the rate when the tidal stress is zero, and ΔCFS = $\chi\sigma_v$, $\sigma_v$ being the vertical tidal stress. The control parameters for the rate-state friction version are the effective normal stress $\sigma$ and the viscous friction term A. For the stress corrosion version, they are the stress corrosion index n and the earthquake stress drop $\Delta\tau$.

The fit of these equations to the data is shown in Fig. 5, where the solid blue curve and the dashed red curves are Eqs. (1) and (2), respectively. These two formulations cannot be distinguished and fit the data equally well. There is no detectable phase shift between the seismicity and the tides (Fig. 2), nor is there any hysteresis observed—data for rising and falling stresses fit the triggering curves equally well (Supplementary Figs. 4 and 5). We,

<figure><img src="image_6.png" /><figcaption>**Fig. 5** Normalized seismicity rate change vs. change in Coulomb stress. Coulomb stress is converted from tidal vertical stress using $\chi = 0.176$ (from the state indicated by the plus sign in Fig. 4). Blue curve is the rate-state friction version and the red curve is the stress corrosion version</figcaption></figure>

therefore, conclude that poroelastic relaxation is negligible in the response to the semi-diurnal tides.

The agreement with the theories is excellent, and extends them to far smaller stresses than previously seen⁵, even into the negative stress regime. In the case shown in Fig. 5, the value of $\chi$ used was 0.176, from the illustrative example. The goodness of fit to the theories does not depend on the value of $\chi$ obtained from the deformation model: that merely determines the scale of the stress axis. The various implications of this result will be deferred to the discussion section.

**Applications to other areas.** Wilcock⁶ searched for tidal triggering on the mid-ocean ridge systems of the NE Pacific, using mainly land-based networks. He found a 15% increase in seismicity within 15° of the lowest tides. The focal mechanisms of the earthquakes, however, were undetermined. With an OBS deployment on the Endeavour segment of the Juan de Fuca ridge, some 2° NE of Axial Volcano, the correlation of seismicity with low tides became much better defined⁹. Most of the triggered seismicity there was near the ridge axis, where the focal mechanisms indicate normal faulting²⁵. This situation is therefore quite similar to Axial Volcano and the same effect of the magma chamber seems necessary to explain these observations.

At the hydrothermal field at 9°50'N on the East Pacific rise, an OBS deployment also showed evidence for tidal triggering⁸. There the ocean tides are much smaller than at Axial Volcano and a significant contribution to tidal stresses is from the solid earth tides. The seismicity maximum correlates with the maximum extensional tidal stress, which can reach 1.3 kPa. The dependence of the seismicity on stress is similar to that observed at Axial Seamount (compare Fig. 3c in ref. ⁸ to our Fig. 5). Evidence for the mechanism of the earthquakes is equivocal: scant focal mechanism data has indicated strike-slip, normal faulting, and reverse faulting²⁶,²⁷, and others have proposed that the seismicity is due to hydrothermally induced extension cracking²⁸. There is also a variation in the tidal phase angle of earthquakes along the strike of the ridge axis. This indicates the earthquake triggering is also modulated by pore pressure changes brought about by hydrothermal circulation²⁹. With this degree of ambiguity and complexity, we cannot assess how deformation of the magma chamber may be related to the tidal triggering in this location.

The unloading model used here was initially tested at Katla volcano (Iceland), where earthquakes show an annual cycle with the maximum seismicity rate occurring in the late summer³⁰ when the snow cover of the glacier above the volcano is minimum (annual fluctuation—6 m). The model³¹ showed that this was also the period of maximum Coulomb stresses in the area above the magma chamber. However, in this case the focal mechanisms were not known³⁰ so it was not possible to determine if the system was in the red or blue regions of Fig. 4.

## Discussion

Our observations of seismicity rate change as a function of stress, shown in Fig. 5, show a much stronger agreement with the models and extends to much smaller stresses than previously observed⁵. In the rate-state friction version, the representative value $\chi = 0.176$ yields $A\sigma = 2$ kPa, about an order of magnitude smaller than found in earlier studies of earthquake triggering⁵,³²,³³. In the earlier studies the earthquakes were deeper (8–20 km), so the difference could be from that factor alone. In those papers, to accommodate lab values for A of 0.003–0.007, near-lithostatic pore pressures were assumed to get low enough values of $\sigma$ to match the observed $A\sigma$. At Axial Volcano the normal stress is 7.2 MPa at the average earthquake depth of 1.2 km, assuming $\mu = 0.8$, a hydrostatic pore pressure gradient,
and a dip of 67°. It is highly unlikely that overpressures can be maintained in the top 1 km of very young oceanic crust where there is no sediment cover and there is vigorous hydrothermal circulation throughout the caldera³⁴,³⁵. Using 7.2 MPa for σ we conclude that $A = 0.0003$, much smaller than lab values. Considering the entire spread of the solution space for $\mu' \le 0.4$, $\chi$ ranges from 0.2 to zero, and with $\sigma = 7.2$ MPa, the equivalent range of $A$ is $0.0004 \ge A > 0$. If hydrostatic pore pressure was assumed in the other studies, estimates for $A$ in that range would also be obtained. The higher sensitivity to triggering at Axial Volcano is due to the shallow depths of earthquakes there and its observation made clear by the high rate of background seismicity during the inflation stage. Considering the Vidale et al.³ study, which failed to detect a tidal correlation for Southern California earthquakes, if we use 8 km depth with a hydrostatic pore pressure gradient, their typical tidal stress level of 1 kPa, and our $A$ value of 0.0003, we calculate³⁶ that they would need 22k events to detect a correlation, whereas their catalog contained only 13k. Thus, the other studies are consistent with our finding that the $A$ parameter at geologic rates must be considerably smaller than lab values. The few laboratory studies that explore if the friction rate parameters $A$ and $B$ depend on sliding rate³⁷,³⁸ find that they do, and experiments at plate tectonic slip rates³⁹ indicate that the friction parameters at those rates may differ significantly from those measured at the much higher rates usually employed in laboratory experiments.

Beeler and Lockner³⁶ noted that there are two triggering regimes: a threshold regime, in which the earthquake nucleation time $t_n$ is shorter than the tidal period and a nucleation regime, in which it is longer. In the former, maximum seismicity rate would correlate with the maximum stressing rate, in the latter, with maximum stress amplitude. Our data clearly confirm the latter (Fig. 2), and the latter is also implicit in the fit in Fig. 5. The uplift rate prior to the 2015 eruption was 61 cm/y¹³. From our inflation model (e.g., Fig. 3a) we found that the corresponding fault stressing rate $\dot{\tau}$ is 5 MPa/y. Using³⁶ $t_n = \frac{2\pi\Delta\sigma}{\dot{\tau}}$, we get $t_n = 24$ h confirming that the system is indeed in the nucleation regime for the semi-diurnal tides examined.

For the stress corrosion version of the triggering equation, if we take the stress corrosion index to be the laboratory values for basalt, $22 < n < 44$⁴⁰, then the best fitting stress drop would be $0.04 < \Delta\tau < 0.09$ MPa. This is a bit lower than the $0.18 < \Delta\tau < 2.8$ range⁴¹ for earthquakes at 1 km depth in Southern California, although these estimates are from mainly strike-slip earthquakes, which have systematically higher stress drops than normal faults⁴². If we take the rule that stress drop is about 3% of the shear strength⁴², then for strength $\tau = \mu\sigma = 5.7$ MPa we get $\Delta\tau = 0.17$ MPa, not much greater than our fit values. Thus, for this version of the triggering law, we do not have any serious conflict with independent estimates.

Thresholds for static or dynamic triggering have been much discussed⁴³⁻⁴⁵. Van der Elst and Brodsky¹⁶ showed that dynamic triggering could be detected at very small strains, and suggested that the lower limit may simply be a matter of detectability. Our results (Fig. 5) show that seismicity rate falls smoothly as the tidal stress falls to zero, indicating that there is no threshold for triggering. Seismicity rate continues to fall when the tidal Coulomb stress becomes negative, indicating that what is often called stress shadowing is a continuous quantifiable function of stress reduction.

It has often been remarked that hydrothermal areas seem particularly susceptible to dynamic triggering from distant earthquakes⁴⁷⁻⁴⁹. Attempts to explain this have invoked various effects of dynamic stresses on the permeability and/or pressure of the pore fluid⁵⁰⁻⁵². The excellent agreement of our data with the dry triggering models indicates that additional mechanisms are

not required to explain the tidal triggering at Axial Volcano. In the case of tidal triggering, some of those proposed mechanisms, such as unclogging of fluid pathways, are less likely because the tides are continually jostling the faults so that clogs, such as from mineralization as suggested in the Yellowstone case⁴⁸, will not have time to form.

Low-frequency earthquakes in the tremor deep in subduction zones also show high sensitivity to tidal stresses⁵³,⁵⁴. This undoubtedly implies very low effective normal stresses, although if our $A$ value is used instead of the laboratory values, effective stresses would be about an order of magnitude larger than those reported. This is more in line with the values of 2–4 MPa obtained from the scaling of the amplitude and frequency of the periodic slow slip events which the tremor accompanies⁵⁵.

The value of $A$ we obtained is the first estimate of a rate-state friction parameter for earthquakes at tectonic loading rates. As shown in the examples given above, using this value provides much more plausible estimates of effective normal stresses in those cases. This finding shows that laboratory values of friction parameters should not be blithely adopted in theoretical studies of earthquake phenomena.

## Methods

**Coulomb stress modeling.** Coulomb stress calculation is performed with the commercial Finite Element Modelling software COMSOL Multiphysics® (https://www.comsol.com). We use a 100 × 100 × 50 km domain designed to limit boundary effects. Boundary conditions are zero-displacement for the bottom and lateral boundaries and free-displacement for the top boundary corresponding to the Earth's surface. For the host rock, we assume an isotropic and homogeneous elastic medium with a bulk modulus $K_r$ of 55 GPa and a Poisson's ratio $\nu_r$ of 0.25, which is in accordance with seismic velocities recorded on the East Pacific Rise⁵⁶. At Axial Seamount, multichannel seismic-reflection has inferred a 14-km long by 3-km-wide shallow magma reservoir located at 1.1–2.3 km depth¹⁵,¹⁶. We, therefore, model the magma reservoir as a 3D ellipsoid cavity with semi-axis: $a = 7$ km, $b = 1.5$ km, and $c = 0.5$ km, and top depth located at 2 km below the surface.

In our modeling, the initial stress field is lithostatic and stress perturbations are calculated considering two scenarios: first, from the pressurization of the magma reservoir and second, for the effect of ocean tides. For the first scenario, the overpressure inside the reservoir is modeled by applying a constant normal stress applied at the boundary of the ellipsoid. For the second scenario, the stress changes due to ocean low tides are modeled by applying a boundary load at the surface corresponding to a 1 m decrease in the water level. Surface unloading causes the reservoir expansion resulting in a magma pressure change, which depends on the reservoir volume, the bulk modulus of the magma and the elastic properties of the host rock. The pressure change is applied on the reservoir's wall considering different bulk modulus $K_m$ from 0 to 20 GPa. For each model, the Coulomb failure stress change is calculated on specific fault planes using $\Delta CFS = \Delta\tau + \mu\Delta\sigma$, where $\Delta\sigma$ is the normal stress change, $\Delta\tau$ the tangential stress changes and $\mu$ the friction coefficient.

**Seismicity catalog.** The cabled seismic network started streaming time-corrected data in late January 2015, 3 months before Axial Volcano erupted on April 24th, 2015 (Wilcock et al., 2016)⁷. In this study, we examine the ~60,000 earthquakes located between January 22nd and April 23rd, 2015. The earthquake catalog is from Wilcock et al.⁷ and is archived in the Interdisciplinary Earth Data Alliance Marine Geoscience Data System (DOI: 10.1594/IEDA/323843)⁵⁷ Tan et al.¹⁴ estimated the $M_c$ of the catalog and find $M_c = 0.1$ when using the goodness-of-fit method⁵⁸ and $M_c = 0.3$ based on b-value stability⁵⁹.

**Tidal stress calculations.** The tidal stresses were modeled following Tan et al.¹⁴. Briefly, we estimate the horizontal strains at the earthquake region (−130.009, 45.955) due to solid earth tide and regional ocean tidal loading with the SPOTL package⁶⁰ using the predicted tidal height in the Oregon State University regional ocean tidal model⁶¹. We then convert strain to stress using p-wave velocity of 5.55 km/s, s-wave velocity of 3.20 km/s, and a density of 2800 kg/m³, assuming plane stress. We also calculate the vertical stress variations due to near-field direct ocean tidal loading from the predicted tidal height taking seawater density of 1030 kg/m³ and gravitational acceleration of 9.8 m/s². We then estimate the horizontal stresses from the vertical stress assuming uniaxial strain. An example is shown in Supplementary Fig. 1.

The ratio between the horizontal and vertical tidal stresses is symmetrically distributed with a median value of 0.058 (Supplementary Fig. 4). Therefore, we use this ratio in the numerical modeling to obtain the $\chi$ value which relates the tidal vertical stress to the average Coulomb stress change on the normal fault (Fig. 4). This $\chi$ value (0.176 for $\mu' = 0.4$ and $K_m = 1$ GPa) is then used to convert the
tidal-vertical-stress time series to the Coulomb-stress time series for tidal triggering analysis (Fig. 5). This approximation gives almost identical results as when we use two separate $\chi$ values, $\chi_1 = 0.133$ relating the tidal vertical stress to the Coulomb stress change on the fault (with horizontal stress = 0) and $\chi_2 = 0.740$ relating the tidal horizontal stress to the Coulomb stress change on the fault (with vertical stress = 0), to convert the different component tidal-stress time series to the Coulomb-stress time series for tidal triggering analysis (Supplementary Fig. 5).

## Data availability

The earthquake catalog is from Wilcock et al.⁷ and is archived in the Interdisciplinary Earth Data Alliance Marine Geoscience Data System (DOI: 10.1594/IEDA/323843).

## Code availability

Numerical models for Coulomb stress calculation have been performed using the software COMSOL Multiphysics® software [https://www.comsol.com] and source files are available from the co-author [F.A.] upon reasonable request.

Received: 18 December 2018 Accepted: 17 May 2019

Published online: 07 June 2019

## References

1. Schuster, A. On lunar and solar periodicities of earthquakes. *Proc. R. Soc. Lond.* **61**, 455–465 (1897).

2. Emter, C. in *Tidal Phenomena, Lecture Notes in Earth Sciences*, Vol. 66 (eds H. Wilhelm, W. Zurn, & H. G. Wenzel) 293–310 (Springer, Berlin, 1997).

3. Vidale, J., Agnew, D., Johnston, M. & Oppenheimer, D. Absence of earthquake correlation with earth tides: an indication of high preseismic fault stress rate. *J. Geophys. Res.* **103**, 567–524, 572 (1998).

4. Metivier, L. et al. Evidence of earthquake triggering by the solid earth tides. *Earth Planet. Sci. Lett.* **278**, 370–375 (2009).

5. Cochran, E. S., Vidale, J. E. & Tanaka, S. Earth tides can trigger shallow thrust fault earthquakes. *Science* **306**, 1164–1166 (2004).

6. Wilcock, W. S. D. Tidal triggering of earthquakes in the Northeast Pacific Ocean. *Geophys. J. Int.* **179**, 1055–1070 (2009).

7. Wilcock, W. S. D. et al. Seismic constraints on caldera dynamics from the 2015 Axial Seamount eruption. *Science* **354**, 1395–1399 (2016).

8. Stroup, D. F., Bohnenstiehl, D. R., Tolstoy, M., Waldhauser, F. & Weekly, R. T. Pulse of the seafloor: Tidal triggering of microearthquakes at 9 degrees 50'N East Pacific Rise. *Geophys. Res. Lett.* **34**, L15301 https://doi.org/10.1029/2007gl030088 (2007).

9. Wilcock, W. S. D. Tidal triggering of micro earthquakes on the Juan de Fuca Ridge. *Geophys. Res. Lett.* **28**, 3999–4002 (2001).

10. Tolstoy, M., Vernon, F. L., Orcutt, J. A. & Wyatt, F. K. Breathing of the seafloor: tidal correlations of seismicity at axial volcano. *Geology* **30**, 503–506 (2002).

11. Chadwick, W. W. Jr., Nooner, S. L., Butterfield, D. A. & Lilley, M. D. Seafloor deformation and forecasts of the April 2011 eruption at Axial Seamount. *Nat. Geosci.* **5**, 474–477 (2012).

12. Levy, S. et al. Mechanics of fault reactivation before, during, and after the 2015 eruption of Axial Seamount. *Geology* **46**, 447–450 (2018).

13. Nooner, S. L. & Chadwick, W. W. Inflation-predictable behavior and co-eruption deformation at Axial Seamount. *Science* **354**, 1399–1403 (2016).

14. Tan, Y. J., Waldhauser, F., Tolstoy, M. & Wilcock, W. S. D. Axial Seamount: periodic tidal loading reveals stress dependence of the earthquake size distribution (b value). *Earth Planet. Sci. Lett.* **512**, 39–45 (2019).

15. Arnulf, A. F. et al. Anatomy of an active submarine volcano. *Geology* **42**, 655–658 (2014).

16. Arnulf, A. F., Harding, A. J., Kent, G. M. & Wilcock, W. S. D. Structure, seismicity, and accretionary processes at the hot spot-influenced Axial Seamount on the Juan de Fuca Ridge. *J. Geophys. Res.* **123**, 4618–4646 (2018).

17. Cocco, M. & Rice, J. R. Pore pressure and poroelasticity effects in Coulomb stress analysis of earthquake interactions. *J. Geophys. Res.* **107** https://doi.org/10.1029/2000JB000138 (2002).

18. Johnson, D. J., Sigmundsson, F. & Delaney, P. T. Comment on “Volume of magma accumulation or withdrawal estimated from surface uplift or subsidence with application to the 1960 collapse of Kilauea volcano” by P. T. Delaney and D. F. McTigue. *Bull. Volcanol.* **61**, 491–493 (2000).

19. Huppert, H. E. & Woods, A. W. The role of volatiles in magma chamber dynamics. *Nature* **420**, 493–495 (2002).

20. Le Voyer, M., Kelley, K. A., Cottrell, E. & Hauri, E. H. Heterogeneity in mantle carbon content from CO²-undersaturated basalts. *Nat. Commun.* **8**, 14062 https://doi.org/10.1038/ncomms14062 (2017).

21. Kanamori, H. & Brodsky, E. E. The physics of earthquakes. *Rep. Prog. Phys.* **67**, 1429–1496 (2004).

22. Dieterich, J. H. A constitutive law for rate of earthquake production and its application to earthquake clustering. *J. Geophys. Res.* **99**, 2601–2618 (1994).

23. Das, S. & Scholz, C. Theory of time-dependent rupture in the earth. *J. Geophys. Res.* **86**, 6039–6051 (1981).

24. Main, I. Applicability of time to failure analysis to accelerated strain before earthquakes and volcanic eruptions. *Geophys. J. Int.* **139**, F1–F6 (1999).

25. Wilcock, W. S. D., Archer, S. D. & Purdy, G. M. Microearthquakes on the endeavour segment of the Juan de Fuca Ridge. *J. Geophys. Res.* **107**, EPM 4-1-EPM 4-21 https://doi.org/10.1029/2001jb000505 (2002).

26. Sohn, R. A., Hildebrand, J. A. & Webb, S. C. A microearthquake survey of the high-temperature vent fields on the volcanically active East Pacific Rise (9 degrees 50'N). *J. Geophys. Res.* **104**, 25367–25377 (1999).

27. Waldhauser, F. & Tolstoy, M. Seismogenic structure and processes associated with magma inflation and hydrothermal circulation beneath the East Pacific Rise at 9 degrees 50' N. *Geochem. Geophys. Geosyst.* **12**, ISSN 1525-2027 https://doi.org/10.1029/2011gc003568 (2011).

28. Tolstoy, M., Waldhauser, F., Bohnenstiehl, D. R., Weekly, R. T. & Kim, W. Y. Seismic identification of along-axis hydrothermal flow on the East Pacific Rise. *Nature* **451**, 181–U187 (2008).

29. Stroup, D. F. et al. Systematic along-axis tidal triggering of microearthquakes observed at 9 degrees 50' N East Pacific Rise. *Geophys. Res. Lett.* **36**, L18302 https://doi.org/10.1029/2009gl039493 (2009).

30. Einarsson, P. & Brandsdottir, B. Earthquakes in the Myrdalsjokull area, Iceland, 1978–1985: Seasonal correlation and connection with volcanoes. *Jokull* **49**, 59–73 (2000).

31. Albino, F., Pinel, V. & Sigmundsson, F. Influence of surface load variations on eruption likelihood: application to two Icelandic subglacial volcanoes, Grimsvotn and Katla. *Geophys. J. Int.* **181**, 1510–1524 (2010).

32. Toda, S., Stein, R. S., Reasenberg, P. A., Dieterich, J. H. & Yoshida, A. Stress transferred by the 1995 M-w = 6.9 Kobe, Japan: shock: effect on aftershocks and future earthquake probabilities. *J. Geophys. Res.* **103**, 24543–24565 (1998).

33. Belardinelli, M. E., Cocco, M., Coutant, O. & Cotton, F. Redistribution of dynamic stress during coseismic ruptures: Evidence for fault interaction and earthquake triggering. *J. Geophys. Res.* **104**, 14925–14945 (1999).

34. Embley, R. W., Murphy, K. M. & Fox, C. G. High-resolution studies of the summit of Axial Volcano. *J. Geophys. Res.* **95**, 12785–12812 (1990).

35. Chase, R. et al. Hydrothermal vents on an axis seamount of the Juan de Fuca ridge. *Nature* **313**, 212–214 (1985).

36. Beeler, N. M. & Lockner, D. A. Why earthquakes correlate weakly with the solid Earth tides: Effects of periodic stress on the rate and probability of earthquake occurrence. *J. Geophys. Res.* **108** https://doi.org/10.1029/2001jb001518 (2003).

37. Carpenter, B. M., Collettini, C., Viti, C. & Cavallo, A. The influence of normal stress and sliding velocity on the frictional behaviour of calcite at room temperature: insights from laboratory experiments and microstructural observations. *Geophys. J. Int.* **205**, 548–561 (2016).

38. Chen, J., Verberne, B. A. & Spiers, C. J. Effects of healing on the seismogenic potential of carbonate fault rocks: Experiments on samples from the Longmenshan Fault, Sichuan, China. *J. Geophys. Res.* **120**, 5479–5506 (2015).

39. Ikari, M. J., Ito, Y., Ujiie, K. & Kopf, A. J. Spectrum of slip behaviour in Tohoku fault zone samples at plate tectonic slip rates. *Nat. Geosci.* **8**, 870–87 (2015).

40. Atkinson, B. K. Subcritical crack growth in geological materials. *J. Geophys. Res.* **89**, 4077–4114 (1984).

41. Shearer, P. M., Prieto, G. A. & Hauksson, E. Comprehensive analysis of earthquake source spectra in southern California. *J. Geophys. Res.* **111**, B06303 https://doi.org/10.1029/2005jb003979 (2006).

42. Boyd, O. S., McNamara, D. E., Hartzell, S. & Choy, G. Influence of lithostatic stress on earthquake stress drops in North America. *Bull. Seismol. Soc. Am.* **107**, 856–868 (2017).

43. Dieterich, J. In *Earthquake Source Mechanics*. AGU Geophys. Mono. (eds Das, S., J. Boatwright, & C. Scholz) 37–49 (American Geophysical Union, Washington, D.C.).

44. Roy, M. & Marone, C. Earthquake nucleation on model faults with rate- and state-dependent friction: effects of inertia. *J. Geophys. Res.* **101**, 13919–13932 (1996).

45. Stein, R. S. The role of stress transfer in earthquake occurrence. *Nature* **402**, 605–609 (1999).

46. van der Elst, N. J. & Brodsky, E. E. Connecting near-field and far-field earthquake triggering to dynamic strain. *J. Geophys. Res.* **115**, B07311 https://doi.org/10.1029/2009jb006681 (2010).

47. Hill, D. P. Seismicity remotely triggered by the magnitude 7.3 Landers, California, Earthquake. *Science* **260**, 1617–1623 (1993).
48. Husen, S., Taylor, R., Smith, R. B. & Healser, H. Changes in geyser eruption behavior and remotely triggered seismicity in Yellowstone National Park produced by the 2002 M 7.9 Denali fault earthquake, Alaska. *Geology* **32**, 537–540 (2004).

49. Prejean, S. G. et al. Remotely triggered seismicity on the United States west coast following the M-W 7.9 Denali fault earthquake. *Bull. Seismol. Soc. Am.* **94**, S348–S359 (2004).

50. Brodsky, E. E., Roeloffs, E., Woodcock, D., Gall, I. & Manga, M. A mechanism for sustained groundwater pressure changes induced by distant earthquakes. *J. Geophys. Res.* **108**, 2390 https://doi.org/10.1029/2002jb002321 (2003).

51. Linde, A. T., Sacks, I. S., Johnston, M. J. S., Hill, D. P. & Bilham, R. G. Increased pressure from rising bubbles as a mechanism for remotely triggered seismicity. *Nature* **371**, 408–410 (1994).

52. Manga, M. et al. Changes in permeability caused by dynamic stresses: field observations, experiments, and mechanisms. *Rev. Geophys.* **50**, RG2004 https://doi.org/10.1029/2011rg000382 (2012).

53. Yabe, S., Tanaka, Y., Houston, H. & Ide, S. Tidal sensitivity of tectonic tremors in Nankai and Cascadia subduction zones. *J. Geophys. Res.* **120**, 7587–7605 (2015).

54. Beeler, N. M., Thomas, A., Burgmann, R. & Shelly, D. Constraints on friction, dilatancy, diffusivity, and effective stress from low-frequency earthquake rates on the Deep San Andreas Fault. *J. Geophys. Res.* **123**, 583–605 (2018).

55. Liu, Y. & Rice, J. R. Spontaneous and triggered aseismic deformation transients in a subduction fault model. *J. Geophys. Res.* **112**, B0904 (2007).

56. Vera, E. et al. The structure of 0-to 0.2-my-old oceanic crust at 9° N on the East Pacific Rise from expanded spread profiles. *J. Geophys. Res.* **95**, 15529–15556 (1990).

57. Wilcock, W., Waldhauser, F. and Tolstoy, M. Catalogs of earthquake recorded on Axial Seamount from January, 2015 through November, 2015. doi:10.1594/IEDA/323843 (Interdisciplinary Earth Data Alliance, 2016).

58. Wiemer, S. & Wyss, M. Minimum magnitude of completeness in earthquake catalogs: examples from Alaska, the western United States, and Japan. *Bull. Seismol. Soc. Am.* **90**, 859–869 (2000).

59. Cao, A. M. & Gao, S. S. Temporal variations of seismic b-values beneath northeastern japan island arc. *Geophys. Res. Lett.* **29**, 48-1-48-3 https://doi.org/10.1029/2001GL013775 (2002).

60. Agnew, D. C. NLOADF: a program for computing ocean-tide loading. *J. Geophys. Res.* **102**, 5109–5110 (1997).

61. Egbert, G. D. & Erofeeva, S. Y. Efficient inverse modeling of barotropic ocean tides. *J. Atmos. Ocean Tech.* **19**, 183–204 (2002).

62. Chadwick, W., Embley, R. W., Claque, D. & Butterfield, D. A. Processed swath bathymetry grids (NetCDF-GMT format) derived from Multibean Sonar Data from the Juan de Fuca -Axial Seamount assembled as part of the Jdf_Gorda NOAA_VENTS Data Compilation Interdisciplinary Earth Data Alliance (IEDA) https://www.iedadata.org/ (2015).

## Acknowledgements

Mike Burton is thanked for instructing us on the properties of magmas at mid-ocean ridges. The earthquake catalog is the same as in Wilcock (2016) and is archived in the Interdisciplinary Earth Data Alliance Marine Geoscience Data System (DOI: 10.1594/IEDA/323843). That data was collected with funding from NSF under grant OCE-1536320. F. Albino is supported by the NERC Center for the Observation and Modelling of Earthquakes, Volcanoes and Tectonics.

## Author contributions

C.H.S. conceived the main ideas with discussions with Y.J.T., led the project, and wrote the initial draft of the paper. Y.J.T. performed the analysis of the tides and earthquake data, fit the triggering equations, contributed to their interpretation, and prepared Figs. 1, 2, and 5. F.A. developed the deformation model, carried out all deformation modeling, and prepared Figs. 3 and 4. All authors contributed to the writing of the final paper.

## Additional information

**Supplementary Information** accompanies this paper at https://doi.org/10.1038/s41467-019-10605-2.

**Competing interests:** The authors declare no competing interests.

**Reprints and permission** information is available online at http://npg.nature.com/reprintsandpermissions/

**Journal Peer Review Information:** *Nature Communications* thanks Elizabeth Cochran, Andrew Delorey and John Vidale for their contribution to the peer review of this work. Peer reviewer reports are available.

**Publisher's note:** Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

<figure><img src="image_7.png" /></figure>

**Open Access** This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons license, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons license and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this license, visit http://creativecommons.org/licenses/by/4.0/.

© The Author(s) 2019