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

Eric Laloy; Romain Hérault; John Lee; Diederik Jacques; Niklas Linde

Source-linked AI summary

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

Eric Laloy, Romain Hérault, John Lee, Diederik Jacques, Niklas Linde

arXiv:1710.09196v1Oct 25, 2017stat.MLphysics.geo-ph

TL;DR

Complex geological inversion requires efficient prior sampling that preserves high-order spatial structure. The paper uses a deep variational autoencoder to construct a low-dimensional representation for binary geological media, then applies it in probabilistic inversion. The representation outperforms PCA, OPCA, and DCT for the tested channelized prior, while achieving 200–500 compression ratios and better 2D results than SGR.

Problem
Efficient, high-fidelity prior sampling and inversion for complex geological media remains a largely unsolved challenge.
Method
A deep variational autoencoder constructs a low-dimensional parametric representation whose decoder generates training-image-consistent binary geological realizations for Bayesian inversion.
Results
The representation outperforms PCA, OPCA, and DCT for unconditional simulation of a channelized prior and gives better 2D inversion results than SGR, with encouraging 3D results.
Takeaways & Limitations
The approach supports efficient probabilistic inversion of complex binary geological media using a low-dimensional model space.
Takeaways & Limitations
Constructing the parameterization requires several tens of thousands of training models, and the global representation can hinder convergence for categorical-field inversions.

Abstract

from arXiv · show

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.

1 Introduction

The paper addresses the difficulty of representing complex, non-Gaussian geological priors for inversion by proposing a parametric low-dimensional approach based on a deep variational autoencoder. The approach is intended to generate TI-consistent realizations and enable global parameter-based inverse methods.

Inverse modeling seeks subsurface property fields whose simulated responses agree with measured data while accounting for prior information and uncertainty.The desired outcome is an ensemble spanning geologic models that honor both inversion data and the prior model.
Training images provide higher-order spatial information for complex or non-multi-Gaussian geology, but MPS lacks an explicit parameterization for classical inversion.This limitation prevents direct use of gradient-based parameter-based methods such as randomized maximum likelihood.
The paper proposes a parametric dimensionality-reduction approach that transforms a lower-dimensional space into realizations consistent with complex binary prior models.This representation is designed to make global parameter-based inverse methods applicable to models informed by training images.
Existing compressed inversion representations include PCA, kernel-PCA, level-set, wavelet, DCT, SVD, and K-SVD methods.Optimization-PCA was identified as an approach designed to build a lower-dimensional representation consistent with a training image.
A convolutional deep variational autoencoder generates TI-consistent geological realizations from a low-dimensional standard Gaussian noise vector and is used to sample posterior binary 2D and 3D fields.The convolutional architecture is motivated by the suitability of adaptive convolutional filters for image processing.
The study evaluates the representation for unconditional and conditionally informed geostatistical simulation before applying it to inverse problems.This evaluation checks whether randomly sampled realizations reproduce the spatial statistics of the selected training image.

2.1 Deep neural network architecture

The architecture uses autoencoders to learn nonlinear low-dimensional parameterizations, with variational sampling and convolutional layers enabling generation of spatially structured geological realizations.

Generalities: A neuron combines weighted inputs and a bias through a nonlinear activation function, with weights and biases learned for the task.Differentiable activation functions allow learning by gradient descent.
Autoencoders: An autoencoder compresses higher-dimensional input data into a central low-dimensional code and reconstructs the input through encoder and decoder parts.The decoder's output has the same dimensionality as the input.
Variational autoencoders: A variational autoencoder represents the code with mean and standard-deviation vectors, rescales sampled standard-normal variables, and decodes them into reconstructed outputs.Training minimizes reconstruction differences while encouraging the latent vector to follow a multivariate standard normal distribution.
Variational autoencoders: After training, sampling a d-dimensional standard-normal vector and passing it through the decoder generates new model realizations.This combines a standard-normal low-dimensional parameterization with a mapping to complex MPS geostatistical realizations.
Convolutional layers: Convolutional layers use small filters over local receptive fields, preserving spatial structure while sharing parameters across feature-map elements.Stacked layers encode increasingly larger patterns and require fewer parameters than fully connected layers.

2.2 Training the deep variational autoencoder

The variational autoencoder uses probabilistic encoding and decoding to learn a low-dimensional representation of binary images, trained with a differentiable objective and a standard-normal latent prior. Training combines reconstruction and regularization terms, with substantial computational requirements.

For a 100×100 binary image generated from a 50-dimensional code, the network contains 34,039,237 optimized parameters.
The probabilistic decoder maps latent codes to distributions over binary images, while the probabilistic encoder maps images to distributions over latent codes.
The latent code prior is an uncorrelated standard normal distribution, and training encourages encoded distributions to remain close to it through Kullback-Leibler divergence.
For binary images, the reconstruction objective can be expressed as binary cross-entropy minimization.
Training used a GPU for 100 epochs and required 5–13 hours across model domains, training-set sizes, and network configurations.

2.3 Assessing geostatistical simulation quality

Geostatistical simulation quality is assessed by comparing reproduced spatial patterns and variability across realizations. The evaluation uses connectivity functions for pattern closeness and Jensen-Shannon distances between multiple-point histograms for uncertainty-space characterization.

Simulation quality is evaluated through pattern reproduction and between-realization variability, called the space of uncertainty.
The two-point cluster function measures the probability of a continuous same-facies path between points separated by a given lag distance.
The 2D space of uncertainty is estimated from average Jensen-Shannon distances between multiple-point histograms across 100 realizations.
Each multiple-point histogram vector contains 216 = 65,536 components when computed with a 4 × 4 template.
The space of uncertainty is defined as the average distance dJS across either training images or dimensionality-reduction realizations.

2.4 Bayesian inversion

The Bayesian inversion framework represents measurements as a deterministic forward-model response plus error and infers parameters through an unnormalized posterior density. Because the nonlinear posterior lacks an exact analytical solution, DREAM(ZS) MCMC is used for approximation.

The forward problem represents measurement data as d = F(θ) + e, where F(θ) is deterministic and e aggregates all error sources.
Bayesian inference treats parameters θ as random variables with posterior density p(θ|d).
The likelihood L(θ|d) is defined as p(d|θ), linking parameter values to the observed data.
For fixed parameter dimensionality, inference uses the unnormalized posterior p(θ|d) ∝ p(θ)L(θ|d), avoiding computation of the normalization factor.
With normally distributed, uncorrelated errors of known constant variance, the log-likelihood is used to avoid numerical over- or underflow.
DREAM(ZS) MCMC approximates the posterior because no exact analytical solution is available for the considered nonlinear inverse problems.

3 Geostatistical simulation results

The deep-neural-network dimensionality-reduction approach generates unconditional and conditioned binary channelized aquifer realizations that broadly reproduce training-set spatial statistics, while outperforming PCA, OPCA, and DCT generation in the unconditional case. It also achieves high compression in 3D, though conditioning accuracy and artifacts involve trade-offs.

2D unconditional simulation: The DR-based unconditional 2D realizations resemble the training models and have similar cluster or connectivity-function curves.Their dJS ratio is 1.03, indicating slightly more between-realization variability than the training set.
2D unconditional simulation: PCA, OPCA, and DCT fail to sample a proper low-dimensional parameter space, whereas the proposed DR approach generates substantially better channelized realizations.PCA and OPCA used 70 random variables; DCT used 250 coefficients, and varying these settings did not improve generation quality.
Computational setup: Constructing conditioned training sets required substantial computation, including 40,000 2D models and 19,000 3D models generated with MPS procedures.The 2D conditioned training set required 7–8 hours, while the 3D conditioned set required approximately 14 hours.
2D conditional simulation: Conditioned DR realizations show similar visual generation performance and connectivity-function curves to the conditioned training images.The conditioned example uses nine known facies locations.
2D conditional simulation: 68% of conditioned DR realizations honor all nine prescribed facies, while 97% contain at most one mismatching conditioning datum.Without sequential relooping, these figures increase to 92% and 100%, respectively, but broken channels and isolated pixels occur more frequently.
3D unconditional simulation: The 3D DR approach produces realizations similar to the training set, reproduces its connectivity metrics well, and reaches a compression ratio of approximately 518.The 3D domain contains 25,920 voxels represented with 50 latent dimensions, although isolated voxels of each facies are over-represented.

4 Inverse problems

The DR-based inversion is evaluated on 2D steady-state and 3D transient hydraulic problems using DREAM(ZS), with comparisons against SGR in the 2D case. The approach retrieves plausible posterior models, but convergence remains incomplete within the tested iteration budgets.

2D steady-state flow: The first 2D test used a 100 × 100 channelized aquifer, 49 head observations, and Gaussian noise with standard deviation 0.02 m.The reference model was generated using the proposed DR approach, with 50 DR variables and eight DREAM(ZS) chains.
2D steady-state flow: After 120,000 serial iterations, sampled models closely resembled the true model, with fPO = 0.86 compared with the prior value fPR = 0.60.The resulting fPO/fPR ratio was 1.43.
Convergence: Formal convergence was incomplete: the Gelman–Rubin criterion was satisfied for 37 of 50 parameters in one test and 29 of 50 dimensions in another.The authors therefore state that longer chains are needed for full posterior exploration.
2D steady-state flow: In the second 2D test, using an MPS-generated true model, posterior realizations remained visually close to the reference, with fPO = 0.78.Posterior variability was somewhat larger than in the first test.
3D transient hydraulic tomography: In the 3D transient case, 16 DREAM(ZS) chains reached data misfits of 0.0202 m–0.0204 m after 60,000 iterations per chain, near the 0.0200 m target.The posterior mode had not yet been sampled within the computational budget.

5 Discussion

The discussion identifies useful inversion capabilities alongside limitations of the global DR representation, including non-local parameter effects, training-data demands, and restricted field types. It also notes potential extensions to gradient-based inversion and prediction-focused modeling.

Performance: The DR approach performed well for the considered 2D steady-state inversion and provided useful solutions for the 3D transient problem, but further improvements are needed.This conclusion is limited to the considered inversion applications.
Global representation: Changing one low-dimensional parameter influences the entire generated realization, and small realization changes can cause large RMSE increases.The discussion identifies this non-locality as a problem for categorical-field inversion, especially with rich, low-error datasets.
Limitations: The parameterization is not totally accurate for direct conditioning and requires several tens of thousands of MPS-generated training models.This may be impractical when MPS generation is computationally intensive.
Limitations: The approach currently produces continuous realizations requiring thresholding for categorical fields and is not yet extended to multicategorical or continuous training images.The authors also identify 3D convolutional layers and alternative reconstruction losses as possible improvements.
Potential extensions: Without hard thresholding, automatic differentiation provides the Jacobian dX/dz, which may support gradient-based inversion with an adjoint forward solver.The derivatives are available because the generator is differentiable.
Potential extensions: Deep generative compression may also support prediction-focused modeling and could potentially apply beyond MPS to priors with discrete structure.Prediction-focused modeling links past data variables with prediction data variables without inverting for material properties.

6 Conclusion

The paper presents a VAE-based low-dimensional representation for efficient sampling and probabilistic inversion of complex binary geological media. It outperforms PCA, OPCA, DCT, and SGR in the reported comparisons, while requiring substantial training data and further development for broader applicability.

Contribution: A VAE-based low-dimensional manifold generates complex binary geologic realizations consistent with the training set and supports fast prior sampling for probabilistic inversion.The representation uses a low-dimensional parametric base model.
Dimensionality reduction: The proposed representation outperformed PCA, OPCA, and DCT for unconditional simulation of a channelized prior model.The comparison concerns parametric dimensionality-reduction methods.
Scope: The approach was designed for probabilistic or deterministic inversion rather than unconditional or point-conditioned geostatistical simulation because construction requires several tens of thousands of training realizations.The training realizations are obtained using multiple-point statistics simulation.
Inversion results: Synthetic 2D steady-flow and 3D transient hydraulic-tomography studies demonstrated the DR-based probabilistic inversion approach.The 2D case yielded better results than SGR, while the 3D results were described as encouraging.
Future work: Future work targets reducing the required training-model set, extending the method to multicategorical and continuous variables, and developing a local rather than global representation.The authors expect localization to improve inversion performance.