Optimal importance sampling for Lévy processes

We develop importance sampling estimators for Monte Carlo pricing of European and path-dependent options in models driven by Lévy processes. Using results from the theory of large deviations for processes with independent increments, we compute an explicit asymptotic approximation for the variance of the pay-off under a time-dependent Esscher-style change of measure. Minimizing this asymptotic variance using convex duality, we then obtain an importance sampling estimator of the option price. We show that our estimator is logarithmically optimal among all importance sampling estimators. Numerical tests in the variance gamma model show consistent variance reduction with a small computational overhead.     

Multilevel Monte Carlo in approximate Bayesian computation

In the following article, we consider approximate Bayesian computation (ABC) inference. We introduce a method for numerically approximating ABC posteriors using the multilevel Monte Carlo (MLMC). A sequential Monte Carlo version of the approach is developed and it is shown under some assumptions that for a given level of mean square error, this method for ABC has a lower cost than i.i.d. sampling from the most accurate ABC approximation. Several numerical examples are given.     

Weak antithetic MLMC estimation of SDEs with the Milstein scheme for low-dimensional Wiener processes

In this paper, we implement a weak Milstein Scheme to simulate low-dimensional stochastic differential equations (SDEs). We prove that combining the antithetic multilevel Monte-Carlo (MLMC) estimator introduced by Giles and Szpruch with the MLMC approach for weak SDE approximation methods by Belomestny and Nagapetyan, we can achieve a quadratic computational complexity in the inverse of the Root Mean Square Error (RMSE) when estimating expected values of smooth functionals of SDE solutions, without simulating Lévy areas and without requiring any strong convergence of the underlying SDE approximation method. By using appropriate discrete variables this approach allows us to calculate the expectation on the coarsest level of resolution by enumeration, which, for low-dimensional problems, results in a reduced computational effort compared to standard MLMC sampling. These theoretical results are also confirmed by a numerical experiment.     

Statistical Inference for the Parameter of Rayleigh Distribution Based on Progressively Typ e-I Interval Censored Sample

In this paper, the estimation of parameters based on a progressively type-I interval censored sample from a Rayleigh distribution is studied. Different methods of estimation are discussed. They include...     

Honest Importance Sampling With Multiple Markov Chains

Importance sampling is a classical Monte Carlo technique in which a random sample from one probability density, π₁, is used to estimate an expectation with respect to another, π. The importance sampling estimator is strongly consistent and, as long as two simple moment conditions are satisfied, it obeys a central limit theorem (CLT). Moreover, there is a simple consistent estimator for the asymptotic variance in the CLT, which makes for routine computation of standard errors. Importance sampling can also be used in the Markov chain Monte Carlo (MCMC) context. Indeed, if the random sample from π₁ is replaced by a Harris ergodic Markov chain with invariant density π₁, then the resulting estimator remains strongly consistent. There is a price to be paid, however, as the computation of standard errors becomes more complicated. First, the two simple moment conditions that guarantee a CLT in the iid case are not enough in the MCMC context. Second, even when a CLT does hold, the asymptotic variance has a complex form and is difficult to estimate consistently. In this article, we explain how to use regenerative simulation to overcome these problems. Actually, we consider a more general setup, where we assume that Markov chain samples from several probability densities, π₁, …, π_k, are available. We construct multiple-chain importance sampling estimators for which we obtain a CLT based on regeneration. We show that if the Markov chains converge to their respective target distributions at a geometric rate, then under moment conditions similar to those required in the iid case, the MCMC-based importance sampling estimator obeys a CLT. Furthermore, because the CLT is based on a regenerative process, there is a simple consistent estimator of the asymptotic variance. We illustrate the method with two applications in Bayesian sensitivity analysis. The first concerns one-way random effect models under different priors. The second involves Bayesian variable selection in linear regression, and for this application, importance sampling based on multiple chains enables an empirical Bayes approach to variable selection.     

Rare Events in Random Geometric Graphs

This work introduces and compares approaches for estimating rare-event probabilities related to the number of edges in the random geometric graph on a Poisson point process. In the one-dimensional setting, we derive closed-form expressions for a variety of conditional probabilities related to the number of edges in the random geometric graph and develop conditional Monte Carlo algorithms for estimating rare-event probabilities on this basis. We prove rigorously a reduction in variance when compared to the crude Monte Carlo estimators and illustrate the magnitude of the improvements in a simulation study. In higher dimensions, we use conditional Monte Carlo to remove the fluctuations in the estimator coming from the randomness in the Poisson number of nodes. Finally, building on conceptual insights from large-deviations theory, we illustrate that importance sampling using a Gibbsian point process can further substantially reduce the estimation variance.

     

Approximating zero-variance importance sampling in a reliability setting

We consider a class of Markov chain models that includes the highly reliable Markovian systems (HRMS) often used to represent the evolution of multicomponent systems in reliability settings. We are interested in the design of efficient importance sampling (IS) schemes to estimate the reliability of such systems by simulation. For these models, there is in fact a zero-variance IS scheme that can be written exactly in terms of a value function that gives the expected cost-to-go (the exact reliability, in our case) from any state of the chain. This IS scheme is impractical to implement exactly, but it can be approximated by approximating this value function. We examine how this can be effectively used to estimate the reliability of a highly-reliable multicomponent system with Markovian behavior. In our implementation, we start with a simple crude approximation of the value function, we use it in a first-order IS scheme to obtain a better approximation at a few selected states, then we interpolate in between and use this interpolation in our final (second-order) IS scheme. In numerical illustrations, our approach outperforms the popular IS heuristics previously proposed for this class of problems. We also perform an asymptotic analysis in which the HRMS model is parameterized in a standard way by a rarity parameter ε, so that the relative error (or relative variance) of the crude Monte Carlo estimator is unbounded when ε→0. We show that with our approximation, the IS estimator has bounded relative error (BRE) under very mild conditions, and vanishing relative error (VRE), which means that the relative error converges to 0 when ε→0, under slightly stronger conditions.     

Monte Carlo Integration With Acceptance-Rejection

This article considers Monte Carlo integration under rejection sampling or Metropolis-Hastings sampling. Each algorithm involves accepting or rejecting observations from proposal distributions other than a target distribution. While taking a likelihood approach, we basically treat the sampling scheme as a random design, and define a stratified estimator of the baseline measure. We establish that the likelihood estimator has no greater asymptotic variance than the crude Monte Carlo estimator under rejection sampling or independence Metropolis-Hastings sampling. We employ a subsampling technique to reduce the computational cost, and illustrate with three examples the computational effectiveness of the likelihood method under general Metropolis-Hastings sampling.     

A Markov Chain Monte Carlo comparison of variance estimators for the sampling of particulate mixtures

During the sampling of particulate mixtures, samples taken are analyzed for their mass concentration, which generally has non‐zero sample‐to‐sample variance. Bias, variance, and mean squared error (MSE) of a number of variance estimators, derived by Geelhoed, were studied in this article. The Monte Carlo simulation was applied using an observable first‐order Markov Chain with transition probabilities that served as a model for the sample drawing process. Because the bias and variance of a variance estimator could depend on the specific circumstances under which it is applied, Monte Carlo simulation was performed for a wide range of practically relevant scenarios. Using the ‘smallest mean squared error’ as a criterion, an adaptation of an estimator based on a first‐order Taylor linearization of the sample concentration is the best. An estimator based on the Horvitz–Thompson estimator is not practically applicable because of the potentially high MSE for the cases studied. The results indicate that the Poisson estimator leads to a biased estimator for the variance of fundamental sampling error (up to 428% absolute value of relative bias) in case of low levels of grouping and segregation. The uncertainty of the results obtained by the simulations was also addressed and it was found that the results were not significantly affected. The potentials of a recently described other approach are discussed for extending the first‐order Markov Chain described here to account also for higher levels of grouping and segregation. Copyright © 2013 John Wiley & Sons, Ltd.     

Numerical approximation of conditional asymptotic variances using Monte Carlo simulation

We consider in this paper the use of Monte Carlo simulation to numerically approximate the asymptotic variance of an estimator of a population parameter. When the variance of an estimator does not exist in finite samples, the variance of its limiting distribution is often used for inferences. However, in this case, the numerical approximation of asymptotic variances is less straightforward, unless their analytical derivation is mathematically tractable. The method proposed does not assume the existence of variance in finite samples. If finite sample variance does exist, it provides a more efficient approximation than the one based on the convergence of finite sample variances. Furthermore, the results obtained will be potentially useful in evaluating and comparing different estimation procedures based on their asymptotic variances for various types of distributions. The method is also applicable in surveys where the sample size required to achieve a fixed margin of error is based on the asymptotic variance of the estimator. The proposed method can be routinely applied and alleviates the complex theoretical treatment usually associated with the analytical derivation of the asymptotic variance of an estimator which is often managed on a case by case basis. This is particularly appealing in view of the advance of modern computing technology. The proposed numerical approximation is based on the variances of a certain truncated statistic for two selected sample sizes, using a Richardson extrapolation type formulation. The variances of the truncated statistic for the two sample sizes are computed based on Monte Carlo simulations, and the theory for optimizing the computing resources is also given. The accuracy of the proposed method is numerically demonstrated in a classical errors-in-variables model where analytical results are available for the purpose of comparisons.     

Bayesian inference and life testing plans for generalized exponential distribution

Recently generalized exponential distribution has received considerable attentions. In this paper, we deal with the Bayesian inference of the unknown parameters of the progressively censored generalized exponential distribution. It is assumed that the scale and the shape parameters have independent gamma priors. The Bayes estimates of the unknown parameters cannot be obtained in the closed form. Lindley’s approximation and importance sampling technique have been suggested to compute the approximate Bayes estimates. Markov Chain Monte Carlo method has been used to compute the approximate Bayes estimates and also to construct the highest posterior density credible intervals. We also provide different criteria to compare two different sampling schemes and hence to find the optimal sampling schemes. It is observed that finding the optimum censoring procedure is a computationally expensive process. And we have recommended to use the sub-optimal censoring procedure, which can be obtained very easily. Monte Carlo simulations are performed to compare the performances of the different methods and one data analysis has been performed for illustrative purposes. This work was partially supported by a grant from the Department of Science and Technology, Government of India     

The cross-entropy method with patching for rare-event simulation of large Markov chains

There are various importance sampling schemes to estimate rare event probabilities in Markovian systems such as Markovian reliability models and Jackson networks. In this work, we present a general state-dependent importance sampling method which partitions the state space and applies the cross-entropy method to each partition. We investigate two versions of our algorithm and apply them to several examples of reliability and queueing models. In all these examples we compare our method with other importance sampling schemes. The performance of the importance sampling schemes is measured by the relative error of the estimator and by the efficiency of the algorithm. The results from experiments show considerable improvements both in running time of the algorithm and the variance of the estimator.     

Efficiency of Monte Carlo EM and Simulated Maximum Likelihood in Two-Stage Hierarchical Models

Likelihood estimation in hierarchical models is often complicated by the fact that the likelihood function involves an analytically intractable integral. Numerical approximation to this integral is an option but it is generally not recommended when the integral dimension is high. An alternative approach is based on the ideas of Monte Carlo integration, which approximates the intractable integral by an empirical average based on simulations. This article investigates the efficiency of two Monte Carlo estimation methods, the Monte Carlo EM (MCEM) algorithm and simulated maximum likelihood (SML). We derive the asymptotic Monte Carlo errors of both methods and show that, even under the optimal SML importance sampling distribution, the efficiency of SML decreases rapidly (relative to that of MCEM) as the missing information about the unknown parameter increases. We illustrate our results in a simple mixed model example and perform a simulation study which shows that, compared to MCEM, SML can be extremely inefficient in practical applications.     

Computing the variance of a conditional expectation via non-nested Monte Carlo

Computing the variance of a conditional expectation has often been of importance in uncertainty quantification. Sun et al. has introduced an unbiased nested Monte Carlo estimator, which they call $1\frac{1}{2}$-level simulation since the optimal inner-level sample size is bounded as the computational budget increases. In this letter, we construct unbiased non-nested Monte Carlo estimators based on the so-called pick-freeze scheme due to Sobol’. An extension of our approach to compute higher order moments of a conditional expectation is also discussed.     

Laplace Variational Approximation for Semiparametric Regression in the Presence of Heteroscedastic Errors

Variational approximations provide fast, deterministic alternatives to Markov chain Monte Carlo for Bayesian inference on the parameters of complex, hierarchical models. Variational approximations are often limited in practicality in the absence of conjugate posterior distributions. Recent work has focused on the application of variational methods to models with only partial conjugacy, such as in semiparametric regression with heteroscedastic errors. Here, both the mean and log variance functions are modeled as smooth functions of covariates. For this problem, we derive a mean field variational approximation with an embedded Laplace approximation to account for the nonconjugate structure. Empirical results with simulated and real data show that our approximate method has significant computational advantages over traditional Markov chain Monte Carlo; in this case, a delayed rejection adaptive Metropolis algorithm. The variational approximation is much faster and eliminates the need for tuning parameter selection, achieves good fits for both the mean and log variance functions, and reasonably reflects the posterior uncertainty. We apply the methods to log-intensity data from a small angle X-ray scattering experiment, in which properly accounting for the smooth heteroscedasticity leads to significant improvements in posterior inference for key physical characteristics of an organic molecule.     

Population Monte Carlo Algorithm in High Dimensions

The population Monte Carlo algorithm is an iterative importance sampling scheme for solving static problems. We examine the population Monte Carlo algorithm in a simplified setting, a single step of the general algorithm, and study a fundamental problem that occurs in applying importance sampling to high-dimensional problem. The precision of the computed estimate from the simplified setting is measured by the asymptotic variance of estimate under conditions on the importance function. We demonstrate the exponential growth of the asymptotic variance with the dimension and show that the optimal covariance matrix for the importance function can be estimated in special cases.     

The Use of Dimensional and Similarity Analyses in the Propagation of Uncertainties: A Physical Example

Dimensional and similarity analyses are used in physics and engineering, specially in fluid mechanics, to reduce the dimension of the input variable space with no loss of information. Here, we apply these techniques to the propagation of uncertainties for computer codes by the Monte Carlo method, in order to reduce the variance of the estimators of the parameters of the output variable distribution. In the physics and engineering literature, dimensional analysis is often formulated intuitively in terms of physical quantities or dimensions such as time, longitude, or mass; here we use the more rigorous and more abstract generalized dimensional analysis of Moran and Marshek. The reduction of dimensionality is only successful in reducing estimator variance when applying variance-reduction techniques and not when using ordinary random sampling. In this article we use stratified sampling, and the key point of the success of the reduction in dimensionality in improving the precision of the estimates is a better measurement of the distances betwen the outputs, for given inputs. We illustrate the methodology with an application to a physical problem, a radioactive contaminant transport code. A substantial variance reduction is achieved for the estimators of the mean, variance, and distribution function of the output. Last, we present a discussion on which conditions are necessary for the method to be successful.     

Estimating Orthant Probabilities of High-Dimensional Gaussian Vectors with An Application to Set Estimation

The computation of Gaussian orthant probabilities has been extensively studied for low-dimensional vectors. Here, we focus on the high-dimensional case and we present a two-step procedure relying on both deterministic and stochastic techniques. The proposed estimator relies indeed on splitting the probability into a low-dimensional term and a remainder. While the low-dimensional probability can be estimated by fast and accurate quadrature, the remainder requires Monte Carlo sampling. We further refine the estimation by using a novel asymmetric nested Monte Carlo (anMC) algorithm for the remainder and we highlight cases where this approximation brings substantial efficiency gains. The proposed methods are compared against state-of-the-art techniques in a numerical study, which also calls attention to the advantages and drawbacks of the procedure. Finally, the proposed method is applied to derive conservative estimates of excursion sets of expensive to evaluate deterministic functions under a Gaussian random field prior, without requiring a Markov assumption. Supplementary material for this article is available online.     

Computational Aspects of Bayesian Spectral Density Estimation

Gaussian time-series models are often specified through their spectral density. Such models present several computational challenges, in particular because of the nonsparse nature of the covariance matrix. We derive a fast approximation of the likelihood for such models. We propose to sample from the approximate posterior (i.e., the prior times the approximate likelihood), and then to recover the exact posterior through importance sampling. We show that the variance of the importance sampling weights vanishes as the sample size goes to infinity. We explain why the approximate posterior may typically be multimodal, and we derive a Sequential Monte Carlo sampler based on an annealing sequence to sample from that target distribution. Performance of the overall approach is evaluated on simulated and real datasets. In addition, for one real-world dataset, we provide some numerical evidence that a Bayesian approach to semiparametric estimation of spectral density may provide more reasonable results than its frequentist counterparts. The article comes with supplementary materials, available online, that contain an Appendix with a proof of our main Theorem, a Python package that implements the proposed procedure, and the Ethernet dataset.     

Adaptive Monte Carlo Variance Reduction for Lévy Processes with Two-Time-Scale Stochastic Approximation

We propose an approach to a twofold optimal parameter search for a combined variance reduction technique of the control variates and the important sampling in a suitable pure-jump Lévy process framework. The parameter search procedure is based on the two-time-scale stochastic approximation algorithm with equilibrated control variates component and with quasi-static importance sampling one. We prove the almost sure convergence of the algorithm to a unique optimum. The parameter search algorithm is further embedded in adaptive Monte Carlo simulations in the case of the gamma distribution and process. Numerical examples of the CDO tranche pricing with the Gamma copula model and the intensity Gamma model are provided to illustrate the effectiveness of our method.