Adaptive Markov Chain Monte Carlo Research Articles

Bayesian adaptive Markov chain Monte Carlo estimation of genetic parameters

Accurate and fast estimation of genetic parameters that underlie quantitative traits using mixed linear models with additive and dominance effects is of great importance in both natural and breeding populations. Here, we propose a new fast adaptive Markov chain Monte Carlo (MCMC) sampling algorithm for the estimation of genetic parameters in the linear mixed model with several random effects. In the learning phase of our algorithm, we use the hybrid Gibbs sampler to learn the covariance structure of the variance components. In the second phase of the algorithm, we use this covariance structure to formulate an effective proposal distribution for a Metropolis-Hastings algorithm, which uses a likelihood function in which the random effects have been integrated out. Compared with the hybrid Gibbs sampler, the new algorithm had better mixing properties and was approximately twice as fast to run. Our new algorithm was able to detect different modes in the posterior distribution. In addition, the posterior mode estimates from the adaptive MCMC method were close to the REML (residual maximum likelihood) estimates. Moreover, our exponential prior for inverse variance components was vague and enabled the estimated mode of the posterior variance to be practically zero, which was in agreement with the support from the likelihood (in the case of no dominance). The method performance is illustrated using simulated data sets with replicates and field data in barley.

Uncertainty in dual permeability model parameters for structured soils

Successful application of dual permeability models (DPM) to predict contaminant transport is contingent upon measured or inversely estimated soil hydraulic and solute transport parameters. The difficulty in unique identification of parameters for the additional macropore- and matrix-macropore interface regions, and knowledge about requisite experimental data for DPM has not been resolved to date. Therefore, this study quantifies uncertainty in dual permeability model parameters of experimental soil columns with different macropore distributions (single macropore, and low- and high-density multiple macropores). Uncertainty evaluation is conducted using adaptive Markov chain Monte Carlo (AMCMC) and conventional Metropolis-Hastings (MH) algorithms while assuming 10 out of 17 parameters to be uncertain or random. Results indicate that AMCMC resolves parameter correlations and exhibits fast convergence for all DPM parameters while MH displays large posterior correlations for various parameters. This study demonstrates that the choice of parameter sampling algorithms is paramount in obtaining unique DPM parameters when information on covariance structure is lacking, or else additional information on parameter correlations must be supplied to resolve the problem of equifinality of DPM parameters. This study also highlights the placement and significance of matrix-macropore interface in flow experiments of soil columns with different macropore densities. Histograms for certain soil hydraulic parameters display tri-modal characteristics implying that macropores are drained first followed by the interface region and then by pores of the matrix domain in drainage experiments. Results indicate that hydraulic properties and behavior of the matrix-macropore interface is not only a function of saturated hydraulic conductivity of the macroporematrix interface (Ksa ) and macropore tortuosity (lf ) but also of other parameters of the matrix and macropore domains.

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

Spatially 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).

Bayesian calibration of a large‐scale geothermal reservoir model by a new adaptive delayed acceptance Metropolis Hastings algorithm

The aim of this research is to estimate the parameters of a large‐scale numerical model of a geothermal reservoir using Markov chain Monte Carlo (MCMC) sampling, within the framework of Bayesian inference. All feasible parameters that are consistent with the measured data are summarized by the posterior distribution, and hence parameter estimation and uncertainty quantification are both given by calculating expected values of statistics of interest over the posterior distribution. It appears to be computationally infeasible to use the standard Metropolis‐Hastings algorithm (MH) to sample the high dimensional computationally expensive posterior distribution. To improve the sampling efficiency, a new adaptive delayed‐acceptance MH algorithm (ADAMH) is implemented to adaptively build a stochastic model of the error introduced by the use of a reduced‐order model. This use of adaptivity differs from existing adaptive MCMC algorithms that tune proposal distributions of the Metropolis‐Hastings algorithm (MH), though ADAMH also implements that technique. For the 3‐D geothermal reservoir model we present here, ADAMH shows a great improvement in the computational efficiency of the MCMC sampling, and promising results for parameter estimation and uncertainty quantification are obtained. This algorithm could offer significant improvement in computational efficiency when implementing sample‐based inference in other large‐scale inverse problems.

Bayesian calibration and uncertainty analysis of hydrological models: A comparison of adaptive Metropolis and sequential Monte Carlo samplers

Bayesian statistical inference implemented by stochastic algorithms such as Markov chain Monte Carlo (MCMC) provides a flexible probabilistic framework for model calibration that accounts for both model and parameter uncertainties. The effectiveness of such Monte Carlo algorithms depends strongly on the user‐specified proposal or sampling distribution. In this article, a sequential Monte Carlo (SMC) approach is used to obtain posterior parameter estimates of a conceptual hydrologic model using data from selected catchments in eastern Australia. The results are evaluated against the popular adaptive Metropolis MCMC sampling approach. Both methods display robustness and convergence, but the SMC displays greater efficiency in exploring the parameter space in catchments where the optimal solutions lie in the tails of the prescribed prior distribution. The SMC method is also able to identify a different set of parameters with an overall improvement in likelihood and Nash‐Sutcliffe efficiency for selected catchments. As a result of its population‐based sampling mechanism, the SMC method is shown to offer improved efficiency in identifying parameter optimization and to provide sampling robustness, in particular in identifying global posterior modes.

Estimation of ECHAM5 climate model closure parameters with adaptive MCMC

Abstract. Climate models contain closure parameters to which the model climate is sensitive. These parameters appear in physical parameterization schemes where some unresolved variables are expressed by predefined parameters rather than being explicitly modeled. Currently, best expert knowledge is used to define the optimal closure parameter values, based on observations, process studies, large eddy simulations, etc. Here, parameter estimation, based on the adaptive Markov chain Monte Carlo (MCMC) method, is applied for estimation of joint posterior probability density of a small number (n=4) of closure parameters appearing in the ECHAM5 climate model. The parameters considered are related to clouds and precipitation and they are sampled by an adaptive random walk process of the MCMC. The parameter probability densities are estimated simultaneously for all parameters, subject to an objective function. Five alternative formulations of the objective function are tested, all related to the net radiative flux at the top of the atmosphere. Conclusions of the closure parameter estimation tests with a low-resolution ECHAM5 climate model indicate that (i) adaptive MCMC is a viable option for parameter estimation in large-scale computational models, and (ii) choice of the objective function is crucial for the identifiability of the parameter distributions.

Model selection and adaptive Markov chain Monte Carlo for Bayesian cointegrated {VAR} models

This paper develops a matrix-variate adaptive Markov chain Monte Carlo (MCMC) methodology for Bayesian Cointegrated Vector Auto Regressions (CVAR). We replace the popular approach to sampling Bayesian CVAR models, involving griddy Gibbs, with an automated efficient alternative, based on the Adaptive Metropolis algorithm of Roberts and Rosenthal (2009). Developing the adaptive MCMC framework for Bayesian CVAR models allows for efficient estimation of posterior parameters in significantly higher dimensional CVAR series than previously possible with existing griddy Gibbs samplers. For a n-dimensional CVAR series, the matrix-variate posterior is in dimension $3n^2 + n$, with significant correlation present between the blocks of matrix random variables. Hence, utilizing a griddy Gibbs sampler for large n becomes computationally impractical as it involves approximating an $n \times n$ full conditional posterior using a spline over a high dimensional $n \times n$ grid. The adaptive MCMC approach is demonstrated to be ideally suited to learning on-line a proposal to reflect the posterior correlation structure, therefore improving the computational efficiency of the sampler. We also treat the rank of the CVAR model as a random variable and perform joint inference on the rank and model parameters. This is achieved with a Bayesian posterior distribution defined over both the rank and the CVAR model parameters, and inference is made via Bayes Factor analysis of rank. Practically the adaptive sampler also aids in the development of automated Bayesian cointegration models for algorithmic trading systems considering instruments made up of several assets, such as currency baskets. Previously the literature on financial applications of CVAR trading models typically only considers pairs trading (n=2) due to the computational cost of the griddy Gibbs. We are able to extend under our adaptive framework to $n >> 2$ and demonstrate an example with n = 10, resulting in a posterior distribution with parameters up to dimension 310. By also considering the rank as a random quantity we can ensure our resulting trading models are able to adjust to potentially time varying market conditions in a coherent statistical framework.

Model selection and adaptive Markov chain Monte Carlo for Bayesian cointegrated {VAR} models

This paper develops a matrix-variate adaptive Markov chain Monte Carlo (MCMC) methodology for Bayesian Cointegrated Vector Auto Regressions (CVAR). We replace the popular approach to sampling Bayesian CVAR models, involving griddy Gibbs, with an automated efficient alternative, based on the Adaptive Metropolis algorithm of Roberts and Rosenthal (2009). Developing the adaptive MCMC framework for Bayesian CVAR models allows for efficient estimation of posterior parameters in significantly higher dimensional CVAR series than previously possible with existing griddy Gibbs samplers. For a n-dimensional CVAR series, the matrix-variate posterior is in dimension $3n^2 + n$, with significant correlation present between the blocks of matrix random variables. Hence, utilizing a griddy Gibbs sampler for large n becomes computationally impractical as it involves approximating an $n \times n$ full conditional posterior using a spline over a high dimensional $n \times n$ grid. The adaptive MCMC approach is demonstrated to be ideally suited to learning on-line a proposal to reflect the posterior correlation structure, therefore improving the computational efficiency of the sampler. We also treat the rank of the CVAR model as a random variable and perform joint inference on the rank and model parameters. This is achieved with a Bayesian posterior distribution defined over both the rank and the CVAR model parameters, and inference is made via Bayes Factor analysis of rank. Practically the adaptive sampler also aids in the development of automated Bayesian cointegration models for algorithmic trading systems considering instruments made up of several assets, such as currency baskets. Previously the literature on financial applications of CVAR trading models typically only considers pairs trading (n=2) due to the computational cost of the griddy Gibbs. We are able to extend under our adaptive framework to $n >> 2$ and demonstrate an example with n = 10, resulting in a posterior distribution with parameters up to dimension 310. By also considering the rank as a random quantity we can ensure our resulting trading models are able to adjust to potentially time varying market conditions in a coherent statistical framework.

A computational framework for empirical Bayes inference

In empirical Bayes inference one is typically interested in sampling from the posterior distribution of a parameter with a hyper-parameter set to its maximum likelihood estimate. This is often problematic particularly when the likelihood function of the hyper-parameter is not available in closed form and the posterior distribution is intractable. Previous works have dealt with this problem using a multi-step approach based on the EM algorithm and Markov Chain Monte Carlo (MCMC). We propose a framework based on recent developments in adaptive MCMC, where this problem is addressed more efficiently using a single Monte Carlo run. We discuss the convergence of the algorithm and its connection with the EM algorithm. We apply our algorithm to the Bayesian Lasso of Park and Casella (J. Am. Stat. Assoc. 103:681---686, 2008) and on the empirical Bayes variable selection of George and Foster (J. Am. Stat. Assoc. 87:731---747, 2000).

Divide and Conquer: A Mixture-Based Approach to Regional Adaptation for MCMC

The efficiency of Markov chain Monte Carlo (MCMC) algorithms can vary dramatically with the choice of simulation parameters. Adaptive MCMC (AMCMC) algorithms allow the automatic tuning of the parameters while the simulation is in progress. A multimodal target distribution may call for regional adaptation of Metropolis–Hastings samplers so that the proposal distribution varies across regions in the sample space. Establishing such a partition is not straightforward and, in many instances, the learning required for its specification takes place gradually, as the simulation proceeds. In the case in which the target distribution is approximated by a mixture of Gaussians, we propose an adaptation process for the partition. It involves fitting the mixture using the available samples via an online EM algorithm and, based on the current mixture parameters, constructing the regional adaptive algorithm with online recursion (RAPTOR). The method is compared with other regional AMCMC samplers and is tested on simulated as well as real data examples.Relevant theoretical proofs, code and datasets are posted as an online supplement.

Aerosol model selection and uncertainty modelling by adaptive MCMC technique

Abstract. We present a new technique for model selection problem in atmospheric remote sensing. The technique is based on Monte Carlo sampling and it allows model selection, calculation of model posterior probabilities and model averaging in Bayesian way. The algorithm developed here is called Adaptive Automatic Reversible Jump Markov chain Monte Carlo method (AARJ). It uses Markov chain Monte Carlo (MCMC) technique and its extension called Reversible Jump MCMC. Both of these techniques have been used extensively in statistical parameter estimation problems in wide area of applications since late 1990's. The novel feature in our algorithm is the fact that it is fully automatic and easy to use. We show how the AARJ algorithm can be implemented and used for model selection and averaging, and to directly incorporate the model uncertainty. We demonstrate the technique by applying it to the statistical inversion problem of gas profile retrieval of GOMOS instrument on board the ENVISAT satellite. Four simple models are used simultaneously to describe the dependence of the aerosol cross-sections on wavelength. During the AARJ estimation all the models are used and we obtain a probability distribution characterizing how probable each model is. By using model averaging, the uncertainty related to selecting the aerosol model can be taken into account in assessing the uncertainty of the estimates.

Adaptive evolutionary Monte Carlo algorithm for optimization with applications to sensor placement problems

In this paper, we present an adaptive evolutionary Monte Carlo algorithm (AEMC), which combines a tree-based predictive model with an evolutionary Monte Carlo sampling procedure for the purpose of global optimization. Our development is motivated by sensor placement applications in engineering, which requires optimizing certain complicated "black-box" objective function. The proposed method is able to enhance the optimization efficiency and effectiveness as compared to a few alternative strategies. AEMC falls into the category of adaptive Markov chain Monte Carlo (MCMC) algorithms and is the first adaptive MCMC algorithm that simulates multiple Markov chains in parallel. A theorem about the ergodicity property of the AEMC algorithm is stated and proven. We demonstrate the advantages of the proposed method by applying it to a sensor placement problem in a manufacturing process, as well as to a standard Griewank test function.

A low complexity user scheduling algorithm for uplink multiuser MIMO systems

A low complexity user scheduling algorithm based on a novel adaptive Markov chain Monte Carlo (AMCMC) method is proposed to achieve the maximal sum capacity in an uplink multiple-input multiple-output (MIMO) multiuser system. Compared with the existing scheduling algorithms, our algorithm is not only more efficient but also converges to within 99% of the optimal capacity obtained by exhaustive search. We demonstrate the convergence of the proposed scheduling algorithm and study the tradeoff between its complexity and performance.

Adaptive independence samplers

Markov chain Monte Carlo (MCMC) is an important computational technique for generating samples from non-standard probability distributions. A major challenge in the design of practical MCMC samplers is to achieve efficient convergence and mixing properties. One way to accelerate convergence and mixing is to adapt the proposal distribution in light of previously sampled points, thus increasing the probability of acceptance. In this paper, we propose two new adaptive MCMC algorithms based on the Independent Metropolis---Hastings algorithm. In the first, we adjust the proposal to minimize an estimate of the cross-entropy between the target and proposal distributions, using the experience of pre-runs. This approach provides a general technique for deriving natural adaptive formulae. The second approach uses multiple parallel chains, and involves updating chains individually, then updating a proposal density by fitting a Bayesian model to the population. An important feature of this approach is that adapting the proposal does not change the limiting distributions of the chains. Consequently, the adaptive phase of the sampler can be continued indefinitely. We include results of numerical experiments indicating that the new algorithms compete well with traditional Metropolis---Hastings algorithms. We also demonstrate the method for a realistic problem arising in Comparative Genomics.

Higher Accuracy for Bayesian and Frequentist Inference: Large Sample Theory for Small Sample Likelihood

Recent likelihood theory produces p-values that have remarkable accuracy and wide applicability. The calculations use familiar tools such as maximum likelihood values (MLEs), observed information and parameter rescaling. The usual evaluation of such p-values is by simulations, and such simulations do verify that the global distribution of the p-values is uniform(0, 1), to high accuracy in repeated sampling. The derivation of the p-values, however, asserts a stronger statement, that they have a uniform(0, 1) distribution conditionally, given identified precision information provided by the data. We take a simple regression example that involves exact precision information and use large sample techniques to extract highly accurate information as to the statistical position of the data point with respect to the parameter: specifically, we examine various p-values and Bayesian posterior survivor s-values for validity. With observed data we numerically evaluate the various p-values and s-values, and we also record the related general formulas. We then assess the numerical values for accuracy using Markov chain Monte Carlo (McMC) methods. We also propose some third-order likelihood-based procedures for obtaining means and variances of Bayesian posterior distributions, again followed by McMC assessment. Finally we propose some adaptive McMC methods to improve the simulation acceptance rates. All these methods are based on asymptotic analysis that derives from the effect of additional data. And the methods use simple calculations based on familiar maximizing values and related informations. The example illustrates the general formulas and the ease of calculations, while the McMC assessments demonstrate the numerical validity of the p-values as percentage position of a data point. The example, however, is very simple and transparent, and thus gives little indication that in a wide generality of models the formulas do accurately separate information for almost any parameter of interest, and then do give accurate p-value determinations from that information. As illustration an enigmatic problem in the literature is discussed and simulations are recorded; various examples in the literature are cited.

Bayesian inversion methods in ocean geoacoustics

This paper describes a complete approach to the inversion of ocean acoustic data for environmental model parameters, which is also applicable to other inverse problems. Within a Bayesian formulation, the general solution is given by the posterior probability density (PPD) of the model parameters, which includes both data and prior information. Properties of the PPD, such as optimal parameter estimates, variance/covariance, inter-parameter correlations, and marginal probability distributions, are computed numerically for nonlinear inverse problems using Markov-chain Monte Carlo (MCMC) importance sampling methods. Since the data uncertainty distribution (including measurement and theory errors) is generally not known a priori, a simple, physically-reasonable form, such as a Gaussian or double-exponential distribution, is assumed, with statistical properties estimated from data residual analysis. In many cases, the full error covariance matrix (including off-diagonal terms) is required, and in some cases effects of nonstationary errors must be included. If biased data errors are suspected, additional unknown parameters representing the biases are included explicitly in the inversion. The validity/applicability of the above assumptions and estimates is examined a posteriori by applying both qualitative and quantitative statistical tests. New advances in efficient and adaptive MCMC sampling for nonlinear inversion will be presented.

On the efficiency of adaptive MCMC algorithms

We study a class of adaptive Markov Chain Monte Carlo (MCMC) processes which aim at behaving as an ``optimal'' target process via a learning procedure. We show, under appropriate conditions, that the adaptive MCMC chain and the ``optimal'' (nonadaptive) MCMC process share many asymptotic properties. The special case of adaptive MCMC algorithms governed by stochastic approximation is considered in details and we apply our results to the adaptive Metropolis algorithm of [Haario, Saksman, Tamminen].

Bayesian modelling of algal mass occurrences—using adaptive MCMC methods with a lake water quality model

Our study aims to estimate confounded effects of nutrients and grazing zooplankton (Crustacea) on phytoplankton groups—specifically on nitrogen-fixing Cyanobacteria—in the shallow, mesotrophic Lake Pyhäjärvi in the northern hemisphere (Finland, northern Europe, lat. 60°54′–61°06′, long. 22°09′–22°22′). Phytoplankton is modelled with a non-linear dynamic model which describes the succession of three dominant algae groups (Diatomophyceae, Chrysophyceae, nitrogen-fixing Cyanobacteria) and minor groups summed together as a function of total phosphorus, total nitrogen, temperature, global irradiance and crustacean zooplankton grazing. The model is fitted using 8 years of in situ observations and adaptive Markov chain Monte Carlo (MCMC) methods for estimation of model parameters. The approach offers a way to deal with noisy data and a large number of weakly identifiable parameters in a model. From our posterior simulations we calculate the lower limit for zooplankton carbon mass concentration (45 μgC L −1) and the upper limit for total phosphorus concentration (16 μg L −1) that satisfy with 0.95 probability our predefined water quality criteria (Cyanobacteria concentration during late summer period does not exceed the value 0.86 mg L −1). Within the observational range total phosphorus has marginal effect on Cyanobacteria compared to the zooplankton grazing effect, which is temperature-dependent. Extensive fishing efforts are needed to attain the criteria.

Adaptive MCMC methods for inference on affine stochastic volatility models with jumps

Summary In this paper we propose an efficient Markov chain Monte Carlo (MCMC) algorithm to estimate stochastic volatility models with jumps and affine structure. Our idea relies on the use of adaptive methods that aim at reducing the asymptotic variance of the estimates. We focus on the Delayed Rejection algorithm in order to find accurate proposals and to efficiently simulate the volatility path. Furthermore, Bayesian model selection is addressed through the use of reduced runs of the MCMC together with an auxiliary particle filter necessary to evaluate the likelihood function. An empirical application based on the study of the Dow Jones Composite 65 and of the FTSE 100 financial indexes is presented to study some empirical properties of the algorithm implemented.