Confronting uncertainty in model-based geostatistics using Markov Chain Monte Carlo simulation

This paper demonstrates for the first time the use of Markov Chain Monte Carlo (MCMC) simulation for parameter inference in model-based soil geostatistics. We implemented the recently developed DiffeRential Evolution Adaptive Metropolis (DREAM) algorithm to jointly summarize the posterior distribution of variogram parameters and the coefficients of a linear spatial model, and derive estimates of predictive uncertainty. The DREAM method runs multiple different Markov chains in parallel and jumps in each chain are generated from a discrete proposal distribution containing a fixed multiple of the difference of the states of randomly chosen pairs of other chains. This approach automatically scales the orientation and scale of the proposal distribution, and is especially designed to maintain detailed balance and ergodicity, thereby generating an exact approximation of the posterior probability density function (pdf) of the parameters of the linear model and variogram. This approach is tested using three different data sets from Australia involving variogram estimation of soil thickness, kriging of soil pH, and spatial prediction of soil organic carbon content. The results showed some advantages of MCMC over the conventional method of moments and residual maximum likelihood (REML) estimation. The posterior pdf derived with MCMC conveys important information about parameter uncertainty, multi-dimensional parameter correlation, and thus how many significant parameters are warranted by the calibration data. Parameter uncertainty constitutes only a small part of total prediction uncertainty for the case studies considered here. The prediction accuracies using MCMC and REML are similar. The variogram estimated using conventional approaches (method of moments, and without simulation) lies within the 95% prediction uncertainty interval of the posterior distribution derived with DREAM. Altogether our results show that conventional kriging and regression-kriging still remain a viable option for production mapping.

[1]  H. Wackernagel,et al.  Mapping temperature using kriging with external drift: Theory and an example from scotland , 1994 .

[2]  Alex B. McBratney,et al.  Comparison of several spatial prediction methods for soil pH , 1987 .

[3]  Jasper A. Vrugt,et al.  High‐dimensional posterior exploration of hydrologic models using multiple‐try DREAM(ZS) and high‐performance computing , 2012 .

[4]  R. Reese Geostatistics for Environmental Scientists , 2001 .

[5]  J. R. Wallis,et al.  An Approach to Statistical Spatial-Temporal Modeling of Meteorological Fields , 1994 .

[6]  David B. Dunson,et al.  Bayesian Data Analysis , 2010 .

[7]  Johan Alexander Huisman,et al.  Hydraulic properties of a model dike from coupled Bayesian and multi-criteria hydrogeophysical inversion , 2010 .

[8]  George Kuczera,et al.  Bayesian analysis of input uncertainty in hydrological modeling: 1. Theory , 2006 .

[9]  A. McBratney,et al.  Further results on prediction of soil properties from terrain attributes: heterotopic cokriging and regression-kriging , 1995 .

[10]  Budiman Minasny,et al.  Mapping continuous depth functions of soil carbon storage and available water capacity , 2009 .

[11]  Anthony N. Pettitt,et al.  Sampling Designs for Estimating Spatial Variance Components , 1993 .

[12]  W. Bouten,et al.  Validity of First-Order Approximations to Describe Parameter Uncertainty in Soil Hydrologic Models , 2002 .

[13]  R. Lark,et al.  Geostatistics for Environmental Scientists , 2001 .

[14]  H. Haario,et al.  An adaptive Metropolis algorithm , 2001 .

[15]  J. Gallant,et al.  A multiresolution index of valley bottom flatness for mapping depositional areas , 2003 .

[16]  W. K. Hastings,et al.  Monte Carlo Sampling Methods Using Markov Chains and Their Applications , 1970 .

[17]  D. Rubin,et al.  Inference from Iterative Simulation Using Multiple Sequences , 1992 .

[18]  D. Higdon,et al.  Accelerating Markov Chain Monte Carlo Simulation by Differential Evolution with Self-Adaptive Randomized Subspace Sampling , 2009 .

[19]  A. Lawson,et al.  Validation of Bayesian kriging of arsenic, chromium, lead, and mercury surface soil concentrations based on internode sampling. , 2009, Environmental science & technology.

[20]  R. Larka,et al.  Fitting a linear model of coregionalization for soil properties using simulated annealing , 2002 .

[21]  Heikki Haario,et al.  Adaptive proposal distribution for random walk Metropolis algorithm , 1999, Comput. Stat..

[22]  R. M. Lark,et al.  Estimating Variogram Uncertainty , 2004 .

[23]  S. Sorooshian,et al.  A Shuffled Complex Evolution Metropolis algorithm for optimization and uncertainty assessment of hydrologic model parameters , 2002 .

[24]  Cajo J. F. ter Braak,et al.  Differential Evolution Markov Chain with snooker updater and fewer chains , 2008, Stat. Comput..

[25]  R. Lark,et al.  On spatial prediction of soil properties in the presence of a spatial trend: the empirical best linear unbiased predictor (E‐BLUP) with REML , 2006 .

[26]  Cajo J. F. ter Braak,et al.  A Markov Chain Monte Carlo version of the genetic algorithm Differential Evolution: easy Bayesian computing for real parameter spaces , 2006, Stat. Comput..

[27]  Spectral tempering to model non-stationary covariance of nitrous oxide emissions from soil using continuous or categorical explanatory variables at a landscape scale. , 2010 .

[28]  B. Minasny,et al.  The Matérn function as a general model for soil variograms , 2005 .

[29]  B. Minasny,et al.  Spatial prediction of soil properties using EBLUP with the Matérn covariance function , 2007 .

[30]  A. Raftery,et al.  How Many Iterations in the Gibbs Sampler , 1991 .

[31]  N. Cressie,et al.  Mean squared prediction error in the spatial linear model with estimated covariance parameters , 1992 .

[32]  R. Lark,et al.  Model‐based analysis using REML for inference from systematically sampled data on soil , 2004 .

[33]  Cajo J. F. ter Braak,et al.  Treatment of input uncertainty in hydrologic modeling: Doing hydrology backward with Markov chain Monte Carlo simulation , 2008 .

[34]  O. F. Christensen Monte Carlo Maximum Likelihood in Model-Based Geostatistics , 2004 .

[35]  P. Diggle,et al.  Model‐based geostatistics , 2007 .

[36]  Hao Zhang On Estimation and Prediction for Spatial Generalized Linear Mixed Models , 2002, Biometrics.

[37]  R. M. Lark,et al.  The Matérn variogram model: Implications for uncertainty propagation and sampling in geostatistical surveys , 2007 .

[38]  Peter J. Diggle,et al.  geoR and geoRglm: Software for Model-Based Geostatistics , 2003 .

[39]  R. M. Lark,et al.  A comparison of some robust estimators of the variogram for use in soil survey , 2000 .

[40]  J. Vrugt,et al.  A formal likelihood function for parameter and predictive inference of hydrologic models with correlated, heteroscedastic, and non‐Gaussian errors , 2010 .

[41]  A. Papritz,et al.  An Empirical Comparison of Kriging Methods for Nonlinear Spatial Point Prediction , 2002 .

[42]  M. Stein,et al.  A Bayesian analysis of kriging , 1993 .

[43]  K. Mardia,et al.  Maximum likelihood estimation of models for residual covariance in spatial regression , 1984 .