# Scalable Gaussian Process Computations Using Hierarchical Matrices

arXiv: https://arxiv.org/abs/1808.03215 · HTML: https://mlsanity.com/paper/1808.03215/

## Key facts

- Cited by 8 papers in the map (data through 14 Oct 2026).
- First posted on arXiv 9 Aug 2018.

## Built on by (most cited first, 8 of 8)

- [Sparse Cholesky Factorization for Solving Nonlinear PDEs via Gaussian Processes](https://mlsanity.com/paper/2304.01294/) (2023)
- [Generative Modeling via Hierarchical Tensor Sketching](https://mlsanity.com/paper/2304.05305/) (2023)
- [A Unifying Perspective on Non-Stationary Kernels for Deeper Gaussian Processes](https://mlsanity.com/paper/2309.10068/) (2023)
- [Iterative Methods for Full-Scale Gaussian Process Approximations for Large Spatial Data](https://mlsanity.com/paper/2405.14492/) (2024)
- [Non-stationary Gaussian process discriminant analysis with variable selection for high-dimensional functional data](https://mlsanity.com/paper/2109.14171/) (2021)
- [Gaussian Process Regression under Computational and Epistemic Misspecification](https://mlsanity.com/paper/2312.09225/) (2023)
- [Variance Reduction for Matrix Computations with Applications to Gaussian Processes](https://mlsanity.com/paper/2106.14565/) (2021)
- [Deterministic Bounds and Random Estimates of Metric Tensors on Neuromanifolds](https://mlsanity.com/paper/2505.13614/) (2025)

## Abstract (arXiv)

We present a kernel-independent method that applies hierarchical matrices to the problem of maximum likelihood estimation for Gaussian processes. The proposed approximation provides natural and scalable stochastic estimators for its gradient and Hessian, as well as the expected Fisher information matrix, that are computable in quasilinear O(n log2 n) complexity for a large range of models. To accomplish this, we (i) choose a specific hierarchical approximation for covariance matrices that enables the computation of their exact derivatives and (ii) use a stabilized form of the Hutchinson stochastic trace estimator. Since both the observed and expected information matrices can be computed in quasilinear complexity, covariance matrices for MLEs can also be estimated efficiently. After discussing the associated mathematics, we demonstrate the scalability of the method, discuss details of its implementation, and validate that the resulting MLEs and confidence intervals based on the inverse Fisher information matrix faithfully approach those obtained by the exact likelihood.