Despite considerable progress in seismology, mineral physics, geodynamics, paleomagnetism, and mathematical geophysics, Earth’s inner core structure and evolution remain enigmatic. One of the most significant issues is its thermal history and the current thermal state. Several hypotheses involving a thermally-convecting inner core have been proposed: a simple, high-viscosity, translational mode, or a classical, lower-viscosity, plume-style convection. Here, we use state-of-the-art seismic imaging to probe the outermost shell of the inner core for its isotropic compressional speed and compare it with recently developed attenuation maps. The pattern emerging in the resulting tomograms is interpreted with recent data on the viscosity of iron as the inner core surface manifestation of a thermally-driven flow, with a positive correlation among compressional speed and attenuation and temperature. Although the outer-core convection controls the heat flux across the inner core boundary, the internally driven inner-core convection is a plausible model that explains a range of observations for the inner core, including distinct anisotropy in the innermost inner core.
Welcome to the wonderful world of AxiSEM3D, a flexible, high-performance computational method for solving the elastodynamic wave equation in three-dimensional, hetereogeneous media. It has been used in seismology for the Earth, the Sun, the Moon, and many of other rocky planets and moons of the Solar System. In this manual, we will explore the practicalities of using this code, touching on the mathematical aspects as relevant to the casual user, the input and output files, and pre- and post-processing. This guide is not a replacement for any of the AxiSEM3D papers, which you will find listed in the next section. Nor does it go into any detail about how to install, build, configure, or make the code as these are dealt with in a separate installation guide. However, it does contain everything to do with the scientific aspects of the code and the user interface post-installation.
Peatlands are a major store of soil carbon, due to their high concentration of carbon-rich decayed plant material. Consequently, accurate assessment of peat volumes is important for determining land-use carbon budgets, especially in the Northern hemisphere. Determination of carbon stocks at the scale of individual peat sites has principally relied on either mechanical probing or electromagnetic geophysical methods. In this study, we investigated the use of seismic nodal instrumentation for quantifying peat depth. We used Stryde™ nodes for a deployment at the Whixall moss in Shropshire, England. We measured seismic arrival times from peat-bottom reflections, as well as dispersive surface waves to invert for a model of variable peat depth along a linear cross-section. The use of very small seismic nodes (micronodes) allows for particularly rapid deployment on challenging terrain.
Eikonal tomography has become a popular methodology for deriving phase velocity maps from surface wave phase delay measurements. Its high efficiency makes it popular for handling datasets deriving from large-N arrays, in particular in the ambient-noise tomography setting. However, the results of eikonal tomography are crucially dependent on the way in which phase delay measurements are predicted from data, a point which has not been thoroughly investigated. In this work, I provide a rigorous formulation for eikonal tomography using Gaussian processes (GPs) to smooth observed phase delay measurements, including uncertainties. GPs allow the posterior phase delay gradient to be analytically derived. From the phase delay gradient, an excellent approximate solution for phase velocities can be obtained using the saddlepoint method. The result is a fully Bayesian result for phase velocities of surface waves, incorporating the nonlinear wavefront bending inherent in eikonal tomography, with no sampling required. The results of this analysis imply that the uncertainties reported for eikonal tomography are often underestimated.
<p>The Macquarie Ridge Complex, located at the boundary between Indo-Australian and Pacific plates in the southwest Pacific Ocean, hosts the largest sub-marine earthquakes in the 20<sup>th</sup> century, not associated with ongoing subduction. We deployed 27 ocean-bottom seismometers, of which 15 have been recovered successfully, to understand the origin of the sub-marine earthquakes and their potential earthquake and tsunami hazards to Australia and New Zealand. Additionally, we deployed five land-based seismometers on Macquarie Island.</p><p>We explore state-of-the-art processing methods to analyze the new seismic dataset from the retrieved seismic stations. One of the goals is to image the tectonic settings beneath the MRC. Here, we present a first-order tomographic model and its relevant uncertainty estimate of the region constructed from ambient noise surface waves using a probabilistic inversion framework. The tomographic image will be complemented with receiver-based imaging results such as those from P-wave coda autocorrelations and receiver functions to confirm the existence of possible geometries. The results are expected to supply a fresh understanding of the tectonic settings under the MRC and unpuzzle the origin of the significant underwater earthquakes in the 20th century.</p>
SUMMARYThe spatio-temporal properties of seismicity give us incisive insight into the stress state evolution and fault structures of the crust. Empirical models based on self-exciting point processes continue to provide an important tool for analysing seismicity, given the epistemic uncertainty associated with physical models. In particular, the epidemic-type aftershock sequence (ETAS) model acts as a reference model for studying seismicity catalogues. The traditional ETAS model uses simple parametric definitions for the background rate of triggering-independent seismicity. This reduces the effectiveness of the basic ETAS model in modelling the temporally complex seismicity patterns seen in seismic swarms that are dominated by aseismic tectonic processes such as fluid injection rather than aftershock triggering. In order to robustly capture time-varying seismicity rates, we introduce a deep Gaussian process (GP) formulation for the background rate as an extension to ETAS. GPs are a robust non-parametric model for function spaces with covariance structure. By conditioning the length-scale structure of a GP with another GP, we have a deep-GP: a probabilistic, hierarchical model that automatically tunes its structure to match data constraints. We show how the deep-GP-ETAS model can be efficiently sampled by making use of a Metropolis-within-Gibbs scheme, taking advantage of the branching process formulation of ETAS and a stochastic partial differential equation (SPDE) approximation for Matérn GPs. We illustrate our method using synthetic examples, and show that the deep-GP-ETAS model successfully captures multiscale temporal behaviour in the background forcing rate of seismicity. We then apply the results to two real-data catalogues: the Ridgecrest, CA 2019 July 5 Mw 7.1 event catalogue, showing that deep-GP-ETAS can successfully characterize a classical aftershock sequence; and the 2016–2019 Cahuilla, CA earthquake swarm, which shows two distinct phases of aseismic forcing concordant with a fluid injection-driven initial sequence, arrest of the fluid along a physical barrier and release following the largest Mw 4.4 event of the sequence.
Abstract The structure of the lowermost mantle and the core‐mantle boundary (CMB) has profound implications for Earth's evolution and current‐day dynamics. Whilst tomographic studies of VS show good agreement in the lowermost mantle, consensus as to VP and especially CMB radius has not yet been reached. We perform a hierarchical Bayesian inversion for VP in the lowermost 300 km of the mantle and the radius of the CMB using differential travel time data. Concurrent with finding VP perturbations of 0.56% RMS amplitude that spatially agree with previous studies in areas of low posterior variance, we find 4.5 km RMS amplitude CMB radius perturbations with a broadly north‐south hemispherical character, with spherical harmonic power evenly distributed between degrees 1–3. These results suggest that CMB radial processes are set by a longer scale process than the VP perturbations.
The proliferation of dense arrays promises to improve our ability to image geological structures at the scales necessary for accurate assessment of seismic hazard. However, combining the resulting local high-resolution tomography with existing regional models presents an ongoing challenge. We developed a framework based on the level-set method that infers where local data provide meaningful constraints beyond those found in regional models - e.g. the Community Velocity Models (CVMs) of southern California. This technique defines a volume within which updates are made to a reference CVM, with the boundary of the volume being part of the inversion rather than explicitly defined. By penalizing the complexity of the boundary, a minimal update that sufficiently explains the data is achieved. To test this framework, we use data from the Community Seismic Network, a dense permanent urban deployment. We inverted Love wave dispersion and amplification data, from the Mw 6.4 and 7.1 2019 Ridgecrest earthquakes. We invert for an update to CVM-S4.26 using the Tikhonov Ensemble Sampling scheme, a highly efficient derivative-free approximate Bayesian method. We find the data are best explained by a deepening of the Los Angeles Basin with its deepest part south of downtown Los Angeles, along with a steeper northeastern basin wall. This result offers new progress towards the parsimonious incorporation of detailed local basin models within regional reference models utilizing an objective framework and highlights the importance of accurate basin models when accounting for the amplification of surface waves in the high-rise building response band.
Earthquake ground motion depends strongly on near‐surface structure, which is challenging to image in urban areas at high resolution. Distributed acoustic sensing (DAS) is an emerging technique that provides a scalable solution by converting preexisting fiber‐optic cables into dense seismic arrays. After the July 2019 M7.1 Ridgecrest earthquake, we converted an underground dark fiber across the city of Ridgecrest, CA, into a DAS array. The recorded aftershocks show substantial lateral variability in site amplification over only 8‐km in distance. To understand the cause of such variability, we used three months of continuous data, dominated by traffic‐generated seismic noise, to image near‐surface structure along the fiber path. We find that the lateral variations of earthquake shaking correlate well with the shallow shear velocity model at sub‐kilometer scales, in particular micro‐basins filled with soft sediments. These results highlight the great potential of DAS for high‐resolution seismic hazard mapping in urban areas.
SUMMARY Distributed acoustic sensing (DAS) networks promise to revolutionize observational seismology by providing cost-effective, highly dense spatial sampling of the seismic wavefield, especially by utilizing pre-deployed telecomm fibre in urban settings for which dense seismic network deployments are difficult to construct. However, each DAS channel is sensitive only to one projection of the horizontal strain tensor and therefore gives an incomplete picture of the horizontal seismic wavefield, limiting our ability to make a holistic analysis of instrument response. This analysis has therefore been largely restricted to pointwise comparisons where a fortuitious coincidence of reference three-component seismometers and colocated DAS cable allows. We evaluate DAS instrument response by comparing DAS measurements from the PoroTomo experiment with strain-rate wavefield reconstructed from the nodal seismic array deployed in the same experiment, allowing us to treat the entire DAS array in a systematic fashion irrespective of cable geometry relative to the location of nodes. We found that, while the phase differences are in general small, the amplitude differences between predicted and observed DAS strain rates average a factor of 2 across the array and correlate with near-surface geology, suggesting that careful assessment of DAS deployments is essential for applications that require reliable assessments of amplitude. We further discuss strategies for empirical gain corrections and optimal placement of point sensor deployments to generate the best combined sensitivity with an already deployed DAS cable, from a wavefield reconstruction perspective.
Datasets and associated code for reading the data (and generating tomographic models) for the CMB topography / lowermost mantle tomographic model of Muir et al.
High resolution earthquake hypocentral locations are of critical importance for understanding the regional context driving seismicity. We introduce a scheme to reliably approximate a hypocenter posterior in a continuous domain that relies on recent advances in deep learning. Our method relies on a differentiable forward model in the form of a deep neural network, which is trained to solve the Eikonal equation (EikoNet). EikoNet can rapidly determine the travel-time between any source-receiver pair for a non-gridded solution. We demonstrate the robustness of these travel-time solutions are for a series of complex velocity models. For the inverse problem, we utilize Stein Variational Inference, which is a recent approximate inference procedure that iteratively updates a configuration of particles to approximate a target posterior by minimizing the so-called Stein discrepancy. The gradients of this objective function can be rapidly calculated due to the differentiability of the EikoNet. The particle locations are updated until convergence, after which we utilize clustering techniques and kernel density methods to determine the optimal hypocenter and its uncertainty. The inversion procedure outlined in this work is validated using a series of synthetic tests to determine the parameter optimisation and the validity for large observational datasets, which can locate earthquakes in 439s per event for 2039 observations. In addition, we apply this technique to a case study of seismicity in the Southern California region for earthquakes from 2019.
Output model for Muir, Tanaka and Tkalčić Description of included data files: corrmat.dat - correlation matrix between slowness & radius perturbation coefficients; order of coefficients is (0,0), (1,-1), (1,0)....(8,8) for slowness, and then the same for radiusdrpowers.dat - power per degree l for radius perturbation, column 1 = l, column 2 = mean, column 3 = 5%ile, column 4 = 95%ilevppowers.dat - power per degree l for Vp perturbation, column 1 = l, column 2 = mean, column 3 = 5%ile, column 4 = 95%ile (using 13.61 km/s reference velocity)m_dr.dat - summary statistics for radius perturbation coefficients, column 1 = l, column 2 = m, column 3 = mean, column 4 = 5%ile, column 5 = 95%ilem_ds.dat - summary statistics for slowness perturbation coefficients, column 1 = l, column 2 = m, column 3 = mean, column 4 = 5%ile, column 5 = 95%ilespatialcorr.dat - spatial correlation between slowness and radius, column 1 = latitude, column 2 = longitude, column 3 = correlationtomodata.dat - summary statistics for Vp, column 1 = latitude, column 2 = longitude, column 3 = mean, column 4 = std devtopodata.dat - summary statistics for radius, column 1 = latitude, column 2 = longitude, column 3 = mean, column 4 = std dev Spherical harmonics are given byY^m_l(phi, theta) = sqrt(2)*sqrt(((2l+1)(l-m)!)/(4pi(l+m)!)) cos(m phi) P^m_l(cos(theta)); m>0Y^m_l(phi, theta) = sqrt(2)*sqrt(((2l+1)(l-m)!)/(4pi(l+m)!)) sin(-m phi) P^(-m)_l(cos(theta)); m<0Y^m_l(phi, theta) = sqrt(((2l+1)(l-m)!)/(4pi(l+m)!)) P^(-m)_l(cos(theta)); m=0 where P^m_l is the associated legendre function including the Condon-Shortley phase
The proliferation of large seismic arrays have opened many new avenues of geophysical research; however, most techniques still fundamentally treat regional and global scale seismic networks as a collection of individual time-series rather than as a single unified data product. Wavefield reconstruction allows us to turn a collection of individual records into a single structured form that treats the seismic wavefield as a coherent 3-D or 4-D entity. We propose a split processing scheme based on a wavelet transform in time and pre-conditioned curvelet-based compressive sensing in space to create a sparse representation of the continuous seismic wavefield with smooth second-order derivatives. Using this representation, we illustrate several applications, including surface wave gradiometry, Helmholtz–Hodge decomposition of the wavefield into irrotational and solenoidal components, and compression and denoising of seismic records.
A combined understanding of both the velocity structures of the Lowermost Mantle (LM) and the underlying topography of the Core-Mantle Boundary (CMB) is an essential input into models of the dynamics of the deep Earth. While long-period tomographic Vs models of the lowermost mantle have begun to agree in recent years, there is poorer agreement between Vp models and even less agreement between models of CMB topography, which reflects the sparser sampling and poor sensitivity of most seismic observations to CMB topography, as well as the tradeoff between velocity in the LM and CMB topography. In order to better resolve both of these features, we utilize a joint inversion of LM velocity and CMB topography using a meticulously curated hand-picked dataset of PcP-P, PKPab-PKPbc and P4KP-PcP differential travel time phases. This collection features short-period phases both in CMB reflection and transmission, allowing both velocity and topography to be independently controlled. We utilize hierarchical Hamiltonian Monte Carlo sampling to marginalize over the unknown data error and a priori characteristic size of velocity and topography perturbations, and posterior-predictive cross-validation to determine the highest model resolution that may be supported by the data. We find that the perturbation spectrum of the velocity spectrum is broadband, supporting observations of multiscale thermally dominated features in the LM, while the topography is roughly hemispherical with the spectrum concentrated at lower spherical harmonic degrees. The disconnect between these two spectra suggest that CMB topography is not heavily influenced by isostasy but is instead driven by dynamic topography, with the amplitude of the observed topography not requiring a low viscosity channel in the LM.
Historical seismic data are essential to fill in the gaps in geophysical knowledge caused by the low rate of significant seismic events. Handling historical data in the context of geophysical inverse problems requires special care, due to the large errors in the data collection process. Using Oldham's data for the discovery of Earth's core as a case study, we illustrate how a hierarchical Bayesian model selection methodology using leave-one-out cross validation can robustly and efficiently answer quantitative questions using even poor-quality geophysical data. We find that there is statistically significant evidence for the existence of the core using only the P-wave data that Oldham effectively discarded in his discussion.
Bayesian methods, powered by Markov Chain Monte Carlo estimates of posterior densities, have become a cornerstone of geophysical inverse theory. These methods have special relevance to the deep Earth, where data are sparse and uncertainties are large. We present a strategy for efficiently solving hierarchical Bayesian geophysical inverse problems for fixed parametrizations using Hamiltonian Monte Carlo sampling, and highlight an effective methodology for determining optimal parametrizations from a set of candidates by using efficient approximations to leave-one-out cross-validation for model complexity. To illustrate these methods, we use a case study of differential traveltime tomography of the lowermost mantle, using short period P-wave data carefully selected to minimize the contributions of the upper mantle and inner core. The resulting tomographic image of the lowermost mantle has a relatively weak degree 2—instead there is substantial heterogeneity at all low spherical harmonic degrees less than 15. This result further reinforces the dichotomy in the lowermost mantle between relatively simple degree 2 dominated long-period S-wave tomographic models, and more complex short-period P-wave tomographic models.
SUMMARY Tomography is one of the cornerstones of geophysics, enabling detailed spatial descriptions of otherwise invisible processes. However, due to the fundamental ill-posedness of tomography problems, the choice of parametrizations and regularizations for inversion significantly affect the result. Parametrizations for geophysical tomography typically reflect the mathematical structure of the inverse problem. We propose, instead, to parametrize the tomographic inverse problem using a geologically motivated approach. We build a model from explicit geological units that reflect the a priori knowledge of the problem. To solve the resulting large-scale nonlinear inverse problem, we employ the efficient Ensemble Kalman Inversion scheme, a highly parallelizable, iteratively regularizing optimizer that uses the ensemble Kalman filter to perform a derivative-free approximation of the general iteratively regularized Levenberg–Marquardt method. The combination of a model specification framework that explicitly encodes geological structure and a robust, derivative-free optimizer enables the solution of complex inverse problems involving non-differentiable forward solvers and significant a priori knowledge. We illustrate the model specification framework using synthetic and real data examples of near-surface seismic tomography using the factored eikonal fast marching method as a forward solver for first arrival traveltimes. The geometrical and level set framework allows us to describe geophysical hypotheses in concrete terms, and then optimize and test these hypotheses, helping us to answer targeted geophysical questions.