## Summary

Controversy and uncertainty about the physics of earthquake triggering mean that stress interactions are rarely incorporated into formal probabilistic forecasts. Statistical methods capture and predict complex features in short-term aftershock forecasts, but we ultimately want to understand the physics behind earthquake triggering. In this paper, we show two fully prospective static stress forecasts that have failed to reproduce spatial patterns of microseismicity. We demonstrate that these failures are not the result of complex main shock ruptures, but are instead caused in part by secondary triggering as deduced from an epidemic type aftershock sequence (ETAS) based stochastic declustering calculation. Prospective testing highlights difficulties in validating physics-based forecasts using microseismicity that can evolve rapidly in time and space. We therefore turn to a global catalogue of larger potentially triggered earthquakes. Prior study with this database found that *M* ≥ 5 earthquakes after main shocks had a ratio of stress-increased to total number of events of 61 per cent, a barely significant result relative to the null value of 50 ± 4.6 per cent. The initial study included every catalogue event; our conclusions from the prospective tests cause us to revisit this choice and to conduct a systematic study of test-event selection. Free parameters include main shock and aftershock magnitude thresholds, as well as calculated probability that aftershocks are background events. We find a mean ratio of stress-increased to total number of events of ~70 per cent across the test parameter range, with high values greater than 80 per cent. The most important parameter is triggered-event magnitude. We therefore conclude that the static stress change hypothesis is significantly more consistent with observation of large earthquake clustering than random chance.

## Introduction

The idea that slip on a fault rearranges stress in the crust and triggers earthquakes (Yamashina 1978; Das & Scholz 1981; Stein & Lisowski 1983; King *et al.* 1994; Stein 1999) provides an appealingly simple physical explanation for the phenomenon. Incorporating static stress changes into earthquake forecasts might account for some of the observed recurrence irregularity, thus improving long-term probability calculations (e.g. Stein *et al.* 1997; Toda *et al.* 1998; Console *et al.* 2008). However, difficulties in identifying statistical correlations between aftershock occurrence and static stress change calculations have been noted (Felzer & Brodsky 2005; Mallman & Zoback 2007; Hainzl *et al.* 2009). Controversies (Felzer & Brodsky 2006; Richards-Dinger *et al.* 2010) about the relative roles of dynamic triggering by seismic waves versus static stress raise further questions about their application to forecasting, as does parameter uncertainty (e.g. Parsons 2005; Field *et al.* 2009).

The purpose of this paper is to evaluate the performance of static stress change calculations made in rapid forecast mode. Forecasting very large aftershocks is of primary importance for hazard assessment, however, the methods are usually tested against the occurrence of small aftershocks. Retrospective evaluations of static stress based forecasts are abundant and many take a variety of calculation uncertainties into account (e.g. Steacy *et al.* 2004; Parsons 2005; Hainzl *et al.* 2009; Cocco *et al.* 2010; Woessner *et al.* 2011).

In this paper, we evaluate static stress change calculations by reporting on the results of two intermediate-term (1-4 yr) prospective forecast efforts. We then attempt to explain some of the observations from prospective testing by considering more complex, stochastic methods of static stress transfer, as well as secondary triggering or epidemic-type aftershocks. Finally, we revisit static stress change calculations made previously using large earthquakes from the global centroid moment tensor catalogue to study the efficacy of including static stress transfer in long-term (decades) earthquake forecasts.

Static stress change calculations involve parameter choices; in forecast mode these choices must be made before they can be fully validated by the occurrence of aftershocks. In this study, we evaluate the success or failure of prospective and blind calculations and therefore hold the parameter choices used to make static stress change calculations as fixed and unchangeable. The point of this paper is to focus on a completely different set of choices, the free parameters that govern how events used to test static stress change calculations are selected. These include magnitude thresholds of main shocks and triggered events, as well as statistical assessments on the likelihood that events are triggered by secondary sources.

## Results from intermediate-term prospective tests

We define a prospective test as follows. After a major earthquake, stress change calculations are made as rapidly as possible that consider a wide array of parameters. These parameters include main shock slip-distribution/geometry, friction coefficient on receiver faults and either receiver fault geometry/rake or an estimate of regional tectonic stress orientations. A preferred solution is developed, which is submitted to a peer review journal. On the day the report is accepted for publication, no further technical changes can be made and the results are fixed. This day marks the beginning of the test period and all subsequent earthquakes that occur can be considered for testing whether the calculated patterns of stress increase and decrease are associated with corresponding seismicity rate increases or decreases.

Here, we present two prospective tests. In the first, the spatial pattern of Coulomb stress was mapped and in the second, stress changes were resolved on individual fault surfaces.

### The 2005 October 8 *M* = 7.6, Kashmir earthquake

Parsons *et al.* (2006) calculated the static stress change expected from this thrust earthquake that struck the Himalayan front causing almost 100 000 deaths (Fig. 1). The Coulomb failure stress is given by

_{n}is stress change normal to optimal fault planes and Δ

*p*is pore pressure change. We applied the concept of an effective coefficient of friction, where Skempton′s coefficient

*B*

_{k}(varies from 0 to 1) is used to incorporate pore fluid effects, in which friction becomes µ′=µ (1 -

*B*

_{k}) and the Coulomb failure criterion becomes after Rice (1992). A range of friction coefficient from 0.0 to 0.8 was assumed to encompass this simplified role of pore fluid pressure.

At the time these calculations were made (using the methods of Toda *et al.* 1998), results from relatively new methods (e.g. Yagi *et al.* 2004) for rapid calculation of large earthquake slip distributions were just starting to be routinely made available. This in turn made rapid stress change calculations possible. In the report, it was demonstrated that the broadest mapped stress change calculations were stable through much of the seismogenic crust and were not very sensitive to parameter choices (Figs 2 and 3). The report was formally accepted on 2006 February 8, marking the beginning of the test period.

At the time of acceptance, most aftershocks [source: Advanced National Seismic System (ANSS) joint catalogue] from the Kashmir earthquake were located close to the main shock plane and were in relatively good agreement with the preferred calculation of spatial Coulomb stress distribution (Fig. 1a) with 94 per cent of the 204 initial aftershocks falling into the stress-increased region of the preferred solution of Parsons *et al.* (2006). However, the forecast did not fare so well during the test period. Aftershocks that occurred through the end of 2009 are plotted in Fig. 1(b) and it seems that the forecast map did a poor job of identifying locations of seismicity rate increases and decreases. A simple count shows that the stress shadow area and the stress-increased region have increased rates during this period, with a ~16-fold increase (11.3 yr^{-1}/0.7 yr^{-1}) and a ~9-fold (4.5 yr^{-1}/0.5 yr^{-1}) jump, respectively. We make formal rate-change calculations that take stress-change-parameter and magnitude-completeness uncertainties into account later.

Mapped stress decreases-often referred to as stress shadows (Harris & Simpson 1996)-are especially helpful in assessing the performance of static stress calculations because they are expected to cause seismicity rate reductions. A rate reduction is diagnostic of a static stress change process because dynamic stressing is expected only to cause rate increases. We quantify the performance of the static stress change calculation by examining regions of mapped stress increase and decrease (Fig. 4). The stress-shadowed zone southwest of the rupture is preferred for study because the shadow zone northeast of the rupture has continued aftershocks from a series of *M* = 5.5 events on the Raikot fault that predate the 2005 *M* = 7.6 Kashmir earthquake (Parsons *et al.* 2006). The stress shadowed and stress increased zones we examine in Fig. 4 are located at similar minimum (30-50 km) distance from the main shock plane and are stable with respect to different boundaries calculated across the parameter ranges. The general observation is an earthquake rate increase in the studied stress shadow area as well as in the stress-increased zone. We note significant seismicity rate increases in both zones when using the *Z* statistic of Habermann (1981) and the β statistic of Matthews & Reasenberg (1988). This is contrary to the static-stress change forecast model in which the seismicity rate changes should correlate with the sign of calculated stress changes.

There are issues associated with making this comparison related to the quality of the earthquake catalogue from the Kashmir region. We calculate a range of completeness thresholds for the catalogue (*M*_{c} = 3.7-4.0) using the maximum curvature method of Wiemer & Wyss (2000), the goodness-of-fit test (Wiemer & Wyss 2000) and the *b*-value stability test of Cao & Gao (2002) as modified by Woessner & Wiemer (2005). The significance of the observed seismicity rate increases has inverse dependence on the magnitude of completeness (*M*_{c}) threshold and decreases with increasing *M*_{c} (Fig. 4).

In the intermediate term (~4 yr), it seems that the post-Kashmir earthquake forecast can be regarded as a failure in terms of its ability to anticipate the spatial distribution of earthquakes in the *M* =~ 3.5 -~5.5 range. Later, we briefly discuss a second prospective test.

### The 2008 May 12, *M* = 7.9 Wenchuan earthquake

The 2008 Wenchuan earthquake struck the Sichuan basin along the eastern edge of the Tibetan Plateau, killing ~70 000 people. Stress change calculations were submitted for peer-reviewed publication 4 d after the event and the final version was accepted on 2008 June 19, marking the beginning of the test period (Parsons *et al.* 2008). A similar calculation was made by Toda *et al.* (2008) using Coulomb stress mapping. Calculations were made on individual faults by Parsons *et al.* (2008) because it seemed that there was detailed knowledge about specific structures in the regions, with many of them underlying highly populated regions. In Fig. 5, we show faults that we calculated were mostly stress-increased in red and stress-decreased faults are shaded blue. Calculations were tested for consistency across a wide parameter range (Figs 6 and 7).

Similarly to the Kashmir case, the early-occurring aftershocks were located nearby the main shock rupture and seemed consistent with the calculated stress change distribution, with 89 per cent of the 244 aftershocks in the catalogue most closely associated with stress-increased faults. Also like the Kashmir case, the aftershocks (source: China Earthquake Administration) that occurred during the test period are inconsistent with our spatial forecast, with many faults shaded in blue showing greater spatial association with aftershocks than faults with calculated stress increases (assuming that earthquakes located adjacent to the mapped faults are indeed occurring on them). We were unable to secure a comprehensive earthquake catalogue for a significant period before the earthquake, so quantitative rate change measurements for *M* = 2 are not possible. Cocco *et al.* (2010) found better spatial matches to the 1992 *M* = 7.4 Landers earthquake aftershocks using optimal-planes mapping as in the Kashmir example of Section 2.1 and Toda *et al.* (2008) applied that method (including spatially variable background stress) to the Wenchuan earthquake. However, many events have occurred in mapped shadow zones. In an attempt to quantify the effects, we calculate before and after rates in mapped stable stress shadow zones (Fig. 5b) using the ANSS joint catalogue for *M* = 4.0 events and find that the rates increased in both areas, by a factor of 1.9 southeast of the rupture and 1.7 to the northwest when rates from 10 yr before and 3 yr after the main shock are normalized and compared.

In summarizing results from the prospective tests, we find initial accordance between calculated stress changes and aftershock occurrence with 94 per cent of the initial aftershocks at Kashmir in the stress increased zone and 89 per cent of Wenchuan aftershocks associated with stress-increased faults. However, over a span of 1-4 yr, we find that forecast spatial distributions of earthquakes are violated, with clear rate increases in stable stress-shadowed zones at Kashmir (Fig. 4) and at Wenchuan (Fig. 5). One might conclude from this that either some delayed, near-source dynamic triggering processes could be important (e.g. Felzer & Brodsky 2006), that the static stress change approach to forecasting is flawed or that the evaluation is flawed. In the next two sections, we study alternative approaches to making and evaluating static stress change calculations.

## Role of main shock slip complexity

Static stress change calculations are usually made by simulating faults as slipping dislocations embedded in an elastic half-space (Okada 1992). Typically a dislocation model of a main shock is necessarily simplified to lie on a single plane. Estimates of main shock slip distribution from teleseismic waves and/or geodetic observations employ damping or smoothing to get a stable solution. Therefore, it is likely that resulting static stress change calculations are also a smoothed, simplified representation of the actual more heterogeneous pattern (e.g. Marsan 2006; Dieterich & Smith 2009).

To investigate the potential role of main shock slip and/or structural heterogeneity, we apply a stochastic method (Herrero & Bernard 1994; Geist 2002), in which the main shock cumulative moment and rake are preserved, but multiscale spatial slip variability is introduced. We examined 100 different realizations of the 2005 October 8 *M* = 7.6 Kashmir earthquake (nine examples shown in Fig. 8) to see if slip or structural heterogeneity during main shock rupture might introduce more complex static stress patterns. For each main shock realization, an accompanying fine-meshed dislocation model is developed and slipped with the same preferred parameter set as put forward in the prospective test of Parsons *et al.* (2006). Stress changes were calculated on optimal planes at 10 km depth and friction coefficient of µ= 0.4.

The resulting stress change mapping has strong variability, but it is spatially limited to the zone adjacent to the main shock plane. This effect might explain aftershock occurrence and spatial variability near the rupture, but regional-scale stress changes are virtually unaffected even by strong introduced slip heterogeneity (Fig. 9). We thus have difficulty explaining the failure of the prospective static stress change forecast in Kashmir because of oversimplification of the main shock source model. In the next section, we apply stochastic declustering (Zhuang *et al.* 2002, 2004; Zhuang 2006) using an epidemic type aftershock sequence (ETAS) model (Ogata 1988, 1992; Ogata & Zhuang 2006) to the Kashmir example to investigate the potential role of secondary triggering.

## Role of secondary triggering determined from stochastic declustering

It has long been recognized that aftershock clusters can have secondary triggering and clustering that extends their reach in time and space (e.g. Ogata 1988, 1992; Guo & Ogata 1997; Marsan *et al.* 2000; Marsan & Lengliné 2008). The static stress change forecast methods we test only include the effects of individual main shocks, thus any secondary triggering or cascading is not accounted for in the initial calculations. Here, we employ stochastic declustering (Zhuang *et al.* 2002, 2004) to estimate the relative probability that each event after the main shock is related to every other event and its probability of being part of the background (spontaneous) seismicity. The idea is that if we limit the aftershock catalogue to those most likely to be related to the main shock, then, we can minimize secondary triggering not accounted for in the initial stress change calculations, with the result being a better test catalogue.

We calculate the probability that events are related to each other (*P _{ij}*) by fitting aftershock characteristic parameters with space-time ETAS (Ogata & Zhuang 2006) and using the algorithm of Zhuang

*et al.*(2002, 2004) and Zhuang (2006) (variables defined in Table 1) as

*f*(

*x*-

*x*,

_{j}*y - y*) defines a short-range Gaussian decay and a long-range inverse power decay (Ogata 1998). We make these calculations for the 832 Kashmir aftershocks above the completeness threshold of

_{j}*M*= 3.8 in the ANSS catalogue through the end of 2010 and developed a matrix of relative probabilities that each event is linked to all subsequent events, as well as probability that an event is part of the background seismicity (

*P*

_{background}; e.g. Wang

*et al.*2010) as

The results are shown in Fig. 10, where aftershocks most associated with the *M* = 7.6 main shock are identified (Fig. 10a), along with an example of secondary triggering by *M* = 6.4 aftershock (Fig. 10b). When the aftershocks are limited to those most probably triggered by the main shock, we find slightly better qualitative agreement with the spatial pattern calculated stress change, though there are still main shock-associated aftershocks located within stress-decreased areas. The spatial mismatch can be quantified by counting the number of events found in the stress-reduced region examined in detail in Section 2.1. There are 10 *M* = 3.8 aftershocks in the stable stress shadow, of which four are linked directly to the main shock (Fig. 11). An additional five events more probably occurred as cascade events and a single aftershock is most associated with the background. This result suggests that, because the stress field likely evolves quickly, a physics-based forecast for all event sizes must also evolve by continuously updating stress changes caused by each aftershock.

We so far have focused primarily on individual main shocks and the small to moderate (*M* = 3.8-6.4) aftershocks that surround them. About half of the aftershocks that occurred in the calculated stress shadow region southwest of the main shock rupture plane can be ascribed to secondary triggering (Fig. 11). However, the other half seem most associated with the main shock and remain unexplained. These events are problematic with respect to the static stress change hypothesis and might be dynamically triggered. These complications suggest to us that the use of small earthquakes to validate stress transfer is not ideal, especially because a key goal for applying static stress change forecasts is to consider long-term interactions amongst the largest, most hazardous earthquakes. Another reason to test static stress calculations on larger earthquakes is that they do not seem to be triggered dynamically (Parsons & Velasco 2011), which reduces the number of physical processes operating. In the next section, we expand our analysis to a global database of main shocks and larger (*M* = 5) aftershocks.

## Assessing long-term static stress change forecasting using the global centroid moment tensor earthquake catalogue

The best way to assess a forecast method is to produce it prospectively (e.g. Woessner *et al.* 2011), as described in Section 2. However, the quantity of available prospective static stress change forecasts is small and their duration is short relative to long-term (30-50 yr) forecast periods. Their validation necessarily involves smaller earthquakes, which we find to have many physical complications with regard to the initial stress change associated with the main shocks. We thus turn to a 20-yr database of more than 100 *M*_{s} = 7 main shocks with known nodal planes that was developed by Parsons (2002) to resolve shear stress changes on subsequent *M*_{s} = 5 earthquake fault planes from the global centroid moment tensor catalogue. A population of 2077 earthquakes selected from within a 2° radius around main shocks that also had an absolute value of shear stress change range of 1 MPa > = 0.1 MPa showed 61 per cent of these were associated with shear stress increase, whereas the rest occurred despite a calculated shear stress decrease. Stress change thresholds were designed to exclude events located too far or too close to the main shock to be reliably associated with static stress calculations. Nodal plane ambiguity prohibited Coulomb stress change calculations (eq. 1) because normal stress changes are not symmetric like shear stresses. These calculations are retrospective, but blind to the calculator because the only parameters involved are the centroid moment tensor planes, magnitude-determined rupture areas and average slip that are drawn directly from empirical relationships (Wells & Coppersmith 1994) and the catalogues. Further, calculations were fully automated because of their large number.

The testing circumstances in the global catalogue are different than those involved with individual main shock forecasts because only larger events are present in the centroid moment tensor catalogue and most importantly, the majority of early aftershocks are undetected (Iwata 2008). Here, then, we are interested in testing the idea that static stress change from one large earthquake can trigger another large event after significant delay (e.g. Stein *et al.* 1997). The Parsons (2002) study indicated that the static stress change model performed slightly better than random chance (61 per cent versus 50 ± 4.6 per cent). We test this ratio against effects of main shock and triggered event magnitudes, as well as probability that potential events are part of the background and are unrelated to main shocks.

We assume that a higher ratio of positively stressed events to the total number implies increasing odds of a static stress change process operating, whereas a reduced ratio suggests decreasing odds. This is not an exercise in tuning parameters to find higher or lower ratios, but rather an examination of the effects of varying them. Tests with subcatalogues drawn from the centroid moment tensor catalogue with randomized times show that the expected ratio is 50 ± 4.6 per cent if there is no static stressing occurring (Parsons 2002). We therefore systematically vary different parameters and observe the effect on the ratio of positively stressed fault planes to total. We point out at the outset that in most cases, varying parameters can ultimately reduce the sample size such that the ratio of stress increased events to total can become unstable. For example, setting a very high minimum magnitude threshold would reduce the number of events significantly. Thus, sample size becomes another parameter to consider.

We study the effects of four parameters: main shock magnitudes, aftershock magnitudes, the probability of aftershocks being part of the background such that they should be excluded from the analysis and sample size used in calculating stress-increase ratios. Magnitudes are an issue because the accuracy of focal plane determination depends on them, with the poorest resolution occurring at the lowest and highest ranges (e.g. Helffrich 1997; Kagan 2003; Bernardi *et al.* 2004). Global centroid moment tensor catalogue solutions are thought to have focal solution axis orientations that are accurate to 15° (Helffrich 1997); stress change calculations are very sensitive to fault orientation and centroid moment tensor uncertainties (Steacy *et al.* 2004; Hainzl *et al.* 2009).

### Influence of high background-event probability

To get background probability, we make calculations using stochastic declustering for 117 individual main shock events of *M* = 7 that were studied previously. Unique ETAS parameters are fit to each main shock and surrounding events/aftershocks for use in declustering. This exercise yields linking probabilities amongst aftershocks, which shows a significant degree of secondary triggering. The global centroid moment tensor catalogue is known to be missing a majority of aftershocks during the first day after main shocks (Iwata 2008). This is a problem with virtually all catalogues (e.g. Kagan 2003), but is particularly pronounced with moment tensor solutions that depend on waveform inversion rather than phase arrival times used for hypocentre detection. We thus did not use the criterion we applied to the Kashmir earthquakes in Section 4, because too many global events went unrecorded immediately after and nearby main shocks to get an accurate probability of their links and because the calculations excluded very near source events (Parsons 2002). We instead focus on the other temporal extreme, when aftershock rates are approaching background rates, which is most appropriate for long-term forecasting.

We calculate the probability that events are related to each other by fitting aftershock characteristic parameters with standard ETAS (Ogata 1988, 1992), where *k*_{0}, *c* and *p* are Omori-law constants (e.g. Utsu 1961). The probability *P _{ij}* that event

*i*is descended from event

*j*is given by the ratio

We note a significant increase in the ratio only if events with 90 per cent or greater odds of being background earthquakes are excluded from the calculations. We see an apparent increase in the ratio if lower background probabilities are excluded (Fig. 12a). If high-probability background events are excluded from the analysis, then the stress-increased ratio is found to be 70-73 per cent for *M* = 5.0 aftershocks in the global centroid moment tensor catalogue.

### Influence of aftershock magnitude

We systematically calculate the ratio of calculated stress increased events to the total with different thresholds of aftershock magnitudes included in the analysis and note an increase with increasing magnitude cut-off (Fig. 12b). Combination of excluding higher background probabilities as well as lower aftershock magnitudes raises the stress-increased ratio to 75 per cent levels, which is a significant increase over the 61 per cent value calculated by Parsons (2002) that included all events. If we consider the extreme case of foreshocks (where aftershock magnitudes exceed the main shock magnitudes), we find 82-100 per cent stress-increased ratios depending on the time interval between foreshock and main shock (Fig. 13). This analysis has a very small sampling of only 11 events. We conclude that magnitude of aftershocks is an important influence on the stress-increased ratio.

### Influence of main shock magnitude

We test one further variable, main shock magnitude. As we have already noted influences of background probability and aftershock magnitude, we now examine simultaneous effects of three variables while tracking sample size. We systematically used different lower main shock-magnitude thresholds while calculating the ratio of calculated stress-increased events to the total. We repeat the calculations for three aftershock-magnitude cut-offs, along with three background probability inclusion levels. The results are shown in Fig. 14; the lower threshold of main shock magnitude has no influence on the stress-increased ratios, except that when higher thresholds (*M* = 7.5) are used, there are too few main shocks included in the analysis and the calculated ratios become unstable. Generally, the only impact of limiting main shock magnitudes to larger values is reduction of sample size.

To summarize, our analysis of the global centroid moment tensor stress-change database finds two influences on its consistency with the static stress change forecast hypothesis. Inclusion of all aftershocks and all main shocks results in a stress change ratio of 61 per cent stress increase to total. We show that this ratio is most strongly affected by the magnitudes of the apparently triggered events and it increases with their increasing magnitude. We further find by stochastic declustering that the stress-change ratio increases when the events with highest calculated probability of being background seismicity are excluded. We finally note that varying main shock magnitudes has no influence.

## Discussion and conclusions

A desire to ultimately use physical models of earthquake interactions, be they dynamic stressing by seismic waves or static stressing by fault slip, motivates our assessment of static stress change forecasting. We show two fully prospective static-stress-change forecasts using available catalogues of microseismicity and find that they fail over intermediate periods of 1-4 yr from the beginning of the test periods. We further note that introducing main shock rupture complexity does little to explain prospective forecasting failures.

We find a different result using an ETAS-based stochastic declustering model, wherein each event in an aftershock series can be classified by its relative probability of being linked to subsequent events in the overall sequence based on their temporal and spatial location and fits to Omori-law decay curves. When microseisms were classified as secondary events (likely triggered by other aftershocks instead of the main shock) and excluded from the prospective test, we find increased consistency with the static stress hypothesis. However, aftershocks directly associated with the main shock are still observed in stress shadow zones, contravening the static stress change concept. We suggest that to properly use microseismicity to validate or invalidate a physics-based forecast likely requires modelling every event that occurs before and during the test period and to have complete understanding of the role that dynamic triggering plays in the near-source region. Such a model is difficult to perform with static stress changes because focal mechanism observations are usually lacking for small earthquakes and has not been shown to improve forecasts (e.g. Rhoades *et al.* 2010), perhaps because of the dynamic effect.

Issues related to using small earthquakes to test forecasts are echoed in the second part of our analysis. We find, when we examine a pre-calculated global database of static stress change on potentially triggered earthquake planes, that there is clearly a much stronger influence on the consistency of earthquake occurrence with the static stress change hypothesis when only *M* > 6.0 rupture planes are considered. Exclusion of *M* < 6.5 events leads to stress-increased ratios >80 per cent (Figs 12 and 15). The ratio of calculated stress-increased events to the total is also influenced by the probability that a potentially triggered earthquake is actually a background event as determined by stochastic declustering, but this is a far milder effect (Figs 12 and 15).

Earthquakes are clustered in space and time and those classified as aftershocks (especially first generation aftershocks) by the ETAS model or by other criteria tend to be spatially-temporally close to preceding earthquakes. Static stress changes are also in these locations. Therefore one might suggest that spatial concurrence with static stress changes is coincidental. We observe *M* > 6 aftershocks preferentially (~70 per cent of cases) occurring where static stress changes are calculated to be positive versus 30 per cent of them happening in stress shadows. We think this ratio is not coincidental for two reasons: (1) stress change calculations were not made in the very near source regions where most aftershocks happen, but in fact the majority were made for events 50-100 km away (Fig. 16). (2) Monte Carlo simulations of the global earthquake catalogue show the maximum stress increased/decreased ratios that happen by chance never exceed 54.6 per cent (Fig. 16).

We conclude from our testing that if static stress change forecasting is applied, it is likely to be most successful for identifying locations and calculating temporal advances of future larger (*M* > 6) earthquakes because the mean values of the stress-increased ratio can be greater than 80 per cent for these rupture sizes and are ~70 per cent across the entire parameter range (Fig. 17). In addition, 96-98 per cent of the calculated results across the parameter range (depending on sample size) lie above the maximum threshold expected by random chance (54.6 per cent). A challenge is of course that, unlike with the centroid moment tensor catalogue of past events, we don′t know where future *M* > 6 earthquake focal planes will be. However, we generally do have structural knowledge about the major faults that can sustain high-magnitude earthquakes. Therefore, if a static stress forecasting technique is to be used, then a fault-based forecast aimed at *M* > 6 earthquakes seems to be the best approach.

## Acknowledgments

Parsons thanks the Institute of Statistical Mathematics for generously hosting a visit during which the bulk of this work was completed. We thank David Jackson, Jeremy Zechar and Daniel Schorlemmer for chairing an AGU special session on ‘Predictive Power of Coulomb Stress’, which helped inspire this manuscript. We are also indebted to Yehuda Ben-Zion, Massimo Cocco, Sandy Steacy and two anonymous reviewers, who all provided extraordinary reviews.

## References

**116,**B05305, doi:10.1029/2010JB007846, in press.