Learning from Uncertain Data: From Possible Worlds to Possible Models
Authors: Jiongli Zhu
Su Feng
Boris Glavic
Babak Salimi
Venue: NeurIPS 2024 poster (2024)
Keywords: data uncertainty
robustness verification
predictive multiplicity
abstract interpretation
linear regression
Primary area: machine_learning_for_other_sciences_and_fields

Abstract: We introduce an efficient method for learning linear models from uncertain data, where uncertainty is represented as a set of possible variations in the data, leading to predictive multiplicity. Our approach leverages abstract interpretation and zonotopes, a type of convex polytope, to compactly represent these dataset variations, enabling the symbolic execution of gradient descent on all possible worlds simultaneously. We develop techniques to ensure that this process converges to a fixed point and derive closed-form solutions for this fixed point. Our method provides sound over-approximations of all possible optimal models and viable prediction ranges. We demonstrate the effectiveness of our approach through theoretical and empirical analysis, highlighting its potential to reason about model and prediction uncertainty due to data quality issues in training data.

Learning from Uncertain Data: From Possible Worlds to Possible Models
Jiongli Zhu
Su Feng"
Boris Glavic
Babak Salimi*
University of California, San Diego Nanjjngg Tech University
UUniversity of Illinois, Chicago
Abstract
We introduce an efficient method for learning linear models from uncertain data, where uncertainty is represented set of possible variations in the data, leading to predictive multiplicity. Our approach leverages abstract interpretation and zono- topes, a type of convex polytope, to compactly represent these dataset variations, enabling the symbolic execution of gradient descent on all possible worlds simul- taneously. We develop techniques to ensure that this process converges to a fixed point and derive closed- form- solutions for this fixed point. Our method provides sound per-approximaaionn of all possible optimal models and viable prediction ranges. We demonstrate the effectiveness of our approach through theoretical and empirical analysis, highlighting its potential to reason about model and prediction uncertainty due to data quality issues in training data.
Introduction
This paper addresses the challenges of learning from uncertain datasets by employing the framework of possible world semantics, a well-established concept in Al and database theory [59, 23, 41, 28). In this approach, uncertainty in a dataset D is conceptualized through collection of possible datasets {Di, D)....,, each representing a potential state of the real world, reflecting variations due missing entries, errors, inconsistencies. and biases. Given this framework and a learning algorithm, our objective is to construct a set of models {f1, f2.ff....., where each model fi is trained on corresponding potential dataset D,. This method, that we implement in system called ZORRO (Z0notope-basse Robustness Analysis for Regression with aver-apprxximations) allows for thorough evaluation of how data uncertainties affect the robustness, reliability, and fairness of models in predictive modeling and statistical inference, particularly in scenarios where the ground truth unidentifiable, necessitating consideration all possible dataset variations.
While the framework of possible world semantics is essential for modeling dataset uncertainties, poses significant challenges due the potentially infinite number of scenarios each dataset might represent. Exploring every possibility and training model for each is impractical. The concept of model multiplicity, which highlights situations where models with similar accuracy differ in individual predictions, has gained traction [7, 42], yet it primarily focuses on competing models without addressing the full range of dataset variations. Similarly, dataset multiplicity introduced in [47] recognizes data variations due to uncertainty but [47] only proposed a solution for linear models with label errors. Our approach expands these ideas by using possible world semanties to systematically manage uncertainty across all features and labels, thus creating a comprehensive framework for evaluating model robustness amid data uncertainties.
To address the challenges of learning from uncertain datasets, we employ the method of abstract interpretation [13]. Utilizing zonotopes- type of convex polytope well-suited for compactly representing high-dimensional data spaces [74, 4]-we over-approximate the set of all possible dataset variations. This framework allows for the simultaneous symbolic execution of gradient descent
38th Conference on Neural Information Processing Systems (NcurIPS 2024).

---

across all possible datasets, compactly over aee-apppoximaiing all possible optimal model parameters as zonotope. We demonstrate that for linear regression with & regularization (Ridge), our method admits a non-trivial closed-form solution. The zonotope representation of model parameters enables efficient inference, facilitating reasoning about the range of possible predictions or of specific parameters.
Contributions. The key contributions of this research are:
1. We introduce an abstract gradient descent algorithm for learning linear regression models from uncertain data. This method over-approximates data variations using zonotopes and symbolically executes gradient descent on all possible datasets concurrently. We define and prove the existence of fixed point that soundly over-approximates all potential models.
2. Symbolic execution generates intractable polynomial zonotopes for gradients due to non- lonelinear terms and monomial growth that is exponential in the number of iterations. We use lincarization and order reduction compactly over-approximate these polynomial zonotopes using linear zonotopes at each step, introducing an efficient version of abstract gradient descent.
3. The efficient version, however, does not guarantee a fixed point for arbitrary order reduction techniques. To address this, we develop advanced order reduction techniques that ensure fixed points and provide non-trivial closed- form solution for these fixed points in ridge regression. 4. We implement our approach in system called ZORRO and use it to evaluate the impact of data uncertainty on linear regression models. Our empirical results and analytical solutions validate the effectiveness of our approach, demonstrating its efficacy in computing prediction ranges and verifying robustness.
Related Work. Predictive multiplicity has shown that single dataset can produce multiple op- timally fitted models due to variations in training processes [8, 11, 44, 66, 69, 9, 15, 71, 42] or modifications to training parameters such as random seeds, data ordering, and hyperparameters. For predictive multiplicity due to missing data, Khosravi et al.[34] address the issue using a probabilistic method that computes expected predictions for all possible imputations. Our approach can be seen as an extreme case of multiple imputation [62, 61], where we consider all possible data variations rather than just a few plausible scenarios. Meyer et al. [47] recently introduced the concept of dataset multiplicity, using possible world semantics to model how uncertain, biased, or noisy training data can lead to predictive multiplicity. However, their focus is on uncertainty in training labels, and they use interval arithmetic for over-approximation of prediction intervals for linear regression. In contrast, our approach handles arbitrary uncertainty in features and labels during both training and testing using zonotope- based learning for approxxmmmtnonnn of prediction ranges and model parameters. We show that interval arithmetic fails to provide tight prediction ranges even for uncertainty in labels.
Our work is broadly related to robust model learning, which ensures robustness against data quality issues such as attacks [31, 76, 70, 60, 55, 52, 29]. Distributional robustness [6, 53, 65] studies model reliability against varying data distributions, while robust statisties [19, 18] examines model performance under outliers or data errors. Our approach provides exact provable robustness guarantees by exploring the entire range of models under extreme dataset variations, which is crucial for individual- levell predictions and reasoning about the robustness of specific parameters.
Our work is also related to robustness certification, which certifies ML models' robustness against data perturbations and uncertainties [27, 68, 50]. These efforts mainly focus on test-time robustness, validating predictions for inputs in the vicinity of test sample. In contrast, we address training- taining-time robustness, considering the effects of possible datasets on training models. Closest to our approach is the work by Meyer et al. [46] for decision trees and Karlas et al. [33] for nearest neighbor classifiers. We use zonotopes to over-approximate prediction ranges for linear regression, generating robustness certificates. While zonotopes have been used for test-time robustness [49, 20, 22, 51], our work is the first to apply zonotopes for training-time robustness for an iterative learning algorithm.
Approaches for uncertainty quantification (UQ) aim to understand the range of outcomes a model may produce using Bayesian methods, ensembling, conformal prediction, and bootstrapping [54, 43, 35, 17, 16]. UQ focuses on epistemic and aleatoric uncertainty, stemming from insufficient data, noisy data, or uncertainty about the model parameters, and does not account for uncertainty due to systematic data quality issues, such as non-random data errors or missing values [25], which induce a multiplicity of possible datasets. In this case, UQ methods might underestimate the uncertainty as they rely on eritical assumptions. Bayesian methods, for instance, require correctly specified priors to accurately model uncertainty, often failing under conditions with unknown erroncous
2

---

priors [73, 67, 72], while conformal prediction (CP) assumes data exchangeability- condition conitoo that breaks down when data errors are systematic [21, 75]. In contrast, our approach addresses this distinct challenge by computing sound over-approximations that guarantee complete coverage of potential predictions across all variations of the dataset. This sound coverage is essential in high-stakes settings such evaluating the robustness of predictive models for medical use.
2 Notation, Problem Formulation and Background
Denote a training dataset D (X, y) with X
x= [
zn," J" Rnxdd as the matrix of features,
Un] E R the corresponding ground truth labels. Let f(;; w) be model parameterized by w E Rº that maps an input data point to a label. learning algorithm A maps training dataset D to the parameters of the trained model, w* A(D). Given test dataset Xeett E Rexd for any test sample from X the function ff computes prediction ÿ = f (z; w*).
and y [y
2.1 Learning Possible Models from Possible Worlds
We possible world semantics to represent the uncertainty in dataset D.
Definition 2.1 (Possible Datasets). Given an uncertain dataset D, the uncertainty in D can be represented by set of possible datasets D D = DD........
Each dataset D; E D is "possible world", i.e., hypothetical variation of the dataset D that could potentially exist in the real world based on our knowledge about the uncertainty in D.
Example 2.2. Consider an e-commerce dataset D where some product price is missing, meaning the exact price unknown. Using possible world semantics, we represent this uncertainty with a set of possible datasets D each containing possible clean price, which could be obtained from prices of the same items on the market.
In App. D, we discuss construction methods for common data quality issues. Our goal is to efficiently construct the of all possible models from uncertain data and understand their behavior in making predictions.
Definition 2.3 (Possible Models and Prediction Range). Given set of possible datasets D associated with an uncertain dataset D, the possible models, denoted fo, are obtained by applying the learning algorithm A to each training dataset D; within D to obtain the set of all possible optimal model parameters w i.e.,
8 {w:
w
D Do}
For test data point the viable prediction range V(x) is defined as the interval between the least upper bound and the greatest lower bound of the outputs produced by all models f i.e.,
V(x)
inf - Ew
[(x. w²),
sup Ewo- fax w²)
This prediction range quantifies the minimum and maximum predictions that can be expected for highlighting the variability in model outputs due to differences the training data. Our framework supports uncertainty in training and test data. We discuss test data uncertainty in App. F.3.
2.2 Sound Approximation of Possible Models with Abstract Interpretation
The set of all possible datasets associated with uncertain data be intractable. We use abstract
interpretation [12] to over-approximate sets of elements of a Alstract concrete domain D (the training data and model weights) with elements from an abstract domain D'd Specifically, we use the abstract domain of zonotopes, a type of convex polytope, to appr-ippro the possible datasets D using a zonotope D that has compact symbolic representation. Instead of applying the learning algorithm A to each possible dataset to compute possible optimal model parameters @h we develop an abstract learning algorithm AF that operates directly on the abstract domain of zonotopes. Given Dd. AF generates a zonotope -- = (FFFF) that over-approximates o= (demonstrated in the graph on the right). Intuitively, this represents the symbolic execution of the learning algorithm across all possible datasets simultancously.
Abstract Learning
Abstract
Dataset
At
Weightt
Dx
W
Possible Dutasets
Possilble Weightss
D
D2
Concrete Learning
3

---

Definition 2.4 (Abstract Domain). Let ID be a concrete domain. An abstract domain for D is a ser ID4 paired with two functions:
Abstraction a P(D) + TDF
Concretization y IDA - P(T))
which satisfy the following condition for any subset S _ D: y (a 2 S. Two abstract elements di and d2 are equivalent, written d1 =; d2, i (d1) = 7 (d2).
Def. 2.4 ensures that the abstract element a{S) associated with set 5 through application of the abstraction function encodes an over-approximation of S. We will use an abstract element DF = a (D)) to over-approximate the possible worlds of an uncertain training dataset D We discuss abstraction functions for specific types of training data uncertainty in App. D.
Definition 2.5 (Abstract Transformer). Consider a function F : D1 + D2 on concrete domains D1 and D. An abstract transformer F ID > ID2 over-approximaats F the abstract domain:
VS E P(D) (L*()()))) 2 F(S)
An abstract transformer is exact (does not loose precision) if E DA :7:77 (F (d))
Importantly, (exact) abstract transformers compose (see App. B.2, Prop. B.1) and, thus, we can construct abstract transformer for complex functions from simpler parts.
To over-approximate the set of possible model parameters w we will develop an abstract trans- former A for the learning algorithm A to get 7 (AF(D))) 2 w
Symbolic Abstract Domains and Zonotopes. We consider a symbolic abstract domain of vectors and matrices (marked with -9 with elements that are polynomials v over variables E { The concretization of polynomial
the result of evaluating to on all assignments E - -,,1. encoded as vectors [-1, fl; 7 (+) {v(e)
[-1, EI.. We lift concretization to vectors and matrices through point-wise application. Such an object 2 is typically referred to as a polynomial zonotope or zonotope if all symbolic expressions are linear (see App. C.2). The concretization of ø² is: y (x+) = (=) E []1,1
Exact Abstract Transformers for Learning Linear Models
Given uncertain training dataset D we aim to over-approximate the set of possible optimal linear models 0* {wi w..... where w E IP represents the optimal parameters of linear model trained on D; E c These optimal parameters are the fixed point of the sequence {w)) ) generated by: w +++1 = (ww where the operator & R² + Hº captures one step of gradient descent, i.e., (bw) = - L(w),, for learning rate 77 and loss function L(w) [57].
In the abstract domain, we use the zonotope representations DE and w to abstract the possible datasets Do and the set of possible model weights @ While is theoretically possible compute symbolic expressions for the standard closed form solution for linear regression, this can result in large expressions that contain fractions with polynomial numerators and denominators and, thus, computing prediction intervals based on such expressions is computationally infeasible (see App. L). Instead, we over-approximate the optimal parameters using an abstract operator & exrect w - w that generates a sequence of abstract elements { [ } Jk=o-
2 k+1 = Φb -biract(uut (wik)
where the abstract operator & Pxzet given by Φt & (wx) w L/V")( captures one gradient descent step in the abstract domain. Specifically, for any loss function L whose gradient V consists of linear polynomial expressions, such as the mean squared error (MSE) loss, & Φkznet is an exact abstract transformer. This follows from the existence of exact abstract transformers for addition and multiplication over polynomial zonotopes [37] and the fact that abstract transformers compose (see App. E, Prop. B.I and Prop. E)
Propusition 3.1. The abstract gradient descent operator & Dexact is an exact abstract transformer for the concrete gradient descent operator Formally, for any abstract w²a
7 { (w))))((0))) $D(7 (wo*))

---

Unfortunately, the sequence {w 1 Th=0 does not have fixed point as the highest-order symbolic terms in w are multiplied with symbolic terms the gradient. Thus, w k++ & Lzanet (w has terms of higher orders than any term in w k To identify when the abstract model weights represent all fixed points the concrete domain, we define abstract fixed points using concretization.
Definition 3.2 (Fixed Point of Abstract Gradient Descent), An abstract model weight w *²* is a fixed point for an abstract transformer Få for the gradient descent operator & if.
7 (w Dt) 7 (1"("=")
Here, Y (:-)) denotes the joint concretization of w² and Då (see App. for more details).
Any abstract fixed point w** over-approximates all possible optimal model weights ox Proposition 3.3. Consider an abstract transformer F² for the gradient descent operator tb Let w be a fixed point according to Def. 3.2. Then holds that y (w x) 2 u
As we demonstrate in App. F.3, w can also be used to compute prediction ranges. In general, an abstract fixed point is not guaranteed to exist for every abstract transformer for 4. While o D as defined above has fixed point (see Prop. F.1 in the appendix), two significant barriers remain. First, the representation size of the polynomial zonotope w for abstract model parameters grows exponentially i, the number of the steps of abstract gradient descent. Additionally, testing for convergence is challenging, checking for containment is already NP-hard for linear zonotopes [39].
4 Efficient Sound Approximation for Learning Linear Models
In this section, we present am efficient abstract gradient descent method for ridge regression that guarantees fixed point and addresses the exponential growth of generated abstract model parameters. The core idea is to use the linearization and order reduction techniques introduced in the following to deal with intractable polynomial zonotopes. We develop specific order reduction technique that admits a fixed point and enables us to find a closed-fone solution for the fixed point.
4.1 Gradient Descent With Order Reduction
A linearization operator L maps polynomial zonotope z² to linear zonotope l that ever-approximates z²: y (L(z+)) 2 (z*). An order reduction operator R takes a linear zonotope é as input and returns another linear zonotope e of reduced order (representation size) such that: > (R(E)) 2 (EF). Specific linearization and order reduction operators are discussed in App. G. The graph on the right shows the progression of the zonotope of parameters toward a fixed point from random FP initialization by applying 200 iterations of the following abstract gradient descent operator Φ² The number inside each zonotope indicates corresponding iteration. After 200 iterations the resulting zonotope is close to fixed point shown in red and computed using techniques we introduce in the following. This operator first linearizes the gradient zonotope using L, subtracts it from the current abstract model parameters w and then reduces the order of the resulting zonotope using R:
30
(w") R. ( -((VV((w"))
This operator is an abstract transformer for the gradient descent operator (Prop. H.1). While each step can computed efficiently by bounding the order of the resulting zono- tope, di may not always converge to a fixed point. The graph on the left illustrates a real example where, despite the existence of a zonotope contain- ing all optimal possible parameters (shown again in red), abstract gradient descent with linearization and order reduction keeps diverging, generating larger and larger zonotopes due to over-approximation error. In the following, we develop an order reduction operator that ensures abstract gradient descent with linearization converges to fixed point for linear regression with & regularization. Additionally, we derive closed-form solution for this fixed point. The red dotted zonotope in the graphs represents the fixed point generated by our method.
189
100
5

---

Decompasition of Gradients. A core idea in our approach is to decompose the abstract dataset into real and symbolic parts and use an order reduction operator that allows for the decomposition of the abstract gradient descent operator into components processing the real and symbolic parts in specific evaluation order. Specifically, we observe that the components of the abstract training dataset Di = (X*, y²) can be decomposed into sum of a real part (without error symbols) and symbolic part (containing only symbolic terms). Furthermore, the abstract model weight w² produced by $4 can be decomposed as shown below. wa is the real part of the zonotope and does not contain any emoor symbols, w/ contains only symbolic terms with error symbols that also occur in DF, and w contains only symbolic terms introduced by linearization and order reduction, hence does not contain any error symbols from D since these methods introduce fresh symbols.
X X R + X
y² + y
w² = wR + w + w
The key observation that enables our approach is that one step of abstract gradient descent can now be decomposed into three components:
(d(w"") 2 (n(wn) 4º D(wn, wh)) + 48 DoNN(u (wwn wh wix)
(2)
where
2 [xx" Xnwn Xu" unl 22wR
(n(wx)= wg 11'
4 D(wgg wx)) = n. (2Nl 2x "xxx" Xx)"x xxxx"X x3 + X3 Xudwa 2 _xtux 2 xxx" vv:)
2 xxx X XR)wn
N (WWRR wh, w w - 7. (2Xl +
X3" xx'xn) w/)) X3 X)(wn wh wxy) x L 2 (Xn X/.
4 (wR)), the real part updater, only relies on and updates wR, the real part of the abstract weight, and coincides with the abstract gradient
wn
(x(wx)
wb (wRR wh wxy
D (wR wfi
operator & on the real part abstract weights. & D (WR, wh), the symbolic data-dependent updater, which is a function of w itself and wn. updates wh, based some linear terms and consists of symbols that come from the abstract dataset DF and does not include any non-linear terms; hence applying this does not generate higher- order terms, so the output is a linear zonotope, and no lincarization is needed. Additionally, no order reduction is needed since the number of eгror symbols remains constant. Φ N (WRR wh wix)
w
the symbolic data-independent updater, a function of w itself, wf,, and w. It consists of the non-linear terms in the gradient, hence performs linearization to over-approximate higher-order terms and performs order reduction reduce the number of error symbols, which would otherwise grow exponentially with the number of iterations. Therefore, generates the updated wh, as a linear zonotope that does not share any symbol with DF, i.e., only consists of fresh symbols introduced by linearization and order reduction. See Appendix I for formal statements and proofs for this section.
Proposition 4.1. Any abstract model weights wi w + w + w is a fixed point for the abstract gradient descent operator 4i if the following conditions are satisfied:
w n = ()(w)),
:h = xb pww wh")
w me &² N (w w w
(3)
We summarize the fixed point construction process in Alg. 1. The algorithm takes as input an abstract dataset D learning rate 17, regularization coefficient X, and a transformation matrix A. The construction or- der is: first w then w and finally w each step building on the previous results. Closed-form solu- tions for wi, and w exist, as detailed in Lem. I.4. The fixed point of R matches the fixed point of gra- dient descent on the real part of the abstract data, in- dependent of other parts. Given wh H- is obtained
Algorithm 1: Abstract Learning
Input: abstract dataset DJ, learning rate 1).
regularization cocff. X, transformation A
w f- closedFormReal (D'
" closedFormSymbl [Dd, wh)
genNonDataaE(DD X, we. wh")
w < solveFixedPointta
wd + -wh + +w w
return u

---

by solving system of linear equations from Eq. (3), ensuring zonotope containment by equalizing the coefficients of the same error symbols, which is sufficient condition for zonotope equivalence. This involves solving system of m equations, where m is the number of error symbols in DF.
To construct " in Eq. (3), which involves zonotope equiv- alence, we encounter challenge because the RHS and LHS have different error symbols, making symbolic enforcement of equivalence inapplicable. We exploit the properties of the order reduction operator to create system of linear equations whose solution yields a fixed point for t N. Re-
ProjectAN
Interiaal
Transffomaation Baseal Houll
Aogmonntal
Innere Rosectiina!
Zhootope
call that this component uses linearization to handle non-linear terms and applies order reduction to efficiently manage the order of the resulting linear zonotope. As shown on the right side, the order reduction operator in Ø n utilizes a transformation matrix A to project the higher-order linear zonotope into new space and then over-approximates it using an interval hull-a zonotope that is a box enclosing the input higher-order zonotope in the projected space and finally projects this box back into the original space using A A-1 The core idea of transformation-basse order reduction [38] is that interval hulls provide better over-approximation in the projected space, where the shape of the input zonotope is closer to a box than in the original space, as shown in the graph on the right.
At a fixed point, the input and output interval hulls of Φ N should be equivalent in the projected space. This equiv- alence can be effectively enforced by checking for equal edge lengths along each dimension, which is sufficient for two interval hulls to be equivalent. This translates system of k equations, where A is the number of model parameters, hence can be done efficiently. We show that this system of equations is guaranteed to have solution. We are now ready to state our main technical result.
Wa
W
(known)
(known)
(unknown)
Theorem 4.2 (Correctness of Alg. 1). Given set of possible training datasets D associated with uncertain data D, and an appropriate abstraction function in the zonotope abstract domain, and given regularization coefficient X, Alg. I, when provided with D = ()) as inpul, computes abstract model parameters w such that:
7 (w-
w 0
where o* denotes the set of all optimal model parameters corresponding to D for linear regression with &2 regularization with regularization coefficient X.
5 Experiments
We implement ZORRO using SymPy [45], a Python library for symbolic computations and evaluate the system two key applications: (1) computing prediction ranges and robustness certification for linear models trained on uncertain data, and (2) robustness of model weights for causal inference using linear models case study. We also measured the performance of ZORRO under varying conditions, including varying the degree of training data uncertainty. All our experiments are performed on a single machine with Apple MI chip, cores, and 16 GB RAM. Experiments are repeated 5 times with different random seeds, and we report the mean (error bars denote 3r). The code is shared at htpps//gituub.com//oonn////
Datasets, Baselines, and Metrics. For robustness verification we use regression tasks: for MPG [58] (392 instances) predict fuel consumption based on car features (cylinders, horse- power, weight); for Insurance [30] (1338 instances) we predict medical insurance charges based on demographies (age, gender, BMJ, habits (smoking), and geographical features. We use 80:20 train-test split and inject random errors to the training data varying (i) the Uncertain Data Percentage, the percentage of instances that have uncertain features labels, and (ii) the Uncertainty Radius, the difference between the minimum and maximum possible value of an uncertain feature expressed as fraction of the feature's domain. There no direct baseline capturing prediction ranges for learning linear regression models from uncertain data. Most existing works in robustness certification focuses on test-time robustness (cf. Sec. 1). The exception is [47], which only supports uncertainty in labels using interval arithmetic. We compare ZORRO against this approach for robustness certification, re- ferred to as MEYER in the following, only for training label uncertainty. For robustness verification, prediction is robust if the size of the prediction interval is smaller than given threshold. We use
7

---

Uncertainly Radius: 5 %
Uncertainty Radius: 10
Uncertainly Radius:
Uncertainty Radius: 20 %
1.0
10
Robustness
0.5
0.5
0.5
0.5
0.0
0.0
0.0
0.0
0.0
5.0
75
5.0
2.5
75
10.0
0.0
7.5
Uncertain Data Percentage (%)
20RRO
-8- Meyer et
(a) Robustness verification with uncertain labels (MPG data).
Robustness Ratio
Uncertainty Radius: 4 %
Uncertainty Radius: 6 %
Uncertainty Radius: 8
Uncertainty Radius: 10 %
1.0
1.0
0.5
0.5
0.5
0.5
0.0
0.0
0.0
0.0
7.5
5.0
2.5
10.0
7.5
Uncertain Data Percentage
ZORRO
-0
Meyer
(b) Robustness verification with uncertain labels (Insurance data).
Uncertainty Radius: 6 %
Uncertainty Radius: 8 %
Uncertainty Radius: 10 %
Uncertainty Radius: 12 %
1.0
Robustness
0.5
0.0
0.0
2.5
5.0
7.5
10.0
0.0
5.0
7.5
10.0
0.0
2.5
5.0
10.0
2.5
5.0
10.0
Uncertain Percentage (%)
(c) Robustness verification with uncertain features on MPG data.
Figure 1: Robustness verification on using intervals (MEYER [47]) and zonotopes (ZORRO).
the robustness ratio which is the fraction of the test data receiving robust predictions as metric in all robustness verification experiments. The robustness threshold is set 5% of the label range for the MPG data, and 0.8% of the label range for the Insurance data. Additionally, we assess the worst-case test loss using certain test data and uncertain model weights trained from uncertain training data.
5.1 Robustness Verification
Prediction Robustness (Uncertain Labels). Fig. (a)((b compares ZORRO with the baseline MEYER [47] on setting where only training labels suffer from uncertainty. We vary the uncertainty radius and uncertain data percentage. As both systems provide sound over-approximation of prediction ranges, they may underestimate model's robustness. As shown, ZORRO consistently certifies significantly higher robustness ratios than MEYER. This is due to the fact that MEYER uses -arthmei-iiitmmeticc which ignores the correlation between model weights in different dimensions and between the training labels and the weights, leading to overly conservative prediction ranges. Specifically, for higher uncertainty radius values, MEYER fails certify robustness for most of the data while ZORRO still can certify robustness for 100% of the instances.
Prediction Robustness (Uncertain Features). We also evaluate the impact of training feature uncertainty (not supported by MEYER). Specifically, introduce uncertainty into the vehicle weight column for the MPG dataset. As shown in Fig. I(c), uncertainty in the features results in relatively less robust predictions compared to uncertain labels for similar uncertainty radius (Fig. l(a)). This primarily because uncertain features result in more high-order terms in the closed form solution than uncertain labels which in turn leads larger over-approximation emrors during linearization.
In all experiments the standard deviation of the robustness ratio, calculated by repeating experiments with different random seeds, is large when the average robustness ratio is close to 0.5. This is because when the uncertain data percentage low, the model will be robust no matter which training instances are selected to be uncertain. Likewise, when the uncertain data percentage is high, then most predictions will be uncertain no matter which training data points are uncertain.
Parameter Robustness. Next, we apply ZORRO for robustness certification of parameters in linear regression models, crucial for statistical estimation and causal analysis. We compare the ground truth coefficients for a treatment variable with the results obtained by ZORRO and through KNN imputation dataset with injected missing data. Fig. 2 shows the treatment variable coefficient and

---

Intercept
0.03
0.03
ZORRO
0.05
0.00
GT
0.00
-0.03
KNN
0.50
0.55
0.055
0.00
0.05
0.50
0.45
Treatment Effect
Treatment Effect
Treatment Effect
(a) Consistent Positive Effects
(b) Contradictory Effects
(c) Consistent Negative Effects
Figure 2: Applying ZORRO to causal inference. The intercept (y-axis) is the model's bias term, the
treatment effect (X-axis) is the coefficient for the treatment variable.
0.20
Zorro
0.12
Zorro
GT
0.12
0.10
0.09
0.04
0.08
0.06
0,12 0.19 0.25 0.31
0.38
0.44
10
12
13
uncertainty radius
dimension
(a) Varying uncertainty radius
(h) Varying data dimension
Figure 3: Range of the loss, through enumeration of all possible worlds (GT) and ZORRO.
the regression model's intercept. The intercept captures the baseline level of the outcome variable when all predictors are zero, highlighting how baseline values can shift under uncertainty. While the model trained after KNN imputation sometimes correctly identifies the directionality of the treatment effect, this is not always the case shown in Fig. 2(b). This highlights the needs for techniques like ZORRO which guarantee that the true treatment effect within certain bounds.
5.2 Solution Quality with Varying Uncertainty and Hyperparameters
We evaluate ZORRO'S effectiveness by testing the tightness of the over-approximation and the accuracy of possible models. Specifically, we examine how the over-upproximation quality and worst-case loss are influenced by the level of data uncertainty and the regularization coefficient.
Varying Data Uncertainty. To evaluate how specific characteristics of the data affect the effec- tiveness of ZORRO, we injected errors into real datasets, varying uncertain data percentage and uncertainty range. We compare ZORRO with the ground truth range of the loss computed by enu- merating all possible worlds (GT). The results shown in Fig. 3(a) demonstrate that ZORRO tightly over-approximates the ground truth loss range, especially for smaller uncertainty radius values. As uncertainty increases, the aver-approximation gap widens due to the increased coefficient of higher- order terms in the gradient, which are linearized, leading to higher linearization errors. As shown in Fig. 3(b), the tightness of ZORRO's over-approximation is not affected by the dimension of the data.
investigate the impact of the regularization coefficient on the ro- bustness of predictions and the worst-case loss of possible models. Following similar approach to Sec. 5.1, we introduce uncertainty in both features and labels in the MPG dataset and use zonotopes with varying levels of uncertainty to over-approximate the training data uncertainty. The results, shown in Fig. 4, indicate that higher regularization coefficient leads to more robust predictions, as regularization tends to "compress" all possible model weights towards the origin. Interestingly, the worst-case loss shows that À = 0 is not optimal across all scenarios, especially when the fraction uncertain instances is high. Instead, small, positive (c.g.. 0.02 or 0.025) generally yields the best worst-case losses. Combining these results, the optimal regularization coefficient should enhance robustness (i.e., higher robustness ratio) while maintaining an acceptable worst-case loss. Therefore, the regularization coefficient should be tuned based on validation dataset to achieve a small range of accurate possible models.
Effect of Regularization.
Conclusions, limitations and broader impacts
We introduce an approach for propagating uncertainty through model training and inference for linear models. Given an abstract uncertain training dataset that over-approximates the possible worlds of training dataset, we develop abstract interpretation techniques to over-approximate the set of possible models and inference results for this set of models using zonotopes. This challenging, we need to compute fixed points of gradient descent in the abstract domain. Our main technical

---

Worst-cass Loss
2.60
Uncertain Data: 8%
1.00 60
Uncertain Data: 10%
1.0
3.00
Uncertain Data: 12%
1.0
Uncertain Data: 14%
1.0
Robustnes Ratio
2.40
0.98-40
0.98
40
0.5
2.20
0.8
2.50
3.00
2.20
0.0
0.00
0.04
0.00
0.04
0.00
0.02
0.06
Figure 4: Varying regularization coefficient A: Robustness ratio (green) and worst-case test loss (red).
contribution is the development of closed-form solutions for such fixed points that can be solved efficiently. Our techniques efficiently over-approximate models and inference for several use cases, including robustness verification, uncertainty management in causal reasoning, and improving the interpretability and reliability of predictions and inferences. This framework can be particularly valuable in critical applications where data quality and robustmess are paramount. While we propose an effective method for abstract learning linear models, extending our approach to more complex models is challenging. Non-linear models, such as neural networks, would require more advanced linearization and order reduction techniques, as well as parallelization, to manage the increased complexity of the involved symbolic operations. In addition, adapting our method to broader classes of models through efficiently approximating the fixed points remains a challenging and promising future direction.
7
Acknowledgment
This research supported by NSF awards I[S-2340124, IIS-2420691, and IIS-2420577, well as by NIH grant U54HG012510. The views, opinions, and findings presented are those of the authors and do not necessarily represent those of the NSF NIH.
References
[1] Amr Alanwar, Anne Koch, Frank Allgöwer, and Karl Henrik Johansson. Data-driven reachability analysis using matrix zonotopes. In Learning for Dynamics and Control, pages 163-175. PMIR, 2021.
Matthias Althoff. Reachability analysis and application the safety assessment of autonomous cars. PhD thesis, Technische Universität München, 2010.
[3] Matthias Althoff. Reachability analysis of nonlinear systems using conservative polynomialization and non-convex sets. In Proceedings of the 16th international conference Hybrid systems: computation and control, pages 173 182, 2013.
[4] Matthias Althoff and Bruce H Krogh. Zonotope bundles for the efficient computation of reachable sets. In 2011 50th IEEE conference on decision and control and European control conference, pages 6814-6821. TEFF, 2011.
[5] Stanley Bak, Sergiy Bogomolov, Brandon Hencey, Niklas Kochdumper. Ethan Lew, and Kostiantyn Potomkin. Reachability of koopman lincarized systems using random fourier feature observables and polynomial zonotope refinement. In International Conference on Computer Aided Verification, pages 490-510. Springer, 2022.
[6] Aharon Ben- Tal, Dick Hertog, Anja De Waegenaere, Bertrand Melenberg, and Gijs Rennen. bust solutions of optimization problems affected by uncertain probabilities. Advanced Risk & Porifolio Management Research Paper Series, 2011.
[7] Emily Black, Manish Raghavan, and Solon Barocas. Model multiplicity: Opportunities, concerns, and solutions. Proceedings of the 2022 ACM Conference Fairness, Accountability, and Transparency, pages 850-863, 2022.
[8] Xavier Bouthillier, Pierre Delaunay, Mirko Bronzi, Assya Trofimov, Brennan Nichyporuk, Justin Szeto, Nazanin Mohammadi Sepahvand, Edward Raff, Kanika Madan, Vikram Voleti, et al. Accounting for variance machine learning benchmarks. Proceedings of Machine Learning and Systems, 3:747-766 2021.
[9] Leo Breiman. Heuristics of instability and stabilization in model selection. The annals of statistics. 24(6) -383, 1996.
10

---

[10] Christophe Combastel. Functional sets with typed symbols: Mixed zonotopes and polynotopes for hybrid nonlinear reachability and filtering. Automatica, 143:110457, 2022.
[11] A Feder Cooper, Yucheng Lu, Jessica Forde, and Christopher M De Sa. Hyperparameter optimization is deceiving us, and how to stop it. Advances in Neural Information Processing Systems, 34:3081 2021.
[12] P. Cousot. Abstract interpretation. CSUR, 28(2):324-33 328, 1996.
[13] Patrick Cousot and Radhia Cousot. Abstract interpretation: a unified lattice model for static analysis of programs by construction or approximation of fixpoints. In Proceedings of ACM SIGACT-SIGPLAN symposium on Principles of progrumming languages, pages 238-252, 1977.
[14] Patrick Cousot and Nicolas Halbwachs. Automatic discovery of linear restraints among variables of a program. In Proceedings of the 5th ACM SIGACT-SIGPLAN symposium on Principles of programming languages, pages 84-96, 1978.
[15] Alexander D' Amour, Katherine Heller, Dan Moldovan, Ben Adlam, Babak Alipanahi, Alex Beutel, Christina Chen, Jonathan Deaton, Jacob Eisenstein. Matthew D Hoffman, et al. Underspecification presents challenges for credibility in modern machine learning. Journal of Machine Learning Research. 23(226): 2022.
[16] Fabricio Olivetti de Franca and Gabriel Kronberger. Prediction intervals and confidence regions for symbolic regression models based likelihood profiles, 2022.
[17] Nicolas Dewolf, Bernard De Bacts, and Willem Waegeman. Valid prediction intervals for regression problems. Artificial Intelligence Review, 56(1):577-61) 613,, April 2022.
[18] Ilias Diakonikolas, Gautam Kamath, Daniel M. Kane, Jerry Li, Ankur Moitra, and Alistair Stewart. Robustness meets algorithms. Commun. ACM, 6465):107-115, apr 2021.
[19] Tlas Diakonikolas and Daniel M. Kane. Recent advances in algorithmic high-dimensional robust statistics. CoRR, abs/191 1. 1.091 1, 2019.
[20] Tianyu Du, Shouling Lujia Shen, Zhang, Jinfeng Li, Shi, Chengfang Fang, Jianwei Yin, Raheem Beyah, and Ting Wang. Cert-rnn: Towards certifying the robustness of recurrent neural networks. CCS, 21(2021):15 19, 2021.
[21] Rina Foygel Barber, Emmanuel Candes, Aaditya Ramdas, and Ryan Tibshirani. The limits of distribution-free conditional predictive inference. Information and Inference: A Journal of the IMA, 10(2):455-82, 2021.
[22] Timon Gehr, Matthew Mirman, Dana Drachsler-Cohen Petar Tsankov, Swarat Chaudhuri, and Martin Vechev. Ai2: Safety and robustness certification of neural networks with abstract interprctation. In 2018 IEEE symposium on security and privacy (SP), pages 3-18. IEEE, 2018.
[23] Matthew L Ginsberg. Counterfactuals. Artificial intelligence. 30(1):35 79, 1986.
[24] Eric Goubault, Tristan Gall, and Sylvie Putot. An accurate join for zonotopes, preserving affine input/output relations. In Jan Midtgaard and Matthew Might, editors, Proceedings of International Workshop on Numericca and Symbolic Abstract Domains, volume 287, pages 65-76, 2012.
[25] Coreela Gruber, Patrick Oliver Schenk, Malte Schierholz, Frauke Kreuter, and Göran Kauermann. Sources of Uncertainty Machine Learning A Statisticians' View, arXiv, 2023.
[26] Arthur Hoerl and Robert W. Kennard, Ridge regression: Biased estimation for nonorthogonal problems. Technometrics, 42(1):80-86 2000.
[27] Chengqiang Huang, Zheng Hu, Xiaowei Huang, and Ke Pei. Statistical certification of acceptable robustness for neural networks. In Artificial Neural Networks and Machine Learning, pages 79-90. Springer, 2021.
[28] Tomasz Imieliński and Witold Lipski Jr. Incomplete information in relational databases. Journal of the ACM (JACM), 31(4)7761-79, 1984.
[29] Matthew Jagielski, Alina Oprea, Battista Biggio, Chang Liu, Cristina Nita-Rotaru, and Bo Li. Manipulating machine lcarning: Poisoning attacks and countermcasures for regression learning. In 2018 IEEE symposium security and privacy (SP), pages 19-35. IEEE, 2018.
[30] Arun Jangir and Willian Oliveira. Healthcare insurance, 2023. This dataset is licensed under CCO: Public Domain.

---

[31] Jinyuan Jia, Xiaoyu Cao, and Neil Zhenqiang Gong, Intrinsic certified robustness of bagging against data poisoning attacks. In Proceedings of the AAAI conference artificial intelligence, volume 35, pages 7961-7969, 2021.
[32] Matt Jordan, Jonathan Hayase, Alex Dimakis, and Sewoong Oh. zonotope domains for lagrangian neural network verification. In NIPS, 2022.
[33] Bojan Karlas, Peng Li, Renzhi Wu, Nezihe Merve Gürel, Xu Chu, Wentao Wu, and Zhang. Nearest neighbor classifiers over incomplete information: From certain answers to certain predictions. Proc. VLDB Endow., 14(3):255-22) 267, 2020.
[34] Pasha Khosravi, Yitao Liang, YooJung Choi, and Guy Van den Braeck. What to Expect of Classifiers? Reasoning about Logistic Regression with Missing Features. arXiv, 2019.
[35] Byol Kim, Chen Xu, and Rina Foygel Barber. Predictive inference is free with the jackknifet-after bootstrap. In International Conference Neural Information Processing Systems. Curran Associates Inc., 2020.
[36] Niklas Kochdumper and Matthias Althoff. Sparse polynomial zonotopes: A novel representation for reachability analysis. IEEE Transactions on Automatic Control, 66(9):4043-4058, 2020.
[37] Niklas Kochdumper and Matthias Althoff. Constrained polynomial zonotopes, Acta Informatica, 60(3):279- 316, 2023.
[38] Anna- Kath-ia Kopetzki, Bastian Schürmann, and Matthias Althoff. Methods for order reduction of zonotopes. In 2017 JEEE 56th Annual Conference on Decision and Control (CDC), pages 5626-5633. IEEE, 2017.
[39] Adrian Kulmburg and Matthias Althoff. On the co-np-completeness of the zonotope containment problem. Eur. J. Control, 62:84 91, 2021.
[40] Sungyoon Lee, Hoki Kim, and Jacwook Lec. Graddiv: Adversarial robustness of randomized neural networks via gradient diversity regularization. IEEE Transactions on Pattern Analysis and Machine Intelligence, 2022.
[41] David Lewis. Counterfactuals and comparative possibility. In IFS: Conditionals, Belief, Decision, Chance and Time, pages 57-85. Springer, 1973.
[42] Charles Marx, Flavio Calmon, and Berk Ustun. Predictive multiplicity in classification. In International Conference on Machine Learning, pages 6765-6774. PMLR, 2020.
[43] Peter McCullagh, Vladimir Vovk, Ilia Nouretdinov, Dmitry Devetyarov, and Alex Gammerman. Condi- tional prediction intervals for linear regression. In International Conference on Machine Learning and Applicutions, pages 1-1-138, 2009.
[44] Johannes Mehrer, Courtney J Spoerer, Nikolaus Kriegeskorte, and Tim C Kietzmann. Individual differences among deep neural network models. Nature communications, 11(1):5725, 2020.
[45] Aaron Meurer, Christopher Smith, Mateusz Paprocki, Ondřej Čertík, Sergey B Kirpichev, Matthew Racklin, AMIT Kumar, Sergiu Ivanov, Jason K Moore, Sartaj Singh, et al. Sympy: symbolic computing in python. PeerJ Computer Science. 3:e 103, 2017.
[46] Anna P. Meyer, Aws Albarghouthi, and Loris D' A'Atoni. Certifying data-bias robustness in linear regression. CoRR, abs/2206.03575. 2022.
[47] Anna P Meyer, Aws Albarghouthi, and Loris D' Antoni The dataset multiplicity problem: How unreliable data impacts predictions. In Proceedings of the 2023 ACM Conference on Fairness, Accountability. and Transparency, pages 193-204, 2023.
[48] Antoine Miné. The octagon abstract domain. Higher-order and symbolic computation, 19:31-100 100, 2006.
[49] Matthew Mirman, Timon Gehr, and Martin Vechev. Differentiable abstract interpretation for provably robust neural networks. International Conference on Machine Learning, pages 3578-3586. PMLR, 2018.
[50] Matthew Mirman, Alexander Hägele, Pavol Bielik, Timon Gehr, and Martin Vechev. Robustness certifi- cation with generative models. In ACM SIGPLAN International Conference on Programming Language Design and Implementation. pages 1141-1154, 2021.
12

---

[51] Mark Niklas Müller, Mare Fischer, Robin Staab, and Martin Vechev. Abstract interpretation of fixpoint iterators with applications to neural networks. Proceedings of the ACM Programming Languages. 7(PLDI)778-881) 2023.
[52] Nicolas Müller, Daniel Kowatsch, and Konstantin Böttinger. Data poisoning attacks regression learning and corresponding defenses. Pucific Rim International Symposium on Dependable Computing (PRDC), pages 80-89. IEEE, 2020.
[53] Hongseok Namkoong and John Duchi. Stochastic gradient methods for distributionally robust optimiza- tion with f-divergences. In D. Lee, M. Sugiyama, U. Luxburg. I. Guyon, and R. Garnett, editors, Advances in Neural Information Processing Systems, volume 29, Curran Associates, Inc., 2016.
[54] David J. Olive. Prediction intervals for regression models. Computational Statistics & Data Analysis, 3166:3115-3122 2007.
[55] Andrea Paudice. Luis Muñoz-Gonzáled and Emil C Lupu. Label sanitization against label flipping poisoning attacks. In ECML PKDD Workshops, pages 15, Springer. 2019.
[56] Robert Plemmons. M-matrix characterizations. i-nonsingular m-matrices. Linear Algebra and its applications, 18(2):175-188, 1977.
[57] Boris T Polyak. Some methods of speeding up the convergence of iteration methods. Ussr computational mathematics and mathematical physicx, 1(5) 1964.
[58] R. Quinlan. Auto MPG. UCI Machine Learning Repository, 1993. DOL: https://bi.org// wttg//00.2432//588500 This dataset is licensed under Creative Commons Attribution 4.0 International (CC BY 4.0) license.
[59] Raymond Reiter. A theory of diagnosis from first principles. Artificial intelligence, 32(1 5557955,,, 1987.
[60] Flan Rosenfeld, Ezra Winston, Pradeep Ravikumar, and Zico Kolter. Certified robustness to label- ahhel-flipping attacks via randomized smoothing. In International Conference on Machine Learning, pages 8230-8241. PMIR, 2020.
[61] Donald B Rubin. Multiple imputations in sample surveys-a phenomenological bayesian approach to nonresponse. In Proceedings of the survey research methods section of the American Statistical Association, volume 1, pages 20-34. American Statistical Association Alexandria, VA, USA, 1978.
[62] Donald B Rubin. Multiple imputation after 18+ years. Journal of the American statistical Association, 974444)473-489 1996.
[63] Anian Ruoss, Mislav Balunovic, Mare Fischer, and Martin Vechev. Learning certified individually fair representations. In Advances Neural Information Processing Systems, pages 7584-7596, 2020.
[64] Christian Schilling, Marcelo Forcts, and Scbastián Guadalupe. Verification of neural-network control systems by integrating taylor models and zonotopes. Proceedings of the AAAI Conference on Artificial Intelligence, volume 36, pages 8169-8177, 2022.
[65] Soroosh Ahaaee, Peyman Mohajerin Esfahani, and Daniel Kuhn. Distributionally robust logistic regression. In Neural Information Processing Systems, 2015.
[66] Ili Shumailov, Zakhar Shumaylov, Dmitry Kazhdan, Yiren Zhao, Nicolas Papernot, Murat A Erdogdu, and Ross J Anderson. Manipulating sgd with data ordering attacks. Advances in Neural Information Processing Systems, 34: 1 11108 2021.
[67] Danicle Silvestro and Tobias Andermann. Prior choice affects ability of bayesian neural networks to identify unknowns. arXiv preprint arXiv:2005. 04987, 2020.
[68] Gagandeep Singh, Timon Gehr, Matthew Mirman, Markus Püschel, and Martin Vechev, Fast and effective robustness certification. Advances in newral information processing systems, 31, 2018.
[69] Jasper Snock, Hugo Larochelle, and Ryan P Adams. Practical bayesian optimization of machine learning algorithms. Advances in neural information processing systems, 25, 2012.
[70] Jacob Steinhardt, Pang Wei Koh, and Percy Liang. Certified defenses for data poisoning attacks. In NewrIPS, 2017.
[71] Damien Teney, Maxime Peyrard, and Ehsan Abbasnejad. Predicting is not understanding: Recognizing and addressing underspecification in machine learning. In European Conference on Computer Vision, pages 458 476. Springer, 2022.
13

---

[72] Ba-Hien Tran. Simone Rossi, Dimitrios Milios, and Maurizio Filippone. All you need is a good functional prior for bayesian deep learning. Journal Machine Learning Research, 23(74):1-55, 2022.
[73] Rens van de Schoot, Sarah Depaoli, Ruth King. Bianca Kramer, Kaspar Märtens, Mahlet G Tadesse, Marina Vannucci, Andrew Gelman, Duco Veen, Joukje Willemsen, al. Bayesian statistics and modelling. Nature Reviews Methods Primers, 1(1):1, 2021.
[74] Marianne Southall Winslett. Reasoning about action using a possible models approach. Number 1427-1432. Department Computer Science, University Illinois at Urbana Chamaaagn,,,,,n 1988.
[75] Margaux Zaffran, Aymeric Dieuleveut, Julie Josse, and Yaniv Romano. Conformal prediction with missing values. In International Conference on Machine Learning, pages 40578 40604, PMLR, 2023.
[76] Xuezhou Zhang, Xiaojin Zhu, and Stephen Wright. Training set dehugging using trusted items. In AAAI, 2018.
14

---

A Linear Regression and Ridge Regression
In this section, we review standard loss functions for linear regression, specifically mean squared error (MSE) and the loss function used in ridge regression. can be trained via either the closed-fom solution or gradient descent. Specifically, the closed-form solution requires computing the inversion of the covariance matrix, while the gradient descent only involves matrix addition and multiplication. Since the linear regression model is convex, gradient descent is guaranteed to converge to the global optimum.
Suppose we have training data D (X. y) with n i.id. samples where the feature matrix X [2) In. E Rnxd, and the labels [,,,.... [M]] E RR. Given an input and the model weight w. the prediction of the linear regression model ÿ w x (the bias term can be integrated into w, corresponding to an added column in with constant 1's).
Loss Functions
The mean squared error (MSE) loss on D is defined as shown below.
1 (w" v1)² -X (X - y).
L(X, y. w)
(1/1 y:)²
(4)
In practice, regularization terms, e.g., based the 4--norm of the model parameters, often added to the original MSE loss to prevent overfitting by penalizing large weights. Using Iz-reullaizzation with a regularization coefficient À which determines the strength of regularization is often called ridge regression. The loss function for ridge regression is:
L(X, y. w) = (Xw (X y) + A w w.
(5)
to the convexity of linear models, the locally optimal point, which can be obtained by gradient descent with an appropriate learning rate n, is globally optimal. In gradient descent, the model weights w are iteratively updated (with some learning rate 77)
Gradient Descent for Linear Regression
towards the reverse direction of the gradient aL(X.y.0) Om
Thus, for ridge regression:
O w)
2 2
Xw " Xx-XXy) 2Xw.
(6)
Ow
Thus, one step of gradient descent + is:
2 2 (x " XT 2Xw
w w²
Closed-form Solution for Linear Regression The convex nature of linear models ensures that the optimal weight w = w", which minimizes the loss L(X, y, w), can be computed by establishing BL, X m.) 0:
OL(X, y, w²)
2 =²(XT Xm" -x + 2\ =
Ow
11
X (XT T X + =
w" (X x + Anl) T
B Background on Abstract Interpretation
Abstract interpretation [13] is a technique aver-appro the results computations over a set of inputs. This is achieved by associating sets of elements of concrete domain TB with elements from an abstract domain DF. this context, various abstract domains are employed: interval domains represent variables ranges of possible values [13], octagon domains allow for constraints between
15

---

pairs of variables within specific bound [48]. and polynomial zonotopes which represent convex polytopes using polynomial constraints [14]. Abstract interpretation, while originally designed for static program analysis such as strictness analysis, has also found applications in wide range of other domains including reachability analysis [1, 5, 2, 10, 36, 4], robustness verification for neural networks [32, 64], learning robust models by providing bounds on the loss for set of inputs [51, 63, 49, 22], and many others.
We now state some important, but well-known, facts about abstract transformers that we utilize in our derivations.
B.1 Abstract Transformers
In Def. 2.5, we presented the standard definition of abstract transformers as functions over abstract domains that over-approximate the application of functions in the concrete domains to sets of elements. To clarify the connection between abstract transformers and possible world semantics observe that both take a concrete function F D + D' and lift it to sets of inputs S P(B)) through point-wise applications:
F(S) = (F(c) e
Thus, abstract interpretation is natural fit for over-approximating PWS by over-approximating sets of possible worlds using abstract elements and then over-approximates computations with PWS for a function f using abstract transformers for such a function.
B.2 Abstract Transformers Compose
Importantly, abstract transformers compose. This enables us decompose a complex computation into simpler operations and build an abstract transformer for the computation by composing abstract transformations for these simpler operations.
Proposition B.1 (Abstract Transformers Compose). Consider (exact) abstract transformers f and ge for functions f and g. then g² ft is an (exact) abstract transformer for o f.
C Symbolic and Standard Representation of Zonotopes
We now provide a more detailed account of the correspondence between the symbolic and standard geometric representation of zonotopes and polynomial zonotopes and discuss why matrices and sets of symbolic matrices as used in our abstract domain can equivalently be thought of as (polynomial) zonotopes.
C.1 Symbolic vs. Geometric Representation
We defined zonotopes z vectors ddr where Vi denotes polynomials over variables (called error symbols) E. The concretization of such symbolic representation of zonotope is the set of vectors in gd that can be derived from z² by assigning values from [-1, 1] each variable E E. We encode such variable assignments as vectors [-1, 1,
Definition C.1 (Polynomial zonotopes Symbolic Representation), A d-dimensional polynomial zonotope is vector z E T [10]. Let E). The concretization of z is defined as:
(z (e) ee 1-1,1]181]
common measure of the representation size of a zonotope is its order. The order of a d-dimensional (polynomial) z² in symbolic representation is the total number of distinct monomials in z divided by
ORD(2*) #M(=*)
where #MLb) denotes the cardinality of the set of distinct monomials in polynomial vb
more common way to represent d-dimensional zonotopes is by fixing a set S of monomials over and writing the zonotope as central point E Rd and sum of generator vectors 9i E Bdd multiplied with the monomials from S. The generator vector gi assigns coefficient to monomial S(i] for each of the dimensions.
16

---

Definition C.2 (Polynomial zonotopes Geometric Representation). The geometric representation of d-dimensional polynomial zonotope sum of center point с R and the monomials over error symbols in & from set S multiplied by coefficients encoded in sei of generator vectors E IR:
1S]
2F = c+
9
For the geometric representation, the order of a d-dimensional (polynomial) z² over monomials S is defined as:
ORD(=*) S
Note that the order of the symbolic and geometric representation of zonotope are the same. As an example consider the 4-dimensional polynomial zonotope z² shown in symbolic representation (left) and standard representation (right).
I 10p + lc1 C263
4+ 262 + 2c1 E2263
62c
0 262 3 3e) F223
0 2
2
2 + 162 + 461 2263
1
S[I]
S(2]
C.2 Symbolic Matrices and Sets Symbolic Matrices
We matrices over symbolic expressions and sets of heterogeneous abstract matrices represent the state of computation and allow matrices from such set to share variables to encode relationships between the elements of such matrices. The semanties we associate with such set is that of a joint concretization. For 8 fM,', we define:
(st)= = (M) M S AcE
zonotopes encode convex sets of points in Rd. We can think of symbolic matrix M va as m-dimmenionall polynomial zonotope. Similarly, set of symbolic matrices S {M," Ra me4xm;; can be thought of as I-dimensional polynomial zonotope R' where E, n4 m,, If every symbolic expression in a matrix or set of matrices linear combination of error symbols (an affine form), then such objects can equivalently be represented as zonotopes.
For instance, below the left show zonotope matrix W (all expressions are linear) and two possible worlds in its concretization (for assignments [-1, -1] and [o, 0.5)).
c1 + 62
(
2 -1
0 0.5
[3
0.5
W =
W
WF
15
2261
15
22
15
11
C.3 Abstract Training Data and Model Weights
Using possible world semantics to compute all possible model weights @ given a set of possible training datasets D there will be natural correspondence between model weight w E w and the dataset D e from which is was derived. When training in the abstract domain such correlations can be preserved by sharing error symbols between the abstract training data DF and corresponding model weights w When such sharing occurs, then it is critical to reason about the joint concretization when determining whether an abstract fixed point as been achieved. Testing equivalence (equal concretization) of w alone can lead to false positives, as w"1 2t w does not in general imply (D², wi >; (D², w²2). This is due to the fact that a particular concrete model weight w E (wii) = (wiz) may be associated with different datasets in (D², w²1) and D w²r) because of shared error symbols between DF and the abstract model weights. Thus,
17

---

concretization equivalence between model weights only does not imply equivalent results after application of gradient decent step.
To further illustrate this consider, two pairs of abstract model weights w²1 and wF2.
w
2 26
D
61
3+
5
w²2 3 62
2
D
E1 5
The only difference between _ and w²z is that w² shares a variable (ei) with DF while w² does not. Observe that 7 (wf, y (wf2). However, applying step of gradient decent w² and w²2 may lead to abstract models that are not equivalent. For sake of this example, consider a hypothetical gradient operator & dummy that subtracts the data from the current model weight.
T dammy (uo*1)
2 61
DE
61
-2 62
sb duummyy (w² 2)
2 253 61 -2 -2
D
c1
5
Note that &f duummy (w-)) 2 L2 dummy concretization this does not guarantee that fixed point has been reached when also taking the data into account. This is important as otherwise it is not possible that the abstract model weights contain all possible optimal model weights w Instead, need to consider the joint concretization of model weight and data, s.t. if 1St = (5t.2) then concretization equivalence guaranteed to hold for all subsequent iterations.
(w=²). Thus, even if two abstract model weights have equal
D Examples Of Training Data Uncertainty and Abstraction
In this section, discuss several causes of training (and test) data uncertainty, how to encode them as possible worlds, and abstraction functions for approximating such uncertainty in our symbolic model.
D.1 Measurement Uncertainty
Sensors typically have some measurement uncertainty. If the measurement error is know, e.g., provided by the instrument manufacturer or estimated through repeated measurements and calibration, then for a set of sensor readings used in training D (( y.). then each possible world in D is derived from by replacing values az[]] with values xz +
D.2 Missing Values and Imputation
Consider a training dataset D where some of the features missing for some of the datapoints in D. If each feature X;'s domain is an interval [4 ug] and assuming that missing values can be represented as independent random variables, then the set of possible worlds of the training data e are all training datasets that can be derived from D by replacing each missing value in feature X; with value from [4, u:).
We can use an abstraction function Ommsssing that represents each missing value as an interval [4: u.]. In the symbolic representation, this is encoded as the central point of the interval of an error symbol with coefficient half of the interval's length. We associate a separate error symbol with each missing value. For D pu
we define Amisssing for i [1, n]] and E
[1, m] as:
unthy + u un as 2
if D =1 otherwise
Amissing (Di,)
Dis
18

---

For instance, consider the training dataset D E R shown below where - denotes a missing value.
3
D =
1
6
9
Assuming that both features have domain [0, 10], we get:
5 5
5
Omissing (D)
6
+ 563,1
9
If the independence assumption on missing values holds, then y (Amisssng (D) = D If the independence assumption does not hold, then D would be subset of the worlds described above and Amiesinnn (D) is still valid abstraction function, albeit an vveaapproximating one:
? (Omissiig (D)) D
Imputation If we make the stronger assumption the unknown ground truth value corresponding to missing value is in set of estimations {ai,.... am | returned by a set of imputation methods, then we can use (min((a;))) max(fa;)) instead of [4, u:]-
Training Label Uncertainty The abstract transformer Amissing can also be used training data labels are missing.
E Abstract Transformers for Zonotopes
In this section we introduce (exact) abstract transformers for (polynomial) zonotopes that have been introduced related work and discuss their computation complexity and the space requirements for the output zonotope zf in terms of its order ORD(z (2) (see Sec. 2.2). Kochdumper et al. [37, Table 1] shows an overview of which operations are exact for linear and polynomial zonotopes (and other set representations). Relevant for our purpose is that exact transformers exist for polynomial zonotopes for all operations used in the learning algorithms for linear models we consider.
Proposition E.1 (Exact Transformers for Polynomial Zonotopes). There exist exuct abstract trans- formers for scalur addition and multiplication well for matrix addition and multiplication for polynomial zonotopes [3]. There exist exact transformers for scalar addition and matrix addition for zonotopes. Abstract transformers for multiplication and matrix multiplication for zonotopes exist, but are not exact.
We present the details of these operations in the following.
E.1 Arithmetic Operations
Addition. Scalar and matrix addition are exact in both zonotopes and polynomial zonotopes. Given two matrix (polynomial) zonotopes V and W in maxm their addition Zt V W is defined by adding entries. For each i [1, n] and j E m.::
Z デ
V ig W
The order of Z is the sum of the orders of V and Wt:
(RD(Z*) = ORD(V=) + ORD(W"))
Scalar multiplication.
Multiplying a (polynomial) zonotope matrix Wil with a scalar c exact using the abstract transformer "poty defined below. We simply multiply each entry in the matrix by c.For each i [I,m]] and j E [1,m]:
(c pply wt) W
19

---

Multiplication with scalar abstract value d E Li is exact for polynomial zonotope matrices, but increases the order of the input zonotope:
ORD(d- W) = (RD(d) ORD(W)
Here define the order of scalar abstract value d v to be the number of monomials in the representation of d. For a linear zonotope W, d i is in general not linear as contains higher- order terms. Thus, matrix multiplication for linear zonotopes requires application of linearization (see App. G):
d -ii W L(d poly W)
Matrix multiplication. Matrix multiplication is defined using scalar multiplication and addition. As discussed above, addition is exact for both linear and polynomial zonotopes. However, scalar multiplication of symbolic expressions is only exact for polynomial zonotopes, but requires lineariza- tion and, thus, over-approximation for linear zonotopes. That is, matrix multiplication is exact for polynomial zonotopes only. Consider two matrices V E vmnxm and V E frmxk then matrix mul- tiplication is defined as usual, but using abstract transformers for scalar addition and multiplication. Each symbolic entry in the matrix V W is sum of m elements, each the multiplication of one entry of V with one entry of W Thus,
ORD(V* W m ORD(V*) ORD(WW) For linear zonotopes, we have to again apply linearization to make sure that the output a linear zonotope.
F Abstract Transformers for Gradient Descent
F1 Abstract Fixed Points Over-Approximate Possible Fixed Points
Proof of Prop. 3.3. Initially, will assume that wis is computed through repeated application of F#. Let n be the smallest number such that w = w = (w DF) for iteration with abstract transformer F². As F is an abstract transformer for & and abstract transformers compose (Prop. I), we know that for every D; D and E N, we have (w D;) E 7 (wi .Di). To prove the claim it is sufficient to show that for every such D; we have (w, D.) E 7 (w D*). WLOG consider some D; D and (w." D.) E 7 (w D²). Let m E N be the smallest number such that : = wT m S n, then based on the fact that F is an abstract transformer (over-approxinates 4), the result holds as (w", D2) E y (WF. Dt). Now consider the case where m > n. We will show through induction that for all j [n + 1, 1, (w,, D,) E (wi. Df) and, thus (w) D:) E (wi" D
n-I,
Induction start. For j = n the result trivially holds based on the definition of abstract fixed points.
Induction step. Assume that (w? D.) E (wf" DF) for j E [n, m - 1]. we have to show that this implies that (w? Di) E y (whn, Di). By definition, have
w +1 = ((ww))
As Ft is an abstract transformer for Φ. we have, (w!+1.D. D.) E (wi | .D). Now based on the fact that win is an abstract fixed point according to Def. 3.2 and, thus, (win, D¹) 2 7 (ww++t) D²). it follows that (wi+ D.) E (D.)).
So far we have demonstrated that fixed point wi that appears in the iteration sequence {w { contains all optimal model weights w @+ We now prove the stronger result that as long as wi fulfills the condition of Def. 3.2, no matter not, its concretization encloses w Consider one D, Y (D#) and as above let w denote its optimal model weight. We will demonstrate that (W,, D:) E 7 (wi x D²). First note that given that gradient descent for linear models is convex, for any initial model weight wo. the sequence {w/ h converges to w. Specifically, let denote the model weight such that (w, D.) 7 (wi Dt)
the result of an iteration sequence using F#t or
20

---

and let denote the smallest integer such that w." = " for the sequence generated starting from
w. w. However, now we can apply the same proof by induction shown above demonstrate that (w D.) E (w² D*). This concludes the proof.
F.2 Exact Abstract Transformer for Gradient Descent
In this section we show that the abstract transformer for gradient descent introduced sec. 3 has fixed point.
Proposition F.1. The abstract gradient descent operator & has a fixed point w
Proof of Prop. F.I. The existence of a fixed point is implied by the fact that & is an exact abstract transformer of the concrete gradient descent operator 4. Consider y (w DD) {ww,. D.)). Let n; be the smallest integer such that w." w and let = max; n4 i.e., at iteration the concrete model weights have converged in every possible world in the concretization of (w .DF). ). As o is an exact abstract transformer, can show by induction that y (w* , {(w D.)} } for any As all computations in the concretization have converged at n. we know that w/ ((w()) w
and, thus, [w., D.)] {(w," +1 D.)). As Φ exact, we get the desired result: 7 (wf DF) ++ D.)] = y (w n+1 ID+).
{(wF. D.)] = {[w,
F,3 Prediction with Abstract Model Weight Fixed Points
We now discuss how to use the abstract model weights wi* returned by our abstract transformer for learning linear models during inference to over-approximate the possible set of predictions for test data point. We start by discussing test data that is not uncertain and then extend the discussion to the case where the data also uncertain.
Deterministic Test Data For now let assume that the test data is not uncertain. The following corollary then enables to use abstract gradient descent for inference and for over-approximating the prediction ranges.
Corollary F.2. Let (x) denote the linear model for parameters w. Given an incomplete training dataset D associated with set of possible worlds D the prediction range V V(x) for test data point x can be over-approximated by:
V²(x)
fw( (wEn(w>))
min
max
fw)
"EET(
where w x is the fixed point of the abstract gradient descent operator tiF applied to DI = a(D9),
the abstract representation of D the zonotope domain.
The prediction returned by linear model with parameters w for data point w x. As does not contain any symbolic term this is the sum of linear terms multiplied by constants, i.c. the result is a 1-dimensional linear zonotope which is a linear expression of the form:
m
Co
i=1
where E R and E E. The minimum and maximum value of 1-d linear zonotope can be determined efficiently as shown below:
Σ 1=1
c...c
+
c
Uncertain Test Data In Sec. 2 we modeled uncertainty in the training data sets of possible worlds D As mentioned in that section, our techniques also supports uncertain test data, i.e., both the training and test data may be uncertain.
21

---

In the most general case, the training and test data may be correlated.' ! We model this as set of possible worlds (D の X (est )) where each world is pair of a training dataset D, and test dataset
Xiests
(D, (iet)) = ( Xest)...... (Dms X
If training and test data are independent, then can specify their worlds separately ( n and X ) and assume the worlds of (est) to be their cross product. Uncertainty propagation for inference then requires us to compute the set of possible predictions:
= (e(())))))) (; E (D,X = A(Di)
For inference in the abstract domain first have select appropriate abstraction function for the test data, e.g., using the same abstraction function a we use for the training data. The only difference the case discussed above is that test data point is now also abstract element of and the prediction wFx* is a polynomial zonotope as it contains higher order terms that are the result of multiplying linear terms. To efficiently determine the minimum and maximum we can, e.g., employ lincarization to map the polynomial zonotope into a linear one and apply the solution described above for finding the minimum and maximum of linear zonotope.
G Linearization and Order Reduction Techniques
G.1 Linearization
The purpose of linearization is to over-approximate an input polynomial zonotope with a linear zonotope.
Definition G.1 (Lincarization Operator L). A linearization operator L maps a polynomial zonotope z² linear zonotope L It replaces high-order polynomial terms with new ervor symbols, ensuring that all expressions are linear while maintaining an aprrxpppro
Y (L(z+)) 2 (2f).
During inference, computing the prediction intervals using the linear zonotope representation can be done efficiently with linear programming. However, more importantly, in construction of abstract fixed points we use a specific order reduction technique that requires prior lincarization in each gradient descent step to enforce the existence of a fixed point. For a d-dimensional polynomial zonotope z with a set of monomials S:
S
#
we have its linearization:
L(2*) c+
Σ S
Σ g&c 100
Linear monomials
Replaced with lincar monomials
where denotes the set of indices of all linear terms z Here, the third part over-approximmtes each high-order term in z by replacing it with a new error symbol.
G.2 Order Reduction Operators
Order reduction operators are used to reduce the representation size of zonotope. That is, an order reduction operator takes as input linear zonotope é and return linear zonotope of smaller order
NNote that this does necessarily imply violation of the i.i.d. assumption. For instance, consider dataset with a textual feature race that first translated into categorical feature that then one-hot encoded into multiple binary attributes. If we are uncertain about the meaning of a particular value of the original attribute, then this leads to correlation between the uncertainty of both datasets as the interpretation of this value affects all data points with that particular value in the race feature before preprocessing.
22

---

(smaller representation size) that -vr-approximates lt, For linear zonotopes this means that order reduction operators reduce the number of distinct eπror symbols that occur a zonotope by merging emrrr symbols.
Definition G.2 (Order Reduction Operator R). An order reduction operator R takes linear zonotope e as input and returns another linear zonotope of reduced order:
Y (R((F) 2 2 (et)
ORD(R((1) < ORD(E)
We now present details about two commonly adopted order reduction techniques: Interval Hull (IH) and transformation-bee baasdd Interval Hull (TIH) in App. G. In this section, we use the geometric representation of zonotopes (see Def. C.2). Note that for a linear zonotope, each Sli] consists of single error term e;- Thus, we write such zonotope & as:
S
e c
g:
IH merges selected subset 8g of the error symbols of zonotope and their corresponding generator vectors. For d-dimensional zonotope f::
(Eies, [9.(1]]) é
[8|
Rin( (F) = RIH 8,, c + Σ gic;
=c Σ
g:c.
+
iml
185,
(Eies, [d]] c
Retained error symbols
Overaapproximated with d-dimensional box
The selected terms 8, are merged into a d-dimensional box described by new error symbols {ti,, ..g). We will drop 6, if all error symbols of the input zonotope are selected.
The error symbols getting merged (6.) are often determined based some heuristic, e.g., symbols with lowest coefficients. TIH RITH first projects the zonotope to another space using an invertible linear transformation matrix A Rdd then applies IH, and finally projects the resulting zonotope back with A- IH is special case of TIH where A =
RTIH (6,, A, xf) A Ri (A=f).
Intuitively, the purpose of the linear transformation A is project the zonotope space where its shape is closer to box, to reduce the loss of precision brought by order reduction [2, 38). One of the most widely used TIH is PCA-based order reduction, whose transformation is obtained from the PCA of the set of all generator vectors of the input zonotope [38, 51].
Order Reduction for PTIME Zonotope Training With linearization, the number of error symbols in the model weights zonotope w still grows exponentially with rate (pp)), leading to exponential time complexity for gradient descent. We can overcome this challenge through order reduction, which
enforces the maximum number of terms in the symbolic expressions of model weights zonotopes [38]. Note that for a linear zonotope, each S[i] consists of a single error term G, Thus, we can write such zonotope l as:
[ Σ g&
é c
Order reduction R. reduces the order (representation size) of a linear zonotope é., This is achieved by merging error symbols in S, while ensuring that the result over-approximation the input zonotope. i.e.,
(R(E)) 27 (et)
23

---

Two commonly adopted order reduction approaches [38] are Interval Hull (IH), denoted as RIH- and Transformation- -aasa Interval Hull (TIH), denoted R 12, 38, 64). Specifically, IH merges a set of error symbols 8 S {Gi and their corresponding generator vectors (here d%.
{() d.):
[S] c +
RIH(2*) RIH
=c
Eg.5 diagl xEAx
Σ :E4,
[1
(d) E
The selected terms 8g are merged into a d-dimensional box described by d new error symbols {ei, [''). The error symbols getting merged (6.) are often determined based on some heuris- tic [38], e.g., the symbols with lowest coefficients. Similar to IH, TIH also merges enror terms using IH. but some projected space. TIH first projects the zonotope to another space using an invertible linear transformation matrix A Rdxd then conducts IH, and finally projects the resulting zonotope back into the original space with 1-1:
R11= A- RR(n
Intuitively, the linear transformation aims to project the zonotope to a space where its shape is closer to a box, to reduce the loss of precision brought by order reduction [2, 38). One of the most widely used TIH is PCA-based order reduction, whose transformation is obtained from the PCA of all generator vectors [38, 51].
H Efficient Abstract Gradient Descent with Order Reduction
As mentioned in sec. 4.1, address the tractability issues with abstract gradient descent, employ two key techniques: linearization and order reduction. We employ linearization at each step of gradient descent to ensure that the resulting abstract representation of model parameters remains linear zonotope.
Given linearization operator L and order reduction operator R, construct abstract gradient descent operator, $2:
'' R ( L(VV((w²))
which ensures that the abstract representation size remains bounded while providing efficient over- approximation. This operator is abstract transformer for the concrete gradient descent operator. Proposition H.1. For any linearization and order reduction R. the abstract gradient descent &F an abstract transformer for the concrete gradient descent operator t. Formally, for any abstract wi,
7 (fF()) Φ(> (w=)),
Proof. Given that abstract transformers compose (Prop. 1. we can decompose into separate steps and construct abstract transformer for & by composing abstract transformers for the individual steps. We can inject identity functions anywhere into the computation without changing its result. Let ident represent identify function of an appropriate type and 8 be the abstract operator resulting from injecting identity functions as shown below, then:
Fd ("w") ident w - ident (1 ((w))( = =w (nVL(w"))
4b birrect (wo)) (w²)
Since $² exract is an abstract transformer for Φ,, so FF, Now observe that both and are abstract transformers for ident as
(L(w")) 2 (w- ident(7 (wf))
Thus, T is abstract transformer for t.
Note that Prop. 3.3 then implies that a fixed point wi for oFis an overaapproximation of all possible model weights "
(wi I w o+

---

H.1 Abstract Gradient Descent With Order Reduction And Fixed Points
While $ ensures that every step of gradient descent can be computed efficiently as we bound the order of the resulting zonotope in each step, this operator typically does not have fixed point and even if it does, we still would have to solve an NP-hard problem [39] to detected that we have converged. The reason for the lack of fixed point is that both linearization and order reduction results in over-approximation and the over-approximation error may grow each iteration. Recall that we did show real example of where this operator diverges in Sec. 4.1.
I An Efficient Approximate Abstract Transformer for Ridge Regression
In this section and the next section we present the proof of our main technical result: Thm. 4.2. Recall that Thm. 4.2 states that Alg. computes abstract model parameters w that are a fixed point for the abstract transformer @Fusing a closed form solution. Because oF shown to be an abstract transformer for gradient decent this then implies that w² over-approximates all possible model weights w We start by presenting additional details of the decomposition we employ to force fixed points, formally prove that the condition on w * from Prop. 4. is a sufficient for w being fixed point for ΦF, and then develop closed form solution that requires solving system of linear equations. For y the parts of the abstract model weights that exclusively contains symbols that do not appear in the abstract training dataset D the equation system only has a solution if the regularization coefficient À is larger then or equal to a constant 3 that depends on D Based on our experience and extensive experimental evaluation, we typically have / 0, i.e., the closed form solution exists for any regularization coefficient X. Nonetheless, we present a technique for achieving any desired B 2 0 by splitting zonotopes into smaller parts with lower B3 computing fixed points for each split individually, and merging the final result. These techniques will be presented in app. J.
L.1 Decomposing Fixed Point Equations For Abstract Gradient Descent
We now present additional details about our decomposition of the gradient from Sec. 4.1 for linear regression with & regularization (ridge regression). The loss function I. for ridge regression an & penalty, is given by:
LXX, w) - y) (Xw 9) + x w0" w,
where À is the regularization coefficient (see App. A for details).
Recall that observed that an abstract model weight zonotope w" can be decomposed into parts that can be dealt with separately in the sense that we will show that a sufficient condition for achieving a fixed point w" is that each component has a fixed point. Given an abstract training dataset D use X and yu to denote its feature matrix and labels. We decompose them into real (concrete) and symbolic components:
x = X + X y yR + V
where XR Raxd and yR E 1- represent the real centers of zonotopes X and y while X Andd and y A contain the symbolic terms.
Similarly, we decompose the abstract model weights at iteration into real and symbolic components:
w i w/ + w
where w/ E df represents the real center, and w E Ad contains the symbolic terms. The symbolic terms are further decomposed into those containing data symbols " i.e., symbols that are shared with D;, and those introduced via linearization and order reduction in linear abstract gradient descent (w ) i.e., that do not appear in DE
#i w
w + WN-
25

---

Accordingly, we decompose the abstract gradient presented in Eq. (9) into several distinct components: real numbers nt),. linear symbolic expressions that share symbols with D (GP). linear symbolic expressions that do not share symbols with Di (G ). and high-order symbolic expressions (9#).
Using the loss function for ridge regression (see Eq. (5) in App. A):
T5 (w= R. uw ((]VI(u'))) R ww nI ( (x* X'w x²² y²) 22xur)
R wn w w 7 2 (Xa + x3) (Xn X3) (wn + wi wix)
(10)
(XR + x3) (yx v + 2X(wx + w - w ) R wR - 1 GR + w/ - 17 0o + w - n. 2. - -- L(9H))
(11)
where
2 (XR2 RRWR XR² yn) + 2A
GR
X
(12)
n
9. (2X1 +
2 R X x + [
2
R
T
x + X T
XR)WR 2 n
2
X ys
(13)
n
n
n
2
G (2X +
n
T X )
(14)
n
Gn 2 (x R
x
( + wh + wxx) x T x (u) + wxx + X3 T
yss
(15)
In Eq. (10) we substitute definitions which are further expanded in Eq. (11) using new notation for rearranged parts of the gradient (Eq. (12) to (15)). Furthermore, in this step we also use the fact that linearization only affects higher order terms and Gи is the only component of the gradient with non-linear terms. Now consider the relationship of Eq. (12) to (15)
Enforcing Fixed Points With Parameterized Order Reduction. We now introduce an order reduction operator A that is specific type of TIH (see App. G) which only merges non-data symbols and is parameterized by its transformation matrix A, identify a sufficient condition for achieving an abstract fixed point for gradient descent with this order reduction operator, and demonstrate that the fixed point exists for any coefficient À of lg-regularization that is larger than a data-dependent constant 3. Furthermore, we show how such fixed point wn can be computed using closed form solutions for w R. w and wxx
Recall from App. G that the TIH order reduction operator RITH (ō., A) is parameterized by d, (the set of error symbols to merge) and a transformation matrix A. Consider the input to order reduction as shown in Eq. (10), We use Ep to denote the error symbols that this input shares with the data (that occur D and En to denote those that only occur the input to order reduction. Then we define RA- = RTTH (EN. A). i.e., we only merge symbols that are not shared with the data to ensure that the correspondence between model weights and possible worlds is preserved by letting w share symbols in linear expressions with D The linear part that does not share symbols with DF has (ppag) monomials (where p- = EEE and is the number of error symbols in w²). Now expanding the definition of abstract gradient descent with order reduction using RA the order reduction operator (for some given transformation matrix A), we get:
T (L. ( = RR - n9 + w/ P.P.. + wiy - n + L(GH)))
(16)
= un 7 + w ng + R (w 19 L(Gn)))
(17)
WR - 1/9 + wb noP. A (RIH (A(w)) - nO + L(Gn))))
(18)
(n(wn)
n (w
Eq. (18) shows the components of & AL, ( (wf and related them to the notation ΦR(WR).
('D(wR, w and (" wh wx) ) we introduced in Eq. (2). Specifically, ΦR(WR) = (wn - n)n) does not contain any error symbols and, thus, corresponds to the real number part of the updated model weight. Similarly, T' (WWR, wh) = (w) - ngP)) only contains linear symbolic expression with symbols that are from Dd. The remaining part, 4h x (WRR wh wxx
26

---

A (A (wfy 7 - L(GH)) M)) does not share symbols with DE, as it consists of new symbols generated from order reduction.
I.2 Constructing Fixed Points For Abstract Gradient Descent
Having defined of (A, which was denoted by & (L.R) in app. H, we are now ready to state the first part of our main technical result: under some mild assumptions on À being larger than some data- dependent constant B, (L. has fixed point and, importantly this fixed point can be computed efficiently. In practice, we often have 3 O. However, if this is not the case, we demonstrate later that at the cost of reduced computational efficiency we can reduce / by splitting DF into several zonotopes and solving the problem independently for each zonotope we will details in App.
Theorem I.1 (Existence of Abstract Fixed Points). Consider an abstract training dalaset D There exists a constant specific to D such that 4 (L, has a fixed point for any 2 B.
In the remainder this section, we prove Thm. I.1 by demonstrating how to construct such a fixed point efficiently. We start by identifying sufficient condition for w to be a fixed point, then demonstrate that the problem of finding w that fulfills this sufficient condition can be decomposed based on our decomposition w w "h + w presented earlier. Specifically, certain components are independent of other components suggesting an evaluation order where we find a fixed point for one component treating the previously computed fixed points for other components as constants. Then we proceed to prove closed form solutions for wi and for wh given w.. Finally, given w and wh" we construct system of equations whose solution gives fixed point for w and demonstrate that this system of equations has closed form solution for any À 2
Sufficient Condition for Abstract Fixed Points. Recall from Def. 3.2 that the fixed point w must satisfy 7 (w+ DD) 2 2 (Ø) (A) (wd We now present sufficient condition for this to hold. Intuitively this condition requires the invariance of different components of Φ1 (I. (wF) when doing joint concretization with D
Proposition 1.2. If an abstract model weight w = = + wh + w satisfies the following three conditions, then an abstract fixed point of the abstract gradient descent operator T
(A.)
wi w 19R
(19)
DD = w -ng?.
(20)
i+ =; RA (w - n6 L(Gn))) W
(21)
Proof. Consider the definition of abstract fixed point w (Def. 3.2 applied to of (L.RA) ):
(wF+,D*) ,D)) 2
2 (wF Dt)
This can certainly be achieved if
w Dd 24 T (L. (I..Ra)(="))) (w DF
This requires equal joint concretization of the fixed point the pairs (w* (w*) and (L,R,) (( DF). To understand the need for joint concretization here note that certain model weights upaairdd with certain datasets D in the joint concretization. If we would only consider equivalence of the abstract model weight, two zonotopes wi and wt' may be equivalent w =; w"" but pair particular model weight E (w-) = (wt") with different datasets D2 E (Dt) and D2 E y (Dt). Thus, it is possible that w" me 43 (L. (WF) holds but y (w=) ( (wot) does not contain all model weights from St
27

---

First observe that requiring equality of abstract elements is stricter condition than equivalence wrt. 2 as
wi = w²2 + w" = wi2
Also note that w is the only component of w-* that contains symbols. Thus, requiring equality in Eq. (20) is sufficient for ensuring equivalent join concretization with DF as long as Eq. (19) and Eq. (21) also hold. Furthermore, Eq. (19) and Eq. (20) imply:
GR = 0
GP =
Expanding w w + w/ + " and substituting into Eq. (18) we get
ΦF (L. (wi+)) wi ngn + why n RA (w n9( L(Gn)))
= w R (w "99 L(Gn)))
Thus, using Eq. (21) and the fact that w does not share any symbols with Dř, we get the desired result:
(19) A (20) A w = RA (w - n99 + L(GH))) (w"r,D*) = (# (wf"),
Independence of Fixed Point Components. Recall that we observed that some components of the gradient depend only on some components of wn Specifically, GR only depends on w and G only depends on wi and wh This implies that we can find a fixed point by first computing wi.- then finding wh," treating wR constant, and finally determine w treating both w and wh)" as a constant.
Lemma 1.3 (Independence of Fixed Point Components). The following independence relationships hold:
(i) GR does not depend on wh nor w
(ii) gP, does not depend on w
Proof. The definition of GR is shown in Eq. (12). The only component of w that appears in the definition of Gn is wT. Thus, the equations defining wi only depend on w and are independent of wh and w Using a similar argument we can show that GL does not depend on w
Closed Form Solutions for Real And Linear Data Components of Abstract Model Weights. Before explaining how determine w utilizing the lemma above we derive closed form solutions for w and wh.
Lemma 1.4 (Closed Form Solutions for w and w Given Dł, wie as defined below fulfills Eg. (19)
w/ (X T Xn + AnI) T yn
(22)
Given Dd and wi wh as defined below fulfills Eq. (20).
wh = (XR T Xx R AnI) U Xx R T
(X T X X XR)WĀ)
(23)
Proof. Fixed Point for w Substituting the definition of GR from Eq. (12) into Eq. (21):
Wn w - (2 (2 [xR" Xnw Xx yn) Tnn
2xw
21
(nn' T X Xnw - Xnn yn + \nank) 0
R T XRwp XxR yR Anw =0
(n>
* (XR,² T Xx + >n = Xnn yn
S (X R
X R + XnI) -1 X yR
(24)
28

---

The last step relies on the fact that XR XR"XR + Anl is invertible which is guaranteed for À > 0 [26]. Note that, as expected since the fixed point equation for wi does not contain any symbolic expressions, the closed form for wxR is equal to the well-known closed form for ridge regression. Fixed Point for wh Note that given the independence of wi from wh we can assume that win has been determined using the closed form solution shown above. Thus, we can treat wij constant in the following derivation. We start by substituting the definition or Eq. (13) into Eq. (20).
wh = =w n (2N
2 XxR" R T X KR)W +
2 (X R X n
X XR)W3
n
(x SUR + X y'
n
<0 2n (X T X + Aml " + (Xri T X + X T
Xn)wf - X - Xn' v's)
+0 (AnI+ + Xx X x + (X f1 X + x T
T
n| xn)wh XLyn Xn T y
(n) 0)
(Anl + XR X R w
(Xx' X + x X R ) - x X - XR T n/s
+ = R T + Xnl)
X
X syRR + xxxx y (X x - x
X R xn)wn
(25)
Note that matrix inversion is applied to (X T x + XnI) which does not contain any symbolic terms and as discussed above this matrix guaranteed to be invertible.
Analysis of the Fixed Point Equations for w Using the closed form solutions established in Lem. 1.4, we can compute wa and w fulfilling Eq. (19) and Eq. (20). In the remainder of this subsection we show how to compute a solution w to Eq. (21) for a given À 2 B (regular-
ization coefficient), i..., prove Thm. I.I that postulated the existence of such fixed points wi constructively.
To find a fixed point w we are searching for an abstract model weight that fulfills Eq. (21):
w = Ra w n + L(Gn)))
(26)
This implies that w needs to be equivalent to the result of RA on some input zonotope w "9 + L(G Recall the definition of R for linear input zonotope t!:
( A RTHH)
(27)
Note that R applies RIH to all error symbols in its input. Thus, (-H(-) is d'dimensional box in the projected space into which the input l is mapped into by A. Furthermore, as w +L(GH) does only contain symbolic terms, the interval hull (the box computed by RTH) is centered at O.
Lemma I.5. Any solution w'y of Eq. (26) is equivalent to the image of a box ba centered at 0 produced by A i.e.,
w = A b
Proof. IH does not change the center of a zonotope and A is a linear map. Now observe that all terms in Eq. (14) and (15) are symbolic (contain error symbols). Thus, the RHS of Eq. (26) that is the input to order reduction is a zonotope with center 0.
A linear zonotope that is box in d-dimensional space with center O is. up to equivalence, uniquely
determined by the diameter of the box in each dimension. That is, if we use K; 2 0 to denote the diameter of the box in the tth dimension and & e Ie to denote the vector of these elements, then we can write any such zonotope b² as:
k161
h1
G1
ba =
:
= diag
kaed
0
kdd
Ed
Ed
29

---

Combining Eq. (26) with Eq. (27), we know that w is a zonotope that is equivalent to the image of A-1 on box bi parameterized by k. We will only consider solutions of this form:
61
w A-'diag kt)
kg
(28)
where ",,".., add are the column vectors of A-1 As the naming of error symbols is irrelevant, finding a solution to Eq. (26) now amounts to finding vector ls as the error symbols are now treated fixed.
Now we can substitute the formula for w from Eq. (28) into the RHS Eq. (26) as shown below. Recall from app. G.2 that IH merges different error symbols by adding up the absolute values of their coefficients. For the 7th dimension, each kÚ will appear with a certain coefficient. For some terms, the symbolic expression for the 7th dimension may contain error symbols that do not have any ky in their coefficients. Thus, without considering for now the precise values for each coefficient we know that the box produced by IH the RHS is of the form shown below where C E RRdd is a matrix with non-negative entries Ci.j 0 (at position (i,j)), and Co E R is vector with non-negative elements C,, (at dimension i). Intuitively, Ct,j is the coefficient of kj in the 7th dimension and C,o0 is the sum of the coefficients of error symbols that do not have any kj in their coefficients (e.g., the fresh symbols introduced by linearization). Note since IH does calculate coefficients of the error symbols in the output as a sum of absolute values, we know that all entries of C' and Co are positive. Furthermore, recall that assumed that kj 2 0.
Ra(wN 799 + L(GH))) A RI (A(w)y n(9 L(Gn))))
(C1.0 + Ej= ki)k
(29)
A-
A diagteo + Ck}
(ca.0 + Ej= kz)
Substituting the LHS and RHS of Eq. (21) with Eq. (28) and (29), we get:
E1
A A }
2g A diagteu + Ck}
cd
[e]]
- (kag
24 diag(co + Ck)
(30)
6d
* Co + Ck
(31)
(I - C² = cu
s.t. k >> 0
(32)
+11 - (111,,)kk
(Grkk)) C.U s.t. k. 2 1 [1,4]
(33)
sf1
Intuitively, the equivalence of boxes from Eq. (30) is equivalent to the equality of their diameters which yields Eq. (31). Therefore, any & E Rd satisfying Eq. (31) gives us a fixed point w In Eq. (32) make explicit our assumption that k; 2 0 and finally in Eq. (33) write Eq. (32) as system of linear equations.
"Thi does not restrict possible concretization any such zonotope we can achieve. As invertible, A has full rank and thus, its column vectors ...... as form a basis of a vector space. Given A 1 the concretization of a zonotope of this form is uniquely determined by & (recall that we assume that all k. are positive).
30

---

Existence of a Closed Form Solution for the Fixed Point Equations for w What remains to be shown to prove Thm. 1.1 that Eq. (33) is guaranteed to have solution. For that we will investigate the structure of C' and co- Recall that the matrix C and co are obtained from IH order reduction on input (see Eq. (29)
A(w - n99 + L(GH))).
w and 3. share the same set of symbols, while they do not share any error symbols with L(G ) as linearization always generates new enror symbols. This implies that the terms of A 7 ) and -r)AL((H)) will not cancel out. In other words, they contribute separately to the diameter & of the merged box produced by IH. Therefore, we will look into their contributions separately.
For A(w ), we substitute the definitions of G., and w to get:
Alwix - 1922 - A((1 222AA1
27 X
X R
w
(34)
KIEi
= A (1 (1-21)2) 27 XR XRR
A
(expand wxx)
kaEa
k1c
27
A(1 - 111-22AJA-1-
A(1
A-1 AXR" "XRA-1
=
n
kiged
(I= AA -1 and X²
Xn = AAn
XAAA
kier
= -(AN - 2MA)A -1 2n
AXR XRA-1)
(35)
Kaed
This contributes (11- 2m/ 2n
% 1kg
29 Σ5 =1
16.5 lkt)
to the diameter of the merged box along
dimension i, where qb is the element in AX XRA- at position (f,j). By choosing a small
learning rate 7) to let
2n/2- 21 qr 2 0 the contribution is equal to ((1 2nd - 2n 4.,);
27 E (9.ss11) Ik;) From Eq. (34) is obvious that this component has a shared term w containing jfi
k and as w/ is centered at O, we know that this component has only symbolic terms. This in turn implies that this component only contributes to the coefficient matrix C' and does not contribute to the constant part eo of the system of linear equations.
Now consider the other component
- A9((H)
(36)
- nAL
7 (X
2
X +
X
(R) + w/x + x
X/(WRR + w/ + w//)
x3² n/s
AL
n
x
R
x x )
T
contributes to C
(37)
+ x T x (S + why) x
T
+ X
x
+ X
X XR)WS,
s ys
contributes
Note that the first part of the equation is the only part that contains w which in turn contains the variables & for which we want to solve. As every element of this part of the equation contains some & it only contributes to C. In contrast the second part of the equation does not contain any variable k and, thus, is the only part of both components that contribute to co. Let 2n C'k to denote the matrix of coefficients that is the sum of absolute values of coefficients of
31

---

(Xx + x Xx + X w that will be projected by A stemming from the first part n of the equation in the expanded version of L( We use di 2 to denote the element of C" at (i.j). As mentioned before this is the only part of this component that contains terms containing unknowns k;- We use do to denote the contribution from the second part of Eq. (37) towards the output of linearization L. Then based on the fact that the second part of this component is the only part of both components that contributes to Co, we know that the whole contribution of the component from Eq. (36) towards C + Co is:
nd + 20 C'k:
(38)
In summary, the system of linear equations from Eq. (33) can be written as:
2n
20
20
(2// + 11,s
Ki
(1/6ol + dggDisg = 114,0 s.l. ko 20
Vi Ei [1.d]
f=1 jデt
(Am + 9o.e Ds
(1as) + c2.)):g
n s.t. kg 20 Vi [1.]]
2 c2.
(39)
In the expanded system of linear equations, the diagonal elements of the coefficient matrix are (An + qú,: c,,)), while the afi-diagonal elements are -1g.,s d..)- Lemma 1.6. for all i E [1,d]
(An + q8 - d,)) 2
(gul + d,y))
then the system of linear equations Eq. (39) has a solution.
Proof. When this condition holds, the coefficient matrix of this system of linear equations is diag- onally dominant. Also considering all diagonal elements non-negative ((n + qú, <,,)) 2
EE= (195,) 4.4) 2 0) and all off-diagonal elements are non-positive, the coefficient matrix is 331 an M-matrix For any M-matrix M, the inverse M exists and is has all non-negative entries. Furthermore, for any system of linear equations represented as matrix equation Mæ = b where æ the variables of the equation system, M b is solution. As our matrix I - C is an M-matrix, we know that its inverse exists and has only positive entries and we can compute solution & as:
= (1 - C co
Given that all elements Cro0 are positive and that (7 - C) positive, the solution k for Eq. (39) also positive.
The correctness of Thm. I.1 follows immediately based on the sufficient condition for the existence of solution stated in Lem. I.6..
Corollary 1.7 (Proof of Thm. I.1). Given input abstract training dataset D let B R max, Eh (1q.aa) + c2 +c,, -(.)) For any regularization coefficient A 2 B, we can compute
141
an abstruct fixed point u²* for the model weights according to Def. 3.2.
Proof. To find an abstract fixed point of the model weights
w wi + w x + w
'There are over 40 sufficient conditions for a matrix to an matrix [56]. The one we use here requires that the matrix fultillss two conditions: (1) all off ofi-diagonal elements are non- non-positive, and (2) each diagonal element is no less than the sum of absolute values of the off -ii-iiaa elements at the same row (also known as diagonally dominant).
32

---

that fulfills the sufficient condition for an abstract fixed point from Prop. 1.2, we first use the closed form solutions from Lem. 1.4 to compute wi and then using wi compute wh. These are guaranteed to fulfill the first two conditions of Prop. I.2. As the regularization coefficient 2 3 = # max; E== (19s5) + (42)) + C, qt.i), the condition from Lem. 1.6 holds. Thus, iti
the linear system of equations for w that fulfillss the last condition of Prop. I.2 has a solution that we can compute using any techniques for solving a system of linear equation. It follows that w wh + wh + w fulfills all conditions of Prop. I.2 and is an abstract fixed point according to Def. 3.2.
Weakening the Requirements on Regularization
The existence of a fixed point wi as postulated in Thm. I.1 only holds for sufficiently large regularization À 2 B where 3 depends on D.. The majority of real world examples without multicolinearity in the real part of the training data that we have investigated have 3 = 0. However, we still provide a technique to reduce 3 if the need should arise at the cost of increasing the runtime of our approach. Specifically, we will make use of splitring [2] which divides an input zonotope into multiple zonotopes of smaller expand. We apply splitting to Dd, then compute the fixed point for each zonotope in the result of splitting, and finally combine these fixed point using a join operation for zonotopes that wper-approximates the union of the concretizations of two zonotopes (allows us to merge multiple zonotope such that the merged zonotope has a concretization that over-approximmtes the concretizations of all of the inputs). We will demonstrate that this approach enables us to reduce to any value 2 O. Intuitively, this works because we can reduce the extend of the input zonotope for the fixed point computation to the extend necessary to ensure 3 = 0. Albeit in worst-case this can require a splitting operator that generates a number of zonotopes that exponential in the dimension.
J.1 Splitting and Join for zonotopes
Intuitively, a splitting operator divides an input zonotope into two or more zonotopes such that the union of their concretization is the same as the coneretization of the input zonotope. Definition J.1 (Split Operator). An operator S that maps a zonotope et to set of zonotopes [et,.. et Y is for any input zonotope E we have:
S(E=) te,,
e 3
s.t. UT 7
As an example of a split operator consider the binary split (2-split) dimension denoted by S2, that splits the input zonotope into two parts by scaling the jth generator gi by 1 and shifting the center of the zonotope by 29: (- 2-)::
(2) {ei { fot Σgse .c Σgse }
1.0
10.1
ather E|([-111]
other cyh [-1]
= =fc- fc 2
+ 2 g9iG +
c -+7g: + g +
Σ [o J
(40)
As 7 (<4) U 2-split to an (m, i)-split Sm,, that divides the zonotope evenly into m parts along the 7th dimension. We will call the composition of (m, i)-spiits across all dimensions as a p-split. Given μ for m E N as input, the effect of the y-split S,, is the generators of each zonotope in the result of splitting the generators of the input zonotope scaled by μ. Thus, j--slittinng allows us downscale the generators of a zonotope. For a dl d-dimensional zonotope and μ . Sp returns d zonotopes.
(0))= (25) > (et), (-)(()) is split operator according to Def. J.I. We can generalize
Definition J.2 (μ-Split), For d-dimensional zonotopes and l m for E N, the p-split Sp is defined as:
S Sm,d ....... Sm,,1
33

---

where the (m. sp-st Sm, is defined as shown below and the application of split operator to set of zonotopes defined as applying the split operator every element the set.
(f,) =1E%
for
E = -c m + 2j
9: + gi + Σ9is
1
371
For example. S2 + (1) {0.5 + 0.5c1, 1.5 0.5e1). It is easy that S is indeed split according to Def. J.1.
Proposition J.3 (S,, Split Operator). Sy is a split operator for any p. were N.
Next we introduce join operators for zonotopes and then demonstrate that for any abstract transformer F for function F, we can construct another abstract transformer for F by (i) splitting the input zonotope using some split operator S, (ii) applying r separately on each zonotope returedd by the split, and (iii) merge the results using a join that introduce next.*
Definition J.4 (Join), An operator Li is a join for zonotopes if it over-approximates union of con- cretization. That is, for any set S of zonolopes:
U (2") sfEst
Several join operators have been proposed for zonotopes, e.g., [24]. Next we establish that new abstract transformers be constructed by wrapping existing transformers with split and join. This ensures that our approach of combining our abstract fixed point calculation with splitting yields valid fixed points.
Lemma J.S (Abstract Transformers Are Sound on Splits). Consider zonotope z and let F be an abstract transformer for a function Furthermore, let be a join operator, be a split operator, and let ((tt) 1....... m. Then, Fs. 2 Ft S is an abstract transformer for
Proof. We have to show that
(FSL)((Et) 2 F(y (E)).
We have
F(y (e))
=F( U (Et))
(Def. J.1 for S)
(EE(()F)
m
= 7 (E)))
m =Ü F(r (E.))
(function application sets of concrete elements is point-wise application)
1m1
S
U (F+(et)
(Fi is an abstract transformer)
1=1
(r ( uzb) 27 (zt)u> u (2t)
S7
=y (Fs.u((Et)
'Technically, join is normally detined as an over-approximation of the least upper bound of two abstract elements wrt. partial order of the abstract domain. We order abstract elements based on set inclusion of their concretizations here and only define join regarding to this partial order.
34

---

J.2 Weakening Regularization Requirements
Equipped with the y-split and a join operator, we are ready analyze how the introduction of a y--plii can be used to achieve an arbitrarily small 3. Given the input data DF, the learning rate 1) 2 B, the regularization coefficient X, the transformation matrix A, we use I(D4, 17. X, A) denote the function that computes the fixed point w for these inputs using the process outlined in Corollary 1.7. From Lem. J.5 follows that ss,,.u is abstract transformer for gradient descent and, thus, according
to Prop. 3.3 werr-pppooxmmaes w
Next, we show that given an abstract training dataset Då and a desired regularization coefficient Atarget, we can find a value of μ such that the abstract fixed point construction with jt-splitting computes a fixed point. Recall that I' can construct fixed point if Z where 3 is constant that depends on D We will show that by choosing μ carefully, we can achieve 3 S Atarget for each zonotope in the result of the split and, thus, will return a valid fixed point for each zonotope the split result.
Lemma J.6. For any abstract training dataset DF and desired regularization coefficient X, there exists m E N such that for μ m [ss..u returns a fixed point w²*
Proof. Consider an abstract training dataset D = =X*, with X = X X WLOG consider single abstract training dataset DA; in the result of f--plit on D
Dd. E (D²))
We will show that using the regularization coefficient A, we can find j such that the precondition À for T to compute fixed point holds. Then applying Lem. J.5 we get the desired result.
First observe that since D is in the result of p-splitting, we know that has the following shape where the symbolic component X is scaled by j and the real component XRi of D typically differs from X as j-splitting changes the center of the zonotope.
X = XRi + 4
Recall that 2n C'k denotes the matrix of coefficients that is the sum of absolute values of coefficients of 2 AAx n T x + X Xn + X! T ) Plugging in the definition of X from above, we get n
a matrix C'k that is the sum of absolute values of coefficients of
2 A(xRi Ri T "XX 2
XRi uxx T px3)wi
(41)
n
K16
2 A(pX Ri TX3 uxx XRi + 12x3 A
:
(42)
n
kigea
k1
2p A(X R.
X X T X Rc 1-1 2n² A(X T
X3)A-1
;
(43)
n
Kaea
contributes to c/.g
contributes oZ,,
Recall that ce denotes entry of C' at row i and column j. Here, x T X has higher order than (XR² T Tx + xt XR). thus their terms will not cancel out. Therefore, we can decompose the contribution to G, into the above two parts, denoted as
4, Ghs + eZ.
Given scalar μt,
c max Cu,y u,
i.c., cy denotes the largest entry <. in DA,.
35

---

Let i and such that c = c'.j we have:
4 =et + + S pags + pm = pémaz
where Smaz a constant that is equal to the maximum c'y before scaling X through u-splitting.
Next, we prove the lemma for a specific transformation A V during order reduction derived through singular value decomposition (SVD). Using SVD, we can decompose the covariance matrix X XR, I X nu of the real number part of features X Rn as shown below
X R X Ri VΣVT.
For A = VT, and using the fact that V is rotation in SVD which implies V-1 = VT, we get:
AX A A
=((VVV A
(VE(V"+"))"
-(V-V))XV~-V-11 =
Thus, the matrix AX T X RR Σ is diagonal matrix, meaning that all off-diagonal elements of it (19.,) for all (i+)) j) in Corollary 1.7) are zero. In addition, the diagonal element Eit, i] (gs, in Corollary L.7) is the 1th eigenvalue of the covariance matrix X XRe X Ri. which must be positive when assuming no multicolinearity.
To sum up, 3 from Corollary 1.7 can re-written as:
1 8 n
max (4, +
(lass) + +chg)) (1.s)
(44)
1
max (cp +
Σ
(i.s)
(45)
d
1
max (pc +
ucmax 9t.s)
(46)
n2
3#1
1 max (dpemar a1,a) n
djccmaa min(q
(47)
n
n
where 1 min;( is a positive constant. Therefore, when setting μ n
min,(m.e)
0, we have
de
S 1 n duccmaz min; (q.44) O.
K An Approximate Abstract Transformer for Ridge Regression
Alg. 2 is the detailed version of Alg. we presented in the main paper. Like Alg. 1, Alg. 2 takes as input an abstract dataset DF that over-approximmtes D learning rate 7), regularization coefficient X, and transformation matrix A used for order reduction (e.g., the SVD-based transformation discussed in App. J.2). In Line 1, we first use function constructFPEEuatioon points for wi and w using the closed form solution from Lem. 1.4 and construct the equation system E\t for w (Eq. (39)). Furthermore, this function also computes the threshold 3 on the regularization coefficient for this system of equations to have solution. If 3 is smaller or equal to the desired regularization coefficient À (Line 10), then solve the linear equations 3 (Line 11) and then returns the fixpoint wi * Otherwise, we determine a sufficiently small splitting factor j that ensures that 3 / X for each unt in the result of the split (Line 3), compute the fixed point for each such w (Lines 6 and 8), and merge these fixed points using join to compute the final result (Line 9).
to compute fixed
36

---

Algorithm 2: Constructing Fixed Point for Abstract Gradient Descent
Input:
abstract dataset D (X*, , )²), learning rate 7), Regularization coefficient X, transformation matrix
Output: abstract fixed point w²
wk, wh" E, B constructtPoouaations(> 1))
1 > then
// Is splitting necessary?
determineelumSSlittt......
[D', th < (D²) for D [m do 1
// Apply p-splitting // Construct FP for each D in the split
wh: w 5 costraccePPPPeeettooon1( X,1)
w'; wi wh +
return w
// Merge FPs using Join
10 else
// No splitting required
11
sdleeFixcicocnn@taaatiooSS... u wh wh + u
return w
14 def ConstructPPEqutiionsDD/..
W eralClosddFormheeal(DD....
// Eq- (22)
16
+ eraCCCoseeFormmyyyoooolcC..Aaaaa
Eq. (23)
caearracccoaaaaairiiaeeetttttttt n, w w
// Eq. (39)
18
return wt wh E, 8
19
L Infeasibility of Symbolically Evaluating the Closed-form Solution for Ridge Regression
In this section, demonstrate why evaluating the closed-form solution for linear regression with MSE is not feasible. It produces symbolic expressions with many variables and large symbolic terms (their representation size) that include the fractions with denominators and enumerators that are polynomial expressions. Apart from the size of these expressions, computing viable prediction intervals from such symbolic representation of model weights computationally hard. We use MSE loss for simplicity, but the same arguments also apply to ridge regression.
We illustrate the issues using mandomly-generatte toy dataset with 10 samples and 3 uncertain data cells the third column (corresponding to 3 error symbols c0, 61, c2)
1.0
1.6
1.3cg + 1.4
4.6
1.0
0.8
1.5
- -1.5
1.0
-1.7
-1.9
0.39
1.0
-0.57
1.1
- 5.7
1.0
-0.39
0.36
3.0
X
1.0
0.035
1.2
3.6
1.0
-0.34
-0.73
00.4
1.0 0.038
0.3e, 0.44
0.46
1.0
1.5
1.2G2+14
2.2
1.0
0.98
-0.66
3.0
'The first column of X consisting of I's corresponds to the bias term the intercept in the weights.
37

---

Then, symbolically evaluating the closed-fom solution for linear regression with MSE loss, we obtain the symbolic expression representing all possible model weights:
w + (X:²" x y
1.5-10
0|0000/0000001100000000000000100 - -007777700770 000
Next, we consider one test sample æt = 1, and apply this closed-form model weight expres- sion to infer its prediction range (*+) (cf. Corollary F.2). With merely 3 uncertain data points in the training data, the prediction of this test data point
yt- = - xt
1.2 1.3 4f
3.7 3.7
8.0-10 +rg.
-1.1.00
22
1.0.10
Lo:110 2 8 - 10000000000000000000000000000000 20t.:/00.00....00000000000201110
21
1.8.10
11111100 f..4 11 1o0 8.3 1.3 17..1f1177....1111111...... 2.0.10000000000 12.6.1002
10
(48)
already consists of fractions of complex polynomial expressions. Finding the viable prediction range for this data point, i.c. the minimum and maximum possible prediction for the point across all models in the concretization of the symbolic expression, using this expression is infeasible, as it is harder than finding an extrema for multivariate polynomials under linear constraints which is known to be NP-hard. In practice, expressions will be significantly larger as they will involve more error symbols. Assume each column has p uncertain cells, and is the number of uncertain labels, the number of distinct monomials is O(p" of q), where d is the number of dimensions, and the d in exponents mainly comes from the matrix inversion when computing the determinant.
In summary, the symbolic expressions obtained from symbolically evaluating the closed- form solution for linear regression are not suitable for representing the space of possible model weights.
M Additional Experiments
M.1 Robustness Verification Additional Results
For a better comparison with MEYER [47], we heatmaps visualize the robustness ratios averaged over repeated experiments with different random seeds. In App. M.1, ZORRO exhibits much less yellow regions (representing higher uncertainty) compared MEYER. Especially when the uncertainty is high (bottom-right part). there are many cases that ZORRO returns high robustness ratio but MEYER does not. Recall that both approach gives sound over-approaximation of prediction robustness, thus the robustness ratios returned by both are lower bounds of the ground truth prediction robustness ratio. Therefore, the aforementioned cases correspond to highly robust ground truth (no less than the result of ZORRO), where ZORRO gives tight over- ape-xaproximations but MEYER does not.
M.2 Micro Benchmark
We use range of loss the metric measure the quality (tightness) of uncertain training result. range of loss for samples are calculated by max loss min loss across all sampled results. Loss range for symbolic representations are calculated my measuring the range of the interval concretization of the loss function result. Notice that samplings return a subset of all possible result which has no soundness guarantees, as a result, produces an under-approximation of the loss range. The symbolic fixed point approach produces an over-approximated loss range. We take 1,000 (Ik) and 10,000 (10k) samples from all possible worlds as under approximations the concrete range when number of possible worlds are to large to compute the ground truth range.
'The expression is a sum of fractions in which both numerator and denominators are multivariate polynomial expressions. The range 1,, 1] for each error symbol can be formulated as linear constraints over the input space.
38

---

Meyer al.
ZORRO
100
Uncertain Radius (86)
100
100
97.2
80.8 67.5
60.2
48.8
40.5
Uncortain Radius (86
0.05
100
100
100
100
Robustries Ratio (%)
0.1
98.2
89.0
64.5
45.0
26.5 17.2
11.8
72
4.8
0.1
100
100
100
100
1.00
1.00
100
100
100
60
100
08.5
57.5
32.0
18.0
6.8
3.5
1.2
0.15
100
100
100
100
100
1.00
100
100
100
0.2
07.5
68.2
35.2
15.5
6.0
1.0
0.2
100
100
100
100
100
1000
97.7
20
94.8
49.8
20.8
5.8
1.5
0.25
100
1.00
100
100
100
98.2
93.5
87.2
74.0
66.5
0.01
0.03
0.04
6.05 0.06 0.07
0.06
0.09
0.1
0.01
0.02
0.03
0.04
0.05
0.06
0.07
0.08
0.09
Percentage of Uncertain Data
Percentage of Uncertain Data
(a) Mpg
Meyer al.
ZORRO
100
Uncertain Radius (86)
0.02
Uncortain Radlus (56)
0.02
100
1.00
Robustriss Rablo (%)
100
100
100
99.1
86.0 62.0
36.3
21.7
14.3
9.7
0.04
100
100
100
100
100
99.9
99.6
98.4
94.9
86.3
60
0.06
300
100
95.9
64.7
27.5
14.3
8.0
4.0
1.5
0.4
0.06
100
100
100
99.9
99.3
93.8
63.7
06.7
46.0
20.3
-40
20
0.4
100
100
Percentage of Uncertain Data
Percentage Uncertain Data
(b) Insurance
Figure 5: Robustness verification under label errors using intervals (MEYER) or zonotopes (ZORRO).
Zono
0:30
0.35
Zono
0.15
0.16
Zoro
U22
20.5
8012
0.11
0.10
0.10
0.02 OLdk 0.07 0.00 0.22 uncertain ralio
0.19 0.25
unceriainly radius
(a) varying amount of uncertainty
(b) varying uncertainty radius
(c) varying dimension
Figure 6: Fixed point versus samplings
Varying ratio of tuples contains uncertain data Figure 6a shows the loss range results while updating the ratio of uncertain rows from 0.01 0.12. ZORRO has a tight range considering all samples are under approximations of ground truth range especially in As uncertainty increases, the over-approximation gap widens due to the increased coefficient of higher order terms in the gradient, which are linearized, leading to higher linearization errors.
Increasing uncertainty radius Figure 6b shows the result by increasing the radius of uncertain data (imputation results) by multiplier of 0.06 to 0.44 where 0.06 means the range for each uncertain data is increased by 6%. Similar to amount of uncertainty.
Varying dimensions Figure 6c Shows the result by change number of dimensions to the training data. Result indicates tightness of ZORRO 's over-approximation is not affected by the dimension of the data.
Figure 7 added sampling result to Figure 3 shows additionally that sampling result is under approximation of the ground truth result.
M.3 Additional Experiments on Varying Regularization Coefficients
Fig. 8 shows the effect of varying regularization coefficients on the worst-case loss and prediction robustness ratio. To avoid much overlapping in the plots, we used one standard deviation as the error bar. Similar to the conclusions in Sec. 5.2, X = O is often not the optimal regularization coefficient in fact, a small, positive could result in higher accuracy and better
terms of accuracy or robustness. robustness, and this optimal À varies across different settings. In general, higher data uncertainty requires higher X, which coincides with the intuition of using regularization in traditional settings of linear regressions, which aims to mitigate the effect of data noises errors.
39

---

0.24
0.13
0.20
0.12
Zorro
g 00.30
0.11
0.12
GT
0.08
0.09
0.04
0.08
0.00
0.07
0.06
0.12
0.31 0.38 range multiplier
0.19
0.25
0.44
6
11
13
dimension
(a) varying range size
(b) varying data dimension
Figure 7: The range of losses obtained by enumerating all possible worlds (GT), sampling 1000 possible worlds (1k), and ZORRO.
Uncertain
Uncertain
Uncertain
Uncertain
Rebustness Ratio
Worst case
2.70
0.90
0.95
90
0.8
0.15
2.60
0.85
0.6
3:10
0.10
2.70
3.20
0.05
Figure 8: Results with uncertain training labels.
Figure 9: Robustness ratio (red y-axis) and worst-case test loss (green y-axis) vs. regularization coefficient Ả (x-axis), with varying percentages of uncertain labels.
M.4 Comparing with Bayesian Regression
We ran empirical evaluations to demonstrate how data quality issues pose challenges for Bayesian linear regression (implemented with torchbnn [40]), making them inapplicable the setting from the third plot in Fig Ic, where the uncertain data percentage set to 10%, we tested Bayesian linear regression on different possible worlds using two methods: impute-and-pe -nd-eeddd---dd and sampling from possible worlds. The results show that the prediction intervals generated by Bayesian methods do not cover the ground truth prediction, i.e., the prediction by the model trained on the ground truth training data. In contrast, our approach guarantees 100% coverage across all cases.
setting. Using
Imputation
Imputation
Imputation
Radius (%)
Coverage (%)
Avg. Intv.
Coverage (5)
Avg. Inty.
Coverage (4)
Avg. Inv.
Coverage (%)
Avg. lnts.
2.5
48.8
0.157
42.5
0.153
46.3
0.176
100
0.183
0.176
100
38.5
0.157
37.5
35
0.176
100
10
36.3
0.157
31,3
0.154
32,5
0.176
100
2,387
Table Bayesian approach on imputed data.
Uncertain Radius (%)
Possible
Possible World
Possible World
Coverage (%)
Avg. Intv.
Coverage (4)
Avg. Intv.
Coverage (%)
Avg. Ins.
Coverage (%6)
Avg. loty.
2,5
58.75
0.151
0.154
55
0,176
190
0.183
30
0.158
0.176
100
42.5
0.155
100
10
57.5
0.351
31.3
0.169
21.5
0.168
100
2.587
Table 2: Bayesian approach on sampled possible worlds.
M.5
Varying Robustness Threshold
Different practical applications may differ in how much uncertainty they are willing to tolerate, thus leading to varying choices of robustness thresholds. To account for this in the experiments, we selected both a low threshold (0.5% for the insurance data) and high threshold (5% for the MPG data). As presented in Fig 10, we also explored other thresholds, which did not impact the trends significantly: Zorro can consistently certify a larger fraction of the test data points than the baseline due to its use of the more expressive zonotope domain.
40

---

Robustness Threshold: 0.4 Rabustness Threshold: 0.8 Robustness Threshold: 1.6 % Rebustness Threshald:
10
1,0
1.0
1.0
Robustness
0.
0.5
0.5
0.5
0.0
0.0
0.0
0.0
0.0
5.0
7.5
5.0
2.5
10.0
0.0
7.5
Uncertain Data Percentage (%)
20RRO
-8- Meyer et
Figure 10: Varying robustness threshold from 0.4% to 2.4% with fixed uncertainty radius 6%
M.6
Varying Uncertain Feature
The variability of trained models will depend on the correlation that features with uncertainty have with the label, which then also affects over-approximation of this set. As rule of thumb, model variability will increase with the correlation between the uncertain feature and the labels, leading to less robust predictions. We conducted an additional experiment using the MPG dataset, focusing on the feature "acceleration," which has relatively low correlation with the label. Unlike the feature "weight" (Fig. Ic), where the robustness drops when the uncertainty radius is 10%, the robustness for "acceleration" starts to drop only when the uncertainty radius increased to 16%. The result shown in Fig 11 indicates that features with lower predictive power have less impact on the robustness of model compared to features that are highly correlated with the label.
Robustness Ratio
Uncertaing Radius:
Uncertainty Radius:
Uncertainty Radius:
Uncertainty Radius:
1.0
1.0
0.5
0.0
23
5.0
7.5
10.0
0.0
2.5
5.0
7.5
10.0
0.0
2.5
5.0
7.5
10.0
0.0
25
5.0
7.5
10.0
Uncertain Data
Percentage (%)
ZORRO
Uncertainty
Uncertainty Radius:
Uncertainty Radius:
Uncertainty
1.0
1.0
1.0
1.0
Robustness
0.5
0.5
0.5
0.5
0.0
0.0
0.0
0.0
25
5.0
7.5
10.0
0.0
2.5
5.0
7.5
10.0
0.0
2.5
5.0
7.5
10.0
0.0
25
5.0
7.5
10.0
Uncertain Data Percentage (%)
ZORRO
Figure 11: Robustness verification with errors in feature "acceleration".
M.7 Runtime of ZORRO
We evaluated the runtimes for computing the closed-form solution on the MPG dataset, varying the numbers of uncertain data points, in Fig 12. With the same number of uncertainty data points, uncertain features lead to more complex computations covariance matrices, and result in higher runtimes compared to label uncertainty, where the covariance matrix remains real-valued.
Runtime (s)
Runtime (s)
0.02
0.04
0.06
0.08
0.10
0.02
0.04
0.06
0.08
0.10
Uncertainty Percentage
Uncertainty Percentage
(a) Uncertain feature
(b) Uncertain label
Figure 12: ZORRO runtimes MPG data.
41

---

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 state our contributions and scope in Sec. 1,
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 in the paper and important assumptions and limitations. A No or
NA answer to this question will not perceived well by the reviewers.
The claims made should match theoretical and experimental results, and reflect how be expected to generalize other settings.
much the results
is fine to include aspirational goals as motivation as long is clear that these goals are not attained by the paper.
2. Limitations
Question: Does the paper discuss the limitations of the work performed by the authors?
Answer: [Yes]
Justification:
We discuss the limitations of our approach in Sec. and app. and J.
Guidelines:
The answer NA means that the paper has no limitation while the answer No means that the paper has limitations, but those 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 how these assumptions might be violated in practice and what the implications would be.
The authors should reflect on the scope of the claims made, e.g., if the approach was only tested on few datasets or with few runs. In general, empirical results often depend on 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 images are taken in low lighting. Or speech-to-text system might be used reliably to provide closed captions for online lectures because it fails to handle technical jargon.
The authors should discuss the computational efficiency of the proposed algorithms and how they scale with dataset size.
If applicable, the authors should discuss possible limitations of their approach address problems of privacy and faimess.
While the authors might fear that complete honesty about limitations might be used by reviewers as grounds for rejection, a 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 an 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 Proofs
Question:
For each theoretical result, does the paper provide the full set of assumptions and complete (and correct) proof?
Answer: [Yes]
42

---

We provide intuitions in the main body of the paper (Sec. 3 and 4). We include the complete proofs in the appendix, and refer to them in the main body of the paper (App. and H to J)
Justification:
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 the statement of any theorems.
The proofs can either appear in the main paper or the supplemental material, but if they appear in the supplemental material, the authors are encouraged provide short proof sketch to provide intuition.
Inversely, any informal proof provided in the core of the paper should be complemented by formal proof's provided in appendix or supplemental material.
Theorems and Lemmas that the proof relies upon should be properly referenced.
4. Experimental Result Reproducibility
Question: Does the paper fully disclose all the information needed to reproduce the main ex- perimental results of the paper to the extent that it affects the main claims and/or conclusions of the paper (regardless whether the code and data are provided or not)?
Answer: [Yes]
Justification: Descriptions in the main body of the paper (Sec. 5), and the formulas and technical details in App. I and J are sufficient for reproducing our approach and experimental results.
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 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 a dataset and/or model, the authors should describe the steps taken 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, describing the architecture fully might suffice, or if the contribution a specific model and empirical evaluation, may be necessary to either make it possible for others to replicate the model with the same dataset, or provide access to the model. In general. releasing code and data 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., in the case of large language model), releasing of a model checkpoint, or other means that are appropriate to the research performed.
While NeurlPS does not require releasing code, the conference does require all submis- sions to provide some reasonable avenue for reproducibility, which may depend on the nature of the contribution. For example
(a) If the contribution is primarily a new algorithm, the paper should make it clear how to reproduce that algorithm.
(b) the contribution is primarily a new model architecture, the paper should describe the architecture clearly and fully.
(c) If 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- open-source dataset or instructions for how to construct the dataset).
(d) We recognize that reproducibility may be tricky in some cases, in which authors are welcome to describe the particular way they provide for reproducibility. In the case of closed-source models, may be that access to the model is limited in some way (e.g., to registered users), but it should be possible for other researchers have some path reproducing or verifying the results.
43

---

5. Open access to data and code
Does the paper provide open access to the data and code, with sufficient instruc- tions to faithfully reproduce the main experimental results, as described in supplemental material?
Question:
Answer: [Yes]
We share the anonymize Github repository containing the data and code in the paper.
Justification:
Guidelines:
The answer NA means that paper does not include experiments requiring code.
Please see the NeurlPS code and data submission guidelines (https://nips. cc/
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 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: / //niss. cccccc/pp.o wwwpccpwpcc//eeeecceeeceeeeeetle wwtcc//gggg//eeeeeeeeeeee Sotps//ooeeeuesiioiooonnoo ssioolllicc 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. 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 as much information as possible in supplemental material (appended to the paper) recommended, but including URLs to data and code permitted.
6. Experimental Setting/Details
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?
Question:
Answer: [Yes]
Justification:
We describe the experiment settings in detail in Sec.
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 level of detail that is necessary to appreciate the results and make sense of them.
The full details can be provided either with the code, in appendix, 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: We repeat experiments with different random seeds and show the mean and standard deviation the results.
Guidelines:
The answer NA means that the paper does not include experiments.
The authors should answer "Yes" if the results are accompanied by emor bars, confi- dence intervals, or statistical significance tests, least for the experiments that support the main claims of the paper.
44

---

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 run with given experimental conditions).
The method for calculating the error bars should be explained (closed form formula, call to 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 the standard error of the mean.
It is OK to report 1-sigma error bars, but one should state it. The authors should preferably report a 2-sigma error 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 show in tables figures symmetric error bars that would yield results that are out range (e.g. negative error rates).
If error bars are reported in tables or plots, The authors should explain the text how they were calculated and reference the corresponding figures or tables in 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 describe the sufficient information on the computer resources in Sec. 5.
Guidelines:
The answer NA means that the paper does not include experiments.
The paper should indicate the type of compute workers CPU or GPU, internal cluster, or cloud provider, including relevant memory and storage.
The paper should provide the amount compute required for each of the individual experimental runs as well 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 into the paper).
9. Code Of Ethics
Question: Does the research conducted the paper conform, in every respect, with the NeurlPS Code of Ethics httpss//neurs: hww.rpuuu..pp/.e cc/ cctppblic///eccccecccecccceccllt ceuid ines?
Answer: [Yes]
Justification: We reviewed the NeurIPS Code of Ethics and make sure the paper conforms with it.
Guidelines:
The answer NA means that the authors have reviewed the NeurIPS Code of Ethics.
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 (c.g.. if there 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: [Yes]
Justification:
We discuss the broader impacts of our technique in Sec. 6.
Guidelines:
The answer NA means that there is no societal impact of the work performed.
45

---

If the authors answer NA or No, they should explain why their work has no societal impact or why the paper does not address societal impact.
Examples of negative societal impacts include potential malicious or unintended uses (e.g., disinformation, generating fake profiles, surveillance), faimess 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 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 an 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 a 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 to monitor how system learns from feedback over time, improving the efficiency and accessibility of ML).
11. Safeguards
Does the paper describe safeguards that have been put in place for responsible release of data models that have high risk for misuse (e.g., pretrained language models, image generators, or scraped datasets)?
Question:
Answer: [NA]
Justification: This paper has no such risk.
Guidelines:
The answer NA means that the paper poses 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 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 challenging, and many papers do not require this, but encourage authors to take this into account and make best. faith effort.
12. Licenses for existing assets
Question: 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?
Answer: [Yes]
Justification: We credit all existing assests used in the paper, and include the license information in the reference.
Guidelines:
The answer NA means that the paper does existing assets.
The authors should cite the original paper that produced the code package or dataset.
The authors should state which version of the asset used and, possible, include a 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.
46

---

If assets are released, the license, copyright information, and terms of use in the package should be provided. For popular datasets, paperswi thcode. com datasttss has curated licenses for some datasets. Their licensing guide can help determine the license of dataset.
For existing datasets that are re-packaged, both the original license and the license of the derived asset (if it has changed) should be provided.
If this information 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: [Yes]
Justification:
We provide documentations for the released anonymized code (Sec. 5),
Guidelines:
The answer NA means that the paper does not release new assets.
Researchers should communicate the details of the dataset/code/moodl as part of their submissions via structured templates. This includes details about training, license, limitations, etc.
The paper should discuss whether and how consent obtained from people whose asset is used.
At submission time, remember to anonymize your assets (if applicable). You can either create an anonymized URL or include an anonymized zip file.
14. Crowdsoureing 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 details about compensation (if any)?
Answer: [NA]
Justification: This paper does not involve crowdsourcing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing nor research with human subjects.
Including this information 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 to the NeurIPS Code of Ethics, 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 institution) were obtained?
Answer: [NA]
Justification: This paper does involve crowdsoureing nor research with human subjects.
Guidelines:
The answer NA means that the paper does not involve crowdsourcing research with human subjects.
Depending 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.
47

---

We recognize that the procedures for this may vary significantly between institutions and locations, and we expect authors adhere the NeurIPS Code of Ethics and the guidelines for their institution.
For initial submissions, do not include any information that would break anonymity (if applicable), such as the institution conducting the review.
48


Reviews:
[{"review_id": "E7UquVlSDH", "reviewer": "Reviewer_h3Xr", "summary": "This paper attempts to represent the uncertainty from dataset variations using a type of convex polytope that can compactly represent high-dimensional spaces, as dataset multiplicity from imperfections like missing data may scale exponentially. The paper develops techniques to ensure the tractability of gradient descent and obtain a fixed solution for linear models in this regime. Over-approximations of possible model parameters are provided by performing suitable order reductions and linearizations to ensure that the gradient descent is tractable. Meyer et al. (2023) introduce the dataset multiplicity problem but only account for uncertainty in output labels. This paper also allows uncertainty in input features, thus providing a more general framework to represent and learn under the regime of dataset multiplicity. They also provide a SymPy implementation of their method and perform an experimental evaluation of their techniques.", "questions": "Are there scenarios where the overapproximation of dataset variations and corresponding model parameter variations results in too broad outcome intervals? I would like to understand the limitations of the proposed method apart from the linearity assumption. \n\nIn the experimental evaluation, the dataset uncertainty is varied by changing the uncertain data percentage and uncertainty radius. Did the authors perform any experiments with missing data as the source of uncertainty? (I imagine that you could represent missing data by using a high uncertainty radius). How was the robustness threshold set for the experiments?", "limitations": "At the moment, the last section of the paper (Conclusions/broader impact/limitations) does not include any explicit limitations of their methods. It would be helpful to have a separate discussion on limitations in the main paper instead of appendix. Broader societal impacts may not be relevant for this work.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "This work provides a general framework to represent dataset multiplicity in both features and labels, allowing us to represent the potentially exponential number of dataset variations with a suitable convex polytope. The paper carefully performs several order reductions and linearizations, overapproximating the dataset variations and fixed point model parameter variations. Although limited to linear models, this work addresses an important source of uncertainty stemming from dataset variations and non-trivially extends existing work in the area by incorporating uncertainty in input features in the dataset. \n\nThe SymPy implementation of the techniques introduced in the paper helps in evaluating this method extensively.", "weaknesses": "The techniques introduced in the paper are only applicable to linear models.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their thoughtful review. Below are the response to reviewer's questions/comments.\n\n**Weaknesses: The techniques introduced in the paper are only applicable to linear models.**\n\nPlease see our response in the main rebuttal.\n\n**Are there scenarios where the overapproximation of dataset variations and corresponding model parameter variations results in too broad outcome intervals? I would like to understand the limitations of the proposed method apart from the linearity assumption.**\n\nOver-approximation of a set of possible worlds in the training data may occur, for instance, when an uncertain value with few discrete possible values is represented through a range, leading to a superset of the actual possible worlds (all values in the range as opposed to all the discrete set of values). For instance, consider the following extreme example: there are two possible worlds in the training data where feature $A$ is $-c$ in one world for all data points while in the other world it is $c$. Using ranges leads to a large over-approximation as the ranges will contain all values between the two extremes $[-c,c]$. This leads to a large increase of the number of possible worlds which could possibly lead to a large over-approximation of prediction ranges. However, in our experiments, we did not encounter significant issues with over-approximation. The over-approximation error in our setup mostly resulted from over-approximation during learning. We will enhance the discussion on these limitations in the main paper to provide a more comprehensive understanding.\n\n\n**In the experimental evaluation, the dataset uncertainty is varied by changing the uncertain data percentage and uncertainty radius. Did the authors perform any experiments with missing data as the source of uncertainty? (I imagine that you could represent missing data by using a high uncertainty radius). How was the robustness threshold set for the experiments?**\n\nCorrect, missing values can be modeled as ranges with a large uncertainty radius. Regarding experiments with missing data: our experiments with ground truth results are based on missing values. We randomly selected some data points to have missing values and then for each missing value determined a range by running several imputation techniques to get a set of possible values $G = \\{c_1, \\ldots, c_n\\}$. We then use the smallest range that covers all values returned by the imputers: $[ min_{c \\in G} c, max_{c \\in G} c ]$. For these experiments (Figures 3 and 6), we did assume that each ground truth value $c$ is equal to one of the guesses returned by the imputers ($c \\in G$). For experiments where we controlled for the *uncertainty radius* (Figure 3(a) and 6(b)), we did further over-approximate the ranges to achieve the desired uncertainty radius. We will add clarifications to the main paper and add a detailed description of how uncertainty was generated in the experiments to the appendix.\n\nRegarding the robustness threshold, the rationale for this choice was that different practical applications may differ in how much uncertainty they are willing to tolerate. To account for this in the experiments, we selected both a low threshold (0.5% for the insurance data) and a high threshold (5% for the MPG data). We also explored other thresholds which did not impact the trends significantly: Zorro consistently can successfully certify a larger fraction of the test data points than the baseline due to its use of the more expressive Zonotope domain. We have included additional results in the main rebuttal varying the robustness threshold (Figure 1 in the attached pdf), and we will add these plots to the appendix.\n\n**Limitations: At the moment, the last section of the paper (Conclusions/broader impact/limitations) does not include any explicit limitations of their methods. It would be helpful to have a separate discussion on limitations in the main paper instead of appendix. Broader societal impacts may not be relevant for this work.**\n\nWe will include a summarized discussion on the limitations of our techniques into the main paper and expand on this discussion in the appendix. As mentioned above, there are potential scenarios when the over-approximation can be large because of over-approximation of the uncertainty in a training dataset when there are few possibilities for an uncertain value induce large ranges or in general when the uncertainty in the training data is highly nonlinear.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "Obr6KW7EyF", "reviewer": "Reviewer_7DaV", "summary": "The paper introduces a method called ZORRO that allows for the over-approximation of the set of ridge regression models learned over a set of possible datasets. The authors use abstract interpretation, propose abstract gradient descent algorithms, and develop efficient techniques to learn the set of models. The empirical results illustrate that ZORRO is effective in robustness verification and other analyses related to linear regression.", "questions": "1. Could you please comment on why the robustness threshold was set to different values for the two datasets?\n2. For modeling uncertainty in features, the distribution was built around one feature. Does the choice of the feature, such as based on the correlation value between the feature and the label, affect the results\n3. How much time does it take to compute the set of possible linear models for the discussed datasets?", "limitations": "The authors discuss the limitations of the paper.", "rating": 7, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 3, "strengths": "1. The paper proposes a framework that can potentially be used for various methods of defining uncertainty in datasets.\n2. The method seems to be efficient in computing prediction ranges for the overapproximation of set of models learned from uncertain data and discusses how this set can be useful for causal inference applications.\n3. ZORRO's results remain consistent regardless of the number of dimensions for ridge regression. The variance of the robustness ratio and worst case loss decreases with larger regularization coefficients.\n4. The authors provide an extensive overview of related literature from related fields. Literature on the Rashomon Effect, which generalizes predictive multiplicity, might also be of additional interest.\n5. Overall, the paper is well-written, and the problem of learning from uncertain datasets is interesting.", "weaknesses": "1. As the authors acknowledge in the Limitation section, it is unclear how to scale the proposed approach to other, more complex hypothesis spaces or problems.\n2. Despite ZORRO being a more general approach, sampling seems to provide results closer to the ground truth even with 1k samples (Figure 7).\n\nMinor:\n1. Please specify that the numbers on figures on page 5 represent gradient descent iterations.\n2. Definition 2.1 and the surrounding text around it is very abstract. While Appendix D is helpful, examples become clearer after reading Section 2.2. Providing one sentence with specific examples of uncertain data around Definition 2.1 might be helpful.", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their thoughtful review and suggestions. We will mention the numbers on figures on page 5 represent the iterations, and we will add a brief example of uncertain data around Definition 2.1 to make this more concrete. In response to the reviewer's other questions/comments:\n\n**1. Could you please comment on why the robustness threshold was set to different values for the two datasets?**\n\nDifferent practical applications may differ in how much uncertainty they are willing to tolerate. To account for this in the experiments, we wanted to test an application with a very tight threshold (0.5% for the insurance data) and a use case with a looser requirements (5% for the MPG data). We also explored other thresholds which did not impact the trends significantly: Zorro can consistently certify a larger fraction of the test datapoints than the baseline due to its use of the more expressive Zonotope domain. We have included the results of one such experiment in Figure 1 in the main rebuttal pdf document and will include them in the appendix.\n\n**2. For modeling uncertainty in features, the distribution was built around one feature. Does the choice of the feature, such as based on the correlation value between the feature and the label, affect the results**\n\nYes, the variability of models will depend on the correlation that features with uncertainty have with the label, which then also affects our over-approximation of this set. As a rule of thumb, model variability will increase with correlation between the uncertain feature and the labels, leading in turn to less robust predictions. We conducted an additional experiment using the MPG dataset, focusing on the feature “acceleration,” which has relatively low correlation with the label. Unlike the feature “weight” (Figure 1(c) of the paper), where the robustness drops when the uncertainty radius is 10%, the robustness for “acceleration” starts to drop only when the uncertainty radius is increased to 16%. As expected, this result indicates that features with lower predictive power have less impact on the robustness of a model compared to features that are highly correlated with the label. Plese see Figure 2 in the main rebuttal pdf document that will include in the appendix. We will also convey this intuition in the main paper.\n\n**3. How much time does it take to compute the set of possible linear models for the discussed datasets?**\n\nWe show the runtimes for computing the closed form solution on the MPG dataset, varying the numbers of uncertain data points, in Figure 3 in the main rebuttal pdf. \n\nWe will include these results in the appendix and add a brief discussion to the main paper. To evaluate what can be achieved by moving to a faster library for symbolic expressions, we implemented some components of the algorithm (symbolic matrix multiplication) with the generic symbolic computation C++ library GiNaC, which gives > 40x speedup on symbolic matrix multiplications. Therefore, we can anticipate a runtime of at most 0.2 seconds for 10% uncertain MPG data if we switch to an implementation using GiNaC or a custom C++ implementation.\n\nNote that SymPy supports more general expressions than needed for our approach. It internally represents symbolic expressions as expression trees which are python objects. This is relatively inefficient for representing linear or even polynomial expressions. We plan to add runtime results using GiNaC. \n\n**Weaknesses: 1. As the authors acknowledge in the Limitation section, it is unclear how to scale the proposed approach to other, more complex hypothesis spaces or problems.**\n\nIt is true that for linear models we benefit from the fact that we were able to derive a closed form solution which has the advantage of avoiding the repeated injection of an over-approximation that happens in each step of gradient descent. We are currently developing methods to (1) further improve the over-approximation quality by considering more complex abstract elements such as (constrained) polynomial zonotopes, and (2) scale our approach to gradient descent for larger models, e.g., deep neural networks, which involves parallelization, optimized data structures, etc. While challenging, we believe that it will be worthwhile to explore these future research directions.\n\n**Weaknesses: 2. Despite ZORRO being a more general approach, sampling seems to provide results closer to the ground truth even with 1k samples (Figure 7).**\n\nWhile sampling provides in this experiment an approximation that is closer to the ground truth, note that sampling inherently under-estimates the possible outcomes and, thus, cannot be used for high-stake decision-making applications where guaranteeing robustness is crucial for ensuring dependability and safety of models and systems. For such applications, falsely claiming robustness based on sampling is unacceptable. In Figure 7, sampling bounds are close to the ground truth ones because we have a small amount of uncertain data points (8 data points). This was necessary to ensure that we could compute the ground truth by enumerating all worlds (even with just k possible options per uncertain value, there are $n^k$ possible worlds if there are $n$ uncertain values in the training dataset). Thus, we did not conduct a corresponding experiment with more uncertainty as it is only computationally feasible to obtain the ground truth ranges for low uncertainty settings with small $n$ and $k$ . We will clarify this in the paper. In summary, (i) sampling cannot provide guarantees on covering the ground truth and (ii) in these experiments sampling did behave relatively well as we had to choose scenarios with a small number of possible worlds.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}, {"review_id": "pyQ7ARSyTM", "reviewer": "Reviewer_DX9v", "summary": "The paper introduces a new approach for uncertainty quantification of linear models. The proposed approach is interesting, but the empirical comparison to other means of uncertainty quantification is quite limited.", "questions": "* The paragraph starting from line 74 appears to argue that one cares more about possible models (training-time) than possible predictions (test-time). Why is that? This is not apparent to me in most practical applications (including the experimental settings of this work).\n* The paragraph starting from line 24 argues that it is impractical to train a separate model for each of the potentially infinite scenarios a dataset might represent. Quite commonly this is addressed by biasing towards structural simplicity, e.g., selecting the model(s) with the \"shortest description length\" (see Occam's razor). I am missing a reference to this very common approach. Why is this approach not sufficient?\n* In Section 2.1 I am wondering where the set of datasets come from. What is the benefit of incorporating, e.g., missing values as suggested in App. D as opposed to Bayesian approaches? I believe it would help to give one practical motivating example in the main text.\n* In Appendix D I did not find the definition of $\\epsilon$. What is it? And why isn't the maximum entropy distribution used in this case of missing data? Also, I am wondering what the relevant concretizations would look like. Perhaps this section can be extended.\n* Why is not the size of the prediction set (which is common, e.g., in conformal prediction) used in place of the \"robustness ratio\"? Also, perhaps something like \"efficiency ratio\" captures the meaning of this ratio more accurately.\n* How are the robustness thresholds chosen?\n* Why are the abovementioned means of uncertainty quantification / confidence estimation not applicable in the experimental settings of the paper?\n* For Figure 3: How exactly were the errors injected? How was the multidimensional data generated?", "limitations": "An extension of the proposed approach beyond linear models seems difficult.", "rating": 5, "confidence": 1, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "* The proposed approach for propagating the uncertainty about the correct linear model seems interesting.\n* The proposed approach empirically constructs tighter confidence sets than a prior work.", "weaknesses": "My main concern is that this work is lacking a proper empirical comparison to most confidence set / uncertainty quantification methods such as Bayesian linear regression, conformal prediction, and likelihood ratios (e.g., [1]).\nFor example, the manuscript could benefit from a more concrete argument for why (and when) the proposed approach would be advantageous over probabilistic (Bayesian) inference.\nI believe, provided a more thorough empirical validation, this work could add value to a conference.\n\nIn my opinion, the paragraph starting from line 82 dismisses such approaches too easily. From what I can see, reference 26 does not provide conclusive evidence that noise models cannot estimate effects of \"data quality issues\" or that estimating \"uncertainty about the model parameters\" is inherently insufficient (which is in effect what the presented approach does by returning a set of candidate models).\n\n[1]: Emmenegger, Nicolas, Mojmir Mutny, and Andreas Krause. \"Likelihood ratio confidence sets for sequential decision making.\" Advances in Neural Information Processing Systems 36 (2024).\n\nTypos: lines 285, 331", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "We thank the reviewer for their thoughtful review. Due to space limitations, we have summarized each question and numbered them from Q1 to Q9.\n\n**Q1. Are possible models more important than possible predictions?**\n\nWe did not mean to argue that understanding possible models is more important than understanding possible predictions. We wanted to convey that our approach and adversarial robustness differ in which sources of uncertainty they support: uncertain test data in the case of adversarial robustness and uncertain training (and optionally also test) data for our approach.  In our setting, we have to reason about possible models even if we are only interested in possible predictions as the uncertainty in the training data leads to model uncertainty which in turn also affects prediction uncertainty. We will rephrase this paragraph to make this distinction more clear. There are applications, like the causal inference use case we cover in the experiments, where the outcome of interest are the parameters of the model.\n\n**Q2. (line 24) why not use structurally simple models?**\n\nWhile selecting models with structural simplicity can be effective under certain assumptions—specifically where a model exists that performs near optimal across all possible worlds—it is not applicable in situations where data quality issues result in significant variability across the possible worlds. In such cases, this approach may yield a highly suboptimal model. However, the primary application of our approach is to evaluate or certify robustness for a given learning procedure based on worst-case analysis with prediction intervals which is crucial in high-stakes decision-making where the ground truth is not recoverable.\n\n**Q3. Sec 2.1 where do sets of datasets come from?**\n\nWe take an erroneous dataset as input, which might have missing values, detected outliers, or other data quality issues or sources of uncertainty. When these errors appear not at random, the ground truth values are typically not identifiable. Therefore, we consider a range of all possible clean data values in place of erroneous values. As mentioned in Appendix D, the ranges could be obtained from background knowledge, or reliable imputation techniques. To the best of our knowledge, there is no direct application of the Bayesian approach that derives the full distribution of models for capturing uncertainty due to general data multiplicity. We will clarify in the paper that our approach does not introduce any data quality issues, but rather enables modeling of data quality issues and evaluating their impact on model and prediction uncertainty. Regarding the difference with Bayesian approaches, please see Q7.\n\n**Q4. App. D missing definition of ε. Use maximum entropy distribution for missing data?**\n\nIn Appendix D, $\\epsilon$ refers to the uncertainty range radius. For instance, a temperature sensor may be accurate within 2 degrees. When encoding uncertainty, we do not make any assumptions about the distribution within the ranges. This allows for worst-case analysis even if probabilities are unknown. Intuitively, our over-approximation guarantees on prediction ranges and model parameters hold independent of how values are distributed within a range. We will add clarifications.\n\n**Q5. Terminology: prediction set, or \"efficiency ratio\" vs.\"robustness ratio\"?**\n\nWe agree that showing the distribution of prediction ranges will provide more fine-grained information. We will add additional plots to show distributions of prediction interval sizes. Regarding the nomenclature, while prediction set / efficiency is the common terminology used in conformal prediction (CP), robustness has been used in more closely related work [a]. We will add clarifications to relate our terminology to CP. Please also see our response to Q7.\n\n[a] Meyer et al. Certifying robustness to programmable data bias in decision trees. NeurIPS 2021\n\n**Q6. How are the robustness thresholds chosen?**\n\nPlease see our response in the main rebuttal section.\n\n**Q7. Experimental comparison with uncertainty quantification?**\n\nWhile many uncertainty quantification techniques output prediction intervals including Bayesian neural networks and conformal prediction (CP), we did not compare against these techniques as they are either not applicable to our setting (e.g., because their assumptions are violated or they are designed for a different type of uncertainty) and / or they provide only probabilistic guarantees, underestimating the true prediction ranges. For instance, CP’s assumption of exchangeability is violated when data errors are non-random and CP only guarantees that the ground truth label is within the prediction set in expectation and provides no guarantees for individual datapoints. Contrast this with our guarantees: we cover the prediction made by a model trained on the unknown ground truth training dataset with 100% probability. To also confirm this experimentally, we did evaluate prediction intervals produced by BNNs. We provide a detailed discussion in the general rebuttal. \n\n**Q8. For Figure 3: How were errors injected? How was multidimensional data generated?**\n\nIn Figure 3, we uniformly at random picked a set of 8 data points to inject errors by running several missing value imputation techniques to generate 4 possible options for each of the 8 uncertain values, which is small enough such that we could enumerate all possible worlds and calculate the ground truth prediction ranges. In Figure 3(a), we varied the uncertainty radius by creating ranges directly that over-approximate the possible worlds (4 possible imputed values). In Figure 3(b), we appended additional random columns with normal distributions to the original data, and used different imputation techniques to derive the uncertainty radiuses, i.e., ranges.\n\n**Q9. An extension of the proposed approach beyond linear models seems difficult**\n\nPlease see our response in the main rebuttal.", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "extra_scores": {}}]

Author rebuttal:
We thank the reviewers for their thoughtful reviews. We have responded to each reviewer in the individual responses and use the general response to present experimental results to justify some of our claims made in the responses and to respond to common concerns and concerns that require more extensive discussion.

**How were robustness thresholds chosen?**

Regarding the robustness threshold, the rationale for this choice was that different practical applications may differ in how much uncertainty they are willing to tolerate. To account for this in the experiments, we selected both a low threshold (0.5% for the insurance data) and a high threshold (5% for the MPG data). We also explored other thresholds which did not impact the trends significantly: Zorro can consistently certify a larger fraction of the test data points than the baseline due to its use of the more expressive Zonotope domain. We have included additional results in the main rebuttal pdf varying the robustness threshold (Figure 1), and we will add these plots to the appendix.

**Possibility to extend the approach beyond linear models in future work**

It is true that for linear models we benefit from the fact that we were able to derive a closed form solution which has the advantage of avoiding the repeated injection of an over-approximation that happens in each step of gradient descent. We are currently developing methods to (1) further improve the over-approximation quality by considering more complex abstract elements such as polynomial zonotopes, and (2) scale our approach to gradient descent for larger models, e.g., deep neural networks, which involves parallelization, optimized data structures, etc. While challenging, we believe that it will be worthwhile to explore these future research directions. 

**Experimental comparison against uncertainty quantification techniques like conformal prediction or BNNs**

There are many methods that produce prediction intervals, such as BNNs, Conformal prediction, bootstrap, and quantile regression. We did not compare experimentally against these approaches because they are generally not applicable to our setting or practical tools for producing reliable prediction intervals in the presence of data quality issues. 
- **Bayesian approaches**: In principle, we can use Bayesian approaches to generate competing models and then use them to derive viable prediction ranges. However, the intervals that we get will differ from those we obtain since Bayesian methods require (1) a correctly specified prior to provide any guarantees, and (2) a method to fit the optimal posterior distribution (which can typically only be approximated). In practice, applying this would fit competing models but typically these methods will not be able to consider all possible worlds and models. As a result, misspecification of prior and the approximation during training will lead to biased estimates and sampling will only give an under-estimation of the actual prediction ranges. 
- **Conformal Prediction (CP)**: Conceptually, Conformal Prediction (CP) is non-parametric and usually requires training a standalone model rather than a set of models. It is designed to capture a different source of uncertainty than what we consider here. The conformal guarantee ensures that the ground truth outcome is contained within the prediction interval probabilistically, meaning it covers the true outcome for test datapoints with a specified probability. In contrast, our approach provides guaranteed individual predictions that contain the ground truth prediction, i.e., the one obtained by a model trained on the ground truth training data. CP does not capture uncertainty due to missing data or data multiplicity, as this violates two key assumptions: (i) exchangeability and (ii) known ground truth labels for the validation set, both of which fail in our setting with prevalent data quality issues. While recent work [a-c] has investigated CP under uncertainty in test or training data, these approaches still rely on strong assumptions about data quality and build only a single model. An interesting future direction is to integrate our method with CP to establish prediction ranges that not only contain the ground truth prediction but also the ground truth outcome with rigorous theoretical guarantees.

We will add an extended discussion to the appendix to contrast our proposed techniques with other uncertainty quantification methods like CP, covering references [a-d]. Additionally, we ran empirical evaluations to demonstrate how data quality issues pose challenges for Bayesian approaches, making them inapplicable to our setting. Using the setting from the third plot in Figure 1(c) in the paper, where the uncertain data percentage is set to 10%, we tested Bayesian approaches on different possible worlds using two methods: impute-and-predict and sampling from possible worlds. The results show that the prediction intervals generated by Bayesian methods do not cover the ground truth prediction, i.e., the prediction by the model trained on the ground truth training data. In contrast, our approach guarantees 100% coverage across all cases (please see the tables in the attached PDF).

[a] Jeary et al. Verifiably Robust Conformal Prediction. arXiv 2024.

[b] Javanmardi et al. Conformal Prediction with Partially Labeled Data. Conformal and Probabilistic Prediction with Applications 2023

[c] Zaffran et al., Conformal Prediction with Missing Values. ICML 2023

[d] Barber et al. The limits of distribution-free conditional predictive inference. Information and Inference: A Journal of the IMA 2021.

[e] Meyer, et al. Certifying data-bias robustness in linear regression. FAccT 2023

**Additional experimental results**
- Table 1 and 2: Comparison between Bayesian methods and Zorro
- Figure 1: Effect of varying robustness ratio
- Figure 2: Robustness verification with for uncertainty in non-predictive features
- Figure 3: Runtime of Zorro

Meta-review: The submitted paper was reviewed by three knowledgeable reviewers who unanimously recommended its acceptance. While it is unclear how to scale the proposed approach to more complicated hypothesis spaces, the proposed approach is sensible and non-trivial and provides interesting and promising experimental results. In line with the reviewers, I am thus recommending acceptance of the paper. Nevertheless, I would encourage the authors to include clarifications along the lines of the discussions with the reviewers in the final paper, e.g., clarify differences to other "approaches" (e.g. Bayesian, conformal production) to help readers not that familiar with the topic understand the conceptual differences. Last but not least, the additional results provided during the discussion phase could make helpful additions to the appendix of the paper.
Decision: Accept (poster)