On the Relationship Between Offshore Geodetic Coverage and Slip Model Uncertainty: Analog Megathrust Earthquake Case Studies

We apply a geodetic slip inversion technique to analog subduction megathrust earthquakes to demonstrate how limited offshore geodetic coverage affects coseismic slip models. We analyzed two archetypical megathrust earthquakes: trench‐breaking and non‐trench‐breaking earthquakes. Slip inversion models of analog earthquakes show quantitative and qualitative changes as a function of offshore coverage. Shallow slip cannot be resolved if the observation coverage of the offshore segment is <50%. Moreover, the slip pattern of shallow events flips from landward to trenchward skewed as offshore coverage reduces to <40%. The estimated slip for both event types converges to a similar unimodal pattern when there is no offshore coverage. We infer 5–20% slip overestimation when the observations are above the high slipping zone during trench‐breaking events versus 5–10% underestimation during non‐trench‐breaking events if observations are land limited. The moment magnitude derived for trench‐breaking ruptures might be significantly affected (ΔMw ~ 0.5).


Introduction
Geodetic observations deliver the key data for megathrust earthquakes through inversion of surface deformation for slip at depth (e.g., Bedford et al., 2013;Hoffmann et al., 2018;Konca et al., 2008;Lin et al., 2013;Loveless, 2017;Ozawa et al., 2011).While the Earth surface above the shallow source of megathrust earthquakes can be observed directly in continental settings (e.g., Cheloni et al., 2014;Stevens & Avouac, 2015), subduction megathrust earthquakes occur generally offshore.Simultaneous use of the global navigation satellite system (GNSS) and offshore/near-source observations such as GPS-acoustic and/or bottom pressure sensors (e.g., Iinuma et al., 2012;Maeda et al., 2011) has provided a powerful tool to better monitor slip behavior in subduction zones.The effects of network geometries and measurement density with respect to the location of the megathrust fault and the resulting uncertainties are generally acknowledged, and their impact on the magnitude of slip, size, and pattern of slip area has been shown in case studies using optimization (e.g., Romano et al., 2012).However, a case study demonstrating how limited coverage affects slip models qualitatively and quantitatively is missing.To close this gap, we here use a new set of analog subduction megathrust earthquakes generated by a seismotectonic scale modeling approach (e.g., Corbi et al., 2013Corbi et al., , 2017Corbi et al., , 2019;;Rosenau et al., 2009Rosenau et al., , 2010Rosenau et al., , 2017Rosenau et al., , 2019;;van Rijsingen et al., 2019).Such a laboratory approach provides experimental observations of well-monitored analog earthquakes to study the sensitivity of slip inversion models to geodetic coverage more systematically in an emergent setting.In the laboratory-scale

RESEARCH LETTER
10.1029/2020GL088266 Key Points: • Analog earthquake experiments provide case studies to evaluate the linkage between offshore geodetic observations and coseismic slip models • Slip inversion models of analog earthquakes show quantitative and qualitative changes as a function of offshore geodetic coverage • Inverted slip for trench-breaking and non-trench-breaking events converges to a similar pattern when there is no offshore geodetic coverage

Supporting Information:
• Supporting Information S1 experiment presented, physically self-consistent ruptures (see review in Rosenau et al., 2017) nucleate under well-defined boundary conditions and monitored geodesy-like with optimal surface coverage using digital image correlation techniques.With respect to the later, we introduce AGSIT-the analog geodetic slip inversion technique-to resolve the full complexity of the coseismic slip pattern emerging experimentally instead of employing either purely synthetic data (e.g., checkerboard tests) or natural cases with limited coverage.This way, we maximize realism and at the same time minimize the unknown parameters in the inversion method (e.g., fault geometry).
In this study, the established analog model generates two archetypical subduction megathrust earthquakes: those rupturing to the trench (here referred to as "trench-breaking" or A-B type events following Lay et al., 2012) and those confined to intermediate depths ("non-trench-breaking" or B-C type events) on which we base our case studies alongside synthetic tests.

Seismotectonic Scale Modeling
Seismotectonic scale models of subduction zones generating physically self-consistent analog megathrust earthquake ruptures at laboratory scale have been introduced a decade ago (Rosenau et al., 2009(Rosenau et al., , 2017) ) to study such phenomena as the interplay between short-term elastic (seismic) and long-term permanent deformation (Rosenau & Oncken, 2009), intrinsic variability of subduction earthquake ruptures (Rosenau et al., 2010), and earthquake recurrence behavior (Corbi et al., 2017(Corbi et al., , 2019(Corbi et al., , 2020;;Rosenau et al., 2019).In this study, we use this approach to generate a catalog of analog megathrust earthquakes and analyze the corresponding coseismic surface displacement.Based on a set of dimensionless numbers, the experiment dimensions mass, length, and time are scaled to nature ensuring geometric, kinematic, and dynamic similarity of the scale model and its natural prototype (King Hubbert, 1937;Rosenau et al., 2017).We note, however, that we report and analyzed all observations generally at the laboratory scale in this study.
The presented experimental setup is modified from the 3-D setup used in Rosenau et al. (2019).The subduction forearc model wedge is set up in a glass-sided box (1,000 mm across strike, 800 mm along strike, and 300 mm deep) with a 15°dipping, elastic basal conveyor belt, and a rigid backwall.A flat-topped wedge made of an elastoplastic sand-rubber mixture (50 vol.% quartz sand G12: 50 vol.%EPDM-rubber) is sieved into the setup representing a 240 km long forearc segment from the trench to the volcanic arc (Figure S1 in the supporting information).The shallow part of the wedge base that represents the seismogenic zone of the subduction plate interface incorporates a basal layer of rice grains (magenta zone in Figure S1) characterized by unstable stick-slip sliding.Stick-slip sliding in rice is governed by a rate-and-state dependent friction law similar to natural rocks (Rosenau et al., 2009(Rosenau et al., , 2017;;Scholz, 1998).Continuous compression of the wedge via the basal conveyor belt (velocity = 0.05 mm/s) simulates plate convergence and results in the quasiperiodic nucleation of stick-slip events (analog earthquakes) within the rice layer.The wedge itself responds elastically to these basal slip events similar to crustal rebound during natural subduction megathrust earthquakes.

Surface Deformation Monitoring
To monitor the surface deformation of the analog model, a stereoscopic set of two CCD (charge-coupled device) cameras (LaVision Imager pro X 11MPx, 14 bit) images the wedge surface continuously at 2.5 Hz.
To derive observational data similar to those from geodetic techniques, that is, velocities (or incremental displacements) at locations on the model surface, we use digital image correlation (DIC) (e.g., Adam et al., 2005;Rubino et al., 2015;Sutton et al., 2009) and derive the 3-D incremental surface displacements at high spatial resolution (<0.1 mm).The result is an evenly spaced and narrow grid of vectors per 0.4 s time step where each vector represents a virtual GNSS station (vGNSSs) located on the model Earth surface.

Inversion Procedure
Here, we introduce "AGSIT," which estimates slip distribution from surface displacement observations in a customized fashion tailored to analog earthquakes.All three components of the vectors of the static coseismic surface displacement derived from DIC are used as input data for the inversion procedure.We tentatively associate a virtual model coastline with the surface projection of the downdip limit of the model seismogenic zone.As a result, the number of equally spaced vGNSSs covering the model forearc up to the 10.1029/2020GL088266

Geophysical Research Letters
trench exceeds 8,000 (ca.6 km trench-parallel and 3 km trench-normal spacing), while the onshore and offshore regions contain approximately 5,000 and 3,000 vectors, respectively (Figure S1).In nature, the exact position of the seismic source is often imprecisely known and causes uncertainty in the final solution of the source inversion.In contrast, in the scale model, the seismogenic zone is predefined, and therefore the interface geometry and mechanical properties are well constrained.For the inversion, the full model megathrust interface including seismogenic and aseismic areas is discretized by evenly sized dislocations (1,344 rectangular patches).
We tackle the inverse problem using damped bound-constrained least squares (DBCLS), and the problem is designed and solved based on convex optimization (Boyd & Vandenberghe, 2004) (Text S1).To tie coseismic slip in individual fault patches to the observed surface displacement vectors at individual surface points, Green's functions for rectangular dislocations in an elastic half-space (Okada, 1985) are computed and applied and the dip-slip vector is solved for each patch.The Green's functions change from each analog earthquake to another to account for subtle geometric changes of the model wedge (mainly trench position) over multiple seismic cycles.To stabilize the inverse problem and produce a more physically feasible slip inversion model, we use a Laplacian regularization (e.g., Kositsky & Avouac, 2010).
To systematically study the effect of reduced offshore geodetic coverage on the recovery of fault kinematics, we use both synthetic and experimental slip distributions.First, a common procedure of noise-free synthetic three-component observations has been used to establish a checkerboard test on both types of events and differentiate and highlight the well-recovered versus unreliable zones on the fault surface.Afterward, 13 analog earthquakes including A-B and B-C type events are picked from different stages of the model evolution and used as case studies including a reasonable intrinsic variability (Figures S2 and S3).The resulting sets provide us a distribution of experimental slip patterns to reduce the slip pattern similarity effect.We generate more than 20 slip inversion models for each event by sequentially disregarding trenchward rows of equally weighted vGNSSs, thereby systematically reducing simulated offshore coverage (magenta markers in Figures 1 and 2).Because the basal surface geometry of the model does not vary through the experiment, a fixed megathrust geometry has been implemented, which enables a fair comparison between the inversion models.

Results of the Inversions With Limited Offshore Geodetic Coverage
In the following, we document the sensitivity of fault recoverability (Figures 1, S4, and S5) as well as of the slip pattern (Figures 2, 3, and S6-S9) and predicted tsunamigenic uplift (Figure S11) to the offshore coverage in terms of qualitative (i.e., shape) change and quantitative (i.e., location of updip/downdip limit and centroid, seismic moment) variation (Figure 4).After demonstrating the impact of offshore coverage reduction on fault recoverability based on checkerboard tests, we document the variation of the general slip pattern.
Note that all values used are at a laboratory scale or are nondimensional.Effect of the onshore observations weighting on the slip model is presented in Text S2.

Effect on Fault Recoverability
Checkerboard tests are the conventional way to assess fault recoverability using synthetic slip distribution.
Here the checkerboard tests of two sets of Green's functions, one for an A-B and one for a B-C type analog earthquake, have been calculated (Figures 1, S4, and S5).The slip value alternates between 0 and 1 mm of the dip-slip component (strike of 0°and rake of 90°), and the input checker pattern has been alternated in the fault grid while every four adjacent patches have the same input slip values (336 alternating rectangles; Figure S4).The results of the checkerboard tests for the inversion of noise-free synthetic three-component vGNSSs are visualized as the difference between the target (synthetic) slip and the inverted slip in Figure 1 for 100%, 50%, and 10% offshore coverage models.From 100% to 80% coverage, the highest recovered slip zone is in the trench domain.Progressively, reducing offshore coverage results in a corresponding loss of offshore slip recoverability.The outputs suggest that 80% (in A-B case) and 90% (in B-C case) coverages could recover the shallow near-trench zone robustly and reproduce the checker pattern.The slip recoverability reduces landward with increasing depth of the detachment due to distance offsetting of source and receiver.In the case of the land-limited coverage, the increasing offset between parameter and receiver in the trench-normal direction results in a loss of model recoverability (Williamson & Newman, 2018).
10.1029/2020GL088266   S9) exemplify the principle systematics of slip pattern variation as a function of offshore geodetic coverage for two individual events.Accordingly, the slip patterns of A-B type events show a high sensitivity to reducing offshore coverage as changes in both patterns and maximum slip have been observed.More specifically, the slip pattern flips from landward skewed to trenchward skewed (Figures 2, S6, and S7) as offshore coverage is reduced below 40%.Furthermore, the peak slip location retreat reaches a maximum of 25% while the shift in the skewness can increase the peak slip value by up to 20% (Figures 3a, S8, and S9).Along strike, as the coverage decreases, the rather uniform slip model develops into a bimodal distribution at 50% offshore coverage in some cases and eventually to a unimodal slip distribution for no offshore coverage (Figures 3c and S8).In general, we find that onshore information focuses slip in the center of the model, while offshore information tends to diffuse slip across the model interface (Figure S10).As a result, slip models constrained by intermediate coverages are hybrid and serve to 10.1029/2020GL088266 satisfy these ambiguous constraints.Coseismic seafloor uplift predicted from this pattern follows the systematics of slip variation and shows a landward shift of peak uplift along with both overestimation and underestimation of the uplift (Figure S11).Implications for tsunamigenesis are discussed in Text S3.

Effects on Slip Pattern and Peak Slip
Slip inversion models of B-C type events (Figures 2d-2f), unlike for A-B type events, are less sensitive to offshore coverage.They generally keep their unimodal character both across and along strike but show minor systematic changes.Accordingly, the location of maximum slip shifts trenchward then retreats as offshore coverage reduces such that the skew of the slip distribution flips from landward to neutral (at 20% offshore coverage) to slightly trenchward.Along with the change in skew, peak slip increases to a maximum of 110-125% for the neutrally skewed slip distribution and then decreases 5-15% from the maximum for 0% offshore coverage (Trench-parallel transect, Figures 3d and S9) while the trench-normal transect shows up to 10% peak slip underestimation for 0% offshore coverage (Figures 3b and S9).In fact, both type events converge to a very similar unimodal slip distribution if there is no offshore coverage (Figures 2c and 2f).

Effects on Rupture Limits and Centroid Location
Figures 4a and 4b show the spatial variation of the key rupture parameters as a function of offshore coverage and in relation to the intrinsic variability of all events.We use the across-strike cumulative slip distributions   10.1029/2020GL088266

Geophysical Research Letters
(Figures S13a and S13b) to further quantify the impact of offshore coverage on source parameters of the analog earthquakes.Accordingly, the locations of the 90th and 10th percentiles have been used to track the updip and downdip limit of coseismic slip, respectively, and the location of the median (50th percentile) of the cumulative slip distribution serves as a proxy for the depth of the centroid location (Figures 4a, 4b, and S13).Figures 4a and 4b visualize the sensitivity of these locations to offshore coverage.Accordingly, for A-B type events, the updip limit and centroid of the slip models retreat landward by a couple of kilometers when scaled to nature as the offshore coverage reduces to <60%.In contrast, for B-C type events, the updip limit and centroid of the slip models shift trenchward by a few kilometers as coverage is reduced to <70% but returns to the initial locations when offshore coverage is further reduced to <20%.The downdip limit is rather insensitive to the coverage for both types of events.

Effects on Seismic Moment and Magnitude
Figure 4c shows the effect of offshore coverage on the seismic moment of the analog earthquakes.Accordingly, we found no significant sensitivity of the seismic moments for the B-C type events.In contrast, for A-B type events, a seismic moment underestimation can be observed as coverage is reduced to <70%.All our slip inversion models for A-B type events arrive at seismic moments of 1.86 to 8.59 Nm corresponding to moment magnitudes M w of −5.88 to −5.44.For B-C type events, the magnitude range is notably narrower with seismic moments between 2.70 and 4.84 Nm and corresponding moment magnitudes between M w −5.77 and −5.60.

Discussion
Resolving an admissible finite-fault model requires sufficient data in the appropriate locations; otherwise, any change in the inversion process can lead to a significantly underdetermined inverse problem and consequently nonunique solutions (e.g., Minson et al., 2014, and references therein).The problem becomes more intricate in subduction megathrust settings where the observed data are typically limited to onshore regions and where resolvability decreases offshore.This problem has been highlighted by the coseismic slip models of the Tohoku-Oki 2011 earthquake derived from geodetic data sets with different offshore coverage (e.g., Ozawa et al., 2011;Romano et al., 2012;Simons et al., 2011) (Figure S15).Accordingly, substantial qualitative and quantitative differences exist among the 46 published slip models at a distance less than 150 km from the trench, and the differences become generally greater trenchward (Sun et al., 2017;Wang et al., 2018).Lack of offshore coverage may lead to both an underestimation or overestimation of slip and therefore tsunamigenic uplift during megathrust earthquakes (e.g., Williamson & Newman, 2018).Equivalently during the interseismic, limited offshore coverage leads to a limited resolution of locking near the trench (Loveless & Meade, 2010, 2011).Yohler et al. (2019) demonstrated recently how including seafloor absolute pressure gauges to the land-based GNSS stations allowed them to capture not only an early phase of shallow slow slip but how the inclusion decreased uncertainty and increased peak slip by ~23% and ~13%, respectively.Similarly, in our slip models for B-C type events, the slip amount for 0% offshore coverage could also be underestimated by 10% compared to the full coverage model when the geodetic observations only cover onshore region (0% offshore coverage) (Figure 3b).
For B-C events, where the slipping patches mostly occurred in the deeper parts of the plate interface and if vGNSSs are only located far from the source of slip, the slip inversion models could underestimate the amount of slip in a range of 5-15%.Therefore, according to our study and to Williams and Wallace (2015), the coseismic slip models which derive from either offshore observation data directly above the maximum slip zone (e.g., Iinuma et al., 2012;Romano et al., 2012) (Figure S15) or from stations which are located far from the source of slip (e.g., Ozawa et al., 2011;Romano et al., 2012) could involve the potential of overestimation and underestimation slip, respectively.We point out that only for A-B type events the combined effect of coverage-related uncertainty in rupture area and slip may be large enough to result in a significant difference in moment magnitude (ΔM w ~0.5, section 3.4) to be validated independently from seismological constraints.
Besides coseismic slip models, interseismic locking models are potentially sensitive to geodetic coverage.For instance, continental collision settings where the geodetic observations cover the shallow megathrust up to 10.1029/2020GL088266
Locking models of the Cascadia subduction zone exemplify that changing model assumptions can result in substantially different locking distribution from fully locked to trench to a shallow creeping area for the same land-based data set (McCaffrey et al., 2013;Schmalzle et al., 2014).

Conclusions
We demonstrate how coseismic slip inversion models change quantitatively and qualitatively as a function of offshore geodetic coverage via applying an original AGSIT to experimental subduction megathrust earthquakes generated by a seismotectonic scale model.For shallow, trench-breaking (A-B-type) events, changes in the slip estimations become significant when the observations cover less than 70% of the offshore distance to the trench, and a significant seismic moment underestimation (~M w = 0.5), unlike for B-C type event, is suggested as coverage is reduced to <70%.Additionally, a potential trench normal bimodality of slip distributions cannot be resolved if the offshore coverage is less than approximately 50%.An initially trench-parallel slip distribution continuity might break down into apparent patches while the maximum slip retreats landward systemically concurrent with reducing the offshore coverage.Both overestimation and underestimation were inferred in the case of B-C and A-B type events, respectively, when the observations are directly located above the high slip zone.Accordingly, the amount of slip is estimated to be 5-20% larger than the model derived from 90-100% offshore coverage.If vGNSSs are located far from a deep source of slip, such as for B-C events, the slip inversion models underestimate the slip in a range of 5-10%.

Figure 2 (
Figure 2 (Figures S6 and S7) and Figure 3 (Figures S8 andS9) exemplify the principle systematics of slip pattern variation as a function of offshore geodetic coverage for two individual events.Accordingly, the slip patterns of A-B type events show a high sensitivity to reducing offshore coverage as changes in both patterns and maximum slip have been observed.More specifically, the slip pattern flips from landward skewed to trenchward skewed (Figures 2, S6, and S7) as offshore coverage is reduced below 40%.Furthermore, the peak slip location retreat reaches a maximum of 25% while the shift in the skewness can increase the peak slip value by up to 20% (Figures3a, S8, and S9).Along strike, as the coverage decreases, the rather uniform slip model develops into a bimodal distribution at 50% offshore coverage in some cases and eventually to a unimodal slip distribution for no offshore coverage (Figures3c and S8).In general, we find that onshore information focuses slip in the center of the model, while offshore information tends to diffuse slip across the model interface (FigureS10).As a result, slip models constrained by intermediate coverages are hybrid and serve to Figure 2 (Figures S6 and S7) and Figure 3 (Figures S8 andS9) exemplify the principle systematics of slip pattern variation as a function of offshore geodetic coverage for two individual events.Accordingly, the slip patterns of A-B type events show a high sensitivity to reducing offshore coverage as changes in both patterns and maximum slip have been observed.More specifically, the slip pattern flips from landward skewed to trenchward skewed (Figures 2, S6, and S7) as offshore coverage is reduced below 40%.Furthermore, the peak slip location retreat reaches a maximum of 25% while the shift in the skewness can increase the peak slip value by up to 20% (Figures3a, S8, and S9).Along strike, as the coverage decreases, the rather uniform slip model develops into a bimodal distribution at 50% offshore coverage in some cases and eventually to a unimodal slip distribution for no offshore coverage (Figures3c and S8).In general, we find that onshore information focuses slip in the center of the model, while offshore information tends to diffuse slip across the model interface (FigureS10).As a result, slip models constrained by intermediate coverages are hybrid and serve to

Figure 1 .
Figure 1.Checkerboard test for trench-breaking A-B (a-c) and non-trench-breaking B-C (d-f) type events.The magenta and yellow markers indicate the offshore and onshore vGNSSs stations, respectively.The results are visualized as the difference between the target (synthetic) slip and the inverted slip.Note that the number of markers is smaller than the number of digital image correlation derived measurements realized in the experiments for reasons of visualization.The target slip can be found in Figure S4.The uniform colormap is adapted from Crameri (2018).

Figure 2 .
Figure 2. Examples of slip patterns and their sensitivity to offshore geodetic coverage: trench-breaking A-B (a-c) versus non-trench-breaking B-C type events (d-f).The gray and magenta markers represent the vGNSSs coverage of onshore and offshore, separately.

Figure 4 .
Figure 4. Sensitivity to offshore geodetic coverage of the locations of the analog earthquake centroid (Cd, 50th percentile of the cumulative slip distribution; see Figures S8, S9, and S13), its downdip (DD, ~10th percentile), and updip (UD, ~90th percentile) rupture limits for A-B (a) and B-C (b) type events and of the estimated seismic moments of both types of events (c).Shown are mean values (solid line) and 1 sigma standard-deviations (shaded area) for all events analyzed in this study.

Figure 3 .
Figure 3. Examples of slip distributions across and along strike of the subduction zone and their sensitivity to offshore geodetic coverage: trench-normal transect (a, b) and trench-parallel transects (c, d); trench-breaking (A-B) type event versus non-trench-breaking (B-C) type event.