Quadratic covariance estimation and equivalence of predictions

This article illustrates the use of linear and nonlinear regression models to obtain quadratic estimates of covariance parameters. These models lead to new insights into the motivation behind estimation methods, the relationships between different methods, and the relationship of covariance estimation to prediction. In particular, we derive the standard estimating equations for minimum norm quadratic unbiased translation invariant estimates (MINQUEs) from an appropriate linear model. Connections between the linear model, minimum variance quadratic unbiased translation invariant estimates (MIVQUEs), and MINQUEs are examined and we provide a minimum norm justification for the use of one-step normal theory maximum likelihood estimates. A nonlinear regression model is used to define MINQUEs for nonlinear covariance structures and obtain REML estimates. Finally, the equivalence of predictions under various models is examined when covariance parameters are estimated. In particular, we establish that when using MINQUE, iterative MINQUE, or restricted maximum likelihood (REML) estimates, the choice between a stationary covariance function and an intrinsically stationary semivariogram is irrelevant to predictions and estimated prediction variances.     

Minimum norm quadratic estimation of components of spatial covariance

Quadratic estimators of components of a nested spatial covariance function are presented. Estimators are unbiased and possess a minimum norm property. Inversion of a covariance matrix is required but, by assuming that spatial correlation is absent, a priori, matrix inversion can be avoided. The loss of efficiency that results from this assumption is discussed. Methods can be generalized to include estimation of components of a generalized polynomial covariance assuming the underlying process to be an intrinsic random function. Particular attention is given to the special case where just two components of spatial covariance exist, one of which represents a nugget effect.     

The Correlation Structure of Matheron's Classical Variogram Estimator Under Elliptically Contoured Distributions

The classical variogram estimator proposed by Matheron can be written as a quadratic form of the observations. When data have an elliptically contoured distribution with constant mean, the correlation between the classical variogram estimator at two different lags is a function of the spatial design matrix, the covariance matrix, and the kurtosis. Several specific cases are studied closely. A subclass of elliptically contoured distributions with a particular family of covariance matrices is shown to possess exactly the same correlation structure for the classical variogram estimator as the multivariate independent Gaussian distribution. The consequences on variogram fitting by generalized least squares are discussed.     

The Correlation Structure of Matheron's Classical Variogram Estimator Under Elliptically Contoured Distributions

The classical variogram estimator proposed by Matheron can be written as a quadratic form of the observations. When data have an elliptically contoured distribution with constant mean, the correlation between the classical variogram estimator at two different lags is a function of the spatial design matrix, the covariance matrix, and the kurtosis. Several specific cases are studied closely. A subclass of elliptically contoured distributions with a particular family of covariance matrices is shown to possess exactly the same correlation structure for the classical variogram estimator as the multivariate independent Gaussian distribution. The consequences on variogram fitting by generalized least squares are discussed.     

Minimum-variance unbiased quadratic estimation of covariances of regionalized variables

The parameters of covariance functions (or variograms) of regionalized variables must be determined before linear unbiased estimation can be applied. This work examines the problem of minimum-variance unbiased quadratic estimation of the parameters of ordinary or generalized covariance functions of regionalized variables. Attention is limited to covariance functions that are linear in the parameters and the normality assumption is invoked when fourth moments of the data need to be calculated. The main contributions of this work are (1) it shows when and in what sense minimum-variance unbiased quadratic estimation can be achieved, and (2) it yields a well-founded, practicable, and easy-to-automate methodology for the estimation of parameters of covariance functions. Results of simulation studies are very encouraging.     

Variance–Covariance Matrix of the Experimental Variogram: Assessing Variogram Uncertainty

Assessment of the sampling variance of the experimental variogram is an important topic in geostatistics as it gives the uncertainty of the variogram estimates. This assessment, however, is repeatedly overlooked in most applications mainly, perhaps, because a general approach has not been implemented in the most commonly used software packages for variogram analysis. In this paper the authors propose a solution that can be implemented easily in a computer program, and which, subject to certain assumptions, is exact. These assumptions are not very restrictive: second-order stationarity (the process has a finite variance and the variogram has a sill) and, solely for the purpose of evaluating fourth-order moments, a Gaussian distribution for the random function. The approach described here gives the variance–covariance matrix of the experimental variogram, which takes into account not only the correlation among the experiemental values but also the multiple use of data in the variogram computation. Among other applications, standard errors may be attached to the variogram estimates and the variance–covariance matrix may be used for fitting a theoretical model by weighted, or by generalized, least squares. Confidence regions that hold a given confidence level for all the variogram lag estimates simultaneously have been calculated using the Bonferroni method for rectangular intervals, and using the multivariate Gaussian assumption for K-dimensional elliptical intervals (where K is the number of experimental variogram estimates). A general approach for incorporating the uncertainty of the experimental variogram into the uncertainty of the variogram model parameters is also shown. A case study with rainfall data is used to illustrate the proposed approach.     

Updating of Population Parameters and Credibility of Discriminant Analysis

The uncertainty of classification in discriminant analysis may result from the original characteristics of the phenomena studied, the approach of inferring population parameters, and the credibility of the parameters which are estimated by geologist or statistician. A credibility function and a significance function are proposed. Both can be used to appraise the uncertainty of classification. The former is involved with the uncertainty resulting from the errors in the reward-penalty matrix, while the latter may be involved with the uncertainty resulting from the original characteristics of the phenomena studied and the statistical approach. Inappropriate classified results may be originated from the bias estimates of population parameters (mean vector and covariance matrix), which are estimated by bias samples. These bias estimates can be updated by constraining the varying region of the mean vector. The equations for updating Bayesian estimates of the mean vector and the covariance matrix are demonstrated if the mean vector is restricted to a subregion of the entire real space. Results for a gas reservoir indicate that the discriminant rules based on the updated equations are more efficient than the traditional discriminant rules.     

Variogram Fitting by Generalized Least Squares Using an Explicit Formula for the Covariance Structure

In the context of spatial statistics, the classical variogram estimator proposed by Matheron can be written as a quadratic form of the observations. If data are Gaussian with constant mean, then the correlation between the classical variogram estimator at two different lags is a function of the spatial design matrix and the variance matrix. When data are independent with unidimensional and regular support, an explicit formula for this correlation is available. The same is true for a multidimensional and regular support as can be shown by using Kronecker products of matrices. As variogram fitting is a crucial stage for correct spatial prediction, it is proposed to use a generalized least squares method with an explicit formula for the covariance structure (GLSE). A good approximation of the covariance structure is achieved by taking account of the explicit formula for the correlation in the independent situation. Simulations are carried out with several types of underlying variograms, as well as with outliers in the data. Results show that this technique (GLSE), combined with a robust estimator of the variogram, improves the fit significantly.     

The Second-Order Stationary Universal Kriging Model Revisited

Universal kriging originally was developed for problems of spatial interpolation if a drift seemed to be justified to model the experimental data. But its use has been questioned in relation to the bias of the estimated underlying variogram (variogram of the residuals), and furthermore universal kriging came to be considered an old-fashioned method after the theory of intrinsic random functions was developed. In this paper the model is reexamined together with methods for handling problems in the inference of parameters. The efficiency of the inference of covariance parameters is shown in terms of bias, variance, and mean square error of the sampling distribution obtained by Monte Carlo simulation for three different estimators (maximum likelihood, bias corrected maximum likelihood, and restricted maximum likelihood). It is shown that unbiased estimates for the covariance parameters may be obtained but if the number of samples is small there can be no guarantee of good estimates (estimates close to the true value) because the sampling variance usually is large. This problem is not specific to the universal kriging model but rather arises in any model where parameters are inferred from experimental data. The validity of the estimates may be evaluated statistically as a risk function as is shown in this paper.     

考虑量测误差的反分析方法

由于量测误差的存在,使反分析参数的可靠性降低。分析了路堤工程量测误差的来源以及确定方法,以及采用极大似然法导出了包含协方差矩阵的最小化目标函数表达式,并推导了测斜仪的协方差矩阵。将最小化目标函数值作为遗传算法的适应度值可以优化土性参数。给出一个实例,得出考虑测斜仪误差比未考虑测斜仪误差的计算结果更符合实际情况。     

Estimation of a multidimensional covariance function in case of anisotropy

A model of a multivariate covariance function with an ellipsoidal directional correlation scale has been developed. The axes of the ellipsoidal scale are related to the eigenvalues and eigenvectors of a matrix B which characterizes the ellipsoid of the range of influence. The matrix B is found to be related to a matrix T which can be estimated directly from sparse sampling data and can be used to determine estimates of the matrix B. The method has been applied to both two-dimensional and three-dimensional cases. The numerical results show that the satisfactory accuracy is obtained with sparse sampling data from an anisotropic random function.     

增长曲线模型均值矩阵的一种新的相对效率

对于一般增长曲线模型,当设计阵和协方差阵满足一定的条件时,均值矩阵(XBZ)的最小二乘估计(LSE)与最佳线性无偏估计(BLUE)相等,可用LSE代替BLUE。当这些条件不满足时,用LSE代替BLUE就要蒙受一些损失,有时这种损失很大,因而研究这种损失的大小显得相当重要,通常用它们协方差阵的差异来度量,从不同的准则出发,可以定义不同的相对效率,利用广义行列式定义了一种新的相对效率,并给出了它们的下界。     

Combining sensitivities and prior information for covariance localization in the ensemble Kalman filter for petroleum reservoir applications

Sampling errors can severely degrade the reliability of estimates of conditional means and uncertainty quantification obtained by the application of the ensemble Kalman filter (EnKF) for data assimilation. A standard recommendation for reducing the spurious correlations and loss of variance due to sampling errors is to use covariance localization. In distance-based localization, the prior (forecast) covariance matrix at each data assimilation step is replaced with the Schur product of a correlation matrix with compact support and the forecast covariance matrix. The most important decision to be made in this localization procedure is the choice of the critical length(s) used to generate this correlation matrix. Here, we give a simple argument that the appropriate choice of critical length(s) should be based both on the underlying principal correlation length(s) of the geological model and the range of the sensitivity matrices. Based on this result, we implement a procedure for covariance localization and demonstrate with a set of distinctive reservoir history-matching examples that this procedure yields improved results over the standard EnKF implementation and over covariance localization with other choices of critical length.     

Bayesian updating revisited

Bayesian updating methods provide an alternate philosophy to the characterization of the input variables of a stochastic mathematical model. Here, a priori values of statistical parameters are assumed on subjective grounds or by analysis of a data base from a geologically similar area. As measurements become available during site investigations, updated estimates of parameters characterizing spatial variability are generated. However, in solving the traditional updating equations, an updated covariance matrix may be generated that is not positive-definite, particularly when observed data errors are small. In addition, measurements may indicate that initial estimates of the statistical parameters are poor. The traditional procedure does not have a facility to revise the parameter estimates before the update is carried out. alternatively, Bayesian updating can be viewed as a linear inverse problem that minimizes a weighted combination of solution simplicity and data misfit. Depending on the weight given to the a priori information, a different solution is generated. A Bayesian updating procedure for log-conductivity interpolation that uses a singular value decomposition (SVD) is presented. An efficient and stable algorithm is outlined that computes the updated log-conductivity field and the a posteriori covariance of the estimated values (estimation errors). In addition, an information density matrix is constructed that indicates how well predicted data match observations. Analysis of this matrix indicates the relative importance of the observed data. The SVD updating procedure is used to interpolate the log-conductivity fields of a series of hypothetical aquifers to demonstrate pitfalls and possibilities of the method.     

Time-dependent Inversion of Surface Subsidence due to Dynamic Reservoir Compaction

We introduce a novel, time-dependent inversion scheme for resolving temporal reservoir pressure drop from surface subsidence observations (from leveling or GPS data, InSAR, tiltmeter monitoring) in a single procedure. The theory is able to accommodate both the absence of surface subsidence estimates at sites at one or more epochs as well as the introduction of new sites at any arbitrary epoch. Thus, all observation sites with measurements from at least two epochs are utilized. The method uses both the prior model covariance matrix and the data covariance matrix, which incorporates the spatial and temporal correlations between model parameters and data, respectively. The incorporation of the model covariance implicitly guarantees smoothness of the model estimate, while maintaining specific geological features like sharp boundaries. Taking these relations into account through the model covariance matrix enhances the influence of the data on the inverted model estimate. This leads to a better defined and interpretable model estimate. The time-dependent aspect of the method yields a better constrained model estimate and makes it possible to identify non-linear acceleration or delay in reservoir compaction. The method is validated by a synthetic case study based on an existing gas reservoir with a highly variable transmissibility at the free water level. The prior model covariance matrix is based on a Monte Carlo simulation of the geological uncertainty in the transmissibility.     

Bayesian updating revisited

Bayesian updating methods provide an alternate philosophy to the characterization of the input variables of a stochastic mathematical model. Here, a priori values of statistical parameters are assumed on subjective grounds or by analysis of a data base from a geologically similar area. As measurements become available during site investigations, updated estimates of parameters characterizing spatial variability are generated. However, in solving the traditional updating equations, an updated covariance matrix may be generated that is not positive-definite, particularly when observed data errors are small. In addition, measurements may indicate that initial estimates of the statistical parameters are poor. The traditional procedure does not have a facility to revise the parameter estimates before the update is carried out. alternatively, Bayesian updating can be viewed as a linear inverse problem that minimizes a weighted combination of solution simplicity and data misfit. Depending on the weight given to the a priori information, a different solution is generated. A Bayesian updating procedure for log-conductivity interpolation that uses a singular value decomposition (SVD) is presented. An efficient and stable algorithm is outlined that computes the updated log-conductivity field and the a posteriori covariance of the estimated values (estimation errors). In addition, an information density matrix is constructed that indicates how well predicted data match observations. Analysis of this matrix indicates the relative importance of the observed data. The SVD updating procedure is used to interpolate the log-conductivity fields of a series of hypothetical aquifers to demonstrate pitfalls and possibilities of the method.     

Computationally efficient generation of gaussian conditional simulations over regular sample grids

The generation over two-dimensional grids of normally distributed random fields conditioned on available data is often required in reservoir modeling and mining investigations. Such fields can be obtained from application of turning band or spectral methods. However, both methods have limitations. First, they are only asymptotically exact in that the ensemble of realizations has the correlation structure required only if enough harmonics are used in the spectral method, or enough lines are generated in the turning bands approach. Moreover, the spectral method requires fine tuning of process parameters. As for the turning bands method, it is essentially restricted to processes with stationary and radially symmetric correlation functions. Another approach, which has the advantage of being general and exact, is to use a Cholesky factorization of the covariance matrix representing grid points correlation. For fields of large size, however, the Cholesky factorization can be computationally prohibitive. In this paper, we show that if the data are stationary and generated over a grid with regular mesh, the structure of the data covariance matrix can be exploited to significantly reduce the overall computational burden of conditional simulations based on matrix factorization techniques. A feature of this approach is its computational simplicity and suitability to parallel implementation.     

A nonparametric view of generalized covariances for kriging

Fitting trend and error covariance structure iteratively leads to bias in the estimated error variogram. Use of generalized increments overcomes this bias. Certain generalized increments yield difference equations in the variogram which permit graphical checking of the model. These equations extend to the case where errors are intrinsic random functions of order k, k=1, 2, ..., and an unbiased nonparametric graphical approach for investigating the generalized covariance function is developed. Hence, parametric models for the generalized covariance produced by BLUEPACK-3D or other methods may be assessed. Methods are illustrated on a set of coal ash data and a set of soil pH data.     

Generalized covariance functions in estimation

I discuss the role of generalized covariance functions in best linear unbiased estimation and methods for their selection. It is shown that the experimental variogram (or covariance function) of the detrended data can be used to obtain a preliminary estimate of the generalized covariance function without iterations and I discuss the advantages of other parameter estimation methods.