Spherical Slepian functions and the polar gap in geodesy

The estimation of potential fields such as the gravitational or magnetic potential at the surface of a spherical planet from noisy observations taken at an altitude over an incomplete portion of the globe is a classic example of an ill-posed inverse problem. We show that this potential-field estimation problem has deep-seated connections to Slepian's spatiospectral localization problem which seeks bandlimited spherical functions whose energy is optimally concentrated in some closed portion of the unit sphere. This allows us to formulate an alternative solution to the traditional damped least-squares spherical harmonic approach in geodesy, whereby the source field is now expanded in a truncated Slepian function basis set. We discuss the relative performance of both methods with regard to standard statistical measures such as bias, variance and mean squared error, and pay special attention to the algorithmic efficiency of computing the Slepian functions on the region complementary to the axisymmetric polar gap characteristic of satellite surveys. The ease, speed, and accuracy of our method make the use of spherical Slepian functions in earth and planetary geodesy practical.     

A six-parameter statistical model of the Earth's magnetic field

A six-parameter statistical model of the non-dipole geomagnetic field is fitted to 2597 harmonic coefficients determined by Cain, Holter & Sandee (1990) from MAGSAT data. The model includes sources in the core, sources in the crust, and instrument errors. External fields are included with instrument errors. The core and instrument statistics are invariant under rotation about the centre of the Earth, and one of the six parameters describes the deviation of the crustal statistics from rotational invariance. The model treats the harmonic coefficients as independent random samples drawn from a Gaussian distribution. The statistical model of the core field has a correlation length of about 500 km at the core-mantle boundary, too long to be attributed to a white noise source just below the boundary layers at the top of the core. The estimate of instrument errors obtained from the statistical model is in good agreement with an independent estimate based on tests of the instruments (Langel, Ousley & Berbert 1982).     

Modelling of mantle postglacial relaxation in axisymmetric geometry with a composite rheology and a glacial load interpolated by adjusted spherical harmonics analysis

Although studies on glacial isostatic adjustment usually assume a purely linear rheology, we have previously shown that mantle relaxation after the melting of Laurentide ice sheet is better described by a composite rheology including a non-linear term. This modelling is, however, based on axially symmetric geometry and glacial forcing derived from ICE-3G and suffers from a certain amount of arbitrariness in the definition of the ice load. In this work we apply adjusted spherical harmonics analysis to interpolate the ice thicknesses of ICE-3G and ICE-1 glaciological models. This filters out the non-axisymmetric components of the ice load by considering only the zonal terms in the spherical harmonics expansion. The resulting load function is used in finite-element simulation of postglacial rebound to compare composite versus purely linear rheology. Our results confirm that composite rheology can explain relative sea level (RSL) data in North America significantly better than a purely linear rheology. The performance of composite rheology suggests that in future investigations, it may be better to use this more physically realistic creep law for modelling mantle deformation induced by glacial forcing.     

Taking into account truncation problems and geomagnetic model accuracy in assessing computed flows at the core-mantle boundary

SUMMARY
Since the time Roberts & Scott (1965) first expressed the key 'frozen flux' hypothesis relating the secular variation of the geomagnetic field (SV) to the flow at the core surface, a large number of studies have been devoted to building maps of the flow and inferring its fundamental properties from magnetic observations at the Earth's surface. There are some well-known difficulties in carrying out these studies, such as the one linked to the non-uniqueness of the flow solution [if no additional constraint is imposed on the flow (Backus 1968)] which has been thoroughly investigated. In contrast little investigation has been made up to now to estimate the exact importance of other difficulties, although the different authors are usually well aware of their existence. In this paper we intend to make as systematic as possible a study of the limitations linked to the use of truncated spherical harmonic expansions in the computation of the flow. Our approach does not rely on other assumptions than the frozen flux, the insulating mantle and the large-scale flow assumptions along with some simple statistical assumptions concerning the flow and the Main Field. Our conclusions therefore apply to any (toroidal, steady or tangentially geostrophic) of the flow models that have already been produced; they can be summarized in the following way: first, because of the unavoidable truncation of the spherical harmonic expansion of the Main Field to degree 13, no information will ever be derived for the components of the flow with degree larger than 12; second, one may truncate the spherical harmonic expansion of the flow to degree 12 with only a small impact on the first degrees of the flow. Third, with the data available at the present day, the components of the flow with degree less than 5 are fairly well known whereas those with degree greater than 8 are absolutely unconstrained.     

The statistical distribution of magnetotelluric apparent resistivity and phase

The marginal distributions for the magnetotelluric (MT) magnitude squared response function (and hence apparent resistivity) and phase are derived from the bivariate complex normal distribution that describes the distribution of response function estimates when the Gauss–Markov theorem is satisfied and the regression random errors are normally distributed. The distribution of the magnitude squared response function is shown to be non-central chi-squared with 2 degrees of freedom, with the non-centrality parameter given by the squared magnitude of the true MT response. The standard estimate for the magnitude squared response function is biased, with the bias proportional to the variance and hence important when the uncertainty is large. The distribution reduces to the exponential when the expected value of the MT response function is zero. The distribution for the phase is also obtained in closed form. It reduces to the uniform distribution when the squared magnitude of the true MT response function is zero or its variance is very large. The phase distribution is symmetric and becomes increasingly concentrated as the variance decreases, although it is shorter-tailed than the Gaussian. The standard estimate for phase is unbiased. Confidence limits are derived from the distributions for magnitude squared response function and phase. Using a data set taken from the 2003 Kaapvaal transect, it is shown that the bias in the apparent resistivity is small and that confidence intervals obtained using the non-parametric delta method are very close to the true values obtained from the distributions. Thus, it appears that the computationally simple delta approximation provides accurate estimates for the confidence intervals, provided that the MT response function is obtained using an estimator that bounds the influence of extreme data.     

Statistical analysis of strong-motion accelerograms and its application to earthquake early-warning systems

In view of increasing damage due to earthquakes, and the current problems of earthquake prediction, real-time warning of strong ground motion is attracting more interest. In principle, it allows short-term warning of earthquakes while they are occurring. With warning times of up to tens of seconds it is possible to send alerts to potential areas of strong shaking before the arrival of the seismic waves and to mitigate the damage, but only if the seismic source parameters are determined rapidly. The major problem of an early-warning system is the real-time estimation of the earthquake's size.
We investigated digitized strong-motion accelerograms from 244 earthquakes that occurred in North and Central America between 1940 and 1986 to find out whether their initial portions reflected the size of the ongoing earthquake. Applying conventional methods of time-series analyses we calculate appropriate signal parameters and describe their uncertainties in relation to the magnitude and epicentral distance. The study reveals that the magnitude of an earthquake can be predicted from the first second of a single accelerogram within ±1.36 magnitude units. The uncertainty can be reduced to about ±0.5 magnitude units if a larger number (≥8) of accelerograms are available, which requires a dense network of seismic stations in areas of high seismic risk.     

Explicit, approximate expressions for the resolution and a posteriori covariance of massive tomographic systems

We present an approximate method to estimate the resolution, covariance and correlation matrix for linear tomographic systems Ax = b that are too large to be solved by singular value decomposition. An explicit expression for the approximate inverse matrix A ⁻ is found using one-step backprojections on the Penrose condition AA ⁻≈ I , from which we calculate the statistical properties of the solution. The computation of A ⁻ can easily be parallelized, each column being constructed independently.
The method is validated on small systems for which the exact covariance can still be computed with singular value decomposition. Though A ⁻ is not accurate enough to actually compute the solution x , the qualitative agreement obtained for resolution and covariance is sufficient for many purposes, such as rough assessment of model precision or the reparametrization of the model by the grouping of correlating parameters. We present an example for the computation of the complete covariance matrix of a very large (69 043 × 9610) system with 5.9 × 10⁶ non-zero elements in A . Computation time is proportional to the number of non-zero elements in A . If the correlation matrix is computed for the purpose of reparametrization by combining highly correlating unknowns x _i , a further gain in efficiency can be obtained by neglecting the small elements in A , but a more accurate estimation of the correlation requires a full treatment of even the smaller A _ij . We finally develop a formalism to compute a damped version of A ⁻.     

The detection of single-domain greigite (Fe3S4) using rotational remanent magnetization (RRM) and the effective gyro field (Bg): mineral magnetic and palaeomagnetic applications

The intensity of rotational remanent magnetization (RRM) acquired by single-domain greigite at a rotation frequency of 5 rps was combined with measurements of anhysteretic remanent magnetization (ARM) to calculate the effective biasing field ( Bg ) that produced the RRM. Samples of single-domain greigite had Bg values between -137 and -84 μT, and a MDF_RRM of c. 80 mT. By contrast, a suite of natural and synthetic ferrimagnetic iron oxide samples, including single-domain magnetite and y Fe₃O₄ tape particles, acquired Bg values between -3 and -14 μT, and MDF_RRM ranged between 43 and 68 mT (when RRM was acquired). Multidomain magnetite did not acquire a RRM. Bg values at 5 rps were calculated from previously published data for magnetite particles of different grain sizes, which revealed a minimum Bg value of -24 μT and a MDF_RRM of 57 mT for the finest fraction (0.2-0.8 μm in diameter). In a geological example, measurements of Bg and MDF_RRM were used to detect the presence of greigite in a 4 m long Late Weichselian sediment core. Variations in inclination, declination and the intensity of the natural remanent magnetization (NRM) correlate with changes in magnetic mineralogy.