Introduction

Earthquakes contribute to dissipating the energy accumulated in the brittle lithosphere due to the tectonic stress arising from the motion of contiguous crustal volumes with respect to each other. Thrust faulting, usually featured by angles of dip ranging in between 5°−30°, mostly occurs along the margins of plates, where their motions induce elastic strain accumulation, which is released by multifaceted fault slip dynamics ranging from almost periodic silent events to megathrust earthquakes1. Strike-slip-faulting earthquakes are localized along steeply dipping faults (70°−90°) or transcurrent plate boundaries and transfer zones, while normal faults develop along rift zones in extensional regimes having intermediate dip (45°−65°). Structural, morphological, and geophysical differences have been highlighted among the three main tectonic settings2,3,4. Normal faults cause fracturing mainly concentrated in the hanging walls and spaced clusters of parallel faults across rifting areas5,6,7. Intricated geometries are also typical of transcurrent regions, often accompanied by releasing and restraining bends or step-overs and other geological structures shedding light on complex spatial stress patterns8,9; experiments in the lab support geological observations10. Conversely, thrust-faulting earthquakes usually occur along gently dipping subductions characterized by a unique, longitudinally extended front which slip is localized at11,12.

The dynamics of fluids has been also noticed to vary with the faulting type13. In extensional tectonics, fluids percolate into fractures during the interseismic period, whereas they are expelled during the coseismic phase, while the opposite is observed in compressive geodynamic settings14,15. Such a fan of geological manifestations also mirrors a seismological counterpart16: normal fault earthquakes are featured by a steeper power-law frequency-size distribution than transcurrent and reverse seismicity respectively17,18,19, with lower cut-off magnitudes20 and higher fractal dimension of the hypocenters’ time series of aftershocks21,22. The latter property clearly states that the post-seismic relaxation is spread over several fault patches.

Recently, sound evidence has been provided that fault structure plays a key role in driving almost all the crucial large-scale processes characterizing seismic dynamics23 such as the localization of earthquake nucleation, propagation and arrest of rupture and aftershock occurrence. The rheological properties of fault rocks are related to the dynamics of coseismic slip, which, in turn, is connected to the topological features of faulting. The roughness of the dislocation surface is proven to affect stress accumulation, strain accommodation and release24. In this regard, geometrical complexities, and heterogeneity, i.e., state of fracturing and fault topology, represent the hidden factor shaping the seismic dynamics from coseismic to spatial and temporal scale of tectonics25,26,27. Fault topology is in fact suggested to control the statistical properties of seismicity, producing a wide range of behaviors with different varying recurrence times, periodicity, and consistency over time in seismic time series, its magnitude and radiation pattern28.

The moment tensor is the most complete quantitative information that can be extracted from seismic recordings29. The moment tensor is symmetric with six independent components. It can be decomposed into an isotropic part (ISO), a double-couple (DC) and a compensated linear vector dipole (CLVD). The isotropic part provides a measure of the volume change, while the deviatoric contribution, i.e., the second and the third terms, have null trace.

In the simplest case, an earthquake can be thought as a unidirectional slip on a single fault plane, so that its moment tensor can be represented by a double-couple of forces acting with null net torque30. However, real faults are not planar surfaces, nor have they well-defined simple geometries at all. They are fractals31,32,33 which stem from spatial self-organization. Hence, large-scale heterogeneities34, not only friction, are likely to control the development of fault systems producing complex ruptures which can lead to low DC components in the moment tensors35. Non-DC events can have different physical origins: high ISO terms are usually caused by varying fluid flow, landslides, and volcanic eruptions; explosions are classical examples of artificial ISO-earthquakes. Non-planar or multi-patch ruptures are instead characterized by significant CLVD components, e.g.36,37. They arise from the sum of DC-like moment tensors produced by shear faulting on locally planar fractures with different spatial orientations38. Other possible origins are shear fractures in heterogeneous and anisotropic media (e.g., facies transitions). Therefore, the CLVD contribution can be considered a suggestive marker for topological complexity of the seismic source and thus, moment tensors may be suitable to delve into the dynamics of the coseismic fracture. Nevertheless, above all for shallow seismic events, some moment tensor components may not be accurately determined, producing spurious non-DC39. Wrong hypocentre localizations, centroid mis-location, and inaccurate velocity models of seismic waves in the crust and in the mantle can also be responsible for spurious non-DC moment tensors40. For these reasons, non-DC components have been considered to be artifacts in most of the cases41,42.

The previous conclusion was also supported by the large scatter of data and errors in different catalogues43 suggesting that better detection procedures should be required before a minimal accuracy could be get for a reliable analysis of CLVDs.

Nowadays, much larger and detailed moment tensor catalogues are available than in the past, which allow us to perform a new analysis.

Results

The statistical analysis of two global catalogues44,45,46,47 of moment tensors clearly shows that the DC contributions are not uniformly distributed as a function of the angle of rake. The type of event is classified according to the Aki-Richards convention29: an earthquake is a “thrust-faulting event” if its rake is in the range 90° ± 30°, while normal-faulting earthquakes fall in −90° ± 30° and strike-slip events are in 0° ± 30°, ±180°  30°. Thrust-faulting events are characterized by higher DC values, as already reported in42, followed by normal, transpressive and strike-slip ones respectively. Transtensional tectonic settings host seismicity featured by the largest CLVD contribution (Fig. 1a, c). Data confirms the statistical significance of the variability in the composition of moment tensors in different faulting types (Table 1) also in regional moment tensor catalogues48,49, even though with less evident outputs, likely because of asymmetric statistics of focal mechanisms (e.g., the NIED database mostly lists reverse fault earthquakes and only few normal-faulting ones). Different behaviours are also observed regarding the correlation between CLVD values and the size of earthquakes. A clear increase of the DC with the magnitude is measured for compressive and transpressive data sets, while the size of normal and transtensional events do not seem to affect their non-DC components. Transcurrent seismicity displays a positive correlation at moderate magnitudes while it turns negative at larger (> 6) ones (Fig. 1b, d). The average DC contributions are found to be positively correlated with the magnitude of the largest event in the catalogue (Fig. 2a, c), i.e., the corner magnitude of the frequency size-distribution, and negatively correlated with the b-value of the Gutenberg-Richter law (Fig. 2b, d). This output is coherent with the DC distribution as a function of the rake angle and with its correlation with the size of earthquakes. Our results for the Global CMT and the ISC reviewed catalogues are compatible with each other. A plot of the normalized frequency of the rake angle as a function of the DC clarifies that the differences in the mean DC values in Table 1 are only partially produced by displaced statistical modes, while the misfit is mostly due to larger high-order moments of the distributions. Thrust earthquakes are usually featured by elevated plane shear components of the moment tensors with limited dispersion. The variance of the distribution increases for normal and strike-slip earthquakes. It reaches its maximum value in the case of transtensional seismicity, whose distribution is markedly spread over a large interval (Fig. 3a). At last, the geographical localization of seismic events as a function of their moment tensor contents seems to be not spatially homogeneous (Fig. 3b): non-DC earthquakes more likely occur along transfer zones, intra-plate settings and slow deforming continental regions. Differently, almost pure DC events are clustered along the subduction zones. This qualitative observation is also coherent with previous results. While the first part of our analysis is focused on shallow seismicity (hypocentral depth lower than 50 km), we also investigate how the composition of moment tensors is affected by depth (Fig. 4). Our results show that the DC components tend to be more and more uniform as a function of the angle of rake increasing depth.

Fig. 1: Trends of double-couple component of the seismic moment tensors.
figure 1

Results for the GCMT and ISC catalogues are shown in panels (a, b and c, d) respectively, also with respect to different magnitude bins. Seismicity is classified into normal (green), strike-slip (red), thrust (blue) according to the Aki-Richards definition29, while transpressive and transtensional events are respectively in yellow and in cyan, as represented in the legends on the right panels. Error bars in (a, c) stand for 2σ.

Table 1 Double-couple components of moment tensors in different seismic catalogues.
Fig. 2: Correlations between double-couple component and maximum magnitude and b-value of the Gutenberg-Richter distribution.
figure 2

Double-couple trends as a function of the magnitude of the largest event in the sets for the GCMT and ISC catalogues (a, c). b-value as a function of the average double-couple of the earthquakes of each set (b, d). (b, d) suggest that normal-faulting earthquakes have lower double-couple components than strike-slip and thrust-faulting events at a global scale. Error bars stand for σ.

Fig. 3: Normalized DC frequency distributions for different classes of earthquakes.
figure 3

Thrust-faulting events showcase a distribution peaked at elevated DC values, while normal, transpressive, strike-slip and transtensional seismic events tend to have increasingly left heavy-tailed distributions with non-negligible percentages of earthquakes featured by low DC (a). A geographic map with points representing earthquakes coloured according to their DC (b, data from Global CMT, 2000–2021, Mw ≥ 5.3).

Fig. 4: Trends of double-couple component of the seismic moment tensors in the GCMT catalogue (Mw ≥ 5.3, 1990–2021) for different ranges of hypocentral depth.
figure 4

The red line represents DC values for seismicity shallower than 30 km, while the orange one stands for events in between 30 and 60 km, the salmon one is the trend of 60–100 km and the yellow line for deep earthquakes (>100 km). Error bars are not plotted for the sake of readability; however, their size has the same order of magnitude of those in Fig. 1, usually larger for deeper events.

Discussion

Moment tensors are outputs of an inverse problem which can be calculated using different techniques with different strengths and weaknesses. In our analysis we considered catalogues obtained using different approaches: for instance, the inversion procedure of the GCMT is mainly based on surface waveforms, while the ISC-GEM also derives focal mechanisms via body waves polarities and amplitudes. The first method usually guarantees higher quality datasets since it takes advantage of the low-frequency information of the entire seismogram; on the other hand, it is rather sensitive to the adopted seismic velocity model, which can introduce severe errors in the estimation of Mw < 5 earthquakes moment tensors50. In our analysis, we only considered reviewed seismic events occurred in the last two decades with magnitudes Mw > 5 (compare with Table 1 for details) to guarantee a satisfactory data quality. This choice allows us to reject the hypothesis of a uniform distribution of the DC component as a function of the rake angle, as proven by the significant p-values of the χ2 test for both the GCMT and the ISC-GEM catalogues (p 0.001, Table 1).

Prudence is required because even though heterogeneities in source geometry become smaller and smaller with magnitude so that moment tensor inversions are more consistent for large earthquakes, Kagan’s angles51 are still large for Mw 6 (up to 30°) in global catalogues and significant values (5°–15°) also affect the largest events. For this reason, large numbers of events with significant magnitude are needed for keeping uncertainties as small as possible. Despite efforts, the issue of possible systematic errors which may affect both the DC percentage and its uncertainty remains. In this regard, since a direct inspection is extremely complicated for global and large-scale regional catalogues in which uniform criteria are adopted for the inversion algorithms, we focused on the strongest arguments against a physical origin of large percentages (>30%) of CLVD component in moment tensor solutions of moderate and large earthquakes showing that, even though present, they are not sufficient to explain the variability of the non-DC components in moment tensors.

The first argument is the decrease of the non-DC component for larger earthquakes51. The plots in Fig. 1b, d prove something different; the decrease is limited to thrust and transpressive earthquakes, while it is not observed in other faulting types, at least above magnitude 5.5. This observation suggests a physical reason for this phenomenon. For instance, the finite-source ruptures with variations in rake and slip amplitude over the fault plane, or multiple rupture planes are known to lead to non-DC moment tensors38.

A second argument regards the moment tensor inversion procedure, which tend to use longer period seismic signals with increasing magnitudes46 reducing the resolution of the inversion producing an apparent smoothing of the earthquake fault rupture. This effect may contribute to the progressive decrease of the non-DC components in thrust-faulting earthquakes; however, it does not explain why this phenomenon is not observed in other fault types in the same range of magnitudes. Once again, it is possible that large seismic events, along subduction zones, may involve simpler seismogenic sources than those occurring in other tectonic settings.

Another observation concerns the lack of correlation between non-DC terms between moment tensor catalogues41. We calculated the correlation between the DC values of Mw > 5.5 earthquakes of the GCMT and ISC catalogues grouped as a function of the rake angle (intervals of 3°). The Pearson correlation coefficient ρ 0.81 (compare with Fig. S1 in the Supplementary Material) clearly states that the accuracy of the measurements of the components of the moment tensor is sufficient to achieve coherent results between different catalogues.

An additional non-physical effect on the inversion should be taken into account before exploring possible mechanisms affecting the composition of the moment tensors: the dip of faults produces a bias in the estimation of the DC components due to the different radiated wave patterns.

In fact, strike-slip events tend to very strongly radiate SH (Love) waves, but weak body waves, while Rayleigh waves are more dominant in dip-slip mechanisms.

For this reason, we also perform an investigation focusing on large (Mw ≥ 7.0) shallow (depth lower than 50 km) events listed in the GCMT catalog. It shows a weak negative correlation between the values of the DC component of earthquakes and angle of dip of faults (Fig. 5).

Fig. 5: Distribution of DC components as a function of the angle of dip, δ, for large (Mw ≥ 7.0) shallow (hypocentral depth within 50 km) earthquakes in the GCMT catalogue (1990–2021).
figure 5

More than three hundred events are considered and a clear trend, although weak and with large scattering of data, indicates the DC values are negatively correlated to the angle of dip. The best linear fit is given by DC = −0.19δ  + 96.8. This effect has origin in the different composition of the seismic pattern used for moment tensor inversion and must be considered for a correct interpretation of results shown in the previous figures.

This effect introduces an underestimation of DC values in strike-slip events; however, the residual distribution of double-couple percentages is not uniform as a function of the rake angle. Considering together the results shown in Figs. 13, 5, normal-faulting and transtensional events have lower DC components, while thrust-faulting events are those with higher values in agreement with the outcomes reported in42. Therefore, our results suggest that the variations in the composition of moment tensors may also have physical origins. Their interpretation is rather simple: the complexity of earthquake ruptures statistically decreases from normal faulting, strike-slip to thrust-faulting earthquakes, being the last ones usually associated to sharper, more continuous faults with respect to the previous tectonic settings characterized by multiple anastomosed faults and larger damage zones. This structural variation could explain why the coseismic slip along thrust faults shows higher double couple.

However, even inside the same faulting type, a certain degree of variability is observed. Figure 6 shows that thrust-faulting intraplate earthquakes have lower DC components than subduction events, while there is not a statistically significant difference between mid-ocean and continental rifting events (Fig. 7).

Fig. 6: Double-couple component of the seismic moment tensors in the GCMT catalogue for 5.3 ≤ Mw ≤ 7.0 shallow (depth lower than 50 km) thrust-faulting (30° ≤ rake angle ≤ 150°) events (1990–2021).
figure 6

The maximum magnitude is limited such that the mean values of magnitudes of the sets of intraplate and interplate events are not statistically different. This choice is needed to avoid a bias in the estimation of average DC values since subduction events reach larger magnitudes than intraplate thrust-faulting events. Intraplate thrusting events (green points, a) are analysed separately from subduction earthquakes (red points, a). Error bars stand for 2σ. Seismic events are selected within the regions highlighted in the map, intraplate earthquakes are those inside the green contours, while interplate are within the red ones (b). Thrust-faulting intraplate earthquakes show lower DC components than the events occurring along the subduction zones.

Fig. 7: Double-couple component of the seismic moment tensors in the GCMT catalogue for 5.3 ≤ Mw ≤ 7.5 shallow (depth lower than 50 km) normal-faulting (−150° ≤ rake angle ≤ −30°) events (1990–2021).
figure 7

The maximum magnitude is limited such that the mean values of magnitudes of the sets of intraplate and interplate events are not statistically different. This choice is needed to avoid a bias in the estimation of average DC values since mid-oceanic events reach larger magnitudes than those at continental rifting zones (usually showing an upper bound at about Mw 7.5. Intraplate thrusting events (green points, a) are analysed separately from subduction earthquakes (red points, a). Error bars stand for 2σ. Seismic events are selected within the regions highlighted in the map, intraplate earthquakes are those inside the green contours, while interplate events are within the red ones (b). No statistical difference is observed in the DC components of normal-faulting intraplate and interplate earthquakes.

This explanation is coherent with structural and physical observations of a broad fan of dynamic behaviors of seismicity and faulting in different tectonic regimes52 (Fig. 8). Continental normal faults can generate earthquakes with maximum Mw around 7.5, lower than mid-ocean normal-faulting events (~ Mw 8.0); analogously, intraplate thrusts can nucleate events up to Mw 8.0–8.5, whereas subductions may reach at least Mw 9.5 because of larger and more uniform active tectonic settings. Therefore, it is counterintuitive that thrusts can produce thinner damage zones and dislocation widths (Fig. 8).

Fig. 8: The width of the cataclastic band and the sharpness of faults tend to decrease moving from normal faults (NF), strike-slip (SS), to thrust faults (TF).
figure 8

a, b, c: three examples from the Mesozoic carbonate rocks of the central Apennines (a) and the Dolomites in the eastern Alps (b, c) in central and northern Italy. d: Thrusts are generally associated with sharper, continuous faults with respect to strike-slip and normal faults, rather characterized by multiple and anastomosed faults and generating wider cataclastic damage zones. This structural variation could explain why the coseismic slip along thrust faults shows in average higher double couple.

This paradox can be explained by the different degree of complexity of the geological structures and seismogenic sources statistically affected by the local tectonic setting because of various stress patterns53, different balances of forces generating long-term crustal motion along normal faults versus thrusts54 (i.e., mainly gravitational energy and elastic energy respectively) and physical properties of rocks55.

The wide range of possible rupture dynamics also raises concerns about the appropriateness of using a unique, average moment tensor for representing earthquakes involving complex ruptures (e.g., simultaneous or cascade fault activations, with large, >30%, CLVD percentages). The most important reason is related to the estimation of the size of earthquakes56. Because of the triangle inequality,

$${M}^{{TH}} = {\int _{{t}_{0}}^{{t}_{f}}}\left[\begin{array}{ccc}{M}_{{rr}}\left(t\right) & {M}_{r\theta }\left(t\right) & {M}_{r\varphi }\left(t\right)\\ {M}_{\theta r}\left(t\right) & {M}_{\theta \theta }\left(t\right) & {M}_{\theta \varphi }\left(t\right)\\ {M}_{\varphi r}\left(t\right) & {M}_{\varphi \theta }\left(t\right) & {M}_{\varphi \varphi }\left(t\right)\end{array}\right]{dt}\,\leftrightarrow \,{M}_{0}^{{TH}}\\ = {\int _{{t}_{0}}^{{t}_{f}}}{\int _{\varSigma }^{}}\,\mu \left(x,t\right)u\left(x,t\right)\,{d\varSigma}\,{dt}\,{\ge \,M}_{0}^{{OB}}\\ = {\int _{\varSigma }^{}}\,\mu \left(x\right)u\left(x\right)\,{d\varSigma}\leftrightarrow {M}^{{OB}}=\left[\begin{array}{ccc}{M}_{{rr}} & {M}_{r\theta } & {M}_{r\varphi }\\ {M}_{\theta r} & {M}_{\theta \theta } & {M}_{\theta \varphi }\\ {M}_{\varphi r} & {M}_{\varphi \theta } & {M}_{\varphi \varphi }\end{array}\right]$$
(1)

the operational value of the seismic moment is an appropriate estimation of the energy nucleated during coseismic slip only if the rupture occurs along a perfect fault plane, otherwise it is underestimated.

MTH is the theoretical moment tensor associated with the theoretical seismic moment \({M}_{0}^{{TH}}\), considering the activation of complex seismogenic structures, MOB is the output of moment tensor inversion (average), \({M}_{0}^{{OB}}\) is the operational seismic moment, μ(x) is the local shear modulus, u(x) is the local slip in the position x, Σ is the fault surface obtained via moment tensor inversion and interpretation of geophysical data, while t0 and tf represent the times at which nucleation begins and arrest occurs respectively.

Therefore, average radiation coefficients may differ significantly from the pattern produced by the real high-resolution shear slip in the case of roughly non-planar faulting, i.e., non-DC earthquakes (Fig. 9), underestimating magnitudes. Our conclusions agree with the analysis reported in57. It shows that about ten thousand global event solutions updated from the GCMT catalogue to account for the effects of Earth’s heterogeneity are featured by larger scalar moments and double-couple components than previously thought. Since a non-uniform distribution of the DC as a function of magnitude has been highlighted in compressive and transcurrent tectonic settings, further analysis should be done in order to evaluate its impact on the scaling exponent of the Gutenberg-Richter distribution, which we also proved to be negatively correlated to the percentage of DC. This effect might be an additional source of bias for the estimation of the b-value, which is already subject to several other pitfalls58,59,60, reducing its potential reliability for seismic hazard assessment.

Fig. 9: A moment tensor and a focal mechanism are usually linked to each seismic event of large and moderate magnitude assuming the slip occurs on a planar fault.
figure 9

If the fracture involves several faults, the DC contribution decreases producing a large CLVD term, which is then a good proxy for the complexity of faulting. However, the attribution of a single moment tensor may imply an underestimation of the seismic moment of the earthquake if the fault slip involves several fault patches with different orientations (a). The degree of complexity of the fault rupture can be simply characterized looking at the DC contribution of the moment tensor, which seems to differ as a function of the tectonic setting (b, data from ISC Catalogue) accordingly with the well-known positive correlation between the fractal dimension of the hypocentres time series and b-value (c, data from32, 65,66,67,68, error bars from32, 65, 66, estimated using the ranges of values in67 assuming uniform distributions and computed from68 following methods described in Section 4 and references therein).

Conclusions

The increasing quality and completeness of global moment tensor catalogues allow to enhance our knowledge of seismic processes delving into the connection between seismicity, tectonics, and faulting. In our analysis, we draw attention to different compositions of the moment tensors of moderate and large seismic events as a function of the tectonic setting: thrusts host earthquakes with more elevated DC percentages with respect to strike-slip and normal faults. The CLVD component decreases as the size of earthquakes increases in reverse faulting, while this trend is weaker or absent in other classes of seismicity, with also an upstream behavior noticed in transcurrent earthquakes, likely due to noise in the inversion procedure. An apparent significant departure from planar shear, even though of debated origin, is found to be positively correlated to the b-value and negatively related to the corner magnitude of the frequency-size distribution which is compatible with a systematic magnitude underestimation in low DC earthquakes. Our results suggest that, at least for large seismic events featured by suspiciously high non-DC components (e.g., 30/10/2016 Mw 6.5 Norcia61 and 13/11/2016 Mw 7.8 Kaikoura62 earthquakes) should be considered to better assess their size accurately also because of possible impact on seismic forecasting.

Methods

The DC contributions are calculated starting from the moment tensor components63, so that, given the eigenvalues M1 > M2 > M3, the fractional components are given by

$$\left\{\begin{array}{c}{C}_{{ISO}}=\,\frac{{M}_{1}+\,{M}_{2}+{M}_{3}}{3\,M} \hfill \\ {C}_{{CLVD}}=\,\frac{2({M}_{1}+\,{M}_{3}-{2M}_{2})}{3\,M}\hfill\\ {C}_{{DC}}=\frac{{M}_{1}-{M}_{3}-\left|{M}_{1}+{M}_{3}-{2M}_{2}\right|}{2\,M}\end{array}\right.$$
(2)

where \(M=|{C}_{{CLVD}}|+|{C}_{{ISO}}|+{C}_{{DC}}\).

Only earthquakes with size above the completeness magnitude are included in our analysis. Volcanic, anthropic, and collapse-related events are removed (ISO ≈ 0 in shear faulting). The error bars in Figs. 1, 2, 7, 8 represent the standard deviations for the mean of the DC percentages in each rake interval. The compositions of the moment tensors are assumed to be independent and identically distributed within each rake interval, so that the uncertainties of the DC components are dominated by the statistical fluctuations at least for large catalogues. The classification of earthquakes according to their focal mechanisms follows the classical definition by Aki and Richards29. An earthquake is a thrust-faulting event if its rake is in the range 90° ± 30°, while normal-faulting earthquakes fall in −90° ± 30° and strike-slip events are in 0° ± 30°, ±180°  30°. The b-value is calculated according to the Tinti-Mulargia method64.