Progress in Earth andPlanetary Science

Kubota et al.
Progress in Earth and Planetary Science(2022) 9:68
https://doi.org/10.1186/s40645-022-00524-0

Open Access

RESEARCH ARTICLE

A new mechanical perspective
on a shallow megathrust near‑trenchslip from the high‑resolution fault modelof the 2011 Tohoku‑Oki earthquake
Tatsuya Kubota1* , Tatsuhiko Saito1   and Ryota Hino2   

Abstract
The 2011 Tohoku-Oki earthquake generated a surprisingly large near-trench slip, and earth scientists have devoted
significant attention to understanding why. Some studies proposed special rupture mechanisms, such as extensive
dynamic frictional weakening; others simulated this near-trench slip behavior without supposing the extensive
dynamic weakening. However, we have not reached a decisive conclusion for this question due to limited spatial
near-trench slip resolution. Hence, in this study we use new tsunami data recorded just above the large slip area in
addition to offshore and onshore geodetic data to improve the spatial resolution of stress release in the Tohoku-Oki
earthquake and quantitatively examine the mechanical state of the plate interface. A maximum slip of 53 m reaching
the trench and an insignificant stress drop (< 3 MPa) at the shallowest portion of the fault were estimated. Based on
our modeling results and the past experimental studies, it is suggested that friction at the shallow near-trench portion
should be inherently low both before and during the earthquake. This result provides perspectives on the shallow slip
behavior along the plate boundary, in which the strain energy accumulation at the deep portion of the fault accounts
for the anomalous large shallow slip, but shallow mechanical coupling does not. A large shallow slip has been considered as a result of the release of sufficiently large strain energy at the shallow portion of the plate interface, but wesuggest that shallow slips similar to that during the 2011 Tohoku-Oki earthquake may occur in any subduction zones
where the energy sufficiently accumulates only in the deeper portion.
Keywords: The 2011 Tohoku-Oki earthquake, Ocean-bottom pressure gauge, Tsunami, Stress drop, Frictional
strength of megathrust, Plate mechanical coupling
1 Introduction
The devastating tsunami generated by the 2011 TohokuOki earthquake caused severe damage to the coastof Japan. The most surprising feature of this megathrust earthquake is the occurrence of a very large slip(> 50 m) reaching the trench (Kodaira et al. 2012; 2020;
Sun et al. 2017). Before this earthquake, it was widely*Correspondence: kubotatsu@bosai.go.jp
1
National Research Institute for Earth Science and Disaster Resilience, 3‑1
Tennodai, Tsukuba, Ibaraki 305‑0006, Japan
Full list of author information is available at the end of the article

believed that the interface between the two plates at theshallow portion, corresponding to a depth range shallower than ~ 10 km from the seafloor, was stably sliding(creeping) during the interseismic period, resulting in nostress accumulation (i.e., no mechanical coupling) (Byrne
et al. 1988; Scholtz 1998; Bilek and Lay 2002). Only the
deep portion of the plate interface (deeper than ~ 10 km)
was believed to be locked to accumulate the stress (i.e.,mechanically coupled), which exhibited unstable stick–slip (i.e., earthquake) behavior (Byrne et al. 1988; Scholtz
1998; Bilek and Lay 2002). The stable sliding behavior at
the shallow portion was also considered to be consistent© The Author(s) 2022, corrected publication [2022]. Open Access This article is licensed under a Creative Commons Attribution 4.0
International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you
give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes
were made. The images or other third party material in this article are included in the article’s Creative Commons licence, unless indicatedotherwise in a credit line to the material. If material is not included in the article’s Creative Commons licence and your intended use is not
permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To
view a copy of this licence, visit http://​creat​iveco​mmons.​org/​licen​ses/​by/4.​0/.

Kubota et al. Progress in Earth and Planetary Science(2022) 9:68

with almost no interseismic slip deficit near the trench
estimated before the Tohoku-Oki earthquake (Suwa et al.
2006; Hashimoto et al. 2009; Loveless and Meade 2010;2011) that was unlikely to have caused the large coseismic slip. To better understand shallow slip behavior andtsunami generation in subduction zones, it is crucial toclarify the reason for this unexpected extremely large
slip occurring during the Tohoku-Oki earthquake. This
study attempts to address this scientific issue, particularly
focusing on the perspective of fault mechanics (i.e., stressaccumulation and release), whereas many previous studies have discussed from the kinematic perspective (i.e.,
slip deficit and coseismic slip). Although the kinematicand mechanical perspectives are closely related to eachother, we should note that these two concepts are essentially different but are often confusingly used to discuss
the interseismic process (Wang and Dixon 2004).
To explain this shallow slip behavior, some studies haveconsidered that the frictional strength along the shallow
plate interface is high enough to accumulate shear stressduring the interseismic period. It has been proposed thatwhen the Tohoku-Oki earthquake occurred, the shallow
frictional stress extensively reduced with an increase in
the slip rate (dynamic weakening) (Di Toro et al. 2011;
Noda and Lapusta 2013), leading to a large shallow
slip. An example of the dynamic frictional weakening
mechanism related to frictional heating along the plateinterface, called “thermal pressurization,” has also been
proposed to explain the shallow slip behavior (Hirono
et al. 2019; Shibazaki et al. 2019). The extensive reduction in frictional strength basically requires significant
shear stress reduction at the shallow portion of the plateinterface.
Other studies have attempted to explain the shallow
slip behavior without considering these special mechanisms requiring a shallow large stress drop (e.g., Layet al. 2012). The basic standpoint of the plate coupling
and rupture mechanics in these studies is that the shallow portion of the plate boundary does not accumulate
the shear stress but only the deep portion does, as has
been conventionally thought (Sholtz 1998; Bilek and Lay2002). From the observational facts, a seafloor drilling
survey (Chester et al. 2013) showed the existence of alow-friction material along the Tohoku shallow fault zone
(Fulton et al. 2013; Ujiie et al. 2013). A laboratory experiment of this shallow fault-zone material retrieved fromthe Tohoku subduction zone, under a water-dampened
condition mimicking a plate boundary suggested thatthe low-friction material was insensitive to the slip rate(i.e., no extreme dynamic weakening behavior) (Remitti
et al. 2015). Fukuyama and Hok (2015) conducted the
dynamic earthquake rupture simulation of an interplate
earthquake mimicking the Tohoku-Oki earthquake,

Page 2 of 19

who assumed a large stress release region only at the
deep portion of the plate interface corresponding to astrong mechanical plate locking (i.e., no stress release at
the shallow portion). This dynamic rupture simulation
showed the surface-reaching slip can occur without stressdrop at the shallow portion, due to the stress-free boundary condition at the ground surface. In addition, basedon fault mechanics, some recent modeling studies havepointed out that, a large shallow slip can occur becauseof the “pinning” effect of the deep frictionally locked area
which causes slip delay (or slip deficit) in the shallow portion (Almeida et al. 2018; Herman and Govers 2020; Herman et al. 2018; Lindsey et al. 2021).
It is difficult to identify the shallow stress release process in detail because of the large uncertainty in the neartrench slip, owing to the lack of a dataset with sufficientquality to resolve it. Although seafloor geodetic data
(Fujiwara et al. 2011; Kido et al. 2011; Sato et al. 2011)
have made significant contributions toward elucidating
the slip distribution (Iinuma et al. 2012; Sun et al. 2017),
the spatial resolution of the slip inferred from such data,
which contains information that is only relevant at thelocation of the observation point, is limited. The onshore
seismic and geodetic datasets recorded more than
200 km from the trench axis also have difficulty in resolving the near-trench slip (e.g., Loveless and Meade 2015).
The analyses using the seismograms also require theassumption of the rupture propagation velocity acrossthe fault, but this assumption may cause an uncertainty
in estimating the spatial extent of the entire slip region
because of a substantial trade-off between the rupture
velocity and the spatial extent (Kubota et al. 2018). The
distributions of the coseismic slip and stress release proposed in the past vary from model to model, particularly
at the shallowest portion (Brown et al. 2015; Sun et al.
2017; Lay 2018), due to the issues raised above. Therefore, the cause of the large shallow slip has been investigated for the past ten years, but no decisive conclusion
has been reached.
In this study, we attempt to constrain the slip and stressdrop distributions with the highest precision to reveal
the cause of the near-trench large slip during the 2011
Tohoku-Oki earthquake and to understand the mechanics of the shallow trench-reaching rupture. To achievethis, we use the tsunami data recorded by several oceanbottom pressure gauges installed directly above the fault
area (Kubota et al. 2021), which have not before been utilized in a slip inversion for the event (Fig. 1). The tsunami
data are advantageous to overcome this issue because
shallow near-trench slips excite tsunamis more efficiently
than the deeper slip and hence the tsunami waveform
data contain unique and robust spatial information about

Kubota et al. Progress in Earth and Planetary Science(2022) 9:68

Page 3 of 19

Fig. 1 a Station locations for the tsunami observation. Inverted triangles, squares, and triangles are the ocean-bottom pressure gauges, theoffshore GPS buoys, and near-shore seafloor wave gauges, respectively. The colors of the symbols indicate agencies owning each dataset. b
Close-up around the Off-Tohoku region. The white star denotes the epicenter of the Tohoku-Oki earthquake (Suzuki et al. 2012). c The location of
the onshore and offshore geodetic stations used in the present study. Black and red squares denote inland GPS and seafloor geodetic observatories,
respectively. Small triangles represent the configuration of the triangular subfaults used for the analysis

the distribution of the fault slip, in particular at the shallow portion near the trench.

2 Data and methods
2.1 Data

We use the tsunami data from ocean-bottom pressuregauges installed just above the large slip area (near-fielddata), obtained by Tohoku University (Kubota et al. 2021)(green inverted triangles in Fig. 1a and b). We also use
the ocean-bottom pressure gauge data of the Universityof Tokyo (Maeda et al. 2011) off Iwate (dark red inverted
triangles in Fig. 1b). The station information is summarized in Table 1.
We also use the offshore tsunami data obtained far
from the source area (far-field data) from the oceanbottom pressure gauges of the Japan Agency for MarineEarth Science and Technology (JAMSTEC, light blue
inverted triangles in Fig. 1a), Japan Meteorological
Agency (JMA, red inverted triangles), National ResearchInstitute for Earth Science and Disaster Resilience (NIED,
black inverted triangles), and National Oceanic andAtmospheric Administration (NOAA)’s DART (Deepocean Assessment and Reporting of Tsunamis) system(a blue inverted triangle). We also use the tsunami waveforms from the GPS buoys (yellow squares in Fig. 1a andb) and the wave gauges (orange triangles) of NOWPHAS(Nationwide Ocean Wave information network for Portsand HarbourS) of the Port and Airport Research Institute(PARI).
In the data processing of the tsunami data, weremove the tidal fluctuation using the theoretical model
NAO.99Jb (Matsumoto et al. 2000). We then apply a lowpass filter with a cutoff period of 100 s to the near-fieldocean-bottom pressure records and a band-pass filter
with a passband of 100–3600 s to the other far-field tsunami data.
In addition to these tsunami dataset, we also use offshore geodetic observation data (Kido et al. 2011; Satoet al. 2011) and onshore GPS data obtained by the Geospatial Information Authority of Japan (GSI) (Fig. 1c).
2.2 Methods
2.2.1 Fault geometry and crustal deformation calculation

In contrast to most of the past tsunami modeling studies using the simplified fault configuration (e.g., Satake
et al. 2013; Yamazaki et al. 2018), such as implementing a planar fault geometry, large-sized and discontinuous subfaults, and/or a buried fault with a top that does
not reach the seafloor, we incorporate the nonplanar
fault geometry of the Japan Integrated Velocity Structure Model (JIVSM) (Koketsu et al. 2012) for the analysis (Fig. 2a). This is because we intend to obtain thecorrect slip and stress drop distribution by suppressing spurious slip and stress drop and artificial stressKubota et al. Progress in Earth and Planetary Science(2022) 9:68

Page 4 of 19

Table 1 List of tsunami stations used in this study
StationLatitude [°N]Longitude [°E]

Depth [m]Inversion time
window [s]AgencySampling rate oforiginal data ­[s]aTM2

39.2459

142.4526

997

0–1800

ERI

0.1

TM1

39.2283142.7720

1618

0–1800

ERI

0.1

P06

38.6340

142.58381254

0–3600

Tohoku University

1

P02

38.5002142.5016

1104

0–3600

Tohoku University

1

P03

38.1834

142.39981052

0–3600

Tohoku University

1

P07

38.0016142.4495

1059

0–3600

Tohoku University

1

P08

38.2829142.83201418

0–3600

Tohoku University

1

P09

38.2650

143.0002

1556

0–3600

Tohoku University

1

GJT3

38.2945

143.48143293

0–3600

Tohoku University

1

21,41838.7180148.6980

5500

1200–4200

NOAA

15

KPG2

42.2365

144.8454

2210

600–3600

JAMSTEC

0.1

KPG1

41.7040

144.4375

2218

600–3600

JAMSTEC

0.1

KCTD

41.6675

144.3409

2540

600–3600

JAMSTEC

10

NMS09

42.3692145.9167

3316

600–3600

Tohoku University

1

NMS05

42.1667145.82354548

600–3600

Tohoku University

1

BOSO2

34.7550140.7517

2098

600–4200

JMA

1

BOSO3

34.8050140.5067

1912

600–4200

JMA

1

HPG1

35.0031139.2247

1176

1800–5400

JAMSTEC

1

VCM3

35.0712139.3906

1225

1800–5400

NIED

0.1

VCM1

34.5954139.91982125

1800–5400

NIED

0.1

807

40.1167

142.0667

125

0–3600NOWPHAS

5804

39.6272

142.1867200

0–3600NOWPHAS580239.2586

142.0969

204

0–3600NOWPHAS

5803

38.8578141.8944160

0–3600NOWPHAS5

801

38.2325141.6836144

0–3600NOWPHAS

5806

36.9714141.1856137

0–3600NOWPHAS

5613

42.9106

144.3972

50.0

1800–5400

NOWPHAS

5602

42.5439141.4458

50.7

1800–5400

NOWPHAS5

202

40.9250

141.4242

43.8

1800–4800

NOWPHAS5203

40.5608

141.5683

27.7

Not Used

NOWPHAS

5219

40.2178141.8600

49.5

Not Used

NOWPHAS5

205

38.2500

141.0661

21.3

Not Used

NOWPHAS

5aAll observed records were resampled to 1 s in the inversion analyses

discontinuity around the boundaries of the subfaults(e.g., Hu and Wang 2008; Ma 2012; Brown et al. 2015;
Wang et al. 2018; Du et al. 2021). The detail is describedin Texts S1 and S2 in Additional file 1.
To derive the fault slip distribution through Green’s
functions, the nonplanar fault plane is divided into
small triangular subfault elements, in which the lengthof one side of the triangle is approximately 10 km (triangles in Figs. 1 and 2). In this fault triangulation, thelongitude and latitude of the subfault vertices are identical to those of Iinuma et al. (2016) but our triangulation used the plate depths from JIVSM (Koketsu et al.
2012), while those of Nakajima and Hasegawa (2006)
were used in Iinuma et al. (2016). The total number of

the triangular subfaults used in this study is Nsub = 434.
Using these triangular subfault elements, we calculate
the seafloor displacement using a uniform half-spacestructure model (Meade 2007). To correctly simulate
the trench-reaching rupture (Wang et al. 2018) in thehalf-space, we slightly modify the fault configuration
as follows. First, the depths of all the triangular vertices are systematically moved vertically upward by
8 km to discard the seawater layer. Then, the depthsof the uppermost row of triangular subfaults (i.e., thecomputational free surface) are set to 0 km from thefree surface so that the shallowest fault surface lines up
with the trench axis. We show the strike and dip angles
derived from the JIVSM in Fig. 2b and c, which areKubota et al. Progress in Earth and Planetary Science(2022) 9:68

Page 5 of 19

Fig. 2 Configuration of the subfaults used in this study. a Center depths of the subfaults measured from the sea surface. Depth contour of theplate boundary of JIVSM (Koketsu et al. 2012) is also shown, with an interval of 5 km for thin lines and 10 km for thick lines. b Strike angles of thetriangular subfaults. Strike direction is shown by small arrows. c Dip angles of the triangular subfaults. Dipping direction is shown by small arrows. d
Rake angles of the triangular subfaults. Direction of the plate convergence (DeMets et al. 2010) is also shown by thick arrows

2.2.2 Observation equation to be solved

Fig. 3 Schematic illustration of the subfault geometry. The unit
fault-normal vector n (blue arrow) is derived from the locations of
the subfault vertices (gray dots). The unit slip vector s (red arrow) is
derived from the plate convergence vector v (thick black arrow) andthe dot product of the fault-normal and the slip vectors ( s • n = 0)

calculated from a vector normal to the subfault surface (i.e., fault-normal vector) n (e.g., Aki and Richards(2002); the schematic illustration is shown in Fig. 3).
The rake angles at each subfault are defined basedon the direction of the plate convergence (azimuth of
290°, DeMets et al. 2010, Fig. 2d), so that the horizontal component of the slip vector ( s ) has an opposite
direction to the relative plate convergence direction(Fig. 3).

We here describe the procedure to estimate the slip andstress drop distributions along the plate interface (Fig. 4).
The procedure is similar to that used in our previous
study (Kubota et al. 2021) (refer to Text S1 in Additional
file 1 for more detailed information). Considering thatthe observed data are expressed by the linear superposition of Green’s function, we estimate the amount of slip
at each of the triangular subfaults by solving the following observation equation: 
wtsun Gtsun
wtsun dtsun
 wgeod dgeod   wgeod Ggeod 
=
m.

(1)
0
αS
0
βE

The vector d is the data vector, consisting of the observed
data, and G is the matrix consisting of Green’s functions(the synthetics by the unit slip of the subfault). The subscripts denote the types of datasets. The scalar value w
denotes the weight of each dataset (see Sect. 2.2.4 for the
choice of the weight value). Vector m consists of the slip
amount of the triangular subfaults, which is to be estimated. To stabilize the solution, we impose the smoothing constraint (Maerten et al. 2005) expressed by matrix
S and the damping constraint using the identity matrix
E. Parameters α and β are the weights of each constraint.
We set α = 2 and β = 0.1 for the main results, but weKubota et al. Progress in Earth and Planetary Science(2022) 9:68

Page 6 of 19

Fig. 4 Slip and stress distributions of the 2011 Tohoku-Oki earthquake. a Slip distribution obtained in the present study. Contour lines are drawn at
10 m intervals. b Stress drop distribution calculated from the slip distribution. Contour lines are drawn at 5 MPa intervals. Green dashed lines denotethe iso-depth contour of the plate interface at 10 km intervals

will discuss its appropriateness in Discussion section(Sect. 4.1.1).
2.2.3 Calculation of the Green’s Functions

We present the procedure to calculate Green’s functions
for tsunamis (matrix Gtsun in Eq. (1)) and geodetic data
(Ggeod) (See Additional file 1: Text S2 for more detail).
First, the seafloor vertical displacement distribution
from each of the triangular subfaults is simulated using
the elastic dislocation theory with the half-space (Meade
2007). The slip direction of each subfault is confined tobe the direction of the slip vector schematically shownby red arrow in Fig. 3 (i.e., rake-parallel). We also incorporate the effect of the apparent seafloor vertical movement due to the horizontally moving seafloor (Tanioka
and Satake 1996). We then calculate the sea-surface
height change from the seafloor deformation using a spatial smoothing filter related to the seawater as the initialcondition of the tsunami simulation (e.g., Kajiura 1963;
Saito 2019). This filter cuts off the small-wavenumberor spatially short-wavelength component, leading toa smoother initial sea-surface height distribution. We
assume a constant depth of H0 = 6 km for this filter, basedon the average water depth around the main rupture

region. Finally, we simulate tsunamis using the linear dispersive equation (e.g., Saito 2019), by numerically solving the equation on the discretized staggered grid. The
cosine-shaped source time function with the duration of
Tr = 40 s is assumed for the rupture time history of eachsubfault, as used in Kubota et al. (2021). In this simulation, the bathymetry data of GEBCO 2020 are interpolated to 2 km spatial intervals. The time interval for thissimulation is 1 s. In the calculation of the ocean-bottompressure gauge waveforms, we incorporate not only thepressure change due to the tsunami (sea-surface heightchange) but also the static pressure changes related tothe seafloor permanent displacement (Tsushima et al.
2012). We finally apply the same low-pass or band-passfilters to the simulated waveforms, as those applied to theobserved data.
Green’s functions for the geodetic data (Ggeod) are also
calculated from each of the triangular subfaults, in a similar way to the calculation of the seafloor displacement(Meade 2007).
2.2.4 Solving the Observation Equation

In the linear inversion analysis, we consider the spatial
and temporal evolution of the rupture by distributing theKubota et al. Progress in Earth and Planetary Science(2022) 9:68

Page 7 of 19

Fig. 5 Temporal evolution of the slip. a Configuration of the temporal elements of the inversion analysis. b Slip amount of the subfaults in eachtemporal element. b Temporal evolution of the total slip amount, calculated by the summation of the temporal elements

Green’s function in the time domain in addition to thespace domain (Fig. 5a). We distribute the Green’s functions in the time domain for each subfault with the temporal interval of Δt = 0.5 × Tr = 20 s. The slip of the k-th
temporal element begins at t = tkbeg = (k − 1) × Δt andends at t = tkend = tkbeg + Tr. We assume Nt = 9 Green’s
functions in the time domain for each subfault. The total
number of unknown parameters, which compose thevector m in Eq. (1), is N = Nsub × Nt = 3906. Considering the rupture front propagation, the slips of the k-th
temporal element at the i-th subfault are constrained
to be zero when the rupture front does not arrive there(i.e., we allow the slip only when satisfying the condition
Vr × tkend ≥ ai, where ai is the distance between the hypocenter and the center of the i-th subfault and Vr = 4 km/s
is the rupture velocity).
By solving the least squares problem expressed as
Eq. (1), using the nonnegative least squares inversion
method (Lawson and Hanson 1995), we estimate the spatiotemporal evolution of the slip of the subfaults (vector
m in Eq. (1)). The time windows used for this inversion
analysis are shown in the blue traces in Fig. 6a–d. We use
the time series data with 1-s sampling (i.e., approximately
1800 to 3600 samples per station), and hence, the total
number of the tsunami data sample (dtsun) is ~ ­105. The
offshore and onshore geodetic stations used for the inversion are shown in Fig. 1c, which have three components
per station (two horizontal and one vertical), resulting
in the total geodetic data samples (dgeod) of ~ ­103. We set

the relative weight of the tsunami and geodetic dataset
as wtsun = 1 and wgeod = 10, although the actual relative
weight is much larger in the tsunami data than the geodetic data, due to the difference of the total data samples.
The weights of the spatial smoothing and damping constraints (parameters α and β) are set as α = 2 and β = 0.1,
based on the previous modeling results (Kubota et al.
2021) and based on additional inversion analyses assuming various α and β values (discussed later, in Sect. 4).
2.2.5 Variance Reduction

To evaluate the goodness of the waveform fitting between
the observation and simulation, we introduce the variance reduction (VR) expressed as follows:�2 � � obs
cal
− di
i di


VR = 1 − × 100(%),(2)
� � obs �2diiwhere diobs is the i-th data of the observed data (vector
dtsun in Eq. (1)) and dical is the calculated one. To calculate
the VR, the same time window as that used for the inversion analysis is used (marked by blue lines in Fig. 6a–d).
2.2.6 Calculation of Stress Drop Distribution

After estimating the slip distribution of the triangularsubfaults, we calculate the distribution of the shear stresschange along the fault (i.e., stress drop, Fig. 4b) by computing the shear stress change along the slip direction atKubota et al. Progress in Earth and Planetary Science(2022) 9:68

Page 8 of 19

Fig. 6 Comparison between the observed and synthetic data. a–d Comparisons of the seafloor pressure gauge waveforms above the focal area,
the far-field seafloor pressure gauge waveforms, the offshore GPS buoy waveforms, and the nearshore wave gauge waveforms near the coast. Grayand red waveforms are the observed and synthetic waveforms, respectively. Blue lines indicate the time window used for the analysis. Locations ofeach station are shown in Fig. 1. e A comparison of the observed and synthetic data of the geodetic dataset. Black and gray bars and arrows denotethe observed movement of the offshore and onshore geodetic stations. Light and dark red bars and arrows are the synthetic ones. Note that thelength of the bars and arrows differs between the onshore and offshore stations. The distribution of the seafloor vertical displacement calculated
from the fault model is also shown

the center of each subfault (Saito and Noda 2022). The
shear stress change �σi (xi ), at the i-th triangular subfault, in which its center is located at xi, is expressed by asuperposition of the shear stress change contributions by
the fault slip of all the subfaults:�σi (xi ) =�S xi ; xj mj ,(3)
jwhere mj is the slip amount at the j-th subfault (locatedat xj), and ΔS(xi; xj) is the shear stress change at thecenter of the i-th subfault caused by the unit slip of thej-th subfault. We define the direction of the shear stressof each subfault based on the fault-normal and slip vectors shown in Fig. 3. Using the unit fault-normal vectorn(xi ) and the unit slip vector s(xi ) at the i-th subfault,
S xi ; xj is calculated as (e.g., Saito and Noda 2022):�S xi ; xj = s(xi ) · τ xi ; xj n(xi ),(4)where |s(xi )| = 1 and |n(xi )| = 1, and τ xi ; xj is the stresstensor at the center of the i-th subfault element caused by

the unit slip of the j-th subfault. We calculate the stresstensor τ xi ; xj using the half-space media (Meade 2007)
assuming the rigidity of 40 GPa.

3 Results
3.1 Fault slip distribution

We obtain a large slip of 53 m to the trench axis in theregion off Miyagi (approximately at 38.0°N, 144.0°E,
Fig. 4a). The seismic moment is Mo = 5.07 × ­1022 Nm(Mw 9.07, assuming a constant rigidity of 40 GPa). In
Fig. 5b and c, we show the subfault slip amounts ineach temporal element (k = 1, 2, …, Nt). At the initialstage, t ≤ 80 s from the rupture initiation, we observe
relatively large slip around the hypocenter. Then, therupture extends to the trench axis until t = 140 s. At
the later part of the rupture (t > 140 s), we observe therupture at the near-trench area of the northern part,off-Sanriku, Iwate, region (~ 39.5°N). This spatiotemporal evolution of the rupture is consistent with that estimated from the seismogram analyses in the past studies(e.g., Ide et al. 2011). The synthetic tsunami (Fig. 6a–d)

Kubota et al. Progress in Earth and Planetary Science(2022) 9:68

and onshore and offshore displacement (Fig. 6e) fromthis slip distribution model are in good agreement withthe observations (VR = 97.7%). The main slip area,
defined as the region where the slip exceeded 10 m, isconsistent with that reported previously (e.g., Iinuma
et al. 2012; Satake et al. 2013; Yamazaki et al. 2018).
We evaluate the uncertainty of the slip estimation
by additional inversion tests based on a jack-knife
approach (Additional file 1: Figure S1a and red shaded
area in Fig. 7a). In this test, we exclude one of the nine
near-field pressure gauge stations (Fig. 6a) and use theremaining eight stations to estimate the slip distribution. The inversion setting is identical to the original
setting which used all the data. We define the possible
range of the slip amount using the maximum and minimum slip amounts from all tested models. As a result,the possible range of the shallowest slip is between 49
and 55 m at the Off-Miyagi region (red shaded area in
Fig. 7a).
Our model has a peak slip at the trench axis (red line in
Fig. 7a), which is consistent with the study by Sun et al.
(2017) who estimated the near-trench slip profile of the
Off-Miyagi region using the bathymetry change based on
the surveys conducted before and after the Tohoku-Oki
earthquake (blue line in Fig. 7a). In contrast, other previous models (e.g., Yamazaki et al. 2018) locate the peakslip a few tens of kilometers from the trench axis (green
line in Fig. 7a). The gradient of the slip along the dip
direction is also consistent with that of Sun et al. (2017).

Page 9 of 19

The slip amount at the trench axis in our study isslightly smaller (~ 53 m) than that of Sun et al. (2017)(~ 65 m), by ~ 12 m. Therefore, we examine whether theslip profile reaching up to 65 m, estimated by Sun et al.
(2017), can explain the tsunami waveform observedby the near-field pressure gauges (Fig. 8). Because thecoseismic slip estimated by Sun et al. (2017) was available only at the near-trench region, our slip distribution
model is used as a reference to emulate the slip model ofSun et al. (2017) (see Text S3 for the detailed procedure).
As a result, the simulated tsunami waveforms are almost
consistent with the observation (VR = 95.7%), but thepeak amplitudes at the stations near the epicenter (P08and P09) are larger than those observed (Fig. 8b). This
may indicate that the maximum slip near the trench is
not as large as 65 m.
We also examine how the near-field ocean-bottompressure data that we newly used contribute to theimprovement of the modeling resolution, by conducting tsunami simulations using the previous slip distributions of Satake et al. (2013), Iinuma et al. (2012),
and Yamazaki et al. (2018) (Fig. 9, see Text S4 for thedetailed simulation procedure). The model of Iinuma
et al. (2012) has an extremely large slip at the trench(light blue line in Fig. 7a), which produces a shortwavelength tsunami inconsistent with the observations (stations P03, P07). The model of Yamazaki
et al. (2018) has the maximum slip located at ~ 40 km
inboard of the trench (green line), which produces

Fig. 7 Cross-sectional view of the fault slip and stress drop of the Tohoku-Oki earthquake. a Comparison of the fault slip distributions between thepresent study (red) and the ones in the previous studies (Iinuma et al. 2012; Satake et al. 2013; Sun et al. 2017; Yamazaki et al. 2018) along the line
normal to the trench axis shown in the inset map. The error range of the present study’s model is shown by the red shaded area. The blue shaded
area denotes the error range of Sun et al. (2017). The small triangles denote the location of the ocean-bottom pressure gauges. The region of thedifferential bathymetry survey conducted in Sun et al. (2017) is shown by the cyan line. b Comparison of the stress drop distributions. The error
range of the present study’s model is shown by the red shaded area

Kubota et al. Progress in Earth and Planetary Science(2022) 9:68

Page 10 of 19

Fig. 8 Tsunami simulation using a previous slip model from Sun et al. (2017). a Comparison of the along-dip slip profile of the present model (red),the model in which the slip amount is modified (dark blue), and the slip profile estimated by Sun et al. (2017) (blue). The amount of modificationof the slip is indicated by the gray line. The colors, shading, triangles, and bathymetry description are the same as in the legend of Fig. 7a. bComparison between the observed (gray) and simulated waveforms (colored) at representative stations

the larger short-wavelength tsunami than the observation at GJT3. We quantify the goodness of the slip
models based on the VR value (Eq. (2)) using the time
window of 0 to 1800s from the near-field seven pressure gauge waveforms (P02, P03, P06, P07, P08, P09,
and GJT3, hereafter referred to as V
­ Rnear). We obtain­VRnear = 65.4% for the model of Satake et al. (2013),­VRnear = 5.5% for Iinuma et al. (2012), and 35.7% forYamazaki et al. (2018), while our fault model obtains­VRnear = 99.2%. The fault models of these previous studies cannot explain the near-field ocean-bottom pressuregauges waveforms, which indicates that the use of thenear-field tsunami data obtained by the pressure gauges
contributes to revealing the detailed shallow slip profileat the trench.
3.2 Stress drop

We show the shear stress change along the fault plane(i.e., the stress drop) in Fig. 4b. The stress release at the
deep portion (> 10 km) is large (> 5 MPa), where the slip
amount is smaller than ~ 40 m, whereas the stress releaseat the shallowest portion (> ~ 40 m slip) is much smaller.
Similar to the slip uncertainty evaluation, we evaluate the stress drop uncertainty based on the jack-knife
approach (Additional file 1: Figure S1b, red shaded
area in Fig. 7b). We calculate the stress drop distribution based on the slip distribution models estimated by
excluding one of the nine near-field pressure gauges. The
possible range of the stress drop is then defined by the

maximum and minimum slip amounts from all tested
models. Taking the uncertainty of the estimation into
account, the stress drop at the shallow portion is smaller
than 3 MPa (Additional file 1: Figure S1b and 7b). Considerable stress drop at the deeper portion suggests thata strong mechanical coupling at the deeper portion accumulated the significant amount of shear stress before theearthquake, while insignificant coseismic stress releaseat the shallow portion suggests much weaker shallow
mechanical coupling than the deep portion before theearthquake (discussed later in Sect. 4.2).

4 Discussion
4.1 Evaluation of estimated fault model
4.1.1 Constraint weights in the inversion analysis

In our inversion analysis, we use the weights of thesmoothing and damping constraints of α = 2 and β = 0.1.
However, as the slip distribution may be sensitive to theseparameters, shear stress change calculation from the slip
distribution will be also highly sensitive to the choicesof these parameters. Therefore, we conduct additional
inversions assuming various constraint weight parameters (Additional file 1: Figures S2–S10).
Our tests of smoothing and damping parameters (Additional file 1: Figures S2–S10) show that the inversions
are more sensitive to the smoothing weight α than thedamping weight β. When we assume the smaller smoothing constraint (α = 0.2, Additional file 1 Figures S2–S4),
although the observed waveforms used in the inversion

Kubota et al. Progress in Earth and Planetary Science(2022) 9:68

Page 11 of 19

Fig. 9 Tsunami simulation results of previous fault models of the Tohoku earthquake. a–c Slip distribution and initial sea-surface height distribution
assumed from the fault model proposed by Satake et al. (2013), Iinuma et al. (2012), and Yamazaki et al. (2018), respectively. d Comparison between
the observed (gray) and simulated waveforms (colored) at representative stations

are reproduced well (VR > 98%), the calculation fails to
reproduce the later part of the tsunamis observed GPS
buoys not used in the inversion analysis (Additional file 1:
Figures S2e, S3e, and S4e). In addition, the slip distributions significantly include short-wavelength heterogeneities and this leads to the extremely large stress dropwith finer-scale spatial heterogeneity. This overly roughdistributions of slip and stress drop distributions must
be artifacts due to the inappropriately smaller smoothing constraint, considering that the tsunami data, whichhave the dominant period of larger than 100 s, cannot
resolve such small-scale heterogeneities. The heterogeneous slip and stress drop distributions related to thesmaller smoothing constraint should be arisen from themodeling error, mainly due to the overfitting of the data

which contain various noises such as the observational
noises (e.g., Yabuki and Matsu’ura 1992; Fukahata andWright 2008).
When we impose smoothing constraints equal to or
larger than the original value (α = 2 or 20, Additional
file 1: Figures S5–S10), we basically observe the trenchreaching slip. When supposing the stronger smoothing
constraint (α = 20, Additional file 1: Figures S8–S10),
the amount of the near-trench slip is smaller than theoriginal results (~ 40 m) and the impulsive tsunami components at some stations (e.g., TM1, TM2, 21,418, 807,
and 804) could not be reproduced. This means that thesmoothing weight of α = 20 is inappropriately large. Even
though the α = 20 models are inappropriate, we notethat the stress drop distributions calculated from theseKubota et al. Progress in Earth and Planetary Science(2022) 9:68

models are similar to the original model, in which the
maximum stress drop is located at the deeper portion
but no stress release at the shallower portion. Our main
results, that the stress release occurred only at the deeper
portion, do not change in a possible range of the smoothing constraint.
4.1.2 Robustness of shallow stress drop estimation

Our result suggests that the stress release at the shallowest portion corresponding to the near-trench large slip(> ~ 40 m slip) is lower in magnitude than the deep portion (Fig. 4b). Brown et al. (2015) suggested that somepast fault models had larger stress drop in the shallow
portion than that in the deeper portion (> 30 MPa, e.g.,
Yamazaki et al. 2018). We here evaluate the exact location

Page 12 of 19

of the main stress drop and examine how insignificant
the shallow stress is, based on tsunami forward simulations assuming virtual earthquakes.
We assume large stress drop regions at shallow, deep,
and deeper portions of the plate interface, respectively(Fig. 10a–c, shown by thick black lines). We assign astress drop of 5 MPa in each of the portions, based on the
stress drop amount estimated in our inversions (Fig. 4b).
Then, assuming the given stress drop amount as the virtual observed data (left-hand side of Eq. (3)), we estimate
the slip amount of the j-th subfault (mj) by solving the linear inversion problem (Noda et al. 2021; Saito and Noda2022) (contour lines in Fig. 10a–c). Finally, using this virtual fault slip model, we forwardly simulate tsunamis tocompare with the pressure gauge waveforms (Fig. 10d).

Fig. 10 Numerical evaluation of plate coupling area location. a–c Slip distribution assuming coupling at the shallowest, middle, and deeper
portions of the plate boundary, respectively. Polygons drawn with thick black lines denote the regions where the stress drop was assigned. dComparison between the observed (gray) and simulated waveforms (colored) at representative stations

Kubota et al. Progress in Earth and Planetary Science(2022) 9:68

The temporal evolution of the fault is assumed so that all
the slip across the fault starts simultaneously at t = 0 andends after the duration of Tr = 40 s.
As a result, when assuming a large stress drop area atthe deep portion of the plate boundary where the large
stress drop was estimated in our modeling (Fig. 10b), we
obtain the near-trench slip up to ~ 50 m and a moderate spatial slip gradient consistent with our fault model(red line in Fig. 11). This model explains the featuresof the observed tsunamis well (red traces in Fig. 10d).
On the other hand, when assuming a large stress droparea of 5 MPa at the shallowest portion near the trench(Fig. 10a), a maximum slip of > 80 m and a large spatial
gradient of the slip amount are necessary to cause thelarge stress drop at the shallowest portion (blue line in
Fig. 11). This shallow extremely large slip generates verylarge short-wavelength tsunamis but cannot explain theobservation (blue traces in Fig. 10d). The ~ 20 m shallow
slip based on the assumption of a further deeper stressdrop area (Fig. 10c, green line in Fig. 11) cannot explain
the observation as well (green traces in Fig. 10d). Fromthese forward simulations, we conclude that the large
stress drop area should be located around deeper part ofthe plate boundary near the hypocenter.
We also examine how insignificant the stress was
released in the shallow portion, based on another tsunami forward simulation (Additional file 1: Figure S11).
We modify the stress drop distribution of our modelFig. 11 Cross-sectional view of fault slip for the numerical evaluation
test. The slip distribution profiles assuming coupling at the shallowest(blue), middle (red), and deeper (green) portions of the plateboundary, respectively, are shown. The colors, shading, triangles, andbathymetry description are the same as in the legend of Fig. 7aPage 13 of 19(Fig. 4b) so that the shallow stress drop in the Off-Miyagi region (~ 37–39°N) becomes zero (Additional file 1:
Figure S11a) and construct the virtual slip distribution
with the same procedure as the forward simulations
above (Additional file 1: Figure S11b). Using this fault
slip distribution, we then simulate tsunamis (Additional
file 1: Figure S11c). As a result, the simulation reasonably
explains the features of the observed tsunamis well evenif we suppose no stress drop in the shallow portion (darkblue traces in Additional file 1: Figure S11c). This result
supports that the shallow stress drop is much smaller inmagnitude than the deeper portion.
4.1.3 Effects of heterogeneous structure

To check the validity of the use of the JIVSM plate boundary model for the fault modeling, we compare the JIVSM
with a plate boundary model of Iwasaki et al. (2015), constructed from the precise seismicity and active survey
data (Fig. 12a). We also compare the JIVSM with the Vp
structure along the MY102 survey line in the Off-Miyagi
region estimated by an air-gun survey (Miura et al. 2005,
Fig. 12b). There are some differences between the models, but we confirm the plate depth is reasonably similar
to each other at the shallowest portion where the large
near-trench slip occurred. This indicates that the use of
the plate geometry from the JIVSM is reasonable to estimate the realistic stress drop distribution.
It is possible that the heterogeneity of the crustal
structure may affect the calculation of the stress drop(Figs. 12b and c, e.g., Yagi and Fukahata 2011; Fukahata
et al. 2012). Because the rigidity at the shallowest portion of the plate boundary should be smaller than thatassumed in this study, the stress change at the shallow
portion is expected to be smaller than that obtained
from our study. Even if the heterogeneity of the structure
is considered, we would be able to reasonably state thatthe shallow stress drop is much less significant than the
deeper portion. In the next step, we further consider thisheterogeneous crustal structure using the numerical simulation approach (e.g., Ma 2012; Maeda et al. 2017; Sunet al. 2017; Herman et al. 2018; Du et al. 2021). In addition, it will be necessary to consider the seafloor bathymetry, which cannot be incorporated in the half-spacecalculation (Wang et al. 2018). It will be also important toconsider the inelastic deformation effect due to dynamic
coseismic stress change, leading to an enhancement of
the tsunami (e.g., Ma 2012; Ma and Nie 2019; Du et al.
2021).
4.2 Mechanism of near‑trench large slip: kinematicand dynamic perspectives of the fault boundary

The use of the near-field ocean-bottom pressure gauge
made it possible to resolve the shallow slip and stressKubota et al. Progress in Earth and Planetary Science(2022) 9:68

Page 14 of 19

Fig. 12 Comparison of the plate depth and velocity structure of the previous studies. a Comparison of the plate depth of JIVSM (Koketsu et al.
2012, dark red contours) and that of Iwasaki et al. (2015) (dark blue contours). b The vertical profile along the MY102 line (Miura et al. 2005, thickblack line in Fig. 10a). Colors indicate the Vp structure estimated by Miura et al. (2005). The plate depth contours of both JIVSM and Iwasaki et al.
(2015) are also drawn. c The Vp structure of JIVSM along the MY102 line

drop distributions of the Tohoku-Oki earthquake, whichshows the shallow stress drop was insignificantly small.
In addition, the seafloor drilling survey conducted after
the Tohoku-Oki earthquake suggested that the absolute shear stress level after the earthquake at the shallow
plate interface was almost zero (Lin et al. 2013; Brodsky
et al. 2017; 2020). According to the experimental studies using the shallow fault-zone material of the Tohoku
subduction zone (Ujiie et al. 2013; Remitti et al. 2015),

the temperature measurements of the interplate frictional heating due to the earthquake (Fulton et al. 2013),
and the seismic active survey study (Hondori and Park2022), it is suggested that the friction coefficient in theshallow plate interface is considerably low, resulting invery low absolute level of the shear stress. Therefore, it
is suggested that the shear stress at the shallow portion
should be also nearly zero during, and even before, theTohoku-Oki earthquake, and the shallow portion did notKubota et al. Progress in Earth and Planetary Science(2022) 9:68

accumulate the shear stress during the long-term interseismic period (i.e., no shallow mechanical locking).
Based on these results, we propose that the main reason
for the large shallow coseismic slip is a stress release of adeeper locked zone, where the stress accumulates during
the interseismic period. In other words, the Tohoku-Oki
earthquake slip occurred to compensate for the interseismic slip deficit at both deep and shallow parts of the plateinterface, which was provoked only by deep mechanical coupling (Fig. 13) (Almeida et al. 2018; Herman et al.
2018). Before the Tohoku-Oki earthquake, it was widelybelieved that the shallow plate interface does neither
accumulate the shear stress (Scholtz 1998; Bilek and Lay2002) nor the slip deficit (Suwa et al. 2006; Hashimoto
et al. 2009; Loveless and Meade 2010; 2011). Our study
indicates that the shallow stress accumulation behavior
is consistent with that expected before the Tohoku-Oki
earthquake, but the slip deficit accumulation behavior isnot.
Most previous studies have conventionally understoodthe occurrence of large earthquakes based on a kinematic perspective (e.g., Nishikawa et al. 2019; Uchida

Page 15 of 19

and Bürgmann 2021) (Fig. 13b). As has been seen inother large earthquakes, the main slip of the Tohoku-Oki
earthquake was located in the slip deficit area during theinterseismic period (Lindsey et al. 2021), indicating thatthe earthquake released the interseismic slip deficit. Seismic waves were radiated mainly at deeper depths, butthere was almost no radiation at shallow depths (Ide et al.
2011). Afterslips occurred in areas without coseismic slip
(Watanabe et al. 2021). Other typical kinematic picturesare shown in Fig. 13b. To explain the extremely large
shallow slip which was unusual within this kinematic
perspective, some studies considered the possibility ofan additional mechanism causing an extensive dynamic
reduction of friction, such as thermal pressurization,
which results in an extremely large stress drop process
at the shallow part (Hirono et al. 2019; Shibazaki et al.
2019).
Here, based on the data analysis focusing on the stressdrop distribution, we expanded this kinematic view to
a new mechanical picture (Fig. 13a). From a mechanical perspective, the rupture area of the Tohoku-Oki
earthquake can be divided into deeper (> ~ 10 km) and

Fig. 13 Interpretation of the mechanical and kinematic properties along the Tohoku plate boundary. a Mechanical perspective associated withthe megathrust Tohoku-Oki earthquake. The regions where the plate boundary is mechanically locked are shown in red, and the regions where themechanical coupling is weak are shown in blue. The regions surrounded by dashed lines are less certain than those with the solid lines. The regions
where the mechanical property is unclear are shown in gray. b Kinematic perspective of the seismic activities along the plate boundary. The spatial
relationship of the coseismic (magenta) and postseismic (cyan) slips are shown. Aftershock areas are also shown in gray, and tremor areas are shownin greenKubota et al. Progress in Earth and Planetary Science(2022) 9:68

shallower (< ~ 10 km) portions. We propose that thedriving force of the entire slip was the accumulated
strain energy at the deep mechanically coupled area. The
amount of shallow slip seemed incredibly large, but it canbe reasonably interpreted by considering the effect of the
deep stress release and its interaction with the free surface (Almeida et al. 2018; Herman et al. 2018). In addition, a large seismic wave radiation area is well correlatedwith the area of the large stress drop, whereas the shallow
weak seismic wave radiation area corresponds to the low
stress drop area.
Although the kinematic concept of afterslips is simple and complementary to the main shock (Watanabe
et al. 2021), we propose two different mechanisms ofafterslip: The afterslips located just north (~ 39.5°N) andsouth (~ 37°N) of the rupture area were driven to dissipate stress concentration due to the mainshock, while theslip which occurred ~ 200 km south (~ 35.5°N) was notdirectly driven by the stress concentration, but it spontaneously started (i.e., triggered aseismic slip).
4.3 Toward understanding of the Tohoku‑Oki earthquake
and megathrust earthquake physics

In addition to the conventional kinematic perspectiveof the megathrust earthquake, the stress drop distribution provided us with additional information about themegathrust earthquake, including the cause of the driving force that triggered the shallow large slip, the sourceof seismic wave excitation, and the existence of different
types of afterslip generation mechanisms (Fig. 13). Theseobservations are consistent with the basic mechanical model of faulting (Kostrov 1974), in which the strainenergy stored in the lithosphere between the interseismic
period excites the fault slip and seismic wave radiation.
Without assuming any special mechanism requiring an
extremely large shallow stress drop, the anomalous shallow slip can be explained by a combination of free surface and deep stress release (Fukuyama and Hok 2015;
Almeida et al. 2018; Herman et al. 2018).
The large shallow slip of the Tohoku-Oki earthquake
was mainly due to the effect of the free surface and the
deep stress release (Almeida et al. 2018; Herman et al.
2018). This indicates that the shallow slip behavior
depends largely on the location and amount of the stressdrop. More specifically, the earthquake slip behavior
relies on the amount of strain energy accumulated aroundthe mechanically coupled area. The rupture of the Mw 7.2
foreshock on 9 March 2011, which preceded the TohokuOki earthquake by two days, had a much smaller stressdrop of 1.0 MPa on average (Kubota et al. 2017), resulting in a smaller amount of the released strain energy
as compared to the March 11 mainshock. Although we

Page 16 of 19

need further efforts to deepen our understanding of themechanics of the foreshock and the mainshock, the difference in the amounts of the strain energy released bythese earthquakes would be essentially important tounderstand the diversity of the earthquake dynamics.
The past studies using the onshore geodetic data suggested that there were no interplate slip deficits in theshallowest portion of the plate interface (Suwa et al.
2006; Hashimoto et al. 2009; Loveless and Meade 2010;2011), leading to no shallow large seismic slip. However,recent studies suggest interplate slip deficits have been
geodetically detected even in the shallowest portion of
the plate in many subductions (e.g., Loveless and Meade
2015; Noda et al. 2018; Herman and Govers 2020; Lindsey et al. 2021). Our results show that unusually large
shallow slips and giant tsunamis such as those occurring due to the Tohoku-Oki earthquake can occur in any
subduction zones even without a shallow stress release,
if enough strain energy is accumulated to generate earthquakes around the locked portion which is commonly
thought to be at the seismogenic depth (> 10 km). Thisidea requires methods to evaluate whether enough shear
stress or strain energy is accumulated or not around the
deeper locked areas which can generate large coseismic
slips during the anticipated megathrust earthquakes in
the future. For example, some studies assess the potentials of future megathrust earthquakes by evaluating theenergy budget associated with the earthquake faulting
based on the detected mechanically coupled distributions along the plate interface (Noda et al. 2021; Saito
and Noda 2022). Our results showed that observational
earthquake science is steadily progressing from kinematic modeling toward mechanical modeling to achieve
quantitative evaluation.

5 Conclusions
To understand the reason for the large near-trench slip
during the 2011 Tohoku-Oki earthquake, this study estimated the slip and stress drop distributions with the highspatial resolution using the tsunami data recorded byocean-bottom pressure gauges installed above the fault
area, which had not before been used in the past. The
estimated model had a large slip of > 40 m at the shallowest portion (z < 10 km) in the Off-Miyagi region and theslip peaked at 53 m at the Japan Trench. However, thestress release at the shallowest portion was insignificantly
small (< 3 MPa). The main stress drop region (> 5 MPa)
was located at the deep portion (> 10 km) where the slip
amount was smaller than ~ 40 m. The results suggested
the deep mechanical plate locking corresponding to thelarge stress drop provoked the interseismic slip deficit
in both shallow and deep regions of the plate boundary.
Although a large shallow slip had been considered as aKubota et al. Progress in Earth and Planetary Science(2022) 9:68

result of the release of large strain energy in the shallow
portion of the plate boundary in the past, our analyses
provided us with a new mechanical perspective along theplate boundary, in which shallow slips can occur without the shallow energy accumulation but only with theenergy accumulation in the deeper portion.
Abbreviations
DART​: Deep-ocean Assessment and Reporting of Tsunamis; GSI: Geospatial
Information Authority of Japan; JAMSTEC: Japan Agency for Marine-Earth
Science and Technology; JIVSM: Japan Integrated Velocity Structure Model;
JMA: Japan Meteorological Agency; JPL: Jet Propulsion Laboratory; NASA:
National Aeronautics and Space Administration; NIED: National Research Institute for Earth Science and Disaster Resilience; NOAA: National Oceanic andAtmospheric Administration; NOWPHAS: Nationwide Ocean Wave information
network for Ports and HarbourS; PARI: Port and Airport Research Institute.

Supplementary Information
The online version contains supplementary material available at https://​doi.​
org/​10.​1186/​s40645-​022-​00524-0.
Additional file 1. Supplementary Texts and Figures. This file includes
Texts S1 to S4, Figures S1 to S11, and the detailed captions for Datasets S1
to S6.
Additional file 2. Dataset S1. Comma-separated CSV file for the locationsof tsunami stations used in this study.
Additional file 3. Dataset S2. Comma-separated CSV file for the locationsof the geodetic stations used in this study.
Additional file 4. Dataset S3. Text file for slip amounts of each subfault
estimated by the present study.
Additional file 5. Dataset S4. NetCDF grid file for final slip distribution
interpolated from the slip amounts of each subfault.
Additional file 6. Dataset S5. Text file for stress drops of each subfault
calculated from the slip distribution model.
Additional file 7. Dataset S6. NetCDF grid file for stress drop distribution
interpolated from the stress drop amount of each subfault.

Page 17 of 19

Tohoku Earthquake” (https://​www.​jma.​go.​jp/​jma/​kishou/​books/​gizyu​tu/​133/​
gizyu​tu_​133.​html, accessed on 1 September 2022, available only in Japanese).
The NIED pressure gauge data were provided on request. The tsunami dataof the nearshore GPS buoy and wave gauges were downloaded from theNOWPAHS webpage (https://​nowph​as.​mlit.​go.​jp/​pastd​ata/, only available
in Japanese, accessed on 1 September 2022). The coseismic displacements
at offshore geodetic stations (Kido et al. 2011; Sato et al. 2011) are listed inIinuma et al. (2012). The coseismic displacement data at the onshore geodetic
stations were downloaded from the website of Jet Propulsion Laboratory(JPL), National Aeronautics and Space Administration (NASA) (https://​gipsy-​oasis.​jpl.​nasa.​gov/​index.​php?​page=​pppda​ta, accessed on 1 September 2022),
which were originally acquired by the GSI. The Japan Integrated Velocity Structure Model (Koketsu et al. 2012) was downloaded from https://​www.​jishin.​go.​jp/​evalu​ation/​seism​ic_​hazard_​map/​lpshm/​12_​chosh​uki_​dat/ (accessedon 1 September 2022, available only in Japanese). Plate boundary model of
Iwasaki et al. (2015) (Fig. 12) was downloaded from http://​evrrss.​eri.u-​tokyo.​ac.​jp/​datab​ase/​PLATE​model/​PLMDL_​2016/, and the readme file is available at
http://​evrrss.​eri.u-​tokyo.​ac.​jp/​datab​ase/​PLATE​model/​PLMDL_​2016/​README_​
en.​pdf (accessed on 1 September 2022). The two-dimensional Vp structure
model along the MY102 line (Fig. 12, Miura et al. 2005, https://​www.​jamst​ec.​
go.​jp/​obsmcs_​db/e/​survey/​data_​area.​html?​cruise=​KY99-​05) was provided
by JAMSTEC upon request via JAMSTEC Seismic Survey Database (JAMSTEC
2004, https://​doi.​org/​10.​17596/​00020​69; https://​www.​jamst​ec.​go.​jp/​obsmcs_​
db/e/). The GEBCO 2020 bathymetry data were downloaded from https://​www.​gebco.​net/​data_​and_​produ​cts/​histo​rical_​data_​sets/#​gebco_​2020(accessed on 1 September 2022). We used a triangular dislocation element(tde) program (Meade 2007, https://​github.​com/​brend​anjme​ade/​tde) to calculate seafloor deformation (accessed on 1 September 2022). The digital dataof the slip distribution and stress drop estimated by this study are available in
Additional files 2–7: Datasets S1 to S6 and the detailed caption of the dataset
are shown in Additional file 1. The dataset is also available at the external datarepository Zenodo (https://​doi.​org/​10.​5281/​zenodo.​71042​97).

Declarations
Competing interests
The authors declare that they have no competing interest.
Author details
1
National Research Institute for Earth Science and Disaster Resilience, 3‑1 Tennodai, Tsukuba, Ibaraki 305‑0006, Japan. 2 Graduate School of Science, Tohoku
University, 6‑6 Aza‑Aoba, Aramaki, Aoba‑Ku, Sendai, Miyagi 980‑8578, Japan.
Received: 2 June 2022 Accepted: 17 November 2022

Acknowledgements
We thank two anonymous reviewers and the Editor Takeshi Iinuma for theirreviews and insightful comments. The figures in this manuscript were prepared using Generic Mapping Tools (GMT) version 6 (Wessel et al. 2019). We
also thank Editage for the English language review.
Author contributions
TK conducted the analyses and numerical experiments described in thispaper. TS and RH interpreted the results. All authors drafted the manuscript. All
authors read and approved the final manuscript.
Funding
This work was financially supported by JSPS KAKENHI Grant Numbers
JP19H02409 (TK, TS), JP19H05596 (RH), JP19K04021 (TS), JP19K14818 (TK), andJP22K22K14126 (TK).
Availability of data and materials
The ocean-bottom pressure gauge data installed by Tohoku University areavailable in the supplementary data of Kubota et al. (2021), at https://​doi.​
org/​10.​5281/​zenodo.​44203​93. DART tsunami data were downloaded fromhttps://​www.​ngdc.​noaa.​gov/​hazard/​dart/​2011h​onshu_​dart.​html (accessedon 1 September 2022). The tsunami data of JAMSTEC were downloaded fromhttp://​www.​jamst​ec.​go.​jp/​scdc/​top_e.​html (accessed on 1 December 2019
but currently unavailable on 1 September 2022). The pressure data of the
JMA are available in the supplementary data of Technical Report of the Japan
Meteorological Agency Vol. 1 “3 “Report on the 2011 Off the Pacific Coast ofReferences
Aki K, Richards PG (2002) Quantitative seismology, 2nd edn. University Science
Books, Mill Valley, CA
Almeida R, Lindsey EO, Bradley K, Hubbard J, Mallick R, Hill EM (2018) Can the
updip limit of frictional locking on megathrusts be detected geodetically? Quantifying the effect of stress shadows on near-trench coupling.
Geophys Res Lett 45:4754–4763. https://​doi.​org/​10.​1029/​2018G​L0777​85
Bilek SL, Lay T (2002) Tsunami earthquakes possibly widespread manifestationsof frictional conditional stability. Geophys Res Lett 29:1673. https://​doi.​
org/​10.​1029/​2002G​L0152​15
Brodsky EE, Saffer FP, Chester CM, Huffman K, Moore JC, Wu H-Y (2017) The
postearthquake stress state on the Tohoku megathrust as constrained byreanalysis of the JFAST breakout data. Geophys Res Lett 44:8294–8302.
https://​doi.​org/​10.​1002/​2017G​L0740​27
Brodsky EE, Mori JJ, Anderson L, Chester FM, Conin M, Dunham EM, Eguchi N,
Fulton PM, Hino R, Hirose T, Ikari MJ, Ishikawa T, Jeppson T, Kano Y, Kirkpatrick J, Kodaira S, Lin W, Nakamura Y, Rabinowitz HS, Regalla C, Remitti F,
Rowe C, Saffer DM, Saito S, Sample J, Sanada Y, Savage HM, Sun T, Toczko
S, Ujiie K, Wolfson-Schwehr M, Yang T (2020) The state of stress on thefault before, during, and after a major earthquake. Annu Rev Earth PlanetSci 48:49–74. https://​doi.​org/​10.​1146/​annur​ev-​earth-​053018-​060507

Kubota et al. Progress in Earth and Planetary Science(2022) 9:68

Brown L, Wang K, Sun T (2015) Static stress drop in the Mw 9 Tohoku-oki earthquake: heterogeneous distribution and low average value. Geophys Res
Lett 42:10595–10600. https://​doi.​org/​10.​1002/​2015G​L0663​61
Byrne DE, Davis DM, Sykes LR (1988) Loci and maximum size of thrust earthquakes and the mechanics of the shallow region of subduction zones.
Tectonics 7:833–857. https://​doi.​org/​10.​1029/​TC007​i004p​00833
Chester FM, Rowe C, Ujiie K, Kirkpatrick J, Regalla C, Remitti F, Moore, JC, Toy
V, Wolfson-Schwehr M, Bose S, Kameda J, Mori JJ, Brodsky EE, Eguchi
N, Toczo S, Expedition 343 and 343T Scientists (2013) Structure andcomposition of the plate-boundary slip zone for the 2011 Tohoku-Oki
earthquake. Science 342:1208–1211. https://​doi.​org/​10.​1126/​scien​ce.​12437​19
DeMets C, Gordon RG, Argus DF (2010) Geologically current plate motions.
Geophys J Int 181:1–80. https://​doi.​org/​10.​1111/j.​1365-​246X.​2009.​04491.x
Di Toro G, Han R, Hirose T, De Paola N, Nielsen S, Mizoguchi K, Ferri F, Cocco
M, Shimamoto T (2011) Fault lubrication during earthquakes. Nature
471:494–499. https://​doi.​org/​10.​1038/​natur​e09838
Du Y, Ma S, Kubota T, Saito T (2021) Impulsive tsunami and large runup along
the Sanriku coast of Japan produced by an inelastic wedge deformation
model. J Geophys Res Solid Earth 126:e2021JB022098. https://​doi.​org/​10.​1029/​2021J​B0220​98
Fujiwara T, Kodaira S, No T, Kaiho Y, Takahashi N, Kaneda Y (2011) The 2011
Tohoku-Oki earthquake: displacement reaching the trench axis. Science
334:1240. https://​doi.​org/​10.​1126/​scien​ce.​12115​54
Fukahata Y, Wright TJ (2008) A non-linear geodetic data inversion using ABICfor slip distribution on a fault with an unknown dip angle. Geophys J Int
173:353–364. https://​doi.​org/​10.​1111/j.​1365-​246X.​2007.​03713.x
Fukahata Y, Yagi Y, Mitui Y (2012) Absolute strain release in the 2011 Tohoku-oki
Earthquake: waveform inversion and dynamic fault weakening. J Geol
Soc Japan 118:396–409 (in Japanese with English abstract). https://​doi.​
org/​10.​5575/​geosoc.​2012.​0040
Fukuyama E, Hok S (2015) Dynamic overshoot near trench caused by large
asperity break at depth. Pure Appl Geophys 172:2157–2165. https://​doi.​
org/​10.​1007/​s00024-​013-​0745-z
Fulton P, Brodsky E, Kano Y, Mori J, Chester F, Ishikawa T, Harris R, Lin W, Eguchi
N, Toczko S (2013) Low coseismic friction on the Tohoku-Oki fault determined from temperature measurements. Science 342:1214–1217. https://​doi.​org/​10.​1126/​scien​ce.​12436​41
Hashimoto C, Noda A, Sagiya T, Matsu’ura M (2009) Interplate seismogenic
zones along the Kuril-Japan trench inferred from GPS data inversion. Nat
Geosci 2:141–144. https://​doi.​org/​10.​1038/​NGEO4​21
Herman MW, Govers R (2020) Locating fully locked asperities along the SouthAmerica subduction megathrust: a new physical interseismic inversion approach in a Bayesian framework. Geochem Geophys Geosyst
21:e2020GC009063. https://​doi.​org/​10.​1029/​2020G​C0090​63
Herman MW, Furlong KP, Govers R (2018) The accumulation of slip deficit insubduction zones in the absence of mechanical coupling: implications
for the behavior of megathrust earthquakes. J Geophys Res Solid Earth
123:8260–8278. https://​doi.​org/​10.​1029/​2018J​B0163​36
Hirono T, Tsuda K, Kaneki S (2019) Role of weak materials in earthquake rupture
dynamics. Sci Rep 9:6604. https://​doi.​org/​10.​1038/​s41598-​019-​43118-5
Hondori EJ, Park JO (2022) Connection between high pore-fluid pressureand frictional instability at tsunamigenic plate boundary fault of 2011
Tohoku-Oki earthquake. Sci Rep 12:12556. https://​doi.​org/​10.​1038/​
s41598-​022-​16578-5
Hu Y, Wang K (2008) Coseismic strengthening of the shallow portion of
the subduction fault and its effects on wedge taper. J Geophys Res
113:B12411. https://​doi.​org/​10.​1029/​2008J​B0057​24
Ide S, Baltay A, Beroza GC (2011) Shallow dynamic overshoot and energetic
deep rupture in the 2011 Mw 9.0 Tohoku-Oki earthquake. Science
332:1426–1429. https://​doi.​org/​10.​1126/​scien​ce.​12070​20
Iinuma T, Hino R, Kido M, Inazu D, Osada Y, Ito Y, Ohzono M, Tsushima H, Suzuki
S, Fujimoto H, Miura S (2012) Coseismic slip distribution of the 2011
off the Pacific Coast of Tohoku Earthquake (M9.0) refined by means ofseafloor geodetic data. J Geophys Res 117:B07409. https://​doi.​org/​10.​1029/​2012J​B0091​86
Iinuma T, Hino R, Uchida N, Nakamura W, Kido M, Osada Y, Miura S (2016)
Seafloor observations indicate spatial separation of coseismic andpostseismic slips in the 2011 Tohoku earthquake. Nat Commun 7:13506.
https://​doi.​org/​10.​1038/​ncomm​s13506

Page 18 of 19

Iwasaki T, Sato H, Shinohara M, Ishiyama T, Hashima A (2015) Fundamental
structure model of island arcs and subducted plates in and around
Japan. Abstract T31B-2878 presented at the AGU Fall Meeting 2015, San
Francisco, CA.
Japan Agency for Marine-Earth Science and Technology [JAMSTEC] (2004)
JAMSTEC seismic survey database. https://​doi.​org/​10.​17596/​00020​69
Kajiura K (1963) The leading wave of a tsunami. Bull Earthq Res Inst 41:535–
571. https://​doi.​org/​10.​15083/​00000​33711
Kido M, Osada Y, Fujimoto H, Hino R, Ito Y (2011) Trench-normal variation in
observed seafloor displacements associated with the 2011 Tohoku-Oki
earthquake. Geophys Res Lett 38:L24303. https://​doi.​org/​10.​1029/​2011G​
L0500​57
Kodaira S, No T, Nakamura Y, Fujiwara T, Kaiho Y, Miura S, Takahashi N, Kaneda Y,
Taira A (2012) Coseismic fault rupture at the trench axis during the 2011
Tohoku-oki earthquake. Nat Geosci 5:646–650. https://​doi.​org/​10.​1038/​
ngeo1​547
Kodaira S, Fujiwara T, Fujie G, Nakamura Y, Kanamatsu T (2020) Largecoseismic slip to the trench during the 2011 Tohoku-Oki earthquake.
Annu Rev Earth Plane Sci 48:321–343. https://​doi.​org/​10.​1146/​annur​
ev-​earth-​071719-​055216
Koketsu K, Miyake H, Suzuki H (2012) Japan integrated velocity structure
model version 1. In: Proceedings of the 15th world conference on earthquake engineering. https://​www.​iitk.​ac.​in/​nicee/​wcee/​artic​le/​WCEE2​
012_​1773.​pdf
Kostrov VV (1974) Seismic moment and energy of earthquakes, and seismic
flow of rock. Izv Earth Phys 1:23–40
Kubota T, Hino R, Inazu D, Ito Y, Iinuma T, Ohta Y, Suzuki S, Suzuki K (2017)
Coseismic slip model of offshore moderate interplate earthquakes on
March 9, 2011 in Tohoku using tsunami waveforms. Earth Planet Sci Lett
458:241–251. https://​doi.​org/​10.​1016/j.​epsl.​2016.​10.​047
Kubota T, Saito T, Ito Y, Kaneko Y, Wallace LM, Suzuki S, Hino R, Henrys S(2018) Using tsunami waves reflected at the coast to improve offshore
earthquake source parameters: application to the 2016 Mw 7.1 Te Araroa
earthquake, New Zealand. J Geophys Res Solid Earth 123:8767–8779.
https://​doi.​org/​10.​1029/​2018J​B0158​32
Kubota T, Saito T, Tsushima H, Hino R, Ohta Y, Suzuki S, Inazu D (2021) Extracting near-field seismograms from ocean-bottom pressure gauge inside
the focal area: application to the 2011 Mw 9.1 Tohoku-Oki earthquake.
Geophys Res Lett 48:e2020GL091664. https://​doi.​org/​10.​1029/​2020G​
L0916​64
Lawson CL, Hanson BJ (1995) Solving least squares problems. Soc Ind Appl
Math Phila Pennsylvania. https://​doi.​org/​10.​1137/1.​97816​11971​217
Lay T (2018) A review of the rupture characteristics of the 2011 Tohoku-oki Mw
9.1 earthquake. Tectonophysics 733:4–36. https://​doi.​org/​10.​1016/j.​tecto.​2017.​09.​022
Lay T, Kanamori H, Ammon CJ, Koper KD, Hutko AR, Ye L, Yue H, Rushing TM(2012) Depth-varying rupture properties of subduction zone megathrustfaults. J Geophys Res 117:B04311. https://​doi.​org/​10.​1029/​2011J​B0091​33
Lin W, Conin M, Moore JC, Chester FM, Nakamura Y, Mori JJ, Anderson L,
Brodsky EE, Eguchi N, Expedition 343 Scientists (2013) Stress state in thelargest displacement area of the 2011 Tohoku-Oki earthquake. Science
339:687–690. https://​doi.​org/​10.​1126/​scien​ce.​12293​79
Lindsey EO, Mallick R, Hubbard JA, Bradley KE, Almeida RV, Moore JDP,
Bürgmann R, Hill EM (2021) Slip rate deficit and earthquake potential on
shallow megathrusts. Nat Geosci 14:321–326. https://​doi.​org/​10.​1038/​
s41561-​021-​00736-x
Loveless JP, Meade BJ (2010) Geodetic imaging of plate motions, slip rates, andpartitioning of deformation in Japan. J Geophys Res 115:B02410. https://​doi.​org/​10.​1029/​2008J​B0062​48
Loveless JP, Meade BJ (2011) Spatial correlation of interseismic coupling andcoseismic rupture extent of the 2011 ­MW = 9.0 Tohoku-oki earthquake.
Geophys Res Lett 38:L17306. https://​doi.​org/​10.​1029/​2011G​L0485​61
Loveless JP, Meade BJ (2015) Kinematic barrier constraints on the magnitudesof additional great earthquakes off the east coast of Japan. Seismol Res
Lett 86:202–209. https://​doi.​org/​10.​1785/​02201​40083
Ma S (2012) A self-consistent mechanism for slow dynamic deformation andlarge tsunami generation for earthquakes in the shallow subduction
zone. Geophys Res Lett 39:L11310. https://​doi.​org/​10.​1029/​2012G​L0518​
54

Kubota et al. Progress in Earth and Planetary Science(2022) 9:68

Ma S, Nie S (2019) Dynamic wedge failure and along-arc variations of tsunamigenesis in the Japan Trench margin. Geophys Res Lett 46:8782–8790.
https://​doi.​org/​10.​1029/​2019G​L0831​48
Maeda T, Furumura T, Sakai S, Shinohara M (2011) Significant tsunami observed
at ocean-bottom pressure gauges during the 2011 off the Pacific coast of
Tohoku Earthquake. Earth Planets Space 63:803–808. https://​doi.​org/​10.​5047/​eps.​2011.​06.​005
Maeda T, Takemura S, Furumura T (2017) OpenSWPC: An open-source integrated parallel simulation code for modeling seismic wave propagationin 3D heterogeneous viscoelastic media. Earth Planets Space 69:102.
https://​doi.​org/​10.​1186/​s40623-​017-​0687-2
Maerten F, Resor P, Pollard D, Maerten L (2005) Inverting for slip on threedimensional fault surfaces using angular dislocations. Bull Seismol Soc
Am 95:1654–1665. https://​doi.​org/​10.​1785/​01200​30181
Matsumoto K, Takanezawa T, Ooe M (2000) Ocean tide models developedby assimilating TOPEX/POSEIDON altimeter data into hydrodynamical
model: a global model and a regional model around Japan. J Oceanogr
56:567–581. https://​doi.​org/​10.​1023/A:​10111​57212​596
Meade BJ (2007) Algorithms for the calculation of exact displacements, strains,and stresses for triangular dislocation elements in a uniform elastic half
space. Comput Geosci 33:1064–1075. https://​doi.​org/​10.​1016/j.​cageo.​
2006.​12.​003
Miura S, Takahashi N, Nakanishi A, Tsuru T, Kodaira S, Kaneda Y (2005) Structuralcharacteristics off Miyagi forearc region, the Japan Trench seismogenic
zone, deduced from a wide-angle reflection and refraction study. Tectonophysics 407:165–188. https://​doi.​org/​10.​1016/j.​tecto.​2005.​08.​001
Nakajima J, Hasegawa A (2006) Anomalous low-velocity zone and linearalignment of seismicity along it in the subducted Pacific slab beneath
Kanto, Japan: reactivation of subducted fracture zone? Geophys Res Lett
33:L16309. https://​doi.​org/​10.​1029/​2006G​L0267​73
Nishikawa T, Matsuzawa T, Ohta K, Uchida N, Nishimura T, Ide S (2019) The
slow earthquake spectrum in the Japan Trench illuminated by the S-net
seafloor observatories. Science 365:808–813. https://​doi.​org/​10.​1126/​
scien​ce.​aax56​18
Noda H, Lapusta N (2013) Stable creeping fault segments can become
destructive as a result of dynamic weakening. Nature 493:518–521.
https://​doi.​org/​10.​1038/​natur​e11703
Noda A, Saito T, Fukuyama E (2018) Slip-deficit rate distribution along the
Nankai Trough, southwest Japan, with elastic lithosphere and viscoelastic
asthenosphere. J Geophys Res Solid Earth 123:8125–8142. https://​doi.​
org/​10.​1029/​2018J​B0155​15
Noda A, Saito T, Fukuyama E, Urata Y (2021) Energy-based scenarios for great
thrust-type earthquakes in the Nankai trough subduction zone, southwest Japan, using an interseismic slip-deficit model. J Geophys Res Solid
Earth 126:e2020JB020417. https://​doi.​org/​10.​1029/​2020J​B0204​17
Remitti F, Smith SAF, Mittempergher S, Gualtieri AF, Di Toro G (2015) Frictional
properties of fault zone gouges from the J-FAST drilling project (Mw 9.0
2011 Tohoku-Oki earthquake). Geophys Res Lett 42:2691–2699. https://​doi.​org/​10.​1002/​2015G​L0635​07
Saito T (2019) Tsunami generation and propagation. Springer Japan, Tokyo.
https://​doi.​org/​10.​1007/​978-4-​431-​56850-6
Saito T, Noda A (2022) Mechanically coupled areas on the plate interface in
the Nankai trough, Japan and a possible seismic and aseismic rupture
scenario for megathrust earthquakes. J Geophys Res Solid Earth
127:e2022JB023992. https://​doi.​org/​10.​1029/​2022J​B0239​92
Satake K, Fujii Y, Harada T, Namegaya Y (2013) Time and space distribution of
coseismic slip of the 2011 Tohoku earthquake as inferred from tsunami
waveform data. Bull Seismol Soc Am 103:1473–1492. https://​doi.​org/​10.​
1785/​01201​20122
Sato M, Ishikawa T, Ujihara N, Yoshida S, Fujita M, Mochizuki M, Asada A (2011)
Displacement above the hypocenter of the 2011 Tohoku-Oki earthquake.
Science 332:1395. https://​doi.​org/​10.​1126/​scien​ce.​12074​01
Scholtz CH (1998) Earthquakes and friction laws. Nature 391:37–42. https://​doi.​
org/​10.​1038/​34097
Shibazaki B, Noda H, Ikari MJ (2019) Quasi-dynamic 3D modeling of thegeneration and afterslip of a Tohoku-oki earthquake considering
thermal pressurization and frictional properties of the shallow plateboundary. Pure Appl Geophys 176:3951–3973. https://​doi.​org/​10.​1007/​
s00024-​018-​02089-w

Page 19 of 19

Sun T, Wang K, Fujiwara T, Kodaira S, He J (2017) Large fault slip peaking at
trench in the 2011 Tohoku-oki earthquake. Nat Commun 8:14044. https://​doi.​org/​10.​1038/​ncomm​s14044
Suwa Y, Miura S, Hasegawa A, Sato T, Tachibana K (2006) Interplate coupling
beneath NE Japan inferred from three-dimensional displacement field. J
Geophys Res 111:B04402. https://​doi.​org/​10.​1029/​2004J​B0032​03
Suzuki K, Hino R, Ito Y, Yamamoto Y, Suzuki S, Fujimoto H, Shinohara M, Abe M,
Kawaharada Y, Hasegawa Y, Kaneda Y (2012) Seismicity near the hypocenter of the 2011 off the Pacific coast of Tohoku earthquake deduced byusing ocean bottom seismographic data. Earth Planets Space 64:1125–
1135. https://​doi.​org/​10.​5047/​eps.​2012.​04.​010
Tanioka Y, Satake K (1996) Tsunami generation by horizontal displacement ofocean bottom. Geophys Res Lett 23:861–864. https://​doi.​org/​10.​1029/​
96GL0​0736
Tsushima H, Hino R, Tanioka Y, Imamura F, Fujimoto H (2012) Tsunami waveform inversion incorporating permanent seafloor deformation and its
application to tsunami forecasting. J Geophys Res 117:B03311. https://​doi.​org/​10.​1029/​2011J​B0088​77
Uchida N, Bürgmann R (2021) A decade of lessons learned from the 2011
Tohoku-Oki earthquake. Rev Geophys 59:e2020RG000713. https://​doi.​
org/​10.​1029/​2020R​G0007​13
Ujiie K, Tanaka H, Saito T, Tsutsumi A, Mori J, Toczko S (2013) Low coseismic
shear stress on the Tohoku-Oki megathrust determined from laboratory
experiments. Science 342:1211–1214. https://​doi.​org/​10.​1126/​scien​ce.​12434​85
Wang K, Dixon T (2004) “Coupling” semantics and science in earthquake
research. Eos Trans Am Geophys Union 85:180. https://​doi.​org/​10.​1029/​
2004E​O1800​05
Wang K, Sun T, Brown L, Hino R, Tomita F, Kido M, Iinuma T, Kodaira S, Fujiwara
T (2018) Learning from crustal deformation associated with the M9 2011
Tohoku-oki earthquake. Geosphere 14:2. https://​doi.​org/​10.​1130/​GES01​
531.1
Watanabe S, Ishikawa T, Nakamura Y, Yokota Y (2021) Co- and postseismic
slip behaviors extracted from decadal seafloor geodesy after the 2011
Tohoku-oki earthquake. Earth Planets Space 73:162. https://​doi.​org/​10.​
1186/​s40623-​021-​01487-0
Wessel P, Luis JF, Uieda L, Scharroo R, Wobbe F, Smith WHF, Tian D (2019) The
generic mapping tools version 6. Geochem Geophys Geosys 20:5556–
5564. https://​doi.​org/​10.​1029/​2019G​C0085​15
Yabuki T, Matsu’ura M (1992) Geodetic data inversion using a Bayesian
information criterion for spatial distribution of fault slip. Geophys J Int
109:363–375. https://​doi.​org/​10.​1111/j.​1365-​246X.​1992.​tb001​02.x
Yagi Y, Fukahata Y (2011) Rupture process of the 2011 Tohoku-oki earthquake
and absolute elastic strain release. Geophys Res Lett 38:L19307. https://​doi.​org/​10.​1029/​2011G​L0487​01
Yamazaki Y, Cheung KF, Lay T (2018) A self-consistent fault slip model forthe 2011 Tohoku earthquake and tsunami. J Geophys Res Solid Earth
123:1435–1458. https://​doi.​org/​10.​1002/​2017J​B0147​49

Publisher’s NoteSpringer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.