Limited Shallow Slip for the 2020 Simeonof Earthquake, Alaska, Constrained by GNSS‐Acoustic

The 22 July 2020 Mw7.8 Simeonof earthquake was a deep megathrust event that ruptured along the Shumagin segment of the Alaska‐Aleutian subduction zone. This earthquake occurred ∼250 km from a seafloor geodetic GNSS‐Acoustic site IVB1, where we observed a velocity of 3.78 ± 1.15 cm/yr with the down‐going slab prior to the earthquake followed by 0.6 ± 0.7 eastward and −15.5 ± 0.8 cm northward coseismic offset. We computed a slip model of the coseismic rupture using the static offset at IVB1 alongside regional continuous GNSS and strong motion stations. The small static horizontal offset at the site precludes significantly shallower rupture than previously inferred from terrestrial observations, confirming that the Simeonof earthquake was a deep megathrust earthquake. The observed site velocity implies partial locking prior to the earthquake, implying significant shallow strain accumulation such that the small coseismic offset is unlikely to have relieved all of the accumulated strain since the last coseismic rupture.

• We observed seafloor offset at a GNSS-Acoustic site along the Alaska subduction zone following the 22 July 2020 Mw7.8 Simeonof earthquake • Additional slip updip of the main earthquake rupture area is not required to generate the observed seafloor offset • The shallow plate interface was partially locked prior to the Simeonof earthquake and may still host unrelieved strain

Supporting Information:
Supporting Information may be found in the online version of this article. and a Mw 7.5 in 1948 (Davies et al., 1981). Ryan et al. (2012) concluded that the Shumagin gap could potentially rupture as part of a multi-segment rupture with Mw > 9.0 regardless of whether the shallow plate interface was locked or creeping, a conclusion supported by paleoseismic studies that identified tsunami deposits along the eastern Alaska-Aleutian Subduction Zone resulting from large megathrust earthquakes with a recurrence interval of ∼560 years (Carver & Plafker, 2008).
The fault locking models of Li and Freymueller (2018) and Drooff and Freymueller (2021) resolve the Shumagin gap as a region with low locking fraction at depth and variable locking fraction along-strike at depths shallower than 25 km, transitioning from a primarily locked shallow plate interface in the Kodiak segment to the east to a primarily creeping shallow plate interface in the Sanak segment to the west. The Shumagin segment was inferred to be partially locked (locking fractions of <0.4). These models resolve the diminishing plate locking as a function of increasing fault depth but the fraction of the shallow plate locking is poorly constrained by the available land-based GNSS data due to the absence of offshore measurements. Despite low locking at depth, sufficient strain to enable a Mw 8.0 earthquake could have accumulated since the previous rupture along this segment of the plate interface, a Mw 7.4 that occurred in 1917 (Estabrook & Boyd, 1992).
Previous studies of the Simeonof earthquake fail to uniquely resolve the updip boundary of coseismic slip or interseismic locking fraction using data only collected at land stations, far from the shallow subduction zone interface. This ambiguity is not unique to the Alaska-Aleutian Subduction Zone but rather a persistent limitation of models that attempt to constrain fault parameters such as locking in offshore locales using exclusively terrestrial data (see analogous locking models for the Cascadia subduction zone from ). As a result, we cannot accurately estimate the strain accumulated in the shallow subduction zone, the rate at which it accumulated, or how much may have been relieved without independent seafloor geodetic measurements. The state of locking on the updip portion of the megathrust remains an important question because of the implications for tsunamigenesis, as demonstrated by the huge tsunami from the rupture near the Japan Trench during the Tohoku-Oki earthquake (Fujiwara et al., 2011).
In this paper we present offshore Global Navigation Satellite System-Acoustic (GNSS-Acoustic) data collected at a site, IVB1, ∼250 km from the epicenter of the Simeonof earthquake ( Figure 1). GNSS-Acoustic is a seafloor geodetic technique used to position points on the seafloor with centimeter-level horizontal accuracy (Chadwell & Spiess, 2008;Fujita et al., 2006;Spiess et al., 1998). This technique has previously been used to capture the offshore surface deformation from past earthquakes, including two M7 earthquakes off the Kii Peninsula in 2004 (Kido et al., 2006;Tadokoro et al., 2006), the 2005 Mw 7.2 earthquake off Miyagi prefecture (Matsumoto et al., 2006;Sato et al., 2011b), and the 2011 Mw 9.0 Tohoku-Oki earthquake (Kido et al., 2011a;Sato et al., 2011a). The data we present are part of an array ( Figure 1) that directly captured the horizontal co-seismic and early postseismic surface deformation of the Simeonof, Sand Point, and Chignik  earthquakes along the trench of the Alaska-Aleutian Subduction Zone, and are the first seafloor displacement observations of an earthquake outside of Japan coastal waters. With these new measurements we are able to further resolve the updip extent of the seismic rupture and its implications for strain in the shallow subduction zone.

GNSS-Acoustic Data
GNSS-Acoustic combines two types of instrumentation, an array of acoustic transponders resting on the seafloor and a surface platform to determine the position of a seafloor benchmark in a global reference frame (Bürgmann & Chadwell, 2014). An autonomous waveglider, which uses wave energy for propulsion, with solar powered electronics is used as the surface platform. The waveglider is outfitted with a transducer on its underside, and two GNSS antennas on its top side. During a GNSS-Acoustic survey, the waveglider is remotely piloted to hold station above the center of an array of three transponders deployed in a concentric circle with radius approximately equal to the seafloor depth, from which it interrogates all of the seafloor transponders in the array ( Figure  S1 in Supporting Information S1). The waveglider collects acoustic two-way travel times between itself and the seafloor transponders to determine its position relative to the seafloor while also collecting data from multiple GNSS constellations (i.e., GPS, GLONASS). If sufficient information is known about the sound velocity profile of the water column, these data are used with a one-dimensional acoustic ray-trace inversion to locate the position of the transponder array center in a global reference frame. A further overview of the GNSS-Acoustic technique is provided by Bürgmann and Chadwell (2014).
We established three GNSS-Acoustic sites along the trench of the Alaska-Aleutian Subduction Zone in 2018 with the goal of estimating the plate locking of the shallow megathrust of the Shumagin gap and nearby segments ( Figure 1 We deployed wavegliders at IVB1 and SPT1 in the fall of 2020 to measure the co-seismic and initial post-seismic deformation from the Simeonof earthquake. (These surveys concluded before the 19 October 2020 Sand Point earthquake.) Thus, we shall only consider the data collected at IVB1 for this study since that is the only site that has position measurements both before and after the Simeonof earthquake as well as the closest site to the rupture area.
The IVB1 site was occupied during three separate campaigns in 2018, 2019, and 2020. GNSS data were processed using GipsyX (Bertiger et al., 2020), and a ray-trace inversion was performed on the combined GNSS positions and acoustic two-way travel times using the methods described in Chadwell and Sweeney (2010) to solve for the GNSS-Acoustic array positions. Each occupation has a time series of residuals that capture the difference between the modeled and measured two-way travel times for each transponder; these residuals are used to assess the data quality of each epoch. We remove epochs during which the waveglider was >150 m from the center of the GNSS-Acoustic array or the array center derived from the GNSS-Acoustic residuals was >5 cm removed from the mean, as these correspond to GNSS cycle slips or breakdowns of the horizontal layered sound velocity assumption. We prefer to collect data for extended occupations, ideally long enough to capture 3-5 "equivalent days" of data to average out fluctuations due to sound speed variations.
The uncertainty due to sound velocity variations are unconstrained since we only have one Conductivity Temperature Depth profile collected in 2018 constraining the sound velocity at this site. However, because the waveglider path during the survey adheres closely to the transponder array center the unknown variations due to sound velocity are mostly common mode along the three ray paths. If the sound velocity profile approximates a horizontal layered structure, the sound velocity variations primarily affect the vertical position estimate and have little effect on the estimate of the horizontal array position. This is verified in Figure 2a, which plots the GNSS-Acoustic residuals between measured and modeled two-way travel times converted into "apparent" horizontal positions of the IVB1 array center by spatial averaging and multiplying by the harmonic mean sound velocity. With only three transponders in the array, we cannot decouple the vertical component of the transponder positions from changes in the Nadir Total Delay (NTD) and will not report these values.
The IVB1 site moved 3.2 ± 0.6 North and 2.4 ± 0.6 cm West between 2018 and 2019 at an apparent velocity of 3.78 ± 1.15 cm/yr parallel to the downgoing slab (Figures 2b and 2c). Assuming this is the interseismic velocity, we estimate that the coseismic offset is 0.6 ± 0.7 east and −15.5 ± 0.8 cm  north between the projected position immediately prior to the earthquake and the 2020 position observation. The reported coseismic uncertainties include the propagated uncertainties from the positions, but does not reflect uncertainty in the proportion of coseismic versus postseismic displacement. Furthermore, since the 2020 survey was conducted before the 19 October 2020 Mw7.6 strike-slip earthquake, we do not consider that event in this analysis.
The primary source of uncertainty in the GNSS-Acoustic measurements is from variation in the sound velocity profile due to surface tides, eddies, internal tides, and internal waves. Errors are reduced by collecting and averaging over many days of data, ideally enough that the time series would be 3-5 days long if the data gaps were removed. In Figure 2a, the 2018 survey appears more precise, but is actually more poorly constrained than the other surveys due to a shorter survey duration that cannot properly constrain the oceanographic signal in the residuals (See also Table S1 in Supporting Information S1). Due to the limited solar energy offshore Alaska in the early fall we had to stop data collection to recharge the waveglider batteries often. The waveglider remained on station for up to a month to collect enough data and average out the effects of ocean variability. As a result, the positions must be considered as an average array position during that time window. This has implications on how best to estimate the displacement following the earthquake since there was likely 1-2 cm of postseismic deformation at the IVB1 site during this time period if it was deforming at a similar rate to the nearby land stations (see Supporting Information S1). With only one data point collected after the earthquake it is difficult to estimate the postseismic deformation, but this may be improved in future field seasons as these stations are re-surveyed. In addition, which dates are chosen in a processing for the interseismic site velocity can subtly alter the estimate. However, these data easily capture the surface displacements due to the earthquake, showing horizontal offsets in excess of 15 cm compared to prior surveys with an uncertainty on the order of a centimeter.

Slip Model
We repeated the analysis of Crowell and Melgar (2020) with the addition of the offshore GNSS-A data, and generated a kinematic rupture model of the Simeonof earthquake using seven regional continuous 1-Hz GNSS stations, six strong motion sensors, and the static coseismic offset at IVB1. We used this data as input for an inverse model using the MudPy software developed by Melgar and Bock (2015) and a fault interface generated by SLAB 2.0 (Hayes et al., 2018). The hypocenter, rupture speed (3.0 km/s), origin time, filtering, 1-D velocity structure, and source time windows are all identical to those used in Crowell and Melgar (2020). Green's functions are computed using the frequency-wavenumber approach of Zhu and Rivera (2002). The input data were weighted by their formal uncertainty with the exception of the vertical component of IVB1, which was left unconstrained since we do not have the requisite information to decouple the vertical position from the NTD.
A comparison between rupture models published by Crowell and Melgar (2020) and this study is shown in Figure 3. MudPy generates a suite of models with different smoothness constraints. The models shown have the same inversion parameters in order to facilitate a direct comparison, resulting in similar rupture areas along the deep subduction zone interface. However, the moment of the model calculated including the static GNSS-A offset is slightly smaller, being 95% as large as the moment calculated in the Crowell and Melgar (2020) model. The models' fit to the static offset at IVB1 reveals another subtle difference between the two models. Although the Crowell and Melgar (2020) inversion does not consider the static offset at IVB1, we may still calculate the static offset predicted by the inversion's rupture patch by multiplying it by the Green's functions calculated for the IVB1 site in our inversion. Doing so reveals a predicted offset at IVB1 of 20.6 cm oriented normal to the subduction zone, consistent with the primarily dip-slip rupture reported in the Crowell and Melgar (2020) study. In comparison, the offset we directly measured at IVB1 is smaller and oriented almost directly South. Taken with the smaller moment, the observed offset implies that the rupture patch may not have extended as far up dip as previously thought and may have a smaller oblique component, at least in the shallow subduction zone. A checkerboard test for resolution performed by Crowell and Melgar (2020) showed weak resolution in the region around IVB1, so the result shown here is more reliable in the up-dip region.
Since this model only differs from the previous study with the addition of the IVB1 offset, its resolution has already been well-documented (see Crowell and Melgar (2020)). However, it is still important to clarify the additional resolution provided by IVB1. We performed a sensitivity analysis of the IVB1 site by predicting the horizontal static offset resulting from a number of synthetic rupture patches comparable to the Simeonof rupture ( Figure S2 in Supporting Information S1). As we shifted the rupture patch toward IVB1, we found that the predicted seafloor offset would have exceeded the measured coseismic offset if the rupture model was <60 km from IVB1. This agrees with the resolution matrix of the site ( Figure S3 in Supporting Information S1) and places a bound on the updip extent of coseismic rupture during the Simeonof earthquake comparable to that predicted by Crowell and Melgar (2020), although that study considered the predicted vertical offset from a fault slip patch rather than the predicted horizontal offset.

Discussion
There are two important observations to make regarding the displacement of seafloor geodetic station IVB1: its apparent velocity prior to the Simeonof earthquake and the offset after the earthquake. The IVB1 site moved 3.2 ± 0.6 North and 2.4 ± 0.6 cm West between 2018 and 2019 at an apparent velocity of 3.78 ± 1.15 cm/yr in the North American reference frame, parallel to the downgoing slab and a significant fraction (62%) of the 6.1 cm/ yr plate convergence rate between the Pacific and North American plates at this section of the subduction zone. Because this velocity is only constrained by two observations, it has a large uncertainty that makes it difficult to formally compare with previous predictions. However, this initial estimate is comparable to that predicted by the locking models of Li and Freymueller (2018) and Drooff and Freymueller (2021) and larger than that predicted by Xiao et al. (2021). If we assume that the seafloor geodetic station IVB1 moved at this velocity prior to the Simeonof earthquake, then we observed an offset 0.90 ± 0.70 in the east component and −18.62 ± 0.81 cm in the north component that we may attribute to the Simeonof earthquake, which includes some postseismic deformation. Using an estimate of the postseismic deformation to correct these data yields an estimate of 0.6 ± 0.7 east and −15.5 ± 0.8 cm north coseismic offset. These measurements were collected before 19 October 2020, so this offset estimate has no contamination from the Mw 7.6 strike-slip earthquake that occurred near Sand Point on that date. This coseismic offset is much larger than expected based on a linear extrapolation of the 2018-2020 velocity constrained prior to the earthquake. However, adding this static offset to the previous slip model by Crowell and Melgar (2020) did not significantly alter the earthquake rupture and in fact decreased the moment of the rupture by ∼5%, consistent with previous observations that the Simeonof rupture occurred primarily in the deeper subduction zone (Crowell & Melgar, 2020;Liu et al., 2020;Xiao et al., 2021;Ye et al., 2021). The minimal shallow offset is also consistent with the tsunami observations, which were <30 cm on the nearby Shumagin Islands (Larson et al., 2020). Our synthetic models shed some further context to the result by showing that the IVB1 station is not sensitive to the deeper subduction zone interface where the primary earthquake rupture occurred but rather to the fault interface beneath an area within ∼120 km from the seafloor station. Thus, we can infer that there may have been some slip along the shallow plate interface near the site, but it must have been smaller than the primary rupture by about an order of magnitude.  Crowell and Melgar (2020). Slip patches with <0.3 m of slip are omitted. Orange vectors denote the measured co-seismic offset at regional geodetic sites and gray vectors denote the offsets predicted by the slip inversions. For the left plot, the gray vector for IVB1 plots on top of the orange vector.
10.1029/2023GL105045 6 of 7 The coseismic offset we observed at station IVB1 indicates the Simeonof earthquake is unlikely to have relieved all of the strain that may have accumulated along a partially locked shallow plate interface as constrained by the site velocity from 2018 to 2019. Furthermore, the full postseismic response to the Simeonof earthquake has not yet occurred and will not resolve for many years. Thus, we do not yet know how much strain may ultimately be relieved in the shallow subduction zone and how much of a seismic and tsunami hazard it may pose in the future, although a baseline hypothesis is that it remains partially locked and able to rupture seismically.
An alternative hypothesis is that the offset we observe at station IVB1 is not due to coseismic slip along the plate interface but rather from another source. For instance, the seafloor could have moved as part of a slump or landslide following the earthquake, but this seems unlikely because the offset at IVB1 is consistent with deep slip and smaller in magnitude than we would expect from a slump or landslide covering such a large area. Furthermore, there was not a large tsunami that such a displacement would likely have produced. Another hypothesis is that the offset at IVB1 is primarily postseismic rather than coseismic, as was inferred at the SEM1 array following the neighboring Chignik earthquake in 2021 . However, this hypothesis would also yield a rupture model with little shallow slip similar to the one we present in this study. We will need to collect further measurements at IVB1 so that we can properly assess its postseismic response.

Conclusions
We observed seafloor offset at a seafloor geodetic GNSS-A site IVB1 near the 2020 Mw7.8 Simeonof earthquake. This site was moving 3.78 ± 1.15 cm/yr parallel to the downgoing slab prior to the earthquake, after which it was offset 0.90 ± 0.70 to the east and −18.62 ± 0.81 cm to the north. This offset is small enough in magnitude to confirm that the Simeonof earthquake must have been a primarily deep rupture as previously observed but large enough to require some minor slip along the shallow plate interface in the vicinity of the seafloor site. Synthetic tests confirmed that the station IVB1 is capable of capturing meter-scale co-seismic slip approximately 120 km away from the site. In this study the earthquake we observed had a primarily deep rupture that was verified by the small offset at the seafloor GNSS-A site. In comparison, other studies have identified cases where coseismic slip greater than 1 m did extend to shallow depths, such as has been observed for the Chignik  and Tohoku (Kido et al., 2011;Sato et al., 2011a) earthquakes. Thus, even though this study confirms previously published models of the Simeonof earthquake (Crowell & Melgar, 2020;Liu et al., 2020;Ye et al., 2021), the seafloor geodetic measurements are necessary to understand the shallow slip behavior in subduction zone earthquakes.
The interseismic site velocity we observed at station IVB1 implies partial locking along the shallow plate interface prior to the Simeonof earthquake as predicted by previous studies Li and Freymueller (2018) and Drooff and Freymueller (2021). Because of the small magnitude of shallow slip inferred following the Simeonof earthquake, it is likely that the accumulated strain along the shallow subduction zone interface was not entirely relieved and the potential for further rupture in a shallow megathrust earthquake remains. Future observations will help to resolve the postseismic slip and afterslip at station IVB1 in order to further quantify the amount of strain release along the shallow plate interface and the potential for a future rupture.

Data Availability Statement
The GNSS-Acoustic data used in this study , including acoustic two-way travel times, INS orientations, RINEX, CTDs, a priori site information, and final positions, are available at the GAGE Facility operated by the EarthScope Consortium and may be accessed at the following link: https://doi.org/10.7283/ C8GY-SG81.