Sequential sparse Gaussian process quantile regression

Abstract: Quantile regression aims to estimate the conditional quantiles of a response variable from observed data. In a Bayesian setting, Gaussian process quantile regression provides uncertainty quantification but faces significant computational challenges due to the nonconjugacy of the asymmetric Laplace likelihood and the cost of posterior inference. We develop a sparse Gaussian process framework in which the quantile function is represented through a reduced set of inducing variables and posterior inference is performed using a Laplace approximation. A decomposition of the predictive uncertainty into conditional-prior and posterior-induced variance components is then exploited to drive two complementary adaptive mechanisms: inducing-input infilling and data acquisition. These mechanisms are combined within a sequential algorithm that allocates computational effort toward the dominant source of predictive uncertainty and adaptively controls model complexity. Numerical experiments on benchmark problems demonstrate the accuracy of the Laplace approximation, the benefits of variance-based inducing-input placement, and the effectiveness of the proposed sequential enrichment strategy compared with predefined data-acquisition strategies.

• The paper proposes a Laplace approximation method to efficiently perform Bayesian sparse GP quantile regression by reducing the optimization space.
• The paper presents a variance-based adaptive inducing-point selection strategy to minimize conditional-prior variance and control model complexity.
• The paper demonstrates that sequential enrichment through adaptive data acquisition outperforms uniform sampling in reducing IMSE and improving predictive accuracy.

Sequential Sparse Gaussian Process Quantile Regression

Introduction and Motivation

Quantile regression provides direct estimation of conditional quantiles, offering richer characterization of predictive distributions than mean regression, particularly relevant for risk assessment and uncertainty quantification. In Bayesian frameworks, Gaussian Process (GP) quantile regression enables principled uncertainty quantification of conditional quantiles by marrying a GP prior on the quantile function with an asymmetric Laplace likelihood. While attractive, this Bayesian formulation presents two challenges: (1) nonconjugacy of the asymmetric Laplace likelihood renders exact inference analytically intractable, and (2) scaling to large datasets and complex function spaces in GPs mandates sparse approaches, where selecting and adapting the set and placement of inducing points is critical.

The paper proposes a methodology addressing both challenges with three primary contributions:

1. Laplace Approximation in Sparse GP Quantile Regression: Posterior inference over inducing variables is achieved with a Laplace approximation, reducing the optimization space from $\mathcal{O}(M^2)$ to $\mathcal{O}(M)$ by directly solving for the MAP and the Hessian, bypassing full variational parameterization.
2. Variance-based Adaptive Inducing-Point Selection: Inducing points are adaptively placed to maximize the integrated reduction in conditional-prior variance, systematically controlling model complexity based on the global reduction of predictive uncertainty.
3. Sequential Uncertainty Decomposition for Enrichment: The total predictive variance decomposes into conditional-prior and posterior-induced components. Each is targeted by either inducing-point placement or training data acquisition, interconnected by a sequential strategy driven by a variance-based switching rule.

Bayesian Sparse GP Quantile Regression Framework

Sparse GP Formulation

The core setting models the unknown conditional quantile $q_\tau(\mathbf{x})$ as a GP with $M \ll N$ inducing points $\mathbf{Z}$, yielding inducing variables $\mathbf{u} = q_\tau(\mathbf{Z})$. The likelihood for each observation is asymmetric Laplace; given non-conjugacy, the marginal likelihood is intractable. The core inference proceeds by:

• Approximating the posterior $$\pi(\mathbf{u} \mid \mathbf{y}, \mathbf{X}, \mathbf{Z})$$ with a Gaussian via Laplace's method centered at the MAP $\hat{\mathbf{u}}$, with covariance from the negative Hessian.

This contrasts variational approaches where mean and covariance of $\mathbf{u}$ are both optimized; the present approach dramatically lowers computational costs and optimization complexity, as supported in experimental sections.

Predictive Uncertainty Decomposition

The predictive posterior GP for $q_\tau(\mathbf{x})$ conditionally decomposes its variance into two terms:

• Conditional-prior variance: Reflects uncertainty due to the sparse representation, depending on width and placement of inducing points.
• Posterior-induced variance: Propagates uncertainty in $\mathcal{O}(M)$0 to $\mathcal{O}(M)$1, dependent on uncertainty in $\mathcal{O}(M)$2's posterior, i.e., reflecting limited data.

This motivates targeting each component in a dedicated manner.

Adaptive Sequential Enrichment Strategies

Inducing-Point Infilling

New inducing points are placed by maximizing the integrated reduction in conditional-prior variance across $\mathcal{O}(M)$3. This criterion is computationally tractable due to closed-form expressions and admits gradient-based optimization owing to the smoothness conferred by the Laplace approximation.

Figure 1: Integrated conditional-prior variance.

Data Acquisition

Training data acquisition is carried out via rejection sampling, with locations selected according to high posterior-induced variance, thus focusing sampling where uncertainty due to limited data dominates.

Variance-based Switching Criterion

A key heuristic alternates between the two enrichment types: if the integrated conditional-prior variance exceeds (scaled by a factor) the posterior-induced variance, another inducing point is added; otherwise, additional data is acquired. This adaptively allocates computational budget to address the bottleneck to predictive accuracy.

Figure 2: Evolution of normalized integrated conditional-prior and posterior-induced variance as the number of inducing points increases, showing the inflection where enrichment strategy should switch.

Empirical Analysis

Validity and Accuracy of Laplace Approximation

The Laplace posterior approximation is compared to MCMC sampling. Q--Q plots between MCMC and Laplace-approximated posteriors for $\mathcal{O}(M)$4 show convergence as dataset size grows, supporting the large-sample validity of the Laplace approximation.

Figure 3: Q--Q plots comparing ordered statistics of samples from the true posterior versus the Laplace approximation for several dataset sizes.

Similarly, predictive posterior-induced variance from Laplace and MCMC match closely, particularly beyond moderate sample sizes.

Figure 4: Comparison of posterior-induced variance from MCMC and Laplace approximations for different dataset sizes.

Effectiveness of Adaptive Inducing-Point Allocation

Adaptive placement (integrated variance reduction, IVR) is compared to non-adaptive (QMC, e.g., Halton sequence) placement. IVR yields exponential decay in conditional-prior variance and achieves lower IMSE over a range of $\mathcal{O}(M)$5, outperforming non-adaptive placement until overfitting occurs due to insufficient data.

Figure 5: Integrated conditional-prior variance reduction for adaptive (IVR) and QMC allocation strategies.

Additionally, tracking both variance terms shows that once conditional-prior variance falls below posterior-induced variance, further addition of inducing points offers diminishing returns or even worsens IMSE, highlighting the efficacy of the switching rule.

Global Model Accuracy and Asymptotic Behavior

Experiments varying both $\mathcal{O}(M)$6 and $\mathcal{O}(M)$7 indicate that optimal $\mathcal{O}(M)$8 increases with $\mathcal{O}(M)$9, and that global IMSE decreases asymptotically as both grow—consistent with theoretical expectations for sparse GP models.

Figure 6: IMSE versus number of inducing points for varying training dataset sizes, showing optimal complexity increases with available data.

Sequential Enrichment and Data Acquisition

Applying the full sequential enrichment—interleaving adaptively-allocated inducing points and data—demonstrates that adaptive data acquisition offers the most impact for highly nonstationary/heteroskedastic settings. On functions with localized, high-variance regions (e.g., Michalewicz), adaptive data acquisition via the posterior-induced variance yields faster and more robust reduction in IMSE than uniform sampling.

Figure 7: Sequential enrichment on the Sabater 2D function, tracking IMSE under different data acquisition strategies.

Practical and Theoretical Implications

Practically, the methodology supports fully adaptive, end-to-end GP quantile regression where model complexity is not fixed a priori but grows only as warranted by the information content of the data. The computational cost is dominated by $q_\tau(\mathbf{x})$0 per iteration, but is mitigated via the Laplace approach and closed-form derivatives, yielding practical gains over variational methods. Theoretically, the decomposition of uncertainty and adaptive allocation of model capacity and data addresses overfitting and inefficiencies endemic to both under- and over-specified sparse models. The architecture is directly extensible to risk-averse Bayesian optimization and can be adapted for heteroskedastic or multi-quantile settings.

Conclusion

This work provides an efficient, theoretically principled framework for Bayesian quantile regression that integrates sparse GP modeling, Laplace-approximated Bayesian inference, and sequential adaptive enrichment driven by variance decomposition. The result is a scalable and robust method for high-fidelity, uncertainty-aware quantile estimation that adapts to the inferential bottleneck at each stage—whether data or model capacity. The proposed approach is supported by both strong empirical evidence and sound theoretical arguments, with immediate applications in risk-sensitive Bayesian optimization and uncertainty quantification. Future work should pursue deeper theoretical guarantees for the surrogate objective and explore extensions to fully nonparametric and hierarchical forms.