# Characteristics of secondary slip fronts associated with slow earthquakes in Cascadia

Quentin Bletery, Amanda Thomas, Jessica Hawthorne, Robert Skarbek, Alan Rempel, Randy Krogstad

► To cite this version:

Quentin Bletery, Amanda Thomas, Jessica Hawthorne, Robert Skarbek, Alan Rempel, et al.. Characteristics of secondary slip fronts associated with slow earthquakes in Cascadia. Earth and Planetary Science Letters, 2017, 463, pp.212-220. 10.1016/j.epsl.2017.01.046. hal-02389208

HAL Id: hal-02389208

https://hal.science/hal-02389208v1

Submitted on 2 Jul 2025

HAL is a multi-disciplinary open access
archive for the deposit and dissemination of sci-
entific research documents, whether they are pub-
lished or not. The documents may come from
teaching and research institutions in France or
abroad, or from public or private research centers.

L'archive ouverte pluridisciplinaire HAL, est
destinée au dépôt et à la diffusion de documents
scientifiques de niveau recherche, publiés ou non,
émanant des établissements d'enseignement et de
recherche français ou étrangers, des laboratoires
publics ou privés.

<figure><img src="image_1.png" /></figure>

Distributed under a Creative Commons Attribution 4.0 International License
Accepted Manuscript

Characteristics of secondary slip fronts associated with slow earthquakes in Cascadia

Quentin Bletery, Amanda M. Thomas, Jessica C. Hawthorne, Robert M. Skarbek, Alan
W. Rempel, Randy D. Krogstad

PII: S0012-821X(17)30058-4

DOI: http://dx.doi.org/10.1016/j.epsl.2017.01.046

Reference: EPSL 14310

To appear in: *Earth and Planetary Science Letters*

Received date: 10 November 2016

Revised date: 30 January 2017

Accepted date: 31 January 2017

<figure><img src="image_2.png" /></figure>

Please cite this article in press as: Bletery, Q., et al. Characteristics of secondary slip fronts associated with slow earthquakes in Cascadia. *Earth Planet. Sci. Lett.* (2017), http://dx.doi.org/10.1016/j.epsl.2017.01.046

This is a PDF file of an unedited manuscript that has been accepted for publication. As a service to our customers we are providing this early version of the manuscript. The manuscript will undergo copyediting, typesetting, and review of the resulting proof before it is published in its final form. Please note that during the production process errors may be discovered which could affect the content, and all legal disclaimers that apply to the journal pertain.
## Highlights

* We implemented an algorithm to automatically detect secondary fronts (SSFs).

* We provide a catalog of >1,000 SSFs in Cascadia including some of their physical characteristics.

* Estimated SSF slip rates and stress drops vary over 3 orders of magnitude.

* SSF duration seems related to SSF direction.

* SSFs fill a major gap in the Moment - Duration scaling law of slow slip processes.
Characteristics of secondary slip fronts associated with
slow earthquakes in Cascadia

Quentin Bletery<sup>a</sup>, Amanda M. Thomas<sup>a</sup>, Jessica C. Hawthorne<sup>b</sup>, Robert M. Skarbek<sup>c</sup>, Alan W. Rempel<sup>a</sup>, Randy D. Krogstad<sup>a</sup>

<sup>a</sup>Department of Earth Sciences, University of Oregon, 1272 University of Oregon, Eugene, Oregon 97403, USA.

<sup>b</sup>University of Leeds, School of Earth and Environment, Leeds, United Kingdom.

<sup>c</sup>Lamont-Doherty Earth Observatory, Columbia University, Palisades, New York, USA.

Abstract

We implement an algorithm to automatically detect migrations of low fre-
quency earthquakes at time scales between 30 min and 32 h during the 2003,
2004 and 2005 slow slip events in Cascadia. We interpret these migrations
of seismicity as a passive manifestation of secondary slip fronts (SSFs) that
propagate faster than the main front. We identify the dominant features of
383 SSFs, including time, location, duration, area, propagation velocity and
estimate: their moment, stress drop, slip, and slip rate. We apply the same
algorithm to continuous tremor detection in Cascadia between 2009 and 2015
and characterize 693 SSFs at time scales between 4 h and 32 h. We identify
— to our knowledge for the first time — numerous 10-24 h long SSFs that
propagate at velocities intermediate between slow slip events and previously
reported SSFs. The systematic detection of SSFs fills a gap between seismi-
cally and geodetically detectable slow earthquake processes. Analyses of SF
basic features indicates a wide range of stress drops and slip rates (with me-

Email address: qbletery@uoregon.edu (Quentin Bletery)

Preprint submitted to Earth and Planetary Science Letters
dians of 5.8 kPa and 1.1 mm/h) as well as an intriguing relationship between SSF direction and duration that was observed in other contexts and could potentially help discriminate between the different physical models proposed to explain slow slip phenomena.

*Keywords:* Slow slip, low-frequency earthquakes, nonvolcanic tremor, secondary slip fronts, Cascadia

# 1. Introduction

Slow slip events (SSEs) in subduction zones propagate along the plate interface at velocities of 5-10 km/day and are largely confined to the region known as the transition zone, located down-dip of the seismogenically locked zone (Dragert et al., 2001; Bartlow et al., 2011). As SSEs propagate, small on-fault asperities capable of generating seismic radiation fail in earthquake-like events known as low-frequency earthquakes (LFEs) and nonvolcanic tremor (Obara, 2002; Kao et al., 2006; Ide, 2010; Ueno et al., 2010; Houston et al., 2011). Multiple lines of evidence indicate that LFEs and tremor — which has been interpreted as a superposition of LFEs (Shelly et al., 2007b) — are a passive manifestation of surrounding aseismic slow slip (Shelly et al., 2007b). The occurrence of LFEs and tremor during SSEs are not random. Instead, spatially and temporally coherent migrations of LFEs and tremor can often be identified within the actively slipping portion of the fault after the main slip front associated with the SSE has passed. These migration patterns are thought to represent secondary slip fronts (SSFs) occurring over several different length and timescales. For example, Houston et al. (2011) and Royer et al. (2015) used tremor and LFEs to image approximately 40 km long slip
19 fronts — known as rapid tremor reversals (RTRs) — that propagate back
20 in the opposite direction but 10 to 40 times faster than the main slip front.

21 Ghosh et al. (2010) used tremor in Cascadia to image seismic migrations
22 — known as streaks — that propagate back and forth along the main slip
23 front (i.e., perpendicular to its propagation direction) but approximately 100
24 times faster. Rubin and Armbruster (2013) reported tremor migrations at the
25 kilometer scale that propagate along strike either towards or away from the
26 main front, but approximately ten times faster. Finally, Peng et al. (2015a)
27 and Peng and Rubin (2016) found evidence for complex tremor migration
28 patterns beneath that evolve as the main front propagates. While all of
29 the aforementioned studies were conducted in Cascadia, similar observations
30 have been reported in Japan (Shelly et al., 2007a; Ueno et al., 2010; Obara
31 et al., 2012), Mexico (Radiguet et al., 2012; Frank et al., 2014; Peng and
32 Rubin, 2017) and along the San Andreas fault (Ryberg et al., 2010; Shelly,
33 2009, 2010; Peng et al., 2015b), suggesting that SSFs may be a ubiquitous
34 feature of slow slip.

35 Slow-slip phenomena require the fault slip rate to increase to observed
36 speeds, typically 10 to 100 times the plate rate, but abstain from accelerating
37 to typical earthquake slip rates. The specific physical mechanisms responsible
38 for imposing this speed limit are still debated. SSFs generally propagate at
39 least an order of magnitude faster than the main front associated with the
40 SSE, and estimates of the ratio of the SSF stress drop to that of the main front
41 range between one tenth and one (Rubin and Armbruster, 2013; Royer et al.,
42 2015; Hawthorne et al., 2016). These stress drop estimates combined with the
43 propagation speeds mentioned above imply that SSFs must be slipping either
as fast as or faster than the main front. Hence, reproducing stress drops, slip and propagation speeds of both the main and secondary fronts requires a rheology that limits the speed of the main front while allowing the SSFs to slip one to two orders of magnitude faster. Such behavior is difficult to reproduce with a single state variable friction law (Rubin, 2011). Additionally, SSFs appear to be sensitive to very small (≈ 1 kPa) tidal stress perturbations (Thomas et al., 2013; Royer et al., 2015; Houston, 2015) allowing them to be used to probe in situ conditions in deep faults. Observations of SSFs, in addition to SSEs, may provide useful constraints on physical models of slow slip processes (Rubin, 2011; Hawthorne and Rubin, 2013).

To date, most studies that have documented SSFs have relied on subjective identification methods such as visual identification, often using only one dimension of space and time (Shelly et al., 2007a; Ueno et al., 2010; Ghosh et al., 2010; Houston et al., 2011; Rubin and Armbruster, 2013; Royer et al., 2015; Peng et al., 2015a; Savard and Bostock, 2015). While these approaches successfully find some SSFs, such subjective identification makes it possible and even likely that many SSFs go undetected. Here we develop an algorithm that systematically identifies SSFs based on spatial and temporal migration of tremor/LFEs and discuss its advantages compared to the similar method developed by Obara et al. (2012). We apply our algorithm to both tremor and LFE datasets resulting from SSEs that occurred between 2003 and 2015 in Cascadia. We detect over 1000 SSFs, determine their locations, duration, directions, and estimate some of their physical properties, such as moment, average slip, stress drop and slip rate.
## 2. Data and algorithm

We use two types of data: a catalog of LFEs located near Vancouver Island during the 2003, 2004 and 2005 SSEs (Savard and Bostock, 2015) and a catalog of continuous tremor detections from August 2009 to December 2015 that occur throughout the Cascadia subduction zone (Wech and Creager, 2008; Wech, 2010). The LFE catalog was obtained using cross-station waveform correlation; thanks to a dense campaign seismic network, it has very good spatial (location errors of < 1 km) and temporal resolution (Savard and Bostock, 2015). The tremor catalog contains locations and detection times of tremor events determined using waveform envelope correlation and clustering (Wech and Creager, 2008; Wech, 2010). Detections correspond to 5-minute seismogram windows that likely represent a superposition of LFEs (Shelly et al., 2007b). Therefore, resolution of both time (each detection is separated by at least 2.5 min) and location (5-10 km) is worse than for the LFE catalog. Because the tremor detection algorithm is independently implemented on several adjacent seismic networks, some tremor detections appear twice in the final catalog (Wech, 2010). To eliminate duplicate detections, when two tremors are detected at the same time within 25 km, we exclude one (Wech, 2010).

These two catalogs provide a large number of identified LFE and tremor events. To identify possible SSFs, we analyze LFE and tremor migration patterns in 3 dimensions: distance along strike, distance along dip, and time. The algorithm is illustrated in Figure 1. As in the approach of Obara et al. (2012), we look for SSFs in time windows of various durations. The different time windows considered in our analysis are 30 min, 1 h, 2 h, 4 h, 8 h, 16 h,
32 h for the LFE catalog and 4 h, 8 h, 16 h, 32 h for the tremor catalog. As individual detections are all separated by at least 2.5 min (Wech and Crea-ger, 2008), the number of potential detections in a ≤ 2 h time window too small to analyze their potential migration. The analysis is performed inde-pendently on the different time windows. Therefore, tremor/LFE detections may theoretically be part of more than one SSF in different time windows; thus a short duration SSF may be a part of a longer duration SSF.

We first identify clusters of LFE and tremor detections in the individ-
ual catalogs. We use the subclust matlab routine to identify LFE/tremor
clusters of different duration. The subclust routine estimates the centers
of subclusters of a given duration in a data set using subtractive clustering
(Chiu, 1994; Yager and Filev, 1994). We set the duration of the subclusters
to correspond to the different time windows. The subclust routine then re-
turns the tremor/LFE clusters of 30 min to 32 h duration with the highest
concentration of events. For each identified cluster, we then calculate the lin-
ear regression of the distance along strike as a function of time and exclude
points farther than two standard deviations (2σ) from the regression line.
After removing those points, we do the same linear regression but for the
distance along dip as a function of time and remove events farther than 2σ
from the linear regression. Finally, we project the remaining points along the
best-fitting migration axis and remove again events that lie outside 2σ from
the linear regression of the distance along the migration axis as a function
of time. The remaining points are then tested: if the cluster still contains
more than 20 LFE or tremor detections and if the root mean square (RMS)
of their distance to the last regression is smaller than 15% of the migration
length, then the cluster is identified as an SSF. For a given time window, if
the SSF length is small, then the front is not moving quickly (the signal is
low). Thus we require the deviation from the propagation axis (the noise)
to be very low (to ensure a decent signal to noise ratio). If the SSF length
is larger, then the front is clearly moving (the signal is larger) and we can
be more tolerant with the deviation from the propagation axis (as we have
a clear signal). We repeat the process for all the clusters. We then remove
the LFE or tremor detections identified as part of an SSF from the catalog
and rerun the clustering routine on the new catalog with identified SSFs re-
moved. We iterate until the catalog no longer contains any LFE or tremor.
If no SSF is identified on one iteration, we remove all the tested clusters from
the catalog and iterate once more in the truncated catalog. All the detected
SSFs are then tested to ensure the continuity of the migration. We subdivide
the successfully identified SSFs into 4 equal portions along their best-fitting
migration axis and 3 equal periods and if one of the parts or periods contains
less than 1/20th of the number of LFE/tremor detections that compose the
SSF, the SSF is rejected in order to ensure that each SSF is continuous rather
than composed of two distinct SSFs.

Our approach is inspired by Obara et al. (2012) (who implemented an
automatic SSF detection based on principal component analysis of tremor in
Japan), but has some significant differences. Obara et al. (2012) identified
SSFs using moving time windows shifted by discrete increments. As it is
unlikely that SSFs all start and end at the specific times of arbitrary time
windows, this may result in detecting SSFs that either are just a part of a
larger undetected SSFs or include LFE/tremor that should not be associated
with the SSF, potentially affecting the apparent migration velocities. The
use of a clustering subroutine allows for flexibility in the start and end time
of the explored SSFs — as tested time windows are not shifted by discrete
increments but determined in respect of high concentrations of tremor/LFEs
— substantially increasing the number of SSF detections. Another significant
difference is the detection thresholds. In the work of Obara et al. (2012),
SSF acceptance is based on linearity, angular difference and in the number
of tremor detections, all of which vary with the given time window. We
apply only two thresholds, which are the same for all time windows: the
number of LFE/tremor detections (20 for all time windows) and a deviation
criteria based on the ratio of the misfit to the SSF length (required to be
$\le 15\%$), which can be interpreted as a signal to noise ratio. Our choices of
20 LFE/tremor detections and 15% are arbitrary — as are the thresholds
used by Obara et al. (2012) — but these values do provide visually coherent
patterns in all the detected SSFs (Figures 2, 3).

To ensure that our algorithm does not interpret noise as SSF, we test
it on synthetic catalogs. To produce the synthetic catalogs, we use the ac-
tual catalog locations and randomize the time of occurrence. We run the
algorithm for all the time windows of our analysis and detect no SSF in the
synthetic catalogs.

3. Spatial distribution of secondary slip fronts

We identified 383 SSFs in the LFE catalog of the 2003, 2004, and 2005
SSEs beneath Vancouver Island and 693 SSFs in the tremor catalog between
August 2009 and December 2015 in Cascadia. While SSFs lasting from 10
min to 11 h have been documented in multiple subduction zones (Ghosh et al., 2010; Houston et al., 2011; Obara et al., 2012; Rubin and Armbruster, 2013; Royer et al., 2015; Peng et al., 2015a; Savard and Bostock, 2015), we identified — to our knowledge for the first time — numerous SSFs lasting between 11 and 22.5 h in both catalogs (243 in the tremor catalog, 57 in LFEs catalog). They propagate at intermediate speeds (≈ 3 km/h) between those of the main front (≈ 0.2–0.5 km/h) and previously documented SSFs (≈ 40 km/h) (Ghosh et al., 2010; Houston et al., 2011; Obara et al., 2012; Rubin and Armbruster, 2013; Royer et al., 2015; Peng et al., 2015a; Savard and Bostock, 2015).

The spatial distribution of SSFs is of particular interest because it can reveal atypical behaviors of some fault areas and may highlight spatial variations of fault properties. Figure 4 shows the distribution of tremor derived SSFs during the 2009-2015 period along the Cascadia megathrust by time window. Note that some long duration tremor derived SSFs may correspond to the main front of short duration SSEs. In the inset diagrams we plot the number of SSFs as a function of their propagation directions. Some areas show a clear tendency to host numerous SSFs (for instance at latitude 40°N-41.5°N, 42°N-42.5°N, or 44.5°N-45.5°N) whereas other regions had no detectable SSFs. Additionally, there is a clear preference for along-strike propagation for the majority of the tremor derived SSFs. Both the spatial distribution and propagation directions are consistent from SSE to SSE, suggesting that the physical mechanism responsible for generating SSFs does not change over a timescale of several years.

Because both the temporal and spatial resolution are much better in the
LFE catalog, we are able to detect short timescale SSFs beneath Vancouver island (Figure 5). Since the 3 SSEs that drove LFEs in this area all roughly propagated in the same direction (blue arrow in Figure 5), this catalog is also more convenient to analyze the direction of SSFs with respect to the direction of SSEs. SSFs propagate in all directions and encompass phenomena previously described as RTRs or streaks. They also encompass slip fronts that propagate in the direction of the SSE (Figure 5). Such forward-propagating SSFs are hard to distinguish by eye from LFE or tremor migration associated with the main front. They typically initiate after the main front has passed and migrate towards the tip of the main front. Interestingly, we observe a clear relationship between the direction and the duration of SSFs. The vast majority of short duration (i.e. 30 min and 1 h) SSFs propagate along dip. As the SSF duration increases, the proportion of SSFs propagating along strike gradually increases (see Rose diagrams in Figure 5) resulting in the longest duration LFE derived SSFs propagating almost solely along strike. Similar relationships were observed by Ghosh et al. (2010) and Houston et al. (2011) in Cascadia and by Obara et al. (2012) beneath the Kii peninsula in the Nankai subduction zone. Potential physical mechanisms responsible for SSFs will have to explain this enigmatic dependence between duration and direction (Rubin, 2011). The ratio of along-strike SSFs propagating in the SSE direction to along-strike SSFs propagating in the opposite direction also seems to depend on the duration. While reversal SSFs are predominant in the 8 h time window, almost all SSFs propagate forward in the 32 h time window.

The distribution of propagation velocity $V_{prop}$ is shown at the bottom of
each subfigure in Figures 4 and 5. In individual time windows, along-dip SSFs do not appear faster than along-strike SSFs. However, as the length $L$ of SSFs is limited by the extent of the fault, the average propagation velocity decreases with the duration $T$ (because $V_{prop} = L/T$). Accordingly, the limited number of distinct time windows in the analysis of the tremor catalog results in little difference in the propagation velocity (Figure S1.f).

## 4. Physical characteristics of secondary slip fronts

Our detection algorithm provides a number of direct observations characterizing each SSF including start time, location, direction (Figures 4, 5), duration $T$, length $L$, width $W$, pulse (i.e. actively slipping) length $l$ (Figure 6), propagation velocity $V_{prop}$ (defined as $L/T$), and number of LFE or tremor detections $n_e$. The distributions of these SSF properties are shown in Figure S1 (catalogs are available in the supplementary materials). Assessing other physical characteristics of SSFs, such as moment or stress drop, is challenging because these quantities require an estimate of the total slip associated with each SSF.

Individual SSFs are too small to generate a surface deformation signal detectable with GPS and hence most studies that have estimated slip for individual SSFs have relied on the total SSE moment and the fraction of LFE detections or bursts that occured as part of the SSF (Rubin and Armbruster, 2013; Royer et al., 2015). Here, we adopt a similar approach and assume that the total moment $M_e$ associated with each LFE or tremor detection is the same (although different for LFEs and tremor). This assumption is demonstrably false for individual LFE repeats, which have been shown to
have very different moments (Bostock et al., 2015), and for individual tremor events, as their definition (i.e. the 2.5 minute timescale used to separate individual events) is arbitrary and includes no contribution from waveform amplitudes (Wech and Creager, 2008). However, since moment estimates for each SSF include a minimum of 20 LFEs or tremors, it is likely that the average moment will be a more accurate measure of total SSF moment. Nevertheless, our estimates of the physical characteristic of SSFs (such as moment, slip, stress drop, slip rate) are tentative and caution should be exercised when interpreting them.

We obtain $M_e$, the average moment associated with each LFE or tremor, by dividing the moment $M_0$ (assessed by geodetic studies (Schmidt and Gao, 2010; Krogstad and Schmidt, 2015)) associated with an SSE by the number of events $N_e$ detected during the event:

$$M_e = \frac{M_0}{N_e}. \qquad (1)$$

An implicit assumption behind equation (1) is that moment associated with aseismic slip that does not generate LFE or tremor detections is negligible. The moment, $M_{SSF}$, associated with one SSF, is defined as the sum of the moment associated with all the LFE or tremor detections that are part of that SSF,

$$M_{SSF} = \sum_{i=1}^{n_e} M_e(i), \qquad (2)$$

where $n_e$ is the number of LFE/tremor detections part of the given SSF. As we assumed that for all $i$, $M_e(i) = \frac{M_0}{N_e}$,

$$M_{SSF} = \frac{n_e}{N_e} M_0. \qquad (3)$$
Note that possible deviations from constant $M_e(i)$ would tend to be averaged out over the number of events in equation (3) because $M_e(i)$ is sometimes larger and sometimes smaller than the assumed average $M_e = \frac{M_0}{N_e}$. We further estimate an average slip $d$ for each SSF,

$$d = \frac{M_{SSF}}{\mu A}, \qquad (4)$$

where $\mu$ is the shear modulus and $A = WL$ is the SSF slip area, with the SSF width $W$ and length $L$ both estimated from our analysis (Figure 6). To be consistent with geodetic inversions, we use $\mu = 40$ GPa (Schmidt and Gao, 2010; Krogstad and Schmidt, 2015). The stress drop $\Delta\sigma$ is then given by

$$\Delta\sigma = \frac{4(\lambda + \mu)}{\pi(\lambda + 2\mu)}\mu\frac{d}{W}, \qquad (5)$$

where $\lambda$ is Lamé's parameter (Lay and Wallace, 1995), which we set to $\lambda = 40$ GPa. The total slip $d$ accumulates in the time the front propagates the distance $l$. Thus, we can deduce the slip rate $V_{slip}$ from the propagation velocity $V_{prop}$ as

$$V_{slip} = \frac{d}{l}V_{prop}. \qquad (6)$$

The distributions of assessed physical characteristics (moment, slip, stress drop and slip rate) are shown in Figures S2 and S3 (catalogs in the supplementary materials). Because we did not estimate these parameters independently, their distributions are correlated with each other (Figure S4). However, most of the correlation coefficients between the different parameters are relatively weak (i.e. below 0.6 in magnitude, see Figure S4). As expected, high correlation coefficients are found for slip and stress drop, since the stress drop is estimated from slip (eq. 5), and between the different spatial dimensions ($L$, $W$ and $l$), implying that the aspect ratios of SSFs do not vary
appreciably. Significant anticorrelations are observed between propagation
rate $V_{prop}$ and duration $T$, as well as between the spatial dimensions and both
slip $d$ and slip rate $V_{slip}$. It is also noteworthy that stress drop $\Delta\sigma$ appears
to be independent of moment $M$ in both catalogs (Figure S4).

The median stress drop among all SSFs is 5.8 kPa, which is very con-
sistent with stress drop estimates obtained for RTRs using strainmeters (8
kPa) (Hawthorne et al., 2016). These values of stress drop are several times
smaller than typical stress drops in SSEs (10 - 100 kPa), and over 2 or-
ders of magnitude smaller than regular earthquakes (Gao et al., 2012). The
estimated stress drops are also all consistently smaller than effective nor-
mal stress (0.2-3 MPa) inferred in LFE source regions (Liu and Rice, 2007;
Hawthorne and Rubin, 2013; Royer et al., 2015). The median slip is 3 mm
and the median slip rate is 1.1 mm/h, which is much faster than slip rates in-
ferred for SSEs (on the order of 1-2 mm/day (Wech and Bartlow, 2014)). The
observation that SSF slip rates are significantly faster than SSE slip rates is
consistent with previous studies (Rubin and Armbruster, 2013; Royer et al.,
2015; Hawthorne et al., 2016) and implies that the fault rheology limits the
main front speed and then allows secondary fronts to slip faster. As ac-
knowledged in the introduction, such a time dependent behavior is not easy
to reproduce with a single state variable friction law (Rubin, 2011). SSF slip
rates may then be used to constrain physical models of slow slip phenomena.

Though variations in the moments and slip distributions of SSFs range
within a single order of magnitude (4.10<sup>16</sup> to 4.10<sup>17</sup> N.m and 0.1 to 40 mm)
for all SSFs, stress drop varies by over 3 orders of magnitude (from 0.1 to
300 kPa), as does slip rate (between 0.02 and 30 mm/h) (Figures S2 and
S3). We recall here that our estimates of stress drop are based on the original assumption that SSF moment scales with the number of LFE/tremor detections, which could be inaccurate for some SSFs. Large scatter in slip, stress drop and slip rate is also expected from equations (4), (5) and (6) as a result of multiplying quantities ($M_{SSF}$, $L$, $W$, $l$ and $V_{prop}$) that have large variability. Therefore a large part of the scatter in slip rate and stress drop might be attributed to assumption inaccuracy and propagation of uncertainties, but some of this scatter could also be sourced to changes in fault properties. As SSF areas widely overlap — this is particularly clear in the LFE catalog (Figure 5) — the varying SSF properties might suggest that, in addition to being able to slip at different rates, the same fault area is also able to release very different amounts of stress.

## 5. Moment – duration scaling

Ide et al. (2007) suggested that slow slip phenomena — including LFEs, very low frequency earthquakes (VLFEs) and SSEs — follow a moment – duration scaling law distinct from regular earthquakes. While regular (fast) earthquakes have moment proportional to the cube of their duration, slow slip phenomena would have moment proportional to duration (with $M/T = 10^{12}$–$10^{13}$ N.m/s) (Figure 7). In this representation, important gaps separate the various slow slip phenomena (Ide et al., 2007). The question of whether these gaps result from detection limitations or from the inherent behavior of fault slip is of importance because it could provide us with physical constraints on earthquake processes (Ide et al., 2007; Ide, 2008; Ide et al., 2008; Peng and Gomberg, 2010; Gao et al., 2012; Ide and Yabe, 2014; Liu, 2014; Gomberg
et al., 2016).

In the far field, seismic amplitudes scale with moment rate (Aki and Richards, 2002). Thus, the threshold for seismically detectable signals follows the noise level in the different frequency bands. We use the USGS low noise model (Peterson et al., 1993; Ide, 2014) and empirically calibrate it on the smallest moment rate events that have been seismically detected (VLFEs from Ide and Yabe (2014)) to determine a threshold for seismically detectable events (red dashed line in Figure 7). Existing geodetic instrumentation has proven capable of detecting magnitude 5.7–6.0 events and above regardless of their duration (Peng and Gomberg, 2010) (blue dashed line in Figure 7). The gap between LFEs and VLFEs (Figure 7) may be attributed to particularly large noise in the frequency band corresponding to the gap (Ide et al., 2008) while the gap between the VLFEs and SSEs corresponds to moments too low to be detected by geodesy (Peng and Gomberg, 2010) and characteristic periods too long to be detected in seismograms (Ide et al., 2008).

The large range of duration and moment among the SSFs we detected follows the slow slip scaling law and fills the gap between VLFEs and SSEs, where signal is not directly detectable either from seismic or geodetic records (Figure 7). The vertical (moment) threshold in SSFs corresponds to the minimum number of LFE/tremor detections (20) in one SSF times $M_e$ and the distinct horizontal lines in SSF duration correspond to our different time windows. These two noncontinuous characteristics are consequences of our approach. Interestingly, SSFs appear to follow the same threshold of constant moment rate — on the order of $10^{13}$ N.m/s — as the other slow slip phenomena (Figure 7). The fact that SSFs follow the slow slip phenomena
scaling law may be a consequence of estimating SSF moment from moment of SSEs (eq. 3), as SSEs follow this law. Nevertheless, it is worth noting that our moment estimates are in agreement with SSF moments estimated from strain meters (Hawthorne et al., 2016) (green square in Figure 7).

While individual slow slip phenomena do not appear to have moment proportional to duration (e.g. Ide et al. (2008); Bostock et al. (2015); Gomberg et al. (2016)), they all fit a constant moment rate scaling law (within less than 2 orders of magnitude variation) when considered together on a 10 orders of magnitude scale (Figure 7). More generally, fast and slow slip phenomena appear to be well distinguished in Figure 7. Fault slip processes that would range between the two scaling laws are, in theory, detectable — and should therefore have been detected — by seismology, geodesy or both (Figure 7). The large gap in detection of fault slip processes between the two scaling laws suggests that earthquakes and slow slip phenomena are two very distinct fault slip processes.

## 6. Conclusion

We identified a total of 1076 SSFs in Cascadia at timescales between 30 min and 32 h. We cataloged their duration, width, length, pulse length and time of occurrence and provided estimates for some of their physical properties such as moment, slip, slip rate, and stress drop. The median moment is $5.3 \times 10^{16}$ N.m (corresponding to a moment magnitude of 5.1), the median slip is 3 mm, the median slip rate is 1.1 mm/h and the median stress drop is 5.8 kPa. Moment estimates suggest that SSFs bridge a major gap between seismically and geodetically detectable events in the slow slip phenomena.
nomena scaling law (Ide et al., 2007). As previously observed with different approaches (Ghosh et al., 2010; Houston et al., 2011; Obara et al., 2012), we found that short duration SSFs dominantly propagate along dip while long duration SSFs dominantly propagate along strike. We also found that inferred SSF stresses drops and slip rates range over 3 orders of magnitude, suggesting that the same fault region could be able to slip at very different rates and to release very different stresses. These findings may be used to constrain models of slow earthquake processes.

## Acknowledgments

This work was funded by the NSF grant EAR-1520238. The tremor catalog was made available online by A. Wech and K. Creager (tunk.ess.washington.edu/map_display/). We thank G. Savard and M. Bostock for providing the LFE catalog.

## References

Aki, K., Richards, P. G., 2002. Quantitative seismology. Vol. 1.

Bartlow, N. M., Miyazaki, S., Bradley, A. M., Segall, P., 2011. Space-time correlation of slip and tremor during the 2009 Cascadia slow slip event. Geophysical Research Letters 38 (18).

Bletery, Q., Sladen, A., Delouis, B., Vallée, M., Nocquet, J.-M., Rolland, L., Jiang, J., 2014. A detailed source model for the mw9.0 Tohoku-oki earthquake reconciling geodesy, seismology, and tsunami records. Journal of Geophysical Research: Solid Earth 119 (10), 7636–7653.
Bletery, Q., Sladen, A., Jiang, J., Simons, M., 2016. A bayesian source model for the 2004 great sumatra-andaman earthquake. Journal of Geophysical Research: Solid Earth 121 (7), 5116–5135.

Bostock, M., Thomas, A., Savard, G., Chuang, L., Rubin, A., 2015. Magnitudes and moment-duration scaling of low-frequency earthquakes beneath southern vancouver island. Journal of Geophysical Research: Solid Earth 120 (9), 6329–6350.

Chiu, S. L., 1994. Fuzzy model identification based on cluster estimation. Journal of Intelligent & Fuzzy Systems 2 (3), 267–278.

Delouis, B., Nocquet, J.-M., Vallée, M., 2010. Slip distribution of the febru-
ary 27, 2010 mw= 8.8 maule earthquake, central chile, from static and
high-rate gps, insar, and broadband teleseismic data. Geophysical Research
Letters 37 (17).

Dragert, H., Wang, K., James, T. S., 2001. A silent slip event on the deeper
cascadia subduction interface. Science 292 (5521), 1525-1528.

Frank, W. B., Shapiro, N. M., Husker, A. L., Kostoglodov, V., Romanenko, A., Campillo, M., 2014. Using systematically characterized low-frequency earthquakes as a fault probe in guerrero, mexico. Journal of Geophysical Research: Solid Earth 119 (10), 7686–7700.

Gao, H., Schmidt, D. A., Weldon, R. J., 2012. Scaling relationships of source
parameters for slow slip events. Bulletin of the Seismological Society of
America 102 (1), 352–360.
Ghosh, A., Vidale, J. E., Sweet, J. R., Creager, K. C., Wech, A. G., Houston, H., Brodsky, E. E., 2010. Rapid, continuous streaking of tremor in cascadia. *Geochemistry, Geophysics, Geosystems* 11 (12).

Gomberg, J., Wech, A., Creager, K., Obara, K., Agnew, D., 2016. Reconsidering earthquake scaling. *Geophysical Research Letters* 43 (12), 6243–6251.

Hawthorne, J., Rubin, A., 2013. Tidal modulation and back-propagating fronts in slow slip events simulated with a velocity-weakening to velocity-strengthening friction law. *Journal of Geophysical Research: Solid Earth* 118 (3), 1216–1239.

Hawthorne, J. C., Bostock, M. G., Royer, A. A., Thomas, A. M., 2016. Variations in slow slip moment rate associated with rapid tremor reversals in cascadia. *Geochemistry, Geophysics, Geosystems*, in review.

Houston, H., 2015. Low friction and fault weakening revealed by rising sensitivity of tremor to tidal stress. *Nature Geoscience* 8 (5), 409–415.

Houston, H., Delbridge, B. G., Wech, A. G., Creager, K. C., 2011. Rapid tremor reversals in cascadia generated by a weakened plate interface. *Nature Geoscience* 4 (6), 404–409.

Ide, S., 2008. A brownian walk model for slow earthquakes. *Geophysical Research Letters* 35 (17).

Ide, S., 2010. Striations, duration, migration and tidal response in deep tremor. *Nature* 466 (7304), 356–359.
Ide, S., 2014. Modeling fast and slow earthquakes at various scales. Proceedings of the Japan Academy, Series B 90 (8), 259-277.

Ide, S., Beroza, G. C., Shelly, D. R., Uchide, T., 2007. A scaling law for slow earthquakes. Nature 447 (7140), 76-79.

Ide, S., Imanishi, K., Yoshida, Y., Beroza, G. C., Shelly, D. R., 2008. Bridging the gap between seismically and geodetically detected slow earthquakes. Geophysical Research Letters 35 (10).

Ide, S., Yabe, S., 2014. Universality of slow earthquakes in the very low frequency band. Geophysical Research Letters 41 (8), 2786-2793.

Ito, Y., Obara, K., 2006. Very low frequency earthquakes within accretionary prisms are very low stress-drop earthquakes. Geophysical Research Letters 33 (9).

Kao, H., Shan, S.-J., Dragert, H., Rogers, G., Cassidy, J. F., Wang, K., James, T. S., Ramachandran, K., 2006. Spatial-temporal patterns of seismic tremors in northern cascadia. Journal of Geophysical Research: Solid Earth 111 (B3).

Krogstad, R. D., Schmidt, D. A., 2015. Assessing the updip spatial offset of tremor and slip during ets events in cascadia. In: 2015 AGU Fall Meeting. Agu.

Lay, T., Wallace, T. C., 1995. Modern global seismology. Vol. 58. Academic press.
Liu, Y., 2014. Source scaling relations and along-strike segmentation of slow slip events in a 3-d subduction fault model. Journal of Geophysical Research: Solid Earth 119 (8), 6512-6533.

Liu, Y., Rice, J. R., 2007. Spontaneous and triggered aseismic deformation transients in a subduction fault model. Journal of Geophysical Research: Solid Earth 112 (B9).

Matsuzawa, T., Obara, K., Maeda, T., 2009. Source duration of deep very low frequency earthquakes in western shikoku, japan. Journal of Geophysical Research: Solid Earth 114 (B11).

Obara, K., 2002. Nonvolcanic deep tremor associated with subduction in southwest japan. Science 296 (5573), 1679-1681.

Obara, K., Matsuzawa, T., Tanaka, S., Maeda, T., 2012. Depth-dependent mode of tremor migration beneath kii peninsula, nankai subduction zone. Geophysical Research Letters 39 (10).

Peng, Y., Rubin, A. M., 2016. High-resolution images of tremor migrations beneath the olympic peninsula from stacked array of arrays seismic data. Geochemistry, Geophysics, Geosystems 17, 587-601.

Peng, Y., Rubin, A. M., 2017. Intermittent tremor migrations beneath guer-
rero, mexico and implications for fault healing within the slow slip zone.
Geophysical Research Letters.

Peng, Y., Rubin, A. M., Bostock, M. G., Armbruster, J. G., 2015a. High-resolution imaging of rapid tremor migrations beneath southern vancou-
ver island using cross-station cross correlations. Journal of Geophysical Research: Solid Earth 120 (6), 4317-4332.

Peng, Z., Gomberg, J., 2010. An integrated perspective of the continuum between earthquakes and slow-slip phenomena. Nature Geoscience 3 (9), 599-607.

Peng, Z., Shelly, D. R., Ellsworth, W. L., 2015b. Delayed dynamic triggering of deep tremor along the parkfield-cholame section of the san andreas fault following the 2014 m6. 0 south napa earthquake. Geophysical Research Letters 42 (19), 7916-7922.

Peterson, J., et al., 1993. Observations and modeling of seismic background noise. US Geological Survey Albuquerque, New Mexico.

Radiguet, M., Cotton, F., Vergnolle, M., Campillo, M., Walpersdorf, A., Cotte, N., Kostoglodov, V., 2012. Slow slip events and strain accumulation in the guerrero gap, mexico. Journal of Geophysical Research: Solid Earth 117 (B4).

Royer, A., Thomas, A., Bostock, M., 2015. Tidal modulation and triggering of low-frequency earthquakes in northern cascadia. Journal of Geophysical Research: Solid Earth 120 (1), 384-405.

Rubin, A. M., 2011. Designer friction laws for bimodal slow slip propagation speeds. Geochemistry, Geophysics, Geosystems 12 (4).

Rubin, A. M., Armbruster, J. G., 2013. Imaging slow slip fronts in cascadia with high precision cross-station tremor locations. Geochemistry, Geophysics, Geosystems 14 (12), 5371-5392.
Ryberg, T., Haberland, C., Fuis, G., Ellsworth, W., Shelly, D., 2010. Locating non-volcanic tremor along the san andreas fault using a multiple array source imaging technique. Geophysical Journal International 183 (3), 1485-1500.

Savard, G., Bostock, M. G., 2015. Detection and location of low-frequency earthquakes using cross-station correlation. Bulletin of the Seismological Society of America.

Schmidt, D., Gao, H., 2010. Source parameters and time-dependent slip distributions of slow slip events on the cascadia subduction zone from 1998 to 2008. Journal of Geophysical Research: Solid Earth 115 (B4).

Shelly, D. R., 2009. Possible deep fault slip preceding the 2004 parkfield earthquake, inferred from detailed observations of tectonic tremor. Geophysical Research Letters 36 (17).

Shelly, D. R., 2010. Migrating tremors illuminate complex deformation beneath the seismogenic san andreas fault. Nature 463 (7281), 648-652.

Shelly, D. R., Beroza, G. C., Ide, S., 2007a. Complex evolution of transient slip derived from precise tremor locations in western shikoku, japan. Geochemistry, Geophysics, Geosystems 8 (10).

Shelly, D. R., Beroza, G. C., Ide, S., 2007b. Non-volcanic tremor and low-frequency earthquake swarms. Nature 446 (7133), 305-307.

Thomas, A. M., Beroza, G. C., Shelly, D. R., 2016. Constraints on the source parameters of low-frequency earthquakes on the san andreas fault. Geophysical Research Letters.
Thomas, T. W., Vidale, J. E., Houston, H., Creager, K. C., Sweet, J. R., Ghosh, A., 2013. Evidence for tidal triggering of high-amplitude rapid tremor reversals and tremor streaks in northern cascadia. Geophysical Research Letters 40 (16), 4254-4259.

Ueno, T., Maeda, T., Obara, K., Asano, Y., Takeda, T., 2010. Migration of low-frequency tremors revealed from multiple-array analyses in western shikoku, Japan. Journal of Geophysical Research: Solid Earth 115 (B9).

Wech, A. G., 2010. Interactive tremor monitoring. Seismological Research Letters 81 (4), 664-669.

Wech, A. G., Bartlow, N. M., 2014. Slip rate and tremor genesis in cascadia. Geophysical Research Letters 41 (2), 392-398.

Wech, A. G., Creager, K. C., 2008. Automated detection and location of cascadia tremor. Geophysical Research Letters 35 (20).

Yager, R. R., Filev, D. P., 1994. Generation of fuzzy rules by mountain clustering. Journal of Intelligent & Fuzzy Systems 2 (3), 209-219.
<figure><img src="image_3.png" /><figcaption>Figure 1: Algorithm scheme. See main text for description.</figcaption></figure>
<figure><img src="image_4.png" /><figcaption>Figure 2: Examples of SSFs detected in the LFE catalog for the different time windows of our analysis. Grey arrows show SSF direction determined by the algorithm.</figcaption></figure>
<figure><img src="image_5.png" /></figure>

<figure><img src="image_6.png" /></figure>

<figure><img src="image_7.png" /></figure>

<figure><img src="image_8.png" /></figure>

Figure 3: Examples of SSFs detected in the tremor catalog for the different time windows of our analysis. Grey arrows show SSF direction determined by the algorithm.
<figure><img src="image_9.png" /><figcaption>Figure 4: Detected secondary slip fronts in Cascadia from August 2009 to December 2015 for time windows of 4 h, 8 h, 16 h and 32 h (tremor catalog). White to gray color scale indicates the density of tremor detections. Arrows show lengths and directions of the SSFs propagating along strike (red) and along dip (black). Rose diagrams show the number of SSFs in the different directions. Histograms show velocity distributions for SSFs propagating along strike (red) and along dip (black). Red and black vertical lines superimposed on histograms show average velocity for SSFs propagating along strike and along dip, respectively.</figcaption></figure>
<figure><img src="image_10.png" /></figure>

Figure 5: Detected secondary slip fronts during the 2003, 2004 and 2005 SSEs for time windows of 30 min, 1 h, 2, h, 4 h, 8 h, 16 h and 32 h (LFE catalog). White to gray color scale indicates the relative density of LFEs. Blue thick arrow show the approximate direction of the three SSEs (45° west from north). Red arrows show lengths and directions of SSFs propagating in the direction of SSE propagation ±45° (along strike). Black arrows show lengths and directions of SSFs propagating perpendicular to the direction of SSE propagation ±45° (along dip). Rose diagrams show the number of SSFs in the different directions. Histograms show velocity distributions for SSFs propagating along strike (red) and along dip (black). Red and black vertical lines superimposed on histograms show average velocity for SSFs propagating along strike and along dip, respectively.
<figure><img src="image_11.png" /><figcaption>Figure 6: Scheme illustrating how SSF length (*L*), width (*W*), and pulse length (*l*) are deduced from the analysis. *L* is defined as the length of the best projection axis (parallel to the propagation direction) (a and b). *l* is defined as the 2σ dispersion from the linear regression of the distance along the best projection axis as a function of time (a). *W* is defined as the 2σ spatial dispersion (perpendicular to the propagation direction) from the best projection axis (b).</figcaption></figure>
Moment - Duration scaling laws for regular earthquakes and slow slip phenomena

<figure><img src="image_12.png" /></figure>

Figure 7: Moment – duration scaling law for regular earthquakes and slow slip phenomena (gray shaded areas). Colored solid lines are isolines of equal moment rate. Black stars indicate regular earthquakes over a large range of magnitude (compilation of Peng and Gomberg (2010) augmented with M4.8-5.8 events from the GCMT catalog and recent mega-earthquakes (Delouis et al., 2010; Bletery et al., 2014, 2016)). Gray stars indicate tsunami earthquakes (compilation from Peng and Gomberg (2010)). Blue dots are SSEs (compilation of Peng and Gomberg (2010) augmented with events from Gao et al. (2012) and Liu (2014)). Orange dots are VLFEs (Ito and Obara, 2006; Matsuzawa et al., 2009; Ide and Yabe, 2014). Orange circles are unnamed 20-200 s long duration earthquakes assimilated to VLFEs (Ide et al., 2008). Red circles are LFE moment and duration averaged over 100 (Thomas et al., 2016) and >100,000 events (Bostock et al., 2015). Cyan dashed line shows an approximate moment threshold ($M_w \approx 5.7 - 6.0$) for geodetically detectable events. Red dashed line shows the limit of seismically detectable events (see text). In agreement with moment of large SSFs determined from strain meters (green square) (Hawthorne et al., 2016), SSFs detected in this study (green dots) span a moment – duration region that cannot be directly recorded either by common seismic or geodetic instrumentation.