Migration diffusivity as a controlling factor in the duration of earthquake swarms

Earthquake swarms exhibit highly uncertain temporal behavior. We investigated the relationship between the swarm duration and the diffusivity of hypocenter migration for triggered earthquake swarms in northeastern Japan. These parameters were systematically estimated by applying a diffusion model and using a unified definition of time windows for the initial and final stages of swarm activity. This approach detected a clear negative correlation between the diffusivity and swarm durations. The relation follows a power-law with an exponent of -0.5\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$-\,0.5$$\end{document} to -1.0\documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$-\,1.0$$\end{document}. Examination of published data confirmed that this relationship globally holds under various localities and tectonic environments. These results suggest that diffusivity, and by extension, crustal permeability and fluid viscosity play a key role in controlling the duration of the fluid-driven swarms.


Introduction
Earthquake swarms are a type of seismic activity that is not associated with a clear mainshock or a specific period. Swarms often occur in volcanic and geothermal fields (e.g., Chen and Shearer 2011; Shelly et al. 2016), in the vicinity of anthropogenic fluid injection wells (e.g., Keranen et al. 2014), and even in the areas where none of these causal mechanisms apply (e.g., Hauksson et al. 2019). Swarms show considerable variation in terms of their duration, but no theories exist to explain this variation. This study reports on a novel relationship observed between the hypocenter migration speed in the early stage and the duration of swarms.
Earthquake swarms differ significantly from ordinary earthquake sequences (e.g., mainshock-aftershock sequences) in their temporal behavior. Typical earthquake sequences consist of a large mainshock and subsequent aftershocks whose number decay obeys the Omori-Utsu law. Instead, earthquake swarms do not obey the Omori-Utsu law, making it difficult to evaluate or predict their temporal behavior. Durations of earthquake swarms exhibit considerable variation from a few days to indefinite periods of months and years. The Matsushiro earthquake swarm, for example, occurred from 1965 to 1967. Studies of the swarm found that subsurface migration of water and/or magma played a key role in driving the swarm (e.g., Hagiwara and Iwata 1968;Cappa et al. 2009). The volume and behavior of geofluids may thus affect earthquake swarm duration. Observing factors that influence swarm duration can help further elucidate physical mechanisms for the swarms. Better constraints on the mechanisms or refined expectations on durations of swarms can also help assess earthquake risks when the swarm area is close to human activity.
Hypocenter migration or crustal deformation often observed from earthquake swarms has been utilized for studying driving mechanisms for earthquake swarms. Several candidates have been proposed to explain their origins and behavior. Shapiro et al. (1997) proposed a diffusion model that explains the spatiotemporal evolution of hypocenter migration by pore-pressure diffusion from a point pressure source. Following their convention, we define the triggering front as the outermost hypocenter area surrounding the initiation point of the earthquake swarm. Previous studies have found that hypocenter migration associated with volcanic, geothermal, and fluid injection earthquake swarms adhere to this diffusion model (e.g., Parotidis et al. 2003;Shelly et al. 2016).
The triggering front's diffusivity thus represents an essential observable factor in interpreting earthquake swarms. However, previous studies employing the diffusion model lack a unified approach for selecting the triggering front and an end time for model fitting. This study addresses a set of criteria used to interpret several earthquake swarms occurring in northeastern Japan. The analyzed swarms were triggered by the 2011 Off the Pacific Coast of Tohoku earthquake (M w = 9.0; hereafter referred to as the Tohoku earthquake). The analysis identifies a relationship between the diffusivity in the early stage of hypocenter migration and the swarm duration.

Earthquake swarms in northeastern Japan triggered by the Tohoku earthquake
We investigated earthquake swarms in northeastern Japan triggered by the Tohoku earthquake (Fig. 1). Swarms triggered around Moriyoshi (Fig. 2a), Kakunodate (Fig. 2b), Gassan (Fig. 2c), Sendai (Fig. 2d), and Yonezawa-Kitakata (Fig. 2e) also showed hypocenter migration. Previous studies have suggested that these swarms located in the stress shadow of the Tohoku earthquake (Terakawa et al. 2013) were driven by pressurized crustal fluids upwelling from the mid-crust in response to the stress change by the Tohoku earthquake (e.g., Okada et al. 2015). Though the swarms have a similar spatial extent of about 4-6 km, the swarm duration is highly variable. Swarms around Kakunodate, Gassan, and Sendai ceased within a year while swarms in Moriyoshi and Yonezawa-Kitakata continue even 10 years after the triggering event. The latter swarms represent some of the most long-lived ones ever documented.
This study analyzed hypocenter data from the unified catalog of the Japan Meteorological Agency (hereafter referred to as the JMA catalog) from March 2011 to December 2019. To ensure completeness of the detectability of earthquakes, we analyzed earthquakes with a magnitude greater than or equal to 1.0. We constructed independent datasets for five target areas. For the Yonezawa-Kitakata area, we divided earthquakes into subclusters labeled A to D (Fig. 2e). Following an approach in Yoshida and Hasegawa (2018b), we subsequently divided the clusters of Yonezawa-Kitakata A and B into two sequences of earthquakes occurring in the first ~ 50 days after the Tohoku earthquake and over an indefinite period afterward (Additional file 1: Figure S1).

The duration of earthquake swarms
Earthquake swarm activity usually shows no clear termination, but instead dies off gradually. To evaluate the duration of earthquake swarms, we defined apparent duration EVT90 [days] as the elapsed time between the start time of the swarm and when the cumulative number of events reached 90% of the total number. This definition works for ongoing swarms as well due to their considerable decrease in recent seismic activity. Similarly, as the end time for fitting the diffusion model, we also defined EVT30 [days] that is the elapsed time from the start of the swarm to when the cumulative number of events reached 30% of the total number. The start time of event counting is the origin time of the first earthquake that occurred in each swarm.

Estimation of the diffusivity of hypocenter migration
We estimated the diffusivity of hypocenter migration by fitting the diffusion model proposed by Shapiro et al. (1997): where r [m] is the distance from the pressure source to the triggering front, t [s] is the elapsed time from the start of diffusion, and D [m 2 /s] is the diffusivity. By squaring Eq. (1), the diffusivity D is estimated by linear regression between squared distance and elapsed time. We prepared the distance data as the 90th percentile distance for the moving time bins that contain 20 earthquakes with overlapped 10 events.
Estimates of time and distance require the definition of spatial and temporal origin. Previous studies have typically set the origin as the origin time and hypocenter of the first earthquake. Because the true diffusion origin is unknown, this assumption may bias results and avoid suitable comparisons among different swarms. This study used a grid search to identify a spatial origin that gives the minimum RMS misfit between the observed and predicted triggering front. We set the grid interval of 0.5 km. We also tried to search the start time of the swarm in the temporal direction; the results were almost the same as those obtained by searching only in the spatial direction. However, in some cases, we observed unreasonable results in which the time origin became exceptionally far in the past. Thus, we decided to fix the start time of the swarm as the origin time of the first event in each swarm.
The selection of the end time for model fitting affects the estimate of the diffusivity. Sometimes the distancetime variation within EVT90 cannot be explained by Eq. (1) because the hypocenter migration seems to stop at a certain distance (e.g., Fig. 3c). Therefore, we employed EVT30 [days] as the end time for fitting.   Additional file 1: Figure S1 in supplementary information gives results for other sequences. The best-fit synthetic curve for the triggering front well explains the initial part of each swarm. It shows almost good fitness for the entire sequence as well. For comparison, we also calculated results using the entire period and assigning the spatiotemporal origin of migration to the first earthquake (Additional file 1: Figure S2). The approach described in the former section, using EVT30 and grid-searched spatial origin, clearly provided superior results (Fig. 3). For some swarms, spatiotemporal evolution for the triggering front seems to stop for a certain distance and elapsed time (for example, Fig. 3f ). Estimated diffusivities ranged approximately from 0.01 to 2 m 2 /s. These values did not differ in the order of magnitude from those reported by previous studies for the same swarms (Kosuga 2014; Okada et al. 2015;Yoshida and Hasegawa 2018a, b). Figure 4a summarizes the relation between the diffusivity (D) and EVT90 (apparent duration). Horizontal bars attached to the symbols show the 2 σ errors of the diffusivity estimate. A clear negative correlation, with a correlation coefficient of − 0.85, exists between the diffusivity and the duration, even considering the 2 σ errors. This relationship adheres to a power-law of EVT90 ∼ D −0.5 to D −1.0 . Thus, swarm duration decreases with increasing diffusivity coefficients of hypocenter migration.
We assessed the uncertainty of the diffusivity differently by examining their range that gives similar RMS values of model fitting. The fitting RMS for each spatial origin of the diffusion is shown in Additional file 1: Figure S3. The first event is generally outside the range of small RMS values, showing the necessity of searching the spatial origin. We additionally assessed the effect of the cut-off percentage N that defines the swarm duration EVT-N (Additional file 1: Figure S4). The horizontal bars in Additional file 1: Figure S4 indicate the range of the diffusivity that gives less than 1.5, 2.0, and 3.0 times the minimum RMS. We found that the negative correlation holds even the uncertainty of the above parameters is considered.

Hypocenter accuracy and catalog completeness
The diffusivity estimated from Eq. (1) depends on the hypocenter catalog from which the distance is calculated. The JMA catalog represents the complete data available. On the other hand, the Japan Unified high-resolution relocated Catalog for Earthquakes (JUICE) by Yano et al. (2017) provides better location accuracy. JUICE lists locations for shallow earthquakes relocated by the double-difference earthquake location method (Waldhauser and Ellsworth 2000) and gives tightly centered epicentral distributions for the target swarms than those from the JMA catalog (Additional file 1: Figure S5). However, the JUICE catalog covers a shorter period up to 2012.
We performed the same analysis described above using data from the JUICE catalog to examine the dependence of results on the hypocenter catalog (Additional file 1: Figure S6). Again, a clear negative correlation appears between the diffusivity and the duration (Additional file 1: Figure S7). The diffusivities from the JUICE catalog are slightly different from those of the JMA catalog; however, they do not influence the overall correlation. This suggests that the duration-diffusion relationship may be biased to a smaller exponent. Note that the duration of some earthquake swarms was shortened in the JUICE catalog because of its abridged time coverage.
We investigated whether catalog completeness influenced the observed correlation. In this study, we set the minimum magnitude to 1.0, assuming the earthquakes larger than or equal to this threshold were collected entirely in the catalog. Re-analysis with a more conservative threshold of 1.5 results in a similar negative correlation between the diffusivity and the duration (Additional file 1: Figure S8).

Comparison with other swarms from various localities and tectonic environments
The relationship between the diffusivity of hypocenter migration and swarm duration shown in Fig. 4a is for earthquake swarms in northeastern Japan. We additionally investigated the relation under various tectonic environments from the result of previous studies (Saccorotti et al. 2002;Parotidis et al. 2003; Hill and Prejean  . We used swarm durations either cited or estimated from the time series in the literature. If more than one diffusivity was reported, we treated them as independent data. Figure 4b shows the diffusivities (D) and durations reported by previous studies. We found a similar negative correlation as shown in Fig. 4a for swarms in northeastern Japan. Even if we group the data into different tectonic categories such as volcanic, non-volcanic, subduction zone, and non-subduction zone, the negative correlation holds for all groups (Additional file 1: Figure  S9). This fact implies the negative correlation is a global relationship regardless of the tectonic environment. The distribution is best fitted by either a straight line with a single power exponent (Fig. 4a) or two lines that intersect at around (D, EVT90) = (0.6, 30) (Fig. 4b). Above this point, lines exhibit power exponents greater than 1.0. The population adhering to these parameters may represent faster temporal decay, but the elevated parameters may also arise as an artifact of uncertainties in the swarm duration.

Interpretation of the negative correlation between diffusivity and swarm duration
The negative correlation between the diffusivity and swarm duration suggests that rapidly migrating swarms tend to cease within a short period. According to Shapiro et al. (1997), the diffusivity is proportional to the hydraulic permeability and inversely proportional to the dynamic viscosity of pore fluid. These two parameters thus influence swarm duration. Hauksson et al. (2019) reported diffusive hypocenter migration with a very low diffusivity of 0.006-0.01 m 2 /s for the long-lived swarms (greater than ~ 2 years) around Cahuilla Valley pluton, California. They cited the low permeability geologic environment of the Mesozoic plutonic rocks as a reason for longer swarm durations. Literature sources have also described cases of short-lived swarms in high-permeability environments. For example, Yukutake et al. (2011) reported a complex spatiotemporal pattern outlined by several vertical, thin, planar hypocenter distributions in Hakone caldera, Japan. They estimated the diffusivities of ~ 0.5-1.0 m 2 /s and inferred the existence of high-permeability fault zones occupied by highly pressurized pore fluid to account for the relatively short swarm duration of a few days. Shelly et al. (2016) reported hypocenter migrations with a high diffusivity (~ 1.0 m 2 /s) for a short earthquake swarm that ceased within a few days at Long Valley caldera, California. They suggested the contribution of low-viscosity fluids likely composed of water or CO 2 and supplied by underlying magma chambers.
The above reports show that swarm duration declines with low viscosity of fluid or high permeability of the crust. These two parameters, fluid viscosity and crustal permeability, may influence pore fluid pressure, enhancing seismicity. A quantitative estimate is difficult for these parameters. Instead, we can evaluate their effect through the diffusivity from the migration speed of the triggering front of seismicity.
The relationship between the diffusivity and swarm duration can be described by a power-law with an exponent of − 0.5 to − 1.0. In the following argument, we need to assume that the duration of the triggering front diffusion (t in Eq. (1)) equals the swarm duration. According to Eq. (1), t ∝ D −0.5 to D −1.0 derives r ∝ D 0−0.25 . This suggests that even if an earthquake swarm continues for a long time, its source area does not expand unlimitedly. For example, in the Moriyoshi area (Fig. 2a), the swarm duration is more than 10 3 days (Figs. 3a, 4), but the source area is within a few kilometers. Other analyzed swarms also showed a similar spatial size. These facts are consistent with the above discussion. For an infinitely large reservoir volume, the swarm persists regardless of high or low diffusivities. Thus, the limited source size implies either a relatively consistent volume of fluid reservoirs or a roughly constant fluid supply rate in the crust, regardless of swarm locality. Some swarms last for a long time with no apparent hypocenter migration (Fig. 3). Here, we show a possible reason to prevent the hypocenter migration. If the change in pore fluid pressure follows the diffusion equation, the pore fluid pressure decreases with the distance from the diffusion origin. Thus, the effective normal stress becomes insufficient to slip faults beyond a certain distance. Hainzl (2004) simulated an earthquake swarm by the chained earthquake occurrence driven by the stress fluctuation due to the afterslip of each earthquake fault under high pore fluid pressure. In this case, seismic activity continues even after the halt of hypocenter migration. However, since the simulation shows the temporal decay of seismicity, the above mechanism might not be enough to explain the very long swarm sequence. Possible additional factors of the long-lived swarms may be the rapid fault healing and the multiple fluid intrusion to the source area.
This study found a negative correlation between the diffusivity of the hypocenter migration and the swarm duration. It has been pointed out that the diffusivity relates to crustal permeability and fluid viscosity (Shapiro et al. 1997). Some long-lived swarms that last after the halt of hypocenter migration require other factors such as stress fluctuation by afterslip of each earthquake, fault healing, and continued fluid supply. For a detailed discussion, it is also necessary to consider other factors such as aseismic slip that initiate swarm activity (e.g., Vidale and Shearer 2006;Lohman and McGuire 2007), the temporal change in the crustal permeability (e.g., Shelly et al. 2015), and complex 3D fault architecture as a barrier to fluid intrusion (Ross et al. 2020). While many factors make earthquake swarms complex phenomena, some studies have distinguished multiple mechanisms in swarm sequences (e.g., Hatch et al. 2020). Quantitate estimates of the above factors, and the evaluation of their interaction are the next step of this research.

Conclusions
We estimated the diffusivity of hypocenter migration and the duration of earthquake swarms triggered by the Tohoku earthquake. We employed a consistent method to estimate elapsed time and distance from a spatiotemporal origin for the appreciation of the diffusion model and comparison of the diffusivity among different localities. We defined the distance as the 90th percentile distance from the spatial origin determined by a grid search. End times were selected as the point at which the cumulative number of earthquakes reached 30% of the total events (EVT30). Similarly, we defined EVT90 as the duration of swarm activity. This approach delineated a clear negative correlation between the diffusivity and swarm duration, which adhered to an approximate power-law of EVT90 ∼ D −0.5 to D −1.0 . This relationship also holds for other swarms under various tectonic settings throughout the world. These results suggest that the product of the diffusivity and swarm duration approaches constant values, which implies a similar volume of fluid reservoirs in swarm areas regardless of locality. These results also indicate that diffusivity, related to crustal permeability and fluid viscosity, represents a critical factor in controlling swarm duration.
Abbreviations EVT: Event number; JMA: Japan Meteorological Agency; JUICE: Japan Unified hIgh-resolution relocated Catalog for Earthquake; RMS: Root mean square.

Supplementary Information
The online version contains supplementary material available at https:// doi. org/ 10. 1186/ s40623-021-01480-7.  Figure S3. Spatial distribution (plan view and E-W cross-section) of RMS misfit to the diffusion model with spatial origin assigned to the 3D grids in each panel. Panels (a) to (i) show the result for all the areas shown in Fig. 1 Figure S4. Scatter plot between the diffusivity and swarm duration to assess their uncertainty for the areas shown in Fig. 4a. The color of diamonds represents N value for EVT-N [days] as the duration of each swarm, which are common for all panels. Horizontal bars attached to the diamonds in the panel (a), (b), and (c) indicate the range of diffusivity with allowed normalized RMS by the minimum RMS of 1.5, 2.0, 3.0, respectively. Two lines proportional to D −0.5 and D −1.0 are shown for reference. Figure S5. Hypocenter distribution plotted by the JUICE catalog (Yano et al., 2017). Colored circles denote epicenters of earthquakes (M ≥ 1.0) that occurred from 11 March 2011 to 31 December 2012. The color and size of each circle represent depth and magnitude, respectively. Blue stars represent the first event in each area. Figure S6. Same as Figure S2 but using the JUICE catalog and the spatial origin determined by a grid search. Gray crosses denote distance and time calculated from the spatial origin and the origin time of the first event. Orange circles denote the data used for the fitting to the diffusion model.  Table S1. The duration (EVT90) and diffusivity (D) with its range ( D − 2σ and D + 2σ ) for the areas shown in Fig. 2. These values were estimated for earthquakes with M ≥ 1.0 in the JMA catalog. Table S2. Same as Table S1, but for earthquakes in the JUICE catalog. Table S3. Same as Table S1, but for earthquakes with M ≥ 1.5 in the JMA catalog. Table S4. Same as Table S2, but for earthquakes with M ≥ 1.5 in the JUICE catalog. Table S5. The durations and diffusivity (D) from previous studies. For the same areas we investigated, we employed our estimates of duration (colored by blue) using earthquakes with M ≥ 1.0 in the JMA catalog.