Geophysical Journal International

SummaryFull waveform inversion (FWI) is a powerful tool in seismic imaging, capable of producing high-resolution models of the subsurface. However, the method remains computationally intensive and sensitive to initial models due to its nonlinearity and ill-posed nature. To quantify uncertainty in FWI results, variational inference (VI) methods, such as Stein Variational Gradient Descent (SVGD), have been increasingly explored. These approaches approximate the posterior distribution by evolving a set of particles using gradient information from the log-posterior. Despite their promise, their effectiveness heavily depends on the quality of the prior used for initialization. In this work, we propose a hybrid framework that improves the efficiency and robustness of VI-based FWI by initializing SVGD with samples drawn from a reconstruction-guided diffusion model. Rather than replacing SVGD with a generative sampler, our approach preserves the theoretical foundations of VI while leveraging the expressive capacity of deep generative models. The diffusion model is trained to generate geologically plausible models conditioned on seismic images, thereby guiding the SVGD initialization toward regions of high posterior support. This initialization significantly reduces the number of required SVGD updates and improves convergence, while keeping the core VI formulation intact. Our results show enhanced posterior approximation and more geologically consistent solutions, with an order of magnitude lower computational cost compared to naïvely initialized SVGD. However, challenges remain, such as the computational demands of likelihood evaluations, the formation of a training set that encompasses all plausible realizations, and sensitivity to reconstruction-guidance weights during sampling. Overall, this method provides a principled and efficient approach to uncertainty-aware FWI, integrating physics-informed inference with data-driven generative modeling for practical applications in full waveform inversion.

SummaryThe degree-2 spherical harmonic coefficients of Earth’s time-variable gravity field are highly sensitive to large-scale mass redistribution within the hydrosphere and cryosphere. Under contemporary global warming, climate-driven mass changes in these reservoirs are a dominant source, yet their individual contributions remain incompletely quantified. Traditional estimates based on hydrospheric models, filtered GRACE spherical harmonic solutions, or GRACE mascon products are limited by incomplete cryospheric representation, spatial leakage, and regularization biases. Here, we apply the Fingerprint Approach that solves the sea-level equation on an elastic Earth to generate geoid fingerprints for four barystatic processes: terrestrial water storage, the Greenland Ice Sheet, the Antarctic Ice Sheet, and mountain glaciers. Using these fingerprints and unfiltered GRACE/GRACE-FO Stokes coefficients for 2003–2024, we reconstruct the individual degree-2 coefficients C20, C21, and S21, along with the associated time series of Earth’s dynamic oblateness (J2) and the mass terms of polar motion excitation ($\chi _{\rm{1}}^{{\rm{mass}}}$, $\chi _{\rm{2}}^{{\rm{mass}}}$). We evaluate the reconstructed contributions against residual geodetic observations from satellite laser ranging and Earth orientation parameters, after removing atmospheric, oceanic, and glacial isostatic adjustment effects. The combined hydrospheric and cryospheric reconstructions reproduce both the secular trends and annual cycles of the residual observed J2, $\chi _{\rm{1}}^{{\rm{mass}}}$, and $\chi _{\rm{2}}^{{\rm{mass}}}$ series. Terrestrial water storage dominates the seasonal variability of J2, $\chi _{\rm{1}}^{{\rm{mass}}}$, and $\chi _{\rm{2}}^{{\rm{mass}}}$, whereas accelerated ice-mass loss from Greenland and Antarctica controls the secular trend, with mountain-glacier mass loss also contributing to the long-term trend in J2. The resulting polar motion excitation drift closely matches the residual geodetic estimate in both magnitude and direction, indicating that contemporary climate-driven mass redistribution can largely account for recent changes in residual geodetic observations, and demonstrating the value of fingerprint-based reconstructions for monitoring climate impacts on the Earth system.

SummaryWe present a forward waveform modeling study to investigate the regional crustal structure of the central-southern Apennines, along a NNE-SSW profile. The profile cross-cuts the Apenninic chain axis and extends from the eastern Adriatic domain, characterized by a thick crust and thick seismogenic layer, to the western Tyrrhenian domain, dominated by tectonic thinning, distributed CO2 gas emissions at the surface and volcanic structures. This region hosted the largest earthquakes in recent history, making precise knowledge of the crustal structure crucial for a comprehensive understanding of seismogenesis and seismic hazard assessment. We analyzed and modelled seismic data from two lower-crustal strike-slip earthquakes in the eastern segment of the profile (2018 Mw 5.1 and 2023 Mw 4.6), recorded by the Italian National Seismic Network (IV). The two events, located along the target NNE-SSW linear profile, provide a unique opportunity to study the westward propagation and evolution of seismic phases. Using a 2D numerical modeling approach, we modelled direct (Sg) and Moho-reflected (SmS) phases on transverse component seismograms, comparing the synthetic to the observed waveforms in terms of arrival times and waveform shapes. A faster Adriatic lower crust with an average shear-wave velocity of 3.85 km/s supports the hypothesis of distributed crustal mafic intrusions at the margin between the Apenninic chain and Adriatic foreland. We estimate a local Adriatic Moho depth of 38 km, in agreement with previous investigations. Furthermore, we identify a strong attenuation zone across the Apenninic chain axis, extending directly from the surface down to 10 km depth, significantly impacting both seismogenic processes and waveform propagation. This first, regional-scale waveform study highlights the significance of waveform analysis for constraining seismic velocities and interface velocity contrasts in the Southern Apennines mountain range.

SummarySeismoelectromagnetic (SE) signals are created in wet poroelastic media by electrokinetic conversions occurring at the pore scale. Two of these signals are frequently investigated: a co-seismic wavefield, bounded to the propagating seismic waves, and an electromagnetic (EM) wave created when a seismic wave passes through the interface between two porous media. Driven by the possibility offered by the SE signals in terms of vertical resolution, compared with seismic reflection methods, a restricted number of authors studied how a thin layer, or a combination of thin layers, can modify SE signals. In this context, the present work presents original experimental and numerical data as a means to further investigate the relation between the thickness of a layer and the EM interface-generated response (IR). In addition, the fact that the IR is very sensitive to contrasts in fluid conductivity, whereas seismic waves are not, is an appealing characteristic of SE exploration. Consequently, we explored the influence of the pore-fluid electric conductivity, combining experimental and numerical methodologies. We identified a strong similarity in how layer thinning affects both seismic waves and IR signals. Moreover, we observed that the maximum enhancement of an IR signal occurs when the layer thickness is approximately half the P-wave wavelength (λP) in the layer. When analysing the influence of fluid conductivity, we showed that the electrokinetic theory adopted in this study provides satisfactory predictions for the waveforms and amplitudes observed experimentally. Finally, we extrapolated our experimental results with field-scale simulations in order to understand how the effects observed experimentally translate to georesources applications.

SummaryOcean-bottom seismometers (OBSs) are reliable instruments to record ground motions and acoustic signals on the sea floor. Precise timing of the data is essential for most seismological analyses. The internal clocks of the OBSs are not GNSS-controlled, so the clock drift mainly caused by ageing of the quartz crystal and temperature effects must be corrected. As part of the BRAVOSEIS experiment, eight OBSs were deployed in the Antarctic Bransfield Strait for 13 months. All OBSs suffered from a large (-15.3 to 5.4 s) and non-linear (-1.2 to -0.6 s residual to linear drift) clock drift. We used noise cross-correlations to determine the clock drift. The parameters for data pre- and post-processing such as filtering and normalisation had to be selected carefully. Overlapping correlation windows were stacked to derive daily Green’s functions. The time shift between consecutive days was calculated and distributed linearly over the data to obtain a continuous data set without gaps or overlaps. Airgun shots were used to constrain the cumulative clock drift of one station without initial synchronisation. Two onshore stations served as GNSS-controlled reference for four OBSs in the northern part of the Central Bransfield Basin. A different noise regime prevailed in the southern part of the basin; therefore, two already corrected OBSs from the northern part were used as reference stations for the southern OBSs. In this way, the clock drift of all OBSs could be corrected accurately.

SummarySpectral induced polarization (SIP) is a promising technique for detecting microbial activity in porous media, yet its interpretation remains limited by the absence of mechanistic models that account for microbial cell structure. In this study, we present a new semi-analytical model for the electrical polarization of microbial cells that treats both the cell plasma and the surrounding medium as electrolytes, and accounts for the cell membrane as well as the influence of the charged surface structures of the cell. We validate our model through numerical simulations based on the Poisson–Nernst–Planck equations. The model builds upon the membrane capacitance model by Sun and Morgan and integrates surface conductivity effects via the O’Konski model and low-frequency polarization using an adapted Dukhin–Shilov approach. The agreement between the numerical results and our new semi-analytical model is good. The model accounts for three dominant polarization mechanisms: (1) diffuse layer polarization at low frequencies (102–104 Hz), (2) membrane-related capacitive effects at intermediate frequencies (105–107 Hz), and (3) Maxwell-Wagner-type polarization at high frequencies (107–109 Hz). In experimental studies, polarization of bacteria typically appears at frequencies around 0.05 and 20 Hz. As the characteristic frequency of polarization processes usually decreases with increasing polarization length scales, the remaining discrepancy between model and experimental observations suggests that measured signals may be influenced by cell aggregates, biofilms, or metabolic byproducts. Our findings provide a foundation for a mechanistic understanding of microbial polarization and highlight the need for future work to extend the model to conglomerates of microbial cells.

AbstractHigh-frequency induced polarisation (HFIP) measurements enable quantification of ground ice content in frozen media by capturing ice relaxation within the frequency range of 1 to 100 kHz. Existing parameterised inversion approaches may bias results by imposing an ice relaxation signature where none exists, assuming a Cole-Cole-type response that may not reflect the true dielectric behaviour of ice, and neglecting low-frequency polarisation. These limitations can lead to high data misfits and ambiguities in interpretation. This study presents an alternative approach that applies independent frequency inversion to directly derive complex resistivity spectra from field measurements, avoiding reliance on predefined models. The resulting inverted spectra provide a representation that more closely captures the true subsurface response. A second, petrophysical, inversion is then performed by fitting a two-component mixture model to the inverted spectra, weighted by the volumetric fractions of its components. One of these components is ice, allowing for the estimation of the volumetric ice content. The approach was applied at Heliport Mire (Abisko, Sweden), a permafrost peatland site, using two complementary profiles: a 50-m 2D profile that captured broad lateral variations of frozen to unfrozen conditions, and an 8-m high-resolution 2D profile that resolved the vertical transition between the upper unfrozen and underlying frozen layers. Independent frequency inversion, across 1 Hz to 57 kHz, successfully produced smooth, coherent spectral responses of true resistivity and phase shift across both profiles. Petrophysical inversion results show diverse conditions along the profile, identifying three distinct zones: ice-rich frozen peat (40-77% ice content), a thawed or degraded peat region (<10% ice content), and unfrozen forest (<5% ice content, effectively representing ice-free conditions). HFIP-derived ice content values were consistent with those derived from laboratory measurements on a permafrost core extracted along the profile. The high-resolution profile distinctly identified the boundary between unfrozen and frozen ground, as confirmed by direct probing measurements. Additionally, the petrophysical model resolves parameters such as shape factor and matrix permittivity, offering further insight into subsurface properties. This methodology advances ground ice characterisation by providing robust quantitative estimates of ice content while retaining spectral information with broader interpretative potential.

SummaryIn the last two decades, the improvement of both instruments and theory, as well as the broadened scope of applications, led to a spectacular development of the use of induced polarization. In particular, the richness of complex conductivity spectra is driving the scientific community towards vast deployment of this measurement method often referred to as Spectral Induced Polarization (SIP). In this contribution, we describe an innovative multichannel instrument that we develop for fast monitoring of critical zone processes. The spectral content of a signal with line spectrum resulting from square-wave current is exploited by injecting successively three square-wave currents with periods of 1, 10 and 100 s, covering the frequency range of 10−2 to 102 Hz in less than four minutes. One dataset consists of eight successive current injections at different depths. For each current injection, the electrical potential is simultaneously measured at seven dipoles. The time-series are recorded with a 2 kHz sampling rate, allowing to calculate by Fourier transform the amplitude and phase spectra up to 1 kHz for each quadrupole. The complex conductivity data was validated by a comparison with the commercial SIP-Fuchs instrument, despite a significant discrepancy below 0.1 Hz which may be due to a worse signal-to-noise ratio at low frequencies. The prototype version of the instrument has been installed in 2018 at a wetland at Ploemeur-Guidel hydrogeological observatory to monitor reactive processes with high spatial resolution across the top meter of soil. The instrumental device, controlled by a Gantner data acquisition system connected to a solar panel, is fully autonomous and consumes little energy. Acquisitions are made several times a day and recorded on a SD card. Seven-year continuous monitoring highlights significant temporal variations of both subsurface resistivity and phase angle. The absence of correlation between resistivity and phase variations in the continuously saturated soil thickness highlights the potential of the system to monitor and separate different types of dynamics processes, such as groundwater/surface water mixing and mineral precipitation/dissolution.

SummaryFrom April until the end of June 2025, we deployed a dense seismic network of 271 three-component stations within an 8 km radius around Lavey-les-Bains, Switzerland, to investigate the structure of the country’s hottest known natural geothermal system. The site hosts a 3 km-deep exploration well (Lavey-1), drilled in 2022, that revealed unexpectedly low flow rates despite temperatures exceeding 120°C, prompting the suspension of the project. The site lies within the narrow Rhône Valley, characterized by steep topography, strong lateral structural heterogeneity, and elevated anthropogenic noise, complicating seismic imaging. The dense nodal array was complemented by a distributed acoustic sensing (DAS) system along a buried telecommunication cable, providing a hybrid dataset suited for passive seismic imaging. We describe the network geometry, instrumentation and deployment logistics; assess data completeness and noise characteristics and present first examples of ambient noise and earthquake recordings. Preliminary analyses demonstrate a high data quality and spatial coverage. This experiment establishes a benchmark dataset for developing advanced passive imaging techniques in complex Alpine environments.

SummaryThe joint use of data from GRACE-like gravity missions and various ocean altimetry missions in a global inversion approach allows to quantify the individual contributions to global and regional sea level budgets. However, the contribution from the Antarctic Ice Sheet (AIS) is subject to large uncertainties mainly depending on the applied strategy to account for effects due to glacial isostatic adjustment (GIA). The large uncertainty of GIA affects estimates of AIS contributions as well as other elements of sea level budgets. Here, we investigate strategies to improve the representation of AIS mass changes within an existing global inversion framework. The framework employs pre-defined, time-invariant spatial patterns, so-called fingerprints, for representing the individual sea-level budget components, including AIS contributions. We improve this inversion method by including additional observations of satellite altimetry over ice sheets, and by further developing the parameterization of AIS ice mass changes. We extend from a basin-wise spatial resolution to a parameterization that resolves time-variable ice mass changes at about 50 km, enabling a better localization of the AIS contributions to global and regional sea level change. From real-data experiments, we obtain ice mass balance estimates that are well within the uncertainty bounds of published reconciled estimates utilizing similar datasets. In particular, inclusion of ice altimetry improves the spatial resolution and at the same time keeps the global inversion results in line with those from regional GRACE analyses. We find differences between inversion results with and without including ice altimetry as an additional observation. These differences are smaller for the time period after 2010 with the availability of CryoSat-2 altimetry having improved sensor technology and high-latitude coverage. This indicates that these differences are caused by ice altimetry errors, whose further characterization and consideration within the estimation remains a future task. Furthermore, the spatial distribution of the differences suggests that they are also related to GIA errors. The improved representation of ice sheets in the global framework developed here provides a prerequisite for working towards minimizing GIA-related errors while assessing the ice sheets’ mass balance.

SummaryGiven the scarcity of seismometers in marine environments, traditional seismology has limited effectiveness in oceanic regions. Submarine Distributed Acoustic Sensing (DAS) systems offer a promising alternative for seismic monitoring in these areas. However, the existing machine learning model trained on land-based DAS data does not perform well with submarine DAS due to differences in noise characteristics, deployment conditions, and environmental factors. This study presents a machine learning approach tailored specifically to submarine DAS data to enable automated seismic event detection and P and S wave identification. Leveraging DeepLab v3, a neural network architecture optimized for semantic segmentation, we developed a specialized model to handle the unique challenges of submarine DAS data. Our model was trained and validated on a dataset comprising nearly 57 million manually and semi-automatically labeled seismic records from multiple globally distributed submarine sites, providing a robust basis for accurate seismic detection. The model adapts to a variety of deployment scenarios and can process DAS data from cables with different lengths, configurations, and channel spacings, making it versatile for various ocean environments. We thus provide an adaptable and efficient tool for automated earthquake analysis of DAS data, which has the potential to enhance real-time earthquake monitoring and tsunami early warning in submarine environments.

SummaryWe investigate upper crustal seismic anisotropy in the southeastern termination of the Zagros Mountains and Qeshm Island through shear wave splitting analysis on the aftershock records of 3 local earthquakes: the 2006 Tiab earthquake (Mw = 6) in the north of the Main Zagros Reverse Fault (MZRF), the 2006 Fin event (Mw = 5.9) in SE Zagros, and the 2005 Qeshm Island event (Mw = 5.8). The results show that north of the MZRF in the Faryab region, the local anisotropy as measured by fast axis orientations of the shear waves, is mainly controlled by crustal structures such as fault shear fabrics or plane fractures around a set of orthogonal active strike-slip faults of a dominantly NE-SW strike. Stress-induced anisotropy caused by fluid-filled microcracks aligned with local maximum horizontal compression seems to be of minor importance. In the vicinity of the MZRF and the Zendan-Minab-Palami Fault, another major fault of the region, local anisotropy is controlled by the structural fabric developed by those faults. Far from the faults, the fast orientations rotate to become aligned with the regional compressional stress. The MZRF is located on an old subduction suture, and it seems that the metamorphic or mylonite fabrics that prevail in this part have shaped the anisotropic structure in the upper crust. In the Simply-Folded Belt of the Zagros and in Qeshm Island, where the deformation of the upper crust is younger, local anisotropy is mostly stress-induced, and the role of fault structures, even though the region is affected by extensive and complex sets of active faulting, is relatively minor. The intensity of anisotropy, as witnessed by splitting delay times, decreases from the metamorphic belt north of the MZRF, to the deformation front in Qeshm Island. This observation is in agreement with the decrease in the accumulated deformation from north to south everywhere in the Zagros. Crustal anisotropy in the Zagros does not show large-scale uniform patterns, rather, it varies over relatively small distances. This is partly due to the fact that both stress-induced and structurally-controlled anisotropies are at work, and partly because of the short-distance variations of the modern stress field in the Zagros.

SummaryWe present a systematic approach to optimise distributed acoustic sensing (DAS) fibre-optic cable layouts using global optimisation techniques. Our method represents cable geometries using splines, enabling efficient exploration of layouts while respecting physical deployment constraints. The use of evolutionary algorithms enables single and multi-objective optimisation, taking into account complex design constraints such as terrain, accessibility, exclusion zones, cable length, and coupling-related or local site effects, while allowing efficient parallelisation of the optimisation process. We demonstrate the approach on a real-world case study, optimising the layout of a DAS cable for monitoring slope stability in the Cuolm da Vi area of Switzerland. We adapt design criteria for seismic source location problems, and for ambient noise surface wave tomography, to account for the unique characteristics of DAS, such as directional sensitivity patterns. The results show significant potential for improvements in source location accuracy and surface wave tomographic resolution by optimising cable layouts, highlighting the potential of this approach for optimising DAS deployments in various geophysical applications.

SummaryIn this study, we present a new integrated experimental approach to investigate simultaneously the electrical spectral induced polarisation (SIP), mechanical, hydraulic and chemical properties of synthetic clayey soil mixed with different types and quantities of organic matter. It addresses knowledge gaps that aim to advance SIP as a non-destructive analysis tool for soils. We used an inorganic clay as a proxy for clayey soil with a moderate cation exchange capacity to achieve more realistic test conditions since most studies use sand. Three organic matter (OM) types with contrasting properties, biosolids, peat and sugar cane residue, broaden the range of organic carbon materials that have been tested. Our study demonstrates a strong relationship between the imaginary part of the complex conductivity and the total organic carbon content of the soil-OM mixtures. It indicates that the relationships depend on the degree of aromaticity, with the slope angle increasing as the degree of aromaticity of the OM increases. Hence, the quantity of OM, as well as its chemical structure, plays a key role in SIP response. Interestingly, these relationships are independent of soil water saturation and bulk density. These findings are of paramount importance for enabling field-scale applications and confirm the potential of SIP as a non-invasive tool for monitoring and characterising soil in situ.

SummaryThis paper is the second part of a series examining the effects of ground polarization in airborne electromagnetic (AEM) data collected with fixed-wing platforms. Induced polarization (IP) effects can be detected using airborne electromagnetic methods; however, most geophysical studies have focused on helicopter-borne systems whose sensitivity to subsurface polarizable features is well established. In contrast, the potential of fixed-wing AEM systems for IP detection remains largely unexplored, and their effects have not yet been modelled. Building on Part A of this series, which examined the sensitivity of TEMPEST™ system to ground chargeability with numerical analysis and dataspace inspection, we extend the study using field survey data to model subsurface IP effects in inversion. This study is defined at three different exploration scales: deposit scale, survey-line and regional scale. The first experiment focuses on a comparative modelling analysis between the TEMPEST™ and SkyTEM312FAST helicopter-borne system along two overlapping survey lines. The results show highly comparable chargeability and resistivity distributions, with consistent outcomes across the TEMPEST™ measured components (X and Z) and with geological interpretation of the area. These findings demonstrate that fixed-wing AEM can effectively resolve IP anomalies with resolution and depth penetration similar to helicopter-borne systems, despite differences in acquisition geometry and system design. Then, to assess regional-scale applicability, the entire Musgrave Province in South Australia was inverted incorporating IP effects and comparing the results with the non-IP modelling of the area. The IP modelling shown a systematically reduction of inversion misfit, when compared with non-AIP modelling with differences between the resistivity models higher than 100%. To conclude, the ground truthing of regional modelling has been carried over the well-characterized Nemo-Babel mineralization. This confirmed that TEMPEST™ derived chargeability anomalies align closely with known mineralized zones, validating both spatial accuracy and correspondence with mineralization of the modelled resistivity and chargeability. Overall, this study demonstrates that fixed-wing AEM platforms, such as TEMPEST™, can detect and quantify ground chargeability from regional to deposit scale, providing a valuable tool to target exploration and to characterize mineralized bodies.

SummaryLandslide disasters are typically triggered by various environmental factors, making it crucial to understand the interaction between subtle internal changes and these factors for accurate risk assessment. Noise-based velocity change measurement offers a promising tool, yet its widespread application is limited by the inherent instability of noise sources, constraining temporal resolution. Here, we employ an wave-packet-based nine-component spatial stacking approach with a dense seismic array deployed at the Xishan Village landslide. This advancement allows for the extraction of extraction of high temporal velocity change (20-minute) at different frequencies, enabling four-dimensional dynamic analysis of landslide internal changes. Our findings reveal complex spatial distributions of velocity changes influenced by solar thermal radiation and rainfall at different locations and depths. Notably, during rainfall of approximately 20 mm, the observed maximum velocity reduction correlates closely with a fracture zone at ∼8 m depth, suggesting that pre-existing deformation structures significantly enhance local permeability, and promote the now deeper rainwater infiltration. This infiltration leads to increased pore pressure and velocity reduction. These results highlight the ambient noise method potential for urban landslide monitoring, providing technical support for early warning and risk assessment.

SummaryTeleseismic P-wave receiver function (RF) analysis is a powerful tool for probing deep structures; however, its application to regions with sedimentary cover remains challenging. The interference by converted and reflected waves related to the sediments can form strong reverberations and render individual phases unidentifiable, complicating the investigation of sedimentary structures. Moreover, the strong, long-lasting sediment-induced waves can mask seismic phases from deeper layers, making it difficult to investigate the underlying layers. We investigated all the phases in RFs and establish simple relationships between the underlying layer structure and the phases and phase groups. By measuring the times of all phases and phase groups or stacking the amplitudes of all phases, we can investigate the structure more reliably; however, when interference is so strong that the peaks and troughs do not align with individual phases, the models resulting from using these phases may contain significant uncertainty. In these cases, simple RF waveform fitting constrained by the one-way vertical travel time (the product of wave vertical slowness and layer thickness) or the period of the phase groups associated with the sedimentary layer can resolve the sedimentary structure reliably. Deconvolving the P-wave RFs predicted by the resulting sedimentary model from the original RF can remove the sedimentary waveforms effectively and correct the sediment-induced delays in deeper-layer phases. We demonstrate the effectiveness of our proposed approaches for inverting the structure and removing the sedimentary response using both synthetic and real data.

SummaryGlobal Navigation Satellite System (GNSS) plays a fundamental role in monitoring time-dependent ground displacement. However, GNSS daily position time series can often contain significant outliers, reaching up to several centimeters. These are likely of non-tectonic origin, and, if not properly accounted for, they can significantly impact the accuracy and dependability of the estimation of key parameters for geophysical analyses, such as long-term velocities and transient deformations. Characterizing these outliers can provide information about their possible sources and help us implementing mitigation strategies. Asymmetric outliers, i.e., those characterized by a primary direction, therefore occurring on one side of the mean time series, are of particular interest since they could point to the presence of recurring or repeatable sources of error. Their key features and potential causes are, however, still not fully analyzed and understood. We analyze asymmetric outliers in thousands of GNSS time series across three regions – Central-Southern Italy, New Zealand, and the Western U.S. – using data from multiple processing centers, and we reveal some persistent features among all datasets. Tens of the analyzed sites exhibit hundreds of large outliers (10-50 mm), far exceeding typical position uncertainties (∼1-6 mm). Remarkably, the outliers are numerous in the horizontal component, and tend to occur near mountainous regions, with preferred direction roughly orthogonal to the local topography. The results consistency across different datasets and instrumental features suggest a physical origin for these outliers rather than a specific processing approach or instrumental configuration. Further analyses at local scales align with previous studies linking skewed position errors to uncorrected tropospheric delays driven by the coupling between atmospheric conditions and local terrain (e.g., trapped lee waves). However, other factors – such as multipath, snow accumulation on GNSS antennas or obstructed sky visibility – could also contribute to the observed asymmetric outliers. We explore mitigation strategies at both processing and post-processing stages, but further analyses and more sophisticated approaches, such as high-resolution tropospheric modeling, are needed to better understand the involved processes and achieve meaningful improvements.

SummaryThe Gravity Recovery and Climate Experiment (GRACE) and its Follow-On mission provide essential observations of Earth’s surface mass redistribution. However, inherent north-south striping noise in the GRACE spherical harmonic (SH) products limits their application at sub-basin scales. To address this, we introduce a novel spatial domain decorrelation filter, the Physical-Informed Spatial Pattern (PISP) filter, which leverages the structured physical characteristics of the noise for its precise identification and removal. Comprehensive numerical experiments validated that PISP effectively eliminates striping noise globally and yields a consistent noise background across latitudes, with noise reduced to a uniform level in more than 90% of the months examined and with stable performance under strong-noise conditions. In a case study of water storage variations in Lake Victoria, PISP preserves the primary signal amplitude and reduces the root-mean-square error relative to reference data to 5.84 cm after spatial smoothing, outperforming the 6.81 cm achieved by the MVMDS + DDK6. Furthermore, for three earthquakes with magnitudes exceeding 8.8, PISP effectively removes striping noise using regional masking, successfully recovering the co-seismic signal morphology. By further verifying the method’s stability across various noise scenarios, the results demonstrate PISP’s potential for future global research integrating multi-satellite gravity data.

SummarySlowslip events (SSEs) modulate the earthquake cycle in subduction zones, yet understanding their physics remains challenging due to sparse observations and high computational cost of physicsbased simulations. We present a scientific machine-learning approach using a data-driven reduced-order modelling (ROM) framework to efficiently simulate the SSE cycle governed by rate-and-state friction in a Cascadia-like 2D subduction setting. Our approach projects fault slip, sliprate, and state-variable trajectories onto a splinebased latent space, which is subsequently emulated using properorthogonal decomposition and radialbasisfunction interpolation. Achieving a speedup of ∼360, 000 × compared to volumetric simulations, the ROMs enable comprehensive parameter exploration and Bayesian Markov chain Monte Carlo (MCMC) inversion. By smoothly interpolating between the physics-based simulations, the ROMs reveal complex dependencies that might be overlooked with coarser parameter space sampling. Our analysis reveals complex, non-linear dependencies of SSE characteristics on the width and magnitude of the deep, low-effective-normal-stress region while holding frictional parameters constant. We show that, to first order, the recurrence time of SSEs is linearly dependent on the normalized fault width, defined as the SSE zone width divided by the characteristic nucleation length, in agreement with previous work. However, at second order, the recurrence interval can change more rapidly with small variations in the SSE zone width. We identify a region of steep, non-linear dependence of the recurrence interval on the normalized fault width, which we attribute to the extent of the velocity-weakening portion of the subduction interface that produces SSEs. Our MCMC inversion constrained by Northern Cascadia SSEs observations indicates near-lithostatic pore fluid pressure (99.6 ± 0.17% lithostatic) and positions the upper frictional transition zone at 30.4 ± 2.8 km depth, consistent with geophysical observations. The inversion resolves the deep SSE-portion of the slab spanning 45 ± 16 km with low effective normal stress of 3.8 ± 1.4 MPa. We discuss how varying the fault frictional parameters affects the MCMC-inverted parameter values. This framework provides a new tool for advancing the physics-based understanding of SSEs and subduction zone faulting mechanics. By systematically linking megathrust properties such as fluid pressure and fault strength to rate-and-state friction governed slow slip cycle characteristics, such as recurrence interval, our approach helps to constrain the first and second-order physics-based controls and the uncertainties of how subduction zones slip.