Close
About
FAQ
Home
Collections
Login
USC Login
Register
0
Selected
Invert selection
Deselect all
Deselect all
Click here to refresh results
Click here to refresh results
USC
/
Digital Library
/
University of Southern California Dissertations and Theses
/
Advancing probabilistic seismic hazard analysis through short-term time-dependent forecasts and deterministic simulations
(USC Thesis Other)
Advancing probabilistic seismic hazard analysis through short-term time-dependent forecasts and deterministic simulations
PDF
Download
Share
Open document
Flip pages
Contact Us
Contact Us
Copy asset link
Request this asset
Transcript (if available)
Content
Advancing Probabilistic Seismic Hazard Analysis Through Short-
Term Time-Dependent Forecasts and Deterministic Simulations
By
Kevin Ross Milner
A Dissertation Presented to the
FACULTY OF THE USC GRADUATE SCHOOL
UNIVERSITY OF SOUTHERN CALIFORNIA
In Partial Fulfillment of the
Requirements for the Degree
DOCTOR OF PHILOSOPHY
GEOLOGICAL SCIENCES
May 2020
Copyright 2020 Kevin Ross Milner
ii
Acknowledgements
I have been extremely fortunate to have the patient and unwavering support of
many loved ones, colleagues, and mentors, all of whom have shaped my life and career
tremendously. My family has always encouraged me to believe in myself and supported
me at every step along the way. My mom, Christine, is my biggest cheerleader and
advocate, encouraging me from day one (and even editing much of my undergraduate
writing). My dad, Brian, is the most patient man that I know and a gifted engineer;
without his help I may not have passed a number of undergraduate math classes. My
sister, Amy, has always been my best friend and the person who I aspire to be. My
partner, Jordan, has loved and encouraged me through the most stressful times of
trying to finish a dissertation, and improved it as a gifted editor; every piece of my
writing that she touched has benefited tremendously from it. Jordan, my “d-word” is
almost behind us, and I stand ready to support you in your next ventures.
Professionally, a number of mentors and colleagues repeatedly encouraged me to
expand beyond my initial role as a computer programmer at the Southern California
Earthquake Center (SCEC) and to produce original research. Tom Jordan hired me and
called me a “scientist” long before I believed the label applied. Every time that I visited
his office and he started scribbling on the white board I felt as though I was part of
something important, even when I didn’t fully understand the details. He encouraged me
to advance my career by becoming his student, and bait-and-switched me into pursuing
a doctoral degree under the guise of a masters. Ned Field has always been my closest
collaborator and gave me way more respect and responsibility than I deserved as a
iii
fresh-out-of-college programmer. As I learned the ropes and gained more domain
expertise, he continued to trust me and delegate increasingly complex research tasks,
and has always provided me with encouragement and friendship. Bruce Shaw has been
an incredibly generous collaborator, sharing his excitement, vision, and expertise with
me as we worked together to build the model described in Chapter 4. My supervisor at
SCEC, Christine Goulet, has always supported my doctoral research while encouraging
me to finish, and continues to give me the space to pursue a diverse and interesting
range of research projects. Scott Callaghan has provided me with years of great
friendship, and also taught me to run CyberShake (and spent many hours
troubleshooting); I could not have done this without his help. I would like to thank the
rest of my co-authors for their valuable contributions, especially Keith Richards-Dinger
for his clear explanations, Morgan Page for encouraging me to write up the UCERF3-
ETAS Ridgecrest results and helping along the way, and Bill Savran for being a great
sounding board to test ideas off of. Phil Maechling has always supported my
professional development and works tirelessly and selflessly in support of SCEC
computational science. The SCEC professional staff, including Tran Huynh, Deborah
Gormley, and the recently-retired John McRaney, have always been on my team. I
would like to thank my graduate advisor, Cindy Waite, who has helped me navigate the
complexities of degree progress and paperwork with a consistent friendly smile, and my
committee, Professors Sylvain Barbot, Carl Kesselman, and John Vidale, who have
been so generous with their time. I would like to thank the National Science Foundation
(NSF), whose Research Experiences for Undergraduates program, which funded my
SCEC UseIT internships, turned me from a probable career in the software industry to
iv
one in science. NSF, the U.S. Geological Survey, and the W. M. Keck Foundation each
continue to invest in important research projects (see each chapter’s
acknowledgements section for a listing of these generous grants). Finally, I would like to
thank the University of Southern California for not only providing a wonderful work
environment, but supporting the professional development of their staff through the
extremely generous tuition assistance benefit program which has allowed me to pursue
this degree; Fight on! I count myself truly lucky to have such an incredible network of
support.
v
Table of Contents
Acknowledgements ..........................................................................................................ii
List of Tables .................................................................................................................. vii
List of Figures ................................................................................................................ viii
Abstract ...........................................................................................................................xi
1. Introduction ................................................................................................................. 1
1.1 Probabilistic Seismic Hazard Analysis ................................................................... 1
1.2 Empirical Earthquake Rupture Forecasts for California ......................................... 2
1.3 Empirical Ground Motion Models for California...................................................... 6
1.4 The Ergodic Assumption and Ground Motion Simulation ...................................... 9
1.5 Multi-Cycle Earthquake Sequence Simulators ..................................................... 11
1.6 Tables .................................................................................................................. 13
2. Operational Earthquake Forecasting during the 2019 Ridgecrest, California,
Earthquake Sequence with the UCERF3-ETAS Model ................................................. 16
2.1 Abstract ............................................................................................................... 16
2.2 Overview and Operational Status of UCERF3-ETAS .......................................... 17
2.3 UCERF3-ETAS Ridgecrest Response and Subsequent Model Development ..... 18
2.4 Impact of Modeled Faults and Characteristicness on M 6.4 Aftershock
Probabilities ............................................................................................................... 21
2.5 Spatial Distribution of M 7.1 Aftershock Probabilities ........................................... 24
2.6 Sensitivity to Finite Fault Surfaces ....................................................................... 26
2.7 Sensitivity to Model Parameters .......................................................................... 29
2.8 Discussion and Conclusions ................................................................................ 32
2.9 Data and Resources ............................................................................................ 34
2.10 Acknowledgements ............................................................................................ 34
2.11 Tables ................................................................................................................ 36
2.12 Figures ............................................................................................................... 38
3. Short-Term Forecasting of Earthquake Ground Motions Using Conditional
Hypocenter Distributions ............................................................................................... 47
3.1 Abstract ............................................................................................................... 47
vi
3.2 Introduction .......................................................................................................... 49
3.3 Scenario Definition ............................................................................................... 51
3.4 Methodology ........................................................................................................ 52
3.5 Results ................................................................................................................. 55
3.6 Discussion ........................................................................................................... 58
3.7 Data and Resources ............................................................................................ 60
3.8 Acknowledgements .............................................................................................. 60
3.9 Figures ................................................................................................................. 61
4. A Prototype Probabilistic Seismic Hazard Model for California Constructed with
Fully-Deterministic Physical models .............................................................................. 72
4.1 Abstract ............................................................................................................... 72
4.2 Introduction .......................................................................................................... 73
4.3 Source Model: RSQSim ....................................................................................... 74
4.4 Ground Motion Model: CyberShake ..................................................................... 80
4.5 Ground Motion Variability in 1D: SCEC BroadBand Platform .............................. 81
4.6 Ground Motion Variability In 3D: CyberShake ..................................................... 84
4.7 Procedure to Decompose Variance Through Rotation and Translation ............... 85
4.8 Variance Decomposition Results ......................................................................... 91
4.9 Fully Physics-Based Hazard Curves .................................................................... 93
4.10 Discussion ......................................................................................................... 96
4.11 Data and Resources ........................................................................................ 100
4.12 Acknowledgements .......................................................................................... 101
4.13 Tables .............................................................................................................. 102
4.14 Figures ............................................................................................................. 105
4.15 Supplemental Figures ...................................................................................... 120
5. Discussion and Prospectus ..................................................................................... 125
5.1 Operational Earthquake Forecasting ................................................................. 125
5.2 Non-Ergodic Simulation-Based PSHA ............................................................... 128
References .................................................................................................................. 132
vii
List of Tables
Table 1.1: Co-authored peer-reviewed publications from 2008-2020 ............................ 13
Table 2.1: Aftershock probabilities from UCERF3-ETAS simulations performed with
the same model parameters, varying only the rupture geometries used ....................... 36
Table 2.2: Aftershock probabilities from UCERF3-ETAS simulations performed with
the same input rupture geometries (preferred model, ShakeMap source), varying
only model parameters .................................................................................................. 37
Table 4.1: CyberShake site name, locations, and parameters considered in this
study ............................................................................................................................ 102
Table 4.2: Between-event standard deviations estimated for scenarios defined in
the text ........................................................................................................................ 103
Table 4.3: Within-event standard deviations estimated for scenarios defined in the
text .............................................................................................................................. 103
Table 4.4: Risk-targeted ground motions (RTGM; Luco et al., 2007) computed from
hazard curves from CyberShake and the Abrahamson et al. (2014; henceforth
ASK2014) empirical GMM for each of the sites listed in Table 4.1 .............................. 104
viii
List of Figures
Figure 2.1: Percentile scenarios: Map view of UCERF3-ETAS synthetic catalogs of
aftershocks to the Ridgecrest sequence (preferred model) ........................................... 39
Figure 2.2: Map view of model input events from ComCat from one week preceding
to one week following the M 7.1 earthquake ................................................................. 40
Figure 2.3: UCERF3-ETAS forecasted magnitude-number distribution for one week
following the 4 July M 6.4 earthquake (point source) .................................................... 41
Figure 2.4: UCERF3-ETAS forecasted spatial probability distribution of M ≥ 3.5
aftershocks for 30 days immediately following the M 7.1 earthquake ........................... 42
Figure 2.5: Same as Figure 2.4, but for 30 day probabilities beginning seven days
after the M 7.1 earthquake ............................................................................................ 43
Figure 2.6: Same as Figures 2.4b and 2.5b but with the no-faults UCERF3-ETAS
model ............................................................................................................................ 44
Figure 2.7: Finite rupture surfaces used as input to UCERF3-ETAS simulations .......... 45
Figure 2.8: Cumulative number of M ≥ 3.5 events over time for the region mapped in
Figures 2.4-2.6, starting immediately following the M 7.1 event .................................... 46
Figure 3.1: Map showing the geographic location of the two scenario ruptures (stars)
and the STNI and PDU CyberShake sites (triangles) .................................................... 61
Figure 3.2: Time-dependent probabilities of M ≥ 7 aftershocks within 100 km of a M 6
scenario earthquake near Bombay Beach, as a function of the forecast horizon .......... 62
Figure 3.3: Same as Figure 3.2, except for aftershocks of a Mojave M 6 scenario
earthquake .................................................................................................................... 63
Figure 3.4: Conditional hypocenter distribution (CHD) for one week following the
Bombay Beach M 6 scenario on the Coachella section of the San Andreas fault .......... 65
Figure 3.5: Same as Figure 3.4, but with the CHD for the Mojave M 6 scenario
plotted along-strike and down-dip of the Mojave (South) section of the San Andreas
fault ............................................................................................................................... 66
Figure 3.6: One-week hazard curves with the empirical ASK 2014 ground motion
model for site STNI following the (a) Bombay Beach M 6 and (b) Mojave M 6
scenarios ....................................................................................................................... 67
Figure 3.7: Maps of one-week probability gain 0.2 g 5 s RotD50 spectral acceleration
following the Bombay Beach M 6 (a, b) and Mojave M 6 (c, d) scenarios ..................... 69
Figure 3.8: One-week hazard curves with deterministic CyberShake ground motions
for site STNI following the (a) Bombay Beach M 6 and (b) Mojave M 6 scenarios ........ 69
Figure 3.9: Maps of one-week probability gain 0.2 g 5 s RotD50 spectral acceleration
following the Bombay Beach M 6 (a, b) and Mojave M 6 (c, d) scenarios ...................... 71
ix
Figure 3.10: One-week hazard curves with the (a) CyberShake model and (b)
empirical ASK 2014 ground motion model for site PDU following the Mojave M 6
scenario ......................................................................................................................... 71
Figure 4.1: Ground motion variability controls hazard at large intensities (e.g.,
>0.2 g). This example is computed for the USC site with the Third Uniform
California Earthquake Rupture Forecast (UCERF3) model and Abrahamson et al.
(2014; henceforth ASK2014) ground motion model (GMM), modified with three
different fixed total sigma values for illustration purposes ........................................... 105
Figure 4.2: Three-dimensional (3D) perspective view looking north of faults
considered in Southern California, highlighting an M 7.45 simulated Rate-State
Earthquake Simulator (RSQSim) rupture on the Mojave section of the San Andreas
Fault ............................................................................................................................ 106
Figure 4.3: Two-dimensional side view of the M 7.45 simulated RSQSim rupture
from Figure 4.2 ............................................................................................................ 107
Figure 4.4: 3D perspective view looking north of a complex synthetic M 7.69
RSQSim rupture in the Los Angeles basin .................................................................. 108
Figure 4.5: Mean propagation velocity as a function of patch hypocentral distance
for four different RSQSim parameterizations, each of which incorporates a new
feature over the previous model .................................................................................. 109
Figure 4.6: Slip-time history of a single patch from the rupture depicted in Figures
4.2 and 4.3 .................................................................................................................. 110
Figure 4.7: RotD50 spectra for site USC from ruptures on the Mojave section of the
San Andreas fault, computed with a one-dimensional (1D) velocity structure with
Vs=500 m/s in the Southern California Earthquake Center (SCEC) BroadBand
Platform (BBP) ............................................................................................................ 111
Figure 4.8: 3-second RotD50 ground motion histograms of accelerations for all
7 ≤ M ≤ 7.5 San Andreas (Mojave) RSQSim events (black histogram) for a site at
USC ............................................................................................................................. 112
Figure 4.9: Map view of sites considered in this study, with site locations marked
with triangles next to their names ................................................................................ 113
Figure 4.10: z-score (equation 4.3) distributions (gray histograms) of all ruptures
for the sites from Figure 4.9, compared against the ASK2014 GMM .......................... 114
Figure 4.11: Stacked normal distributions, with the total sum represented with a
thick black line and the cumulative sum after each of N=50 individual distributions
with thin gray lines ....................................................................................................... 115
Figure 4.12: Map view schematic plot of rupture rotation and translation procedure
to emulate empirical records within CyberShake’s reciprocity framework ................... 116
Figure 4.13: Simulated τ(k,n) as a function of the number of synthetic recordings
per rupture (Nrec) for the M 6.6 vertical strike-slip scenario at 50 km distance ........... 117
Figure 4.14: RSQSim simulation hazard curves at USC ............................................. 118
x
Figure 4.15: Source contributions for RSQSim hazard curves computed at USC ....... 119
Figure 4.16: Same as Figure 4.14a, except for (a) site SBSM and (b) site OSI .......... 119
Figure 4.S1: z-score (equation 4.3) distributions (gray histograms) of all RSQSim
ruptures for each individual site (locations mapped in Figure 4.9 and listed in Table
4.1), simulated with the Southern California Earthquake Center (SCEC)
CyberShake platform and compared against the ASK2014) empirical ground
motion model (GMM) ................................................................................................... 121
Figure 4.S2: RSQSim simulation hazard curves for each individual site (locations
mapped in Figure 9 and listed in Table 1) ................................................................... 124
xi
Abstract
Ground motions radiated by earthquakes pose a significant hazard to human life
and the built environment. I present three papers which contribute to the goal of more
accurately forecasting earthquake hazard on multiple timescales. The first gives results
and lessons learned from short-term operational earthquake forecasting during the 2019
Ridgecrest, CA, earthquake sequence. The second demonstrates how consideration of
short-term perturbations to the conditional hypocenter distribution can achieve a ground
motion probability gain of an order of magnitude relative to short-term models which
ignore rupture directivity. The third introduces a first-of-its-kind probabilistic seismic
hazard model that is constructed with fully-deterministic physics-based models. That
model is non-ergodic and forecasts lower long-term hazard for some sites in the Los
Angeles region than empirical ergodic models due to its characterization of repeated
source, site, and path effects, while still containing a sufficient amount of total ground
motion variability.
1
1. Introduction
1.1 Probabilistic Seismic Hazard Analysis
Ground motions radiated by earthquakes pose a significant hazard to human life
and the built environment, especially in active tectonic regions including California.
Despite strong societal interest, the scientific community has yet to identify any
precursory signals that reliably predict the timing, location, and magnitude of
earthquakes (Jordan et al., 2011). Instead, probabilistic seismic hazard analysis (PSHA)
formulated by Cornell (1968) has become the primary method for quantifying seismic
hazard. The fundamental product of PSHA, the hazard curve, expresses the probability
of exceeding an intensity measure level (IML) at a site for some period of time,
combined across all possible earthquakes in the surrounding region and their expected
ground motions. Hazard curves are especially useful to structural engineers, who aim to
design buildings that will withstand a wide variety of potential earthquakes. For
example, the city of Los Angeles faces hazards from many nearby low or moderate slip
rate faults (e.g., the Puente Hills, Hollywood, and Santa Monica faults), each with a
relatively low probability on short timespans, as well as from the San Andreas, a high
slip rate plate boundary fault 50 km northeast of downtown generally considered more
probable. PSHA aims to aggregate the probabilities from all sources of damaging
earthquakes into useful products for seismic design.
2
PSHA employs two major categories of models: Earthquake Rupture Forecasts
(ERFs) and Ground Motion Models (GMMs). ERFs aim to enumerate the magnitude,
geographic extent, and probability (over some timespan) of each of the possible
damaging earthquakes in a region. Some ERFs are time-dependent, for example by
incorporating Reid (1910)’s elastic-rebound theory, but most are time-independent
Poissonian models (the latter are used to establish building codes, as buildings are
designed to survive multiple seismic cycles). For each earthquake described by an
ERF, GMMs predict the distribution of expected ground motions at a site of interest.
Various intensity measures (IMs) exist for quantifying the intensity of an earthquake
ground motion; the most common in current practice is 5% damped Rotd50 spectral
acceleration, the median (across all possible horizontal orientations) peak acceleration
of a damped, single degree of freedom oscillator with a given resonant frequency.
1.2 Empirical Earthquake Rupture Forecasts for California
The most comprehensive and complex ERFs developed to date are for
California, developed by the Working Group on California Earthquake Probabilities
(WGEP). These models assimilate data from paleoseismology, geologic and geodetic
slip rate estimates, and seismology in an effort to describe all possible ruptures in
California and assign time-dependent probabilities. The first WGCEP (1988) model
assumed strict fault segmentation (i.e., that ruptures start and stop only at predefined
fault segment boundaries) and focused only on major, well-studied faults. This strict
segmentation approach did not allow multiple segments of the San Andreas Fault to
rupture together coseismically, thus excluding the major 1906 and 1857 observed
3
earthquakes (Field, 2007). This shortcoming was known to model developers, but was
included because it allowed the modeler to compute earthquake rates directly from a
slip rate or paleoseismic event rate observations. Multi-segment ruptures introduce non-
uniqueness to earthquake rate calculations; WGCEP (1995) included them through a
moment-balanced multi-step procedure which approximated the Gutenberg-Richter (G-
R; 1944) magnitude-frequency distribution. Instead, WGCEP (2003) and the first
published UCERF model, UCERF2 (Field et al., 2009), set multi-segment rupture rates
through aggregated expert opinion. Both of these models forecasted a higher rate of
moderate (𝗠 ∈ [6.5,7]) earthquakes than that predicted by a statewide G-R distribution
(referred to colloquially as “the bulge”), though unlike WGCEP (2003), the UCERF2
rates fell within the 95% confidence bounds of the empirical target magnitude-frequency
distribution.
The Third Uniform California Earthquake Rupture Forecast (UCERF3; Field et
al., 2014; Field, Jordan, et al., 2017) is the most recent and sophisticated of its kind.
UCERF3 is composed of three submodels, each corresponding to different timescales:
long-term time-independent (UCERF3-TI), long-term time-dependent (UCERF3-TD),
and a short-term spatiotemporal clustering model (UCERF3-ETAS). I was the lead
computer programmer and a core participant throughout the UCERF3 project (see
Table 1.1 for a list of co-authored work related to UCERF3). UCERF3-TI addressed
many of the shortcomings of prior models, largely through a new “grand inversion”
methodology (Page et al., 2014) that can solve for the rate of each possible earthquake
while better satisfying all available data constraints. A key difference between UCERF3
4
and prior models is that UCERF3 relaxes the fault segmentation requirement by
breaking each fault into a set of “subsections,” each with length equal to approximately
half of the seismogenic fault width. UCERF3 considers two types of ruptures:
supraseismogenic on-fault ruptures, unique sets of two or more subsections, and
gridded seismicity, a proxy for both smaller subseismogenic on-fault ruptures and
background seismicity on unmodeled faults. Milner et al. (2013) describes the
plausibility filter process by which individual supraseismogenic ruptures were selected
from the set of all possible subsection combinations. This process was designed to be
inclusive of all ruptures from prior models, with the addition of multi-fault ruptures that
passed various criteria including a maximum jump distance of 5 km, and a fault-to-fault
azimuthal change of 60 degrees. Inclusion of complex multi-fault ruptures was a point of
controversy during model development despite prior examples in nature and
subsequent validations including the 2016 Kaikoura M 7.8. Once the set of possible
ruptures was defined, a simulated annealing solver inverted for the rate of each rupture
(the “grand inversion”) minimizing misfit to data constraints such as slip rates, regional
magnitude frequency distributions, and paleoseismic event rates. This new
methodology fit all available data better than did prior models (including removing the
glut of extra 𝗠 ∈ [6.5,7] events), though the inversion was underdetermined, nonunique,
and model trade offs exist.
Through Reid (1910) renewal models, all WGCEPs have incorporated long-term
time-dependence on well-studied faults where the date of the last event is known.
UCERF3-Time-Dependent (UCERF3-TD; Field et al., 2015) extended this to faults for
5
which the date of last event is unknown, through the addition of an historical open
interval (the last event must have been prior to ~1875), integrating over probabilities
from a range of possible prehistoric dates of last event. The largest probability gains
from UCERF3-TD over UCERF3-TI are relatively modest, a factor of two on the
southern San Andreas Fault, but it is possible to forecast much larger short-term
probability gains during active aftershock sequences. To that end, we developed
UCERF3-ETAS (Field, Milner, et al., 2017), the first short-term, time-dependent
WGCEP model. Epidemic-Type Aftershock Sequence (ETAS) models have been in
development since the 1980s (Ogata, 1988), but historically have only been applied to
point source earthquake representations, ignoring the presence of mapped faults. This
is contrary to our intuition that aftershock probabilities may be higher near active faults.
UCERF3-ETAS formalizes that intuition into a quantitative forecast, merging the ETAS
point process model and UCERF3-TD into a hybrid model where on-fault
supraseismogenic ruptures can trigger point source events, and vise versa. The model
outputs ensembles of stochastic aftershock sequence catalogs, from which various
statistics including fault triggering probabilities can be calculated. A major conclusion of
Field, Milner, et al. (2017) is that incorporation of elastic rebound is essential for a fault-
based ETAS model; without elastic rebound, ruptures on faults preferentially trigger
themselves leading to runaway sequences.
UCERF3-ETAS is designed to be operationalizable and, although not yet running
automatically in response to recent seismicity, is effective and informative as an on-
demand resource. Chapter 2 describes how the model was run, improved upon, and its
6
results disseminated during the July 2019 Ridgecrest earthquake sequence (Milner et
al., 2020). I had the first UCERF3-ETAS aftershock simulations running on a local high
performance computing cluster 33 minutes after the July 4 2019 M 6.4, and updated the
model results throughout the sequence. The scientific community and public at large
were interested in the probability that the sequence would trigger a large rupture on the
nearby Garlock fault; this is a question that UCERF3-ETAS is uniquely able and
explicitly designed to quantify. Throughout the sequence, I added many new features to
the model including the ability to fetch finite fault surface estimates directly from either
U.S. Geological Survey (USGS) ShakeMap or finite-fault inversion data products. We
found that estimated triggering probabilities were sensitive to rupture surface details in
expected and surprising manners, which are detailed in Chapter 2. Despite these
sensitivities, we believe that the model can be useful to stakeholders and the public at
large.
1.3 Empirical Ground Motion Models for California
GMMs describe the probability of exceeding an IML at a site, conditioned on the
occurrence of an earthquake (be it an observed earthquake, scenario event, or, in
PSHA, one supplied by the ERF). The most common GMMs are empirical models which
regress a curated dataset of observed earthquake ground motion records against
explanatory variables (including distance, magnitude, fault style, and site conditions) to
produce a log-normal probability density distribution of expected IMLs. Historically,
these models have been referred to as Ground Motion Prediction Equations or
Attenuation Relationships by the literature, though more recently those terms have
7
fallen out of use in favor of “GMM”, a more general term which encompasses multiple
types of models.
The Next Generation Attenuation Models (NGA) project collected and curated
data to build a comprehensive dataset and facilitated development of empirical GMMs
for use in PSHA. The first NGA-West project (Choiu et al., 2006; Abrahamson et al.,
2008) produced five GMMs that were used in the 2008 update to the USGS National
Seismic Hazard Model (NSHM). This effort was a major step forward, as
standardization of both input data and model applicability allowed all five models to be
directly compared and judged against the same data, and represented as an epistemic
uncertainty in PSHA calculations. Each model used a different subset of available
explanatory variables to describe a range of physical processes; one example of this is
basin amplification, in which seismic waves reverberate in deep sedimentary basins
resulting in constructive interference and long duration of shaking. Basin response is an
inherently three-dimensional (3D) problem, but Field (2000) demonstrated that a scalar
proxy, the depth to a fixed shear wave velocity (Vs), can account for average
amplification. Campbell and Bozorgnia (2008) implemented this using the depth to Vs =
2.5 km/s isosurface proxy (Z2.5), Abrahamson and Silva (2008) and Choiu (2008) the
depth to Vs = 1.0 km/s isosurface proxy (Z1.0), and Boore and Atkinson (2008) and Idriss
(2008) did not include basin terms.
The NGA-West2 project (Borzognia et al., 2014) produced the current state-of-
the-art GMMs for active tectonic regions used in the 2014 update to the NSHM. It
8
improved upon NGA-West by expanding the dataset to include significantly more
observed earthquake data, and switched to the RotD50 horizontal component (Boore,
2010), the median orientation-independent IML at a site (previous approaches used the
geometric-mean of each horizontal seismogram component). Basin modeling was also
expanded in NGA-West2; Campbell and Bozorgnia (2014) uses the Z2.5 isosurface
proxy, Abrahamson et al. (2014), Boore et al. (2014), and Choiu and Youngs (2014) the
Z1.0 isosurface proxy, and only Idriss (2014) omits a basin response model.
NGA-West2 also produced models of rupture directivity (Spudich et al., 2014), a
phenomena where ground motions are larger for sites in the direction of rupture
propagation. These models can be used in conjunction with the NGA-West2 GMMs for
earthquakes with known hypocenters and finite fault extents. That said, directivity
models are rarely used in PSHA, because ERFs (including UCERF3-TI/TD) typically
assume a uniform probability density distribution of hypocenter location along-strike of a
potential rupture surface. Following Wang and Jordan (2014), we call that distribution
the “conditional hypocenter distribution” (CHD). NGA-West2 directivity models are zero-
centered, which means that directivity has no net effect when averaged over a uniform
CHD. Chapter 3 discusses how UCERF3-ETAS can be used to forecast the CHD in
short-term ground-motion forecasts that achieve an order of magnitude of gain over
models that ignore the CHD.
9
1.4 The Ergodic Assumption and Ground Motion Simulation
Observational data for the large earthquakes which control seismic hazard at
long return periods are often sparse for a specific region of interest, so modelers
exercise the assumption that the distribution of ground motion recordings from global
earthquakes can forecast the distribution of expected ground motions at a site (“the
ergodic assumption”, see Anderson and Brune, 1999). One consequence of ergodic
modeling is large aleatory variability terms in GMMs, typically represented with a total
standard deviation (𝜎 ) greater than 0.7 in natural-log g units. This variability is--by
construction--irreducible in ergodic models, despite its significant contribution to design-
level hazard (see Figure 4.1 and Anderson and Brune, 1999). Variance is further
decomposed into two categories: within-event variability (with standard deviation 𝜑 )
which describes the expected spatial variability in ground motions for a particular
rupture, and between-events variability (with standard deviation 𝜏 ) in overall spatially
averaged intensity between different ruptures of similar description (e.g., magnitude,
faulting style) (Al Atik et al., 2010). 𝜑 is typically the larger component, and could be
reduced in a non-ergodic model which accounts for repeated path and site amplification
effects for individual site-source pairs. Non-ergodic empirical GMMs are an area of
active research (e.g., Landwehr et al., 2016), but have not yet become standard
practice for PSHA calculations.
The Southern California Earthquake Center (SCEC) CyberShake platform
(Graves et al., 2011; Jordan et al., 2018) is a non-ergodic PSHA model that replaces
parameterized distributions from empirical GMMs with ensembles of 3D ground motion
10
simulations through a realistic 3D seismic velocity field. In cybershake, each UCERF2
rupture is extended with a suite of kinematic rupture variations, each of which has a
unique hypocenter and slip-time distribution that is generated with the Graves and
Pitarka (2014; henceforth GP2014) method. GP2014 does not support automated
rupture generation for the complex, multi-fault ruptures abundant in UCERF3,
necessitating CyberShake’s dependence on the prior UCERF2 model. CyberShake
studies compute synthetic seismograms for each rupture variation at a site, which are
simulated through the latest SCEC community velocity models (e.g., Lee et al., 2014,
and Small et al., 2017), and extract IMs (e.g., RotD50) from each synthetic seismogram
in order to compute seismic hazard curves. CyberShake utilizes seismic reciprocity to
efficiently simulate ground motions from hundreds of thousands of rupture variations at
each site rather than forward simulation, which would be intractable on even the largest
supercomputers available today. Despite this efficiency, CyberShake studies require
significant computing resources--currently around 1000 compute node-hours per site on
Oak Ridge Leadership Computing Center’s Summit supercomputer, which was listed as
the most powerful system on the November 2019 Top 500 list.
CyberShake simulations capture 3D basin excitation, which Olsen et al. (2006)
found can be magnified by rupture directivity, a dynamic not captured by current
empirical approaches (which include the Z1.0 and Z2.5 proxies). In Chapter 3, we
compare short-term ground motion forecasts computed with empirical directivity models
with those from CyberShake, and find that CyberShake’s strong directivity-basin
coupling yields higher hazard in the Los Angeles basin. Although CyberShake simulates
11
the physics of wave propagation, it currently relies on the UCERF2 and GP2014
statistical models, making the oft-applied “physics-based model” label not entirely
accurate.
1.5 Multi-Cycle Earthquake Sequence Simulators
Multi-cycle earthquake sequence simulators have become a useful tool for
quantifying earthquake recurrence and for use in PSHA studies. They simulate
hundreds to millions of years of the seismic cycle on either a single fault or for a
complicated fault system, including rupture nucleation, propagation, coseismic stress
change, and interseismic strain accumulation. Faults are typically broken up into small
(~1 km
2
) triangular or rectangular patches, and model outputs include synthetic catalogs
of ruptures on those patches, each with a magnitude, time, hypocenter, and slip
distribution. Development of UCERF3-TD was informed by interevent time distributions
from multiple earthquake simulators, including the Rate-State Earthquake Simulator
(RSQSim), developed by Dieterich & Richards-Dinger (2010, 2012). RSQSim builds
upon the rate- and state-dependent friction law proposed by Dieterich (1992).
RSQSim’s marked efficiency is largely due to use of various analytic approximations,
notably a quasi-static approximation to compute coseismic stress changes. Although a
fully-dynamic code that simulates dynamic stress changes during seismic wave
propagation would be preferable, these efficiencies allow for simulation of hundreds of
thousands of years of earthquakes on the complex California fault system. Catalogs of
seismicity produced by RSQSim also include many complex multi-fault ruptures,
consistent with recent observations and similar to those in the UCERF3 model.
12
Shaw et al. (2018), which I co-authored, found strong agreement in hazard maps
and fault recurrence intervals between UCERF3 and an RSQSim-based ERF. That
study used a long RSQSim catalog simulated on the UCERF3 fault system, and each
RSQSim rupture was used to parameterize an empirical GMM for PSHA calculations.
RSQSim simulations also produce full slip-time histories for each rupture, and can thus
be used as an extended ERF for ground motion simulation in CyberShake, taking the
place of both UCERF2 and GP2014. Chapter 4 introduces a prototype RSQSim-
CyberShake model--the first of its kind composed entirely of deterministic, physics-
based models. The RSQSim-CyberShake model is non-ergodic, and Chapter 4
demonstrates how hazard curves for individual sites can differ significantly from ergodic
models due to lower variability for repeated source-site paths. Importantly, the model
still matches the total ground motion variability that has been observed empirically when
aggregated across multiple sites. As RSQSim or similar models evolve to incorporate
more earthquake physics, they might play an increased role in future PSHA studies, and
ultimately help reduce the uncertainties which drive hazard.
13
1.6 Tables
Table 1.1
Co-authored peer-reviewed publications from 2008-2020. There are 26 papers in total
with a total of 1168 citations. Source: Google Scholar
(https://scholar.google.com/citations?hl=en&user=hWS0lpgAAAAJ, accessed 4 March
2020).
Year Reference
2008
Field, E. H., & Milner, K. R. (2008). Forecasting California's earthquakes: What can we
expect in the next 30 years? (No. 2008-3027). US Geological Survey.
2008
Callaghan, S., Maechling, P., Deelman, E., Vahi, K., Mehta, G., Juve, G., Milner, K.,
Graves, R., Field, E., Okaya, D. and Gunter, D., 2008, December. Reducing time-to-
solution using distributed high-throughput mega-workflows-experiences from SCEC
CyberShake. In 2008 IEEE Fourth International Conference on eScience (pp. 151-158).
IEEE.
2009
Milner, K., Becker, T. W., Boschi, L., Sain, J., Schorlemmer, D., & Waterhouse, H. (2009).
New software framework to share research tools. Eos, Transactions American Geophysical
Union, 90(12), 104-104.
2010
Callaghan, S., Deelman, E., Gunter, D., Juve, G., Maechling, P., Brooks, C., Vahi, K.,
Milner, K., Graves, R., Field, E. and Okaya, D., 2010. Scaling up workflow-based
applications. Journal of Computer and System Sciences, 76(6), pp.428-446.
2011
Graves, R., Jordan, T.H., Callaghan, S., Deelman, E., Field, E., Juve, G., Kesselman, C.,
Maechling, P., Mehta, G., Milner, K. and Okaya, D., 2011. CyberShake: A physics-based
seismic hazard model for southern California. Pure and Applied Geophysics, 168(3-4),
pp.367-381.
2011
Callaghan, S., Maechling, P., Small, P., Milner, K., Juve, G., Jordan, T.H., Deelman, E.,
Mehta, G., Vahi, K., Gunter, D. and Beattie, K., 2011. Metrics for heterogeneous scientific
workflows: A case study of an earthquake science application. The International Journal of
High Performance Computing Applications, 25(3), pp.274-285.
2012
Porter, K. A., Field, E. H., & Milner, K. (2012). Trimming the UCERF2 hazard logic tree.
Seismological Research Letters, 83(5), 815-828.
2012
Parsons, T., Field, E. H., Page, M. T., & Milner, K. (2012). Possible earthquake rupture
connections on mapped California faults ranked by calculated Coulomb linking stresses.
Bulletin of the Seismological Society of America, 102(6), 2667-2676.
2013
Field, E.H., Biasi, G.P., Bird, P., Dawson, T.E., Felzer, K.R., Jackson, D.D., Johnson, K.M.,
Jordan, T.H., Madden, C., Michael, A.J., Milner, K.R., Page, M.T., Parsons, T., Powers,
P.M., Shaw, B.E., Thatcher, W.R., Weldon, R.J., II, and Zeng, Y., 2013, Uniform California
earthquake rupture forecast, version 3 (UCERF3)—The time-independent model: U.S.
Geological Survey Open-File Report 2013–1165, 97 p., California Geological Survey
14
Special Report 228, and Southern California Earthquake Center Publication 1792,
http://pubs.usgs.gov/of/2013/1165/.
2013
Milner, K. R., Page, M. T., Field, E. H., Parsons, T., Biasi, G. P., & Shaw, B. E. (2013).
Appendix T—Defining the inversion rupture set using plausibility filters. US Geol. Surv.
Open‐File Rept. 2013‐1165.
2014
Edward H. Field, Ramon J. Arrowsmith, Glenn P. Biasi, Peter Bird, Timothy E. Dawson,
Karen R. Felzer, David D. Jackson, Kaj M. Johnson, Thomas H. Jordan, Christopher
Madden, Andrew J. Michael, Kevin R. Milner, Morgan T. Page, Tom Parsons, Peter M.
Powers, Bruce E. Shaw, Wayne R. Thatcher, Ray J. Weldon, Yuehua Zeng; Uniform
California Earthquake Rupture Forecast, Version 3 (UCERF3)—The Time‐Independent
Model. Bulletin of the Seismological Society of America ; 104 (3): 1122–1180. doi:
https://doi.org/10.1785/0120130164
2014
Page, M. T., Field, E. H., Milner, K. R., & Powers, P. M. (2014). The UCERF3 grand
inversion: Solving for the long‐term rate of ruptures in a fault system. Bulletin of the
Seismological Society of America, 104(3), 1181-1204.
2014
Parsons, T., Segou, M., Sevilgen, V., Milner, K., Field, E., Toda, S., & Stein, R. S. (2014).
Stress‐based aftershock forecasts made within 24 h postmain shock: Expected north San
Francisco Bay area seismicity changes after the 2014 M= 6.0 West Napa earthquake.
Geophysical Research Letters, 41(24), 8792-8799.
2014
Edward H. Field, Glenn P. Biasi, Peter Bird, Timothy E. Dawson, Karen R. Felzer, David D.
Jackson, Kaj M. Johnson, Thomas H. Jordan, Christopher Madden, Andrew J. Michael,
Kevin R. Milner, Morgan T. Page, Tom Parsons, Peter M. Powers, Bruce E. Shaw, Wayne
R. Thatcher, Ray J. Weldon, Yuehua Zeng; Long‐Term Time‐Dependent Probabilities for
the Third Uniform California Earthquake Rupture Forecast (UCERF3). Bulletin of the
Seismological Society of America ; 105 (2A): 511–543. doi:
https://doi.org/10.1785/0120140093
2017
Field, E. H., Milner, K. R., Hardebeck, J. L., Page, M. T., van der Elst, N., Jordan, T. H., ...
& Werner, M. J. (2017). A spatiotemporal clustering model for the third Uniform California
Earthquake Rupture Forecast (UCERF3‐ETAS): Toward an operational earthquake
forecast. Bulletin of the Seismological Society of America, 107(3), 1049-1081.
2017
Porter, K., Field, E., & Milner, K. (2017). Trimming a hazard logic tree with a new model-
order-reduction technique. Earthquake spectra, 33(3), 857-874.
2017
Field, E. H., Jordan, T. H., Page, M. T., Milner, K. R., Shaw, B. E., Dawson, T. E., ... &
Weldon, R. J. (2017). A synoptic view of the third Uniform California Earthquake Rupture
Forecast (UCERF3). Seismological Research Letters, 88(5), 1259-1267.
2017
Field, E., Porter, K., & Milner, K. (2017). A Prototype Operational Earthquake Loss Model
for California Based on UCERF3-ETAS–A First Look at Valuation. Earthquake spectra,
33(4), 1279-1299.
2018
Jordan, T. H., Callaghan, S., Graves, R. W., Wang, F., Milner, K. R., Goulet, C. A.,
Maechling, P. J., Olsen, K. B., Cui, Y., Juve, G., Vahi, K., Yu, J., Deelman, E., and Gill, D.
(2018, January). CyberShake models of seismic hazards in Southern and Central
California. In Proceedings of the US National Conference on Earthquake Engineering.
2018
Crouse, C. B., Jordan, T. H., Milner, K. R., Goulet, C. A., Callaghan, S., & Graves, R. W.
(2018, June). Site‐specific MCER response spectra for Los Angeles region based on 3‐D
15
numerical simulations and the NGA West2 equations. In Oral Presentation at 11th National
Conference in Earthquake Engineering. Paper (Vol. 518).
2018
Moschetti, M. P., Chang, S., Crouse, C. B., Frankel, A., Graves, R., Puangnak, H., Luco,
N., Goulet, C. A., Rezaeian, S., Shumway, A., Powers, P. M., Petersen, M. D., Callaghan,
S., Jordan, T. H., & Milner, K. R. (2018, June). The science, engineering applications, and
policy implications of simulation‐based PSHA. In Proceedings of the 11th National
Conference in Earthquake Engineering (11NCEE), June (pp. 25-19).
2018
Field, E. H., & Milner, K. R. (2018). Candidate products for operational earthquake
forecasting illustrated using the HayWired planning scenario, including one very quick (and
not‐so‐dirty) hazard‐map option. Seismological Research Letters, 89(4), 1420-1434.
2018
Shaw, B. E., Milner, K. R., Field, E. H., Richards-Dinger, K., Gilchrist, J. J., Dieterich, J. H.,
& Jordan, T. H. (2018). A physics-based earthquake simulator replicates seismic hazard
statistics across California. Science advances, 4(8), eaau0688.
2019
Minson, S.E., Baltay, A.S., Cochran, E.S., Hanks, T.C., Page, M.T., McBride, S.K., Milner,
K.R. and Meier, M.A., 2019. The limits of earthquake early warning accuracy and best
alerting strategy. Scientific reports, 9(1), pp.1-13.
2019
Michael, A.J., McBride, S.K., Hardebeck, J.L., Barall, M., Martinez, E., Page, M.T., van der
Elst, N., Field, E.H., Milner, K.R. and Wein, A.M., 2020. Statistical Seismology and
Communication of the USGS Operational Aftershock Forecasts for the 30 November 2018
M w 7.1 Anchorage, Alaska, Earthquake. Seismological Research Letters, 91(1), pp.153-
173.
2020
Milner, K. R., Field, E. H., Savran, W. H., Page, M. T., & Jordan, T. H. (2020). Operational
Earthquake Forecasting during the 2019 Ridgecrest, California, Earthquake Sequence with
the UCERF3‐ETAS Model. Seismological Research Letters.
16
2. Operational Earthquake Forecasting during the 2019
Ridgecrest, California, Earthquake Sequence with the UCERF3-
ETAS Model
Kevin R. Milner, Edward H. Field, William H. Savran, Morgan T. Page, Thomas H.
Jordan
2.1 Abstract
The first Uniform California Earthquake Rupture Forecast, Version 3-epidemic-
type aftershock sequence (UCERF3-ETAS) aftershock simulations were running on a
high-performance computing cluster within 33 minutes of the 4 July 2019, M 6.4 Searles
Valley earthquake. UCERF3-ETAS, an extension of the third Uniform California
Earthquake Rupture Forecast (UCERF3), is the first comprehensive, fault-based,
epidemic-type aftershock sequence model. It produces ensembles of synthetic
aftershock sequences both on and off explicitly modeled UCERF3 faults to answer a
key question repeatedly asked during the Ridgecrest sequence: What are the chances
that the earthquake that just occurred will turn out to be the foreshock of an even bigger
event?
As the sequence unfolded--including one such larger event, the 5 July 2019 M
7.1 Ridgecrest earthquake almost 34 hours later--we updated the model with observed
aftershocks, finite rupture estimates, sequence-specific parameters, and alternative
UCERF3-ETAS variants. Although configuring and running UCERF3-ETAS at the time
17
of the earthquake was not fully automated, considerable effort had been focused in
2018 on improving model documentation and ease of use with a public GitHub
repository, command line tools, and flexible configuration files. These efforts allowed us
to quickly respond and efficiently configure new simulations as the sequence evolved.
Here, we discuss lessons learned during the Ridgecrest sequence, including
sensitivities of fault triggering probabilities to poorly constrained finite-rupture estimates
and model assumptions, as well as implications for UCERF3-ETAS operationalization.
2.2 Overview and Operational Status of UCERF3-ETAS
The epidemic-type aftershock sequence (ETAS; Ogata, 1988) extension to the
Third Uniform California Earthquake Rupture Forecast (UCERF3; Field, Jordan, et al.,
2017), referred to as UCERF3-ETAS (Field, Milner, et al., 2017), is the first
comprehensive, fault-based, ETAS model. It produces ensembles of synthetic
aftershock sequences (Figure 2.1) both on and between modeled UCERF3 faults to
answer a key question in earthquake forecasting: What are the chances that an
earthquake that just occurred will turn out to be the foreshock of an even bigger event?
Standard short-term forecasting models (e.g., Reasenberg and Jones, 1989) currently
employed by the U.S. Geological Survey (USGS) and evaluated by the California
Earthquake Prediction Evaluation Council (CEPEC) do not explicitly consider the
proximity of seismic activity to mapped faults. However, Jordan and Jones (2010) report
that CEPEC does respond differently to seismicity near certain faults, for example, their
“go-to-war” protocol in the event of an M ≥ 6 within 3 km of the southern San Andreas
18
fault. UCERF3-ETAS was explicitly designed to formalize this intuition that large
aftershocks are more likely near active faults into an operational forecast model.
Although UCERF3-ETAS is not yet fully operationalized (i.e., running
automatically in response to recent seismicity), largely due to the significant investment
in computational infrastructure required, it is available on demand. Over the summer of
2018, we significantly improved documentation and developed tools to make
configuration and deployment of the model more intuitive, culminating in a training
session at the Southern California Earthquake Center (SCEC) 2018 annual meeting.
Documentation and instructions are posted in a public GitHub repository (see Data and
Resources), although high-performance computing resources are often required to
generate a suite of synthetic catalogs large enough (e.g., 10,000-100,000 catalogs) for
statistically stable results in short time frames. Each simulation ensemble discussed
here took approximately 8 hours to generate 100,000 synthetic catalogs on 3600
processors, which would take 120 days on a single processor.
2.3 UCERF3-ETAS Ridgecrest Response and Subsequent Model
Development
By 11:07 a.m. (PDT) on 4 July 2019 (33 minutes after the M 6.4 Searles Valley
earthquake), the first UCERF3-ETAS aftershock simulations were running at the
University of Southern California’s high-performance computing center. Although the full
suite of 100,000 simulations would not finish until 04:02 p.m. that day, the first results
from an initial subset of 3500 simulations were posted to SCEC’s post-earthquake
19
response page (see Data and Resources) at 11:39 a.m., just over an hour after the M
6.4, and were updated regularly as the simulation progressed. The model forecasted a
2.8% chance of an M ≥ 6.4 aftershock in the next week.
A larger M 7.1 event occurred the next day, 5 July 2019 at 08:19 p.m., when the
primary UCERF3-ETAS developer (Milner) was unavailable to run the model. Based on
online documentation and information presented during the 2018 training session,
another researcher (Savran) was able to submit simulations of the larger earthquake
(then reported as M 6.9) less than an hour later, and run follow-up simulations with the
updated magnitude (M 7.1) and depth (17 km, although later it would be further revised
to 8 km) at 04:25 a.m. the following morning. Those revised simulations, which used a
point-source representation of the M 7.1, were completed and posted for review just
before CEPEC, a committee that advises the governor of California about earthquake
predictions, forecasts, and hazards, convened at 10:30 a.m. The probability of this
sequence triggering a large event on the Garlock fault was of particular concern to
CEPEC, which included in its scientific review the 0.65% probability of such an M ≥ 7
event forecasted by UCERF3-ETAS at that time.
Fault trigger probabilities are sensitive to the distance between the input ruptures
and modeled faults, so, as the finite extent of the rupture became apparent (initially
through the spatial distribution of aftershock locations), we updated our simulations with
finite surface estimations. Prior to the Ridgecrest sequence, all finite rupture surfaces
considered in UCERF3-ETAS were on previously modeled UCERF3 faults, and no
20
capability existed to configure simulations with arbitrarily defined finite surfaces. We
implemented this capability later on 6 July and resubmitted simulations with a fault
surface drawn based on observed seismicity. Over the following week, we updated the
model to fetch earthquake data directly from the USGS Comprehensive Catalog
(ComCat), as well as the finite rupture surfaces used to generate ShakeMaps (Wald et
al., 1999). We also improved post-processing routines to return probabilities for multiple
time frames (e.g., day, week, month, year) and added comparisons with the catalog of
actual aftershocks from ComCat, which automatically update periodically. Simulation
results from these and future simulations are publicly available, posted on GitHub (see
Data and Resources).
We also added another potentially useful visualization, percentile scenarios,
shown in Figure 2.1. Here, we plot a subset of catalogs from the simulated ensemble
with total simulated aftershock counts that fall at various percentile markers. This view
allows users to quickly contextualize typical, extreme, and worst-case scenarios. Figure
2.1a shows a “typical” simulation with a median number of total simulated aftershocks
and no aftershocks above M 6. Figure 2.1b shows a scenario catalog at the 97.5th
percentile in which the M 7.1 is a foreshock to an even larger earthquake on the Garlock
fault. Figure 2.1c shows the most productive simulated catalog from the suite of 100,000
simulations in which an M ≥ 7 Garlock aftershock triggered an M ≥ 8 earthquake on the
San Andreas fault. Although other catalogs around each percentile will not necessarily
have the same history, or even the same fault(s) triggered, these scenarios can offer a
more intuitive snapshot of the scale of aftershock activity than would be possible from
21
just a probability table. Figure 2.1d shows observed aftershocks for the same region
and time period. The count of observed M ≥ 2.5 events in this region, 1924, is lower
than that from the least productive simulated catalog, 1967, though this comparison
does not account for time-dependent magnitude-of-completeness. The count for M ≥ 3.5
events, 304, lies at the 81st percentile of the simulated catalog distribution.
2.4 Impact of Modeled Faults and Characteristicness on M 6.4 Aftershock
Probabilities
UCERF3 models earthquakes both on and between explicitly defined faults.
Long-term, time-independent, supraseismogenic event rates on modeled faults are
determined by the UCERF3 “grand inversion” (Page et al., 2014) with constraints that
include geologically or geodetically estimated fault slip rates, paleoseismic event rates
(where available), and the total regional magnitude-frequency distribution (MFD). In the
context of UCERF3, supraseismogenic describes earthquakes with length at least as
long as the seismogenic fault width. UCERF3 imposes a statewide Gutenberg-Richter
(GR) distribution (Gutenberg and Richter, 1944) but has large earthquakes occurring
preferentially as characteristic earthquakes on modeled faults. The aggregate MFD
across all modeled faults follows a characteristic distribution, meaning a surplus of large
(M ≥ 7) events compared to the GR distribution. Consequently, the off-fault MFD is
"anti-characteristic," with a deficit of large earthquakes compared to the GR distribution
(see Figure 2 of Field, Milner, et al., 2017).
22
This on- and off-fault dichotomy affects UCERF3-ETAS results for the 4 July M
6.4, with a hypocentral distance of 12.2 km to the nearest UCERF3 fault, Airport Lake
(Figure 2.2). Faults in UCERF3 are each treated as proxies for a more complex but
poorly understood fault system, and are assigned polygons with geologic constraints
and buffers, as described in Powers and Field (2013). The epicentral distance to the
Airport Lake fault polygon is 2.7 km, which means that the M 6.4 nucleated in a purely
off-fault area of the UCERF3 model. Furthermore, the nearby Airport Lake and Little
Lake faults are both anti-characteristic in the UCERF3 model, as shown in Figure 4 of
Field, Milner, et al. (2017). The resulting one-week forecasted magnitude-number
distribution is shown in Figure 2.3a, which gives a 2.8% probability of an M ≥ 6.4
aftershock, and a 0.29% chance of an M ≥ 7.1 event (which, given the occurrence of the
latter, had a 59% chance of occurring on an explicitly modeled fault versus 41% off-
fault). These probabilities are lower than those predicted by the no-faults variant of
UCERF3-ETAS (4.5% and 0.90% respectively) which ignores the presence of mapped
faults and imposes a uniform GR MFD everywhere (Figure 2.3b).
The USGS initially posted a 2.3% chance of an M ≥ 7 event in the first week on
its aftershock forecast product page. This was based on a version of the Reasenberg
and Jones (1989) model updated in Page et al. (2016). The corresponding probabilities
from the UCERF3-ETAS and no-faults UCERF3-ETAS models were 0.4% and 1.2%,
respectively. The USGS model includes regionally localized parameters from
Hardebeck et al. (2018) that take into account higher aftershock productivity observed in
geothermal areas in California, as well as epistemic uncertainty in sequence
23
productivity, neither of which are implemented in UCERF3-ETAS. In addition, the USGS
forecast assumes a GR distribution everywhere. All of these factors likely contribute to
the higher probabilities reported in the USGS forecast. The USGS forecast was updated
four times to reflect changing probabilities due to observed aftershock rates and the
passage of time and was subsequently removed from the USGS event page after the M
7.1 superseded it, but it is still available through the ComCat application programmer
interface. The lower probabilities implied by UCERF3-ETAS are due largely to the
model’s anti-characteristic MFD in this region, which highlights the influence of MFD
shape on conditional triggering probabilities. Clearly this deserves more attention, as
does the representation with respect to epistemic uncertainties.
The 5 July M 7.1 also occurred on a structure that was not explicitly modeled in
UCERF3. It nucleated within the UCERF3 polygon for Airport Lake (Figure 2.2) and, as
such, could be considered an Airport Lake event according to UCERF3 rules. However,
the rupture is not a good match for the surface defined for Airport Lake, and much of the
rupture was not associated with any fault polygons. Thus, according to UCERF3, this
event was partially on and partially off explicitly modeled faults, but any sober
assessment would clearly call it an off-fault event. Although this ambiguity has minimal
practical implications for this particular earthquake (we simply used the actual estimated
rupture surface), it does point out problems with the strict on- versus off-fault dichotomy
in UCERF3, as well as a need for something better than rough polygons with respect to
distinguishing between the two. The deeper question is how to best approximate the
fault system when and where it is poorly known.
24
2.5 Spatial Distribution of M 7.1 Aftershock Probabilities
The spatial distribution of expected aftershocks is a potentially useful product of
operational earthquake forecasts (Field et al., 2016). We ran the UCERF3-ETAS model
with finite-rupture and point-source representations of the Ridgecrest sequence to
compare the predicted spatial distribution of aftershocks with the observed sequence
(obtained from ComCat). Although the analysis presented here is qualitative in nature,
Savran et al. (unpublished manuscript) will provide a quantitative analysis of these
models through new tests developed within the Collaboratory for the Study of
Earthquake Predictability.
Figure 2.4 maps the simulated probability of one or more M ≥ 3.5 aftershocks
occurring in the region within 30 days of the M 7.1 event. That probability is calculated
as the fraction of UCERF3-ETAS simulated catalogs (of which there are 100,000 for
each model variation) that contain at least one such event with an epicenter in each
0.02° by 0.02° cell within 30 days of the mainshock. Although simulated catalogs
include aftershocks throughout California, these maps plot probabilities within the
immediate vicinity of the Ridgecrest sequence. Figure 2.4a shows results from the initial
point-source model, and Figure 2.4b shows results from our preferred model, which
uses the USGS ShakeMap finite-rupture source. The point-source model predicts a
circular distribution of aftershocks, with high-probability areas to the northeast and
southwest of the epicenter that lack observed aftershocks. In contrast, the finite-rupture
source model predicts a more compact distribution, which quickly falls off with distance
25
from the fault plane. This is generally consistent with the observed aftershock
distribution, depicted with circles, which mostly occurred in areas of high probability.
As time from the mainshock passes, point-source models can more accurately
forecast the observed spatial distribution through inclusion of triggering from observed
aftershocks. Figure 2.5 shows results from simulations beginning 7 days after the M 7.1
(including all M ≥ 2.5 aftershocks from ComCat during that period). Here, the point-
source spatial probability distribution predicts highest probabilities in the vicinity of
aftershocks which occurred in the first week, but is still broader than the finite-rupture
surface model, which more closely matches the observed seismicity.
Figures 2.4 and 2.5 show higher probabilities along mapped faults relative to
surrounding areas. These increased probabilities reflect the UCERF3-ETAS hypothesis
that large, supraseismogenic aftershocks are more likely to occur in an area if it
contains a mapped active fault. Because no supraseismogenic aftershocks of the M 7.1
have occurred as of 9 December 2019, spatial probabilities generated with the no-faults
UCERF3-ETAS model appear more consistent with observed seismicity (Figure 2.6).
Advantages of the full UCERF3-ETAS model to forecasting spatial probabilities would
only be realized in rare sequences which trigger supraseismogenic rupture on
previously mapped faults, or when knowledge of triggering probabilities on known faults
is important.
26
2.6 Sensitivity to Finite Fault Surfaces
The probability of triggering nearby faults is of particular interest to the public, as
evident in a 25 July 2019 article in the Los Angeles Times (Lin, 2019). In the article,
USGS scientists Morgan Page and Susan Hough were quoted discussing short-term
probability gains on the Garlock and San Andreas faults, respectively. Both of their
statements were informed by UCERF3-ETAS simulations because this model is
uniquely able to quantify the probability of triggering these and other nearby faults.
However, UCERF3-ETAS had not been used prior to Ridgecrest to address these
questions during a real, supraseismogenic earthquake sequence; the model was
primarily run either on scenario ruptures on UCERF3 faults or smaller earthquakes
modeled as point-sources. As such, the model was never tested for sensitivity to poorly
constrained real-time information on finite-rupture surfaces outside of mapped UCERF3
faults.
The capability to specify custom rupture surfaces not on UCERF3 faults was not
developed until 6 July 2019 (the day after the M 7.1), and the first finite-rupture surface
used was quickly drawn through the observed aftershock sequence. Initial simulations
with that surface, which is shown in Figure 2.7a, gave a 1.7% and 0.21% chance of
triggering M ≥ 7 aftershocks within 30 days of the M 7.1 on the Central Garlock and San
Andreas (South Mojave) fault sections, respectively. After viewing implied surface
rupture from early Interferometric Synthetic Aperture Radar (InSAR) data posted to the
SCEC post-earthquake response coordination page, we revised the southern end of the
rupture surface in our simulations to be nearer to the Garlock fault, which more than
27
doubled those probabilities to 4.4% and 0.43%, respectively. In search of a more
authoritative source (e.g., not hand-drawn by the authors), we subsequently added a
new UCERF3-ETAS configuration capability to download and use the finite fault source
used in ShakeMap. This surface, version 14 (v.14, the most recent version as of 1
October 2019), accessed on 16 July 2019 and first published in ComCat on 11 July
2019, was derived from InSAR data and includes a more detailed, curved surface with
two small strands parallel to the main rupture surface (also shown in Figure 2.7a).
Triggering probabilities calculated with the ShakeMap surface for the Central Garlock
and San Andreas (South Mojave) faults are 3.2% and 0.32%, respectively, which lie
between the probabilities generated using the original planar rupture surfaces. We also
computed probabilities from the first finite ShakeMap surface, v.10, a planar source
which was derived from teleseismic inversion and published in ComCat on 9 July 2019.
That source yields the highest probabilities of 5.1% and 0.45%, respectively. All finite-
source probabilities computed for Ridgecrest are larger than when only point sources
are used, which yields Central Garlock and San Andreas (Mojave South) probabilities of
0.65% and 0.13%, respectively.
Because of these probability differences, the choice of finite-rupture model is
clearly an important source of epistemic uncertainty. The differences in fault trigger
probabilities between finite-rupture models (summarized in Table 2.1) are generally
larger than those between model parameterizations all using the same finite rupture
model (summarized in Table 2.2 and discussed in the next section). This sensitivity to
28
the proximity of a rupture to a mapped UCERF3 fault was anticipated and was explored
during model development, although only with point source ruptures.
We also observed an unexpected sensitivity of fault probabilities to the
complexity of the rupture surfaces themselves. UCERF3-ETAS parameterizes each
rupture surface as an evenly discretized point surface in which each point on the
surface has an equal probability of triggering subsequent events (except in certain
circumstances with the no elastic-rebound-triggering, NoERT, model variant, which
does not apply to this event). As the size of a surface increases, either in overall
geographic extent or through increased complexity (e.g., curves and parallel strands),
the number of points on that surface also increases, which decreases the triggering
probability assigned to each point. This can result in changes to fault trigger
probabilities despite no change in distance between the rupture surface and fault. For
example, if we cull the short parallel strands from the ShakeMap rupture surface (Figure
2.7b) without modifying the end point nearest the Garlock fault, the 30 day Garlock M ≥
7 probability increases from 3.2% to 3.8%. If we instead draw a planar source between
the extents of the ShakeMap surface (also shown in Figure 2.7b), the Garlock M ≥ 7
probability increases further to 5.0%, again holding distance between the rupture and
the Garlock fault fixed. These sensitivities to rupture complexity require careful
consideration, as resulting trigger probabilities would be higher near portions of rupture
surfaces which are more carefully mapped (e.g., due to ease of field access) than areas
of a rupture that are crudely drawn. Again, this sensitivity is due to the fact that each
point on the fault surface is assumed to have equal triggering power. Future model
29
development could relax this assumption so as to negate this effect with a coarse-
graining of local triggering potential or by using other local finite-fault information (e.g.,
slip on different fault patches).
2.7 Sensitivity to Model Parameters
UCERF3-ETAS has a number of model parameters, which are described in detail
in Field, Milner, et al. (2017). We ran alternate simulations of the Ridgecrest sequence,
using the preferred ShakeMap (v.14) finite surface representation and varying only
these parameters from their default values. The results are summarized in Table 2.2.
The largest change in aftershock probabilities we observed came from the UCERF3-
ETAS no-faults model, which imposes a GR MFD everywhere and ignores the presence
of faults. This model gives an 8.0% probability of triggering an M ≥ 7 aftershock
anywhere in the system (compared to 4.6% from the default model).
By default, UCERF3-ETAS imposes elastic rebound to prevent larger
supraseismogenic ruptures from nucleating within the rupture surface of a previous,
supraseismogenic event. An alternative probability model, NoERT, allows immediate re-
rupture. Although this was identified in Field, Jordan, et al. (2017) as the primary
UCERF3-ETAS epistemic uncertainty, it does not have a significant impact here
because the Ridgecrest sequence did not occur on previously mapped UCERF3 faults;
thus the same structure cannot re-rupture. The NoERT model gives a 3.9% chance of
an M ≥ 7 aftershock in the first 30 days on the Garlock fault and a 0.42% chance on the
30
Mojave South section of the San Andreas fault (as opposed to 3.2% and 0.32%,
respectively, with the default elastic-rebound model).
UCERF3-ETAS applies ETAS parameters from Hardebeck (2013) to every
rupture by default, but those parameters can be overridden globally and/or for individual
ruptures (the latter is a new capability added during the Ridgecrest sequence). In the
ETAS model (Ogata, 1988), the instantaneous rate of aftershocks at time 𝑡 due to
earthquakes with magnitudes 𝑀 𝑖 occurring at times 𝑡 𝑖 is given by 𝜆 (𝑡 ) =
∑
𝑖 10
𝑎 +𝛼 (𝑀 𝑖 −𝑀 𝑚𝑖𝑛 )
(𝑡 − 𝑡 𝑖 + 𝑐 )
−𝑝 . The direct Omori parameters 𝑎 , 𝑝 , and 𝑐 can be
estimated from historical earthquake catalogs or a particular sequence. This is typically
done with a maximum-likelihood ETAS parameter solver that assumes aftershock
magnitudes are drawn from a GR distribution. This is, however, incompatible with the
forward model in UCERF3-ETAS, which has deviations from GR magnitude scaling at
small spatial scales. This difference affects more than just aftershock magnitudes in the
UCERF3-ETAS model; it also affects secondary triggering productivity because
aftershocks themselves trigger further aftershocks and that triggering is a function of
aftershock magnitude. In the case of the Ridgecrest sequence, aftershock productivity is
lower in UCERF3-ETAS than would be predicted by a typical ETAS model (or the no-
faults UCERF3-ETAS variant) with the same parameters. We thus estimated sequence-
specific temporal parameters by making manual adjustments to the maximum-likelihood
ETAS parameters (through trial and error, we found that we were able to produce a
reasonable fit to the cumulative number of M ≥ 3.5 events by adjusting only the c-value).
In addition, to match the complex behavior for the Ridgecrest sequence, we estimated a
31
unique productivity (a-value) for primary triggering from the M 7.1 mainshock. The
parameters we estimated (on 30 August, 56 days post-mainshock, using all available
data at the time) are 𝑎 = −2.3 for triggering from the M 7.1 (primary triggering), 𝑎 =
−3.03 for triggering from other earthquakes (secondary triggering) in the sequence, 𝑝 =
1.15, and 𝑐 = 0.002 days. In addition, UCERF3-ETAS assumes 𝛼 = 1, and all
simulations are run with a minimum magnitude of 𝑀 𝑚𝑖𝑛 = 2.5. The high primary
productivity of the mainshock, coupled with low secondary aftershock productivity,
results in more early and fewer late aftershocks than predicted using default
parameters. Although the sequence-specific parameters retrospectively improve the
performance of the UCERF3-ETAS model for the Ridgecrest earthquake, it is worth
noting that, with the exception of the first day, the cumulative number of aftershocks lies
well within the 95% confidence limits of the default model as shown in Figure 2.8. The
lower 𝑎 -value for secondary triggering with these sequence-specific parameters results
in a lower upper extremum and 97.5th percentile, shown in Figure 2.8b, than with
generic parameters, shown in Figure 2.8a, which extends above the chart almost
immediately. This is because the upper extremum is controlled by secondary triggering
from large aftershocks in the most productive sequences (e.g., aftershocks of the M≥8
in the simulation plotted in Figure 2.1c). Calculated 30 day fault probabilities are slightly
higher, 3.7% and 0.35% on the Central Garlock and San Andreas (Mojave South) faults,
respectively. In future work, we may develop an automated maximum-likelihood
parameter solver consistent with the UCERF3-ETAS forward model so that sequence-
specific parameters can be calculated early in a sequence without manual intervention.
32
2.8 Discussion and Conclusions
Development of the UCERF3-ETAS code base during the past year allowed us
to rapidly prepare simulations in response to the Ridgecrest sequence. Although a more
fully automated system is the ultimate goal, on-demand manual submission worked for
this sequence, especially because multiple researchers were trained to run the model
and were able to take over when the primary developer was unavailable. This sequence
also highlighted potential obstacles to deployment of an automated system. For
example, the system might erroneously associate an event like the M 7.1 with the
Airport Lake fault due to its nucleation within the polygon, when it clearly ruptured a
different structure. Another issue is the immediate availability of finite-rupture surfaces.
A sensible option could be to use the finite surface from ShakeMap preferentially,
though such a surface might not be available until days after significant events. In this
case, the first finite ShakeMap source (v.10) for the M 7.1 was posted 4 days after the
event, with the detailed source (v.14) posted after six days. Even if we had a finite
surface available immediately (such as through the finite-fault rupture detector algorithm
described in Böse et al., 2015, and in use by ShakeAlert), or more quickly by fitting early
aftershocks, another algorithm is needed to determine if any UCERF3 fault sections
with polygons that intersect the finite surface should be reset in an elastic-rebound
sense.
Moreover, sensitivities to finite rupture surfaces suggest possible model
improvements. Table 2.1 suggests that even a slight change in rupture geometry can
have a significant effect on triggering possibilities. UCERF3-ETAS could be improved by
33
accounting for uncertainties in input rupture surfaces, similar to how the model spreads
fault probabilities across their respective polygons. The resulting source polygons would
smear out rupture surfaces, potentially also mitigating the effects of the observed
sensitivity to source complexity. Ultimately, we believe that inclusion of finite-fault
sources is important despite these sensitivities. The differences between triggering
probabilities from all of the finite sources are much smaller than the difference between
those calculated using our preferred finite source and a point source. The finite source
differences are, furthermore, small enough to be unlikely to change user actions in
response to the forecast. The Los Angeles Times article (Lin, 2019) highlights a desire
by the public to understand how probabilities on familiar faults might have changed after
a significant earthquake, and even if the numbers provided by UCERF3-ETAS are
uncertain, they are likely better than a naive guess.
The primary lesson learned through this experience is that real events, with real,
uncertain, and evolving data, are vital to better understand the limitations and guide
future development of earthquake forecasting models. Prior to this sequence we had
only run UCERF3-ETAS with ideal scenarios on UCERF3 faults, such as the
“HayWired” scenario, discussed in Field and Milner (2018), or for smaller real events
represented as point sources. This sequence has led us to better understand the model
and its limitations and to develop tools to more quickly and robustly respond to future
large California earthquakes.
34
2.9 Data and Resources
This article uses earthquake data and ShakeMap products from the U.S.
Geological Survey (USGS) Comprehensive Catalog (ComCat) described at
https://earthquake.usgs.gov/data/comcat/. The Uniform California Earthquake Rupture
Forecast, Version 3-epidemic-type aftershock sequence (UCERF3-ETAS) model,
including all scripts used to run the model, are posted at
https://github.com/opensha/ucerf3-etas-launcher. UCERF3-ETAS is built on the
OpenSHA platform, available at https://github.com/opensha. Plots were made with
JFreeChart (http://www.jfree.org/jfreechart/), with perceptually uniform color maps from
Peter Kovesi (https://peterkovesi.com/projects/colourmaps/). Southern California
Earthquake Center’s (SCEC) post-earthquake response page is available at
https://response.scec.org. Simulation results posted on GitHub are available at
https://github.com/kevinmilner/ucerf3-etas-results/tree/master/#etas-simulations-list.
2.10 Acknowledgements
The third Uniform California Earthquake Rupture Forecast (UCERF3)
development was supported by the California Earthquake Authority, U.S. Geological
Survey (USGS), USGS-Southern California Earthquake Center (SCEC) Cooperative
Agreement G17AC00047, and National Science Foundation (NSF)-SCEC Cooperative
Agreement EAR-1600087. Calculations were performed at the University of Southern
California (USC) Center for High-Performance Computing and Communications. The
authors thank the USGS Powell Center for Analysis and Synthesis for supporting
workshops on UCERF3–epidemic-type aftershock sequence (ETAS) development, and
35
the W. M. Keck Foundation for supporting the integration of UCERF3 into the
Collaboratory for Interseismic Modeling and Simulation. Further UCERF3-ETAS
development in support of operationalization was supported by SCEC awards 18142
and 19227. The SCEC Contribution Number for this article is 9937. The authors thank
Andrew Michael and Maximilian Werner who provided valuable feedback on early
UCERF3-ETAS simulation results, Nicholas van der Elst for use of his ETAS parameter
estimation tool, and Jordan Landers, Andrea Llenos, and two anonymous reviewers
whose constructive criticisms improved this article.
36
2.11 Tables
Table 2.1
Aftershock probabilities from UCERF3-ETAS simulations performed with the same
model parameters, varying only the rupture geometries used (sources described in
Sensitivities to Finite-Fault Surfaces). Probabilities are listed for a time span starting
immediately after the M 7.1 event. Distances listed are to the Garlock fault and its
UCERF3 polygon respectively. A distance of 0.0 km to the polygon indicates that the
source extends into the Garlock polygon. Results for the default model are italicized.
Source
Description
30 Day Any M
≥ 7 Probability
Distance to
Garlock Fault,
Polygon (km)
30 Day
Garlock Fault
M ≥ 7
Probability
30 day San
Andreas
(South
Mojave) Fault
M ≥ 7
Probability
Point source 1.8% 31.9, 13.5 0.65% 0.13%
First finite
source
3.1% 11.4, 0.0 1.7% 0.21%
Second finite
source
5.6% 5.7, 0.0 4.4% 0.43%
Previous
ShakeMap
(v.10) source
6.3% 3.0, 0.0 5.1% 0.45%
ShakeMap
(v.14) source
4.6% 5.6, 0.0 3.2% 0.32%
ShakeMap
(v.14) culled
source
5.1% 5.6, 0.0 3.8% 0.36%
ShakeMap
(v.14) source
extents
6.3% 5.6, 0.0 5.0% 0.43%
37
Table 2.2
Aftershock probabilities from UCERF3-ETAS simulations performed with the same input
rupture geometries (preferred model, ShakeMap source), varying only model
parameters. Results with the default model parameters are given in the first row and
italicized. Variations in which only the listed parameters are changed from default, are
listed below (parameters are described in the Sensitivity to Model Parameters).
Probabilities are listed for a time span starting immediately after the M 7.1 event.
Parameter
Description
30 Day Any M ≥ 7
Probability
30 Day Garlock
Fault M ≥ 7
probability
30 day San
Andreas (South
Mojave) Fault M ≥
7 probability
Default parameters 4.6% 3.2% 0.32%
No-faults model 8.0% n/a n/a
NoERT probability
model
5.3% 3.9% 0.41%
Sequence-specific
parameters
5.1% 3.7% 0.35%
38
2.12 Figures
(a) (b)
(c) (d)
39
Figure 2.1
Percentile scenarios: Map view of UCERF3-ETAS synthetic catalogs of aftershocks to
the Ridgecrest sequence (preferred model). Shown are 30 days of simulated
seismicity immediately following the M 7.1 event. UCERF3 faults in the map region
are plotted with gray outlines, except for faults which participate in events in the
synthetic catalogs, which are plotted in the same color as the event hypocenter. (a) A
typical catalog, defined as the catalog that lies at the 50th percentile with respect to
the total number of simulated aftershocks across all 100,000 simulations. That catalog
has 2584 M ≥ 2.5 events and multiple M ≥ 5 events, but no supraseismogenic events
on UCERF3 faults. (b) A less probable catalog (97.5th percentile, 7574 M ≥ 2.5
events) in which the Garlock fault was triggered in an M ≥ 7 supraseismogenic
aftershock. (c) The most extreme catalog from the suite of 100,000 simulations with
39,550 M ≥ 2.5 events in which an M ≥ 7 Garlock aftershock triggered an M ≥ 8
earthquake on the San Andreas fault. (d) The observed earthquakes with M ≥ 2.5
within 30 days of the M 7.1, accessed from Comprehensive Catalog (ComCat) on 6
December 2019.
40
Figure 2.2
Map view of model input events from ComCat (M ≥ 2.5, accessed 4 September 2019)
from one week preceding to one week following the M 7.1 earthquake (circles).
Events that occurred up to and including the M 7.1 are plotted as dark circles, and
those that occurred in the week following are light circles. UCERF3 faults traces are
drawn as thick black lines, and the extent of the dipping Airport Lake fault surface is
outlined in dotted gray lines. UCERF3 polygons for each fault are drawn with thin
dashed dark lines. The M 7.1 epicenter (largest dark circle) occurred inside the
polygon for the Airport Lake fault, while the M 6.4 (second largest dark circle)
occurred outside of all UCERF3 fault polygons.
41
(a) (b)
Figure 2.3
UCERF3-ETAS forecasted magnitude-number distribution for one week following the
4 July M 6.4 earthquake (point source). (a) Results from the default UCERF-ETAS
model, which is anti-characteristic in this region, and (b) the no-faults model which
follows a Gutenberg-Richter (GR) distribution to large magnitudes. The mean
expectation (computed from 100,000 simulations) is plotted with a thick black line and
95% confidence bounds with thin black lines. The probability of one or more
aftershocks at or above each magnitude is given with a solid gray line, and the
probability of such an event on a supraseismogenic UCERF3 fault rupture is with a
dashed gray line (not applicable to the no-faults model). The light gray shaded region
at high magnitudes represents sampling uncertainties on probabilities calculated
assuming a binomial distribution.
42
(a) (b)
Figure 2.4
UCERF3-ETAS forecasted spatial probability distribution of M ≥ 3.5 aftershocks for 30
days immediately following the M 7.1 earthquake. Each 0.02° x 0.02° cell is colored
on a logarithmic scale according to the probability of at least one such aftershock
occurring within that cell, with darker colors representing higher probability. Actual
epicenters which occurred during this forecast period (obtained from ComCat on 24
September 2019) are overlaid as circles. (a) Results with point source representations
of all input events, and (b) results for the preferred model, with a finite-rupture surface
obtained from U.S. Geological Survey (USGS) ShakeMap. UCERF3 faults are shown
with gray outlines.
43
(a) (b)
Figure 2.5
Same as Figure 2.4, but for 30 day probabilities beginning seven days after the M 7.1
earthquake.
44
(a) (b)
Figure 2.6
Same as Figures 2.4b and 2.5b but with the no-faults UCERF3-ETAS model. (a)
Probabilities for 30 days immediately following the M 7.1 earthquake and (b) for 30
days starting seven days after the M 7.1 earthquake.
45
(a) (b)
Figure 2.7
Finite rupture surfaces used as input to UCERF3-ETAS simulations. UCERF3
mapped fault traces are drawn as thin black lines (with the Garlock, Airport Lake, and
Little Lake faults annotated), and observed M ≥ 2 seismicity that occurred between the
5 July 2019 M 7.1 and 28 October 2019 are plotted as small light gray circles
(accessed from ComCat on 28 October 2019). (a) The first hand-drawn finite rupture
surface as a thick dark gray dashed line, the second finite rupture surface, which dips
to the southwest, as a thick gray solid line, and the preferred ShakeMap surface as
thick black lines. (b) The primary ShakeMap strand as thick black lines, parallel
strands (which were removed for the “Culled” model in sensitivity tests) as thick gray
lines, and the surface extents as a thick dashed dark gray line.
46
(a) (b)
Figure 2.8
Cumulative number of M ≥ 3.5 events over time for the region mapped in Figures 2.4-
2.6, starting immediately following the M 7.1 event. The thick gray line shows catalog
data from ComCat, accessed on 28 October 2019. The thick black line shows the
median from UCERF3-ETAS. Extrema from 100,000 simulations are plotted with thin
dotted black lines, and 2.5%, 16%, 84%, and 97.5% percentiles are plotted with thin
black lines. (a) UCERF3-ETAS results for the default model, and (b) with the same
inputs except sequence-specific parameters. The time at which the parameters for
panel (b) were estimated is annotated with a dashed black vertical line.
47
3. Short-Term Forecasting of Earthquake Ground Motions Using
Conditional Hypocenter Distributions
Kevin R. Milner, Thomas H. Jordan
3.1 Abstract
Rupture directivity is typically absent from long-term probabilistic seismic hazard
estimates which assume a uniform conditional hypocenter distribution (CHD), but
should be an important component of future operational earthquake forecasting (OEF)
efforts. OEF aims to provide authoritative and actionable information about the time
dependence of seismic hazards to help communities prepare for earthquakes, and such
hazards are best expressed in terms of earthquake ground motion probabilities. In
2017, the Working Group on California Earthquake Probabilities released the first
comprehensive fault-based epidemic-type aftershock sequence (ETAS) model, an
extension of the Third Uniform California Earthquake Rupture Forecast (UCERF3) long-
term model. A key advantage of UCERF3-ETAS is the opportunity to more accurately
forecast ground motions of triggered aftershocks over prior point-process models. One
such advantage relates to rupture directivity. Although the long-term UCERF3 models
assume uniform CHDs along-strike of forecasted ruptures, UCERF3-ETAS forecasts
short-term increases in localized hypocenter probabilities in the vicinity of recent
seismicity, which can have a significant effect on expected ground motions when
properly accounting for rupture directivity.
48
We explore the effect of the CHD in short-term ground motion forecasts following
two magnitude 6 scenario earthquakes in Southern California. Calculations are
performed with both empirical estimates of ground motion intensity and directivity, as
well as 3D deterministic ground motion simulations with the CyberShake platform. We
find significant gains from the inclusion of the short-term CHD, on the order of a factor of
ten with respect to short-term ground motion models which ignore it. This suggests that
the CHD, and resultant ground motions which include rupture directivity, should be
included in future OEF models.
49
3.2 Introduction
The directivity of seismic radiation depends primarily on the location of the
hypocenter relative to the main area of fault slip. In fault-based parameterizations of
empirical ground motion models (GMMs), the shaking intensities depend on the location
of the site relative to the fault source, but usually not on the hypocenter. Consequently,
the mean shaking intensities determined by empirical GMMs include only mean
directivity effects; i.e., averages taken over the distributions of hypocenters and space-
time slip fields sampled by the GMM datasets. The empirical sampling of hypocenters
and slip fields thus underlies the GMM regressions, though the statistics of the samples
are usually not explicitly characterized (Spudich et al., 2014).
There is considerable interest in improving the GMMs by accounting more
explicitly for directivity effects (Bozorgnia et al., 2014). The empirical approach is to
extend the parameterizations to hypocenter-dependent models by defining directivity
relative to a hypocenter-independent model (Spudich & Chiou, 2008) or, better yet,
through a joint regression involving all parameters (Spudich et al., 2014).
In this paper, we take a simulation-based approach, assessing directivity effects
using intensities measured from three-dimensional (3D) wavefield simulations computed
on the CyberShake platform (Graves et al. 2011), and we compare the results with an
empirical directivity model. CyberShake calculates large suites of synthetic
50
seismograms using fault ruptures sampled from the Second Uniform California
Earthquake Rupture Forecast (UCERF2; Field et al., 2009). The UCERF2 rupture
probabilities are augmented with conditional probability distributions that describe the
hypocenter locations and space-time slip fields given a pre-chosen fault area. Following
Wang & Jordan (2014), we call these two model components the “conditional
hypocenter distribution” (CHD) and the “conditional slip distribution” (CSD).
Previous CyberShake calculations have assumed CHDs that are spatially
uniform and CSDs given by the Graves-Pitarka stochastic slip model (Graves & Pitarka,
2010, 2014). Post-processing of the CyberShake intensity data using averaging-based
factorization has quantified the mean and variance of the directivity effects and their
dependence on the CSD (Wang & Jordan, 2014). These calculations have shown that
the overly strong directivity effects predicted by simple Haskell-type source models
(Haskell, 1969) are considerably reduced in slip models with more realistic space-time
complexity (Graves & Pitarka, 2014; Wang & Jordan, 2014). In particular, the Graves &
Pitarka (2014) stochastic model used here predicts directivity effects that are in general
agreement with the observations (Jordan et al., 2018).
Our focus in this paper is to extend the time-independent CyberShake forecasts
to time-dependent models in which the CHD conforms to the short-term probabilities of
the Third Uniform California Earthquake Rupture Forecast-Epidemic-Type Aftershock
Sequence (UCERF3-ETAS; Field, Milner, et al., 2017; Field, Jordan, et al., 2017)
51
model. This time-dependent CHD is conditioned on the current state of seismicity using
a fault-based version of the ETAS model (Ogata, 1988). We investigate the ground-
motion probabilities during aftershock sequences of M 6 scenario events near the San
Andreas fault and show that the directivity effects in the fault-based model can affect the
probabilities of strong shaking by an order of magnitude relative to uniform CHD
models. We conclude that the combination of UCERF3-ETAS and CyberShake could be
a useful tool for operational earthquake forecasting (OEF).
3.3 Scenario Definition
We compute short-term hazard from aftershocks of two scenarios, both M 6
point-source earthquakes. One scenario is for Bombay Beach near the Coachella
section of the San Andreas Fault (SAF). During the 2009 and 2016 swarms in this
region, it was thought that not only has the probability of a large rupture on the SAF
Coachella section increased, but the probability of such a rupture with a southernmost
hypocenter (near the swarms) has increased. This is the type of northward-propagating
SAF rupture identified in Olsen et al. (2006) as likely to cause large ground motions in
Los Angeles due to waveguide and directivity-basin coupling effects. We use the
hypocenter location of the largest earthquake in the 2009 swarm (latitude N 33.3172
degrees, longitude W 115.7280 degrees, depth 5.96 km). The 2009 earthquake was an
M 4.8, however we increase the magnitude in our simulations in order to highlight the
gains possible in higher probability environments when the public is more likely to take
preventative actions based on short-term forecasts.
52
The second scenario is for a M 6 near the Mojave section of the SAF, the closest
section of the SAF to Los Angeles. SAF ruptures originating due North of Los Angeles
on the Mojave section would either rupture Northwest toward Bakersfield, Southeast
toward the Coachella Valley, or bilaterally in both directions. All three cases are typically
predicted to produce lower ground motions in Los Angeles than ruptures of the same
fault which propagated toward Mojave. We place this scenario in the middle of the
Mojave section (latitude N 34.42295 degrees, longitude W 117.80177 degrees, depth
5.8 km). Locations for both scenarios are mapped in Figure 3.1.
3.4 Methodology
We ran the UCERF3-ETAS model to compute 500,000 stochastic realizations of
ten years of seismicity following each scenario on the UCERF2 (Field et al., 2009) fault
system. We use the prior UCERF2 fault model in order to also utilize the results of the
Southern California Earthquake Center (SCEC)’s CyberShake three-dimensional (3D)
deterministic physics-based probabilistic seismic hazard analysis platform (Graves et
al., 2011; Jordan et al., 2018), which does not yet support complex multi-fault ruptures
included in UCERF3 and has already been computed for UCERF2 in the Los Angeles
region. Each UCERF3-ETAS synthetic catalog consists of a list of simulated
earthquakes with origin time, hypocenter, magnitude, and, for on-fault supra-
seismogenic ruptures, finite-rupture extent. We consider only these on-fault supra-
seismogenic ruptures to compute hazard for this study, as only these types of ruptures
are presently supported by CyberShake.
53
Figures 3.2 and 3.3 show the time-dependent increase in probability of large M ≥
7 aftershocks within 100 km of the M 6 Bombay Beach and Mojave scenarios
respectively. Gains relative to long-term time-dependent forecasts are highest for the
first hour immediately following the scenario events, a factor of 426 for Bombay Beach
and 769 for Mojave. The higher gain for Mojave is due to a higher background rate of
events on that section in UCERF3 than the Coachella section near Bombay Beach, and
the placement of the scenario event directly on the SAF Mojave as opposed to 4 km
away from SAF Coachella for Bombay Beach. While gains are more dramatic for the
shortest time horizons (e.g., the first hour), we focus on one week following the scenario
event as gains are still high (>10 for each scenario), and the longer duration allows
enough time for both the scientific community to generate a forecast and the user
community to react. UCERF3-ETAS simulations were performed at the Center for High-
Performance Computing of the University of Southern California, requiring
approximately 25,000 CPU-hours for each scenario (500,000 simulations each). Figure
3.4 shows the short term aftershock CHD on the SAF Coachella section following the
Bombay Beach scenario, with most hypocenters at the Southeast end of the fault, and
Figure 3.5 shows the CHD on the SAF Mojave following the Mojave scenario, with most
hypocenters near the scenario event.
Empirical ground motion predictions from the UCERF3-ETAS simulations are
made with the Abrahamson et al. (2014; henceforth ASK2014) ground motion model,
modified to account for rupture directivity with the Bayless-Sommerville directivity model
(Bayless, 2013). Hazard calculations are done with the OpenSHA platform (Field et al.
54
2003). ASK2014 uses depth to 𝑉 𝑠 = 1.0 km/s (Z1.0) and VS30 as proxies for 3D basin
effects and site amplification. We set VS30 from the Wills et al. (2015) geologically-based
map, and Z1.0 from the tomographically inverted SCEC Community Velocity Model,
version 4.26-M01 (CVM-S4.26-M01; Lee et al. 2014; Small et al., 2017).
Additionally, we compute 3D deterministic ground motions with the SCEC
CyberShake platform. CyberShake Study 15.4 (CS15.4) extends the UCERF2
earthquake rupture forecast with large ensembles of kinematic ruptures generated with
the Graves and Pitarka (2014) method, sampling many different hypocenters and slip
distributions for each UCERF2 rupture. These ensembles are large enough (>400,000
total kinematic rupture variations for all 14,400 UCERF2 on-fault ruptures) to sample the
statistical variability of ground motions, as demonstrated by Wang and Jordan (2014).
CyberShake utilizes seismic reciprocity to efficiently compute long-period (𝑇 = [2,10] s)
synthetic seismograms for all 400,000 kinematic rupture variations of the UCERF2
model through the CVM-S4.26-M01 3D velocity model for 336 sites in the Los Angeles
region. Even with seismic reciprocity, this calculation requires significant computational
resources. CS15.4 used over 30 million CPU-hours and generated over a petabyte of
data on two of the nation's largest supercomputers at the time: Titan at Oak Ridge
National Laboratories and Blue Waters at the National Center for Supercomputing
Applications.
A key point is that while CyberShake’s computational requirements are large, all
work described here uses previously computed CyberShake RotD50 spectra. This study
55
modifies the probabilities of each UCERF2 rupture variation in order to compute time-
dependent hazard curves, which are then computed rapidly from the CS15.4 database
of pre-computed ground motions. For each UCERF3-ETAS aftershock (in each of the
500,000 synthetic catalogs for each scenario), we identify the corresponding UCERF2
rupture for which CyberShake ground motions have already been computed. Those
UCERF2 ruptures are represented by hundreds to thousands of kinematic rupture
variations in CyberShake, each with a different hypocenter and slip distribution. We
select just those rupture variations with similar hypocenters (within 10km) to the original
UCERF3-ETAS hypocenter. This results in an average of 6.8 mapped kinematic rupture
variations for each UCERF3-ETAS aftershock with the Bombay Beach M 6 scenario,
and 8.4 variations for each Mojave M 6 aftershock. This mapping is demonstrated in
Figures 3.4 and 3.5, which shows the original UCERF3-ETAS hypocenters and the
mapped CyberShake rupture variation hypocenters.
3.5 Results
Probabilistic seismic hazard curves of 5 s RotD50 spectral acceleration for one
week following each scenario show large gains relative to long-term hazard estimates.
This is especially true in the thick sedimentary basin CyberShake site STNI (location
mapped in Figure 3.1), which is located at the deepest part of the Los Angeles basin
according to CVM-S4.26-M01. UCERF3-ETAS hazard curves computed at STNI with
the ASK2014 ground motion model and Bayless-Sommerville directivity model are
shown as solid red lines in Figure 3.6. We find a maximum gain factor of 13.8 (at 0.63 g)
56
relative to the long term model for the Bombay Beach M 6 scenario, and 37.4 (at 0.05 g)
for the Mojave M 6 scenario.
In order to isolate the effects of the directivity model, we also compute UCERF3-
ETAS curves with only the ASK2014 model (no directivity modifier), shown as dashed
red lines in Figure 3.6. For the Bombay Beach scenario, consideration of the CHD and
resultant directivity effects increases hazard at all intensity levels, and contributes most
to the probability of exceeding larger intensities (>0.2 g) with a gain factor of 4.8 relative
to ETAS curve without directivity. The opposite is true for the Mojave case, where
inclusion of directivity lowers the hazard estimate, an effect which is also greatest at
large (>0.2 g) intensities. The primary difference between these cases is that the
Bombay Beach scenario tends to trigger ruptures with southeastern hypocenters that
propagate toward Los Angeles, and the Mojave scenario triggers ruptures which
propagate away northwest or southeast. Even though the Mojave scenario is both
closer to and more productive (in terms of the rate that of large aftershocks are
triggered) than Bombay Beach, its probability of exceeding intensities larger than 0.5 g
is lower when considering directivity; Mojave exceedance probabilities are always larger
without the directivity model. This is illustrated spatially in Figure 3.7, which plots gain
from inclusion of directivity for a probability of exceeding 0.2 g. The effect of inclusion of
the CHD and directivity for the Bombay Beach scenario is a relatively spatially uniform
increase in hazard near Los Angeles (Figure 3.7b), and for the Mojave scenario, a
pattern of decreased hazard extending fault-normal from SAF Mojave, with increased
fault-parallel hazard (Figure 3.7d).
57
Gains are more dramatic with the CS15.4 3D deterministic ground motion
simulation method, especially in the Los Angeles Basin. CyberShake hazard curves for
STNI are shown in Figure 3.8. We find a maximum gain factor of 69.4 (at 0.38 g)
relative to the long term model for the Bombay Beach M 6 scenario, and 51.5 (at 0.07 g)
for the Mojave M 6 scenario. We isolate directivity effects in CyberShake by also
computing curves with a uniform CHD, instead of only selecting those rupture variations
with hypocenters nearest to the ETAS aftershock (red dashed curves in Figure 3.8). The
largest gain from consideration of directivity at site STNI is 11.3 (at 0.38 g) for Bombay
Beach, as opposed to a gain of 4.8 with the empirical model. This is consistent with the
findings of Graves et al. (2011), which noted stronger directivity effects in a prior
CyberShake model in the Los Angeles basin than predicted by empirical models. The
downward effect of inclusion of directivity (relative to the no-directivity CyberShake
ETAS model) at site STNI is also stronger for the Mojave scenario than with the
empirical model. Gain maps comparing the two CyberShake ETAS models (Figure 3.9)
show a wider and less uniform pattern of decreased hazard for Mojave M 6, and a
variable pattern with largest gains in the sedimentary basins for Bombay Beach M 6.
Unlike the empirical model, gains of 0.2 g 5 s RotD50 spectral acceleration for
the Bombay Beach M 6 scenario with CyberShake are not always greater than one.
This is apparent in the Inland Empire region of Los Angeles County near Rancho
Cucamonga at the site PDU. CyberShake and empirical hazard curves for site PDU,
whose location is mapped in Figure 3.1, are shown in Figure 3.10. Below 0.2 g, the
58
CyberShake ETAS probabilities are larger when including directivity, but at 0.2 g they
fall off rapidly to below the no-directivity CyberShake ETAS model. The empirical
hazard curve for PDU which includes directivity remains above the model without
directivity at all intensity levels. We believe this falloff at 0.2 g with CyberShake to be a
sampling issue, introduced by the reduced counts of rupture variations when
considering only specific hypocenters relative to the complete CS15.4 model. It is likely
that with larger rupture ensembles the CyberShake model with directivity would remain
above the no-directivity model, as is the case with empirical models. CyberShake
represents ground motion variability through large ensembles of hypocenters and slip
distributions; although this may be adequate for long term models which weight each of
those rupture variations equally, more variations of rupture slip are likely needed for
robust short term hazard estimation.
3.6 Discussion
Although consideration of rupture directivity is not yet standard practice in
seismic hazard assessments, it can and should play an important role in future OEF
models (including the operationalization of the UCERF3-ETAS model). This is
especially important for deep sedimentary basins such as Los Angeles, as
demonstrated by the larger gains due to inclusion of directivity with CyberShake 3D
deterministic ground motion simulations than with the empirical estimates of Bayless-
Sommerville. In the case of the Bombay Beach and Mojave M 6 scenarios, inclusion of
directivity results in the probability of large (e.g. >0.5 g RotD50 5 s SA) ground motions
in Los Angeles from Bombay Beach to surpass Mojave, even though the latter produced
59
more productive aftershock sequences in the ETAS model and is nearer to Los
Angeles. Jordan et al. (2018) noted that as slip complexity in kinematic rupture
generators increases, the amplitude of directivity effects decreases through reduced
constructive interference of resultant waves. They also found that the increased
complexity in the Graves and Pitarka (2014) model used in this study, relative to their
2010 model, is primarily at scales smaller than the wavelength of the 5 s period ground
motions we consider. This suggests that future refinements to the rupture generation
method are unlikely to strongly affect the conclusions of the study.
We believe that increased attention should be given to modeling of directivity in
future empirical ground motion model development efforts (e.g. the Next Generation
Attenuation project), rather than an add-on model approach as in NGA-West2.
Deterministic models such as CyberShake can inform this model development, and
even be used directly in short-term long-period hazard assessment, so long as the
ensemble of rupture variations is large enough to adequately sample expected ground
motion variability with a fixed hypocenter location.
This highlights a unique advantage of the UCERF3-ETAS model, relative to
simpler OEF candidate models (e.g. STEP) which ignore finite fault effects. Significant
gains in ground motion probabilities are possible relative to a STEP-type model, due not
only to the inclusion of finite fault ruptures, but also forecasting rupture directivity. With
the addition of CyberShake, even further probability gains are possible through
directivity-basin coupling effects which are not accounted for in empirical models.
60
3.7 Data and Resources
The Uniform California Earthquake Rupture Forecast, Version 3-epidemic-type
aftershock sequence (UCERF3-ETAS) model, including all scripts used to run the
model, are posted at https://github.com/opensha/ucerf3-etas-launcher. UCERF3-ETAS
is built on the OpenSHA platform, available at https://github.com/opensha. Plots were
made with JFreeChart (http://www.jfree.org/jfreechart/) and the Generic Mapping Tools
(https://www.generic-mapping-tools.org/). More information on CyberShake Study 15.4
is available at https://scec.usc.edu/scecpedia/CyberShake_Study_15.4.
3.8 Acknowledgements
The third Uniform California Earthquake Rupture Forecast (UCERF3)
development was supported by the California Earthquake Authority, U.S. Geological
Survey (USGS), USGS-Southern California Earthquake Center (SCEC) Cooperative
Agreement G17AC00047, and National Science Foundation (NSF)-SCEC Cooperative
Agreement EAR-1600087. UCERF3-epidemic-type aftershock sequence (ETAS)
calculations were performed at the University of Southern California (USC) Center for
High-Performance Computing and Communications. The authors thank the USGS
Powell Center for Analysis and Synthesis for supporting workshops on UCERF3-ETAS
development, and the W. M. Keck Foundation for supporting the integration of UCERF3
into the Collaboratory for Interseismic Modeling and Simulation, which allowed for
CyberShake integration. Further UCERF3-ETAS development was supported by SCEC
awards 18142 and 19227. The authors thank Edward Field for his leadership of
UCERF3-ETAS development.
61
3.9 Figures
Figure 3.1
Map showing the geographic location of the two scenario ruptures (stars) and the
STNI and PDU CyberShake sites (triangles). Cities of Los Angeles and Palm Springs
are annotated with circles.
62
Figure 3.2
Time-dependent probabilities of M ≥ 7 aftershocks within 100 km of a M 6 scenario
earthquake near Bombay Beach, as a function of the forecast horizon. The thick red
curve gives the Third Uniform California Earthquake Rupture Forecast-Epidemic-Type
Aftershock Sequence (UCERF3-ETAS) probabilities from 500,000 simulations, the
blue curve long-term time-dependent probabilities from UCERF3-Time-Dependent
(UCERF3-TD), and the black curve time-independent probabilities from UCERF3-
Time-Independent (UCERF3-TI). UCERF-ETAS aftershock probabilities plotted are
the sum of raw UCERF3-ETAS simulation rates and UCERF3-TD rates, probabilities
of aftershocks from UCERF3-ETAS only are shown in with the thin red line. The gain
for the first hour following the scenario is 426 relative to UCERF3-TD. Gains for one
day and one week are 50 and 11 respectively. As forecast horizon increases, the
background probability from UCERF3-TI and UCERF3-TD increase, causing ETAS
forecasted gains to decrease.
63
Figure 3.3
Same as Figure 3.2, except for aftershocks of a Mojave M 6 scenario earthquake. The
gain for the first hour following the scenario is 769 relative to UCERF3-TD. Gains for
one day and one week are 129 and 27 respectively.
64
65
Figure 3.4
Conditional hypocenter distribution (CHD) for one week following the Bombay Beach
M 6 scenario on the Coachella section of the San Andreas fault. The top panel shows
histograms of the probability of a hypocenter at each point along-strike of Coachella,
conditioned on the occurrence of a rupture which initiates on Coachella in the
UCERF3-ETAS simulated aftershock catalogs. The distribution of original UCERF3-
ETAS hypocenters is shown in black, those hypocenters remapped to the nearby
(within 10 km) precomputed CyberShake hypocenters in blue, and the complete set of
hypocenters for all CyberShake rupture variations in red. The bottom panel shows the
distribution of hypocenters along-strike and down-dip of the Coachella section, with
the same colors. The scenario rupture hypocenter is annotated with a green triangle.
66
Figure 3.5
Same as Figure 3.4, but with the CHD for the Mojave M 6 scenario plotted along-strike
and down-dip of the Mojave (South) section of the San Andreas fault.
67
(a) (b)
Figure 3.6
One-week hazard curves with the empirical ASK 2014 ground motion model for site
STNI following the (a) Bombay Beach M 6 and (b) Mojave M 6 scenarios. Long-term
time-independent and time-dependent hazard curves are shown in solid blue and
green respectively. UCERF3-ETAS hazard curves are shown in red. Calculations
made with UCERF3-ETAS hypocenters and the Bayless-Sommerville directivity
model are shown in solid red, and those with only the ASK 2014 ground motion model
(without considering rupture hypocenter or directivity) in dashed red.
68
(a) (b)
(c) (d)
69
Figure 3.7
Maps of one-week probability gain 0.2 g 5 s RotD50 spectral acceleration following
the Bombay Beach M 6 (a, b) and Mojave M 6 (c, d) scenarios. Gain from the full
time-dependent model including directivity relative to the long-term time-dependent
model are shown in (a) and (c). The component of that gain from the inclusion of the
directivity model only is shown in (b) and (d). Sites between which hazard was
interpolated are annotated as black triangles, and the Mojave M 6 scenario
hypocenter is annotated with a red star in (c) and (d).
(a) (b)
Figure 3.8
One-week hazard curves with deterministic CyberShake ground motions for site STNI
following the (a) Bombay Beach M 6 and (b) Mojave M 6 scenarios. Long-term time-
independent and time-dependent hazard curves are shown in solid blue and green
respectively. UCERF3-ETAS hazard curves are shown in red. Calculations made with
the UCERF3-ETAS CHD (mapped to CyberShake rupture variations within 10km of
corresponding UCERF3-ETAS hypocenters) solid red, and with a uniform CHD in
dashed red.
70
(a) (b)
(c) (d)
71
Figure 3.9
Maps of one-week probability gain 0.2 g 5 s RotD50 spectral acceleration following
the Bombay Beach M 6 (a, b) and Mojave M 6 (c, d) scenarios. Gain from the full time-
dependent model including directivity relative to the long-term time-dependent model
are shown in (a) and (c). The component of that gain from the inclusion of the
UCERF3-ETAS CHD is shown in (b) and (d). Sites between which hazard was
interpolated are annotated as black triangles, and the Mojave M 6 scenario
hypocenter is annotated with a red star in (c) and (d).
(a) (b)
Figure 3.10
One-week hazard curves with the (a) CyberShake model and (b) empirical ASK 2014
ground motion model for site PDU following the Mojave M 6 scenario. Long-term time-
independent and time-dependent hazard curves are shown in solid blue and green
respectively. UCERF3-ETAS hazard curves are shown in red. Calculations made with
UCERF3-ETAS CHD, which include directivity are shown in solid red, and those with
a uniform CHD in dashed red.
72
4. A Prototype Probabilistic Seismic Hazard Model for California
Constructed with Fully-Deterministic Physical models
Kevin R. Milner, Bruce E. Shaw, Christine A. Goulet, Keith B. Richards-Dinger, Scott
Callaghan, Thomas H. Jordan, Jim H. Dieterich
4.1 Abstract
We investigate the efficacy of a multi-cycle deterministic earthquake simulator as
an extended earthquake rupture forecast (ERF) for use in generating simulated ground
motions for probabilistic seismic hazard analyses (PSHA). While use of deterministic
ground motion simulations in PSHA calculations is not new (e.g. CyberShake), prior
studies relied on kinematic rupture generators to extend empirical ERFs. Fully-dynamic
models, which simulate rupture nucleation and propagation of static and dynamic
stresses, are still computationally intractable for the large simulation domains and many
seismic cycles required to perform PSHA. Instead, we employ the Rate-State
Earthquake Simulator (RSQSim) to efficiently simulate hundreds of thousands of years
of M > 6.5 earthquake sequences on the California fault system. RSQSim produces full
slip-time histories for each rupture, which, unlike kinematic models, emerge from
frictional properties, fault geometry, and stress transfer; all intrinsic variability is
deterministic. We use these slip-time histories directly as input to wave propagation
codes with the Southern California Earthquake Center (SCEC) BroadBand Platform for
one-dimensional models of the Earth and SCEC CyberShake for three-dimensional
models to obtain simulated deterministic ground motions.
73
Resultant ground motions match empirical ground motion model (GMM)
estimates of median and variability of shaking well. When computed over a range of
sources and sites, the variability is similar to that of ergodic GMMs. Variability is
reduced for individual pairs of sources and sites, which repeatedly sample a single path,
which is expected for a non-ergodic model. This results in increased exceedance
probabilities for certain characteristic ground motions for a source-site pair, while
decreasing probabilities at the extreme tails of the ergodic GMM predictions. We
present these comparisons and preliminary fully deterministic physics-based RSQSim-
CyberShake hazard curves, as well as a new technique for estimating within- and
between-event variability through simulation.
4.2 Introduction
Probabilistic seismic hazard analysis (PSHA), first formalized by Cornell (1968),
is typically performed by combining an earthquake rupture forecast (ERF) with a set of
empirical ground motion models (GMMs). ERFs have traditionally relied on observed
fault slip rates, scaling relationships, and regional magnitude-frequency distributions to
estimate the rate of large earthquakes on pre-defined fault segments (e.g., Field et al.,
2014). Empirical GMMs are developed through regression of recorded ground motions
using equations that are consistent with basic physical phenomena of source scaling,
wave propagation and site effects. Even as recorded ground motions databases are
expanding, there is still a lack of data for short site-rupture distances and for large,
complex ruptures that often drive the risk to human infrastructure. This lack of data,
along with the ergodic assumption which applies global ground motion statistics from
74
recordings across many sites and ruptures to each individual site and rupture pair, leads
to large aleatory variability terms in GMMs which Anderson and Brune (1999) noted
controls hazard at large intensities (Figure 4.1), and are by definition irreducible within
an ergodic modeling framework. We aim to reduce this uncertainty with physical models
of rupture occurrence, propagation, and resultant ground motions. Although current
physical models rely on many assumptions and imperfect seismic velocity models, their
associated uncertainties are epistemic, and can therefore be reduced as knowledge is
gained.
We present results from one set of physical models. First, we compare total
variance to an ergodic GMM, and then examine the variance structure through a new
technique developed to estimate within- and between-event variability through
simulation. Then, we compute hazard curves for sites in the Los Angeles region. We
conclude that a non-ergodic physics-based PSHA model can reduce design-level
ground motions through quantifying repeated path, source, and site effects.
4.3 Source Model: RSQSim
We replace traditional ERFs with a multi-cycle physics-based earthquake
simulator, the Rate-State Earthquake Simulator (RSQSim), developed by Dieterich &
Richards-Dinger (2010). Shaw et al. (2018) found that RSQSim simulations on the Third
Uniform California Earthquake Rupture Forecast (UCERF3; Field et al., 2014) fault
system produce seismicity catalogs that match long term rates on major faults and
agrees with UCERF3 when carried through to empirical GMM-based PSHA
75
calculations. RSQSim is very computationally efficient due to numerous analytic
approximations and a boundary-element event-driven three-state algorithm, which
allows for the generation of long synthetic catalogs (hundreds of thousands of years) on
the complex California fault system. Ruptures nucleate and propage in RSQSim through
rate- and state-dependent friction, first proposed by Dieterich (1992). In states 0
(healing) and 1 (nucleating slip), RSQSim employs a quasi-static approximation, where
the shear stresses applied to each fault patch are balanced by the frictional shear stress
(Richards-Dinger and Dieterich, 2012). During earthquake slip (state 2), RSQSim uses
a first order quasi-dynamic approximation with a stepwise constant sliding speed and
dynamic overshoot. Unlike traditional ERFs, RSQSim produces full slip-time histories for
all simulated ruptures that can be used directly as input to deterministic wave
propagation simulations. Cumulative slip from an M 7.45 rupture on the Mojave section
of the San Andreas fault is shown in Figure 4.2, along with the rest of the fault model for
Southern California. Rupture stress drops, hypocenters, and roughness are fully
deterministic and dependent only on global frictional parameters and the state of stress
at nucleation. Figure 4.3 shows a side view of the rupture from Figure 4.2, highlighting
the roughness in total displacement present in RSQSim ruptures. Other notable
features of this rupture include tapering of slip at the bottom of the seismogenic zone,
as well as back propagation apparent through late time-of-last-slip values for patches
near the hypocenter. Resultant catalogs of seismicity also include many complex multi-
fault ruptures consistent with recent observations (e.g., the 2016 Kaikoura M 7.8
earthquake) and the UCERF3 model. Figure 4.4 shows one such particularly complex
multi-fault rupture, an M 7.69 which nucleates on the Compton fault before spreading to
76
the Newport-Inglewood, San Pedro Escarpment, Palos Verdes, and ten other UCERF3
fault sections.
We use a synthetic catalog simulated on the full statewide UCERF3 fault system
with the hybrid-loading approach discussed in Shaw (2019), discretized into 533,887
individual triangular fault patches each with an average area of 0.67 km
2
. This catalog
incorporates three RSQSim improvements implemented since the catalog used in Shaw
et al. (2018), each of which is aimed to improve rupture propagation velocity, 𝑣 𝑝𝑟𝑜𝑝 , and
rupture ground motions: finite receiver patch geometry, variable slip speed, and static-
elastic time delay. Figure 4.5 plots mean 𝑣 𝑝𝑟𝑜𝑝 as a function of patch hypocentral
distance, 𝑅 ℎ𝑦𝑝
, for each model iteration. Mean propagation velocities for the catalog
used in Shaw et al. (2018) are 𝑣 𝑝𝑟𝑜𝑝 ∈ [0.9,1.6] km/s for 𝑅 ℎ𝑦𝑝
> 10 km (dashed line in
Figure 4.5), which is significantly lower than typical values of 𝑣 𝑝𝑟𝑜𝑝 ∈ [2.1,2.4] km/s if we
assume a shear-wave velocity, 𝛽 = 3 km/s, and the relations established in Andrews
(1976) and Geller (1976) of 𝑣 𝑝𝑟𝑜𝑝 ∈ [0.7𝛽 , 0.8𝛽 ]. This propagation velocity deficit
resulted in low fault-parallel rupture directivity in numerical ground motion simulations.
The first modification improves the accuracy of the stiffness matrices, 𝐾 𝜏 and 𝐾 𝜎
for shear and normal stresses respectively. Previously, this calculation considered the
finite geometry of the source patch which slips, but employed a single point
representation of the receiving patch. Now, the calculated stiffness is averaged over the
finite geometry of the receiver patch. This results in increased 𝑣 𝑝𝑟𝑜𝑝 ∈ [1.5,2.25] km/s for
𝑅 ℎ𝑦𝑝
> 10 km, and greatly increased 𝑣 𝑝𝑟𝑜𝑝 of patches neighboring the hypocenter to
77
supershear speeds on average, in this case, 𝑣 𝑝𝑟𝑜𝑝 = 4.4 km/s for 𝑅 ℎ𝑦𝑝
≃ 1 km (dotted
line in Figure 4.5). This increases ground motions in aggregate, but fault-parallel
directivity with this model is still less than that predicted empirically.
The second modification eliminates the fixed sliding speed approximation during
earthquake slip (state 2). That slip velocity is determined by the shear impedance
relationship (Brune, 1970),
𝑣 𝑒𝑞
=
2𝛽𝛥𝜏 𝐺 ,
(4.1)
where 𝛥𝜏 is the difference between the shear stress at initiation of slip and and sliding
friction, and G is the shear modulus. Previously, that value has been estimated a priori
for all patches using typical values of 𝛥𝜏 , with a uniform 𝑣 𝑒𝑞
= 1 m/s used in Shaw et al.
(2018). Here we introduce a variable slip speed version of RSQSim, where 𝑣 𝑒𝑞
for each
patch is set according to equation (4.1) from 𝛥 𝜏 at the moment it enters state 2, and
then updated stepwise during a slip episode. This stepwise discretization is
necessitated by the efficient state-driven RSQSim computational scheme, which cannot
accommodate continuously variable 𝑣 𝑒𝑞
. Updates occur whenever the instantaneous
calculated velocity leaves the range [𝑣 𝑒𝑞
∗ 𝑣 𝑓𝑎𝑐𝑡 , 𝑣 𝑒𝑞
/𝑣 𝑓𝑎𝑐𝑡 ], where 𝑣 𝑓𝑎𝑐𝑡 ∈ (0,1) is a
dimensionless constant and 𝑣 𝑓𝑎𝑐𝑡 = 0.8 in this study. A velocity floor is also applied as a
multiplicative factor of the initially calculated slip speed in each slip episode, 0.3 for this
study. Figure 4.6 shows the slip-time evolution on a single patch with seven slip
episodes from the rupture depicted in Figures 4.2 and 4.3. Variable slip speed models
78
exhibit greatly increased propagation velocities, with mean 𝑣 𝑝𝑟𝑜𝑝 ∈ [2,4] km/s for 𝑅 ℎ𝑦𝑝
>
10 km and 𝑣 𝑝 𝑟 𝑜𝑝
= 17.5 km/s for 𝑅 ℎ𝑦𝑝
≃ 1 km (dotted and dashed line in Figure 4.5).
Thirdly, we add a time-delay to the static elastic interaction in an attempt to
suppress the supershear (and sometimes super-P) rupture velocities prevalent in the
variable slip speed model (especially when 𝑅 ℎ𝑦𝑝
< 10 km). Previously, all stress
changes were felt instantaneously on all patches in RSQSim. Ideally, the stress change
experienced on the 𝑖 𝑡 ℎ
patch due to slip on the 𝑗 𝑡 ℎ
patch would be delayed by
𝛥 𝑡 𝑐 𝑎 𝑢 𝑠𝑎𝑙
(𝑖 , 𝑗 ) =
𝑑 (𝑖 ,𝑗 )
𝛽 , where 𝑑 (𝑖 , 𝑗 ) is the distance between the centers of the two
patches. This would necessitate a unique time-delay for each source and receiver patch
pair, breaking the RSQSim paradigm of fixed state-driven updates. Instead, we simplify
this to a single fixed-time delay applied for all 𝑖 ≠ 𝑗 , tuned to roughly match the expected
delay for each patch's immediate neighbors (as it is those interactions that control the
rupture velocity). We calculate this time delay from the patch's shear self-stiffness,
𝐾 𝜏 (𝑗 , 𝑗 ), valid if neighboring elements are of similar size and shape:
𝛥 𝑡 𝑐𝑎𝑢𝑠𝑎𝑙 (𝑖 , 𝑗 ) =
𝜂𝐺
𝛽 |𝐾 𝜏 (𝑗 ,𝑗 )|
,
(4.2)
where 𝜂 is a tuneable dimensionless constant. We use value of 𝜂 = 0.7, which leads to
reasonable average propagation velocities for our final model in the range 𝑣 𝑝𝑟𝑜𝑝 ∈
[1.5,3] km/s for all 𝑅 ℎ𝑦𝑝
(solid line in Figure 4.5).
79
The total synthetic catalog length 269,521 years, from which we consider all
95,803 ruptures with M ≥ 6.5. Considerable runtime is required to generate catalogs
which are sufficiently long for robust PSHA calculations, with ~32,000 synthetic years
possible in a single 48 wall-clock hour job (the maximum allowed) on 64 Frontera
compute nodes (3584 processors) at the Texas Advanced Computing Center. To
reduce the time to solution, we configured ten jobs in parallel with different sets of small
spatially uncorrelated additive perturbations to the initial shear stress of each patch in
the range [0,1) MPa. We threw out the first 5000 simulated years of each calculation to
avoid spin-up effects and ensure that each of the ten catalogs diverged uniquely from
initial conditions, stitching the remaining synthetic time from each catalog together (with
time gaps sampled randomly from the distribution of prior M ≥ 6.5 recurrence intervals)
to form the final long catalog. Model output includes the time at which each triangular
fault patch transitions between states (i.e., starts or stops sliding, or updates sliding
velocity), which we convert to the Standard Rupture Format (SRF; Graves, 2014)
source description for use as input to our wave propagation codes. The SRF description
requires that slip velocity for each patch be evenly discretized in time. We achieve this
by computing the total displacement that occurred during each SRF time step, and then
setting the velocity as that displacement divided by the time step duration. This ensures
that the total displacement (and thus seismic moment released) is preserved for each
rupture, but results in time steps where the velocity differs from that in the RSQSim
model output, notably smoothing out short-lived velocity peaks. Figure 4.6 shows the
slip-time history for a single patch before (black lines) and after (gray lines) this
adjustment.
80
4.4 Ground Motion Model: CyberShake
The CyberShake platform (Jordan et al., 2018; Graves el al., 2011) replaces
empirical GMMs with deterministic three-dimensional (3D) ground motion simulations
through 3D seismic velocity models, characterizing the effects of basin response and
other path effects which are either parameterized or treated as aleatory variability in
empirical GMMs. In prior studies, CyberShake has relied on traditional ERFs to define
rupture properties (magnitude, geometry, probability) and a stochastic kinematic rupture
generator (Graves and Pitarka, 2014) to estimate a suite of slip-time histories for each
rupture with prescribed rupture properties. Ground motion variability in CyberShake has
previously been represented through ensembles of rupture variations, which vary the
hypocenter location and kinematic slip distributions of each input ERF rupture, as well
as aleatory magnitude variability for a given source surface to produce ruptures with a
variety of stress drops. The process of automated kinematic rupture generation is
currently limited to contiguous planar or ribbon-like faults, which has restricted
CyberShake to the older UCERF2 (Field et al., 2009) model which lacks multi-fault
ruptures. As such, previous CyberShake studies do not represent complex multi-fault
ruptures such as the one depicted in Figure 4.4, which are also prevalent in the
UCERF3 model.
For this study, we instead couple the RSQSim model with CyberShake to create the first
PSHA model for a complex fault system composed entirely of fully-deterministic
physical models. We use the Southern California Earthquake Center (SCEC) CVM-
S4.26.M01 tomographically inverted velocity model with geotechnical layer (Lee et al.,
81
2014; Small et al., 2017) to perform 0.5 hz simulations for the entire RSQSim catalog for
a set of 10 sites in the Los Angeles region. We consider three-second RotD50 (Boore,
2010) pseudo-spectral accelerations (PSAs), as this is the shortest spectral period
which can be reliably represented by CyberShake 0.5 hz simulations. CyberShake
utilizes seismic reciprocity (Zhao et al., 2006) to efficiently compute synthetic
seismograms from hundreds of thousands of sources at a single site through just two
3D strain Green tensor simulations (one for each horizontal component of resultant
seismograms). This efficiency makes simulation of all ruptures within 200 kilometers of
a site in a 269,521 year RSQSim catalog tractable on modern supercomputers.
4.5 Ground Motion Variability in 1D: SCEC BroadBand Platform
We initially prototyped this calculation with 3D simulations in a one-dimensional
(1D) layered earth structure in the SCEC BroadBand Platform (BBP; Maechling et al.,
2014) version 19.4.0, calculating RotD50 spectra (𝑇 ∈ [1,10] s) for each RSQSim
rupture at a number of sites in Los Angeles. We use the Los Angeles basin velocity
model with 𝑉 𝑆 30
= 500 m/s and the deterministic low-frequency component of the
Graves and Pitarka (2016) simulation method to simulate ground motions for RSQSim
slip-time histories (after conversion to the SRF representation). Figure 4.7a plots
spectra at the USC site for a single M 7.45 rupture on the Mojave section of the San
Andreas Fault (SAF) calculated with the BBP and compared with the four leading
empirical GMMs from the Enhancement of Next Generation Attenuation Relationships
for Western US project (henceforth NGA-West2; Bozorgnia et al. 2014; Abrahamson et
al., 2014; Boore et al., 2014; Campbell and Bozorgnia, 2014; Chiou and Youngs, 2014).
82
GMM comparisons with BBP are parameterized with the same 𝑉 𝑆 30
= 500 m/s site
condition, and basin depth (𝑍 1.0
= 0.2 or 𝑍 2.5
= 2.5 (km) depending on the model) set
according to the 1D velocity profile. We use the methodology described in Shaw et al.
(2018) to parametrize RSQSim ruptures for empirical GMM estimation, first mapping the
arbitrarily complex RSQSim ruptures to their corresponding UCERF3 fault subsections.
This mapping step is included in order to correct rupture-site distances (𝑅 𝑅 𝑢𝑝
and 𝑅 𝐽𝐵
)
for complex ruptures which may include one or more patches from faults other than the
main rupture surface that might be close to the site of interest, but not significantly
contribute to the ground motion at that site. In that case, distances are instead
calculated to the nearest UCERF3 fault subsection surface where at least 20% of the
subsection (by area) participates in the RSQSim rupture.
Although no individual rupture is expected to match any given mean GMM
spectra exactly, the distribution of spectra from all similar ruptures (Figure 4.7b) can be
statistically compared to GMMs. In aggregate, 1D BBP results for this fault are similar to
and within the uncertainties of GMM predictions, though the variability is less than
predicted by ergodic GMMs. We plot a histogram of simulated amplitudes for a single
spectral period of 3 s in Figure 4.8a, along with comparisons with empirical GMM
predicted log-normal distributions. Here, reduced variability of the one-dimensional BBP
results is more apparent (especially at the tails of the GMM distributions). The lower
variability, however, is expected when contrasting simulations of a single seismic
source, repeatedly sampling the same path from source to site (here in a 1D layered
velocity structure) and with the same rupture incidence angle, with ergodic GMMs.
83
Expanding this analysis to ten sites (listed in Table 4.1 and mapped in Figure
4.9) and across all RSQSim ruptures that within 200 km of each site in our catalog, we
turn to the z-score statistical measure to explore total variability of RSQSim simulated
ground motions. The z-score, 𝑧 (𝑖 , 𝑗 ), for synthetic ground motion 𝑌 𝑠𝑖𝑚 (𝑖 , 𝑗 ) of the 𝑖 𝑡 ℎ
rupture simulated at the 𝑗 𝑡 ℎ
site, in terms of log-mean, 〈𝑙𝑛 𝑌 𝐺𝑀𝑀
(𝑖 , 𝑗 ) 〉, and total
standard deviation, 𝜎 𝐺𝑀𝑀
(𝑖 , 𝑗 ), GMM predictions is
𝑧 (𝑖 , 𝑗 ) =
𝑙𝑛 𝑌 𝑠 𝑖 𝑚 (𝑖 ,𝑗 ) − 〈𝑙𝑛 𝑌 𝐺𝑀𝑀 (𝑖 ,𝑗 ) 〉
𝜎 𝐺𝑀𝑀 (𝑖 ,𝑗 )
.
(4.3)
z-scores measure, for a single site and rupture, how many standard deviations away the
simulated ground motion, 𝑌 𝑠𝑖𝑚 (𝑖 , 𝑗 ), is from the log-mean GMM prediction,
〈𝑙𝑛 𝑌 𝐺𝑀𝑀
(𝑖 , 𝑗 ) 〉. This allows us to statistically compare the variability of RSQSim-BBP
ground motions with those predicted of GMMs to quantify bias and to compare the total
variability across many site locations, magnitudes and distances. A histogram of z-
scores from all RSQSim-BBP simulations compared against the Abrahamson et al.
(2014; henceforth ASK2014) GMM is plotted as a gray histogram in Figure 4.10a. The
mean of -0.11 (in units of 𝜎 𝐺𝑀𝑀
) means that this RSQSim-BBP model has a slight
underprediction bias relative to the GMM predictions. The narrower distribution of the z-
score histogram (relative to the standard normal plotted as a black line) with a fractional
standard deviation (𝜎 -fract) of 0.63 illustrates the lower variability from 1D RSQSim-BBP
relative to the ASK2014 GMM. This is expected with a 1D velocity structure lacking
basin structures and heterogeneities to focus and scatter seismic energy.
84
4.6 Ground Motion Variability In 3D: CyberShake
We next moved to the SCEC CyberShake platform to perform full 3D
deterministic simulations. Heterogeneities in the CVM-S4.26.M01 velocity model
introduce additional ground motion variability (relative to the 1D case), and deep
sedimentary basins amplify shaking in many areas. z-score histograms for the same set
of sites and ruptures as before, but instead computed with CyberShake (Figure 4.10b)
illustrate this increased variability, along with larger ground motions in aggregate than
those from both the GMM and BBP. Each CyberShake site considered has the same
surface 𝑉 𝑆 = 500 (m/s), which we consider equivalent to the 𝑉 𝑆 30
= 500 (m/s) site proxy
in the ASK2014 GMM, as the mesh gridpoint spacing in 0.5 hz CyberShake calculations
is larger than 30 m. We set the 𝑍 1.0
basin depth proxy value in ASK2014 from the CVM-
S4.26.M01 velocity model, with values listed in Table 4.1. The mean z-score increases
to 0.40, meaning that RSQSim-CyberShake 3 s RotD50 ground motions are on average
0.40 ASK2014 standard deviations above from the mean ASK2014 prediction, and the
𝜎 -fract of 0.83 indicates that the RSQSim-CyberShake model for these ten sites
contains 83% of the total variability in ASK2014. z-score histograms for each individual
site are given in supplemental material. The bias toward higher ground motions for soft-
soil and basin sites has been observed in previous CyberShake studies, though the bias
in RSQSim-CyberShake is lower than in the prior CyberShake Study 15.4 which used
the UCERF2 ERF extended with the Graves and Pitarka (2014) kinematic rupture
generator. z-scores from Study 15.4 are plotted in Figure 4.10c, with a mean z-score of
1.00 for this same set of ten sites, and a similar amount of total variability (𝜎 -fract of
0.85) to the RSQSim-CyberShake model.
85
For the single source-site pair of SAF Mojave with 7 ≤ M ≤ 7.5 ruptures simulated
at USC (Figure 4.8b), variability remains lower than the ergodic GMM prediction, though
for this site, the distribution shifts to the right due to deeper soft soils in the CVM-
S4.26.M01 model at this site than in the one-dimensional BBP model plotted in Figure
4.8a. This is expected due to the repeated sampling of a single path, this time through a
three-dimensional medium, from the same source to the USC site.
4.7 Procedure to Decompose Variance Through Rotation and Translation
Figure 4.10b suggests that RSQSim-CyberShake ground motions contain a
similar (though slightly lower) amount of total variability as empirical GMMs, but that
alone is not sufficient to validate the variance structure of the model. The z-score
histogram of Figure 4.10b is convenient in that it summarizes variability across many
sites, magnitudes, and distances, but it can hide model trade-offs and biases. To
illustrate this, consider two alternative candidate models: one ergodic which matches
the empirical GMM median and total standard deviation perfectly across all unique
combinations of explanatory variable values (e.g., magnitude and distance), and
another non-ergodic whose predicted standard deviation is smaller than the empirical
GMM for each unique combination and mean residuals (relative to the empirical GMM)
are sampled from a normal distribution. Both models could have identical z-score
histograms (with zero mean), even with very different variance structures. This scenario
is depicted in Figure 4.11 with a cartoon of z-scores from an ergodic model in Figure
4.11a and a non-ergodic model with the sum of many narrow distributions in Figure
86
4.11b. We must thus decompose the model variance to verify that its individual
components are in line with GMM predictions, not only its total variance.
Al Atik et al. (2010) explains how empirical GMMs typically decompose total
variance into two sources, both zero-mean and independent normally distributed
random variables: between-event (denoted 𝛿𝐵 with standard deviation 𝜏 ) and within-
event (denoted 𝛿𝑊 with standard deviation 𝜙 ). 𝜏 and 𝜙 are typically estimated through
mixed-effects regression (Abrahamson and Youngs, 1992; Al Atik and Abrahamson,
2010), which accounts for irregularly sized and sparse datasets. Wang and Jordan
(2014) demonstrated an averaging-based factorization approach to estimate 𝜏 and 𝜙 for
simulation-based PSHA using a prior Southern California CyberShake study, though
their methodology is best suited for larger datasets (they used 235 CyberShake sites).
Results from both approaches have been shown to be similar for the CyberShake
datasets (X. Meng and C. A. Goulet, personal communication 9 March 2020).
We developed a technique to estimate 𝜏 and 𝜙 𝑠𝑠
(the latter being the single-site
component of 𝜙 , without unexplained site-to-site variance) from RSQSim ruptures using
a smaller set of 𝑁 𝑠𝑖𝑡𝑒 = 10 CyberShake sites, all chosen with the same approximate soil
classification (in this case modeled surface 𝑉 𝑠 = 500 𝑚 /𝑠 ) and distributed throughout
Southern California (Figure 4.9). Instead of performing hundreds of computationally
expensive forward ground motion simulations for individual RSQSim ruptures in order to
estimate 𝜏 and 𝜙 𝑠 𝑠 , we take advantage of the seismic reciprocity assumption employed
in CyberShake by rotating and translating ruptures around each CyberShake site. This
87
allows us to perform only two CyberShake Green’s function simulations (one for each
horizontal component) at each site, and recover ground motions sampling hundreds of
thousands of different combinations of ruptures, rupture orientations, paths, and
distances. The simulation scheme we are about to describe is summarized in Figure
4.12 and defines a scenario as a type of earthquake with explicitly enumerated criteria
(e.g., M 6.6 strike-slip), a rupture as a unique event (slip-time history, magnitude,
geographic extent) which occurred in the RSQSim simulation, and a rotated rupture as
a modified rupture slip-time history which has been rotated and translated horizontally in
space. We first selected three sets of𝑁 𝑟𝑢𝑝 = 50 unique ruptures which satisfy the criteria
for the following scenarios:
1) M 6.6 vertical strike-slip: 𝗠 ∈ [6.55,6.65], 𝑧 𝑇𝑂𝑅
∈ [0,1], 𝑟𝑎 𝑘 𝑒 ∈ {−180,0,180},𝑑𝑖𝑝 =
90, linear (maximum 0.5 km deviation from ideal linear source)
2) M 6.6 reverse: 𝗠 ∈ [6.55,6.65], 𝑧 𝑇𝑂𝑅
∈ [1,5], 𝑟𝑎𝑘𝑒 ∈ [80,100],𝑑𝑖𝑝 ∈ [35,55]
3) M 7.2 vertical strike-slip: 𝗠 ∈ [7.15,7.25], 𝑧 𝑇𝑂𝑅
∈ [0,1], 𝑟𝑎𝑘𝑒 ∈ {−180,0,180},𝑑𝑖𝑝 =
90, linear (maximum 5% deviation from ideal linear source)
More than 50 candidate ruptures exist for each scenario; we narrowed the set
down to 𝑁 𝑟𝑢𝑝 = 50 ruptures for each scenario by choosing a random set which
minimizes the number of repeat ruptures on individual fault sections. This results in
ruptures from a diverse set of faults, with 52, 14, and 51 different fault sections
participating in ruptures for scenarios one, two, and three respectively. The counts for
scenarios one and three are greater than 𝑁 𝑟𝑢𝑝 because multiple sections can participate
88
in a single rupture. The two M 6.6 scenarios were selected to approximate the larger
magnitude scenarios proposed in the Goulet et al. (2014) BBP “Part B” validation
exercise, which were in turn chosen due to an abundance of data at those magnitudes
to constrain empirical GMMs. We consider 𝑁 𝑑𝑖𝑠𝑡 = 3 fixed 3D distances (𝑅 𝑅𝑢𝑝 ) between
each rupture and site: {20,50,100} km (the first two distances were also taken from the
BBP “Part B” exercise). First, we rotate and translate each rupture such that its scalar
moment centroid is due north of our site, its strike is zero (due north, following the Aki &
Richards, 1980, convention), and the minimum distance between the site and the
rotated rupture is as prescribed. We then repeat this for a total of 𝑁 𝑎𝑧
= 18 different
strike azimuths (each 20 degrees apart) by first rotating each rupture in place about its
centroid, then translating the rotated rupture toward or away from the site to correct for
any changes to 𝑅 𝑅𝑢𝑝 as a result of this rotation. We repeat this procedure sampling
𝑁 𝑝𝑎𝑡 ℎ
= 18 different paths from the site to the rotated rupture (each 20 degrees apart),
rotating the rupture about the site in addition to rotating the rupture in place for each
path. This procedure is more clearly illustrated through a schematic with just 5 rupture
strike azimuth and path rotations each in Figure 4.12. We simulate the ground motions
in CyberShake for each unique combination of rupture 𝑖 , site 𝑗 , distance 𝑘 , strike
azimuth 𝑙 , and path 𝑚 , for scenario 𝑛 , 𝑌 𝑠𝑖𝑚 (𝑖 , 𝑗 , 𝑘 , 𝑙 , 𝑚 , 𝑛 ). Weighting of irregularly
sampled data is a primary motivation for the mixed-effects regression technique
employed in empirical studies; we avoid this by prescribing the same number of
simulations for each unique RSQSim rupture, 𝑁 𝑎𝑧
× 𝑁 𝑝𝑎𝑡 ℎ
= 324, with a total of
𝑁 𝑎𝑧
× 𝑁 𝑝𝑎𝑡 ℎ
× 𝑁 𝑑𝑖𝑠𝑡 × 𝑁 𝑟𝑢𝑝 = 48,600 simulated ground motions for each scenario
computed at a single site.
89
We estimate 𝜏 and 𝜙 𝑠𝑠
separately for each scenario 𝑛 , 𝜏 (𝑛 ) and 𝜙 𝑠𝑠
(𝑛 ). Values
are reported in Table 4.2 for 𝜏 and Table 4.3 for 𝜙 𝑠𝑠
, including each individual distance
in order to test distance-dependence, the range of distance-independent values across
each site, and also total site- and distance-independent 𝜏 and 𝜙 𝑠𝑠
. First we compute the
event term (natural-log median ground motion) for each rupture simulated at site 𝑗 ,
𝐵 (𝑖 , 𝑗 , 𝑘 , 𝑛 ) = 𝑚𝑑𝑛 . {𝑙𝑛 (𝑌 𝑠𝑖𝑚 (𝑖 , 𝑗 , 𝑘 , 𝑙 , 𝑚 , 𝑛 )) | 𝑙 ∈ 1. . 𝑁 𝑎𝑧
, 𝑚 ∈ 1. . 𝑁 𝑝𝑎𝑡 ℎ
}.
(4.4)
This is analogous to 𝛿 𝐵 𝑒 in Al Atik et al. (2010) except that here it is computed for single
site at a fixed distance (but still sampling many paths and rupture orientations, because
of the experiment design), and not expressed as a residual relative to an empirical
model (we’re after the variability component only). We take the between-event standard
deviation computed for site 𝑗 , for ruptures matching scenario 𝑛 , and at a fixed distance
𝑘 , to be the standard deviation of the set of all 𝐵 (𝑖 , 𝑗 , 𝑘 , 𝑛 ),
𝜏 (𝑗 , 𝑘 , 𝑛 ) = 𝑠𝑡𝑑 . 𝑑𝑒𝑣 . {𝐵 (𝑖 , 𝑗 , 𝑘 , 𝑛 ) | 𝑖 ∈ 1. . 𝑁 𝑟𝑢𝑝 }.
(4.5)
ASK2014 and the other of the NGA-West2 relations outlined in Bozorgnia et al. (2014)
assume distance-independence of 𝜏 for 3-second RotD50 (which is consistent with our
results), so we compute distance-independent 𝜏 (𝑗 , 𝑛 ) as the mean of 𝜏 (𝑗 , 𝑘 , 𝑛 ) across all
distances,
90
𝜏 (𝑗 , 𝑛 ) = 〈{𝜏 (𝑗 , 𝑘 , 𝑛 ) | 𝑘 ∈ 1. . 𝑁 𝑑𝑖𝑠𝑡 } 〉. (4.6)
We extend this to multiple sites by computing event terms across all sites, again
ensuring that each site has the same approximate soil classification (in this case
surface 𝑉 𝑠 = 500 𝑚 /𝑠 for all sites),
𝐵 (𝑖 , 𝑘 , 𝑛 ) = 𝑚𝑑𝑛 . {𝑙 𝑛 (𝑌 𝑠𝑖𝑚 (𝑖 , 𝑗 , 𝑘 , 𝑙 , 𝑚 , 𝑛 )) | 𝑗 ∈ 1. . 𝑁 𝑠𝑖𝑡𝑒 , 𝑙 ∈ 1. . 𝑁 𝑎𝑧
, 𝑚 ∈ 1. . 𝑁 𝑝𝑎𝑡 ℎ
}.
(4.7)
Similar to (5) and (6), we take the between-event standard deviation across all sites for
scenario 𝑛 at distance 𝑘 as the standard deviation of the set of all 𝐵 (𝑖 , 𝑘 , 𝑛 ),
𝜏 (𝑘 , 𝑛 ) = 𝑠𝑡𝑑 . 𝑑𝑒𝑣 . {𝐵 (𝑖 , 𝑘 , 𝑛 ) | 𝑖 ∈ 1. . 𝑁 𝑟𝑢𝑝 },
(4.8)
and the total between-event standard deviation for all ruptures matching scenario 𝑛 as
the mean of all 𝜏 (𝑘 , 𝑛 ),
𝜏 (𝑛 ) = 〈{𝜏 (𝑗 , 𝑘 , 𝑛 )} 〉. (4.9)
We also compute residuals for each synthetic ground motion with respect to 𝐵 (𝑖 , 𝑗 , 𝑘 , 𝑛 ),
to compute the remaining within-event variability,
𝛿𝑊 (𝑖 , 𝑗 , 𝑘 , 𝑙 , 𝑚 , 𝑛 ) = 𝑙𝑛 (𝑌 𝑠𝑖𝑚 (𝑖 , 𝑗 , 𝑘 , 𝑙 , 𝑚 , 𝑛 )) − 𝐵 (𝑖 , 𝑗 , 𝑘 , 𝑛 ), (4.10)
which is analogous to 𝛿 𝑊 𝑒𝑠
in Al Atik et al. (2010). We take the single-site within-event
standard deviation, computed at site 𝑗 , for ruptures matching scenario 𝑛 , and at a fixed
distance 𝑘 to be
91
𝜙 𝑠𝑠
(𝑗 , 𝑘 , 𝑛 ) = 𝑠𝑡𝑑 . 𝑑𝑒𝑣 . {𝛿𝑊 (𝑖 , 𝑗 , 𝑘 , 𝑙 , 𝑚 , 𝑛 ) | 𝑖 ∈ 1. . 𝑁 𝑟𝑢𝑝 , 𝑙 ∈ 1. . 𝑁 𝑎 𝑧 , 𝑚 ∈ 1. . 𝑁 𝑝𝑎𝑡 ℎ
}.
(4.11)
We estimate single-site within-event standard deviation computed for site 𝑗 and across
all distances as
𝜙 𝑠𝑠
(𝑗 , 𝑛 ) = 𝑠𝑡𝑑 . 𝑑𝑒𝑣 .
{𝛿𝑊 (𝑖 , 𝑗 , 𝑘 , 𝑙 , 𝑚 , 𝑛 ) | 𝑖 ∈ 1. . 𝑁 𝑟𝑢𝑝 , 𝑘 ∈ 1. . 𝑁 𝑑𝑖𝑠𝑡 , 𝑙 ∈ 1. . 𝑁 𝑎𝑧
, 𝑚 ∈ 1. . 𝑁 𝑝𝑎𝑡 ℎ
},
(4.12)
across all sites for a single distance 𝑘 as,
𝜙 𝑠𝑠
(𝑘 , 𝑛 ) = 𝑠𝑡𝑑 . 𝑑𝑒𝑣 .
{𝛿𝑊 (𝑖 , 𝑗 , 𝑘 , 𝑙 , 𝑚 , 𝑛 ) | 𝑖 ∈ 1. . 𝑁 𝑟𝑢𝑝 , 𝑗 ∈ 1. . 𝑁 𝑠𝑖𝑡𝑒 , 𝑙 ∈ 1. . 𝑁 𝑎𝑧
, 𝑚 ∈ 1. . 𝑁 𝑝𝑎𝑡 ℎ
},
(4.13)
and for all ruptures matching scenario 𝑛 as the standard deviation of all residuals for 𝑛 ,
𝜙 𝑠𝑠
(𝑛 ) = 𝑠𝑡𝑑 . 𝑑𝑒𝑣 . {𝛿𝑊 (𝑖 , 𝑗 , 𝑘 , 𝑙 , 𝑚 , 𝑛 ) | 𝑖 ∈ 1. . 𝑁 𝑟𝑢𝑝 , 𝑗 ∈ 1. . 𝑁 𝑠𝑖𝑡𝑒 , 𝑘 ∈ 1. . 𝑁 𝑑𝑖𝑠𝑡 , 𝑙 ∈
1. . 𝑁 𝑎𝑧
, 𝑚 ∈ 1. . 𝑁 𝑝𝑎𝑡 ℎ
}.
(4.14)
4.8 Variance Decomposition Results
Table 4.2 summarizes computed between-event standard deviations for each
scenario, 𝜏 (𝑛 ), along with the range of values computed at a single site, 𝜏 (𝑗 , 𝑛 ), and
values for each distance (computed across all sites), 𝜏 (𝑘 , 𝑛 ). The simulated 𝜏 (𝑛 ) values
of 0.17 for both M 6.6 scenarios and 0.20 for the M 7.2 scenario are significantly lower
than the values of 0.38 and 0.36 determined in the ASK2014 regressions for M 6.6 and
M 7.2 respectively, indicating that the model lacks sufficient between-event variability by
92
this metric. One potential explanation for this discrepancy is that uncertainty on 𝛿 𝐵 𝑒 in
GMM regressions (e.g., due to earthquakes with few recordings per event) may
artificially inflate 𝜏 . Figure 4.13 plots 𝜏 (𝑘 , 𝑛 ) computed from suites of downsampled
synthetic data as a function of the number of simulated recordings per event, 𝑁 𝑟𝑒𝑐 . For
each 𝑁 𝑟𝑒𝑐 value, we estimate downsampled 𝐵 (𝑖 , 𝑘 , 𝑛 ) from a randomly sampled subset
of 𝑁 𝑟𝑒𝑐
ground motions for each rupture, then compute downsampled 𝜏 (𝑘 , 𝑛 ) following
equation 4.8. We do this 100 times for each 𝑁 𝑟𝑒 𝑐 in order to compute a distribution of
downsampled 𝜏 (𝑘 , 𝑛 ), the median of which decreases with increasing 𝑁 𝑟𝑒𝑐
and
asymptotes to the full simulated 𝜏 (𝑘 , 𝑛 ) value. This indicates that our modeled 𝜏 (𝑘 , 𝑛 ),
where 𝐵 (𝑖 , 𝑘 , 𝑛 ) is well constrained by a large set of 𝑁 𝑠𝑖𝑡𝑒 × 𝑁 𝑎𝑧
× 𝑁 𝑝𝑎𝑡 ℎ
= 3240
simulations, should be lower than 𝜏 estimated from sparse data; however, this effect
only explains part of the 𝜏 discrepancy.
Another potential explanation is that the magnitudes in the GMM may be
associated with a wider range of stress drop and magnitude-area estimates than in our
simulations, which are coming directly out of the RSQSim model. More generally,
magnitudes and distances are perfectly known in our experiment, whereas those source
parameters of real earthquakes are uncertain. We chose to select ruptures from the
RSQSim catalogs with magnitudes in a tight range (e.g., 𝗠 = [7.15,7.25] for the M 7.2
scenario), and distances are fixed by construction. These differences, along with the
previously discussed recordings-per-event dependence, highlights how the paucity of
data used in empirical models can inflate 𝜏 , and thus our simulated 𝜏 (𝑛 ) cannot be
considered a perfect analog to empirical 𝜏 . That said, this analysis still reveals a
93
potential lack of between-event variability that warrants consideration in the
development of future RSQSim models. Table 4.2 also indicates that, as in the empirical
models, we do not see a significant dependence of 𝜏 (𝑘 , 𝑛 ) on distance.
Table 4.3 summarizes computed within-event standard deviations for each
scenario, 𝜙 𝑠𝑠
(𝑛 ), along with the range of values computed at a single site, 𝜙 𝑠𝑠
(𝑗 , 𝑛 ), and
values for each distance (computed across all sites), 𝜙 𝑠𝑠
(𝑘 , 𝑛 ). The simulated 𝜙 𝑠𝑠
(𝑛 )
values are 0.49, 0.45, and 0.51 for the M 6.6 strike-slip, M 6.6 reverse, and M 7.2 strike-
slip scenarios respectively. The NGA-West2 relations don’t explicitly define 𝜙 𝑠𝑠
, but Al
Atik (2015) estimates a value of 0.37 (their Table 5.11) and Lin et al. (2011) a value of
0.42 (their Table 3, column 𝜎 𝑟 ). These estimates are both below our simulated 𝜙 𝑠𝑠
(𝑛 ),
indicating that our model produces a sufficient amount of within-event ground motion
variability. One feature of note from Table 4.3 is that our 𝜙 𝑠𝑠
(𝑛 , 𝑘 ) increases significantly
with distance, which is not included in ASK2014, Al Atik (2015), or Lin et al. (2011).
Boore et al. (2014) is the only NGA-West2 relation that includes distance dependence
of 𝜙 , but that dependence is only for distances greater than 130 km.
4.9 Fully Physics-Based Hazard Curves
We compute fully non-ergodic PSHA hazard curves from the suite of simulated
ground motions for the complete synthetic catalog of 𝑁 𝑟𝑢𝑝 = 95,803 RSQSim ruptures.
This process is similar to traditional GMM-based hazard curves, except that each 𝑖 𝑡 ℎ
rupture has a single simulated ground motion intensity at the 𝑗 𝑡 ℎ
site, 𝑌 𝑠𝑖𝑚 (𝑖 , 𝑗 ), rather
than an empirically predicted log-normal distribution, and the same annual rate of
94
occurrence, 𝑅 , equal to the reciprocal of the simulated catalog duration, 𝛥𝑇 , expressed
in years:
𝑅 = 𝛥𝑇
−1
. (4.15)
The number of ruptures whose simulated intensity, 𝑌 𝑠𝑖𝑚 (𝑖 , 𝑗 ), exceeds intensity
measure level 𝐿 at site 𝑗 , modified from equation 11 of Wang and Jordan (2014) where
𝐻 is the Heaviside step function, is
𝑁 𝑒𝑥𝑐𝑒𝑒𝑑 (𝑗 , 𝐿 ) = ∑
𝑁 𝑟𝑢𝑝
𝑖 =1
𝐻 [𝑌 𝑠𝑖𝑚 (𝑖 , 𝑗 ) − 𝐿 ].
(4.16)
The time-independent annual probability of exceeding 𝐿 at site 𝑗 can then be
expressed, assuming Poissonian behavior, as
𝑃 (𝑌 𝑠𝑖𝑚 > 𝐿 |𝑗 , 𝑅 ) = 1 − 𝑒 −𝑅 ∗𝑁 𝑒𝑥𝑐𝑒𝑒𝑑 (𝑗 ,𝐿 )
. (4.17)
This results in a minimum possible nonzero probability for a single exceedance,
𝑁 𝑒𝑥𝑐𝑒𝑒𝑑 = 1, of
𝑃 𝑚𝑖𝑛 = 1 − 𝑒 −𝑅 (4.18)
for any L. Longer simulations give more opportunities to exceed large intensity
measures, and thus can probe lower probability regions, but typical hazard levels of
interest (e.g., 1,000-10,000 year return periods) are well resolved with our 269,521 year
RSQSim catalog. Hazard curves expressing the probability of exceeding 3 s RotD50
PSA at the USC site are shown in Figure 4.14a. Figure 4.14b illustrates the effect of
simulated catalog length on Pmin and resultant hazard curves, which show little
95
variability for return periods less than 10,000 years. At low exceedance probabilities (<
10
−4
), the slope of the simulated hazard curve is most similar to truncated GMM curves
(depicted with dashed and dotted lines in Figure 4.13a).
The effect of narrower source-site distributions in the non-ergodic RSQSim-
CyberShake model can be seen around 0.01-0.05 g. Here, both the CyberShake and
ASK2014 models are controlled by contributions from the San Andreas, but the tight
distribution centered about 0.02 g for M 7-7.5 events (previously observed in Figure
4.8a) leads to a bend in the slope of the CyberShake curve. This is better illustrated
through hazard curves for individual faults in Figure 4.15, where each individual source
curve represents the hazard only considering ruptures for which the named fault
participates. Here we can see that the San Andreas controls the hazard at higher
exceedance probabilities (> 10
−3
), and the shape of the full simulated hazard curve in
Figure 4.15 follows the shape of the San Andreas curve up to 0.1 g before transitioning
to track the lower probability nearby faults which produce the largest ground motions.
While the San Andreas also controls the GMM curve up to 0.1 g in Figure 4.15b, the
ergodic model predicts a smoother and wider exceedance distribution for each fault
which results in smooth hazard curves.
Similar effects can be seen at the rest of the sites, which are plotted in
supplemental material. The most extreme hazard curves are plotted in Figures 4.16a
and 4.16b for sites SBSM and OSI, respectively. Site SBSM is located on a bed of soft
sediments in the San Bernardino basin, nearby and directly in-between both the San
96
Andreas and San Jacinto faults. Exceedance probabilities for CyberShake hazard curve
up to 0.8 g are significantly larger than those from ASK2014, and the mean CyberShake
z-score at this site is 0.99. Site OSI is located in a mountainous area with just a thin
layer of soft sediments overlaying high-velocity hard rock. In this case, the CyberShake
hazard curve is below even the 1-sigma truncated ASK2014 curve at all ground motion
levels. These two sites are examples where physics-based approaches might be able to
better capture the local site-conditions than empirical models.
Values of risk-targeted ground motion (RTGM; Luco et al., 2007), the
probabilistic design level specified in the American Society of Civil Engineers (ASCE) 7-
10 and 7-16 provisions, are listed in Table 4.4 for each site. There, we see that the
RTGMs computed from CyberShake are lower than or equal to those from ASK2014 for
eight of ten sites. SBSM is the only site in our study with a significantly larger
CyberShake RTGM, in this case 50% larger than the ASK2014 RTGM.
4.10 Discussion
The fully physics-based PSHA hazard curves presented in Figures 4.14-4.16
show the potential of future fully non-ergodic PSHA studies. Resultant hazard curves
change dramatically when individual source and path effects are uniquely characterized.
Individual sources, such as the San Andreas Fault in the Los Angeles region, can
produce characteristic ground motions with significantly less variability than predicted in
ergodic GMMs. The slope at low probabilities may also be much steeper for non-ergodic
curves; empirical GMMs extend to infinite ground motion at infinitely low probability, but
97
simulation curves are bound by a minimum possible probability, Pmin, at the level of the
largest simulated ground motion. Here, simulated curves appear more similar to GMM
curves computed with truncated distributions above Pmin. Even though CyberShake
ground motions in this study are larger on average than those from ASK2014 (Figure
4.10b), design-level ground motions (Table 4.4) for the CyberShake model are lower
than or equal to those from ASK2014 at a majority of sites. This reduction is precisely
the goal of non-ergodic PSHA: to characterize sites most likely to experience higher
than average ground motions due to local site or path effects (e.g., at site SBSM in this
model), without extending that higher hazard to all sites uniformly.
Although consistency of mean ground motion predictions in simulation-based
PSHA with data and empirical models is important, ground motion variability controls
hazard at levels of engineering interest (Anderson and Brune, 1999). A viable non-
ergodic PSHA model should contain similar total variability to ergodic models after
sampling a sufficiently diverse set of sites and sources. We outlined a method for
evaluating the total variability of such a model through z-score histograms, which
indicates that the RSQSim-CyberShake model has a comparable amount total
variability to empirical models. This is a necessary validation step, but it is not sufficient,
as the individual components of variability, both between- and within-events, are well
studied in the literature and constrained with data. We outlined a method to estimate
these components through rotation and translation of ruptures about various sites in a
3D medium. Calculated within-event variability estimates indicate that the RSQSim
model at present produces ruptures which are sufficiently heterogeneous as to produce
98
realistic ground motion fields at 𝑇 = 3 s. These calculations also revealed a potential
deficit of between-event variability, which will be a focus of future RSQSim-CyberShake
studies, as well as possible inflation of empirically estimated between-event variability
due to data uncertainties.
A benefit to the use of RSQSim ruptures directly in CyberShake is the presence
of multi-fault ruptures, which are a key component of the UCERF3 model. Even though
UCERF3 was first released in 2013, CyberShake (with the exception of this study) still
employs the prior UCERF2 model. This is largely due to limitations with the current
kinematic rupture generator, which has not yet been configured to simulate complex
multi-fault ruptures (such as the one depicted in Figure 4.4) in a fully automatic fashion
to account for moment partitioning and timing. If CyberShake is to continue to rely on
empirical ERFs, then as they increase in complexity, so must the rupture generators in
order to accommodate more complex geometries. At some point, it becomes simpler to
use a simulated extended-ERF like RSQSim directly, which eliminates the kinematic
rupture generator entirely, though that approach comes with its own set of uncertainties.
The similarity of RSQSim to UCERF3 observed by Shaw et al. (2018) suggests that an
RSQSim-derived ERF could be used by CyberShake to capture some of the UCERF2
to UCERF3 improvements without the computational challenges. We are focused on
calibrating the model in California, a region for which we have well-characterized
earthquake sources and a relatively well populated ground motion dataset. However,
because RSQSim is based on defensible seismological and physical processes, we
expect that the calibrated version will be portable to other tectonic regions as well.
99
Extensive validation, model development, and inclusion of epistemic
uncertainties are required before fully physics-based PSHA can replace current ergodic
approaches. This study illustrates initial results and validations for a single RSQSim
model, whereas a robust hazard estimate should consider epistemic uncertainties of
fault geometry and frictional parameters. Inclusion of the latter might help with the deficit
of between-event variability noted, for example by sampling a wider range of stress
drops. An ideal physics-based PSHA model would include fully-dynamic multi-cycle
(e.g., thousands to millions of years of seismicity) simulation of rupture nucleation and
resultant ground motions in a realistic 3D nonlinear medium, but such a calculation on a
large simulation domain with a complex fault system is intractable on even the largest
supercomputers. Although Richards-Dinger and Dieterich (2012) found that RSQSim
compares favorably to fully-dynamic models, it includes quasi-static approximations in
the name of computational efficiency, and other competing models with different
assumptions should be considered. Additionally (and as with regular CyberShake
studies), characterization of ground motion path effects is highly dependent on the
underlying velocity model, which is uncertain and should be represented as an
epistemic uncertainty. Even if the underlying models were thoroughly validated, current
computational methods are intractable at high frequencies (e.g. 10 Hz) for large
ensembles of ruptures required for PSHA, as higher frequencies require consideration
of nonlinearities. Still, as models evolve and computing power increases, these
deficiencies can be addressed and epistemic uncertainties reduced. Until then,
simulation-based PSHA models like RSQSim-CyberShake might inform hybrid models
100
that combine empirical datasets (or ergodic GMMs) with non-ergodic source, path and
site effects from simulations, and further exploration of the RSQSim-CyberShake model
(e.g., a hazard map for the Los Angeles region) is warranted.
4.11 Data and Resources
Hazard and empirical GMM calculations were performed with the OpenSHA
platform, available at https://github.com/opensha. The Southern California Earthquake
Center (SCEC) BroadBand Platform is available at https://github.com/SCECcode/bbp.
Plots were made with JFreeChart (http://www.jfree.org/jfreechart/) and the Generic
Mapping Tools (https://www.generic-mapping-tools.org/). More information on
CyberShake Study 15.4 is available at
https://scec.usc.edu/scecpedia/CyberShake_Study_15.4. All three-second RotD50
ground motions required to reproduce the z-scores histograms shown in Figures 4.10b
and 4.S1, and hazard curves shown in Figures 4.14, 4.16, and 4.S2 are posted at
https://github.com/kevinmilner/cybershake-analysis/tree/2020_03_09-dissertation-
supplement/study_20_2_rsqsim_4860_10x/supplemental_material. This includes both
simulated CyberShake values, and Arbrahamson et al. (2014) median and standard
deviations. All three-second RotD50 ground motions for the rotation-variability study,
required to reproduce Tables 4.2 and 4.3, are posted at
https://github.com/kevinmilner/cybershake-analysis/tree/2020_03_09-dissertation-
supplement/study_20_2_rsqsim_rot_4860_10x/supplemental_material.
101
4.12 Acknowledgements
Work here was performed within the Collaboratory for Interseismic Simulation
and Modeling (CISM) through a grant from the W. M. Keck Foundation. Significant
support also came from the U.S. Geological Center-Southern California Earthquake
Center (SCEC) Cooperative Agreement G17AC00047 and National Science Foundation
(NSF)-SCEC Cooperative Agreement EAR-1600087. SCEC BroadBand Platform (BBP)
Calculations were performed at the University of Southern California Center for High-
Performance Computing and Communications. Rate-State Earthquake Simulator
(RSQSim) catalogs were generated on the Frontera supercomputer at Texas Advanced
Computing Center (TACC) at The University of Texas at Austin, and the Blue Waters
sustained-petascale computing project, which is supported by the National Science
Foundation (awards OCI-0725070 and ACI-1238993) the State of Illinois. Blue Waters
is a joint effort of the University of Illinois at Urbana-Champaign and its National Center
for Supercomputing Applications. This work is also part of the "Improving Earthquake
Forecasting and Seismic Hazard Analysis Through Extreme-Scale Simulations"
Leadership-Class Computing Resource Allocations (PRAC) allocation supported by
NSF (award number OAC-1713792). RSQSim-CyberShake calculations were
performed at the Oak Ridge Leadership Computing Facility, which is a Department of
Energy Office of Science User Facility supported under Contract DE-AC05-00OR22725.
The authors thank Jacquelyn Gilchrist for her significant early contributions to CISM,
and Fabio Silva for his excellent documentation and help integrating RSQSim ruptures
with the SCEC BBP.
102
4.13 Tables
Table 4.1
CyberShake site name, locations, and parameters considered in this study. Site
locations are also plotted in Figure 4.9. Surface VS refers to the VS value of the
uppermost mesh grid point used in CyberShake calculations, which is taken from the
CVM-S4.26.M01 velocity model after application of a minimum VS threshold of 500 m/s.
All sites listed had modeled VS values below that threshold, and were thus set to the
threshold value. We use this surface VS value as an equivalence to the VS30 site proxy
used in empirical GMMs because the uppermost grid point has a width of more than 30
meters. Z1.0 and Z2.5 refers to the depth to the VS = 1000 m/s and 2500 m/s isosurfaces,
respectively, as reported by the CVM-S4.26.M01 model.
Site Name Latitude Longitude Surface VS
(m/s)
Z1.0 (km) Z2.5 (km)
LAF 33.86889 -118.33143 500 0.71 3.42
OSI 34.6145 -118.7235 500 0.31 0.31
PDE 34.44199 -118.58215 500 0.15 1.84
s022 34.24505 -119.18086 500 0.51 4.24
SBSM 34.064987 -117.29201 500 0.33 1.82
SMCA 34.009090 -118.48939 500 0.59 2.45
STNI 33.930880 -118.17881 500 0.88 5.57
USC 34.019200 -118.28600 500 0.58 4.10
WNGC 34.041824 -118.06530 500 0.5 3.49
WSS 34.1717 -118.64971 500 0.09 2.42
103
Table 4.2
Between-event standard deviations estimated for scenarios defined in the text. Columns
include the total simulated between-event standard deviations, 𝜏 (𝑛 ) (equation 4.9), the
range of simulated 𝜏 (𝑗 , 𝑛 ) from each site (equation 4.6), and values for individual
distances computed across all sites, 𝜏 (𝑘 , 𝑛 ) (equation 4.8).
Scenario Simulated
𝜏 (𝑛 )
Range
of𝜏 (𝑗 , 𝑛 )
𝜏 (𝑘 , 𝑛 ), 𝑘 =
20 km
𝜏 (𝑘 , 𝑛 ), 𝑘 =
50 km
𝜏 (𝑘 , 𝑛 ), 𝑘 =
100 km
M 6.6,
vertical
strike-slip
0.17 [0.15 0.21] 0.17 0.17 0.16
M 6.6,
reverse
0.17 [0.16 0.26] 0.18 0.18 0.16
M 7.2,
vertical
strike-slip
0.20 [0.18 0.22] 0.18 0.20 0.21
Table 4.3
Within-event standard deviations estimated for scenarios defined in the text. Columns
include the total simulated within-event standard deviations, 𝜙 𝑠𝑠
(𝑛 ) (equation 4.14), the
range of simulated 𝜙 𝑠𝑠
(𝑗 , 𝑛 ) from each site (equation 4.12), and values for individual
distances computed across all sites, 𝜙 𝑠𝑠
(𝑘 , 𝑛 ) (equation 4.13).
Scenario Simulated
𝜙 𝑠𝑠
(𝑛 )
Range
of𝜙 𝑠𝑠
(𝑗 , 𝑛 )
𝜙 𝑠𝑠
(𝑘 , 𝑛 ),
𝑘 = 20 km
𝜙 𝑠𝑠
(𝑘 , 𝑛 ),
𝑘 = 50 km
𝜙 𝑠𝑠
(𝑘 , 𝑛 ),
𝑘 = 100 km
M 6.6,
vertical
strike-slip
0.49 [0.40 0.65] 0.44 0.49 0.54
M 6.6,
reverse
0.51 [0.44 0.62] 0.46 0.51 0.56
M 7.2,
vertical
strike-slip
0.48 [0.38 0.63] 0.41 0.47 0.53
104
Table 4.4
Risk-targeted ground motions (RTGM; Luco et al., 2007) computed from hazard curves
from CyberShake and the Abrahamson et al. (2014; henceforth ASK2014) empirical
GMM for each of the sites listed in Table 4.1.
Site Name CyberShake RTGM ASK2014 RTGM
LAF 0.22 0.26
OSI 0.16 0.34
PDE 0.26 0.26
s022 0.24 0.28
SBSM 0.69 0.46
SMCA 0.26 0.26
STNI 0.25 0.27
USC 0.26 0.26
WNGC 0.30 0.27
WSS 0.14 0.15
105
4.14 Figures
Figure 4.1
Ground motion variability controls hazard at large intensities (e.g. >0.2 g). This
example is computed for the USC site with the Third Uniform California Earthquake
Rupture Forecast (UCERF3) model and Abrahamson et al. (2014; henceforth
ASK2014) ground motion model (GMM), modified with three different fixed total sigma
values for illustration purposes. A typical GMM total sigma value of 0.7 is plotted in
black, and reduced values of 0.5 and 0.3 are plotted in gray and light gray,
respectively. As sigma decreases, which may be possible with non-ergodic modeling,
the probability of exceeding the largest ground motions (>1 g 3-second spectral
acceleration) decreases dramatically.
106
Figure 4.2
Three-dimensional (3D) perspective view looking north of faults considered in
Southern California, highlighting an M 7.45 simulated Rate-State Earthquake
Simulator (RSQSim) rupture on the Mojave section of the San Andreas Fault. Darker
colors represent higher patches of total cumulative slip using the color scale from the
top panel of Figure 4.3. All other fault patches which did not participate in the rupture
are shown in gray.
107
Figure 4.3
Two-dimensional side view of the M 7.45 simulated RSQSim rupture from Figure 4.2.
The common x-axis is the along-strike distance of the rupture in kilometers, with zero
at the southeast end of the rupture and the maxim at the northwest end. The y-axis of
each panel is the depth of each triangular RSQSim patch from the free surface, which
is indicated with a dashed line. The top panel gives total cumulative slip with darker
colors indicating areas of greater slip. The middle panel plots the time from rupture
nucleation (in seconds) that each patch first slipped, and the bottom panel the time
that each patch last slipped in the event, contoured in 10 second intervals. The
rupture hypocenter is noted with a star.
108
Figure 4.4
3D perspective view looking north of a complex synthetic M 7.69 RSQSim rupture in
the Los Angeles basin. It nucleates on the Compton fault at the north of the rupture
surface and then spreads to the Newport-Inglewood, San Pedro Escarpment, Palos
Verdes, and other nearby faults. 14 different UCERF3 fault sections participate in
total, but the four previously named account for greater than 99% of the total seismic
moment. Darker colors represent higher patches of total cumulative slip using the
color scale from Figure 4.3.
109
Figure 4.5
Mean propagation velocity as a function of patch hypocentral distance for four
different RSQSim parameterizations, each of which incorporates a new feature over
the previous model. The base model is the catalog used in Shaw et al. (2018), plotted
with a dashed line. The first modification, plotted with a dotted line, adds a new finite
receiver patch capability to the stiffness matrix calculations. The second modification,
plotted with a dotted and dashed line, adds variable slip speed capabilities to RSQSim
with stepwise updating of sliding velocity on a patch during earthquake slip. The final
model, plotted with a solid line and used for PSHA calculations in this study, also
includes a time-delay to the static-elastic interaction.
110
Figure 4.6
Slip-time history of a single patch from the rupture depicted in Figures 4.2 and 4.3.
The top panel shows cumulative slip and the bottom instantaneous velocity, both as a
function of time from rupture nucleation. The actual slip-time history of the patch from
the RSQSim transitions file, which alternates between locked, slipping, and slip speed
updates, is plotted as black lines with circles at each transition point. The gray lines
show a discretized (at 0.1 s) slip-time history which adjusts velocities in order to
preserve total displacement (used as input to CyberShake).
111
(a) (b)
Figure 4.7
RotD50 spectra for site USC from ruptures on the Mojave section of the San Andreas
fault, computed with a one-dimensional (1D) velocity structure with Vs=500 m/s in the
Southern California Earthquake Center (SCEC) BroadBand Platform (BBP). (a)
Spectrum for the M 7.45 rupture on the Mojave section of the San Andreas Fault in
Figures 4.2 and 4.3 plotted as a thick black line. (b) Spectra for 185 different 7.0 ≤ M ≤
7.5 RSQSim ruptures on the Mojave section of the San Andreas Fault simulated at
USC plotted with thin gray lines, the mean of all 185 ruptures as a thick black line, and
the mean plus and minus one standard deviation with dashed black lines. GMM
comparisons (with plus and minus one standard deviation bounds marked with
dashed lines) are plotted with colored lines. GMM predictions are slightly different for
(b) because distributions are averaged across those predicted for each of the 185
RSQSim ruptures (rather than for a single M 7.45 rupture in (a)).
112
(a) (b)
Figure 4.8
3-second RotD50 ground motion histograms of accelerations for all 7 ≤ M ≤ 7.5 San
Andreas (Mojave) RSQSim events (black histogram) for a site at USC. (a) 1D velocity
structure in the SCEC BroadBand platform (BBP). (b) 3D velocity structure for a site at
USC. GMM predicted log-normal distributions are plotted with colored lines.
113
Figure 4.9
Map view of sites considered in this study, with site locations marked with triangles
next to their names. Fault traces from UCERF3 are drawn with thin, gray lines.
114
(a) (b)
(c)
Figure 4.10
z-score (equation 4.3) distributions (gray histograms) of all ruptures for the sites from
Figure 4.9, compared against the ASK2014 GMM. A standard normal distribution is
overlaid with a black line, and the mean z-score value is indicated with a thick dashed
vertical line. (a) RSQSim rupture ground motions simulated with BBP and a 1D
velocity structure. (b) RSQSim rupture ground motions simulated with the SCEC
CyberShake platform and a 3D velocity structure. (c) UCERF2 ruptures extended with
the Graves & Pitarka (2014) kinematic rupture generator simulated with the SCEC
CyberShake platform (Study 15.4).
115
(a) (b)
Figure 4.11
Stacked normal distributions, with the total sum represented with a thick black line and
the cumulative sum after each of 𝑁 = 50 individual distributions with thin gray lines.
(a) plots the sum where each distribution is a standard normal, and (b) where the
mean values are uniformly sampled from a standard normal and the standard
deviation given by 𝜎 = √
1
𝑁 . As 𝑁 approaches infinity, the sum of all distributions in (b)
is normally distributed and matches that from (a).
116
Figure 4.12
Map view schematic plot of rupture rotation and translation procedure to emulate
empirical records within CyberShake’s reciprocity framework. The surface of the initial
rotated rupture, in this case one which matches the M 6.6 vertical strike-slip scenario,
is depicted with a thick black line, its epicenter a star, and its scalar moment centroid
a small green square. Each rotated rupture is initially translated such that it is a fixed
distance (in this case 𝑅 𝑅𝑢𝑝 = 20 km) due north of the site (USC in this example,
depicted as a blue square in the center of the map). The rotated rupture is then
rotated multiple times both in place about its centroid and about the site, holding 𝑅 𝑅𝑢𝑝
constant. These additional rotated rupture surfaces are depicted with thin gray lines,
and the Los Angeles coastline is drawn in black at the bottom left of the map.
117
Figure 4.13
Simulated 𝜏 (𝑘 , 𝑛 ) as a function of the number of synthetic recordings per rupture
(𝑁 𝑟 𝑒𝑐
) for the M 6.6 vertical strike-slip scenario at 50 km distance. For each 𝑁 𝑟𝑒𝑐 value,
we randomly sample 𝑁 𝑟𝑒𝑐
simulated ground motions from the full set, from which we
compute downsampled 𝜏 (𝑘 , 𝑛 ). We repeat this 100 times for each 𝑁 𝑟𝑒𝑐
to estimate a
distribution of downsampled 𝜏 (𝑘 , 𝑛 ) values, the median of which is plotted with a thick
black line, and 95% range a light gray shaded region. The full simulated 𝜏 (𝑘 , 𝑛 ) value
(using all available simulations) is plotted with a horizontal dashed line, and the 𝜏
value estimated from the ASK2014 GMM is plotted with a dotted horizontal line. The
average number of recordings of 3s RotD50 spectral acceleration per event for 𝗠 ∈
[6.4,6.8] ruptures (at any distance) in the database used in the ASK2014 regressions
is indicated with a vertical dashed line, and for those ruptures within a distance range
of 𝑅 𝑅𝑢𝑝 ∈ [40,60] km is indicated with a vertical dotted line.
118
(a) (b)
Figure 4.14
RSQSim simulation hazard curves at USC. CyberShake (3D) is plotted with thick,
black lines. (a) ASK2014 GMM comparisons curves in blue, with the complete hazard
curve plotted as a thick solid line. GMM curves computed from truncated log-normal
distributions at three-, two-, and one-sigma are plotted with dashed, dotted, and
dotted and dashed lines respectively. The 1D BBP hazard curve is included in yellow,
and 95% confidence bounds assuming a binomial distribution (representing sampling
uncertainty from a finite catalog duration) on the 3D simulated curve as a gray shaded
region. (b) A zoomed in view of CyberShake hazard curves, including curves
computed with different RSQSim catalog lengths. The complete catalog is shown with
a thick black line, and subsets of the catalog, starting with the first 25,000 simulated
years in light gray, are shown in thin and increasingly dark lines with increasing
duration.
119
(a) (b)
Figure 4.15
Source contributions for RSQSim hazard curves computed at USC. Individual fault
contributions are shown with thin lines. (a) Results computed with CyberShake, with
the total CyberShake curve on top as a thick black line. (b) Results computed with the
empirical ASK2014 GMM, with the total GMM curve on top as thick blue line.
(a) (b)
Figure 4.16
Same as Figure 4.14a, except for (a) site SBSM and (b) site OSI. Site locations are
listed in Table 4.1 and plotted in Figure 4.9.
120
4.15 Supplemental Figures
(a) (b)
(c) (d)
(e) (f)
121
(g) (h)
(i) (j)
Figure 4.S1
z-score (equation 4.3) distributions (gray histograms) of all RSQSim ruptures for each
individual site (locations mapped in Figure 4.9 and listed in Table 4.1), simulated with
the Southern California Earthquake Center (SCEC) CyberShake platform and
compared against the Abrahamson et al. (2014; henceforth ASK2014) empirical
ground motion model (GMM). A standard normal is overlaid with a black line, and the
mean value is indicated with a thick dashed vertical line. (a) Site LAF. (b) Site OSI. (c)
Site PDE. (d) Site s022. (e) Site SBSM. (f) Site SMCA. (g) Site STNI. (h) Site USC. (i)
Site WNGC. (j) Site WSS.
122
(a) (b)
(c) (d)
123
(e) (f)
(g) (h)
(i) (j)
124
Figure 4.S2
RSQSim simulation hazard curves for each individual site (locations mapped in Figure
4.9 and listed in Table 4.1). CyberShake (3D) is plotted with thick, black lines. (a)
ASK2014 GMM comparisons curves in blue, with the complete hazard curve plotted
as a thick solid line. GMM curves computed from truncated log-normal distributions at
three-, two-, and one-sigma are plotted with dashed, dotted, and dotted and dashed
lines respectively. The 1D BBP hazard curve is included in yellow, and 95%
confidence bounds assuming a binomial distribution (representing sampling
uncertainty from a finite catalog duration) on the 3D simulated curve as a gray shaded
region. (a) Site LAF. (b) Site OSI. (c) Site PDE. (d) Site s022. (e) Site SBSM. (f) Site
SMCA. (g) Site STNI. (h) Site USC. (i) Site WNGC. (j) Site WSS.
125
5. Discussion and Prospectus
The field of earthquake forecasting is not new, and best practices in probabilistic
seismic hazard analysis (PSHA) have seen mostly incremental change in recent
decades, but that pace is increasing. I have outlined new approaches that will contribute
to a burgeoning field of short-term forecasting, and demonstrated how a physics-based,
deterministic, non-ergodic PSHA model can reduce aleatory variability terms that drive
long-term hazard estimates. In this chapter, I will offer my vision for the future of
simulation-based earthquake forecasting and how the work discussed previously might
contribute.
5.1 Operational Earthquake Forecasting
Now is the time to embrace and expand operational earthquake forecasting
(OEF); mature models, quality real-time data streams, and on-demand cloud computing
resources exist, and the public and stakeholders have an appetite for data. A major
challenge in earthquake forecasting is that we generally operate in a low probability
environment, which is to say, one where the chance of a damaging earthquake is very
low in an ordinary day, week, or month. We may never be able to reliably detect
precursory earthquake signals, but once a sequence is underway, we can produce
forecasts that are informative--if not actionable--for the public at large. In the case of the
2019 Ridgecrest M 6.4, the U.S. Geological Survey (USGS) publicly released a forecast
that gave a 2.3% chance of a M ≥ 7 within the following week--a factor of 460 times
126
greater than the ambient probability of such an M ≥ 7 for one week within 50 km of the
M 6.4 epicenter, which was just 0.005% according to the Third Uniform California
Earthquake Rupture Forecast (UCERF3) long-term time-dependent model. Short-term
forecasts also evolve during sequences: the probability forecasted by the USGS
increased after the M 7.1, and then decreased as time elapsed giving a 0.1% chance at
the time of this writing of another M ≥ 7 for one week beginning 29 February 2020.
The current USGS OEF effort, based on Reasenberg and Jones (1989) and
updated in Page et al. (2016), which runs automatically in response to M ≥ 5 seismicity
within the United States and its territories (e.g., Michael et al., 2019), is an important
step toward quantitatively informed decision making. I am involved as a programmer on
that project, as well as an effort to implement an efficient two-dimensional epidemic-type
aftershock sequence (ETAS) model to replace it. That ETAS model will add an
important spatial component to published aftershock forecasts, but it does not
incorporate our knowledge of active faults, or our intuition that aftershock probabilities
may be higher near those faults. The UCERF3-ETAS model was explicitly designed to
quantify that intuition into a usable forecast, and I believe that it, or a subsequent
version (e.g., a possible UCERF4), should be run operationally for California. That said,
real-time imperfect data pose a number of challenges for operationalizing UCERF3-
ETAS (as discussed in Chapter 2), but even in its pre-operationalized state, the model
provides useful products of great potential interest to the public at large and user
communities not possible with simpler models. For example, Lin’s July 2019 Los
Angeles Times article posed a question asked by many at the time: what is the
127
probability that the Ridgecrest sequence could trigger the nearby Garlock and/or San
Andreas faults, impacting millions of people? UCERF3-ETAS, despite its sensitivities
and uncertainties, is able to provide quantitative insight on that question. Additionally,
visualization of realistic scenario catalogs (Figure 2.1) may be a better tool for
presenting the range of possible outcomes than a table of probabilities from a simpler
model (such as those currently published on USGS event pages).
Ground motion forecasts may be the best tool for communicating aftershock
likelihood, as they relate directly to the level of shaking that persons, buildings, and
infrastructure will experience. In chapter 3, I demonstrated that short-term perturbations
to the conditional hypocenter distribution (CHD) can cause ground motion forecasts to
differ by up to an order of magnitude from those that don’t account for rupture directivity.
This is only possible with fault-based OEF models, highlighting another strength of
UCERF3-ETAS which is uniquely able to forecast changes to the CHD on short time
horizons. Although initial OEF ground motion models are likely to employ empirical
ground motion models (GMMs), those models should include directivity. Further
probability gains are possible through numerical simulation which better captures the
large directivity-basin coupling found by Olsen et al. (2006). To that end, I described
how the Southern California Earthquake Center (SCEC) CyberShake platform can be
used to quantify directivity-basin coupling in short-term forecasts through use of
precomputed CyberShake studies and combined with UCERF3-ETAS.
128
5.2 Non-Ergodic Simulation-Based PSHA
Despite substantial advances in seismic data collection and assimilation of that
data into increasingly complex empirical models over the past forty-plus years, Strasser
et al. (2009) observed that aleatory variability terms have remained remarkably
constant. It is unlikely that there will be any meaningful reductions in total sigma, and
thus hazard at building design levels, within an ergodic modeling framework. Chapter 4
presents a prototype non-ergodic PSHA approach, coupling CyberShake with the Rate-
State Earthquake Simulator (RSQSim), a deterministic physics-based multi-cycle
earthquake simulator. In this first comprehensive fully physics-based PSHA study,
reduced variability across repeated source-site paths yields design-level ground
motions less than or equal to those from an ergodic GMM for eight of ten sites. This
pattern of lower design-level ground motions occurred despite larger average ground
motions than from the ergodic GMM, while still maintaining an appropriate level of total
variability. However, some sites will show increased hazard in a nonergodic model; this
is most apparent in our model at site SBSM (Figure 4.16a), which exhibits large
directivity and basin amplification due to its location in soft sediments near the junction
of the San Andreas and San Jacinto faults in San Bernardino. This is precisely the goal
of non-ergodic PSHA: to accurately characterize higher hazard where it exists (e.g., at
site SBSM in this model), without extending it to all sites uniformly. Ultimately, I believe
that there is great potential in non-ergodic modelling, and that iterations of this and
future simulation-based models will play an important role in PSHA, especially as
computational capacity increases and simulations can tractably incorporate more
129
physics of earthquake processes (e.g., dynamic stress transfer and nonlinear plastic
deformation and site response).
Stewart et al. (2017) showed the potential to reduce unexplained site variance
(e.g., with 𝜙 𝑆𝑆
instead of 𝜙 ), but their empirical methodology requires expensive
instrumentation to accurately characterize local site-response, and thus is costly to
generalize across a large geographic area. More comprehensive empirical non-ergodic
simulations aim to further reduce total variance by quantifying repeated path effects
from observational data, but data are sparse for large earthquakes, and it is not clear
that path effects from large finite-fault rupture can be reliably extrapolated small
earthquakes. Simulation-based approaches, when calibrated against empirical data,
can estimate path effects for not-yet observed earthquake sources, as well as site-
response without costly instrumentation. Hybrid empirical-simulation models may be a
key step before the adoption of a purely-deterministic simulation-based model; the
USGS is currently working on one such approach to integrate basin response from
CyberShake with empirical GMMs (Moschetti et al., 2018).
Attitudes toward use of simulations in building design are improving. CyberShake
simulations are currently under consideration for tall-building design in Los Angeles
(Crouse et al., 2018), and the most recent American Society of Civil Engineers (ASCE)
7-16 code provisions explicitly permit use of ground motion simulations (ASCE, 2016).
However, those CyberShake results were calculated with the older UCERF2 model,
which might curtail adoption now that UCERF3 is the building code standard.
130
CyberShake has not yet adopted UCERF3 because the kinematic rupture generator
does not support automatic generation of the hundreds of thousands of complex multi-
fault ruptures included in UCERF3. The most recent (2016) update to the Graves and
Pitarka method adds limited support for multi-fault ruptures if the location and timing of
jumps and moment partitioning between faults are prescribed, but this is not yet an
automated process for large rupture ensembles. Multi-fault ruptures, which are
prevalent in nature (e.g., the 2016 Kaikoura M 7.8 and 2002 Denali M 7.9 earthquakes),
will surely be a significant component of the next UCERF model, and CyberShake must
address this limitation or risk falling further behind. RSQSim is a potential path forward
for CyberShake as a UCERF3 surrogate, or it can provide a rich dataset of multi-fault
rupture slip-time histories to help develop an automated kinematic approach.
Multi-fault ruptures also pose a challenge for empirical GMMs, which are
regressed against simple explanatory variables that are not appropriately descriptive of
complex ruptures. Consider the RSQSim rupture depicted in Figure 4.4: the STNI site
(location shown in Figure 4.9) in the middle of the Los Angeles basin is located on the
hanging wall of the Compton and San Pedro Escarpment faults and the footwall of the
Elysian Park fault, and is closest in 3D to the strike-slip Newport-Inglewood fault. How
should a practitioner set the faulting style, hanging wall term, and distances 𝑅 𝑋 or 𝑅 𝑌 in
an empirical GMM for such a scenario? These factors make ground motion simulation
using an extended earthquake rupture forecast (ERF) like RSQSim particularly
attractive, and more straightforward than layers of empirical models.
131
RSQSim is unlikely to replace empirical ERFs in the short term, but it will
influence the next UCERF model. The set of considered ruptures is an important input
to the UCERF rate inversion; UCERF3 employed a rule-based criteria (Milner et al.,
2013), which could be replaced by a physically-plausible set emergent from RSQSim.
RSQSim could be used as a starting model for UCERF4 inversions, or even occupy its
own rate-model branch on the logic tree. Other branching levels of the UCERF logic
tree could also be informed by RSQSim simulations, including the slip-along-rupture
models for multi-fault ruptures and poorly constrained slip-length relationships. Although
the timescale for a transition to fully-physics based PSHA is unclear, I believe that this
transition is inevitable and that simulation methods can already contribute significantly
to current practice.
132
References
Abrahamson, N. A., and Youngs, R. R., 1992. A stable algorithm for regression
analyses using the random effects model, Bulletin of the Seismological Society of
America, 82, 505–510.
Abrahamson, N., & Silva, W. (2008). Summary of the Abrahamson & Silva NGA ground-
motion relations. Earthquake spectra, 24(1), 67-97.
Abrahamson, N., Atkinson, G., Boore, D., Bozorgnia, Y., Campbell, K., Chiou, B., ... &
Youngs, R. (2008). Comparisons of the NGA ground-motion relations.
Earthquake Spectra, 24(1), 45-66.
Abrahamson, N. A., Silva, W. J., & Kamai, R. (2014). Summary of the ASK14 ground
motion relation for active crustal regions. Earthquake Spectra, 30(3), 1025-1055.
Aki, K., & Richards, P. G. (1980). Quantitative seismology: theory and methods (p. 13).
Freeman.
American Society of Civil Engineers, 2016. Minimum Design Loads for Buildings and
Other Structures, ASCE 7-16, Reston, VA.
Anderson, J. G., & Brune, J. N. (1999). Probabilistic seismic hazard analysis without the
ergodic assumption. Seismological Research Letters, 70(1), 19-28.
Andrews, D. J. (1976). Rupture velocity of plane strain shear cracks. Journal of
Geophysical Research, 81(32), 5679-5687.
Al Atik, L., and Abrahamson, N. A., 2010. Nonlinear site response effects on the
standard devia- tions of predicted ground motions, Bull. Seismol. Soc. Am. 100,
1288–1292.
Atik, L. A., Abrahamson, N., Bommer, J. J., Scherbaum, F., Cotton, F., & Kuehn, N.
(2010). The variability of ground-motion prediction models and its components.
Seismological Research Letters, 81(5), 794-801.
Al Atik, L. (2015). NGA-East: Ground-Motion Standard Deviation Models for Central and
Eastern North America. PEER Report 2015-07, Berkeley, CA.
Bayless, J. (2013). Bayless-Somerville directivity model. Chapter 3 of Pacific
Earthquake Engineering Research Center Report PEER-2013, 9.
Boore, D. M., & Atkinson, G. M. (2008). Ground-motion prediction equations for the
average horizontal component of PGA, PGV, and 5%-damped PSA at spectral
periods between 0.01 s and 10.0 s. Earthquake Spectra, 24(1), 99-138.
133
Boore, D. M. (2010). Orientation-independent, nongeometric-mean measures of seismic
intensity from two horizontal components of motion. Bulletin of the Seismological
Society of America, 100(4), 1830-1835.
Boore, D. M., Stewart, J. P., Seyhan, E., & Atkinson, G. M. (2014). NGA-West2
equations for predicting PGA, PGV, and 5% damped PSA for shallow crustal
earthquakes. Earthquake Spectra, 30(3), 1057-1085.
Böse, M., Felizardo, C., & Heaton, T. H. (2015). Finite-fault rupture detector (FinDer):
Going real-time in Californian ShakeAlert warning system. Seismological
Research Letters, 86(6), 1692-1704.
Bozorgnia, Y., Abrahamson, N. A., Atik, L. A., Ancheta, T. D., Atkinson, G. M., Baker, J.
W., ... & Darragh, R. (2014). NGA-West2 research project. Earthquake Spectra,
30(3), 973-987.
Brune, J. N. (1970). Tectonic stress and the spectra of seismic shear waves from
earthquakes. Journal of geophysical research, 75(26), 4997-5009.
Campbell, K. W., & Bozorgnia, Y. (2008). NGA ground motion model for the geometric
mean horizontal component of PGA, PGV, PGD and 5% damped linear elastic
response spectra for periods ranging from 0.01 to 10 s. Earthquake Spectra,
24(1), 139-171.
Campbell, K. W., & Bozorgnia, Y. (2014). NGA-West2 ground motion model for the
average horizontal components of PGA, PGV, and 5% damped linear
acceleration response spectra. Earthquake Spectra, 30(3), 1087-1115.
Chiou, B., Power, M., Abrahamson, N., & Roblee, C. (2006, September). An overview of
the project of next generation of ground motion attenuation models for shallow
crustal earthquakes in active tectonic regions. In Proceedings of the fifth national
seismic conference on bridges and highways, San Mateo, California.
Chiou, B. S. J. (2008). Chiou-Youngs NGA ground motion relations for the geometric
mean horizontal component of peak and spectral ground motion parameters.
Earthquake Spectra, 24(1), 173-215.
Chiou, B. S. J., & Youngs, R. R. (2014). Update of the Chiou and Youngs NGA model
for the average horizontal component of peak ground motion and response
spectra. Earthquake Spectra, 30(3), 1117-1153.
Cornell, C. A. (1968). Engineering seismic risk analysis. Bulletin of the seismological
society of America, 58(5), 1583-1606.
Crouse, C. B., Jordan, T. H., Milner, K. R., Goulet, C. A., Callaghan, S., & Graves, R.
W. (2018, June). Site‐specific MCER response spectra for Los Angeles region
based on 3‐D numerical simulations and the NGA West2 equations. In
134
Proceedings of the 11th National Conference in Earthquake Engineering. Los
Angeles, CA: Earthquake Engineering Research Institute.
Dieterich, J. H. (1992). Earthquake nucleation on faults with rate-and state-dependent
strength. Tectonophysics, 211(1-4), 115-134.
Dieterich, J. H., & Richards-Dinger, K. B. (2010). Earthquake recurrence in simulated
fault systems. In Seismogenesis and Earthquake Forecasting: The Frank Evison
Volume II (pp. 233-250). Springer, Basel.
Field, E. H. (2000). A modified ground-motion attenuation relationship for southern
California that accounts for detailed site classification and a basin-depth effect.
Bulletin of the Seismological Society of America, 90(6B), S209-S221.
Field, E. H., Jordan, T. H., & Cornell, C. A. (2003). OpenSHA: A developing community-
modeling environment for seismic hazard analysis. Seismological Research
Letters, 74(4), 406-419.
Field, E. H. (2007). A summary of previous working groups on California earthquake
probabilities. Bulletin of the Seismological Society of America, 97(4), 1033-1053.
Field, E. H., Dawson, T. E., Felzer, K. R., Frankel, A. D., Gupta, V., Jordan, T. H., ... &
Wills, C. J. (2009). Uniform California earthquake rupture forecast, version 2
(UCERF 2). Bulletin of the Seismological Society of America, 99(4), 2053-2107.
Field, E. H., Arrowsmith, R. J., Biasi, G. P., Bird, P., Dawson, T. E., Felzer, K. R., ... &
Michael, A. J. (2014). Uniform California earthquake rupture forecast, version 3
(UCERF3)—The time‐independent model. Bulletin of the Seismological Society
of America, 104(3), 1122-1180.
Field, E. H., Biasi, G. P., Bird, P., Dawson, T. E., Felzer, K. R., Jackson, D. D., ... &
Milner, K. R. (2015). Long‐term time‐dependent probabilities for the third Uniform
California Earthquake Rupture Forecast (UCERF3). Bulletin of the Seismological
Society of America, 105(2A), 511-543.
Field, E. H., Jordan, T. H., Jones, L. M., Michael, A. J., Blanpied, M. L., & Workshop
Participants. (2016). The potential uses of operational earthquake forecasting.
Seismological research letters, 87(2A), 313-322.
Field, E. H., Milner, K. R., Hardebeck, J. L., Page, M. T., van der Elst, N., Jordan, T. H.,
... & Werner, M. J. (2017). A spatiotemporal clustering model for the third Uniform
California Earthquake Rupture Forecast (UCERF3‐ETAS): Toward an
operational earthquake forecast. Bulletin of the Seismological Society of
America, 107(3), 1049-1081.
Field, E. H., Jordan, T. H., Page, M. T., Milner, K. R., Shaw, B. E., Dawson, T. E., ... &
Weldon, R. J. (2017). A synoptic view of the third Uniform California Earthquake
Rupture Forecast (UCERF3). Seismological Research Letters, 88(5), 1259-1267.
135
Field, E. H., & Milner, K. R. (2018). Candidate products for operational earthquake
forecasting illustrated using the HayWired planning scenario, including one very
quick (and not‐so‐dirty) hazard‐map option. Seismological Research Letters,
89(4), 1420-1434.
Geller, R. J. (1976). Scaling relations for earthquake source parameters and
magnitudes. Bulletin of the Seismological Society of America, 66(5), 1501-1523.
Goulet, C. A., Abrahamson, N. A., Somerville, P. G., & Wooddell, K. E. (2014). The
SCEC broadband platform validation exercise: Methodology for code validation in
the context of seismic‐hazard analyses. Seismological Research Letters, 86(1),
17-26.
Graves, R. W., & Pitarka, A. (2010). Broadband ground-motion simulation using a
hybrid approach. Bulletin of the Seismological Society of America, 100(5A),
2095-2123.
Graves, R., Jordan, T. H., Callaghan, S., Deelman, E., Field, E., Juve, G., ... & Okaya,
D. (2011). CyberShake: A physics-based seismic hazard model for southern
California. Pure and Applied Geophysics, 168(3-4), 367-381.
Graves, R. (2014). Standard Rupture Format, Version 2.0. Retrieved December 2,
2019, from http://equake-rc.info/static/paper/SRF-Description-Graves_2.0.pdf.
Graves, R., & Pitarka, A. (2014). Refinements to the Graves and Pitarka (2010)
broadband ground‐motion simulation method. Seismological Research Letters,
86(1), 75-80.
Graves, R., & Pitarka, A. (2016). Kinematic ground‐motion simulations on rough faults
including effects of 3D stochastic velocity perturbations. Bulletin of the
Seismological Society of America, 106(5), 2136-2153.
Gutenberg, B., & Richter, C. F. (1944). Frequency of earthquakes in California. Bulletin
of the Seismological Society of America, 34(4), 185-188.
Hardebeck, J. L. (2013). Appendix S: Constraining epidemic type aftershock sequence
(ETAS) parameters from the Uniform California Earthquake Rupture Forecast,
Version 3 catalog and validating the ETAS model for magnitude 6.5 or greater
earthquakes. U.S. Geol. Surv. Open‐File Rept. 2013‐1165‐S, and California
Geol. Surv. Special Rept. 228‐S.
Hardebeck, J. L., Llenos, A. L., Michael, A. J., Page, M. T., & Van Der Elst, N. (2018).
Updated California Aftershock Parameters. Seismological Research Letters,
90(1), 262-270.
Haskell, N. A. (1969). Elastic displacements in the near-field of a propagating fault.
Bulletin of the Seismological Society of America, 59(2), 865-908.
136
Idriss, I. M. (2008). An NGA empirical model for estimating the horizontal spectral
values generated by shallow crustal earthquakes. Earthquake Spectra, 24(1),
217-242.
Jordan, T. H., & Jones, L. M. (2010). Operational earthquake forecasting: Some
thoughts on why and how. Seismological Research Letters, 81(4), 571-574.
Jordan, T. H., Chen, Y. T., Gasparini, P., Madariaga, R., Main, I., Marzocchi, W., ... &
Zschau, J. (2011). Operational earthquake forecasting. State of knowledge and
guidelines for utilization. Annals of Geophysics, 54(4).
Jordan, T. H., Callaghan, S., Graves, R. W., Wang, F., Milner, K. R., Goulet, C. A., ... &
Vahi, K. (2018, June). CyberShake models of seismic hazards in Southern and
Central California. In Proceedings of the 11th National Conference in Earthquake
Engineering. Los Angeles, CA: Earthquake Engineering Research Institute.
Landwehr, N., Kuehn, N. M., Scheffer, T., & Abrahamson, N. (2016). A nonergodic
ground‐motion model for California with spatially varying coefficients. Bulletin of
the Seismological Society of America, 106(6), 2574-2583.
Lee, E. J., Chen, P., Jordan, T. H., Maechling, P. B., Denolle, M. A., & Beroza, G. C.
(2014). Full‐3‐D tomography for crustal structure in southern California based on
the scattering‐integral and the adjoint‐wavefield methods. Journal of Geophysical
Research: Solid Earth, 119(8), 6421-6451.
Lin, P. S., Chiou, B., Abrahamson, N., Walling, M., Lee, C. T., & Cheng, C. T. (2011).
Repeatable source, site, and path effects on the standard deviation for empirical
ground-motion prediction models. Bulletin of the Seismological Society of
America, 101(5), 2281-2295.
Lin II, R. (2019, July 25). San Andreas fault is a 730-mile monster. Ridgecrest
earthquake was a tiny taste of the possible destruction. Los Angeles Times,
available at http://www.latimes.com, accessed September 2019.
Luco, N., Ellingwood, B. R., Hamburger, R. O., Hooper, J. D., Kimball, J. K., & Kircher,
C. A. (2007). Risk-targeted versus current seismic design maps for the
conterminous United States.
Maechling, P. J., Silva, F., Callaghan, S., & Jordan, T. H. (2014). SCEC Broadband
Platform: System architecture and software implementation. Seismological
Research Letters, 86(1), 27-38.
Michael, A.J., McBride, S.K., Hardebeck, J.L., Barall, M., Martinez, E., Page, M.T., van
der Elst, N., Field, E.H., Milner, K.R. and Wein, A.M., 2020. Statistical
Seismology and Communication of the USGS Operational Aftershock Forecasts
for the 30 November 2018 M w 7.1 Anchorage, Alaska, Earthquake.
Seismological Research Letters, 91(1), pp.153-173.
137
Milner, K. R., Page, M. T., Field, E. H., Parsons, T., Biasi, G. P., & Shaw, B. E. (2013).
Appendix T—Defining the inversion rupture set using plausibility filters. US Geol.
Surv. Open‐File Rept. 2013‐1165.
Milner, K. R., Field, E. H., Savran, W. H., Page, M. T., & Jordan, T. H. (2020).
Operational Earthquake Forecasting during the 2019 Ridgecrest, California,
Earthquake Sequence with the UCERF3‐ETAS Model. Seismological Research
Letters.
Moschetti, M. P., Chang, S., Crouse, C. B., Frankel, A., Graves, R., Puangnak, H.,
Luco, N., Goulet, C. A., Rezaeian, S., Shumway, A., Powers, P. M., Petersen, M.
D., Callaghan, S., Jordan, T. H., & Milner, K. R. (2018, June). The science,
engineering applications, and policy implications of simulation‐based PSHA. In
Proceedings of the 11th National Conference in Earthquake Engineering
(11NCEE), June (pp. 25-19).
November 2019 Top500 List. (2019, November 18). Retrieved from
https://www.top500.org/lists/2019/11/
Ogata, Y. (1988). Statistical models for earthquake occurrences and residual analysis
for point processes. Journal of the American Statistical association, 83(401), 9-
27.
Olsen, K. B., Day, S. M., Minster, J. B., Cui, Y., Chourasia, A., Faerman, M., ... &
Jordan, T. (2006). Strong shaking in Los Angeles expected from southern San
Andreas earthquake. Geophysical Research Letters, 33(7).
Page, M. T., Field, E. H., Milner, K. R., & Powers, P. M. (2014). The UCERF3 Grand
Inversion: Solving for the Long‐Term Rate of Ruptures in a Fault SystemThe
UCERF3 Grand Inversion: Solving for the Long‐Term Rate of Ruptures in a Fault
System. Bulletin of the Seismological Society of America, 104(3), 1181-1204.
Powers, P. M, & Field, E. H. (2013). Appendix O: Gridded Seismicity Sources, U.S.
Geol. Surv. Open‐File Rept. 2013‐1165‐O, and California Geol. Surv. Special
Rept. 228‐O.
Reasenberg, P. A., & Jones, L. M. (1989). Earthquake hazard after a mainshock in
California. Science, 243(4895), 1173-1176.
Reid, H. F. (1910). The mechanics of the earthquake. The California Earthquake of April
18, 1906, Report of the State Earthquake Investigation Commission.
Richards‐Dinger, K., & Dieterich, J. H. (2012). RSQSim earthquake simulator.
Seismological Research Letters, 83(6), 983-990.
Shaw, B. E., Milner, K. R., Field, E. H., Richards-Dinger, K., Gilchrist, J. J., Dieterich, J.
H., & Jordan, T. H. (2018). A physics-based earthquake simulator replicates
seismic hazard statistics across California. Science advances, 4(8), eaau0688.
138
Shaw, B. E. (2019). Beyond Backslip: Improvement of Earthquake Simulators from New
Hybrid Loading Conditions. Bulletin of the Seismological Society of America,
109(6), 2159-2167.
Small, P., Gill, D. M., Maechling, P. J., Taborda, R., Callaghan, S., Jordan, T. H., …
Goulet, C. A. (2017). The SCEC unified community velocity model software
framework. Seismological Research Letters, 88(6), 1539–1552.
https://doi.org/10.1785/0220170082
Spudich, P., & Chiou, B. S. (2008). Directivity in NGA earthquake ground motions:
Analysis using isochrone theory. Earthquake Spectra, 24(1), 279-298.
Spudich, P., Rowshandel, B., Shahi, S. K., Baker, J. W., & Chiou, B. S. J. (2014).
Comparison of NGA-West2 directivity models. Earthquake Spectra, 30(3), 1199-
1221.
Stewart, J. P., Afshari, K., & Goulet, C. A. (2017). Non-ergodic site response in seismic
hazard analysis. Earthquake Spectra, 33(4), 1385-1414.
Strasser, F. O., Abrahamson, N. A., & Bommer, J. J. (2009). Sigma: Issues, insights,
and challenges. Seismological Research Letters, 80(1), 40-56.
Wald, D. J., Quitoriano, V., Heaton, T. H., Kanamori, H., Scrivner, C. W., & Worden, C.
B. (1999). TriNet “ShakeMaps”: Rapid generation of peak ground motion and
intensity maps for earthquakes in southern California. Earthquake Spectra, 15(3),
537-555.
Wang, F., & Jordan, T. H. (2014). Comparison of probabilistic seismic‐hazard models
using averaging‐based factorization. Bulletin of the Seismological Society of
America, 104(3), 1230-1257.
Wills, C. J., Gutierrez, C. I., Perez, F. G., & Branum, D. M. (2015). A Next Generation
VS30 Map for California Based on Geology and TopographyA Next Generation
VS30 Map for California Based on Geology and Topography. Bulletin of the
Seismological Society of America, 105(6), 3083-3091.
Working Group on California Earthquake Probabilities (wgcep) (1988). Probabilities of
large earthquakes occurring in California on the San Andreas fault, U.S. Geol.
Surv. Open-File Rept. 88-398 , p. 62.
Working Group on California Earthquake Probabilities. (1995). Seismic hazards in
southern California: probable earthquakes, 1994 to 2024. Bulletin of the
Seismological Society of America, 85(2), 379-439.
Working Group on California Earthquake Probabilities (wgcep) (2003). Earthquake
probabilities in the San Francisco Bay Region: 2002– 2031, U.S. Geol. Surv.
Open-File Rept. 03-214.
139
Zhao, L., Chen, P., & Jordan, T. H. (2006). Strain Green’s Tensors, Reciprocity, and
Their Applications to Seismic Source and Structure Studies. Bulletin of the
Seismological Society of America, 96(5), 1753-1763. doi: 10.1785/0120050253.
Abstract (if available)
Abstract
Ground motions radiated by earthquakes pose a significant hazard to human life and the built environment. I present three papers which contribute to the goal of more accurately forecasting earthquake hazard on multiple timescales. The first gives results and lessons learned from short-term operational earthquake forecasting during the 2019 Ridgecrest, CA, earthquake sequence. The second demonstrates how consideration of short-term perturbations to the conditional hypocenter distribution can achieve a ground motion probability gain of an order of magnitude relative to short-term models which ignore rupture directivity. The third introduces a first-of-its-kind probabilistic seismic hazard model that is constructed with fully-deterministic physics-based models. That model is non-ergodic and forecasts lower long-term hazard for some sites in the Los Angeles region than empirical ergodic models due to its characterization of repeated source, site, and path effects, while still containing a sufficient amount of total ground motion variability.
Linked assets
University of Southern California Dissertations and Theses
Conceptually similar
PDF
Integration and validation of deterministic earthquake simulations in probabilistic seismic hazard analysis
PDF
Rupture synchronicity in complex fault systems
PDF
Forecasting directivity in large earthquakes in terms of the conditional hypocenter distribution
PDF
Preparing for the next major southern California earthquake: utilizing HAZUS with soils maps and ShakeMaps to predict regional bridge damage and closures
PDF
The impact of faulting and folding on earthquake cycles in collisional orogens
PDF
Stress-strain characterization of complex seismic sources
Asset Metadata
Creator
Milner, Kevin Ross
(author)
Core Title
Advancing probabilistic seismic hazard analysis through short-term time-dependent forecasts and deterministic simulations
School
College of Letters, Arts and Sciences
Degree
Doctor of Philosophy
Degree Program
Geological Sciences
Publication Date
04/08/2020
Defense Date
03/24/2020
Publisher
University of Southern California
(original),
University of Southern California. Libraries
(digital)
Tag
conditional hypocenter distribution,CyberShake,earthquake,earthquake forecasting,earthquake simulators,Earthquakes,OAI-PMH Harvest,OEF,operational earthquake forecasting,probabilistic seismic hazard analysis,PSHA,Ridgecrest,RSQSim,UCERF3
Language
English
Contributor
Electronically uploaded by the author
(provenance)
Advisor
Jordan, Thomas (
committee chair
), Barbot, Sylvain (
committee member
), Kesselman, Carl (
committee member
), Vidale, John (
committee member
)
Creator Email
kmilner@usc.edu,steel@kevinmilner.net
Permanent Link (DOI)
https://doi.org/10.25549/usctheses-c89-280775
Unique identifier
UC11674763
Identifier
etd-MilnerKevi-8255.pdf (filename),usctheses-c89-280775 (legacy record id)
Legacy Identifier
etd-MilnerKevi-8255.pdf
Dmrecord
280775
Document Type
Dissertation
Rights
Milner, Kevin Ross
Type
texts
Source
University of Southern California
(contributing entity),
University of Southern California Dissertations and Theses
(collection)
Access Conditions
The author retains rights to his/her dissertation, thesis or other graduate work according to U.S. copyright law. Electronic access is being provided by the USC Libraries in agreement with the a...
Repository Name
University of Southern California Digital Library
Repository Location
USC Digital Library, University of Southern California, University Park Campus MC 2810, 3434 South Grand Avenue, 2nd Floor, Los Angeles, California 90089-2810, USA
Tags
conditional hypocenter distribution
CyberShake
earthquake forecasting
earthquake simulators
OEF
operational earthquake forecasting
probabilistic seismic hazard analysis
PSHA
Ridgecrest
RSQSim
UCERF3