Compression, inversion, and approximate PCA of dense kernel matrices at near-linear computational complexity

DOI10.1137/19M129526XMaRDI QIDQ92247zbMATH OpenOpenAlexFDO

Authors Florian Schäfer, T. J. Sullivan, Houman Owhadi, Florian Schäfer, T. J. Sullivan, Houman Owhadi

Publication date 7 June 2017

Published in Multiscale Modeling & Simulation (Search for Journal in Brave)

Full work available at URL https://arxiv.org/abs/1706.02205

zbMATH Keywords

principal component analysisCholesky factorizationcovariance functiongamblet transformkernel matrixsparsity

Mathematics Subject Classification ID

Factor analysis and principal components; correspondence analysis (62H25) Nontrigonometric harmonic analysis involving wavelets and other special systems (42C40) Computational methods for sparse matrices (65F50) Martingales with discrete parameter (60G42) Numerical linear algebra (65F99) Multigrid methods; domain decomposition for boundary value problems involving PDEs (65N55) Probabilistic methods, particle methods, etc. for boundary value problems involving PDEs (65N75) Analysis of algorithms and problem complexity (68Q25) Analysis of algorithms (68W40)

Abstract: Dense kernel matrices $Thetainmathbb{R}^{NimesN}$ obtained from point evaluations of a covariance function $G$ at locations ${x_i}_{1leqileqN}subsetmathbb{R}^d$ arise in statistics, machine learning, and numerical analysis. For covariance functions that are Green's functions of elliptic boundary value problems and homogeneously-distributed sampling points, we show how to identify a subset $Ssubset{1,dots,N}^2$, with $$, such that the zero fill-in incomplete Cholesky factorisation of the sparse matrix $Thet{a}_{ij}{1}_{(i,j)inS}$ is an $epsilon$-approximation of $Theta$. This factorisation can provably be obtained in complexity $O(Nlog(N)lo{g}^d(N/epsilon))$ in space and $O(Nlo{g}^2(N)lo{g}^{2d}(N/epsilon))$ in time, improving upon the state of the art for general elliptic operators; we further present numerical evidence that $d$ can be taken to be the intrinsic dimension of the data set rather than that of the ambient space. The algorithm only needs to know the spatial configuration of the $x_i$ and does not require an analytic representation of $G$. Furthermore, this factorization straightforwardly provides an approximate sparse PCA with optimal rate of convergence in the operator norm. Hence, by using only subsampling and the incomplete Cholesky factorization, we obtain, at nearly linear complexity, the compression, inversion and approximate PCA of a large class of covariance matrices. By inverting the order of the Cholesky factorization we also obtain a solver for elliptic PDE with complexity $O(Nlo{g}^d(N/epsilon))$ in space and $O(Nlo{g}^{2d}(N/epsilon))$ in time, improving upon the state of the art for general elliptic operators.

Recommendations

• On the Nyström method for approximating a gram matrix for improved kernel-based learning
• Learning Theory
• Low-Rank Eigenvector Compression of Posterior Covariance Matrices for Linear Gaussian Inverse Problems
• An \(O(N \log N)\) hierarchical random compression method for kernel matrices by sampling partial matrix entries
• Approximation of high-dimensional kernel matrices by multilevel circulant matrices

Cites work

• scientific article; zbMATH DE number 3958638 (Why is no real title available?)
• scientific article; zbMATH DE number 3576661 (Why is no real title available?)
• scientific article; zbMATH DE number 3581563 (Why is no real title available?)
• scientific article; zbMATH DE number 3633702 (Why is no real title available?)
• scientific article; zbMATH DE number 635657 (Why is no real title available?)
• scientific article; zbMATH DE number 2107836 (Why is no real title available?)
• scientific article; zbMATH DE number 3210432 (Why is no real title available?)
• scientific article; zbMATH DE number 3273551 (Why is no real title available?)
• 10.1162/15324430260185619
• 10.1162/153244303768966085
• 2010 Rietz lecture: When does the screening effect hold?
• A Fast $ULV$ Decomposition Solver for Hierarchically Semiseparable Representations
• A Full Scale Approximation of Covariance Functions for Large Spatial Data Sets
• A class of multi-resolution approximations for large spatial datasets
• A fast algorithm for particle simulations
• A sparse \({\mathcal H}\)-matrix arithmetic. II: Application to multi-dimensional problems
• A sparse matrix arithmetic based on \({\mathfrak H}\)-matrices. I: Introduction to \({\mathfrak H}\)-matrices
• A unifying view of sparse approximate Gaussian process regression
• Adaptive low-rank approximation of collocation matrices
• An Efficient Block-Oriented Approach to Parallel Sparse Cholesky Factorization
• An Iterative Solution Method for Linear Systems of Which the Coefficient Matrix is a Symmetric M-Matrix
• An \(\mathcal O(N\log N)\) fast direct solver for partial hierarchically semi-separable matrices. With application to radial basis function interpolation
• An analysis of a class of variational multiscale methods based on subspace decomposition
• An explicit link between Gaussian fields and Gaussian Markov random fields: the stochastic partial differential equation approach
• Approximation of solution operators of elliptic partial differential equations by \({\mathcal H}\)- and \({\mathcal H}^2\)-matrices
• Block red-black ordering: A new ordering strategy for parallelization of ICCG method
• Compressing Rank-Structured Matrices via Randomized Sampling
• Constructing continuous stationary covariances as limits of the second-order stochastic difference equations
• Correlation priors
• Data-sparse approximation by adaptive \({\mathcal H}^2\)-matrices
• Decay Rates for Inverses of Band Matrices
• Decay bounds for functions of Hermitian matrices with banded or Kronecker structure
• Efficient numerical methods for non-local operators. \(\mathcal H^2\)-matrix compression, algorithms and analysis.
• Existence of \(\mathcal H\)-matrix approximants to the inverse FE-matrix of elliptic operators with \(L^\infty\)-coefficients
• Fast algorithms for hierarchically semiseparable matrices
• Fast wavelet transforms and numerical algorithms I
• Gaussian Predictive Process Models for Large Spatial Data Sets
• Gaussian processes for machine learning.
• Generalized Nested Dissection
• Graph colouring algorithms
• Hierarchical interpolative factorization for elliptic operators: integral equations
• Hierarchical matrices. A means to efficiently solve elliptic boundary value problems
• Kernel methods in machine learning
• LU factorization of non-standard forms and direct multiresolution solvers
• Localization in matrix computations: theory and applications
• Localization of elliptic multiscale problems
• Localization of matrix factorizations
• Mitigating the influence of the boundary on PDE-based covariance operators
• Multigrid with Rough Coefficients and Multiresolution Operator Decomposition from Hierarchical Information Games
• Nested Dissection of a Regular Finite Element Mesh
• New efficient and robust HSS Cholesky factorization of SPD matrices
• Numerical homogenization of heterogeneous fractional Laplacians
• ON STATIONARY PROCESSES IN THE PLANE
• On Finding Supernodes for Sparse Matrix Computations
• Operator-adapted wavelets, fast solvers, and numerical homogenization. From a game theoretic approach to numerical approximation and algorithm design
• Orderings for Factorized Sparse Approximate Inverse Preconditioners
• Polyharmonic homogenization, rough polyharmonic splines and sparse super-localization
• Propriétés des matrices ``bien localisées près de leur diagonale et quelques applications. (Properties of matrices ``well localized near the diagonal and some applications)
• Randomized matrix-free trace and log-determinant estimators
• SelInv---An Algorithm for Selected Inversion of a Sparse Symmetric Matrix
• Sparse Compression of Expected Solution Operators
• Sparse operator compression of higher-order elliptic operators with rough coefficients
• Stochastic Models That Separate Fractal Dimension and the Hurst Effect
• Studies in the history of probability and statistics XLIX: On the Matérn correlation family
• The Bramble--Hilbert Lemma for Convex Domains
• The Evolution of the Minimum Degree Ordering Algorithm
• The Schur complement and its applications
• The analysis of a nested dissection algorithm
• Whittle-Matérn priors for Bayesian statistical inversion with applications in electrical impedance tomography

Cited in

(34)

• Gamblets for opening the complexity-bottleneck of implicit schemes for hyperbolic and parabolic ODEs/PDEs with rough coefficients
• Multi-scale Vecchia approximations of Gaussian processes
• Sum of Kronecker products representation and its Cholesky factorization for spatial covariance matrices from large grids
• GParareal: a time-parallel ODE solver using Gaussian process emulation
• A modern retrospective on probabilistic numerics
• Do ideas have shape? Idea registration as the continuous limit of artificial neural networks
• Scalable Bayesian Transport Maps for High-Dimensional Non-Gaussian Spatial Fields
• Scaled Vecchia approximation for fast computer-model emulation
• An adaptive factorized Nyström preconditioner for regularized kernel matrices
• Hierarchical sparse Cholesky decomposition with applications to high-dimensional spatio-temporal filtering
• Sparse Recovery of Elliptic Solvers from Matrix-Vector Products
• A fast hierarchically preconditioned eigensolver based on multiresolution matrix decomposition
• Bayesian numerical methods for nonlinear partial differential equations
• Fast eigenpairs computation with operator adapted wavelets and hierarchical subspace correction
• Approximation of high-dimensional kernel matrices by multilevel circulant matrices
• Continuous relaxations for the traveling salesman problem
• Gaussian process hydrodynamics
• Coarsening in algebraic multigrid using Gaussian processes
• GPvecchia
• Fast macroscopic forcing method
• Analytical low-rank compression via proxy point selection
• Sparse Cholesky Factorization by Kullback--Leibler Minimization
• Kernel methods are competitive for operator learning
• Symmetry exploits for Bayesian cubature methods
• Error analysis of kernel/GP methods for nonlinear and parametric PDEs
• A Bayesian conjugate gradient method (with discussion)
• Bayesian nonparametric generative modeling of large multivariate non-Gaussian spatial fields
• Numerical homogenization beyond scale separation
• Multiresolution operator decomposition for flow simulation in fractured porous media
• De-noising by thresholding operator adapted wavelets
• Learning elliptic partial differential equations with randomized linear algebra
• Samplets: construction and scattered data compression
• Solving and learning nonlinear PDEs with Gaussian processes
• Bayesian nonstationary and nonparametric covariance estimation for large spatial data (with discussion)

Describes a project that uses

Uses Software

• George
• SelInv
• hlib

Summary This article provides an overview of numerical methods for approximating and solving problems related to Gaussian processes, kernel matrices, and elliptic PDEs. It discusses various approaches for efficiency and accuracy, including sparse linear solvers (e.g., \multigrid, sparse Cholesky factorization) and low-rank approximations (Nyström, rank-revealing Cholesky factorization). Hierarchical methods (hierarchical matrices, HODLR, HSS matrices) and wavelet-based techniques are also highlighted. A notable contribution mentioned involves a modified Cholesky factorization for kernel matrices from elliptic PDE Green's functions, aiming for scalability and accuracy. Overall, the article showcases a range of solutions for computational challenges in these domains.

Summary_simple Scientists use various methods to solve complex math problems related to Gaussian processes and partial differential equations (PDEs). These problems often involve large amounts of data, making them hard to compute. To tackle this, researchers employ techniques like sparse linear solvers, low-rank approximations, hierarchical methods, and wavelet-based approaches. A new method, modifying an existing math tool called Cholesky factorization, shows promise in solving certain PDEs efficiently without needing overly complex structures. This variety of approaches helps scientists find the best way to solve their specific computational challenges.

This page was built for publication: Compression, inversion, and approximate PCA of dense kernel matrices at near-linear computational complexity

Report a bug (only for logged in users!)Click here to report a bug for this page (MaRDI item Q92247)