Gaussian processes (GPs) provide principled probabilistic predictions while encoding prior knowledge, including equivariances. Yet, their use in large-scale scientific problems is limited by computational cost. Equivariant neural networks are common but typically lack the uncertainty quantification offered by GPs, which is valuable in applications such as molecular research. High-dimensional inputs and large symmetry groups further demand scalability. We establish results pertaining to the interplay of GP equivariance and conditioning and leverage them to obtain equivariant sparse GPs through suitable mean functions and covariance kernels. We instantiate this framework with a flexible class of integration-free equivariant kernels, yielding scalable and data-efficient GP inference. In particular, we introduce triply scalable equivariant Gaussian processes. We employ equivariant sparse variational Gaussian processes for SO(2)-equivariant vector fields and molecular property prediction. Alongside the SVGP, we develop a matrix-free equivariant full-GP implementation that combines an exact Kronecker reduction with preconditioned conjugate-gradient solves, enabling fast and scalable evaluation of the full joint predictive density. We further compare different approaches for selecting inducing points in the equivariant sparse GP models. Our test cases include synthetic SO(2)-equivariant fields as well as the prediction of electric dipole moments of N-methylformamide based on quantum chemistry simulations, achieving accurate, uncertainty-aware predictions at a fraction of the computational cost of classical GP inference.
Figures & tables
Figure 1: SO(2) -equivariant GP regression with the integration-free equivariant kernel KΠ ; see Experiment 3.3.2 . Left: a stochastically equivariant GP realization serving as ground truth, split into 40% training and 60% test data. Middle and right: posterior means from a full GP and a sparse variational GP (SVGP). The SVGP uses inducing points in the fundamental region and achieves substantial speedups while preserving equivariance.
Symbol
Meaning
D⊂Rd
input domain
x,x′
input locations in D
f=(f(x))x∈D
Rp -valued Gaussian process
f(x)∈Rp
(GP) response vector at input x
f(r)(x)
r th component of f(x)
K(x,x′)∈Rp×p
matrix-valued covariance kernel
Table 1: Main notation used throughout the paper.
Figure 2: Synthetic SO(2) experiment with n=10000 . Top row: realized equivariant field, full GP posterior mean error with respect to the realization, and SVGP posterior mean error. Bottom row: predictive uncertainty of the full GP and SVGP models.
Figure 3: Wedge equivariance check. Top left: training observations in S and rotated copy S′ . Top right: visualizes 21logdet(ΣSVGP)−21logdet(ΣFull) on test points in S and S′ . Bottom: full GP posterior mean on S and S′ (left) and SVGP posterior mean on S and S′ with inducing points in the fundamental region indicated by green triangles.
d(x,x′)2=tr(K(x,x))+tr(K(x′,x′))−2tr(K(x,x′)),
Algorithm 1 Greedy inducing point selection
Figure 4: RMSE and LogS learning curves comparing SVGP variants, PCG-based matrix-free full-GP inference, and dense full-GP baselines across increasing training-set sizes.
Figure 5: Learning curves for training sizes n≤2000 . Shaded areas represent the range of the mean score plus or minus one standard deviation over all inducing-point allocation schemes for SVGP inducing-set ratios r∈{0.1,0.2,0.3} .
Figure 6: Wall-clock computation times and peak resident memory for the different GP models.
Appendix figures & tables3 assets
Supplementary material from the paper’s appendix.
Appendix
Figure 7: Boxplots comparing SVGP predictive performance at training size n=20000 for several minibatch ratios and inducing-point selection strategies.
Figure 8: Extended SVGP learning curves in the small-data regime, showing three rows of metrics across inducing-point allocation methods and inducing-set ratios.
Figure 9: Further prior samples of GP(0,KΠ) , with hyperparameters sampled uniformly over [0,20]4 .
Gaussian processes (GPs) offer appealing properties but are costly to train at scale. Sparse variational GP (SVGP) approximations reduce cost yet still rely on Cholesky decompositions of kernel matrices, ill-suited to low-precision, massively parallel hardware. While one can construct valid variational bounds that rely only on matrix multiplications (matmuls) via an auxiliary matrix parameter, optimising them with off-the-shelf first-order methods is challenging. We make the inverse-free approach practical by proposing a better-conditioned bound and deriving a matmul-only natural-gradient update for the auxiliary parameter, markedly improving stability and convergence. We further provide simple heuristics, such as step-size schedules and stopping criteria, that make the overall optimisation routine fit seamlessly into existing workflows. Across regression and classification benchmarks, we demonstrate that our method 1) serves as a drop-in replacement in SVGP-based models (e.g., deep GPs), 2) recovers similar performance to traditional methods, and 3) can be faster than baselines when well tuned.
Gaussian processes constitute a cornerstone of probabilistic machine learning, yet scaling them to large datasets typically forces a trade-off between computational efficiency and model fidelity. This work bridges this gap by presenting an adaptive multi-resolution Gaussian process framework that is both scalable and exact. Our key innovation is constructing a naturally data-sparse covariance matrix with adaptive multi-resolution basis functions. These basis functions are directly anchored to samples, eliminating the need for auxiliary points. By shrinking the support domains of multi-resolution basis, the matrix block sizes are limited, guaranteeing sparsity. The inverse of the data-sparse covariance matrix is computed exactly and efficiently via the sparse Cholesky inverse algorithm. To further improve predictive uncertainties, we construct an augmented basis function. Theoretical analysis and numerical experiments demonstrate that our model achieves exact inference with O(nlog2n) training cost and O(logdn) prediction cost, establishing a principled framework for scalable and high-fidelity Gaussian process regression.
Yanchuang Cao, Jun Liu, Tengchao Yu +1
Institute of Applied Physics and Computational Mathematics, Beijing 100094, China
Gaussian processes (GPs) provide a principled Bayesian framework for uncertainty estimation, but their computational complexity severely limits scalability to large datasets. We propose SIKA-GP, which accelerates GP inference using sparse inducing kernel approximations based on a dyadic ordered template basis, incurring only O(logM) complexity dependence on the number of inducing points. Our approach constructs compact and expressive kernel representations from sparsely activated bases, enabling efficient tensorized GPU computation and seamless integration with modern large-scale models. SIKA-GP can be naturally embedded into Bayesian neural networks (BNNs) with sparse activations, yielding significant speedups in both training and inference without sacrificing predictive performance. The method naturally extends to deep feature learning, addressing the scalability challenges introduced by deep architectures and high-dimensional feature representations. Empirical results on vision and transformer-based language benchmarks demonstrate that our approach consistently delivers fast and accurate GP models, providing a principled path toward scalable kernel learning.
Wenyuan Zhao, Rui Tuo, Chao Tian
Department of Electrical and Computer Engineering, Texas A&M University, College Station, US · Department of Industrial and Systems Engineering, Texas A&M University, College Station, US