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.
Figures & tables
Figure 1: Loss traces on snelson and banana datasets. The N and P suffixes refer to the use of NG updates ( Section 3.1 ) and inducing mean preconditioning ( Section 3.2 ), respectively. R-SVGP (NP) is the only inverse-free variant that matches the performance of W-SVGP and L-SVGP (P), which use Cholesky decompositions.
Figure 2: NLPD/runtime on elevators and kin40k datasets for different choices of M and batch size B=100 . Lines show the mean over 5 seeds; shaded regions indicate ±1 standard error (often smaller than line width).
L
W-SVGP
L-SVGP
R-SVGP
1
-0.38 (0.001)
-0.39 (0.002)
-0.39 (0.002)
2
-1.86 (0.03)
-1.88 (0.02)
-1.88 (0.02)
3
-2.39 (0.03)
-2.64 (0.01)
-2.64 (0.01)
Table 1: NLPD on kin40k for different choices of L , M=500 and batch size B=1000 , after 200k iterations. Entries are mean (standard error) over 5 seeds.
Kernel
W-SVGP
L-SVGP
R-SVGP
RBF
1.96 (0.02)
1.95 (0.01)
1.94 (0.02)
Conv
1.31 (0.02)
1.53 (0.02)
1.53 (0.03)
Table 2: Test error (%) on mnist with RBF and weighted convolutional kernels, M=750 , B=200 . Entries are mean (standard error) over 5 seeds.
Appendix figures & tables1 asset
Supplementary material from the paper’s appendix.
Appendix
Figure S1: NLPD/runtime on elevators and kin40k datasets for different choices of the number of probes K , with M=2000 and B=100 . Lines show the mean over 5 seeds; shaded regions indicate ±1 standard error.
Department of Information Management, National Taiwan University · 2H. Milton Stewart School of Industrial and Systems Engineering, Georgia Institute of Technology