Scalable and adaptive prediction bands with kernel sum-of-squares
Authors: Louis Allain
Sébastien Da Veiga
Brian Staber
Venue: NeurIPS 2025 poster ()
Keywords: kernel sum-of-squares
uncertainty quantification
conformal prediction
kernel methods
Primary area: probabilistic_methods

Abstract: Conformal Prediction (CP) is a popular framework for constructing prediction bands with valid coverage in finite samples, while being free of any distributional assumption. A well-known limitation of conformal prediction is the lack of adaptivity, although several works introduced practically efficient alternate procedures. In this work, we build upon recent ideas that rely on recasting the CP problem as a statistical learning problem, directly targeting coverage and adaptivity. This statistical learning problem is based on reproducible kernel Hilbert spaces (RKHS) and kernel sum-of-squares (SoS) methods. First, we extend previous results with a general representer theorem and exhibit the dual formulation of the learning problem. Crucially, such dual formulation can be solved efficiently by accelerated gradient methods with several hundreds or thousands of samples, unlike previous strategies based on off-the-shelf semidefinite programming algorithms. Second, we introduce a new hyperparameter tuning strategy tailored specifically to target adaptivity through bounds on test-conditional coverage. This strategy, based on the Hilbert-Schmidt Independence Criterion (HSIC), is introduced here to tune kernel lengthscales in our framework, but has broader applicability since it could be used in any CP algorithm where the score function is learned. Finally, extensive experiments are conducted to show how our method compares to related work. All figures can be reproduced with the accompanying code at https://www.gitlab.com/drti/ksos-bands.

Scalable and adaptive prediction bands with kernel sum-of-squares
Louis Allain 1.2
Sébastien Da Veiga²
Brian Staber
Safran Tech, Digital Sciences & Technologies,78 14 Vagny-Les-Hameaux. France Univ Rennes, Ensai, CNRS, CREST UMR 9194, F-35000 Rennes, France (louis. briinniiinnn.. @ a sebastien da ta-veiga@ensai...
Abstract
Conformal Prediction (CP) is a popular framework for constructing prediction bands with valid coverage in finite samples, while being free of any distributional assumption. A well-known limitation of conformal prediction is the lack of adap- tivity, although several works introduced practically efficient alternate procedures. this work, we build upon recent ideas that rely on recasting the CP problem statistical learning problem, directly targeting coverage and adaptivity. This statistical learning problem is based on reproducible kernel Hilbert spaces (RKHS) and kernel sum-of-squares (SoS) methods. First, extend previous results with a general representer theorem and exhibit the dual formulation of the learning problem. Crucially, such dual formulation can be solved efficiently by accelerated gradient methods with several hundreds or thousands of samples, unlike previous strategies based on off-the-shelf semidefinite programming algorithms. Second, we introduce a new hyperparameter tuning strategy tailored specifically to target adaptivity through bounds on test-conditional coverage. This strategy, based on the Hilbert-Schmme Independence Criterion (HSIC), is introduced here to tune kernel lengthscales our framework, but has broader applicability since could be used in any CP algorithm where the score function is learned. Finally, extensive experi- ments are conducted show how our method compares to related work. All figures can be reproduced with the accompanying code at colllabcom/drt/kkoossssss
Introduction
In many applications, machine learning regression models require trustworthy uncertainty quan- tification their predictions. This is especially true for high-stakes applications such design optimization, non-destructive testing, medical diagnostics, autonomous vehicles, or financial fore- casting, where decisions based on model predictions can have significant impacts. Having reliable uncertainty quantification is thus fundamental. Several machine learning models come with un- certainty quantification their predictions, such as Gaussian Processes [Rasmussen and Williams/ 2005), Random Forests [Breiman] 2001] Bayesian Neural Networks [Wang and Yeung. [2020], among many others, but these models generally provide inaccurate prediction bands: coverage guar- antees typically hold asymptotically or with strong distributional assumptions. In practice, however, especially for high-stakes decisions, we should at least provide marginal coverage guarantees that hold in finite sample and without making any distributional assumptions desirable feature is adaptivity: we would like prediction bands to be wide when either the model lacks confidence if the variability in the data is high, and narrow when both the model is confident and the variability is low. Having adaptive prediction bands means having an uncertainty quantification that is informative either the performance of the model or the variability of the data, which is key in crucial applications.
the data. In addition,
39th Conference on Neural Information Processing Systems (NcurIPS 2025).

---

Conformal Prediction (CP) (see, e.g., Gammerman et al. 1998, Papadopoulos et al., 2002 Shafer and Vovk 2008 or [Angelopoulos and Bates 2023 for modern introduction) has been designed from the ground up to be an uncertainty quantification statistical framework that provides marginal coverage guarantees in finite sample while being distribution-free In particular, the split conformal procedure proposed by [Papadopoulos et al. [2002] is especially easy to implement. CP is becoming widely used in many different applications (see [alasubraaaniann et al. [014 Vazquez and Facelli, 2022) and references therein), but unfortunately, by construction, standard CP does not provide adaptive prediction bands. A lot of research has been done in this direction, which we will review later.
In parallel, for specific statistical learning problems, Marteau-Ferey et al. [020] introduced a new kernel framework known as kernel sum-of-square (SoS), tailored specifically to estimate non-negative functions. Their key idea to characterize such functions of interest by a linear positive semidefinite Hermitian operator, which admits a fiite-dimensional representation through a representer theorem. Since then, it has been leveraged for non-convex optimization (Rudi et al. 2025]. estimation of optimal transport distances [Vacher et al., 2021]. modeling of probability densities Rudi and Ciliberto 2021] and PSD-constrained functions [Muzellec et al. 2022]. Very recently, kernel SoS was also identified as a powerful framework for constructing more adaptive prediction bands by Liang 2022 and Fan et al. [2024). They were the first to propose to build prediction bands solutions of such a learning problem, where adaptive coverage is targeted with additional constraints. This point view is the one we adopt and generalize in this paper.
Outline and contributions. We start by presenting CP in Section as well recent variants developed to improve adaptivity. We also introduce the kernel SoS framework, since our proposed method extensively relies on it. In Section 3 we introduce our approach that learns CP score function by solving a statistical optimization problem with several ingredients: an objective function controlling both the width and the regularity of the prediction bands, and constraints for coverage. We also discuss new criterion dedicated to the tuning of kernel lengthscales, with local coverage as an objective. In Section 4 we finally conduct extensive experiments to compare our method to other conformal prediction methods that provide adaptive bands.
Our contributions are as follows:
We generalize the previously introduced kernel SoS point of view for prediction bands and precisely analyze the contribution and practical effect of each term in the objective function,
We provide a representer theorem that makes the problem numerically tractable,
We derive a dual formulation of this problem and propose an accelerated gradient algorithm to enable faster computation on large datasets, unlike previous work that was limited to small datasets,
We introduce new criterion to tune kernel hyperparameters based on the Hilbert-Schmin Independence Criterion (HSIC), which is also applicable to any other CP method. In particular, we provide both theoretical and empirical evidence of the effectiveness of this metric to achieve better adaptivity.
2 Conformal prediction and kernel sum-of-squares
2.1 Conformal prediction
Split conformal prediction. CP was introduced by Gammerman et al. [1998]. with the so-called full-variant, but focus here the split variant, introduced by Papadopoulos et al. 2002 Suppose we have a training dataset D = {(,,, ) from a pair (X,Y) ~ Pxy where X E C R and Y y C R.. This dataset is split in two parts: pre-training dataset D,, {(x, , Y:) j and calibration one D,,. = ((X,, Y.). with N = n + m.
The pre- training dataset D, is used to first fit predictive model (.).). which can be any machine learning algorithm. The second step consists in computing performance scores associated to â on the hold-out calibration dataset Dnn- The usual score in the literature is defined as the absolute emoos S S(X,, Yi) == S, = Y m. Xx) for i E Dm, which are used to compute the quantile To of the set (S.he... with an adjusted level [[1 (m + 1)1/m, where is the desired error rate. Finally,

---

for a new observation +1+ the split CP prediction bands are Čx = [mn + Qul. which satisfy the marginal coverage
P YN++ E (x(Xx+1)) 2 21-a
(I)
for any N as long as (X,, Y1) (XN, Yn), (XN+1, YN++1, are exchangeable. Unfortunately, as underlined by Romano et al. (2019), these prediction bands cannot be adaptive since they have constant width 24... The research direction in recent years has been to adjust CP procedures to target more adaptive bands.
The quest for adaptivity. Historically the first idea was change the score function by rescaling the scores with an estimate of the variability (.) 2 O. The new scores are thus defined as S; = |Y;
[m., (XN+1) = T. (X++1)) Lei and Wasserman (014)) first proposed to use an estimate of the conditional mean absolute de-
(X.)) T (X+). i E Dm with prediction bands CN (x)+1)
viation for (-). Another sensible choice is to scale the scores by an estimate of the standard deviation, as was suggested for several machine learning models like Gaussian Processes, Ran- dom Forests and Bayesian Neural Networks [Johansson et al.) 2014. Papadopoulos] [0224] Jaber ct al. 2025]. However, such scaling functions are rarely estimated in a goal corennede way, without quantitative and explicit consideration for adaptive coverage. This often leads in practice to poorly adaptive prediction bands. In a parallel line of work, the popular Conformalized Quantile Regres- sion (CQR) [Romano et al. 2019] proposes to change the score function by leveraging quantile regression. Instead of using an interval built around an estimate m.)..) of the regression function, they rely estimates Tntower (.) and In appper (.) of the conditional quantiles, and build an interval
(x(XN+1) = [fftoowrr (XN+1) 11) - Co: In Gupper (XN+1) + [o] where now Z. is the adjusted quantile of the set [max (ga (X,) - Y,,Y Y; In
Qupper (x.)) i E D. In other words, the score function is chosen as S(X,Y) = max (g (X) Y.Y qm gounver (X)). By design, the CQR model is adaptive and generally provides sensible prediction bands. However, suffers from two practical limitations: (a) decision maaonnn people usually prefer a point estimate with an interval around this estimate and (b) quantile regression in small data regime can be quite challenging. Also note that regularized version of CQR tailored to target test-conditional coverage was recently proposed by Feldman et al. 2021].
Another line of work consists in modifying the calibration step. For example, Guan 2023] and Hore and Barber (2024) propose to weight the scores when computing the adjusted level quantile, where the weights depend on the test point X N1 this directly implies that the quantile changes with XN+1 and as a result the prediction bands are adaptive. More precisely, given a kernel H...)) .) that defines a density H (z. :) for all x E X, sample X N441 from (XNN+1. The quantile qe. (x1,, XN+1) 11) is computed the empirical distribution > Dx=1 w;ds, + WNN+10+ac where the weights are computed as W; = (X,, X / === HX,,, XN+1)= + ((N++1, =((X)=1.N++1). Although appealing, such modifications to the score function have some practical shortcomings. First, computing different weights for all test points can be computationally demanding during inference. Second, the method suffers in practice from the randomization induced by the sampling of XN+1. To overcome this issue Hore and Barber (0024) propose the m-RLCP methods that averages the predictions bands over m sampling of X +1 +1 However this leads to a marginal coverage equal to 1 - 2a and the computational cost increases significantly. Finally, and perhaps more importantly, the kernel (:.-)) involved in the definition of the weights depends on bandwidth hyperparameter that must be tuned. This choice has a strong impact on the shape of the prediction bands as shown in Hore and Barber [2024). We will come back to this point later since our framework shares the same characteristic. Finally, another reweighting method based instead on Jackknife+ was proposed by |eutschmann et al. 2024).
2.2 Kernel sum-of-squares
Since our proposed method is based on kernel SoS for positive functions, give a brief overview following Marteau- Ferey et al. [2020).
Let H be a RKHS with associated kernel h and Φ: X - H one of its feature map such that k(x,,')) = o(X) (xX). Let S(H) be the set of bounded Hermitian linear operator from H to H. For A E S(H), we write A 2 0 when A positive semi-definite (PSD) operator, and S (HH the
3

---

set of such PSD operators, For all X. we define f((X) = 6 (6(X), A (+H which is non-negative by construction. Kernel SoS refers statistical learning problem where the unknown nonparametric function is constrained to be non-negative and obtained as the solution of
inf f((X,))))).. + A|I+ + + A2 A=0
where || and || F denote the nuclear norm and the Frobenius norm operator A. respectively, Interestingly, Martcau-Ferey et al. (2020) show a representer theorem for Equation (2) which makes
kernel SoS computationally tractable finite-dimension and recalled below.
Theorem (Marteau-Ferey et al. [2020] Assume L: IR + R U {+00} to be lower semi- continuous and bounded below loss function. Equation 2 admits a solution A which can be written A = Ei,==1 (x)(()
for some matrix B E R B > 0. Furthermore
A is unique if L is convex and X2 0. The corresponding non-negative function is given by
(X) = Ei.=1 Bi,k(X, x),,
Theorem provides dinee--immnsionall equivalent problem which involves an unknown PSD matrix B. But Marteau-Ferge et al. [2020] also propose an equivalent formulation: considering K the kernel matrix with elements Ki (<X,, X,) and V its upper Cholesky decomposition, we can define ((X) = V V- _x with kx = X)))= and (X) ( (X). With these notations, the following proposition shows that obtain the same solution if we optimize the PSD) matrix A instead of B.
Proposition 1 (Marteau-Ferey al. 2020) Under the assumptions of Theorem the following problem has at least one solution, which is unique if X2 >> O and I. is convex:
inf ( AllAAI- -
+
+ A2
(3)
A
For any given solution A R of Equation 3 the function fA. is also solution of Equation 2
Although such a result may appear of minor impact, the A formulation actually yields significant computational savings in practice, as we illustrate in our numerical experiments (see Appendix B.5).
Remark 1 Operator A is PSD, hence it admits an eigendecompositiion A = Σizo dna & with A 2 and u, E H. By the reproducing property, we have fA(X) = -E10 ( (X) Hence, fa(X) is an infinite sum of squared functions in H. Equivalently in finite dimension, for a PSD matrix B UDUT we have [B( = k UDU kx = -i= > uak X, x.)) Thus JB function defined as linear combination of squared functions from H.
3 Regularized kernel SoS for adaptive prediction bands
We now come back to our initial problem of building adaptive prediction bands. To do so, we focus on the split CP setting with two i.i.d. datasets D, {X,,, k2) and D,, {(x,, )m 1- and consider estimating a score function in specific supervised learning problem to achieve better adaptivity.
3.1 Learning the scores through an optimization problem
A general framework for learning score functions was recently introduced by Xie et al. 2024a given task-specific loss (c.g. conditional coverage or minimum interval width), optimized score function is obtained via boosting algorithm. In this work, the score function S (xXY) Y) = max (X) Y, Y M2 (X)) /x(X) is inspired by CQR and is parameterized by three unknown functions (14,, P2, or) such that (+(+) < p2(-) and o(.) 2 O, which are iteratively optimized during boosting rounds. From practical viewpoint however, the constraints these functions not inherently accounted for in the boosting algorithm.
We advocated before the use of prediction intervals built around point estimate, which corresponds to the particular case m = f11 = 1/2 and yields a score function S(X,Y) Y m(()l//((X) or equivalently S(X, Y) (Y -m(X)) (XX with f(-) a-)) since only the score quantiles

---

involved in the final prediction interval: we thus recover the rescaled conformal score setting where we learn the rescaling function f(.), with an additional non-negativity constraint. The kernel SoS framework is thus a natural candidate for this learning task.
Kernel SoS formulation for prediction intervals. We introduce two RKHSs 7 and H associ- ated to kernels km with lengthscales 8t and & with lengthscales er, respectively. The regression function m will estimated in the RKHS 7thh while the non-negative scaling function will be estimated using the kernel SoS framework in the RKHS H. The first step is to derive our proposed infiniteddimensionl learning problem from the properties that we impose on prediction bands:
inf mEH", (HL)
Σ (Y m(X,))
Σ fA(X)) + X1lLAII, + 22|
(4)
s.t. fa(x,) 2 (Y, mxx)) ie [n],
(5)
|| 2 8.
(6)
In this problem, we propose include several key components:
1. Accurate mean estimation (first term in (4) with regularity penalty (6) as in standard kernel ridge regression [Schölkopf and Smola 2002].
2. Prediction intervals with minimum mean width (second term in and an additional regularity penalty (third and fourth terms in 4.
3. 100% coverage on pre-training data (5) later calibrated on the calibration dataset.
Minimizing the nuclear norm with coverage constraints was originally proposed by Liang 2022 and minimizing the mean width was later proposed by Fan et al. (2024). Our proposition essentially differs from their work for broader and more efficient practical applicability: (a) we place ourselves in the split CP setting, whereas Liang [0022] and Fan et al. (2024) propose different calibration procedures that are harder implement in practice, with theoretical coverage guarantees that depend on hyperparameters (see Appendix A.3), (b) we rely on a dual formulation which can handle several hundreds of samples, (c) we propose goal-oriented tuning strategy for 0f and (d) our proposal is more general and we give better understanding of why this targets adaptivity"] In particular, adaptivity can actually be controlled through the complexity of the scaling function f. three different ways:
1. The sparsity in the linear combination. This can be controlled by the nuclear norm |||||| which acts as a lasso-type penalty (see Recht et al. (2010)).
2. The & norm of the coefficients in the linear combination, which is equal the Frobenius norm || and is similar to ridge penalty,
3. The RKHS H which impacts the functions uz(-), and notably the lengthscales ff.
To further illustrate point 3., if we take kf (2) (z, x²) (RHKS of linear functions), the rescaling function will be second-order polynomial (thus not complex), while with Gaussian kernel gf + +oo the prediction intervals will be constant (not adaptive enough) and if 0 is small they will be very wiggly (too adaptive). We thus see that in between, we can target better adaptivity: we propose in Section 3.2 new goal-oriented criterion related to local coverage to tune the lengthscales 01. Note also that a instead of 1 100% coverage can be considered the constraints, but this unfortunately leads to a non-convex optimization problem [Braun et al. 22225.
Before diseussing detail the choice of our problem hyperparameters, we first derive representer theorem which makes the optimization numerically tractable.
Theorem (Representer theorem) Let (a, b, s, A1) T and X2 > 0. Then Equation admits
a unique solution (m* ffA)
) of the form m² (X) = Ei= Yi kmh (X, X) 7 k2 and fa (X)
4 (X)T A (X) for some vector x* E R and matrix A E Ru
O.
The detailed proof, based on Marteau- et al. Appendix A. This representer theorem leads to the following tractable semi definite programming
2020 and Muzellec et al.
2022 can be found in
'Tis lack of understanding is for example pointed out by Liang al. 2024
5

---

(SDP) problem:
b
inf
(Yi 7 kx,
n Σ ( A2/A
YER", AEST n i=1
(X, 2 (Y, 7 k?, iE
Km S
In practice, such SDP problem can be solved using off- off-the-shelf solvers like SCS [0 D'Donogh al., 2016], as was advocated by Liang [2022] and Fan al. [024). However, our numerical experiments show that this strategy does not scale past few hundreds pre- pre-training samples: this is
thus severe practical limitation which, in our opinion, heavily weakens the kernel SoS point of view.
To circumvent this major issue, we rely instead on dual formulation of Equation (7)
Dual formulation.
First note that this formulation is possible only when X2 (thus excluding
Liang [022] framework) and that it consists of an optimization problem over (n + 1) variables rather
than (1 + n) variables.
Proposition 2 (Dual formulation) Let (a, b, 8, A1) E TR4, A2 > and A = RR+11 Equation admits dual formulation of the form
sup r(I,0) Diag(Ta)r (TT, 0) ((T(T, 0) K (T, 5) 2- (VDiaag(T-bVV) (8) (T))
where r(T, Y Y--K-" (I,
(T, 0) TI 0) C (T, 0) Diag T. Y. C [T 0) Diag (Ta) K 1 B - I and Vx E R, Diag (T,) == Diag (T) + Ln. Moreover, if
AL, Ω2" (B)
= 4r2
(f, 0) is solution of Equation p. a solution of Equation (4) can be retrieved as
5 (Diag(f.)) + OI,, Diag(f.)) and À 2X2
VDiag P AIIn
where [A]+ denotes the positive part of A
detailed proof is given in Appendix A.2 as well as the dual analytical gradient. Interestingly, this dual formulation be efficiently optimized using accelerated gradient-ascent algorithms [Ruder 2017 Xie et al. 2024b), see Figure 3 for an illustration of the numerical speed-up.
Final prediction bands. Solving the dual formulation yields estimates f(() = 7 kT and (x) = (X) T (xX), from which we derive the estimated score function S(X, Y) IV m(XI/V (X). Following standard split CP, our prediction interval is given by
(x(4))1) = (((N+1)) aa fA(NN+1): (mm.:m11.mm.mmmmmmmmmm ((Xx+1) + Za V (A(*4)))
(9)
where T. is the adjusted quantile of the estimated score function on the calibration set. By design (see Appendix A.3), they satisfy the same marginal coverage guarantee as split CP, as formalized below.
Proposition 3 If the samples (X,, Yi). (XN, YN) (XN+1) YN41) are exchangeable, then the prediction interval (9) satisfies the marginal coverage 7
Extensions. Similarly to Liang 2022], if a point estimate (-)) is available beforehand from any machine learning model trained on separate dataset, our dual formulation can be easily modified. More importantly, it is also straightforward to consider non non-syymetric intervals by estimating two functions fAw (-) and fA (.), at the cost of increasing the dual problem dimension to (2n + 1). see Appendix B.5. We also strongly believe that incorporating constraints of the form μ1 (-) < (22())
proposed by Xie al. [2024a] is possible with kernel SoS, by taking inspiration from Aubin Frankowski and Rudi 20244. Kernel SoS is also complementary to Hore and Barber 20244) as our learned score can be post-processed with their approach. Finally, since we aim for 100% coverage, our framework may lack robustness with respect to outliers: however, we can take advantage of our preliminary estimate mGP(.) (see next section) filter out samples with large residuals.
For PSD matrix A with cigendecomposition A UDUT, its positive part is defined as [A] U max(0, D)UT..

---

3.2 Hyperparameter tuning
First, the hyperparameters 9 and related to the regression function m(-) are fixed with a preliminary estimate mar(? (-) obtained from Gaussian Process (GP) regression with kernel km 4 is optimized by maximum likelihood and s = || lm Second, our numerical experiments in Appendix B.2 show that X2 has very small impact on the shape of the intervals, as long as X2 > 0. Not only does it make the initial problem strongly convex and makes the dual formulation possible, but also facilitates numerical optimization. In the following, we thus fix it at X2 1.
Hyperparameters a, b, A and Of necessitate more attention. We observe numerically that a has very little influence on mean squared-errors when the noise is symmetric (see Appendix B.2). This phenomenon was commented in |an et al. [2024], where they argue that the coverage constraints have the additional effect to reduce mean squared-errow
but here we insist that this is true only for symmetric noise. Our experiments confirm this intuition, we simply set a 0 for the symmetric case. Furthermore, they also show that b and Å1 have highly different influence. This may seem surprising. given how they interact the dual formulation through the last term in (8) potential explanation is that b controls data-driven term, unlike Åι and the latter does not have a strong impact on the optimal solution provided is positive. As a result, we choose to set À 1. On the contrary, has the expected behavior, in the sense that higher values yield narrower intervals, at the cost of increasing the nuclear norm, that is the function complexity. But the last hyperparameter of also controls such complexity: we thus expect compensation between these two hyperparameters through an interaction. Numerical experiments in Section f actually confirm this phenomenon and show that as long as b is sufficiently large (e.g. 10 or 100, recommend to test these two values), it is possible to tune only o to reach the same level of adaptivity. of thus requires specific consideration, since has the largest impact when b is fixed. This is the one we focus on below.
Goal-oriented criterion. Adaptivity is tightly related to local coverage P(YN+1 E Č(X ((N++ N+1
X x) 2 1 Or which is impossible to achieve exactly in a distribution- setting [Vovk] 2012, Barber et al. (02T]. We focus ourselves on weaker version, where we condition on X being in a small neighborhood wi E Fx from the event space Fx such that for all x E X, P(x E wx)= 2
P(YN+1 C(XN+i) 1 XN+11 wx) 2 1 Recently, Deutschmann et al. [2244] showed that
such coverage with split CP can be controlled with the mutual information (Ml) between the inputs and the score function, namely
P(Y E (X)INN, X x1 1 cxp- (Y,Y)))
(10)
Note that contrary to Deutschmann al. 20024] explicitly write the coverage conditionally on the training set D Inherently this bound is uninformative for low- low-probability sets, but interestingly can counterbalance this effect by choosing a score function which is as independent possible of the inputs in order to get MI(X, S(X, Y)) close to O, see Deutschmann et al. [2024). However, this implies computing MI for random vector in dimension d: from practical perspective, MI suffers from the curse of dimensionality and rapidly becomes numerically unstable. Instead, for a re-scaled score function, we show that similar bound holds with MI between one-dimensional variables, which is more robust. Furthermore, using recent inequalities between the total variation distance and the maximum mean discrepancy [Wang and Tay 2023), we also extend our result to replace MI with HSIC (Gretton et al. 2005), for which we observe in practice much improved numerical stability.
Proposition 4 Let ÖDy be the prediction intervals built from score function S(X,Y) -Y-- m(X)// ((XX)///(X) ((() through split CP with D = Dn U D'- Then for any wx Fx such that P(X wx) 2 8, denoting PDx = P(YN-1 E ČDY (XX+)))DD XN+1 E wx) we have:
PDx 21-a
/ eeteeeepppMMMp... (XN (XN+1, YN+1), Fo.n),FP,.CC Fo. (XN+1)))
(11)
where On does not depend on f(.) and rD, (X,Y) =0000 - mip (X), J. In addition, we have
PD* 2 1 - 5 1 azHSIC(PD. (XNSII YN+1), Fo.. (XN+1))
(12)
where 02 only depends the kernel used for HSIC, which must be characteristic.
7

---

To target local coverage, we can then maximize ISI in our kernel SoS procedure, this means that of can be tuned efficiently according to this criterion. For HSIC estimation, we need samples N independent of D,, and thus rely on cross-validation procedure. Finally, handle cases where the noise is homoscedastic, we implement as a last step HSIC test of independence: the p-value is large, we 0 to an arbitrary large value. See Appendix A.3 for the proof and B.3 for implementation details.
Remark 2 We believe that Proposition is also of interest beyond kernel SoS, in the sense that provides principled way to tune hyperparameters in all score functions, such in Hore and Barber (2024] or Braun et al. [2025] for example. However, its true potential would lie in generalizing Equation I2 to other score functions, which is left as future work. It could also be directly used as loss function in Xie et al.) 2024a], instead of estimating local coverage.
4 Experiments
4.1 Small data experiments
For all experiments, we consider a Matérn 5/2 kernel for kt and kf Our proposed kernel SoS algorithm is first compared in the small data regime to standard competitors for split CP: CQR (with random forests, following experiments from Romano al. (2019, Hore and Barber [2024) and rescaled scores. For the latter, we only consider two variants of GP for fair comparison, since we also place ourselves in the RKHS setting. We focus on homoscedastic and heteroscedastic [Binois et al.) 2018) GP model, and consider here Y = m(X) + (X )E with:
m (X),, 7 10X +1
(XX v0.1 2X² [snn(X]
X
E N 1) N( 1)
Case 1 Case 2
sin( 28
0.2 1x2
19 + ( 1)1x996
u[-1,1] (O
5
X/2
We begin by illustrating the interaction between and 0f in Figure left. For all values of observe a consistent HSIC behaviour: it first increases with 0f, thus improving adaptivity, until reaches peak and then decreases. Interestingly, we also note that the higher b, the higher the optimal 01: this clearly shows that both hyperparameters have opposite effects on adaptivity. Furthermore, we see that for 2 10 we reach a plateau for the optimal HSIC. In practice is thus sufficient to only optimize o as soon as b is fixed at large enough value. Figure I also shows that small value for 0 leads to overly adaptive bands, while the HSIC-optimized of produces smooth and adaptive ones.
0.014
18RO' IE
BEC LIX ID
0.012
RSIC (10
0.0000
1000
00008
10000
Siood & 100 Moss prodiction Farak Preckition
1.74 and Mose prodi
000044
Fardb Preekttoo
0.0002
Callentiee Duta Cou Data
Cullearae
Tes Data
0.0000
3g33
6'
Figure 1: Test case 2 with 100. Left: HSIC criterion between r(X, Y) and f(X) as a function of b and of (confidence intervals obtained by bootstrap and optimal values of of in dashed lines). Middle Right: optimal prediction bands with too small and optimized lengthscale, respectively.
In the following we now fix b = 10 and compare kernel SoS with HSIC-optimized o to CQR and rescaled GPs on test case 1 with replications. Figure we investigate several metrics related to adaptivity: mean width, MI and local coverage, see also Appendix B.5 for complementary discussion. We note that CQR and standard GP yield larger intervals and local coverage with great variability around the target. Heteroscedastic GP produces intervals similar to ours, but with higher mean width. In contrast, kernel SoS gives prediction intervals with both small width and satisfying local coverage properties. Additional experiments confirming this behavior to be found in Appendix B.5.

---

MIIXX
KSoS
CQR with RF
Homoscedastic GP
Heteroscedastic GP
Demstry
= 0.9
5
Wadth
8.4
0.5
0.6
0.7
0..5 1.0
Looal Coverage
Figure 2: Test case with n = 100. Adaptivity metrics and density of local coverage.
4.2 Real-world and large seale experiments
Finally, we demonstrate that our dual formulation for solving the kernel SoS problem can scale to several hundreds or thousands training points, unlike previous work that relies on SDP solvers. Figure 33 left, shows that a SDP solver can only handle up 200 training samples, while our dual solver scales casily to = 1000 (see Appendix B.4). When we optimize of with HSIC, we retrieve the optimal solution in Figure right.
Draladlo
Bunds Prediction
40
Prinsoaaber (SCS)
prediction
Oroele
Test Data
Titee (.)
20
10
-2
11
2001
01X
600
MCD
20101
1.0
Sumplesize size (m)
X
Figure 3: Test case 1. Left: time for SDP and dual formulation as a function n (a h 0, 05 0.3, max iter = 104) Right: optimal solution of dual formulation for n 2000.
At last, we consider six real-world datasets commonly used for regression. which are detailed in Appendix B.5. Table the mean width of prediction intervals on test set for all methods. We consider two variants of kSoS: one where of is chosen to minimize the mean width, and the other where we optimize with HSIC. Note that for Bio, which exhibits strong asymmetric noise, use our asymmetric kSoS variant and do not use the HSIC criterion to tune ff, since at the moment it only applies to symmetric intervals. Such asymmetric noise heavily favors CQR, which outperforms all methods on this example. For all other cases except for Concrete which we comment below, kernel SoS targeting mean width achieves better adaptivity when measured by mean width. other competitor, although they generally have strong performance, achieves this robustness. For completeness, we also provide marginal coverage for all datasets and methods in Appendix B.5.
Dataset
CQR
Het GP
Hom GP
kSoS
kSoS Opt. HSIC
Best mean width
Concrete
0.586 i 0.032
0.508 1 0.052
0.543 1 0.044
0.556 1: 0.044
0.568 : 0.06 0.803 1 0032
Bike
1.114 1 0.062
1.000 + 0.079
0.809 0.024
0.803 1 0.031
Bio
1.879 0.046
2.194 + 0.119
2.03 0.07
Diabctes
188.62 9.33
191.24 11.95
190.58 + 11.19
185.83 + 14.47
187.6 + 16.18
MPG
9.89 1 0.82
9.70 + 1.06
9.71 1 0.73
9.15 0.8
9.36 1 0.82
Housing
1.816 0.045
1.585 0.099
1.453 0.068
1.468 1 0.094
1.586 1 0.104
Table 1: Mean width of prediction intervals for six real-world datasets (seantisd on 10 repetitions).

---

At first sight, the kSoS variant based on HSIC produces larger bands and appears to be less adaptive only mean width considered. However, a more meaningful criterion evaluate adaptivity is local coverage. Since it is unfortunately out of reach for real datasets, we rely instead on the worst-set coverage min;=1,....L P(Y E Cpy (XN+1)]XN+++ E R4) where (Reht=1..... is a partition of the input space, see Thurin et al. [2025] We observe in Figure 4 that kSoS with HSIC stays closer to the target 0.9 than all other methods. For Concrete, heteroscedastic GP has the smaller mean width, but at the cost of producing intervals with lower worst-set coverage. Similarly, on Housing width exhibit the lowest mean width, but further from the target. This illustrates once again that our HSIC criterion focuses on local adaptivity.
both homoscedastic GP and kSoS focusing
KSo8 Best meon width
1808 Opt. IISIC
Heteroseddsstic GP
Homosrodastic GP
CQR with IP
0.825
0.850
0.875
0.875
0.900
Housing
Conerete
Figure 4: Housing (left, with b = 100) and Concrete (right, with b = 10) datasets. Mean and standard deviation of worst-set coverage on 10 regions with 100 test points each, 10 repetitions.
5 Limitations
Although our experiments show that kernel SoS a competitor worth trying to build adaptive bands in moderate sample size regime, we can still identify limitations that would necessitate extensions:
Kernel SoS performs well for dimensions up to 10 or 15, a limitation shared by all kernel methods without specific strategies. For high-dimensional objects with structure (time series, probability distributions, molecules, strings, graphs.....). is however straightforward to plug specifically designed kernels in our framework. On the other hand, for nigh-dimensional tabular data, a common strategy is to rely on additive kernels, which are left as future work.
Kernel SoS can scale up to 2000 samples thanks to the dual formulation, but handling much more samples than that would be highly time consuming. This is inherently a limitation that comes from kernels rather than from the dual formulation. Indeed, each objective function and gradient computation involves the usual kernel regression formulas, which scale as O(n²) ). A usual workaround is to use Nystrom-type approximations or Random Fourier Features, which could be casily integrated in problem.
Due to the 100% coverage constraint, kernel SoS is less robust to outliers than CQR. In such cases, a pre-processing step using another model (e.g. the initial GP model we use in our procedure or an heteroscedastic GP model) is recommended.
6 Conclusion
We introduce a generalized kernel sum-of-squares framework for building scalable and adaptive prediction bands for conformal prediction. Scalability is achieved through new representer theorem together with a dual formulation which can be solved efficiently with accelerated gradient algorithms. Unlike previous work, this makes the kernel SoS paradigm for CP able to scale up to thousands of training points, as we illustrate in our experiments. On the other hand, adaptivity and local coverage targeted by optimizing hyperparameters with new HSIC-based criterion, which is numerically robust and has excellent practical performance. Since such a criterion appears promising, as a perspective, plan to investigate its extension more general score functions.
10

---

References
Carl Edward Rasmussen and Christopher K. Williams. Gaussian Processes for Machine Learning.
The MIT Press, 11 2005. ISBN 9780262256834. doi: 04775//mippress332200000..00 URL
https://doi. ptttp:://do.oor//100 http///.7751mmmppreeet.... httppprwws//......eeee//// /00///20000000 1.1000000011
Leo Breiman. Random forests.
Machine Learning, 45(1):5-32, 2001.
doi: 10.1023/A:
1010933404324.
Hao Wang and Dit-Yan Yeung. A survey on bayesian deep learning. Association for Computing Machinery, 53(5), 2020. ISSN 0360-0300, doi: 10.1145/3409383, URL URLhttps://dda org/10 1145/3409383
A Gammerman, V Vovk, and V Vapnik. Learning by transduction. In Conference on Uncertainty in Artificial Intelligence, 1998.
Harris Papadopoulos, Kostas Proedrou, Vladimir Vovk, and Alexander Gammerman. Inductive confidence machines for regression. In European Conference on Machine Learning. 2002. URL https://api. semanticscholes org/ CorpusID:42084299
Glenn Shafer and Vladimir Vovk. A tutorial on conformal prediction. Journal of Machine Learn- ing Reseurch, 9:371-421, mar 2008.
shafer08a/
shafee088. pdf Submitted 8/07; Published 3/08.
Anastasios N. Angelopoulos and Stephen Bates. Conformal prediction: A gentle introduction. Found. Trends Mach. Learn., 16(4):494 March 2023. ISSN 1935-8237. doi: 10.1561/2200000101.
URL https://doi. org/10. 1561/2200000011
Vineeth Balasubramanian, Shen-Shyang Ho, and Vladimir Vovk. Conformal prediction for reliable machine learning: theory, adaptations and applications. Newnes, 2014.
Janette Vazquez and Julio C Facelli. Conformal prediction in clinical medical sciences. Journal of Healthcare Informatics Research, 6(3):241-252, 2022.
Ulysse Marteau- Ferey, Francis Bach, and Alessandro Rudi. Non-parametric models for non-negative functions. In H. Larochelle, M. Ranzato. R. Hadsell, M.F. Balcan, and H. Lin, editors, Ad- vances in Neural Information Processing Systems. volume 33, pages 12816-12826. Curran As- sociates, Inc., 2020. URL UtL./tttpss//prrccedddagg www.gis............ cc/ http:/ppaeeeeeeeeeeeeeaaaeeeeeee www.eeeeppaaaa pttpeeee//pppp/22200 (F11:668887741117777777888888882 .dff
Alessandro Rudi, Ulysse Marteau-Ferey, and Francis Bach. Finding global minima via kernel approximations. Mathematical Programming, 209(1) jan 2025. ISSN 1436-4646. doi: 10.10077510107-0044-00000
Adrien Vacher, Boris Muzellec, Alessandro Rudi, Francis Bach, and Francois-Xavin Vialard. A dimension-frree computational upper-bound for smooth optimal transport estimation. In Mikhail Belkin and Samory Kpotufe, editors, Proceedings of Thirty Fourth Conference on Learning Theory, volume 134 of Proceedings of Machine Learning Research, pages 4143-4173. PMLR, 15-19 Aug 2021. URL htRpsttss:/// htttps::/ppocceediigggg.o momdiiinnnnmmm.......mmmmm.mm..m wwwpss Titess/v134/vacher222 html
Alessandro Rudi and Carlo Ciliberto. Psd representations for effective probability models. In M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang and J. Wortman Vaughan, editors, Ad- vances in Neural Information Processing Systems. volume 34, pages 19411-19422. Curran As- sociates, Inc., 2021. URL httpstttpp///procccinng www.iiiinnnnnnniinneiee....e httpp//pcrpppppppppprroreeeeaaae http://pee..iileeppp//////00///0 filcaalccscaaa11222227aaaaaaaaas
Boris Muzellec, Francis Bach, and Alessandro Rudi. Learning pad-valued functions using kernel sums-of-squares, 2022. URL
Tengyuan Liang. Universal prediction band via semi-definite programming. Journal of the Royal Statistical Society Series Statistical Methodology, 84(44:155--1588)
August 2022. ISSN
1467-9868. doi: 10.1111/rssb. 12554 URL //tp:///xxx.a htctp://dx.xdi.ooo/
/ttpppp:/dd..ooo htt.:/ooo.org//..11///ee.... 1111/rssb.....e 12542.
11

---

Jianqing Fan, Jiawei Ge, and Debarghya Mukherjee. Utopia: Universally trainable optimal prediction intervals aggregation, 2024. BRKLittss///rxii.orrg
16549
16549
Yaniv Romano, Evan Patterson, and Emmanuel Candes. Conformalized quantile regression. In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc. Fox, and R. Garnett, ed- itors, Advances in Neural Information Processing Systems, volume 32. Curran Associates, Inc., 2019. URL htttps://proceediinng http://ooooeddnn....i........... wttp:/ss.erriiiccpppppeeeee LL0G654000cc$1l101866666666----3
files/ hatpapaefflleeppppr///////////// ptteee//pppee////ffieelel/
pdf
Jing Lei and Larry Wasserman. Distribution-free prediction bands for non- -an-ammmeiiicc regression. Journal of the Royal Statistical Society Series B: Statistical Methodology, 76(1):11-96 jan 2014. doi: 10.1 1/rssb. 12021. URL httph: org/ 1111/rssb. 12021
U. Johansson, Henrik Boström, Tuwe Löfström, and Henrik Linusson. Regression conformal prediction with random forests. Machine Learning, 97:155 176, 2014. URL https://api
semanticscholar wwwpilaaa.rrroooooeeee wtgg/CorpusID... ID:01510 14015369
Harris Papadopoulos. Guaranteed coverage prediction intervals with gaussian process regression. IEEE Transactions on Pattern Analysis and Machine Intelligence, 46(12):9072-9083 December 2024. ISSN 1939-3539, doi: 10. 10.11099ppami22044444122222 2024.3418244,
URL http: http//xxxxooooo/o//ooo/// www. org/10
1109/TPAMI 3418214
Edgar Jaber, Vincent Blot, Nicolas Brunel, Vincent Chabridon, Emmanuel Remy, Bertrand Iooss, Didier Lucor, Mathilde Mougeot, and Alessandro Leite. Conformal approach to gaussian process surrogate evaluation with coverage guarantees. Journal of Machine Learning for Modeling and Computing, 6(3), 2025.
Shai Feldman, Stephen Bates, and Yaniv Romano. Improving conditional coverage via orthogo- nal quantile regression. In M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang, and J. Wort- man Vaughan, editors, Advances in Neural Information Processing Systems, volume 34, pages 2060-2071. Curran Associates, Inc., 2021. URL URLAatps: http::///ppproooddddoooooooossss neurips..ccmmm cttppppccccppp...ppppppppp wtte://pprooooffrooooooff22222a2
pdf
Leying Guan. Localized conformal prediction: A generalized inference framework for conformal prediction. Biometrika, 100())33-50, mar 2023. doi: 10. 1...0///ioome//sscc0000... URL https
1tg//00.1093/iiomettt asac000
Rohan Hore and Rina Foygel Barber. Conformal prediction with local weights: randomization enables robust guarantees. Journal of the Royal Statistical Society Series B: Statisticul Methodology, 2024.
doi: 10. 10.10999jssssbb/kaee1000
URL h/tp://do..ooooo org/oo..r.... 1093/ h rkae13
Nicolas Deutschmano, Mattia Rigotti, and Maria Rodriguez Martinez. Adaptive conformal regression with split-jackknife+ scores. Transactions on Machine Learning Research, 2024. ISSN 2835-8856. URL
Ran Xie, Rina Foygel Barber, and Emmanuel J. Candès. Boosted conformal prediction intervals. In A. Globerson, I.. Mackey, D. Belgrave, A. Fan, U. Paquet, J. Tomezak, and C. Zhang, editors, Advances in Neural Information CrciionnPPoeeeee Processing Systems, volume 37, pages 11868-71899 Curran Associates, Inc., 2024a. URL httptttt/////eeeedingg mewpeeeddddinnnnenii.eeeeeeeeeee cc/paper.ppp....pppp.. http:/pprrrffeeeeeee fflees/pper/ 2002:/5lee2227777780000000000222
hat.ssaaaaaapppppeeooooeeeeeeeee Conference..dttee
Bernhard Schölkopf and Alexander J Smola. Learning with kernels: support vector machines, regularization, optimization, and beyond. MIT press, 2002.
Ruiting Liang, Wanrong Zhu, and Rina Foygel Barber. Conformal prediction after efficiency-oriented model selection, 2024. Utt.https://art/ org/ abss/208..0770 07066
Benjamin Recht, Maryam Fazel, and Pablo A Parrilo. Guaranteed minimum-rani solutions of linear matrix equations via nuclear norm minimization. SIAM review, 52(3):471-50) 2010.
Sacha Braun, Liviu Aolaritei, Michael Jordan, and Francis Bach. Minimum volume conformal sets for multivariate regression. arXiv preprint 210322503..906.. 2025.
12

---

Brendan O' D'Doogguee, Eric Chu, Neal Parikh, and Stephen Boyd. Conic optimization via operator splitting and homogeneous self-dual embedding. Journal of Optimization Theory and Applications, 16963):1042-1068, June 2016.
wwtp/"boyddpappers//eaaaee html
Sebastian Ruder. An overview of gradient descent optimization algorithms, 2017. URL https
//arxiv. wwwaaiiwrg/aaaa.aae..aa.aa.oooae aaaaabas11000 14747
Xingyu Xie, Pan Zhou, Huan Li, Zhouchen Lin, and Shuicheng Yan. Adan: Adaptive nesterov momentum algorithm for faster optimizing deep models. IEEE Transactions Pattern Analysis and Machine Intelligence, 46( (0(12)9958--5222 2024b. doi: 10. 10.11999TAAMM2224444223333
Pierre-Cyri Subin-Frankowskk and Alessandro Rudi. Approximation of optimization problems with constraints through kernel sum-of-squares. Optimization, pages 1-26, 2024.
Vladimir Vovk. Conditional validity of inductive conformal predictors. In Steven C. H. Hoi and Wray Buntine, editors, Proceedings of the Asian Conference on Machine Learning, volume 25 of Proceedings of Machine Learning Research, pages 475 490, Singapore Management University, Singapore, 04-06 Nov 2012. PMLR.
mlr press/v25/vovkks
html
Rina Foygel Barber, Emmanuel Candès, Aaditya Ramdas, and Ryan Tibshirani. The limits of distribution-fres conditional predictive inference. Information and Inference: Journal of the IMA, 10(4):455-45) 482, aug 2021. doi: 10. 10.10////maia//aaa011. URL http:/tps:///oi org/10
1093/imainnn iaa1017
Chong Xiao Wang and Wee Peng Tay. Semi-nonparametric estimation distribution divergence in non-euclidean spaces, 2023. URL htttt:///axxxi org/ at./a/2...20331111111
Arthur Gretton, Olivier Bousquet, Alex Smola, and Bernhard Schölkopf. Measuring statistical dependence with hilbert-schmidt norms. International conference on algorithmic learning theory, pages 63-77. Springer, 2005.
Mickael Binois, Robert B Gramacy, and Mike Ludkovski. Practical heteroscedastic gaussian process modeling for large simulation experiments. Journal of Computational and Graphical Statistics, 27 (4):808-82) 2018.
Gauthier Thurin, Kimia Nadjahi, and Claire Boyer. Optimal transporttbased conformal prediction. PMLR, 267:59509-59527, 13-19 Jul 2025. URL http:p::/prrrceeeinsslllll/ mlr. / thurin25a. htmll
13

---

NeurIPS Paper Checklist
1. Claims
Question:
Do the main claims made in the abstract and introduction accurately reflect the paper's contributions and scope?
Answer: [Yes]
Justification: We rely on numerical experiments in the main paper to support our claims, and also provide detailed proofs and additional experiments in the Appendix.
Guidelines:
The answer NA means that the abstract and introduction do not include the claims made in the paper.
The abstract and/or introduction should clearly state the claims made, including the contributions made the paper and important assumptions and limitations. A No or NA answer to this question will not be perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how much the results can be expected to generalize to other settings.
It is fine to include aspirational goals as motivation as long as is clear that these goals attained by the paper.
2. Limitations
Question:
Does the paper discuss the limitations the work performed by the authors?
Answer: [Yes]
Justification: We comment the limitations in several remarks and in the final section that guide towards future work. We also provide a discussion on computational efficiency and how it scales with sample size in the main paper and the Appendix.
Guidelines:
The answer NA means that the paper has limitation while the answer No means that the paper has limitations, but those are not discussed in the paper.
The authors are encouraged to create separate "Limitations" section in their paper.
The paper should point out any strong assumptions and how robust the results are to violations of these assumptions (e.g., independence assumptions, noiseless settings, model well-specification, asymptotic approximations only holding locally). The authors should reflect on how these assumptions might be violated in practice and what the implications would be.
The authors should reflect the scope of the claims made, e.g., if the approach was only tested on a few datasets or with few runs. In general, empirical results often depend implicit assumptions, which should be articulated.
The authors should reflect the factors that influence the performance of the approach. For example, facial recognition algorithm may perform poorly when image resolution low or images are taken in low lighting. Or speech-to-text system might not be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach to address problems of privacy and fairness.
While the authors might fear that complete honesty about limitations might used by reviewers as grounds for rejection, worse outcome might be that reviewers discover limitations that aren't 't acknowledged in the paper. The authors should use their best judgment and recognize that individual actions in favor of transparency play impor- tant role in developing norms that preserve the integrity of the community. Reviewers will be specifically instructed to not penalize honesty concerning limitations.
3. Theory assumptions and prool's
Question: For each theoretical result, does the paper provide the full of assumptions and complete (and correct) proof?
14

---

Answer: [Yes]
Justification: Complete proofs and assumptions are available in the Appendix.
Guidelines:
The answer NA means that the paper does not include theoretical results.
All the theorems, formulas, and proofs in the paper should be numbered and cross- referenced.
All assumptions should be clearly stated or referenced in the statement of any theorems.
The proof's can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged to provide a short proof sketch provide intuition.
Inversely, any informal proof provided in the core of the paper should complemented by formal proofs provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental result reproducibility
Does the paper fully disclose all the information needed to reproduce the main ex- perimental results the paper to the extent that affects the main claims and/or conclusions of the paper (regardless of whether the code and data are provided not)?
Question:
Answer: [Yes]
Justification: We always explicitly mention the hyperparameter values used in all experi- ments, and also provide reproducible code the supplemental material.
Guidelines:
The answer NA means that the paper does not include experiments.
If the paper includes experiments, a No answer to this question will not be perceived well by the reviewers: Making the paper reproducible is important, regardless of whether the code and data are provided or not.
If the contribution is dataset and/or model, the authors should describe the steps taken to make their results reproducible or verifiable.
Depending on the contribution, reproducibility can be accomplished in various ways. For example, if the contribution is a novel architecture, descríbing the architecture fully might suffice, or if the contribution is specific model and empirical evaluation, it may be necessary to either make possible for others to replicate the model with the same dataset, provide access to the model. In general. releasing code and data is often one good way to accomplish this, but reproducibility can also be provided via detailed instructions for how to replicate the results, access to a hosted model (e.g.. g., in the case of a large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurIPS does not require releasing code, the conference does require submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) the contribution is primarily a new algorithm, the paper should make clear how reproduce that algorithm. (b) If the contribution is primarily new model architecture, the paper should describe the architecture clearly and fully.
(c) the contribution is a new model (e-g., a large language model), then there should either be way to access this model for reproducing the results or way to reproduce the model (e.g., with an open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which case authors are welcome describe the particular way they provide for reproducibility. In the case of closed-source models, it may be that access to the model is limited in some way (e.g., to registered users), but it should possible for other researchers to have some path to reproducing verifying the results.
5. Open access to data and code
15

---

Does the paper provide open access to the data and code, with sufficient instruc- tions faithfully reproduce the main experimental results, as described in supplemental material?
Question:
Answer: [Yes]
Justification:
The accompanying code can used to reproduce all our results.
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips.co
for more details.
While we encourage the release of code and data, we understand that this might not be possible, "No" an acceptable answer. Papers cannot be rejected simply for including code, unless this is central to the contribution (e.g., for a new open-source benchmark).
The instructions should contain the exact command and environment needed to run to reproduce the results. See the NeurlPS code and data submission guidelines (https
/ praips.cc/ubblicgguuuee CoodeSumeeieeeeee ssi iessonPolll mollicyy for more details.
The authors should provide instructions on data access and preparation, including how to access the raw data, preprocessed data, intermediate data, and generated data, etc.
The authors should provide scripts to reproduce all experimental results for the new proposed method and baselines. If only subset of experiments are reproducible, they should state which ones are omitted from the script and why.
At submission time, to preserve anonymity, the authors should release anonymized versions (if applicable).
Providing much information possible in supplemental material (appended to the paper) recommended, but including URLs to data and code is permitted.
6. Experimental setting/details
Question: Does the paper specify all the training and test details (e.g., data splits, hyper- parameters, how they were chosen, type of optimizer, etc.) necessary to understand the results?
Answer: [Yes]
Justification: We provide the Appendix an entire section dedicated to all the experimental setting.
Guidelines:
The answer NA means that the paper does not include experiments.
The experimental setting should be presented in the core of the paper to a level of detail that is necessary appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, or as supplemental material.
7. Experiment statistical significance
Question: Does the paper report error bars suitably and correctly defined or other appropriate information about the statistical significance of the experiments?
Answer: [Yes]
Justification: When applicable, we always add error bars or boxplots on our experiments. We also describe the factors of variability.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by error bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
The factors of variability that the error bars are capturing should be clearly stated (for example, train/test split, initialization, random drawing of some parameter, or overall with given experimental conditions).
16

---

The method for calculating the error bars should be explained (closed form formula, call a library function, bootstrap, etc.)
The assumptions made should be given (e.g.. Normally distributed errors).
It should be clear whether the error bar is the standard deviation or the standard error of the mean.
is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma enror bar than state that they have 96% CI, if the hypothesis of Normality of errors is not verified.
For asymmetric distributions, the authors should be careful not to show in tables or figures symmetric error bars that would yield results that are range (e.g. negative enror rates).
If error bars reported in tables plots, The authors should explain in the text how they were calculated and reference the corresponding figures or tables the text.
8. Experiments compute resources
For each experiment, does the paper provide sufficient information on the com- puter resources (type of compute workers, memory, time of execution) needed reproduce the experiments?
Question:
Answer: [Yes]
Justification: We detail the Appendix the computer resources we used for all our experi- ments.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount of compute required for each of the individual experimental runs as well as estimate the total compute.
The paper should disclose whether the full research project required more compute than the experiments reported in the paper (e.g., preliminary or failed experiments that didn't make it into the paper).
9. Code of ethics
Question: Does the research conducted the paper conform, in every respect, with the NeurlPS Code of Ethics htttps://neut wttpiippp..pppp//pp httpp///bc.//...ccccsseeeeeeeeee
Answer: [Yes]
Justification: Our work has potential harmful consequence, and as far as data concerned, we only used synthetic and publicly available datasets.
Guidelines:
The answer NA means that the authors have not reviewed the NeurlPS Code of Ethics.
If the authors answer No, they should explain the special circumstances that require deviation from the Code of Ethics.
The authors should make sure to preserve anonymity (e.g., there a special consid- eration due to laws or regulations in their jurisdiction).
10. Broader impacts
Question: Does the paper discuss both potential positive societal impacts and negative societal impacts of the work performed?
Answer: [NA]
Justification: NA.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
If the authors answer NA or No, they should explain why their work has no societal impact why the paper does not address societal impact.
17

---

Examples of negative societal impacts include potential malicious unintended uses (e.g., disinformation, generating fake profiles, surveillance), faireess considerations (e.g., deployment of technologies that could make decisions that unfairly impact specific groups), privacy considerations, and security considerations.
The conference expects that many papers will be foundational research and not tied to particular applications, let alone deployments. However, if there is direct path to any negative applications, the authors should point it out. For example, it is legitimate to point out that improvement in the quality of generative models could be used to generate deepfakes for disinformation. On the other hand, it is not needed to point out that generic algorithm for optimizing neural networks could enable people to train models that generate Deepfakes faster.
The authors should consider possible harms that could arise when the technology is being used intended and functioning correctly, harms that could arise when the technology is being used as intended but gives incorrect results, and harms following from (intentional unintentional) misuse of the technology.
If there are negative societal impacts, the authors could also discuss possible mitigation strategies (e.g., gated release of models, providing defenses in addition to attacks, mechanisms for monitoring misuse, mechanisms monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper deseribe safeguards that have been put in place for responsible release of data or models that have a high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: We only use synthetic datasets.
Guidelines:
The answer NA means that the paper poses no such risks.
Released models that have high risk for misuse or dual-use should be released with necessary safeguards to allow for controlled use of the model, for example by requiring that users adhere to usage guidelines or restrictions to access the model or implementing safety filters.
Datasets that have been scraped from the Internet could pose safety risks. The authors should describe how they avoided releasing unsafe images.
We recognize that providing effective safeguards is challenging, and many papers do not require this, but we encourage authors to take this into account and make best faith effort.
12. Licenses for existing assets
Are the creators or original owners of assets (e.g., code, data, models), used in the paper, properly credited and are the license and terms of use explicitly mentioned and properly respected?
Question:
Answer: [NA]
Justification:
We do not rely on existing assets.
Guidelines:
The answer NA means that the paper does use existing assets.
The authors should cite the original paper that produced the code package dataset.
The authors should state which version of the asset used and, possible, include URL.
The name of the license (e.g., CC-BY 4.0) should be included for each asset.
For scraped data from a particular source (e.g., website). the copyright and terms of service of that source should be provided.
If assets are released, the license, copyright information, and terms use in package should be provided. For popular datasets, paperswithcone com/ datasets has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
18

---

For existing datasets that are re-packaged, both the original license and the license of the derived asset (if has changed) should be provided.
If this information is not available online, the authors are encouraged to reach out to the asset's s creators.
13. New assets
Question: Are new assets introduced in the paper well documented and is the documentation provided alongside the assets?
Answer: [NA]
Justification: NA.
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/model as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent was obtained from people whose asset used.
At submission time, remember to anonymize your assets (if applicable). You can either create anonymized URL or include anonymized zip file.
14. Crowdsourcing and research with human subjects
Question: For crowdsourcing experiments and research with human subjects, does the paper include the full text of instructions given to participants and screenshots, if applicable, as well as details about compensation (if any)?
Answer: [NA]
Justification: NA.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information in the supplemental material is fine, but if the main contribu- tion of the paper involves human subjects, then as much detail as possible should be included in the main paper.
According the NeurIPS Code of Bthies, workers involved in data collection, curation, or other labor should be paid at least the minimum wage in the country of the data collector.
15. Institutional review board (IRB) approvals or equivalent for research with human subjects
Question: Does the paper describe potential risks incurred by study participants, whether such risks were disclosed to the subjects, and whether Institutional Review Board (IRB) approvals (or equivalent approval/review based the requirements of your country or institution) were obtained?
Answer: [NA]
Justification: NA.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Depending on the country in which research is conducted, IRB approval (or equivalent) may be required for any human subjects research. If you obtained IRB approval, you should clearly state this in the paper.
We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors to adhere to the NeurlPS Code of Ethics and the guidelines for their institution.
For initial submissions, not include any information that would break anonymity (if applicable), such the institution conducting the review.
19

---

16. Declaration of LLM usage
Question: Does the paper describe the usage of LLMs if it is an important, original, or non-standard component of the core methods in this research? Note that if the LLM is used only for writing, editing, or formatting purposes and does not impact the core methodology. scientific rigorousness, originality of the research, declaration is not required.
Answer: [NA]
Justification: NA.
Guidelines:
The answer NA means that the core method development in this research does not involve LLMs as any important, original, non-standard components.
Please refer to L.MM policy (https: hate:p///eerrippc.cc////feeerers for what should or should not be described.
20


Reviews:
[{"review_id": "2nTxol6nXr", "reviewer": "Reviewer_3PxP", "summary": "Conformal prediction provides a framework to construct prediction bands with coverage guarantees. This paper builds on adaptive prediction band methods for regression tasks. Specifically, the paper targets coverage and adaptivity using reproducible kernel Hilbert spaces and kernel sum-of-squares methods. It derives a dual formulation for efficient optimization. Further, it tunes hyperparameters to provide better adaptivity.", "questions": "__Questions and Comments__\n\n1.\tHow is the local coverage metric computed?\n\n2.\tWhat is the reason for opting to have a 100% coverage constraint on the training data (Eq. 5)? Is it for convexity reasons only?\n\n3.\tIt is not easy to read the legends of the figures in the main paper.\n\n4.\tIt is worth highlighting that the width in lines 91-92 is constant for a fixed train and calibration data, not otherwise.\n\n5.\tIt is worth highlighting that the conformal prediction procedure detailed in Section 2.1 assumes $\\alpha \\geq 1/(m+1)$.\n\n6.\tThe paper uses the nuclear and Frobenius norms for the first time in Eq. 2. It would be helpful to define the two norms there.\n\n__Typos__\n\n1.\t\"Then Equation (4) admits a unique solution...\" $\\rightarrow$ \"Then the problem defined in Equations (4), (5), and (6) admits a unique solution...\" (lines 219-220).", "limitations": "A more thorough discussion on limitations would benefit the paper.", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 3, "originality": 3, "strengths_and_weaknesses": "__Strengths__\n\n1.\tThe paper studies the problem of adaptive prediction bands for regression problems, which is of practical relevance.\n\n2.\tThe paper is well written and is easy to understand.\n\n3.\tThe theoretical results are promising. They are also well explained in the paper, providing a good intuitive understanding of the same.\n\n__Weaknesses__\n\n1.\tThe experiments include results on synthetic datasets, good for analyzing the behavior of the proposed algorithm. However, experiments on real-world datasets are highly encouraged and will add to the comparison between the proposed algorithm and baseline methods.", "was_revised": true, "final_justification": "The authors have included experimental results on real-world datasets. I would encourage them to explore a few more UCI datasets (potentially with larger dimensions) for the camera-ready version.", "rebuttal": "Thank you very much for your careful reading and insightful feedback. In what follows, we address each of the points you raised.\n\n#### **Experiments on real-world datasets**\n\nWe greatly appreciate your positive assessment as well as the constructive criticism regarding the lack of experiments on real-world datasets, as was pointed out by all reviewers.\nFollowing Reviewer bubu’s suggestion, we have performed additional numerical experiments on 3 datasets from Romano et al. 2019: concrete (dimension 8), bio (dimension 9, with 1000 training points) and bike (dimension 18, with 1000 training points). The results of kSoS versus CQR, homoscedastic and heteroscedastic GPs is given below, in terms of mean width following Romano et al. 2019 (mean $\\pm$ sd):\n\n| **Case**   | CQR                 | Het GP                  | Hom GP                 | kSoS                 |\n|------------|---------------------|--------------------------|-------------------------|----------------------|\n| Concrete   | $0.586 \\pm 0.032$   | $\\mathbf{0.508} \\pm 0.052$ | $0.543 \\pm 0.044$     | $0.556 \\pm 0.044$    |\n| Bike       | $1.114 \\pm 0.062$   | $1.000 \\pm 0.079$         | $\\mathbf{0.809} \\pm 0.024$ | $\\mathbf{0.804} \\pm 0.032$ |\n| Bio        | $\\mathbf{1.879} \\pm 0.046$ | $2.21 \\pm 0.100$       | $2.194 \\pm 0.119$     | $2.017 \\pm 0.067$         |\n\nOn the Concrete dataset, heteroscedastic GP clearly outperforms all other methods, and kSoS is in second place with homoscedastic GPs. On the contrary, for the Bike dataset, homoscedastic GP and kSoS provide much thinner bands than competitors. Finally, the Bio dataset is the only case study where CQR performs much better than other methods. This is due to the asymmetric nature of the noise. With our asymmetric extension of kSoS, we are able to build better prediction bands than both GPs, even though we do not reach the performance of CQR. All in one, kSoS is systematically among the best methods, whatever the dataset. We hope that these new results demonstrate that our method is both effective and generalizes well across diverse settings.\n\n#### **Discussion on limitations**\n\nFollowing your comment on adding a more thorough discussion on limitations (which is also shared by all reviewers), we completely agree, and propose to add a new paragraph before the conclusion. To summarize, we will mention the following limitations:\n1. Kernel SoS performs well for dimensions up to 10 or 15, a limitation shared by all kernel methods without specific strategies. For high-dimensional objects with structure (time series, probability distributions, molecules, strings, graphs,….), it is however possible and straightforward to plug specifically designed kernels in our framework. On the other hand, for high-dimensional tabular data, a common strategy is to rely on additive kernels: we plan to investigate how they could be incorporated in the kernel SoS framework in future work.\n2. We showed in our numerical experiments that kernel SoS can scale up to 2000 samples thanks to the dual formulation, but handling much more samples than that would be highly time consuming. This is inherently a limitation which comes more from the kernel point of view than from the dual formulation. Indeed, each objective function and gradient computation in the dual formulation involves the usual kernel regression formulas, which scale as $O(n^3$). A usual workaround is to use Nystrom-type approximations or Random Fourier Features, which could be easily integrated in our problem. We plan to study how they perform on datasets with large number of samples in future work.\n3. Due to the 100% coverage constraint, kernel SoS is less robust to outliers than CQR. In such cases, a pre-processing step using another model (e.g. the initial GP model we use in our procedure) is recommended.\n4. Our findings on the hyperparameter influence only hold for the symmetric noise case (for $a$, as demonstrated in the supplemental, but also for the other ones). Additional intensive benchmarks should be conducted on asymmetric datasets.\n\n#### **Answers to specific questions**\n\n1. The local coverage metric computation is detailed in Section B.4 of the supplemental, but we recap the details in what follows. For a given seed (training/calibration/test dataset fixed), we generate $n_X=100$ additional independent samples $X_{n+1}$. For each of these samples, we generate $n_Y=1000$ samples from $Y_{n+1}$, and compute the empirical frequency of the event $Y_{n+1}\\in\\hat{C}(X_{n+1})$. This serves as an estimator of the local coverage at $X_{n+1}$. Finally, the empirical frequencies for all samples $X_{n+1}$ and all seeds are aggregated to draw the density of local coverage.\n2. You are completely right: 100\\% coverage is chosen only for convexity reasons. On a side note, a few months ago we would have also said that a smaller value would have led to a combinatorial problem which is much harder to optimize, but the very recent paper of Braun et al. (2025) seems to indicate that this is not necessarily the case.\n3. Thank you for pointing this out, we will adjust the legend size.\n4. Thank you for this remark, we will modify the sentence accordingly.\n5. This is true, we will modify the sentence accordingly.\n6. Sorry about this mistake, we will introduce the notations for the norms right after Eq. 2."}, {"review_id": "IBEtoisS23", "reviewer": "Reviewer_4Zih", "summary": "Conformal Prediction (CP) is a popular framework for uncertainty quantification in machine learning, valued for its ability to provide valid, non-asymptotic coverage guarantees regardless of the data distribution or underlying model. However, standard CP methods do not adapt prediction intervals to individual inputs, which has motivated recent work on adaptive extensions. This paper introduces a new adaptive CP method based on kernel sum-of-squares, a technique for approximating nonnegative functions in a reproducing kernel Hilbert space. The key idea is to learn the CP score function by solving an optimization problem regularized with a kernel sum-of-squares term (eqs. 4 to 6), which enables the construction of adaptive prediction bands. The authors generalize the representer theorem from Marteau-Ferey et al. (2020) and Muzellec et al. (2022) to their CP setting (Theorem 2), and derive a dual formulation that can be efficiently optimized with gradient-based methods, making their approach scalable to large datasets. They also study the influence of each hyperparameter and propose a novel tuning practical method using the Hilbert-Schmidt Independence Criterion (HSIC) to target improved marginal coverage in local regions of the test point, based on the theoretical bound in Proposition 4. Finally, they empirically evaluate their method in synthetic regression problems and compare against a few methods (Conformalized Quantile regression with random forests, and rescaled scores with Gaussian Processes) in terms of mean width, mutual information and local coverage.", "questions": "Given my comments in \"Strengths And Weaknesses\", I would be glad if the authors could address the following.\n\n**Q1.** On clarifying Proposition 3: why is the exchangeability assumption preserved when applying a standard split conformal scheme to their adaptive learned score function? How does this justify bypassing more elaborate calibration procedures used in prior work?\n\n**Q2.** On the validity of the HSIC-based bound (Proposition 4): how does the dependence of this bound on the kernel and the parameter $\\alpha_2$ influence its applicability? Are there cases where the square root constraint limits the practical use of this result?\n\n**Q3.** What kernel function was used in the experiments, and how sensitive are the results to this choice?\n\n**Q4.** Could the authors discuss and illustrate the relevance of their method to more realistic or complex regression tasks?", "limitations": "While the authors offer a detailed discussion of their method, the paper would benefit from a more explicit acknowledgment of its limitations, maybe through a short dedicated paragraph toward the end.\n\nConcerning potential negative societal impacts, I did not find a specific discussion on this topic. However, I think the content of this paper does not raise any major concerns in this regard.", "rating": 4, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 2, "originality": 2, "strengths_and_weaknesses": "### Strengths\n\n- The paper is very well-written, with a clear motivation that aligns well with the interests of the NeurIPS community\n- The theoretical contributions and the derived methodology appear sound and well-justified\n- I particularly appreciate the authors' efforts in providing an in-depth discussion of their framework, including Section 3.2, where they examine the influence and interactions of the multiple hyperparameters. This reflects their commitment to understand the method's behavior and improving its interpretability. Additionally, their analysis leads to novel theoretical results with potential relevance beyond the scope of this work. For example, the discussion on hyperparameter tuning directly motivates Proposition 4, which establishes weak conditional validity based on a lower bound involving HSIC, and could be of interest for other adaptive CP methods\n- The empirical results show that the proposed CP method is promising\n\n### Weaknesses\n\n1) Some aspects of the theoretical coverage guarantees would benefit from further clarification:\n\na) The proof of Proposition 3 is missing. While this suggests that the result is a trivial extension of standard CP arguments to their adaptive setting, it would be helpful for the authors to explicitly explain why this is the case. In particular, they should clarify how the exchangeability of the learned scores is preserved, and why their approach based on a simple split conformal scheme suffices, whereas prior related works have relied on more elaborate calibration procedures [Liang, 2022; Fan et al., 2024].\n\nb) In Proposition 4, HSIC is actually subject to a constraint to ensure that the square root expression remains defined (i.e., its argument is nonnegative). Since this constraint depends on the value of $\\alpha_2$, and hence on the kernel choice, it would be useful to discuss whether this restricts the applicability or effectiveness of the theoretical bound. If I am not mistaken, this aspect is not addressed in the paper, and it raises the question of whether there exist scenarios in which the proposed methodology (based on maximizing HSIC) may be affected.\n\n\n2) The experimental evaluation is promising but limited in scope:\n\na) Some implementation details could be clarified. For instance, the choice of kernel used in their experiments is not given, and the sensitivity of the results to this choice is not discussed. Additionally, while Figure 3 claims scalability to \"thousands of training points\", Figure 3 left goes up to $n=1000$, which limits the strength of this claim.\n\nb) More importantly, the experiments are conducted only on low-dimensional (mostly $d=1$ and up to $d=5$) synthetic regression tasks with relatively small to moderate sample sizes. While the main text reports figures showing that the proposed method performs the best, the supplementary material questions this conclusion. For example, in Cases 3 and 4, the kernel SoS CP does not always outperform competing methods. In particular, in Case 3, heteroscedastic GP performs better than kernel SoS; in Case 4 (or Case 2 with larger dimensions, in light of Remark 2) where the ideal prediction intervals are constant, kernel SoS appears to struggle in capturing this, which seems a bit counterintuitive as I would expect its adaptivity to naturally account for such case. This makes me question the relevance of kernel SoS for more realistic or complex benchmarks.", "was_revised": true, "final_justification": "The authors have addressed my main concerns by providing additional clarifications, empirical results, and discussions. Their explanations resolved my initial misunderstandings (e.g., 1.b), while the added details and experiments (e.g., implementation details, low-dimensional benchmarks, real-world datasets) improved the overall evaluation. The proposed addition of a section on limitations and future directions also strengthens the paper's exposition. \n\nOverall, I find the response satisfactory, the contributions clearer, and the empirical evidence more convincing. However, while the proposed framework is promising and flexible enough to foster future improvements and extensions (e.g., the independence test introduced during the rebuttal), its practical relevance in realistic scenarios remains more questionable to me. This justifies maintaining my current score: I vote for acceptance.", "rebuttal": "We appreciate your positive remarks and helpful feedback. We address your concerns and questions point by point below.\n\n1. \"Some aspects of the theoretical coverage guarantees would benefit from further clarification\".\n   \t\n\ta) The proof of Proposition 3 is indeed missing because this is a trivial extension of standard split CP: our kernel SoS procedure is conducted on training data only, and thus training conditional coverage is straightforward, yielding marginal coverage (see for example Lei et al. 2018). We actually followed the presentation of the boosting paper of Xie et al. (2024). We will explicitly mention this fact in the paper, thank you for highlighting that — it is important to clarify to avoid any misunderstanding.\n\n\tRegarding the differences between our split CP procedure and Fan et al. (2023), actually they are very slim. Indeed, the calibration procedure in Fan et al. (2023) is also based on a splitting strategy as for split CP, but without adjusting the quantile level. Their calibration method is thus almost equivalent to ours, and from a theoretical perspective they recover the target coverage level but with different tools (which, as a side note, allows them to derive explicit coverage bounds for their bands before calibration). The main advantage of using split CP as in our framework, is that we can benefit from all previous and future extensions of CP to handle specific situations (e.g. covariate shift, jackknife+, ...) at almost no cost, unlike Fan et al.'s strategy. On the other hand, Liang (2022) also provides pre-calibration coverage bounds, while relying on sample splitting. The main variation with respect to split CP or Fan et al. (2023) is that the hyperparameter which must be optimized on the calibration set no longer appears as a multiplicative constant in front of the band width, but as an additive constant inside a square root, thus implying a more complex dyadic search (in contrast to our framework or Fan et al.'s one, where it is simply given as a sample quantile). This simple difference partly explains why they require other proof strategies to get the target coverage. Once again, we believe that the traditional split CP is more straightforward and easier to apprehend, with no obvious disadvantage (to the best of our knowledge) with respect to Fan et al. (2023) or Liang (2022) procedures.\n\n\tb) Concerning the HSIC term appearing in the square-root, looking closely at the proof in Section A.3 of the supplemental, the inside term is always nonnegative. Indeed,\n\t$$0 \\leq 1-\\exp(-\\textrm{MI}(X,S)) \\leq 1 - \\alpha_1\\exp(\\textrm{MI}(V,R)) \\leq 1 - \\frac{\\alpha_1}{\\alpha_2 \\textrm{HSIC}(V,R)}$$\n\twhatever $\\alpha_2$ and the kernel choice for HSIC (as long as it is characteristic for the rightmost inequality to hold). Note that for this task, we choose here the energy distance kernel: it is characteristic, it does not depend on additional hyperparameters (unlike the Gaussian or Matérn one), and has been used with great success for MMD and HSIC (Mak & Joseph 2018, Da Veiga 2021).\n\n2. \"The experimental evaluation is promising but limited in scope\".\n   \n\ta) Although we detailed the kernel used for HSIC, we forgot to mention the one used for $K^m$ and $K^f$, we are very sorry about that. We actually use a Matérn 5/2 kernel since this is the default choice in most GP software, we will add this information in the paper. We did not consider other kernels in our study, because our main point was to thoroughly investigate all other hyperparameters in the procedure, with a particular focus on $\\theta^f$.\n\n\tb) We agree that we only considered here synthetic examples in low dimension, but this was a necessary intensive benchmark to really grasp how hyperparameters interact with each other in our new kernel SoS proposal. Regarding the results in the supplemental, actually for Case 3 kSoS and heteroscedastic GP perform similarly, with equivalent metrics. It may seem that the density of local coverage is better for heteroscedastic GP because its mode coincides with 0.9, but it tends to undercover while kSoS tends to overcover.\n\n\tHowever, it is true that for Case 4 (and all cases which are almost homoscedastic), our results show that kSoS struggles more. This is due to the fact that, while it is always possible to find the highest value of HSIC, this value may not be significantly different from 0, which indicates that the residuals and $f$ are independent. In this situation, we should enforce an homoscedastic kSoS model, i.e. choose an arbritary large $\\theta^f$. To detect such cases, we now propose to perform a HSIC test of independence on the best $\\theta^f$ at the end of our procedure and select this best value only if the null hypothesis is rejected. If we use now our HSIC test of independence on the best $\\theta^f$ and force an homoscedastic model if the test does not reject independence, we actually obtain on Case 4 new results that are almost identical (as expected) to the homoscedastic GP, which is the best performing method. We will replace Figures 16 and 18 with our new procedure including this final test, and describe the independence test in Section B.2.\n\n\tFinally, regarding the lack of more realistic cases, we answer your comment below.\n\n\n#### **Experiments on real-world datasets**\n\nWe greatly appreciate your positive assessment as well as the constructive criticism regarding the lack of experiments on real-world datasets, as was pointed out by all reviewers.\nFollowing Reviewer bubu’s suggestion, we have performed additional numerical experiments on 3 datasets from Romano et al. 2019: concrete (dimension 8), bio (dimension 9, with 1000 training points) and bike (dimension 18, with 1000 training points). The results of kSoS versus CQR, homoscedastic and heteroscedastic GPs is given below, in terms of mean width following Romano et al. 2019 (mean $\\pm$ sd):\n\n| **Case**   | CQR                 | Het GP                  | Hom GP                 | kSoS                 |\n|------------|---------------------|--------------------------|-------------------------|----------------------|\n| Concrete   | $0.586 \\pm 0.032$   | $\\mathbf{0.508} \\pm 0.052$ | $0.543 \\pm 0.044$     | $0.556 \\pm 0.044$    |\n| Bike       | $1.114 \\pm 0.062$   | $1.000 \\pm 0.079$         | $\\mathbf{0.809} \\pm 0.024$ | $\\mathbf{0.804} \\pm 0.032$ |\n| Bio        | $\\mathbf{1.879} \\pm 0.046$ | $2.21 \\pm 0.100$       | $2.194 \\pm 0.119$     | $2.017 \\pm 0.067$         |\n\nOn the Concrete dataset, heteroscedastic GP clearly outperforms all other methods, and kSoS is in second place with homoscedastic GPs. On the contrary, for the Bike dataset, homoscedastic GP and kSoS provide much thinner bands than competitors. Finally, the Bio dataset is the only case study where CQR performs much better than other methods. This is due to the asymmetric nature of the noise. With our asymmetric extension of kSoS, we are able to build better prediction bands than both GPs, even though we do not reach the performance of CQR. All in one, kSoS is systematically among the best methods, whatever the dataset. We hope that these new results demonstrate that our method not only performs well in the original synthetic examples but also generalizes well on other datasets.\n\n#### **Discussion on limitations**\nAs all other reviewers, you mention that our paper could benefit from a separate limitation section. We completely agree, and propose to add a new paragraph before the conclusion. To summarize, we will mention the following limitations:\n1. Kernel SoS performs well for dimensions up to 10 or 15, a limitation shared by all kernel methods without specific strategies. For high-dimensional objects with structure (time series, probability distributions, molecules, strings, graphs,….), it is however possible and straightforward to plug specifically designed kernels in our framework. On the other hand, for high-dimensional tabular data, a common strategy is to rely on additive kernels: we plan to investigate how they could be incorporated in the kernel SoS framework in future work.\n2. We showed in our numerical experiments that kernel SoS can scale up to 2000 samples thanks to the dual formulation, but handling much more samples than that would be highly time consuming. This is inherently a limitation which comes more from the kernel point of view than from the dual formulation. Indeed, each objective function and gradient computation in the dual formulation involves the usual kernel regression formulas, which scale as $O(n^3$). A usual workaround is to use Nystrom-type approximations or Random Fourier Features, which could be easily integrated in our problem. We plan to study how they perform on datasets with large number of samples in future work.\n3. Due to the 100% coverage constraint, kernel SoS is less robust to outliers than CQR. In such cases, a pre-processing step using another model (e.g. the initial GP model we use in our procedure) is recommended.\n4. Our findings on the hyperparameter influence only hold for the symmetric noise case (for $a$, as demonstrated in the supplemental, but also for the other ones). Additional intensive benchmarks should be conducted on asymmetric datasets.\n\n#### **References**\n\nDa Veiga, S. (2021). Kernel-based ANOVA decomposition and Shapley effects--Application to global sensitivity analysis. arXiv preprint arXiv:2101.05487.\n\nMak, S., & Joseph, V. R. (2018). SUPPORT POINTS. The Annals of Statistics, 46(6A), 2562-2592.\n\nLei, J., G’Sell, M., Rinaldo, A., Tibshirani, R. J., & Wasserman, L. (2018). Distribution-free predictive inference for regression. Journal of the American Statistical Association, 113(523), 1094-1111."}, {"review_id": "gAyrnHax2O", "reviewer": "Reviewer_bubu", "summary": "This paper addresses the lack of adaptivity problem in conformal prediction by reframing the construction of prediction intervals as a constrained learning problem in RKHS. Specifically, the authors introduce adaptivity to the score function by estimating the uncertainty measure ($\\sigma^2$) through kernel sum of squares, which is a natural technique for positive function estimation. The new algorithm is termed, kernel sum of squares (kSoS). The authors derive a new representor theorem for the band prediction problem, which shows the unique solution is supported on the finite samples. Dual formulations are provided, which can be solved through accelerated gradient descent, therefore achieving significant speed up. Furthermore, by linking local coverage to HSIC, the authors provide a principled way for hyperparameter tuning (kernel bandwidth). Synthetic data simulation is provided to show the efficacy of the proposed kSOS algorithm.", "questions": "Please see weakness.", "limitations": "NA", "rating": 4, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 3, "clarity": 3, "significance": 2, "originality": 2, "strengths_and_weaknesses": "Strength\n1. The paper is mostly well written and easy to follow.\n2. Strong theoretical results are provided. The new representor theorem ensures that the algorithm can operate under finite dimensional case. Dual formulation achieves impressive speed up. The HSIC based hyperparameter tunning also intuitively makes sense.\n\nWeakness:\n1. (Major) The paper only have a set of synthetic data simulation and lacks real data simulation. One suggestion is to add a small set of real data simulation comparison just like [1]\n[1] Romano, Yaniv, Evan Patterson, and Emmanuel Candes. \"Conformalized quantile regression.\" Advances in neural information processing systems 32 (2019).\n2. (Major) Following weakness 1, the actual merit of Prop. 4 and the hyperparameter tuning strategy is questionable. Since the HSIC is simply appearing as part of the lose lower bound for local coverage, there is no guarantee that optimizing $\\theta_f$ according to HSIC will not harm the actual performance. With very limited synthetic data simulation, the claim that such a tuning scheme would work is not well supported.\n3. (Minor) The choice of hyperparameters (other than tunable parameters) are discussed in Section 3.2 with claims supported in Section 4. However, with the lack of real data simulations, it would be difficult to judge the claims. (For example, b=10 is sufficiently large)\n4. (Minor) Notation definitions can be improved. For example, nuclear and Frobenius norm are defined after Eq. (4), but used as early as Eq. (2). Hyperparameter a and b are not defined.", "was_revised": true, "final_justification": "I believe the real data simulation response is satisfactory.", "rebuttal": "We appreciate your positive remarks and constructive review. We clarify the points you raised in what follows.\n\n- \"(Major) The paper only have a set of synthetic data simulation and lacks real data simulation\".\n  \n  Following your helpful suggestion, we have performed additional numerical experiments on 3 datasets from Romano et al. 2019, see the details in the dedicated section below.\n- \"(Major) Following weakness 1, the actual merit of Prop. 4 and the hyperparameter tuning strategy is questionable\".\n  \n  We are not sure to understand your comment \"there is no guarantee that optimizing $\\theta^f$ according to HSIC will not harm the actual performance\". First, since to the best of our knowledge there exists no other data-driven tuning strategy for $\\theta^f$ in the literature, there is no baseline performance available for kSoS, therefore optimizing the HSIC cannot « harm » anything. But following your comment, we carried out our benchmark again with two naive values for $\\theta^f$: a) taking $\\theta^f=\\theta^m$ from the GP mean and b) taking $\\theta^f$ equal to the median of the feature distances (a usual rule-of-thumb - ROT below - for kernel methods). The results are given in the table below (mean $\\pm$ sd). In general HSIC yields mean widths similar or better than the naive strategies, but importantly it produces much more adaptive bands since the MI is significantly lower.\n\n|           | Mean-Width  | Mean-Width  | Mean-Width  |  MI| MI | MI |\n|-----------|----------------------------|--------------------------------|-------------------------------------|----------------------------|--------------------------------|-------------------------------------|\n|           | $\\theta^m$ | $\\textrm{ROT}$ | $\\max\\textrm{HSIC}$ |  $\\theta^m$ | $\\textrm{ROT}$ | $\\max\\textrm{HSIC}$ |\n| **Case 1** | $3.097 \\pm 0.583$          | $\\mathbf{2.992} \\pm 0.167$              | $\\mathbf{2.991} \\pm 0.279$                | $0.076 \\pm 0.072$          | $0.065 \\pm 0.031$              | $\\mathbf{0.058} \\pm 0.037$                |\n| **Case 2** | $2.307 \\pm 0.152$          | $2.082 \\pm 0.200$              | $\\mathbf{1.999} \\pm 0.122$               | $0.307 \\pm 0.027$          | $0.228 \\pm 0.057$              | $\\mathbf{0.136} \\pm 0.029$               |\n| **Case 3** | $1.486 \\pm 0.216$          | $1.419 \\pm 0.138$              | $\\mathbf{1.378} \\pm 0.129$                | $0.189 \\pm 0.053$          | $0.093\\pm 0.059$              | $\\mathbf{0.074} \\pm 0.033$                |\n| **Case 4** | $4.127 \\pm 0.248$          | $4.145 \\pm 0.254$              | $\\mathbf{4.044} \\pm 0.225$                  | $0.022 \\pm 0.018$          | $0.031 \\pm 0.023$              | $\\mathbf{0.013} \\pm 0.008$                  |\n\n  However, from the results on Case 4 in the supplemental, it is clear that HSIC optimization struggles for homoscedastic cases. This is due to the fact that, while it is always possible to find the highest value of HSIC, this value may not be significantly different from 0, which indicates that the residuals and $f$ are independent. In this situation, we should enforce an homoscedastic kSoS model, i.e. choose an arbritary large $\\theta^f$. To detect such cases, we now propose to perform a HSIC test of independence on the best $\\theta^f$ at the end of our procedure and select this best value only if the null hypothesis is rejected. We will detail this extra step in Section B.2 and update Figures 16 and 18 accordingly.\n\nFinally, the bound with HSIC may indeed be lose, but in our opinion this is not a practical issue, since we mainly see this bound as a mathematical justification for using HSIC over MI to measure dependence (and consequently ensure as much as possible local coverage), which has a major impact on practical robustness of the procedure. Of course using a lose bound would be detrimental if we were using it to actually compute the weak local coverage, but this is not the case.\n- \"(Minor) The choice of hyperparameters (other than tunable parameters) are discussed in Section 3.2 with claims supported in Section 4\".\n\n  With the addition of new experimental results on real datasets, I hope you will find our choice of hyperparameters more convincing. We will also add a comment about a potential limitation related to hyperparameter choices for the asymmetric case in a new paragraph, see details below.\n- \"(Minor) Notation definitions can be improved\".\n\n  Thank you for pointing out this mistake, we will introduce the norms right after Eq. 2.\n\n#### **Experiments on real-world datasets**\n\nFollowing your suggestion, we have performed additional numerical experiments on 3 datasets from Romano et al. 2019: concrete (dimension 8), bio (dimension 9, with 1000 training points) and bike (dimension 18, with 1000 training points). The results of kSoS versus CQR, homoscedastic and heteroscedastic GPs is given below, in terms of mean width following Romano et al. 2019 (mean $\\pm$ sd):\n\n| **Case**   | CQR                 | Het GP                  | Hom GP                 | kSoS                 |\n|------------|---------------------|--------------------------|-------------------------|----------------------|\n| Concrete   | $0.586 \\pm 0.032$   | $\\mathbf{0.508} \\pm 0.052$ | $0.543 \\pm 0.044$     | $0.556 \\pm 0.044$    |\n| Bike       | $1.114 \\pm 0.062$   | $1.000 \\pm 0.079$         | $\\mathbf{0.809} \\pm 0.024$ | $\\mathbf{0.804} \\pm 0.032$ |\n| Bio        | $\\mathbf{1.879} \\pm 0.046$ | $2.21 \\pm 0.100$       | $2.194 \\pm 0.119$     | $2.017 \\pm 0.067$         |\n\nOn the Concrete dataset, heteroscedastic GP clearly outperforms all other methods, and kSoS is in second place with homoscedastic GPs. On the contrary, for the Bike dataset, homoscedastic GP and kSoS provide much thinner bands than competitors. Finally, the Bio dataset is the only case study where CQR performs much better than other methods. This is due to the asymmetric nature of the noise. With our asymmetric extension of kSoS, we are able to build better prediction bands than both GPs, even though we do not reach the performance of CQR. All in one, kSoS is systematically among the best methods, whatever the dataset. We believe these new results enhance the original findings and provide a more comprehensive view of the method’s generalization ability or limitations. For the latter, we will add a new paragraph as proposed below.\n\n#### **Discussion on limitations**\n\nWe will add a new paragraph before the conclusion, where we will mention the following limitations:\n1. Kernel SoS performs well for dimensions up to 10 or 15, a limitation shared by all kernel methods without specific strategies. For high-dimensional objects with structure (time series, probability distributions, molecules, strings, graphs,….), it is however possible and straightforward to plug specifically designed kernels in our framework. On the other hand, for high-dimensional tabular data, a common strategy is to rely on additive kernels: we plan to investigate how they could be incorporated in the kernel SoS framework in future work.\n2. We showed in our numerical experiments that kernel SoS can scale up to 2000 samples thanks to the dual formulation, but handling much more samples than that would be highly time consuming. This is inherently a limitation which comes more from the kernel point of view than from the dual formulation. Indeed, each objective function and gradient computation in the dual formulation involves the usual kernel regression formulas, which scale as $O(n^3$). A usual workaround is to use Nystrom-type approximations or Random Fourier Features, which could be easily integrated in our problem. We plan to study how they perform on datasets with large number of samples in future work.\n3. Due to the 100\\% coverage constraint, kernel SoS is less robust to outliers than CQR. In such cases, a pre-processing step using another model (e.g. the initial GP model we use in our procedure) is recommended.\n4. Our findings on the hyperparameter influence only hold for the symmetric noise case (for $a$, as demonstrated in the supplemental, but also for the other ones). Additional intensive benchmarks should be conducted on asymmetric datasets."}, {"review_id": "aU9fVrRJRa", "reviewer": "Reviewer_ca1L", "summary": "This paper introduces a new method for building adaptive prediction intervals using conformal prediction. \nThe main idea is to learn both a mean function and a non-negative scaling function using a kernel sum-of-squares framework. \nThe method produces adaptive intervals that adjust to the local uncertainty in the data and still guarantee valid coverage. \nThe authors provide a solid theoretical foundation, a dual formulation that scales well, and a new tuning strategy based on HSIC to guide adaptivity.", "questions": "See Strengths and Weaknesses", "limitations": "The paper is missing a clear limitations section. The authors should briefly mention key assumptions and when the method might not work well. Even a short note would help readers understand the boundaries of the approach.", "rating": 4, "confidence": 5, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": 2, "clarity": 3, "significance": 2, "originality": 2, "strengths_and_weaknesses": "The paper's biggest strength is how well it connects ideas from different areas—conformal prediction, kernel methods, and structured low-rank learning—and turns them into something coherent. \nIt doesn't just throw these tools together. Each part is motivated: kernel SoS is used to learn a non-negative function, the split CP step gives finite-sample guarantees, and HSIC is brought in to tune adaptivity. \nThe method is scalable too, thanks to a dual formulation that avoids solving SDPs. \nAnd the hyperparameter section is thoughtful. The authors actually explain what each parameter does and how they interact.\n\nThat said, the evaluation is limited. \nAll experiments are done on synthetic data. \nThese examples are useful for showing how the method behaves under different types of noise, but they don’t really tell us how it handles the messiness of real-world data—things like outliers, covariate shifts, or higher-dimensional structure. Even a single real regression benchmark would have gone a long way. \nAs it stands, the results mostly validate the theory, not the practical utility.\n\nThere's also no discussion of limitations, as requested.\nSome of these issues are mentioned in passing (like the difficulty of local coverage), but there's no structured effort to acknowledge the boundaries of the approach. \nA short limitations section would have made the paper feel more complete and self-aware.", "was_revised": true, "final_justification": "The rebuttal adds a useful limitations section, which clarifies some boundaries of the proposed method and improves transparency. However, the main concern about the evaluation remains unresolved. The added experiments are limited to three datasets, with kSoS clearly outperforming baselines on only one of them, and results on Concrete are somewhat misrepresented. More importantly, coverage, central to conformal prediction, is not reported for real-world data, leaving a core claim unverified. Overall, while the method is technically solid and the theoretical contributions are interesting, the experimental evidence remains too limited and the conclusions overstated, which constrains the confidence in its practical impact.", "rebuttal": "We are grateful for your positive assessment and insightful comments. Below, we respond to each of your remarks in turn.\n\n#### **Experiments on real-world datasets**\n\nThe lack of experiments on real-world datasets was pointed out by all reviewers.\nFollowing Reviewer bubu’s suggestion, we have performed additional numerical experiments on 3 datasets from Romano et al. 2019: concrete (dimension 8), bio (dimension 9, with 1000 training points) and bike (dimension 18, with 1000 training points). The results of kSoS versus CQR, homoscedastic and heteroscedastic GPs is given below, in terms of mean width following Romano et al. 2019 (mean $\\pm$ sd):\n\n| **Case**   | CQR                 | Het GP                  | Hom GP                 | kSoS                 |\n|------------|---------------------|--------------------------|-------------------------|----------------------|\n| Concrete   | $0.586 \\pm 0.032$   | $\\mathbf{0.508} \\pm 0.052$ | $0.543 \\pm 0.044$     | $0.556 \\pm 0.044$    |\n| Bike       | $1.114 \\pm 0.062$   | $1.000 \\pm 0.079$         | $\\mathbf{0.809} \\pm 0.024$ | $\\mathbf{0.804} \\pm 0.032$ |\n| Bio        | $\\mathbf{1.879} \\pm 0.046$ | $2.21 \\pm 0.100$       | $2.194 \\pm 0.119$     | $2.017 \\pm 0.067$         |\n\nOn the Concrete dataset, heteroscedastic GP clearly outperforms all other methods, and kSoS is in second place with homoscedastic GPs. On the contrary, for the Bike dataset, homoscedastic GP and kSoS provide much thinner bands than competitors. Finally, the Bio dataset is the only case study where CQR performs much better than other methods. This is due to the asymmetric nature of the noise. With our asymmetric extension of kSoS, we are able to build better prediction bands than both GPs, even though we do not reach the performance of CQR. All in one, kSoS is systematically among the best methods, whatever the dataset. We hope that these new results demonstrate that our method not only performs well in the original synthetic examples but also generalizes well on other datasets.\n\nIn particular, we now tackle problems up to dimension 18 with 1000 training points (see our comments related to limitations below for additional discussion on some challenges you mentioned - outliers and high-dimensional problems). Finally, you also bring up covariate shift: since our method only operates at the score construction step of split CP, it can directly be generalized to the common covariate shift setting through weighted quantiles.\n\n#### **Discussion on limitations**\n\nWe agree that readers would profit from a dedicated limitation section so that they can grasp the boundaries of our approach. We will add a new paragraph before the conclusion, where we will mention the following limitations:\n1. Kernel SoS performs well for dimensions up to 10 or 15, a limitation shared by all kernel methods without specific strategies. For high-dimensional objects with structure (time series, probability distributions, molecules, strings, graphs,….), it is however possible and straightforward to plug specifically designed kernels in our framework. On the other hand, for high-dimensional tabular data, a common strategy is to rely on additive kernels: we plan to investigate how they could be incorporated in the kernel SoS framework in future work.\n2. We showed in our numerical experiments that kernel SoS can scale up to 2000 samples thanks to the dual formulation, but handling much more samples than that would be highly time consuming. This is inherently a limitation which comes more from the kernel point of view than from the dual formulation. Indeed, each objective function and gradient computation in the dual formulation involves the usual kernel regression formulas, which scale as $O(n^3$). A usual workaround is to use Nystrom-type approximations or Random Fourier Features, which could be easily integrated in our problem. We plan to study how they perform on datasets with large number of samples in future work.\n3. Due to the 100\\% coverage constraint, kernel SoS is less robust to outliers than CQR. In such cases, a pre-processing step using another model (e.g. the initial GP model we use in our procedure) is recommended.\n4. Our findings on the hyperparameter influence only hold for the symmetric noise case (for $a$, as demonstrated in the supplemental, but also for the other ones). Additional intensive benchmarks should be conducted on asymmetric datasets."}]

Author rebuttal:

Meta-review: The paper develops a kernel sum-of-squares framework for adaptive conformal prediction, introducing a dual optimization formulation for scalability, a new HSIC-based hyperparameter tuning strategy, and theoretical guarantees supported by empirical evaluation. The paper is clearly written, the methodology is well-motivated, and the combination of theory, optimization, and empirical evaluation makes for a solid and timely contribution. The additional experiments on real-world datasets strengthen the empirical case and show the method’s potential across different settings.

For the camera-ready version, the authors should report coverage metrics alongside band width for all datasets to make the empirical results more transparent. They should expand the limitations section, especially discussing scalability beyond a few thousand samples, robustness to outliers under the 100% coverage constraint, and performance under asymmetric noise. Clarifying cases where the method does not outperform alternatives, improving figure readability, and introducing notation earlier would further strengthen the final version.
Decision: Accept (poster)