[Paper Review] Conjugate Nearest Neighbor Gaussian Process Models for Efficient Statistical Interpolation of Large Spatial Data

Shinichiro Shirota, Andrew O. Finley+2 more|arXiv (Cornell University)|Jul 23, 2019

Soil Geostatistics and Mapping29 references4 citations

TL;DR

This paper proposes conjugate Bayesian nearest-neighbor Gaussian process models that combine low-rank Gaussian predictive processes with sparsity-inducing nearest-neighbor Gaussian processes to enable exact, computationally efficient spatial interpolation of massive spatial datasets. The method achieves full Bayesian kriging inference for 17 million LiDAR observations in under a minute, enabling high-precision uncertainty quantification for large-scale remote sensing applications like carbon monitoring in Interior Alaska.

ABSTRACT

A key challenge in spatial statistics is the analysis for massive spatially-referenced data sets. Such analyses often proceed from Gaussian process specifications that can produce rich and robust inference, but involve dense covariance matrices that lack computationally exploitable structures. The matrix computations required for fitting such models involve floating point operations in cubic order of the number of spatial locations and dynamic memory storage in quadratic order. Recent developments in spatial statistics offer a variety of massively scalable approaches. Bayesian inference and hierarchical models, in particular, have gained popularity due to their richness and flexibility in accommodating spatial processes. Our current contribution is to provide computationally efficient exact algorithms for spatial interpolation of massive data sets using scalable spatial processes. We combine low-rank Gaussian processes with efficient sparse approximations. Following recent work by [1], we model the low-rank process using a Gaussian predictive process (GPP) and the residual process as a sparsity-inducing nearest-neighbor Gaussian process (NNGP). A key contribution here is to implement these models using exact conjugate Bayesian modeling to avoid expensive iterative algorithms. Through the simulation studies, we evaluate performance of the proposed approach and the robustness of our models, especially for long range prediction. We implement our approaches for remotely sensed light detection and ranging (LiDAR) data collected over the US Forest Service Tanana Inventory Unit (TIU) in a remote portion of Interior Alaska.

Motivation & Objective

• To address the computational intractability of full-likelihood Gaussian process models on massive spatial datasets with tens of millions of locations.
• To develop a scalable, exact inference framework that avoids expensive MCMC sampling while preserving full uncertainty quantification.
• To enable high-precision spatial prediction and uncertainty estimation for large-scale remote sensing applications, such as forest biomass and carbon monitoring.
• To demonstrate the effectiveness of combining low-rank GPP and sparse NNGP components in a conjugate Bayesian framework for full-scale spatial modeling.
• To provide a practical, computationally efficient alternative to MCMC-based inference for spatial processes in large domains like the US Forest Service Tanana Inventory Unit.

Proposed method

• The model decomposes the spatial process into a low-rank Gaussian predictive process (GPP) for long-range dependence and a sparsity-inducing nearest-neighbor Gaussian process (NNGP) for fine-scale residuals.
• The joint model is formulated as a sparse plus low-rank Gaussian process (SLGP), enabling exact conjugate Bayesian inference by avoiding iterative MCMC algorithms.
• Parameter estimation fixes key covariance parameters (e.g., spatial decay φ and noise-to-signal ratio α) at reasonable values, reducing computational burden without sacrificing predictive accuracy.
• The method leverages parallelized computation and thrifty memory management to scale to datasets with up to 17 million spatial locations.
• Posterior predictive distributions are derived exactly using multivariate normal theory, enabling fast prediction and uncertainty quantification at arbitrary locations.
• Model fitting uses K-fold cross-validation to tune hyperparameters and validate performance across large spatial domains.

Experimental results

Research questions

• RQ1 Can conjugate Bayesian inference be effectively applied to large-scale spatial models to avoid MCMC while maintaining predictive accuracy?
• RQ2 How does the combination of low-rank GPP and sparse NNGP components improve computational efficiency and predictive performance on massive spatial datasets?
• RQ3 To what extent does fixing key covariance parameters (e.g., φ and α) impact predictive accuracy in large-scale spatial interpolation?
• RQ4 Can the proposed SLGP model achieve sub-minute inference times for datasets with 17 million spatial locations while maintaining uncertainty quantification?
• RQ5 How well do the models generalize to long-range spatial predictions in non-stationary regions like Interior Alaska?

Key findings

• The SLGP model achieved full Bayesian kriging inference for 17,357,816 LiDAR observations in approximately 22 minutes using conjugate inference, compared to prohibitive MCMC runtimes.
• The NNGP model completed in just 9 seconds, demonstrating the computational advantage of the conjugate approach for large-scale spatial data.
• Both NNGP and SLGP models produced indistinguishable predictions and uncertainty estimates, with CRPS of 0.38 and RMSPE of 0.69 on the TIU dataset.
• The models successfully captured spatial patterns in canopy height and variance, with non-zero coefficients for tree cover (β_TC = 0.12) and fire occurrence (β_Fire = -0.03), indicating meaningful predictor effects.
• The spatial range parameter φ was estimated at 0.6 (equivalent to ~5 km), and the noise-to-signal ratio α was fixed at 0.13, which supported stable and accurate predictions.
• The results confirm that fixing key parameters like φ and α does not degrade predictive performance, validating the use of conjugate models for massive data.