smashGP: Large-Scale Spatial Modeling via Matrix-Free Gaussian Processes

Abstract

Gaussian processes are essential for spatial data analysis. Not only do they allow the prediction of unknown values, but they also allow for uncertainty quantification. However, in the era of big data, directly using Gaussian processes has become computationally infeasible as cubic run times are required for dense matrix decomposition and inversion. Various alternatives have been proposed to reduce the computational burden of directly fitting Gaussian processes. These alternatives rely on assumptions on the underlying structure of the covariance or precision matrices, such as sparsity or low-rank. In contrast, this article uses hierarchical matrices and matrix-free methods to enable the computation of Gaussian processes for large spatial datasets by exploiting the underlying kernel properties. The proposed framework, smashGP, represents the covariance matrix as an ${\mathcal{H}}^2$ matrix in $\mathcal{O}(n)$ time and is able to estimate the unknown parameters of the model and predict the values of spatial observations at unobserved locations in $\mathcal{O}(n\hspace{0.17em}\mathrm{\text{log}}\hspace{0.17em}n)$ time thanks to fast matrix-vector products. Additionally, it can be parallelized to take full advantage of shared-memory computing environments. With simulations and case studies, we illustrate the advantage of smashGP to model large-scale spatial datasets. Supplementary materials for this article are available online.

Keywords:

• Gaussian processes
• Hierarchical matrices
• Matrix-free methods
• Spatial data analysis

Supplementary Materials

Code: Code used to run smashGP can be found at https://gitlab.com/libsmash_public/smashgp. See the Readme file for detailed instructions.

Supplementary Materials: Document containing: (A) derivations for GPs with nonconstant mean; (B) additional details on SMASH and its matrix-free operations; (C) definitions of the evaluation metrics used in the simulations and case studies; (D) additional details on the computational complexity of smashGP; (E) parameter estimation evaluation via simulations; (F) comparison with state-of-art methods for large datasets; and (G) case study for a million data points. (smashGP_supplementary.pdf,.pdf file)

Acknowledgments

The authors thank the editor, associate editor, and two anonymous reviewers for their constructive comments and suggestions that have considerably improved the article.

Disclosure Statement

The authors report there are no competing interests to declare.

Additional information

Funding

This work was partially supported by NSF grants CMMI-1839591 and OAC-2003683.