Book of abstracts
12th International Conference on Statistical
Seismology (STATSEI 12)
Cargèse, 17th to 21st of October 2022

Physics-based estimates of the maximum magnitudeof induced earthquakes in the Groningen gas field
Jean Paul Ampuero∗1 , Huihui Weng2 , and Loes Buijze3
1Géoazur – Université Côte d’Azur, IRD, CNRS, Observatoire de la Côte d’Azur, Géoazur – France
2
Geoazur – IRD, UCA – France3
Utrecht University – Netherlands

Abstract
Induced seismicity in the Groningen gas field, The Netherlands, has led to public concern
and comprehensive investigations to assess future seismic hazard. One of the main challenges
is to estimate the maximum possible earthquake magnitude (Mmax) that could be induced
by gas extraction. Here, we estimate Mmax in the Groningen gas field through a novel
physics-based approach, based on the dynamic fracture mechanics theory, constrained byfield and laboratory observations. The method provides estimates of the along-dip rupture
width, along-strike rupture length, and slip of potential maximum earthquakes for each of
the more than 1000 mapped faults in the reservoir, at any given time in the production
and post-production histories. The rupture widths are constrained by the along-dip twodimensional fracture mechanics theory, and the rupture lengths by the three-dimensional
fracture mechanics theory of dynamic ruptures with large aspect ratio. These two theoreticalelements are extended, combined and validated in this study. The estimated Mmax is timedependent due to the gas extraction history. It saturates at a value that depends on thefriction parameters of fault rocks and on the pre-production tectonic stresses. Applying thephysics-based rupture model, while accounting for plausible ranges of key parameters, yields
an estimate of Mmax in the Groningen gas field. While there are considerable uncertainties insome parameters, such as dynamic friction coefficient and pre-production stresses below thereservoir, these may be reduced through future focused efforts. The proposed approach is apractical step toward physics-based seismic hazard assessment for other gas fields worldwide.
Ongoing extensions of the model to include ruptures across multiple faults will be discussed.

∗

Speakersciencesconf.org:statsei12:418104

The predictive skill of global and regional earthquake
forecasting models for California, New Zealand, and
Italy
José Bayona∗1 , William Savran2 , and Maximilian Werner1
1

2School of Earth Sciences, University of Bristol – United Kingdom
Southern California Earthquake Center, University of Southern California – United States

Abstract
The increasing availability and quality of geophysical datasets, including earthquake catalogs, geological records, and interseismic strain rates, has enabled the creation of global
and regional earthquake forecasting models that can underpin probabilistic seismic hazard
analyses. Regional models provide detailed seismicity forecasts given the high resolution of
their input data, while global models offer greater testability of the large damaging earthquakes due to the more frequent seismic activity observed around the world. Quantitative
comparisons between global and regional models can be useful to assess the performance of
global models at regional scales, identify regional discrepancies that require further investigation, and better inform seismic hazard programs. Here, we prospectively evaluate theability of the Global Earthquake Activity Rate (GEAR1) model, developed by Bird et al.
(2015), and nineteen regionally calibrated models to forecast earthquakes in California, New
Zealand, and Italy from 2014 through 2021. We project GEAR1 onto these testing regions,previously defined by the Collaboratory for the Study of Earthquake Predictability (CSEP),
and scale its M5.95+ rates to a magnitude threshold of 4.95 assuming a b-value of unity. We
then use comparative and consistency likelihood tests implemented in the CSEP’s pyCSEP
Python toolkit to contrast the performance of each regional model with that of GEAR1 andto evaluate the consistency between forecasts and observations. In addition, we compute Gini
coefficients for each forecasting model to measure the absolute spatial localization of theirpredicted earthquake rates and to interpret the spatial test results. The testing dataset contains thirty-eight M4.95+ earthquakes in California, including the 2019 Mw 7.1 Ridgecrest
earthquake, forty-seven M4.95+ earthquakes in New Zealand, including the 2016 Mw 7.8
Kaikōura quake, and eleven target events in Italy, including the 2016 Mw 6.6 Norcia earthquake. Comparative test results show that, in California, a model that adaptively smooths
the locations of small earthquakes is the most informative, while GEAR1 outperforms allother regional models. In Italy, only amodel based on relative seismic intensities and twomodels based on adaptively smoothed seismicity can be considered statistically more informative than GEAR1. In New Zealand, GEAR1 is the most informative model, obtaining
information gains per earthquake of about 0.5 over a fault-based and a smoothed seismicity
model, respectively. These results provide evidence that GEAR1 is a robust first-order model
to forecast M4.95+ seismicity in California, New Zealand, and Italy. When analyzing the
spatial dimension of the models, we find that those that provide forecasts that are either
too dispersed (low Gini coefficients) or too localized (high Gini coefficients) are the worst∗

Speakersciencesconf.org:statsei12:417232

performing models. Thus, we recommend applying an ”intermediate” smoothing procedure,e.g. adaptive smoothing, to improve the forecasting capacity of these models. Finally, we
note thatthe New Zealand models obtain, overall, the highest spatial likelihood scores perearthquake, followed by the models for California and Italy, which could indicate a possible
link between the tectonic setting and the spatial distribution of earthquakes.

Keywords: Statistical seismology, prospective model evaluation, global and regional models

Monitoring the criticality of the crust over a decadenear Ridgecrest, CA, with temporal clustering andtidal triggering of microearthquakes
Eric Beauce∗1 , Benjamin Holtzman1 , and Felix Waldhauser1
1

Lamont-Doherty Earth Observatory, Columbia University – United States

Abstract
Earthquakes occur in the brittle part of the earth’s crust to release the strain produced
by tectonic loading over long time scales. Observing that earthquake sizes are power-law
distributed, that relatively small stress perturbations (kPa) seem enough to trigger earthquakes, or that earthquakes temporal clustering follows fractal statistics, it appears that thecrust is at or near a critical state. In this study, we built a decade-long catalog of microseismicity before the M7.1 2019-07-06 Ridgecrest earthquake. We systematically study thespatial and temporal evolution of the strength of temporal clustering and the sensitivity ofseismicity rates to tidal stresses. We discuss how these observables relate to the mechanicalstate of the fault zone, and thus we produce new observations of the mechanical evolutionof the crust over long time scales and prior to a large earthquake.

∗

Speakersciencesconf.org:statsei12:411433

Change point analysis to fit time variable parametersof the ETAS model
Amel Benali∗†1 , Jiancang Zhuang2 , and Abdelhak Talbi1
1

Centre de Recherche en Astronomie, Astrophysique et Géophysique – Algeria
2
Institute of Statistical Mathematics – Japan

Abstract
This study uses change points based analysis to detect and fit potential temporal variation observed in the parametrisation of the stationary Epidemic Type Aftershock Sequence(ETAS) model. A series of fit-test steps are applied to account for possible time-varying
background seismicity rate and earthquakes triggering parameters. Specifically, the standard ETAS model is fitted to a series of time periods identified by change points at one
or both of their interval limits, while testing is applied to assess the goodness of fit between different ETAS models obtained by combining different fits at adjacent subperiods.
For the Central Apennines (Italy) seismicity, two change points are detected in the time
period 2005-2017. They depict starts of apparent quiescence occurred before the Amatrice
M6 (2016) earthquake and activation observed around Campotosto M5.5 (2017) earthquake,
respectively. Subsequently, a three-stage ETAS model is proposed to fit seismicity, withtime varying parameters. A new time variable ETAS model fitting procedure is proposed
and tested with an increasing background seismicity rate and significant temporal variationof the earthquakes triggering parameters in the three detected adjacent subperiods.

Keywords: ETAS model. Change point analysis. Central Apennines.

∗†

Speaker
Corresponding author: benali.amel1992@gmail.com

sciencesconf.org:statsei12:416140

A Bayesian model selection approach to estimate the
distribution of earthquake magnitudes (including
below completeness)
Thomas Bodin∗1 , Stephanie Durand2 , Blandine Gardonio1 , and David Marsan3
1

Laboratoire de Géologie de Lyon - Terre, Planètes, Environnement [Lyon] (LGL-TPE) – Centre
National de la Recherche Scientifique : UMR5276 – 2, rue Raphaël Dubois69622 Villeurbanne Cedex,
France2
Laboratoire de Géologie de Lyon - Terre, Planètes, Environnement [Lyon] (LGL-TPE) – Centre
National de la Recherche Scientifique : UMR5276 – 2, rue Raphaël Dubois69622 Villeurbanne Cedex,
France3
Institut des Sciences de la Terre (ISTerre) (ISTerre) – Université Savoie Mont Blanc – France

Abstract
The statistical distribution of seismicity is described by the so-called Gutenberg-Richterlaw where the distribution of magnitudes p(M) follows an exponential law p(M) = b exp(bM). The slope of this law, called ”b-value”, is used to compute seismic hazard maps. A
well-known problem is that small events are often not detected, and the actual GutenbergRichter law thus needs to be corrected for small events. This can be done by multiplying the
exponential law by the probability q(M) of detecting an earthquake of magnitude M. In thisway, the distribution of observed magnitudes writes p(M)= q(M)*b exp(-bM). Characterizing the detection function q(M) is crucial as it allows to estimate the number of earthquakesthat are missed at small magnitudes. This is particularly useful when analyzing aftershock
sequences at early times, when the detection function threshold varies quickly.
Here we propose a Bayesian approach to jointly estimated all parameters that define p(M),i.e. the b-value as well as the parameters of the detection function q(M). Different detection
functions, described by different number of parameters, are tested. Of course, a better fit
is achieved as more parameters are used to describe q(M). We show how a Bayesian model
selection approach can be used to compare different detection functions with different levelsof complexity, to properly estimate q(M) and the b-value, together with associated uncertainties.

Keywords: Gutenberg Richter Law, Bayesian inference, b value, microseismicity∗

Speakersciencesconf.org:statsei12:417289

Towards Network Resilience: a collective effort by
the expert community to define ”good” seismic data
quality using New Zealand specific ambient noisemodelsConrad Burton∗1 , Sepideh J Rastin†‡1 , Sam Taylor-Offord§1 , Mark Chadwick1 , Tim
Mcdougall1 , and Gill Jolly1
1

GNS Science – 1 Fairway drive, Avalon, P.O. Box 30–368, Lower Hutt, New Zealand

Abstract
New Zealand National Seismic Network (NZNSN) totaling some 53 stations is managed
by GeoNet and underpins earthquake, volcano, and tsunami monitoring. GeoNet provides
free and open data in support of national and international research into geophysical research as well as natural hazard risk analysis and loss modeling. Data quality affects every
application/purpose in these domains and yet little is known about whether it meets stakeholder and data user needs. Having a definition of ”good” in a common language and indata-terms offers a practical, quantitative reference for end-users and network operators to
support communication and decision-making. Funded by GNS Science, the ”Towards Network Resilience” project has aimed to investigate this fundamental matter at a national
level in three stages: the first stage was to establish a ”feedback style” consultation model
to define what constitutes ”good” seismic data quality collaboratively with the New Zealandexpert end-user communities. The second stage was to develop New Zealand specific reference models based on the community consensus in data terms and formalize it as a common
standard for seismic data quality in New Zealand. This reference took the form of ambient
noise models like those produced by S. J. Rastin (2012). The final stage was to explore how
well historic seismic station data compares against our reference and how such comparisonscould inform key performance indicators for network management and stakeholder relations.
Ultimately we hope to integrate this work with geological hazard monitoring operations in
New Zealand to build the capability of data-driven decision-making, empower network managers with awareness and understanding of station health, and to engineer processes for ahigh degree of network resilience.

Keywords: New Zealand, Ambient Noise Models, GeoNet, Seismic Networks, Data Quality∗

Corresponding author: c.burton@gns.cri.nz
Speaker‡
Corresponding author: s.rastin@gns.cri.nz§
Corresponding author: s.taylor-offord@gns.cri.nz†

sciencesconf.org:statsei12:417033

Spatio-temporal evolution of foreshock sequences inthe damage zone
Camilla Cattania∗1
1

Massachusetts Institute of Technology – United States

Abstract
Fault zones comprise multiple fractures which can generate seismicity in response to mainshocks, tectonic loading, or other stress perturbations, and in turn modify the stress state on
the main fault. The physical processes controlling seismic swarms, and their implications forthe time-dependent seismic hazard posed by the mainshock fault, are not fully understood.
Here we study the precursory phase of large earthquakes in a fault zone comprised of arough (fractal) fault surrounded by damage, represented as a collection of smaller faults
with a power-law decay of density with distance from the main fault. We model seismic
cycles assuming uniform velocity-weakening rate-state friction, and loading by a uniform
far-field stressing rate.
We find that during the interseismic period between large ruptures on the main fault, seismicity takes place predominantly in the damage zone, and it is clustered in episodic swarms.
Seismicy rates increase prior to the mainshock and exhibit a clear migration towards its
hypocenter, both along strike and across the fault zone; in other words, seismicity becomes
more localized in the final phase of the seismic cycle. This coalescence is controlled by two
processes: 1. recovery of negative stress changes imparted by the previous mainshock, which
are largest in the near field; 2. a triggering cascade in which each successive event is more
likely to occur closer to the main fault, due to the power-law increase in fracture density
towards the main fault.
In addition to foreshocks, simulations exhibit aseismic slip in subparallel secondary faults.
We find that the mainshock is triggered by stress changes from both foreshocks and aseismic
slip in the damage zone, and occurs earlier than it would have on a planar fault without
damage.

∗Speakersciencesconf.org:statsei12:421706

Potential correlation between intraslab
intermediate-depth and shallow earthquakes in Japanand Chile
Audrey Chouli∗1 , David Marsan2 , Sophie Giffard Roisin2 , Michel Bouchon2 , and Anne
Socquet2
1

Institut des Sciences de la Terre (ISTerre) – Université Grenoble Alpes, Université Savoie Mont Blanc,
CNRS, IRD, IFSTTAR, ISTerre, 38000 Grenoble, France – France
2
Institut des Sciences de la Terre (ISTerre) – Université Grenoble Alpes, Université Savoie Mont Blanc,
CNRS, IRD, IFSTTAR, ISTerre, 38000 Grenoble, France – France

Abstract
An increase of both shallow and intraslab intermediate-depth seismicity has been observed days to years before some great subduction earthquakes, as before Tohoku-oki (Mw9.0, 2011), Maule (Mw 8.8, 2010) or Iquique (Mw 8.2, 2014) earthquakes (Bouchon et al.,2016, Jara et al., 2017). These observations suggest that a link exists between these deep andshallow foreshocks, but it is still poorly understood and not characterized in a systematicmanner. Some studies have attempted to address this lack of systematic characterization byusing a statistical approach (Delbridge et al., 2017).
The aim of this study is to systematically and statistically identify and characterize the potential correlations between deep and shallow seismicity. We want to assess whether or not
such interactions exist. If they exist, we plan to characterize when and where they occur, atwhat frequency, their characteristic duration, and with what spatial pattern.
For this purpose, we develop a statistical method to assess the relevance of deep-shallow
interactions, that allows to identify statistically significant correlations between deep andshallow seismicity. We focused on the seismicity of the Japan trench and the Chile trench
subduction zones, during the decades prior to the Tohoku-oki and Iquique earthquakes respectively, because deep-shallow interactions were identified there, and because we can test
the events picked by our method against the correlations highlighted in published papers(Bouchon et al., 2016). The correlation values between the deep and shallow events arecalculated on various sliding-windows with durations from days to months. These correlation values are then compared to the ones obtained using parts of catalogs far in time andtherefore uncorrelated, and the significance of the correlation is calculated.
Some windows show a strong correlation. The dependence of our results to different parameters, such as the completeness magnitude, the length of the window or the smoothing
are evaluated. The spatio-temporal analysis of the seismicity on maps for these windows
is also explored. While the results are still preliminary, we believe that this method has
the potential to systematically and quantitatively assess the current presumptions on thelink between deep and shallow seismicity, that would lead to a better understanding of themechanisms leading to megathrust earthquakes.
∗

Speakersciencesconf.org:statsei12:417184

Regressions of moment magnitude on local
magnitude in the New Zealand earthquake catalogue
for the revision of the National Seismic Hazard model
Annemarie Christophersen∗1 , Sandra Bourguignon1 , David Rhoades1 , Trevor Allen2 ,Jerome Salichon1 , John Ristau1 , Chris Rollins1 , and Matt Gerstenberger1
1

2

GNS Science – New ZealandGeoscience Australia – Australia

Abstract
Earthquake catalogues are by nature heterogeneous over time due to changes in network
configuration and seismic processing. Statistical studies, including probabilistic seismic hazard analysis (PSHA), require homogeneous data to study changes in seismicity rather thanchanges in data processing. PSHA requires earthquake size estimates in the form of moment
magnitude to be consistent with modern ground motion characterisation models. In New
Zealand, routine moment magnitude calculations were introduced in 2007, although someearlier estimates for larger earthquakes are also available. The most prevalent magnitude in
the New Zealand earthquake catalogue is a New-Zealand-specific local magnitude, MLNZ77,
which has been consistently calculated for the earliest instrumental records in 1931 until
2011. The MLNZ77 scale is known to estimate larger values compared to moment magnitude for moderate-sized earthquakes in the range 4.5–6.0.
A new seismic processing software, which allows for different magnitude-type calculations,
was introduced in New Zealand in 2012. Routinely calculated magnitudes include a local
magnitude MLSC3 and a summary magnitude M, which is a mix of MLSC3 and a broadband magnitude mB . Recently, a new local-type magnitude has been derived, MLNZ20, thatis calibrated against New Zealand moment magnitude and is applicable for earthquakes ofmagnitude 5.5 and less. We have calculated MLNZ20 for more than 250,000 earthquakes
using more than 2 million individual waveforms with sufficient quality from before and after
the introduction of the new processing software. We have derived regressions of MLNZ20
on the existing local magnitudes, MLNZ77 and MLSC3. These regressions allow us to derive moment-magnitude-consistent size estimates for more than 99% of earthquakes in theearthquake catalogue.
Here we present the updated earthquake catalogue and some of its features.

Keywords: Magnitude regression, moment magnitude, PSHA, New Zealand∗

Speakersciencesconf.org:statsei12:417194

Insights and lessons from more than a decade of
public earthquake forecasting in New Zealand
Annemarie Christophersen∗1 , Matt Gerstenberger2 , and David Rhoades2
1
2

GNS Science – New Zealand
GNS Science – New Zealand

Abstract
Researchers from New Zealand have engaged in developing and testing earthquake forecasting models for more than four decades. They participated in the establishment of the
Collaboratory for the Study of Earthquake Predictability (CSEP) in the late 2000s and set
up a New Zealand earthquake testing centre. The existing earthquake forecasting models
can be classified into three categories, (1) short-term or aftershock models that describe thedecay of the earthquake occurrence rate following a large event, (2) medium-term models
that aim to forecast upcoming large events, and (3) time-invariant models.
Aftershock models include the Short-Term Earthquake Probability (STEP) model, initiallydeveloped by Matt Gerstenberger for California, and the Epidemic Type Aftershock Sequence (ETAS) model, first developed by Yoshi Ogata and implemented in New Zealand by
David Harte. The medium-term models are variations of the Every Earthquake a Precursor
According to Scale (EEPAS) model that was developed by David Rhoades based on theobservation that the size and frequency of earthquake occurrence tends to increase in thevicinity of an upcoming large earthquake. There is a suite of time-invariant models, mostly
based on smoothed seismicity but also including other factors such as strain-rate and proximity to mapped faults.
In 2010, following the M7.1 Darfield earthquake that initiated the Canterbury earthquake
sequence, GNS Science began to provide earthquake forecasts to the public. Initially, theforecasts were derived from the STEP model only. As the sequence continued with reactivation by a series of M6 plus earthquakes for more than one year, hybrid models including
EEPAS and time-invariant models contributed to the forecast.
There were several further earthquake responses over the years. Social science research
has helped in understanding what different stakeholders expected from forecasts and in finetuning the communications.
By the time of the 2016 M7.8 Kaikoura earthquake, the hybrid had evolved to include boththe STEP and ETAS models, two versions of EEPAS, and a mixture of three time-invariant
models. For automatic earthquake forecasting, a Hybrid Forecasting Tool (HFT) has beendeveloped.
In this presentation, we review the progress of providing earthquake forecast to the New
Zealand public and other stakeholders over the last decade. We conclude with some insights
and lessons.
∗

Speakersciencesconf.org:statsei12:417193

Keywords: Earthquake forecasting, New Zealand, large earthquakes, public communications

Exploring Afterslip and Aftershock Co-migration,
Following no Evidence Afterslip Drives Aftershock
productivity
Robert Churchill∗1 , Maximilian Werner1 , Juliet Biggs1 , and Ake Fagereng2
1

University of Bristol – United Kingdom
2
Cardiff University – United Kingdom

Abstract
Aseismic afterslip (transient, gentle fault readjustment) has been proposed as the principal driver of aftershock sequences. Evidence of this relationship includes observationsof similar Omori-type decays, spatio-temporal co-migration, and explanatory mechanicalmodels (e.g. rate-dependent brittle-creep afterslip). However, as this evidence is typically
case-study based, first-order relationships between sequence-scale aftershock behaviours andafterslip abundance are not well-constrained. Given that relative afterslip moment (Mrel,afterslip moment/coseismic moment) can vary < 1% to > 300% for different earthquakes,
we first look for robust correlations between Mrel and aftershock sequence characteristics,
to support or rebut a global driving relationship. We select aftershock sequences from theglobal PDE catalog using three methods, for 41 mainshocks with available afterslip models,and find that Mrel does not correlate with aftershock productivity (number relative to an
Utsu-Seki prediction, relative cumulative moment, or number relative to the backgroundrate), b-value, Omori decay exponent (p) or the background-rate itself. Furthermore, wefind that adding afterslip moment to the mainshock moment does not improve predictions ofabsolute aftershock number. Whilst data bias and uncertainty mean we cannot refute a weakrelationship, our results do not support afterslip as a major driving mechanism of aftershocks
globally. We then explore whether afterslip might control specific spatio-temporal characteristics of aftershock sequences, and whether we can support evidence that afterslip is an
important aftershock triggering mechanism in specific case-studies (just not ubiquitously).
We investigate the spatio-temporal migrations of afterslip and aftershocks for several casestudies in California, including the 2004 Parkfield and 2014 South Napa earthquakes, using high resolution afterslip distribution models and regional seismic catalogs with better
completeness.

Keywords: Afterslip Aftershocks Triggering∗

Speakersciencesconf.org:statsei12:417252

Linking fault roughness and earthquake behavior
Elizabeth Cochran∗†1 , Morgan Page2 , Nicholas Van Der Elst1 , Zachary Ross3 , andDaniel Trugman4
1

United States Geological Survey (USGS) – United States
United States Geological Survey (USGS) – United States
3
California Institute of Technology (Caltech) – United States
4
University of Nevada, Reno (UNR) – United States
2

Abstract
It remains an open question how earthquake behavior may be influenced by the physicalproperties of a fault. The 2016-2019 Cahuilla, California earthquake swarm allows for highresolution mapping of fault roughness, b-values, and stress drop. The sequence comprises
18,250 earthquakes (M-0.27 to M4.4) that ruptured across a 4-km by 4-km fault plane at
4-9 km depth. We first estimate roughness using 2D fault profiles – 150 m-wide bins in thealong-strike and along-dip directions of the fault plane. Roughness is estimated from theresiduals to a line fit to the profiles. We find that the fault is, on average, 50% rougher inthe slip perpendicular direction than parallel to slip. Next, 3D mapping of fault roughness,
estimated from residuals to a plane fit to moving subsets of events, suggests roughness varies
by a factor of 6 over length scales of ˜500 m. Similarly, b-values determined using the b-pos
method vary by a factor of 2 across the fault. In general, we find that b-values are weakly
positively correlated with fault roughness. However, following the largest earthquake (M4.4),
we observe a distinct population of earthquakes with low b-values occurring in an area ofhigh roughness values. While stress drops also vary across the fault, we find no correlation
between roughness and stress drops. We infer that earthquake sequence behavior is partiallycontrolled by fault roughness, especially at low driving stresses.

Keywords: fault roughness, b value, stress drop∗†

Speaker
Corresponding author: ecochran@usgs.gov

sciencesconf.org:statsei12:417200

Statistical relations between global earthquake
magnitudes and parameters of seismograms
Rebecca Colquhoun∗†1 and Jessica Hawthorne1
1

University of Oxford – United Kingdom

Abstract
It is relatively straightforward to determine the magnitude of an earthquake after it has
happened. However, it is unclear if an earthquake ’knows’ its final magnitude before ruptureends. We are interested in whether earthquakes are deterministic: whether features of theinitial stages of an earthquake make accurate predictions about the earthquakes’ final size.
A number of approaches have been taken to investigate determinism. For example, Olson
and Allen (2005) found a relationship between the final magnitude and predominant periodof the early stages of an earthquake. However, the results remain controversial, partly because Olson and Allen (2005) and subsequent researchers analysed only moderate numbers ofearthquakes: of the order of tens or hundreds of events. The increase in data availability andquality over recent years allow us to re-investigate this issue in a more statistically robust way.
We examine the beginning of thousands of earthquakes from around the world. We quantify
predominant period, average period, integral of velocity squared (IV2) and delay time foreach of these earthquakes, and investigate these parameters’ relationship with magnitude.
We find that IV2 and predominant period scale better with final magnitude than average
period.
The thousands of earthquakes used are now allowing us to also robustly analyse smallergroups of events, such as those grouped by depth or tectonic setting. We are also attempting
to relate differences in the fits of these parameters to their physical basis, allowing us tobetter understand earthquake rupture processes.

Keywords: earthquake determinism∗†

Speaker
Corresponding author: rebecca.colquhoun@univ.ox.ac.uk

sciencesconf.org:statsei12:417156

Detection of small slow slip events by means of deeplearning
Giuseppe Costantino∗†1 , Anne Socquet1 , Sophie Giffard-Roisin1 , Mathilde Radiguet1 ,David Marsan1 , and Mauro Dalla Mura2,3
1

Univ. Grenoble Alpes, Univ. Savoie Mont Blanc, CNRS, IRD, Univ. Gustave Eiffel, ISTerre –Université Grenoble Alpes, Université Savoie Mont Blanc, CNRS, IRD, IFSTTAR, ISTerre, 38000
Grenoble, France – France
2
Univ. Grenoble Alpes, CNRS, Grenoble INP, GIPSA-lab – GIPSA-lab CNRS, Université Grenoble
Alpes, Grenoble INP – France
3
Institut Universitaire de France (IUF) – Institut Universitaire de France, Paris, France – France

Abstract
Detecting small Slow Slip Events (SSEs) is still an open challenge. It is a crucial issue
for the characterization of the slip spectrum and for the understanding of the mechanics ofactive faults. The difficulty in revealing low magnitude events is related to their detectionin the geodetic data, notably subduction zones where tremors cannot serve as a proxy forthe slow slip detection. It is therefore necessary to provide new observations by employing
more powerful equipment or to develop novel methods to perceive potential bursts of slow
slip. Unlike seismic data, geodetic observations provide the sole direct access to slow slip
deformation, which need to be analyzed within a large time span (weeks to years). Moreover, they contain spatiotemporal information, which needs to be considered to consistently
constrain the slip evolution. However, the detection of small potential events is still debatedand hard to automatize. Traditionally employed methods rely on the visual inspection of thedata and dedicated modelling methods with a fine-tuning of the parameters (e.g. Ozawa etal., 2001, 2013; Hirose & Obara, 2005; Wallace & Beavan, 2010; Radiguet et al., 2011, 2016;; Wallace et al., 2016; Socquet et al., 2017; Bletery & Nocquet, 2020; Wallace, 2020; Itoh etal., 2022). Here we develop a method based on deep learning for the systematic detection ofSSEs, which may outperform traditional techniques. We focus our analysis on the Cascadia
subduction zone, where a link between slow slip and bursts of tremor activity has been established (Rogers & Dragert, 2003). In this direction, tremor catalogues can be used to validate
potential SSEs detections against the spatiotemporal distribution of tremors. Moreover, acatalog of SSEs has been recently assessed by (Michel et al., 2019), providing an additional
benchmark to our analyses. Nevertheless, the number of available labelled events is still
scarce and does not allow for an adequate training of deep learning models. Therefore, we
make use of synthetic data by generating synthetic slow slip events obtained from synthetic
dislocations (Okada, 1985) and by means of the slab2 model. Each SSE template, assumed
as a sigmoidal-shaped transient, is further added to a window of noise obtained from real
GNSS data.
∗†

Speaker
Corresponding author: giuseppe.costantino@univ-grenoble-alpes.fr

sciencesconf.org:statsei12:417311

We develop a method based on deep learning, using a Convolutional Neural Network (CNN),
followed by a Transformer model, to detect small slow slip events in GNSS data. The CNN
extracts features from the input GNSS time series by learning optimal convolution filters,
while the Transformer focuses on specific parts of the input window thanks to a self-attention
mechanism. We test our method both on synthetic and real position time series from Nevada
Geodetic Laboratory (NGL). Results on synthetics are promising and show a detection tradeoff between the SSEs location, magnitude and the density of the GNSS network. The model
is further compared against real data and tremor catalogues, proving that it exists a magnitude threshold, which depends on the source parameters, down to which the model canreliably (90% confidence) detect slow slip events. The spatio-temporal information is thus
taken into account, showing that our model is able to correctly learn the network geometry
and exploit the temporal dependencies in the data, which is indeed the key feature to explore
the slow slip variability and to address the detection of very small slow slip events.

Keywords: geodesy, seismology, machine learningDeclustering an earthquake catalogue in moderate
seismicity context : strengths and weaknesses ofexisting methods
Guillaume Daniel∗1 and Pierre Arroucau1
1Electricité de France (EDF) – Electricité de France - EDF – 905 avenue du camp de menthe 13097
Aix-en-Provence cedex 02 France, France

Abstract
We present an analysis of the declustering and of the completeness of an updated earthquake catalogue for France. This catalogue was prepared by merging the former FCAT-17
version (Manchuel et al., 2017) with more recent catalogues from RESIF and from neighboring countries. It covers the 250 – 2020 AD period and is composed of events with magnitude
greater than 2.0. This dataset was produced in a seismic hazard assessment perspective andwill be used in a future update of the seismic hazard map previously produced by Drouet etal (2020) for France. This update will mobilize more recent and state-of-the -art methods
and models. Our objective is to provide seismic hazard-oriented end-users with a turn-key
dataset for the estimation of (Poissonian) seismicity rates in any region of France.
We compare the outcomes of several declustering techniques in the aftershock removal process for this long-term catalogue in a low-seismicity context. We discuss the efficiency of
two approaches (Gardner and Knopoff, 1974; Zaliapin and Ben-Zion, 2020) for this removal
operation and discuss the preservation of the Poissonian character of (declustered) seismic
activity. We also revisit completeness periods for France in distinct magnitude bins from theoutputs of several methods (Stepp, 1972; the slope method; Daniel et al., 2011).
These results and this updated catalogue for France (both original and declustered versions)will be made available to the community.

Keywords: France, Seismicity, Declustering, Catalogue Completeness∗Speakersciencesconf.org:statsei12:417317

What injection-induced earthquake sequences tell us
about natural earthquake swarms ?
Philippe Danre∗†1 , Louis De Barros1 , and Frédéric Cappa11Géoazur, UCA, CNRS, IRD, Sophia Antipolis, France – Université Côte d’Azur, CNRS, IRD,
Observatoire de la Côte d’Azur, Géoazur – France

Abstract
Philippe Danré, Louis De Barros, Frédéric Cappa
Université Côte d’Azur, CNRS, Observatoire de la Côte d’Azur, IRD, Geoazur, Sophia Antipolis, France
Whether they are due to natural processes or to fluid injection at depth (wastewater storage,
geothermal activities, etc.), earthquake swarms are singular sequences that present insightful
similarities. For instance, it has been observed for both types of sequences that seismicity migrates with time or that the seismic moment release is weak compared to the spatial extent of
the seismic cluster. If aseismic slip has been observed in a few earthquake swarms, its observation is limited by the depth, low deformation and long duration of those sequences. Moment,
duration and migration velocity scale in a similar way for both types of sequences, indicatingthat dynamics of injection-induced and natural swarms result from the same processes. As
injection-induced and natural swarms seem to follow the same physics, we take advantageof what we know on the former to better understand the latter. Indeed, for decades, studies
have tried to relate seismicity to injection operational parameters, especially seismic momentwith injected fluid volume. We develop two new methods to relate fluid volume to seismic
moment for injection-induced sequences considering the aseismic contribution and validatethem statistically on a set of injection-induced catalogs. We can then propose, for the firsttime, a precise estimation of the fluid volume circulating during different natural earthquake
swarms. This sheds a new light over the processes occurring at depth during those episodes,
and opens promising perspectives to better understand their evolution.

Keywords: Earthquake swarms, Fluid injection, Aseismic slip∗†

Speaker
Corresponding author: danre@geoazur.unice.fr

sciencesconf.org:statsei12:412368

Flexible and Scalable Earthquake Forecasting
Kelian Dascher-Cousineau∗1,2 , Oleksandr Shchur3 , Emily Brodsky4 , and Stephan
Gunnemann3
1

University of California Berkeley (UC Berkeley) – 120 Haviland Path, Berkeley, CA, United States
2
University of California, Santa Cruz (UCSC) – 1156 High St, United States
3
TUM School of Computation, Information and Technology (TUM) – 85748 Garching, Germany
4
University of California, Santa Cruz – 1156 High St, United States

Abstract
Seismology is witnessing an explosive growth in the diversity and scale of earthquake
catalogs owing to improved seismic networks and increasingly automated data augmentation
techniques. A key assumption in this community effort is that more detailed observations
should translate into improved earthquake forecasts. Current operational earthquake forecasts build on seminal work designed for sparse earthquake records that combines the canonical statistical laws of seismology. This parsimonious approach is remarkably robust andstubbornly difficult to improve upon. Advances in the past decades have mainly focused on
the regionalization of the models, the recognition of catalog peculiarities, and the extension
to spatial forecasts; but have failed to leverage the wealth of new geophysical data. Here, wedevelop a neural-network based earthquake forecasting model that leverages the new data
in an adaptable forecasting framework: the Recurrent Earthquake foreCAST (RECAST).
We benchmark temporal forecasts generated by RECAST against the widely used Epidemic
Type Aftershock Sequence (ETAS) model using both synthetic and observed earthquake
catalogs. We consistently find improved model fit and forecast accuracy for Southern California earthquake catalogs with more than 104 events. The approach provides a flexible andscalable path forward to incorporate additional data into the earthquake forecast.

Keywords: Earthquake forecasting, machine learning, temporal point processes∗

Speakersciencesconf.org:statsei12:417346

Analysis of seismic migrations highlight that naturalswarms and injection-induced seismicity are bothdriven by fluid-induced aseismic slip
Louis De Barros∗1 , Philippe Danré2 , and Frédéric Cappa11Géoazur - Université Côte d’Azur – Université Côte d’Azur, CNRS, IRD, Observatoire de la Côted’Azur, Géoazur – France2Géoazur - Université Côte d’Azur – Université Côte d’Azur, CNRS, IRD, Observatoire de la Côted’Azur, Géoazur – France

Abstract
Seismic swarms are characterized by a high rate of small to moderate earthquakes clustered in time and space and they occur in a large variety of natural contexts. They arealso induced by anthropogenic fluid injections in geological reservoirs, such as in enhanced
geothermal systems or for shale gas extraction. In both cases, the driving mechanism required to maintain a seismic activity for days or years is still an open question. Moreover,despite similarities in the behaviors of anthropogenic and natural swarms, it is still not clear
if they may share common driving mechanisms.
Classically, fluid pressure diffusion or slow-slips were proposed as driving mechanisms, but
recent models suggest that swarm seismicity is rather triggered by fluid-induced aseismic
fault slip. Here, using a catalog of natural and anthropogenic swarms, we observe that duration, migration velocity and total moment scale similarly for all swarms. This underlines acommon aseismic slip driving process for both types of swarms, and allow a new possibility
to discriminate fluid-induced swarms from aseismic slip triggered seismicity.
Moreover, one specific characteristic of seismic swarms is the migration of the earthquake
hypocenters. Instead of being directly induced by the fluid pressure diffusion, the seismic
front may map the shear stress concentration at the tips of propagating aseismic slip thatis primarily induced by the fluid pressure increase. In this case, seismic migration depends
on the hydromechanical properties that control the dynamics of aseismic slip propagation
rather than the hydraulic diffusivity.
Using synthetic seismic catalogs computed from the hydromechanical modeling of a slipweakening fault response to a fluid injection, we show that the shape of the seismic front
in a distance-time plot depends on the initial fault criticality to failure, and not on the hydraulic properties. We then extrapolate this numerical result to injection-induced earthquake
swarms. The seismic front has a diffusive behavior when the aseismic slip is directly driven
by the pressure increase, either at the swarm beginning or later if the faults are not initiallycritically stressed. On the contrary, an accelerating seismic front indicates a high-criticality
stressed fault, on which the aseismic slip runs away from the fluid pressurized area. Therefore, the migration behaviors of seismicity prove to be a useful indicator of the swarm driving
process as well as of the stress state of faults over which the swarms occur.
∗

Speakersciencesconf.org:statsei12:417090

The asymmetric seismic moment tensor in
micropolar media
Stephanie Durand∗†1 and Rafael Abreu2
1

Laboratoire de Géologie de Lyon - Terre, Planètes, Environnement [Lyon] (LGL-TPE) – Centre
National de la Recherche Scientifique : UMR5276, Institut national des sciences de lÚnivers, UniversitéClaude Bernard Lyon 1, École Normale Supérieure - Lyon, Institut national des sciences de lÚnivers,
Institut national des sciences de lÚnivers, Institut National des Sciences de l’Univers, Institut nationaldes sciences de lÚnivers, Institut national des sciences de lÚnivers, Institut national des sciences delÚnivers, Institut national des sciences de lÚnivers, Institut national des sciences de lÚnivers, Institut
national des sciences de lÚnivers – 2, rue Raphaël Dubois69622 Villeurbanne Cedex, France
2
Instituet fuer Geophysik, WWU – Germany

Abstract
Seismic events produced by block rotations about vertical axis occur in many geodynamic
contexts. In this study, we show that these rotations can be accounted for using the proper
theory, namely micropolar theory, and a new asymmetric moment tensor can be derived.
We then apply this new theory to the Kaikoura earthquake (2016/11/14), Mw 7.8, one of
the most complex earthquakes ever recorded with modern instrumental techniques. Using
advanced numerical techniques, we compute synthetic seismograms including a full asymmetric moment tensor and we show that it induces measurable differences in the waveforms
proving that seismic data can record the effects of the block rotations observed in the field.
Therefore, the theory developed in this work provides a full framework for future dynamic
source inversions of asymmetric moment tensors.

Keywords: seismic moment tensor, asymmetric, micropolar media∗†

Speaker
Corresponding author: stephanie.durand@ens-lyon.fr

sciencesconf.org:statsei12:416841

Joint Bayesian estimation of b-values and Mc :application to spatio-temporal variations
Stephanie Durand∗†1 , Thomas Bodin2 , Blandine Gardonio2 , and David Marsan3
1

Laboratoire de Géologie de Lyon - Terre, Planètes, Environnement [Lyon] (LGL-TPE) – Centre
National de la Recherche Scientifique : UMR5276, Institut national des sciences de lÚnivers, UniversitéClaude Bernard Lyon 1, École Normale Supérieure - Lyon, Institut national des sciences de lÚnivers,
Institut national des sciences de lÚnivers, Institut National des Sciences de l’Univers, Institut nationaldes sciences de lÚnivers, Institut national des sciences de lÚnivers, Institut national des sciences delÚnivers, Institut national des sciences de lÚnivers, Institut national des sciences de lÚnivers, Institut
national des sciences de lÚnivers, Institut national des sciences de lÚnivers, Institut national dessciences de lÚnivers, Institut national des sciences de lÚnivers – 2, rue Raphaël Dubois69622
Villeurbanne Cedex, France
2
Laboratoire de Géologie de Lyon - Terre, Planètes, Environnement [Lyon] (LGL-TPE) – Centre
National de la Recherche Scientifique : UMR5276, Institut national des sciences de lÚnivers, UniversitéClaude Bernard Lyon 1, École Normale Supérieure - Lyon, Institut national des sciences de lÚnivers,
Institut national des sciences de lÚnivers, Institut National des Sciences de l’Univers – 2, rue RaphaëlDubois69622 Villeurbanne Cedex, France
3
Institut des Sciences de la Terre (ISTerre) – Institut National des Sciences de l’Univers, Institut derecherche pour le développement [IRD] : UR219, Université Grenoble Alpes, Université Gustave Eiffel,
Centre National de la Recherche Scientifique : UMR5275, Université Savoie Mont Blanc – ISTerre,
OSUG-C (Maison des Géosciences), 1381, rue de la Piscine, 38610 GIERES ISTerre, Université Savoie
Mont Blanc, Campus Scientifique, 73376 Le Bourget-du-Lac Cedex, France

Abstract
The statistical distribution of seismicity is described by the so-called GutenbergRichter law which tells us that the logarithm of the probability of having an earthquake ofenergy larger than M0, P(M > M0), decreases linearly with increasing energy, M0, suchthat log(P (M > M0)) = bM0 + a. The slope of this law, the so-called ”b-value”, is investigated to map different seismic regions and is used to compute seismic hazard maps soit has strong society implications. The way the b-value is currently estimated relies on thedetermination of a magnitude of completness Mc which defines the magnitude above which
every earthquake is detected. The estimation of Mc is done independantly of the estimationof b while it can strongly affect the estimation of b. Moreover, as a result any study of seismicity requires to first remove all the earthquakes with a magnitude below Mc which implies
loosing a lot of information. Instead we propose here to revisit this problem of determining
b-values by jointly inverting for b and Mc using a Bayesian formalism and including thestatistical modeling of the events below Mc. We also go further by extending the method to∗†

Speaker
Corresponding author: stephanie.durand@ens-lyon.fr

sciencesconf.org:statsei12:417088

analyze spatio-temporal variations of b and Mc with a Bayesian formalism which provides
the probability of any temporal or spatial changes in these parameters. The algorithm istested on synthetic data, and then applied to a 20 years long catalogue of Japan including
the Tohoku earthquake for which we detect temporal variations of b before and after Tohoku.

Keywords: bayesian, seismicity statistics, b value, Mc, spatio, temporal variations

Overview of Earthquake Rupture Forecasts for the2023 USGS NSHM update
Edward Field∗1
1USGS, Golden – United States

Abstract
On behalf of a very large working group, this presentation provides a summary of anticipated improvements for the time-dependent Earthquake Rupture Forecasts (also known asseismic source characterizations) for the 2023 update of the USGS National Seismic Hazard
Model. The main goals are a more uniform application of methodologies across the country, amore complete representation of epistemic uncertainties, and establishing a basis for an eventual operational earthquake forecasting capability nationwide. Efforts also include a better
representation of multi-fault ruptures and the ability to apply variable degrees of segmentation throughout each fault system (e.g., a wider range of models than applied previously
in California). The presentation also summarizes updated geologic constraints, deformation
models, statistical seismology components, efforts to operationalize various computer codes,
hazard and risk sensitivity analyses (e.g., for logic tree refinement), and any loose ends that
continue to pose a challenge.

Keywords: System Level Earthquake Forecast∗

Speakersciencesconf.org:statsei12:417176

Seismicity under a dormant volcano: identification ofan active crustal fault below Piton des Neiges
Lise Firode∗1 , Zacharie Duputel1 , Valérie Ferrazzini1 , and Olivier Lengliné2
1

Observatoire Volcanologique du Piton de la Fournaise (OVPF) – IPGP, 14 RN3-Km27, 97418, La
Plaine des Cafres, La Réunion, France, Réunion
2
Institut Terre Environnement Strasbourg (ITES) – Ecole Nationale du Génie de l’Eau et del’Environnement de Strasbourg, université de Strasbourg, Institut National des Sciences de l’Univers :UMR7063, Centre National de la Recherche Scientifique : UMR7063 – 5 rue René Descartes, 67084
Strasbourg, France

Abstract
Volcanic environments frequently exhibit seismic phenomena. This is the case for the two
major volcanic edifices on La Reunion Island (Piton des Neiges and Piton de la Fournaise),
where significant seismic activity is recorded. If the seismicity at Piton de la Fournaise
can be easily explained by its volcanic activity (more than three eruptions per year since2014), the seismic activity of Piton des Neiges (inactive for about 27000 years) is still poorly
understood. Following several earthquakes that have been detected and felt by the population
under the north flank of Piton des Neiges, several seismic stations have been deployed in thearea. Here, we use template matching and double-difference relocation to build a highresolution catalog of the region since 2013. Our results reveal that the seismicity is densely
clustered on a northeast dipping fault located in the oceanic crust, under the volcanic edifice.
The orientation of this planar structure is consistent with the analysis of focal mechanisms
indicating reverse fault activity. We also investigate variation in seismicity after correcting
for apparent variations due to changes in the completeness of the catalog (following thedensification of the seismic network since 2016). Preliminary analysis suggests that there
is no significant variation of the b-value over time. The absence of intense seismic swarms
suggests an origin that is not directly related to transient magmatic transfer (contrarily towhat is observed at Piton de la Fournaise). However, some moderate changes in seismicity
rates are observed over the studied time-period and remain to be explained. These results
are intended to contribute to the understanding of deep seismicity (i.e., at depth between 5
and 30 km bsl) under dormant volcanoes. As far as La Reunion island is concerned, they
constitute a first step towards an explanation of the seismicity at the scale of the island.

Keywords: Dormant volcanoe, crustal fault, seismicity rate∗

Speakersciencesconf.org:statsei12:417041

Seismic Clusters Growth and Underlying Triggering
Mechanisms
Tomas Fischer∗1 , Sebastian Hainzl2 , and Josef Vlcek3
1

Charles University, Faculty of Science (CUNI) – Albertov 6, Prague, Czech Republic2
Helmholtz Centre Potsdam, GFZ (GFZ) – Potsdam, Germany
3Charles University, Faculty of Science (CUNI) – Czech Republic

Abstract
Migration of hypocenters is a common attribute of induced injection seismicity and ofearthquake swarms, which usually distinguishes them from aftershock sequences. Spreadingof the triggering front is usually examined by fitting the square-root time dependence of
hypocenter distances from the origin by the pore pressure diffusion model. The earthquake
migration patterns however often exhibit not only spreading envelopes, but also fast-growing
streaks embedded in the overall migration trends. It also turns out that when distance-event
index plots are used instead of distance-time plots, the originally non- linear and interrupted
migration front becomes linear and continuous. This is consistent with a front-rupture
model where the new ruptures are triggered at the edge of previous ruptures, due to stresstransfer supported by the fluid flow in the newly opened pore space. By simulation of theseismicity caused by pore-pressure diffusion and by the front-rupture growth we succeeded
to reproduce the distance-time and distance-event index plots constructed for the observed
earthquake swarms. We apply the new model to selected swarms from Iceland, California
and West Bohemia and to injection induced seismicity. The seismic clusters are identified
semi-automatically in the migration plots and the dependence of the activated area and total
seismic moment are compared with the models of crack growth and front-rupture growth.
The results allow distinguish between different types of seismicity growth and suggest thebackground triggering mechanisms of the observed seismicity.

Keywords: Statistical seismology, fluid induced seismicity, earthquake interaction, seismicity migration∗

Speakersciencesconf.org:statsei12:417341

The symptomatic spatiotemporal clustering of
low-frequency earthquakes and tectonic tremor
William Frank∗1
1

Department of Earth, Atmospheric and Planetary Sciences [MIT, Cambridge] (EAPS) – Massachusetts
Institute of Technology 77 Massachusetts Avenue Cambridge, MA 02139-4307, United States

Abstract
Tectonic tremor and their constituent low-frequency earthquakes often accompany slow
fault slip on plate boundaries. While tectonic tremor are difficult to pin down due to emergent
waveforms without identifiable seismic phases, low-frequency earthquakes can be precisely
located in both time and space with impulsive P- and S-waves. The characteristic feature of
low-frequency earthquakes is that each seismic source emits many repeating events throughtime. We will show evidence that the evolution of these repeating sources follow a well-defined
power-law distribution (e.g. clusters of activity) and are powerful indicators of aseismic fault
motion along plate boundaries. The amplitude distribution of low-frequency earthquakes
during spatiotemporal clustering suggests that clustered activity is driven by an underlying
process, rather than event-to-event interactions as is the case for classic mainshock-aftershock
sequences. We will also highlight similar clustering behavior driven by underlying processes
in volcanic long-period seismicity. Systematically generated earthquake catalogs with a stationary detection capacity are what make the above observations possible and interpretable.
We will thus conclude with recent work on developing a new approach to detecting tremor
with a ”single-station array”. We implement an AI-based framework to learn how to detectand locate seismic sources from the records of a single seismic station based on the knowledge
previously inferred from a dense seismic array. This allows us to detect and locate tectonic
tremors with a uniform detection capacity over the entire lifespan of a single station. We
demonstrate a proof-of-concept of this method in northern Cascadia, allowing us to explore
the spatiotemporal distribution of tremor activity over more than two decades and 20 major
slow slip events.

∗

Speakersciencesconf.org:statsei12:420234

From the Izu-Bonin to the north of Hokkaido : how
did the M9.0 Tohoku earthquake affect the Pacificplate seismicity ?
Blandine Gardonio∗1 , David Marsan2 , Thomas Bodin1 , Stephanie Durand1 , andAlexandre Schubnel3
1

Laboratoire de Géologie de Lyon - Terre, Planètes, Environnement [Lyon] (LGL-TPE) – Centre
National de la Recherche Scientifique : UMR5276, Institut national des sciences de lÚnivers, UniversitéClaude Bernard Lyon 1, École Normale Supérieure - Lyon, Institut national des sciences de lÚnivers,
Institut national des sciences de lÚnivers, Institut National des Sciences de l’Univers, Institut nationaldes sciences de lÚnivers, Institut national des sciences de lÚnivers, Institut national des sciences delÚnivers, Institut national des sciences de lÚnivers, Institut national des sciences de lÚnivers, Institut
national des sciences de lÚnivers, Institut national des sciences de lÚnivers, Institut national dessciences de lÚnivers, Institut national des sciences de lÚnivers – 2, rue Raphaël Dubois69622
Villeurbanne Cedex, France
2
Institut des Sciences de la Terre (ISTerre) – Institut National des Sciences de l’Univers, Institut derecherche pour le développement [IRD] : UR219, Université Grenoble Alpes, Université Gustave Eiffel,
Centre National de la Recherche Scientifique : UMR5275, Université Savoie Mont Blanc – ISTerre,
OSUG-C (Maison des Géosciences), 1381, rue de la Piscine, 38610 GIERES ISTerre, Université Savoie
Mont Blanc, Campus Scientifique, 73376 Le Bourget-du-Lac Cedex, France3Département des Géosciences - ENS Paris – École normale supérieure - Paris – 24 rue Lhomond, 75005
Paris, France, France

Abstract
The last twenty years have seen a number of large, devastating earthquakes on subduction
zones. In many ways, the M9.0 Tohoku-oki earthquake was bewildering for the seismologicalcommunity. It occurred on a previously identified coupled area but ruptured a larger zone
than expected and, above all, the large amount of near-trench coseismic slip was a surprise.
Because Japan is one of the best area instrumented in the world, the 2011 Mw 9.0 Tohoku-oki
earthquake is one of the world’s best-recorded ruptures. Many studies have analyzed withgreat details the pre-seismic, co-seismic and post-seismic phases of the Tohoku earthquake.
Researchers also focused on the triggering of on-land seismicity following the mega-thrustearthquake. However, no study zoom out and considered the consequences of this earthquakeon
the Pacific plate in this area.In this study, we analyzed the Japanese Meteorological Agency
seismic catalog over ten years of data to assess the consequences of such large mega-thrustearthquake over the Pacific plate from the Izu-Bonin area to the north of Hokkaido island.
∗

Speakersciencesconf.org:statsei12:415156

We studied the seismicity from 0 to 700km depth, taking advantage of one of the most complete subduction zone catalogue.
Our results show that the seismic rate south of Japan experienced a decrease at the time of
Tohoku about 30% and an increase of 20% underneath the Hokkaido island. The subduction
zone that is downdip Tohoku doesn’t seem affected by the megathrust earthquake. While itis difficult to understand and to model such large scale effects of the Tohoku earthquake on
the Pacific plate, we think it is primordial to observe and detail them with precision.

The machine learning-based algorithm NESTORE
for strong aftershocks forecasting in seismic clusters
becomes a free available software: application to Italy
Stefania Gentili∗1 , Piero Brondi2 , and Rita Di Giovambattista3
12

National Institute of Oceanography and Applied Geophysics - OGS (OGS) – Italy
National Institute of Oceanography and Applied Geophysics - OGS (OGS) – Italy
3
Istituto Nazionale di Geofisica e Vulcanologia (INGV) – Italy

Abstract
After an intense earthquake, the forecasting of strong subsequent events is strategic for
civil defense purposes because already weakened structures may suffer further damage, increasing the risk of collapse and casualties. NESTORE (Next STrOng Related Earthquake) isan algorithm for probabilistic forecasting of clusters in which a strong mainshock is followed
by at least one aftershock of comparable magnitude. More specifically, NESTORE labels
clusters as type A if the magnitude difference between the mainshock and its strongest aftershock is less than or equal to 1, otherwise as type B. The goal of NESTORE is providing
a near real-time estimation of the probability that the ongoing cluster is type A.
Some earlier versions of the algorithm have been successfully applied to California and Italian
seismicity (Gentili & Di Giovambattista 2017, 2020, 2022). The new version, due to accurate
statistical analysis of the results and the possibility of being applied to regions of higher completeness magnitude, is more robust and provides more reliable results. After the occurrenceof a strong mainshock, the latest version of the NESTORE software, called NESTOREv1.0,
analyzes nine selected seismicity features at increasing time intervals after the occurrence ofthe mainshock, through a multiparameter pattern recognition approach based on one-nodedecision trees. The main issues concerning this type of application are essentially two. The
first is the number of examples (clusters) available in good quality catalogs, which is generally
one order of magnitude lower than the one generally needed in machine learning applications. The second is class imbalance, because a subsequent strong earthquake is recorded
in a percentage lower than 50% of the observed clusters. In recent years, the algorithm has
been improved to increase its generalization capability, despite the small number of clusters
available in seismological databases and the presence of some outliers. A renovation of theold version of the software, both in terms of algorithm and usability, was developed as partof the project ” Analysis of seismic sequences for strong aftershock forecasting” funded bya grant from the Ministry of Foreign Affairs and International Cooperation as part of thescientific and technological collaboration between Italy and Japan.
The current version of the software, NESTOREv1.0, is now mature enough to be distributed
to the scientific community for application and testing in new areas, and will be available
on GitHub in the coming months. NESTOREv1.0 is composed of four modules that respectively: (1) identify clusters from a seismic catalog, (2) find appropriate thresholds for features∗

Speakersciencesconf.org:statsei12:416579

to distinguish cluster types on a training set, (3) use the result of the previous module toverify the performance of the trained algorithm on an independent test set, and (4) allow
near real-time application to clusters as they occur.
As a test of the algorithm, we considered for the Italian territory the seismicity recordedin the last 40 years by ISIDE and OGS catalogs, the latter benefiting from a dense network
in Northeast Italy. Several regions with different seismicity characteristics emerge from theanalysis. The area of northeastern Italy and western Slovenia shows lower productivity bothin terms of number of events and earthquake energy than the rest of the territory. For
this reason, a separate training of NESTORE was applied, which showed results in good
agreement with those obtained with previous versions of NESTORE (Gentili and Di Giovambattista, 2020). The second area is located in a smaller region between northwestern
Tuscany and the southwestern part of Emilia-Romagna, with seismicity characterized byshort-duration bursts with high productivity. Since the available training set for this region is almost exclusively populated by B-type clusters, a local training for that area is notpossible and we excluded that region from our analysis. The remaining larger area, which
includes most of the Italian territory, is characterized by type A clusters mainly in the central Apennines and type B clusters in the remaining territory. We trained NESTOREv1.0
using the ISIDE catalog on the time interval 1980-2009 (corresponding to 22 clusters) andtested the forecasting capability on clusters occurring between 2010 and 2021 (14 cases). The
software test provides a correct cluster type forecasting for 86% of the cases, encouraging
the application of NESTOREv1.0 in the Italian territory.
Acknowledgements
Funded in part by a grant from the Italian Ministry of Foreign Affairs and International
Cooperation

References
Gentili S. and Di Giovambattista R.; 2017: Pattern recognition approach to the subsequent
event of damaging earthquakes in Italy. Physics of the Earth and Planetary Interiors, 266,
pp. 1-17.
Gentili S. and Di Giovambattista R.; 2020: Forecasting strong aftershocks in earthquake
clusters from northeastern Italy and western Slovenia. Physics of the Earth and Planetary
Interiors, 303, 106483.
Gentili, S. and Di Giovambattista, R.; 2022: Forecasting strong subsequent earthquakes in
California clusters by machine learning, Physics of the Earth and Planetary Interiors, 327,
106879.

Keywords: Largest aftershock, Seismicity, Machine learning, Bath’s law, Aftershock sequence,
Italian seismicity

The statistical distribution of time intervals between
consecutive earthquakes worldwideÁlvaro González∗1 , Isabel Serra2,3 , and Alvaro Corral4
1

CRM Centre for Mathematical Research (CRM) – Bellaterra (Barcelona), Spain
2
Autonomous University of Barcelona (UAB) – Bellaterra (Barcelona), Spain
3
Barcelona Supercomputing Center (BSC) – Barcelona, Spain
4
Centre for Mathematical Research (CRM) – Bellaterra (Barcelona), Spain

Abstract
Earthquakes cluster in time in a tighter way than in a Poisson process, in which eventswould be independent from each other and from when each one occurred. This tight clustering must be considered for forecasting the probability of earthquake occurrence in a future
time period, or when testing hypothesis of short- or long-range earthquake triggering.
Nevertheless, the standard analysis of temporal earthquake occurrence usually proceeds by”declustering” the earthquake time series, trying to identify aftershocks or other triggered
events and then pruning them from the sample, leaving only the supposedly independent
events. This procedure attempts to artificially make the process Poissonian-like (so that theprobability of occurrence of the next earthquake is forced to be constant in time).
Since there is not a unique way of identifying triggered earthquakes, this removal is subjective to some degree (it involves lack of knowledge about the process, that is, epistemic
uncertainty). Such a method also reduces the sample itself, reducing the power of any
statistical inference made with it (in other words, with fewer events it is more difficult todistinguish which model best fits the data).
An example of this issue is the debate on whether the recent surge of great earthquakes(magnitude 8 or larger) since 2004 is random or not. If they were Poissonian, the distribution of time intervals between them should be exponential. The answer may depend onwhether triggered events are artificially removed from the sample or not.
In this research, we explore in a comprehensive way the statistical distribution of time intervals between consecutive earthquakes occurred worldwide. As a previous step, we analyse
in detail the completeness of global earthquake catalogues, considering periods separated
by milestones in global seismic monitoring. We use a complete earthquake sample for eachanalysis, and do not attempt to separate triggered from independent events.
We consider different magnitude thresholds, and for each of them test which statistical distribution (such as Weibull, gamma or exponential) best fits the data. This enables us toquantitatively assess whether there is a universal distribution or, on the contrary, if it depends on the magnitude range considered. Also, we test whether Poissonian occurrence canbe rejected for the whole series of the largest earthquakes.
∗

Speakersciencesconf.org:statsei12:417204

Keywords: interevent times, point processes, renewal processes, earthquake probabilities, earthquake forecastingBuilding an Earthquake Forecasting Tool
Kenny Graham∗†1 , Annemarie Christophersen1 , David Rhoades1 , Matt Gerstenberger1 ,Katie Jacobs1 , Rand Huso1 , and Silvia Canessa1
1

GNS Science [Lower Hutt] – 1 Fairway Dr, Avalon, Lower Hutt 5011, New Zealand

Abstract
Operational earthquake forecasting involves the circulation of authoritative information
to help communities understand the likely way sequences will evolve and aid recovery, andfuture resilience planning. Earthquake forecasting informs major decisions aimed at reducing
seismic risk and improving earthquake resilience. In New Zealand, GNS Science through the
GeoNet programme is the official provider of geological hazard information and has been
publishing operational earthquake forecasts since the September 2010 Darfield earthquake.
However, the requirement for human input in the production of forecasts has limited the
speed at which initial forecasts are produced and the frequency with which they can be
updated. Hence, there is a need for an automated tool that can produce forecasts quickly,
regularly and systematically for the whole country or a selected region. Here we discuss a software tool (the Hybrid Forecast Tool, HFT) that integrates different forecasting models andhas the capability to routinely generate forecasts. The incorporated forecast models cover
three different timescales; short-term (from hours to few years, which are constrained byearthquake-clustering statistics and mainly describe aftershock decay), medium-term (fromyears to decades, which captures the increase of seismicity prior to large earthquakes) andlong-term (from decades to centuries, which smooth out the spatial distribution of earthquake
occurrence and can include information on faults and strain-rate). The individual models
were originally developed in different software programming languages (Fortran, Java, andR) and on different operating systems, and they have been tested and published in various forms. The HFT combines all models in a single platform, using a Docker container,
which circumnavigates technical challenges such as the incompatibility of different software
libraries. The HFT can be operated through a graphic user interface (GUI) and commandline that allow for setting up automatic and regular forecast runs. The HFT enables GNS
Science to reliably provide forecasts for future responses to large earthquakes and regularforecast updates.

Keywords: Operational earthquake forecasting, earthquake, clustering statistics∗†

Speaker
Corresponding author: k.graham@gns.cri.nz

sciencesconf.org:statsei12:417268

State of the Art of Time-Dependent Probabilistic
Seismic Hazard Assessment
Kenny Graham∗†1 , Matt Gerstenberger1 , Annemarie Christophersen1 , and David
Rhoades1
1

GNS Science [Lower Hutt] – 1 Fairway Dr, Avalon, Lower Hutt 5011, New Zealand

Abstract
The ability to successfully forecast the future behaviour of a system is a strong indication that the system is well understood. Certainly, many details of the earthquake systemremain obscure, but several hypotheses related to earthquake occurrence and seismic hazard
have been offered. Probabilistic seismic hazard analysis (PSHA) provides the conceptual
framework to quantify the rate (or likelihood) of exceeding various ground-motion levels fora spatial grid within a given time span given all possible earthquakes. For decades, PHSA
has been utilised by governments and industry for making official national hazard maps,
developing building code requirements, determining earthquake insurance rates and even for
deciding safety criteria for large infrastructure such as nuclear power plants and hydroelectric
dams. PSHA in its classical formulation adopts the removal of foreshocks and aftershocks
from earthquake catalogs to separate independent seismic activity from background seismicity. This treatment has been shown to likely cause an underestimation of the seismic hazard,
hence the need for time-dependent PSHA. Time-dependent applications can be classified into
two groups; long-term models based on renewal statistics and short- and medium-term models based on clustering statistics. The latter honours the fact that earthquake occurrence is
not stationary by including spatiotemporal clustering (e.g., aftershock sequences and precursory seismicity). Combining seismic hazard with different spatiotemporal scales could be apractical way to achieve a reduction in the exposure (e.g., the medium-to-long term planningof land use) or the vulnerability (e.g., short-term structural reinforcement or assessment).
Thus, in terms of reducing seismic disaster risk and answering the increasing needs of social
sustainability (towards earthquake resilience), integrating a fully time dependent PSHA iswarranted. However, the current state of this art is under-reported in the literature. Here
we review the concept and recent progress on the study of time-dependent seismic hazard
assessment, with emphasis on short- to medium-term time dependence.

Keywords: Time Dependent, Probabilistic Seismic Hazard Assessment, aftershock sequences, precursory seismicity∗†

Speaker
Corresponding author: k.graham@gns.cri.nz

sciencesconf.org:statsei12:417270

Stochastic Chaos and Predictability in Laboratory
Earthquakes
Adriano Gualandi∗1 , Davide Faranda2,3,4 , Chris Marone5,6 , Massimo Cocco7 , andGianmarco Mengaldo81Istituto Nazionale di Geofisica e Vulcanologia, Osservatorio Nazionale Terremoti (INGV) – Italy
Laboratoire des Sciences du Climat et de l’Environnement, University of Paris-Saclay (LSCE) –LSCE-IPSL (CNRS-CEA-UVSQ), Université Paris-Saclay, Gif-sur-Yvette, France – France
3
London Mathematical Laboratory (LML) – United Kingdom
4
Laboratoire de Meteorologie Dynamique, Ecole Normale Superieure, Paris (ENS) – Laboratoire deMétéorologie Dynamique – France
5”La Sapienza” University of Rome – Italy6
Pennsylvania State University, Geosciences (PSU) – United States
7
Istituto Nazionale di Geofisica e Vulcanologia, Sezione Roma 1 (INGV) – Italy
8
National University of Singapore, Department of Mechanical Engineering, College of Design andEngineering (NUS) – Singapore
2

Abstract
Laboratory earthquakes exhibit characteristics of a low-dimensional random attractor
with a dimension similar to that of natural slow earthquakes. A model of stochastic differential equations based on rate- and state-dependent friction explains the laboratory observations. We study the transition from stable sliding to stick-slip events and find that aperiodic
behavior can be explained by small perturbations (< 0.1%) in the stress state. Friction’s
nonlinear nature amplifies small scale perturbations, reducing the predictability of the otherwise periodic macroscopic dynamics. Similar to weather forecast, we think that ensemble
forecast is the most reasonable way to assess the future state of the system because of theunavoidable stochastic terms that affect the dynamics.

∗

Speakersciencesconf.org:statsei12:417171

A Coulomb Stress response model fortime-dependent earthquake forecasts
Sebastian Hainzl∗1 and Torsten Dahm1
1

German Research Centre for Geosciences (GFZ) – Potsdam, Germany

Abstract
Physics-based models allow extrapolating previously unsampled parameter ranges andenable conclusions on underlying tectonic or human-induced processes. Two widely-usedphysics-based seismicity models are the Coulomb Failure (CF) and Rate-State (RS) models,
which assume pre-existing populations of faults that respond to Coulomb stress changes. The
CF model depends on the absolute Coulomb stress and assumes instantaneous triggering
if stress exceeds a threshold, while the RS model only depends on stress changes. Bothmodels can predict background earthquake rates and time-dependent stress effects, but the
RS model, with its three independent parameters, can additionally explain delayed aftershock
triggering. This study introduces a modified CF model in which instantaneous triggering is
replaced by a mean time-to-failure that depends exponentially on the absolute stress value.
For critical initial stresses, we show that the model leads to identical forecasts as the RS modeland reproduces the Omori-Utsu relation for aftershock decays and stress-shadowing effects.
Thus, both CF and RS forecasts can be seen as special cases of the new model. However, the
new stress response model can also account for subcritical initial stress conditions, which isparticularly relevant for induced seismicity in intraplate regions.

Keywords: stress triggering, earthquake forecast, seismicity model∗

Speakersciencesconf.org:statsei12:416452

Repeating low frequency icequakes at the base ofcold glaciers in the Mont-Blanc massif triggered bysnowfalls
Agnes Helmstetter∗1
1Institut des Sciences de la Terre (ISTerre) – Université Grenoble Alpes – BP 53 - 38041 Grenoble
cedex 9, France

Abstract
We have detected clusters of repeating low frequency icequakes located all over the MontBlanc massif between 2017 and 2021.These signals show a narrow frequency content with acharacteristic frequency of about 5 Hz at all stations, whatever the distance from the source.
These events occur as bursts of tens to several thousand events, with typical inter-event
times of about 2-3 minutes, and last for days or weeks. Bursts of events usually start as
small events with irregular occurrence times and progressively evolve toward more regular
recurrence times. Seismic amplitudes and seismic waveforms also progressively evolve withtime within each cluster, suggesting changes in either rupture length or rupture velocity.
Most events are detected in autumn, winter and spring during and after snowfalls. We used
all available temporary and permanent stations within or around the Mont-Blanc massif tolocate these events. Most clusters are located above 3000 m above sea level, mainly on steep
and North-facing glaciated slopes. Some clusters are clearly associated with cold ice (nearMont-Blanc summit) while others below 2700 m a.s.l. are likely located under temperateglaciers. Other clusters are located near the limit between cold-based and temperate glaciers.
Cold-based glacier mainly move by creep, with little or no basal sliding. This suggests thatthe loading mechanism for these low-frequency icequakes may be viscous deformation within
the ice, while high-frequency icequakes under Argentère glacier have been explained by therepeated failure of asperities surrounded by stable sliding (Helmstetter al al. JGR 2015).

Keywords: low, frequency icequakes repeaters glaciers triggering∗

Speakersciencesconf.org:statsei12:417312

Informative modes in the nearest-neighbor
earthquake diagram
Yu-Fang Hsu∗1 , Ilya Zaliapin2 , and Yehuda Ben-Zion1
1

Department of Earth Sciences, University of Southern California, Los Angeles – United States
2
Department of Mathematics and Statistics. University of Nevada Reno – United States

Abstract
The nearest-neighbor analysis of seismicity in space-time-magnitude domains (e.g., Zaliapin and Ben-Zion, J. Geophys. Res., 2013) can be used for multiple problems of statisticalseismology, including detection of clusters, catalog declustering, and characterization of induced seismicity. Most of these problems use a representation of a proximity measure h as
a product of rescaled time T and rescaled distance R between events, and are facilitatedby the existence of two prominent modes (background and clustered) in the 2-dimensional
joint distribution (T, R). The goal of this work is to utilize additional informative modesin the nearest-neighbor diagram. We analyze earthquakes in California and observe several
distinguishing modes of seismicity, including induced seismicity, repeating earthquakes, andlarge vs. small aftershock sequences. In particular, most of the induced seismicity observedin the Geyser and Coso geothermal fields occurs in the diagram’s background domain. This
induced seismicity mode has similar overall values of h as the tectonic background seismicity,
but is characterized by larger rescaled time T and smaller rescaled distances R. The modeof repeating earthquakes observed prominently in the Parkfield area is characterized by verylarge rescaled time T and small rescaled distances R, with overall h values comparable tothose in the cluster mode of tectonic earthquakes. Aftershock sequences of large and smallto-medium magnitude mainshocks occupy different parts of the cluster mode. We also findmodes that likely represent artifacts generated during the construction of the analyzed earthquake catalog. Furthermore, we illustrate these features of the (T, R) diagram and discuss
how the results provide important information on spatiotemporal distribution features indifferent types of seismicity.

Keywords: nearest, neighbor analysis of seismicity, earthquake clustering, aftershock sequences,
declustering, induced seismicity, repeating earthquakes∗

Speakersciencesconf.org:statsei12:417032

Numerical Investigation of Cascading Foreshock andAftershock in Discrete Fault Network
Kyungjae Im∗†1 and Jean-Philippe Avouac1
1

California Institute of Technology – United States

Abstract
Earthquakes often show a foreshock-mainshock-aftershock sequence. Analyzing the sequence is crucial to understanding the clustered seismicity nature and eventually improving
the predictability of large earthquakes. Cascading has long been considered to play a crucialrole in generating both foreshock and aftershock. The widely adopted earthquake forecast
models, e.g., CRS and ETAS models, are based on the cascading process. In this work, we
numerically investigated the foreshock and aftershock sequence in a discrete fault networkwith rate and state friction law and compared the result with the ETAS model. We set a
fault zone consisting of dense fault segments and an off-fault area consisting of sparsely distributed smaller faults in the simulation domain. The simulations are conducted 100 times
with randomly generated fault location, initial velocity, and fault length within the weighted
distribution, yielding a Gutenberg-Richter law. The simulation results reproduce the properties of foreshocks and aftershocks in nature. Simulated aftershocks occur in the area ofincreased Coulomb stress and decay following Omori law as observed in nature. We observe
immediate foreshocks that cascade toward mainshock within a day and early foreshocks that
do not immediately trigger mainshock but still significantly destabilize the fault zone andsubsequently advance the time to mainshock rupture. Individual foreshock sequences do notshow a clear trend. However, once stacked, it shows apparent inverse-Omori law acceleration. In general, our model similarly behaves to the ETAS model since it incorporates thetwo statistical frameworks of the ETAS model (Omori law and Gutenberg-Richter law), andall seismic events have magnitude-dependent productivity to trigger another event. Comparison to the ETAS model shows that the triggering productivity is lower in the aftershock
sequence than in the foreshocks due to the depletion of critically stressed faults in our simulations. Our model considers physics at an individual fault level and allows us to incorporate
the actual geometry of the fault zone, enhancing the predictability of earthquake behavior.

∗†

Speaker
Corresponding author: kjim@caltech.edu

sciencesconf.org:statsei12:417327

Accounting for earthquake rates’ temporal andspatial uncertainties through least-information
forecasts. Applications to the New Zealand Seismic
Hazard Model.
Pablo Iturrieta∗1 , Gerstenberger Matt†2 , Russ Van Dissen2 , Chris Rollins2 , Ting Wang3 ,and Danijel Schorlemmer1
1

GFZ German Research Center For Geosciences, Potsdam (GFZ) – Germany
2
GNS Science – New Zealand3
University of Otago [Dunedin, Nouvelle-Zélande] – 362 Leith St, North Dunedin, Dunedin 9016,
New-Zeland, New Zealand

Abstract
Earthquake rates are seldom stationary, precluding the elaboration of seismicity forecasts
for low-seismicity regions, where the scarcity of data limits a significant statistical analysis.
We investigate the spatial and temporal performance limits of stationary Poisson seismicity
forecasts, with applications to Probabilistic Seismic Hazard Assessment (PSHA). We design abootstrapping experiment using fast-forward catalogues, i.e. catalogues from high-seismicity
regions as proxy of long-term low-seismicity region catalogues. First, we propose a temporal process that relaxes the condition of stationarity, but rather assumes local-stationarity
between intervals and power-law self-similarity in earthquake number distributions. Results
show variability of earthquake rates to be up to 10 times higher/lower as predicted by Poisson models, where also least-information forecasts in low-seismicity areas, such as uniform
area models (UA), are preferred to smoothed-seismicity-models (SSM). Second, we devise adata-driven method based on strain rates to constrain uniform areas that encompasses thespatial uncertainty, from which a rate distribution can be inferred from the earthquake counts
within. We provide a set of forecasts for the update of the New Zealand National Seismic
Hazard model, which have an increased mean rate values up to 10 times in low-seismicity
regions compared to SSM, whereas also presenting a decrease close to 20% in regions exhibiting high seismicity. Rate variability is represented as Mixed-Poisson processes, such as
the Negative Binomial or Lognormal-Poisson distributions. We further analyse the impactof the forecasts in the seismic hazard space. For 10% probabilities of exceedance in 50 years
in low-seismicity regions, using UA provides an increase of ˜0.1 g of the expected peakground acceleration (PGA) compared to SSM. Including both temporal and spatial uncertainties increases the expected PGA to 0.2 g, whereas the effect on high-seismicity regions
is minimal ( ˜0.01–0.05 g). These results suggest the need to address the non-stationarity
found in training data, which would provide a more realistic forecast than any Poisson-based
method. It also encourages the use of non-Poisson analytical closed forms in forecasts tosolve the seismic hazard integral in PSHA, which can retain the non-stationary features of
the training data.
∗†

Speaker
Corresponding author: m.gerstenberger@gns.cri.nz

sciencesconf.org:statsei12:417122

A size distribution of deep low-frequency tremors in
the Nankai Trough zone: Modeling with a mixture
distribution
Takaki Iwata∗†1 and Yoshihiro Hiramatsu2
1

Prefectural University of Hiroshima – Japan
2
Kanazawa University – Japan

Abstract
As a size distribution of slow earthquakes, a power-law distribution (e.g. Bostock etal., 2015, JGR, doi:10.1002/2015JB012195) and an exponential distribution (e.g., Chestler& Creager, 2017, JGR, doi:10.1002/2016JB012717) have been proposed. Nakano (2019,
GRL, doi:10.1029/2019GL083029) suggested that a tapered Gutenberg-Richter (GR) distribution (Vere-Jones et al., 2001, GJI, doi:10.1046/j.1365-246x.2001.01348.x), which contains
the characteristics of the above two distributions, fits well the size distribution of slow events.
In this study, we propose a mixture distribution of a power-low distribution and an exponential distribution for the size distribution of slow earthquakes. Let π p and π e, whereπ p+π e=1, be the mixture coefficients (or mixture ratio) of power-low and exponential distributions, respectively. We suppose that these coefficients have a dependency on the sizeof events. We also suppose that the value of π e changes from 0 to 1 (and the value of π pchanges from 1 to 0) as the size of events is larger and it is represented by the cumulative
distribution of a normal distribution. Consequently, the mixture distribution approaches an
exponential distribution if the size of events is sufficiently large, and it does a power-low
distribution if the size is sufficiently small.
Then, we consider the product of the mixture distribution and the cumulative distribution of another normal distribution. The second cumulative distribution corresponds tothe detection capability of events (e.g., Ogata & Katsura, 1993, GJI, doi:10.1111/j.1365246X.1993.tb04663.x; Iwata, 2013, GJI, doi:10.1093/gji/ggt208). Thus, this product (ormodel) can be used to analyze the size distribution of observed events without any cutoffsize.
We applied this model to the deep low-frequency (DLF) tremors in the Nankai Troughzone, southwestern Japan. The dataset of the size of DLF tremors is based on Daiku etal. (2018, Tectonophys., doi:10.1016/j.tecno.2017.11.016)(D2018). D2018 constructed five
regions along the Nankai Trough (see Figs. 2 and 3 of D2018 for the details.) and calculated the apparent moment A of a few hundreds of DLF tremors in each of the five
regions from April 2002 to July 2013. In this study, following Nakamoto et al. (2021, JGR,
doi:10.1029/2020JB021138), we used the seismic moment M 0=4πρV sˆ3 A/F s as the sizeof the DLF tremor, where ρ, V S, and F S are the density of material, the S-wave velocity
and the average of S-wave radiation pattern coefficients, respectively. These values are takenfrom Maeda and Obara (2009, JGR, doi:10.1029/ 2008JB006043).
∗†

Speaker
Corresponding author: tiwata@pu-hiroshima.ac.jp

sciencesconf.org:statsei12:417217

We attempted to fit the proposed model to the observed size distributions of the tremorsin the five regions. We also applied an exponential distribution, a power-low distribution,
and a tapered GR distribution to the five size distributions. The parameters of each of thefour models (or distributions) are estimated with the maximum likelihood method. Then,the goodness-of-fit of the models is compared through Akaike Information Criterion (AIC)(Akaike, IEEE Trans. Auto. Control, doi: 10.1109/TAC.1974.1100705). As the value of AIC
is smaller, the goodness-of-fit is better.
For all of the five regions, out of the examined four models, the model proposed in this study
is selected as the best model; the value of AIC is the smallest. The results suggest that thesize distribution of DLF tremors changes from a power-law distribution to an exponential
distribution and that the proposed model is efficient in representing the property of thischange.

Keywords: size distribution, deep low, frequency tremor, slow earthquake, mixture distribution

Testing subduction-zone earthquake models
David Jackson∗1
1David D. Jackson, UCLA, Los Angeles California USA – United States

Abstract
Several global forecast models are under testing to explore earthquake magnitudes andrates in subduction zones. GEAR1, published in 2015, assumed that all subduction zones
had identical dip angle, seismogenic thickness, shear modulus, corner magnitude, coupling,
and all other properties except length and convergence rate. Another recent model, SMERF,
assumed separate values of dip, coupling, corner magnitude, and other parameter for 37 separate subduction zones, thus depending on 6 x 37 fitted parameters. To protect against
over-parameterization, Bayona-Viveros et al. combined the 37 zones into 14 super-zones,
assuming uniform parameters in each. That modification allowed estimation of parameters
over shorter time periods, and resulted in separate forecast models SMERF2, TEAM, andWHEEL. The aggregation was based on a 2016 paper by Kagan and Jackson, who aggregated earthquake data over Flinn-Engdahl zones. While the aggregation helped stabilize
the estimation procedure, the Flinn-Engdahl zones were drawn before plate tectonics was
accepted and they were based almost entirely on earthquake data. Here I’ll present resultswith a similar approach, but with aggregation of zones by tectonic features such as dip angle,
convergence rate etc.. This allows tests of the influence of those properties, without biasing
the tests by selecting earthquake data. I’ll start with the popular hypothesis that dip angles
best explain the largest magnitudes

Keywords: testing, zones, dip angle, regionalization, corner, magnitude∗

Speakersciencesconf.org:statsei12:417221

Coda-Based Estimation of Source Parameters ofLaboratory Acoustic-Emission Events
Tatiana Kartseva∗1,2 , Nikolai Shapiro3 , Andrey Patonin4 , Vladimir Smirnov2,5 , andAlexander Ponomarev2
1

2

CNRS, ISTerre, UGA – Isterre University of Grenoble – France
Schmidt Institute of Physics of the Earth RAS, Moscow, Russia – Russia
3
CNRS, ISTerre, UGA – Isterre University of Grenoble – France
4
CPGR Geophysical observatory ”Borok”, Borok, Russia – Russia
5
Lomonosov Moscow State University, Moscow, Russia – Russia

Abstract
We propose an approach that is aimed to enrich the catalogs of acoustic emission events
recorded in laboratory experiments with such parameters as seismic moment and corner frequency. Because of the difficulty of separation of direct waves in experiments performed on
small samples, we use the coda waves that are composed of the reverberation of the acoustic
field in the tested sample. After multiple reverberations, the resulting wavefield can be approximated as nearly homogeneously distributed over the sample and with signal amplitudes
decaying exponentially in time (linearly in a logarithmic scale).
Within the framework of this model, the frequency-dependent coda amplitude at any moment of time is described as combination of a source spectra, of a decay rate combining
internal attenuation with reverberation losses, and of a sensor response. One of the maindifficulties with the laboratory experiments is that acoustic sensors are very difficult to calibrate and their absolute response function in most of cases remains unknown. With thesimple reverberation model, the logarithms of coda amplitudes at different times and sensors and for multiple events are described by a system of linear equations that we solve ina least-square sense to find frequency-dependent coda-decay rates, relative source spectraand sensor responses. In a next step, we compute spectral ratios between different events
to eliminate the sensor responses and to estimate main source parameters such as cornerfrequencies and relative seismic moments.
We provide details of our data analyses technique and present time-dependent corner frequency vs relative moment diagrams for two experiments on granite of the Voronezh massif
and Berea sandstone under pseudo-triaxial loading. The main trend of inverse proportionality of corner frequencies and relative moments is maintained throughout both experiments.
The dependence close to the cubic that is frequently estimated for tectonic earthquakes observed on the first stages of both experiments when confining pressure steps applied to theintact rock. After applying axial load changes in stress-drop is being observed.
Experiments carried out in the RESC of IPE RAS ”Petrophysics, Geomechanics and Paleomagnetism” on a controlled hydraulic press INOVA-1000 of the Geophysical Observatory”Borok”, IPE RAS.
∗

Speakersciencesconf.org:statsei12:417120

Keywords: acoustic emission, rock physics, coda of a signal, spectral ratio, corner frequency, relative
seismic moment

Modelling the probabilities of anomalous seismicity
rate increases prior to mainshocks
Russell Azad Khan∗†1 , Max Werner1 , Juliet Biggs1 , and Ake Fagereng2
1

University of Bristol (UoB) – University of Bristol, Beacon House, Queens Road, Bristol, BS8 1QU,
UK, United Kingdom
2
Cardiff University – Cardiff University Cardiff Wales CF10 3AT, United Kingdom

Abstract
The idea of using foreshocks to forecast larger earthquakes is not new, however
this idea has never been realised. Regional micro-seismicity catalogues are
relatively new, and research shows you only observe foreshocks if you analyse
micro-seismicity catalogues. Exploiting these new catalogues may reveal key
insights about foreshocks, and how they relate to the earthquake nucleation
phase. Recent work has defined foreshocks as statistically significant increases
in the seismicity rate around where the mainshock occurred, compared to thenormal rates. Using this definition and searching within 20 days and 10km ofmainshocks in Southern California (the Quake Template matching catalogue),around 20% of mainshocks have foreshocks. But once you take into account
how often anomalous seismicity rate increases occur and no mainshock follows,
only 6% of mainshocks have foreshocks. How the normal seismicity rates aremodelled determines how many mainshocks have foreshocks. There is no currently
agreed upon best method, and past models have not been tested to see
how well they represent these rates. We test the following models from theliterature: A Poisson distribution based on the background rate; a Gamma distribution
based on the inter-event time distribution; and a Poisson distribution∗†

Speaker
Corresponding author: russell.azadkhan@bristol.ac.uk

sciencesconf.org:statsei12:417151

based on epidemic type aftershock sequence model conditional intensities. We
define a simple model of normal earthquake behaviour by performing MonteCarlo sampling of the local catalogue around mainshocks. We count how many
earthquakes occur in random 20 day windows. This gives a distribution of thenumbers of events within 20 days. This simple method retains the clusters ofevents that are typical for that area. We then use the Kolmogorov-Smirnov
test to see if the above models could have produced the distribution. We also
propose new models. We fit a negative binomial model to this distribution. We
also produce another distribution in the same manner based on seismic moment
rather than the seismicity rate. Preliminary results suggest the probability amainshock has foreshocks strongly depends on the chosen model, and even negative
binomial models don’t fit the empirical data well. As a result, foreshockprobabilities stated in the literature should be treated with caution.

Keywords: Foreshocks, modelling, seismicity rates

Statistical Power of Tests for Evaluation ofEarthquake Forecast Models
M. Asim Khawaja∗1 , Sebastian Hainzl1 , Behnam Maleki Asayesh1 , Pablo Iturrieta1 , andDanijel Schorlemmer1
1

GFZ German Research Center For Geosciences, Potsdam (GFZ) – Germany

Abstract
The Collaboratory for the Study of Earthquake Predictability (CSEP) is an international effort to evaluate earthquake forecasting models independently and provide the cyberinfrastructure together with a suite of testing methods. For testing forecasts, CSEP definesa grid-based format to represent the forecast in terms of the expected number of earthquakes
for each cell in the grid. In CSEP, the spatial consistency of a forecast is evaluated using
the Spatial test (S-test), based on joint log-likelihood evaluations. The high-resolution gridcombined with inhomogeneous earthquake distributions leads to many empty cells that may
never experience an earthquake and cause a huge disparity in the number of cells and thenumber of observed earthquakes. This imbalance can lead to the lack of statistical powerof the S-test and forecast testing in general. The power of the S-test defines the ability ofthe test to correctly identify whether two distributions of data (e.g., earthquake forecast andactual seismicity) are the same or different.
To explore this issue, we take global earthquake forecasts as a case study. We find that for avery fine-meshed single-resolution grid, an unrealistic high number of approximately 32000
earthquakes is required in the test catalog to achieve the highest power of the S-test. The
required number of earthquakes in the test catalog decreases for decreasing grid resolution,
but uniformly decreasing the resolution leads to the loss of spatial information offered by theforecast model in seismically active regions. Thus, we use data-driven multi-resolution grids
to test the forecast models in which the resolution is determined by the availability of data.
The data-driven multi-resolution grids offer powerful testing with as low as 16 earthquakes inthe test catalog. Therefore, we propose using multi-resolution grids in future CSEP forecast
experiments and further studying its application in regional and local experiments.
We further explore multi-resolution grids for spatial forecasts of aftershock activity. For that
purpose, we study stress-based forecast models with a binary classification approach using
measures, i.e., ROC curve and MCC-F1, based on the contingency matrix. We observe thatthe choice of grid resolution also affects the binary test results. Currently, we are striving to
explore the quantitative effect of the grid resolution on those performance metrics.

Keywords: Earthquake forecasts evaluation, Multi resolution grid∗

Speakersciencesconf.org:statsei12:416655

Machine learning DAS signal denoising without clean
ground-truth signals
Sacha Lapins∗†1 , Antony Butcher1 , J-Michael Kendall2 , Thomas Hudson2 , Maximilian
Werner1 , and Anna Stork3
1

University of Bristol – United Kingdom
2
Oxford University – United Kingdom
3
Silixa – United Kingdom

Abstract
We present a weakly supervised machine learning (ML) method for suppressing strong
random noise in distributed acoustic sensing (DAS) recordings. The method aims to map
random noise processes to a chosen summary statistic, such as the distribution mean, median or mode, whilst retaining the true underlying signal. The dataset used is a DAS array
deployed at the surface of the Rutford Ice Stream in Antarctica. Despite being a very low
anthropogenic noise environment, strong random noise processes, such as weather-induced
instrument noise, heavily dominate the signal from microseismic icequake events. Here, we
demonstrate that the proposed method greatly suppresses incoherent noise and enhances
the signal-to-noise ratios (SNR) of microseismic events, enhancing the performance of subsequent processing steps, such as event detection. Our best performing model for this task isextremely lightweight by deep learning standards (three hidden layers, 47,330 model parameters), processing 30 secs data recorded at a sample rate of 1000 Hz over 985 channels ( ˜ 1
km of fibre) in < 0.3 secs. Furthermore, it is trained in a ‘weakly supervised’ manner, suchthat it requires no manually-produced labels (i.e., pre-determined examples of clean event
signals or sections of noise) for training. Lastly, our GPU-optimised Tensorflow implementation is considerably faster than standard Python filtering routines, such as bandpass filtering
using the open-source obspy framework, and doesn’t require any prior assumptions on the
distribution of the noise or event signals, unlike traditional spectral methods. We argue that
efficient data-driven denoising methods, such as the one presented, will prove essential totime-critical DAS detection and early warning processing workflows, particularly in the caseof microseismic monitoring.

∗†

Speaker
Corresponding author: sacha.lapins@bristol.ac.uk

sciencesconf.org:statsei12:415708

Source characterization with 3 seconds of records ona single station
Pablo Lara1,2 , Quentin Bletery∗1 , Jean-Paul Ampuero1 , Adolfo Inza2 , and Hernando
Tavera2
1

Université Côte d’Azur, IRD, CNRS, Observatoire de la Côte d’Azur, Géoazur – Institut de recherche
pour le développement [IRD] : UR082 – France
2
Instituto Geofı́sico del Perú (IGP) – Peru

Abstract
How fast can we produce an estimation of the magnitude and source location of anearthquake? Present early-warning systems require at least 4 stations to provide a firstestimation. This introduces inherent delays which depend on the instrumentation density.
We present here a Machine Learning (ML) approach for earthquake detection and source
characterization running on a single station. We developed two separate algorithms for
detection and source characterization. Our detector is trained on noise and earthquake
waveforms from the STEAD global database, plus a regional dataset from Peru. It achieves
over 99.9 % accuracy. The second algorithm is trained on a database including STEAD,
Peruvian, Chilean and Japanese waveforms. We extract 164 attributes (in the time, spectral
and cepstral domains) from seismic traces including only 3 seconds of records after the arrivalof the P wave. We then train the ML algorithm to predict the magnitude, epicentral distance,
back-azimuth and depth of the earthquake from these 164 attributes. We repeat the processwith 4 sec of data, then 5, up until 50 seconds, allowing us to provide updated estimates
every second. The algorithm is currently running in beta version in the early-warning system
of Peru and will go online in July 2022.

Keywords: Machine Learning, early warning, single, station∗

Speakersciencesconf.org:statsei12:415093

Instantaneous Tracking of Earthquake Growth With
Elasto-Gravity Signals
Andrea Licciardi∗1 , Quentin Bletery1 , Bertrand Rouet-Leduc2,3 , Jean-Paul Ampuero1 ,and Kévin Juhel1,41Géoazur – Université Côte d’Azur, IRD, CNRS, Observatoire de la Côte d’Azur, Géoazur, Sophia
Antipolis, France – France
2
Disaster Prevention Research Institute, Kyoto University, Kyoto. Japan – Japan3
Los Alamos National Laboratory, Geophysics Group, Los Alamos, NM, USA – United States
4
Laboratoire de Planétologie et Géodynamique, UMR 6112, Nantes University, CNRS, Nantes –Laboratoire de Planétologie et Géodynamique, UMR 6112, Nantes University, CNRS, Nantes – France

Abstract
Rapid and reliable magnitude (M w ) estimation for large earthquakes (M w > = 8) is
key tomitigate the risk associated with strong shaking and tsunamis. Standard early warning
systems based on seismic waves fail to rapidly estimate the size of such large earthquakes.
GNSS-based approaches provide better estimations but are also subjected to large
uncertainties and latency associated with the slowness of seismic waves. The recently
discovered speed-of-light Prompt Elasto-Gravity Signals (PEGS) has raised hopes toovercome these limitations, but has never been tested for operational early warning. Here,
we show that PEGS can be used in real-time to track earthquake growth instantaneously
after
the event reached a certain magnitude. We develop a deep learning model that leverages theinformation carried by PEGS recorded by regional broadband seismometers in Japan before
P-waves. After training on a database of synthetic waveforms augmented with empiricalnoise, we present the first example of instantaneous tracking of an earthquake source timefunction on real data. Our model unlocks ”true real-time” access to the rupture evolutionof large earthquakes using a portion of seismograms that is routinely treated as noise, andcan immediately be transformative for tsunami early warning.

∗

Speakersciencesconf.org:statsei12:418270

Comparing multiple deep learning models to detect
hydraulic-fracturing induced seismicity
Cindy Lim∗1 , Sacha Lapins1 , Margarita Segou2 , and Maximilian Werner1
1

2School of Earth Sciences [Bristol] – United Kingdom
British Geological Survey [Edinburgh] – United Kingdom

Abstract
The monitoring of hydraulic fracturing-induced seismicity (HFIS) is crucial to assess andmitigate risks to protect infrastructure and livelihood. Deep learning phase picking models
have been proposed as a solution for efficiently processing large volumes of data typically
produced during microseismic monitoring. These models have detected regional earthquakes
not catalogued by existing traditional methods (e.g., STA/LTA, autocorrelation), which canlead to changes in the statistical characteristics (e.g., Gutenberg-Richter b-value) of the catalogue which can be used to predict the maximum magnitude that can occur at a site aswell as to predict the rate of earthquake in relation to injection rate. It also presents theopportunity to observe new changes in spatiotemporal evolution of induced seismicity anddiscern the relative importance of different driving mechanisms of HFIS.
We chose four existing phase picking models (GPD, U-GPD, PhaseNet and EQTransformer)
trained on large volumes of regional earthquakes recorded on surface station datasets (100 Hz)
and investigated how well they generalise and pick microseismic phases in higher frequency(2000 Hz) borehole data.The PNR-1z dataset comprises continuously recorded injection operations on the PNR-1z well at a hydraulic fracturing site in Preston New Road, England
where the operators catalogued > 38,000 events using the Coalescence Microseismic Mapping(CMM) method. We examined phase classification abilities of each model and subsequently
generated earthquake origin time catalogues for the PNR-1z dataset to compare against theinitial catalogue.
Results show that some of the models (particularly PhaseNet) recognise seismic phases wellwithin our data: they recover up to 94.6% of the initial catalogue and also found a large
proportion of undetected events (> 32,800, an 80% increase in events). PhaseNet’s successful
application on our dataset might be attributed to its exposure to different instrument data
during training, as well as its comparatively small model size which likely reduces overfitting
to its initial training set. In addition, because PhaseNet’s architecture is fully convolutional,
the model produces a phase probability for each sample in a continuous seismic trace, thus
resulting in precise phase picks. We conclude that in its current state, PhaseNet can beapplied off-the-shelf on high frequency borehole data to detect HFIS. The many events thatwent undetected previously might reveal new insights into the spatio-temporal evolution ofseismicity during fluid injections.

∗

Speakersciencesconf.org:statsei12:417487

An extension of the ETAS model to 3-D
Behnam Malekiasayesh∗1,2 , Sebastian Hainzl2 , Gert Zöller3 , and Matthias Holschneider3
1

Institute of Geosciences, University of Potsdam, Potsdam – Germany
GFZ German Research Centre for Geosciences, Potsdam – Germany
3Institute of Mathematics, University of Potsdam, Potsdam – Germany
2

Abstract
The Epidemic Type Aftershock Sequence (ETAS) model is a widely used stochastic point
process model for describing the distribution and clustering of earthquakes in space andtime. For historical reasons, the ETAS model is two-dimensional and uses only epicentral
information because of the large depth uncertainties in the past. Nowadays, however, earthquake catalogs also provide high-precision depth values that contain valuable information.
Therefore, 3-D versions of the ETAS model are desirable, exploiting the full hypocentral
information of earthquake catalogs. We develop a new 3-D ETAS model that incorporates a
depth kernel into the trigger function. To explore its most appropriate parametric form, wefirst examine the empirical depth distribution of aftershocks as a function of the depth andmagnitude of the mainshock. To this end, we apply the nearest-neighbor distance (NND)
method to select direct aftershocks observed in California and other regions and fit truncated power-law kernels that account for the free surface and finite depth of the seismogenic
zone. Finally, we show initial applications of the developed 3D ETAS model for several casestudies.

Keywords: 3D ETAS model, Nearest Neighbor Distance (NND), statistical seismology, point processmodeling.

∗

Speakersciencesconf.org:statsei12:417078

Lab Earthquake Prediction, Precursors, and the
Evolution of Size-Frequency Distributions During the
Lab Seismic Cycle
Chris Marone∗1,2
1Dipartimento di Scienze della Terra, La Sapienza Università di Roma, Italia – Italy2
Department of Geosciences, Pennsylvania State University, USA – United States

Abstract
Recent work shows that frictional stick-slip failure events, the lab equivalent of earthquakes, can be produced under controlled conditions for the full range of slip behaviors observed on tectonic faults. The spectrum of lab fault slip modes includes a continuous range
from aseismic slip and creep events to aperiodic slow slip and fast, elastodynamic failure
events. These data show that labquakes are routinely preceded by foreshocks in the form of
acoustic emissions. They also reveal systematic changes in fault zone elastic properties that
represent robust precursors to failure for the full spectrum of slip modes. Foreshocks include
a cascade of microfailure events that show systematic evolution of the frequency magnitude(FM) relations, with an increase in event rate and a greater proportion of larger events astime to failure approaches. Machine learning (ML) techniques based on these changes arecapable of predicting the failure time of labquakes and also the fault zone shear stress andother aspects of labquakes.. These works provide an opportunity to use ML to interrogate
the physics of impending failure. I summarize the physical processes of ML-based failure
predictions, precursors, and changes in FM relations during the lab seismic cycle.

∗

Speakersciencesconf.org:statsei12:417490

Earthquake swarms along the Chilean subduction
zoneDavid Marsan∗1 , Anne Socquet2 , and Thomas Reverso3
1

David Marsan (ISTerre) – Université Savoie Mont Blanc, Université Grenoble Alpes, CNRS, IRD,IFFSTAR – France
2
ISTerre – UGA-CNRS – France
3
ISTerre – ISTerre USMB – France

Abstract
We search for earthquake swarms along the Chilean subduction, from $-18ˆo$ to $38ˆo$ of latitude, 2003-2020, as a way to provide better insights on where and when slow
slip episodes occur. We develop an objective, statistical, model-based method that detectswarms as anomalous changes in earthquake rate. A set of 16 swarms is obtained after
careful inspection of the candidate swarms; this assessment allows to reject cases for whichmodel errors, rather than ’true’ anomalous rate changes, are likely the cause of the detection.
Averaging the activity over these 16 episodes, we find indirect evidence for a mostly aseismic
driving mechanism, and a mean aseismic to seismic ratio estimated to 89 when using theseismicity rate as a proxy for slip. All the swarms are found in the 20-50 km depth range
with the notable exception of one 60-100 km-deep swarm that occurs several days after the2010 Maule earthquake and downdip of it. The dominant depth range (20-50 km) is inagreement with previous studies that suggest this range to be a transition zone from theshallower, locked part of the subduction, to the freely slipping interface at greater depth andintraslab earthquake activity. The swarms can be separated into three spatial groups, thatare all related to a subducting oceanic ridge. This structural control by fluid-rich geologicalfeatures is modulated by stress control, i.e., swarms cluster in time with intermediate tolarge ruptures, both prior and following them, pointing to a close interplay between seismic
and aseismic slip in relatively well segmented parts of the Chilean subduction.

∗

Speakersciencesconf.org:statsei12:417916

Influence of sea level changes on seismicity rates
from a hydrothermal system in the Marmara region,
TurkeyPatricia Martı́nez-Garzón∗1 , Gregory Beroza2 , Gian Maria Bocchini3 , Marco Bohnhoff4 ,and Grzegorz Kwiatek5
1Helmholtz-Centre Potsdam GFZ German Research Centre for Geosciences – Potsdam, Germany
2
Department of Geophysics, Stanford University – United States
3
Ruhr-University Bochum – Germany
4Helmholtz Centre GFZ Potsdam – Germany
5
Helholtz Centre GFZ Potsdam – Germany

Abstract
For decades, periodic Earth phenomena, such as tidal loading have been studied to investigate whether their associated stress changes are enough to bring a fault to slip. Triggeringof tectonic tremor by solid earth tides has been observed almost everywhere that tremor occurs, in regions typically displaying lithostatic pore pressures and low effective normal stress.
Observations of tidal triggering of earthquakes are far less common, and are mostly related toregions where the Earth’s crust is close to failure, such as geothermal and magmatic systems,
where small stress changes can induce fault slip. As the stress changes are small, typicallyseismicity of small magnitude displays a larger response to tidal loading. Therefore, highresolution seismicity catalogs with low detection thresholds are needed to capture it. In the
Sea of Marmara region of the North Anatolian Fault Zone, Turkey, a fault segment of about
150 km is late in its seismic cycle. The Armutlu Peninsula, directly south of Istanbul, is ahydrothermal region rich in crustal fluids and with abundant hot springs. We document astrong effect of sea level changes in the Sea of Marmara on the seismicity of the Armutlu
Peninsula. During a 90-day time period corresponding to the winter months, changes in thesea level up to 1 m were recorded. The sea level fluctuations display an overall period of 13.6
days, corresponding to the fortnightly tides, but are also strongly affected by storms. We
generated new high-resolution seismicity catalogs for this region using both AI-based andtemplate matching techniques, and compared the results with those obtained from catalogs
employing traditional techniques for earthquake detection and association. Both templatematching and AI-based catalog show that local seismicity rates are significantly larger during
time periods when the sea level is rising. Available data from local strainmeters suggested
that the associated strain changes on the order of 100-400 nstrain are sufficient to bring the
local faults to failure.

Keywords: earthquake triggering, tidal loading, hydrothermal∗

Speakersciencesconf.org:statsei12:417082

Beyond cascade and pre-slip: complexities involvedin the generation of earthquake sequences
Patricia Martı́nez-Garzón∗1 and Piero Poli2
1

Helmholtz-Centre Potsdam GFZ German Research Centre for Geosciences – Potsdam, Germany
2
CNRS, Université de Grenoble – ISTerre lab., University Savoie Mont-Blanc - CNRS, Le
Bourget-Du-Lac, France – Université Grenoble Alpes ISTerre CS 40700 38058 GRENOBLE Cedex 9,FranceAbstract
Detection and analysis of earthquake precursory processes is a challenging topic thathave concerned geoscientist from different disciplines over almost a century. The classicalview on how earthquakes nucleate evolved for decades around two widely accepted models
where the mainshock is promoted by either the stress change from previous earthquakes
or the occurrence of an aseismic transient. The improvements in earthquake monitoring,
integration of geodetic techniques to capture slow deformation, and the incorporation of
artificial intelligence for seismic data processing recently allowed a more complete picture on
how earthquake sequences evolve before the occurrence of the mainshock. Here, we review
and compare the available analyses describing observations on precursory processes of 32
earthquake sequences covering a magnitude range MW (3.2, 9.1). We combine and present
both seismological and geodetic analysis and discuss the patterns that emerge collectively,
with a special emphasis on relating the observations with the physical processes driving theevolution of earthquake sequences. Our analysis highlights various structural, tectonic andboundary conditions that play a role in the dynamics of the earthquake sequences, as well as
various physical processes that may be acting superimposed. From it, we propose the ideaof earthquake nucleation as a complex imbricated process, involving several temporal andspatial scales.

Keywords: earthquake initiation, precursory processes, review∗

Speakersciencesconf.org:statsei12:416775

The space-time variability of theMagnitude-Frequency distribution
Warner Marzocchi∗1 , Marcus Herrmann1 , Ilaria Spassiani2 , and Ilaria Spassiani1
1

University of Naples, Federico II, Dept. Earth, Environmental, and Resources Sciences, Complesso diMonte Sant’Angelo, Via Cupa Nuova Cintia, 21, 80126 Napoli – Italy2
Istituto Nazionale di Geofisica e Vulcanologia - Sezione di Roma (INGV) – Via di Vigna Murata, 605 00143 Roma, Italy

Abstract
The Magnitude-Frequency Distribution (MFD) is one of the most prominent aspects ofearthquake occurrence and, as commonly used, it represents a major limitation in improving
earthquake predictability.
According to the first round of experiments carried out by the Collaboratory for the studyof earthquake predictability (CSEP), the best earthquake forecasting models use a common
space-time invariant Gutenberg-Richter (GR) law, performing significantly better than models with space varying b-values (the parameter of the GR law). However, despite a GR-likelaw is found to describe well magnitude distribution on large space-time windows, many independent pieces of evidence – stemming from the analysis of seismicity, model requirements(e.g., UCERF3), and physical constraints – raise doubts on its reliability when zooming onmuch smaller space-time windows.
Here, we show the space-time complexity of MFD in some recent seismic sequences, where
high-resolution earthquake catalogs are available. The results of these analyses highlightsome important new scientific results and new challenges. First, the importance of defining
carefully a priori the space-time window in which the MFD is calculated, for example using geological constraints and/or machine learning techniques. Second, a signature of theelastic rebound in faults that have just slipped, which we have modeled using an energy(time) dependent corner magnitude of the tapered GR distribution. Third, some important drawbacks in calculating the MFD in high-resolution catalogs, which derives from astill common inhomogeneity of the magnitudes; this represents a major problem in modeling unbiased MFDs and comparing them when they come from different earthquake catalogs.
Finally, we discuss the implications of these results in terms of the statistical and physical behavior of seismic sequences, and of the improvement of earthquake predictability andforecasting.

Keywords: magnitude, frequency distribution, earthquake predictability, high, resolution earthquake catalogs∗

Speakersciencesconf.org:statsei12:417210

A deep catalogue of focal mechanisms for the 2016
Amatrice, Italy earthquake sequence
Men-Andrin Meier∗1
1

ETH Zurich – Switzerland

Abstract
The 2016 Amatrice, Italy earthquake sequence occurred on a normal fault system in theCentral Apennines. With ˜140 permanent or temporary seismic stations directly above theseismic activity, the sequence has been exceptionally well recorded. Tan et al., 2021, The Seismic Record, have used PhaseNet (Zhu and Beroza, 2019, GJI) to detect ˜15 million seismicphases, and compiled a high-precision hypocentre catalogue with ˜900’000 events with the
HypoDD method of Waldhauser and Ellsworth, 2000, BSSA. Starting with this catalogue,
we use a convolutional neural network classifier to predict P-wave first motion polarities,
and the HASH algorithm of Hardebeck and Shearer, 2002, BSSA to compute a deep focal
mechanism catalogue. The catalogue consists of > 56’000 focal mechanisms, about 8’000 ofwhich have nodal plane uncertainties below 25 degrees. In contrast to existing, conventional
focal mechanism catalogues, the deep catalogue samples almost the entire study region, andalmost the entire magnitude range ( ˜M0-M5), although nodal plane uncertainties generally
tend to increase with decreasing magnitude. We use the focal mechanism catalogue to study
the kinematics of the Amatrice earthquake sequence, to test the hypothesis of a coseismic
rotation of the static stress field by large and small events, and to analyse the complexity ofthe stress field before, during and after the earthquake sequence.

∗

Speakersciencesconf.org:statsei12:425781

Updates to the Fault-System Inversion Approach foruse in the 2023 U.S. National Seismic Hazard Model
Kevin Milner∗1 and Edward Field2
1Southern California Earthquake Center, University of Southern California (SCEC/USC) – UnitedStates2
Geologic Hazards Science Center, U.S. Geological Survey – United States

Abstract
A 50-state update to the U.S. National Seismic Hazard Model will be released in 2023
(NSHM23). NSHM23 will utilize an inversion-based methodology for active fault systems in
the Western U.S., building upon the approach used in the 3rd Uniform California Earthquake
Rupture Forecast (UCERF3). The fault-system inversion approach enables relaxation of
fault segmentation assumptions, inclusion of multi-fault ruptures, and mediation between andbetter fitting of competing data constraints (e.g., slip rates, paleoseismic recurrence intervals,
and regional seismicity). We will describe a number of model improvements, including: avariable segmentation constraint to dial the propensity of multi-fault ruptures up and down;
a decoupling from the gridded seismicity model; a simpler methodology for determining theregional target magnitude frequency distribution (MFD); a wider range of fault section MFD
nucleation constraint options; better constraint weighting to more evenly fit all available data
relative to their uncertainties; and several improvements to the simulated annealing algorithm
that solves the inversion. In aggregate, these changes lead to a suite of models that better
fit the available data constraints, encompass a broader range of uncertainties, and are morecomputationally efficient than UCERF3. We will then present preliminary NSHM23 fault
system inversion solutions and compare them with prior models.

Keywords: NSHM23, UCERF3, inversion, simulated annealing, ERF, hazard∗

Speakersciencesconf.org:statsei12:417314

Evolution of the seismicity rate during and after
hydraulic fracturing operations in Preston New
Road, UK
Riccardo Minetto∗1 and Agnès Helmstetter2
1

Institut des Sciences de la Terre (ISTerre) – Institut National des Sciences de l’Univers, Institut derecherche pour le développement [IRD] : UR219, Université Savoie Mont Blanc, Centre National de laRecherche Scientifique, Université Gustave Eiffel, Université Grenoble Alpes – France
2
Institut des Sciences de la Terre (ISTerre) – Institut National des Sciences de l’Univers, Institut derecherche pour le développement [IRD] : UR219, Université Savoie Mont Blanc, Centre National de laRecherche Scientifique, Université Gustave Eiffel, Université Grenoble Alpes – France

Abstract
In August 2019 an hydraulic fracturing operation was carried out at the PNR-2 well in
Preston New Road, UK. Hydraulic fracturing caused abundant seismic activity that culminated with a ML2.9 event. The seismic activity was recorded by a downhole array of 12
sensors located in a nearby monitoring well. 55555 events were detected and located in realtime during the operations. After improving the number of detections with template matching, distinct seismicity patterns can be identified by separately analysing the periods during
the injections and the periods after the injections. More specifically, peaks of seismicity canbe observed during the injections. These peaks gradually occur at later times, while theirintensity gradually decreases. After the injections, the seismicity decay can be described by
the Omori-Utsu law. The decay does not remain constant, but tends to slow down with each
injection. The observed variations in seismicity rate can be explained in terms of rate- andstate- theory (RTS) and critical pressure theory (CPT).

∗

Speakersciencesconf.org:statsei12:417240

Towards an Operational Earthquake Forecasting
Model for Europe
Leila Mizrahi1 , Shyam Nandan1 , Marta Han∗1 , and Stefan Wiemer1
1

Swiss Seismological Service, ETH Zurich – Switzerland

Abstract
The go-to models for developing time-dependent earthquake forecasts are EpidemicType Aftershock Sequence (ETAS) models. They model earthquake occurrence as a spatiotemporal self-exciting point process, using basic empirical laws such as the Omori-Utsu law
for the temporal evolution of aftershock rate, the Gutenberg-Richter law to describe the size
distribution of earthquakes, the exponential productivity law and so on. The main focus andcore strength of ETAS lie in modelling aftershock occurrence. An important aspect whichholds great potential for improvement is the modelling of background seismicity. In mostformulations of ETAS, background seismicity is either treated to be constant in time andspace during parameter inversion or only considered spatially varying. The spatial variationof the background seismicity is generally inferred by smoothing the location of past declustered earthquakes. Here, we develop an algorithm based on expectation maximization, whichallows the inversion of a covariate-dependent background seismicity rate and the remaining
ETAS parameters to be self-consistent. We test a wide range of covariates ranging fromstrain-rate maps, fault slip rates, and relative seismic intensity in expert defined zonations.
We develop these ETAS models for continental Europe leveraging the rich data sets and expert solicitations acquired for building the European Seismic Hazard Model (ESHM) 2020.
We rigorously test these suites of covariates pseudo-prospectively to identify the ones thatlead to the best operational forecasting model for Europe.

∗

Speakersciencesconf.org:statsei12:417106

Component-driven ensemble modelling withsemi-parametric ETAS models
Leila Mizrahi∗1 , Shyam Nandan1 , William Savran2 , Stefan Wiemer1 , and Yehuda
Ben-Zion2
1

Swiss Seismological Service, ETH Zurich – Switzerland
2
University of Southern California – United States

Abstract
Epidemic-Type Aftershock Sequence (ETAS) models represent the state-of-the-art benchmark against which to test any newly developed earthquake forecasting models. These newmodels can come in all shapes and flavours, and the idea of ensemble modelling has received
attention in recent years within the statistical seismological community. Ensemble modellingaims to produce more informative forecasts by capturing the beneficial properties of different
individual models within an overarching model ensemble. Traditionally, such an ensemble
is a linear (or multiplicative) combination of ingredient models, and the challenge is to optimize the weights given to each ingredient model. We here propose a novel approach forcomponent-driven ensemble modelling, where a forecast is created by marrying componentsthat are derived from different models.
The components are driven by the following questions: How many background events do we
expect? Where do we expect them? When do we expect them? How many aftershocks do we
expect? Where do we expect them? When do we expect them?
By answering each of these questions with different ingredient models, we can easily createa suite of ensembles. Using pseudo-prospective forecasting experiments, we then compare
the performance of these ensembles to those of the ingredient models, and that of a standard ETAS as a benchmark. Preliminary results combining semi-parametric ETAS variants
suggest that component-driven ensembles can outperform all ingredient models, with similar
trends observed in Southern California and Italy. The novel model combination techniqueenables these results to reveal the strengths and weaknesses of the ingredient models, which
accelerates the development of better and more robust earthquake forecasts and aids theunderstanding of earthquake dynamics.

Keywords: ensemble modeling, forecasting, ETAS∗

Speakersciencesconf.org:statsei12:415952

Bayesian time-dependent earthquake forecasting
Christian Molkenthin∗†1 , Christian Donner2 , Gert Zöller3 , Sebastian Hainzl4 , Matthias
Holschneider3 , and Manfred Opper5
1

University of Potsdam, Institute of Mathematics – University of Potsdam, Karl-Liebknecht-Str. 24-25
D-14476 Potsdam, Germany
2
ETH Zürich, Swiss Data Science Center – ETH Zürich, Swiss Data Science Center, CH-8006 Zürich,
Switzerland
3
University of Potsdam, Institute of Mathematics – University of Potsdam, Institute of Mathematics,
Karl-Liebknecht-Str. 24-25 D-14476 Potsdam, Germany
4
GFZ German Research Centre for Geosciences, Germany – GFZ German Research Centre forGeosciences, Germany, Germany
5
TU Berlin, Artificial Intelligence Group – TU Berlin, Artificial Intelligence Group, D-10587 Berlin,
Germany

Abstract
The spatio-temporal Epidemic Type Aftershock Sequence (ETAS) model is one of themost successful models for time-dependent statistical modeling of earthquake occurrence.
Bayesian inference of the ETAS model allows a natural way to include prior knowledge andto quantify the uncertainty of the model estimates. However, using a Bayesian approach inthe context of the ETAS model is more complicated than using powerful standard maximum
likelihood methods for the model inference, and can be limited to small or intermediate
data sets. We show that a recently introduced Bayesian semiparametric ETAS version using
Gaussian process modeling can produce improved forecasts including uncertainty quantification compared with results obtained by the standard ETAS model. Examples of test sites
in California illustrate the advantages of the suggested Bayesian approach. In addition, we
present how the Bayesian inference can be accelerated which makes the approach suitablealso for larger data sets.

∗†

Speaker
Corresponding author: molkenth@uni-potsdam.desciencesconf.org:statsei12:417037

Seismic evidence of an aseismic slip event preceding
and following the 2017, Valparaiso, Mw 6.9
earthquake
Luc Moutote∗†1 , Zacharie Duputel2 , and Olivier Lengliné3
1

Institut Terre Environnement Strasbourg (ITES) – université de Strasbourg, Centre National de laRecherche Scientifique : UMR7063 – 5 rue René Descartes, 67084 Strasbourg, France2
Observatoire Volcanologique du Piton de la Fournaise (OVPF) – IPG PARIS – France
3
Institut Terre Environnement Strasbourg (ITES) – ITES, UMR 7063 CNRS, Université de Strasbourg– France

Abstract
Understanding the processes that lead to large earthquakes remains one of the key goals
in seismology. This is central to our understanding of seismic hazard, including the longcontroversial issue of earthquake predictability. Recently, Slow Slip Events (SSEs) have been
identified before a few large earthquakes, which might reveal the nucleation process of themain rupture. However, the current space and time resolution of geodetic measurement
makes detection and monitoring of SSEs difficult. SSEs are believed to drive their own seismicity. A fine monitoring of such seismicity may help monitor SSEs and understand theirlink with the rupture process. Here, we study the case of the 2017 Valparaiso earthquakethat was preceded both by an SSE and an intense seismicity. We build a high-resolution
catalog (Mc=2) of the region using a deep-learning phase picking model and by locating
detected earthquakes in a local 3D Earth model. The resulting catalog consists of more than
100 000 events from 2015 to 2021 (compared with the ˜8000 events reported by the Centro
Sismológico Nacional over the same time-period). The temporal evolution of the seismicity
is tested against a modified Epidemic Type Aftershock Sequences model that accounts for
short-term incompleteness (Hainzl 2021). We show that an over-productive earthquake rate
emerges within the foreshock sequences, continuing at least for several days after the mainshock. This transient background activity seems unaffected by large magnitudes observed
during the sequence. In addition, we identify several repeater events between 2016 and 2021
using waveform cross-correlation and double difference relocation. The strongest repeater activity is observed within the foreshock seismicity and seems to decelerate over a few months.
These observations support the idea that a part of the seismicity surrounding the mainshockdiffers from the usual earthquake triggering processes. The joint observation of an SSE, atransient background seismicity and a strong repeater activity initiating during foreshock
times suggest a close connection between the SSE and the seismicity. Such connection seems
to persist momentarily after the sequence, enhancing the aftershock activities. Therefore, we
believe that analyzing high resolution seismicity helps to monitor the activity of the rupture
process preceding and following large earthquakes, and can provide new insight for real-time
seismic hazard assessments.
∗†

Speaker
Corresponding author: lmoutote@unistra.fr

sciencesconf.org:statsei12:415167

Keywords: Foreshock, Slow Slip, Repeaters, Nuclation, ETAS, I, Short term Incompletness, High
resolution catalog, Rupture process

SeisBench: Benchmarking and applying deeplearning based phase pickers
Jannes Münchmeyer∗†11

Université Grenoble Alpes – Université Grenoble Alpes, Université Grenoble Alpes – France

Abstract
Seismic event detection and phase picking are the base of many seismological workflows.
In recent years, several publications demonstrated that deep learning approaches significantly
outperform classical approaches, even achieving human-like performance under certain circumstances. However, there are two shortcomings. First, given the differences in datasets
and evaluation metrics, it is not possible to compare the performance of different picking
models to each other. In particular, it is unclear how models perform when applied to data
outside the training region. Second, the models do not offer a standardised API and are therefore difficult to apply in production. Here, we present SeisBench - A toolbox for machinelearning in seismology. SeisBench is an open-source toolbox written in Python, providing
unified access to benchmark datasets, machine learning models, and training pipelines. Using SeisBench, we conducted a large scale benchmark of deep learning based phase pickers.
We evaluated six picking and detection algorithms on eight datasets. We’ll present selected
results of this benchmark here. In addition to the benchmark, we’ll show spotlights how
SeisBench can be applied for seismological analysis. This includes the application to ocean
bottom data, the creation of highly complete regional earthquake catalogs, and the usage inrouting teleseismic monitoring at GEOFON.

∗†

Speaker
Corresponding author: jannes@munchmeyer.desciencesconf.org:statsei12:416987

A Deep Gaussian Process Model for Seismicity
Background Rates
Jack Muir∗1 and Zachary Ross2
1

2

University of Oxford – United Kingdom
California Institute of Technology – United States

Abstract
The spatiotemporal properties of seismicity give us incisive insight into the stress stateevolution and fault structures of the crust. Empirical models based on self-exciting pointprocesses continue to provide an important tool for analysing seismicity, given the epistimic uncertainty associated with physical models. In particular, the epidemic-type aftershock sequence (ETAS) model acts as a reference model for studying seimicity catalogs.
The traditional ETAS model uses simple parametric definitions for the background rate oftriggering-independent seismicity, which limits its applicability to the temporally complex
seimicity of seismic swarms that are not dominated by triggering. In order to robustly
capture time-varying seismicity rates, we introduce a deep Gaussian process formulation
for the background rate as an extension to ETAS. Gaussian processes (GPs) are a robust
non-parametric model for function spaces with covariance structure. By conditioning thelengthscale structure of a GP with another GP, we have a deep-GP: a probabilistic, hierarchical model that automatically tunes its structure to match data constraints. We illustrate
our method using synthetic examples and real data from the Cahuilla swarm and Ridgecrest
aftershock sequence.

∗

Speakersciencesconf.org:statsei12:412122

Understanding the seismic behavior of a complex
fault structure through the use of a novel 3D
stochastic declustering algorithm
Maura Murru∗1 , Rodolfo Console1 , Giuseppe Falcone1 , Matteo Taroni1 , and CaterinaMontuori1
1

Istituto Nazionale di Geofisica e Vulcanologia - Sezione di Roma (INGV) – Via di Vigna Murata, 605 00143 Roma, Italy

Abstract
The Alto Tiberina Fault (ATF) system in Central Italy (Northern Apennine) is a complex
fault structure composed of a principal low-angle fault and other minor synthetic/antithetic
faults located in the hanging wall of the principal fault. In this work, we developed a novel 3D
stochastic declustering algorithm based on the ETAS model, in order to study the clustering
properties of the seismicity by considering the earthquakes’ depth. This algorithm wasapplied to a new high-quality instrumental catalog (TABOO catalog) spanning from 2010
to 2015. We found two distinct behavior for the seismicity in the principal fault and thesynthetic/antithetic faults, both in terms of clustering properties and magnitude frequency
distribution. The principal ATF shows a high b-value and low clustering (i.e. a majorityof spontaneous events), while the synthetic/antithetic faults show a low b-value and highclustering (i.e. a majority of triggered events). These results also imply a different magnitude
frequency distribution for spontaneous and triggered events, despite a similar rake of the focal
mechanisms (normal events). Such a complex fault structure, therefore, shows peculiar andunusual seismic properties.

Keywords: Alto Tiberina Fault (ATF) system, Central Italy, Stochastic declustering, ETAS model∗

Speakersciencesconf.org:statsei12:416515

A new Bayesian approach to modelling of short andlong term earthquake forecasting using inlabru
Mark Naylor∗1 , Kirsty Bayliss2 , Francesco Serafini2 , Farnaz Kamranzad2 , Finn
Lindgren3 , and Ian Main2
1

School of Geosciences [Edinburgh] – The Kingś Buildings, James Hutton Road, Edinburgh EH9 3FE,
United Kingdom
2
School of Geosciences [Edinburgh] – United Kingdom
3
School of Mathematics [Edinburgh] – United Kingdom

Abstract
Earthquake hazard maps in space and/or time provide a rational basis for long-term
hazard mitigation and shorter term operational forecasting respectively. Both require thedetermination of optimal combinations of model input parameters, and a realistic estimate of
the uncertainties involved. While Bayesian methods are well suited to the problem, the construction of posteriors using random sampling (e.g. MCMC) is prohibitively slow for fitting
the complex spatio-temporal point process models required, particularly if the parameters
are strongly correlated and we wish to include many relevant spatial covariates, such as maps
fault maps, or satellite-based deformation rates. The Integrated Nested Laplace Approximation (INLA) offers a potential solution to this problem, because it was designed exactly tohandle large spatial models with a range of spatially varying covariates and strongly correlated parameters. It is based on a deterministic method which is also very fast (compared to
MCMC methods), and has been applied successfully to analyse complex ecological models
involving point process modelling.
Here, we present examples of how this approach can be used for time-independent andtime-dependent seismicity forecasting. We show how diverse spatial covariates can be incorporated into forecasts, and how forecasts can be generated from synthetic catalogues basedon different posterior samples of the parameters to capture the full range of uncertainty. For
example, we show how the approach can be used to test the utility of maps of the Coulomb
stress redistribution after an earthquake for potential application in real time seismicity forecasting. We also show how surface geology maps can be used to distinguish between regions
where we have an ‘evidence of absence’ of faults from regions where there is an ‘absence ofevidence’ – and consequently treat these differently.
In a related abstract, Serafini builds on this work through the implementation of a Hawkes
process within inlabru that allows us to better model the self-exciting clustering of earthquakes in space and time.

Keywords: seismicity forecasting Bayesian patio, temporal∗

Speakersciencesconf.org:statsei12:416174

Imaging a Complex Earthquake Sequence in Sparta,
North Carolina, Eastern United States
Miguel Neves∗1 , Lindsay Chuang1 , Wei Li2 , Zhigang Peng1 , and Sidao Ni3
1

Georgia Institute of Technology (GT) – 311 Ferst Drive Atlanta, GA 30332-0340 USA, United States
2
University of Science and Technology of China (USTC) – University of Science and Technology ofChina No.96, JinZhai Road Baohe District Hefei, Anhui, 230026 P.R.China, China
3
State Key Laboratory of Geodesy and Earth’s Dynamics, Institute of Geodesy and Geophysics,
Chinese Academy of Sciences – China

Abstract
On August 9, 2020, a Mw5.1 earthquake ruptured the uppermost Crust near the townof Sparta, North Carolina. The earthquake ruptured a 2-km-long surface and it is likely the
first surface-rupturing event in the Eastern United States since the New Madrid earthquake
sequence (Figueiredo et al., 2022). Following the mainshock, 5 temporary seismic stations
were also deployed by scientists from CERI@Univ. Memphis (Horton et al., 2021) in thearea. Together with state-of-the-art earthquake detection and relocation techniques, thissequence provides an opportunity to better understand the local fault structures and thesource processes of intraplate earthquakes.
Here we combine a deep learning earthquake phase picker (EQTransformer, Mousavi et al.,2020) and a matched filter technique to compile a more complete earthquake catalog forthe sequence. We perform detections for the period of August 1 to November 18, 2020.
We first detect P and S waves and associate them using EQTransformer. Then we perform a matched filter detection to search for additional events using the combined USGS
and EQTransformer catalog as templates (498 templates). Using cross-correlation derived
differential-travel times, we relocate all the detected events to obtain a final catalog of 1339
earthquakes with high-resolution locations and a magnitude of completeness of about -1.5.
Our relocated catalog shows two major lineaments that appear to illuminate two fault
branches with a slight kink. These two lineaments are consistent with our focal mechanism solution using the generalized cut-and-paste (gCAP) method that shows a thrust eventwith some strike-slip component and a moment magnitude of 5.13. In addition, the non-DC
components of about 12%, implying the mainshock may have ruptured along complex fault
geometries. The two double-couple nodal planes are 349◦/69◦/116◦ and 115◦/32◦/41◦ for an
optimal centroid depth of 1.3 km. One of our identified lineaments is consistent with surface
observations of thrust faulting by Figueiredo et al. (2022) and the dip of our preferred nodal
plane. Analysis of first motion polarities of all events M> 1.5 shows variable focal mechanisms, consistent with weak and near-surface structures, but with mostly strike-slip/oblique
strike-slip mechanisms in our second lineament. Our results suggest that the mainshock rupture occurred primarily along a thrust fault and a secondary blind strike-slip fault. Theseobservations have potential implications on seismic hazard estimation in intraplate settings.
∗

Speakersciencesconf.org:statsei12:417199

Keywords: Earthquake detection and location, Intraplate earthquake sequence, Complex ruptureSecond-order smoothing method over the Delaunay
Tessellation
Yuanyuan Niu∗1 and Jiancang Zhuang2
1

The Graduate University for Advanced Studies – Japan
2
The Institute of Statistical Mathematics – Japan

Abstract
Smoothness prior is usually used in seismicity analysis and geophysical or geodetic Bayesian
inversion. Calculating smoothness over regular grids greatly simplifies the numerical calculation. When observed data are irregularly distributed, interpolation over the regular grids
increases the uncertainty of the inversion results, whereas smoothing over the Delaunay
Tessellation (DT) avoids the additional interpolation. However, the implementation of thenumerical calculation of the second-order derivate (smoothness) of a function, give its values on the DT nodes, is much difficult than implementing the calculation of the first-orderderivative (flatness). We proposed a Bayesian smoothing method over the DT, in which
not only the flatness prior but also smoothness prior is implemented. For calculating thesmoothness over the DT, we provide two methods: 1-layer method and 2-layers method.
Numerical experiments and an application on gravity inversion problem show the effectivityof our method.

Keywords: Delaunay Tessellation, smoothness prior, Bayesian analysis.

∗

Speakersciencesconf.org:statsei12:417189

Aftershocks Preferentially Occur in Previously
Active Areas
Morgan Page∗1 and Nicholas Van Der Elst1
1

U.S. Geological Survey – Pasadena, CA, United States

Abstract
The clearest statistical signal in aftershock locations is that most aftershocks occur
close to their mainshocks. More precisely, aftershocks are triggered at distances following a power-law decay in distance (Felzer & Brodsky, 2006). This distance decay kernel isused in Epidemic-type Aftershock Sequence (ETAS) modeling and is typically assumed to
be isotropic, even though individual sequences show more clustered aftershock occurrence.
The assumption of spatially isotropic triggering kernels can impact the estimation of ETAS
parameters themselves, such as biasing the magnitude-productivity term, alpha, and assigning too much weight to secondary rather than primary (direct) triggering. Here we showthat aftershock locations in Southern California, at all mainshock-aftershock distances, preferentially occur in areas of previous seismicity. For a given sequence, the scaling between
aftershock rates and the previous seismicity rate is approximately linear. However, the total
number of aftershocks observed for a given sequence is independent of background rate. We
explain both of these observations within the framework of rate-and-state friction (Dieterich,1994).

∗

Speakersciencesconf.org:statsei12:417015

A unified scaling law of seismic to aseismic moment
release from laboratory to nature.
Luigi Passarelli∗1 , Antonio Piio Rinaldi2 , Louis De Barros3 , Philippe Danre4 , Fréderic
Cappa5 , Paul Selvadurai2 , and Stefan Wiemer2
1

Istituto Nazionale di Geofisica e Vulcanologia - Sezione di Bologna (INGV) – Via Donato Creti, 12,40100 Bologna, Italy
2
Swiss Seismological Service, ETH Zurich (SED, ETH Zurich) – Switzerland
3Géoazur (GEOAZUR) – Université Nice Sophia Antipolis (UNS), INSU, Observatoire de la Coted’Azur, Institut de recherche pour le développement [IRD], CNRS : UMR7329 – Bât 1, 250 rue Albert
Einstein Les Lucioles 1, Sophia Antipolis 06560 VALBONNE, France
4
Geoazur – Université Côte d’Azur, CNRS, IRD, Observatoire de la Côte d’Azur, Géoazur – France
5
University of Nice Sophia-Antipolis – University of Nice Sophia-Antipolis – France

Abstract
Aseismic and slow slip transients are known mechanisms by which faults release theaccumulated strain. They occur in all tectonic settings but are abundantly detected andinvestigated at subduction zones. However, a growing number of new observations give theopportunity to study the slow slip transients at shallow crustal depth in extensional andtransform settings. Swarm-like earthquakes and non-volcanic tremor very often accompany
aseismic/slow slip transients. The investigation of the release of seismic and aseismic energy
has been successful to improve our understanding of the physical mechanics at play during
aseismic slip ruptures. Global dataset of slow slip events (SSEs) accompanied by tremor or
earthquake swarms indicate a power-law scaling between the seismic and aseismic moment
over several order of magnitude. The proposed power-law scaling depends on the tectonic
setting and characteristics depth of the events. In particular, shallow (< 10km) SSEs triggering swarm-like seismicity and outside subduction zones form a separating power-law scaling.
A physical interpretation invokes a fault pressurization mechanism that creates the conditionof nucleating both aseismic and seismic slip where the seismicity is driven by the stressing
rate imparted at the leading edges of the aseismic rupture front. A similar mechanism isinvoked to explain aseismic and seismic slip during injection operation on fault for induced
seismicity during in-situ experiments. We here extended the scaling of crustal SSEs to thescale of induced seismicity, in-situ, and laboratory earthquake. We find that the power-law
scaling remarkably extends to lower magnitude slow slip transients up to the laboratory scale,although an observational gap remains for SSEs between Mw 0-4. In order, to fill in the gap
we use numerical simulations of SSEs and associated seismicity coupling a poro-elastic modelwith a stochastic earthquake simulator. The simulated distributions of aseismic and seismic
moments fit well with the natural and induced cases preserving the scaling. Our findings
indicates that the scaling is compatible with the hypothesis that stress changes due to theaseismic slip itself is the main driving mechanism, while fluid fault pressurization remains
the required condition for nucleation of aseismic slip.
∗

Speakersciencesconf.org:statsei12:417329

Temporal variations of seismicity rates andGutenberg-Richter b-values: an example from ahigh-definition Italian catalog
Anna Eliana Pastoressa∗1 , Maura Murru1 , Matteo Taroni1 , Rodolfo Console1 , CaterinaMontuori1 , Giuseppe Falcone1 , Raffaele Di Stefano1 , and Lauro Chiaraluce11Istituto Nazionale di Geofisica e Vulcanologia (INGV) – Via di Vigna Murata 605 - 00143 Rome, Italy

Abstract
The high-resolution earthquake catalog related to the seismic activity of the Alto Tiberina fault system (Northern Apennine, Central Italy), constructed from the seismic networkof the Alto Tiberina Near Fault Observatory (TABOO), were analysed by statistical methods
from a seismological perspective.
With the aim to bring new information on the study of earthquakes preparatory phases andin order to identify possible tectonic stress indicators, we focused on the reliable time seriescomputation of the seismicity rates and of the Gutenberg-Richter b-value.
In this work, the time series were evaluated starting from the generation of the stochastic declustered catalog (Console, Jackson & Kagan, Pageoph 2010) obtained using the 2D
ETAS Model.
The capability of the algorithm based on ETAS model to recognize and distinguish thetriggered seismicity from the background (e.g. Console & Murru 2001; Console, Murru &Lombardi 2003), allowed to evaluate the daily seismicity rates and b-value variations for boththe spontaneous and triggered earthquakes.
Moreover, in order to overcome the problems due to build b-value time series for a nonuniform temporal distribution of seismicity, we considered the weighted likelihood approach
for a more robust b-value estimation (Taroni, Vocalelli & De Polis 2021).
The obtained time series, in particular the b-value time series, are different with respect
to the ones obtained using the classical rolling window approach for the b-value computation, indicating that some of the b-value fluctuations are related to the natural statisticalvariability of the parameter.

Keywords: time, series analysis, b, value, statistical seismology, Alto Tiberina fault system.

∗

Speakersciencesconf.org:statsei12:416534

An automatic polarity determining method of theseismic recordings based on Information theory
Weilai Pei∗1,2 , Zhou Shiyong2 , and Jiancang Zhuang3
1

2

Sinopec Petroleum Exploration and Production Research Insititute – China
School of Earth and Space Sciences [Beijing] – Peking University, Haidian District, Beijing, China
3
The Institute of Statistical Mathematics – Japan

Abstract

Arrival detection and polarity determination are crucial and fundamental problems in seismic data
analysis. We develop a new method that calculates probability distribution funtion (pdf) of arrival times
based on Information theory. The initial arrival time t of the seismic wave is decided by independent parameter thatrepresentsamplitudethresholdinwavef ormsegment.Af terapplingthresholdselectioncriterion, webuildama
motionpolarityof theseismicwavef orm.W etestthismethodwithsyntheticwavef ormsaswellasrecordingsf romverif iedda

Keywords: Focal mechanism inversion, Polarities determination, Automatic Picking, Markov chain,
Information theory∗

Speakersciencesconf.org:statsei12:417165

Earthquake productivity within general ETAS
modelsAntonella Peresan∗1 , George Molchan2 , and Elisa Varini†31

Istituto Nazionale di Geofisica e di Oceanografia Sperimentale (OGS) – via Treviso 55 - 33100 - Udine,
Italy
2Institute of Earthquake Prediction Theory and Mathematical Geophysics (IEPT) – Profsoyuznaya str.
84/32, Moscow 117997 Russia, Russia
3
Consiglio Nazionale delle Ricerche [Milano] (CNR) – Via Bassini, 15 - Via Corti 12 - 20133 Milano,
Italy

Abstract
Earthquake clustering is a prominent feature of seismic catalogs. Several models havebeen proposed to describe the triggering mechanism of earthquakes; nonetheless, some theoretical aspects have not yet been thoroughly explored. In this study we investigate the
distribution of the number of triggered events as described by a branching process, whichmight provide useful practical constraints for empirical analysis of seismic catalogs.
According to recent literature (e.g. Shebalin et al., 2020, and references therein), the productivity of a magnitude m event can be defined as the number of triggered events of magnitude
above m-D, where D is a positive default value. In particular we distinguish between thenumber of direct descendants and the number of all descendants, denoted respectively by
the random variables v(D) and V(D).
In the standard Epidemic Type Aftershock Sequence (ETAS) model, the distribution of v(D)is Poissonian. However, evidence has recently emerged in favor of the discrete exponential
distribution for v(D) and for V(D) with a dominant initial magnitude m (the case of a simple
aftershock cluster). Therefore we consider the general ETAS model adapted to any distribution of v(D) and prove that the branching structure of the model excludes the possibilityof having a common distribution type for both v(D) and V(D) at once (Molchan et al., 2022).
We then investigate the features of the V(D) distribution within a wide class of ETAS models. We show that there is a fundamental difference in tail behavior of the V(D)-distributions
for general-type clusters and for clusters with a dominant initial magnitude: the tail is heavy
in the former case and light in the latter. The real data display similar behavior. Theoreticalconclusions are also illustrated through the analysis of a synthetic earthquake catalog.
References
Molchan, G., Varini, E., Peresan, A. (2022). Productivity within the ETAS seismicity model.
arXiv preprint arXiv:2203.05253∗†

Speaker
Corresponding author: elisa.varini@cnr.it

sciencesconf.org:statsei12:416873

Ogata, Y. (1998). Space-time point-process models for earthquake occurrences. Ann. Inst.
Stat. Math., 50, 379-402.
Shebalin, P.N., Narteau, C., Baranov, S.V. (2020). Earthquake productivity law. Geophys.
J. Int., 222, 1264-1269, doi: 10.1093/gji/ggaa252

Defining and testing the predictive efficacy of arealistic spring block model
Giuseppe Petrillo∗1 , Alberto Rosso2 , and Eugenio Lippiello3
1

2

The Institute of Statistical Mathematics (ISM) – JapanLPTMS, CNRS, Univ. Paris-Sud – Centre National de la Recherche Scientifique - CNRS – France
3
Università della Campania ”Luigi Vanvitelli” – Italy

Abstract
The occurrence of earthquakes is a singular characteristic of every part of the globe, but
is characterized by statistical laws that are universal, regardless of the geographical location.
The temporal and spatial correlations are a striking example.
Here we propose a spring block-like model coupling an elastic layer, subject to a velocity
weakening heterogeneous friction, to a second layer where friction is velocity strengthening.
We show that these ingredients allow us to recover all the most relevant earthquake statisticalfeatures, such as the Gutenberg-Richter law, the Omori law, the productivity law and thespatial clustering of aftershocks. The extensive simulations that can be obtained with thismodel, allow us to test the seismic gap hypothesis. Indeed, the estimate of the occurrence
probability of the next big shock on the basis of the time delay from the last earthquake still
represents a big challenge. The problem is that this issue cannot be addressed only on thebasis of historical catalogs which contain too few well documented big shocks, and decades of
future observations appear necessary. Our results support recent findings which have shownthat important insights can be obtained from the spatial organization of aftershocks andits relationship to the mainshock slip profile. In particular, we show that large earthquakes
do not regularly repeat in time, but it is possible to achieve insights on the time until thenext big shock from the percentage of aftershocks occurring inside the mainshock slip contour.

∗

Speakersciencesconf.org:statsei12:414005

Integrated monitoring of fault systems
Piero Poli∗1
1

ISTerre, UMR 5275 CNRS – ISTerre lab., University Savoie Mont-Blanc - CNRS, Le Bourget-Du-Lac,
France – Université Grenoble Alpes ISTerre CS 40700 38058 GRENOBLE Cedex 9, France

Abstract
Faults, or more in general fault systems, are complex geological objects. This complexity
controls the style of seismicity, and is often convoluted with external stress perturbations (e.g.
static and dynamic stress, variation of pore pressure, tidal stress perturbation). Taken together, this stress and rheological complexity make it hard to untangle the physical processes
controlling the seismicity evolution in time and space style, and limits our understanding of
how large earthquakes or seismic swarms generate.
In this work we integrate different geophysical observables to get insights about the different
physical processes occurring in the Alto Tiberian fault (ATF) system in the central Apennines (ITALY). We first present a catalog containing up to half-millions microearthquakes
over a period of ˜3yrs. The analysis (statistics interaction, repeating earthquakes, stressdrop) of this seismicity reveals the existence of various style of deformation in the ATF.
Combining the seismicity evolution with geodetic observations, we infer the rheological parameters of the shallow fault system, which hosts aseismic deformation and large swarms.
We further use observations of velocity variations estimated from ambient noise correlation,
to reveal a significant stress sensitivity of the shallow crust. We finally put some thresholdson the stress perturbation required to initiate seismicity in the ATF, by analysis of dynamic
triggering.
This integrated analysis, permits to derive significant new insights about the physical stateof the ATF, and highlight the importance of integrated analysis of geophysical data.

∗

Speakersciencesconf.org:statsei12:415853

An Overview of medium-term earthquake forecasting
using the EEPAS model: Where are we at?
Sepideh J Rastin∗†1 , David Rhoades1 , and Annemarie Christophersen1
1

GNS Science – 1 Fairway drive, Avalon, P.O. Box 30–368, Lower Hutt, New Zealand

Abstract
EEPAS stands for ”Every Earthquake a Precursor According to Scale” and is a New
Zealand born earthquake forecasting model. It works based on the precursory scale increase
phenomenon (Ψ) to forecast the future earthquake occurrence rates. Ψ refers to the observed
increase in the rate of occurrence and magnitude of minor earthquakes before most major
earthquakes. Each example of Ψ provides a value of mainshock magnitude Mm, precursor
time TP, and precursory area AP. Accumulated data of Ψ from many earthquakes shows thesequantities to be correlated. Predictive regressions of Mm, log TP, and log AP on MP arewhat EEPAS uses to forecast the contribution of any earthquake to the expectation of future
larger earthquakes. EEPAS is a medium-term forecasting model which can provide forecastswith time horizons ranging from months to decades. Since its introduction in 2004, EEPAS
has been successfully applied to well-catalogued regions including New Zealand, California,
Mexico, Japan and Greece. It has also been applied in studies of physics-based earthquake
simulators. EEPAS has been one of the core elements of NZ hybrid forecasting tool, andcontinues to be developed. A milestone in EEPAS improvement was compensation of theforecasts for missing precursory earthquakes in the time-lag between the end of the catalogue
and the forecasting time-horizon. This enables the forecasting horizon to be extended froma few months to a decade or more without much loss of forecasting quality (informationgain). We have also learned how to compensate EEPAS forecasts for the limited recordof precursory information before any target earthquake. Overall, EEPAS is better able toadjust for limitations of any earthquake catalogue than previously. However, there are still
significant challenges to move toward a global EEPAS model. One challenge is the existenceof a space-time trade-off in Ψ. This trade-off gets reflected in a corresponding trade-off of the
EEPAS scaling parameters for area and time. We have not yet mastered exploiting this tradeoff to optimise forecasting performance. There is also work to be done to address surprises
in low-seismicity regions, where recorded precursory seismicity can be sparse or seem to
be insignificant. Studies on synthetic catalogues, in which strain rates were systematicallyvaried, suggest much longer precursor times for such regions. Incorporating strain rate into
EEPAS would be another step forward towards a global EEPAS model.

Keywords: earthquake forecasting, precursory seismicity, New Zealand∗†

Speaker
Corresponding author: s.rastin@gns.cri.nz

sciencesconf.org:statsei12:417026

The rates of hazardous earthquakes in the New
Zealand region, and their uncertainties andvariabilities
Chris Rollins∗1 , Matt Gerstenberger1 , Annemarie Christophersen1 , David Rhoades1 ,Kirankumar Thingbaijam1 , Sepideh Rastin1 , and Russ Van Dissen1
1

GNS Science – New Zealand

Abstract
The 2022 New Zealand National Seismic Hazard Model (NSHM) consists of a Seismicity
Rate Model (SRM), a Ground Motion Characterization Model, and an uncertainty model.
The SRM itself has several components. For earthquakes on known active crustal faults,
the SRM is a ”grand inversion” (Field et al., 2014) that solves for the rates of earthquake
ruptures subject to several constraints, including a prescribed total magnitude-frequency distribution (MFD) and its uncertainty. This is blended with a spatially continuous distributed
seismicity model, which is also normalized to sum to a prescribed MFD at each magnitude.
The Hikurangi-Kermadec and Puysegur subduction zones are described by separate grandinversions, each constrained by its own total MFD among other inputs. To compute reliable MFDs for these uses, we make several modifications and additions to the New Zealandearthquake catalogue. We start from the Christophersen et al. (2022) catalogue, which homogenized the magnitudes of the GeoNet catalogue into an Mw-like magnitude. We perform
event-matching between this catalogue and several relocated and global catalogues, fromwhich we import improved earthquake depths, focal mechanisms, and improved magnitudes
and locations and missing events in the offshore subduction regions. We then use the improved earthquake depths, focal mechanisms, subduction interface geometries, and relative
plate-motion directions to classify earthquakes as upper-plate, interface or intraslab. Then,
we compute magnitude-frequency distributions in the 3D volumes covering the crustal grandinversions, the distributed seismicity models, and the Hikurangi-Kermadec and Puysegur
interfaces. While we use the more recent part of the catalogue (generally 1965-present)
to compute b-values, our estimate of earthquake rates includes events back to 1840 andincorporates magnitude uncertainties and rounding errors, depth uncertainties, and changing completeness thresholds through time. Perhaps more importantly, it also includes atreatment of the ambiguity in earthquake rates over decadal to centenary timescales thatarises due to alternating periods of activity and quiescence in the instrumental and historic
catalogue.

Keywords: catalogue, magnitude, frequency distribution, earthquake rate∗

Speakersciencesconf.org:statsei12:415549

Geometrical Properties of Seismicity in California
Zachary Ross∗†1 , Yehuda Ben-Zion2 , and Ilya Zaliapin3
1

California Institute of Technology – United States
University of Southern California – United States
3
University of Nevada, Reno – United States

2

Abstract
The spatial geometry of seismicity encodes information about loading and failure processes, as well as properties of the underlying fault structure. Traditional approaches tocharacterizing geometrical attributes of seismicity rely on assumed locations and geometryof fault surfaces, particularly at depth, where resolution is overall quite poor. In this study,
we develop an alternative approach to quantifying geometrical properties of seismicity using techniques from anisotropic point process theory. Our approach does not require priorknowledge about the underlying fault properties. We characterize the geometrical attributesof 32 distinct seismicity regions in California and introduce a simple four class classification
scheme that covers the range of geometrical properties observed. Most of the regions classified as having localized seismicity are within northern California, while nearly all of theregions classified as having distributed seismicity are within southern California. In addition,
we find that roughly 1 out of 4 regions exhibit orthogonal seismicity structures. The resultsof this study provide a foundation for future analyses of geometrical properties of seismicity
and new observables to compare with numerical modeling studies.

∗†

Speaker
Corresponding author: zross@caltech.edu

sciencesconf.org:statsei12:415782

Bayesian ETAS for Improved Earthquake Rate
Models for the Pacific Northwest
Max Schneider∗1,2 and Peter Guttorp2,3
1

United States Geological Survey – United States
2
University of Washington – United States
3
Norwegian Computing Center – Norway

Abstract
The Pacific Northwest (PNW) has substantial earthquake risk, both due to the offshore
Cascadia megathrust fault and other fault systems that produce earthquakes under theregion’s population centers. Forecasts of aftershocks following large earthquakes are thus
highly desirable and require statistical models calibrated to a catalog of the PNW’s pastearthquakes and aftershock sequences. This is complicated by the fact that the PNW contains multiple tectonic regimes hypothesized to have different aftershock dynamics as well as
two types of earthquake clustering (aftershock sequences and swarms). We use the EpidemicType Aftershock Sequence (ETAS) model to describe the characteristics of earthquakes andaftershocks for the PNW, accounting for these different types of seismicity. Typically, maximum likelihood estimation is used to fit ETAS to an earthquake catalog; however, the ETAS
likelihood suffers from flatness near its optima, parameter correlation and numerical instability, making likelihood-based estimates less reliable. We present a Bayesian procedure for
ETAS estimation, such that parameter estimates and uncertainty can be robustly quantified, even for small and complex catalogs like the PNW. The procedure is conditional on
knowing which earthquakes triggered which aftershocks; this latent structure and the ETAS
parameters are estimated iteratively. The procedure uses a Gibbs sampler to conditionallyestimate the posterior distributions of each part of the model. We simulate several synthetic
catalogs and test the modelling procedure, which produce well-mixed posterior distributions
centered on true parameter values. We also use the procedure to model the earthquakesof the continental PNW, using a new catalog formed by algorithmically combining US andCanadian data sources and then identifying and removing earthquake swarms. We perform
a completeness analysis that supports two complete subcatalogs split by latitude, and withdiffering start years and magnitudes of completeness. While MLEs are unstable and dependon both the optimization procedure and its initial values, Bayesian estimates are insensitive
to these choices. Bayesian estimates also fit the subcatalogs better than do MLEs. We use
the Bayesian method to rigorously estimate ETAS parameters and their uncertainty when
including swarms in the model, modelling across different tectonic regimes and complete
subcatalogs, as well as from catalog measurement error. Many parameter estimates change
substantially when considering these catalog issues, indicating their importance for seismicity
rate modelling and aftershock forecasting in the PNW.

Keywords: ETAS, Bayesian statistics, aftershock forecasting, Pacific Northwest∗

Speakersciencesconf.org:statsei12:417339

Periodic modulation of seismicity - a laboratory andmodeling perspective
Alexandre Schubnel∗1 , Martin Colledge2 , Jérôme Aubry3 , Kristel Chanard4 , François
Petrelis5 , Clara Duverger3 , and Laurent Bollinger3
1

Laboratoire de géologie de l’ENS (LGE) – INSU, CNRS : UMR8538, École normale supérieure [ENS] Paris – 24 Rue Lhomond 75231 PARIS CEDEX 05, France
2
LGENS – CNRS, UMR.8538 – Laboratoire de Géologie, ENS-PSL Research University, 24 rue
Lhomond, 75231 Paris, France – France
3
CEA DASE – CEA DASE – France
4
IPGP – IPGP Paris – France
5
LP ENS – Laboratoire de Physique de l’ENS, Ecole Normale Supérieure, PSL, Université Paris
Diderot, Sorbonne Université, 75005, Paris – France

Abstract
Small transient stress perturbations are prone to trigger earthquakes and a better understanding of the dynamic of earthquake triggering by transient stress perturbations is essential
in order to improve our understanding of earthquake physics and our consideration of seismic
hazard. In the Earth’s crust, these transient stress changes can be caused by various sources(passing of seismic waves, forcing by tides, hydrological load, and other seasonal climatic
loads). Recent observations suggest that at first order, seismicity increases with the rate of
tectonic loading while the modulation of the seismicity, increases both with the period andthe amplitude of the disturbances. Slow earthquakes and micro-seismicity also appear to
be more sensitive to tides, so the modulation seems to depend on an effective normal stressdisturbance rather than on the absolute value of the disturbance itself.
Here, we first study an experimental fault system, which mimics the natural seismic cycle,that is to say periodic perturbations of the stresses (of variable amplitude and frequency)superimposed to a far field loading, at crustal pressure conditions. The experimental system
is instrumented so that micro-seismicity (acoustic emissions), deformations and stresses, arecontinuously recorded during the experiments, in order to study the response of a systemproducing micro-seismicity as a function of loading rate, amplitude and frequency of a periodic stress perturbation. The susceptibility of the system to mechanical parameters (e.g.
shear stress oscillations, friction oscillations, ...) is derived from the amplitude of the distribution of acoustic emissions normalized by the amplitudes of the parameters in question.
Susceptibility to the confinement pressure oscillations, to the shear stress oscillations, and to
the friction oscillations all exhibit trends that suggest more complex mechanical behaviours
than a simple Mohr-Coulomb-based criterion.
Second, we investigate theoretically the effects of periodic-in-time modulations on the properties of earthquakes by considering successively the one dimensional Burridge–Knopoff (BK)model and the two dimensional Olami–Feder–Christensen (OFC) model. Each model is modified to take into account either a modulation of normal stress or of shear stress acting on∗

Speakersciencesconf.org:statsei12:421379

a fault. Despite the differences between the BK and the OFC model, several results areobserved in both models. In particular, we observe that earthquake occurrences correlatewith stress modulation. The correlation is strongly dependent on parameters such as thetype of modulation, its frequency and amplitude, and in some cases on the magnitude of theconsidered earthquakes.

The physics of earthquake forecasting: challenges
and recent advances
Margarita Segou∗1
1

British Geological Survey – United Kingdom

Abstract
The development of short-term earthquake forecasts is critical for the design of risk mitigation actions and raising social awareness during evolving earthquake sequences. In the last
decade statistical and physics-based models, or their combination, were proven successful intracking the spatio-temporal evolution of different sequences around the world. The enhancement of earthquake catalogs and the on-the-fly availability of otherwise time-intensive data
products, such as seismic source models, provide the necessary input for the improvement ofphysics-based models using continuum mechanic principles. Although the golden standard
for short-term earthquake forecasting remains the statistical ETAS (Epidemic-Type Aftershock Sequence) model, physics-based, Coulomb Rate-and-State (CRS) forecasts have now
reached, and sometimes even exceed, the performance of their statistical counterparts.
The ETAS and CRS forecast models face different challenges: all-positive statistical forecasts predict seismicity at no-occurrence locations while physics-based models underpredict
observed seismicity at the near-source area around the onset of earthquake sequences. Interestingly, both models improve their performance when high-resolution datasets are incorporated and enhanced implementations are pursued. Recent forecasting experiments in
Italy and California bring us closer to the limits of our predictability but they also challenge
us to evaluate our long-standing assumptions. The static stress transfer theory has known
different implementations driven by varying assumptions but remains puzzling how stressredistribution works at a granular, sub-km, spatial level with combined static and dynamic
interactions. Additionally, well-standing empirical relationships can be now tested more extensively, e.g. is event magnitude the controlling factor of its productivity and, what is thenature of frequency-magnitude distributions at small-magnitude ranges? Outside of the governing principles and assumptions behind our forecasts, our experimental design strategy has
also changed these last few years. The initial, 5-year-long forecast experiments over statewide testing regions have provided scientists with the critical infrastructure required andmotivated the implementation of sequence-specific experiments. It remains open whether
more specific performance metrics should be developed to investigate specific aspects behindthe performance of our current parametric models.
From the above, it is clear that the earthquake forecasting discipline of seismology has
gone a long way in the last decade. Recent advances in machine-learning-based earthquake
catalog development suggest that the capability of short-term earthquake forecasting moves
now closer to operational time-frames worldwide, but at the same time the underlying mechanisms of localized earthquake triggering for large-magnitude occurrences remain a mystery.
We can anticipate that machine learning will revolutionize earthquake forecasting by both∗Speakersciencesconf.org:statsei12:421372

providing very rich ML-based datasets that spark the scientific curiosity toward discovery,
and by enabling the spontaneous discovery itself through the development of non-parametricmachine-learning-based forecast expressions.

Unsupervised Artificial Intelligence and Seismicity
Declustering : Application to the Corinth Rift area
(central Greece)
Antoine Septier∗1 , Alexandra Renouard2 , Jacques Deverchère1 , and Julie Perrot1
1

UBO – UMR 6538 LGO, IUEM-UBO – France
2
Imperial College London – United Kingdom

Abstract
The intracontinental Corinth Rift is one of the fastest extending and most seismically
active regions worldwide (Bernard et al., 2006). It represents an outstanding land-accessible
natural laboratory to better understand normal fault mechanics and continental tectonics.
The western Gulf of Corinth is the most active part, hosting an extensive microseismic
activity with intense seismic swarming. In this context, high-quality catalogues of earthquake
clusters are critical to further study the tectonic complexity of such swarm-prone regions.
Here, we introduce a new declustering technique using a self-organising map (SOM), an
unsupervised machine learning clustering approach based on an artificial neural network(Vettigliet et al, 2018). The SOM neural network performs a non-linear mapping of large
input spaces into a two-dimensional grid, which hopefully preserves the topological andmetric relationships of the data and gives a better representation of the resulting clusters. It
therefore provides easier visualisation and interpretation of the results. The clusters identified
by this method can then be studied to allow a probabilistic classification of the seismic events
as background event, aftershock or swarm element, based on the features we previously used
to train the SOM network (inter-event space-time distances and waveform similarity). To
estimate the classification confidence level, we use the neighbourhood structure of the SOM
grid to train a second SOM network that clusters data points according to their topologicallocation. Preliminary results show classification confidence of about 80% for the ten-year
seismic catalogue of the National Hellenic Seismological Network. This method has thepotential to be applied in different tectonic settings and in large geographical areas withcomplex spatio-temporal earthquake distributions and behaviours.

∗

Speakersciencesconf.org:statsei12:417149

Bayesian time-independent and time-dependent
modelling of seismicity with Inlabru
Francesco Serafini∗1 , Finn Lindgren2 , Mark Naylor1 , and Ian Main1
1

University of Edinburgh, School of Geosciences – United Kingdom
2
University of Edinburgh - School of Maths – United Kingdom

Abstract
Earthquake hazard maps in space and/or time provide a rational basis for long-term
hazard mitigation and shorter-term operational forecasting respectively. Both require thedetermination of optimal combinations of model input parameters, and a realistic estimate of
the uncertainties involved. While Bayesian methods are well suited to the problem, the construction of posteriors using random sampling (e.g. MCMC) is prohibitively slow for fitting
the complex Spatio-temporal point process models required, particularly if the parameters
are strongly correlated and we wish to include many relevant spatial covariates, such as fault
maps, or satellite-based deformation rates. The Integrated Nested Laplace Approximation(INLA) offers a potential solution to this problem because it was designed to handle large
spatial models with a range of spatially varying covariates and strongly correlated parameters. It is based on a deterministic method which is also very fast (compared to MCMC
methods) and has been applied successfully to analyse complex ecological models involving
point process modelling. Here, using the R-package Inlabru, we present examples of how thisapproach can be used for time-independent and time-dependent seismicity forecasting. We
show how diverse spatial covariates can be incorporated into forecasts, and how forecasts canbe generated from synthetic catalogues based on different posterior samples of the parameters to capture the full range of uncertainty. In particular, we show how the approach can beused to test the utility of maps of the Coulomb stress redistribution after an earthquake forpotential application in real-time seismicity forecasting. We also show how surface geology
maps can be used to distinguish between regions where we have an ‘evidence of absence’ offaults from regions where there is an ‘absence of evidence – and consequently, treat these differently. In a related abstract, Serafini builds on this work through the implementation of a
Hawkes process within the R-package Inlabru that allows us to better model the self-exciting
clustering of earthquakes in space and time.

Keywords: Bayesian computing, Hawkes process, INLA, time, dependent seismicity forecasts, time,
independent seismicity forecasts∗

Speakersciencesconf.org:statsei12:416466

Ranking earthquake forecasts using proper scoring
rules: binary events in a low probability environment
Francesco Serafini∗1 , Mark Naylor1 , Finn Lindgren2 , Maximilian J. Werner3 , and Ian
Main1
1

University of Edinburgh, School of Geosciences – United Kingdom
2
University of Edinburgh - School of Maths – United Kingdom
3
University of Bristol - School of Earth Sciences – United Kingdom

Abstract
Operational earthquake forecasting for risk management and communication during seismic sequences depends on our ability to select an optimal forecasting model. To do this,we need to compare the performance of competing models in prospective experiments, andto rank their performance according to the outcome using a fair, reproducible and reliable method, usually in a low-probability environment. The Collaboratory for the Study ofEarthquake Predictability conducts prospective earthquake forecasting experiments aroundthe globe. In this framework, it is crucial that the metrics used to rank the competing
forecasts are ‘proper’, meaning that, on average, they prefer the data generating model. We
prove that the Parimutuel Gambling score, proposed, and in some cases applied, as a metricfor comparing probabilistic seismicity forecasts, is in general ‘improper’. In the special casewhere it is proper, we show it can still be used improperly. We demonstrate the conclusions
both analytically and graphically providing a set of simulation based techniques that canbe used to assess if a score is proper or not. They only require a data generating model
and, at least two forecasts to be compared. We compare the Parimutuel Gambling score’s
performance with two commonly used proper scores (the Brier and logarithmic scores) using
confidence intervals to account for the uncertainty around the observed score difference. We
suggest that using confidence intervals enables a rigorous approach to distinguish between
the predictive skills of candidate forecasts, in addition to their rankings. Our analysis showsthat the Parimutuel Gambling score is biased, and the direction of the bias depends on theforecasts taking part in the experiment. Our findings suggest the Parimutuel Gambling score
should not be used to distinguishing between multiple competing forecasts, and for care to
be taken in the case where only two are being compared.

Keywords: Probabilistic forecasting, Proper scores, Forecasts ranking∗

Speakersciencesconf.org:statsei12:416469

Deep-learning-based single-station-array tremors
detection
Leonard Seydoux∗1 and William Frank1
1

Department of Earth, Atmospheric and Planetary Sciences [MIT, Cambridge] (EAPS) – Massachusetts
Institute of Technology 77 Massachusetts Avenue Cambridge, MA 02139-4307, United States

Abstract
Relative to a single seismic station, seismic arrays (or networks of seismic stations) drastically improve our ability to detect and locate seismic sources in space. Traditional location methods use the relative arrival times across the array to locate the seismic source.
Waveform-based methods can simultaneously detect and locate seismic events by projecting
the coherent waveforms back towards the suspected source location. The seismic array’s
detection and location capacity depend mainly on the array density and the waveforms coherency across the array, which depends on the subsurface structure. The signature of seismic
sources consists of the convolution of the source time function and the impulse response of
the Earth’s structure between the source and the receivers. This implies that the source
location is encoded in the waveforms even if we are not yet able to explicitly disentangle thesignature of the source-receiver path from the seismic record. We implement an AI-based
framework to learn how to detect and locate seismic sources from the records of a single
seismic station based on the knowledge previously inferred from a dense seismic array. We
detect and locate tectonic tremors within the continuous seismic record of a single station in
the northern Cascadia subduction area. This proof-of-concept work shows how and underwhich assumptions the information from a seismic array of multiple stations can be decodedfrom the continuous seismic record of a single station. Our approach allows us to extend thetemporal coverage backward and forward of existing tremor catalogs without relying on theseismic arrays used to build the catalogs.

Keywords: Tectonic tremors, Cascadia, detection and location, deep learning∗

Speakersciencesconf.org:statsei12:417347

Magnitude and Slip Scaling Relations for
Fault-Based Seismic Hazard
Bruce Shaw∗1
1Lamont Doherty Earth Observatory, Columbia University, Palisades, New York, U.S.A – United States

Abstract
Scaling relations play an important role in fault-based seismic haard, where slip rates on
faults are balanced by cumulative slip in an earthquake rupture forecast. Here we discuss aseries of issues which arise in constructing sets of logic tree branches for use in fault-based
hazard estimates. These issues include:
1) depth of rupture; 2) transient deepening following large earthquakes; 3) uncertainty inmean; 4) difficulties in using classical continuum exponent fits; 5) multifault compatability;
6) and minimal models.
Recommendations for minimal branch models are discussed.

∗

Speakersciencesconf.org:statsei12:417163

Stochastic models for induced seismicity
Robert Shcherbakov∗1,2
1

Department of Earth Sciences, University of Western Ontario (UWO) – London, Ontario, N6A 5B7,
Canada2
Department of Physics and Astronomy, University of Western Ontario (UWO) – London, Ontario,
N6A 3K7, Canada

Abstract
Induced seismicity has become a significant factor in evaluating earthquake hazard associated with the recent development of unconventional energy resources. It is typically driven
by anthropogenic activities such as hydraulic fracturing, mining excavation, subsurface water injection, and others. Therefore, it is crucial for earthquake hazard mitigation and risk
assessment to develop stochastic models that can describe the observed induced seismicity
rate and forecast its evolution. In this work, I consider three such stochastic models to quantify the rate of the occurrence of induced earthquakes. For this, I generalize the Epidemic
Type Aftershock Sequence (ETAS) model, where the background rate is driven by the anthropogenic activities, for example, by the subsurface fluid injection. In the first model, thebackground rate is proportional to the fluid injection rate and also accounts for the time lag
associated with the fluid diffusion. The second model uses a convolution operation, where
the injection rate is convolved with a specific kernel to describe the influence of the fluid
propagation on the seismicity rate. Finally, the last model assumes that the effects of thefluid injection operations can be described by the rate-and-state mechanism. For all threemodels, the maximum likelihood method is developed to estimate the model parameters.
Subsequently, the models are incorporated into a Bayesian predictive framework, where theprobabilities for the occurrence of the largest expected earthquakes are computed during
specified forecasting time intervals. To compute the Bayesian predictive distribution theMarkov Chain Monte Carlo sampling of the posterior distribution is performed. This allows
the quantification of the uncertainties of the model parameters within the Bayesian framework. The three stochastic models are applied to several instances of induced seismicity inAlberta, Canada, and Oklahoma, USA. The main result of this work is the development of
a Bayesian approach for the estimation of the probabilities for the largest expected induced
earthquakes governed by the modified ETAS process. The proposed methodology can beused in the probabilistic assessment of the earthquake hazard associated with fluid injection
operations.

Keywords: ETAS, induced seismicity, Bayesian analysis∗

Speakersciencesconf.org:statsei12:416981

A Bayesian Framework for Earthquake/Aftershock
Forecasting and Testing
Robert Shcherbakov∗1,2
1

Department of Earth Sciences, University of Western Ontario (UWO) – London, Ontario, N6A 5B7,
Canada2
Department of Physics and Astronomy, University of Western Ontario (UWO) – London, Ontario,
N6A 3K7, Canada

Abstract
Earthquakes are driven by tectonic loading and typically trigger subsequent events. They
form clusters and swarms in space and in time. This is a direct manifestation of the nonPoisson behavior in the occurrence of earthquakes, where earthquake magnitudes and timeintervals between successive events are not independent and are influenced by past seismicity.
As a result, the distribution of the number of earthquakes is no longer strictly Poisson andthe statistics of the largest events deviate from the GEV distribution. In statistical seismology, the occurrence of earthquakes is typically approximated by a stochastic marked point
process. Among different models, the ETAS model is the most successful in reproducing
several key aspects of seismicity. Recent analysis suggests that the ETAS model generates
sequences of events which are not Poisson. This becomes important when the ETAS basedmodels are used for earthquake forecasting (Shcherbakov et al., 2019). In this work, I consider the Bayesian framework combined with the spatio-temporal ETAS model to constrain
the magnitudes of the largest expected earthquakes/aftershocks during a future forecastingtime interval. The formulated approach includes the MCMC sampling of the posterior distribution of the ETAS parameters and computation of the Bayesian predictive distribution
for the magnitudes of the largest expected events. It also allows full incorporation of the parameter uncertainties in computation of the relevant probabilities. To validate the forecasts,
the statistical tests developed by the CSEP are reformulated for the Bayesian framework.
In addition, I define and compute the Bayesian p-value to evaluate the consistency of theforecasted extreme earthquakes during each forecasting time interval. The Bayesian p-value
gives the probability that the largest forecasted earthquake can be more extreme than theobserved one. The suggested approach is applied to the recent 2019 Ridgecrest earthquake
sequence to forecast retrospectively the occurrence of the largest aftershocks (Shcherbakov,
2021). The results indicate that the Bayesian approach combined with the ETAS model
outperformed the approach based on the Poisson assumption, which uses the extreme value
distribution and the Omori law.
R. Shcherbakov, ”Statistics and forecasting of aftershocks during the 2019 Ridgecrest, California, earthquake sequence”, J. Geophys. Res. 126, (2021), doi: 10.1029/2020JB020887.
R. Shcherbakov, J. Zhuang, G. Zoeller, and Y. Ogata, ”Forecasting the magnitude of thelargest expected earthquake”, Nature Communications 10, (2019), doi: 10.1038/s41467-01911958-4.
∗

Speakersciencesconf.org:statsei12:415792

Keywords: ETAS model, Bayesian predictive distribution, extreme earthquakes

The 2020 Westmorland, California earthquake swarm
as aftershocks of a slow slip event sustained by fluid
flow
Krittanon Sirorattanakul∗†1 , Zachary Ross1 , Mostafa Khoshmanesh1 , Elizabeth
Cochran2 , Mateo Acosta1 , and Jean-Philippe Avouac1
1

California Institute of Technology (Caltech) – United States
2
United States Geological Survey (USGS) – United States

Abstract
Swarms are bursts of small magnitude earthquakes without an obvious mainshock. Somehave been observed to be associated with aseismic transients, while others are thought tobe related to fluids. However, the association is rarely quantitative due to the difficulty
measuring small geodetic signals related to swarms. We use high-quality GPS/GNSS, radar
interferometry (InSAR), relocated seismicity, and novel data processing techniques to study
a swarm of > 2,000 earthquakes that occurred between September 30 and October 6, 2020,
near Westmorland, California. Using 5-min sampled GPS supplemented with InSAR, wedocument a spontaneous shallow Mw 5.19 aseismic transient that preceded the swarm by 0.5
days. Seismicity in the early phase can be modeled as an inhomogeneous Poisson processsuggesting that they were non-interacting and driven primarily by the aseismic transient
rather than a product of cascade triggering. A stress-driven model based on rate-and-state
friction successfully explains the overall spatial and temporal evolution of seismicity, including the time lag between the aseismic transient and the swarm and the non-linear expansionof seismicity. Later in the sequence, a distinct back front and a square root of time expansion of clustered seismicity on the en-echelon structures suggest that fluids helped sustain
the swarm. Coulomb stress changes and statistics of interevent times suggest that 45 – 65%of seismicity was driven by the aseismic transient, 20 – 35% by inter-earthquake static stresstriggering, and 10 – 30% by fluids. Our model also provides constraints on the friction
parameter and the pore pressure and suggests that this swarm behaved like an aftershock
sequence but with the mainshock replaced by the aseismic transient. With geodetic data
available in real-time, the model could allow for forecasting the seismicity rate evolution ofswarms before they start and as they evolve.

Keywords: swarm, aftershocks, aseismic, fluid∗†

Speaker
Corresponding author: ksirorat@caltech.edu

sciencesconf.org:statsei12:417185

Statistical analysis to assess the skill of the
Operational Earthquake Forecasting system in Italy
Ilaria Spassiani∗†1 , Warner Marzocchi2 , Giuseppe Falcone1 , and Maura Murru1
1

Istituto Nazionale di Geofisica e Vulcanologia - Sezione di Roma (INGV) – Via di Vigna Murata, 605 00143 Roma, Italy
2
University of Naples, Federico II, Dept. Earth, Environmental, and Resources Sciences, Complesso diMonte Sant’Angelo, Via Cupa Nuova Cintia, 21, 80126 Napoli – Italy

Abstract
Statistical seismology represents a very active and productive research field, which integrates modelling, data analysis and forecasting in the study of the earthquake process. A
relevant application in this field is the Operational Earthquake Forecasting (OEF, Jordan
et al. (2011)), an experiment started by several Countries worldwide (Field et al. (2017),
Gerstenberger et al. (2014), Falcone et al. (2021)) to collect authoritative information about
the spatiotemporal dependence of seismic hazard. On the basis of stochastic models, thisexperiment provides the probabilistic forecast for the occurrence of a seismic event in a specific space-time-magnitude domain. Such information can both help stakeholders to establish
seismic risk reduction strategies, and enhance communities’ resilience to potentially destructive earthquakes.
The idea of the OEF experiment in Italy was carried out after the Mw 6.2 event that struck
on April 6th, 2009, beneath L’Aquila (central Italy), on proposal of the International Commission for Earthquake Forecasting (ICEF) for Civil Protection. In particular, the reportemphasized some important features not yet considered at that time, such as the need toprovide continuous information and the importance of testing the forecasts. According to
the standards proposed by the Collaboratory for the Study of Earthquake Predictability(CSEP, Schorlemmer et al. (2018)), the OEF-Italy experiment is set on a specific 0.1◦ x 0.1◦
grid lattice covering the whole Italian territory, and is based on a stochastic ensemble model
which consists in a combination of ETAS, ETES and STEP models (Lombardi & Marzocchi(2010), Falcone et al. (2010), Gerstenberger et al. (2005)), opportunely weighted according
to the relative past performances. Since 2013, the OEF-Italy system is operative in real-time,and produces time-dependent seismic maps of the weekly number of earthquakes, expectedby the ensemble model, occurring in each cell of the spatial grid above a fixed threshold
magnitude or macroseismic intensity.
As highlighted by the ICEF report in 2009, proper stochastic modeling and testing areindeed essential to deliver reliable seismic rates for societal purposes and, we add, to quantify our real level of knowledge about earthquake predictability. From this perspective, inthis work we perform some statistical analyses and tests to assess the forecasting skill of
the OEF-Italy system, by means of both classical methodologies (e.g., Zechar et al. (2010),
Taroni et al. (2018)) and new performance measures borrowed from other research fields,∗†

Speaker
Corresponding author: ilaria.spassiani@ingv.it

sciencesconf.org:statsei12:417224

like meteorology, specific to validate alarm-based systems by a binary criterion (forecast =yes or no, occurrence = yes or no; Casati et al. (2008), Roberts & Lean (2008)). Our final
aim is to: i) investigate possible weaknesses and room for improvements in the OEF-Italy
stochastic modelling, ii) get insights on potential proper probability thresholds to activate
state of alert procedures, iii) highlight possible features in the Italian tectonic seismic activity.

Main references
Casati, B., Wilson, L.J., Stephenson, D.B., Nurmi, P., Ghelli, A., Pocernich, M., Damrath,U., Ebert, E.E., Brown, B.G., and Mason, S. (2008). Forecast verification: current status andfuture directions. Meteorological Applications, 15(1):3-18, doi: https://doi.org/10.1002/met.52.
Falcone, G., Console, R. and Murru, M. (2010). Short-term and long-term earthquake occurrence models for Italy: ETES, ERS and LTST. Annals of Geophysics, 53(3):41-50. ISSN
2037-416X. doi: 10.4401/ag-4760
Falcone, G., Spassiani, I., Ashkenazy, Y., Shapira, A., Hofstetter, R., Havlin, S. and Marzocchi, W. (2021). An Operational Earthquake Forecasting Experiment for Israel: Preliminary
Results. Frontiers in Earth Science, 791
Field, E. H., Milner, K. R., Hardebeck, J. L., Page, M. T., van der Elst, N., Jordan, T.
H., Michael, A. J., Shaw, B. E. and Werner, M. J. (2017). A spatiotemporal clusteringmodel 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
Gerstenberger, M. C., Wiemer, S., Jones, L. M. and Reasenberg, P. A. (2005). Real-time
forecasts of tomorrow’s earthquakes in California. Nature, 435(7040): 328-331
Gerstenberger, M., McVerry, G., Rhoades, D. and Stirling, M. (2014). Seismic hazard modeling for the recovery of Christchurch. Earthquake Spectra, 30 (1), 17-29
Jordan, T., Chen Y., Gasparini P., Madariaga R., Main, I., Marzocchi, W., Papadopoulos,
G., Sobolev, G., Yamaoka, K. and Zschau, J. (2011). Operational Earthquake Forecasting.
State of Knowledge and Guidelines for Utilization. Ann. Geophys., 54(4). ISSN 2037-416X.
doi: 10.4401/ag-5350. URL: https://www.annalsofgeophysics.eu/ index.php/annals/article/view/5350
Roberts, N.M., and Lean, H.W. (2008). Scale-selective verification of rainfall accumulations from high-resolution forecasts of convective events. Monthly Weather Review, 136(1):78–97
Schorlemmer, D., Werner, M. J., Marzocchi, W., Jordan, T. H., Ogata, Y., Jackson, D.
D., Mak, S., Rhoades, D. A., Gerstenberger, M. C., Hirata, N., Liukis, M., Maechling, P. J.,
Strader, A., Taroni, M., Wiemer, S., Zechar, J. D. and Zhuang, J. (2018). The collaboratory for the study of earthquake predictability: Achievements and priorities. SeismologicalResearch Letters, 89 (4), 1305-1313.
Lombardi, A.M., Marzocchi, M. (2010). The ETAS model for daily forecasting of Italian seismicity in the CSEP experiment. Annals of Geophysics, 53(3): 155-164. doi:http://dx.doi.org/10.4401/ag4848
Taroni, M., Marzocchi, W., Schorlemmer, D., Werner, M. J., Wiemer, S., Zechar, J. D.,
Heiniger, L. and Euchner, F. (2018). Prospective CSEP evaluation of 1-day, 3-month, and5-yr earthquake forecasts for Italy. Seismological Research Letters, 89 (4), 1251-1261
Zechar, J. D., Gerstenberger, M. C. and Rhoades, D. A. (2010). Likelihood-based tests forevaluating space–rate–magnitude earthquake forecasts. Bulletin of the Seismological Societyof America, 100 (3), 1184-1195

Keywords: Operational earthquake forecasting, statistical analysis, performance skill

Targeting Temporal Seismicity with a Neural Point
Process
Sam Stockman∗1 , Dan Lawson1 , and Max Werner1
1

University of Bristol – United Kingdom

Abstract
Point processes have been the dominant way of describing the evolution of seismicity fora while, with the Epidemic Type Aftershock Sequence (ETAS) model being the most widelyused. Recent advances in the field of machine learning have constructed highly flexible
point process models using neural networks to seek improvements on existing parametricmodels. We conduct an investigation into the whether these flexible point process models
can be applied to short term seismicity forecasting by extending an existing temporal neuralmodel to the magnitude domain. We first find that the neural model can learn to capture
synthetic ETAS data, however, requiring less computational time due to the fact that it
is not dependent on the full history of the sequence. On data from a new high-resolution
earthquake sequence in the Central Apennines we investigate the sensitivity of both ETAS
and the neural model with respect to the lowest magnitude passed into each. We find that as
the magnitude is lowered the performance of ETAS reduces unlike that of the neural model.
We argue that some of the gains made by the neural model are due to its ability to handle
incomplete data. Since neither model has an overall forecasting advantage over the other, ourwork suggests that there is not enough signal in 2-dimensional data to give a neural modelan advantage over ETAS, but does offer a more efficient way to model data incompleteness.
Our work also stresses the need for considering temporal incompleteness when using new
machine learning generated earthquake catalogs.

∗

Speakersciencesconf.org:statsei12:415273

Using Hidden Earthquakes to Forecast Future Events
in Southern California
Matteo Taroni∗11

Istituto Nazionale di Geofisica e Vulcanologia (INGV) – Via di Vigna Murata 605 - 00143 Rome, Italy

Abstract
The spatial distribution of the past seismicity is often used to feed earthquake models,in order to forecast the zones that can be affected by earthquakes in the future. A prospective earthquake forecasting experiment in California showed that in the last 15 years thebest-performing model is the one that used epicenters of small earthquakes (M≥2) to define the earthquakes’ spatial distribution. In this work we performed a pseudo-prospective
experiment, focused on the seismicity of the southern California region from 2008 to 2021,
comparing the performance of two spatial seismicity models, based on a classical instrumental catalog and a modern template matching catalog. The results of such an experimentclearly demonstrate that by using ”hidden earthquakes”, i.e. very small events (M≥0.3)
from the modern catalog, it is possible to increase the forecasting performance of the spatialmodels. These results also suggest a profound similarity between the spatial distribution of
the earthquakes, from very small events to larger ones.

∗

Speakersciencesconf.org:statsei12:416969

A Testable Worldwide Earthquake Faulting
Mechanism Model
Matteo Taroni1 and Jacopo Selva∗†11Istituto Nazionale di Geofisica e Vulcanologia (INGV) – Italy

Abstract
We present a simple model to forecast global focal mechanisms. This model is based ona simple discrete counting distribution of the global centroid moment tensor catalog, and it
also includes, using a Bayesian scheme, the a priori informationfrom the Anderson theory of faulting. Our model is tested in hindcasting mode against independent data of global large earthquakes with M s ≥ 7. We obtained statistically significant
good agreement between model and data using consistency test, demonstrating that thissimple model can satisfactorily forecast focal mechanisms at the global scale. The defined
testing procedure can be used to test the model in prospective mode against future events.
These forecasts may inform short- to long-term hazard quantifications that require a finite
source characterization, as well as real-time source inversion algorithms. This model has
been recently published in SRL (https://doi.org/10.1785/0220200445).

Keywords: Focal Mechanism, Forecast, CMT∗†

Speaker
Corresponding author: jacopo.selva@ingv.it

sciencesconf.org:statsei12:417074

Resolving Differences in the Rupture Properties ofProminent California Earthquakes Using Bayesian
Source Spectral Analysis
Daniel Trugman∗1
1

Nevada Seismological Laboratory, University of Nevada, Reno (UNR) – United States

Abstract
The spectra of earthquake waveforms can provide useful insights into rupture processes,
yet the analysis and interpretation of earthquake spectra remains a challenging problem.
Here we develop a Bayesian framework that embraces the inherent data and modeling uncertainties of spectral analysis to infer key source properties. The method uses a spectral ratio
approach to correct the observed S-wave spectra of nearby earthquakes for path and siteattenuation. The objective then is to solve for a joint posterior probability distribution of
three source parameters – seismic moment, corner frequency, and high-frequency falloff rate
– for each earthquake in the sequence (both EGF and target), as well as a measure of rupture
directivity for all target events with good azimuthal station coverage. While computationallyintensive, this technique provides a quantitative understanding of parameter tradeoffs anduncertainties and allows one to impose physical constraints through prior distributions on
all source parameters, which guide the inversion when data is limited. We demonstrate themethod by analyzing in detail the source properties of 14 different target events of magnitude M5 in southern California that span a wide range of tectonic regimes and fault systems.
These prominent earthquakes, while comparable in size, exhibit marked diversity in theirsource properties and directivity, with clear spatial patterns, depth-dependent trends, and apreference for unilateral directivity. Coherent spatial variations in source properties suggestthat regional differences in tectonic setting, hypocentral depth or fault zone characteristics
may drive variability in rupture processes. This finding may have important implications forour understanding of earthquake physics and its relation to hazard.

∗

Speakersciencesconf.org:statsei12:415538

Spatial correlation of strain rate and crustal
seismicity in Japan
Taku Ueda∗1 and Takuya Nishimura1
1

Disaster Prevention Research Institute, Kyoto University – Japan

Abstract
Earthquakes are caused by the stress release in the earth’s interior due to the relative
plate motion. Therefore, it is considered that strain rate estimated from geodetic data may
correlate with seismicity, and this relationship has been discussed in several areas (Zeng etal., 2018; Stevens and Avouac, 2021; Xiong et al., 2021; Wu et al., 2021). Zeng et al., (2018)
revealed that the maximum shear strain rate along the San Andreas Fault and Eastern California Shear Zone positively correlated with the distribution of the background seismicity
rate in the shallow-crustal environment.
In Japan, there is a tendency that the higher strain rate, the more frequent shallow M > 6
earthquakes occur (Nishimura, 2017) and time-independent forecast model for large crustal
earthquakes based on geodetic data is suggested (Nishimura, 2022). In this study, we investigate the spatial correlation between background seismicity rate (smaller earthquakes) andstrain rate.
I used the JMA hypocenter catalog (the Preliminary Determination of Epicenters) as theearthquake catalog. I applied the HIST-ETAS model (e.g., Ogata, 2004) to the M > 3.0
crustal earthquakes which were located at < 25 km depth from 1980 through 2010 to estimate spatial variations in background seismicity rate. Then, I compare the spatial variations
in background seismicity rate estimated via the HIST-ETAS model with the maximum shear
strain rate estimated using the GNSS data from 2006 through 2009 by Okazaki et al. (2021).
Spatial variations in background seismicity rate positively correlate with that in the maximum shear strain rate. In the region including active volcanoes or whose lower crust has
fluid concentrated zone, background seismicity rate is higher at the same maximum shear
strain rate compared with other regions, suggesting that fluids weaken the fault strength.

Keywords: Background seismicity rate, Strain rate, HIST, ETAS model∗

Speakersciencesconf.org:statsei12:417230

Mapping spatio-temporal changes in the magnitude
frequency distribution of aftershocks with b-positive
Nicholas Van Der Elst∗1
1

United States Geological Survey (USGS) – United States

Abstract
The earthquake magnitude-frequency distribution is characterized by the b-value. Someobservations suggest that the b-value for aftershocks is larger than for spontaneous events(mainshocks), and that exceptions to this rule occur exclusively in foreshock sequences. Unfortunately, measurements of b-value are biased during aftershock sequences by short-term
incompleteness of the earthquake catalog. I present an alternative estimator of the b-value,‘b-positive’ that is insensitive to transient changes in catalog completeness and that allows
for the real-time, unbiased measurement of b-value during an ongoing sequence. The estimator is based on the differences in magnitude between successive earthquakes, which aredescribed by a double-exponential distribution with the same b-value as the magnitude distribution itself. The new estimator mostly confirms the findings of Gulia and others (Nature,
2019; SRL 2020 ), showing a decrease in b-value between the foreshock and mainshock inseveral prominent sequences over the past decade. More generally, sequences may show an
evolving magnitude-frequency distribution with time, with relatively low b-values early inthe sequence, and a gradual shift to higher b-values as the sequence progresses. Mapping the
spatio-temporal evolution of the 2020 SW Puerto Rico sequence reveals an interesting patternof lower b-values at the margins of the aftershock zone and in the vicinity of bursts of newactivity, followed by gradually increasing b-values as activity subsides. The non-stationary
distribution of magnitudes within a sequence may hold information about source processes
and has important implications for intermediate and long term aftershock forecasting.

∗Speakersciencesconf.org:statsei12:419799

Bayesian analysis of temporal variations of seismicity
based on non-extensive statistical mechanics
Elisa Varini∗†1 and Renata Rotondi1
1Consiglio Nazionale delle Ricerche [Milano] (CNR) – Via Bassini, 15 - Via Corti 12 - 20133 Milano,
Italy

Abstract
In 1988, Tsallis introduced a non-additive entropy that has found applications in various fields such as physics, chemistry, medicine, informatics, linguistics, and economics. The
concept also gained attention in statistical seismology due to the evidence that Earth’s crust
behaves as a complex, non-linear dynamic system characterized by long-range correlations.
By addressing the problem in the context of non-extensive statistical mechanics, the classicalmethods of analysis are unsuitable: the Boltzmann-Gibbs (BG) entropy is appropriate to
study systems with short range interactions, while new entropies that generalize some of theproperties of BG entropy are better suited to describe systems with long range interactions.
In this respect a new probability distribution of the magnitude was derived by maximizing
a non-extensive generalization of the Boltzmann–Gibbs entropy, namely the Tsallis entropy(Tsallis, 2009).
The result is the q-exponential distribution, whose shape parameter q is called the entropic
index and is found to characterize the subadditive (q> 1 ) and superadditive (q< 1 ) regimes;
the exponential distribution is recovered as q tends to 1.
Two active seismogenic areas in central Italy are analyzed, the former including the L’Aquila
sequence that occurred in 2009 and the latter the Amatrice-Norcia sequence that occurredin 2016. Parameter estimation is performed by following the Bayesian approach in order toexploit the prior knowledge on the phenomenon (Vallianatos et al., 2016) and to assess theuncertainties.
A detailed analysis of the variations of both q index and Tsallis entropy is performed byestimating them over time windows of a fixed number of events that shift at each new event.
Both the q index and the Tsallis entropy show significant and lasting decreases before thefirst strong earthquake in the sequences, and sudden increases after them. This indicatesthat these quantities offer clues on the activation state of the systems (Rotondi et al., 2022).
In the literature the q-exponential distribution has been used to study other seismic parameters like seismic moment, inter-event times, Euclidean distance between successive events.
Following the example of Schoenberg et al. (2009), we have thought that the q-exponential
distribution could describe, in addition to earthquake magnitude, also the spatial propertiesof seismicity such as the area of the Voronoi cells generated by the epicentral locations. We∗†

Speaker
Corresponding author: elisa@mi.imati.cnr.it

sciencesconf.org:statsei12:416985

extend our analysis to the areas of Voronoi cells generated by the recent seismicity recordedin the study regions. Then we compare the performance of the estimated q-distribution withthat of some promising probability distributions of the Voronoi cell area, namely the tapered
Pareto distribution and the generalized Gamma distribution. Also in this case the analysisof the time variations of the probability model that provides the best fitting to the datasets
provides us information on the phases of the seismic cycle.
Future work aims to analyze different seismic areas in order to test and consolidate theresults.
References
Rotondi R., Bressan G. and Varini E.; 2022: Analysis of temporal variations of seismicity through non-extensive statistical physics, Geophys. J. Int., 230, 1318-1337.
Schoenberg F.P., Barr C., and Seo J.; 2009: The distribution of Voronoi cells generated
by Southern California earthquake epicenters. Environmetrics, 20, 159-171.
Tsallis C.; 1988: Possible generalization of Boltzmann-Gibbs mechanics. J. Stat. Phys.,52, 479-487.
Tsallis C.; 2009: Introduction to Nonextensive Statistical Mechanics. Approaching a Complex World, Springer, doi:10.1007/978-0-387-85359-8.
Vallianatos F., Papadakis G. and Michas G.; 2016: Generalized statistical mechanics approaches to earthquakes and tectonics, Proc. R. Soc. A, 472, 20160497.

Bayesian Assimilation of Multicycle Earthquake
Simulations into Probabilistic Forecasting Models
Luis Vazquez∗1 and Thomas Jordan1
1

University of Southern California (USC) – United States

Abstract
A general problem in earthquake forecasting is how to assimilate deterministic physicalsimulations into probabilistic forecasting models. Here we focus on recalibrating the timeindependent Uniform California Earthquake Rupture Forecast Version 3 (UCERF3-TI) of
Field et al. (2014) against long earthquake catalogs ( ˜10ˆ6 years) generated by the multicycle Rate-State Quake Simulator (RSQSim) of Dieterich & Richards-Dinger (2010). We map
RSQSim ruptures from the Shaw et al. (2018) catalog onto equivalent UCERF3 rupturesby maximizing the mapping efficiency while preserving the seismic moment. We assume thesequence of equivalent UCERF3 ruptures is Poisson distributed, i.e., each rupture occurs at
a time-independent rate, our knowledge of which is uncertain. We use the full UCERF3 logic
tree to construct a joint prior distribution of rupture rates, which we represent by independent
Gamma distributions. Updating the UCERF3 gamma priors with the empirical RSQSim
Poissonian rate yields Gamma posterior distributions that can be calculated analytically.
Our results show that participation rates are decreased in the northern section of the San
Andreas Fault System (SAF) with respect to the UCERF3 model, which we attribute to theslip rate difference between RSQSim and UCERF3. Furthermore, our posterior model resultsin substantial lower rates in the Coachella, San Gorgonio, and San Bernardino SAF sections,reflecting the fact that RSQSim does not propagate ruptures through the San Gorgonio
knot. We also briefly address the more challenging problem of constructing time-dependent
forecasts that have been conditioned on knowledge of the previous rupture history over timeintervals of the last 100 years or so.

∗

Speakersciencesconf.org:statsei12:414435

Can we reproduce fault mechanics in the lab?
Jérôme Weiss∗1
1

Institut des Sciences de la Terre (ISTerre) – Université Grenoble Alpes, Université Savoie Mont Blanc,
CNRS, IRD, IFSTTAR, ISTerre, 38000 Grenoble, France – France

Abstract
Faults and earthquakes are generally considered as a paradigm for slowly driven physicalsystems that release the incoming energy through intermittent bursts spanning a wide rangeof sizes and energies. Owing to the (so far elusive) search for robust prediction tools, it has
been dreamed for a long time to reproduce faults and earthquake physics under controlled
conditions and detailed monitoring in the lab, in order to better understand the underlying
physics (1).
Earthquakes and faults dynamics are characterized by several scale-free statistics (2), suchas the distribution of their seismic moments (Gutenberg-Richter’s law) as well as temporaland spatial correlations, e.g. expressed through a decaying rate of aftershocks (Omori’s law)(3) and scale-free sub-diffusion (4). In the last decades, it has been also realized that tectonic
loading can be released, besides classical earthquakes, through a variety of mechanisms suchas aseismic slip, tremors, repeating earthquakes, ect..
Before to be proposed as an experimental model allowing to decipher these empirical observations, a laboratory analog should therefore reproduce this complexity. In other words,
a power law distribution of burst energies is a necessary, but not a sufficient condition toclaim a formal analogy. As an example, compressive rupture under monotonic loading,characterized by very clear precursory phenomena but scale-free statistics only at failure(5),cannot be considered as a good candidate, as previously claimed (1).
Recently, several experimental set-ups, based on the propagation of a quasi-brittle crack
within a disordered landscape (∗Speakersciencesconf.org:statsei12:416396

Analysis of the 2016 Central Italy earthquake
sequence by using a refined earthquake catalog
Louise Xiang∗1 and David Marsan∗†2
1

Xiang Louise Claire (ISTerre) – ISTerre USMB – Campus scientifique, 73376 Le Bourget du Lac cedex,
France2Marsan David (ISTerre) – ISTerre USMB – Campus scientifique, 73376 Le Bourget du Lac cedex,
FranceAbstract
By using an outstanding earthquake catalog obtained by machine learning method (Tanet al., 2021), that contains more than 900 000 relocated earthquakes, we aim at modeling the seismicity of the 2016 Central Italy earthquake sequence, and to find hints of thehigh-pressurized fluids located under the studied area by identifying small swarms occurring
during the sequence.
The magnitude of completeness mC is found to fluctuate significantly over the course of thesequence, and down to sub-daily time scales. Rather than defining one single value for mC,
we here propose to keep all earthquakes, and to account for changes in detection with time.
We define an earthquake stochastic model based on the ” Epidemic-Type-Aftershock-Sequence” (ETAS) model, by making it an 4-dimensional problem (x, y, z, t — m) and by introducing an earthquake detection probability. We assume that there is no background seismicity,
an earthquake is considered as a point source and each large earthquake has its own specific productivity prefactor. The parameters optimization of this model is done by using an”Expectation-Maximization ” (EM) algorithm that requires a splitting of triggering earthquakes according to the threshold magnitude mP, i.e. the minimum magnitude that an
earthquake need to reach in order to be modeled as potentially triggering other earthquakes.
By fixing ETAS parameters and exploring varying mP values, the results show that the accuracy of our models improves with decreasing mP value, as expected, but that this increasein the quality comes with over-fitting the data. We analyze this by assessing the prediction
capacity of the models: we separate the sequence in two time intervals (a ” training ” interval to fit the parameters, and a ” validation ” interval to measure the predictive valueof the model). This approach improves on the classical method that simply fits the model
parameters over the whole available dataset without preventing overfitting. We then analyze
how the observed seismicity rate locally departs from the best-predictive model so to isolate
possible swarm activity during the sequence.

Keywords: Earthquake modelling, aftershock sequence, earthquake nucleation, aseismic deformation, swarms∗†

Speaker
Corresponding author: david.marsan@univ-smb.fr

sciencesconf.org:statsei12:417169

The research on the spherical space-time ETAS
model
Ziyao Xiong∗1 and Jiancang Zhuang1
1

The Institute of Statistical Mathematics – Japan

Abstract
The widely used space-time ETAS (epidemic-type aftershock sequence) model was developed by Ogata (1998). This model successfully explains foreshocks and high order aftershocks
in earthquake sequences and provides a very effective tool for seismic activity analysis. However, this space-time ETAS model is only suitable for studying within a small space range.
When the study range is wide or in a high latitude, since the earth is a sphere, the simulation results of the model will produce errors. This study reformulates the space-time
ETAS model from its planar version to a spherical version by modifying the widely usedplanar space-time ETAS model (Zhuang et al., 2002; Zhuang, 2012). The new version isapplicable to analyze and forecast seismicity in high latitude regions or on a global scale.
And the model is verified by applying it to the global CMT catalog. The results show thatthe new model can well simulate the global seismicity variation. It supports better modelingof seismic activities and seismic interactions in global regions. Moreover, our model has highcomputational efficiency, even for a large range of global data.

∗

Speakersciencesconf.org:statsei12:416840

The role of three-dimensional fault interactions increating complex seismic sequences
Yifan Yin∗1 , Elı́as Heimisson1 , Percy Galvez2 , and Stefan Wiemer1
1

2

Swiss Seismological Service, ETH Zurich (SED, ETH Zurich) – Switzerland
King Abdullah University of Science and Technology (KAUST) – King Abdullah University of Scienceand Technology P.O.Box 55455 Jeddah 21534, Saudi Arabia, Saudi Arabia

Abstract

A physics-based earthquake simulator should reproduce first-order empirical power-law behaviors ofmagnitudes and inter-event time clustering. However, sequences exhibiting these laws have only been produced in discrete and low-dimension continuum simulations. We show that the same emergence also occurs
in 3-D continuum simulations with fault interactions. Our model approximates a strike-slip fault systemslipping under rate-and-state friction. We produce spatio-temporally clustered earthquake sequences exhibiting characteristic Gutenberg-Richter scaling as well as empirical inter-event time distribution. With
fault interaction, partial ruptures emerge when seismogenic width W over characteristic nucleation length
L ∞islargerthan16.24, butnoneoccurswithoutf aultinteraction.T hemainshockrecurrencetimesof individualf aultsrema
periodicandf itaBrownianpassagetimedistribution.T hesystemmainshockrecurrencetimehasashort−termOmori−typedecay, indicatinga22%chanceof mainshockclustering.W eextendthesimulationtoa10−f aultsystemmimickingthecen
basedmulti−cyclemodelsadequatelyref lectobservedstatisticalsignaturesandhaspotentialf orlong−termhazardassessm
termf orecasting.

Keywords: earthquake cycle, fault interaction, earthquake rupture∗

Speakersciencesconf.org:statsei12:416950

Keywords: modeling seismicity, branching processes, Galton, Watson processEvaluating earthquake forecasts with likelihood
based marginal and conditional scores
Jiancang Zhuang∗1
1

The Institute of Statistical Mathematics – Japan

Abstract
Earthquake forecasting consists of the forecasts of the numbers, occurrence times, locations, and magnitudes of seismic events, and even the correlations among them. In order toevaluate each component in the forecasting, CSEP provides varieties of testing criterions but
still far from enough. In this study, starting from the essential meanings of the full likelihoodof a point process, I will explain how to evaluate the marginal likelihoods for each component
in the forecasts and the conditional likelihoods for some components under the condition thatsome other components are given. These marginal and conditional likelihoods are the bases
for the marginal and conditional scores in performance evaluation when multiple models areused in the forecasting. I will also illustrate how to use Monte Carlo simulations to accelerate
the computation.

Keywords: scoring forecast, point process, CSEP, Monte Carlo, simulation∗

Speakersciencesconf.org:statsei12:416860

Neural-network based models for earthquake
prediction
Oleg Zlydenko∗†1 , Yohai Bar-Sinai∗2 , Neri Berman2 , Gal Elidan1 , Avinatan Hassidim1 ,Doron Kukliansky1 , Yossi Matias1 , Brendan Meade1,3 , and Sella Nevo‡11

Google Research – Israel
Tel Aviv University (TAU) – Israel
3
Harvard University – United States
2

Abstract
State-of-the-art models in statistical seismology, like ETAS and its variants, combine wellestablished statistical properties of earthquakes (e.g. Gutenberg-Richter law, and Omori’slaw) into a simple functional parametric form, whose parameters are learned from the data by
likelihood optimization. A natural way to improve on such models is using neural networks,
due to their ability to model complicated functional dependence, the relative ease withwhich different measurement modalities can be included, and the ability to simultaneously
model multiple aspects of the prediction problem. In this work we present such a neuralarchitecture, combining recent developments in neural point-process modeling and standard
ML techniques. Our models are trained to predict various statistical properties of future
seismicity, outperforming or competitive with ETAS models in various common metrics.
The architecture makes it easy to incorporate other information sources and modalities, andpose different prediction tasks, and we will discuss preliminary results in this direction.

Keywords: Rate prediction, neural networks, machine learning, earthquake forecasting∗

Speaker
Corresponding author: olegzl@google.com‡
Corresponding author: sellanevo@google.com†

sciencesconf.org:statsei12:415824

New findings on maximum magnitudes and future
seismicity scenarios for the gas field in Groningen,
The Netherlands
Gert Zöller∗1 and Sebastian Hainzl2
1

Institute for Mathematics / University of Potsdam – Karl-Liebknecht-Str. 24-25, 14476
Potsdam-Golm, Germany, Germany
2
GFZ German Research Centre for Geosciences, Potsdam – Germany

Abstract
The Groningen gas field is one of the most frequently considered sites for studying induced seismicity. Zöller and Holschneider (Bull. Seismol. Soc. Am. 106, 2917-2921, 2016)
calculated confidence intervals for the maximum possible earthquake magnitude and the
maximum expected magnitude in a future time interval. In recent years, gas production
has been reduced significantly and is planned to be terminated in about ten more years.
During this time interval, earthquakes are still expected to occur. Even later, a peteringout of seismicity is likely due to memory effects. We show that the confidence interval of
the estimated maximum possible magnitude decreases compared to 2016. This observation
may be explained by decreasing uncertainties arising from the growing size of the earthquake
catalog. Furthermore, we present seismicity scenarios for the next decades based on various
data-driven models. In particular, we employ physical models footing on Coulomb-FailureStress (CFS) and Rate-and-State-dependent friction (RS), as well as statistical models basedon inhomogeneous Poisson processes. The results indicate that earthquakes with a magnitude higher than 4.0 can be expected in the next 30 years, even if the gas field is shut downin about ten years.

Keywords: induced seismicity, maximum earthquake magnitude, statistics∗

Speakersciencesconf.org:statsei12:417111

Stress-Strain Characterization of Complex
Earthquakes and Seismic Sequences Using
Probabilistic Hierarchical Models
Alan Juarez Zuñiga∗†1 and Thomas H. Jordan1
1

University of South California (USC) – United States

Abstract
Earthquake ruptures and seismic sequences are the inelastic strain response to tectonic
stress fields. We develop a physics-based probabilistic model to characterize stress-strain
fields using moment measures of mechanism complexity. Distributions of focal mechanisms
in seismic sequences and earthquake ruptures are used to estimate the ratio of the Aki seismic moment to the total moment and the partitioning of the moment tensor field over anorthonormal basis of five deviatoric mechanisms. These moment measures of complexity
are then used to estimate the principal-stress directions, a differential stress ratio (R), thecomplexity factor (F’+), and a stress-strain sensitivity parameter (k). Large values of F’+
indicate high mechanism complexity, and high k indicates high-stress concentration. We
quantify the effect of focal mechanism uncertainty on the moment measures using a hierarchical probabilistic model. We apply the model to characterize the complexity of stress-strain
fields of seismicity in California. We show that the observed complexity is related to evolutionary geologic variables such as the cumulative slip in the fault. For example, the seismicity
along the San Andreas Fault in the Parkfield and Loma Prieta sections, which have over 150
km of cumulative slip, shows significant stress-strain sensitivity (k=12), i.e., low complexity.
In contrast, faults in the Mojave Block with less than 5 km of cumulative slip show a low
sensitivity (k=4.5). The stress field before the Mw7.2 Ridgecrest earthquake of 2019 showshigh complexity (k=7), while the aftershocks show less complexity with k=7.25, indicatingthat the aftershocks are more aligned with the mainshocks than the background seismicity.
These preliminary results suggest that seismic monitoring of mechanism complexity may
provide information that is useful in time-dependent earthquake forecasting.

Keywords: Mechanism Complexity, Probabilistic Forecasting, Stress Inversion∗†

Speaker
Corresponding author: alanjuar@usc.edu

sciencesconf.org:statsei12:416284