Adaptive Inducing-Point Selection in Sparse GPs

• The paper demonstrates a variance-driven algorithm that adds inducing points where integrated conditional-prior variance is maximally reduced.
• Variance-based methods decompose predictive uncertainty into conditional-prior and posterior-induced components, guiding precise model refinement.
• Empirical results show up to 90% reduction in conditional-prior variance and notable improvements in IMSE, validating adaptive inducing-point selection.

Variance-based adaptive inducing-point selection comprises a suite of algorithmic strategies for constructing and refining the set of inducing variables within sparse Gaussian process (GP) models, using explicit decompositions or estimates of predictive variance to guide where, when, and how additional inducing points should be introduced. These methods aim to balance model fidelity and computational tractability by allocating inducing points adaptively to regions with high residual uncertainty, as quantified by explicit variance components in the GP posterior or variational approximation. Recent developments integrate these variance-driven mechanisms with sequential data acquisition, yielding self-tuning procedures capable of both precision placement and automated stopping.

1. Fundamentals of Sparse Gaussian Processes and Variance Decomposition

Sparse GP models introduce a reduced set of inducing variables, indexed by inducing inputs $Z = [z_1,\dots,z_M]$, to approximate the full GP posterior at a computational cost that is sublinear in the number of training data $N$. In quantile regression under a GP prior, as described by Nicolas & Le Maître (2024), the Laplace posterior over the inducing variables $u = q_\tau(Z)$ is approximated as $p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$ (Nicolas et al., 30 Jun 2026).

The approximate predictive variance at any test location $x$ is rigorously decomposed as

$V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$

where:

• $\Sigma_z^2(x) = k(x, x) - K_{xZ} Z^{-1} K_{Zx}$ is the conditional-prior variance, quantifying uncertainty due to the finite inducing set $Z$ alone,
• $$\Sigma_u^2(x) = K_{xZ} Z^{-1} \hat{C} Z^{-1} K_{Zx}$$ is the posterior-induced variance, measuring uncertainty originating from the finite data $y$ as propagated through $N$0.

This decomposition underpins variance-based adaptive schemes by separately identifying model uncertainty reducible by refining the inducing set versus that reducible only by acquiring more data.

2. Variance-Driven Inducing-Point Infilling Algorithms

Variance-based infilling algorithms add inducing points to locations that most effectively reduce integrated conditional-prior variance. The canonical objective, as formalized in "Sequential sparse Gaussian process quantile regression" (Nicolas et al., 30 Jun 2026), is the integrated reduction

$N$1

which, by the block-inverse identity, admits an efficient closed form for the reduction at each $N$2.

A typical infilling step selects $N$3 as the maximizer of $N$4 over candidate locations, evaluated (and optionally optimized) via Monte Carlo integration. This process, termed Information-based Variance Reduction (IVR) by Nicolas & Le Maître, ensures that each inducing point is maximally informative with respect to predictive uncertainty.

In streaming or online settings, as in Galy-Fajou & Opper (2020), additions may instead rely on direct evaluation of the maximum kernel correlation $N$5, with new points added where $N$6 drops below a threshold $N$7, or where the marginal predictive variance $N$8 exceeds a fixed tolerance (Galy-Fajou et al., 2021).

3. Practical Evaluation and Computational Aspects

Efficiently evaluating $N$9 and $u = q_\tau(Z)$0 for many candidates proceeds via precomputed Cholesky or inverse factorizations of the kernel Gram matrix for $u = q_\tau(Z)$1, resulting in $u = q_\tau(Z)$2 compute per candidate. The full posterior-covariance term $u = q_\tau(Z)$3 (from the Laplace Hessian inversion) is updated once per re-training with $u = q_\tau(Z)$4 cost (Nicolas et al., 30 Jun 2026).

A summary of principal computational steps is as follows:

Step	Expression/Operation	Complexity per $u = q_\tau(Z)$5
$u = q_\tau(Z)$6	$u = q_\tau(Z)$7	$u = q_\tau(Z)$8
$u = q_\tau(Z)$9	$p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$0	$p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$1
IVR $p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$2	MC integral over $p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$3 samples	$p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$4 per $p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$5

Adaptive point selection maintains model efficiency by stopping when the marginal benefit of additional points vanishes (i.e., when integrated conditional-prior variance, $p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$6, reaches parity with the integrated posterior-induced variance, $p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$7).

4. Integrated Sequential Enrichment: Data Acquisition and Point Placement

Variance-based selection operates not only on the inducing set but also interleaves with targeted data acquisition. The complementary sources of uncertainty, $p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$8 and $p(u|Z, y, X) \approx \mathcal{N}(\hat{u}, \hat{C})$9, enable sequential algorithms to switch between infilling (inducing-point addition) and acquiring new labels. The selection is governed by a switching rule: if $x$0, infilling is performed; otherwise, sampling data is preferred ($x$1, by default).

When data acquisition is triggered, input locations $x$2 are drawn proportionally to $x$3, prioritizing observation in regions where model uncertainty is attributed to insufficient labeled data (Nicolas et al., 30 Jun 2026).

This unified framework enables model complexity (the number and location of inducing points) and sample complexity (the number and location of new training points) to be allocated adaptively, focusing resources where estimated gains in predictive confidence are greatest.

5. Empirical Performance and Comparative Analysis

Empirical results on benchmark problems (Sabater and Michalewicz synthetic functions) demonstrate rapid, exponential decay of conditional-prior variance and integrated mean-squared error (IMSE) under variance-driven infilling, with up to a $x$4 reduction in $x$5 and $x$6--$x$7 improvements in IMSE over nonadaptive QMC/halton placement, particularly for moderate $x$8 ($x$9 IVR vs $V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$0 QMC for comparable $V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$1) (Nicolas et al., 30 Jun 2026). The intersection of the curves $V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$2 and $V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$3 in experiments occurs near the IMSE plateau, validating the switching criterion.

In data-limited regimes, acquisition guided by $V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$4, rather than uniform sampling, accelerates IMSE reduction by $V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$5--$V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$6 on the Michalewicz task. When training data are abundant, performance plateaus as infilling becomes the rate-limiting step. This confirms that variance-based strategies automatically adjust to the dominant source of remaining model uncertainty.

6. Comparison with Probabilistic and Fixed-Set Selection Approaches

Variance-driven adaptive selection schemes contrast sharply with probabilistic inducing-point selection methods, such as those based on point-process priors (Uhrenholt et al., 2020). In the latter, a latent random set $V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$7 is learned via variational inference, trading off data fit and a smooth cardinality penalty, with uncertainty over $V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$8 (the inducing-point count) controlled by a single hyperparameter.

Variance-based methods, in contrast, typically operate in a greedy or streaming fashion, requiring recomputation of predictive variance and selection for each candidate point, and treating $V(x) = \Sigma_z^2(x) + \Sigma_u^2(x)$9 as a dynamic control variable. These procedures directly exploit model-specific uncertainty quantification without explicit regularization over set size, and can be interpreted as optimizing a geometric coverage (or information reduction) criterion in input space.

Empirical observations suggest that probabilistic point-process methods closely track the performance of greedy variance-based algorithms for the same $\Sigma_z^2(x) = k(x, x) - K_{xZ} Z^{-1} K_{Zx}$0, but differ in convenience, computational requirements, and interpretability of the model selection process (Uhrenholt et al., 2020).

7. Implementation and Practical Considerations

In implementation, initialization with minimal covering sets (Halton or initial data points) ensures numerical stability. Hyperparameters governing the variance thresholds or kernel behavior (e.g., lengthscale) must be adapted or cross-validated to match model coverage to data and predictive tolerance (Galy-Fajou et al., 2021). Efficient maintenance of the kernel matrix via Cholesky and rank-one updates is crucial for scalability as $\Sigma_z^2(x) = k(x, x) - K_{xZ} Z^{-1} K_{Zx}$1 grows.

For streaming or online data, the addition of new inducing points is performed in tandem with expected improvement in coverage, and model retraining is interleaved to synchronize the variational parameters and hyperparameter estimates. The interplay between kernel configuration and the addition criterion leads to adaptive tuning—shrinking the kernel lengthscale or reducing the addition threshold yields finer coverage with increased $\Sigma_z^2(x) = k(x, x) - K_{xZ} Z^{-1} K_{Zx}$2, while more diffuse kernels admit leaner models.

---

The variance-based adaptive inducing-point selection paradigm thus leverages a principled, uncertainty-driven mechanism to dynamically allocate model complexity in sparse GP frameworks. Through explicit variance decomposition, efficient sequential algorithms, and integrated data-driven enrichment, these methods produce models that are empirically and theoretically superior to fixed design or randomized selection, particularly in quantile regression and other GP applications demanding robust uncertainty quantification (Nicolas et al., 30 Jun 2026, Galy-Fajou et al., 2021, Uhrenholt et al., 2020).