of 7
Reply to
Comment on
Models of Stochastic, Spatially Varying Stress in
the Crust Compatible with Focal-Mechanism Data, and How Stress
Inversions Can Be Biased toward the Stress Rate
by Deborah
Elaine Smith and Thomas H. Heaton
by Jeanne L. Hardebeck
by Deborah Elaine Smith
*
and Thomas H. Heaton
Introduction
In her comment (
Hardebeck, 2015
) on our stress hetero-
geneity article (
Smith and Heaton, 2011
), Hardebeck sug-
gests a different focal-mechanism error distribution than
what we used in our 2011 article and suggests that this new
error distribution will reduce our estimates of stress hetero-
geneity. In response to this, we have rerun our calculations
three ways: (1) with the original mechanism error distribu-
tion from
Smith and Heaton (2011)
, (2) with a mechanism
error distribution similar to the one presented by
Hardebeck
(2015)
, and (3) with a mechanism error distribution derived
from repeating earthquake statistics. We find the two new
mechanism error models, relative to the original mechanism
error distribution, reduce the heterogeneity ratio (
HR
) esti-
mates by approximately 35%
40% (using Hardebeck
s sug-
gested distribution) and by approximately 8%
10% (using
the repeating earthquake based error distribution).
Applying these two new mechanism error distribution
models helps parameterize the estimates of stress hetero-
geneity amplitude but does not change the main novel points
of the
Smith and Heaton (2011)
article. Namely, we find that
focal-mechanism data are still compatible with a hetero-
geneous stress that is more dissimilar at large interevent dis-
tances and more correlated at small interevent distances and
that a heterogeneous stress can bias traditional stress inver-
sions toward the stressing rate function.
Last, we demonstrate that as the size of the stress inver-
sion region decreases and as the maximum variability of the
heterogeneous stress decreases, the normalized stress inver-
sion bias also decreases. This is consistent with taking the
model of
Smith and Heaton (2011)
to the limit where region
size decreases to a point source; however, most stress inver-
sions may require dimensions closer to the outer scale of the
stress (
60
km for southern California) and hence experience
significant stress inversion biasing toward the stressing rate.
Hardebeck (2015)
also refers to her 2010 article (
Har-
debeck, 2010
) as a basis for refuting biasing of stress inver-
sions. We disagree with Hardebeck
s stated refutation of
stress inversion orientations in
Smith and Heaton (2011)
and
Smith and Dieterich (2010)
; however, addressing the
Harde-
beck (2010)
article is best accomplished via a direct model-
ing test in a separate comment that takes into account the
complexities of the nucleation process, including rapidly
evolving slip over the period of days.
Background
Smith and Heaton (2011)
introduced a 3D stochastic
model of stress heterogeneity with power-law-like scaling
for heterogeneity amplitude. When coupled with a stressing
rate and a plastic yield criterion, this model was used to gen-
erate synthetic earthquake focal mechanisms to simulate back-
ground seismicity. To compare the synthetic focal mechanisms
with real data, they added focal-mechanism error with a nor-
mally distributed rotation amplitude centered on zero and uni-
form randomly distributed rotation axes. The statistics of these
synthetic focal mechanisms with modeled measurement error
were compared with the southern California data set by
Har-
debeck and Shearer (2003)
to invert for estimated stress
heterogeneity parameters in southern California. The param-
eter of particular interest is
HR
, as defined in equation (2) of
Smith and Heaton (2011)
,inwhich
HR
measures the ampli-
tude of the stress variation relative to the mean stress level.
The parameter of particular interest is
HR
, which mea-
sures the amplitude of the stress variation relative to the mean
stress level.
The key conclusions in
Smith and Heaton (2011)
are as
follows. First, similar to one of the
Hardebeck (2006)
obser-
vations,
Smith and Heaton (2011)
found the average angular
difference (
AAD
) decreased with decreasing interevent dis-
tance when comparing pairs of synthetic focal mechanisms
generated by their stochastic stress model. This implies the
power-law filtering generates stress that is fairly uncorrelated
at large distances and is increasingly correlated or similar at
smaller distances.
Second, they found that when one applies a failure cri-
terion to generate events from a stress grid with spatially var-
iable initial stress and a spatially uniform loading rate, the
events tend to cluster toward the loading rate direction. In
other words, by treating earthquakes as a critical phenome-
non, the points that preferentially fail already have a prestress
*Now at U.S. Geological Survey, 345 Middlefield Road, Menlo Park, Cal-
ifornia 94025.
452
Bulletin of the Seismological Society of America, Vol. 105, No. 1, pp. 452
458, February 2015, doi: 10.1785/0120140132
orientation somewhat aligned with the loading stress. This
produced a biased sampling, in contrast to a simple uniform
random sampling of the stress in the grid.
This second finding has two implications. The events
tend to cluster more tightly than one would calculate for a
uniform random sampling of the stress grid. This can lead
to underestimates of the stress heterogeneity parameter,
HR
,
if the tight clustering is not properly taken into account. If a
spatially uniform loading stress
and the spatial mean of prestress
are misaligned, there can be a bias in the focal-mechanism
event orientations that are more consistent with the loading
stress orientation.
Error Distributions for Modeling
Focal-Mechanism Uncertainty
As pointed out by
Hardebeck (2015)
, taking into ac-
count the geometry of the problem generates an estimated
error distribution for modeled focal-mechanism uncertainty
that has a peak offset from zero, in contrast with our normal
distribution with a peak centered at zero. The only questions
are: What does that distribution look like and to what degree
is it offset from zero? Then to what degree does this new dis-
tribution affect our estimates of stress heterogeneity,
HR
?
To answer these questions, we use three models of focal
mechanism error to generate noisy synthetic focal mecha-
nisms. These three error models include: (1) our original nor-
mal error distribution, (2) a distribution similar to what
Hardebeck (2015)
presents, and (3) a distribution based
on Shearer and Hardebeck
s Northridge repeating earthquake
aftershock data (
Shearer
et al.
, 2003
). Then by comparing
noisy synthetic mechanisms with mechanisms from Shearer
and Hardebeck (see
Data and Resources
) we recalculate our
estimates of stress heterogeneity (
HR
) and percent bias in
light of these three mechanism error models.
The first new error distribution we explore, based on the
Hardebeck (2015)
suggestion, may be representative of the
uncertainty for individual focal-mechanism measurements;
therefore, we choose it as the upper bound for our uncertainty
estimates. However, the calculations that estimated
HR
in
figure 9 of
Smith and Heaton (2011)
, based on data from
Hardebeck (2006)
, determined the best fit to
AAD
s between
pairs of focal mechanisms. So the question arises: To what
extent could the relative focal-mechanism error between
pairs of mechanisms be significantly less than this upper
bound, especially for closely spaced mechanisms?
To address this question, we choose the second new er-
ror distribution, based on statistics of repeating earthquake
clusters, as a means of estimating relative focal-mechanism
error. The idea is fairly simple. Because repeating earth-
quakes have similar waveforms and occur essentially at the
same location, any variation in the mechanisms for a particu-
lar cluster can be approximated to first order as relative meas-
urement uncertainty. To test this, Shearer and Hardebeck
were willing to share with us their Northridge aftershock re-
peating cluster catalog from
Shearer
et al.
(2003)
. From this
catalog, we extract A and B quality events (the same quality
level used in our 2011 study and in Hardebeck
s study),
choose clusters with a minimum of three events, and calcu-
late the mean angular spread for each cluster. From a total of
294 events, we have the following distribution of mechanism
error as shown in Figure
1a
. It has a mean angular spread of
approximately 19°, which is 5° more than the normal distri-
bution studied by Smith and Heaton, and an
AAD
of approx-
imately 27°, which is very similar to the minimum
AAD
seen
by
Hardebeck (2006)
at distances less than 100 m.
In modeling the focal-mechanism error, we represent
each rotation error as a rotation amplitude
ω
about a rotation
axis

θ
;
φ

. For the normal distribution error in
Smith (2006)
0
10
20
30
40
50
60
70
80
0
0.02
0.04
Mean Angular Spread from Northridge Repeating Earthquakes,
Synthetic Data (Red Line) Mean = 18
°
Mean Angular Spread About the Average Orientation (°)
0
10
20
30
40
50
60
70
80
0
0.01
0.02
0.03
0.04
Average Angular Difference (°)
Average Angular Difference from Northridge Repeating Earthquakes,
Synthetic Data (Red Line) Mean = 24
°
0
10
20
30
40
50
60
70
80
0
0.02
0.04
Mean Angular Spread About the Average Orientation (°)
Mean Angular Spread for Hardebeck’s Error Estimate
Synthetic Data (Red Line) Mean = 31
°
0
10
20
30
40
50
60
70
80
0
0.01
0.02
0.03
0.04
Average Angular Difference (°)
Average Angular Distance for Hardebeck’s Error Estimate
Synthetic Data (Red Line) Mean = 31
°
(a)
(b)
Figure
1.
(a) Mean angular spread mechanism distribution and
average angular difference (
AAD
) from repeating earthquake clus-
ters in the Northridge aftershock sequence (
Shearer
et al.
, 2003
).
Only A and B quality events were selected, with a minimum cluster
size of three events each. Modeled mechanism error statistics are
plotted with dashed lines using the parameterization explained in
the text. (b) Similar histograms are plotted from the
Hardebeck
(2015)
data. Modeled mechanism error statistics are again plotted
with dashed lines using a parameterization explained in the text. The
color version of this figure is available only in the electronic edition.
Reply
453
and
Smith and Heaton (2011)
, we had a normally distributed
ω
and uniform random

θ
;
φ

. For the repeating earthquake
distribution statistics in Figure
1a
and for the
Hardebeck
(2015)
estimate of focal-mechanism error in Figure
1b
,we
apply lognormal statistics for the rotation amplitude
ω
and
Von-Mises distributions for the rotation axes

θ
;
φ

. This has
the effect of offsetting the rotation amplitude from zero, us-
ing the lognormal statistics and clustering the error rotations
due to the circular statistics on the rotation axes. In Figure
1a
,
when simulating repeating earthquake statistics, we applied a
μ

2
:
7
and
σ

0
:
67
for the
ω
lognormal statistics and a
θ

0
and
κ

1
:
0
for each angle in the rotation axis

θ
;
φ

.
In Figure
2a
, when simulating the
Hardebeck (2015)
error
statistics, we applied a
μ

3
:
375
and
σ

0
:
41
for the
ω
lognormal statistics and a
θ

0
and
κ

2
:
25
for each angle
in the rotation axis

θ
;
φ

. Note that the lognormal random
numbers were drawn via the standard MATLAB toolbox, and
the circular random numbers were generated by Beren
s
circ_vmrnd
MATLAB toolbox (2006; see
Data and Resources
).
How the Error Distributions Affect Stress
Heterogeneity Amplitude,
HR
,
α
, and Percent
Biasing toward the Stressing Rate Orientation
In Figure
2
, we examine the synthetic
AAD
as a function
of focal-mechanism pair separation distance to invert for
HR
and
α
as in figure 9 of
Smith and Heaton (2011)
. This plot is
based on
AAD
statistics from
Hardebeck (2006)
for A and B
quality southern California focal-mechanism data. We do this
for each error distribution by showing the best-fitting synthetic
solutions, using a grid-search technique similar to but simpler
than
Smith and Heaton (2011)
. This grid search focuses on the
parameter space close to the solution and ignores location
error for 2D grids. The best-fitting synthetic with normal dis-
tribution error is shown in Figure
2a
,with
α

0
:
8
and
HR

2
:
25
; the best-fitting synthetic with repeating earth-
quake error is shown in Figure
2b
,with
α

0
:
9
and
HR

2
:
0
; and the best-fitting synthetic with Hardebeck
s
error is shown in Figure
2c
,with
α

1
:
0
and
HR

1
:
5
.
These plots show three things. (1) The new error distri-
butions generate smaller
HR
s. The best-fit synthetic with
repeating earthquake error predicts an
HR
11% smaller than
the estimate with normal error and 16% less than
HR

2
:
375
, which was originally estimated by
Smith and
Heaton (2011)
. The best-fit synthetic with
Hardebeck (2015)
error predicts an
HR
33% smaller than the estimate with
normal error and 37% smaller than the
Smith and Heaton
(2011)
estimate. (2) The more error that is added to the focal-
mechanism synthetics, the more spatial smoothing is needed
to generate similar
AAD
curves, which results in higher
α
.
(3) It is impossible to generate the minimum of 26°
AAD
with
Hardebeck
s error as seen in Figure
2a
, in which the minimum
AAD
for mechanism error alone is 31°
32°. Consequently,
both distributions predict a small-to-moderate reduction in
the estimate of
HR
, and the statistics of the repeating events
10
−1
10
0
10
1
10
2
0
10
20
30
40
50
60
70
Distance Between Events (km)
Average Angular Difference, AAD (°)
AAD vs. Distance for Normal Error, HR = 2.25,
α
= 0.8
10
−1
10
0
10
1
10
2
0
10
20
30
40
50
60
70
AAD vs. Distance for Repeating Quakes, HR = 2.0,
α
= 0.9
Distance Between Events (km)
Average Angular Difference, AAD (°)
10
−1
10
0
10
1
10
2
0
10
20
30
40
50
60
70
AAD vs. Distance for Hardebeck Error, HR = 1.5,
α
= 1.0
Distance Between Events (km)
Average Angular Difference, AAD (°)
(a)
(b)
(c)
Figure
2.
Approximations of the
Hardebeck (2006)
calculation
of
AAD
between pairs of focal mechanisms as a function of inter-
event distance for southern California A and B quality mechanisms
(black lines) and the best-fitting synthetic focal-mechanism data
from 2D grids given a particular mechanism uncertainty model
(dashed gray lines). (a) Using a normal distribution centered on zero
for mechanism uncertainty like
Smith and Heaton (2011)
, we esti-
mate the best-fitting parameters to be approximately
α

0
:
8
and
HR

2
:
25
. (b) Using mechanism uncertainty similar to repeating
earthquakes, we estimate the best-fitting parameters to be approx-
imately
α

0
:
9
and
HR

2
:
0
. (c) Using mechanism uncertainty
similar to
Hardebeck (2015)
, we estimate the best-fitting parameters
to be approximately
α

1
:
0
and
HR

1
:
5
; however, because the
minimum
AAD
for this distribution is approximately 31°, it is
impossible to precisely fit the data at small interevent distances.
(a)
(c) Figures were modified from
Smith and Heaton (2011)
,
which was originally modified from
Hardebeck (2006)
.
454
Reply
may yield the most appropriate relative mechanism error dis-
tribution, especially for short interevent distances.
Smith and Heaton (2011)
estimated the stress hetero-
geneity amplitude,
HR
, and percent biasing for 12 subregions
in southern California, so we redo those estimates with our new
focal-mechanism error distributions. One method that
Smith
and Heaton (2011)
appliedtoestimatethe
HR
value was the
AAD
between pairs of focal mechanisms. Generally, the noisier
the data the larger the
AAD
, which encompasses both mecha-
nism measurement error and underlying stress heterogeneity
hence, the need for good mechanism error estimates when us-
ing this statistic to estimate stress heterogeneity. In other words,
the wrong model of error can cause either an under- or over-
estimate of the stress heterogeneity.
Smith and Heaton (2011)
developed a scaling plot that
relates the
AAD
to the
HR
. To create this, suites of synthetic
focal-mechanism catalogs with different values of
HR
and an
assumed mechanism error distribution were employed. The
scaling plot then allowed them to overlay the data for each
region as horizontal lines in their figure 12a. This mapped
out a range of
HR
, the values of which are compatible with
southern California.
In this reply to Hardebeck
s comment, we calculate
HR
scaling plots for three distributions. In Figure
3a
,weshowthe
scaling plot for our original normal distribution with zero
mean; in Figure
3b
, we show the scaling plot for the distribu-
tion suggested from the statistics of repeating earthquakes; and
in Figure
3c
, we show the scaling plot for a distribution similar
to the one suggested by
Hardebeck (2015)
. The intersection of
the data (shown as horizontal lines) maps out the range of
stress heterogeneity values observed for each assumed mecha-
nism error distribution. The results are shown in Table
1
.
We find that by updating the mechanism error model,
the estimates of
HR
are reduced but by varying amounts.
The mechanism error distribution suggested by repeating
earthquake statistics generates
HR
values in the 0.96
2.62
range, approximately 8%
10% less than the original calcu-
lation. The mechanism error distribution suggested by
Har-
debeck (2015)
generates
HR
values in the 0.44
2.20 range,
approximately 35%
45% less than the original calculation.
When we compare these new
HR
estimates with the per-
cent bias scaling plotted in figure 17 of
Smith and Heaton
(2011)
, we have the following estimates of percent bias in
stress inversions. The distribution for mechanism error based
on repeating earthquakes generates our upper bound for
HR
and biasing, where biasing is estimated in the 30%
60%
range, with an average of about 51%. The distribution similar
to the one suggested by
Hardebeck (2015)
generates our lower
bound for
HR
and biasing, where biasing is estimated in the
15%
50% range, with an average of about 38%.
The Scale Dependence of the Inversion Region for
Stress Orientation Biasing
The other effect we investigate is how the power-law
filtering of the stress affects stress inversion biasing. We note
10
−1
10
0
10
1
10
2
0
10
20
30
40
50
60
70
80
90
Heterogeneity Ratio (HR)
Average Angular Distance,
AAD
(°)
HR
Scaling for Normal Error Using
AAD
10
−1
10
0
10
1
10
2
0
10
20
30
40
50
60
70
80
90
Heterogeneity Ratio (HR)
Average Angular Distance,
AAD
(°)
HR
Scaling for Repeating Earthquake Error Using
AAD
10
−1
10
0
10
1
10
2
0
10
20
30
40
50
60
70
80
90
Hetero
g
eneit
y
Ratio
(
HR
)
Average Angular Distance,
AAD
(°)
HR
Scaling for Hardebeck’s Error Using
AAD
(a)
(b)
(c)
Figure
3.
How different mechanism error distributions affect the
estimates of heterogeneity ratio (
HR
). The solid black line shows the
scaling relationship derived from calculating the
AAD
of synthetic fo-
cal mechanisms with varying
HR
and an assumed mechanism error
distribution. The horizontal lines show the
AAD
of 12 regions in
southern California, as listed in Table
1
. The intersection of the hori-
zontal lines with the solid black scaling relationship helps define the
span of
HR
estimates suggested by the data and assumed mechanism
error. In (a), the scaling relationship for the normal distribution with
zero mean suggested by
Smith and Heaton (2011)
is recalculated
and plotted with
HR
in the 1.08
3.14 range. In (b), a similar scaling
relationship is calculated for the mechanism error distribution sug-
gested by the repeating earthquake data with
HR
in the 0.92
2.62
range, an approximately 8%
10% decrease. In (c), a scaling relation-
ship with a mechanism error distribution similar to
Hardebeck (2015)
yields an
HR
in the 0.44
2.20 range, an approximately 35%
45% de-
crease. The Apple Valley region experiences a slightly higher decrease
of 60%. The color version of this figure is available only in the elec-
tronic edition.
Reply
455
that both the data analyzed by
Hardebeck (2006)
and the
power-law filtered stress heterogeneity generated by
Smith
and Heaton (2011)
suggest that stress is more correlated at
small interevent distances up to an outer scale. The grid-
search inversion of the synthetic focal-mechanism data sug-
gests a spatial smoothing parameter of
α
in the 0.8
1.0 range.
If
α

1
:
5
, there would be no biasing at all because the
stress would be almost homogeneous. On the other hand,
if
α

0
:
0
, then every point is uncorrelated with the other
and biasing would occur uniformly everywhere regardless
of the length scale. However, for
α
ranging from 0.8 to
1.0, there is heterogeneity but at a small enough length scale
that it looks relatively homogeneous. In this case, some re-
gions actively fail and others do not, depending upon their
initial stress amplitude and orientation. We therefore hypoth-
esize that shrinking the dimension of the stress inversion re-
gion will result in less stress inversion biasing because there is
less stress variation within the smaller regions. In other words,
as the size of the stress inversion region shrinks to zero, the
variation in the stress becomes zero, and the stress inversion
biasing disappears within individual regions. Interestingly, as
we go toward this limit, we now have a biased sampling of
which regions fail where the regions are points, which is ex-
actly the biased sampling of individual events as in the original
model of
Smith and Heaton (2011)
.
To test this, we select 2D synthetic focal-mechanism cat-
alogs from
Smith and Heaton (2011)
and calculate normal-
ized stress inversion biasing for different sizes of subregions.
We use the best grid-search result for the repeating earthquake
error distribution (
α

0
:
9
;
HR

2
:
0
). Using catalogs with
10 different random seeds, we randomly subdivide the simu-
lation space for each catalog into smaller subregions. Then we
randomly sample the events within each subregion and add
random mechanism error consistent with repeating earthquake
statistics. This random subdividing and resampling is repeated
at least 50 times. Last, we calculate a normalized stress inver-
sion biasing similar to the percent biasing statistic presented in
Smith and Heaton (2011)
, where it represents the percent bias-
ing toward the stressing rate tensor. In this case, it is normal-
ized because, for a small number of events in the subregions,
sometimes the percent biasing is larger due to random mecha-
nism error. Therefore, we use the summation of two statistics
from
Smith and Heaton (2011)
,themean
ψ
BI
=
ψ
BT
and the
mean
ψ
IT
=
ψ
BT
(the sum of which should equal one in the ab-
sence of noise), to normalize our estimates of biasing. Then
we plot the results in Figure
4
.
As defined in equation (23) of
Smith and Heaton (2011)
,
note that
ψ
BI
is the four component angular difference be-
tween two deviatoric stress tensors, the spatial mean back-
ground stress tensor and the stress tensor from inversion
of focal-mechanism data. The stress inversion is calculated
using Michael
s program
slick
(see
Data and References
)
based on
Michael (1984
,
1987)
.
ψ
BT
is the angular difference
between the background stress tensor and the tectonic load-
ing stress tensor. Last,
ψ
IT
, is the angular difference between
the tectonic stress tensor and the inverted stress tensor.
We find the following: (1) Normalized stress inversion
biasing decreases as the size of the subregion decreases.
(2) The percent of subregions with sufficient events to fail
dramatically decreases; therefore, at the smaller dimensions,
only a small fraction of the subregions fail. (3) The fall-off of
normalized stress inversion biasing with decreasing spatial
scale follows the fall-off of
AAD
with decreasing interevent
distance simulated by
Smith and Heaton (2011)
and calcu-
lated by
Hardebeck (2006)
. At the same time, due to limita-
tions in our simulation, it is difficult to take the region size
below 2 km. Regional stress inversion studies may face sim-
ilar limitations; therefore, our estimates of stress inversion
biasing in the previous section should hold for most regional
stress inversions.
Conclusions
We appreciate the work
Hardebeck (2015)
has done in
demonstrating that spherical geometry requires a mechanism
error distribution with a nonzero mean. Her suggestion of a
more realistic mechanism error distribution has been a
Table 1
The Effect of Mechanism Error Distributions on the Estimates of Heterogeneity Ratio (
HR
)
Average Angular
Difference
HR
for Normal
Distribution with
Mean of Zero
HR
for Distribution
from Repeating
Earthquakes
Percent Change
Decrease from Normal
Distribution
HR
for Distribution
Similar to
Hardebeck (2015)
Percent Change
Decrease from
Normal Distribution
Ventura basin
57.7
1.90
1.74
8.1
1.26
33.7
San Gabriel Mountains
61.2
2.42
2.23
8.0
1.69
30.2
Los Angeles basin
62.1
2.62
2.40
8.3
1.84
29.8
Apple Valley
46.2
1.08
0.96
11.0
0.44
59.8
Landers
63.8
3.14
2.86
9.0
2.20
29.9
Banning
58.7
2.02
1.86
7.9
1.36
32.4
Palm Springs
62.9
2.83
2.59
8.4
1.99
29.6
Coachella
56.9
1.81
1.66
8.1
1.18
34.6
Northern Elsinore
53.3
1.49
1.36
8.4
0.89
40.0
Central Elsinore
53.0
1.47
1.34
8.6
0.87
40.6
Anza
54.0
1.55
1.42
8.7
0.94
39.0
Borrego
54.4
1.58
1.45
8.5
0.98
38.3
456
Reply
helpful contribution. We thank her and Peter Shearer for
sharing the Northridge repeating earthquake catalog (
Shearer
et al.
, 2003
).
We suggest that the focal-mechanism error derived from
repeating earthquake data (our upper bound for
HR
) is the
best error for angular differences between pairs of events,
especially for small interevent distances. Indeed, the relative
error between pairs of mechanisms close to one another may
be significantly less than the absolute error for individual
mechanisms. This relative uncertainty most likely increases as
the interevent distance increases; hence, we consider the abso-
lute uncertainty estimated by
Hardebeck (2015)
for individual
mechanisms to be the upper limit on mechanism error and
more applicable for the outer-scale limit of the simulations.
When applying a grid search in Figure
2
with our pre-
ferred repeating earthquake error distribution,
HR
is reduced
by only 16% (from 2.375 to 2.0) from the
Smith and Heaton
(2011)
value, and
α
is slightly increased from 0.8 to 0.9. The
error distribution suggested by
Hardebeck (2015)
suggests
an
HR
of approximately 1.5 and an
α

1
:
0
, but it has trou-
ble fitting the small interevent data well. When using these
new error statistics to estimate
HR
for the 12 subregions, we
find that
HR
is decreased anywhere from 8%
11% (for our
preferred statistics from repeating earthquake clusters) to
35%
45% (for the mechanism error distribution suggested
by
Hardebeck, 2015
). We consider these estimates to provide
upper and lower bounds on the stress heterogeneity. In our
updated estimates, the biasing toward the stressing rate as in
Smith and Heaton (2011)
can be anywhere from 15% to 50%
for the lower
HR
bound and anywhere from 30% to 60% for
the upper
HR
bound, which we prefer.
Although this has helped refine the estimated values, we
also find the following major points of
Smith and Heaton
(2011)
to still be valid: (1) focal-mechanism statistics are
consistent with stress heterogeneity that is more highly cor-
related at small length scales and less correlated at longer
length scales, and (2) stress inversions applied to regions of
10 km or greater width can exhibit significant biasing toward
the stressing rate.
The stress heterogeneity model used by
Smith and Hea-
ton (2011)
and in this reply is a very simple heterogeneity
model. It would be surprising if the real Earth can be repre-
sented by so few parameters. In fact, there could easily be
time dependence to the heterogeneity amplitude and a spatial
dependence with distance from major fault traces. Ulti-
mately, the best way to measure heterogeneous stress is to
look at borehole breakout data, which is the actual measure-
ment of stress orientation variability over the length scale of
meters. Recent data, such as figure 2 from
Day-Lewis (2010)
,
show maximum horizontal stress rotations of 45° over length
scales of order 20 m in the San Andreas Fault Observatory at
Depth (SAFOD) Pilot Hole, indicative of larger stress variabil-
ity over short length scales close to the fault. This may show
even greater variability of stress than what our model contains;
therefore, how much stress variability there is in the real
Earth remains an outstanding question to be further refined
and studied.
Data and Resources
The focal-mechanism data used were A and B quality
events from the 1984 to 2003 southern California data set
(
Hardebeck and Shearer, 2003
). This catalog is documented
at the website
www.data.scec.org/research/altcatalogs.html
(last accessed May 2014). Repeating earthquake data for the-
Northridge aftershocks were obtained from Shearer and Har-
debeck (
Shearer
et al.
, 2003
). Some fitting of the hyperbolic
functions was aided by Kaleidagraph. The
circm_vrnd.m
Circular Statistics Toolbox by Philip Berens was obtained
from
http://www.mathworks.com/matlabcentral/fileexchange/
10676-circular-statistics-toolbox
directional-statistics-/content/
circ_vmrnd.m
(last accessed November 2014). All other code
was developed using MATLAB (
www.mathworks.com/
products/matlab
; last accessed November 2014) and run on
MacPro computers with OS X. The program
slick
by Michael
was obtained from
http://earthquake.usgs.gov/research/
software/
(last accessed December 2014).
10
−1
10
0
10
1
10
2
Distance Between Events/Inversion Region (km)
Average Focal Mechanism Angular Difference
26-27° minimum average angular difference
for repeating earthquakes mechanism error
inner-scale of = 60 m
outer-scale of = 60 km
10°
20°
30°
40°
50°
60°
70°
0
0.5
0.4
0.3
Normalized % Biasing Toward Tectonic Orientation
0.2
0.1
0.0
% Biasing
Comparison of AAD and % Biasing with Distance
Figure
4.
Approximation of the
Hardebeck (2006)
calculation
of
AAD
between pairs of focal mechanisms as a function of inter-
event distance for southern California A and B quality mechanisms
(thin line), and the calculated normalized stress inversion biasing for
α

0
:
9
and
HR

2
:
0
and focal-mechanism error consistent with
repeating earthquake observations (thick gray line). The modeled
normalized stress inversion biasing toward the tectonic stressing
rate depends on the outer size of the inversion region. We assume
that no biasing occurs when there is little to no heterogeneity at the
inner scale and that the maximum biasing occurs at the outer scale.
If the outer scale of the stress is 60 km, as suggested by
Smith and
Heaton (2011)
, then the stress inversion biasing is about 50% for
subregions
60
km and decreases as the inversion subregions de-
crease. At a subregion size of 3 km, the normalized stress inversion
biasing is approximately 35%. However, at the same time, as the
size of the subregions decreases, the percentage of the subregions
that have sufficient number of events to fail also decreases. Conse-
quently, the seismicity becomes a biased sampling of which subre-
gions fail where these failing subregions represent local stress
variations. The figure was modified from
Smith and Heaton (2011)
,
which was originally modified from
Hardebeck (2006)
.
Reply
457
Acknowledgments
We thank Jeanne Hardebeck and Peter Shearer for sharing the
Northridge aftershock repeating earthquake cluster catalog with us. We also
thank Jeanne Hardebeck and the Southern California Earthquake Center
(SCEC) for the focal-mechanism catalog of southern California events.
References
Day-Lewis, A., M. Zoback, and S. Hickman (2010). Scale-invariant stress
orientations and seismicity rates near the San Andreas fault,
Geophys.
Res. Lett.
37,
L24304, doi:
10.1029/2010GL045025
.
Hardebeck, J. L. (2006). Homogeneity of small-scale earthquake faulting,
stress and fault strength,
Bull. Seismol. Soc. Am.
96,
1675
1688.
Hardebeck, J. L. (2010). Aftershocks are well aligned with the background
stress field, contradicting the hypothesis of highly heterogeneous
crustal stress,
J. Geophys. Res.
115,
doi:
10.1029/2010JB007586
.
Hardebeck, J. L. (2015). Comment on
Models of stochastic, spatially vary-
ing stress in the crust compatible with focal-mechanism data, and how
stress inversions can be biased toward the stress rate
by Deborah
Elaine Smith and Thomas H. Heaton,
Bull. Seismol. Soc. Am.
105,
no. 1, doi:
10.1785/0120130127
.
Hardebeck, J. L., and P. M. Shearer (2003). Using
S=P
amplitude ratios to
constrain the focal mechanisms of small earthquakes,
Bull. Seismol.
Soc. Am.
93,
2434
2444.
Michael, A. J. (1984). Determination of stress from slip data: Faults and
folds,
J. Geophys. Res.
89,
11,517
11,526.
Michael, A. J. (1987). Use of focal mechanisms to determine stress: A con-
trol study,
J. Geophys. Res.
92,
357
368.
Shearer, P. M., J. L. Hardebeck, L. Astiz, and K. B. Richards-Dinger
(2003). Analysis of similar event clusters in the aftershocks of the
1994 Northridge, California, earthquake,
J. Geophys. Res.
108,
doi:
10.1029/2001JB000685
.
Smith, D. E. (2006). A new paradigm for interpreting stress inversions from
focal mechanisms: How 3D stress heterogeneity biases the inversions
toward the stress rate, California Institute of Technology, Pasadena,
California, available at
http://thesis.library.caltech.edu/2060/
(last ac-
cessed December 2014).
Smith, D. E., and J. H. Dieterich (2010). Aftershock sequences modeled with
3D stress heterogeneity and rate-state seismicity equations: Implications
for crustal stress estimation,
Pure Appl. Geophys.
167,
1067
1085.
Smith, D. E., and T. H. Heaton (2011). Models of stochastic, spatially vary-
ing stress in the crust compatible with focal-mechanism data, and how
stress inversions can be biased toward the stress rate,
Bull. Seismol.
Soc. Am.
101,
1396
1421.
Department of Terrestrial Magnetism
Carnegie Institution of Washington
5241 Broad Branch Road
Washington, D.C. 20015
desmith144@gmail.com
(D.E.S.)
California Institute of Technology
252-21 Seismolab
Pasadena, California 91125
heaton@caltech.edu
(T.H.H.)
Manuscript received 14 May 2014;
Published Online 13 January 2015
458
Reply