{'fraction_non_alphanumeric': 0.03288324066719619, 'fraction_numerical': 0.020135027799841144, 'mean_word_length': 4.86969696969697, 'pattern_counts': {'":': 0, '<': 2, '<?xml version=': 0, '>': 6, 'https://': 0, 'lorem ipsum': 0, 'www.': 0, 'xml': 0}, 'pii_count': 0, 'substrings_counts': 0, 'word_list_counts': {'cursed_substrings.json': 0, 'profanity_word_list.json': 0, 'sexual_word_list.json': 0, 'zh_pornsignals.json': 0}}
Injection induced seismicity size distribution dependent on shear stress

Like natural seismicity, induced seismicity caused by fluid injection also shows a power law size distribution, and its gradient b-value (ratio of small to large earthquakes) is often used for seismic hazard analysis. Despite well-known relationship that b-value is negatively correlated with differential stress for natural earthquakes, there is no understanding of the physical causes for b-value variations in injection-induced seismicity in the scale where the differential is nearly constant. We investigate a b-value dependence on the relative shear stress acting on existing fractures and show that the seismicity occurring along existing fractures with high shear stress have significantly lower b-values than does that associated with lower shear stress fractures. The b-value for injection induced seismicity is dependent on relative shear stress on faults. Our results provide a novel physical explanation for the b-value variations of induced seismicity.


Introduction
Earthquake size distribution follows a power law (Gutenberg and Richter, 1944) known as the Gutenberg-Richter (GR) relationship over a large range of earthquake scales from laboratory to megathrust, in addition to anthropogenic earthquakes. The gradient of the power law, the b-value or slope of that power law is obtained from the relation log10(N)=a-bM (Gutenberg and Richter, 1944), where N is the number of events having magnitude at M, and a is a scaling factor. The low b-value means high relative number of larger earthquakes and vice versa. The b-value is often monitored for earthquake prediction and hazard analysis (Nanjo and Yoshida, 2018). Recent interest on this topic has led to several physical explanations for b-value variations (Ide et al., 2016;Nishikawa and Ide, 2014;Schorlemmer et al., 2005). Insights from various studies ranging from laboratory experiments (Amitrano, 2003;Goebel et al., 2013;Scholz, 1968) to natural earthquakes (Mori and Abercrombie, 1997;Spada et al., 2013;Wiemer and Wyss, 1997) have led to the assertion that b-value negatively correlates with differential stress (Scholz, 2015(Scholz, , 1968. Thus, spatial and temporal variations of b-value are generally understood to reflect spatial and temporal variations in the stress loading level. Recently, the importance of induced seismicity to earthquake science and to society has increased along with the demand for subsurface development using fluid injection (Ellsworth, 2013;Evans et al., 2012;Majer et al., 2007). In the case of induced seismicity, bvalue and its variation have also often been observed and used in various assessments of seismic risk (Bachmann et al., 2012(Bachmann et al., , 2011Langenbruch et al., 2018;Langenbruch and Zoback, 2016;McGarr, 2014;Shapiro et al., 2010). Though variations in b-value have often been observed in injection induced seismicity data (Bachmann et al., 2012), unlike natural seismicity, there is no general understanding of the physical causes for these variations and we have no basis for evaluating the importance of such variations (Gaucher et al., 2015). We cannot explain the b-value variations during induced seismicity as being caused by variations in differential stress as we do for natural earthquakes because the time scale, tectonic setting and the stress state of the fractures are different from natural earthquakes. The region activated by fluid injection is on the scale of several km and its depth rarely exceeds 5 km. A significant variation in differential stress is not expected on these scales. The static stress change caused by earlier induced seismicity can change stress state locally, but in the induced seismicity case, these changes are small so their effect on induced seismicity is much less significant than pore pressure (Catalli et al., 2013), which does not alter differential stress.
Therefore, the question we address is what is the physical explanation of b-value variation for injection induced seismicity?
During hydraulic stimulation, shear slip occurs along fractures having a range of orientations. The geometrical relationship between an existing fault and in-situ stress determine the stress state of the fault, which is the balance of shear stress and normal stress. Therefore, there is some variation of shear stress acting on the existing faults which are candidate locations for induced seismicity according to the geometrical relationship to the insitu stress and fault orientation (see methods). The required pore pressure increases to destabilize existing fault varies by its stress state.
In many b-value studies in natural earthquake settings (Mori and Abercrombie, 1997;Nanjo and Yoshida, 2018;Scholz, 2015;Schorlemmer et al., 2005) and laboratory experiments (Amitrano, 2003;Goebel et al., 2013;Scholz, 1968), b-value variation has been found to correlate with the strength, differential stress, and shear stress. First, we might consider the correlation to the stress relative to shear strength (Scholz, 1968) in the form of Coulomb failure stress (CFS) or pore pressure increase. However, if we assumed that all shear failure is caused by pore pressure and that the static stress change is negligible condition, CFS changes only due to the pore pressure increase in injection-induced setting.
Pore pressure promotes shear slip by reducing the effective normal stress but it does not affect shear stress. In other words, pore pressure moves the Mohr stress circle closer to the failure line, but pore pressure does not change the differential stress, which is the diameter Mohr stress circle. Therefore, for induced seismicity studies, we cannot consider CFS or pore pressure change as a substitute for differential stress change as used in natural earthquake studies.
Instead of using differential stress or CFS for a natural earthquake case, we investigate the relationship between b-value and the shear stress on existing faults that slip.
We consider that the variation of shear stress acts as driving potential controlling the size of the slip. We use an integrated analysis of microseismicity data and measured in-situ stress information to infer the shear stress along the existing faults that slip. Also, this is the first case study where the b-value dependency on the stress is confirmed for field seismicity and directly measured in-situ stress data beyond the laboratory scale.

Study field: Basel, Swizterland EGS project
We use the induced seismicity data recorded during the Basel Enhanced Geothermal System  hydraulic stimulation since high-quality microseismic data and in-situ stress information are available. Hydraulic stimulation was performed in the Basel urban area, in December 2006. Basel is located near the border of Switzerland, Germany and France, where is the southern end of the Rhine Graben. This area is thought to be the area of highest geothermal potential in northern Europe ; where other classic EGS projects were operated (Evans et al., 2005). Injected water penetrated into the granite basement throughout the permeable zone in an open-hole section (4603-5000 m) of the injection well; Basel-1 (true vertical depth; TVD 5000 m) . The flow rate was increased step by step and finally reached to 3300 L/min on the 5 th day of the stimulation at which time the wellhead pressure (WHP) was 29.6 MPa. Hydraulic stimulation successfully reactivated a number of existing fractures causing microseismicity from the start of injection. The locations of the microseismicity indicated that an artificial geothermal reservoir was created successfully by enhancing the permeability in the vicinity of the injection well.

Microseismic data
The operating company deployed a microseismic monitoring network that consisted of six borehole stations and one temporary station in the injection well. All stations had 3component seismometers. Microseismic monitoring station locations are shown in Fig. 1(d).
Using the waveforms from 6 downhole stations, we determined the hypocenters of 2800 microseismic events from manual picks having high confidence (Asanuma et al., 2008). Our target events occurred from the start of the injection until after 6 months from the start of the injection Mukuhira et al., 2013). We then performed a cluster analysis based on waveform similarity (Moriya et al., 2003) and consequently relocated hypocenters of 60 % of the microseismic events (Asanuma et al., 2008) with the double difference method (Waldhauser and Ellsworth, 2000) (Fig. 1). The seismically activated region is roughly a cube 1,000 m on a side and the reservoir depth is 4,000~5,000 m in the granite basement ( Fig. 1).

In-situ stress data
The orientation of SHmax (maximum horizontal stress) has been well constrained by borehole breakout analysis for well Basel-1. The mean orientation of SHmax (maximum horizontal stress) in the basement granite is N144°E±14° Evans, 2015, 2009).
Stress magnitude was also estimated from borehole breakout analysis using Schlumberger ultrasonic borehole televiewer (UBI) logging run data taken from the granite section (2569-4992 m) of Basel-1 (Valley and Evans, 2015). The consistency of the stress information results was investigated by comparing with other geophysical observations. The consensus estimates of in-situ stress magnitude as a function of depth z in meters from the free surface are (Valley and Evans, 2015): (1) %"&' = 0.01990 − 17.78 (2) where stress is in MPa. The stress magnitudes given by equations (1-3) are shown in Fig. 2(a).
SHmax has a very low variation with depth, and this causes a transition from a strike slip to a normal faulting stress state regime at around 4800 m depth. Six depth dependent Mohr stress circles at 200 m intervals in the 4,000-5,000 m depth range are presented in Fig. 2(b). The Mohr stress circles become smaller with depth meaning that differential stress decreases with depth. However, it varies from 60 MPa at 4000 m to 40 MPa at 5000m, which is quite small and not sufficient to cause a significant change in differential stress variation since the differential stress discussed in natural earthquake is much larger than this value (Scholz, 2015).

Fault orientation data
Orientations of fault planes are obtained from the focal mechanisms estimated for about 100 of the relatively large events using data from the Swiss natural seismicity network (Deichmann and Giardini, 2009;Terakawa et al., 2012). We select the slip plane using geomechanical information and confirmed them with the hypocentral distribution of other microseismicity in the vicinity of those lager events (Mukuhira et al., 2017). The poles to the selected fault planes are shown in Fig. 3(a).
To add more fault orientation information, we then extracted the fault orientation information from the result of the seismic cluster analysis. Clustering analysis of microseismic events is based on waveform similarity obtained from coherence. We expect that the hypocenters of events having similar waveforms are located close to each other and that they have the same source mechanism (Moriya et al., 2003;Shelly et al., 2016). The locations of events in many clusters delineated a planer or streak shape that are visualized in

Geomechanical analysis
Using the in-situ stress model and orientation of an existing fault, we estimate the shear stress on the fault; shear stress τ and normal stress σn working on a fault plane with a known orientation can be determined from the magnitude of the principal stress Sn (n=1, 2, 3; S1 > S2 > S3) and their orientations, and geometry of the fault plane as: where an (n=1, … ,3) are direction cosines from each principal stress to the vector normal to the fault plane.
Due to the state of stress transitions and the differential stress perturbation with depth, the absolute value of shear stress can be influenced by the differential stress and state of

b-value analysis
We employ standard and well accepted method to estimate b-value. In following section, we divide the catalog by NSS value and then estimate Mc for each sub-catalog using the entire magnitude range method (EMR) (Woessner and Wiemer, 2005). We estimate bvalue with maximum likelihood method (Aki, 1965;Utsu, 1999) and binning size of magnitude in b-value estimation is 0.1. All computations were performed using ZMAP (Wiemer, 2011). To test the discrepancy of the b-values in sub-catalogs divided by NSS, we employed the Utsu's test (Utsu, 1999). The equation is in Appendix.

Correlation between NSS and magnitude
First, we simply examine the correlation between the magnitude of seismic events and NSS (Fig. 4). We observe a correlation between the magnitude and NSS with the large events showing a tendency to occur along the faults with higher NSS. Smaller events also occur along faults with high NSS, so the correlation between event size and NSS is not a simple one. In order to investigate the relation in a statistical sense, we introduce the concepts of bvalue and investigate the relationship between b-value and shear stress.
We divide the catalog by NSS threshold value into those with NSS higher than the threshold as the higher group: HG and those with lower NSS as the lower group: LG (see Fig.   4). We plot the GR magnitude frequency distribution of the events and estimate the Mc with the EMR method (Woessner and Wiemer, 2005), and the b-value with the maximum likelihood method (Aki, 1965;Utsu, 1999). Then b-values for both groups are compared in We first look the b-value dependency of HG using NSS thresholds of 0.71, 0.76 and 0.87, which correspond to the NSS of well oriented fractures having friction coefficients of 0.6, 0.85, and 1.0 respectively (Fig. 5(a)-(c)). We observe that b-value for NSS > 0.71 is higher than that for NSS > 0.76 and 0.87, but b-value for NSS > 0.76 is lower than that for NSS > 0.87. Meanwhile, for the LGs (Fig. 3d-f), b-value systematically decreases with increasing NSS threshold from 1.31 for NSS < 0.71 to 1.12 for NSS < 0.87. The other observation made from the comparison of b-values for the same NSS threshold (e.g. Fig. 5(a) and (d)) is that the HGs show significantly lower b-values than those of LGs.

b-value dependency to normalized shear stress
We further investigate the relationship observed in Fig. 5, by determining b-value for a range of NSS thresholds. We estimate Mc for each sub-catalog defined by NSS threshold and estimate corresponding b-values (Fig. 6). Overall dependency of b-values on NSS threshold is clearly observed for every 0.05 increment in the threshold for HGs (Fig. 6(a)) and LGs (Fig. 6(b)). We estimate b-values when a sub-catalog contains more than 100 events LGs become stable as the threshold increases because of the increase in the number of events used for the estimation (Fig. 6(e)). Consequently, both cases show a bvalue dependency on NSS threshold. Meanwhile, we observe that the b-values for HG groups are significantly lower than those for LG groups with the same threshold. The significance of the differences between them is evaluated with statistical test (Utsu, 1999) and they show ( Fig. 6(d)) that differences are highly significant (AIC = 5 is highly significant line).

Uncertainty
We discovered the b-value dependency on NSS based on the b-value and stress analysis using the in-situ stress model, which we assume laterally constant stress condition. We can observe that both b-values for HG and LG decrease with NSS, though the bvalue for HG slightly decreases. In addition, b-values for HG are lower than b-value for LG except for the points more than 0.9, where we cannot estimate b-value for HG due to a small number of the events. We can confirm three main observations we had for the observed stress model, even when we employ the stress model considering the stress heterogeneity in magnitudes and orientation.
We further examined b-value dependency on shear stress with Mc estimated with another method (maximum curvature method), constant Mc, and intentionally overestimated Mc and confirmed our findings (Supplementary Fig. 1 -5). We also performed a statistical simulation of the influence of the approximately 1400 events for which we cannot determine NSS. We found that our observations and conclusions are not changed (Supplementary Fig. 6 -8).

Discussion
Our main observations are that grouping events using both upper and lower bounds on NSS show b-value variation associated with changes in NSS threshold and that the HG groups always show lower b-values than those estimated for LG groups with the same NSS value. Our observations demonstrate a clear b-value dependency on NSS for injection induced seismicity. The discovered b-value dependency seems to be linked to the general interpretation of b-value dependency on differential stress, in terms of the dependency on the stress. However, the interpretation of our discovered dependency is slightly different from the general interpretation of b-value dependency.
We first consider the general interpretation of b-value dependency on differential stress using the Mohr stress circle. The Mohr stress circle shown in Fig. 8(a) shows how the diameter of the circle changes with increasing differential stress. The b-value negatively correlates with differential stress, i.e., the diameter of the Mohr stress circle. Many studies of the b-values of natural earthquakes found indirect indicators of differential stress such as depth (Mori and Abercrombie, 1997), stress regime type (Schorlemmer et al., 2005), and subduction plate age (Nishikawa and Ide, 2014). Other b-value studies from laboratory experiments found a b-value reduction with increasing loading stress (Goebel et al., 2013;Scholz, 1968). These studies implicitly or explicitly consider that most of the earthquakes occurred along the representative fault or fractures that are well oriented for slip under the given stress condition. Thus, we indicate the stress state of well oriented faults with dots on the Mohr stress circle in Fig. 8(a). As differential stress increases, the shear stress along well oriented faults increase as well. So, the general interpretation of b-value dependency on differential stress can be equivalent to the b-value dependency on the shear stress (here absolute shear stress). Some of the b-value studies discussed b-value dependency with shear stress (Ide et al., 2016;Nishikawa and Ide, 2014;Yoshida et al., 2017).
Let us recall our geomechanical setting considered in this study. On the reservoir scale, differential stress does not change significantly with time or space. Therefore, the differential stress is mostly constant, leading to a constant diameter of the Mohr stress circle ( Fig. 8(b)). In the reservoir, various type of the existing fault is there and stress state of them also varies according to the geometrical relationship between faults and in-situ stress. The existing faults can have shear slip according to the pore pressure increase by fluid injection to destabilize it. During the stimulation, Coulomb failure criterion line moves to the right in relation to the Mohr stress circle due to the injection pressure increases (this is equivalent to Mohr stress circle shifting to the left in the case where the horizontal axis is effective normal stress). As this figure shows, pore pressure increase does not affect to differential stress nor shear stress. Pore pressure increase just prompt the failure but the potential of shear slip (differential stress) remains as before. Then, the faults on the left of the failure line reach a critical state and fail. Here, we can see the difference in the process to reach the failure. In natural earthquake case, differential stress increase leads the failure because of increasing loading, so Mohr stress circle getting bigger. Meanwhile, induced seismicity case, pore pressure increase leads the failure because of pore pressure increase which decreases effective stress, so Mohr stress circle shift to the failure line.
As we explained NSS is the relative height of Mohr stress circle for the point that shows the stress state of the fault. Then, stress states of the high NSS faults are situated in the higher portion and critical part of the Mohr stress circle (black shaded) and these are the events that lead a lower b-value ( Fig. 8(b)). Thus, when more events occur within the higher shear stress fault, the b-value decreases. Conceptually, this relationship is shown by the upward pointing yellow arrow in Fig. 8(b) suggesting that the more events from high shear stress fault cause the lower b-value. This is consistent with the conceptual explanation of bvalue dependency on applied shear stress with upward yellow arrow shown in Fig. 8(a) suggesting the (absolute) high shear stress cause the lower b-value. So, we could explain the ubiquitous phenomena of b-value variation but in the case of injection induced seismicity, where tectonic stress condition and triggering process are totally different, with a different theory compared to the general understanding in seismology. At the end, our explanations are quite consistent with the insight of b-value variation. Also, the larger event occurring from relatively higher shear stress fault is reasonable since the higher shear stress fault has much potential to the size of shear slip.
One might consider that pore pressure increase correlates to the b-value variation as the frictional strength does (Scholz, 1968). CFS is also used as the stress relative to the failure in natural earthquake case. Pore pressure and CFS are essentially equivalent in the injection induced seismicity setting since shear and normal stresses are constant, and pore pressure is the only variable parameter in the definition of CFS (Beeler et al., 2000). As we discussed, the pore pressure and CFS cannot change the differential stress and shear stress.
So, b-value variation interpretation with pore pressure is not consistent with general interpretation. This hypothesis should be investigated with data analysis in future.
The largest event in this field and several other large felt events occurred in the deeper part of the reservoir where the differential stress is slightly smaller than at shallower depths and where the state of stress is normal fault slip type. From previous studies, the larger events would thus be expected to have occurred within a strike slip stress regime rather than the normal fault slip stress regime (Schorlemmer et al., 2005). Our results show opposite tendency suggesting that the differential stress interpretation behind of the b-value dependency on stress regime is not the case. In addition, we do not observe a significant preference for larger events to occur shallower in the reservoir where differential stress is larger (Mukuhira et al., 2017). From these observations, there is a tendency for the occurrence of large events to be correlated with normalized shear stress rather than differential stress.

Conclusions
We conclude that, on the reservoir scale, a lower b-value associated with fluid injection can be attributed to the occurrence of events along fractures that are oriented to have high shear stress. Our conclusion provides the first clear physical interpretation of variation in b-values for induced seismicity where earthquakes occur along fractures of various orientations in a constant regional stress condition. Also, this is the first case to show the bvalue dependency on the directly measured stress data beyond the laboratory scale. This helps to fill the gap in our understanding of b-value between the laboratory scale and earthquakes along plate boundaries or in the continental crust, clarifying that even for the microseismicity that occurs within one local region, b-value can be related to shear stress level.
Our findings contribute to our ability to assess seismic risk and to mitigate it when undertaking subsurface development associated with fluid injection since the b-value reduction can to some degree be controlled by the injection pressure. With an increase of injection pressure, the Coulomb failure criterion line shifts to the right (Fig. 5(b)

Significance test of difference in b-values from sub-catalogs (Utsu's test)
For the two groups of earthquakes having N1 and N2 events and b-values b1 and b2 respectively, the significance of the difference in b-value is tested with Akaike's Information Criterion (AIC). ΔAIC is given by ΔAIC=5 means that the difference in b-values between the two data sets is highly significant and ΔAIC=2 means the difference is significant (Aki, 1965).
two data sets is highly significant and ΔAIC=2 means the difference is significant (Aki, 1965).    Red dots correspond the events having FPS by SED and blue dots are the events whose fault geometry was estimated by cluster analysis. As an example, the case of threshold = 0.76 is shown. The events are subdivided into two groups; those having NSS higher than the threshold identified as HG (0.76) and those with NSS lower than the threshold identified as LG (

Supplementary information
for Injection induced seismicity size distribution dependent on shear stress