On the utility of graphics cards to perform massively parallel simulation of advanced Monte Carlo methods

We present a case-study on the utility of graphics cards to perform massively parallel simulation of advanced Monte Carlo methods. Graphics cards, containing multiple Graphics Processing Units (GPUs), are self-contained parallel computational devices that can be housed in conventional desktop and laptop computers and can be thought of as prototypes of the next generation of many-core processors. For certain classes of population-based Monte Carlo algorithms they offer massively parallel simulation, with the added advantage over conventional distributed multi-core processors that they are cheap, easily accessible, easy to maintain, easy to code, dedicated local devices with low power consumption. On a canonical set of stochastic simulation examples including population-based Markov chain Monte Carlo methods and Sequential Monte Carlo methods, we nd speedups from 35 to 500 fold over conventional single-threaded computer code. Our findings suggest that GPUs have the potential to facilitate the growth of statistical modelling into complex data rich domains through the availability of cheap and accessible many-core computation. We believe the speedup we observe should motivate wider use of parallelizable simulation methods and greater methodological attention to their design.     

On the Utility of Graphics Cards to Perform Massively Parallel Simulation of Advanced Monte Carlo Methods

We present a case study on the utility of graphics cards to perform massively parallel simulation of advanced Monte Carlo methods. Graphics cards, containing multiple Graphics Processing Units (GPUs), are self-contained parallel computational devices that can be housed in conventional desktop and laptop computers and can be thought of as prototypes of the next generation of many-core processors. For certain classes of population-based Monte Carlo algorithms they offer massively parallel simulation, with the added advantage over conventional distributed multicore processors that they are cheap, easily accessible, easy to maintain, easy to code, dedicated local devices with low power consumption. On a canonical set of stochastic simulation examples including population-based Markov chain Monte Carlo methods and Sequential Monte Carlo methods, we find speedups from 35- to 500-fold over conventional single-threaded computer code. Our findings suggest that GPUs have the potential to facilitate the growth of statistical modeling into complex data-rich domains through the availability of cheap and accessible many-core computation. We believe the speedup we observe should motivate wider use of parallelizable simulation methods and greater methodological attention to their design. This article has supplementary material online.     

Combining stochastic and deterministic approaches within high efficiency molecular simulations

Generalized Shadow Hybrid Monte Carlo (GSHMC) is a method for molecular simulations that rigorously alternates Monte Carlo sampling from a canonical ensemble with integration of trajectories using Molecular Dynamics (MD). While conventional hybrid Monte Carlo methods completely re-sample particle’s velocities between MD trajectories, our method suggests a partial velocity update procedure which keeps a part of the dynamic information throughout the simulation. We use shadow (modified) Hamiltonians, the asymptotic expansions in powers of the discretization parameter corresponding to timestep, which are conserved by symplectic integrators to higher accuracy than true Hamiltonians. We present the implementation of this method into the highly efficient MD code GROMACS and demonstrate its performance and accuracy on computationally expensive systems like proteins in comparison with the molecular dynamics techniques already available in GROMACS. We take advantage of the state-of-the-art algorithms adopted in the code, leading to an optimal implementation of the method. Our implementation introduces virtually no overhead and can accurately recreate complex biological processes, including rare event dynamics, saving much computational time compared with the conventional simulation methods.     

资产跳跃下CM策略多期收入保证价格模拟

恒定混合策略(CM策略)多期收入保证价格是保本基金发行方采取设置止损的CM\linebreak策略作为投资策略时收取保 本费的理论依据, 其中标的资产由复合泊松过程和维纳过程共同驱动, 这一定价问题内嵌奇异期权, 蒙特卡罗模拟方法擅长处理这种高维数量金融问题. 基于风险中性测度推导出多期收入保证价格的现值表达式, 用条件蒙特卡罗推导出这一现值表达式的模拟公式. 在给定参数下分别用普通蒙特卡罗和条件蒙特卡罗计算CM策略多期收入保证价格的数值解, 结果显示两种蒙特卡罗方法均能有效计算其数值解, 之后通过给定显著性水平下的置信区间长度评价两种方法的精确度, 结果显示条件蒙特卡罗比普通蒙特卡罗有很大改进. 接着运用条件蒙特卡罗模拟研究多期收入保证价格对不同参数范围的变化情况.     

基于光的波粒二象性猜想的数学仿真

从现有的经典物理光学理论和专业实验结果出发,运用数学思维,综合光子理论,建立了基于光的波粒二象性猜想的四种数学模型.针对光微子碰撞猜想,建立了基于光子碰撞后概率分布的模型.针对光子作为电磁场自我旋转的猜想,分别从专业证明和数学模型分析方面建立了电磁场偏转模型和光子旋转模型.最后建立了我们自己的猜想模型——光子蜂窝网络模型.该模型引入了"光子域"、"光子电力"、"光子磁力"、"光子键"等概念,从五个子模型出发,定性解释了四个光学现象,合理回答了题目提出的三大问题,并定量证明了衍射光强分布.     

Efficient implementation of the Barnes-Hut octree algorithm for Monte Carlo simulations of charged systems

Computer simulation with Monte Carlo is an important tool to investigate the function and equilibrium properties of many biological and soft matter materials solvable in solvents.The appropriate treatment of long-range electrostatic interaction is essential for these charged systems,but remains a challenging problem for large-scale simulations.We develop an efficient Barnes-Hut treecode algorithm for electrostatic evaluation in Monte Carlo simulations of Coulomb many-body systems.The algorithm is based on a divide-and-conquer strategy and fast update of the octree data structure in each trial move through a local adjustment procedure.We test the accuracy of the tree algorithm,and use it to perform computer simulations of electric double layer near a spherical interface.It is shown that the computational cost of the Monte Carlo method with treecode acceleration scales as log N in each move.For a typical system with ten thousand particles,by using the new algorithm,the speed has been improved by two orders of magnitude from the direct summation.     

A Nonlinear Deterministic Mode Decomposition Strategy for High-Dimensional Monte Carlo Simulations

We present a new reduced Monte Carlo simulation strategy for nonlinear high-order systems based on an extension of the proper orthogonal decomposition for transient excited structures. This novel approach enables the evaluation of the response statistics in a fractional amount of time compared to the full Monte Carlo simulation procedure. (© 2016 Wiley-VCH Verlag GmbH & Co. KGaA, Weinheim)     

Measuring Technical,Allocative and Overall Efficiencies Using a DEA Algorithm

This study presents an effectively-designed algorithm for measuring technical, allocative and overall efficiencies, using data envelopment analysis (DEA). The specially designed DEA code takes advantage of unique features related to the DEA algorithm. Computational efficiency of the proposed DEA algorithm is confirmed by a Monte Carlo simulation study.     

基于Copula的VsR度量与事后检验

估计VaR的传统方法有三种:协方差矩阵法、历史模拟法和蒙特仁洛模拟法。通常,文献中认为刚蒙特卡洛模拟法度量VaR有很多方面的优点。但是,本文通过实证检验发现,使用传统蒙特卡洛模拟法估计的VaR偏小,事后检验效果很不理想。本文引入Copula函数来改进传统的蒙特卡洛模拟法。Copula函数能将单个边际分布和多元联合分布联系起来,能处理非正态的边际分布,并且它度量的相关性不再局限于线性相关性。实证检验表明,基于Copula的蒙特卡罗模拟法可以更加准确地度量资产组合的VaR。     

MCMC Estimation of Restricted Covariance Matrices

This article is motivated by the difficulty of applying standard simulation techniques when identification constraints or theoretical considerations induce covariance restrictions in multivariate models. To deal with this difficulty, we build upon a decomposition of positive definite matrices and show that it leads to straightforward Markov chain Monte Carlo samplers for restricted covariance matrices. We introduce the approach by reviewing results for multivariate Gaussian models without restrictions, where standard conjugate priors on the elements of the decomposition induce the usual Wishart distribution on the precision matrix and vice versa. The unrestricted case provides guidance for constructing efficient Metropolis–Hastings and accept-reject Metropolis–Hastings samplers in more complex settings, and we describe in detail how simulation can be performed under several important constraints. The proposed approach is illustrated in a simulation study and two applications in economics. Supplemental materials for this article (appendixes, data, and computer code) are available online.     

Method of Moments Using Monte Carlo Simulation

Abstract

We present a computational approach to the method of moments using Monte Carlo simulation. Simple algebraic identities are used so that all computations can be performed directly using simulation draws and computation of the derivative of the log-likelihood. We present a simple implementation using the Newton-Raphson algorithm with the understanding that other optimization methods may be used in more complicated problems. The method can be applied to families of distributions with unknown normalizing constants and can be extended to least squares fitting in the case that the number of moments observed exceeds the number of parameters in the model. The method can be further generalized to allow “moments” that are any function of data and parameters, including as a special case maximum likelihood for models with unknown normalizing constants or missing data. In addition to being used for estimation, our method may be useful for setting the parameters of a Bayes prior distribution by specifying moments of a distribution using prior information. We present two examples—specification of a multivariate prior distribution in a constrained-parameter family and estimation of parameters in an image model. The former example, used for an application in pharmacokinetics, motivated this work. This work is similar to Ruppert's method in stochastic approximation, combines Monte Carlo simulation and the Newton-Raphson algorithm as in Penttinen, uses computational ideas and importance sampling identities of Gelfand and Carlin, Geyer, and Geyer and Thompson developed for Monte Carlo maximum likelihood, and has some similarities to the maximum likelihood methods of Wei and Tanner.     

Risk analysis with contractual default. Does covenant breach matter?

Mergers and acquisitions (M&A), private equity and leveraged buyouts, securitization and project finance are characterized by the presence of contractual clauses (covenants). These covenants trigger the technical default of the borrower even in the absence of insolvency. Therefore, borrowers may default on loans even when they have sufficient available cash to repay outstanding debt. This condition is not captured by the net present value (NPV) distribution obtained through a standard Monte Carlo simulation. In this paper, we present a methodology for including the consequences of covenant breach in a Monte Carlo simulation, extending traditional risk analysis in investment planning. We introduce a conceptual framework for modeling technical and material breaches from the standpoint of both lenders and shareholders. We apply this framework to a real case study concerning the project financing of a 64-million euro biomass power plant. The simulation is carried out on the actual model developed by the financial advisor of the project and made available to the authors. Results show that both technical and material breaches have a statistically significant impact on the net present value distribution, and this impact is more relevant when leverage and cost of debt increase.     

Implementations of the Monte Carlo EM Algorithm

The Monte Carlo EM (MCEM) algorithm is a modification of the EM algorithm where the expectation in the E-step is computed numerically through Monte Carlo simulations. The most exible and generally applicable approach to obtaining a Monte Carlo sample in each iteration of an MCEM algorithm is through Markov chain Monte Carlo (MCMC) routines such as the Gibbs and Metropolis–Hastings samplers. Although MCMC estimation presents a tractable solution to problems where the E-step is not available in closed form, two issues arise when implementing this MCEM routine: (1) how do we minimize the computational cost in obtaining an MCMC sample? and (2) how do we choose the Monte Carlo sample size? We address the first question through an application of importance sampling whereby samples drawn during previous EM iterations are recycled rather than running an MCMC sampler each MCEM iteration. The second question is addressed through an application of regenerative simulation. We obtain approximate independent and identical samples by subsampling the generated MCMC sample during different renewal periods. Standard central limit theorems may thus be used to gauge Monte Carlo error. In particular, we apply an automated rule for increasing the Monte Carlo sample size when the Monte Carlo error overwhelms the EM estimate at any given iteration. We illustrate our MCEM algorithm through analyses of two datasets fit by generalized linear mixed models. As a part of these applications, we demonstrate the improvement in computational cost and efficiency of our routine over alternative MCEM strategies.     

Circulant Embedding of Approximate Covariances for Inference From Gaussian Data on Large Lattices

Recently proposed computationally efficient Markov chain Monte Carlo (MCMC) and Monte Carlo expectation–maximization (EM) methods for estimating covariance parameters from lattice data rely on successive imputations of values on an embedding lattice that is at least two times larger in each dimension. These methods can be considered exact in some sense, but we demonstrate that using such a large number of imputed values leads to slowly converging Markov chains and EM algorithms. We propose instead the use of a discrete spectral approximation to allow for the implementation of these methods on smaller embedding lattices. While our methods are approximate, our examples indicate that the error introduced by this approximation is small compared to the Monte Carlo errors present in long Markov chains or many iterations of Monte Carlo EM algorithms. Our results are demonstrated in simulation studies, as well as in numerical studies that explore both increasing domain and fixed domain asymptotics. We compare the exact methods to our approximate methods on a large satellite dataset, and show that the approximate methods are also faster to compute, especially when the aliased spectral density is modeled directly. Supplementary materials for this article are available online.     

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.     

Numerical experiments on the Helmholtz equation derived from the solar radiation transfer equation in three-dimensional space

A numerical scheme using the finite-difference approach to solve the modified Helmholtz partial differential equation derived from the solar radiative transfer equation is developed and tested along with the method of evaluating the slant-path optical depth. For overhead solar incidence, we obtain good agreement between the finite-difference approach and the semianalytical solution in terms of the local intensity, local flux, average intensity, and average flux. For face-parallel oblique solar incidence in which the semianalytical method is not applicable, we compare the results with those of previous studies utilizing the Monte Carlo method and the approximate semianalytical method. We show that the present numerical scheme can be applied to any incident solar angle which the approximate semianalytical method is incapable of. Comparisons with results from Monte Carlo method reveal reasonable agreement for the averaged intensity and flux density.     

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.

     

Monte Carlo EM加速算法

EM算法是近年来常用的求后验众数的估计的一种数据增广算法, 但由于求出其E步中积分的显示表达式有时很困难, 甚至不可能, 限制了其应用的广泛性. 而Monte Carlo EM算法很好地解决了这个问题, 将EM算法中E步的积分用Monte Carlo模拟来有效实现, 使其适用性大大增强. 但无论是EM算法, 还是Monte Carlo EM算法, 其收敛速度都是线性的, 被缺损信息的倒数所控制, 当缺损数据的比例很高时, 收敛速度就非常缓慢. 而Newton-Raphson算法在后验众数的附近具有二次收敛速率. 本文提出Monte Carlo EM加速算法, 将Monte Carlo EM算法与Newton-Raphson算法结合, 既使得EM算法中的E步用Monte Carlo模拟得以实现, 又证明了该算法在后验众数附近具有二次收敛速度. 从而使其保留了Monte Carlo EM算法的优点, 并改进了Monte Carlo EM算法的收敛速度. 本文通过数值例子, 将Monte Carlo EM加速算法的结果与EM算法、Monte Carlo EM算法的结果进行比较, 进一步说明了Monte Carlo EM加速算法的优良性.     

Participating life insurance policies: an accurate and efficient parallel software for COTS clusters

In this paper we discuss the development of a parallel software for the numerical simulation of Participating Life Insurance Policies in distributed environments. The main computational kernels in the mathematical models for the solution of the problem are multidimensional integrals and stochastic differential equations. The former is solved by means of Monte Carlo method combined with the Antithetic Variates variance reduction technique, while differential equations are approximated via a fully implicit, positivity-preserving, Euler method. The parallelization strategy we adopted relies on the parallelization of Monte Carlo algorithm. We implemented and tested the software on a PC Linux cluster.     

Integrated structural approach to Credit Value Adjustment

This paper proposes an integrated pricing framework for Credit Value Adjustment of equity and commodity products. The given framework, in fact, generates dependence endogenously, allows for calibration and pricing to be based on the same numerical schemes (up to Monte Carlo simulation), and also allows the inclusion of risk mitigation clauses such as netting, collateral and initial margin provisions. The model is based on a structural approach which uses correlated Lévy processes with idiosyncratic and systematic components; the pricing numerical scheme, instead, efficiently combines Monte Carlo simulation and Fourier transform based methods. We illustrate the tractability of the proposed framework and the performance of the proposed numerical scheme by means of a case study on a portfolio of commodity swaps using real market data.