Increased likelihood of induced seismicity in highly overpressured shale formations

S U M M A R Y Fluid-injection processes such as disposal of saltwater or hydraulic fracturing can induce earthquakes by increasing pore pressure and/or shear stress on faults. Natural processes, including transformation of organic material (kerogen) into hydrocarbon and cracking to produce gas, can similarly cause fluid overpressure. Here, we document two examples from the Western Canada Sedimentary Basin where earthquakes induced by hydraulic fracturing are strongly clustered within areas characterized by pore-pressure gradient in excess of 15 kPa m−1. Despite extensive hydraulic-fracturing activity associated with resource development, induced earthquakes are virtually absent in the Montney and Duvernay Formations elsewhere. Statistical analysis suggests a negligible probability that this spatial correlation developed by chance. This implies that, in addition to known factors such as anthropogenic pore-pressure increase and proximity to critically stressed faults, high in situ overpressure of shale formations may also represent a controlling factor for inducing earthquakes by hydraulic fracturing. On a geological timescale, natural pore-pressure generation may lead to fault-slip episodes that regulate the magnitude of formation overpressure.

Large-volume disposal of produced brine (Keranen et al. 2014;Schultz et al. 2014;Weingarten et al. 2015) or hydraulic fracturing (Schultz et al. 2015a(Schultz et al. , 2017;;Atkinson et al. 2016;Bao & Eaton 2016;Lei et al. 2017) can induce earthquakes by increasing pore pressure or stress on faults (Ellsworth 2013;Segall & Lu 2015).Necessary conditions for the occurrence of injection-induced seismicity include a source of elevated pore pressure, a proximal, critically stressed fault and a pathway for fluid pressure to propagate from the injection site to the fault (Ellsworth 2013;Eaton 2018).Extraction of hydrocarbons from unconventional low-permeability reservoirs, such as organic-rich shale, makes extensive use of hydraulic fracturing to increase reservoir permeability.Earthquakes have been observed in association with a small fraction of well completions (Atkinson et al. 2016), yet pre-development assessment is hindered by a paucity of validated predictive models to forecast site-specific seismic hazard.Of particular importance is the need to achieve a better understanding of specific geological factors that impact the likelihood of induced seismicity (Schultz et al. 2016;Ghofrani & Atkinson 2016;Pawley et al. 2018).
Prospective fairways containing pervasive hydrocarbon accumulations in low-permeability (tight) reservoirs occur most frequently within epicontinental basins (O'Connor et al. 2014).Shale reservoir units are often overpressured, a state wherein the fluid pressure within pore spaces exceeds the hydrostatic gradient (∼10 kPa m −1 ).Overpressure may develop during basin evolution as a result of inhibited pore-water expulsion, in the case of undercompacted sediments (i.e.load transfer to the pore fluid) and/or fluid expansion due to hydrocarbon generation and cracking to gas, particularly in cases where rates of thermogenic gas accumulation exceed loss from the formation (Law & Dickinson 1985;Hansom & Lee 2005).From a resource-exploitation perspective, overpressure can be beneficial; more fluid stored in the formation reduces the work required to stimulate hydraulic fractures (Wang & Gale 2009).Consequently, the most prolific parts of a fairway often coincide with areas of highest overpressure (Cander 2012).

DATA A N D M E T H O D O L O G Y
Two major resource fairways in the Western Canada Sedimentary Basin are considered in this study, both prone to induced seismicity from hydraulic fracturing (BCOGC 2014;Bao & Eaton 2016;Mahani et al. 2017;Schultz et al. 2017).The Triassic Montney Formation is an extensive siliciclastic unit deposited in environments ranging from shallow-water shoreface sands to offshore marine muds (Dixon 2000).Reservoir rocks within this trend consist of interbedded shale, siltstone and sandstone layers.Approximately 5500 multistage hydraulically fractured (MSHF) horizontal wells have been completed in the Montney Formation (Fig. 1a).The De- vonian Duvernay Formation consists of mudstones interfingered with tight limestones (Creaney et al. 1994), deposited under anoxic to euxinic conditions that allowed for abundant preservation of organic matter (Stasiuk & Fowler 2004).In this study, approximately 400 MSHF horizontal wells have been used in this emerging unconventional fairway, primarily within two distinct geographic areas (Fig. 1b).Both of these fairways are characterized by anomalous maximum overpressure gradients relative to other shale plays (Supporting Information Table S1).The fairways vary between 1.0 and 5.5 km deep and 20-60 m thick for the Duvernay and 0.5-4.5 km deep and 10-300 m thick in the Montney (Rokosh et al. 2012;Playter et al. 2017).
Fig. 1 shows locations of MSHF horizontal wells, together with earthquakes of magnitude greater than local magnitude (M L ) 2.5 (Schultz et al. 2015b) for the period 2009 January 1 to 2016 January 19.These events are interpreted as primarily induced by hydraulic fracturing (Atkinson et al. 2016).Seismicity and well locations are overlaid on contours of pore-pressure gradient within the middle unit of the Montney Formation (Canadian Discovery 2014) and the Duvernay shale (Canadian Discovery 2015).Throughout this paper, we utilize empirically derived pressure gradients to normalize for the effect of lithostatic loading.The in situ shale formation-pressure gradient maps were contoured using 555 and 77 data points for the Montney and Duvernay, respectively, based primarily on static gradient, bottom-hole buildup and extended leak-off tests in penetrating wellbores (e.g.Nguyen & Cramer 2013).Further details to the input formation pressure and formation depth data set can be found in prior reports (Canadian Discovery 2015).Both fairways are characterized by strong lateral variations in formation pressure, and in both cases the distribution of drilling locations spans areas from near hydrostatic pressure (∼ 10 kPa m −1 ) to overpressure, including regions where the pressure gradient exceeds 15 kPa m −1 .Lateral variations in formation pressure are thought to be related to variabilities in conditions responsible for hydrocarbon maturation (thermal history, overburden and deposited organic content) and ongoing hydrodynamic effects of heterogeneous permeability focusing or impeding fluid flow (e.g.Fox & Soltanzadeh 2015;Seifert et al. 2015).These overpressured regions are typically sealed within the high permeability shales-any communication into adjacent strata, such as high permeability carbonate platforms, readily diffuses through fluid flow (Bachu 1999).Visual correlations are evident between areas of high overpressure and clusters of induced seismic events.Although a causal link between anthropogenic porepressure increase and fault activation is well established (e.g.Segall & Lu 2015), the potential seismogenic influence of formation pressure within sedimentary basins has not been documented.
To explore the statistical robustness of these observations, we have conducted Monte Carlo simulations to characterize the probability that spatial correlations evident in Fig. 1 occurred purely by chance.Fig. 2 shows histograms of formation-pressure gradient, interpolated at earthquake epicentres and well locations.In the case of the Montney fairway, earthquake locations are characterized by a multimodal distribution of formation-pressure gradient with a mean of 15.1 kPa m −1 ; in the case of the Duvernay fairway, the formation-pressure distribution is strongly clustered around a mean of 17.7 kPa m −1 .In terms of formation-pressure gradients extracted at well locations, the Montney fairway has a mean value of 12.2 kPa m −1 with a distribution tailing to high pressure gradients.The Duvernay fairway is characterized by a mean formationpressure gradient of 16.8 kPa m −1 with a tail extending to slightly less overpressured values.Using a bootstrap approach (Efron & Tibshirani 1986), we constructed 10 000 random sets of formationpressure gradient values by choosing a set of well locations equal in number to the earthquake sets.Despite the tendency for drilling of unconventional hydrocarbon resources to concentrate within areas of overpressure (Hansom & Lee 2005;Wang & Gale 2009), none of the bootstrapped samples produced a distribution with a mean value as high as the earthquake distributions.Indeed, the mean value of the bootstrapped distributions is more than one standard deviation less than the earthquakes for the Duvernay and 15 standard deviations less than the earthquakes for the Montney.Consideration of earthquake location uncertainty shows that mean formation pressure at earthquake locations in Fig. 1 is likely underestimated.These results are suggestive of a potential earthquake-overpressure relationship.To further test the validity and robustness of the observation that Montney and Duvernay related induced earthquakes preferentially occur within regions of relatively higher formation overpressure, we consider the Mann-Whitney U test (Gibbons & Chakraborti 2011).This statistical test allows us to determine if overpressures encountered at earthquake and well locations have significantly different median values.In this sense, the Mann-Whitney U test is utilized in a hypothesis testing approach with the p-value indicating the level of statistical confidence that the two distributions are different.We use the standard p-value of 0.05 to assert statistical significance.To do this test, first the earthquake catalogues are declustered to prevent overcounting of repeated events.The declustering process for the Duvernay fairway is based on the centroid averages of 17 uniquely identified spatiotemporal swarms within Schultz et al. (2017).Each centroid has a formation-pressure mean and standard deviation value based on sampling of values within the centroid's error ellipse (∼10 km radius).Montney declustering is spatially based on binning swarms into susceptible regions previously identified including the Caribou, Beg Town, Altares, Septimus and Doe Dawson swarms (BCOGC 2014).These swarms are further subdivided temporally based on well associations which do not overlap within the 3 month association window of Ghofrani & Atkinson (2016).Likely, our declustering process underestimates the numbers of unique swarms, especially due to the large temporal association window used for the Montney.We note that fewer identified swarms will result in a systematically larger p-value in Mann-Whitney U tests, and hence more difficulty to discern significance.Despite this conservative approach, we determine that the spatial distribution of earthquake swarms for both the Montney and Duvernay fairways tends to be encountered in regions of relatively higher formation overpressure (Figs 3a and b).
Lastly, to scrutinize whether this result has been influenced by the quality and uncertainty of the input data we perform bootstrap resampling tests (Efron & Tibshirani 1986).In these bootstrap trials we repeat Mann-Whitney U tests with randomly reduced sample counts in the distribution of well data (90 per cent) and randomly perturbed formation pressures for earthquake centroids based on the magnitude of their standard errors.Centroid location errors are on the order of 5 km, which translates into formation-pressure gradient perturbations on the order of ∼1 kPa m −1 .The 10 000 repetitions of these bootstrap trials indicate p-values with statistical significance in most scenarios (Figs 3c and d).Thus, these results are robust and are unlikely to have been significantly impacted by error in our approach.From the results of these statistical tests, we conclude that there is negligible likelihood that the spatial correlation between higher formation-pressure gradient and locations of earthquakes occurred by chance.

D I S C U S S I O N
The geomechanical principles controlling the non-coincidental spatial association between regions of higher reservoir overpressure and induced seismic activity merit discussion.In general, injectionrelated seismicity is thought to be the result of anthropogenic perturbation to effective stresses acting on a fault, moving the stress condition closer to failure (Segall & Lu 2015;Sibson 2017).Often, incomplete subsurface information precludes a robust understanding of the in situ state of stress and orientations of faults that are in hydraulic communication with injection sites.Instead, injection volumes (and resultant pore-pressure perturbations) become a firstorder proxy for the criticality of pre-existing faults (Sibson 2017).By analogy, natural processes also have the propensity to drive this faulting mechanism, for example, via dynamic processes related to hydrocarbon maturation (Law & Dickinson 1985;Hansom & Lee 2005).Specific to the Duvernay and Montney fairways, we argue that the natural decomposition of kerogen has driven these reservoirs towards failure equilibrium and resulted in a greater likelihood of hydraulic fracturing induced earthquakes.Thus, the natural processes driving overpressure in shale formations are important for understanding the susceptibilities related to hydraulic fracturing induced earthquakes (Pawley et al. 2018;Schultz et al. 2018).
Prior numerical modelling has established the capacity of gas generation to produce formation pressure in excess of disequilibrium compaction, either through late-stage thermal decomposition of kerogen or by in situ cracking of crude oil (Luo & Vasseur 1996).When an organic-rich mudrock is subjected to increasing temperature and pressure during burial, the kerogen within the unit decomposes to progressively lower density liquid and gas phases.This fluid-volume expansion (due to petroleum generation) is accompanied by a formation-pressure increase.The extremely low matrix permeability of petroleum source rocks, typically on the order of 10 −6 -10 −9 Darcy (D) (∼10 −18 -10 −21 m 2 ) under in situ confining stress (Moghadam & Chalaturnyk 2015) inhibits diffusion of generated fluids away from the source bed.This leads to a transient increase in formation pressure that can temporarily retard petroleum generation, although primary migration ultimately overcomes this suppression and allows hydrocarbon generation to proceed (Carr 2000).Crude oil retained in the source rock may crack into lighter C1-C4 gaseous hydrocarbons in the later stages of hydrocarbon generation.Overall, fluid expansion associated with late-stage gas generation is capable of creating high-magnitude overpressure (Hansom & Lee 2005).
However, the maximum magnitude of overpressure can be limited in the presence of local tectonism.Indeed, recently documented feedbacks between hydrocarbon generation and earthquake cycles in the Cordilleran foreland (adjacent to our study regions) demonstrate that fluctuating formation pressure can produce oscillations between stable and unstable fault states (MacKay 2015).This behaviour implies a process of formation pressure buildup punctuated by episodic fluid expulsion (Fig. 4).In this process, hydrocarbon generation serves as both a trigger for fault activation as well as a means of re-establishing fluid pressure during interseismic times.A link between dilatancy, fluid flow and intermittent fault-slip processes has long been recognized in other settings (e.g.Sibson et al. 1975;Sibson 1981;Cox 2016).Specific to the Duvernay Formation, episodic fluid flow behaviour has been observed in fracture injection tests-possibly as a consequence of natural fracture stimulation during hydraulic fracturing (Zanganeh et al. 2018).Such 'fault-valve' behaviour is characterized by sealed-fault pre-seismic pressurization, co/post-seismic fluid discharge and a gradual transition to interseismic fault sealing by remineralization and cementation; after a full cycle this process returns to the pre-/interseismic phase and then repeats (Sibson 1992;Davies & Smith 2006).A salient feature of this model is a co-seismic permeability enhancement along the rupture area of the fault (Zhang et al. 2013;Guglielmi et al. 2015a).In circumstances where a well-developed fault intersects a hydrocarbon source rock, this model implies that the upper limit for pressure buildup in the formation may be limited by fault properties such as strength, architecture and mineralogy rather than fracture-opening pressure (Luo & Vasseur 2002;De Barros et al. 2016).Episodic fault slip could thus be controlled, in part, by fault rheology and residual formation-pressure generation rate.In the case of hydraulic fracturing, such fault-valve behaviour provides a potential conceptual mechanism by which these induced earthquakes may nucleate.A sealed fault which intersects a shale formation would lie seismically dormant during the interseismic period, accumulating formation pressure (without fault slip).A hydraulic fracturing operation rapidly increases the rate of this process via both pressure stimulation, propagation of tensile failures and transmission along natural fractures (Guglielmi et al. 2015a;Lele et al. 2017); in turn, this stimulation creates fluid-pressure pathways and hastens the arrival of the co-seismic phase on a potential nearby fault (Chopra et al. 2017;Corlett et al. 2018).Within this conceptual model for fault reactivation, the increased permeability along the damage zone provides a conduit for fluid discharge into the over/underlying strata (Zhang et al. 2013;Guglielmi et al. 2015a) and thus the propagation of additional induced earthquakes along the fault patch.
During monitoring, we note that adequately resolved earthquake locations have been achieved within these nearby strata (Bao & Eaton 2016;Schultz et al. 2017;Wang et al. 2017).In this sense, regions accumulating high formation overpressure could anticipate susceptibility to induced earthquakes, since they may reflect local regions that are approaching the co-seismic phase.The fault-valve model allows for rapid transmission of fluid overpressure along the damage zone of the fault to seismogenic fault patches, which may be in the nearby strata.For example, faults in the Duvernay have been observed to have a basement rooted connection to a transtensional strike-slip system extending through the Duvernay shale and overlying strata (Corlett et al. 2018;Weir et al. 2018).This anticipation is consistent with the statistical observations in our dataset (Figs 2 and 3).If we extrapolate this reasoning, the anomalous maximum overpressure of the Montney and Duvernay fairways relative to other shale plays (Supporting Information Table S1) could be one factor that is indicative of why other North American shale plays have been relatively quiescent by comparison.
While we are suggesting that formation overpressure is a contributing factor to the occurrence of induced seismicity, we recognize it is not the sole contributor.Geomechanically speaking, increased pore pressure is often utilized as a proxy for fault stability in the unknown subsurface; however, effective stress change on the fault is a more complete picture (Segall & Lu 2015).In fact, fault stability is a function of numerous factors including the friction coefficient, fault orientation and the in situ stress (Walsh & Zoback 2016;Schoenball et al. 2018), of which only one variable is pore pressure.For example, frictional properties of the fault have the potential to allow for aseismic fault movement (and fluid flow) given rate-strengthening or conditionally stable rate-weakening conditions (Scholz 1998;Guglielmi et al. 2015a,b;De Barros et al. 2016).However, in some situations, rate-strengthening considerations may be overcome or play a secondary role (in comparison to effective stress) for incipient failure (Scuderi & Collettini 2016, 2017).Furthermore, the presence of and hydraulic connectivity to seismogenic faults also plays a role in influencing the expression of induced earthquakes (e.g.Hincks et al. 2018;Skoumal et al. 2018).For example, other geological factors have been recognized to influence induced earthquakes within the Duvernay fairway, via a proxy for hydraulically conductive faults (Schultz et al. 2016;Corlett et al. 2018).Lastly, operational factors such as completion volume used during stage stimulation also play a role in the occurrence of hydraulic fracturing related earthquakes (Schultz et al. 2018); this study, however, indicated that geological factors likely played a significant role in the spatial susceptibility of certain regions to induced seismicity.An incorporation of these factors, alongside formation overpressure, contributes to a deeper understanding of where and why these earthquakes occur (Pawley et al. 2018).

C O N C L U S I O N S
Our analysis of two earthquake-prone unconventional resource fairways in the Western Canada Sedimentary Basin implies that high formation overpressure constitutes a potentially, previously unrecognized controlling factor for induced seismicity associated with hydraulic fracturing of organic-rich, fine-grained formations.This association has a practical utility, namely that formation pressure can be independently assessed prior to hydraulic fracturing and could be of value for characterizing site-specific hazard before initiating operations.By analogy with the fault-valve hypothesis, if a seismogenic fault intersects a formation where this process is occurring, episodic fault slip could occur that is controlled by fault architecture, rheology and formation-pressure generation rate.Oscillatory fault activation of this type has been documented during Downloaded from https://academic.oup.com/gji/article-abstract/214/1/751/4995198 by guest on 30 July 2018 the Laramide orogeny in the adjacent Cordillera, and could be analogous to the process by which hydraulic fracturing earthquakes occur.Lastly, these results could be indicative of why some parts of the Montney and Duvernay fairways have been so seismically prone in comparison to other shale plays: their anomalously high maximum overpressures may be one factor contributing to the seismogenic process.

A C K N O W L E D G E M E N T S
This research was supported by the Natural Sciences and Engineering Research Council of Canada (grants IRCSA 485691,.Sponsors of the Microseismic Industry Consortium are thanked, particularly Chevron Canada Resources, Aquaterra Water Management Inc., ConocoPhillips Canada and Nanometrics Ltd. who have provided support at the University of Calgary, and TransAlta for its support through Western University.Canadian Discovery, Ltd. and Amy Fox are thanked for providing access to formation-pressure data used in this work.Burns Cheadle is thanked for fruitful discussions about pore-pressure eration in hydrocarbon source rocks.We would also like to thank two anonymous reviewers for their comments that helped improve this manuscript.Seismicity data were downloaded from the online data archives of Earthquakes Canada.

Figure 1 .
Figure 1.Hydraulically fractured wells (black dots), induced earthquakes (magenta circles) and formation-pressure gradient (coloured area) for (a) the Triassic Montney fairway and (b) the Devonian Duvernay fairway.Earthquakes are from 2009 to 2016 and limited to M L ≥ 2.5.Known or suspected mining blasts and seismicity associated with other causes are omitted.Drilling in the Montney fairway is extensive, but induced seismicity is strongly clustered (BCOGC 2014).Drilling in the Duvernay fairway is concentrated in the Kaybob and Willesden Green areas, with induced seismicity exclusively concentrated within the Kaybob area.

Figure 2 .
Figure 2. Histograms of formation-pressure gradient for the Montney and Duvernay fairways, obtained by interpolation of contours from Fig. 1.The top rows show formation-pressure distribution based on earthquake epicentres and the middle rows show the distribution based on well locations.Value in parenthesis indicates the number of data values and value in upper right in each graph shows mean of the distribution.The lower panels show the distribution of average formation-pressure gradient values based on bootstrap analysis with 10 000 trials; number in the top right shows standard deviation of the bootstrap distribution.All figure panels have the same units on the x-axis.

Figure 3 .
Figure 3. Statistical tests of the overpressure-earthquake association.(a) and (b) Histograms of the distribution of formation pressure at wells (grey bars) are compared against pressures encountered at the declustered earthquake centroids for both Montney (red bars) and Duvernay (green bars) fairways.Results from the Mann-Whitney U tests are reported as p-values in the superimposed text.(c) and (d) Distributions of p-values (blue bars) during sensitivity tests show that statistical significance (dashed line) is still recovered under reasonable accounts of input error.

Figure 4 .
Figure 4. Model for the fault-valve cycle.Formation pressure builds up due to hydrocarbon generation and then is expelled due to co-seismic permeability enhancement.P denotes co-seismic pressure drop, here assumed to be constant for simplicity.Fault seal is gradually restored during interseismic intervals, thus re-initializing the cycle.Modified from MacKay (2015).