A note on the convergence analysis of a sparse grid multivariate probability density estimator

With the recent growth in volume and complexity of available data has come a renewed interest in the problem of estimating multivariate probability density functions. However, traditional methods encounter the curse of dimensionality (complexity grows exponentially with dimension). Here we provide an outline of a convergence analysis of a sparse grid based probability density estimation, which supports the use of the method for moderately complex (up to 15 dimensions) data sets, as has already been demonstrated for sparse grid quadrature and interpolation. References Hans-Joachim Bungartz and Michael Griebel. A note on the complexity of solving {P}oisson's equation for spaces of bounded mixed derivatives. J. Complexity, 15(2):167--199, 1999. http://dx.doi.org/10.1006/jcom.1999.0499 Hans-Joachim Bungartz and Michael Griebel. Sparse grids. Acta Numerica, 13:147--269, 2004. http://dx.doi.org/10.1017/S0962492904000182 J. Garcke. Sparse grid tutorial. http://www.math.tu-berlin.de/ garcke/paper/sparseGridTutorial.pdf [29 January, 2009], 2008. J. Garcke, M. Griebel, and M. Thess. Data mining with sparse grids. Computing, 67:225--253, 2001. http://dx.doi.org/10.1007/s006070170007 D. Scott. Multivariate Density Estimation: Theory, Practice, and Visualization. Wiley, New York, 1992. C. Zenger. Sparse grids. In W. Hackbusch, editor, Parallel Algorithms for Partial Differential Equations: Proceedings of the Sixth GAMM-Seminar, Kiel, 1990, volume 31 of Notes on Num. Fluid Mech., pages 241--251. Vieweg, 1991.

Summary A sparse grid approach to orbit uncertainty propagation is presented. Efficient and accurate uncertainty propagation methods for nonlinear dynamic systems have been of enormous interest to space object tracking. Recent methods include those based on the time evolution of the probability density function, the statistical moments, the random samples, or a sum of Gaussian components. The idea of the sparse grid method for orbit uncertainty propagation is to represent the initial uncertainty by a sparse grid, propagate the sparse grid points individually through the nonlinear orbit dynamics, and compute the statistical moments from the propagated sparse grid points. The Smolyak rule is used to generate the sparse grid through linear combinations of low-level tensor products; the number of the resultant sparse grid points is polynomial in the dimension of the state. The sparse grid based uncertainty propagation may be considered a generalization of the Unscented Transformation (UT) of the well-known Unscented Kalman Filter but can achieve much higher accuracy levels than UT. Two examples of orbit propagation in different dynamic environments are used to show the efficacy of the sparse grid approach.

Nonparametric density estimation is a fundamental problem of statistics and data mining. Even though kernel density estimation is the most widely used method, its performance highly depends on the choice of the kernel bandwidth, and it can become computationally expensive for large data sets. We present an adaptive sparse-grid-based density estimation method which discretizes the estimated density function on basis functions centered at grid points rather than on kernels centered at the data points. Thus, the costs of evaluating the estimated density function are independent from the number of data points. We give details on how to estimate density functions on sparse grids and develop a cross validation technique for the parameter selection. We show numerical results to confirm that our sparse-grid-based method is well-suited for large data sets, and, finally, employ our method for the classification of astronomical objects to demonstrate that it is competitive to current kernel-based density estimation approaches with respect to classification accuracy and runtime.

Multivariate density estimation is a fundamental problem in Applied Statistics and Machine Learning. Given a collection of data sampled from an unknown distribution, the task is to approximately reconstruct the generative distribution. There are two different approaches to the problem, the parametric approach and the non-parametric approach. In the parametric approach, the approximate distribution is represented by a model from a predetermined family. In this thesis, we adopt the parametric approach and investigate the use of a model family called latent tree models for the task of density estimation. Latent tree models are tree-structured Bayesian networks in which leaf nodes represent observed variables, while internal nodes represent hidden variables. Such models can represent complex relationships among observed variables, and in the meantime, admit efficient inference among them. Consequently, they are a desirable tool for density estimation. While latent tree models are studied for the first time in this thesis for the purpose of density estimation, they have been investigated earlier for clustering and latent structure discovery. Several algorithms for learning latent tree models have been proposed. The state-of-the-art is an algorithm called EAST. EAST determines model structures through principled and systematic search, and determines model parameters using the EM algorithm. It has been shown to be capable of achieving good trade-off between fit to data and model complexity. It is also capable of discovering latent structures behind data. Unfortunately, it has a high computational complexity, which limits its applicability to density estimation problems. In this thesis, we propose two latent tree model learning algorithms specifically for density estimation. The two algorithms have distinct characteristics and are suitable for different applications. The first algorithm is called HCL. HCL assumes a predetermined bound on model complexity and restricts to binary model structures. It first builds a binary tree structure based on mutual information and then runs the EM algorithm once on the resulting structure to determine the parameters. As such, it is efficient and can deal with large applications. The second algorithm is called Pyramid. Pyramid does not assume predetermined bounds on model complexity and does not restrict to binary tree structures. It builds model structures using heuristics based on mutual information and local search. It is slower than HCL. However, it is faster than EAST and is only slightly inferior to EAST in terms of the quality of the resulting models. In this thesis, we also study two applications of the density estimation techniques that we develop. The first application is to approximate probabilistic inference in Bayesian networks. A Bayesian network represents a joint distribution over a set of random variables. It often happens that the network structure is very complex and making inference directly on the network is computational intractable. We propose to approximate the joint distribution using a latent tree model and exploit the latent tree model for faster inference. The idea is to sample data from the Bayesian network, learn a latent tree model from the data offline, and when online, make inference with the latent tree model instead of the original Bayesian network. HCL is used here because the sample size needs to be large to produce accurate approximation and it is possible to predetermine a bound on the online running. Empirical evidence shows that this method can achieve good approximation accuracy at low online computational cost. The second application is classification. A common approach to this task is to formulate it as a density estimation problem: One constructs the class-conditional density for each class and then uses the Bayes rule to make classification. We propose to estimate those class-conditional densities using either EAST or Pyramid. Empiricalevidence shows that this method yields good classification performances. Moreover, the latent tree models built for the class-conditional densities are often meaningful, which is conducive to user confidence. A comparison between EAST and Pyramid reveals that Pyramid is significantly more efficient than EAST, while it results in more or less the same classification performance as the latter.

Nonparametric estimation for a probability density function that describes multivariate data has typically been addressed by kernel density estimation (KDE). A novel density estimator recently developed by Farmer and Jacobs offers an alternative high-throughput automated approach to univariate nonparametric density estimation based on maximum entropy and order statistics, improving accuracy over univariate KDE. This article presents an extension of the single variable case to multiple variables. The univariate estimator is used to recursively calculate a product array of one-dimensional conditional probabilities. In combination with interpolation methods, a complete joint probability density estimate is generated for multiple variables. Good accuracy and speed performance in synthetic data are demonstrated by a numerical study using known distributions over a range of sample sizes from 100 to 106 for two to six variables. Performance in terms of speed and accuracy is compared to KDE. The multivariate density estimate developed here tends to perform better as the number of samples and/or variables increases. As an example application, measurements are analyzed over five filters of photometric data from the Sloan Digital Sky Survey Data Release 17. The multivariate estimation is used to form the basis for a binary classifier that distinguishes quasars from galaxies and stars with up to 94% accuracy.

In simulation technology, computationally expensive objective functions are often replaced by cheap surrogates, which can be obtained by interpolation. Full grid interpolation methods suffer from the so-called curse of dimensionality, rendering them infeasible if the parameter domain of the function is higher-dimensional (four or more parameters). Sparse grids constitute a discretization method that drastically eases the curse, while the approximation quality deteriorates only insignificantly. However, conventional basis functions such as piecewise linear functions are not smooth (continuously differentiable). Hence, these basis functions are unsuitable for applications in which gradients are required. One example for such an application is gradient-based optimization, in which the availability of gradients greatly improves the speed of convergence and the accuracy of the results. This thesis demonstrates that hierarchical B-splines on sparse grids are well-suited for obtaining smooth interpolants for higher dimensionalities. The thesis is organized in two main parts: In the first part, we derive new B-spline bases on sparse grids and study their implications on theory and algorithms. In the second part, we consider three real-world applications in optimization: topology optimization, biomechanical continuum-mechanics, and dynamic portfolio choice models in finance. The results reveal that the optimization problems of these applications can be solved accurately and efficiently with hierarchical B-splines on sparse grids.

To recover or approximate smooth multivariate functions, sparse grids are superior to full grids due to a significant reduction of the required support nodes. The order of the convergence rate in the maximum norm is preserved up to a logarithmic factor. We describe three possible piecewise multilinear hierarchical interpolation schemes in detail and conduct a numerical comparison. Furthermore, we document the features of our sparse grid interpolation software package spinterp for MATLAB.

This model-based design of experiments (MBDOE) method determines the input magnitudes of an experimental stimuli to apply and the associated measurements that should be taken to optimally constrain the uncertain dynamics of a biological system under study. The ideal global solution for this experiment design problem is generally computationally intractable because of parametric uncertainties in the mathematical model of the biological system. Others have addressed this issue by limiting the solution to a local estimate of the model parameters. Here we present an approach that is independent of the local parameter constraint. This approach is made computationally efficient and tractable by the use of: (1) sparse grid interpolation that approximates the biological system dynamics, (2) representative parameters that uniformly represent the data-consistent dynamical space, and (3) probability weights of the represented experimentally distinguishable dynamics. Our approach identifies data-consistent representative parameters using sparse grid interpolants, constructs the optimal input sequence from a greedy search, and defines the associated optimal measurements using a scenario tree. We explore the optimality of this MBDOE algorithm using a 3-dimensional Hes1 model and a 19-dimensional T-cell receptor model. The 19-dimensional T-cell model also demonstrates the MBDOE algorithm’s scalability to higher dimensions. In both cases, the dynamical uncertainty region that bounds the trajectories of the target system states were reduced by as much as 86% and 99% respectively after completing the designed experiments in silico. Our results suggest that for resolving dynamical uncertainty, the ability to design an input sequence paired with its associated measurements is particularly important when limited by the number of measurements.

B-spline quasi-interpolation on sparse grids

The integrative analysis of copy number alteration (CNA) and gene expression (GE) is an essential part of cancer research considering the impact of CNAs on cancer progression and prognosis. In this research, an integrative analysis was performed with generalized differentially coexpressed gene sets (gdCoxS), which is a modification of dCoxS. In gdCoxS, set-wise interaction is measured using the correlation of sample-wise distances with Renyi’s relative entropy, which requires an estimation of sample density based on omics profiles. To capture correlations between the variables, multivariate density estimation with covariance was applied. In the simulation study, the power of gdCoxS outperformed dCoxS that did not use the correlations in the density estimation explicitly. In the analysis of the lower-grade glioma of the cancer genome atlas program (TCGA-LGG) data, the gdCoxS identified 577 pathway CNAs and GEs pairs that showed significant changes of interaction between the survival and non-survival group, while other benchmark methods detected lower numbers of such pathways. The biological implications of the significant pathways were well consistent with previous reports of the TCGA-LGG. Taken together, the gdCoxS is a useful method for an integrative analysis of CNAs and GEs.

Problems with uncertainties arise in many practical fields and traditionally are formulated as dynamic systems with interval parameters. Often the complexity of existing methods is exponential in relation to the number of interval parameters. The adaptive interpolation algorithm and approaches directed to reducing the curse of dimensionality are considered. The main assumption on which these approaches are based is that not all interval parameters make a significant contribution to the solution of the problem. The use of tensor train decomposition and sparse grids allows us to take into account these features and expand the scope of the algorithm for the case of a large number of interval parameters. The effectiveness of the considered approaches is confirmed on several model problems.

Sparse grids have already been successfully used in various high-performance computing (HPC) applications, including data mining. In this article, we take a legacy classification kernel previously optimized for the AVX2 instruction set and investigate the benefits of using the newer AVX-512-based multi-and many-core architectures. In particular, the Knights Landing (KNL) processor is used to study the possible performance gains of the code. Not all kernels benefit equally from such architectures, therefore choices in optimization steps and KNL cluster and memory modes need to be filtered through the lens of the code implementation at hand. With a less traditional approach of manual vectorization through instruction-level intrinsics, our kernel provides a differently faceted look into the optimization process. Observations stem from results obtained for node-and cluster-level classification simulations with up to 2 <sup xmlns:mml="http://www.w3.org/1998/Math/MathML" xmlns:xlink="http://www.w3.org/1999/xlink">^</sup> 28 multidimensional training data points, using the CooLMUC-3cluster of the Leibniz Supercomputing Center (LRZ) in Garching, Germany.

Multivariate density estimation (MDE) suggests that mortgage-backed security (MBS) prices can be well described as a function of the level and slope of the term structure. We analyze bow this function varies across MBSs with different coupons. An important finding is that the interest rate level proxies for the moneyness of the option, the expected level of prepayments, and the average life of the cash flows, while the term structure slope controls for the average rate at which these cash flows should be discounted. Though the origination and prepayment behavior of mortgages differ substantially across coupons, there remains an unexplained common factor in MBS prices. This factor does not seem to be related to the usual suspects and therefore presents a puzzle to financial economists. Coauthors are Robert F. Whitelaw, Matthew Richardson, and Richard Stanton. Article published by Oxford University Press on behalf of the Society for Financial Studies in its journal, The Review of Financial Studies.

We present several multivariate histogram density estimates that are universallyL 1-optimal to within a constant factor and an additive term\(O\left( {\sqrt {\log {n \mathord{\left/ {\vphantom {n n}} \right. \kern-\nulldelimiterspace} n}} } \right)\). The bin widths are chosen by the combinatorial method developed by the authors inCombinatorial Methods in Density Estimation (Springer-Verlag, 2001). The present paper solves a problem left open in that book.

: The 24th Symposium on the Interface: Computing Science and Statistics was held March 18-21, 1992 in College Station, Texas. The theme was Graphics and Visualization. Sixty invited papers were presented to more than 300 attendees. James Blinn gave the Keynote Address while Edward Tufte gave the Plenary Address. There were 38 technical sessions, including a poster session. There was also a short course Visualization. This report presents papers based on the presentations. Also included is a list of attendees.... Exploratory data analysis, Computational statistics, Multivariate density estimation, Scientific visualization, Nonparametric regression.