High-dimensional Numerical Integration Research Articles

Generalized Dimension Truncation Error Analysis for High-Dimensional Numerical Integration: Lognormal Setting and Beyond

Generalized Dimension Truncation Error Analysis for High-Dimensional Numerical Integration: Lognormal Setting and Beyond

An Efficient and Fast Sparse Grid Algorithm for High-Dimensional Numerical Integration

This paper is concerned with developing an efficient numerical algorithm for the fast implementation of the sparse grid method for computing the d-dimensional integral of a given function. The new algorithm, called the MDI-SG (multilevel dimension iteration sparse grid) method, implements the sparse grid method based on a dimension iteration/reduction procedure. It does not need to store the integration points, nor does it compute the function values independently at each integration point; instead, it reuses the computation for function evaluations as much as possible by performing the function evaluations at all integration points in a cluster and iteratively along coordinate directions. It is shown numerically that the computational complexity (in terms of CPU time) of the proposed MDI-SG method is of polynomial order O(d3Nb)(b≤2) or better, compared to the exponential order O(N(logN)d−1) for the standard sparse grid method, where N denotes the maximum number of integration points in each coordinate direction. As a result, the proposed MDI-SG method effectively circumvents the curse of dimensionality suffered by the standard sparse grid method for high-dimensional numerical integration.

Transmissibility-based damage detection with hierarchical clustering enhanced by multivariate probabilistic distance accommodating uncertainty and correlation

This paper proposes a new damage detection method by integrating the advantage of transmissibility function (TF) as a health index sensitive to damage but robust to excitation and agglomerative hierarchical clustering (AHC) with intuitive explanation and visualization but avoiding specifying the number of clusters. Different from conventional AHC-based damage detection methods utilizing deterministic distance as a similarity metric and ignoring the distribution of structural features, a multivariate probabilistic distance-based similarity metric is proposed in this study to account for the uncertainty and correlation of multiple TFs following multivariate complex-valued Gaussian ratio distribution. To realize this, an analytically tractable approximation of the multivariate probabilistic distance is derived by Laplace’s asymptotic expansion to avoid high-dimensional numerical integration. To accelerate the computation of probabilistic distances over a wide frequency band that are fused to formulate the similarity metric in AHC, a function vectorization scheme is proposed to avoid the time-consuming loop operation among different frequency points. A threshold is established via bootstrapped Monte Carlo simulation to cut the dendrogram produced by AHC. Two case studies are used to validate the performance of the proposed method, indicating that, compared to the damage detection methods based on the deterministic distance of the TF, the proposed method exhibits better performance due to improving the similarity metric based on multivariate probabilistic distance properly accommodating the correlation of different TFs.

Challenging the Curse of Dimensionality in Multidimensional Numerical Integration by Using a Low-Rank Tensor-Train Format

Numerical integration is a basic step in the implementation of more complex numerical algorithms suitable, for example, to solve ordinary and partial differential equations. The straightforward extension of a one-dimensional integration rule to a multidimensional grid by the tensor product of the spatial directions is deemed to be practically infeasible beyond a relatively small number of dimensions, e.g., three or four. In fact, the computational burden in terms of storage and floating point operations scales exponentially with the number of dimensions. This phenomenon is known as the curse of dimensionality and motivated the development of alternative methods such as the Monte Carlo method. The tensor product approach can be very effective for high-dimensional numerical integration if we can resort to an accurate low-rank tensor-train representation of the integrand function. In this work, we discuss this approach and present numerical evidence showing that it is very competitive with the Monte Carlo method in terms of accuracy and computational costs up to several hundredths of dimensions if the integrand function is regular enough and a sufficiently accurate low-rank approximation is available.

Quantum Ray Tracing with Simulation

Rendering is a computer graphics process that generates images from a three-dimensional scene including models, lights and cameras. Ray tracing is a category of rendering algorithms that calculates pixel colors by simulating the physical movements of light rays such as reflections, refractions and scatterings. The core problem of ray tracing is a high-dimensional numerical integration, and the classical Monte Carlo ray tracing requires a big number of rays to reduce error. By leveraging the inherent parallelism of quantum computing, we propose a framework of quantum ray tracing algorithm, and do simulation experiments to prove that it has a quadratically faster error convergence than classical Monte Carlo ray tracing.

A note on concatenation of quasi-Monte Carlo and plain Monte Carlo rules in high dimensions

In this note, we study a concatenation of quasi-Monte Carlo and plain Monte Carlo rules for high-dimensional numerical integration in weighted function spaces. In particular, we consider approximating the integral of periodic functions defined over the s-dimensional unit cube by using rank-1 lattice point sets only for the first d(<s) coordinates and random points for the remaining s−d coordinates. We prove that, by exploiting a decay of the weights of function spaces, almost the optimal order of the mean squared worst-case error is achieved by such a concatenated quadrature rule as long as d scales at most linearly with the number of points. This result might be useful for numerical integration in extremely high dimensions, such as partial differential equations with random coefficients for which even the standard fast component-by-component algorithm is considered computationally expensive.

Complex survival trial design by the product integration method

Nonproportional hazards (NPHs) are often observed in survival trials such as the immunotherapy cancer trials. Under NPH, the classical log-rank test can be inefficient, and the estimated hazards ratio from the Cox model is difficult to interpret. The weighted log-rank test, and the tests for comparing the restricted mean survival time or the milestone survival become increasingly popular in handling NPH. The sample size calculation for these tests may require high-dimensional numerical integration. We present a sample size determination method for survival trials via product integration on the basis of a continuous-time multistate Markov model. The main challenge of the method lies in the design of the multistate model under a complex NPH pattern, and this is illustrated for NPH induced by delayed effect with individual heterogeneity in the lag duration, cure fractions, and treatment switching due to disease progression or noncompliance. Numerical examples are presented to demonstrate the accuracy of the proposed method. We obtain the following findings. The powers of the tests for milestone survival and RMST depend on both the trial duration and milestone timepoint, and may not increase as the milestone timepoint increases. If the milestone timepoint is appropriately chosen, the RMST test can be more powerful than the conventional log-rank test in the presence of diminishing treatment effect or in the proportional hazards cure model. In general, the RMST test yields lower power than a proper Fleming-Harrington weighted log-rank test.

Uncertainty quantification of offshore wind farms using Monte Carlo and sparse grid

ABSTRACT The power produced by an offshore wind farm is subject to multiple uncertainties, such as volatile wind, turbine performance wear, and availability losses. Knowledge about the propagation of these uncertainties and their effect on the produced power is crucial in the design stage of a wind farm. Due to the multitude of uncertainties, an analysis requires high-dimensional numerical integration to determine these parameter sensitivities. Such an analysis has not been done in the current literature for the full set of parameters. In this work, a thorough analysis of all uncertainties is provided, modeled from several years of collected data from the existing wind farms Horns Rev 1, DanTysk, and Sandbank. The analysis reveals four major parameters, allowing the other parameters to be neglected in future measurement data acquisitions and sensitivity analysis processes. Furthermore, the accuracy of several Uncertainty Quantification techniques is analyzed and a recommendation for future analysis is given.

I- flow: High-dimensional integration and sampling with normalizing flows

In many fields of science, high-dimensional integration is required. Numerical methods have been developed to evaluate these complex integrals. We introduce the code i-flow, a Python package that performs high-dimensional numerical integration utilizing normalizing flows. Normalizing flows are machine-learned, bijective mappings between two distributions. i-flow can also be used to sample random points according to complicated distributions in high dimensions. We compare i-flow to other algorithms for high-dimensional numerical integration and show that i-flow outperforms them for high dimensional correlated integrals. The i-flow code is publicly available on gitlab at https://gitlab.com/i-flow/i-flow.

Maximum Likelihood Estimation of Multilevel Structural Equation Models with Random Slopes for Latent Covariates

A maximum likelihood estimation routine for two-level structural equation models with random slopes for latent covariates is presented. Because the likelihood function does not typically have a closed-form solution, numerical integration over the random effects is required. The routine relies upon a method proposed by du Toit and Cudeck (Psychometrika 74(1):65-82, 2009) for reformulating the likelihood function so that an often large subset of the random effects can be integrated analytically, reducing the computational burden of high-dimensional numerical integration. The method is demonstrated and assessed using a small-scale simulation study and an empirical example.

Smoothing the payoff for efficient computation of Basket option prices

We consider the problem of pricing basket options in a multivariate Black–Scholes or Variance-Gamma model. From a numerical point of view, pricing such options corresponds to moderate and high-dimensional numerical integration problems with non-smooth integrands. Due to this lack of regularity, higher order numerical integration techniques may not be directly available, requiring the use of methods like Monte Carlo specifically designed to work for non-regular problems. We propose to use the inherent smoothing property of the density of the underlying in the above models to mollify the payoff function by means of an exact conditional expectation. The resulting conditional expectation is unbiased and yields a smooth integrand, which is amenable to the efficient use of adaptive sparse-grid cubature. Numerical examples indicate that the high-order method may perform orders of magnitude faster than Monte Carlo or Quasi Monte Carlo methods in dimensions up to 35.

Copula-Based Random Effects Models for Clustered Data

Sorting and spillovers can create correlation in individual outcomes. In this situation, standard discrete choice estimators cannot consistently estimate the probability of joint and conditional events, and alternative estimators can yield incoherent statistical models or intractable estimators. I propose a random effects estimator that models the dependence among the unobserved heterogeneity of individuals in the same cluster using a parametric copula. This estimator makes it possible to compute joint and conditional probabilities of the outcome variable, and is statistically coherent. I describe its properties, establishing its efficiency relative to standard random effects estimators, and propose a specification test for the copula. The likelihood function for each cluster is an integral whose dimension equals the size of the cluster, which may require high-dimensional numerical integration. To overcome the curse of dimensionality from which methods like Monte Carlo integration suffer, I propose an algorithm that works for Archimedean copulas. I illustrate this approach by analysing labour supply in married couples.

Digital inversive vectors can achieve polynomial tractability for the weighted star discrepancy and for multivariate integration

We study high-dimensional numerical integration in the worst-case setting. The subject of tractability is concerned with the dependence of the worst-case integration error on the dimension. Roughly speaking, an integration problem is tractable if the worst-case error does not grow exponentially fast with the dimension. Many classical problems are known to be intractable. However, sometimes tractability can be shown. Often such proofs are based on randomly selected integration nodes. Of course, in applications, true random numbers are not available and hence one mimics them with pseudorandom number generators. This motivates us to propose the use of pseudorandom vectors as underlying integration nodes in order to achieve tractability. In particular, we consider digital inverse vectors and present two examples of problems, the weighted star discrepancy and integration of Hölder continuous, absolute convergent Fourier and cosine series, where the proposed method is successful.

A CENTRAL LIMIT THEOREM FOR LATIN HYPERCUBE SAMPLING WITH DEPENDENCE AND APPLICATION TO EXOTIC BASKET OPTION PRICING

We consider the problem of estimating 𝔼[f(U1, …, Ud)], where (U1, …, Ud) denotes a random vector with uniformly distributed marginals. In general, Latin hypercube sampling (LHS) is a powerful tool for solving this kind of high-dimensional numerical integration problem. In the case of dependent components of the random vector (U1, …, Ud) one can achieve more accurate results by using Latin hypercube sampling with dependence (LHSD). We state a central limit theorem for the d-dimensional LHSD estimator, by this means generalising a result of Packham and Schmidt. Furthermore we give conditions on the function f and the distribution of (U1, …, Ud) under which a reduction of variance can be achieved. Finally we compare the effectiveness of Monte Carlo and LHSD estimators numerically in exotic basket option pricing problems.

Limited-information goodness-of-fit testing of hierarchical item factor models.

In applications of item response theory, assessment of model fit is a critical issue. Recently, limited-information goodness-of-fit testing has received increased attention in the psychometrics literature. In contrast to full-information test statistics such as Pearson's X(2) or the likelihood ratio G(2) , these limited-information tests utilize lower-order marginal tables rather than the full contingency table. A notable example is Maydeu-Olivares and colleagues'M2 family of statistics based on univariate and bivariate margins. When the contingency table is sparse, tests based on M2 retain better Type I error rate control than the full-information tests and can be more powerful. While in principle the M2 statistic can be extended to test hierarchical multidimensional item factor models (e.g., bifactor and testlet models), the computation is non-trivial. To obtain M2 , a researcher often has to obtain (many thousands of) marginal probabilities, derivatives, and weights. Each of these must be approximated with high-dimensional numerical integration. We propose a dimension reduction method that can take advantage of the hierarchical factor structure so that the integrals can be approximated far more efficiently. We also propose a new test statistic that can be substantially better calibrated and more powerful than the original M2 statistic when the test is long and the items are polytomous. We use simulations to demonstrate the performance of our new methods and illustrate their effectiveness with applications to real data.

Toward efficient computation of the expected relative entropy for nonlinear experimental design

The expected relative entropy between prior and posterior model-parameter distributions is a Bayesian objective function in experimental design theory that quantifies the expected gain in information of an experiment relative to a previous state of knowledge. The expected relative entropy is a preferred measure of experimental quality because it can handle nonlinear data-model relationships, an important fact due to the ubiquity of nonlinearity in science and engineering and its effects on post-inversion parameter uncertainty. This objective function does not necessarily yield experiments that mediate well-determined systems, but, being a Bayesian quality measure, it rigorously accounts for prior information which constrains model parameters that may be only weakly constrained by the optimized dataset. Historically, use of the expected relative entropy has been limited by the computing and storage requirements associated with high-dimensional numerical integration. Herein, a bifocal algorithm is developed that makes these computations more efficient. The algorithm is demonstrated on a medium-sized problem of sampling relaxation phenomena and on a large problem of source–receiver selection for a 2D vertical seismic profile. The method is memory intensive but workarounds are discussed.

Dynamic Factor Models for Multivariate Count Data: An Application to Stock-Market Trading Activity

We propose a dynamic factor model for the analysis of multivariate time series count data. Our model allows for idiosyncratic as well as common serially correlated latent factors in order to account for potentially complex dynamic interdependence between series of counts. The model is estimated under alternative count distributions (Poisson and negative binomial). Maximum likelihood estimation requires high-dimensional numerical integration in order to marginalize the joint distribution with respect to the unobserved dynamic factors. We rely upon the Monte Carlo integration procedure known as efficient importance sampling, which produces fast and numerically accurate estimates of the likelihood function. The model is applied to time series data consisting of numbers of trades in 5-min intervals for five New York Stock Exchange (NYSE) stocks from two industrial sectors. The estimated model provides a good parsimonious representation of the contemporaneous correlation across the individual stocks and their serial correlation. It also provides strong evidence of a common factor, which we interpret as reflecting market-wide news.

Dynamic Factor Models for Multivariate Count Data: An Application to Stock-Market Trading Activity

We propose a dynamic factor model for the analysis of multivariate time series count data. Our model allows for idiosyncratic as well as common serially correlated latent factors in order to account for potentially complex dynamic interdependence between series of counts. The model is estimated under alternative count distributions (Poisson and negative binomial). Maximum Likelihood estimation requires high-dimensional numerical integration in order to marginalize the joint distribution with respect to the unobserved dynamic factors. We rely upon the Monte-Carlo integration procedure known as Efficient Importance Sampling which produces fast and numerically accurate stimates of the likelihood function. The model is applied to time series data consisting of numbers of trades in 5 inutes intervals for five NYSE stocks from two industrial sectors. The estimated model accounts for all key dynamic and distributional features of the data. We find strong evidence of a common factor which we interpret as reflecting market-wide news. In contrast, sector-specific factors are found to be statistically insignificant.

Quasi-Monte Carlo for Highly Structured Generalised Response Models

Highly structured generalised response models, such as generalised linear mixed models and generalised linear models for time series regression, have become an indispensable vehicle for data analysis and inference in many areas of application. However, their use in practice is hindered by high-dimensional intractable integrals. Quasi-Monte Carlo (QMC) is a dynamic research area in the general problem of high-dimensional numerical integration, although its potential for statistical applications is yet to be fully explored. We survey recent research in QMC, particularly lattice rules, and report on its application to highly structured generalised response models. New challenges for QMC are identified and new methodologies are developed. QMC methods are seen to provide significant improvements compared with ordinary Monte Carlo methods.

Cubature formulas for function spaces with moderate smoothness

We construct simple algorithms for high-dimensional numerical integration of function classes with moderate smoothness. These classes consist of square-integrable functions over the d-dimensional unit cube whose coefficients with respect to certain multiwavelet expansions decay rapidly. Such a class contains discontinuous functions on the one hand and, for the right choice of parameters, the quite natural d-fold tensor product of a Sobolev space H s [ 0 , 1 ] on the other hand. The algorithms are based on one-dimensional quadrature rules appropriate for the integration of the particular wavelets under consideration and on Smolyak's construction. We provide upper bounds for the worst-case error of our cubature rule in terms of the number of function calls. We additionally prove lower bounds showing that our method is optimal in dimension d = 1 and almost optimal (up to logarithmic factors) in higher dimensions. We perform numerical tests which allow the comparison with other cubature methods.