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

# Double distance dependence in high-frequency ground motion along the plate boundary in Northern Chile

Benjamin Idini<sup>a,*</sup>, Sergio Ruiz<sup>b</sup>, Jean-Paul Ampuero<sup>c</sup>, Felipe Leyton<sup>b</sup>, Efraín Rivera<sup>d</sup>

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

## ARTICLE INFO

### Keywords:

Earthquake ground motion

Anelastic attenuation

Subduction tectonics

Seismic Hazard

Northern Chile

## ABSTRACT

In subduction zones, the forearc crust transitions from a highly fractured wedge near the trench to a less fractured basement near the volcanic arc. Here we study the role of wedge integrity on the frequency content of strong ground motion produced by subduction earthquakes in Northern Chile. Our data includes aftershocks from the 2014 $M_w$ 8.1 Iquique earthquake and intermediate-depth seismicity, including the 2005 $M_w$ 7.6 Tarapacá earthquake. We focus on the S-wave spectral decay parameter $\kappa$ obtained from strong ground motion in the frequency band 2–40 Hz, paying particular attention to how $\kappa$ varies with hypocentral distance. We report two distinct trends of $\kappa$ versus hypocentral distance. One trend applies to events near the trench, where we observe a rapid increase of $\kappa$ with hypocentral distance. Another trend applies to events near the bottom of the megathrust and greater depths, where we observe a relatively slow increase of $\kappa$ with hypocentral distance. We interpret this difference as resulting from a heterogeneous tectonic structure with lateral variations in anelastic attenuation (i.e., attenuation seems stronger in the fractured wedge near the trench). Our results improve our understanding of the role of high-frequency ground motion in seismic hazard analysis near the Chilean coastline. In addition, our results stress the need for a better characterization of anelastic attenuation in the Chilean subduction zone, specially close to the highly fractured wedge near the trench.

## 1. Introduction

The spectrum of far-field strong ground motion (SGM) follows a simple shape determined by earthquake magnitude at low frequency and anelastic attenuation at high frequency. At low frequency, the SGM acceleration spectrum grows in amplitude proportionally to $f^2$ until reaching the corner frequency $f_c$, a frequency that depends on earthquake size (Brune, 1970). The SGM spectrum turns flat at frequencies higher than said corner frequency. Anelastic attenuation imprints an exponential decay with frequency on the SGM spectrum. This effect signals the irreversible loss of energy in seismic waves as they propagate through a dissipative medium from the earthquake focus to the recording station. The parameter $\kappa$ quantifies the spectral decay in the SGM spectrum $a(f)$ at frequencies $f > f_c$, following (Anderson and Hough, 1984)

$$a(f) = a_0 \exp(-\pi \kappa f), \quad (1)$$

where $a_0$ is the SGM spectrum defined by earthquake source parameters in the absence of anelastic attenuation.

The parameter $\kappa$ represents a quantity integrated along the path that depends on the total distance traveled by seismic waves and the dissipative properties of the propagating medium (Hough and Anderson, 1988),

$$\kappa = \int_{0}^{R} (\beta Q)^{-1} dr, \quad (2)$$

where $\beta$ is the S-wave speed, $Q$ is the attenuation quality factor (Knopoff, 1964), and $R$ is the path length covered by seismic waves along a ray. The $Q$ factor represents the fraction of energy lost in a cycle of seismic wave oscillation, with a lower $Q$ indicating more dissipation.

A significant portion of the integral in Eq. (2) involves short distances $\delta r \ll R$ near the recording station (i.e., site effects), where $\beta$ and $Q$ are relatively small compared to the rest of the path. As a consequence, it is often seen in the literature that Eq. (2) is written

$$\kappa = \kappa_0 + \int_{0}^{R-\delta r} (\beta Q)^{-1} dr \approx \kappa_0 + \kappa_r R, \quad (3)$$

<footnote><sup>a</sup> Department of Earth and Planetary Sciences, University of California Santa Cruz, 1156 High St, Santa Cruz, 95064, CA, United States</footnote>

<footnote><sup>b</sup> Departamento de Geofísica, Universidad de Chile, Blanco Encalada 2002, Santiago, 8320000, Region Metropolitana, Chile</footnote>

<footnote><sup>c</sup> Geoazur, Université Cote d Azur, 250 Rue Albert Einstein, Sophia Antipolis, 06560, Valbonne, France</footnote>

<footnote><sup>d</sup> Departamento de Obras Civiles, Universidad Católica del Maule, Av San Miguel 3605, Talca, 06560, Region del Maule, Chile</footnote>

<footnote>* Corresponding author.</footnote>

<footnote>E-mail address: bidini@ucsc.edu (B. Idini).</footnote>
<figure><img src="image_3.png" /><figcaption>Fig. 1. Profile of the Northern Chile seismicity, showing events included in the CSN strong-motion database, and a cartoon of the major tectonic structures discussed in the text. The three event groups shown represent a classification based on distance along dip from the trench. Symbol sizes scale with magnitude.</figcaption></figure>

where $κ_0$ encapsulates the anelastic attenuation near the recording station and $κ_r$ represents the path-averaged anelastic attenuation. Previous calculations of the parameter $κ$ cover various regions around the world (Fernández et al., 2010; Kilb et al., 2012; Ktenidou et al., 2014, e.g.), including the Chilean subduction zone (Neighbors et al., 2015; Cabas et al., 2017; Pozo et al., 2023). Unlike previous work, here we analyze variations in $κ_r$ when seismic waves travel along paths that include distinct tectonic structures in the Northern Chile subduction zone.

The Northern Chile subduction zone hosts morphological structures that vary laterally from the subduction trench to the volcanic arc (Fig. 1). Close to the trench, the marine forearc comprises an eroded and fractured continental wedge resulting from the collision of the Nazca plate and the South American plate (Contreras-Reyes et al., 2012; Geersen et al., 2015; León-Ríos et al., 2016; Ma et al., 2022). Evidence for the existence of fractured material in this wedge comes from seismicity located in crustal faults within the wedge (León-Ríos et al., 2016; Cesca et al., 2016; Pastén-Araya et al., 2021), lateral density gradients (Maksymowicz et al., 2018), the surficial expression of a dense network of marine crustal faults (Reginato et al., 2020; González et al., 2023), and lateral seismic wave velocities that decrease near the wedge (Contreras-Reyes et al., 2012; Pastén-Araya et al., 2021).

Here we evaluate how path effects across these distinct tectonic structures influence high-frequency ground motion recorded in Northern Chile. A similar tectonic structure is visible further south in North-Central Chile (Contreras-Reyes et al., 2014), suggesting that our analysis presented here may carry implications for multiple areas along Chile's convergent margin. A heavily fractured and eroded medium will likely produce more dissipation (i.e., lower $Q$) than a more pristine medium. We pay particular attention to correlations between the location of the earthquake focus along the plate boundary and the path-averaged parameter $κ_r$ calculated over the S-wave window of the ground motion record. Our study focuses in a small number of stations that offer advantageous observation capabilities (i.e., a large number of recorded events and favorable site conditions).

## 2. Data and methods

### 2.1. Strong ground motion in Northern Chile

Seismicity in the Northern Chile subduction zone results from the convergence between the Nazca and South American plates at a rate

of ~7 cm yr$^{-1}$ (Angermann et al., 1999). Our study area is near the city of Iquique, where the Centro Sismológico Nacional (CSN) and the Integrated Plate Boundary Observatory Chile (IPOC) networks capture one of the richest pictures of the Chilean seismicity. All the ground motion records used in this study were compiled from data obtained by the strong-motion instruments in the CSN and IPOC seismic arrays. Our main database consists of earthquakes identified in the CSN earthquake catalog (Figs. 1 and 2), which are mainly aftershocks of the $M_{w}$ 8.1 2014 Iquique (Ruiz et al., 2014; Schurr et al., 2014; Twardzik et al., 2022) and the $M_{w}$ 7.6 2005 Tarapacá (Peyrat et al., 2006; Delouis and Legrand, 2007) earthquakes. The records in our main database contain events with local magnitude $Ml > 4$ and horizontal PGA > 0.005 g.

In addition to our main database, we generate a secondary ground motion database from mining continuous seismic ground motion for events with $Ml > 2$ recorded at station PSGCX during the years 2007–2014 (Fig. 3). We exclude this secondary database from estimations of the parameter $κ$, given that events with $Ml < 4$ most likely overlap the corner frequency with the frequency range where we estimate the anelastic attenuation (2–40 Hz). Despite the presence of confounding effects, this secondary database allows us to qualitatively assess the effect of anelastic attenuation over a larger number of events spanning magnitudes $Ml$ 2–5.5 (Fig. 3). We chose PSGCX station to compile this secondary database because anelastic attenuation localized near the station (or $κ_0$) is expected to be small given the surrounding fast seismic wave speeds ($V_{S30}$ = 1639 m/s Leyton et al., 2018a). We chose the years 2007–2014 because this time window includes enough aftershocks from the largest earthquakes recorded in the area.

We classify the regional seismicity in three groups: up-dip events, down-dip events, and intermediate depth events. Up-dip events consist of earthquakes with hypocenters located as far as 80 km from the trench. Down-dip events consist of earthquakes with hypocenters located farther than 80 km from the trench and depths shallower than 60 km. Intermediate depth events consist of earthquakes with depths ranging from 60 to 350 km (Fig. 1). We explicitly leave events located to the west of the trench out of this study given they are noticeably rare in the CSN earthquake catalog.

### 2.2. The parameter $κ$

We obtain the parameter $κ$ from S-wave Fourier amplitude spectra (FAS) following a simple procedure where we fit a linear model (Eq. (1)) to the smoothed FAS (Anderson and Hough, 1984). We
<figure><img src="image_4.png" /><figcaption>**Fig. 2.** Earthquake hypocenters of the events contained in the CSN strong-motion database. Black triangles represent the stations considered in this study. White contour lines represent the plate interface shown every 20 km depth. The region shown is same included in Fig. 3.</figcaption></figure>

constrain the linear fit to the frequency band $f_1 - f_2$, where $f_1$ is the frequency where the FAS peaks and $f_2$ the frequency where the spectrum decays below 3 times the noise level. We calculate $\kappa$ only on earthquakes with magnitude $M_l > 4$, where the corner frequency $f_c$ is expected to be lower than $f_1$.

We compute the FAS in time windows starting at the S-wave arrival time and ending at the time when the cumulative integral of the squared velocity reached 95% of its final value. Lancieri et al. (2012) used a similar S-wave window selection procedure to study the spectral properties of aftershocks of the $M_{w}$ 7.8 2007 Tocopilla earthquake in Northern Chile. Shorter S-wave windows that include less of the scattered energy present in the coda have been used in other studies (Houtte et al., 2014, e.g.). We found that the spectra were nearly identical using a window length of 5 or 20 sec (Fig. 4), consistent with the notion that even at high frequency the spectrum is dominated by the contribution of the direct S-wave arrival. The selected S-wave windows were detrended and cosine-tapered at 1% of their length, then zero-padded at both sides up to a total window length equal to five times the original one. The spectra were finally smoothed using a **Konno and Ohmachi (1998)** smoothing function with a bandwidth exponent of 40 (Laurendeau et al., 2017, e.g.). Noise windows were extracted from the 20 s preceding the manually-picked P-wave arrival and their spectra were computed following the same procedure used for S-wave windows (Fig. 4).

# 3. Results

## 3.1. High-frequency ground motion down dip the plate boundary

We first qualitatively assess the ground motion spectra (i.e., FAS) of records in our secondary database consisting of events recorded at station PSGCX ($V_{S30}$ = 1639 m/s (Leyton et al., 2018a)). We observe that the ground motion spectrum decays with a different incline for each group of events (Fig. 3). Up-dip events show the steepest decay in high frequency amplitude, whereas down-dip events better preserve high-frequency amplitude out of the three event groups. Given that PSGCX locates directly on top of down-dip events, the path for down-dip events is the shortest of the three groups (Fig. 1). Such relatively short path

naturally produces the lowest anelastic attenuation, represented by the path-integrated $\kappa$ (Eq. (3)), and roughly explains the flatter shape and elevated high-frequency content of down-dip event spectra.

Up-dip and intermediate depth events show a striking difference in decay at high frequency at all observed magnitudes (Fig. 3). This difference cannot be attributed to the length of path because both groups of events locate at roughly similar distance from PSGCX (Fig. 1). A natural alternative explanation comes from the globally-observed increment with depth in seismic wave speed $\beta$ (see Husen et al., 1999 for a regional 1D model of $\beta$). Up-dip events cover a longer horizontal path along the shallow layers of slow $\beta$ compared with the shorter vertical paths of intermediate depth events (Fig. 1). A $Q$ that increases with depth produces a similar effect.

A potential lateral variation of $Q$ introduces a new confounding variable to the interpretation of the differences in high-frequency spectral decay among event groups. An eroded and fractured wedge hypothetically leads to a lower average $Q$ for up-dip events compared with intermediate depth events. In the rest of this paper, we assess the role of an eroded and fractured wedge in the high-frequency spectral decay by comparing the trends of the parameter $\kappa$ versus distance among the three event groups.

## 3.2. Two trends in the parameter $\kappa$ versus hypocentral distance

We now calculate the parameter $\kappa$ for all events in our main database. The magnitude ($M_l > 4$) and PGA (PGA > 0.005 g) thresholds in this database allow us to mostly avoid corner frequency effects and obtain a high signal-to-noise ratio. In addition to PSGCX, we concentrate in stations T01A ($V_{S30}$ = 339 m/s Leyton et al., 2018a) and T07A ($V_{S30}$ = 326 m/s Leyton et al., 2018a), which are among the stations with the largest number of recordings in the area for the considered time window.

We observe two strikingly distinct trends in the parameter $\kappa$ versus hypocentral distance (i.e., $\kappa_r$) when comparing up-dip and intermediate depth events (Fig. 5). For simplicity, we split the $\kappa$ results between events located to the east and west directions from stations PSGCX, T01A, and T07A. In each station, all up-dip events locate to the west and all intermediate depth events locate to the east. Down-dip events
<figure><img src="image_5.png" /><figcaption>Fig. 3. The Fourier amplitude spectrum (FAS) of S-waves recorded at station PSGCX for events located at different distances along dip: (a) up-dip interplate events located close to the trench, (b) down-dip interplate events in the bottom part of the megathrust seismogenic zone, and (c) intermediate-depth intraplate events. Color indicates event magnitude. Black dashed lines indicate the expected spectral shape of the Brune (1970) model for a corner frequency of 3 Hz. (d) The seismicity included in the CSN 2014 catalog that has been recovered from continuous data at PSGCX. The black boxes contain the hypocenters of events shown in a, b, and c.</figcaption></figure>

end up in the east or west groups depending on the station and sometimes they split between east and west. We fit separate linear models (Eq. (3)) to the $κ$ versus distance results for the east and west groups to obtain two pairs of parameters $κ_0$ and $κ_r$.

Up-dip events show a ~10 times higher $κ_r$ than intermediate depth events in all three stations. A previous independent calculation of $κ_r$ at T07A and other stations (Pozo et al., 2023) confirms the result reported here. The anelastic attenuation at zero hypocentral distance (i.e., $κ_0$) is lower at PSGCX than at T01A and T07A. This observation is consistent with PSGCX exhibiting lower site effects given the higher seismic velocities near the recording site (i.e., higher $V_{S30}$) (Leyton et al., 2018b).

For each of the stations analyzed here, we observe a disagreement between $κ_0$ obtained from events to east and west of the station (Fig. 5). We interpret this difference as the result of a tectonic structure with lateral heterogeneity in attenuation, a possibility discussed in further detail in the next section. This observed difference in $κ_0$ agrees with a previous result where station T07A and others exhibit different $κ_0$ for interplate and intermediate depth events (Pozo et al., 2023), an event

classification that roughly follows our own event classification into east/west events from the station. Further interpretation of $κ_0$ results and site effects in the Chilean subduction zone can be found in Pozo et al. (2023).

We have centered our analysis on coastal stations because they favor the observation of ray paths that produce a broad range of hypocentral distances. Given the path geometry in Northern Chile (Fig. 1), the coastal stations used here allow us to access hypocentral distances that broadly range from 20 to 200 km (Fig. 5). Inland stations (i.e., > 200 km away from the trench), on the contrary, only access along-plate events from a relatively long hypocentral distance (Fig. 1) that spans a range from roughly 150 to 200 km. This narrow range in hypocentral distance prevents us from obtaining an accurate enough linear fit that may constrain the slope $κ_r$ of the $κ$ vs. $R$ linear relationship (Eq. (3)).

## 4. Discussion

Our $κ_r$ results confirm previous suggestions of depth-varying high-frequency ground motion in Japan (Ye et al., 2013; Piña-Valdés et al.,
<figure><img src="image_6.png" /><figcaption>Fig. 4. Example of the selection of S-wave and noise windows for the calculation of the parameter κ. (top) Ground motion acceleration record. (bottom) Smooth Fourier amplitude spectrum (FAS) of the windows selected from the recording on top using S-wave windows of varying length.</figcaption></figure>

2018a) and Northern Chile (Piña-Valdés et al., 2018b). In the following, we discuss whether it is possible to attribute the effect to previously known tectonic structures, the role of source effects, and the main implications of these results for regional seismic hazard studies in subduction zones, with particular emphasis in Chile.

## 4.1. A constant Q

Is it possible to explain a ~10 times difference in κr by invoking a constant Q and a depth-varying β? We can consider that

$$ \kappa_r \sim (\bar{Q}\bar{\beta})^{-1}, \qquad (4) $$

with the bar on top of the variables indicating an average taken along the path. According to the 1D regional seismic velocity model by Husen et al. (1999), we consider a conservative estimate for β̃ to be 3 and 5 km/s for up-dip and intermediate depth events, respectively. This difference in β̃ can explain only less than a factor of ~2 difference in κr. Therefore, a medium with a depth-varying β and a constant Q does not explain the difference observed in κr between up-dip and intermediate depth events.

## 4.2. A crust with low Q

We require an additional ~5 times difference between the Q̄ values of up-dip and intermediate depth events to explain the tenfold difference in κr reported here. Can this additional attenuation result from a continental crust that has a significantly lower quality factor Qc compared to that of the upper mantle Qm (Fig. 1)? For waves that travel

<figure><img src="image_7.png" /><figcaption>Fig. 5. The parameter κ as a function of the hypocentral distance. White (black) markers represent events with hypocenters to the east (west) of the station. Only selected stations are shown. The dashed red lines represent least-squares linear regressions for events labeled as west and east. The error bars represent 2-σ uncertainty. The event classification is simplified here for brevity: interplate events correspond to up-dip and down-dip events, and intraslab events correspond to intermediate-depth events.</figcaption></figure>

across a two-layer model, the κr of intermediate depth events can be written as

$$ \kappa_r \approx \frac{x}{\beta_c Q_c} + \frac{1-x}{\beta_m Q_m}, \qquad (5) $$

where x represents the fraction of the path that lies in the continental crust and the subscripts c and m represent medium properties in the crust and mantle, respectively.
<figure><img src="image_8.png" /><figcaption>Fig. 6. A two-layer model with depth-varying seismic speed $\beta$ and $Q$. The crust path fraction $x$ represents the fractional path involving a highly dissipative crustal layer with relatively low $Q_c$. A poorly dissipative upper mantle with high $Q_m$ lies below the crust (Fig. 1). The orange line describes models with an average $Q$ that explains the tenfold difference in $\kappa_r$ observed between up-dip and intermediate depth events. In the blue line, the $\kappa_r$ difference has been arbitrarily relaxed to fivefold. We use the conservative values $\beta_c \sim 3$ km/s and $\beta_m \sim 5$ km/s in Eq. (6).</figcaption></figure>

We evaluate Eq. (5) to show that such a two-layer model requires a relatively thin dissipative crust to explain the reported difference in anelastic attenuation between up-dip and intermediate depth events (Fig. 6). We can write the crustal path fraction that intermediate-depth earthquakes require to fit our observations as

$$x \approx \left( \frac{\beta_m Q_m \kappa_r^\perp}{\beta_c Q_c \kappa_r^\perp} - 1 \right) \left( \frac{\beta_m Q_m}{\beta_c Q_c} - 1 \right)^{-1}. \quad (6)$$

Our observations indicate $\kappa_r^\perp / \kappa_r^\parallel \sim 10$, where $\kappa_r^\perp$ is the $\kappa_r$ of up-dip events. To derive Eq. (6), we considered $\kappa_r^\perp \le (\beta_c Q_c)^{-1}$ because up-dip events may include a path that partially covers the less dissipative layer with $Q_m > Q_c$.

The dissipative crust must be thinner than $\le 12$ km for an intermediate depth event with a path length of $\sim 120$ km, and thinner than $\le 6$ km if we were to limit the vertical increase in $Q$ between layers to the reasonable ratio $Q_m/Q_c \sim 10$ (Fig. 6). We obtain a roughly 2 times thicker crust when we arbitrarily relax the observed difference in $\kappa_r$ between up-dip and intermediate depth events from tenfold to fivefold.

A thin crust represents a limitation of the depth-varying $Q$ model. Up-dip events may locate as deep as 30 km depth, implying that much of the event path would not land within the dissipative crust. The situation worsens when we consider the curvature of seismic rays due to the increasing seismic speed $\beta$ with depth. As a consequence, it becomes challenging to find a two-layer model that satisfies our $\kappa_r$ observations.

Another limitation in the depth-varying $Q$ model is the strong increment in $Q$ with depth. Tomography applied to subduction zones usually reveals a highly heterogeneous $Q$ structure. Instead of a vertical structure where $Q$ increases with depth, regions of low $Q$ seem to correlate with the presence of partial melts or hydrated material (Schurr et al., 2006; Chen and Clayton, 2012). Furthermore, the previously hypothesized global asthenosphere suggests a $Q$ that decreases with depth due to the presence of partial melts (Durek and Ekström, 1996), the opposite trend required by our depth-varying $Q$ model.

### 4.3. The laterally heterogeneous structure of Q

Our analysis above shows that our $\kappa_r$ results cannot be explained by a model with homogeneous $Q$ and depth-varying seismic speed. Our $\kappa_r$ results require a heterogeneous $Q$ structure. In principle, a purely depth-varying $Q$ model explains the double distance dependence of $\kappa$ with hypocentral distance. However, the depth-varying model requires a dissipative thin crust overlying an upper mantle with a $Q$ that is 10 times higher or more.

Instead of a thin dissipative crust and an abrupt increment in $Q$ with depth, we can explain the two different $\kappa_r$ by invoking a lateral variation in $Q$ due to the presence of an eroded and fractured wedge along the path of up-dip events. In this scenario, up-dip events travel across a medium with a low average $Q_w$, whereas intermediate depth events travel across a much more pristine continental basement with a high $Q_b$. The reported $\kappa_r$ can be simply explained if we set $Q_b/Q_w \sim 5$ and maintain the regional 1D model for $\beta$ by Husen et al. (1999). Given the $\kappa_r$ results shown in Fig. 5, this model implies a $Q_b$ ranging from 2000 to 3000 and a $Q_w$ ranging from 500 to 1000.

The presence of an eroded and fractured wedge with fluid-filled cracks could potentially support the value of $Q_b/Q_w \sim 5$ required to explain our $\kappa_r$ results. However, previous studies of attenuation in Northern Chile do not cover the eroded and fractured continental wedge discussed here (Schurr et al., 2006). We require additional field studies to image the attenuation in the eroded and fractured wedge, much of which is underwater. Distributed acoustic sensing on submarine fiber optic cables constitute a new promising technology that may further illuminate an area that has remained hard to resolve to this day (Williams et al., 2019; Sladen et al., 2019; Lior et al., 2022).

### 4.4. New considerations when estimating $\kappa_0$

The laterally heterogeneous $Q$ structure discussed above carries important consequences for the determination of $\kappa_0$, the $R = 0$ extrapolation of a linear $\kappa$ vs. $R$ plot (Eq. (3)). This parameter plays an important role in the modeling of the predicted high-frequency ground motion at a particular site in seismic hazard studies (Houtte et al., 2014, e.g.). Multiple consecutive linear slopes may be required to fit the nonlinearity of a $\kappa$ vs. $R$ plot (Fig. 5) when rays pass through strongly heterogeneous attenuation structures; such is the case for up-dip events with rays that pass through a continental wedge with varying degree of erosion and fractures in Northern Chile (Fig. 1). Alternatively, events with strongly nonlinear $\kappa(R)$ attenuation should be avoided in the determination of $\kappa_0$ due to the risk of extrapolating a linear model into a region where the attenuation properties may be drastically different compared to the trend captured in the data. For example, this is the case of west events in Fig. 5, where $\kappa_0$ is systematically underestimated when compared to the $\kappa_0$ from east events.

### 4.5. Source effects in high-frequency ground motion

In this section, we analyze the possibility that the two strikingly distinct $\kappa_r$ reported here correspond to a misinterpretation of a depth-varying source effect. In the past, a depth-varying high-frequency content has been attributed to systematic changes of fault properties with depth. In these models, the base of the megathrust is rheologically heterogeneous, comprising compact brittle asperities embedded in a ductile fault zone matrix (Meng et al., 2011). Abrupt changes of the rupture speed as these asperities break may explain the stronger high-frequency radiation at depth during large earthquakes (Huang et al., 2012; Galvez et al., 2014).

The lines of evidence for this alternative are rich and diverse, including teleseismic and regional observations. Teleseismic back-projection shows that the sources of radiation at $\sim 1$ Hz during megathrust earthquakes is concentrated in the deep portions of the seismogenic zone (Ide et al., 2011; Meng et al., 2011; Lay et al., 2012; Koper et al., 2012). Analysis of strong-motion data at frequencies between 0.1 Hz and 10 Hz also locate the sources of high-frequency radiation during Japanese events at the bottom of the seismogenic zone (Kurahashi and Irikura, 2011; Lay et al., 2012; Kurahashi and Irikura, 2013; Ye et al., 2013).

Rupture speed and stress drop increase systematically with depth due to increasing wave speeds, elastic moduli, and confining stress (Valée, 2013, e.g.), leading to a higher absolute amplitude in high-frequency content. While these source effects modify the absolute amplitude of high-frequency content, they do not alter the spectral
decay registered in the parameter $κ$. For a given seismic moment $M_0$, an event with higher corner frequency has stronger high-frequency radiation. At high frequencies ($f > f_c$), the amplitude of the far-field acceleration spectrum, assuming the $ω$-squared model (Brune, 1970), is

$$a(f) \propto M_0 f_c^2. \qquad (7)$$

For a circular rupture the corner frequency scales as (Eshelby, 1957; Madariaga, 1976)

$$f_c \propto v_r(\Delta\tau/M_0)^{1/3}, \qquad (8)$$

where $v_r$ is rupture speed and $\Delta\tau$ is stress drop. The amplitude of the far-field acceleration spectrum at $f > f_1 > f_c$, including the spectral decay in Eq. (1), is then

$$a(f) \propto v_r^2 M_0^{1/3} \Delta\tau^{2/3} \exp(-\pi\kappa f). \qquad (9)$$

As suggested in Eq. (9), the $κ_r$ results reported here are independent of additional contributions to high-frequency content from a depth-varying rupture speed and/or stress drop.

The lateral variations in attenuation invoked here may explain the observed high-frequency content in regional data to first order, but are not suitable to explain the teleseismic observations. Direct teleseismic phases are radiated downward from the plate boundary, thus they do not interact with the eroded and fractured wedge invoked here. Consequently, a depth-dependent source effect may still be required to reconcile the overall amplitudes of high-frequency radiation observed at teleseismic and regional distances.

## 4.6. Seismic hazard in Northern Chile coastal zones

The systematically stronger high-frequency ground motion from down-dip and intermediate depth earthquakes identified here has important implications for seismic hazard analyses in subduction zones. This is especially true for Chilean coastal cities overlying the bottom of the megathrust seismogenic zone (e.g., Iquique) or inland cities on top of regions capable of producing large intermediate-depth earthquakes. The intense damage caused by the 1939 Chillán Ms 7.8 intermediate-depth earthquake (Beck et al., 1998; Ruiz and Madariaga, 2018), the deadliest Chilean earthquake in history, may be related to enhanced high-frequency ground motion hitting the stiff adobe buildings that were common at the time. A recent example of the unfavorable combination of stiff structures with unusually elevated high-frequency content was the $M_{w}$ 6.7 2019 Coquimbo earthquake (Ruiz et al., 2019). This event struck the coastal cities of La Serena and Coquimbo causing severe damage despite its moderate magnitude. The most recent Chilean ground-motion prediction equation models (Idini et al., 2017; Montalva et al., 2017, 2022) capture the differences in ground motion between interplate and intermediate-depth earthquakes but they do not consider the remarkably enhanced high-frequency ground motion of down-dip events. When the effects of higher $\bar{Q}$ and higher $\Delta\tau$ (as estimated for down-dip events in Japan) are combined using Eq. (9) (with $r \sim 50$ km, $\bar{\beta} \sim 3$ km/s, $\bar{Q}_{up} \sim 800$, $\bar{Q}_{down} \sim 4000$), down-dip events can show $\sim 3$ times higher spectral amplitude at 8 Hz in a coastal city. Low-story and other stiff structures or structural components are particularly exposed to this high-frequency hazard given that they oscillate with relatively short periods.

## 5. Conclusions

Here we report a double dependency of high-frequency strong ground motion with earthquake location along the plate boundary in the Northern Chile subduction zone. We quantify high-frequency strong ground motion using the spectral decay parameter $κ$. We calculate $κ$ in a frequency band ranging from 2 to 40 Hz for small-to-moderate interplate and intermediate-depth earthquakes with local magnitude $M_l \ge 4$. We interpret our observation as differences in attenuation

properties of the medium encountered by waves from earthquakes along the plate boundary as they reach recording stations located in Northern Chile. Our results suggest that earthquakes located directly below coastal zones in Northern Chile represent additional seismic hazard compared to earthquakes located up-dip on the plate boundary due to the systematically larger amplitude of incoming high-frequency waves.

According to the foregoing interpretation, attenuation in Northern Chile is remarkably low for waves that do not travel through the continental wedge. Such is the case for down-dip interplate events, whose hypocenters are located below the Chilean coastline (Fig. 2). For up-dip events, the strong attenuation in the marine forearc wedge combined with large distances to the Chilean coastline lead to ground motions relatively depleted in high-frequency content. Conversely, for down-dip events, the combined effects of high $\bar{Q}$ and short source-receiver distance lead to ground motions with remarkably elevated high-frequency content. For instance, the $M_l \sim 4$ down-dip events show flat high-frequency spectra up to 40 Hz (Figs. 2b and 3b) and presumably even to higher frequencies. Low-attenuation sites like PSGCX in the Northern Chile subduction zone make this region ideal for conducting earthquake source studies at very high frequencies (several tens of Hz), which can provide unique constraints on the small-scale processes that control earthquake rupture dynamics.

## Data and resources

The CSN strong-motion database is available at evtdb.csn.uchile.cl (last access: September 2016). The continuous acceleration data for station PSGCX can be requested on the German Research Center for Geosciences' website at geofon.gfz-potsdam.de. The CSN earthquake catalog is available at the CSN's webpage www.sismologia.cl. Some figures were made using the Generic Mapping Tools v.5.2.1 (Wessel and Smith, 1998).

## CRediT authorship contribution statement

**Benjamin Idini:** Conceptualization, Data curation, Formal Analysis, Investigation, Methodology, Project administration, Software, Visualization, Writing – original draft, Writing – review & editing. **Sergio Ruiz:** Conceptualization, Funding acquisition, Supervision, Writing – review & editing. **Jean-Paul Ampuero:** Formal Analysis, Validation, Writing – review & editing. **Felipe Leyton:** Data curation, Resources. **Efraín Rivera:** Data curation, Formal analysis.

## Declaration of competing interest

The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.

## Data availability

The data is publicly available in CSN repositories.

## Acknowledgments

We are thankful to the Centro Sismológico Nacional (CSN) and the Integrated Plate Boundary Observatory Chile (IPOC) for acquiring and storing the strong ground motion data used in this study. This work was partially supported by the Programa de Riesgo Sísmico (PRS) of the Universidad de Chile. S. R. acknowledges funding by FONDECYT grant N° 1200779.
## References

Anderson, J.G., Hough, S.E., 1984. A model for the shape of the Fourier amplitude spectrum of acceleration at high frequencies. Bull. Seismol. Soc. Am. 74 (5), 1969-1993.

Angermann, D., Klotz, T., Reigber, C., 1999. Space-geodetic estimation of the nazca-south america Euler vector. Earth Planet. Sci. Lett. 171 (3), 329-334.

Beck, S., Barrientos, S., Kausel, E., Reyes, M., 1998. Source characteristics of historic earthquakes along the central Chile subduction askew et alzone. J. South Amer. Earth Sci. 11 (2), 115-129.

Brune, J.N., 1970. Tectonic stress and the spectra of seismic shear waves from earthquakes. J. Geophys. Res. 75 (26), 4997-5009.

Cabas, A., Rodriguez-Marek, A., Bonilla, L.F., 2017. Estimation of site-specific kappa (κ 0)-consistent damping values at KiK-net sites to assess the discrepancy between laboratory-based damping models and observed attenuation (of seismic waves) in the field. Bull. Seismol. Soc. Am. 107 (5), 2258-2271.

Cesca, S., Grigoli, F., Heimann, S., Dahm, T., Kriegerowski, M., Sobiesiak, M., Tas-sara, C., Olcay, M., 2016. The M w 8.1 2014 Iquique, Chile, seismic sequence: a tale of foreshocks and aftershocks. Geophys. J. Int. 204 (3), 1766-1780.

Chen, T., Clayton, R.W., 2012. Structure of central and southern Mexico from velocity and attenuation tomography. J. Geophys. Res.: Solid Earth 117 (B9).

Contreras-Reyes, E., Becerra, J., Kopp, H., Reichert, C., Díaz-Naveas, J., 2014. Seismic structure of the north-central Chile convergent margin: Subduction erosion of a paleomagmatic arc. Geophys. Res. Lett. 41 (5), 1523-1529.

Contreras-Reyes, E., Jara, J., Grevemeyer, I., Ruiz, S., Carrizo, D., 2012. Abrupt change in the dip of the subducting plate beneath north Chile. Nat. Geosci. 5 (5), 342-345.

Delouis, B., Legrand, D., 2007. Mw 7.8 tarapaca intermediate depth earthquake of 13 June 2005 (northern Chile): Fault plane identification and slip distribution by waveform inversion. Geophys. Res. Lett. 34 (1).

Durek, J.J., Ekström, G., 1996. A radial model of anelasticity consistent with long-period surface-wave attenuation. Bull. Seismol. Soc. Am. 86 (1A), 144-158.

Eshelby, J.D., 1957. The determination of the elastic field of an ellipsoidal inclusion, and related problems. Proc. R. Soc. Lond. Ser. A Math. Phys. Eng. Sci. 241 (1226), 376-396.

Fernández, A.I., Castro, R.R., Huerta, C.I., 2010. The spectral decay parameter kappa in northeastern Sonora, Mexico. Bull. Seismol. Soc. Am. 100 (1), 196-206.

Galvez, P., Ampuero, J.-P., Dalguer, L.A., Somala, S.N., Nissen-Meyer, T., 2014. Dynamic earthquake rupture modelled with an unstructured 3-D spectral element method applied to the 2011 M 9 Tohoku earthquake. Geophys. J. Int. 198 (2), 1222-1240.

Geersen, J., Ranero, C.R., Barckhausen, U., Reichert, C., 2015. Subducting seamounts control interplate coupling and seismic rupture in the 2014 Iquique earthquake area. Nature Commun. 6 (1), 8267.

González, F., Bello-González, J., Contreras-Reyes, E., Tréhu, A., Geersen, J., 2023. Shallow structure of the northern Chilean marine forearc between 19°S-21°S using multichannel seismic reflection and refraction data. J. South Amer. Earth Sci. 123, 104243.

Hough, S., Anderson, J., 1988. High-frequency spectra observed at Anza, California: implications for Q structure. Bull. Seismol. Soc. Am. 78 (2), 692-707.

Houtte, C.V., Ktenidou, O.-J., Larkin, T., Holden, C., 2014. Hard-site κ 0 (kappa) calculations for Christchurch, New Zealand, and comparison with local ground-motion prediction models. Bull. Seismol. Soc. Am. 104 (4), 1899-1913.

Huang, Y., Meng, L., Ampuero, J.-P., 2012. A dynamic model of the frequency-dependent rupture process of the 2011 Tohoku-Oki earthquake. Earth Planets Space 64, 1061-1066.

Husen, S., Kissling, E., Flueh, E., Asch, G., 1999. Accurate hypocentre determination in the seismogenic zone of the subducting nazca plate in northern Chile using a combined on-/offshore network. Geophys. J. Int. 138 (3), 687-701.

Ide, S., Baltay, A., Beroza, G.C., 2011. Shallow dynamic overshoot and energetic deep rupture in the 2011 M w 9.0 Tohoku-Oki earthquake. Science 332 (6036), 1426-1429.

Idini, B., Rojas, F., Ruiz, S., Pastén, C., 2017. Ground motion prediction equations for the Chilean subduction zone. Bull. Earthq. Eng. 15, 1853-1880.

Kilb, D., Biasi, G., Anderson, J., Brune, J., Peng, Z., Vernon, F.L., 2012. A comparison of spectral parameter kappa from small and moderate earthquakes using southern California ANZA seismic network data. Bull. Seismol. Soc. Am. 102 (1), 284-300.

Knopoff, L., 1964. Q. Reviews of geophysics.

Konno, K., Ohmachi, T., 1998. Ground-motion characteristics estimated from spectral ratio between horizontal and vertical components of microtremor. Bull. Seismol. Soc. Am. 88 (1), 228-241.

Koper, K.D., Hutko, A.R., Lay, T., Sufri, O., 2012. Imaging short-period seismic radiation from the 27 february 2010 Chile (Mw 8.8) earthquake by back-projection of P, PP, and PKIKP waves. J. Geophys. Res.: Solid Earth 117 (B2).

Ktenidou, O.-J., Cotton, F., Abrahamson, N.A., Anderson, J.G., 2014. Taxonomy of κ: A review of definitions and estimation approaches targeted to applications. Seismol. Res. Lett. 85 (1), 135-146.

Kurahashi, S., Irikura, K., 2011. Source model for generating strong ground motions during the 2011 off the Pacific coast of tohoku earthquake. Earth Planets Space 63, 571-576.

Kurahashi, S., Irikura, K., 2013. Short-period source model of the 2011 M w 9.0 off the Pacific coast of Tohoku earthquake. Bull. Seismol. Soc. Am. 103 (2B), 1373-1393.

Lancieri, M., Madariaga, R., Bonilla, F., 2012. Spectral scaling of the aftershocks of the Tocopilla 2007 earthquake in northern Chile. Geophys. J. Int. 189 (1), 469-480.

Laurendeau, A., Bonilla, L., Mercerat, D., Courboulex, F., Alvarado, A., Singaucho, J., Guéguen, P., Bertrand, E., 2017. Seismic response of the basin of quito from continuous accelerometric records of RENAC-quito. In: 16th World Conf. on Earthquake Engineering.

Lay, T., Kanamori, H., Ammon, C.J., Koper, K.D., Hutko, A.R., Ye, L., Yue, H., Rushing, T.M., 2012. Depth-varying rupture properties of subduction zone megathrust faults. J. Geophys. Res.: Solid Earth 117 (B4).

León-Ríos, S., Ruiz, S., Maksymowicz, A., Leyton, F., Fuenzalida, A., Madariaga, R., 2016. Diversity of the 2014 Iquique's foreshocks and aftershocks: clues about the complex rupture process of a Mw 8.1 earthquake. J. Seismol. 20, 1059-1073.

Leyton, F., Leopold, A., Hurtado, G., Pastén, C., Ruiz, S., Montalva, G., Saéz, E., 2018a. Geophysical characterization of the Chilean seismological stations: First results. Seismol. Res. Lett. 89 (2A), 519-525.

Leyton, F., Pastén, C., Ruiz, S., Idini, B., Rojas, F., 2018b. Empirical site classification of CSN network using strong-motion records. Seismol. Res. Lett. 89 (2A), 512-518.

Lior, I., Rivet, D., Ampuero, J.P., Sladen, A., Barrientos, S., Sánchez-Olavarría, R., Opazo, G.A.V., Prado, J.A.B., 2022. Harnessing distributed acoustic sensing for earthquake early warning: Magnitude estimation and ground motion prediction.

Ma, B., Geersen, J., Lange, D., Klaeschen, D., Grevemeyer, I., Contreras-Reyes, E., Petersen, F., Riedel, M., Xia, Y., Tréhu, A.M., et al., 2022. Megathrust reflectivity reveals the uplip limit of the 2014 Iquique earthquake rupture. Nature Commun. 13 (1), 3969.

Madariaga, R., 1976. Dynamics of an expanding circular fault. Bull. Seismol. Soc. Am. 66 (3), 639-666.

Maksymowicz, A., Ruiz, J., Vera, E., Contreras-Reyes, E., Ruiz, S., Arraigada, C., Bonvalot, S., Bascuñan, S., 2018. Heterogeneous structure of the northern Chile marine forearc and its implications for megathrust earthquakes. Geophys. J. Int. 215 (2), 1080-1097.

Meng, L., Inbal, A., Ampuero, J.-P., 2011. A window into the complexity of the dynamic rupture of the 2011 Mw 9 Tohoku-Oki earthquake. Geophys. Res. Lett. 38 (7).

Montalva, G.A., Bastías, N., Leyton, F., 2022. Strong ground motion prediction model for PGV and spectral velocity for the Chilean subduction zone. Bull. Seismol. Soc. Am. 112 (1), 348-360.

Montalva, G.A., Bastías, N., Rodríguez-Marek, A., 2017. Ground-motion prediction equation for the Chilean subduction zone. Bull. Seismol. Soc. Am. 107 (2), 901-911.

Neighbors, C., Liao, E., Cochran, E.S., Funning, G., Chung, A., Lawrence, J., Christensen, C., Miller, M., Belmonte, A., Andres Sepulveda, H., 2015. Investigation of the high-frequency attenuation parameter, κ (kappa), from aftershocks of the 2010 M w 8.8 Maule, Chile earthquake. Geophys. J. Int. 200 (1), 200-215.

Pastén-Araya, F., Potin, B., Ruiz, S., Zerbst, L., Aden-Antoniów, F., Azúa, K., Rivera, E., Rietbrock, A., Salazar, P., Fuenzalida, A., 2021. Seismicity in the upper plate of the northern Chilean offshore forearc: Evidence of splay fault south of the Mejillones Peninsula. Tectonophysics 800, 228706.

Peyrat, S., Campos, J., de Chabalier, J.-B., Perez, A., Bonvalot, S., Bouin, M.-P., Legrand, D., Neressian, A., Charade, O., Patau, G., et al., 2006. Tarapacá intermediate-depth earthquake (Mw 7.7, 2005, northern Chile): A slab-pull event with horizontal fault plane constrained from seismologic and geodetic observations. Geophys. Res. Lett. 33 (22).

Piña-Validés, J., Socquet, A., Cotton, F., 2018a. Insights on the Japanese subduction megathrust properties from depth and lateral variability of observed ground motions. J. Geophys. Res.: Solid Earth 123 (10), 8937-8956.

Piña-Validés, J., Socquet, A., Cotton, F., Specht, S., 2018b. Spatiotemporal variations of ground motion in northern Chile before and after the 2014 M w 8.1 Iquique megathrust event. Bull. Seismol. Soc. Am. 108 (2), 801-814.

Pozo, I., Montalva, G., Miller, M., 2023. Assessment of kappa values in the Chilean subduction zone for interface and in-slab events. Seismol. Soc. Amer. 94 (1), 385-398.

Reginato, G., Vera, E., Contreras-Reyes, E., Tréhu, A.M., Maksymowicz, A., Bello-González, J.P., González, F., 2020. Seismic structure and tectonics of the continental wedge overlying the source region of the Iquique mw8.1 2014 earthquake. Tectonophysics 796, 228629.

Ruiz, S., Ammirati, J.-B., Leyton, F., Cabrera, L., Potin, B., Madariaga, R., 2019. The january 2019 (M w 6.7) coquimbo earthquake: insights from a seismic sequence within the Nazca plate. Seismol. Res. Lett. 90 (5), 1836-1843.

Ruiz, S., Madariaga, R., 2018. Historical and recent large megathrust earthquakes in Chile. Tectonophysics 733, 37-56.

Ruiz, S., Metois, M., Fuenzalida, A., Ruiz, J., Leyton, F., Grandin, R., Vigny, C., Madariaga, R., Campos, J., 2014. Intense foreshocks and a slow slip event preceded the 2014 Iquique M w 8.1 earthquake. Science 345 (6201), 1165-1169.

Schurr, B., Asch, G., Hainzl, S., Bedford, J., Hoechner, A., Palo, M., Wang, R., Moreno, M., Bartsch, M., Zhang, Y., et al., 2014. Gradual unlocking of plate boundary controlled initiation of the 2014 Iquique earthquake. Nature 512 (7514), 299-302.

Schurr, B., Rietbrock, A., Asch, G., Kind, R., Oncken, O., 2006. Evidence for lithospheric detachment in the central andes from local earthquake tomography. Tectonophysics 415 (1-4), 203-223.
Sladen, A., Rivet, D., Ampuero, J.P., De Barros, L., Hello, Y., Calbris, G., Lamare, P., 2019. Distributed sensing of earthquakes and ocean-solid earth interactions on seafloor telecom cables. Nat. Commun. 10 (1), 5777.

Twardzik, C., Duputel, Z., Jolivet, R., Klein, E., Rebischung, P., 2022. Bayesian inference on the initiation phase of the 2014 Iquique, Chile, earthquake. Earth Planet. Sci. Lett. 600, 117835.

Vallée, M., 2013. Source time function properties indicate a strain drop independent of earthquake depth and magnitude. Nat. Commun. 4 (1), 2606.

Wessel, P., Smith, W.H., 1998. New, improved version of generic mapping tools released. EOS Trans. Am. Geophys. Union 79 (47), 579.

Williams, E.F., Fernández-Ruiz, M.R., Magalhaes, R., Vanthillo, R., Zhan, Z., González-Herráez, M., Martins, H.F., 2019. Distributed sensing of microseisms and teleseisms with submarine dark fibers. Nat. Commun. 10 (1), 5778.

Ye, L., Lay, T., Kanamori, H., 2013. Ground shaking and seismic source spectra for large earthquakes around the megathrust fault offshore of northeastern Honshu, Japan. Bull. Seismol. Soc. Am. 103 (2B), 1221–1241.