Segmentation of Fault Networks Determined from Spatial Clustering of Earthquakes

We present a new method of data clustering applied to earthquake catalogs, with the goal of reconstructing the seismically active part of fault networks. We first use an original method to separate clustered events from uncorrelated seismicity using the distribution of volumes of tetrahedra defined by closest neighbor events in the original and randomized seismic catalogs. The spatial disorder of the complex geometry of fault networks is then taken into account by defining faults as probabilistic anisotropic kernels, whose structures are motivated by properties of discontinuous tectonic deformation and previous empirical observations of the geometry of faults and of earthquake clusters at many spatial and temporal scales. Combining this a priori knowledge with information theoretical arguments, we propose the Gaussian mixture approach implemented in an Expectation-Maximization (EM) procedure. A cross-validation scheme is then used and allows the determination of the number of kernels that should be used to provide an optimal data clustering of the catalog. This three-steps approach is applied to a high quality relocated catalog of the seismicity following the 1986 Mount Lewis ($M_l=5.7$) event in California and reveals that events cluster along planar patches of about 2 km$^2$, i.e. comparable to the size of the main event. The finite thickness of those clusters (about 290 m) suggests that events do not occur on well-defined euclidean fault core surfaces, but rather that the damage zone surrounding faults may be seismically active at depth. Finally, we propose a connection between our methodology and multi-scale spatial analysis, based on the derivation of spatial fractal dimension of about 1.8 for the set of hypocenters in the Mnt Lewis area, consistent with recent observations on relocated catalogs.

💡 Research Summary

The paper introduces a novel three‑step clustering framework designed to reconstruct the seismically active portion of fault networks using only the spatial locations of earthquakes.
Step 1 – Separation of clustered from background events: The authors compute the volumes of tetrahedra formed by each event’s four nearest neighbors both in the original catalog and in a catalog where event locations have been randomized. In a homogeneous Poisson background the tetrahedron volumes follow a broad distribution, whereas real seismicity exhibits an excess of very small volumes, reflecting spatial clustering. By comparing the two distributions the method flags events that belong to clusters (presumed to lie on faults) and discards the rest as background.
Step 2 – Probabilistic representation of faults: Faults are modeled as anisotropic Gaussian kernels rather than idealized planar rectangles. This choice is motivated by the physical picture that earthquakes nucleate not only on a perfectly defined fault plane but also within a surrounding damage zone of finite thickness. The kernel’s mean vector gives the fault’s central location, while its covariance matrix encodes the orientation, length, and thickness of the fault‑related seismic cloud. Prior geological knowledge (typical fault thickness of a few hundred metres, lengths of kilometres) and information‑theoretic considerations (maximum likelihood, minimum description length) constrain the admissible kernel shapes.
Step 3 – EM mixture estimation and model selection: The catalog is fitted with a mixture of the anisotropic Gaussian kernels using the Expectation‑Maximization (EM) algorithm. EM iteratively updates the probability that each earthquake belongs to each kernel and refines the kernel parameters until convergence. Crucially, the number of kernels K is not preset; a cross‑validation scheme partitions the data into training and validation subsets, evaluates the log‑likelihood penalized for model complexity (e.g., BIC), and selects the K that yields the best predictive performance, thereby avoiding over‑fitting.
Application to the 1986 Mount Lewis sequence: The authors apply the full methodology to a high‑quality relocated catalog of the aftershocks of the M = 5.7 Mount Lewis earthquake in California. The optimal model consists of several planar patches each covering roughly 2 km² with an average thickness of about 290 m. These dimensions are comparable to the size of the mainshock rupture, suggesting that the observed seismicity is not confined to a mathematically sharp fault core but rather occupies a finite‑thickness damage zone that remains seismically active at depth.
Fractal analysis: By estimating the spatial fractal dimension of the hypocenter set, the study finds D ≈ 1.8, a value intermediate between a line (D = 1) and a surface (D = 2). This aligns with recent relocated‑catalog studies and reinforces the view that fault networks possess multi‑scale, non‑Euclidean geometry.
Key contributions and implications:

1. A statistically rigorous, fully automated way to separate clustered seismicity from a Poissonian background using tetrahedral‑volume statistics.
2. A physically motivated anisotropic Gaussian kernel that captures both the planar orientation and finite thickness of fault‑related seismic clouds.
3. Integration of EM mixture modeling with cross‑validation to objectively determine the number of fault‑scale clusters.
4. Empirical evidence that seismicity in the Mount Lewis area is distributed within a damage‑zone thickness of a few hundred metres, challenging the conventional picture of earthquakes occurring strictly on infinitesimally thin fault planes.
5. Quantification of the spatial fractal dimension, supporting the concept of fault networks as multi‑scale structures.
   Overall, the work provides a robust, data‑driven tool for fault network reconstruction that can be applied to other high‑resolution seismic catalogs, offering new insights into fault geometry, damage‑zone activity, and seismic hazard assessment.