Adaptive Monte Carlo Algorithm Research Articles

Exact inference in contingency tables via stochastic approximation Monte Carlo

Monte Carlo methods for the exact inference have received much attention recently in complete or incomplete contingency table analysis. However, conventional Markov chain Monte Carlo, such as the Metropolis–Hastings algorithm, and importance sampling methods sometimes generate the poor performance by failing to produce valid tables. In this paper, we apply an adaptive Monte Carlo algorithm, the stochastic approximation Monte Carlo algorithm (SAMC; Liang, Liu, & Carroll, 2007), to the exact test of the goodness-of-fit of the model in complete or incomplete contingency tables containing some structural zero cells. The numerical results are in favor of our method in terms of quality of estimates.

Marginal Likelihood Estimation with the Cross-Entropy Method

We consider an adaptive importance sampling approach to estimating the marginal likelihood, a quantity that is fundamental in Bayesian model comparison and Bayesian model averaging. This approach is motivated by the difficulty of obtaining an accurate estimate through existing algorithms that use Markov chain Monte Carlo (MCMC) draws, where the draws are typically costly to obtain and highly correlated in high-dimensional settings. In contrast, we use the cross-entropy (CE) method, a versatile adaptive Monte Carlo algorithm originally developed for rare-event simulation. The main advantage of the importance sampling approach is that random samples can be obtained from some convenient density with little additional costs. As we are generating independent draws instead of correlated MCMC draws, the increase in simulation effort is much smaller should one wish to reduce the numerical standard error of the estimator. Moreover, the importance density derived via the CE method is grounded in information theory, and therefore, is in a well-defined sense optimal. We demonstrate the utility of the proposed approach by two empirical applications involving women’s labor market participation and U.S. macroeconomic time series. In both applications the proposed CE method compares favorably to existing estimators.

Monte Carlo sensitivity analysis of an Eulerian large-scale air pollution model

Variance-based approaches for global sensitivity analysis have been applied and analyzed to study the sensitivity of air pollutant concentrations according to variations of rates of chemical reactions. The Unified Danish Eulerian Model has been used as a mathematical model simulating a remote transport of air pollutants. Various Monte Carlo algorithms for numerical integration have been applied to compute Sobol's global sensitivity indices. A newly developed Monte Carlo algorithm based on Sobol's quasi-random points MCA-MSS has been applied for numerical integration. It has been compared with some existing approaches, namely Sobol's ΛΠτ sequences, an adaptive Monte Carlo algorithm, the plain Monte Carlo algorithm, as well as, eFAST and Sobol's sensitivity approaches both implemented in SIMLAB software. The analysis and numerical results show advantages of MCA-MSS for relatively small sensitivity indices in terms of accuracy and efficiency. Practical guidelines on the estimation of Sobol's global sensitivity indices in the presence of computational difficulties have been provided.

Geometric Convergence of Adaptive Monte Carlo Algorithms for Radiative Transport Problems Based on Importance Sampling Methods

Importance sampling is a very well-known variance-reducing technique used in Monte Carlo simulations of radiative transport. It involves a distortion of the physical (analog) transition probabilities with the goal of causing events of interest in the computation to occur more frequently than in the analog process. This distortion is then compensated by a corresponding alteration of the estimating random variable in order to remove any bias from the estimates of quantities of interest. In this paper, we construct several families of estimators based on importance sampling methods to solve general transport problems and prove that the adaptive application of each estimator produces geometric convergence of the approximate solution. We also present numerical results that illustrate important elements of the theory.

Adaptive Monte Carlo Algorithm for Sensitivity Studies of Eulerian Large-scale Air Pollution Models

Sobol’ global sensitivity indices have been used here as main numerical indicators for sensitivity of model output to input parameters. Therefore the basic mathematical problem is computing of multidimensional integrals. Numerical results for an Eulerian large-scale air pollution model have been presented. They have been obtained by using an adaptive Monte Carlo algorithm for numerical integration. The results have been compared with SIMLAB implementation of applied approaches for sensitivity analysis.

Monte Carlo algorithms for evaluating Sobol’ sensitivity indices

Sensitivity analysis is a powerful technique used to determine robustness, reliability and efficiency of a model. The main problem in this procedure is the evaluating total sensitivity indices that measure a parameter’s main effect and all the interactions involving that parameter. From a mathematical point of view this problem is presented by a set of multidimensional integrals. In this work a simple adaptive Monte Carlo technique for evaluating Sobol’ sensitivity indices is developed. A comparison of accuracy and complexity of plain Monte Carlo and adaptive Monte Carlo algorithms is presented. Numerical experiments for evaluating integrals of different dimensions are performed.

Editorial: Special issue on adaptive Monte Carlo methods

Monte Carlo methods have proved indispensible for many modern statistical applications—particularly within Bayesian statistics. When analysing complex models, the intractability of likelihoods and posterior distributions requires the use of numerical approximations, and, particularly in high dimensions, Monte Carlo approximation is often the approach of choice. The most commonly used Monte Carlo method is Markov chain Monte Carlo (MCMC) and reversible jump MCMC (Green 1995), but important alternatives include sequential Monte Carlo (Doucet et al. 2000), population Monte Carlo (Cappe et al. 2004) and Importance Sampling (see Fearnhead 2008, for more details of alternatives to MCMC). In all cases, the efficiency of a given Monte Carlo method will depend on how it is implemented. In many cases, theoretical aspects of the methods are well understood, and can be used to guide implementation. For example much is known about optimal acceptance rates for various MCMC algorithms (Roberts and Rosenthal 2001)—however using such information to tune an algorithm by hand can be very time-consuming. This special issues focusses on a class of Monte Carlo methods that “tune themselves”. The algorithms adapt as they are run, and hence are called adaptive Monte Carlo methods. The papers in this issue review such methods, suggest new adaptive Monte Carlo algorithms, and include case-studies demonstrating their efficiency. One common theme is the simplicity of implementing many adaptive Monte Carlo methods, and one hope of this special issue is that it will help encourage their use.

An adaptive Monte Carlo integration algorithm with general division approach

We propose an adaptive Monte Carlo algorithm for estimating multidimensional integrals over a hyper-rectangular region. The algorithm uses iteratively the idea of separating the domain of integration into 2 s subregions. The proposed algorithm can be applied directly to estimate the integral using an efficient way of storage. We test the algorithm for estimating the value of a 30-dimensional integral using a two-division approach. The numerical results show that the proposed algorithm gives better results than using one-division approach.

PRICING PATH-DEPENDENT OPTIONS ON STATE DEPENDENT VOLATILITY MODELS WITH A BESSEL BRIDGE

This paper develops bridge sampling path integral algorithms for pricing path-dependent options under a new class of nonlinear state dependent volatility models. Path-dependent option pricing is considered within a new (dual) Bessel family of semimartingale diffusion models, as well as the constant elasticity of variance (CEV) diffusion model, arising as a particular case of these models. The transition p.d.f.s or pricing kernels are mapped onto an underlying simpler squared Bessel process and are expressed analytically in terms of modified Bessel functions. We establish precise links between pricing kernels of such models and the randomized gamma distributions, and thereby demonstrate how a squared Bessel bridge process can be used for exact sampling of the Bessel family of paths. A Bessel bridge algorithm is presented which is based on explicit conditional distributions for the Bessel family of volatility models and is similar in spirit to the Brownian bridge algorithm. A special rearrangement and splitting of the path integral variables allows us to combine the Bessel bridge sampling algorithm with either adaptive Monte Carlo algorithms, or quasi-Monte Carlo techniques for significant numerical efficiency improvement. The algorithms are illustrated by pricing Asian-style and lookback options under the Bessel family of volatility models as well as the CEV diffusion model.

Application of an adaptive Monte Carlo algorithm to mixed logit estimation

This paper presents the application of a new algorithm for maximizing the simulated likelihood functions appearing in the estimation of mixed multinomial logit (MMNL) models. The method uses Monte Carlo sampling to produce the approximate likelihood function and dynamically adapts the number of draws on the basis of statistical estimators of the simulation error and simulation bias. Its convergence from distant starting points is ensured by a trust-region technique, in which improvement is ensured by locally maximizing a quadratic model of the objective function. Simulated data are first used to assess the quality of the results obtained and the relative performance of several algorithmic variants. These variants involve, in particular, different techniques for approximating the model’s Hessian and the substitution of the trust-region mechanism by a linesearch. The algorithm is also applied to a real case study arising in the context of a recent Belgian transportation model. The performance of the new Monte Carlo algorithm is shown to be competitive with that of existing tools using low discrepancy sequences.

A novel parallel adaptive Monte Carlo method for nonlinear Poisson equation in semiconductor devices

We present a parallel adaptive Monte Carlo (MC) algorithm for the numerical solution of the nonlinear Poisson equation in semiconductor devices. Based on a fixed random walk MC method, 1-irregular unstructured mesh technique, monotone iterative method, a posterior error estimation method, and dynamic domain decomposition algorithm, this approach is developed and successfully implemented on a 16-processors (16-PCs) Linux-cluster with message-passing interface (MPI) library. To solve the nonlinear problem with MC method, monotone iterative method is applied in each adaptive loop to obtain the final convergent solution. This approach fully exploits the inherent parallelism of the monotone iterative as well as MC methods. Numerical results for p–n diode and MOSFET devices are demonstrated to show the robustness of the method. Furthermore, achieved parallel speedup and related parallel performances are also reported in this work.

UbiquitousCPviolation in a top-inspired left-right model

We explore CP violation in a Left-Right Model that reproduces the quark mass and CKM rotation angle hierarchies in a relatively natural way by fixing the bidoublet Higgs VEVs to be in the ratio m_b:m_t. Our model is quite general and allows for CP to be broken by both the Higgs VEVs and the Yukawa couplings. Despite this generality, CP violation may be parameterized in terms of two basic phases. A very interesting feature of the model is that the mixing angles in the right-handed sector are found to be equal to their left-handed counterparts to a very good approximation. Furthermore, the right-handed analogue of the usual CKM phase delta_L is found to satisfy the relation delta_R \approx delta_L. The parameter space of the model is explored by using an adaptive Monte Carlo algorithm and the allowed regions in parameter space are determined by enforcing experimental constraints from the K and B systems. This method of solution allows us to evaluate the left- and right-handed CKM matrices numerically for various combinations of the two fundamental CP-odd phases in the model. We find that all experimental constraints may be satisfied with right-handed W and Flavour Changing Neutral Higgs masses as low as about 2 TeV and 7 TeV, respectively.

The Metropolis algorithm for on-shell four-momentum phase space

We present several implementations of the Metropolis method, an adaptive Monte Carlo algorithm, which allow for the calculation of multi-dimensional integrals over arbitrary on-shell four-momentum phase space. The Metropolis technique reveals itself very suitable for the treatment of high energy processes in particle physics, particularly when the number of final state objects and of kinematic constraints on the latter gets larger. We compare the performances of the Metropolis algorithm with those of other programs widely used in numerical simulations.

A progressive radiosity algorithm based on piecewise polynomial intensity distribution

A progressive radiosity algorithm based on piecewise polynomial intensity distribution is presented in the paper. Unlike the conventional radiosity method, radiosity across each patch is assumed to vary in a polynomial distribution. A generalized radiosity equation for scenes with non-diffuse surfaces is then set up and solved by adopting the progressive refinement technique. To improve the efficiency of the algorithm, an adaptive Monte Carlo algorithm incorporation with localization technique is employed to evaluate the form-factors from source patches to nodes. Experimental results demonstrate the potential advantage of our method for rendering complex scenes.

Adaptive Search in Quasi-Monte-Carlo Optimization

Motivated by a linear time-complexity result for an adaptive Monte Carlo algorithm, we propose and analyze an adaptive deterministic algorithm. We restrict a grid search to nested subregions that promise to provide improvement of the current solution, and we obtain an exponential rate of convergence in function evaluations. For proving the main result we restrict ourselves to functions on hypercubes. In a final section we outline how to extend the method to the general case and give some numerical examples

Adaptive search in quasi-Monte Carlo optimization

Motivated by a linear time-complexity result for an adaptive Monte Carlo algorithm, we propose and analyze an adaptive deterministic algorithm. We restrict a grid search to nested subregions that promise to provide improvement of the current solution, and we obtain an exponential rate of convergence in function evaluations. For proving the main result we restrict ourselves to functions on hypercubes. In a final section we outline how to extend the method to the general case and give some numerical examples.

The third dielectric and pressure virial coefficients of dipolar hard sphere fluids

We present numerical computations of the third dielectric virial coefficient C*e and the third pressure virial coefficient C* p of a fluid of hard sphere dipoles, over the range 0 < μ*2 < 4. The requisite 9-dimensional quadratures are performed using an adaptive Monte Carlo algorithm. We also compute terms in the high temperature expansions of the virial coefficients up to order μ*20 (C*e) and μ*24 (C* p ). Our work complements and extends earlier investigations by Joslin and by Rushbrooke and Shrubsall. Our results agree well with available experimental data, especially considering the relative crudeness of the molecular model. We also compute the critical point of the fluid.

Monte Carlo calculation of the above-threshold ionization of atomic hydrogen

The authors present a systematic study on the above-threshold ionization of the hydrogen atom using an adaptive Monte Carlo algorithm. The ionizations from both the ground state and the excited states are discussed in the weak-field regime. It is found that the Monte Carlo approach is quite efficient and may be useful to treat the high-order multiphoton ionization of other atoms.

Monte Carlo integration of multiphoton ionization matrix elements in the weak-field regime.

An adaptive Monte Carlo algorithm is applied, as a test case, to the evaluation of multiphoton ionization matrix elements for hydrogen. This study indicates that a Monte Carlo approach may be useful in treating high-order multiphoton ionization of alkali-metal atoms, and possibly other atoms, by weak fields.

Prediction of the native conformation of a polypeptide by a statistical-mechanical procedure. I. Backbone structure of enkephalin.

AbstractA new methodology for theoretically predicting the native, three‐dimensional structure of a polypeptide is presented. Based on equilibrium statistical mechanics, an algorithm has been designed to determine the probable conformation of a polypeptide by calculating conditional free‐energy maps for each residue of the macromolecule. The conditional free‐energy map of each residue is computed from a set of probability integrals, obtained by summing over the interaction energies of all pairs of nonbonded atoms of the whole molecule. By locating the region(s) of lowest free energy for each map, the probable conformation for each residue can be identified. The native structure of the polypeptide is assumed to be the combination of the probable conformations of the individual residues. All multidimensional probability integrals are evaluated by an adaptive Monte Carlo algorithm (SMAPPS—Statistical‐Mechanical Algorithm for Predicting Protein Structure). The Monte Carlo algorithm searches the entire conformational space, adjusting itself automatically to concentrate its sampling in regions where the magnitude of the integrand is largest (“importance sampling”). No assumptions are made about the native conformation. The only prior knowledge necessary for the prediction of the native conformation is the amino acid sequence of the polypeptide. To test the effectiveness of the algorithm, SMAPPS was applied to the prediction of the native conformation of the backbone of Met‐enkephalin, a pentapeptide. In the calculations, only the backbone dihedral angles (ϕ and ψ) were allowed to vary; all side‐chain (χ) and peptide‐bond (ω) dihedral angles were kept fixed at the values corresponding to the alleged global minimum energy previously determined by direct energy minimization. For each conformation generated randomly by the Monte Carlo algorithm, the total conformational energy of the polypeptide was obtained from established empirical potential energy functions. Solvent effects were not included in the computations. With this initial application of SMAPPS, three distinct low‐free‐energy β‐bend structures of Met‐enkephalin were found. In particular, one of the structures has a conformation remarkably similar to the one associated with the previously alleged global minimum energy. The two additional structures of the pentapeptide have conformational energies lower than the previously computed low‐energy structure. However, the Monte Carlo results are in agreement with an improved energy‐minimization procedure. These initial results on the backbone structure of Met‐enkephalin indicate that an equilibrium statistical‐mechanical procedure, coupled with an adaptive Monte Carlo algorithm, can overcome many of the problems associated with the standard methods of direct energy minimization.