Skip to main content

ORIGINAL RESEARCH

Earth Sci. Syst. Soc., 02 November 2023
This article is part of the Special Issue Hazards and Society View all 5 articles

Deformation of the Juan de Fuca Plate Beneath the Central Cascadia Continental Margin (44°-45°N) in Response to an Upper Plate Load

Anne M. Trhu
Anne M. Tréhu1*Kathy Davenport&#x;Kathy Davenport1Christopher B. Kenyon&#x;Christopher B. Kenyon1Suzanne M. CarbotteSuzanne M. Carbotte2John L. NabelekJohn L. Nabelek1Douglas R. ToomeyDouglas R. Toomey3William S. D. WilcockWilliam S. D. Wilcock4
  • 1College of Earth, Ocean and Atmospheric Sciences, Oregon State University, Corvallis, OR, United States
  • 2Lamont-Doherty Earth Observatory, Columbia University, Palisades, NY, United States
  • 3Department of Geosciences, University of Oregon, Eugene, OR, United States
  • 4School of Oceanography, University of Washington, Seattle, WA, United States

A 3D crustal model for the central Cascadia continental shelf and Coast Range between 44°N and 45°N shows that the crystalline crust of the forearc wedge beneath the coastline is characterized by a NW-trending, vertical slab of high-velocity rock interpreted to represent the dike complex that fed the Yachats Basalt, which was intruded into the forearc approximately 37 million years ago. A spatial correlation is observed between downward deflection of the crust of the subducting Juan de Fuca plate, inferred from inversion of PmP arrivals to image the Moho surface, and the high velocity (and consequently high density) anomaly underlying the Yachats Basalt. Apparent subsequent rebound of the subducting plate at greater depth suggests a primarily elastic response of the subducting plate to this load. Calculations for a range of plausible values for the magnitude of the load and the width and depth of the depression indicate that the effective elastic thickness of the subducted Juan de Fuca plate is < 6 km. Although our simple analytical models do not include partial support of the load of the slab by the adjacent upper plate crust or time dependence to account for the motion of the slab beneath the load, incorporation of those effects should decrease rather than increase the apparent strength of the subducted plate. We conclude that the subducted Juan de Fuca plate beneath the central Oregon margin is elastically thin and has the potential to store elastic strain energy before rupturing. Our model of a well-defined, focused and static upper plate load that locally deforms the subducted plate within the nominally seismogenic or transitional part of the Cascadia plate boundary may be unique in providing a relatively straightforward scenario for estimating the mechanical properties of the subducted Juan de Fuca plate. We extrapolate from these results to speculate that elastic deformation of the subducting plate may contribute to the low level of seismicity throughout much of the Cascadia forearc in the inter-seismic period between great earthquakes but note that our local results do not preclude faulting or elasto-plastic deformation of a thin and weak plate as it subducts. These results also suggest that the subducting plate should deform in response to larger scale variations in upper plate thickness and density.

Introduction

Many investigators have attempted to estimate the strength of the subducting plate based on the flexural response of the plate to subduction as measured from topography as the plate descends into the trench (e.g., Caldwell and Turcotte, 1979; Forsyth, 1980; Bry and White, 2007; Contreras-Reyes and Osses, 2010; Craig and Copley, 2014). An expected relationship between plate age and strength at the trench, however, has been elusive. Bry and White (2007) discerned a possible linear increase in apparent elastic thickness from 5–10 km for ages of 0–20 Ma from a global survey of subduction zones, whereas Contreras-Reyes and Osses (2010) argued for an elastic thickness of 13–15 km for a plate age of ∼5 Ma offshore Chile. Craig and Copley (2014) found no consistent relationship between plate age and strength and attributed this to bending beyond the elastic limit, which weakens the plate and results in outer rise seismicity (Forsyth, 1980).

Other studies have highlighted the impact of upper plate structure on subduction zone earthquakes (e.g., Song and Simons, 2003; Wells et al., 2003; Collot et al., 2004; Lotto et al., 2017; Salleres and Ranero, 2019; Egbert et al., 2022). Recently Arnulf et al. (2022) and Kimura et al. (2022) showed that the presence of a large mafic intrusion in the forearc of the Nankai Tough beneath the Kii peninsula deflects the subducting plate. Arnulf et al. (2022) argued that the high-density intrusion, which has an area of ∼8,000 km2, influences upper mantle hydration, earthquake nucleation, and plate boundary segmentation. Kimura et al. (2022) argue that it acts as a barrier to rupture propagation, with large plate boundary earthquakes nucleating near the edges of the intrusion and propagating away from it. We discuss a new three-dimensional (3D) model of the P-wave velocity (Vp) structure of the central Cascadia margin from 44° to 45°N and suggest that similar processes may be occurring on the central Cascadia margin.

The model was constructed from amphibious, controlled-source seismic data obtained through several piggyback programs over the past three decades and includes a high velocity (and by inference density) anomaly in the upper plate that deflects the lower plate. In contrast to the Kii peninsula region of the Nankai Trough, in which an ∼100 km wide pluton loads a 27 million-year-old plate that is subducting at a rate of ∼47 mm/yr, the central Cascadia upper plate load is a narrow (<10 km) vertical slab that loads a 10 million-year-old plate (Wilson, 2002) subducting at a rate of ∼31 mm/yr (Wilson and McCrory, 2022). We use the geometry of the lower plate in response to the load to estimate plate strength. The relict nature and compact geometry of the upper plate load and its direct correlation with a deflection in the Moho of the subducted plate provide a unique opportunity to derive an in situ estimate of the elastic strength of the subducted plate. We discuss the results in the context of models for the thermal structure of the plate (Cozzens and Spinelli, 2012), contemporary interplate locking (DeSanto et al., 2022) and the extent of sediment subduction derived from a 3D electrical resistivity model (Egbert et al., 2022). Because the frictional behavior of a fault and the amount of strain energy that can be stored in the adjacent rock depends on the balance between the strength of the fault surface and the strength of the adjacent rock (Cook, 1981), our model provides new insights into the geodynamic response of the Cascadia megathrust.

Background

The Cascadia Subduction zone is exceptional in its low level of instrumentally recorded seismicity, although the presence of an active volcanic arc and historical and paleoseismic evidence for large plate boundary earthquakes indicate that it remains an active subduction zone and that large earthquakes are likely in the future (Walton, Staisch, et al., 2021). Geologic and historic tsunami data indicate that the most recent subduction-zone-wide great earthquake occurred in 1700 (Atwater et al., 2005), and paleoseismic data indicate an average, but highly variable, recurrence interval of ∼500 years for such events. Figure 1A shows that the offshore forearc between ∼43° and 47°N has been nearly devoid of earthquakes with M >3 since at least 1989, except between 44° and 45°N (Figure 1A).

FIGURE 1
www.frontiersin.org

FIGURE 1. (A) Central Cascadia subduction zone showing earthquakes with magnitude greater than 3 in the ANSS Comprehensive Catalog from June 1989 through April 2023 (red dots). Violet lines offshore show locations of 3 strike-slip faults from Goldfinger et al. (1997): AC, Alvin Canyon fault; DB, Daisy Bank fault; WF, Wecoma fault. PO shows location of Pythia’s Oasis (Philip et al., 2023), an unusually warm and vigorous fluid seep. Black dashed lines offshore show the outline of basins from Wells et al. (2003): WB, Willapa Basin; AB, Astoria Basin; NB, Newport Basin bifurcated by a structural high known as the Stonewall Bank (SB), a NW-trending anticline with Miocene sediments exposed at the seafloor that has been active in the Holocene (Yeats et al., 1998); CB, Coos Bay Basin. White dashed line in (A) and (B) shows the seaward edge of the early Eocene Siletz terrane (Wells et al., 1998). Dashed lines onshore show late Eocene volcanics that penetrated through the underlying Siletz terrane (Wells et al., 2014): YB, Yachats Basalt; CH, Cascade Head; TV, Tillamook Volcanics; GH, Greys River Volcanics. White rectangle shows location of inset. Black rectangle shows location of (B). (B) Aeromagnetic map of the study region with ANSS-catalog earthquakes relocated by Tréhu et al. (2008, 2015, 2018) and Williams et al. (2011) in red. The depth distribution of these events, relocated relative to the M4.7 and M4.8 earthquakes, are shown in the depth section. Additional earthquakes located by Stone et al. (2018) and Morton et al. (2018) in dark and light blue, respectively. Depths of these earthquakes are scattered and not reliable (Tréhu et al., in prep). Regional moment tensor solutions are shown for the four largest events since 2004 (Tréhu et al., 2008, 2015, 2018). Black dashed line shows the extent of the magnetic anomaly indicative of the Yachats basalts. Grey dashed circles labeled “smt” show the location of subducted seamounts inferred from the magnetic anomalies (Tréhu et al., 2012). Cross-section shows the major structural boundaries interpreted from Gerdom et al. (2000) and relocated earthquakes. JdF, Juan de Fuca plate. Moho (solid grey lie) and the contact between the accretionary complex, Newport Basin and Siletz terrane (long grey dashed line) in the upper plate are well constrained by the Vp model; the plate boundary (short dashed grey line) is inferred based on the thickness of oceanic crust prior to subduction.

In 2004, two moderate earthquakes (M4.9, M4.7), occurred approximately 1 month apart and were followed by several aftershocks (Figure 1B). Based on regional moment tensor (RMT) analysis, high frequency pP observations and raytracing of secondary phases through a 2D velocity model of the source region (Gerdom et al., 2000), Tréhu et al. (2008) argued that these events were low-angle thrusts located on or near the plate boundary with slip vectors approximately aligned with the predicted plate motion vector (Figure 1B). Since 2004, smaller earthquakes have been reported by the PNSN (Pacific Northwest Seismic Network) at a rate of ∼6/year in clusters around these larger events (e.g., Tréhu et al., 2015; Tréhu et al., 2018), and inclusion of ocean bottom seismometer data from the Cascadia Initiative (Toomey et al., 2014), which provided offshore data from 2011 and 2015, has resulted in identification of many additional earthquakes in these clusters (Morton et al., 2018; 2023; Stone et al., 2018). Locations of these events using velocity models that approximate the structure beneath the continental shelf rather than the simpler models used by the PNSN suggest that hypocenters in the PNSN catalog are systematically overestimated, and relative relocation of clusters of earthquakes suggest that many occur in a narrow depth range near the plate boundary (Williams et al., 2011; Tréhu et al., 2015).

The relationship between the seismicity and forearc crustal structure is illustrated in Figure 1B based on potential field data calibrated by 2D seismic transects (Tréhu et al., 2012). The seaward edge of the Siletz terrane, an oceanic plateau that was accreted to North America at ∼50 Ma (see Wells et al., 2014 for a comprehensive review), is clearly defined by high-resolution aeromagnetic data (US Geological Survey, 2002), as are several subducted seamounts. A cluster of earthquakes is associated with a subducted seamount that appears to be impinging on the seaward edge of the Siletz terrane (Figure 1B). RMT solutions for 2 M >3.8 earthquakes in this cluster that occurred in 2012 and 2017, show thrust mechanisms that are steeper than those of the 2004 earthquake and may reflect the impact of a subducted seamount on the leading edge of Siletzia (Tréhu et al., 2015; Tréhu et al., 2018). A second cluster is located several kilometers landward of the seaward edge of the Siletz terrane and coincident with a local gravity low that reflects an anomalously deep part of the Newport basin (Tréhu et al., 2012). It was interpreted to indicate erosion of the Siletz terrane from below, resulting in the relatively recent formation of a basin within the larger and older Newport basin (Figure 1A), which was initiated in the late Miocene (McNeill et al., 2000).

The distinctive northwest-trending lineations in the aeromagnetic data (outlined and labeled YB in Figure 1B) coincide with coastal exposures of the Yachats Basalt (Figure 1A insert) and extend ∼17 km offshore (Figure 1B). The Yachats Basalt consists of seaward-dipping basalt flows about 37 m.y. old that were erupted from a west-northwest-oriented dike swarm intruding the accreted Siletz terrane onshore (Snavely and MacLeod, 1974; Snavely et al., 1976; Wells et al., 2014; Wells and Blakely, 2019). Wells et al. (2014) have argued that the Yachats Basalt, as well as several other middle-to-late Eocene volcanic units (Figure 1A), formed in a margin parallel extensional stress field and were subsequently rotated approximately 50° clockwise along with the entire Oregon block of the Siletz terrane. The linear topographic ridges and magnetic anomalies formed by the Yachats Basalt dikes are unique along the Oregon and Washington coastline (Figure 1A), and they project toward the anomaly offshore (Wells and Blakely, 2019).

The forearc between 44° and 45°N is also characterized by several other along-strike changes in geologic and geophysical parameters interpreted as potential indicators of interplate dynamics (Walton et al., 2021). These include: 1) a trio of strike-slip faults that offset the deformation front and can be traced across the outer wedge (Figure 1A; Goldfinger et al., 1997); 2) a lower degree of locking than north of 46°N or south of 43.5°N (Schmalzle et al., 2014); 3) a low slip region between two high slip patches in the 1700 earthquake at 44.4°N (Wang et al., 2013); 4) an offshore gravity high interpreted to be a potential rupture boundary (Wells et al., 2003); 5) segment boundaries at 44.2 and 45°N in the paleo-earthquake record derived from marine turbidites (Goldfinger et al., 2012); and 6) an abrupt southward increase in the slope of the outer wedge at 45°N (Watt and Brothers, 2020).

Data Acquisition

The Cascadia subduction zone between 44° and 45°N has been the focus of several controlled source seismic experiments in the past several decades (Figure 2). The primary dataset for this study was acquired as a piggyback project on the 2012 Ridge-to-Trench project (e.g., Han et al., 2016; Canales et al., 2017; Han et al., 2018). This included, in addition to the primary Ridge-to-Trench data, four coast-parallel lines of airgun shots with shorter connecting lines that were located to complement 2 coast-parallel shot lines acquired in 1996. The airgun sources from the R/V Marcus Langseth (cruise MGL1211) were recorded on seven additional ocean bottom seismometers (OBSs) and 35 temporary short period, continuously-recording seismometers in the Coast Range (Figure 2). No seismometers were deployed in the Willamette Valley because seismometers in the Willamette Valley installed in 1996 had background noise levels in the frequency range of interest that were ∼5 times higher than those in the Coast Range.

FIGURE 2
www.frontiersin.org

FIGURE 2. Map showing location of sources (lines) and receivers (inverted triangles) from the 1989 (orange), 1996 (yellow) and 2012 (red) data acquisition programs that were used for this analysis. The three sites and lines of sources for which data are shown in Figure 3 are highlighted and labeled. Source line T07 is also labeled because it is mentioned in the text. The white rectangle shows the limits of Figure 1B. The dashed blue line shows the location of CASIE21 line PS01.

Data from two prior experiments are also included in the dataset for the model. In 1989, as a piggy-back on an ODP seismic reflection site survey of the deformation front (MacKay, 1995), we acquired a single multichannel seismic reflection (MCS) transect from the Juan de Fuca plate, across the deformation front, to the coast near 44.8°N (Tréhu et al., 1995). The seismic sources for the MCS profile were recorded on a linear array of eight ocean bottom and ten onshore seismometers, which detected energy from diving waves and wide-angle reflections to offsets of up to 160 km (USGS Open-File reports 93-317, 93-318). The data were modeled using MacRay, a 2D forward modeling approach based on user-specified layers and lateral velocity variations with the layers (Tréhu et al., 1994). The result showed that thickened Siletz terrane rocks extended 30 km offshore, forming a nearly vertical backstop to the accretionary prism and restricting sediment subduction to at most a 2 km thick channel.

In 1996, a second amphibious, large-aperture transect was acquired near 44.6°N (Figure 1B). Airgun shots were recorded on 21 ocean bottom seismometers or hydrophones and on 29 onshore stations on a line that extended from the coast to the Cascade foothills (Gerdom et al., 2000). The 1996 experiment also included two north-south profiles on the continental margin (Gerdom et al., 2000), which were used to develop velocity/depth models for locating offshore earthquakes in this region (Williams et al., 2011; Morton et al., 2018; Stone et al., 2018). These margin-parallel shot lines were recorded as “fan” shots on the transect and on five additional seismometers that had been installed north and south of the primary transect to inform on 3D crustal structure (Figure 2). As in the 1989 data set, strong secondary arrivals interpreted to be reflections from the Moho of the subducted plate were observed, indicating that PmP observations could be used to map the Moho of the subducted plate in 3D with a suitable distribution of sources and receivers (as had been done in 1998 beneath the Olympic Peninsula; Tréhu et al., 2002; Preston et al., 2004), providing the impetus for the 2012 onshore experiment. A simplified version of the 2-dimensional Vp model derived from this experiment (Gerdom et al., 2000) served as a starting model for this study. A subset of the 1996 dataset was incorporated into the database for this study (Figure 2).

A much larger network of 755 vertical component nodal seismometers was installed in 2021 with an interstation spacing similar to that in the prior experiments but with a much larger footprint that extended from the California/Oregon border into southwest Washington (Cascadia 2021 Field Team, 2022). The nodal seismometers recorded shots from the 2021 CASIE21 multichannel seismic (Carbotte et al., 2021) and OBS (Canales et al., 2023) experiment and will enable extension of 3D velocity modeling to the margin from the Oregon/California border to the Olympic peninsula. Because only one margin-parallel line of shots was acquired along the forearc and many sites along this segment of the margin were reoccupations sites installed in 2012, the new large aperture data were not used in this study since they provide limited additional information on the structures of interest for this study. The CASIE21 reflection data along line PD01 (Figure 2), however, are important for validating our assumption that the Moho of the subducting plate is a proxy for the position of the plate boundary.

Data Processing and Interpretation

The integrated data set from 1989, 1996, and 2012 used to develop the velocity model presented here includes 57 onshore and 35 ocean bottom recorders. Station coordinates and elevations are given in Table 1 of Kenyon (2016). The continuous data recorded on the onshore and ocean bottom seismometers were merged with the shot instant data to generate SEGY-format data and record sections for each station/line pair. Record sections were examined with different bandpass filters and reduction velocities. Results from 2D transects were used to identify the offset range at which to expect high amplitude secondary arrivals representing reflections from the base of the subducted crust and potentially other distinctive secondary arrivals on fan shots, defined as source to receiver geometries in which the shots and receiver are not aligned. To illustrate the magnitude of this task, we note that each of the stations installed in 2012 resulted in 18 record sections from the vertical component alone; in total, ∼1,000 record sections were considered, with the majority of these showing signals from the offshore shots.

P-wave arrival times were picked from the vertical component of the onshore stations. Either the vertical geophone or hydrophone was used to pick arrivals on ocean bottom seismometers, depending on which component had the greater signal-to-noise ratio. Although some of the stations included horizontal geophones, no horizontal component data were used for this analysis. Pg, PmP, and Pn phases, defined as diving and refracted waves through the crust, wide-angle reflections from the base of the crust, and refracted wave in the uppermost mantle, were picked manually using tlPicker software developed at the University of Washington. Given the large amount of data, we were conservative and picked only arrivals for which the phase identification was unambiguous, resulting in 53,855 Pg, 17,922 PmP, and 13,143 Pn arrival time picks.

Three examples of record sections showing data from a single line of shots recorded at a single receiver are shown in Figure 3, with the location of the shot lines and stations shown in Figure 2. These examples are typical of the data and were chosen to illustrate aspects of the data that control the features of the model discussed in this paper. Figure 3A shows a rapid travel time advance in Pg between offsets of 70 and 80 km, suggesting that crustal P-waves from sources south of this transition are traveling through material with a much faster velocity than waves from sources to the north and that the crust contains a strong horizontal velocity gradient. This record section also displays strong PmP arrivals indicative of reflections from the base of the crust with angles of incidence close to the critical angle. The source-receiver distance range in which PmP is a strong arrival depends on the depth to Moho and the velocity structure and is a∼57–110 km in Figure 3A.

FIGURE 3
www.frontiersin.org

FIGURE 3. Examples of vertical component data recorded onshore that illustrate the impact of 3D crustal structure on the data. The data are displayed as equally spaced traces (an approximately linear function of distance along the shot line given) with the X-axis labeled by the distance between source and receiver, which is a non-linear function. Travel time picks for first arrivals interpreted to be Pg and secondary arrival interpreted to be PmP are overlain on the data. (A) Station S04, Line T03. (B) Station S16, Line T05. (C) Station 21, Line T05.

In Figure 3B, the shortest distance between the line of sources and the receiver is 49 km, near the middle of the record section. A clear difference in the travel time for sources from the north compared to shots from the south is observed, and the shortest time does not correspond to the shortest source/receiver distance. Water depths are nearly constant along most of the line (Figure 2) and cannot explain this asymmetry, which must result from along-strike variations in upper plate velocity. Strong secondary arrivals consistent with PmP are also observed for source/receiver offsets >50 km. No picks were made on the southern end of this section because of uncertainty about how to classify these arrivals.

Figure 3C shows arrivals from the same shots as in Figure 3B recorded on a seismometer located only ∼15 km farther south. As in Figure 3B, for the same source/receiver distance, ray paths from the north are much faster than ray paths from the south. PmP is observed along the entire line, and the time difference between Pg and PmP for a particular source/station distance varies with the azimuth as well as the distance between the source and receiver. Note also strong scattering of Pg on Figure 3B, recorded on Yachats Basalt, compared to Figure 3C.

Inversions were run using Stingray/Tomolab software developed at the University of Oregon (Toomey et al., 1994; Dunn et al., 2005). Crustal Pg arrivals were inverted with a variety of starting models, including a 1D velocity versus depth function fixed to the seafloor surface and an east-west oriented 2D velocity model representing a smoothed version of the model of Gerdom et al. (2000). Within the volume that is well sampled by the ray paths provided by our experimental geometry these two starting models yield similar results, with a slightly higher final misfit for the 1D model (see Figure 14 in Kenyon, 2016). Because of the overall eastward thickening of the crust, however, the 2D starting model results in fewer artifacts around the edges of the well sampled volume and produces a model that can more easily be used for studies for which we must extrapolate from the well-resolved volume (e.g., for earthquake relocations or for depth migration of multichannel seismic reflection data). Different smoothing and regularization approaches were also tested and compared as well as synthetic exercises to test the resolution provided by the experimental geometry. Examples of resolution tests are given in Supplementary Figure S1, with additional tests in Kenyon (2016).

Once the preferred crustal Vp model based on first arrivals had been obtained, PmP and Pn arrival times were included in the arrival time dataset for the inversion, which included a surface corresponding to the base of the subducted crust. Because Pn arrivals were primarily observed from shots located seaward of the deformation, they constrain the upper mantle velocity west of ∼124.7°W. Beneath the continental shelf and Coast Range, the position of the Moho surface is based almost entirely on PmP reflections. The sensitivity of the Moho geometry to different starting models and smoothing parameters is discussed in the next section.

Results

The inversion of travel times converged rapidly and smoothly with most of the decrease in the residuals obtained after 3 iterations and no significant improvement between iterations 5 and 6. The root mean squared misfit between the picked and predicted travel times for the initial 2D starting model was 1.1, 0.73, and 0.36 for Pg, PmP, and Pn, respectively, and decreased to 0.058, 0.072 and 0.04, respectively, after 6 iterations (close to within the estimated pick uncertainty of 0.04 s for Pg and Pn and 0.08 s for PmP). The residual for each iteration and the residuals for individual observations for the sixth iteration as a function of offset are shown in Supplementary Figure S2. Figure 4 shows horizontal and vertical slices spaced at 2 km intervals in depth and at 0.1° (∼8 km) intervals in longitude.

FIGURE 4
www.frontiersin.org

FIGURE 4. Selected equally spaced slices from the Vp model. Left side shows horizontal slices referenced to the seafloor. Right side shows vertical slices along lines of longitude. The thick grey line shows where the Moho surface intersects the slice. Bright green points highlighted by horizontal arrows on the N-S slice along −124.5° are top-of-crust and Moho picks from CASIE21 line PS01. Arrows on slices at 8 km depth and at −124.3° show the high-Vp, near-vertical, NW-trending slab discussed in the text.

The most striking and well-resolved feature of the Vp model is a northwest trending, near vertical, high Vp anomaly in the upper plate. This feature of the model is clearest in horizontal slices at depths of 6–10 km below the seafloor and is shown to extend through the crust in the north-south slices. It is highlighted by arrows in the horizontal slice at 8 km and in the vertical slice at −124.3°. Resolution tests (Supplementary Figures S1A, B) indicate that the shape and internal velocity of this body are well resolved by the experiment geometry, although we are not able to resolve whether its base is in contact with the crust of the subducting plate (possibly even penetrating the upper crust of the subducted plate) or whether a thin, but unresolved, low-velocity zone separates the two. Another well-resolved large anomaly is an oval, low-velocity region located north of the northwest trending high velocity anomaly that extends from the surface to 8 km depth.

The primary feature of the Moho surface model (the thick grey line on the NS slices in Figure 4 and the surface shown in Figure 5A) is an apparent downward deflection of the Moho beneath the high-density crustal slab. Because of the experimental geometry, resolution is best east of −124.7°, which represents the approximate position of reflection points for PmP arrivals from line T7 into stations near the coast. Because of the slab dip, reflection points for PmP are shifted towards the source rather than being midway between source-receiver pairs.

FIGURE 5
www.frontiersin.org

FIGURE 5. (A) Result of inversion of PmP arrivals for the Moho surface corresponding to the model of Figure 4 (model V2). Red line is the coast. Black star shows the apex of the deflection. (B) Cross-section through three 2D starting models used to test the resolution of the Moho surface. (C) Profiles of the Moho surface after six iterations along longitude −124.3°. The result for starting model V2 is shown with two different values of smoothing. The part of the surface sampled by PmP observations is shown as solid lines and the unsampled region is shown by dotted lines. The depth of the starting surface at this longitude is shown by the horizontal dashed lines. Our qualitative estimate of the range of possible values for the amplitude and width of the deflection are also indicated. (D) Profiles of the Moho surface after 6 iterations along latitude 44.35. Starting models V1 and V3 are also shown. (E) Moho surface obtained after six iterations for starting models V3 (i) and V2 with additional smoothing (ii). Red line is the coast; dashed black line indicates the region within which the difference after six iterations for all starting models tested are the same to within 1 km (see Supplementary Figure S4 for difference plots). Discontinuous offsets smaller than ∼1 km may be present but are not resolved by the model.

To test whether this deflection is required by the data, we inverted for the Moho surface for three different starting models. V2 is the starting model that resulted in the surface shown in Figures 4, 5A. V3 was generated by projecting the maximum deflection to the north and south to generate a linear trough rather than a “dimple.” V4 is parallel to V2 but 3 km deeper. Results for the three initial models, including two different smoothing functions for V2, are shown for a north-to-south cross-section at −124.29° in Figure 5C and for a west-to-east cross-section at 44.35° in Figure 5D. The resulting surfaces for V2 with greater smoothing and for V3 are shown in Figure 5E. We conclude that the initial deflection and downdip flattening of the Moho beneath the high velocity slab are robust features of the model that do not depend on the starting model.

Having confirmed that the data require a depression in the Moho beneath the high velocity/density anomaly, we consider how well we can estimate its apex, amplitude and width. The apex of the deflection (shown by a star on Figure 5A) is near 44.35°N, 124.6°W. Estimates for the amplitude and width are correlated, with a width of ∼70 km and amplitude of ∼2 km for the solution for starting models V3 and V4 and a width of ∼120 km and amplitude of ∼4 km for starting model V2 (Figure 5C). Although estimating the width and amplitude of the deflection in the down-dip direction is complicated by an overall dip of ∼14° to the east and limited spatial coverage in the down-dip direction, we estimate a width of ∼50 km and amplitude of 2–3 km, similar in magnitude to the estimate in the strike direction.

Discussion

Geologic Interpretation of the Seismic Model

Figure 6A shows the velocity model at 8 km depth overlain by selected other geological and geophysical features of the continental margin in this region. Figure 6B shows a cross-section through the model that samples the near-vertical, high-velocity anomaly and shows that the contrast in velocity between this feature and the adjacent rocks of the Siletz terrane extends through the upper plate wedge (Figure 6C). The average velocity above the plate boundary (defined as a surface parallel to and 6 km above the Moho surface of Figure 5A for reasons discussed in this section) is shown in Figure 7, highlighting the magnitude and focused nature of the velocity contrast.

FIGURE 6
www.frontiersin.org

FIGURE 6. (A) Velocity slice at 8 km depth with geological and geophysical features discussed in the text overlain. The locations of the vertical slices in (B) and (C) are shown. Thick grey lines labeled in km are contours of depth to the plate boundary for model V2 of in Figure 5 assuming the plate boundary is 6 km above the Moho. Bars above the map show the transition zone from 90 or 100% interplate locking to the zone of episodic tremor and slip (ETS) beneath the Coast Range (from DeSanto et al., 2022). ACF, Alvin Canyon fault; DBF, Daisy Bank Fault; PO, Pythia’s Oasis. (B) Vertical slice through the high Vp slab. Arrows along the top of the cross-section show the seaward edge of Siletz based on magnetic anomaly data and the coastline, respectively. Thin dashed grey line is the interpreted top of Siletzia and a subducted seamount located near the line of this section and associated with seismicity. The outline of the high velocity slab and the Moho and top of crust for model V2 are shown by thick grey lines and labeled. (C) Vertical slice ∼30 km north of (B). Annotations the same as for (B). (D) Temperature 3 km beneath the top of the subducted crust from the model of Cozzens and Spinelli (2012). The red line shows results for a purely conductive model; the blue line shows the high flow (coolest) model of Cozzens and Spinelli (2012).

FIGURE 7
www.frontiersin.org

FIGURE 7. Average velocity above the plate boundary, defined as a surface 6 km above the Moho surface of Figure 5A. Star is the apex of the deflection. Bold dashed line is the region within which the Moho surface inversion is independent of the starting model. Light dashed line is the line along which estimates for the deflection amplitude and width were made. Red line is the coast. Contours are at 0.5 km/s intervals.

The location and orientation with respect to the Yachats Basalt and its offshore extension, as indicated by the outline of their magnetic signature (Figure 6A), suggest that the high-velocity, near vertical anomaly in the upper plate represents the dike complex that fed the near-surface basalt flows, consistent with geologic observations (Snavely and MacLeod, 1974). Geochronologic and paleomagnetic studies indicate that the flows were emplaced within and through the Siletz terrane between 36 and 37.5 million years ago in an extensional forearc environment and have rotated clockwise with Siletz Oregon Block since then (Wells et al., 2014). Since higher Vp implies higher density, the dike complex represents a static load embedded within the upper plate.

The low velocity zone in the upper 8 km north of the dike complex is interpreted to be a basin within the larger Newport Basin and corresponds to topography on a mid-Miocene unconformity mapped throughout the region by McNeill et al. (2000). The model confirms that the earthquakes in the northern cluster are located on or near the plate boundary consistent with prior suggestions that the basin is the results of deep subduction of a seamount beneath the Siletz terrane and erosion of the base of the overlying plate (Wells et al., 2003; Tréhu et al., 2012).

The seaward edge of the Siletz terrane (SES) appears as a north-south ridge of relatively high velocity near −124.5° that clearly bounds the seaward edge of the basin in the slices at 6 and 8 km depth (Figure 4). The cross-sections in Figures 6B, C show an interpretation of the upper and western edges of the Siletz terrane based on the 3D seismic velocity model. This interpretation takes into account the depth-dependent increase in velocity within the terrane (e.g., Tréhu et al., 1994). In Figure 6B, the seaward edge appears to be nearly vertical. The apparent seamount inferred from aeromagnetic anomaly data and associated with the southern cluster of seismicity is not well resolved in this model, although the model suggests that the seismicity is located primarily within the upper plate rather than on the plate boundary, consistent with the conclusion of Tréhu et al. (2012); Tréhu et al. (2015) and Morton et al. (2018) that they result from the interaction between a subducted seamount and the SES. In contrast to the near-vertical SES in Figure 6B, the small velocity inversion along that boundary in Figure 6C suggests the presence of a poorly resolved channel of underthrust sediment (see Supplementary Figure S1).

Figures 6B, C also show the position of the Moho surface from Figure 5A. Based on the Juan de Fuca plate crustal thickness from seismic refraction data from the incoming plate (Tréhu et al., 1994; Gerdom et al., 2000; Canales et al., 2017) and on the multichannel seismic reflection data along line PS01 acquired during cruise MGL2104 (Carbotte et al., 2021) we interpret the plate boundary to be surface located 6 km above the Moho. CASIE21 line PS01 (location shown on Figure 2) shows a strong and nearly continuous reflection, interpreted to be the top of the subducted crust, that overlies an intermittent reflection that is coincident with the Moho surface in the 3D model (Figure 4 on the slice along −124.5°). The agreement between the CASIE21 seismic reflection data and the 3D Vp model validates this assumption. The CASIE21 observations also justify an assumption that the top of the subducted crust and plate boundary in this region are approximately coincident because there is no shallower strong and continuous reflection from the lower crust between 44° and 45°N that would indicate the presence of a layer of subducted sediment between the plate boundary and the top of the subducted crust thick enough to be resolved by our data. The average dip of the subducting plate increases beneath the SES and may flatten toward the coast, and the local deflection beneath the dike complex is superimposed on this regional trend.

We consider it unlikely that the local deflection is a pre-existing feature of the subducting slab that is coincidentally located beneath the crustal slab because we would expect the topography of the overlying plate to sag above such a depression if the depression had not been filled prior to subduction; instead, the Yachats Basalt is associated with topographic highs (Figure 1). On the other hand, if the depression had been filled prior to subduction, we would expect to see a plate boundary reflection above the reflection interpreted to be the top of the subduction crust–a scenario that is not compatible with the CASIE21 reflection data.

Deformation Modeling of the Juan de Fuca Plate

Having ruled out a preexisting depression in the subducted plate, we consider two possible endmember explanations for the deflection: elastic or plastic deformation. If the deformation is purely elastic, the depression should disappear once the slab has subducted beneath the region where the slab is present. If it is plastic, a trough parallel to the convergence vector should remain in the slab at greater depth. The tests of the sensitivity of the Moho surface to the starting model in Figure 5 indicate that the initial depression flattens at ∼30 km depth, consistent with an elastic rebound of a loaded plate. We recognize, however, that this is close to the edge of the resolved model space and cannot resolve whether the rebound is partial (elasto-plastic) or complete (purely elastic).

Assuming that the observed deformation is purely elastic, we can estimate the flexural rigidity and elastic plate thickness of the subducted Juan de Fuca plate by estimating the load resulting from the dike complex and modeling the observed amplitude and wavelength of deflection. We model the deflection as a point load (Brotchie and Silvester, 1969) on an unbroken plate and consider two options for estimating the load. These are depicted in Figure 8A. In the first, the point source is estimated by collapsing the excess load of a column that is 20 km high with a 10 × 10 km cross-section. In the second, the entire mass excess of the dike complex (approximated by a 20 km high prism with a 10 × 40 km cross-section) is collapsed into the point load. This second scenario attempts to qualitatively capture the dynamic effect of a dipping slab moving toward the load. As a third option, we calculate the impact of a line load perpendicular to the deflection (Turcotte and Schubert, 1982).

FIGURE 8
www.frontiersin.org

FIGURE 8. (A) Schematic scenario for flexural modeling. The load was calculated for three scenarios: 1) a point source representing the load of a 10 × 10 × 20 km column; a 10 × 40 × 20 km rectangular slab; and a line source based on a 10 × 20 km cross-section. (B) Average P-wave velocity and corresponding density in the upper 20 km of the model along a north-south transect at −124.29° (light blue line) and a simplified model of a density contrast between the high-velocity dike complex and the adjacent Siletz terrane. (C) Overview of the equations and procedure for calculating the elastic thickness of the subducting plate given a load and the width and amplitude of the depression resulting from the load. Calculations in Table 1 assume E = 80 GPa, n = 0.25 and r of 130, 160 and 190 kg/m3. See Figure 5C for estimates for the amplitude and width of the deflection.

To estimate the density contrast between the slab and the adjacent rock, we convert the average velocity above 20 km depth along a north-south slice of the Vp model at −124.29° and convert Vp to density assuming the Nafe-Drake curve (Ludwig et al., 1970; Figure 8B). This exercise suggests that the average density contrast between the dike complex and the surrounding rock is ∼160 kg/m3.

The equations used for calculating the elastic thickness are shown in Figure 8C. Observables are the width and depth of the observed depression along a north-south profile at longitude −124.29, as shown in Figure 5C, avoiding the eastward increase in dip of the slab. Calculations for the flexural rigidity and elastic plate thickness for the three source scenarios for assumed density contrasts of 130, 160 and 190 kg/m3 and different possible combinations of the width and amplitude of the deflection are given in Table 1 assuming a Young’s modulus of 80 GPa and a Poisson’s ratio of 0.25. Although the resolution tests indicate that observations of the minimum deflection and width are not independent, we included an extreme case of minimum deflection and maximum width to illustrate the sensitivity of the calculation to the observations.

TABLE 1
www.frontiersin.org

TABLE 1. Summary of results of flexural modeling.

All combinations of parameters yield very low elastic plate thicknesses, ranging from 1.3 to 6.4 km, with a thickness >5 km only achieved with the extreme endmember line-load model of 2 km of deflection over a width of 120 km. We recognize that these simple static elastic calculations are an oversimplification of reality, in which the load is partially supported by the adjacent crust of the upper plate, the motion of the slab beneath the load is time dependent, and viscosity of material both below and above the elastic plate may affect the rates of deflection and rebound. Including these effects as well as the superimposed effect of subduction to the east beneath a thickening forearc wedge of increasing and heterogeneous density is beyond the scope of this study.

These estimates are on the low end of the range of observations reported for the flexural strength of young oceanic lithosphere based on seamount loading of oceanic lithosphere (Watts, 2001). In many cases in nature, estimates of the flexural strength of oceanic lithosphere are complicated by time dependence of the load generation and/or plastic deformation. We have argued that the load in the scenario we model here results from an old magmatic event and therefore is independent of time on a scale relevant to this calculation (although off course, it is susceptible to erosion on a longer time scale) and that the apparent rebound of the deflection argues for a primarily elastic response. Moreover, partial support of the load by the upper plate or delayed deflection because of mantle viscosity would imply an even weaker plate given the observed amplitude and width of the deflection.

The lack of seismicity associated with the lines of maximum curvature of the deflection is consistent with our interpretation of dominantly elastic deformation, although, given the limited length of the seismic record and the ∼1 km resolution of the smoothed Moho surface obtained from the inversion, we cannot rule out the presence of small fault offsets. Our conclusion of dominantly elastic deformation in response to the focused load of the dike complex also does not preclude the presence of larger amplitude elastic and plastic deformation of the subducted plate elsewhere in Cascadia as it encounters large scale variations in the upper plate load.

Comparison to Geodetic and Thermal Models

We now evaluate whether elastic deformation is compatible with geodetic and thermal models. Although geodetic models based only on onshore GPS measurements cannot resolve whether interplate locking is restricted to a narrow band of shore or distributed more widely across a zone of partial locking (e.g., Schmalzle et al., 2014), and increasing database of offshore geodetic measurements supports the former model (DeSanto et al., 2022). Bars above the map in Figure 6A illustrate two models that fit the combined onshore/offshore data from this segment of the Cascadia margin (DeSanto et al., 2022; D. Schmidt, personal communication, 2023). Locking degree is shown as a function of distance from the deformation with black bars representing strong locking and the variable grey shade representing an exponential decrease to 0% locking 165 km east of the deformation front. Along this segment of the margin episodic tremor and slip (ETS), which is thought to result from slip on the nearly unlocked portion of the plate boundary (e.g., Walton, Staisch, et al., 2021) is observed in a band ∼25–65 km from the coast (e.g., Boyarko et al., 2015).

Figure 6D shows predicted temperatures 3 km beneath the plate boundary (i.e., at the midline of the subducted oceanic crust) from the model of Cozzens and Spinelli (2012) for a transect across the Cascadia forearc at 44°N. Results are shown for their endmember models of purely conductive heat transport and compared to crust that has been cooled by vigorous fluid circulation. Although Cozzens and Spinelli (2012) were not able to resolve the impact of fluid circulation on the temperature of the plate boundary beneath the central Cascadia margin from heat flow measurements available at the time, recent heat flow measurements over two buried seamounts located ∼25 km west of the deformation front indicate that vigorous hydrothermal circulation in the upper crust of the Juan de Fuca plate persists to the deformation front (Norvell et al., 2023), suggesting that the in situ temperature is between these two endmembers. The temperature in the subducted crust where the deflection is observed is, therefore, likely well below the brittle/ductile transition for basalt (Violay et al., 2012).

We reconcile the apparently elastic behavior of the subducting plate in this region with the GPS-based estimates of interplate locking, which indicate that the plate boundary where we image the depression in the subducting plate is downdip of the zone of strong locking, by calling on different brittle/ductile transitions for the plate boundary and for the subducted crust. A temperature of 350°C is widely cited as marking the transition between locking and transitional behavior of the plate boundary (Hyndman and Wang, 1993), with 450°C representing the boundary between transitional and freely slipping interplate behavior. These temperatures are thought to represent the behavior of a silica-rich sediment and gouge-filled plate boundary shear zone. In contrast, Violay et al. (2012) have argued, based on laboratory measurements extrapolated to in situ conditions, that the brittle/ductile transition for basalt is ∼550°C. The geodetic data suggest that the Cozzens and Spinelli (2012) endmember of rapid fluid circulation is too cold to be consistent with the narrow band of strong locking indicated by the GPS data and that the conductive model is too hot. Stronger constraints on the temperature structure of the incoming plate are needed to better constrain the amount of cooling provided by fluid circulation (Norvell et al., 2023). Regardless of the amount of cooling due to circulation, however, the subducted crust where the deflection is observed is likely cold enough to be dominated by elastic and brittle deformation.

Potential Implications for Seismic Hazard Evaluation

The amount of elastic energy that can be stored in the rocks adjacent to a locked fault, to be released when the locked fault ruptures in an earthquake, depends on the elastic strength of those rocks as well as on the frictional properties of the fault surface. An elastically weak rock can store more elastic energy for a given fault strength than a stiff rock (Cook, 1981). While Salleres and Ranero (2019) recently explored this phenomenon using Vp at the base of the upper plate as a proxy for the elastic modulus of the upper plate wedge, they discounted the impact of the lower plate by concluding that it does not affect potential down-dip variations in strain accumulation across the megathrust.

Our results indicate that, at least in Cascadia, the elastic properties of the lower plate should not be ignored when modeling subduction zone dynamics. The elastic strain implied by the imaged deflection is added to the strain resulting from plate bending due to subduction. We speculate that the transient deformation responsible for this added lower plate strain might play in role in triggering slip on the megathrust. As elastic strain energy accumulates over time in both the upper and lower plates due to locking on the plate boundary, the additional strain in the lower plate due to the ongoing deflection of the plate beneath the upper load may indirectly trigger slip on the plate boundary, and we pose this scenario as one worth exploring with more sophisticated strain accumulation and earthquake triggering models. As Arnulf et al. (2022) and Kimura et al. (2022) noted, the Nankai earthquakes of 1944 and 1946 initiated on either side of large upper plate pluton and propagated away from the pluton. This speculation is also consistent with the observation that great subduction zone earthquakes around the globe more often than not initiate at basement highs on the edge of forearc basins and propagate beneath the basins. (Song and Simons, 2003; Wells et al., 2003; Llenos and McGuire, 2007). The Wang et al. (2013) model for patchy slip during the 1700 Cascadia earthquake may also be consistent with this scenario because two patches of high slip flank the zone of lower plate deformation that we image. By analogy with Nankai, we speculate that rupture may have initiated here and propagated bilaterally in 1700.

We also speculate that our local result can provide an insight into the mechanical properties of the rest of the Cascadia subduction zone. The plate age in our study area is similar to the age of the subducted crust beneath the continental shelf along much of the Cascadia margin (Wilson, 2002). The lack of similar examples of focused lower plate flexure may be due in part to limited data coverage as well as to the fact that other sources of heterogeneous forearc loading have a more complex temporal history or are located farther from the deformation front, further complicating modeling to derive an estimate of effective elastic strength. Ultimately the 2021 onshore nodal (Cascadia 2021 Field Team, 2022), OBS data (Canales et al., 2023), integrated with the CASIE21 seismic reflection images (Carbotte et al., in review), will lead to extension of our 3D model up-dip (Nolan, 2023) and to the north (Ashraf et al., 2022) and south (work in progress).

Our results may also contribute insights into the interpretation of a 3D model for the electrical conductivity structure of Cascadia (Egbert et al., 2022). In Figure 6A, we show the outline of the low conductivity “g” anomaly, which represents a discontinuity in a high conductivity region immediately above the subducting plate that is interpreted to represent sediment subducted beneath the edge of the Siletz terrane (Egbert et al., 2022). We speculate that the load represented by the dike complex has pinched off the subducted sediment channel beneath the margin here. Future analyses can test this hypothesis by increasing the volume imaged seismically to increase overlap with the conductivity model and, potentially, through a joint inversion of seismic and electrical conductivity that uses the different sensitivity of these two imaging approaches to better constrain the spatial resolution of both.

Geodynamic modeling of plate interactions that incorporate the geometric relationship between the upper and lower plates imaged by this study as well as its rheological implications may also provide new insights into the driving forces that have caused brittle deformation of the upper and lower plates farther up-dip in the region characterized by the Alvin Canyon, Daisy Bank and Wecoma faults (Figure 1A; Goldfinger et al., 1997) and on possible connections between these faults, the deeply sourced seafloor seep known as Pythia’s Oasis (Philip et al., 2023), Holocene upper plate deformation that formed Stonewall Bank (Yeats et al., 1998), and clusters of recent seismicity (Morton et al., 2023).

Conclusion

A 3D Vp model of the central Cascadia margin between 44 and 45°N that includes a model of the Moho surface based on wide-angle reflections reveals that the forearc wedge contains a northwest-trending, near-vertical high-velocity anomaly that extends through the upper plate wedge. The anomaly is interpreted to be the dike complex that fed the mid-to-late Eocene Yachats Basalt, which is exposed in the Coast Range. This dike complex acts as a local, static load that deflections the plate boundary. The imaged scenario provides the opportunity to estimate the effective elastic thickness of the subducted Juan de Fuca plate beneath the margin. We conclude that:

• The experiment geometry, in which six parallel lines of offshore sources spaced 6–12 km apart and recorded by a grid of stations onshore spaced ∼8 km apart, complemented by ocean bottom seismometers, provides good resolution of the crustal velocity structure and the geometry of the Moho surface beneath the continental shelf and coastline.

• The velocity anomaly is ∼10 km wide, much narrower than the region of extrusive basalt that was mapped using field samples and aeromagnetic data.

• The Moho is depressed beneath the dike complex. The amplitude and width of the deflection, measured perpendicular to the curvature of the subduction zone, are 2–4 km and 70–120 km, respectively. A down-dip profile suggests at least partial elastic rebound.

• Multichannel seismic reflection data from the 2021 CASIE21 experiment as well as older data from the oceanic plate indicate that the plate boundary is located ∼6 km above Moho in this region.

• Assuming the Nafe-Drake relationship between velocity and density, the contrast between the dike complex and the adjacent mafic rocks of the Siletz terrane is ∼160 kg/m3.

• An apparent elastic thickness of 2–6 km is obtained from simple analytic elastic plate models for three different assumed load configurations, density contrasts of 130–190 kg/m3 and a range of estimates for the amplitude and width of the deflection, implying that the subducted Juan de Fuca plate is very weak.

• Although we assume pure elastic behavior, we cannot resolve whether the observed deformation is purely elastic or elasto-plastic (i.e., only partially recoverable). Elastic behavior is compatible with thermal models and observed seismicity but small offset faults may be present but are not resolved by our analysis.

• Elastic energy stored in the lower plate should be considered in models of interplate dynamics.

Data Availability Statement

The datasets presented in this study can be found in online repositories. The names of the repository/repositories and accession number(s) can be found below: https://ds.iris.edu/ds/nodes/dmc/data/types/time-series-data/.

Author Contributions

AT: conceived of study, obtained funding, organized field work, supervised analysis, drafted manuscript. KD: modeled Moho surface using PmP arrivals, generated figures, discussed interpretation, edited manuscript. CK: modeled Pg arrivals, generated figures, discussed interpretation, edited manuscript. SC: provided source for 2012 study, discussed interpretation, edited manuscript. JN: assisted with field work, discussed interpretation, edited manuscript. DT: developed inversion software, discussed interpretation, edited manuscript. WW: developed analysis software, discussed interpretation, edited manuscript. All authors contributed to the article and approved the submitted version.

Funding

This work was supported by NSF grant EAR1147975 and USGS grant G17AP00046 to Oregon State University.

Conflict of Interest

The authors declare that the research was conducted in the absence of any commercial or financial relationships that could be construed as a potential conflict of interest.

Publisher’s Note

All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article, or claim that may be made by its manufacturer, is not guaranteed or endorsed by the publisher.

Acknowledgments

We thank the many individuals who helped to acquire the field data in 1989, 1996 and 2012 and the public land managers and private landowners who provided access to their land. We thank Kelin Wang for interesting discussions. Ray Wells and Nick Harmon provided thorough and constructive reviews.

Supplementary Material

The Supplementary Material for this article can be found online at: https://www.escubed.org/articles/10.3389/esss.2023.10085/full#supplementary-material

References

Arnulf, A. F., Bassett, D., Harding, A. J., Kodaira, S., Nakanishi, A., and Moore, G. (2022). Upper-Plate Controls on Subduction Zone Geometry, Hydration and Earthquake Behaviour. Nat. Geosci. 15 (2), 143–148. doi:10.1038/s41561-021-00879-x

CrossRef Full Text | Google Scholar

Ashraf, A., Hooft, E., Tréhu, A., Nolan, S., Wirth, E., Ward, K., et al. (2022). A 3-D Seismic Velocity Model Across the South-Central Cascadia Subduction Margin (abs). AGU fall 2022 meeting.

Google Scholar

Atwater, B. F., Musumi-Rokkaku, S., Satake, K., Tsuji, Y., Ueda, K., and Yamaguchi, D. K. (2005). The Orphan Tsunami of 1700 – Japanese Clues to a Parent Earthquake in North America. USGS Prof. Pap. 1707. doi:10.3133/pp1707

CrossRef Full Text | Google Scholar

Boyarko, D. C., Brudzinski, M. R., Porritt, R. W., Allen, R. M., and Tréhu, A. M. (2015). Automated Detection and Location of Tectonic Tremor along the Entire Cascadia Margin From 2005 to 2011. Earth Plan. Sci. Lett. 430, 160–170. doi:10.1016/j.epsl.2015.06.026

CrossRef Full Text | Google Scholar

Brotchie, J. S., and Sylvester, R. (1969). On Crustal Flexure. J. Geophys. Res. 74, 5240–5252. doi:10.1029/JB074i022p05240

CrossRef Full Text | Google Scholar

Bry, M., and White, N. (2007). Reappraising Elastic Thickness Variation at Oceanic Trenches. J. Geophys. Res. 112. doi:10.1029/2005JB004190

CrossRef Full Text | Google Scholar

Caldwell, G., and Turcotte, D. L. (1979). Dependence of the Thickness of the Elastic Oceanic Lithosphere on Age. J. Geophys. Res. 84, 7572–7576. doi:10.1029/jb084ib13p07572

CrossRef Full Text | Google Scholar

Canales, J. P., Carbotte, S. M., Nedimovic, M., and Carton, H. (2017). Dry Juan de Fuca Slab Revealed by Quantification of Water Entering Cascadia Subduction Zone. Nat. Geosci. 10, 864–870. doi:10.1038/NGEO3050

CrossRef Full Text | Google Scholar

Canales, J. P., Miller, N. C., Baldwin, W., Carbotte, S. M., Han, S., Boston, B., et al. (2023). CASIE21-OBS: An Open-Access, OBS Controlled-Source Seismic Dataset for Investigating the Structure and Properties of the Cascadia Accretionary Wedge and the Downgoing Explorer-Juan de Fuca-Gorda Plate System. Seismol. Res. Lett. 94, 2093–2109. doi:10.1785/0220230010

CrossRef Full Text | Google Scholar

Carbotte, S. M., Han, S., and Boston, B. (2021). Cascadia Seismic Imaging Experiment – CASIE21, R/V Marcus Langseth MGL2104, Cruise Report.

Google Scholar

Cascadia2021 Field Team (2022). Field Report for the Cascadia2021 Seismic Node Experiment. Available at: https://blogs.oregonstate.edu/cascadia2021/2022/06/08/field-data-report/ (Accessed March 19, 2022).

Google Scholar

Collot, J.-Y., Marcaillou, B., Sage, F., Michaud, F., Agudelo, W., Charvis, P., et al. (2004). Are Rupture Zone Limits of Great Subduction Earthquakes Controlled by Upper Plate Structures? Evidence From Multichannel Seismic Reflection Data Acquired Across the Northern Ecuador-Southwest Columbian Margin. J. Geophys. Res. 109. doi:10.1029/2004JB003060

CrossRef Full Text | Google Scholar

Contreras-Reyes, E., and Osses, A. (2010). Lithospheric Flexure Modeling Seaward of the Chile Trench: Implications for Oceanic Plate Weakening in the Trench Outer Rise Region. Geophys. J. Int. 182, 97–112. doi:10.1111/j.1365-246X.2010.04629.x

CrossRef Full Text | Google Scholar

Cook, N. G. W. (1981). “Stiff Testing Machines, Stick Slip Sliding, and the Stability of Rock Deformation,” in The Mechanical Behavior of Crustal Rocks, Geophysical Monograph (Am. Geophys. Un), 24.

Google Scholar

Cozzens, B. D., and Spinelli, G. A. (2012). A Wider Seismogenic Zone at Cascadia Due to Fluid Circulation in Subducting Oceanic Crust. Geology 40 (10), 899–902. doi:10.1130/G33019.1

CrossRef Full Text | Google Scholar

Craig, T. J., and Copley, A. (2014). An Explanation for the Age Independence of Oceanic Elastic Thickness Estimates From Flexural Profiles at Subduction Zones, and Implications for Continental Rheology. Earth Plan. Sci. Lett. 492, 207–216. doi:10.1016/j.epsl.2014.02.027

CrossRef Full Text | Google Scholar

DeSanto, J. B., Schmidt, D. A., Chadwell, C. D., Zumberge, M. A., and Sasagawa, G. S. (2022). “Horizontal Deformation Rates Near the Cascadia Subduction Zone Trench Revealed by Offshore GNSS-Acoustic Time Series,” in Fall Meeting 2022 (AGU).

Google Scholar

Dunn, R. A., Leki´cDetrick, R. S., and Toomey, D. R. (2005). Three-Dimensional Seismic Structure of the Mid-Atlantic Ridge (35°N): Evidence for Focused Melt Supply and Lower Crustal Dike Injection. J. Geophys. Res. 119. doi:10.1029/2004JB003473

CrossRef Full Text | Google Scholar

Egbert, G., Yang, B., Bedrosian, P. A., Key, K., Livelybrooks, D. W., Schultz, A., et al. (2022). Fluid Transport and Storage in the Cascadia Forearc Influenced by Overriding Plate Lithology. Nat. Geosci. 15, 677–682. doi:10.1038/s41561-022-00981-8

CrossRef Full Text | Google Scholar

Forsyth, D. W. (1980). Comparison of Mechanical Models of the Oceanic Lithosphere. J. Geophys. Res. 85, 6364–6368. doi:10.1029/jb085ib11p06364

CrossRef Full Text | Google Scholar

Gerdom, M., Tréhu, A. M., Flueh, E. R., and Klaeschen, D. (2000). The Continental Margin off Oregon From Seismic Investigations. Tectonophysics 329, 79–97. doi:10.1016/s0040-1951(00)00190-6

CrossRef Full Text | Google Scholar

Goldfinger, C., Kulm, L. D., Yeats, R. S., McNeill, L., and Hummon, C. (1997). Oblique Strike-Slip Faulting of the Central Cascadia Submarine Forearc. J. Geophys. Res. 102, 8217–8243. doi:10.1029/96jb02655

CrossRef Full Text | Google Scholar

Goldfinger, C., Nelson, C. H., Morey, A. E., Johnson, J. E., Patton, J. R., Karabanov, E., et al. (2012). Turbidite Event History—Methods and Implications for Holocene Paleoseismicity of the Cascadia Subduction Zone. U.S. Geological Survey, 170. Professional Paper 1661–F.

Google Scholar

Han, S., Carbotte, S. M., Canales, J. P., Nedimović, M. R., and Carton, H. (2018). Along-Trench Structural Variations of the Subducting Juan de Fuca Plate From Multichannel Seismic Reflection Imaging. J. Geophys. Res. Solid Earth 123 (4), 3122–3146. doi:10.1002/2017JB015059

CrossRef Full Text | Google Scholar

Han, S., Carbotte, S. M., Canales, J. P., Nedimović, M. R., Carton, H., Gibson, J. C., et al. (2016). Seismic Reflection Imaging of the Juan de Fuca Plate From Ridge to Trench: New Constraints on the Distribution of Faulting and Evolution of the Crust Prior to Subduction. J. Geophys. Res. Solid Earth 121 (3), 1849–1872. doi:10.1002/2015JB012416

CrossRef Full Text | Google Scholar

Hyndman, R. D., and Wang, K. (1993). Thermal Constraints on the Zone of Major Thrust Earthquake Failure: The Cascadia Subduction Zone. J. Geophys. Res. 98, 2039–2060. doi:10.1029/92jb02279

CrossRef Full Text | Google Scholar

Kenyon, C. B. (2016). A 3-D Tomographic Model of the P-Wave Velocity Structure of the Central Cascadia Forearc. MS Thesis. United States: Oregon State University.

Google Scholar

Kimura, G., Nakamura, Y., Shiraishi, K., Fujie, G., Kodaira, S., Tsuji, T., et al. (2022). Nankai Forearc Structural and Seismogenic Segmentation Caused by a Magmatic Intrusion off the Kii Peninsula. Geochem. Geophys. Geosystems 23. doi:10.1029/2022gc010331

CrossRef Full Text | Google Scholar

Llenos, A. L., and McGuire, J. J. (2007). Influence of Fore-Arc Structure on the Extent of Great Subduction Zone Earthquakes. J. Geophys. Res. Solid Earth 112 (B9). doi:10.1029/2007jb004944

CrossRef Full Text | Google Scholar

Lotto, G. C., Dunham, E. M., Jeppson, T. N., and Tobin, H. J. (2017). The Effect of Compliant Prisms on Subduction Zone Earthquakes and Tsunamis. Earth Plan. Sci. Lett. 458, 213–222. doi:10.1016/j.epsl.2016.10.050

CrossRef Full Text | Google Scholar

Ludwig, W. J., Nafe, J. E., and Drake, C. L. (1970). “Seismic Refraction,” in The Sea. Editor A. E. Maxwell (New York: Wiley-Interscience), 53–84. v. 4.

Google Scholar

MacKay, M. E. (1995). Structural Variation and Landward Vergence at the Toe of the Oregon Accretionary Prism. Tectonics 14, 1309–1320. doi:10.1029/95tc02320

CrossRef Full Text | Google Scholar

McNeill, L. C., Goldfinger, C., Kulm, L. D., and Yeats, R. S. (2000). Tectonics of the Neogene Cascadia Forearc Basin: Investigations of a Deformed Late Miocene Unconformity. Geol. Soc. Am. Bull. 112, 1209–1224. doi:10.1130/0016-7606(2000)112<1209:TOTNCF>2.0.CO;2

CrossRef Full Text | Google Scholar

Morton, E. A., Bilek, S. L., and Rowe, C. A. (2023). Cascadia Subduction Zone Fault Heterogeneities From Newly Detected Small Magnitude Earthquakes. Jour. Geophys. Res. 128, e2023JB026607. doi:10.1029/2023JB026607

CrossRef Full Text | Google Scholar

Morton, E. A., Bilek, S. L., and Rowe, C. A. (2018). Newly Detected Earthquakes in the Cascadia Subduction Zone Linked to Seamount Subduction and Deformed Upper Plate. Geology 46 (11), 943–946. doi:10.1130/G45354.1

CrossRef Full Text | Google Scholar

Nolan, S. (2023). MS Thesis at Oregon State University. (defense in May 2023).

Google Scholar

Norvell, B., Kyritz, T., Spinelli, G. A., Harris, R. N., Dickerson, K., Tréhu, A. M., et al. (2023). Thermally Significant Fluid Seepage Through Thick Sediment on the Juan de Fuca Plate Entering the Cascadia Subduction Zone. Geochem. Geophys. Geosystems 24, e2203GC010868. doi:10.1029/2023GC010868

CrossRef Full Text | Google Scholar

Preston, L. A., Creager, K. C., Crosson, R. S., Brocher, T. M., and Tréhu, A. M. (2003). Intraslab Earthquakes: Dehydration of the Cascadia Slab. Science 302, 1197–1200. doi:10.1126/science.1090751

PubMed Abstract | CrossRef Full Text | Google Scholar

Philip, B. T., Solomon, E. A., Kelley, D. S., Tréhu, A. M., Whorley, T. L., Roland, E., et al. (2023). Fluid Sources and Overpressures Within the Central Cascadia Subduction Zone Revealed by a Warm, High-Flux Seafloor Seep. Sci. Adv. 9, eadd6688. doi:10.1126/sciadv.add6688

PubMed Abstract | CrossRef Full Text | Google Scholar

Salleres, V., and Ranero, C. R. (2019). Upper-Plate Rigidity Determines Depth-Varying Rupture Behavior of Megathrust Earthquakes. Nature 576. doi:10.1038/s41586-019-1784-0

CrossRef Full Text | Google Scholar

Schmalzle, G. M., McCaffrey, R., and Creager, K. C. (2014). Central Cascadia Subduction Zone Creep. Geochem. Geophys. Geosystems 15 (4), 1515–1532. doi:10.1002/2013GC005172

CrossRef Full Text | Google Scholar

Snavely, P. D., MacLeod, N. S., Wagner, H. C., and Rau, W. W. (1976). Geologic Map of the Waldport and Tidewater Quadrangles, Lincon, Lane, and Benton Counties, Oregon, U.S. Geological Survey Miscellaneous Investigations Map I-866, Scale 1:62,500. Available at: https://ngmdb.usgs.gov/Prodesc/proddesc_9768.htm.

Google Scholar

Snavely, P. D., and MacLeod, N. S. (1974). Yachats Basalt—An Upper Eocene Differentiated Volcanic Sequence in the Oregon Coast Range: U.S. Geol. Surv. J. Res. 2, 395–403.

Google Scholar

Song, T. A., and Simons, M. (2003). Large Trench-Parallel Gravity Variations Predict Seismogenic Behavior in Subduction Zones. Science 301, 630–633. doi:10.1126/science.1085557

PubMed Abstract | CrossRef Full Text | Google Scholar

Stone, I., Vidale, J. E., Han, S., and Roland, E. (2018). Catalog of Offshore Seismicity in Cascadia: Insights into the Regional Distribution of Microseismicity and its Relation to Subduction Processes. J. Geophys. Res. Solid Earth 123, 641–652. doi:10.1002/2017JB014966

CrossRef Full Text | Google Scholar

Toomey, D. R., Allen, R. M., Barclay, A. H., Bell, S. W., Bromirski, P. D., Carlson, R. L., et al. (2014). The Cascadia Initiative: A Sea Change in Seismological Studies of Subduction Zones. Oceanography 27 (2), 138–150. doi:10.5670/oceanog.2014.49

CrossRef Full Text | Google Scholar

Toomey, D. R., Solomon, S. C., and Purdy, G. M. (1994). Tomographic Imaging of the Shallow Crustal Structure of the East Pacific Rise at 9°30’N. J. Geophys. Res. Solid Earth 99 (12), 24,135–24,157. doi:10.1029/94JB01942

CrossRef Full Text | Google Scholar

Tréhu, A. M., Asudeh, I., Brocher, T. M., Luetgert, J., Mooney, W. D., Nabelek, J. L., et al. (1994). Crustal Architecture of the Cascadia Forearc. Science 265, 237–243. doi:10.1126/science.266.5183.237

PubMed Abstract | CrossRef Full Text | Google Scholar

Tréhu, A. M., Blakely, R. J., and Williams, M. C. (2012). Subducted Seamounts and Recent Earthquakes Beneath the Central Cascadia Forearc. Geology 40 (2), 103–106. doi:10.1130/G32460.1

CrossRef Full Text | Google Scholar

Tréhu, A. M., Braunmiller, J., and Davis, E. (2015). Seismicity of the Central Cascadia Continental Margin Near 44.5°N - a Decadal View. Seis. Res. Lett. 86, 819–829. doi:10.1785/0220140207

CrossRef Full Text | Google Scholar

Tréhu, A. M., Braunmiller, J., and Nabelek, J. L. (2008). Probable Low-Angle Thrust Earthquakes on the Juan de Fuca-North America plate boundary. Geology 36, 127–130. doi:10.1130/g24145a.1

CrossRef Full Text | Google Scholar

Tréhu, A. M., Lin, G., Maxwell, E., and Goldfinger, C. (1995). A Seismic Reflection Profile Across the Cascadia Subduction Zone Offshore Central Oregon: New Constraints on Methane Distribution and Crustal Structure. J. Geophys. Res. 100, 15,101–15,116. doi:10.1029/95JB00240

CrossRef Full Text | Google Scholar

Tréhu, A. M., Brocher, T. M., Creager, K., Fisher, M., Preston, L., Spence, G., et al. (2002). “Geometry of the Subducting Juan de Fuca Plate: New Constraints From SHIPS98,” in The Cascadia Subduction Zone and Related Subduction Systems-Seismic Structure, Intraslab Earthquakes and Processes, and Earthquake Hazards: U.S. Geological Survey Open-File Report 02-328, and Geological Survey of Canada Open File 4350. Editors S. Kirby, K. Wang, and S. Dunlop, 25–32. Available at: http://geopubs.wr.usgs.gov/openfile/of02-328.

Google Scholar

Tréhu, A. M., Wilcock, W. S. D., Hilmo, R., Bodin, P., Connolly, J., Roland, E. C., et al. (2018). The Role of the Ocean Observatories Initiative in Monitoring the Offshore Earthquake Activity of the Cascadia Subduction Zone. Oceanography 31, 104–113. doi:10.5670/oceanog.2018.116

CrossRef Full Text | Google Scholar

Turcotte, D. L., and Schubert, G. (1982). Geodynamics. 1st ed. John Wiley and Sons, Inc.

Google Scholar

US Geological Survey (2002). Digital Aeromagnetic Datasets for the Conterminous United States and Hawaii - A Companion to the North American Magnetic Anomaly Map. U.S. Geological Survey Open-File Report 02-361. Available at: https://pubs.usgs.gov/of/2002/ofr-02-361/mag_home.htm.

Google Scholar

Violay, M., Gibert, B., Mainprice, D., Evans, B., Dautria, J.-M., Azais, P., et al. (2012). An Experimental Study of the Brittle-Ductile Transition of Basalt at Oceanic Crust Pressure and Temperature Conditions. J. Geophys. Res. 117. doi:10.1029/2011JB008884

CrossRef Full Text | Google Scholar

Walton, M. A. L., Staisch, L., Dura, T., Pearl, J. K., Sherrod, B., Gomberg, J., et al. (2021). Toward an Integrative Geological and Geophysical View of Cascadia Subduction Zone Earthquakes. Annu. Rev. Geophys 49, 367–398. doi:10.1146/annurev-earth-071620-065605

CrossRef Full Text | Google Scholar

Wang, P. L., Engelhart, S. E., Wang, K., Hawkes, A. D., Horton, B. P., Nelson, A. R., et al. (2013). Heterogeneous Rupture in the Great Cascadia Earthquake of 1700 Inferred From Coastal Subsidence Estimates. J. Geophys. Res. Solid Earth 118 (5), 2460–2473. doi:10.1002/jgrb.50101

CrossRef Full Text | Google Scholar

Watt, J. T., and Brothers, D. S. (2020). Systematic Characterization of Morphotectonic Variability Along the Cascadia Convergent Margin: Implications for Shallow Megathrust Behavior and Tsunami Hazards. Geosphere 17, 95–117. doi:10.1130/GES02178.1

CrossRef Full Text | Google Scholar

Watts, A. B. (2001). Isostasy and Flexure of the Lithosphere. 1st ed. Cambridge University press, 480.

Google Scholar

Wells, R., Bukry, D., Friedman, R., Pyle, D., Duncan, R., Haeussler, P., et al. (2014). Geologic History of Siletzia, a Large Igneous Province in the Oregon and Washington Coast Range: Correlation to the Geomagnetic Polarity Time Scale and Implications for a Long-Lived Yellowstone Hotspot. Geosphere 10 (4), 692–719. doi:10.1130/ges01018.1

CrossRef Full Text | Google Scholar

Wells, R. E., and Blakely, R. J. (2019). Late Eocene Extension of the Cascadia Forearc and Its Control on Modern Deformation. AGU Fall Meet. Abstr. 2019, T31C–T06.

Google Scholar

Wells, R. E., Blakely, R. J., Sugiyama, Y., Scholl, D. W., and Dinterman, P. A. (2003). Basin-Centered Asperities in Great Subduction Zone Earthquakes: A Link Between Slip, Subsidence, and Subduction Erosion? J. Geophys. Res. Solid Earth 108 (10). doi:10.1029/2002JB002072

CrossRef Full Text | Google Scholar

Wells, R. E., Weaver, C. S., and Blakely, R. J. (1998). Fore-Arc Migration in Cascadia and its Neotectonic Significance. Geology 26 (8), 759. doi:10.1130/0091-7613(1998)026<0759:FAMICA>2.3.CO;2

CrossRef Full Text | Google Scholar

Williams, M. C., Trehu, A. M., and Braunmiller, J. (2011). Seismicity at the Cascadia Plate Boundary Beneath the Oregon Continental Shelf. Bull. Seismol. Soc. Am. 101 (3), 940–950. doi:10.1785/0120100198

CrossRef Full Text | Google Scholar

Wilson, D. S., and McCrory, P. A. (2022). A Late Cenozoic Kinematic Model for Fault Motion within Greater Cascadia. Geochem. Geophys. Geosystems 23. doi:10.1029/2022gc010442

CrossRef Full Text | Google Scholar

Wilson, D. S. (2002). “The Juan de Fuca Plate and Slab: Insochron Structure and Cenozoic Plate Motions,” in The Cascadia Subduction Zone and Related Subduction Systems-Seismic Structure, Intraslab Earthquakes and Processes, and Earthquake Hazards: U.S. Geological Survey Open-File Report 02-328, and Geological Survey of Canada Open File 4350. Editors S. Kirby, K. Wang, and S. Dunlop, 9–12. Available at: http://geopubs.wr.usgs.gov/openfile/of02-328.

Google Scholar

Yeats, R. S., Kulm, L. D., Goldfinger, C., and McNeill, L. C. (1998). Stonewall Anticline: An Active Fold on the Oregon Continental Shelf. Geol. Soc. Am. Bull. 110, 572–587. doi:10.1130/0016-7606(1998)110<0572:saaafo>2.3.co;2

CrossRef Full Text | Google Scholar

Keywords: forearc structure, Cascadia, earthquake hazards, subduction zone structure, controlled source seismic imaging

Citation: Tréhu AM, Davenport K, Kenyon CB, Carbotte SM, Nabelek JL, Toomey DR and Wilcock WSD (2023) Deformation of the Juan de Fuca Plate Beneath the Central Cascadia Continental Margin (44°-45°N) in Response to an Upper Plate Load. Earth Sci. Syst. Soc. 3:10085. doi: 10.3389/esss.2023.10085

Received: 28 April 2023; Accepted: 12 October 2023;
Published: 02 November 2023.

Edited by:

Richard Palin, University of Oxford, United Kingdom

Reviewed by:

Ray Wells, U.S. Geological Survey, United States
Nicholas Hayman, Oklahoma Geological Survey (OGS), United States

Copyright © 2023 Tréhu, Davenport, Kenyon, Carbotte, Nabelek, Toomey and Wilcock. This is an open-access article distributed under the terms of the Creative Commons Attribution License (CC BY). The use, distribution or reproduction in other forums is permitted, provided the original author(s) and the copyright owner(s) are credited and that the original publication in this journal is cited, in accordance with accepted academic practice. No use, distribution or reproduction is permitted which does not comply with these terms.

*Correspondence: Anne M. Tréhu, trehua@oregonstate.edu

Present addresses: Kathy Davenport, Sandia National Laboratory, Albuquerque, NM, United States
Christopher B. Kenyon, Tetra Tech, Bothell, WA, United States

Disclaimer: All claims expressed in this article are solely those of the authors and do not necessarily represent those of their affiliated organizations, or those of the publisher, the editors and the reviewers. Any product that may be evaluated in this article or claim that may be made by its manufacturer is not guaranteed or endorsed by the publisher.