Robert J. Skoumal

Interactions among faulting, earthquakes, and eruptions are fundamental to plate tectonics and hazard forecasting yet rarely observed along mid‐ocean ridges. On Iceland's Reykjanes Peninsula, seismotectonic–volcanic unrest resumed after nearly 800‐year hiatus, providing an opportunity to observe these interactions during 2021–2025 activity. By integrating high‐resolution seismicity, focal mechanisms, satellite geodesy, surface deformation, and eruption data, we document ∼4 m of total extension accommodated through 14 rifting episodes. The largest, in 2023, involved graben reactivation and diking, with seismic swarms and earthquake faulting that matched the surface ruptures, where strike‐slip faulting preceded normal‐faulting earthquakes and extension. The accrued extension was released by extension fractures triggered by magma accumulation. Long‐term observations show no correlation between erupted magma volume, seismicity, and crustal extension. This highlights dynamic relation between rifting, faulting, and magmatism in transtensional settings and their implications for hazard assessment.

We would like to correct a mistake in equation (1) of the Skoumal et al. (2026) article. The equation for parameter ϕ following Angelier (1979) and the

We introduce pySATSI, a Python package for computing earthquake focal mechanism stress inversions. This algorithm can handle a wide variety of types of stress inversion problems with a single script and can duplicate many capabilities of preceding methodologies. We also add new capabilities that include spatiotemporally variable inversion grids, damped stress estimates for clusters with few or no focal mechanisms, and variable fault-plane ambiguities that the user can assign to individual events. In addition, we added the ability to use damped stress inversions with fault-plane ambiguity probabilities that are weighted by fault instabilities. Our algorithm is computationally efficient with faster runtimes than previous algorithms, scales well for large datasets, and can be easily parallelized.

Accurate tectonic models are essential for assessing seismic hazard and fault interactions. However, the plate configuration at the complex Mendocino triple junction, where the San Andreas Fault and the Cascadia subduction zone meet, remains uncertain. We analyzed fault slip associated with a recently identified zone of tectonic tremor and low-frequency earthquakes (LFEs) near the southern edge of the subducting Gorda slab. Based on tidal sensitivity and P-wave first motions, we show that the LFEs are generated by dipping, strike-slip motion. This suggests that a former Farallon slab fragment, now captured by the Pacific plate, is translating northward beneath westernmost North America. This geometry effectively extends the slab interface fault, challenging prevailing interpretations of slab window formation and creating a potential unaccounted earthquake hazard in this region.

Distributed acoustic sensing (DAS), which transforms a fiber optic cable into an array of high frequency strainmeters, has the potential to help us characterize earthquakes with a dense sampling of measurements. While earthquake focal mechanisms are frequently determined using P ‐wave polarities and S / P amplitude ratios with inertial seismometers, the dense sampling of DAS over potentially large portions of the focal sphere can aid our solutions. Here, we consider ∼200 regional earthquakes and thousands of S / P measurements on a DAS cable colocated with a network of inertial seismometers near Arcata, California. We demonstrate the S / P ratio measurements made on DAS are similar to those made on the inertial seismometers and can be used to constrain focal mechanism solutions.

Accurately characterizing 3D fault geometry is vital for improving our understanding of earthquake behavior and informing the development of seismic hazard models. Despite their importance, subsurface fault structures tend to be poorly constrained because of limitations in observational data. Improvements to the seismic networks and earthquake detection algorithms have increased the precision and volume of earthquake catalogs, which help illuminate detailed subsurface fault structure and provide the most direct information available about fault geometries at depth. We present a Python package to automate generating 3D fault geometries directly from hypocentral seismicity patterns. This method begins with clustering events based on their spatial density, identifying coherent patterns. Nearby clusters are then merged based on the similarity of their orientations. We fit nonplanar surfaces using support vector regression to balance surface accuracy with minimal deviations from planarity. The fault models are output as quadrilateral meshes at user-defined resolution. In the process of generating the 3D fault surfaces, we compute the spatial density of seismicity around the surface and the planarity as quantitative metrics of the model outputs. As a proof of concept, we apply this approach to the San Andreas–Calaveras fault junction region and the 2019 Ridgecrest earthquake sequence, both in California, which contain complex subparallel faults well defined at the Earth’s surface and abundant microseismicity. These case studies demonstrate the method’s ability to model complex fault structures, including long continuous fault surfaces, crossing faults, variably dipping segments, and subparallel faults. We test the method on both standard network catalogs and double-difference relocated catalogs. We find that our seismicity-based fault model results align with published 3D models that incorporate additional constraints and interpretations (Plesch et al., 2020; Aagaard and Hirakawa, 2021). This workflow provides a low-user-input solution for estimating fault geometries at depth from earthquake catalogs.

From July to October 2022, a noneruptive volcanic earthquake swarm occurred within ∼15 km of Ta‘ū Island, located in eastern American Samoa. Felt reports from residents were the only available information about the swarm when it started, as American Samoa lacked a seismic monitoring network. We developed a consistent single-station catalog for the entire swarm, using seismic data from the nearest station IU.AFI, ∼250 km away. We applied the EQTransformer deep-learning model (Mousavi et al., 2020), automatically picking Pn and Sn arrivals on IU.AFI continuous data. We retained only events with Sn–Pn times of 22.5–25 s, consistent with the expected locations based on felt reports, then detected smaller swarm events with subsequent template-matching. This single-station catalog characterized the swarm’s onset and escalation to peak activity before a multiagency field response team installed a local seismic network in mid-August 2022. This permanent seismic network captured the swarm’s decline. EQTransformer identified short S–P times on the first two locally deployed seismometers, both Raspberry Shake sensors, to constrain the swarm’s distance from Ta‘ū Island. Modern seismological processing methods, combined with basic observations such as felt reports, can quickly contribute useful information during an earthquake response in a poorly monitored region.

The proximity of wastewater disposal to the Precambrian basement is a critical factor influencing induced earthquake rates in the Central United States, but the impact of reducing injection depths has not been widely demonstrated. Beginning in 2015, state regulatory efforts in Oklahoma and Kansas mandated that wells injecting into the lower Arbuckle Group, a basal sedimentary unit, be backfilled with cement (i.e., “plugged back”) so that they inject into shallower formations. This plug back activity gives us a unique opportunity to investigate the relationship between injection depth and induced seismicity rate. To evaluate the impact that decreased injection rates and plug backs had on the seismicity rates, we created a suite of rate–state earthquake models. Observed seismicity rates are best fit when only lower Arbuckle volumes are considered, suggesting the lower Arbuckle injectors were primarily responsible for the seismicity and that plug backs were effective at isolating the injected volumes to shallower formations. Our models demonstrate that if these wells had not been plugged back, seismicity rates would be multiple times larger than they are today. We find that the combination of well plug backs and injection volume decreases can be an effective strategy for reducing induced seismicity rates.

Due to geological complexities and observational gaps, it is challenging to identify the governing physical processes of geothermal field deformation including ground subsidence and earthquakes. In the west and east regions of the Heber Geothermal Field (HGF), decade-long subsidence was occurring despite injection of heat-depleted brines, along with transient reversals between uplift and subsidence. These observed phenomena contradict current knowledge that injection leads to surface uplift. Here we show that high-yield production wells at the HGF center siphon fluid from surrounding regions, which can cause subsidence at low-rate injection locations. Moreover, the thermal contraction effect by cooling increases with time and eventually overwhelms the pressure effects of pressure fluctuation and poroelastic responses, which keep relatively stable during geothermal operations. The observed subsidence anomalies result from the siphoning effect and thermal contraction. We further demonstrate that thermal contraction dominates long-term trends of surface displacement and seismicity growth, while pressure effects drive near-instantaneous changes.

We introduce a Python package for computing focal mechanism solutions. This algorithm, which we refer to as SKHASH, is largely based on the HASH algorithm originally written in Fortran over 20 yr ago. HASH innovated the use of suites of solutions, spanning the expected errors in polarities and takeoff angles, to estimate focal mechanism uncertainty. SKHASH benefits from new features with flexible input formats and allows users to take advantage of recent advances in constraining focal mechanisms for small magnitude or poorly recorded earthquakes. The 3D locations of earthquakes and the velocity models used are varied when finding acceptable solutions. As a result, source–receiver azimuths are reflective of errors from the earthquake locations and velocity models, in addition to the takeoff angles. Users can consider weighted P-wave first-motion polarities derived from traditional or machine-learning picks, cross-correlation consensus, and/or imputation techniques using SKHASH. Focal mechanism solutions can also be further constrained using traditional, machine learning, and/or cross-correlation consensus S/P amplitude ratios. With improved reporting of individual and collective P polarity and S/P amplitude misfits, users can better evaluate the success of the solutions and the quality of the measurements. The reporting also makes it easier to identify potential issues with metadata, including incorrectly reported station polarity reversals. In addition, by leveraging vectorized operations, taking advantage of an efficient backend Python C Application Programming Interface, and the use of a parallel environment, the Python SKHASH routine may compute mechanisms quicker than the HASH routine.