Density Estimation by Monte Carlo and Quasi-Monte Carlo

Estimating the density of a continuous random variable $X$ has been studied extensively in statistics, in the setting where $n$ independent observations of $X$ are given a priori and one wishes to estimate the density from that. Popular methods include histograms and kernel density estimators. In this review paper, we are interested instead in the situation where the observations are generated by Monte Carlo simulation from a model. Then, one can take advantage of variance reduction methods such as stratification, conditional Monte Carlo, and randomized quasi-Monte Carlo (RQMC), and obtain a more accurate density estimator than with standard Monte Carlo for a given computing budget. We discuss several ways of doing this, proposed in recent papers, with a focus on methods that exploit RQMC. A first idea is to directly combine RQMC with a standard kernel density estimator. Another one is to adapt a simulation-based derivative estimation method such as smoothed perturbation analysis or the likelihood ratio method to obtain a continuous estimator of the cumulative density function (CDF), whose derivative is an unbiased estimator of the density. This can then be combined with RQMC. We summarize recent theoretical results with these approaches and give numerical illustrations of how they improve the convergence of the mean square integrated error.NSERC Discovery Grant, IVADO Grant, ERDF, ESF, EXP. 2019/0043

A Randomized Quasi-Monte Carlo Simulation Method for Markov Chains

We introduce and study a randomized quasi-Monte Carlo method for estimating the state distribution at each step of a Markov chain, under the assumption that the chain has a totally ordered (discrete or continuous) state space. The number of steps in the chain can be random and unbounded. The method simulates $n$ copies of the chain in parallel, using a $(d+1)$-dimensional low-discrepancy point set of cardinality $n$, randomized independently at each step, where $d$ is the number of uniform random numbers required at each transition of the Markov chain. This technique is effective in particular to obtain a low-variance unbiased estimator of the expected total cost up to some random stopping time, when state-dependent costs are paid at each step. We provide numerical illustrations where the variance reduction with respect to standard Monte Carlo is substantial. The variance is reduced by factors of several thousands in some cases. We prove bounds on the convergence rate of the worst-case error and variance for special situations. In line with what is typically observed in RQMC contexts, our empirical results indicate much better convergence than what these bounds guarantee

Variance Reduction with Array-RQMC for Tau-Leaping Simulation of Stochastic Biological and Chemical Reaction Networks

We explore the use of Array-RQMC, a randomized quasi-Monte Carlo method designed for the simulation of Markov chains, to reduce the variance when simulating stochastic biological or chemical reaction networks with $\tau$-leaping. The task is to estimate the expectation of a function of molecule copy numbers at a given future time $T$ by the sample average over $n$ sample paths, and the goal is to reduce the variance of this sample-average estimator. We find that when the method is properly applied, variance reductions by factors in the thousands can be obtained. These factors are much larger than those observed previously by other authors who tried RQMC methods for the same examples. Array-RQMC simulates an array of realizations of the Markov chain and requires a sorting function to reorder these chains according to their states, after each step. The choice of sorting function is a key ingredient for the efficiency of the method, although in our experiments, Array-RQMC was never worse than ordinary Monte Carlo, regardless of the sorting method. The expected number of reactions of each type per step also has an impact on the efficiency gain.ERDF, ESF, EXP. 2019/00432, Canada Research Chair, IVADO Research Grant, NSERC Discovery Grant RGPIN-110050

Monte Carlo and Quasi-Monte Carlo Density Estimation via Conditioning

Estimating the unknown density from which a given independent sample originates is more difficult than estimating the mean, in the sense that for the best popular non-parametric density estimators, the mean integrated square error converges more slowly than at the canonical rate of $\mathcal{O}(1/n)$. When the sample is generated from a simulation model and we have control over how this is done, we can do better. We examine an approach in which conditional Monte Carlo yields, under certain conditions, a random conditional density which is an unbiased estimator of the true density at any point. By averaging independent replications, we obtain a density estimator that converges at a faster rate than the usual ones. Moreover, combining this new type of estimator with randomized quasi-Monte Carlo to generate the samples typically brings a larger improvement on the error and convergence rate than for the usual estimators, because the new estimator is smoother as a function of the underlying uniform random numbers.IVADO Research Grant, NSERC-Canada Discorvery Grant, Canada Research Chair, Inria International Chair, ERDF, ESF, EXP. 2019/0043

A Tool for Custom Construction of QMC and RQMC Point Sets

We present LatNet Builder, a software tool to find good parameters for lattice rules, polynomial lattice rules, and digital nets in base 2, for quasi-Monte Carlo (QMC) and randomized quasi-Monte Carlo (RQMC) sampling over the s-dimensional unit hypercube. The selection criteria are figures of merit that give different weights to different subsets of coordinates. They are upper bounds on the worst-case error (for QMC) or variance (for RQMC) for integrands rescaled to have a norm of at most one in certain Hilbert spaces of functions. We summarize what are the various Hilbert spaces, discrepancies, types of weights, figures of merit, types of constructions, and search methods supported by LatNet Builder. We briefly discuss its organization and we provide simple illustrations of what it can do.NSERC Discovery Grant, IVADO Grant, Corps des Mines Stipend, ERDF, ESF, EXP. 2019/0043

Applying dissipative dynamical systems to pseudorandom number generation: Equidistribution property and statistical independence of bits at distances up to logarithm of mesh size

The behavior of a family of dissipative dynamical systems representing transformations of two-dimensional torus is studied on a discrete lattice and compared with that of conservative hyperbolic automorphisms of the torus. Applying dissipative dynamical systems to generation of pseudorandom numbers is shown to be advantageous and equidistribution of probabilities for the sequences of bits can be achieved. A new algorithm for generating uniform pseudorandom numbers is proposed. The theory of the generator, which includes proofs of periodic properties and of statistical independence of bits at distances up to logarithm of mesh size, is presented. Extensive statistical testing using available test packages demonstrates excellent results, while the speed of the generator is comparable to other modern generators.Comment: 6 pages, 3 figures, 3 table

Physical tests for Random Numbers in Simulations

We propose three physical tests to measure correlations in random numbers used in Monte Carlo simulations. The first test uses autocorrelation times of certain physical quantities when the Ising model is simulated with the Wolff algorithm. The second test is based on random walks, and the third on blocks of n successive numbers. We apply the tests to show that recent errors in high precision simulations using generalized feedback shift register algorithms are due to short range correlations in random number sequences. We also determine the length of these correlations.Comment: 16 pages, Post Script file, HU-TFT-94-