Aftershock Rate Changes at Different Ocean Tide Heights

The differential probability gain approach is used to estimate quantitatively the change in aftershock rate at various levels of ocean tides relative to the average rate model. An aftershock sequences are analyzed from two regions with high ocean tides, Kamchatka and New Zealand. The Omori-Utsu law is used to model the decay over time, hypothesizing an invariable spatial distribution. Ocean tide heights are considered rather than phases. A total of 16 sequences of M ≥6 aftershocks off Kamchatka and 15 sequences of M ≥6 aftershocks off New Zealand are examined. The heights of the ocean tides at various locations were modeled using FES 2004. Vertical stress changes due to ocean tides are here about 10–20 kPa, that is, at least several times greater than the effect due to Earth tides. An increase in aftershock rate is observed by more than two times at high water after main M ≥6 shocks in Kamchatka, with slightly less pronounced effect for the earthquakes of M = 7.8, December 15, 1971 and M = 7.8, December 5, 1997. For those two earthquakes, the maximum of the differential probability gain function is also observed at low water. For New Zealand, we also observed an increase in aftershock rate at high water after thrust type main shocks with M ≥6. After normal-faulting main shocks there was the tendency of the rate increasing at low water. For the aftershocks of the strike-slip main shocks we observed a less evident impact of the ocean tides on their rate. This suggests two main mechanisms of the impact of ocean tides on seismicity rate, an increase in pore pressure at high water, or a decrease in normal stress at low water, both resulting in a decrease of the effective friction in the fault zone.

Laboratory experiments were also conducted to study the effect of tides on seismicity (Lockner and Beeler, 1999;Beeler and Lockner, 2003). A strong correlation was found to exist between periodic stress (tides) and the occurrence of failure (earthquake) at shear stress amplitudes above approximately 0.3 MPa. Shear stress variations between 10 kPa (1 m of ocean tide load) and 0.1 MPa represent a transition region in which correlation with earthquake occurrence may occur. For amplitudes below 10 kPa (the order of the vertical component of the Earth tide), little or no correlation can be detected. These results are consistent with observations. This suggested that tidal stress amplitudes of about 3-10 kPa are required to trigger earthquakes (Hardebeck et al., 1998;Cochran et al., 2004). For example, it was shown that aftershocks after the M 7.4 Landers event were triggered by stress increases greater than 10 kPa (Stein, 1999).
There are relatively few studies where the connection between seismicity and ocean tides was analyzed. This is due to the complexity of numerical calculations of ocean tides. Only in recent years computer programs for the numerical modeling of ocean tides have been developed. These programs were used to study the connection between ocean tides and seismicity (Crockett et al., 2006;Ide and Tanaka, 2014;Baranov et al., 2019). The result was to detect a clear connection between seismicity and ocean tides. Despite a large number of comparative studies of seismicity and tides, the relationship between aftershock rates and tides was studied in few papers (Souriau et al., 1982;Chen et al., 2012b;Datta and Kamal, 2012;Saltykov, 2014;Baranov et al., 2019). A partial dependence of aftershock rate on tidal heights (ocean and Earth tides) for normal-fault and thrust-fault earthquakes was found for the Pyrenees (Souriau et al., 1982). But for strike-slip earthquakes no statistical connection was found. The M >7 earthquakes worldwide since 1900 are more likely to occur during the 0°, 90°, 180°or 2,70°phases (i.e., earthquake-prone phases) of the semidiurnal solid Earth (M 2 ) tidal curve (Chen et al., 2012b). Diurnal and semi-diurnal periodicities in aftershock rate (M ≥4) were found for the aftershock sequence of the Tohoku M 9.1 earthquake of March 11, 2011 in Japan with an apparent weakening of the tidal triggering effect over time (Datta and Kamal, 2012). This suggests that large aftershocks in the fault zone of the Japan 2011 earthquake were strongly influenced by Earth tides.
Partial dependence of aftershock rates on Earth tides was also found for aftershocks of the M 6.8 earthquake of June 21, 1996 off Kamchatka (Saltykov, 2014), and an alternative trigger mechanism was proposed for tidal effects based on an amplitude-dependent dissipation model. Another study of aftershock sequences following Kamchatka earthquakes using the method of differential probability gain, DPG (Baranov et al., 2019), demonstrated a considerable increase in the aftershock rate (by factors of two and more) at low or at high water. An increase in aftershock rate at low water corresponds to unloading of the seafloor, while high water may result in an increase of pore pressure in the fracture zone and therefore decreasing friction forces. The correlation between seismicity and tides in relation to focal mechanisms was studied by several authors. Correlation was mainly found for normal and partly for thrust slip types (Heaton, 1975;Souriau et al., 1982;Tsuruoka, et al., 1995;Tanaka et al., 2002;Cochran et al., 2004). The connection between seismicity and tides at ocean ridges were studied extensively (Wilcock, 2001;Tolstoy et al., 2002;Stroup et al., 2007;Scholz et al., 2019). A correlation between tides and seismicity was found along the Juan de Fuca Ridge (Wilcock, 2001;Tolstoy et al., 2002) and for East Pacific Rise (Stroup et al., 2007). Scholz et al. (2019) suggested pulsation of magma chambers due to vertical stress changes as a mechanism responsible for the tidal triggering of earthquakes for midocean ridges.
The impact of tides on earthquakes is often explained by changes in the stress field in the fault using the Coulomb criterion (Stein, 1999;Cochran et al., 2004). Coulomb failure stress change Δτ c Δτ + μ(Δσ n −ΔP), where τ and σ n are the shear stress change on the fault (assumed to be positive in the direction of slip) and the normal stress change on the fault (positive if the fault is unclamped), respectively, P is the pore pressure change in the fault, and μ is the coefficient of friction (with range 0-1). Failure is encouraged if Δτ c is positive and discouraged if negative; both increased shear and unclamping of faults facilitate failure. Thus friction threshold can be exceeded either when the normal stress is decreased or when pore pressure is increased (Klein, 1976a;Klein, 1976b;Wilcock, 2001;Tolstoy et al., 2002;Cochran et al., 2004;Stroup et al., 2007;Métivier et al., 2009). This can reflect earthquake initiation (Stein, 2004). Some authors suggested alternative mechanisms to explain the effect of tides on seismicity: nonlinear dilatancy-diffusion model (Heaton, 1982), seismic modulation by combination of different tidal waves, impacts of tidal waves occurring in resonance with the selfoscillating system of the focal zone under the influence of tidal waves (Nikolaev, 1996) or an amplitude-dependent dissipation model taking into account tidal variations in physical properties  (Nikolaev, 1994) Hadson`s catalogue from 1962 to 1988, M > 5 Partly for some areas 10 (Tsuruoka, et al., 1995) + + + + Partly for normal fault, maximum seismicity for tidal phase from 0°to 60°− Normal-fault-type earthquakes tend to occur at the time when the cubic tidal stress takes a maximum tensile value or a little bit later 988 globally distributed earthquakes with magnitude of 6.0 or larger from CMT from 1977 to 1992 −For strike-slip or thrust 11 (Nikolaev, 1996)     of the medium (Saltykov, 2014). The researchers used different methods to study the effect of tides on seismicity from simple visual comparison between time series of seismicity and tides to more or less sophisticated statistical methods as binomial distribution test, random distribution test, method of Chapman-Miller (Malin and Chapman, 1970), and Schuster's test using a statistical method of Rayleigh (Schuster, 1897). Schuster's test was described in detail in many papers (e.g., Heaton, 1975;Rydelek et al., 1992). Table 1 summarized studies of tides on seismicity using several key parameters: type of tides (Earth tides, ocean tides), tidal phases, connection with seismicity, tidal amplitude, the mechanism responsible, and the object of study. Why do we study the impact of ocean tides on aftershock rate using a model of tide heights?
(1) Time scales of tides and aftershock occurrence are comparable. The time-dependent distribution of aftershocks without impact of tides can be modeled.
(2) The impact of ocean tides on aftershock rates still remains a challenge (Cochran et al., 2004;Tanaka, 2012;Baranov et al., 2019). (3) The well-known Schuster test using tide phases requires a correct catalogue declustering. Aftershocks that remain in the catalogue may alter the statistics. (4) The impact of ocean tides on seismicity is often related to instabilities or a compliance of fault zones (Cochran et al., 2004). The aftershock zones marked by high stress perturbations perfectly fit those conditions. (5) We concentrate here on ocean tides in areas where the amplitudes are 2 m or more (Kamchatka and New Zealand). The corresponding stress changes (10-20 kPa) are larger compared with Earth tides (Varga and Grafarend, 2017). (6) We study ocean tide heights, not phases. Triggering of earthquakes is caused by stress changes depending on the level of tides. The frequency structure of the ocean tides is complex, and the prevailing periods may vary in time (Lyard et al., 2006).
The purpose of this publication is to assess the quantitative effect of ocean tidal heights on seismic activity and to compare these estimates for different types of faults.

Selection of Aftershock Sequences
The aftershock sequences for analysis were found using the "nearest neighbor" algorithm of Zaliapin and Ben-Zion (Zaliapin and Ben-Zion, 2013). Aftershock sequences of large earthquakes may have complex structure, with significant "splashes" of secondary aftershocks caused by large primary aftershocks. Another component of seismicity is background seismicity, which is often referred to as seismic noise. Using the "nearest neighbor" algorithm, we selected only the direct "offspring" (Zaliapin and Ben-Zion, 2013) of the considered large earthquakes, thus minimizing the presence of secondary aftershocks and background seismicity in the selected sequences. This allowed us to model the aftershocks sequences using Omori's law. One alternative that we used in a previous analysis (Baranov et al., 2019) is the Molchan-Dmitrieva algorithm (Molchan and Dmitrieva, 1992), which selects aftershock sequences together with secondary aftershock subsequences. Large aftershocks often trigger a temporary increase of activity. In such cases the Omori model may be significantly altered by this temporary activation. Direct aftershocks found using the nearest neighbor approach do not contain the secondary aftershock sequences, and thus are correctly modeled by the Omori law. Another alternative could be a stochastic declustering of the earthquake catalogue with much deeper complexity and non-uniqueness of the clusters (Varini et al., 2020). Another disadvantage of the stochastic declustering methods is their basic hypothesis that the number of aftershocks is a function of the magnitude of the corresponding main shock. This hypothesis, as was recently found, is not true (Shebalin et al., 2020). We selected M ≥6 main shocks in the Kamchatka region, and M ≥6 events in the New Zealand region, which have offshore aftershocks. Figures 1A,B show a map of the maximum ocean tide amplitude with resolution 0.125°together with examples of offshore aftershocks for M 7.5 (Kamchatka, June 08, 1993) and M 7.8 (New Zealand, November 13, 2016) earthquakes. All aftershocks (as described above, we consider only direct aftershocks of the large earthquakes considered using the nearest neighbor approach), both onshore and offshore ones, were used to estimate the parameters of the Omori-Utsu law (Utsu, 1961), but only the offshore aftershocks were used to calculate the effects of ocean tides on seismicity. For selection of direct aftershocks in Kamchatka and New Zealand we used the parameters listed in Supplementaty Table S1 of (Shebalin et al., 2020).

The Magnitude of Completeness Magnitude M c and Time of Completeness t c
The catalogue completeness usually decreases after large earthquakes (Helmstetter et al., 2006;Hainzl, 2016;Shebalin and Baranov, 2017). For each aftershock sequence we began by estimating the magnitude of completeness M c using data in the interval (t c , 1 month) after the main shock by the MBS method (Cao and Gao, 2002;Wossner and Wiemer, 2005) and t c 1 h. At this stage we disregarded data within first 6 h after the main shocks, where which the catalogue usually remains incomplete even for relatively large magnitudes. The MBS method allows one to detect the frequently observed effect of the lack of lowmagnitude aftershocks during the first few hours and sometimes days after a large earthquake, not detectable by the Maximum Curvature technique (Wiemer and Wyss, 2000), see the example in Figure 2. Then, with the preliminary value of M c , we found t c using the equation (Helmstetter et al., 2006) M c M m −4.5−0.75 log 10 (t c ), where M m is the mainshock magnitude. The next step was to find the final value M c using the same MBS technique, but with the new value of t c . In the  Frontiers in Earth Science | www.frontiersin.org December 2020 | Volume 8 | Article 559624 following analysis we omit the interval (0, t c ) after the main shocks.

Modelling the Aftershock Rates
We model aftershock rates λ(t) using the Omori-Utsu model (Utsu, 1961) where t is the time after the main shock, and c, p, and K are parameters. We estimated the parameters using aftershocks with M ≥ M c in an interval (t c , 720 h) with the Bayesian approach , assuming that c and p are homogeneous. Figure 3 shows an example of the posterior distribution of (c, p) Bayesian estimates. We subdivide the aftershock zone into 0.2°× 0.2°boxes. The interval (t c , 720 h) is divided in subintervals of 0.2 h. Supposing the c-value and the p-value are equal in all boxes, we calculated the modeled rates of  the aftershocks λ ij , where i is the index of spatial box, and j the time subinterval: , N i is the actual number of aftershocks in the ith bin within the interval (t c , 720 h). We calculated the actual number of aftershocks λ ij in each spatio-temporal volume. Figure 3 shows an example of the p-value diagram for the M 7.8, Kamchatka, December 05, 1997 earthquake.

Modeling Tide Heights
Real ocean tides can be as high as 12-18 m in some bays. The tidal amplitudes can reach 2 m near the Pacific coast of Kamchatka and the coast of New Zealand, producing a pressure contrast of approximately 20 kPa. The elastic response of the solid Earth to the ocean load is obtained by convolution of seawater mass distribution using FES 2004 program (http://www.aviso.altimetry.fr/). FES 2004 gives tide height at a specified point at a specified time instant (Lyard et al., 2006). The program is based on the solution of tidal barotropic equations by finite elements (triangles) on a global element grid (∼1 million elements). It uses numerical models of ocean bottom topography and shoreline (Le Provost et al., 1994;Le Provost et al., 1998). The program can compute 15 main tidal components on a 1/8°grid (amplitudes and phases), as well as 28 additional tidal components. The presence of ice is incorporated for polar regions. The accuracy is within a few centimeters for open ocean and within 10 cm for offshore areas. The grid is not uniform, being denser near the shore and less detailed in the open ocean, according to Le Provost's criterion (Le Provost and Vincent, 1986). The program requires an input file that contains site coordinates and times in hours as measured from January 1, 1985. The program uses the sites as specified in the input file to yield tide heights at a required time instant. If a point is on land, the value is −9999. Using the FES 2004 software we modelled the heights h ij of the ocean tides in each spatio-temporal volume. Next, we built a map of maximum amplitudes of ocean tides for the entire world. This map was calculated using a program of our own, Amplitude. The program iterates over the coordinates of latitude and longitude at steps of 0.125°. For each point of the mesh, the program generates a time series for a year with time step 1 min. Then, for each time series, the FES 2004 program is launched, which calculates the tide height at a given point for the annual interval and finds the maximum tide amplitude at the point. Thus, a grid is obtained where for each point the maximum amplitude of ocean tide is found. Figure 4 shows the map of the maximum ocean tide amplitude with resolution 0.125°for the entire world. New Zealand and Kamchatka that we deal with here are regions with large ocean tide amplitudes.

Differential Probability Gains
The ocean tide height can be considered as a parameter controlling the relative changes of seismicity rates. If a correlation between ocean tides and seismicity rates exists, then it is possible to calculate what is the average change of the rates relative to an average model, at specific values of the control parameter ( Figure 5A). The differential probability gain (DPG) function Shebalin et al., 2014) is defined as the ratio of the actual number of seismic events to the number expected on the average model. It is a function of the control parameter. Here we estimate a smoothed differential probability gain function g(h) using ranges of the control parameter with constant width dh: Here, ω ij and λ ij denote, respectively, the actual numbers of seismic events and the numbers that are expected on the average model in spatial box i and time span j. The ocean tide height h may also be considered as an alarm function of a forecasting model. The error diagram (Molchan, 1991) evaluates retrospectively the performance of the model with respect to the reference model. With the error diagram ( Figure 5A), the differential probability gains are defined as the local slope (derivative) of the diagram ( Figure 5B). Thus we assume that at each point of the considered area and at any moment, a specific tide height corresponds to an increase or a decrease of the aftershock rate with respect to the local rate that would have been observed without the impact of tides. We model this  , 1979-2016 (B). The colour scale shows the relief (m) with respect to sea level; black contours with red labels show maximum amplitude of the ocean tide, m. The mainshock fault-plane solutions were taken from the Global Centroid-Moment-Tensor catalogue (Ekström et al., 2012); epicentres prior to 1976 are not reported in Global Centroid-Moment-Tensor, and they are marked by asterisks.  Table 2); N denotes the number of M ≥ M C aftershocks in the interval (t c , 720 h) measured from the mainshock time; the maximum amplitude of ocean tide in the interval (t c , 720 h) measured from mainshock time was converted to pressure, kPa; relief at the epicentre of main event (m) with respect to sea level; type indicates the focal mechanism of the main event. expected tide-independent aftershock rate by the Omori-Utsu law in time with a spatial distribution obtained by averaging over two months, as defined by Eq. 1. The function g(h) is the corresponding indicator of influence. The larger the difference between g(h) and 1, the more significant is this impact at corresponding values of the ocean tide height h.   In order to check if our results could have been obtained by chance, we estimated the confidence interval of the DPG estimates. For each aftershock sequence we generated 10,000 random synthetic catalogues that satisfy the Omori-Utsu law with parameters determined from the real catalogue, and repeated the procedure for determining the DPG with each synthetic catalogue. Synthetic catalogues were constructed in a standard way for a non-stationary Poisson point process with successive times found using a random number generator: where ξ is a random number uniformly distributed in (0,1), and intensity λ(t) is given by the Omori-Utsu formula. We started at t 0 t c /100 and did not take into account events with t 0 < t c , similarly to the analysis of the real catalogue. To each event we assigned the epicenter coordinates of the j-th event from the real aftershock sequence, randomly choosing an integer j uniformly distributed in (1, N), where N is the number of events in the real sequence. Finally, for each interval (h−dh, h + dh) of the tide height, we determined the mean and standard deviation of DPG from 10,000 values. If the value of the DPG falls outside the limits of the confidence interval determined in this way, the deviation of the DPG based on real data from 1.0 can be considered significant.

SEISMIC DATA
We investigated 16 aftershocks sequences of shallow earthquakes off Kamchatka, 1971-2013 (Table 2, Figure 6A). We used seismic data from the earthquake catalogue produced by the Kamchatka Branch of the Geophysical Survey of the Russian Academy of Sciences (Chebrov et al., 2016) for the period between 1962 and 2016 (http://www.emsd.ru/sdis/earthquake/ catalogue/catalogue.php). Most of the main shocks are thrust events ( Figure 6A;  Figure 6B). Tables 2,3 show the main parameters of these aftershock sequences. Figures 6A,B show a map of maximum ocean tide amplitude with resolution 0.125°. For Kamchatka near the coast, the amplitude of the ocean tides reaches 2 m, whereas for New Zealand it varies between 1.4 and 3 m, occasionally more. In Kamchatka all mainshock epicenters lie in the ocean, while for New Zealand some of the main shocks are on land, having a sufficient number of aftershocks in the ocean. Next, the quality of the catalogue for New Zealand is better, the magnitude of completeness M c is generally smaller, and this generally results in a larger number of aftershocks to study and produces a better statistical significance of results compared with Kamchatka. The quality of hypocenter location in the catalogues is not important for our analysis. The spatial variation of the tide heights in all considered aftershock areas do not exceed 0.01 m in Kamchatka, and 0.02 m except earthquakes of 1993 and November 2016 in New Zealand; for those two earthquakes the spatial variation is about 0.07 m (Figure 1). The accuracy of the determinations of times of events (s) is very high when considered in relation to the time scale of tides (h).

RESULTS
Aftershock sequences have been identified for 16 (M ≥6) events with epicenters near Kamchatka (see Figure 6A; Table 2) and 15 (M ≥6) events with epicenters near New Zealand (see Figure 6B; Table 3) using the methods described above. For these sequences we estimated M c and t c as described in Section "Magnitude of completeness M c and time of completeness t c ", and used these parameters to estimate s and p in 1 as described in Section "Modelling the aftershock rates" using data during the first month after the main shocks (Tables 2, 3). Then for rectangular 0.2°× 0.2°boxes in the ocean we estimated the modeled aftershock rates λ ij in (t c , 720 h) within the aftershock areas and time intervals of 0.2 h assuming no impact of tides (Eq. 1), counted the actual aftershock numbers ω ij , and used the FES 2004 program to find ocean tide heights h ij . Finally, using Eq. 2, for all sequences we calculated the differential probability gain as a function of tide height. Everywhere we used the value dh 0.15 m. We estimate also the confidence intervals for DPG using random catalogues as described in Differential Probability Gains.
The Kamchatka Region Figure 7 shows differential probability gain functions for four large earthquakes (M 6.7, December 28, 1984;M 6.3, July 10, 1987;M 6.8, March 2, 1992;M 6.8, February 28, 2013). Although the differential probability gain function shape may vary from one sequence to another, all examples demonstrate a clear increase of the function at high water (>0.5 m). A similar effect is observed for all sequences, with slightly less pronounced effect for the earthquakes of M 7.8, December 15, 1971 andM 7.8, December 5, 1997 (Figure 8). For those two earthquakes, the maximum of the DPG function is also observed at low water (near −1 m). We have constructed also the integral differential probability gain function for Kamchatka (Figure 9). Eq. 2 was applied to all spatial boxes and time spans corresponding to several earthquakes. This analysis shows a clear maximum impact of the ocean tides at high water (>0.5 m).

New Zealand
Focal mechanisms of large earthquakes in New Zealand are different: mostly normal events in the north, thrusts in the south, and strike-slip earthquakes in the middle. Accordingly, we divided the main shocks into three groups: thrust faults (Figure 10), normal faults ( Figure 11) and strike-slip faults ( Figure 12). For each group we calculated individual and integral differential probability gain functions as we did for the Kamchatka sequences. As can be seen from Figures 10, for thrust faults the most general phenomenon consists in an increasing rate of aftershocks at high water (>0.5 m). A similar effect was found for Kamchatka (Figures 7-9) where almost all large earthquakes considered here were of the thrust type. For normal fault type earthquakes there is no effect at high water. In contrast, we observed a maximum of DPG at low water (<−0.5 m, Figure 11). A similar effect was observed for aftershocks following the earthquakes of M 7.8, December 15, 1971 andM 7.8, December 5, 1997 in Kamchatka (Figure 8). For strikeslip earthquakes the maximum of DPG functions is reached at high water for all sequences (>0.5 m). At low water, high DPG values were clearly observed for aftershocks of the M 6.5 earthquake of July 21, 2013 (<−0.5 m, Figure 12C). These results demonstrate an impact of ocean tides on seismicity in the zones of large earthquakes within a period of 30 days after them. This highlights certain limitations to the analysis. The rupture zone of a large earthquake is in an excited state, so the impact of relatively weak disturbances, such as tides, may be stronger there than in the entire region.

LIMITATIONS OF THE ANALYSIS
In the analysis we used a model of ocean tides. The actual tide heights may be slightly different. However, continuous direct measurements do not exist at the moment. Another shortcoming of our analysis is that we could not ascertain the role of aftershock depth. This was mostly due to lack of data below the magnitude of completeness. An unreliable depth determination for offshore earthquakes is another factor. Accordingly, the physical model of the relationships found is very preliminary.

DISCUSSION AND CONCLUSION
In the present study we used the differential probability gain approach Shebalin et al., 2014) to estimate quantitatively the change in aftershock rate at various levels of ocean tides, relative to an average Omori-Utsu model that supposes no dependence on tides. The differential probability gain function is a numeric factor indicating how much the rate of aftershocks is increased or decreased on average at specific values of the tide heights. The value one indicates no impact. We consider the impact of ocean tides only, because their effect is several times larger than the effect of solid Earth tides in the regions under consideration. Variations of vertical stresses of 10-20 kPa due to ocean tides for the coast of Kamchatka and New Zealand fall within the range of the effect of stresses on seismicity according to laboratory experiments (Lockner and Beeler, 1999). Similar variations due to Earth tides at shallow depths (<50 km) are only approximately 1 kPa, i.e., one order of magnitude less (Varga and Grafarend, 2017).
We considered aftershock sequences of large (M ≥6) earthquakes off Kamchatka, 1971-2013(16 sequences) and New Zealand, 1979-2016. For all sequences we disregarded data within a few first minutes or hours after main shocks, in which the catalogue usually remains incomplete even for relatively large magnitudes. Only aftershocks with epicenters in the ocean were used. In these regions, the amplitudes of ocean tides are large, and their influence is a few times or tens of times larger than the effects of Earth tides. The observed increase in the rate of aftershocks at high and/or low water demonstrated a significant effect of ocean tides on seismicity. One important feature that distinguishes this study from most others is that we studied the heights rather than the phases of ocean tides. Moreover, the goal of this study was to find quantitative estimates of the increase in seismicity rate at specific heights of the ocean tides. The effect varies from place to place and for different focal mechanisms. The main result consists in finding a significant change in seismicity rates at tide heights more than 0.5 m relative to the baseline, positive or negative. For normal fault earthquakes the effect is stronger at low water, for thrust events mostly at high water, and for strike-slip earthquakes at high and/or low water.
Although the differential probability gains function shape may vary from one sequence to another, we observed a clear tendency of increasing aftershock rate by about two times at either low or high water. As an explanation of this effect, we can suggest friction reduction on a fault due to vertical stress decreasing at low waters and due to increased pore pressure at high waters. For normal faults an increase in aftershock rate was observed at low water, and thus the friction reduction mechanism due to the unloading of vertical stress is more likely. For thrust earthquakes, an increase in aftershock rate was observed at high water, and in that case increased pore pressure is the preferable mechanism. For strike-slip events with intermediate stresses both mechanisms can operate.
The main practical result of this study is a quantitative assessment of the effect of ocean tides on earthquake rate. Although the results were obtained for aftershocks, we may suppose that similar dependences are valid for all seismic events. Of course, this needs additional validation, but potentially opens a way to take into account ocean tides for time-dependent seismic hazard assessment.

DATA AVAILABILITY STATEMENT
The raw data supporting the conclusions of this article will be made available by the authors, without undue reservation.

AUTHOR CONTRIBUTIONS
PS provided the design and application of the differential probability gain tests. AB designed and performed calculations of the ocean tide heights.

FUNDING
We acknowledge the support from Russian Science Foundation Grant No. 20-17-00180.