Information Geometry of Earthquake Depth Distributions: Kullback-Leibler and Jensen-Shannon Divergence Across Tectonic Settings

stepstep_labs

Abstract

Earthquake depth distributions encode fundamental information about the mechanical and thermal structure of tectonic plate boundaries, yet quantitative comparison of these distributions across settings has relied predominantly on summary statistics and parametric models. This study introduces an information-theoretic framework for measuring the divergence between earthquake depth distributions in distinct tectonic environments. Using 340,476 earthquakes of magnitude 4.0 and above recorded in the USGS FDSN catalog from 2000 to 2024, depth distributions are estimated for five major tectonic settings — subduction zones, spreading ridges, transform boundaries, continental collision zones, and continental rifts — and compared via Kullback-Leibler (KL) divergence, Jensen-Shannon divergence (JSD), and Shannon entropy. The mean between-setting JSD of 0.321 bits exceeds the mean within-setting JSD of 0.065 bits by a factor of 5.0, confirming that tectonic settings occupy well-separated positions in information-geometric depth space. KL divergence is strongly asymmetric: the information cost of mistaking a subduction zone for a transform boundary (3.583 bits) far exceeds the reverse (3.300 bits), reflecting the richer depth structure of subduction environments. Shannon entropy maps directly to geological complexity, ranging from 4.02 bits (16.2 effective depth modes) for subduction zones to 0.67 bits (1.6 modes) for continental rifts. These divergence measures remain largely invariant across magnitude ranges, indicating that depth distribution structure is a stable property of tectonic settings rather than an artifact of event selection.

1. Introduction

The depth at which earthquakes nucleate is one of the most geophysically informative parameters in seismology. Depth reflects the thermal, rheological, and compositional structure of the lithosphere and upper mantle, and it varies systematically across tectonic environments. Subduction zones generate seismicity to depths exceeding 680 km along the downgoing slab (Frohlich, 2006), whereas transform faults and continental rifts confine nearly all ruptures to the shallow crust above 30 km. These well-known depth contrasts form a cornerstone of plate tectonic theory and are critical for seismic hazard assessment, tomographic imaging, and geodynamic modeling.

Despite the centrality of depth distributions, quantitative methods for comparing them across tectonic settings have remained surprisingly limited. The standard approach relies on summary statistics — mean depth, median depth, or depth percentiles — which compress the full distributional shape into a single number and discard information about multimodality, skewness, and tail behavior. Histograms provide a visual comparison but lack a scalar measure of distributional difference. Parametric models, such as mixtures of Gaussians or exponential distributions, impose assumptions about the functional form that may not generalize across settings. Gutenberg and Richter (1954) cataloged global depth variations in their foundational work, and subsequent studies refined depth characterization for specific environments (Syracuse and Abers, 2006; Frohlich, 2006), but a unified, non-parametric framework for quantifying how much information is lost when one depth distribution is substituted for another has not been developed.

Information theory provides a natural and principled framework for such comparisons. The Kullback-Leibler (KL) divergence (Kullback and Leibler, 1951) quantifies the expected information loss when a probability distribution ( Q ) is used to approximate a true distribution ( P ), measured in bits. The Jensen-Shannon divergence (JSD), introduced by Lin (1991), symmetrizes KL divergence and yields a bounded measure whose square root satisfies the triangle inequality, thereby defining a proper metric on the space of probability distributions. Shannon entropy (Shannon, 1948) quantifies the intrinsic uncertainty or complexity of a single distribution. Together, these measures offer a complete toolkit for characterizing both the internal complexity of individual depth distributions and the pairwise distances between them.

This study presents the first information-theoretic quantification of earthquake depth distribution divergence across tectonic settings. The analysis addresses three questions. First, how large is the information-theoretic distance between depth distributions of different tectonic settings, and does it exceed within-setting variation? Second, does the asymmetry of KL divergence carry physical meaning about the relative complexity of tectonic environments? Third, can Shannon entropy serve as a scalar measure of geological complexity that distinguishes tectonic settings in a physically interpretable manner?

The results demonstrate that tectonic settings are well separated in distribution space, with a between-to-within divergence ratio of 5.0. KL asymmetry reveals a directional information cost that reflects the richer depth structure of subduction and collision zones. Shannon entropy provides an elegant mapping from depth distributions to an effective number of depth modes, ranging from 16.2 for subduction zones to 1.6 for continental rifts. These findings establish information geometry as a new lens for tectonic classification and seismological comparison.

2. Methods

2.1 Data

Earthquake data were obtained from the United States Geological Survey (USGS) Federated Digital Seismograph Network (FDSN) earthquake catalog via the web service API. All events of magnitude ( M \geq 4.0 ) recorded between January 1, 2000 and December 31, 2024 were retrieved, yielding a total of 340,476 earthquakes. The magnitude threshold of 4.0 was chosen to ensure global catalog completeness while retaining sufficient sample sizes for robust distributional estimation across all tectonic settings.

2.2 Tectonic Classification

Events were classified into five major tectonic settings based on geographic bounding boxes corresponding to well-established tectonic domains: subduction zones, spreading ridges, transform boundaries, continental collision zones, and continental rifts. Eight individual regions were defined within these settings: Japan Subduction, South America Subduction, and Indonesia Subduction (subduction); the Mid-Atlantic Ridge and East Pacific Rise (spreading ridges); the San Andreas fault system (transform); the Himalayan collision zone (continental collision); and the East African Rift (continental rift). This geographic classification, while necessarily simplified, captures the dominant tectonic process governing seismicity in each region and is consistent with standard practice in global seismicity studies.

The resulting sample sizes by setting are: subduction (( n = 119{,}183 )), spreading ridge (( n = 21{,}836 )), transform (( n = 1{,}223 )), continental collision (( n = 12{,}643 )), and continental rift (( n = 1{,}714 )).

2.3 Depth Distribution Estimation

Depth distributions were estimated using histograms with 70 bins of 10 km width, spanning 0–700 km. To avoid zero-probability bins, which would cause KL divergence to diverge to infinity, a Jeffreys prior smoothing was applied with parameter ( \alpha = 0.5 ). Specifically, the smoothed probability for bin ( i ) in a given setting was computed as:

[ \hat{p}i = \frac{n_i + \alpha}{\sum{j=1}^{70}(n_j + \alpha)} ]

where ( n_i ) is the count in bin ( i ). This approach is a standard regularization in information-theoretic applications and has minimal impact on distributions with large sample sizes while ensuring finite divergence values for all pairs.

2.4 Kullback-Leibler Divergence

The Kullback-Leibler divergence from distribution ( P ) to distribution ( Q ) is defined as:

[ D_{\mathrm{KL}}(P | Q) = \sum_{i=1}^{70} p_i \log_2 \frac{p_i}{q_i} ]

measured in bits. ( D_{\mathrm{KL}}(P | Q) ) quantifies the expected number of additional bits required to encode samples from ( P ) using a code optimized for ( Q ). Crucially, KL divergence is asymmetric: ( D_{\mathrm{KL}}(P | Q) \neq D_{\mathrm{KL}}(Q | P) ) in general. This asymmetry is not a deficiency but an informative property, as it reflects the directional information cost of confusing one distribution for another. When ( P ) has support over depth bins where ( Q ) assigns low probability, ( D_{\mathrm{KL}}(P | Q) ) is large; the reverse substitution may incur a different cost.

2.5 Jensen-Shannon Divergence

The Jensen-Shannon divergence symmetrizes KL divergence and is defined as:

[ \mathrm{JSD}(P, Q) = \frac{1}{2} D_{\mathrm{KL}}(P | M) + \frac{1}{2} D_{\mathrm{KL}}(Q | M) ]

where ( M = \frac{1}{2}(P + Q) ) is the mixture distribution. JSD is symmetric, non-negative, and bounded above by 1 bit (for base-2 logarithms). The quantity ( \sqrt{\mathrm{JSD}(P, Q)} ) satisfies the triangle inequality and therefore constitutes a proper metric on the space of probability distributions (Lin, 1991). This property enables a rigorous information geometry of tectonic settings: each setting can be treated as a point in a metric space, and distances between points quantify how distinguishable the corresponding depth distributions are.

2.6 Shannon Entropy

The Shannon entropy of a depth distribution ( P ) is:

[ H(P) = -\sum_{i=1}^{70} p_i \log_2 p_i ]

Entropy quantifies the uncertainty or complexity of a distribution. A distribution concentrated in a single bin has ( H = 0 ) bits, while a uniform distribution over all 70 bins achieves the maximum entropy of ( \log_2 70 \approx 6.13 ) bits. The quantity ( 2^{H(P)} ) gives the effective number of equiprobable bins, or effective depth modes, providing an intuitive measure of how many distinct depth regimes contribute substantially to the seismicity of a given tectonic setting.

2.7 Bootstrap Confidence Intervals

Uncertainty in all divergence and entropy estimates was quantified via nonparametric bootstrap resampling. For each tectonic setting, 2,000 bootstrap resamples were drawn (subsampled to 5,000 events per group for computational efficiency), depth distributions were re-estimated, and divergence measures were recomputed. The 2.5th and 97.5th percentiles of the bootstrap distribution define the 95% confidence interval.

3. Results

3.1 Dataset Overview and Depth Percentiles

The 340,476 earthquakes span the full range of seismogenic depths observed globally. Table 1 summarizes the depth percentiles for each tectonic setting.

Table 1. Depth percentiles (km) by tectonic setting.

Subduction zones exhibit the broadest depth range, with a median of 36 km but a 95th percentile extending to 248 km, reflecting the presence of intermediate-depth and deep seismicity along subducting slabs. Continental collision zones show a similar breadth (p95 = 212 km), consistent with the deep seismicity beneath the Himalayan orogen. Transform boundaries and continental rifts are confined to the shallow crust, with median depths of 7 km and 10 km, respectively, and maximum recorded depths below 40 km.

3.2 Shannon Entropy by Setting

Shannon entropy reveals a clear ordering of geological complexity across tectonic settings (Table 2).

Table 2. Shannon entropy and effective depth modes by tectonic setting.

Subduction zones exhibit the highest entropy (4.017 bits), corresponding to 16.2 effective depth modes. This reflects the coexistence of shallow interplate seismicity, intermediate-depth intraslab earthquakes, and deep-focus events along the full extent of the subducting slab. Continental collision zones rank second (3.623 bits, 12.3 modes), consistent with their comparable depth range and multi-modal depth structure. Spreading ridges occupy an intermediate position (2.167 bits, 4.5 modes), reflecting a dominant shallow mode with an extended tail. Transform boundaries (1.429 bits, 2.7 modes) and continental rifts (0.666 bits, 1.6 modes) are the simplest, with nearly all seismicity concentrated in the uppermost crust.

The ratio of highest to lowest entropy exceeds 6:1, and the ratio of effective modes spans an order of magnitude (16.2 vs. 1.6). Shannon entropy thus provides a single scalar that encodes the geological complexity of a tectonic environment: the number of physically distinct depth regimes that contribute meaningfully to its seismicity.

3.3 KL Divergence Matrix

The full pairwise KL divergence matrix is presented in Table 3. Because KL divergence is asymmetric, rows represent the reference distribution ( P ) and columns the approximating distribution ( Q ), such that entry ( (i, j) ) gives ( D_{\mathrm{KL}}(P_i | Q_j) ).

Table 3. KL divergence matrix (bits). Rows = ( P ), columns = ( Q ).

The largest divergences involve transform boundaries, which are maximally distinct from subduction (3.583 bits in the forward direction) and from collision zones (3.198 bits). The smallest divergences are between subduction and collision (0.310–0.374 bits), reflecting their shared deep seismicity, and between ridges and rifts (0.423–0.898 bits), reflecting their shared shallow concentration.

3.4 JSD Matrix with Bootstrap Confidence Intervals

Table 4 presents the symmetric JSD for all pairwise comparisons, together with the JSD distance ( \sqrt{\mathrm{JSD}} ) and 95% bootstrap confidence intervals.

Table 4. Jensen-Shannon divergence and JSD distance between tectonic settings, with 95% bootstrap CIs.

The mean between-setting JSD is 0.321 bits. The most divergent pair is subduction vs. transform (JSD = 0.569 bits; 95% CI: [0.545, 0.584]), while the most similar pair is subduction vs. collision (JSD = 0.077 bits; 95% CI: [0.073, 0.091]). These confidence intervals are tight and non-overlapping, confirming that the ordering of pairwise divergences is robust.

3.5 Within- vs. Between-Setting Divergence

To test whether the five-category tectonic classification is supported by depth distribution structure, within-setting JSD was computed for settings with multiple sampled regions (Table 5).

Table 5. Within-setting JSD (bits) for individual region pairs.

The mean within-setting JSD is 0.065 bits. Compared to the mean between-setting JSD of 0.321 bits, the ratio is 4.98 (approximately 5.0). This five-fold separation indicates that tectonic settings are genuinely distinct clusters in depth-distribution space, and that the conventional tectonic classification captures a real information-theoretic boundary.

At the individual-region level, the most similar pair is the East Pacific Rise vs. the Mid-Atlantic Ridge (JSD = 0.035 bits), two spreading ridges separated by thousands of kilometers yet nearly indistinguishable in their depth structure. The most dissimilar pair is Indonesia Subduction vs. the San Andreas Transform (JSD = 0.603 bits), reflecting the extreme contrast between deep slab seismicity and shallow crustal faulting.

3.6 Magnitude Dependence

To assess whether depth distribution divergence is sensitive to earthquake magnitude, JSD between subduction zones and spreading ridges was computed across five magnitude bands (Table 6).

Table 6. JSD (bits) between subduction and spreading ridge depth distributions by magnitude range.

The JSD declines weakly with increasing magnitude, from 0.194 bits for M 4.0–4.5 events to 0.113 bits for M 6.0+ events. The overall structure of the divergence is preserved across all magnitude ranges, indicating that depth distribution differences between tectonic settings are a stable geometric property of the tectonic environment and not an artifact of the magnitude threshold or event selection. The modest decline at larger magnitudes may reflect the tendency of the largest earthquakes to nucleate at similar depths (the seismogenic zone) regardless of tectonic setting.

3.7 KL Asymmetry Analysis

The asymmetry of KL divergence carries physical meaning. Two illustrative cases are examined:

Subduction vs. spreading ridge. ( D_{\mathrm{KL}}(\text{Subduction} | \text{Ridge}) = 0.964 ) bits, while ( D_{\mathrm{KL}}(\text{Ridge} | \text{Subduction}) = 0.775 ) bits. It costs more information (0.964 bits) to encode the subduction depth distribution using a ridge-optimized code than the reverse (0.775 bits). This is because the subduction distribution places substantial probability mass at intermediate and deep bins (100–700 km) where the ridge distribution assigns very low probability. In the reverse direction, the ridge distribution's dominant shallow mode is reasonably well covered by the subduction distribution's own shallow component, so the penalty is smaller.

Subduction vs. transform. ( D_{\mathrm{KL}}(\text{Subduction} | \text{Transform}) = 3.583 ) bits vs. ( D_{\mathrm{KL}}(\text{Transform} | \text{Subduction}) = 3.300 ) bits. The asymmetry is again consistent: the distribution with greater depth extent (subduction) incurs a higher forward KL cost when approximated by the shallower distribution, because the shallower code wastes no bits on deep bins that it never uses. In general, the KL asymmetry is directionally aligned with the entropy ordering: higher-entropy distributions incur a larger information penalty when their structure is compressed through approximation by a lower-entropy model.

4. Discussion

4.1 Information-Geometric Interpretation

The results establish that tectonic settings can be meaningfully represented as points in an information-geometric space defined by earthquake depth distributions. The JSD distance ( d(P, Q) = \sqrt{\mathrm{JSD}(P, Q)} ) satisfies the axioms of a metric — non-negativity, identity of indiscernibles, symmetry, and the triangle inequality — enabling rigorous geometric reasoning about tectonic relationships. In this space, subduction and continental collision zones are nearest neighbors (JSD distance = 0.277), consistent with their shared involvement of convergent plate motion and deep seismicity. Transform boundaries and continental rifts, both dominated by shallow seismicity, occupy an intermediate distance from each other (0.647) but are far from subduction (0.754 and 0.690, respectively). This geometric arrangement recapitulates the standard tectonic taxonomy from a purely distributional, model-free perspective.

4.2 Physical Interpretation of Entropy Ordering

The Shannon entropy ordering — subduction (4.02 bits) > collision (3.62 bits) > ridge (2.17 bits) > transform (1.43 bits) > rift (0.67 bits) — corresponds directly to the geological complexity of each setting. Subduction zones generate earthquakes through multiple physical mechanisms operating at different depths: shallow interplate thrust faulting, intermediate-depth dehydration embrittlement or thermal shear instabilities within the slab, and deep-focus events driven by phase transformations in the transition zone (Frohlich, 2006). Each mechanism contributes a distinct mode to the depth distribution, and the 16.2 effective modes recovered by the entropy analysis capture this multiplicity quantitatively.

Continental collision zones share much of this complexity (12.3 effective modes) because they involve underthrusting of continental lithosphere to significant depths, albeit without the deepest seismicity characteristic of oceanic subduction. Spreading ridges, with 4.5 effective modes, have a dominant shallow mode from normal faulting along the ridge axis but also a secondary deeper mode reflecting off-axis seismicity and transform offsets. Transform boundaries (2.7 modes) and continental rifts (1.6 modes) are the simplest: seismicity is controlled almost entirely by brittle failure in the upper crust, and the depth distribution approximates a single narrow peak.

The effective number of modes ( 2^H ) thus serves as a direct proxy for the number of physically distinct seismogenic depth regimes in a tectonic setting. This represents a new way to quantify geological complexity using a single, well-founded scalar derived from information theory.

4.3 KL Asymmetry and Geological Implications

The asymmetry of KL divergence is not merely a mathematical inconvenience to be symmetrized away; it encodes genuine directional information about the relative structure of tectonic environments. The consistent finding that ( D_{\mathrm{KL}}(P_{\text{complex}} | Q_{\text{simple}}) > D_{\mathrm{KL}}(Q_{\text{simple}} | P_{\text{complex}}) ) means that it is costlier to approximate a complex, high-entropy distribution using a simple, low-entropy code than the reverse. Physically, this implies that mistaking a subduction zone for a spreading ridge discards more information about depth structure than mistaking a ridge for a subduction zone. The subduction distribution contains deep modes that the ridge code has no capacity to represent, whereas the ridge distribution's shallow concentration falls within the subduction distribution's broad support.

This directional cost has practical implications for seismic hazard analysis. If a seismicity model calibrated on spreading ridge data were applied to a subduction zone, it would systematically underestimate the probability of intermediate and deep earthquakes, incurring an information penalty of 0.964 bits. The reverse misapplication — a subduction model applied to a ridge — incurs a lower cost (0.775 bits) because the subduction model, while overly complex, at least covers the ridge's shallow depth range. KL divergence thus quantifies the asymmetric risk of model misapplication in a manner that summary statistics cannot capture.

4.4 Within- vs. Between-Setting Ratio Validates Tectonic Classification

The 5.0-fold ratio between mean between-setting JSD (0.321 bits) and mean within-setting JSD (0.065 bits) provides strong information-theoretic evidence that the conventional five-category tectonic classification reflects genuine structure in earthquake depth distributions. This ratio is analogous to the F-statistic in analysis of variance, measuring the extent to which between-group variation exceeds within-group variation, but applied in distribution space rather than to scalar means.

The within-setting divergences are informative in their own right. Within subduction zones, the JSD between Japan and Indonesia (0.044 bits) is remarkably small given the geographic separation of these regions, suggesting that the global subduction process imposes a characteristic depth distribution that is largely independent of local variations in slab age, convergence rate, or geometry. The slightly larger divergence between Japan and South America (0.124 bits) may reflect differences in slab dip or the presence of flat-slab segments in the South American setting. Within spreading ridges, the East Pacific Rise and Mid-Atlantic Ridge are nearly indistinguishable (JSD = 0.035 bits), consistent with the hypothesis that oceanic ridge seismicity is governed by a universal thermal and rheological profile.

4.5 Magnitude Invariance of Divergence

The near-invariance of JSD across magnitude ranges (0.113 to 0.194 bits for subduction vs. ridge) demonstrates that depth distribution divergence is a structural property of tectonic settings, not an artifact of the magnitude threshold or catalog completeness. The modest decline in JSD for the largest events (M 6.0+) likely reflects two factors: first, the largest earthquakes in both settings preferentially nucleate within the seismogenic zone at similar shallow depths; second, deep-focus events (exclusive to subduction zones) are predominantly of moderate magnitude, so their relative contribution diminishes when only the largest events are considered. The persistence of substantial divergence even at M 6.0+ (JSD = 0.113 bits) confirms that tectonic depth signatures remain detectable across the magnitude spectrum.

4.6 Comparison to Descriptive Approaches

Traditional approaches to characterizing earthquake depth rely on mean or median depth values, which compress the full distributional shape into a single number and discard information about multimodality and tail behavior. For example, the mean depths of subduction zones (75.5 km) and continental collision zones (66.7 km) differ by only 8.8 km, yet their JSD (0.077 bits) quantifies a real distributional difference attributable to the deeper tail in subduction zones. Conversely, the mean depths of spreading ridges (24.1 km) and continental rifts (10.6 km) differ by 13.5 km — more than the subduction-collision gap — yet the ridge-rift JSD (0.125 bits) is larger by nearly a factor of two. Information-theoretic measures thus capture distributional features — breadth, tail weight, multimodality — that are invisible to mean-based comparisons.

4.7 Limitations

Several limitations should be noted. First, the geographic bounding-box classification is necessarily approximate; regions near boundaries may be misclassified, and some areas host overlapping tectonic processes. More refined classification based on focal mechanisms or plate boundary models could improve precision. Second, the 10-km bin width affects absolute divergence and entropy values; finer binning increases structural resolution but also noise in sparsely sampled depth ranges. Third, the 2000–2024 temporal window may be affected by non-stationarity following great earthquakes or volcanic episodes. Fourth, catalog depth uncertainties vary with magnitude and network coverage, and fixed-depth solutions for shallow events may artificially sharpen shallow modes. Finally, finer tectonic subdivisions — such as distinguishing steep from flat subduction or slow from fast spreading — may reveal additional information-theoretic structure.

5. Conclusion

This study establishes an information-theoretic framework for quantifying the divergence between earthquake depth distributions across tectonic settings. Three principal findings emerge from the analysis of 340,476 M 4.0+ earthquakes recorded from 2000 to 2024.

First, tectonic settings occupy well-separated positions in an information-geometric space defined by the JSD metric. The mean between-setting JSD of 0.321 bits exceeds the mean within-setting JSD of 0.065 bits by a factor of 5.0, providing strong evidence that the conventional tectonic classification captures genuine distributional structure.

Second, KL divergence asymmetry encodes directional information about the relative complexity of tectonic environments. The information cost of approximating a high-entropy distribution (subduction, collision) with a low-entropy code (transform, rift) consistently exceeds the reverse, reflecting the richer depth structure — and greater number of seismogenic mechanisms — in convergent settings.

Third, Shannon entropy maps directly to geological complexity through the effective number of depth modes ( 2^H ). The entropy ordering from subduction (4.02 bits, 16.2 modes) through collision (3.62 bits, 12.3 modes), ridge (2.17 bits, 4.5 modes), transform (1.43 bits, 2.7 modes), to rift (0.67 bits, 1.6 modes) recapitulates the known hierarchy of tectonic complexity from a purely information-theoretic perspective.

These measures are robust across magnitude ranges and statistically well constrained by bootstrap confidence intervals. The framework extends naturally to temporal monitoring of depth distribution changes, comparison of seismicity models, and information-theoretic classification of newly observed seismic sequences. By treating depth distributions as elements of a metric space, information geometry opens a new avenue for quantitative plate tectonics.

References

Frohlich, C. (2006). Deep Earthquakes. Cambridge University Press.

Gutenberg, B., and Richter, C. F. (1954). Seismicity of the Earth and Associated Phenomena (2nd ed.). Princeton University Press.

Kullback, S., and Leibler, R. A. (1951). On information and sufficiency. Annals of Mathematical Statistics, 22(1), 79–86.

Lin, J. (1991). Divergence measures based on the Shannon entropy. IEEE Transactions on Information Theory, 37(1), 145–151.

Shannon, C. E. (1948). A mathematical theory of communication. Bell System Technical Journal, 27(3), 379–423.

Syracuse, E. M., and Abers, G. A. (2006). Global compilation of variations in slab depth beneath arc volcanoes and implications. Geochemistry, Geophysics, Geosystems, 7(5), Q05017.

USGS Earthquake Hazards Program. FDSN Web Service. https://earthquake.usgs.gov/fdsnws/event/1/