Eric Laloy

High‐dimensional posterior exploration of hydrologic models using multiple‐try DREAM(ZS) and high‐performance computing

[1]xa0Spatially distributed hydrologic models are increasingly being used to study and predict soil moisture flow, groundwater recharge, surface runoff, and river discharge. The usefulness and applicability of such complex models is increasingly held back by the potentially many hundreds (thousands) of parameters that require calibration against some historical record of data. The current generation of search and optimization algorithms is typically not powerful enough to deal with a very large number of variables and summarize parameter and predictive uncertainty. We have previously presented a general-purpose Markov chain Monte Carlo (MCMC) algorithm for Bayesian inference of the posterior probability density function of hydrologic model parameters. This method, entitled differential evolution adaptive Metropolis (DREAM), runs multiple different Markov chains in parallel and uses a discrete proposal distribution to evolve the sampler to the posterior distribution. The DREAM approach maintains detailed balance and shows excellent performance on complex, multimodal search problems. Here we present our latest algorithmic developments and introduce MT-DREAM(ZS), which combines the strengths of multiple-try sampling, snooker updating, and sampling from an archive of past states. This new code is especially designed to solve high-dimensional search problems and receives particularly spectacular performance improvement over other adaptive MCMC approaches when using distributed computing. Four different case studies with increasing dimensionality up to 241 parameters are used to illustrate the advantages of MT-DREAM(ZS).

Mass conservative three-dimensional water tracer distribution from Markov chain Monte Carlo inversion of time-lapse ground-penetrating radar data

Received 3 August 2011; revised 15 May 2012; accepted 4 June 2012; published 13 July 2012. [1] Time-lapse geophysical measurements are widely used to monitor the movement of water and solutes through the subsurface. Yet commonly used deterministic least squares inversions typically suffer from relatively poor mass recovery, spread overestimation, and limited ability to appropriately estimate nonlinear model uncertainty. We describe herein a novel inversion methodology designed to reconstruct the three-dimensional distribution of a tracer anomaly from geophysical data and provide consistent uncertainty estimates using Markov chain Monte Carlo simulation. Posterior sampling is made tractable by using a lower-dimensional model space related both to the Legendre moments of the plume and to predefined morphological constraints. Benchmark results using cross-hole groundpenetrating radar travel times measurements during two synthetic water tracer application experiments involving increasingly complex plume geometries show that the proposed method not only conserves mass but also provides better estimates of plume morphology and posterior model uncertainty than deterministic inversion results.

Probabilistic inference of multi-Gaussian fields from indirect hydrological data using circulant embedding and dimensionality reduction

© 2015. American Geophysical Union. All Rights Reserved. We present a Bayesian inversion method for the joint inference of high-dimensional multi-Gaussian hydraulic conductivity fields and associated geostatistical parameters from indirect hydrological data. We combine Gaussian process generation via circulant embedding to decouple the variogram from grid cell specific values, with dimensionality reduction by interpolation to enable Markov chain Monte Carlo (MCMC) simulation. Using the Matern variogram model, this formulation allows inferring the conductivity values simultaneously with the field smoothness (also called Matern shape parameter) and other geostatistical parameters such as the mean, sill, integral scales and anisotropy direction(s) and ratio(s). The proposed dimensionality reduction method systematically honors the underlying variogram and is demonstrated to achieve better performance than the Karhunen-Loeve expansion. We illustrate our inversion approach using synthetic (error corrupted) data from a tracer experiment in a fairly heterogeneous 10,000-dimensional 2-D conductivity field. A 40-times reduction of the size of the parameter space did not prevent the posterior simulations to appropriately fit the measurement data and the posterior parameter distributions to include the true geostatistical parameter values. Overall, the posterior field realizations covered a wide range of geostatistical models, questioning the common practice of assuming a fixed variogram prior to inference of the hydraulic conductivity values. Our method is shown to be more efficient than sequential Gibbs sampling (SGS) for the considered case study, particularly when implemented on a distributed computing cluster. It is also found to outperform the method of anchored distributions (MAD) for the same computational budget. Key Points: Joint Bayesian inference of Gaussian conductivity fields and their variograms A dimensionality reduction that systematically honors the underlying variogram Distributed multiprocessor implementation is straightforward

Inversion using a new low-dimensional representation of complex binary geological media based on a deep neural network

Abstract Efficient and high-fidelity prior sampling and inversion for complex geological media is still a largely unsolved challenge. Here, we use a deep neural network of the variational autoencoder type to construct a parametric low-dimensional base model parameterization of complex binary geological media. For inversion purposes, it has the attractive feature that random draws from an uncorrelated standard normal distribution yield model realizations with spatial characteristics that are in agreement with the training set. In comparison with the most commonly used parametric representations in probabilistic inversion, we find that our dimensionality reduction (DR) approach outperforms principle component analysis (PCA), optimization-PCA (OPCA) and discrete cosine transform (DCT) DR techniques for unconditional geostatistical simulation of a channelized prior model. For the considered examples, important compression ratios (200–500) are achieved. Given that the construction of our parameterization requires a training set of several tens of thousands of prior model realizations, our DR approach is more suited for probabilistic (or deterministic) inversion than for unconditional (or point-conditioned) geostatistical simulation. Probabilistic inversions of 2D steady-state and 3D transient hydraulic tomography data are used to demonstrate the DR-based inversion. For the 2D case study, the performance is superior compared to current state-of-the-art multiple-point statistics inversion by sequential geostatistical resampling (SGR). Inversion results for the 3D application are also encouraging.

Reply to comment by Chu et al. on “High‐dimensional posterior exploration of hydrologic models using multiple‐try DREAM (ZS) and high‐performance computing”

PUBLICATIONS Water Resources Research COMMENTARY 10.1002/2013WR014425 This article is a reply to Chu et al. [2014] doi:10.1002/2012WR013341. Reply to comment by Chu et al. on ‘‘High-dimensional posterior exploration of hydrologic models using multiple-try DREAM (ZS) and high-performance computing’’ Jasper A. Vrugt 1,2,3 and Eric Laloy 4 Correspondence to: J. A. Vrugt, jasper@uci.edu Citation: Vrugt, J. A., and E. Laloy (2014), Reply to comment Chu et al. on ‘‘High- dimensional posterior exploration of hydrologic models using multiple-try DREAM (ZS) and high-performance computing,’’ Water Resour. Res., 50, 2781–2786, doi:10.1002/ 2013WR014425. Received 15 JULY 2013 Accepted 8 FEB 2014 Accepted article online 15 FEB 2014 Published online 21 MAR 2014 Department of Civil and Environmental Engineering, University of California, Irvine, California, USA, 2 Department of Earth System Science, University of California, Irvine, California, USA, 3 Institute for Biodiversity and Ecosystems Dynamics, University of Amsterdam, Amsterdam, Netherlands, 4 Institute for Environment, Health and Safety, Belgian Nuclear Research Centre, Mol, Belgium 1. Introduction The comment of Chu et al. [2014], hereafter referred to as CYG, raises questions about some of the results presented in our paper (Laloy and Vrugt [2012], hereinafter LV, which is not to be confused with Las Vegas, although appropriate concerning the subject of this work on Monte Carlo simulation). As a preamble, we would like to remark that the work presented in Chu et al. [2010] and LV (2012) concern two different fields of study. CYG view the model calibration as an optimization problem, and use common concepts to efficiently find a single realization of the parameter values that minimizes (or maximizes, if appropriate) some user-defined objective function. Our work, on the contrary, uses Bayesian principles coupled with MCMC simulation to derive a distribution of parameter values that honor the observed data. This distribution summarizes parameter and model predictive (simulation) uncertainty, a requirement for probabilistic analysis, operational forecasting, disentangling error sources, and decision making. The maxi- mum a posteriori density (MAP) parameter values derived with MCMC simulation should reside in close vicinity of the ‘‘best’’ solution found with an optimization algorithm, if the exact same data set, prior distribu- tion, and likelihood (objective) function are used. In this reply, we assume that CYG used a correct imple- mentation of the MT-DREAM (ZS) algorithm and similar data set, prior and likelihood function as LV. Otherwise, the comparative analysis is meaningless. We emphasize this for three reasons. First, the results presented herein contradict CYG and are similar to those reported in LV but now with more trials plotted. Second, contributions in physics [Horowitz et al., 2012; Toyli et al., 2012; Yale et al., 2013] and geophysics [Linde and Vrugt, 2013; Laloy et al., 2012, Rosas-Carbajal et al., 2014; T. Lochbuehler et al., Summary statistics from training images as model constraints in probabilistic inversion, Geophysical Journal International, in review, 2014] demonstrate proper convergence behavior of MT-DREAM (ZS) on complex and high- dimensional targets involving hundreds of parameters. Third, to justify their SP-UCI algorithm the original paper of Chu et al. [2010] portrays misleading results of the predecessor of DREAM, called SCEM-UA. Section 4 of this reply will address this latter issue in more detail. We now respond to the comments of CYG. We use different sections with numbering corresponding to CYG. 2. Computational Time Unit CYG find the Computational Time Unit (CTU) diagnostic to be a poor indicator of the performance of MT- DREAM (ZS) . They suggest using the number of function evaluations or clock time instead. The CPU-time (s) scales linearly with CTU, or CPU 5 aCTU, where a (s) denotes the average time it takes to complete a sin- gle function (model) evaluation. As a is dependent on the processor speed (hardware), LV purposely reported the CTU values. Note that we neglect the actual run time of MT-DREAM (ZS) , in the determination of a, which is appropriate given the intended application of this algorithm to CPU-intensive forward models. We purposely do not use the number of function evaluations as performance diagnostic. This metric does not properly convey the CPU-time (CTU) of parallel algorithms such as DREAM (ZS) or MT-DREAM (ZS) . These VRUGT AND LALOY C 2014. American Geophysical Union. All Rights Reserved. V

Training‐Image Based Geostatistical Inversion Using a Spatial Generative Adversarial Neural Network

Probabilistic inversion within a multiple-point statistics framework is still computationally prohibitive for large-scale problems. To partly address this, we introduce and evaluate a new training-image based simulation and inversion approach for complex geologic media. Our approach relies on a deep neural network of the spatial generative adversarial network (SGAN) type. After training using a training image (TI), our proposed SGAN can quickly generate 2D and 3D unconditional realizations. A key feature of our SGAN is that it defines a (very) low-dimensional parameterization, thereby allowing for efficient probabilistic (or deterministic) inversion using state-of-the-art Markov chain Monte Carlo (MCMC) methods. A series of 2D and 3D categorical TIs is first used to analyze the performance of our SGAN for unconditional simulation. The speed at which realizations are generated makes it especially useful for simulating over large grids and/or from a complex multi-categorical TI. Subsequently, synthetic inversion case studies involving 2D steady-state flow and 3D transient hydraulic tomography are used to illustrate the effectiveness of our proposed SGAN-based probabilistic inversion. For the 2D case, the inversion rapidly explores the posterior model distribution. For the 3D case, the inversion recovers model realizations that fit the data close to the target level and visually resemble the true model well. Future work will focus on the inclusion of direct conditioning data and application to continuous TIs.