Can the regional 3D stress field according to the Wallace–Bott hypothesis predict fault slip directions of future large earthquakes?

When evaluating strong ground motions and tsunamis from specified source faults, it is required that the input parameters, such as fault geometry, rake angle, and slip amount, do not deviate from those of a real earthquake. Recently, a regional three-dimensional (3D) tectonic stress field was used to estimate rake angles for mapped sub-marine faults with the Wallace–Bott hypothesis (WBH), the direction of fault slip was parallel to the resolved stress vector on a preexisting fault, and strong ground motions and tsunamis were simulated. However, this modeling technique has not been adequately validated. Additionally, it is necessary to examine how the stress field estimated from seismological data for a limited period (~ 10 years) can be used as a proxy for the long-term tectonic stress field. In this study, to provide such validation, we utilized two catalogs of focal mechanism solutions for earthquakes and compared the observed rake angles with those calculated from the regional 3D tectonic stress field with the WBH by fixing the fault strike and dip angles according to those from the focal mechanism data. The resulting misfit angles between the observed and calculated rake angles are generally small (ranging between − 30° and 30°), excluding several regions (e.g., the source and surrounding regions of the 2011 off the Pacific coast of Tohoku earthquake and swarm-like activities activated after the 2011 quake). We also confirmed that the calculated rake angles and classified fault types are consistent with geomorphologically and geologically evaluated types of faulting for major Quater-nary active faults in the Kyushu district of southwest Japan. These results support the validity and effectiveness of estimating rake angles for a specific fault with known geometry from the above method and data, while also showing that close attention is needed to apply this method to, for example, seismically inactive regions where the inverted stress field includes significant uncertainties and/or near sites of recent and large earthquakes where the stress field has been perturbed.


Introduction
In evaluating strong ground motions and tsunamis from a specified source fault, input parameters such as fault geometry (fault length, fault width, fault strike, and dip angle), rake angle, and slip amount are expected not to deviate from those for a future earthquake.Geometries of faults have been mainly investigated from geomorphological, geological, and seismological studies (e.g., aerial photo interpretation, trenching and coring surveys, seismic reflection and refraction surveys, and gravity anomalies).For expected slip amounts when the fault ruptures from known-geometry faults, empirical scaling relations between fault dimension (e.g., fault length, fault area) and slip amounts have been developed for various types of earthquakes (e.g., Wells and Coppersmith 1994;Takemura et al. 1998;Somerville et al. 1999;Irikura and Miyake 2001;Murotani et al. 2013Murotani et al. , 2015)).In contrast, the rake angle has been conventionally assumed to be a representative value for each fault type (i.e., 90° for reverse-fault type, − 90° for normal fault type, 0° for left-lateral faults, and 180° for right-lateral faults), while setting rake angles is also essential for simulating strong ground motions and tsunamis.In particular, the appropriate setting is an indispensable issue for tsunami hazard assessment (e.g., Annaka et al. 2007;Mulia et al. 2020;Murotani et al. 2022;Satake et al. 2022).
The Japanese Islands are situated under a complicated tectonic setting due to the interaction of four major tectonic plates.The Pacific Plate subducts beneath the Okhotsk Plate along the Kuril and Japan Trenches and subducts beneath the Philippine Sea Plate along the Izu-Bonin Trench.Furthermore, the Philippine Sea Plate subducts beneath the Okhotsk Plate along the Sagami Trough and beneath the Eurasia Plate along the Suruga and Nankai Trough, and Ryukyu Trench.Resulting from this complicated tectonic setting, various types of faulting have occurred, and the resulting tectonic stress field has been investigated over the four past decades (e.g., Yoshii 1979;Huzita 1980;Wesnousky et al. 1982;Nishimura et al. 2004;Imanishi and Kuwahara 2009;Terakawa and Matsu'ura 2010;Matsushita andImanishi 2015, Uchide et al. 2022).Terakawa and Matsu'ura (2008) developed the CMT data inversion method to estimate the three-dimensional (3D) pattern of tectonic stress from the CMT data of seismic events by using Akaike's Bayesian Information Criterion (ABIC; Akaike et al. 1980).The essential difference between the CMT data inversion method and traditional stress inversion methods (e.g., Gephart and Forsyth 1984;Michael 1984Michael , 1987) ) is the use of CMT data, which can be directly related to a tectonic stress field without any knowledge of actual physical processes in a source region and does not directly use the Wallace-Bott hypothesis (WBH; Wallace 1951; Bott 1959) (Fig. 1) to estimate the 3D tectonic stress field.Terakawa and Matsu'ura (2010) (TM2010) applied the CMT data inversion method to more than 12,500 moment tensor solutions (3.0 ≤ M ≤ 5.0) from January 1997 to January 2007 determined by using the Full Range Seismograph Network of Japan (F-net) and estimated the 3D tectonic stress field at depths ranging from 0 to 100 km in and around Japan (Fig. 2).In TM2010, each component of the 3D stress fields is represented by the superposition of 146,848 tricubic B splines with 20 and 10 km equally spaced grid intervals in the horizontal and vertical directions, respectively.This enables us to estimate six components of tectonic stress fields as continuous functions defined in the model region with estimation errors, although only relative values of six components have physical meaning.Figure 2 shows the spatial pattern of tectonic stress fields at a depth of 10 km.The stress pattern is represented with the lower hemisphere projections of focal spheres whose nodal planes indicate maximum shear stress planes.The colors of focal spheres indicate types of faulting according to the classification criteria by Frohlich (1992).TM2010 also verified that the estimated spatial pattern of stress field accorded with the focal mechanisms of large events (M > 5.0), which were not used for the CMT data inversion, and that they were well correlated with the presentday (Quaternary) tectonics of the Japanese islands.
Types and slip directions for inland faults have been mainly evaluated geologically and geomorphologically.Tectonic landforms due to a repetition of faulting include fault scarps and fault flexures for dip-slip faults and offsets for strike-slip faults.However, types and slip directions are not necessarily well known, especially those with low slip rates and/or those in regions with high erosion rates.Furthermore, landforms from strike-slip are usually less distinct than those from dip-slip, and types and slip directions are sometimes insignificant or controversial.Revealing fault types and slip directions for submarine faults generating tsunamis is further arduous due to limited data.
Recently, realistic rake angles for active submarine faults in Japan have been seismologically evaluated from 3D tectonic stress fields with WBH, and those rake angles were utilized for evaluating strong ground motions and tsunamis.For example, Research Committee on Large Earthquakes in the Sea of Japan (2014) and the Integrated Research Project on Seismic and Tsunami Hazards Around the Sea of Japan (Takeda et al. 2014) calculated rake angles expected from the 3D tectonic stress field by TM2010 with WBH, and then strong ground motions and tsunamis were simulated (e.g., Iwata et al. 2018;Satake et al. 2022).
However, the method of estimating rake angles from seismologically estimated stress fields with WBH has two outstanding issues.The first is that it has not been sufficiently validated.Because TM2010 did not directly use WBH, it is necessary to examine whether the tectonic stress field by TM2010 can reproduce the observed rake angles.The second is that the average recurrence interval of a characteristic earthquake in Japan generally ranges from tens to hundreds of years for events on a plate boundary (e.g., Ishibe and Shimazaki 2009) and several thousand years to tens of thousands of years for those on major Quaternary active faults (e.g., Ishibe andShimazaki 2008, 2012).Considering that an earthquake releases the stress accumulated during those periods, it is necessary to examine whether the stress field estimated from seismological data for a limited period (~ 10 years) can be used as a proxy for the long-term tectonic stress fields.
In this study, we mainly focus on the first issue.As the first attempt at such verification, we compare rake angles for focal mechanism solutions of earthquakes from two catalogs: i.e., the F-net mechanism solutions (National Research Institute for Earth Science and Disaster Resilience 2023) and the Japan University Network Earthquake Catalog of First-Motion Focal Mechanisms (JUNEC FM 2 ; Ishibe et al. 2014), with rake angles expected from the 3D tectonic stress field with WBH.We also compare the geomorphologically and geologically evaluated fault types with those estimated from the tectonic stress field with WBH for major Quaternary active faults in the Kyushu district, southwest Japan, including the Futagawa and Hinagu fault zones (FZs) that were ruptured by the 2016 Kumamoto earthquake sequence.

Method and data
As part of this basic method to estimate stress fields, the WBH, in which the fault plane's slip direction is parallel to the direction of shear traction on the fault plane, has been widely utilized (Fig. 1).Given a stress tensor σ, the traction vector at the fault plane whose normal vector is n can be described by the following equation: The normal and tangential components of this vector are the normal and shear traction vectors and are given by: Fault slip occurs to release the shear stress, and hence, the theoretical slip direction can be described by the unit vector t s /|t s | , where |t s | is the length of the shear traction vectors.The rake angles can be obtained Fig. 2 The stress pattern at 10 km in depth after Terakawa and Matsu'ura (2010).The stress pattern is represented by the lower hemisphere projection of focal mechanisms of potential seismic events.The colors indicate types of faulting according to the classification criteria by Frohlich (1992) by taking angles between the fault slip direction and fault strike.The National Research Institute for Earth Science and Disaster Resilience (NIED) began to install a nationwide broadband seismograph network in 1994 under the Fundamental Research on Earthquakes and Earth's Interior Anomaly (FREESIA) project, and after the project terminated in March 2001, the seismograph network was integrated into the network installed as a part of the measure by the Headquarter for Earthquake Research Promotion (HERP) as F-net (Okada et al. 2004).
In the present study, we validate the method of estimating rake angles from the 3D tectonic stress field obtained by TM2010 with WBH for F-net data spanning from January 1997 to December 2021 with the centroid ranging from 125° to 150° E in longitude and 25°-47° N in latitude.To validate the applicability of the method to Quaternary active faults in the shallow crust, we utilize all available focal mechanism solutions of earthquakes with a centroid depth of 30 km or shallower regardless of their variance reductions (VR).The number of validated F-net mechanism solutions is 20,148 with moment magnitudes (M w ) ranging from 3.1 to 8.7.We fix the strike and dip angles of the fault to nodal planes of the focal mechanism solutions and calculate (theoretical) rake angles expected from the 3D tectonic stress field by TM2010 with WBH (hereafter referred to as calculated rake angles).The calculated rake angles are then compared with the actual rake angles of focal mechanisms for earthquakes (hereafter referred to as observed rake angles) by calculating misfit angles.The misfit angles (hereafter, denoted by λ) are defined by the O (observed rake angles by the F-net) − C (calculated rake angles).Here, the smaller absolute misfit angle (hereafter denoted by |λ|) between the first and second nodal planes is adopted as a representative value.We separate the period of the F-net mechanism data into three intervals (Fig. 3): from January 1997 to January 2007 (Period I), from February 2007 to the occurrence of the 2011 off the Pacific coast of Tohoku earthquake (M JMA 9.0 [magnitude determined by the Japan Meteorological Agency (JMA)], M w 8.7 [F-net]; hereafter referred to as the 2011 Tohoku-oki earthquake) (Period II), and the postseismic period of the 2011 Tohoku-oki earthquake (Period III).We set Period I to confirm the applicability of the method because TM2010 does not directly use WBH to invert the 3D tectonic stress field.We further divide the latter period into two periods at the occurrence of the 2011 Tohoku-oki earthquake on March 11 because significant changes in seismicity were reported following the 2011 event (e.g., Hirose et al. 2011;Ishibe et al. 2011aIshibe et al. , 2015Ishibe et al. , 2017;;Toda et al. 2011).The numbers of validated focal mechanisms are 6191, 2623, and 11,334 for Period I, Period II, and Period III, respectively.
We also utilize JUNEC FM 2 (Ishibe et al. 2014; Additional file 1: Fig. S1) for validation.JUNEC FM 2 contains focal mechanism solutions for 14,544 earthquakes that occurred in and around the Japanese Islands from July 1985 to December 1998.This catalog was determined by using first-motion polarities reported by the Japan University Seismic Network and covers small-magnitude earthquakes (M ≥ 2.0) prior to the recent development of seismic observation networks and automated waveform data processing systems.We validate the method for 7221 focal mechanism solutions of earthquakes with qualities of "A", "B", or "C", epicenters ranging from 125° to 150° E in longitude and 25°-47° N in latitude, and depths of 30 km or shallower.All results for JUNEC FM 2 are shown in the supplementary materials (Additional file 1: Fig. S1, Additional file 2: Fig. S2, Additional file 3: Fig. S3, Additional file 4: Fig. S4).
In Japan, the HERP has selected 114 (as of November 2023) major active FZs as basic survey targets, considering the degree of activity and the impact on society when the fault ruptures, to efficiently survey many active faults.HERP also compiled the outcomes of previous surveys, evaluated fault geometry (fault length, width, dip angle, fault strike), fault type, histories of past activities, and average recurrence intervals, and then conducted long-term evaluations (e.g., expected earthquake magnitude, earthquake occurrence probability during the next 30 years) (e.g., HERP 2013).To investigate whether the above geomorphologically and geologically evaluated fault types can be reproduced from the 3D tectonic stress field with WBH, we compare geomorphologically and geologically evaluated fault types with calculated rake angles for major Quaternary active faults in the Kyushu district, southwest Japan.Along the Futagawa and Hinagu FZs in the central part of the Kyushu district, large earthquakes called the 2016 Kumamoto earthquake sequence occurred in April 2016 and caused severe damage near the source region (e.g., Yamada et al. 2017).We use the fault types evaluated by HERP and the fault models obtained from the Japan Seismic Hazard Information Station (J-SHIS) by NIED.We calculate the rake angles at the center of each fault and classify these faults as reverse-fault type (rake angle: 90° ± 45°), normal fault type (rake angle: − 90° ± 45°), left-lateral fault (rake angle: 0° ± 45°), and right-lateral fault (rake angle: 180° ± 45°).These fault types are then compared with geomorphologically and geologically estimated fault types.

Validation of the method for the F-net focal mechanism catalog
Period I (January 1997-January 2007) During Period I, in which the moment tensor solutions of earthquakes were used for the CMT data inversion by TM2010, the calculated rake angles are mostly consistent with the observed rake angles (Figs. 4 and 5).Misfit angles λ mostly range between − 30° and 30° with an average of − 1.81° and a standard deviation of 20.78°.The histogram of λ shows a clear Gaussian distribution with a center of approximately zero.Among 6191 earthquakes, |λ| is ≤ 30° (class A) for 5846 (~ 94.4%) earthquakes, whereas 204 (~ 3.3%), 34 (~ 0.5%) and 107 (~ 1.7%) Earthquakes with large |λ| are concentrated in subduction zones and aftershock regions following large (M w 6.0 or above) earthquakes, such as the M w 7.9 interplate earthquake that occurred along the Kuril Trench on 26th September 2003 in southeastern off Tokachi (called the 2003 Tokachi-oki earthquake by JMA) (No. 13 in Fig. 7a and Table 1), while |λ| is mostly small for the mainshock as discussed in the next paragraph.The close-up figures of the comparison of the calculated rake angles with the observed rake angles in Hokkaido (Additional file 5:  The basic characteristics of faulting [i.e., the reverse faulting dominantly distributed in northeastern Japan (Additional file 9: Fig. S9), the mixture of reverse faulting and strike-slip faulting in central Japan (Additional file 12: Fig. S12), and strike-slip faulting and normal faulting in southwestern Japan (Additional file 15: Fig. S15 and Additional file 18: Fig. S18)] are well reproduced from the 3D tectonic stress field with WBH.
The |λ| values are also small for large (M w 6.0 or above) earthquakes.Among 29 earthquakes with M w 6.0 or above during Period I, |λ| is smaller than 20° for 24 (~ 82.8%) earthquakes and smaller than 30° for 25 (~ 86.2%) earthquakes (Table 1).For example, the calculated rake angle is 82.58° and λ is 4.42° for a reverse-type M w 6.1 earthquake that occurred on 26th July 2003 in northern Miyagi Prefecture (No. 12 in Fig. 7a and Table 1).The observed rake angle of the M w 6.6 earthquake that occurred on 20th March 2005 in the southwestern off Kyushu district (− 177°; named the 2005 Fukuoka-ken Seiho-oki earthquake by JMA) can also be well reproduced by the method (169.49°) with λ = 13.51°(No. 26 in Fig. 7a and Table 1).Furthermore, the M w 7.9 interplate earthquake that occurred along the Kuril Trench on 26th  Frohlich (1992) September 2003 in southeastern off Tokachi (called the 2003 Tokachi-oki earthquake by JMA) is well reproduced by the method with λ = 8.01° (No. 13 in Fig. 7a and Table 1).On the other hand, four earthquakes exhibit |λ| > 30° (Fig. 7a).A typical example is the M w 6.5 earthquake, which occurred on 19th January 2005 in the far southeastern off Boso Peninsula near the triple junction among the Pacific, Philippine, and Okhotsk plates (No. 25 in Fig. 7a and Table 1).The observed focal mechanism solution for this earthquake is a reverse-type with a rake angle of 105°, whereas the calculated one is a left-lateral with a normal-faulting component, and the resulting λ is 142.90°.The possible causes of the relatively large |λ| for these earthquakes are discussed in the next chapter.

Period II (February 2007-the 2011 Tohoku-oki earthquake)
Similarly, the calculated rake angles are mostly consistent with the observed rake angles for Period II (Figs. 8 and  9).|λ| values are also mostly ≤ 30° and the histogram of     Among the 24 earthquakes with M w 6.0 or above that occurred during Period II, |λ| is ≤ 20° for 20 (~ 83.3%) earthquakes and ≤ 30° for 22 (~ 91.7%) earthquakes (Table 2; Fig. 7b).The rake angles for major inland shallow earthquakes such as the 2007 Noto Hanto earthquake (M w 6.7) that occurred on 25th March (λ = − 5.58°, No. 31 in Fig. 7b and Table 2) and the 2008 Iwate-Miyagi Nairiku earthquake (M w 6.9) that occurred on 14th June (λ = 6.28°,No. 40 in Fig. 7b and Table 2) are well reproduced by the 3D tectonic stress field with WBH.Furthermore, observed rake angles for great interplate earthquakes that occurred during Period II, such as the M w 7.2 earthquake that occurred east off Miyagi Prefecture on 9th March, 2 days before the 2011 Tohoku-oki earthquake (λ = − 4.53°, No. 48 in Fig. 7b and Table 2), are also reproduced by the method.On the other hand, only two earthquakes exhibit |λ| > 30°.Both earthquakes occurred near Chichi-jima Island along the Izu-Bonin Trench.One is the M w 6.2 earthquake (No. 37 in Fig. 7b and Table 2) with λ = 41.58°, which occurred on 15th March 2008, and the other is the M w 7.3 earthquake (No. 46 in Fig. 7b and Table 2) with λ = − 38.03°, which occurred on 22nd December 2010.Similar to the result for Period I, these earthquakes occurred far offshore, where the observation stations are sparse.
The earthquakes that occurred in the source and neighboring regions of the major earthquakes (i.e., the 2011 Tohoku-oki earthquake) particularly showed large |λ|.One representative example showing large |λ| values is an earthquake sequence in the prefectural boundary region between Ibaraki and Fukushima, where normal faulting earthquakes abruptly began to occur following the 2011 Tohoku-oki earthquake (Fig. 14; e.g., Kato et al. 2011).The largest shock in the sequence was the M w 6.6 Fukushima Hama-Dori earthquake (λ = − 116.30°,No. 71    and January 2007 in the region, and the uncertainty of inverted stress field was large. Considering that the uncertainties of focal mechanism solutions generally range from 20 to 30° (e.g., Ishibe et al. 2014), our results support the applicability of the WBH method for evaluating the expected rake angles of future large earthquakes from recent tectonic stress fields.It is worth discussing the utilization of rake angles expected from tectonic stress fields with WBH for forecasting strong ground motions and tsunamis.In addition, our study also elucidates the necessity of close attention for applying this method for areas near the occurrences of large earthquakes, seismically inactive areas where the amount of available mechanism data to invert the stress is limited, and/or far offshore areas where the observation stations are sparsely distributed.

Possible causes of large |λ|
The resulting λ between the calculated rake angles and observed rake angles mostly ranged between − 30° and 30°, whereas they were large for the source region and nearby areas of major earthquakes such as the 2011 Tohoku-oki earthquake and regions with large uncertainties where the available focal mechanism data to invert the stress field of TM2010 are limited.Here, we discuss four possible causes of generating large |λ|.

Incomplete understanding of stress field
One possible explanation for the large |λ| is an incomplete understanding of the stress field.The distribution of focal mechanism solutions from January 1997 to January 2007, which was used for the CMT data inversion, is spatially heterogeneous, and the uncertainties of the inverted stress field tend to be larger in the regions with a smaller number of available focal mechanism solutions.As mentioned, there were few available F-net data in the shallow crust (centroid depth ≤ 30 km) between January 1997 and January 2007 in the prefectural boundary region between Ibaraki and Fukushima (Figs. 14 and  15c).In the region, many earthquakes of normal faulting types with E-W tension suddenly increased immediately after the 2011 Tohoku-oki earthquake.Imanishi et al. (2012) found that microearthquakes of a normal faulting type occurred in the region before the 2011 Tohokuoki earthquake and suggested that the stress field in the shallow crust was originally and locally a normal faulting regime, unlike other regions in the Tohoku district which are in E-W compression.This also suggests that  Frohlich (1992) coseismic stress changes by the 2011 Tohoku-oki earthquake revealed the local stress heterogeneity.The large |λ| values in this region are attributed to large uncertainties of stress fields due to the lack of available focal mechanism solutions for the CMT data inversion by TM2010.The recent improvement of the CMT data inversion method to incorporate preceding inversion results into the update analysis as direct prior information about the stress field (Terakawa and Matsu'ura 2023) supports the results of Imanishi et al. (2012) and concludes that the proper stress field characterized by normal faulting appeared after the Tohoku-oki earthquake in the prefectural boundary region between Ibaraki and Fukushima.The M w 6.0 earthquake that occurred on 5th November 2018 near Kunashiri Island (No. 106 in Fig. 7c and Table 3) also occurred in a region where the available mechanism data are very limited (Fig. 15c).
A heterogeneity of the local tectonic stress regime and seismic activation in such regions also generates large |λ|.These regions have typically been highlighted by the stress increase imparted by major earthquakes (e.g., Toda et al. 2011;Terakawa et al. 2013;Terakawa and Matsu'ura 2023).In the Tohoku district, seismicity rates in most areas where the seismic activity drastically increased following the 2011 Tohoku-oki earthquake had been low for the preseismic period.The hypocenter distributions for both the periods are complementary each other (Fig. 16).

Coseismic stress change
Another possible explanation for large |λ| is coseismic stress changes, especially in the source and nearby regions of major earthquakes.Occurrences of earthquakes perturb the stress field and cause changes in seismicity (e.g., Stein et al. 1992;King et al. 1994;Toda et al. 1998;Ishibe et al. 2011b).Temporal changes in λ (Fig. 17) showed that earthquakes with large |λ| explosively increased after the occurrences of major earthquakes (e.g., the 2011 Tohoku-oki earthquake) and then gradually decayed with time.The |λ| values can be large even for regions where seismicity was originally high and the stress field was reliably inverted by abundant focal mechanism data in and near the source region of the 2011 Tohoku-oki earthquake (Figs. 12 and 13).For example, the M w 6.5 earthquake (No. 86 in Fig. 7c and Table 3) that occurred on 12th July 2014 east off Fukushima Prefecture exhibits a large |λ| (λ = 130.77°), while the number of focal mechanism data to invert the stress field is relatively large (Fig. 15c).Terakawa and Matsu'ura (2023) pointed out that the stress orientation change in a region off Fukushima and Ibaraki, the southern margin of the main rupture area of the megathrust event, seems to be real.The concentration of focal mechanism solutions with large |λ| during Period III coincides well with the above region.However, coseismic stress changes in inland regions and coastal regions of the Sea of Japan of the Tohoku district are too small to alter the stress orientation (Terakawa et al. 2013).
The |λ| values can also be larger due to transient stress changes, for example, accompanied by slow slip events.The |λ| value exceeds 30° for the M w 6.0 earthquake which occurred on 9th October 1997 (No. 3 in Fig. 7a and Table 1), approximately 6 years before the 2003 Tokachioki earthquake (No. 13 in Fig. 7a and Table 1).Based on an analysis of GPS data and seismicity, it was suggested that the deeper half of the plate interface within the rupture area of the 2003 Tokachi-oki earthquake was uncoupled and that a slow slip event occurred in several years prior to the mainshock (e.g., Baba and Hori 2006;Ogata 2005).The large |λ| value might be related to transient stress changes imparted by a slow slip event (e.g., Katsumata 2011).

Local deviatoric stress changes caused by the pore-fluid pressure increase
Another possible explanation for the large |λ| is the local deviatoric stress changes caused by the pore-fluid pressure increase, enforced by the intrusion of high-pressure fluid into a fault zone (Matsu'ura and Terakawa 2021).One representative example is the concentration of large |λ| in the source region of the 2011 Tohoku-oki earthquake (Fig. 7c).These events have focal mechanisms of normal faulting with E-W tension, which are largely deviated from the stress pattern in TM2010.Terakawa and Matsu'ura (2023) concluded that strong shaking by the 2011 Tohoku-oki earthquake resulted in a wide damage zone around the main rupture fault and that the rapid intrusion of high-pressure fluid through the suddenly developed fluid-path network triggered clustered events with abnormal focal mechanisms at shallower depths.Matsu'ura and Terakawa (2021) mathematically indicated that the fault orientations of aftershocks are not necessarily consistent with the surrounding deviatoric stress field when the enforced pore-fluid pressure changes, driven by the intrusion of high-pressure fluid into an accidentally chosen preexisting fault from deep reservoirs, are dominant.However, the activity of these events rapidly decreased within a year (Terakawa and Matsu'ura 2023), indicating that they are not the indicators of tectonic stress fields.Therefore, large |λ| values which come from these events do not undermine the applicability of the method estimating realistic rake angles for reliable prediction of strong ground motions and tsunamis caused by future large earthquakes.

Uncertainties of moment tensor solutions
Uncertainties of moment tensor solutions can also be a factor generating large |λ|, and they are basically larger for offshore regions than for inland regions due to the sparse distribution of seismograph observation stations.Variance reductions (VR), which indicate the fit between location and/or moment tensor solutions.For the M w 6.2 earthquake that occurred on 15th March 2008 near Chichi-jima Island along the Izu-Bonin Trench (Period II), which showed a large |λ| (41.58°;No. 37 in Fig. 7b and Table 2), the centroid depth of 5 km is significantly shallower than the hypocentral depth (66 km) determined by the JMA.Relatively large misfit angles are possibly due to both the large estimation errors of the stress field due to the relatively low quality of mechanism solutions and the large uncertainties of CMT solutions themselves that were used for validation.Recently, a large-scale seafloor observation network for earthquakes and tsunamis consisting of 150 observatories (S-net) was established in the Japan Trench area by NIED (Aoi et al. 2020).Furthermore, another seafloor observation network called N-net is now under construction in the Nankai Trough subduction zone by NIED.These networks will contribute to a better understanding of the stress field, reduce |λ| and uncertainties in the hypocentral location and improve early tsunami warning capabilities.

Application to fault zones in the Kyushu District
In the Kyushu district, southwest Japan, 16 FZs are evaluated as major Quaternary active FZs by the HERP.
In addition, 17 faults are evaluated as short faults in the Regional Evaluation of Active Faults in the Kyushu District (1st edition) (HERP 2013).We analyze 16 FZs and two short faults whose fault models are developed by J-SHIS.In the northern part of the Kyushu district, left-lateral faults striking roughly NW-SE directions predominate, and they are distributed at intervals of approximately 10-20 km (Fig. 19a).Conjugated rightlateral faults striking roughly NE-SW directions are also distributed in some parts of the district.In the central and southern parts, right-lateral faults and normal faults reflecting their N-S extension field are mixed.The fault types inferred from the 3D tectonic stress field by TM2010 with WBH generally reproduce the above regional characteristics well, with left-lateral fault types predominating in the northern part and right-lateral and normal fault types predominating in the central and southern parts (Fig. 19 and Table 4).Among 32 fault segments (FSs), fault types predicted by the method coincide with the HERP evaluation for 25 FSs (~ 78.1%).For example, left-lateral faults distributed in the northern part of the Kyushu district (i.e., Nishiyama FZ, Umi F, Hinata-toge-Okasagi-toge FZ and Kego FZ) are well reproduced by the method.The Izumi F is classified as a right-lateral fault by the method in the present study, while the HERP evaluated as a normal fault, including a right-lateral component.However, the calculated rake angle is − 135.65° and is near the classification boundary between the right-lateral and normal faults.The northeastern and southeastern parts of the Koshiki FZ (Koshiki segment) have been evaluated as a north-uplifting normal fault with left-lateral displacement and a northwest-uplifting normal fault, respectively (HERP 2013), while they are estimated to be mainly strike-slip faults from the method.The Ichiki segment, Koshiki-Kaikyo central segment, and Fukiagehama seiho-oki segment of the Ichiki FZ are evaluated as "normal fault with right-lateral component", "normal fault whose strike-slip component is unknown", and "normal fault whose strike-slip component is unknown", respectively (HERP 2013).The fault types inferred from the method in the present study are mainly strike-slip faults.
The 2016 Kumamoto earthquake sequence ruptured the Takano-Shirahata segment of the Hinagu FZ and the Futagawa segment of the Futagawa FZ.The largest foreshock occurred on 14th April 2016 (M w 6.1), with an F-net rake angle of − 164°, and this event was evaluated as the rupture of the northern segment (Takano-Shirahata segment) of the Hinagu FZ, a right-lateral fault striking NE-SW.The mainshock on 16th April 2016 (M w 7.1) ruptured the Futagawa FZ with an F-net rake angle of − 142°.λ is 29.65° for the largest foreshock and − 0.20° for the mainshock.The geomorphologically and geologically evaluated fault types and the observed and calculated rake angles are consistent with each other for both FZs.The fault types for the Futagawa segment of the Futagawa FZ and Takano-Shirahata segment for the Hinagu FZ are evaluated as right-lateral faults.The rake angles calculated from the 3D tectonic stress field with WBH were − 168.267° for the Futagawa segment of the Futagawa FZ and 179.84° for the Takano-Shirahata segment of the Hinagu FZ, and they were also classified as rightlateral faults.These results are consistent with a previous study (Matsumoto et al. 2018) that indicated that the prestate of stress on the fault controls the slip direction of complicated coseismic fault slip for the 2016 Kumamoto earthquake.
The calculated rake angles depend on both the regional stress field and fault geometry, and the robustness of converted fault types may differ from each fault.Furthermore, mismatch between the fault types estimated from the method and geomorphologically and geologically evaluated fault types may be due to uncertainties in the geometry of the targeted faults.In particular, the dip angles are sometimes poorly constrained and are conventionally set to some representative values for forecasting strong ground motions and tsunamis.The comparison (mismatch) of fault types obtained from geomorphology and geology with those from the method in this study would provide a valuable opportunity for the reexamination of fault types and fault geometries.

Concluding remarks and future development
The method of estimating realistic rake angles from a 3D tectonic stress field according to WBH was validated by comparison to focal mechanism data.The calculated rake angles accurately reproduced the observed rake angles; λ, the misfit angles between the observed and calculated rake angles, mostly ranged between − 30° and 30°.During Period I (January 1997 to January 2007), the absolute misfit angles, |λ|, were ≤ 30° for approximately 94.4% (for all M w ) and 86.2% (for M w ≥ 6.0) of earthquakes.|λ| values were ≤ 30° for approximately 75.1% (for all M w ) and 91.7% (for M w ≥ 6.0) of earthquakes in Period II (February 2007 to the occurrence time of the 2011 Tohoku-oki earthquake).After the occurrence of the 2011 Tohokuoki earthquake (Period III), large |λ| values were typically observed in the source and neighboring regions of the 2011 Tohoku-oki earthquake, whereas the calculated rake angles were mostly consistent with the observed rake angles as similar to periods I and II.Considering that the uncertainties of focal mechanism solutions generally range from 20 to 30° (e.g., Ishibe et al. 2014), our study supports the applicability of the WBH method for evaluating the expected rake angles of future large earthquakes from seismologically estimated tectonic stress fields.The |λ| values were large for the focal mechanism solutions of earthquakes in the source and nearby the occurrences of large earthquakes such as the 2011 Tohoku-oki earthquake (e.g., Fig. 12), the seismically inactive areas where the number of available mechanism data to invert the stress is limited (e.g., prefectural boundary region between Ibaraki and Fukushima, Fig. 14), and/or far offshore areas where the observation stations are sparsely distributed.Our study also elucidates the necessity of paying close attention to apply the method to the above areas.
We suggested four possible causes for generating a large |λ|.The coseismic (and postseismic) stress changes imparted by major earthquakes such as the 2011 Tohokuoki earthquake can change the stress pattern and generate a large |λ|.One representative example is the off Fukushima and Ibaraki region, the southern margin of the main rupture area of the 2011 Tohoku-oki earthquake (Terakawa and Matsu'ura 2023).The incomplete understanding of the tectonic stress field because of a lack of available data for short-term periods is also a possible factor for generating a large |λ|.The distribution of focal mechanism solutions of earthquakes used for the CMT data inversion is spatially heterogeneous, and available mechanism data are very limited in several regions, such as the prefectural boundary region between Ibaraki and Fukushima.The local heterogeneity of the stress field that deviates from the regional stress field and the activation of seismicity in such regions would generate a large In the present study, we also confirmed that the fault types obtained from the 3D tectonic stress field with WBH coincided well with the geomorphologically and geologically estimated fault types for FZs in Kyushu district, southwest Japan.In the Futagawa and Hinagu FZs, large earthquakes called the 2016 Kumamoto earthquake sequence occurred in April 2016.The rake angles calculated from the WBH method were consistent with both geomorphologically and geologically estimated fault types and observed rake angles.These results suggest that the rake angles for faults with known geometry can be accurately derived from the 3D tectonic stress field with WBH, whereas there are limitations for applying this method to a region where the stress has been perturbed due to recent large earthquakes (e.g., the 2011 Tohokuoki earthquake) and/or a region where the uncertainty of stress orientations is large, for example, due to a limited number of available focal mechanisms and sparse distribution of observation stations.A comprehensive study for Quaternary active faults in Japan by comparing calculated rake angles from tectonic stress fields according to WBH with geomorphologically and geologically evaluated fault types would be helpful for further understanding the temporal stability of tectonic stress fields from seismological data for ~ 10 years and the availability of the tectonic stress field with WBH to constrain the fault slip directions for future large earthquakes.
Why the tectonic stress field obtained from the CMT data inversion is consistent with the WBH?The moment tensor of a seismic event is mathematically equivalent to the volume integral of coseismic static stress changes over the whole region surrounding the source (Matsu'ura et al. 2019).Based on this relationship, in the CMT data inversion method, they represent the CMT data of a seismic event by the weighted volume integral of the stress field (Terakawa andMatsu'ura 2008, 2023).This formulation is attributed to the idea that seismic events release a part of the stress field, or that seismic events whose moment tensors are consistent with the stress tensor are the most likely.Meanwhile, based on the volume integral representation of the moment tensor, Matsu'ura et al. (2019) further elucidated that seismic slip in the direction of the resolved shear stress maximizes the efficiency of elastic strain energy release under realistic stress conditions.Therefore, WBH assumes that seismic slip occurs (in the direction of the resolved shear stress) at preexisting faults to release elastic strain energies most effectively.A seismic event can release elastic strain energies most effectively when its moment tensor is consistent with the stress tensor (Matsu'ura and Terakawa 2021).Thus, the formulation of the CMT data inversion method as above is supported by the same assumption as the WBH.In other words, the most likely seismic slip is a physical process that releases elastic strain energy most effectively.It is very interesting that the results of the present study substantiate these physical backgrounds of the CMT data inversion and WBH.
The current method for predicting strong ground motion ("Recipe") by HERP (HERP 2020) is based on the characterized source model consisting of several asperities with large slip amounts and background regions with smaller slip amounts (e.g., Irikura and Miyake 2011).The Recipe recommends setting rake angles as 90° for reverse faults, − 90° for normal faults, 0° for left-lateral faults, and 180° for right-lateral faults for active faults for which the rake angle is not specified by the subcommittee of the long-term evaluation.However, the effect of rake angle on tsunami forecasting is not negligible (Satake et al. 2022), and it is worth discussing utilizing rake angles expected from the tectonic stress field with WBH for forecasting strong ground motions and tsunamis.There are several issues in evaluating strong ground motions and/or tsunamis that contribute to disaster prevention and mitigation.Seismicity is spatially and temporally heterogeneous, and the uncertainty of the estimated 3D tectonic stress field depends on the seismicity.Evaluating the uncertainties of the rake angles arising from the uncertainties of the stress field is therefore an important issue (e.g., Terakawa 2017).Information on fault geometry (fault location, strike, and dip angle) is also essential to calculate rake angles expected from tectonic stress fields with WBH, and the improvement of fault imaging techniques and data accumulation through additional surveys are also important issues to improve the reliability of strong ground motion and/or tsunami predictions.

Fig. 1
Fig.1Schematic illustration of the method for estimating rake angles from 3D tectonic stress fields with WBH

Fig
Fig. S25) for each period are shown in Additional files.The basic characteristics of faulting [i.e., the reverse faulting dominantly distributed in northeastern Japan (Additional file 9: Fig. S9), the mixture of reverse faulting and strike-slip faulting in central Japan (Additional file 12: Fig. S12), and strike-slip faulting and normal faulting in southwestern Japan (Additional file 15: Fig. S15 and Additional file 18: Fig. S18)] are well reproduced from the 3D tectonic stress field with WBH.The |λ| values are also small for large (M w 6.0 or above) earthquakes.Among 29 earthquakes with M w 6.0 or above during Period I, |λ| is smaller than 20° for 24 (~ 82.8%) earthquakes and smaller than 30° for 25 (~ 86.2%) earthquakes (Table1).For example, the calculated rake angle is 82.58° and λ is 4.42° for a reverse-type M w 6.1 earthquake that occurred on 26th July 2003 in northern Miyagi Prefecture (No. 12 in Fig.7aand Table1).The observed rake angle of the M w 6.6 earthquake that occurred on 20th March 2005 in the southwestern off Kyushu district (− 177°; named the 2005 Fukuoka-ken Seiho-oki earthquake by JMA) can also be well reproduced by the method (169.49°) with λ = 13.51°(No. 26 in Fig.7aand Table1).Furthermore, the M w 7.9 interplate earthquake that occurred along the Kuril Trench on 26th

Fig. 4 a
Fig. 4 a Distribution of F-net focal mechanism solutions of earthquakes that occurred in Period I. b Distribution of focal mechanism solutions with calculated rake angles from the 3D tectonic stress fields with WBH for Period I.The colors indicate fault types according to the classification criteria byFrohlich (1992)

Fig. 5
Fig. 5 Distribution of λ between the observed and calculated rake angles during Period I

Fig. 7
Fig. 7 Distribution of λ for large (M w ≥ 6.0) earthquakes that occurred during a Period I, b Period II, and c Period III.The numbers in the figures are the IDs in Tables1, 2, 3

Fig. 8 a
Fig. 8 a Distribution of F-net focal mechanism solutions of earthquakes that occurred in Period II.b Distribution of focal mechanism solutions with calculated rake angles from the 3D tectonic stress fields with WBH for Period II.The colors indicate fault types according to the classification criteria byFrohlich (1992)

Fig. 11 a
Fig. 11 a Distribution of F-net focal mechanism solutions of earthquakes that occurred in Period III.b Distribution of focal mechanism solutions with calculated rake angles from the 3D tectonic stress fields with WBH for Period III.The colors indicate fault types according to the classification criteria byFrohlich (1992)

Fig. 12
Fig. 12 Distribution of λ between the observed and calculated rake angles during Period III

Fig. 14
Fig.14Distribution of focal mechanism solutions for earthquakes with a observed and b calculated rake angles in the prefectural boundary region between Fukushima and Ibaraki during Period III.c Distribution of λ between the observed and calculated rake angles during Period III.d Regional tectonic stress field byTerakawa and Matsu'ura (2010) at a depth of 10 km and distribution of F-net mechanism solutions for earthquakes with centroid depths ≤ 30 km spanning from January 1997 to January 2007 that were used for the CMT data inversion

Fig. 15
Fig. 15 Distribution of λ for large (M w ≥ 6.0) earthquakes that occurred during a Period I, b Period II, and c Period III with the number of F-net mechanism solutions used for CMT data inversion (from January 1997 to January 2007) at a spacing of 1.0° × 1.0° for ≥ 10 mechanism data points.The numbers in the figures are the IDs in Tables1, 2, 3

Fig. 16 a
Fig. 16 a Hypocenter distribution of earthquakes during a preseismic (October 1997-occurrence time of the 2011 Tohoku-oki earthquake) and b postseismic (2011 Tohoku-oki earthquake-December 2020) periods.The colors indicate the hypocentral depths.c Distribution of F-net mechanisms and their fault types according to the criteria by Frohlich (1992) for c preseismic and d postseismic periods

Fig. 19 a
Fig.19a Types of faulting for the Quaternary major active faults that were evaluated by the Headquarters for Earthquake Research Promotion (HERP).b Calculated rake angles for the Quaternary major active faults that were evaluated by HERP.A type of faulting for each fault is converted from the calculated rake angle (see the text for the classification criteria).c F-net mechanism solutions for earthquakes that occurred from January 1997 to December 2021.The colors indicate fault types according to the classification criteria byFrohlich (1992)

-
The parameters of the fault models are according to the J-SHIS by NIED a R reverse, N normal, RL right-lateral, LL left-lateral |λ|.Local deviatoric stress changes caused by the porefluid pressure increase enforced by the intrusion of highpressure fluid into a fault zone (Matsu'ura and Terakawa 2021) and the presence of large uncertainties in focal mechanism solutions for earthquakes (e.g., comparatively low VR along the Izu-Bonin Trench and Ryukyu Trench) are other possible factors generating a large |λ|.

Table 1
List of focal mechanism solutions for large (M w ≥ 6.0) earthquakes with calculated rake angles and λ during Period I No.

Year Month Day Hour Minute Second Longitude (deg.) Latitude (deg.) Depth (km) M
No.

Table 3
), occurred 1 month after the 2011 Tohoku-oki earthquake.The |λ| values for this sequence are enormously large because the inverted tectonic stress field is reverse-type with roughly E-W compression.However, there were few available F-net data in the shallow crust (centroid depth ≤ 30 km) between January 1997

Table 2
List of focal mechanism solutions for large (M w ≥ 6.0) earthquakes with calculated rake angles and λ during Period II No.

Year Month Day Hour Minute Second Longitude (deg.) Latitude (deg.) Depth (km) M
No.

Table 3
List of focal mechanism solutions for large (M w ≥ 6.0) earthquakes with calculated rake angles and λ during Period III No.

Year Month Day Hour Minute Second Longitude (deg.) Latitude (deg.) Depth (km) M
No.

Year Month Day Hour Minute Second Longitude (deg.) Latitude (deg.) Depth (km) M
No.

Year Month Day Hour Minute Second Longitude (deg.) Latitude (deg.) Depth (km) M
No.

Year Month Day Hour Minute Second Longitude (deg.) Latitude (deg.) Depth (km) M
No.

Table 4
Quaternary active faults in Kyushu district, southwest Japan, that were evaluated by HERP and the fault models were developed by J-SHIS from NIED, and the geomorphologically and geologically evaluated fault types were compared with those inferred from the regional stress field with WBH