Statistical Features of Foreshocks in Instrumental and ETAS Catalogs

Abstract

We study the spatial distribution of earthquakes in temporal intervals before and after the occurrence of large shocks (mainshocks) in the magnitude range $m \in [2,5]$ for four different regional catalogs. We find that the spatial organization of pre-shock seismicity depends on the mainshock magnitude and is independent of the lower magnitude threshold. These properties are found to be a stable feature of regional catalogs and cannot be reproduced by Epidemic Type Aftershock Sequence models. Our findings suggest that the area fractured during the mainshock is encoded in the foreshock spatial organization and, therefore, enhance the prognostic value of foreshocks.

Appendices

Appendix 1: Model Parameter Selection

In this appendix we present the procedure adopted to select the model parameters $B,\alpha ,p,c$ of the ETAS model (Eq. 2) so that synthetic catalogs reproduce statistical features of the RSCEC. The first observation (Fig. 7) is that for any choice of $\alpha$, the ratio $n_a(m_M)/n_m(m_M)$ for events identified as aftershock by our declustering procedure, follows the productivity law $n_a(m_M)/n_m(m_M) \simeq K_a e^{\alpha _a m_M}$. Numerical results are in agreement with the productivity law for events identified as aftershocks. In particular, for a wide range of values of $\alpha \in [0.8,1] b \log (10)$ implemented in numerical simulations, we always find $\alpha _a \simeq \alpha$. For this reason we set in numerical simulations $\alpha =1.8$ corresponding to the value of $\alpha _a$ obtained as best fit of experimental data from RSCEC. The value $B=0.115$, in numerical simulations, is then fixed to have the same value of $K_a$ in the synthetic and in the RSCEC catalog (Fig. 7). The other parameters p and c affect only weakly the main findings (B) of our study and we implement $p=1.2$ and $c=1$ s. Finally the number of independent background events is chosen to have a similar number of events with $m>m_c$ in the synthetic and instrumental catalogs and their spatial distribution is obtained from the smoothed distribution of instrumental earthquakes.

Fig. 7

The number of $m>2$ identified aftershocks $n_a(m_M)$ divided by the number of identified mainshock $n_M(m_M)$ is plotted versus the magnitude class $m_M$. Black open circles are results obtained for the RSCEC catalogs, and other symbols and colors are used for results in synthetic catalogs with different values of $\alpha$ according to the legend. On the right of each curve we give the value of the best fit parameter $\alpha _a$ from the exponential relation (productivity law) $n_a(m_M)/n_M(m_M) =K_a e^{\alpha _a m_M}$. Open red squares corresponding to $\alpha =1.8$ provide the best agreement with the RSCEC results. Filled symbols are used for the number of $m>2$ identified foreshocks, black circles for the RSCEC and red squares for the ETAS catalog with $\alpha =1.8$. On the right we give the value of the best fit parameter $\alpha _f$ from the foreshock productivity law $n_f(m_M)/n_M(m_M) =K_f e^{\alpha _f m_M}$

Appendix 2: Independence of the Spatial Distribution on Model Parameters

In this appendix we present the aftershock and foreshock spatial densities, $\rho _a(\Delta r,m_M)$ and $\rho _f(\Delta r,m_M)$ measured in synthetic catalogs obtained with different values of the ETAS parameters $\alpha ,p,c$. Results (Figs. 8, 9, 10, 11, 12 and 13) show that both distributions do not present significant differences for different parameters.

Fig. 8

The aftershock linear density distribution in the synthetic catalog for different values of $\alpha$. Different colors correspond to different mainshock magnitude classes $m_M$. Inset the scaling collapse of $\rho _a(\Delta R,m_M)$ according to Eq. (3)

Fig. 9

The foreshock linear density distribution in the synthetic catalog for different values of $\alpha$. Different colors correspond to different mainshock magnitude classes $m_M$

Fig. 10

The aftershock linear density distribution in the synthetic catalog for different values of p. Different colors correspond to different mainshock magnitude classes $m_M$

Fig. 11

The foreshock linear density distribution in the synthetic catalog for different values of p. Different colors correspond to different mainshock magnitude classes $m_M$

Fig. 12

The aftershock linear density distribution in the synthetic catalog for different values of c. Different colors correspond to different mainshock magnitude classes $m_M$

Fig. 13

The foreshock linear density distribution in the synthetic catalog for different values of c. Different colors correspond to different mainshock magnitude classes $m_M$

Appendix 3: The Influence of the Declustering Procedure

In this section we provide additional evidence confirming that the main results of the manuscript are not affected by possible biases due to the declustering criterion. In particular in ref. Richards-Dinger et al. (2010) it has been proposed that couples of events identified as mainshock–aftershock according to the window criterion could be both aftershocks of an exogenous larger event occurred at distances larger than L and/or at times farther than y days. This claim was supported by three main observations: (a) the linear density distributions of events occurring immediately before or after the mainshock roughly coincide; (b) the Omori law decay tends to become less steep if periods of one year just after large shocks ($m>7$) are not included in the analysis; (c) the cumulative number of events identified as aftershocks presents jumps of activity in correspondence to large shocks.

We notice that the observations (b) and (c) apply to aftershocks identified within a spatial range larger than 10 km. In our study we focus on aftershocks on a smaller spatial range, $R \le 10$ km (see Sect. 9.1), where none of the two observations is verified. In Fig. 14a we plot the number of $m>1$ aftershocks in the entire RSCEC catalog and in the catalog where one year after the three $m>7$ earthquakes is removed. The number of identified mainshocks is roughly reduced by half but the direct and inverse Omori law (Fig. 14a) are practically unaffected. In the right panel of Fig. 14 we plot the cumulative fraction of $m>2$ events $C(\tau )$, i.e., the cumulative number of events up to time t divided by the total number of $m>2$ events in the catalog, as function of the normalized time $\tau =(t-t_0)/(t_f-t_0)$. Here $t_0$ and $t_f$ are the initial and final time of the catalog, respectively. As expected, $C(\tau )$ presents jumps soon after the occurrence of $m \ge 6.1$ large shocks plotted in the figure as vertical brown lines. Conversely, the cumulative number of mainshocks with magnitude $m \in [2,4)$, as well as the number of $m>2$ identified aftershocks follow a straight line, corresponding to Poisson occurrence in time, and are totally uncorrelated to jumps in $C(\tau )$.

Fig. 14

The number of aftershocks (Left panel) (b) and foreshocks (a) defined as events with $m>m_{th}=1$, recorded inside acircle of radius of 10km centered at the mainshock epicenter, as function of the time distance $\Delta t=\vert t-t_M\vert$ and divided by the number of mainshocks $n_M$. Different symbols and different colors corresponding to different $m_M=2,3,4$ are results for RSCEC. Continuous green lines are the results for the same data set where periods of one year after each of the three $m>7$ shocks are removed from the catalog. The cumulative fraction $C(\tau )$ of $m \ge 2$ events in RSCEC (Right Panel) (green diamonds) is plotted as function of the rescaled time $\tau$. Black circles and red squares represent the cumulative fraction of $m>2$ mainshocks and aftershocks, respectively, identified according to the window criterion. The continuous blue line is the diagonal $y=x$ of Poisson occurrence in time. Brown vertical line correspond to the occurrence time of $m\ge 6.1$ earthquakes

Furthermore observation (a), as deeply discussed in the paper, is a real physical effect and not a bias of the declustering procedure. Indeed, the linear density distribution of aftershocks and foreshocks are very similar, even if the number of aftershocks is larger than the foreshock one. The latter result cannot be observed if events identified as fore–mainshock couples or mainshock–aftershock couples are both events triggered by the same far-in-time large shock. In addition, this hypothesis cannot explain the following results:

• The aftershock and foreshock productivity law showing that the number of both aftershocks and foreshocks increases with the magnitude class $m_M$ of the events identified as mainshock (features A1, F1) (Fig. 3);

• The increase in the seismicity rate as the occurrence time $t_M$ of the mainshock is approaching and the decrease after $t_M$ (Fig. 3), for all catalogs and for all magnitude classes (features A2,F2);

• The dependence of the average mainshock–aftershock and mainshock–foreshock epicentral distances on the magnitude class $m_M$ (features A3, F3).

All these results clearly show that the organization in time, space and magnitude of events occurring after and before $t_M$ are controlled by the properties of the mainshock and cannot be attributed to exogeneous earthquakes. To further support this conclusion, we explicitly investigate the dependence of the aftershock and foreshock linear density on the parameters L, y and $y_2$ of the declustering procedure. The three parameters represent the extent of the spatial area and the temporal intervals before and after the event, respectively, where larger shocks should not be recorded. Obviously, the larger the three parameters, the more accurate is the declustering procedure, as well as the smaller is the number of events identified as mainshocks. Results, in the main text, are obtained for $L=100$ km, $y=3$ days and $y_2=0.5$ days representing a compromise between an efficient declustering method and a sufficiently large number of mainshocks allowing an accurate statistical analysis. Here we explore several combinations of the set $(L,y,y_2)$ in a wide range of values, $L \in [100,350]$ km, $y \in [1,15]$ days and $y_2 \in [0.25,1.5]$ days and show that, in all cases, no significant differences can be appreciated for $\rho _f(\Delta r,m_M)$ and $\rho _a(\Delta r,m_M)$ with different $m_M$. We have only considered values of $y_2 \le 1.5$ days because it is unlikely that an earthquake is statistically correlated to a subsequent larger shock occurring far in time. In Fig. 15 we compare results for $(L = 100\, km, y = 15 \,days, y_2 = 1.5\, days)$ with three limiting cases, $(L=100\, km,y=15 \,days,y_2=1.5\, days)$, $(L=350\, km,y=3 \,days,y_2=0.5\, days)$ and $(L=100\, km,y=1 \,days,y_2=0.25 \,days)$, for $m_M=2,3$ and different $m_{th}=1,2$. Results show that the aftershock and foreshock spatial linear density are independent of the selection parameters, supporting the efficiency of the declustering procedure.

To better explore this point, in Fig. 15 we plot results obtained implementing a more accurate criterion for the mainshock identification introduced by Baiesi and Paczuski (BP) (Baiesi and Paczuski 2004; Gu et al. 2013). The BP selection criterion defines the time–space–magnitude distance $n_{ij} \propto \tau _{ij}r_{ij}^D 10^{-b m_i}$, for any two earthquakes i and j. Here $\tau _{ij}$ is the temporal distance between the two earthquakes, $r_{ij}$ the epicentral distance, D the (fractal) dimension of earthquake epicenters, and b the parameter of the Gutenberg–Richter law. We use the same parameters of ref. Baiesi and Paczuski (2004); Gu et al. (2013) to obtain clusters of “correlated” earthquakes and we naturally define as a mainshock the event with the largest magnitude in each cluster not belonging to previous clusters. This different mainshock selection procedure does not produce (Fig. 15) significant changes in $\rho _f(\Delta r,m_M)$ and $\rho _a(\Delta r,m_M)$ for different $m_M$, with respect to the window criterion. The same conclusion is also obtained by means of a faster but less accurate selection method based on the correlation coefficient value  (Bottiglieri et al. 2009).

Fig. 15

The aftershock (foreshock) linear density is plotted with open (filled) symbols for different values of $m_M$ in the four geographic regions. Different symbols and colors correspond to different parameters of the window declustering criterion or to the BP declustering method. In particular, black circles are for the window criterion with the parameters $(L=100,y=3,y_2=0.5)$ used in the paper. Green diamonds are for $(L=100,y=15,y_2=2.5)$, blue up triangles are for $(L=350,y=3,y_2=0.5)$ and orange down triangles for $(L=100,y=1,y_2=0.25)$. Red squares are the results obtained with the BP declustering criterion

To further support that the linear density distribution is substantially unaffected by the declustering procedure, we consider the same declustering method with parameters $(L=100,y=3,y_2=0.5)$ but we impose the extra constraint that events of any magnitude $m>m_{th}$ must not be recorded in the temporal interval [1, 5] days before the mainshock within a radius of 10 km from the epicenter. This constraint should further reduce the probability that the mainshock belongs to an aftershock sequence of a previous event. Results plotted in Fig. 16 show that the linear density distribution does not significantly change if this further constraint is applied.

Fig. 16

The same as Fig. 2 with the additional constraint of zero events in the temporal interval [1, 5] days before the mainshock within a radius of 10 km from its epicenter

All these considerations, combined with the tests performed on synthetic catalogs (Figs. 1, 2), provide strong support in favor of our findings being independent of the declustering procedure.

1.1 Appendix 3.1: The temporal interval for aftershock selection

Once mainshocks have been identified, we use the simplest definition of aftershocks (foreshocks) as all events occurring within a time interval of T hours after (before) $t_M$ and within a circle of radius R km centered in the mainshock epicenter. Since the number of background events inside this circle is roughly proportional to $T R^2$, the smaller T and R the larger is the probability that events occurring in the spatio-temporal interval (R, T) are “true” aftershocks or foreshocks. Since we are mainly interested in the aftershock and foreshock spatial organization, we can leave R as a free parameter and fix only the value of T. Results for $\rho _a(\Delta r,m_M)$ and $\rho _f(\Delta r,m_M)$ are obtained for $T=12$ h which ensures both an efficient selection procedure and a sufficiently large number of aftershocks and foreshocks. This value is justified by the analysis performed in Fig. 17 where we plot $\rho _f(\Delta r,m_M=2)$ for $\Delta r \in (0,10]$ km and different values of $T \in [3,48]$ hours. We find that for $T \le 12$ hours, data are substantially T-independent suggesting that our selection procedure is able to capture the correct linear density distribution.

We also notice that the tail of the distribution, for each catalog, becomes flatter and flatter by increasing T which is the clear signature of the influence of background seismicity. For this reason we have considered, for each catalog, epicentral distances only up to a value $R_{max}$ under the condition that, for $\Delta r \,<\,R_{max}$, data for $T=6$ h and $T=12$ h do not present significant differences.

Fig. 17

Foreshocks are defined as all events occurring within an interval T before the mainshock inside a region of radius $\Delta r$. The foreshock linear density is then plotted for different T and instrumental catalogs. The dotted violet vertical line indicates the upper value $R_{max}$ considered for each catalog

Appendix 4: The Role of Catalog Completeness

All experimental catalogs present a completeness level $m_L$ that, in the majority of regions and temporal intervals, is systematically larger than $m_{th}=1$. Nevertheless, incompleteness is expected to influence the number of events but not their spatial distribution and, therefore, it does not affect the aftershock and foreshock linear density distribution. To support this conclusion we have artificially introduced short-term incompleteness in ETAS catalogs according to the following procedure. We use the empirical relation for the short-term completeness level $m_c(t,m_i)$ (Helmstetter et al. 2006) expected at time t after a magnitude $m_i$ earthquake occurred at time $t_i$, $m_c(t,m_i)=m_i-0.76 \log 10(t-t_i)-4.5$, where time is measured in days. At each time t, we define a completeness level $m_L(t)=max_{t_i<t}m_c(t,m_i)$ and remove from the catalog all events with magnitude smaller than $m_L$. This procedure leads to the removal of 41% of events, reflecting in important differences in the aftershock number and their organization in time. Conversely, no relevant differences (Fig. 18) can be observed in $\rho _a$ and $\rho _f$ between the complete and incomplete catalogs.

Fig. 18

The distributions $\rho _a(\Delta r,m_M)$ and $\rho _f(\Delta r,m_M)$ in the ETAS catalog are plotted as continuous lines for $m_M=2,3,4$ and compared with the distributions evaluated in the incomplete ETAS catalog

We further observe that incompleteness decreases by increasing the lower threshold $m_{th}$ and, therefore, can introduce a spurious dependence on $m_{th}$. However, our results are shown to be independent of $m_{th}$ (Fig. 19) with the aftershock distribution $\rho _a$, obtained for $m_{th}=1$, in very good agreement with results for $m_{th}=2$, for all catalogs and different $m_M$ values (feature A4). This behavior is recovered in synthetic catalogs.

On the other hand, the strong similarity (Fig. 4) between $\rho _f$ for $m_{th}=1$ and $m_{th}=2$, obtained for RSCEC, RNCEC and ItEC, can be attributed to incompleteness only if the number of foreshocks, with magnitude in the range $m \in [1,2)$, is too small to significantly affect the $\rho _f$ distribution. This is not the case, since we find that $N_f(1,2)> 1.5 N_f(2,\infty )$, for all values of $m_M$ and for all catalogs. Here $N_f(m_1,m_2)$ is the number of events identified as foreshocks with magnitude m in the range $m \in [m_1,m_2)$. Furthermore, $\rho _f(\Delta r,m_M)$ obtained considering only events with $m \in [1,2)$ roughly collapses on the $\rho _f(\Delta r,m_M)$ for $m\ge 2$. The same argument applies to the JaEC, restricting to foreshocks with $m \in [2,3)$ or considering all $m \ge 3$ foreshocks. Therefore, catalog incompleteness cannot represent an explanation for the different dependence on $m_{th}$ between foreshocks in synthetic and experimental catalogs.

Fig. 19

The aftershock epicentral linear distributions for different instrumental catalogs. Open symbols show results for $m_{th}=2$ as considered in Fig. 4 whereas continuous lines are results for $m_{th}=1$