Numerical Comparison of Ensemble Kalman Filter and Randomized Maximum Likelihood

Summary The objective of this paper is to compare the performance of the ensemble Kalman filter (EnKF) to the performance of a gradient-based minimization method for the problem of estimation of facies boundaries in history matching. The EnKF is a Monte Carlo method for data assimilation that uses an ensemble of reservoir models to represent and update the covariance of variables. In several published studies, it outperforms traditional history-matching algorithms in adaptability and efficiency. Because of the approximate nature of the EnKF, the realizations from one ensemble tend to underestimate uncertainty, especially for problems that are highly nonlinear. In this paper, the distributions of reservoir-model realizations from 20 independent ensembles are compared with the distributions from 20 random-ized-maximum-likelihood (RML) realizations for a 2D waterflood model with one injector and four producers. RML is a gradient-based sampling method that generates one reservoir realization in each minimization of the objective function. It is an approximate sampling method, but its sampling properties are similar to the Markov-chain Monte Carlo (McMC) method on highly nonlinear problems and are relatively more efficient than McMC. Despite the nonlinear relationship between the data (such as production rates and facies observations) and the model variables, the EnKF was effective at history matching the production data. We find that the computational effort to generate 20 independent realizations was similar for the two methods, although the complexity of the code is substantially less for the EnKF.

Four non-iterative parameter estimation methods, applying four different ways to assimilate data, are compared to the (iterative) randomized maximum likelihood (RML) method on problems where the forward model is weakly nonlinear and simplistic. Both original and slightly modified versions of established data assimilation methods used for reservoir history matching are considered. Within the Gauss-linear regime, these methods are equivalent in the sense that they all sample correctly from the posterior pdf. The modification consists of replacing the global Kalman gain by a local gain similar to the one applied with RML. This modification, and the simplicity of the forward model considered, allows for asymptotic calculations, correct to lowest order in the quantity controlling the forward-model nonlinearity, to be performed for the non-iterative methods as well as for RML. Comparison of the resulting asymptotic parameter estimates explicitly reveals interesting relations between different ways to assimilate data, that goes beyond the Gauss-linear regime. Results from extensive numerical experiments with the modified methods and with the original methods support the asymptotic results when considering the mean results over a large number of randomized experimental settings where also the forward models vary. For randomized experiments with fixed forward model settings, trends from the asymptotic results were not visible in the numerical results. Although the nonlinear forward models considered are not closely related to reservoir models, the knowledge obtained should be useful also when the latter type of models are considered.

Recently, the ensemble Kalman filter (EnKF) has gained popularity in atmospheric science for the assimilation of data and the assessment of uncertainty in forecasts for complex, large-scale problems. A handful of papers have discussed reservoir characterization applications of the EnKF, which can easily and quickly be coupled with any reservoir simulator. Neither adjoint code nor specific knowledge of simulator numerics is required for implementation of the EnKF. Moreover, data are assimilated (matched) as they become available; a suite of plausible reservoir models (the set of ensembles) is continuously updated to honor data without rematching data assimilated previously. Because of these features, the method is far more efficient for history matching dynamic data than automatic history matching algorithms based on optimization algorithms. Moreover, the suite of ensembles provides a way to evaluate the uncertainty in reservoir description and performance predictions. Here we establish a firm theoretical relation between randomized maximum likelihood and the ensemble Kalman filter. We also consider examples where the performance of the EnKF does not provide a reliable characterization of uncertainty in the model or performance predictions. Although we have previously generated reservoir characterization examples where the method worked well, here we also provide examples where the performance of EnKF does not provide a reliable characterization of uncertainty.

Assimilation of wheat and soil states for improved yield prediction: The APSIM-EnKF framework

This work focuses on the improvement of an integrated methodology for the automatic history matching of compartmentalised reservoirs using 4D seismic results, stochastic initialization and the Ensemble Kalman Filter method. We show the comparison of two different history matching approaches using the Ensemble Kalman Filter (EnKF) to update the Fault Transmissibility Multipliers (FTM) initially estimated with and without considering the 4D seismic results. In this study, the parameters updated during the history matching are two-phase fault transmissibility multipliers (FTM), absolute permeability and effective porosity of a synthetic realistic 3D reservoir. The true impedance map and the changes in reservoir pressure and saturation were previously computed from 4D seismic results. The systematic estimation of two-phase fault transmissibility multipliers is based on the integration of the collected 4D seismic results and an established method validated in our previous work based on a deterministic model, using the gradient-based History Matching, Levenberg Marquardt method (LM). We present the history matching of a synthetic reservoir using the Ensemble Kalman Filter (EnKF) considering the 4D seismic results to update the models and geostatistical techniques to produce the initial geological models. The stochastic method used is the Sequential Gaussian Simulation (SGS) technique to generate 100 initial models. During history matching using the EnKF, the saturation distributions are computed from the forward modelling of a two-phase system (oil-water). The impedance maps are then estimated using the Gassmann equation and compared with the true impedance map as part of the History Matching process. To validate the results, the cost function consisting of two components is calculated, the first is the structural similarity index of the two reconstructed impedance images to the real impedance image and the second is the RMS cost function value, l2 – norm of the difference between the true (real) and the simulated pressure-production data. The EnKF history matching using the two-phase FTM values considering 4D seismic results produced lower cost function values compared with the model using the initial FTM multiplier without considering 4D seismic results. The EnKF history matching algorithm using 4D seismic presented in this work produced results closer to the true reservoir impedance map compared to our previous 4D gradient based history matching method.

Even though the ensemble Kalman filter (EnKF) is widely used, history matching reservoirs with facies description has proven to be a major challenge. A preferred technique for estimating large-scale facies fields within the petroleum industry is still missing. In this paper we present a new approach to this problem. Instead of applying the EnKF directly on facies realizations, the approach applies a transformation of facies fields to a specific level-set function, representing distances between facies types. This ensures better agreement with the EnKF Gaussianity assumptions, and the method always returns facies realizations with geological authenticity. The method also offers large flexibility in generating the initial ensemble, which can be done using any geostatistical tool. Further, no modifications of the standard EnKF equations are needed. The methodology is evaluated on two synthetic examples with increasing complexity. In both examples we consider reservoirs with channel structure. The results presented show that the updated models give large improvements in matching the measurements, and the uncertainty of the models is decreased. Further, recovery of the true petrophysical parameters is highly dependent on sufficient information in the measurements, but in one of the examples considered we are able to completely recover the true channel structure. Additional improvements in the quality of the updated facies fields are obtained by proper handling of the distances close to the reservoir boundaries, and conditioning on specific statistical measures to better preserve prior information about channel properties. Introduction The ensemble Kalman filter (EnKF) (Evensen [1], Aanonsen et al. [2]) is one of the most promising tools for assisted history matching of reservoir models. The filter has in particular proven to work well on cases where petrophysical properties are modeled using Gaussian random fields with variogram models. However, many operative oil and gas reservoirs have challenging complex geological structures - channelized reservoirs being one example. These structures typically have rapid spatial variations in the petrophysical properties (e.g., permeability and porosity). The effect of applying the EnKF directly on such reservoirs is that the characteristics of the fields are lost, i.e., sharp gradients are smeared out (dissipated). Several approaches for using EnKF on complex reservoirs have been proposed. The goal is to ensure that the updated fields are facies realizations, or at least approximations of such, and that good history match is obtained. The methods involve use of truncated pluri-Gaussian methods (Liu and Oliver [3, 4]), Gaussian mixture models (Dovera and Della Rossa [5], Sun et al. [6]), discrete cosine transformations (Jafarpour and McLaughlin [7]), kernel methods (Sarma and Chen [8]), and level-set methods (Moreno et al. [9]). Many of the above mentioned methods give promising results, but there are still drawbacks related to e.g., flexibility in how the initial ensemble may be generated, or the method involves complex and large modifications of the standard EnKF implementation. Our goal is to develop a methodology which is both simple, flexible and able to preserve facies structure for the updated fields. To achieve this, we propose a method motivated by the work of Moreno et al. [9], in the sense that level-set functions are used to represent the facies fields. However, in Moreno et al. [9], only perturbations of a fixed "base case" are updated, and it is questionable whether the real uncertainty is properly described by the prior ensemble. The method we propose is not limited to updating perturbations, and the initial ensemble can be computed using any geostatistical tool. This means that using our methodology, the initial ensemble can be conditioned on the correct facies types at the well locations. The method relies on a transformation of the facies field to a specific level-set function - a signed distance function (representing distances between opposite facies types). The distances are then updated using the EnKF, and converted to petrophysical properties when the reservoir simulator is run to the next assimilation time. The advantage with this transformation is that the ensemble of parameters to be estimated is closer to a Gaussian distribution, which is in better agreement with the EnKF Gaussianity assumptions. In addition, we are guaranteed that the updated fields are facies realizations.

In this study, the relationship between nonlinear model properties and inverse problem solutions is analyzed using a numerical technique based on the inverse problem theory formulated by Mosegaard and Tarantola. According to this theory, the inverse problem and solution are defined via convolution and conjunction of probability density functions (PDFs) that represent stochastic information obtained from the model, observations, and prior knowledge in a joint multidimensional space. This theory provides an explicit analysis of the nonlinear model function, together with information about uncertainties in the model, observations, and prior knowledge through construction of the joint probability density, from which marginal solution functions can then be evaluated. The numerical analysis technique derived from the theory computes the component PDFs in discretized form via a combination of function mapping on a discrete grid in the model and observation phase space and Monte Carlo sampling from known parametric distributions. The efficacy of the numerical analysis technique is demonstrated through its application to two well-known simplified models of atmospheric physics: damped oscillations and Lorenz’s three-component model of dry cellular convection. The major findings of this study include the following: (i) Use of a nonmonotonic forward model in the inverse problem gives rise to the potential for a multimodal posterior PDF, the realization of which depends on the information content of the observations and on observation and model uncertainties. (ii) The cumulative effect of observations over time, space, or both could render the final posterior PDF unimodal, even with the nonmonotonic forward model. (iii) A greater number of independent observations are needed to constrain the solution in the case of a nonmonotonic nonlinear model than for a monotonic nonlinear or linear forward model for a given number of degrees of freedom in control parameter space. (iv) A nonlinear monotonic forward model gives rise to a skewed unimodal posterior PDF, implying a well-posed maximum likelihood inverse problem. (v) The presence of model error greatly increases the possibility of capturing multiple modes in the posterior PDF with the nonmonotonic nonlinear model. (vi) In the case of a nonlinear forward model, use of a Gaussian approximation for the prior update has a similar effect to an increase in model error, which indicates there is the potential to produce a biased mean central estimate even when observations and model are unbiased.

A localized ensemble Kalman filter (EnKF) method is developed to assimilate transient flow data to calibrate a heterogeneous conductivity field. To update conductivity value at a point in a study domain, instead of assimilating all the measurements in the study domain, only limited measurement data in an area around the point are used for the conductivity updating in the localized EnKF method. The localized EnKF is proposed to solve the problems of the filter divergence usually existing in a data assimilation method without localization. The developed method is applied, in a synthetical two dimensional case, to calibrate a heterogeneous conductivity field by assimilating transient hydraulic head data. The simulations by the data assimilation with and without localized EnKF are compared. The study results indicate that the hydraulic conductivity field can be updated efficiently by the localized EnKF, while it cannot be by the EnKF. The covariance inflation and localization are found to solve the problem of the filter divergence efficiently. In comparison with the EnKF method without localization, the localized EnKF method needs smaller ensemble size to achieve stabilized results. The simulation results by the localized EnKF method are much more sensitive to conductivity correlation length than to the localization radius. The developed localized EnKF method provides an approach to improve EnKF method in conductivity calibration.

Time-lapse seismic data has been widely used for the detailed reservoir characterization, and data assimilation algorithms are commonly used for petroleum reservoir history matching of production data. However, hardly any seismic data has been integrated into the reservoir inverse modeling workflow, due to the large data size, finer gridding, and especially the scarcity in time of the seismic data. Popular ensemble-based reservoir inverse modeling methods such as ensemble Kalman filter (EnKF) also face the problems of high computational cost due to the storage of the intermediate variables, restarting the reservoir simulating process and the inconsistency of the full-step and step-wise simulations. The intrinsic sequential data assimilations characteristics of the EnKF set obstacles to assimilate 4D seismic data. The alternative method can be the ensemble smoother (ES), which is an ensemble based method for data assimilation. The ES is based on a Bayesian updating scheme of the reservoir model to match the production history and improve the production forecast. To improve the algorithm convergence, a multiple data assimilation (MDA) method was proposed by Emerick and Reynolds (2012a). The algorithm we used is called ensemble smoother with multiple data assimilation (ESMDA), we modified the ESMDA to integrate the geophysical data, and created the new workflow to history match both the production and geophysical data. The available production data is well measurements, including oil production rate, well water cut and bottom-hole pressure, as well as time-lapse geophysical data, including P-wave impedance. In this paper, we first formulated the mathematical descriptions of ESMDA, we then down-scaled the seismic data, and modified the data assimilation algorithms, illustrated the history matching workflow for both production and geophysical data, finally, we showed how it works with a case study on a water flooding operation in a synthetic reservoir. In comparison, we also showed the history matching only on the production data, which yields an inferior results than the one matches both production and seismic data.

The ensemble Kalman filter (EnKF) has been shown repeatedly to be an effective method for data assimilation in large-scale problems, including those in petroleum engineering. Data assimilation for multiphase flow in porous media is particularly difficult, however, because the relationships between model variables (e.g., permeability and porosity) and observations (e.g., water cut and gas–oil ratio) are highly nonlinear. Because of the linear approximation in the update step and the use of a limited number of realizations in an ensemble, the EnKF has a tendency to systematically underestimate the variance of the model variables. Various approaches have been suggested to reduce the magnitude of this problem, including the application of ensemble filter methods that do not require perturbations to the observed data. On the other hand, iterative least-squares data assimilation methods with perturbations of the observations have been shown to be fairly robust to nonlinearity in the data relationship. In this paper, we present EnKF with perturbed observations as a square root filter in an enlarged state space. By imposing second-order-exact sampling of the observation errors and independence constraints to eliminate the cross-covariance with predicted observation perturbations, we show that it is possible in linear problems to obtain results from EnKF with observation perturbations that are equivalent to ensemble square-root filter results. Results from a standard EnKF, EnKF with second-order-exact sampling of measurement errors that satisfy independence constraints (EnKF (SIC)), and an ensemble square-root filter (ETKF) are compared on various test problems with varying degrees of nonlinearity and dimensions. The first test problem is a simple one-variable quadratic model in which the nonlinearity of the observation operator is varied over a wide range by adjusting the magnitude of the coefficient of the quadratic term. The second problem has increased observation and model dimensions to test the EnKF (SIC) algorithm. The third test problem is a two-dimensional, two-phase reservoir flow problem in which permeability and porosity of every grid cell (5,000 model parameters) are unknown. The EnKF (SIC) and the mean-preserving ETKF (SRF) give similar results when applied to linear problems, and both are better than the standard EnKF. Although the ensemble methods are expected to handle the forecast step well in nonlinear problems, the estimates of the mean and the variance from the analysis step for all variants of ensemble filters are also surprisingly good, with little difference between ensemble methods when applied to nonlinear problems.

Bayesian estimation theory provides a general approach for the state estimate of linear or nonlinear and Gaussian or non-Gaussian systems. In this study, we first explore two Bayesian-based methods: ensemble adjustment Kalman filter (EAKF) and sequential importance resampling particle filter (SIR-PF), using a well-known nonlinear and non-Gaussian model (Lorenz '63 model). The EAKF, which is a deterministic scheme of the ensemble Kalman filter (EnKF), performs better than the classical (stochastic) EnKF in a general framework. Comparison between the SIR-PF and the EAKF reveals that the former outperforms the latter if ensemble size is so large that can avoid the filter degeneracy, and vice versa. The impact of the probability density functions and effective ensemble sizes on assimilation performances are also explored. On the basis of comparisons between the SIR-PF and the EAKF, a mixture filter, called ensemble adjustment Kalman particle filter (EAKPF), is proposed to combine their both merits. Similar to the ensemble Kalman particle filter, which combines the stochastic EnKF and SIR-PF analysis schemes with a tuning parameter, the new mixture filter essentially provides a continuous interpolation between the EAKF and SIR-PF. The same Lorenz '63 model is used as a testbed, showing that the EAKPF is able to overcome filter degeneracy while maintaining the non-Gaussian nature, and performs better than the EAKF given limited ensemble size.

Most geostatistical methods for generating plausible facies models require prior knowledge of the facies proportions distribution to account for heterogeneity within the reservoir. This parameter affects the flow behavior, predictive performance and fluid volume distribution of the model. In practice, it is hardly known with certainty and it is important to account for this uncertainty in the modeling and history matching phases. The ensemble Kalman filter (EnKF), when modified appropriately, to account for non-linearity and non-Gaussianity, is quite robust for automatic history matching problems. However, because the problem of history matching applied to facies is complex, most attempts at solving this problem typically assume constant petrophysical properties and stationary facies proportions. However, realistic representation of facies distribution requires incorporating the variation in petrophysical properties in the geologic model, in addition to any associated uncertainty. This is critical to correctly defining the flow connectivity within the model and the variations in storativity within facies. In this paper, we investigate a number of practical issues in history matching applied to facies with heterogeneous petrophysical properties and non-stationary proportions: (1) Can we account for heterogeneity in facies petrophysical properties during history matching and subsequently reduce any associated uncertainty in them by sequential calibration to data? To the best of our knowledge, this is the first time that heterogeneity in facies petrophysical properties has ever been applied to problems of this kind; (2) What is the EnKF performance when the facies proportions are stationary but implicitly uncertain? (3) Can EnKF efficiently adjust the facies boundaries and distributions when the prior proportions are assumed known but non-stationary? To answer these questions we designed three different test problems. We show that when the facies proportion is non-stationary the EnKF is a robust in estimating the facies properties and distribution.

The ensemble Kalman filter (EnKF) has become a popular method for history matching production and seismic data in petroleum reservoir models. However, it is known that EnKF may fail to give acceptable data matches especially for highly nonlinear problems. In this paper, we introduce a procedure to improve EnKF data matches based on assimilating the same data multiple times with the covariance matrix of the measurement errors multiplied by the number of data assimilations. We prove the equivalence between single and multiple data assimilations for the linear-Gaussian case and present computational evidence that multiple data assimilations can improve EnKF estimates for the nonlinear case. The proposed procedure was tested by assimilating time-lapse seismic data in two synthetic reservoir problems, and the results show significant improvements compared to the standard EnKF. In addition, we review the inversion schemes used in the EnKF analysis and present a rescaling procedure to avoid loss of information during the truncation of small singular values.

Abstract. Over the years data assimilation methods have been developed to obtain estimations of uncertain model parameters by taking into account a few observations of a model state. The most reliable Markov chain Monte Carlo (MCMC) methods are computationally expensive. Sequential ensemble methods such as ensemble Kalman filters and particle filters provide a favorable alternative. However, ensemble Kalman filter has an assumption of Gaussianity. Ensemble transform particle filter does not have this assumption and has proven to be highly beneficial for an initial condition estimation and a small number of parameter estimations in chaotic dynamical systems with non-Gaussian distributions. In this paper we employ ensemble transform particle filter (ETPF) and ensemble transform Kalman filter (ETKF) for parameter estimation in nonlinear problems with 1, 5, and 2500 uncertain parameters and compare them to importance sampling (IS). The large number of uncertain parameters is of particular interest for subsurface reservoir modeling as it allows us to parameterize permeability on the grid. We prove that the updated parameters obtained by ETPF lie within the range of an initial ensemble, which is not the case for ETKF. We examine the performance of ETPF and ETKF in a twin experiment setup, where observations of pressure are synthetically created based on the known values of parameters. For a small number of uncertain parameters (one and five) ETPF performs comparably to ETKF in terms of the mean estimation. For a large number of uncertain parameters (2500) ETKF is robust with respect to the initial ensemble, while ETPF is sensitive due to sampling error. Moreover, for the high-dimensional test problem ETPF gives an increase in the root mean square error after data assimilation is performed. This is resolved by applying distance-based localization, which however deteriorates a posterior estimation of the leading mode by largely increasing the variance due to a combination of less varying localized weights, not keeping the imposed bounds on the modes via the Karhunen–Loeve expansion, and the main variability explained by the leading mode. A possible remedy is instead of applying localization to use only leading modes that are well estimated by ETPF, which demands knowledge of which mode to truncate.

The Ensemble Kalman Filter (EnKF) has gained popularity over recent years as a Monte-Carlo based technique for assisted history matching and real time updating of reservoir models. The EnKF procedure utilizes an ensemble of model states (e.g. realizations of reservoir properties such as porosity and permeability) to approximate the covariance matrices used in the updating process. EnKF works efficiently with Gaussian variables and linear dynamics, but it often fails to preserve the reference probability distribution of the model parameters and to achieve an acceptable production data match where the system dynamics are strongly nonlinear, especially of the type related to multiphase flow, or if non-Gaussian prior models are used. In order to alleviate these drawbacks, we investigated various weighted averaging techniques for computing the ensemble mean by introducing a weighting factor to each ensemble; two new formulations were implemented. The first weighting factor was calculated based on the mismatch in entropy of the model parameters, a normalized measure of the spread of a given probability distribution. The second weighting factor was computed using the forecast mismatch. In addition, both weights could be applied at a single updating step for reducing the forecast mismatch and maintaining the prior distribution simultaneously. The performance of traditional EnKF and these weighted EnKF methods were evaluated by performing various simulation studies with different reservoir heterogeneity. The qualities of the final matching results were assessed by computing the experimental histogram and variograms of the final ensemble, as well as the Root Mean Square Error (RMSE) of the predicted data mismatch. The results reveal that reasonable improvement in the efficiency of the EnKF is achieved by suggested weighted techniques. The RMSE of the predicted data is improved, and the quantity of spurious model parameters is reduced at each updating step. Taking advantage of the entropy based weighting factor assists the filter to preserve the reference distribution. The improvement indicates that the Entropy weighted EnKF (EWEnKF) has a significant potential to resolve the shortfall of traditional EnKF in reservoir characterization and history matching of challenging reservoirs with non-Gaussian distributions.