Design-based versus model-based sampling strategies: Comment on R. J. Barnes' “bounding the required sample size for geologic site characterization”

Two fundamentally different sources of randomness exist on which design and inference in spatial sampling can be based: (a) variation that would occur on resampling the same spatial population with other sampling configurations generated by the same design, and (b) variation occurring on sampling other populations, hypothetically generated by the same spatial model, using the same sampling configuration. The former leads to the design-based approach, which uses classical sampling theory; the latter leads to the model-based approach and uses geostatistical theory. Failure to recognize these two sources of randomness causes misunderstanding about dependence of variables and the role of randomization in sampling, unwarranted narrowing down the choice of sampling strategies to those that are model-based, and abuse in simulation experiments. This is exemplified in Barnes' publication on the required sample size for geologic site characterization by nonparametric tolerance intervals. A basic design-based strategy like Simple Random Sampling is shown to require smaller sample sizes than the model-based strategy advocated by Barnes. In addition, Simple Random Sampling is completely robust against model errors and less complicated.     

Conditional Latin Hypercube Simulation of (Log)Gaussian Random Fields

In earth and environmental sciences applications, uncertainty analysis regarding the outputs of models whose parameters are spatially varying (or spatially distributed) is often performed in a Monte Carlo framework. In this context, alternative realizations of the spatial distribution of model inputs, typically conditioned to reproduce attribute values at locations where measurements are obtained, are generated via geostatistical simulation using simple random (SR) sampling. The environmental model under consideration is then evaluated using each of these realizations as a plausible input, in order to construct a distribution of plausible model outputs for uncertainty analysis purposes. In hydrogeological investigations, for example, conditional simulations of saturated hydraulic conductivity are used as input to physically-based simulators of flow and transport to evaluate the associated uncertainty in the spatial distribution of solute concentration. Realistic uncertainty analysis via SR sampling, however, requires a large number of simulated attribute realizations for the model inputs in order to yield a representative distribution of model outputs; this often hinders the application of uncertainty analysis due to the computational expense of evaluating complex environmental models. Stratified sampling methods, including variants of Latin hypercube sampling, constitute more efficient sampling aternatives, often resulting in a more representative distribution of model outputs (e.g., solute concentration) with fewer model input realizations (e.g., hydraulic conductivity), thus reducing the computational cost of uncertainty analysis. The application of stratified and Latin hypercube sampling in a geostatistical simulation context, however, is not widespread, and, apart from a few exceptions, has been limited to the unconditional simulation case. This paper proposes methodological modifications for adopting existing methods for stratified sampling (including Latin hypercube sampling), employed to date in an unconditional geostatistical simulation context, for the purpose of efficient conditional simulation of Gaussian random fields. The proposed conditional simulation methods are compared to traditional geostatistical simulation, based on SR sampling, in the context of a hydrogeological flow and transport model via a synthetic case study. The results indicate that stratified sampling methods (including Latin hypercube sampling) are more efficient than SR, overall reproducing to a similar extent statistics of the conductivity (and subsequently concentration) fields, yet with smaller sampling variability. These findings suggest that the proposed efficient conditional sampling methods could contribute to the wider application of uncertainty analysis in spatially distributed environmental models using geostatistical simulation.     

A Fast Approximation for Seismic Inverse Modeling: Adaptive Spatial Resampling

Seismic inverse modeling, which transforms appropriately processed geophysical data into the physical properties of the Earth, is an essential process for reservoir characterization. This paper proposes a work flow based on a Markov chain Monte Carlo method consistent with geology, well-logs, seismic data, and rock-physics information. It uses direct sampling as a multiple-point geostatistical method for generating realizations from the prior distribution, and Metropolis sampling with adaptive spatial resampling to perform an approximate sampling from the posterior distribution, conditioned to the geophysical data. Because it can assess important uncertainties, sampling is a more general approach than just finding the most likely model. However, since rejection sampling requires a large number of evaluations for generating the posterior distribution, it is inefficient and not suitable for reservoir modeling. Metropolis sampling is able to perform an equivalent sampling by forming a Markov chain. The iterative spatial resampling algorithm perturbs realizations of a spatially dependent variable, while preserving its spatial structure by conditioning to subset points. However, in most practical applications, when the subset conditioning points are selected at random, it can get stuck for a very long time in a non-optimal local minimum. In this paper it is demonstrated that adaptive subset sampling improves the efficiency of iterative spatial resampling. Depending on the acceptance/rejection criteria, it is possible to obtain a chain of geostatistical realizations aimed at characterizing the posterior distribution with Metropolis sampling. The validity and applicability of the proposed method are illustrated by results for seismic lithofacies inversion on the Stanford VI synthetic test sets.     

Direct Pattern-Based Simulation of Non-stationary Geostatistical Models

Non-stationary models often capture better spatial variation of real world spatial phenomena than stationary ones. However, the construction of such models can be tedious as it requires modeling both statistical trend and stationary stochastic component. Non-stationary models are an important issue in the recent development of multiple-point geostatistical models. This new modeling paradigm, with its reliance on the training image as the source for spatial statistics or patterns, has had considerable practical appeal. However, the role and construction of the training image in the non-stationary case remains a problematic issue from both a modeling and practical point of view. In this paper, we provide an easy to use, computationally efficient methodology for creating non-stationary multiple-point geostatistical models, for both discrete and continuous variables, based on a distance-based modeling and simulation of patterns. In that regard, the paper builds on pattern-based modeling previously published by the authors, whereby a geostatistical realization is created by laying down patterns as puzzle pieces on the simulation grid, such that the simulated patterns are consistent (in terms of a similarity definition) with any previously simulated ones. In this paper we add the spatial coordinate to the pattern similarity calculation, thereby only borrowing patterns locally from the training image instead of globally. The latter would entail a stationary assumption. Two ways of adding the geographical coordinate are presented, (1) based on a functional that decreases gradually away from the location where the pattern is simulated and (2) based on an automatic segmentation of the training image into stationary regions. Using ample two-dimensional and three-dimensional case studies we study the behavior in terms of spatial and ensemble uncertainty of the generated realizations.     

A geostatistical basis for spatial weighting in multivariate classification

Earth scientists and land managers often wish to group sampling sites that are both similar with respect to their properties and near to one another on the ground. This paper outlines the geostatistical rationale for such spatial grouping and describes a multivariate procedure to implement it. Sample variograms are calculated from the original data or their leading principal components and then the parameters of the underlying functions are estimated. A dissimilarity matrix is computed for all sampling sites, preferably using Gower's general similarity coefficient. Dissimilarities are then modified using the variogram to incorporate the form and extent of spatial variation. A nonhierarchical classification of sampling sites is performed on the leading latent vectors of the modified dissimilarity matrix by dynamic clustering to an optimum. The technique is illustrated with results of its application to soil survey data from two small areas in Britain and from a transect. In the case of the latter results of spatially weighted classifications are compared with those of strict segmentation. An appendix lists a Genstat program for a spatially constrained classification using a spherical variogram as an example.     

Robust Resampling Confidence Intervals for Empirical Variograms

The variogram function is an important measure of the spatial dependencies of a geostatistical or other spatial dataset. It plays a central role in kriging, designing spatial studies, and in understanding the spatial properties of geological and environmental phenomena. It is therefore important to understand the variability attached to estimates of the variogram. Existing methods for constructing confidence intervals around the empirical variogram either rely on strong assumptions, such as normality or known variogram function, or are based on resampling blocks and subject to edge effect biases. This paper proposes two new procedures for addressing these concerns: a quasi-block-bootstrap and a quasi-block-jackknife. The new methods are based on transforming the data to decorrelate it based on a fitted variogram model, resampling blocks from the decorrelated data, and then recorrelating. The coverage properties of the new confidence intervals are compared by simulation to a number of existing resampling-based intervals. The proposed quasi-block-jackknife confidence interval is found to have the best properties of all of the methods considered across a range of scenarios, including normally and lognormally distributed data and misspecification of the variogram function used to decorrelate the data.     

On the stability of the geostatistical method

A thorough geostatistical analysis of spatial data, observed at given spatial locations, includes exploratory data analysis, spatial-model building, diagnosing the model fit, and inference on unknown model parameters or unobserved values (at known locations). Using results from mathematical analysis, exact and asymptotic distribution theory, and simulation studies, we argue that, when used sensibly, the geostatistical method is reassuringly stable.     

Geostatistical Simulation of Regionalized Pore-Size Distributions Using Min/Max Autocorrelation Factors

In many fields of the Earth Sciences, one is interested in the distribution of particle or void sizes within samples. Like many other geological attributes, size distributions exhibit spatial variability, and it is convenient to view them within a geostatistical framework, as regionalized functions or curves. Since they rarely conform to simple parametric models, size distributions are best characterized using their raw spectrum as determined experimentally in the form of a series of abundance measures corresponding to a series of discrete size classes. However, the number of classes may be large and the class abundances may be highly cross-correlated. In order to model the spatial variations of discretized size distributions using current geostatistical simulation methods, it is necessary to reduce the number of variables considered and to render them uncorrelated among one another. This is achieved using a principal components-based approach known as Min/Max Autocorrelation Factors (MAF). For a two-structure linear model of coregionalization, the approach has the attractive feature of producing orthogonal factors ranked in order of increasing spatial correlation. Factors consisting largely of noise and exhibiting pure nugget–effect correlation structures are isolated in the lower rankings, and these need not be simulated. The factors to be simulated are those capturing most of the spatial correlation in the data, and they are isolated in the highest rankings. Following a review of MAF theory, the approach is applied to the modeling of pore-size distributions in partially welded tuff. Results of the case study confirm the usefulness of the MAF approach for the simulation of large numbers of coregionalized variables.     

Optimized Sample Schemes for Geostatistical Surveys

Sample schemes used in geostatistical surveys must be suitable for both variogram estimation and kriging. Previously schemes have been optimized for one of these steps in isolation. Ordinary kriging generally requires the sampling locations to be evenly dispersed over the region. Variogram estimation requires a more irregular pattern of sampling locations since comparisons must be made between measurements separated by all lags up to and beyond the range of spatial correlation. Previous studies have not considered how to combine these optimized schemes into a single survey and how to decide what proportion of sampling effort should be devoted to variogram estimation and what proportion devoted to kriging An expression for the total error in a geostatistical survey accounting for uncertainty due to both ordinary kriging and variogram uncertainty is derived. In the same manner as the kriging variance, this expression is a function of the variogram but not of the sampled response data. If a particular variogram is assumed the total error in a geostatistical survey may be estimated prior to sampling. We can therefore design an optimal sample scheme for the combined processes of variogram estimation and ordinary kriging by minimizing this expression. The minimization is achieved by spatial simulated annealing. The resulting sample schemes ensure that the region is fairly evenly covered but include some close pairs to analyse the spatial correlation over short distances. The form of these optimal sample schemes is sensitive to the assumed variogram. Therefore a Bayesian approach is adopted where, rather than assuming a single variogram, we minimize the expected total error over a distribution of plausible variograms. This is computationally expensive so a strategy is suggested to reduce the number of computations required     

Modelling Short‐Scale Variability and Uncertainty During Mineral Resource Estimation Using a Novel Fuzzy Estimation Technique

For mineral resource assessment, techniques based on fuzzy logic are attractive because they are capable of incorporating uncertainty associated with measured variables and can also quantify the uncertainty of the estimated grade, tonnage etc. The fuzzy grade estimation model is independent of the distribution of data, avoiding assumptions and constraints made during advanced geostatistical simulation, e.g., the turning bands method. Initially, fuzzy modelling classifies the data using all the component variables in the data set. We adopt a novel approach by taking into account the spatial irregularity of mineralisation patterns using the Gustafson–Kessel classification algorithm. The uncertainty at the point of estimation was derived through antecedent memberships in the input space (i.e., spatial coordinates) and transformed onto the output space (i.e., grades) through consequent membership at the point of estimation. Rather than probabilistic confidence intervals, this uncertainty was expressed in terms of fuzzy memberships, which indicated the occurrence of mixtures of different mineralogical phases at the point of estimation. Data from different sources (other than grades) could also be utilised during estimation. Application of the proposed technique on a real data set gave results that were comparable to those obtained from a turning bands simulation.     

Estimation of non-ergodic variograms and their sampling variance by design-based sampling strategies

Design-based sampling strategies based on classical sampling theory offer unprecedented potentials for estimation of non-ergodic variograms. Unbiased and uncorrelated estimates of the semivariance at the selected lags and of its sampling variance can be simply obtained. These estimates are robust against deviations from an assumed spatial autocorrelation model. The same holds for the variogram model parameters and their sampling (co)variances. Moreover, an objective measure for lack of fit of the fitted model can simply be derived. The estimators for two basic sampling designs, simple random sampling and stratified simple random sampling of pairs of points, are presented. The first has been tested in real world for estimating the non-ergodic variograms of three soil properties. The parameters of variogram models and their sampling (co)variances were estimated with 72 pairs of points distributed over six lags.     

地质统计学新进展

地质统计学(空间信息统计学)是数学地质领域中一门发展迅速且有着广泛应用前景的新兴科学。结合地质统计学发展现状,对地质统计学的新进展进行了研究,从地质统计学理论体系、应用及软件开发等方面探讨了地质统计学的发展前缘。并指出时空多元技术、条件模拟、非参数和非线性将是地质统计学今后的发展趋势。     

Sampling design optimization for spatial functions

A new procedure is presented for minimizing the sampling requirements necessary to estimate a mappable spatial function at a specified level of accuracy. The technique is based on universal kriging, an estimation method within the theory of regionalized variables. Neither actual implementation of the sampling nor universal kriging estimations are necessary to make an optimal design. The average standard errorand maximum standard error of estimationover the sampling domain are used as global indices of sampling efficiency. The procedure optimally selects those parameters controlling the magnitude of the indices, including the density and spatial pattern of the sample elements and the number of nearest sample elements used in the estimation. As an illustration, the network of observation wells used to monitor the water table in the Equus Beds of Kansas is analyzed and an improved sampling pattern suggested. This example demonstrates the practical utility of the procedure, which can be applied equally well to other spatial sampling problems, as the procedure is not limited by the nature of the spatial function.     

Variography of submarine morphology: Problems of deregularization, and cartographical implications

During the German Antarctic Expedition VI (leg 3, December 1987 to March 1988), bathymetric surveys were made in the Weddell Sea by the SEABEAM sonar system. For the first time geostatistical methods were applied in the SEABEAM-postprocessing. The investigations of variography that were necessary prior to the cartographical-geomorphological evaluation shed new light on classical geostatistical concerns. SEABEAM data provide a good example of a mean square, differentiable regionalized variable, where data are sampled over a two-dimensional support due to the technique of the sonar device. By deregularizations of the sample variograms, spatial continuity can be shown to be a property of seafloor depth as well as a point variable. The results are discussed in a sedimentological context. As an application of the regional variogram analyses, large-scale kriged bathymetric maps are presented.     

The Bayesian bridge between simple and universal kriging

Kriging techniques are suited well for evaluation of continuous, spatial phenomena. Bayesian statistics are characterized by using prior qualified guesses on the model parameters. By merging kriging techniques and Bayesian theory, prior guesses may be used in a spatial setting. Partial knowledge of model parameters defines a continuum of models between what is named simple and universal kriging in geostatistical terminology. The Bayesian approach to kriging is developed and discussed, and a case study concerning depth conversion of seismic reflection times is presented.     

On the Use of Non-Euclidean Distance Measures in Geostatistics

In many scientific disciplines, straight line, Euclidean distances may not accurately describe proximity relationships among spatial data. However, non-Euclidean distance measures must be used with caution in geostatistical applications. A simple example is provided to demonstrate there are no guarantees that existing covariance and variogram functions remain valid (i.e. positive definite or conditionally negative definite) when used with a non-Euclidean distance measure. There are certain distance measures that when used with existing covariance and variogram functions remain valid, an issue that is explored. The concept of isometric embedding is introduced and linked to the concepts of positive and conditionally negative definiteness to demonstrate classes of valid norm dependent isotropic covariance and variogram functions, results many of which have yet to appear in the mainstream geostatistical literature or application. These classes of functions extend the well known classes by adding a parameter to define the distance norm. In practice, this distance parameter can be set a priori to represent, for example, the Euclidean distance, or kept as a parameter to allow the data to choose the metric. A simulated application of the latter is provided for demonstration. Simulation results are also presented comparing kriged predictions based on Euclidean distance to those based on using a water metric.     

Extrapolating the Fractal Characteristics of an Image Using Scale-Invariant Multiple-Point Statistics

The resolution of measurement devices can be insufficient for certain purposes. We propose to stochastically simulate spatial features at scales smaller than the measurement resolution. This is accomplished using multiple-point geostatistical simulation (direct sampling in the present case) to interpolate values at the target scale. These structures are inferred using hypothesis of scale invariance and stationarity on the spatial patterns found at the coarse scale. The proposed multiple-point super-resolution mapping method is able to deal with “both continuous and categorical variables,” and can be extended to multivariate problems. The advantages and limitations of the approach are illustrated with examples from satellite imaging.     

Measures of Parameter Uncertainty in Geostatistical Estimation and Geostatistical Optimal Design

Studies of site exploration, data assimilation, or geostatistical inversion measure parameter uncertainty in order to assess the optimality of a suggested scheme. This study reviews and discusses measures for parameter uncertainty in spatial estimation. Most measures originate from alphabetic criteria in optimal design and were transferred to geostatistical estimation. Further rather intuitive measures can be found in the geostatistical literature, and some new measures will be suggested in this study. It is shown how these measures relate to the optimality alphabet and to relative entropy. Issues of physical and statistical significance are addressed whenever they arise. Computational feasibility and efficient ways to evaluate the above measures are discussed in this paper, and an illustrative synthetic case study is provided. A major conclusion is that the mean estimation variance and the averaged conditional integral scale are a powerful duo for characterizing conditional parameter uncertainty, with direct correspondence to the well-understood optimality alphabet. This study is based on cokriging generalized to uncertain mean and trends because it is the most general representative of linear spatial estimation within the Bayesian framework. Generalization to kriging and quasi-linear schemes is straightforward. Options for application to non-Gaussian and non-linear problems are discussed.     

Geostatistical methods to identify and map spatial variations of soil salinity

The problem of estimating and predicting spatial distribution of a spatial stochastic process, observed at irregular locations in space, is considered in this paper. Environmental variables usually show spatial dependencies among observations, with lead one to use geostatistical methods to model the spatial distributions of those observations. This is particularly important in the study of soil properties and their spatial variability. In this study geostatistical techniques were used to describe the spatial dependence and to quantify the scale and intensity of spatial variations of soil properties, which provide the essential spatial information for local estimation. In this contribution, we propose a spatial Gaussian linear mixed model that involves (a) a non-parametric term for accounting deterministic trend due to exogenous variables and (b) a parametric component for defining the purely spatial random variation due possibly to latent spatial processes. We focus here on the analysis of the relationship between soil electrical conductivity and Na content to identify spatial variations of soil salinity. This analysis can be useful for agricultural and environmental land management.