Search

We present a theoretical framework that links Fermat's principle of least time to optimal transport theory via a cost function that enforces local transport. The proposed cost function captures the physical constraints inherent in wave propagation; when paired with specific mass distributions, it yields shortest paths in the considered media through the optimal transport plans. In the discrete setting, our formulation results in physically significant optimal couplings, whose off-diagonal entries identify shortest paths in both directed and undirected graphs. For undirected graphs with positive edge weights, commonly used to parameterize seismic media, our method provides solutions to the Eikonal equation consistent with those from the Dijkstra algorithm. For directed negative-weight graphs, corresponding to transportation cost matrices with negative entries, our approach aligns with the Bellman–Ford algorithm but offers considerable computational advantages. We also highlight potential research directions. These include the use of sparse cost matrices to reduce the number of unknowns and constraints in the considered transportation problem, and solving specific classes of optimal transport problems through the Dijkstra algorithm to enhance computational efficiency. [ABSTRACT FROM AUTHOR]

Concerns raised by Okazaki & Ueda (2022) on the paper by Sambridge et al. (2022) are addressed. Two issues are discussed and some new numerical results presented. The first concerns whether the properties of the Wasserstein time-series misfit introduced in our earlier paper will translate to model space non-uniqueness in a seismic waveform inversion setting. It is argued that this is unlikely, given the special conditions, which must exist between all observed/predicted seismic waveform pairs for non-uniqueness to result. The second issue discussed is the efficacy of using the Sliced Wasserstein algorithm of Bonneel et al. (2015) as an alternate to the marginal Wasserstein algorithm, as proposed by Okazaki & Ueda (2022). It is argued that for optimization-based waveform fitting, the Sliced Wasserstein algorithm is a viable alternate provided care is taken to ensure that conditions arise which do invalidate analytical derivative expressions of the resulting Wasserstein misfit. In practice, this would likely mean recasting the 2D Optimal Transport problem posed in our earlier paper onto unstructured grids. [ABSTRACT FROM AUTHOR]

We propose a new approach to measuring the agreement between two oscillatory time-series, such as seismic waveforms, and demonstrate that it can be used effectively in inverse problems. Our approach is based on Optimal Transport theory and the Wasserstein distance, with a novel transformation of the time-series to ensure that necessary normalization and positivity conditions are met. Our measure is differentiable, and can readily be used within an optimization framework. We demonstrate performance with a variety of synthetic examples, including seismic source inversion, and observe substantially better convergence properties than achieved with conventional L 2 misfits. We also briefly discuss the relationship between Optimal Transport and Bayesian inference. [ABSTRACT FROM AUTHOR]

We develop a theoretical framework for framing and solving probabilistic linear(ized) inverse problems in function spaces. This is built on the statistical theory of Gaussian Processes, and allows results to be obtained independent of any basis, avoiding any difficulties associated with the fidelity of representation that can be achieved. We show that the results of Backus–Gilbert theory can be fully understood within our framework, although there is not an exact equivalence due to fundamental differences of philosophy between the two approaches. Nevertheless, our work can be seen to unify several strands of linear inverse theory, and connects it to a large body of work in machine learning. We illustrate the application of our theory using a simple example, involving determination of Earth's radial density structure. [ABSTRACT FROM AUTHOR]

By starting from a general framework for probabilistic continuous inversion (developed in Part I) and introducing discrete basis functions, we obtain the well-known algorithms for probabilistic least-squares inversion set out by Tarantola & Valette. In doing so, we establish a direct equivalence between the spatial covariance function that must be specified in continuous inversion, and the combination of basis functions and prior covariance matrix that must be chosen for discretized inversion. We show that the common choice of Tikhonov regularization (⁠|$\mathbf {C_m^{-1}} = \sigma ^2\mathbf {I}$|⁠) arises from a delta-function spatial covariance, and that this lies behind many of the artefacts commonly associated with discretized inversion. We show that other choices of spatial covariance function can be used to generate regularization matrices yielding substantially better results, and permitting localization of features even if global basis functions are used. We are also able to offer a straightforward explanation for the spectral leakage problem identified by Trampert & Snieder. [ABSTRACT FROM AUTHOR]

Most linear inverse problems require regularization to ensure that robust and meaningful solutions can be found. Typically, Tikhonov-style regularization is used, whereby a preference is expressed for models that are somehow 'small' and/or 'smooth'. The strength of such preferences is expressed through one or more damping parameters, which control the character of the solution, and which must be set by the user. However, identifying appropriate values is often regarded as a matter of art, guided by various heuristics. As a result, such choices have often been the source of controversy and concern. By treating these as hyperparameters within a hierarchical Bayesian framework, we are able to obtain solutions that encompass the range of permissible regularization parameters. Furthermore, we show that these solutions are often well-approximated by those obtained via standard analysis using certain regularization choices which are—in a certain sense—optimal. We obtain algorithms for determining these optimal values in various cases of common interest, and show that they generate solutions with a number of attractive properties. A reference implementation of these algorithms, written in Python, accompanies this paper. [ABSTRACT FROM AUTHOR]

In geophysical inversion, inferences of Earth's properties from sparse data involve a trade-off between model complexity and the spatial resolving power. A recent Markov chain Monte Carlo (McMC) technique formalized by Green, the so-called trans-dimensional samplers, allows us to sample between these trade-offs and to parsimoniously arbitrate between the varying complexity of candidate models. Here we present a novel framework using transdimensional sampling over tree structures. This new class of McMC sampler can be applied to 1-D, 2-D and 3-D Cartesian and spherical geometries. In addition, the basis functions used by the algorithm are flexible and can include more advanced parametrizations such as wavelets, both in Cartesian and Spherical geometries, to permit Bayesian multiscale analysis. This new framework offers greater flexibility, performance and efficiency for geophysical imaging problems than previous sampling algorithms. Thereby increasing the range of applications and in particular allowing extension to trans-dimensional imaging in 3-D. Examples are presented of its application to 2-D seismic and 3-D teleseismic tomography including estimation of uncertainty. [ABSTRACT FROM AUTHOR]

Non-linear inverse problems in the geosciences often involve probabilistic sampling of multimodal density functions or global optimization and sometimes both. Efficient algorithmic tools for carrying out sampling or optimization in challenging cases are of major interest. Here results are presented of some numerical experiments with a technique, known as Parallel Tempering, which originated in the field of computational statistics but is finding increasing numbers of applications in fields ranging from Chemical Physics to Astronomy. To date, experience in use of Parallel Tempering within earth sciences problems is very limited. In this paper, we describe Parallel Tempering and compare it to related methods of Simulated Annealing and Simulated Tempering for optimization and sampling, respectively. A key feature of Parallel Tempering is that it satisfies the detailed balance condition required for convergence of Markov chain Monte Carlo (McMC) algorithms while improving the efficiency of probabilistic sampling. Numerical results are presented on use of Parallel Tempering for trans-dimensional inversion of synthetic seismic receiver functions and also the simultaneous fitting of multiple receiver functions using global optimization. These suggest that its use can significantly accelerate sampling algorithms and improve exploration of parameter space in optimization. Parallel Tempering is a meta-algorithm which may be used together with many existing McMC sampling and direct search optimization techniques. It's generality and demonstrated performance suggests that there is significant potential for applications to both sampling and optimization problems in the geosciences. [ABSTRACT FROM AUTHOR]

SUMMARY A meaningful interpretation of seismic measurements requires a rigorous quantification of the uncertainty. In an inverse problem, the data noise determines how accurately observations should be fit, and ultimately the level of detail contained in the recovered model. A common problem in seismic tomography is the difficulty in quantifying data uncertainties, and thus the required level of data fit. Traditionally, the complexity of the solution model (defined by both the number of basis functions and the regularization) is defined arbitrarily by the user prior to inversion with only limited use of data errors. In the context of multiscale problems, dealing with multiple data sets that are characterized by different noise variances and that span the Earth at different scales is a major challenge. Practitioners are usually required to arbitrarily weigh the contribution of each data type into the final solution. Furthermore, the basis functions are usually spatially uniform across the velocity field and regularization procedures are global, which prevents the solution model from accounting for the uneven spatial distribution of information. In this work we propose to address these issues with a Hierarchical Bayesian inversion. The new algorithm represents an extension of the transdimensional tomography to account for uncertainties in data noise. This approach has the advantage of treating the level of noise in each data set, as well as the number of model parameters, as unknowns in the inversion. It provides a parsimonious solution that fully represents the degree of knowledge one has about seismic structure (i.e. constraints, resolution and trade-offs). Rather than being forced to make decisions on parametrization, level of data fit and weights between data types in advance, as is often the case in an optimization framework, these choices are relaxed and instead constrained by the data themselves. The new methodology is presented in a synthetic example where both the data density and the underlying structure contain multiple length scales. Three ambient noise data sets that span the Australian continent at different scales are then simultaneously inverted to infer a multiscale tomographic image of Rayleigh wave group velocity for the Australian continent. The procedure turns out to be particularly useful when dealing with multiple data types with different unknown levels of noise as the algorithm is able to naturally adjust the fit to the different data sets and provide a velocity map with a spatial resolution adapted to the spatially variable information present in the data. [ABSTRACT FROM AUTHOR]

The reversible jump algorithm is a statistical method for Bayesian inference with a variable number of unknowns. Here, we apply this method to the seismic tomography problem. The approach lets us consider the issue of model parametrization (i.e. the way of discretizing the velocity field) as part of the inversion process. The model is parametrized using Voronoi cells with mobile geometry and number. The size, position and shape of the cells defining the velocity model are directly determined by the data. The inverse problem is tackled within a Bayesian framework and explicit regularization of model parameters is not required. The mobile position and number of cells means that global damping procedures, controlled by an optimal regularization parameter, are avoided. Many velocity models with variable numbers of cells are generated via a transdimensional Markov chain and information is extracted from the ensemble as a whole. As an aid to interpretation we visualize the expected earth model that is obtained via Monte Carlo integration in a straightforward manner. The procedure is particularly adept at imaging rapid changes or discontinuities in wave speed. While each velocity model in the final ensemble consists of many discontinuities at cell boundaries, these are smoothed out in the averaged ensemble solution while those required by the data are reinforced. The ensemble of models can also be used to produce uncertainty estimates and experiments with synthetic data suggest that they represent actual uncertainty surprisingly well. We use the fast marching method in order to iteratively update the ray geometry and account for the non-linearity of the problem. The method is tested here with synthetic data in a 2-D application and compared with a subspace method that is a more standard matrix-based inversion scheme. Preliminary results illustrate the advantages of the reversible jump algorithm. A real data example is also shown where a tomographic image of Rayleigh wave group velocity for the Australian continent is constructed together with uncertainty estimates. [ABSTRACT FROM AUTHOR]

Focuses on the importance of inversion techniques in the Earth Sciences. Provision of indirect observation on the surface; Existence of the geometrical concept of algorithm annealing; Efficacy of inversion techniques.

Describes a numerical method for solving partial differential equations in problems where traditional eulerian and lagrangian techniques fail. Use of the geometrical concept of natural neighbors; `Natural-element method' (NEM); Example of NEM. INSET: A new method to calculate natural-neighbour...

Monte Carlo methods are widespread in geophysics and have proved to be powerful in non-linear inverse problems. However, they are associated with significant practical challenges, including long calculation times, large output ensembles of Earth models, and difficulties in the appraisal of the results. This paper addresses some of these challenges using generative models, a family of tools that have recently attracted much attention in the machine learning literature. Generative models can, in principle, learn a probability distribution from a set of given samples and also provide a means for rapid generation of new samples which follow that approximated distribution. These two features make them well suited for application to the outputs of Monte Carlo algorithms. In particular, training a generative model on the posterior distribution of a Bayesian inference problem provides two main possibilities. First, the number of parameters in the generative model is much smaller than the number of values stored in the ensemble, leading to large compression rates. Secondly, once trained, the generative model can be used to draw any number of samples, thereby eliminating the dependence on an often large and unwieldy ensemble. These advantages pave new pathways for the use of Monte Carlo ensembles, including improved storage and communication of the results, enhanced calculation of numerical integrals, and the potential for convergence assessment of the Monte Carlo procedure. Here, these concepts are initially demonstrated using a simple synthetic example that scales into higher dimensions. They are then applied to a large ensemble of shear wave velocity models of the core–mantle boundary, recently produced in a Monte Carlo study. These examples demonstrate the effectiveness of using generative models to approximate posterior ensembles, and indicate directions to address various challenges in Monte Carlo inversion. [ABSTRACT FROM AUTHOR]

Following the linearized attenuation tomography from our previous study (Pejić et al., 2017, https://doi.org/10.1002/2016JB013692), we perform hierarchical transdimensional Bayesian tomography of the upper ≈400 km of the inner core, using 398 globally distributed t∗ estimates. The results are in good agreement with the ones obtained through linearized tomography: they show more complex attenuation pattern than the purely hemispherical one, and the noise estimated from the hierarchical inversion is in good agreement with estimates obtained from the Discrepancy Principle in the previous study. The attenuation pattern we observe gives more weight to the geodynamical models that couple the thermal anomalies of the lowermost mantle to the inner core boundary. Key Points: Transdimensional sampling of Voronoi cells on a spherical surface applied to a problem in global tomographyNo hemispherical pattern is observed in attenuation structure of the upper inner coreThe study provides a strong seismological hint for a strongly attenuating inner core [ABSTRACT FROM AUTHOR]

Reports on the attempts made by seismologists to study the behavior of seismic properties of the Earth's inner core. Usage of the linearized inversion techniques in geophysics to study the inner core; Description of the approach called ensemble inference; Implications for the nature of iron in the inner core.

We address the problem of unpolarized light spectroscopy of geological materials. Using infrared radiation, the aim of this technique is to learn about the absorbing species, such as hydroxyl. The use of unoriented samples leads to the need to perform a rigorous statistical analysis, so that the three principal absorbances of the crystal can be retrieved. We present here such an analysis based on a derivation of the probability density function for a single random measurement. Previous methods for retrieval of the absorbances are shown to be suboptimal, producing biased results that are sometimes even unphysical (e.g., negative estimates for an inherently positive quantity). The mathematical structure of the problem is developed to use the maximum likelihood estimation method, and we show how to optimize for the three absorbance parameters. This leads to good parameter retrieval on both synthetic and real data sets. [ABSTRACT FROM AUTHOR]

Abstract: Cross correlation of seismograms provides new information on the Earth both through the exploitation of ambient noise and specific components of earthquake records. Here we cross‐correlate recordings of large earthquakes on a planetary scale and identify a range of hitherto unobserved seismic phases in Earth's correlation wavefield. We show that both arrivals with the timing expected for the regular seismic wavefield and previously unexplained phases are produced by interference between seismic paths having the same ray parameter but with only a subset of propagation legs in common. This insight explains the origin and generation mechanism of the features of Earth's correlation wavefield and opens up new ways of addressing issues in global seismology. Strong similarity between observed and synthesized correlation wavefields indicates that the Earth's radial structure is remarkably well constrained in the intermediate period range. [ABSTRACT FROM AUTHOR]

This paper presents the application of a novel trans-dimensional sampling approach to a time domain airborne electromagnetic (AEM) inverse problem to solve for plausible conductivities of the subsurface. Geophysical inverse field problems, such as time domain AEM, are well known to have a large degree of non-uniqueness. Common least-squares optimisation approaches fail to take this into account and provide a single solution with linearised estimates of uncertainty that can result in overly optimistic appraisal of the conductivity of the subsurface. In this new non-linear approach, the spatial complexity of a 2D profile is controlled directly by the data. By examining an ensemble of proposed conductivity profiles it accommodates non-uniqueness and provides more robust estimates of uncertainties. [ABSTRACT FROM AUTHOR]

The mean land-surface temperature represents an important boundary condition for many geothermal studies. This boundary is particularly important to help constrain the models made to analyse resource systems, many of which are shallow in nature and observe relatively small thermal gradients. Consequently, a mean land-surface temperature map of the Australian continent has been produced from 13 years of MODIS satellite imagery, for the period 2003–2015. The map shows good agreement with independent methods of estimating mean land-surface temperature, including borehole surface-temperature extrapolation and long-term, near-surface ground measurements. In comparison to previously used methods of estimating mean land-surface temperature, our new estimates are up to 12 °C warmer. The MODIS-based method presented in this study provides spatially continuous estimates of land-surface temperature that can be incorporated as the surface thermal boundary condition in geothermal studies. The method is also able to provide a quantification of the uncertainties expected in the application of these estimates for the purposes of thermal modelling. [ABSTRACT FROM AUTHOR]

An abstract of the article "Quantitative absorbance spectroscopy with unpolarized light: Part I. Physical and mathematical development," by Malcolm Sambridge, John Fitz Gerald, and Hugh St. C. O'Neill is presented.

This paper outlines a methodology for the estimation of the environmental noise equivalent reflectance in aquatic remote sensing imagery using an object-based segmentation approach. Noise characteristics of remote sensing imagery directly influence the accuracy of estimated environmental variables and provide a framework for a range of sensitivity, sensor specification, and algorithm design studies. The proposed method enables estimation of the noise equivalent reflectance covariance of remote sensing imagery through homogeneity characterization using image segmentation. The method is first tested on a synthetic data set with known noise characteristics and is successful in estimating the noise equivalent reflectance under a range of segmentation structures. Testing on a Portable Hyperspectral Imager for Low-Light Spectroscopy (PHILLS) hyperspectral image in a coral reef environment shows the method to produce comparable noise equivalent reflectance estimates in an optically shallow water environment to those previously derived in optically deep water. This method is of benefit in aquatic studies where homogenous regions of optically deep water were previously required for image noise estimation. The ability of the method to characterize the covariance of an image is of significant benefit when developing probabilistic inversion techniques for remote sensing. [ABSTRACT FROM PUBLISHER]

Abstract: India–Somalia and North America–Eurasia relative motions since Early Miocene ( ) have been recently reconstructed at unprecedented temporal resolution ( ) from magnetic surveys of the Carlsberg and northern Mid-Atlantic Ridges. These new datasets revamped interest in the convergence of India relative to Eurasia, which is obtained from the India–Somalia–Nubia–North America–Eurasia plate circuit. Unless finite rotations are arbitrarily smoothed through time, however, the reconstructed kinematics (i.e. stage Euler vectors) appear to be surprisingly unusual over the past . In fact, the Euler pole for the India–Eurasia rigid motion scattered erratically over a broad region, while the associated angular velocity underwent sudden increases and decreases. Consequently, convergence across the Himalayan front featured significant speed-ups as well as slow-downs with almost no consistent trend. Arguably, this pattern arises from the presence of data-noise, which biases kinematic reconstructions—particularly at high temporal resolution. The rapid and important India–Eurasia plate-motion changes reconstructed since Early Miocene are likely to be of apparent nature, because they cannot result even from the most optimistic estimates of torques associated, for instance, with the descent of the Indian slab into Earth's mantle. Our previous work aimed at reducing noise in finite-rotation datasets via an expanded Bayesian formulation, which offers several advantages over arbitrary smoothing methods. Here we build on this advance and revise the India–Eurasia kinematics since , accounting also for three alternative histories of rifting in Africa. We find that India–Eurasia kinematics are simpler and, most importantly, geodynamically plausible upon noise reduction. Convergence across the Himalayan front overall decreased until , but then systematically increased, albeit moderately, towards the present-day. We test with global dynamic models of the coupled mantle/lithosphere system how a previously proposed uplift of the Tibetan plateau by 1–2km, as well as the break-up of the Indo-Australian plate into the India and Australia units, impact the kinematic record from to . Our modelling results indicate that the convergence slow-down until resulted chiefly from the Indo-Australian break-up. However, the tectonic scenario yielding the best fit to kinematic reconstructions is indeed one where the Tibetan mean elevation increased by 1–2km before the Indo-Australian break-up was accomplished. [Copyright &y& Elsevier]

This paper applies parallel tempering within a Bayesian formulation for strongly nonlinear geoacoustic inverse problems. Bayesian geoacoustic inversion consists of sampling the posterior probability density (PPD) of seabed parameters to estimate integral properties, such as marginal probability distributions, based on ocean acoustic data and prior information. This sampling is usually carried out using the Markov-chain Monte Carlo method of Metropolis-Hastings sampling. However, standard sampling methods can be very inefficient for strongly nonlinear problems involving multi-modal PPDs with the potential to miss important regions of the parameter space and to significantly underestimate parameter uncertainties. Parallel tempering achieves efficient/effective sampling of challenging parameter spaces with the ability to transition freely between multiple PPD modes by running parallel Markov chains at a series of increasing sampling temperatures with probabilistic interchanges between chains. The approach is illustrated for inversion of (simulated) acoustic reverberation data for which the PPD is highly multi-modal. While Metropolis-Hastings sampling gives poor results even with very large sample sizes, parallel tempering provides efficient, convergent sampling of the PPD. Methods to enhance the efficiency of parallel tempering are also considered. [ABSTRACT FROM AUTHOR]

Abstract: We present a method to quantify abrupt changes (or changepoints) in data series, represented as a function of depth or time. These changes are often the result of climatic or environmental variations and can be manifested in multiple datasets as different responses, but all datasets can have the same changepoint locations/timings. The method we present uses transdimensional Markov chain Monte Carlo to infer probability distributions on the number and locations (in depth or time) of changepoints, the mean values between changepoints and, if required, the noise variance associated with each dataset being considered. This latter point is important as we generally will have limited information on the noise, such as estimates only of measurement uncertainty, and in most cases it is not practical to make repeat sampling/measurement to assess other contributions to the variation in the data. We describe the main features of the approach (and describe the mathematical formulation in supplementary material), and demonstrate its validity using synthetic datasets, with known changepoint structure (number and locations of changepoints) and distribution of noise variance for each dataset. We show that when using multiple data, we expect to achieve better resolution of the changepoint structure than when we use each dataset individually. This is conditional on the validity of the assumption of common changepoints between different datasets. We then apply the method to two sets of real geochemical data, both from peat cores, taken from NE Australia and eastern Tibet. Under the assumption that changes occur at the same time for all datasets, we recover solutions consistent with those previously inferred qualitatively from independent data and interpretations. However, our approach provides a quantitative estimate of the relative probability of the inferred changepoints, allowing an objective assessment of the significance of each change. [Copyright &y& Elsevier]

We present a globally distributed data set of ∼400 000 frequency-dependent SH-wave traveltimes. An automated technique is used to measure teleseismic S, ScS and SS traveltimes at several periods ranging from 10 to 51 s. The targeted seismic phases are first extracted from the observed and synthetic seismograms using an automated time window algorithm. Traveltimes are then measured at several periods, by cross-correlation between the selected observed and synthetic filtered waveforms. Frequency-dependent effects due to crustal reverberations beneath each receiver are handled by incorporating crustal phases into WKBJ synthetic waveforms. After correction for physical dispersion due to intrinsic anelastic processes, we observe a residual traveltime dispersion on the order of 1–2 s in the period range of analysis. This dispersion occurs differently for S, ScS and SS, which is presumably related to their differing paths through the Earth. We find that: (1) Wavefront-healing phenomenon is observed for S and to a lesser extent SS waves having passed through very low velocity anomalies. (2) A preferred sampling of high velocity scatterers located at the CMB may explain our observation that ScS waves travel faster at low-frequency than at high-frequency. (3) A frequency-dependent attenuation , with , is compatible with the globally averaged dispersion observed for S waves. [ABSTRACT FROM AUTHOR]

The capillary pressure in a reservoir determines the saturation distribution, and hence the total in situ volumes of fluids (oil/water/gas). The accurate knowledge of the capillary pressure distribution is one of the primary factors that may be decisive in the reliable estimation of hydrocarbon reserves. The centrifuge procedure provides laboratory data, which can be inverted to derive capillary pressure curves in laboratory conditions. The derived laboratory capillary pressure curves are then scaled up for full-field simulation of petroleum reservoirs. The inversion procedure is uncertain and gives nonunique capillary pressure curves. The standard industrial practice, however, has been to derive a single capillary pressure curve, and ignoring the uncertainty. Since the capillary pressure is uncertain, estimates of important reservoir parameters dependent on capillary pressure distribution are uncertain. This paper shows how the uncertainty in centrifuge capillary pressure can be quantified. It also shows how this uncertainty propagates through the scale-up process and impacts on the estimate of the oil recovery potential. We illustrate using a simple, two-phase (oil/water) synthetic reservoir model. [ABSTRACT FROM AUTHOR]

We compare explicit differential operators for unstructured grids and their accuracy with the aim of solving time-dependent partial differential equations in geophysical applications. As many problems suggest the use of staggered grids we investigate different schemes for the calculation of space derivatives on two separate grids. The differential operators are explicit and local in the sense that they use only information of the function in their nearest neighborhood, so that no matrix inversion is necessary. This makes this approach well-suited for parallelization. Differential weights are obtained either with the finite-volume method or using natural neighbor coordinates. Unstructured grids have advantages concerning the simulation of complex geometries and boundaries. Our results show that while in general triangular (hexagonal) grids perform worse than standard finitedifference approaches, the effects of grid irregularities on the accuracy of the space derivatives are comparably small for realistic grids. This suggests that such a finite-difference-like approach to unstructured grids may be an alternative to other irregular grid methods such as the finite-element technique. [ABSTRACT FROM AUTHOR]

The core‐mantle boundary (CMB) is the most extreme boundary within the Earth where the liquid, iron‐rich outer core interacts with the rocky, silicate mantle. The nature of the lowermost mantle atop the CMB, and its role in mantle dynamics, is not completely understood. Various regional studies have documented significant heterogeneities at different spatial scales. While there is a consensus on the long scale length structure of the inferred S‐wave speed tomograms, there are also notable differences stemming from different imaging methods and datasets. Here we aim to overcome over‐smoothing and avoid over‐fitting data for the case where the spatial coverage is sparse and the inverse problem ill‐posed. We present an S‐wave tomography model at a global scale for the Lowermost Mantle (LM) using the Hierarchical Trans‐Dimensional Bayesian Inversion (HTDBI) framework, LM‐HTDBI. Our LM‐HTDBI analysis of ScS‐S travel times includes uncertainty, and the complexity of the model is deduced from the data itself through an implicit parameterization of the model space. Our comprehensive resolution estimates indicate that short‐scale anomalies are significant and resolvable features of the lowermost mantle regardless of the chosen mantle‐model reference to correct the travel times above the D" layer. The recovered morphology of the Large‐Low‐Shear‐wave Velocity Provinces (LLSVPs) is complex, featuring small high‐velocity patches among low‐velocity domains. Instead of two large, unified, and smooth LLSVPs, the newly obtained images suggest that their margins are not uniformly flat. Plain Language Summary: The lowermost mantle sits atop the core‐mantle boundary, the most dramatic boundary within our planet, with contrasts in physical properties that exceed those that exist at the surface. Despite significant progress in recent years, this part of the Earth is not well understood, and various tomographic studies on a global scale, along with regional studies that focus on seismic waveform modeling, pave the path towards higher resolution and new understanding. Important questions to answer are on the distribution, shape, size and composition of inhomogeneities in the lowermost mantle, and their critical role in the mantle and core dynamics. While there is a general consensus on the long‐scale length structures inferred from long‐period shear waves, there are notable differences in details of the tomograms of the lowermost mantle, stemming from the use of different imaging methods and datasets. Here, we utilize a large travel time data set of ScS and S waves with a significant addition of new measurements sensitive to the lowermost mantle to perform a probabilistic shear‐wave tomography, and we retrieve a high‐resolution image of the lowermost mantle. The new shear‐wave speed tomogram and comprehensive resolution‐estimations indicate that short and medium scale inhomogeneities are omnipresent features of the lowermost mantle. Key Points: S‐wave tomography of the lowermost mantle using state‐of‐the‐art Bayesian approach with 2D spherical Voronoi cellsThe inversion technique treats the model complexity and the data noise as free parameters and avoids damping and smoothingThis study provides an important bridge between long‐scale features at a global scale and short‐scale features of regional models [ABSTRACT FROM AUTHOR]

Obtaining slip distributions for earthquakes results in an ill-posed inverse problem. While this implies that only limited and uncertain information can be recovered from the data, inferences are typically made based only on a single regularized model. Here, we develop an inversion approach that can quantify uncertainties in a Bayesian probabilistic framework for the finite fault inversion (FFI) problem. The approach is suitably efficient for rapid source characterization and includes positivity constraints for model parameters, a common practice in FFI, via coordinate transformation to logarithmic space. The resulting inverse problem is nonlinear and the most probable solution can be obtained by iterative linearization. In addition, model uncertainties are quantified by approximating the posterior probability distribution by a Gaussian distribution in logarithmic space. This procedure is straightforward since an analytic expression for the Hessian of the objective function is obtained. In addition to positivity, we apply smoothness regularization to the model in logarithmic space. Simulations based on surface wave data show that smoothing in logarithmic space penalizes abrupt slip changes less than smoothing in linear space. Even so, the main slip features of models that are smooth in linear space are recovered well with logarithmic smoothing. Our synthetic experiments also show that, for the data set we consider, uncertainty is low at the shallow portion of the fault and increases with depth. In addition, a simulation with a large station azimuthal gap of 180° significantly increases the slip uncertainties. Further, the marginal posterior probabilities obtained from our approximate method are compared with numerical Markov Chain Monte Carlo sampling. We conclude that the Gaussian approximation is reasonable and meaningful inferences can be obtained from it. Finally, we apply the new approach to observed surface wave records from the great Illapel earthquake (Chile, 2015, M w = 8.3). The location and amplitude of our inferred peak slip is consistent with other published solutions but the spatial slip distribution is more compact, likely because of the logarithmic regularization. We also find a minor slip patch downdip, mainly in an oblique direction, which is poorly resolved compared to the main slip patch and may be an artefact. We conclude that quantifying uncertainties of finite slip models is crucial for their meaningful interpretation, and therefore rapid uncertainty quantification can be critical if such models are to be used for emergency response. [ABSTRACT FROM AUTHOR]

Antarctica and Greenland hold enough ice to raise sea level by more than 65m if both ice sheets were to melt completely. Predicting future ice sheet mass balance depends on our ability to model these ice sheets, which is limited by our current understanding of several key physical processes, such as iceberg calving. Large-scale ice flow models either ignore this process or represent it crudely. To model fractured zones, an important component of many calving models, continuum damage mechanics as well as linear fracture mechanics are commonly used. However, these methods have a large number of uncertainties when applied across the entire Antarctic continent because the models were typically tuned to match processes seen on particular ice shelves. Here we present an alternative, statistics-based method to model the most probable zones of the location of fractures and demonstrate our approach on all main ice shelf regions in Antarctica, including the Antarctic Peninsula. We can predict the location of observed fractures with an average success rate of 84% for grounded ice and 61% for floating ice and a mean overestimation error rate of 26% and 20 %, respectively. We found that Antarctic ice shelves can be classified into groups based on the factors that control fracture location. [ABSTRACT FROM AUTHOR]

Because of increased emissions of greenhouse gases oceans are warming, causing sea level to rise as the density of seawater falls. Predicting the rates of steric expansion is challenging because of the natural variability of the ocean and because observations are insufficient to adequately cover the ocean basins. Here, we investigate the ability of one ocean reanalysis, two objective analyses, and one combination of satellite geodetic measurements to accommodate data gaps and to reconstruct typical patterns of the steric sea level variability at interannual and multidecadal time scales. Six climate indices are used to identify robust features of the internal variability, using a Least Absolute Shrinkage and Selection Operator (LASSO) regression to select significant predictors of the steric variability. Spatially consistent fingerprints are revealed for all climate indices in the ocean reanalysis dataset, allowing the recovery of most of the steric variability observed in the tropical and North Pacific, as well as large fractions of the Atlantic and Indian Ocean signals. Robust climate mode fingerprints are also identified with high spatial resolution but limited temporal coverage in the geodetic observations. The objective analyses fail to detect many of the patterns expected from climate modes, especially before the Argo era. Climate indices constitute valuable yet underexploited tools to assess the performance of different techniques to reconstruct steric sea levels at interannual and multidecadal scales. Such progress will increase confidence in the historical reconstructions of steric sea levels, which is necessary to improve the closure of regional and global sea level budgets and to validate the predictions of climate models. [ABSTRACT FROM AUTHOR]

Over the past decade, we have applied laser ablation ICP mass spectrometry U-series analysis to a large number of bones and teeth. The method is fast and provides high-resolution data of U-series isotopes, which give insights into the complexity of uranium migration into, within and out of bones. In this paper, we present our laser ablation approach in detail, from the experimental set up, to data reduction and uranium diffusion modelling. Laser ablation analysis can now be applied with minimum damage. Complete diffusion data sets can be obtained by laser drilling which leaves a hole of around 200 μm in diameter and 1–2 mm depth. If the natural U-diffusion follows a simple single-stage process, valuable age information can be obtained. In other cases, highly complex U-migrations have been observed, which make any age assessment impossible. Two examples from Wadjak and Callao illustrate the potential of nearly non-destructive laser ablation U-series analysis of human skeletal tissues. The results demonstrate the physical presence of modern humans in south-east Asia at a time when our species just started to migrate into Europe. [ABSTRACT FROM AUTHOR]

The major cause of sea-level change during ice ages is the exchange of water between ice and ocean and the planet's dynamic response to the changing surface load. Inversion of ~1,000 observations for the past 35,000 y from localities far from former ice margins has provided new constraints on the fluctuation of ice volume in this interval. Key results are: (i) a rapid final fall in global sea level of ~40 m in <2,000 y at the onset of the glacial maximum ~30,000 y before present (30 ka BP); (ii) a slow fall to -134 m from 29 to 21 ka BP with a maximum grounded ice volume of ~52 × 106 km3 greater than today; (iii) after an initial short duration rapid rise and a short interval of near-constant sea level, the main phase of deglaciation occurred from ~16.5 ka BP to ~8.2 ka BP at an average rate of rise of 12 m⋅ka-1 punctuated by periods of greater, particularly at 14.5-14.0 ka BP at ≥40 mm•y-1 (MWP-1A), and lesser, from 12.5 to 11.5 ka BP (Younger Dryas), rates; (iv) no evidence for a global MWP-1B event at ~11.3 ka BP; and (v) a progressive decrease in the rate of rise from 8.2 ka to ~2.5 ka BP, after which ocean volumes remained nearly constant until the renewed sea-level rise at 100-150 y ago, with no evidence of oscillations exceeding ~15-20 cm in time intervals ≥200 y from 6 to 0.15 ka BP. [ABSTRACT FROM AUTHOR]

This paper develops a probabilistic Bayesian approach to the problem of inferring the spatiotemporal evolution of earthquake rupture on a fault surface from seismic data with rigorous uncertainty estimation. To date, uncertainties of rupture parameters are poorly understood, and the effect of choices such as fault discretization on uncertainties has not been studied. We show that model choice is fundamentally linked to uncertainty estimation and can have profound effects on results. The approach developed here is based on a trans-dimensional self-parametrization of the fault, avoids regularization constraints and provides rigorous uncertainty estimation that accounts for model-selection ambiguity associated with the fault discretization. In particular, the fault is parametrized using self-adapting irregular grids which intrinsically match the local resolving power of the data and provide parsimonious solutions requiring few parameters to capture complex rupture characteristics. Rupture causality is ensured by parametrizing rupture-onset time by a rupture-velocity field and obtaining first rupture times from the eikonal equation. The Bayesian sampling of the parameter space is implemented on a computer cluster with a highly efficient parallel tempering algorithm.The inversion is applied to simulated and observed W-phase waveforms from the 2010 Maule (Chile) earthquake. Simulation results show that our approach avoids both over- and underparametrization to ensure unbiased inversion results with uncertainty estimates that are consistent with data information. The simulation results also show the ability of W-phase data to resolve the spatial variability of slip magnitude and rake angles. In addition, sensitivity to spatially dependent rupture velocities exists for kinematic slip models.Application to the observed data indicates that residual errors are highly correlated and likely dominated by theory error, necessitating the iterative estimation of a non-stationary data covariance matrix. The moment magnitude for the Maule earthquake is estimated to be ∼8.9, with slip concentrated in two zones updip of and north and south of the hypocentre, respectively. While this aspect of the slip distribution is similar to previous studies, we show that the slip maximum in the southern zone is poorly resolved compared to the northern zone. Both slip maxima are higher than reported in previous studies, which we speculate may be due to the lack of bias caused by the regularization used in other studies. [ABSTRACT FROM AUTHOR]

We conduct an experiment to investigate whether linearity in the observed velocity gradient in the volume of the inner core sampled by the PKP ray paths beneath Central America is a robust approximation. Instead of solving an optimization problem, we approach it within the Bayesian inference. This is an ensemble approach, where model specification is relaxed so that instead of only one solution, groups of reasonable models are acceptable. Furthermore, in transdimensional Bayesian inference used here, the number of basis functions needed to model observations is by itself an unknown. Our modeling reveals that in the ensemble of models, the most likely are those containing only 2 nodes (linear trend). Thus our result justifies the assumption used for the determination of inner core rotation with respect to the rest of the mantle that the observed gradient is constant in its nature (linear). Recent observations in seismology suggest that it is likely that the spatial variability in elastic parameters is a widespread phenomenon in the inner core. Future array observations will further constrain spatial extent and magnitude of velocity changes and show whether there is a significant difference between these observations in the two quasi-hemispheres of the inner core. [ABSTRACT FROM AUTHOR]

Geodynamical models and seismic observations suggest that the Earth's solid inner core rotates at a different rate than the mantle. However, discrepancies exist in rotation rate estimates based on seismic waves produced by earthquakes. Here we investigate the inherent assumption of a constant rotation rate using earthquake doublets-repeating earthquakes that produce similar waveforms. We detect that the rotation rate of the Earth's inner core with respect to the mantle varies with time. We perform an inverse analysis of 7 doublets observed at the College station, Alaska, as well as 17 previously reported doublets, and reconstruct a history of differential inner-core rotation between 1961 and 2007. We find that the observed doublets are consistent with a model of an inner core with an average differential rotation rate of 0.25-0.48° yr−1 and decadal fluctuations of the order of 1° yr−1 around the mean. The decadal fluctuations explain discrepancies between previous core rotation models and are in concordance with recent geodynamical simulations. [ABSTRACT FROM AUTHOR]

There is significant seismic activity in the region around Australia, largely due to the plate boundaries to the north and to the east of the mainland. This activity results in serious seismic and tsunami hazard in the coastal areas of Australia. Hence seismicity is and will be monitored in real time by Geoscience Australia (GA), which uses a network of permanent broadband seismometers. Seismic moment tensor (MT) solutions are currently determined using 1-D, radially symmetric models of Earth and this requires augmentation by recording stations located outside of Australia. A 3-D model of the Australian continent developed recently using full waveform tomography now offers the opportunity to significantly improve the determination of MT solutions of earthquakes from tectonically active regions. A complete-waveform, time-domain MT inversion method has been developed using a point-source approximation. A series of synthetic tests using first a 1-D and then a 3-D structural model has been performed. The feasibility of deploying 3-D versus 1-D Earth structure for the inversion of seismic data has been studied and the advantages of using the 3-D structural model were illustrated with examples. The 3-D model is superior to the 1-D model, as a number of sensitivity tests show. The ultimate goal of this work is an automated MT inversion system in Australia relying on GA and other international stations, although more work remains to be done before the full implementation of such a scheme in real time. [ABSTRACT FROM AUTHOR]

Abstract: The silicon (Si) isotope (δ 30Si) composition of deep-sea sponges from near Antarctica, subantarctic waters (Tasmania Seamounts) and subtropical waters north of New Zealand vary widely between +0.87‰ and −3.40‰ (vs. NBS28). Depth profiles show that sponge δ 30Si compositions trend to lower values with increasing depth. This is exemplified by sponges from the Tasmania Seamounts where δ 30Si varies from +0.87‰ to −3.13‰ over a depth range from 100m to 1200m. These changes in δ 30Si of sponges are inconsistent with a Rayleigh type isotope fractionation model requiring constant δ 30Si fractionation between sponge and seawater. We conclude that overall Si isotope fractionation Δ 30Si (δ 30Si sponge – δ 30Si seawater) is influenced by seawater Si concentration, with more fractionated (lower) isotope values being associated with sponges collected from waters high in Si. We invoke and fit a model whereby the Δ 30Si fractionation varies as a function Si influx and efflux. Using this model it appears that Δ 30Si fractionation during transport into the sponge is constant at −1.34‰. The model also shows asymptotic behaviour with Δ 30Si trending towards a maximum of −6.02‰ at very high Si concentrations. These results suggest that the δ 30Si composition of fossil spicules may be useful for reconstruction paleo-Si concentrations during the past. [Copyright &y& Elsevier]

Abstract: We present an overview of Markov chain Monte Carlo, a sampling method for model inference and uncertainty quantification. We focus on the Bayesian approach to MCMC, which allows us to estimate the posterior distribution of model parameters, without needing to know the normalising constant in Bayes'' theorem. Given an estimate of the posterior, we can then determine representative models (such as the expected model, and the maximum posterior probability model), the probability distributions for individual parameters, and the uncertainty about the predictions from these models. We also consider variable dimensional problems in which the number of model parameters is unknown and needs to be inferred. Such problems can be addressed with reversible jump (RJ) MCMC. This leads us to model choice, where we may want to discriminate between models or theories of differing complexity. For problems where the models are hierarchical (e.g. similar structure but with a different number of parameters), the Bayesian approach naturally selects the simpler models. More complex problems require an estimate of the normalising constant in Bayes'' theorem (also known as the evidence) and this is difficult to do reliably for high dimensional problems. We illustrate the applications of RJMCMC with 3 examples from our earlier working involving modelling distributions of geochronological age data, inference of sea-level and sediment supply histories from 2D stratigraphic cross-sections, and identification of spatially discontinuous thermal histories from a suite of apatite fission track samples distributed in 3D. [Copyright &y& Elsevier]

• A quantitative method resolves the fine structure of the D″ layer. • Small-scale heterogeneities exist at CMB beneath eastern Alaska & northern Pacific. • There is a lack of correlation between D″ topography and shear velocity variation. • The D″ layer may be a mixed boundary layer comprised of thermo-chemical anomalies. • Inferred accumulated slab material suggests different mantle convection styles. The D″ layer, regarded as a thermal boundary layer and a chemically distinct region above the core-mantle boundary (CMB), has been associated with the phase transition from bridgmanite (Bm) to post-perovskite (pPv) in the lowermost mantle. However, the composition of the lowermost mantle and thermal conditions where Bm-pPv phase-transition occurs is still debatable. The methods typically used to study the fine-scale seismic features in the D″ layer has provided important information. However, trial-and-error seismic waveform modelling cannot uniquely quantify D″-layer properties because of subjective model-parameterization choices and inherent non-uniqueness of solutions and the waveform inversion method has a limited resolution of the velocity gradient and depth of the D″ discontinuity. We develop a grid-search scheme to constrain the detailed 1-D shear-wave velocity structure in the lowermost mantle beneath Alaska and the northern Pacific, accompanied with quantitative assessment of the uncertainty of 1D models. Our results show strong lateral variations of the D″ discontinuity from west to east beneath Alaska, along with the existence of smaller-scale heterogeneities in the east. We find a broad velocity increase, as thick as 240 km, at the top of D″ that indicates this region may involve a composite of downwelling thermo-chemical anomalies at the CMB. There are even smaller scale heterogeneities of approximately 120 km × 120 km in size with larger lateral variations in the lowermost mantle beneath northern Pacific. Both the magnitude and gradient of the velocity at the top of the D″ layer vastly change in adjacent regions, with an increase from 2.8% to 4.5% in magnitude and from 0.08% to 1.2% in gradient, but with a relatively consistent depth of the D″ discontinuity at ∼340 km above the CMB. The weak correlation between D″ topography and velocity variations indicate chemical heterogeneities must be present beneath the northern Pacific, which might come from north-westward subducted Pacific oceanic lithosphere. Our characterisation of the spatial pattern of small-scale heterogeneities in the lowermost mantle supports a hybrid thermo-chemical boundary layer (TCBL) model beneath Alaska and northern Pacific. [ABSTRACT FROM AUTHOR]

We present the results of a new global Bayesian tomographic inversionfor a radial anisotropic shear wave velocity model of the Earth. Thedevelopment of radial anisotropic shear wave models are useful forunderstanding mantle dynamics. Our model improves on existing work byproviding robust estimates of uncertainties.A key aspect of this new model is the careful selection of surfacewave observations produced by earthquakes from the Harvard CMTcatalog. We automatically process a large set of seismic signalsrecorded for these events by simultaneously inverting Love andRayleigh waveforms up to the fifth overtone. After initial processing,only events where high quality Love and Rayleigh surface waves areobserved are retained for the tomographic inversion. This careful dataselection ensures that inference of radial anisotropy is not biaseddue to variations in spatial coverage of Love and Rayleigh observations.Our tomographic inversion uses a novel Bayesian Trans-dimensional Treeapproach that dynamically adapts the complexity of the solution toincrease spatial resolution where observations provide additionalconstraint. This leads to a robust estimation of both the shear wavemodel and its uncertainties. Furthermore, we demonstrate how ourapproach can be extended to a joint inversion to include body shearwaves to further constrain radial anisotropy of the Earths mantle. [ABSTRACT FROM AUTHOR]

Working with a large temporal dataset spanning several decades often represents a challenging task, especially when the record is heterogeneous and incomplete. The use of statistical laws could potentially overcome these problems. Here we apply Benford's Law (also called the 'First-Digit Law') to the traveled distances of tropical cyclones since 1842. The record of tropical cyclones has been extensively impacted by improvements in detection capabilities over the past decades. We have found that, while the first-digit distribution for the entire record follows Benford's Law prediction, specific changes such as satellite detection have had serious impacts on the dataset. The least-square misfit measure is used as a proxy to observe temporal variations, allowing us to assess data quality and homogeneity over the entire record, and at the same time over specific periods. Such information is crucial when running climatic models and Benford's Law could potentially be used to overcome and correct for data heterogeneity and/or to select the most appropriate part of the record for detailed studies. [ABSTRACT FROM AUTHOR]

An abstract of the article "Quantitative absorbance spectroscopy with unpolarized light: Part II. Experimental evaluation and development of a protocol for quantitative analysis of mineral IR spectra," Hugh St. C. O'Neill, John Fitz Gerald, and Malcolm Sambridge is presented.