# Building a stable classifier with the inflated argmax

Jake A. Soloff  Department of Statistics, University of Chicago  Rina Foygel
Barber  Department of Statistics, University of Chicago  Rebecca Willett
Department of Statistics, University of Chicago  Department of Computer
Science, University of Chicago

###### Abstract

We propose a new framework for algorithmic stability in the context of
multiclass classification. In practice, classification algorithms often
operate by first assigning a continuous score (for instance, an estimated
probability) to each possible label, then taking the maximizer—i.e., selecting
the class that has the highest score. A drawback of this type of approach is
that it is inherently unstable, meaning that it is very sensitive to slight
perturbations of the training data, since taking the maximizer is
discontinuous. Motivated by this challenge, we propose a pipeline for
constructing stable classifiers from data, using bagging (i.e., resampling and
averaging) to produce stable continuous scores, and then using a stable
relaxation of argmax, which we call the “inflated argmax”, to convert these
scores to a set of candidate labels. The resulting stability guarantee places
no distributional assumptions on the data, does not depend on the number of
classes or dimensionality of the covariates, and holds for any base
classifier. Using a common benchmark data set, we demonstrate that the
inflated argmax provides necessary protection against unstable classifiers,
without loss of accuracy.

##  1 Introduction

An algorithm that learns from data is considered to be _stable_ if small
perturbations of the training data do not lead to large changes in its output.
Algorithmic stability is an important consideration in many statistical
applications. Within the fairness literature, for instance, stability is one
aspect of reliable decision-making systems (Friedler et al.,, 2019; Huang and
Vishnoi,, 2019). In interpretable machine learning, it similarly serves as a
form of reproducibility (Murdoch et al.,, 2019; Yu and Kumbier,, 2020; Yu and
Barter,, 2024). Caramanis et al., (2011) relates stability to robustness,
where a machine learning algorithm is robust if two samples with similar
feature vectors have similar test error. In the context of generalization
bounds, Bousquet and Elisseeff, (2002) and subsequent authors study stability
of an algorithm’s real-valued output—for instance, values of a regression
function. In the setting of a multiclass classification problem, where the
data consists of features $X\in\mathcal{X}$ and labels
$Y\in[L]=\\{1,\dots,L\\}$, the results of this literature can thus be applied
to analyzing a classifier’s predicted probabilities
$\hat{p}_{\ell}(X)\in[0,1]$—i.e., our learned estimates of the conditional
probabilities $\Pr\\{Y=\ell\mid X\\}$, for each $\ell\in[L]$. However, as we
will see, stability of these predicted probabilities by no means implies
stability of the predicted label itself,
$\hat{y}=\operatorname{argmax}_{\ell\in[L]}\hat{p}_{\ell}(x)$—an arbitrarily
small perturbation in $\hat{p}_{\ell}(x)$ might completely change the
predicted label $\hat{y}$. The distinction matters: for system
trustworthiness, we care about the model’s final output, on which users base
their decisions.

In this paper, we propose a new framework for algorithmic stability in the
context of multiclass classification, to define—and achieve—a meaningful
notion of stability when the output of the algorithm consists of predicted
labels, rather than predicted probabilities. Our work connects to other
approaches to allowing for ambiguous classification, including set-valued
classification (Grycko,, 1993; Del Coz et al.,, 2009; Sadinle et al.,, 2019;
Chzhen et al.,, 2021) and conformal inference (Papadopoulos et al.,, 2002;
Vovk et al.,, 2005; Lei,, 2014; Angelopoulos et al.,, 2023) (we will discuss
related work further, below).

###  1.1 Problem setting

In supervised classification, we take a data set
$\mathcal{D}=\big{(}(X_{i},Y_{i})\big{)}_{i\in[n]}$ consisting of $n$
observations, and train a classifier that maps from the feature space
$\mathcal{X}$ to the set $[L]$ of possible labels.111We remark that all the
results of this paper apply also to the case of countably infinitely many
labels, $L=\infty$, in which case we should take $\mathbb{N}$ to be the label
space instead of $[L]$. Typically, this map is constructed in two stages.
First we run some form of regression to learn a map
$\hat{p}:\mathcal{X}\to\Delta_{L-1}$ from features to predicted probabilities,
with $\hat{p}_{\ell}(x)$ denoting our learned estimate of the conditional
label probability, $\Pr\\{Y=\ell\mid X=x\\}$ (here
$\Delta_{L-1}=\\{w\in\mathbb{R}^{L}:w_{i}\geq 0,\sum_{i}w_{i}=1\\}$ denotes
the probability simplex in $L$ dimensions). We will write
$\hat{p}=\mathcal{A}(\mathcal{D})$, where $\mathcal{A}$ denotes the learning
algorithm mapping a data set $\mathcal{D}$ (of any size) to the map $\hat{p}$.
Then the second step is to convert the predicted probabilities
$\hat{p}_{\ell}(x)$ to a predicted label, which is typically done by taking
the argmax, i.e., $\hat{y}=\operatorname{argmax}_{\ell}\hat{p}_{\ell}(x)$
(with some mechanism for breaking ties). This two-stage procedure can be
represented in the following diagram:

Training data $\mathcal{D}=\big{(}(X_{i},Y_{i})\big{)}_{i\in[n]}$ Predicted
probabilities $\hat{p}(x)\in\Delta_{L-1}$ Test features $x\in\mathcal{X}$
Predicted label $\hat{y}\in[L]$ Learning
---
algorithm $\mathcal{A}$
Selection
---
(via $\operatorname{argmax}$)

When predictions are ambiguous, meaning two or more classes nearly achieve the
maximum predicted probability, the selected label becomes unstable and can
change based on seemingly inconsequential perturbations to the training data.
In other words, the above workflow is fundamentally incompatible with the goal
of algorithmic stability. Consequently, in this paper we instead work with
_set-valued classifiers_ , which return a set of candidate labels—in practice,
this typically leads to a singleton set for examples where we have high
confidence in the label, but allows for a larger set in the case of ambiguous
examples. While the idea of returning a set of candidate labels is not itself
new, we will see that the novelty of our work lies in finding a construction
that offers provable stability guarantees.

In this framework, a feature vector $x\in\mathcal{X}$ is now mapped to a set
of candidate labels $\hat{S}\subseteq[L]$ (rather than a single label
$\hat{y}$), as illustrated here:

Training data $\mathcal{D}=\big{(}(X_{i},Y_{i})\big{)}_{i\in[n]}$ Predicted
probabilities $\hat{p}(x)\in\Delta_{L-1}$ Test features $x\in\mathcal{X}$
Candidate labels $\hat{S}\subseteq[L]$ Learning
---
algorithm $\mathcal{A}$
Selection
---
rule $s$

Formally, given a test point $x\in\mathcal{X}$, this two-stage procedure
returns $\hat{S}=s(\hat{p}(x))\subseteq[L]$, where
$\hat{p}=\mathcal{A}(\mathcal{D})$ is the output of the regression algorithm
$\mathcal{A}$ trained on data $\mathcal{D}$. Here
$s:\Delta_{L-1}\to\raisebox{1.79993pt}{\Large$\wp$}([L])$ denotes a _selection
rule_ , mapping a vector of estimated probabilities to the set of candidate
labels. Of course, the earlier setting—where the procedure returns a single
label $\hat{y}=\operatorname{argmax}_{\ell}\,\hat{p}_{\ell}(x)$, rather than a
subset—can be viewed as a special case by simply taking $s$ to be the
$\operatorname{argmax}$ operator.

If we instead allow for $\hat{S}$ to contain multiple candidate labels, a
commonly used version of this framework is given by a top-5 procedure (or
top-$k$, for any fixed $k$). That is, after running a learning algorithm
$\mathcal{A}$ to estimate probabilities $\hat{p}(x)$, the selection rule $s$
then returns the top 5 labels—the labels $\ell_{1},\dots,\ell_{5}\in[L]$
corresponding to the highest 5 values of $\hat{p}_{\ell}(x)$. This approach is
more stable than a standard argmax. However, the choice of 5 in this setting
is somewhat arbitrary, and the set $\hat{S}$ always contains 5 candidate
labels, regardless of the difficulty of the test sample—while intuitively, we
might expect that it should be possible to return a smaller $\hat{S}$ for test
points $x$ that are easier to classify (and a larger $\hat{S}$ if $x$ is more
ambiguous). In contrast, in our work we seek a more principled approach, where
we provide provable stability guarantees while also aiming for the set
$\hat{S}$ to be as small as possible.

###  1.2 Contributions

In this work, we introduce new definitions of algorithmic stability in the
context of classification, focusing on the stability of the predicted label.
We reduce the problem of stabilizing a classifier to separately stabilizing
the learning and selection stages, described above. For the selection step, we
propose the _inflated argmax_ operation, derive its stabilizing properties and
give an algorithm to compute it efficiently. We combine this operation with
bagging at the learning step, to provably stabilize any classifier. In
particular, our guarantee holds with no assumptions on the data, and no
constraints on the dimensionality of the covariates nor on the number of
classes.

##  2 Framework: stable classification

In this section, we propose a definition of algorithmic stability in the
setting of multiclass classification. To begin, we formally define a
_classification algorithm_ as a map222We will assume without comment that
$\mathcal{C}$ is measurable. Furthermore, in practice classification
algorithms often incorporate randomization (either in the learning stage, such
as via stochastic gradient descent, and/or in the selection stage, such as by
using a random tie-breaking rule for $\operatorname{argmax}$). All of our
definitions and results in this paper can be naturally extended to randomized
algorithms—see Appendix B.

$\mathcal{C}:\cup_{n\geq 0}(\mathcal{X}\times[L])^{n}\ \times\ \mathcal{X}\
\longrightarrow\ \raisebox{1.79993pt}{\Large$\wp$}([L]),$

which maps a training data set $\mathcal{D}$ of any size $n$, together with a
test feature vector $x\in\mathcal{X}$, to a candidate set of labels
$\hat{S}=\mathcal{C}(\mathcal{D},x)$. To relate this notation to our earlier
terminology, the two-stage selection procedure described in Section 1.1
corresponds to the classification algorithm

$\mathcal{C}(\mathcal{D},x)=s(\hat{p}(x))\textnormal{ where
}\hat{p}=\mathcal{A}(\mathcal{D}).$

Abusing notation, we will write this as $\mathcal{C}=s\circ\mathcal{A}$,
indicating that $\mathcal{C}$ is obtained by applying selection rule $s$ to
the output of the learning algorithm $\mathcal{A}$.

We now present our definition of algorithmic stability for a classification
algorithm $\mathcal{C}$. As is common in the algorithmic stability literature,
we focus on stability of the algorithm’s output with respect to randomly
dropping a data point: does the output of $\mathcal{C}$ on a test point $x$
change substantially if we drop a single data point from the training data
$\mathcal{D}$? If the output of the algorithm were a real-valued prediction
$\hat{y}\in\mathbb{R}$, we could assess this by measuring the real-valued
change in $\hat{y}$ when a single data point is dropped. For classification,
however, we will need to take a different approach:

######  Definition 1.

We say a classification algorithm $\mathcal{C}$ has selection stability
$\delta$ at sample size $n$ if, for all datasets
$\mathcal{D}\in(\mathcal{X}\times[L])^{n}$ and all test features
$x\in\mathcal{X}$,

$\frac{1}{n}\sum_{i=1}^{n}\mathbf{1}\left\\{\hat{S}\cap\hat{S}^{\setminus
i}=\varnothing\right\\}\leq\delta,$

where $\hat{S}=\mathcal{C}(\mathcal{D},x)$ and where $\hat{S}^{\setminus
i}=\mathcal{C}(\mathcal{D}^{\setminus i},x)$, for each $i\in[n]$.

Here the notation $\mathcal{D}^{\setminus i}$ denotes the data set
$\mathcal{D}$ with $i$th data point removed—that is, for a data set
$\mathcal{D}=\big{(}(X_{j},Y_{j})\big{)}_{j\in[n]}$, the leave-one-out data
set is given by $\mathcal{D}^{\setminus
i}=\big{(}(X_{j},Y_{j})\big{)}_{j\in[n]\setminus i}$.

###  2.1 Interpreting the definition

Intuitively, selection stability controls the probability that our algorithm
makes wholly inconsistent claims when dropping a single data point at random
from the training set. If we interpret $\hat{S}$ as making the claim that the
true label $Y$ is equal to one of the labels in the candidate set $\hat{S}$,
and similarly for $\hat{S}^{\setminus i}$, then as soon as there is _even one
single value_ that lies in both $\hat{S}$ and $\hat{S}^{\setminus i}$, this
means that the two claims are not contradictory—even if the sets are large and
are mostly nonoverlapping.

At first glance, this condition appears to be too mild, in the sense that we
require the two prediction sets only to have _some_ way of being compatible,
and allows for substantial difference between the sets $\hat{S}$ and
$\hat{S}^{\setminus i}$. However, since standard classification algorithms
always output a single label, they often cannot be said to be stable even in
this basic sense. Thus, we can view this definition as providing a minimal
notion of stability that we should require any interpretable method to
satisfy.

###  2.2 Connection to classical algorithmic stability

Most prior work on algorithmic stability concerns algorithms with continuous
outputs—for example, in our notation above, the algorithm $\mathcal{A}$ that
returns estimated probabilities $\hat{p}(x)$. With such algorithms, there are
more standard tools at our disposal for quantifying and ensuring stability. In
this section, we connect classical algorithmic stability to our notion of
selection stability (Definition 1). We first recall the following definition.

######  Definition 2.

(Soloff et al., 2024b, ) A learning algorithm $\mathcal{A}$ has tail stability
$(\varepsilon,\delta)$ at sample size $n$ if, for all datasets $\mathcal{D}$
of size $n$ and all test features $x\in\mathcal{X}$,

$\frac{1}{n}\sum_{i=1}^{n}\mathbf{1}\left\\{\|\hat{p}(x)-\hat{p}^{\setminus
i}(x)\|\geq\varepsilon\right\\}\leq\delta,$

where $\hat{p}=\mathcal{A}(\mathcal{D})$ and $\hat{p}^{\setminus
i}=\mathcal{A}(\mathcal{D}^{\setminus i})$, and where $\|\cdot\|$ denotes the
usual Euclidean norm on $\mathbb{R}^{L}$.

This intuitive notion of stability is achieved by many well known algorithms
$\mathcal{A}$, such as nearest-neighbor type methods, or methods based on
bagging or ensembling (as established by Soloff et al., 2024a ; Soloff et al.,
2024b —see Section 3.1 below for details).

###  2.3 From stability to selection stability

To construct a classification algorithm using the two-stage procedure outlined
in Section 1.1, we need to apply a selection rule $s$ to the output of our
learning algorithm $\mathcal{A}$. We might expect that choosing a stable
$\mathcal{A}$ will lead to selection stability in the resulting classification
algorithm—but in fact, this is not the case: even if the learning algorithm
$\mathcal{A}$ is itself extremely stable in the sense of Definition 2, the
classification rule obtained by applying $\operatorname{argmax}$ to the output
of $\mathcal{A}$ can still be extremely unstable. The underlying issue is that
$\operatorname{argmax}$ is extremely discontinuous—the perturbation
$\|\hat{p}(x)-\hat{p}^{\setminus i}(x)\|$ in the predicted probabilities can
be arbitrarily small but still lead to different predicted labels, i.e.,
$\operatorname{argmax}_{\ell}\hat{p}_{\ell}(x)\neq\operatorname{argmax}_{\ell}\hat{p}^{\setminus
i}_{\ell}(x)$.

Since combining $\operatorname{argmax}$ with a stable learning algorithm
$\mathcal{A}$ will not suffice, we instead seek a different selection rule
$s$—one that will ensure selection stability (when paired with a stable
learning algorithm $\mathcal{A}$). To formalize this aim, we introduce another
definition:

######  Definition 3.

A selection rule $s:\Delta_{L-1}\to\raisebox{1.79993pt}{\Large$\wp$}([L])$ is
$\varepsilon$-compatible if, for any $v,w\in\Delta_{L-1}$,

$\|v-w\|<\varepsilon\ \Longrightarrow\ s(v)\cap s(w)\neq\varnothing.$

This notion of $\varepsilon$-compatibility allows us to construct
classification algorithms with selection stability. Combining the above
definitions leads immediately to the following result:

######  Proposition 4.

Let $\mathcal{A}$ be a learning algorithm with tail stability
$(\varepsilon,\delta)$ at sample size $n$, and let $s$ be a selection rule
satisfying $\varepsilon$-compatibility. Then the classification algorithm
$\mathcal{C}=s\circ\mathcal{A}$ has selection stability $\delta$ at sample
size $n$.

Therefore, by pairing a stable learning algorithm $\mathcal{A}$ with a
compatible selection rule $s$, we will automatically ensure selection
stability of the resulting classification algorithm
$\mathcal{C}=s\circ\mathcal{A}$.

Of course, $\varepsilon$-compatibility of the selection rule $s$ might be
achieved in a trivial way—for instance, $s$ returns the entire set $[L]$ for
any input. As is common in statistical settings (e.g., a tradeoff between Type
I error and power, in hypothesis testing problems), our goal is to ensure
selection stability while returning an output $\hat{S}$ that is as informative
as possible. In particular, later on we will consider the specific aim of
constructing $\hat{S}$ to be a singleton set as often as possible.

##  3 Methodology: assumption-free stable classification

In this section, we formally define our pipeline for building a stable
classification procedure using any base learning algorithm $\mathcal{A}$. At a
high level, we leverage 4 and separately address the learning and selection
stages described in Section 2.

  1. 1.

In Section 3.1, we construct a bagged (i.e., ensembled) version of the base
learning algorithm $\mathcal{A}$. The recent work of Soloff et al., 2024b
ensures that the resulting bagged algorithm has tail stability
$(\varepsilon,\delta)$, with $\varepsilon\asymp\frac{1}{\sqrt{n\delta}}$.

  2. 2.

In Section 3.2, we introduce a new selection rule, the inflated argmax, and
establish its $\varepsilon$-compatibility. Combined with the bagged algorithm,
then, the resulting classification algorithm will be guaranteed to have
selection stability $\delta$ (by 4).

Before describing these two steps in detail, we first present our main theorem
that characterizes the selection stability guarantee of the resulting
procedure. Given sample size $n$ for the training data set $\mathcal{D}$, the
notation $\widetilde{A}_{m}$ will denote a bagged version of the base
algorithm $\mathcal{A}$, obtained by averaging over subsamples of
$\mathcal{D}$ comprised of $m$ data points sampled either with replacement
(“bootstrapping”, commonly run with $m=n$) or without replacement
(“subbagging”, commonly run with $m=n/2$)—see below for details.

######  Theorem 5.

Fix any sample size $n$, any bag size $m$, and any inflation parameter
$\varepsilon>0$. For any base learning algorithm $\mathcal{A}$, the
classification algorithm
$\mathcal{C}=\operatorname{argmax}^{\varepsilon}\circ\widetilde{A}_{m}$,
obtained by combining the bagged version of $\mathcal{A}$ together with the
inflated argmax, satisfies selection stability $\delta$ where

$\delta=\varepsilon^{-2}\cdot\frac{1-1/L}{n-1}\cdot\frac{p_{n,m}}{1-p_{n,m}},$
(1)

where $p_{n,m}=1-(1-\frac{1}{n})^{m}$ for bootstrapping, and
$p_{n,m}=\frac{m}{n}$ for subbagging.

The guarantee in Theorem 5 holds for any base learning algorithm
$\mathcal{A}$, and applies regardless of the complexity of the feature space
$\mathcal{X}$, and allows the test feature $x\in\mathcal{X}$ to be chosen
adversarially. The dependence on the number of classes $L$ is mild—in fact,
the tail stability parameter $\delta$ in (1) differs only by a factor of two
for $L=2$ versus $L=\infty$. Of course, the guarantee does depend on the
choice of the bag size $m$. In general a smaller value of $m$ leads to a
stronger stability guarantee (since $p_{n,m}$ increases with $m$), but this
comes at the expense of accuracy since we are training the base algorithm
$\mathcal{A}$ on subsampled data sets $\mathcal{D}^{r}$ of size $m$ (rather
than $n$). For the common choices of $m=n$ for bootstrap or $m=n/2$ for
subbagging, we have $\frac{p_{n,m}}{1-p_{n,m}}=\mathcal{O}(1)$ for each.

###  3.1 Bagging classifier weights

In this section, we formally define the construction of a bagged algorithm. We
consider the two most common variants of this ensembling method: bagging
(based on bootstrapping training points) and subbagging (based on
subsampling). For any data set $\mathcal{D}\in(\mathcal{X}\times[L])^{n}$, we
can define a subsampled data set of size $m$ as follows: for a sequence
$r=(i_{1},\dots,i_{m})\in[n]^{m}$ (which is called a _bag_), we define the
corresponding subsampled data set
$\mathcal{D}^{r}=\big{(}(X_{i},Y_{i})\big{)}_{i\in
r}\in(\mathcal{X}\times[L])^{m}$. Note that if the bag $r$ contains repeated
indices (i.e., $i_{k}=i_{\ell}$ for some $k\neq\ell$), then the same data
point from the original data set $\mathcal{D}$ will appear multiple times in
$\mathcal{D}^{r}$.

######  Definition 6 (Bootstrapping or subbagging a base learning algorithm
$\mathcal{A}$).

Given a data set $\mathcal{D}\in(\mathcal{X}\times[L])^{n}$ and some
$x\in\mathcal{X}$, return the output
$\widetilde{\mathcal{A}}_{m}(\mathcal{D})(x)=\mathbb{E}_{r}[\mathcal{A}(\mathcal{D}^{r})(x)]$,
where the expected value is taken with respect to a bag bag $r$ that is
sampled as follows:

  * •

Bootstrapping (sometimes simply referred to as bagging) (Breiman, 1996a, ;
Breiman, 1996b, ) constructs each bag $r$ by sampling $m$ indices
$r=(i_{1},\dots,i_{m})$ uniformly with replacement from $[n]$.

  * •

Subbagging (Andonova et al.,, 2002) constructs each bag $r$ by sampling $m\leq
n$ indices $r=(i_{1},\dots,i_{m})$ uniformly without replacement from $[n]$.

The following result (Soloff et al., 2024b, ) ensures tail stability for any
bootstrapped or subbagged algorithm:

######  Theorem 7 ((Soloff et al., 2024b, )).

For any base learning algorithm $\mathcal{A}$ returning outputs in
$\Delta_{L-1}$, the bagged algorithm $\widetilde{\mathcal{A}}_{m}$ has tail
stability $(\varepsilon,\delta)$ for any $\varepsilon,\delta>0$ satisfying
(1).

##### Computational cost of bagging.

In this section, we have worked with the ideal, derandomized version of
bagging for simplicity—that is, we assume that the expected value
$\mathbb{E}_{r}[\mathcal{A}(\mathcal{D}^{r})(x)]$ is calculated exactly with
respect to the distribution of $r$. In practice, of course, this is
computationally infeasible (since for bootstrapping, say, there are $n^{m}$
possible bags $r$), and so we typically resort to Monte Carlo sampling to
approximate this expected value, defining
$\widetilde{\mathcal{A}}_{m}(\mathcal{D})(x)=\frac{1}{B}\sum_{b=1}^{B}\mathcal{A}(\mathcal{D}^{r_{b}})(x)$,
for $B$ i.i.d. bags $r_{1},\dots,r_{B}$ sampled via bootstrapping or
subbagging. See Appendix C for extensions of our main result, Theorem 5, to
the finite-$B$ version of the method.

###  3.2 The inflated argmax

We are now ready to define the inflated argmax, denoted by
$\operatorname{argmax}^{\varepsilon}$, which will be used as our selection
rule $s$ for building our classification algorithm. First, for each $j\in[L]$,
define

$R_{j}^{\varepsilon}=\left\\{w\in\mathbb{R}^{L}:w_{j}\geq\max_{\ell\neq
j}w_{\ell}+\varepsilon/\sqrt{2}\right\\},$

which is a closed and convex subset of $\mathbb{R}^{L}$. For any
$w\in\mathbb{R}^{L}$, we then define

$\operatorname{argmax}^{\varepsilon}(w):=\left\\{j\in[L]:\textnormal{dist}(w,R_{j}^{\varepsilon})<\varepsilon\right\\},$
(2)

where $\textnormal{dist}(w,R_{j}^{\varepsilon})=\inf_{v\in
R_{j}^{\varepsilon}}\|w-v\|$. Intuitively,
$j\in\operatorname{argmax}^{\varepsilon}(w)$ indicates that a small
perturbation of $w$ would lead to the $j$th entry being largest by some
margin.

While $\operatorname{argmax}^{\varepsilon}$ is defined as an operator on
$\mathbb{R}^{L}$, in the context of classification algorithms, we only apply
it to vectors $w$ lying in the simplex $\Delta_{L-1}$ (in particular, in the
context of a two-stage classification procedure as in Section 1.1, we will
apply the inflated argmax to the vector $w=\hat{p}(x)$ of estimated
probabilities). In Figure 1, we visualize the inflated argmax applied to the
simplex in a setting with three possible labels, $L=3$. Each different shaded
area corresponds to the region of the simplex $\Delta_{L-1}$ where
$\operatorname{argmax}^{\varepsilon}(w)$ returns a particular subset.

Figure 1: The left plot illustrates the inflated argmax (2) over the simplex
$\Delta_{L-1}$ when $L=3$. The numbers in brackets corresponds to the output
of the inflated argmax, $\operatorname{argmax}^{\varepsilon}(w)$, for various
points $w$ in the simplex. The right plot shows the same but for the standard
argmax, which corresponds to the limit of
$\operatorname{argmax}^{\varepsilon}(w)$ as $\varepsilon\rightarrow 0$.

####  3.2.1 Inflated argmax leads to selection stability

Unlike the standard argmax, the inflated argmax allows us to achieve stability
in the context of classification. This is due to the following theorem, which
shows that the inflated argmax is $\varepsilon$-compatible, as introduced in
Definition 3. (The proofs of this result and all properties of the inflated
argmax are deferred to Appendix A.)

######  Theorem 8.

For any $\varepsilon>0$ and any $v,w\in\mathbb{R}^{L}$, if
$\|v-w\|<\varepsilon$ then
$\operatorname{argmax}^{\varepsilon}(v)\cap\operatorname{argmax}^{\varepsilon}(w)\neq\varnothing$.
In particular, viewing $\operatorname{argmax}^{\varepsilon}$ as a map on the
simplex $\Delta_{L-1}$, $\operatorname{argmax}^{\varepsilon}$ is
$\varepsilon$-compatible.

To gain a more concrete understanding of this theorem, we can look again at
Figure 1 to examine the case $L=3$. Theorem 8 ensures that any two regions
corresponding to disjoint subsets—e.g., the top corner region corresonding to
$\\{2\\}$, and the bottom center region corresponding to $\\{1,3\\}$—are
distance at least $\varepsilon$ apart. In particular, the region in the
center, corresponding to vectors $w$ that map to the full set $\\{1,2,3\\}$,
is a curved triangle of constant width $\varepsilon$, also known as a
_Reuleaux triangle_ (Reuleaux,, 1876; Weisstein,, 2005).

##### Proof of Theorem 5.

With the above results in place, we can now see that the inflated argmax
allows us to achieve selection stability, when combined with a bagged
algorithm. In particular, combining Theorem 7, Theorem 8, and 4 immediately
implies our main result, Theorem 5.

####  3.2.2 Additional properties of the inflated argmax

The following result establishes some natural properties obeyed by the
inflated argmax, and also compares to the standard argmax.

######  Proposition 9 (Basic properties of the inflated argmax).

Fix any $\varepsilon>0$. The inflated argmax operator satisfies the following
properties:

  1. 1.

(Including the argmax.) For any $w\in\mathbb{R}^{L}$,
$\operatorname{argmax}(w)\subseteq\operatorname{argmax}^{\varepsilon}(w)$.333In
this result, $\operatorname{argmax}(w)$ should be interpreted as the set of
_all_ maximizing entries of $w$—i.e., in the case of ties,
$\operatorname{argmax}(w)$ may not be a singleton set. Moreover,
$\cap_{\varepsilon>0}\operatorname{argmax}^{\varepsilon}(w)=\operatorname{argmax}(w)$.

  2. 2.

(Monotonicity in $\varepsilon$.) For any $w\in\mathbb{R}^{L}$ and any
$\varepsilon<\varepsilon^{\prime}$,
$\operatorname{argmax}^{\varepsilon}(w)\subseteq\operatorname{argmax}^{\varepsilon^{\prime}}(w)$.

  3. 3.

(Monotonicity in $w$.) For any $w\in\mathbb{R}^{L}$, if $w_{j}\leq w_{k}$,
then $j\in\operatorname{argmax}^{\varepsilon}(w)$ $\Rightarrow$
$k\in\operatorname{argmax}^{\varepsilon}(w)$.

  4. 4.

(Permutation invariance.) For any $v,w\in\mathbb{R}^{L}$, if
$v=(w_{\sigma(1)},\dots,w_{\sigma(L)})$ for some permutation $\sigma$ on
$[L]$, then $j\in\operatorname{argmax}^{\varepsilon}(v)$ $\Leftrightarrow$
$\sigma(j)\in\operatorname{argmax}^{\varepsilon}(w)$.

Next, while we have established that inflated argmax offers favorable
stability properties, we have not yet asked whether it can be efficiently
computed—in particular, it is not immediately clear how to verify the
condition $\textnormal{dist}(w,R_{j}^{\varepsilon})<\varepsilon$ in (2), in
order to determine whether $j\in\operatorname{argmax}^{\varepsilon}(w)$. The
following result offers an efficient algorithm for computing the inflated
argmax set.

######  Proposition 10 (Computing the inflated argmax).

Fix any $w\in\mathbb{R}^{L}$ and $\varepsilon>0$. Let $w_{[1]}\geq\dots\geq
w_{[L]}$ denote the order statistics of $w$, and define

$\hat{k}(w)=\max\bigg{\\{}k\in[L]:\Big{(}\sum_{\ell=1}^{k}(w_{[\ell]}-w_{[k]})\Big{)}^{2}+\sum_{\ell=1}^{k}(w_{[\ell]}-w_{[k]})^{2}\leq\varepsilon^{2}\bigg{\\}}.$

Then

$\operatorname{argmax}^{\varepsilon}(w)=\left\\{j\in[L]:w_{j}>\frac{\varepsilon}{\sqrt{2}}+\hat{A}_{1}(w)-\sqrt{\hat{k}(w)+1}\sqrt{\frac{\varepsilon^{2}}{\hat{k}(w)}+(\hat{A}_{1}(w))^{2}-\hat{A}_{2}(w)}\right\\},$

where $\hat{A}_{1}(w)=\frac{w_{[1]}+\dots+w_{[\hat{k}(w)]}}{\hat{k}(w)}$, and
$\hat{A}_{2}(w)=\frac{w_{[1]}^{2}+\dots+w_{[\hat{k}(w)]}^{2}}{\hat{k}(w)}$.

Finally, thus far we have established that
$\operatorname{argmax}^{\varepsilon}$ enables us to achieve stable
classification with a computationally efficient selection rule—but we do not
yet know whether $\operatorname{argmax}^{\varepsilon}$ is optimal, or whether
some other choice of $s$ might lead to smaller output sets $\hat{S}$ (while
still offering assumption-free stability). For instance, we might consider a
fixed-margin rule,

$s_{\textnormal{margin}}^{\varepsilon}(w)=\\{j:w_{j}>\max_{\ell}w_{\ell}-\varepsilon/\sqrt{2}\\},$
(3)

for which $\varepsilon$-compatibility also holds—might this simpler rule be
better than the inflated argmax? Our next result establishes that—under some
natural constraints—$\operatorname{argmax}^{\varepsilon}$ is in fact the
optimal choice of selection rule, in the sense of returning a singleton set
(i.e., $|\hat{S}|=1$) as often as possible.

######  Proposition 11 (Optimality of the inflated argmax).

Let $s:\Delta_{L-1}\to\raisebox{1.79993pt}{\Large$\wp$}([L])$ be any selection
rule. Suppose $s$ is $\varepsilon$-compatible (Definition 3), permutation
invariant (in the sense of 9), and contains the argmax. Then for any
$w\in\Delta_{L-1}$ and any $j\in[L]$,

$s(w)=\\{j\\}\ \Longrightarrow\
\operatorname{argmax}^{\varepsilon}(w)=\\{j\\}.$

In other words, for any selection rule $s$ satisfying the assumptions of the
proposition (which includes the fixed-margin rule,
$s_{\textnormal{margin}}^{\varepsilon}$), if $s(w)$ is a singleton set then so
is $\operatorname{argmax}^{\varepsilon}(w)$. (See also Appendix D for a closer
comparison between the inflated argmax and the fixed-margin selection rule
given in (3).)

##  4 Experiments

Figure 2: Results on the Fashion MNIST data set. The figure shows the
instability $\delta_{j}$ (defined in (4)) over all test points $j=1,\dots,N$.
The curves display the fraction of $\delta_{j}$’s that exceed $\delta$, for
each value $\delta\in[0,1]$. The vertical axis is on a log scale. See Section
4 for details.

In this section, we evaluate our proposed pipeline, combining subbagging with
the inflated argmax, with deep learning models and on a common benchmark data
set.444Code to fully reproduce the experiment is available at
https://github.com/jake-soloff/stable-argmax-experiments. Training all of the
models for this experiment took a total of four hours on 10 CPUs running in
parallel on a single computing cluster.

##### Data and models.

We use Fashion-MNIST (Xiao et al.,, 2017), which consists of $n=60,000$
training pairs $(X_{i},Y_{i})$, $N=10,000$ test pairs
$(\tilde{X}_{j},\tilde{Y}_{j})$, and $L=10$ classes. For each data point
$(X,Y)$, $X$ is a $28\times 28$ grayscale image that pictures a clothing item,
and $Y\in[L]$ indicates the type of item, e.g., a dress, a coat, etc. The base
model we use is a variant of LeNet-5, implemented in PyTorch (Paszke et al.,,
2019) tutorials as GarmentClassifier(). The base algorithm $\mathcal{A}$
trains this classifier using $5$ epochs of stochastic gradient descent.

##### Methods and evaluation.

We compare three methods:

  1. 1.

Taking the argmax of the base learning algorithm $\mathcal{A}$.

  2. 2.

Taking the argmax of the subbagged algorithm $\widetilde{\mathcal{A}}_{m}$,
with $B=200$ bags of size $m=n/2$.

  3. 3.

Taking the $\varepsilon$-inflated argmax of the subbagged algorithm
$\widetilde{\mathcal{A}}_{m}$, with $B=200$ bags of size $m=n/2$ and tolerance
$\varepsilon=.05$.

We evaluate each method based on several metrics. First, to assess selection
stability, we compute, for each test point $j=1,\ldots,N$,

$\delta_{j}:=\frac{1}{500}\sum_{k=1}^{500}\mathbf{1}\left\\{s(\hat{p}(\tilde{X}_{j}))\cap
s(\hat{p}^{\setminus i_{k}}(\tilde{X}_{j}))\neq\varnothing\right\\},$ (4)

where $i_{1},\dots,i_{500}$ are sampled uniformly without replacement from
$[n]$ (i.e., we are sampling 500 leave-one-out models $\hat{p}^{\setminus i}$
to compare to the model $\hat{p}$ trained on the full training set). Second,
to evaluate how accurate each method is, we compute average precision (how
often the method unambiguously returns the correct label) and set size (the
number of labels in the candidate set):

$\beta_{\textnormal{prec}}=\frac{1}{N}\sum_{j=1}^{N}\mathbf{1}\left\\{s(\hat{p}(\tilde{X}_{j}))=\\{\widetilde{Y}_{j}\\}\right\\},\quad\beta_{\textnormal{size}}=\frac{1}{N}\sum_{j=1}^{N}\big{|}s(\hat{p}(\tilde{X}_{j}))\big{|}.$

Ideally we would want a method to achieve high precision (a value of
$\beta_{\textnormal{prec}}\in[0,1]$ that is close to $1$) and low set size (a
value of $\beta_{\textnormal{size}}\geq 1$ that is close to $1$).

##### Results.

In Figure 2 and Table 1, we present results for the instability and accuracy
of each method. In the figure, we show the survival function of
$(\delta_{j})_{j\in[N]}$. The base algorithm,
$\operatorname{argmax}\circ\mathcal{A}$, has the longest tail, i.e.,
$\delta_{j}$ is large for many test points $j$. The other two methods,
$\operatorname{argmax}\circ\widetilde{\mathcal{A}}_{m}$ and
$\operatorname{argmax}^{\varepsilon}\circ\widetilde{\mathcal{A}}_{m}$, both
improve upon this tail, but the latter method shows the best stability by far:
using the inflated argmax ensures that even the largest $\delta_{j}$ is small,
with $\max_{j}\delta_{j}\approx 0.14$. In the table, we display the average
accuracy measures, $\beta_{\textnormal{prec}}$ and
$\beta_{\textnormal{size}}$. We see that precision is nearly identical for all
three methods. Moreover, while the first two methods always return a set of
size 1 (since we are using the standard argmax), the inflated argmax increases
set size only very slightly in this experiment. In other words, the inflated
argmax allows us to achieve substantially higher stability without sacrificing
accuracy in this high signal data set.

| $\beta_{\textnormal{prec}}$ | $\beta_{\textnormal{size}}$
---|---|---
$\operatorname{argmax}\circ\mathcal{A}$ | 0.879 (0.003) | 1.000 (0.000)
$\operatorname{argmax}\circ\widetilde{\mathcal{A}}_{m}$ | 0.893 (0.003) | 1.000 (0.000)
$\operatorname{argmax}^{\varepsilon}\circ\widetilde{\mathcal{A}}_{m}$ | 0.886 (0.003) | 1.018 (0.001)
Table 1: Results on the Fashion MNIST data set. The table displays the average
precision, $\beta_{\textnormal{prec}}$, and size, $\beta_{\textnormal{size}}$,
as defined in Section 4. Standard errors are shown in parentheses.

##  5 Discussion

##### Related work.

There is an enormous literature on set-valued classification, so we only
reference some of the most relevant works. Much prior work has considered the
possibility of a more relaxed argmax to allow for uncertainty (Bridle,, 1990;
de Brébisson and Vincent,, 2016; Martins and Astudillo,, 2016; Blondel et
al.,, 2020; Zhao et al.,, 2023). The more recent literature has focused on
producing a sparse set of weights, but none of these works offer a formal
stability guarantee for the support of the weights. Our work is the first to
propose a version of the argmax that can provably control a notion of
stability for the classification setting.

Our definition of selection stability relies on a notion of consistency
between sets—two sets of candidate labels are consistent if they have at least
one common element. This is similar in flavor to conformal classification
(Lei,, 2014; Sadinle et al.,, 2019), where a set of candidate values is valid
if it contains the correct label; this method does not rely on any
distributional assumptions, and consequently has been applied to complex and
high-dimensional applications such as generative language models (Quach et
al.,, 2023). These frameworks share the motivation of avoiding ‘overconfidence
in the assignment of a definite label’ (Hechtlinger et al.,, 2018).

##### Overview, limitations and open questions.

We prove a guarantee on the selection stability based on our methodology
combining bagging with the inflated argmax. Theorem 5 does not place any
assumptions on the learning algorithm nor on the data, including any
distributional assumptions. In fact, the training data and test point may be
chosen adversarially based on the algorithm, and the output will still be
stable. Moreover, we do not assume that the sample size is large: the
guarantee holds for any fixed training set size $n$ and improves as $n$
increases. Furthermore, the inflated argmax selection rule ensures that the
returned sets of candidate labels are as small as possible (i.e., that there
is as little ambiguity as possible about the predicted label for any given
test point $x$).

While our theorem does not require assumptions, our method does require
bagging the base learning algorithm. The main limitation of our work is that
bagging is computationally intensive. However, training different bags is
easily parallelizable, which is what allowed us to easily train a
convolutional neural network on $1,000$ total subsets of the training data.

We leave open several avenues for future work. Practitioners may be interested
in other forms of discrete stability, which are more strict than requiring one
common element between $\hat{S}$ and $\hat{S}^{\setminus i}$. For instance,
one popular way to measure set similarity is the Jaccard index, given by
$\frac{|\hat{S}\cap\hat{S}^{\setminus i}|}{|\hat{S}\cup\hat{S}^{\setminus
i}|}$. Another open problem is to extend the stability framework and methods
studied here to more structured discrete outputs, such as in classification.

## Acknowledgements

We gratefully acknowledge the support of the National Science Foundation via
grant DMS-2023109 DMS-2023109. R.F.B. and J.A.S. gratefully acknowledge the
support of the Office of Naval Research via grant N00014-20-1-2337. R.M.W.
gratefully acknowledges the support of NSF DMS-2235451, Simons Foundation MP-
TMPS-00005320, and AFOSR FA9550-18-1-0166. The authors thank Melissa Adrian
for help with running the experiments.

##  Appendix A Proofs of properties of the inflated argmax

In this section, we prove all the results of Section 3.2, establishing the
various properties of the inflated argmax. In fact, all the proofs will rely
on the following theorem, which gives an alternative characterization of the
set $\operatorname{argmax}^{\varepsilon}(w)$.

######  Theorem 12.

Fix $\varepsilon>0$. Define the map
$c_{\varepsilon}:\mathbb{R}^{L}\rightarrow\mathbb{R}$ as

$c_{\varepsilon}(w)\textnormal{ is the unique solution to
}\bigg{(}\sum_{j\in[L]}(w_{j}-c)_{+}\bigg{)}^{2}+\sum_{j\in[L]}(w_{j}-c)^{2}_{+}=\varepsilon^{2},$

and define the map $t_{\varepsilon}:\mathbb{R}^{L}\rightarrow\mathbb{R}$ as

$t_{\varepsilon}(w)=c_{\varepsilon}(w)+\varepsilon/\sqrt{2}-\sum_{j\in[L]}(w_{j}-c_{\varepsilon}(w))_{+}.$

Then for any $w\in\mathbb{R}^{d}$, it holds that

$\operatorname{argmax}^{\varepsilon}(w)=\\{j\in[L]:w_{j}>t_{\varepsilon}(w)\\}.$

Moreover, $t_{\varepsilon}(w)\leq c_{\varepsilon}(w)$ for all $w$ and all
$\varepsilon$, and $t_{\varepsilon}(w)$ and $c_{\varepsilon}(w)$ are both
nondecreasing functions of $\varepsilon$, and both nondecreasing functions of
$w$ (i.e., if $w\leq w^{\prime}$ holds coordinatewise, then
$c_{\varepsilon}(w)\leq c_{\varepsilon}(w^{\prime})$ and
$t_{\varepsilon}(w)\leq t_{\varepsilon}(w^{\prime})$).

With this key result in place, we are now ready to turn to the proofs of the
results stated in Section 3.2. Throughout all the proofs below, the functions
$c_{\varepsilon}(w)$ and $t_{\varepsilon}(w)$ are defined as in the statement
of Theorem 12.

###  A.1 Proof of Theorem 8

Let $w,w^{\prime}\in\mathbb{R}^{L}$ be any vectors with
$\operatorname{argmax}^{\varepsilon}(w)\cap\operatorname{argmax}^{\varepsilon}(w^{\prime})=\varnothing$.
Define

$B=\\{j\in[L]:w_{j}>c_{\varepsilon}(w)\\},\quad
B^{\prime}=\\{j\in[L]:w^{\prime}_{j}>c_{\varepsilon}(w^{\prime})\\},$

and note that we must have $B\subseteq\operatorname{argmax}^{\varepsilon}(w)$
since for $j\not\in\operatorname{argmax}^{\varepsilon}(w)$, $w_{j}\leq
t_{\varepsilon}(w)\leq c_{\varepsilon}(w)$ by Theorem 12. Similarly,
$B^{\prime}\subseteq\operatorname{argmax}^{\varepsilon}(w^{\prime})$. We also
must have $B$ (and similarly $B^{\prime}$) nonempty by definition of
$c_{\varepsilon}(w)$.

For convenience, define $c=c_{\varepsilon}(w)$,
$c^{\prime}=c_{\varepsilon}(w^{\prime})$, $t=t_{\varepsilon}(w)$,
$t^{\prime}=t_{\varepsilon}(w^{\prime})$. We then have

$\|w-w^{\prime}\|^{2}\geq\|w_{B}-w^{\prime}_{B}\|^{2}+||w^{\prime}_{B^{\prime}}-w_{B^{\prime}}\|^{2}\\\
=\|(w_{B}-c\mathbf{1}_{B})+(t^{\prime}\mathbf{1}_{B}-w^{\prime}_{B})+(c-t^{\prime})\mathbf{1}_{B}\|^{2}+\|(w^{\prime}_{B^{\prime}}-c^{\prime}\mathbf{1}_{B^{\prime}})+(t\mathbf{1}_{B^{\prime}}-w_{B^{\prime}})+(c^{\prime}-t)\mathbf{1}_{B^{\prime}}\|^{2},$
(5)

where the first step holds since $B\cap B^{\prime}=\varnothing$ due to
$\operatorname{argmax}^{\varepsilon}(w)\cap\operatorname{argmax}^{\varepsilon}(w^{\prime})=\varnothing$.
Next we will need a technical lemma.

######  Lemma 13.

Let $k,k^{\prime}\geq 1$ be integers and let $\varepsilon>0$. Then for any
$x,y\in\mathbb{R}^{k}_{\geq 0}$,
$x^{\prime},y^{\prime}\in\mathbb{R}^{k^{\prime}}_{\geq 0}$, and
$r,r^{\prime}\in\mathbb{R}$, if

$(x^{\top}\mathbf{1}_{k})^{2}+\|x\|^{2}=(x^{\prime}{}^{\top}\mathbf{1}_{k^{\prime}})^{2}+\|x^{\prime}\|^{2}=\varepsilon^{2}$
(6)

and

$x^{\top}\mathbf{1}_{k}+x^{\prime}{}^{\top}\mathbf{1}_{k^{\prime}}=r+r^{\prime}+\varepsilon\sqrt{2},$
(7)

then

$\|x+y+r\mathbf{1}_{k}\|^{2}+\|x^{\prime}+y^{\prime}+r^{\prime}\mathbf{1}_{k^{\prime}}\|^{2}\geq\varepsilon^{2}.$

To apply this lemma, first we let $k=|B|\geq 1$ and
$k^{\prime}=|B^{\prime}|\geq 1$, and define

$x=w_{B}-c\mathbf{1}_{B},\ y=t^{\prime}\mathbf{1}_{B}-w^{\prime}_{B},\
x^{\prime}=w^{\prime}_{B^{\prime}}-c^{\prime}\mathbf{1}_{B^{\prime}},\
y^{\prime}=t\mathbf{1}_{B^{\prime}}-w_{B^{\prime}},$

and note that all these vectors have nonnegative entries (specifically, $x$
and $x^{\prime}$ are nonnegative by definition of $B$ and $B^{\prime}$, while
$y$ and $y^{\prime}$ are nonnegative by Theorem 12, using the fact that $B\cap
s(w^{\prime})=B^{\prime}\cap s(w)=\varnothing$). Then the condition (6) is
satisfied by definition of $c=c_{\varepsilon}(w)$ and
$c^{\prime}=c_{\varepsilon}(w^{\prime})$. Define also

$r=c-t^{\prime},\ r^{\prime}=c^{\prime}-t.$

and note that

$r+r^{\prime}=(c-t)+(c^{\prime}-t^{\prime})=\sum_{j\in[L]}(w_{j}-c)_{+}-\varepsilon/\sqrt{2}+\sum_{j\in[L]}(w^{\prime}_{j}-c^{\prime})_{+}-\varepsilon/\sqrt{2}\\\
=x^{\top}\mathbf{1}_{k}+x^{\prime}{}^{\top}\mathbf{1}_{k^{\prime}}-\varepsilon\sqrt{2},$

where the second step holds by definition of $t=t_{\varepsilon}(w)$ and
$t^{\prime}=t_{\varepsilon}(w^{\prime})$. Therefore, (7) is also satisfied.
Returning to (5), Lemma 13 then implies that
$\|w-w^{\prime}\|^{2}\geq\varepsilon^{2}$, which completes the proof of the
theorem.

###  A.2 Proof of 9

First we verify that
$\operatorname{argmax}(w)\subseteq\operatorname{argmax}^{\varepsilon}(w)$. Let
$j\in\operatorname{argmax}(w)$, i.e., $w_{j}=\max_{\ell}w_{\ell}$. Let
$v=w+\mathbf{e}_{j}\cdot\varepsilon/\sqrt{2}$, where
$\mathbf{e}_{j}=(0,\dots,0,1,0,\dots,0)$ is the $j$th canonical basis vector.
Then $v\in R_{j}^{\varepsilon}$, and
$\|w-v\|=\varepsilon/\sqrt{2}<\varepsilon$, so we have
$\textnormal{dist}(w,R_{j}^{\varepsilon})<\varepsilon$ and therefore
$j\in\operatorname{argmax}^{\varepsilon}(w)$. We also need to verify that
$\cap_{\varepsilon>0}\operatorname{argmax}^{\varepsilon}(w)=\operatorname{argmax}(w)$,
i.e., for $j\not\in\operatorname{argmax}(w)$, there is some $\varepsilon>0$
with $j\not\in\operatorname{argmax}^{\varepsilon}(w)$. Let $k\in[L]$ be an
index with $w_{k}>w_{j}$, and let $\varepsilon=\sqrt{2}(w_{k}-w_{j})$. Then
for any $v\in R_{j}^{\varepsilon}$, we have

$\|w-v\|^{2}\geq(w_{j}-v_{j})^{2}+(w_{k}-v_{k})^{2}=(w_{k}-\varepsilon/\sqrt{2}-v_{j})^{2}+(w_{k}-v_{k})^{2}\\\
\geq\inf_{t\in\mathbb{R}}\left\\{\left(t-(\varepsilon/\sqrt{2}+(v_{j}-v_{k}))\right)^{2}+t^{2}\right\\}=\frac{1}{2}(\varepsilon/\sqrt{2}+(v_{j}-v_{k}))^{2}\geq\varepsilon^{2},$

where the last step holds since $v_{j}-v_{k}\geq\varepsilon/\sqrt{2}$. This
means that $\textnormal{dist}(w,R_{j}^{\varepsilon})\geq\varepsilon$ and so
$j\not\in\operatorname{argmax}^{\varepsilon}(w)$ at this value of
$\varepsilon>0$.

Next we check monotonicity in $\varepsilon$. Fix
$\varepsilon<\varepsilon^{\prime}$. Then for any $j\in[L]$,

$j\in\operatorname{argmax}^{\varepsilon}(\varepsilon)\ \Longleftrightarrow\
w_{j}>t_{\varepsilon}(w)\ \Longrightarrow\ w_{j}>t_{\varepsilon^{\prime}}(w)\
\Longleftrightarrow\ j\in\operatorname{argmax}^{\varepsilon^{\prime}}(w),$

where each step holds by Theorem 12.

Next we verify monotonicity in $w$. If $w_{j}\leq w_{k}$, then by Theorem 12

$j\in\operatorname{argmax}^{\varepsilon}(w)\ \Longleftrightarrow
w_{j}>t_{\varepsilon}(w)\ \Longrightarrow w_{k}>t_{\varepsilon}(w)\
\Longleftrightarrow\ k\in\operatorname{argmax}^{\varepsilon}(w).$

Finally we check permutation invariance. Suppose
$v=(w_{\sigma(1)},\dots,w_{\sigma(L)})$ is a permutation of $w$. Then by
construction, we have $c_{\varepsilon}(v)=c_{\varepsilon}(w)$ and
$t_{\varepsilon}(v)=t_{\varepsilon}(w)$. In particular,

$j\in\operatorname{argmax}^{\varepsilon}(v)\ \Longleftrightarrow\
v_{j}>t_{\varepsilon}(v)\ \Longleftrightarrow\
w_{\sigma(j)}>t_{\varepsilon}(w)\ \Longleftrightarrow\
\sigma(j)\in\operatorname{argmax}^{\varepsilon}(w).$

###  A.3 Proof of 10

Define a function $f_{w}$ as

$f_{w}(c)=\left(\sum_{j\in[L]}(w_{j}-c)_{+}\right)^{2}+\sum_{j\in[L]}(w_{j}-c)^{2}_{+}.$
(8)

Then $\hat{k}(w)$ can equivalently be defined as

$\hat{k}(w)=\max\\{k\in[L]:f_{w}(w_{[k]})\leq\varepsilon^{2}\\}.$

(Note that this maximum is well defined, since $f_{w}(w_{[1]})=0$ and so the
set is nonempty.) Since $c\mapsto f_{w}(c)$ is strictly decreasing over $c\leq
w_{[1]}$, we see that the solution $c_{\varepsilon}(w)$ to the equation
$f_{w}(c)=\varepsilon^{2}$ must satisfy $c_{\varepsilon}(w)\leq
w_{[\hat{k}(w)]}$, and (if $\hat{k}(w)<L$) also
$c_{\varepsilon}(w)>w_{[\hat{k}(w)+1]}$. In particular, this implies

$\displaystyle\varepsilon^{2}$
$\displaystyle=\left(\sum_{j\in[L]}(w_{j}-c_{\varepsilon}(w))_{+}\right)^{2}+\sum_{j\in[L]}(w_{j}-c_{\varepsilon}(w))_{+}^{2}$
$\displaystyle=\left(\sum_{j=1}^{\hat{k}(w)}(w_{j}-c_{\varepsilon}(w))\right)^{2}+\sum_{j=1}^{\hat{k}(w)}(w_{j}-c_{\varepsilon}(w))^{2}$
$\displaystyle=\left(\hat{k}(w)\cdot(\hat{A}_{1}(w)-c_{\varepsilon}(w))\right)^{2}+\left(\hat{k}(w)\hat{A}_{2}(w)-2\hat{k}(w)\hat{A}_{1}(w)c_{\varepsilon}(w)+\hat{k}(w)c_{\varepsilon}(w)^{2}\right)$
$\displaystyle=\hat{k}(w)(\hat{k}(w)+1)\left(c_{\varepsilon}(w)^{2}-2c_{\varepsilon}(w)\hat{A}_{1}(w)+\frac{\hat{k}(w)(\hat{A}_{1}(w))^{2}+\hat{A}_{2}(w)}{\hat{k}(w)+1}\right).$

This is a quadratic function of $c_{\varepsilon}(w)$, and is solved by

$c_{\varepsilon}(w)=\hat{A}_{1}-\sqrt{\frac{(\hat{A}_{1}(w))^{2}-\hat{A}_{2}(w)+\frac{\varepsilon^{2}}{\hat{k}(w)}}{\hat{k}(w)+1}}.$

We also have

$\displaystyle t_{\varepsilon}(w)$
$\displaystyle=c_{\varepsilon}(w)+\varepsilon/\sqrt{2}-\sum_{j\in[L]}(w_{j}-c_{\varepsilon}(w))_{+}$
$\displaystyle=c_{\varepsilon}(w)+\varepsilon/\sqrt{2}-\sum_{j=1}^{\hat{k}(w)}(w_{[j]}-c_{\varepsilon}(w))$
$\displaystyle=(\hat{k}(w)+1)c_{\varepsilon}(w)+\frac{\varepsilon}{\sqrt{2}}-\hat{k}(w)\hat{A}_{1}(w).$

Plugging in our above expression for $c_{\varepsilon}(w)$, then,

$t_{\varepsilon}(w)=\frac{\varepsilon}{\sqrt{2}}+\hat{A}_{1}(w)-\sqrt{\hat{k}(w)+1}\sqrt{(\hat{A}_{1}(w))^{2}-\hat{A}_{2}(w)+\frac{\varepsilon^{2}}{\hat{k}(w)}}.$

Since
$\operatorname{argmax}^{\varepsilon}(w)=\\{j:w_{j}>t_{\varepsilon}(w)\\}$ by
Theorem 12, this completes the proof.

###  A.4 Proof of 11

First we will need a lemma:

######  Lemma 14.

Fix $w\in\mathbb{R}^{L}$. Then $w\in R_{j}^{\varepsilon}$ if and only if
$\operatorname{argmax}^{\varepsilon}(w)=\\{j\\}$.

For intuition on this result, we can look back at Figure 1—for instance, the
blue region marked by $\\{1\\}$ illustrates the set of vectors
$w\in\Delta_{L-1}$ such that $\operatorname{argmax}^{\varepsilon}(w)=\\{1\\}$,
and according to this lemma, this is equal to the region
$R_{1}^{\varepsilon}$.

Using this lemma, we now need to show that, for any $\varepsilon$-compatible
and permutation invariant selection rule $s$, if $s(w)=\\{j\\}$ for some
$w\in\Delta_{L-1}$ and $j\in[L]$, then $w\in R_{j}^{\varepsilon}$ for some
$j$. First, we must have $j=\operatorname{argmax}(w)$, since $s$ is assumed to
contain the argmax. Fix any $k\neq j$, and define a vector $v\in\Delta_{L-1}$
by permuting $w_{j}$ and $w_{k}$:

$v_{\ell}=\begin{cases}w_{k},&\ell=j,\\\ w_{j},&\ell=k,\\\ w_{\ell},&\ell\neq
j,k.\end{cases}.$

Then by permutation invariance we have $s(v)=\\{k\\}$. Since $s$ is
$\varepsilon$-compatible, we therefore have

$\varepsilon\leq\|w-v\|=\sqrt{(w_{j}-v_{j})^{2}+(w_{k}-v_{k})^{2}+\sum_{\ell\neq
j,k}(w_{\ell}-v_{\ell})^{2}}=\sqrt{2(w_{j}-w_{k})^{2}}.$

Therefore we have $w_{j}-w_{k}\geq\varepsilon/\sqrt{2}$ for all $k\neq j$,
which proves $w\in R_{j}^{\varepsilon}$ and therefore
$\operatorname{argmax}^{\varepsilon}(w)=\\{j\\}$.

###  A.5 Proof of Theorem 12

##### Step 1: verifying that $c_{\varepsilon}(w)$ is well defined.

First we check that the solution $c_{\varepsilon}(w)$ exists and is unique.
Fix $w$ and define the function $f_{w}$ as in (8). Note that $f_{w}(c)\equiv
0$ for $c\geq\max_{i}w_{i}$, and $c\mapsto f_{w}(c)$ is strictly decreasing
over $c\leq\max_{i}w_{i}$, with $\lim_{c\rightarrow-\infty}f_{w}(c)=\infty$;
therefore, a unique solution $c_{\varepsilon}(w)$ to the equation
$f_{w}(c)=\varepsilon^{2}$ must exist.

##### Step 2: verifying $t_{\varepsilon}(w)\leq c_{\varepsilon}(w)$.

Next we check that the claim $t_{\varepsilon}(w)\leq c_{\varepsilon}(w)$ must
hold. First, by definition of $c_{\varepsilon}(w)$, we have

$\varepsilon^{2}=\left(\sum_{j\in[L]}(w_{j}-c_{\varepsilon}(w))_{+}\right)^{2}+\sum_{j\in[L]}(w_{j}-c_{\varepsilon}(w))^{2}_{+}\leq
2\left(\sum_{j\in[L]}(w_{j}-c_{\varepsilon}(w))_{+}\right)^{2},$

where the inequality holds by the properties of the $\ell_{1}$ and $\ell_{2}$
norms. Therefore,

$\sum_{j\in[L]}(w_{j}-c)_{+}\geq\varepsilon/\sqrt{2},$

which verifies that

$t_{\varepsilon}(w)=c_{\varepsilon}(w)+\varepsilon/\sqrt{2}-\sum_{j\in[L]}(w_{j}-c_{\varepsilon}(w))_{+}\leq
c_{\varepsilon}(w).$

##### Step 3: checking monotonicity.

Now we turn to verifying the monotonicity properties of $t_{\varepsilon}(w)$
and $c_{\varepsilon}(w)$. First we check that these functions are
nonincreasing in $\varepsilon$. First, since $f_{w}(c)$ is nonincreasing in
$c$, and $c_{\varepsilon}(w)$ is the solution to $f_{w}(c)=\varepsilon^{2}$,
this immediately implies that $c_{\varepsilon}(w)$ is nonincreasing in
$\varepsilon$. Next we turn to $t_{\varepsilon}(w)$. Define

$t^{\prime}(c)=c+\frac{1}{\sqrt{2}}\sqrt{\left(\sum_{j\in[L]}(w_{j}-c)_{+}\right)^{2}+\sum_{j\in[L]}(w_{j}-c)_{+}^{2}}-\sum_{j\in[L]}(w_{j}-c)_{+},$

so that we have $t_{\varepsilon}(w)=t^{\prime}(c_{\varepsilon}(w))$ by
construction. We can verify that $c\mapsto t^{\prime}(c)$ is nondecreasing,
and therefore, $t_{\varepsilon}(w)=t^{\prime}(c_{\varepsilon}(w))\leq
t^{\prime}(c_{\varepsilon^{\prime}}(w))=t_{\varepsilon^{\prime}}(w)$, where
the inequality holds since $c_{\varepsilon}(w)\leq
c_{\varepsilon^{\prime}}(w)$.

Next we check that $t_{\varepsilon}(w)$ and $c_{\varepsilon}(w)$ are
nondecreasing functions of $w$. Fix any $w\leq w^{\prime}$ (where the
inequality is coordinatewise, i.e., $w_{j}\leq w^{\prime}_{j}$ for all
$j\in[L]$). Since $w\mapsto f_{w}(c)$ is a nondecreasing function, we have
$f_{w}(c)\leq f_{w^{\prime}}(c)$ for all $c$. We therefore have
$f_{w^{\prime}}(c_{\varepsilon}(w))\geq
f_{w}(c_{\varepsilon}(w))=\varepsilon^{2}=f_{w^{\prime}}(c_{\varepsilon}(w^{\prime}))$.
Since $c\mapsto f_{w^{\prime}}(c)$ is nonincreasing, therefore,
$c_{\varepsilon}(w^{\prime})\geq c_{\varepsilon}(w)$. Next we consider
$t_{\varepsilon}$. Let
$t^{\prime\prime}_{\varepsilon}(c)=c+\varepsilon/\sqrt{2}-\sum_{j\in[L]}(w_{j}-c)_{+}$,
which is a nondecreasing function of $c$. Then we have
$t_{\varepsilon}(w)=t^{\prime\prime}_{\varepsilon}(c_{\varepsilon}(w))\leq
t^{\prime\prime}_{\varepsilon}(c_{\varepsilon}(w^{\prime}))=t_{\varepsilon}(w^{\prime})$.

##### Step 4: returning to the inflated argmax.

Finally, fixing any $w\in\mathbb{R}^{L}$, we will prove that
$\operatorname{argmax}^{\varepsilon}(w)=\\{j:w_{j}>t_{\varepsilon}(w)\\}$.
First choose any $j\in[L]$ with $w_{j}>t_{\varepsilon}(w)$. We need to verify
that $j\in\operatorname{argmax}^{\varepsilon}(w)$. Define $v\in\mathbb{R}^{L}$
with entries

$v_{k}=\begin{cases}\max\\{c_{\varepsilon}(w)+\varepsilon/\sqrt{2},w_{j}\\},&k=j,\\\
\min\\{c_{\varepsilon}(w),w_{k}\\},&k\neq j.\end{cases}$

By construction, we have $v\in R_{j}^{\varepsilon}$. We calculate

$\displaystyle\textnormal{dist}(w,R_{j}^{\varepsilon})^{2}$
$\displaystyle\leq\|w-v\|^{2}$
$\displaystyle=\left(w_{j}-\max\\{c_{\varepsilon}(w)+\varepsilon/\sqrt{2},w_{j}\\}\right)^{2}+\sum_{k\neq
j}\left(w_{k}-\min\\{c_{\varepsilon}(w),w_{k}\\}\right)^{2}$
$\displaystyle=\left(c_{\varepsilon}(w)+\varepsilon/\sqrt{2}-w_{j}\right)^{2}_{+}+\sum_{k\neq
j}(w_{k}-c_{\varepsilon}(w))^{2}_{+}$
$\displaystyle<\left(c_{\varepsilon}(w)+\varepsilon/\sqrt{2}-t_{\varepsilon}(w)\right)^{2}+\sum_{k\neq
j}(w_{k}-c_{\varepsilon}(w))^{2}_{+}$
$\displaystyle\leq\left(\sum_{k\in[L]}(w_{k}-c_{\varepsilon}(w))_{+}\right)^{2}+\sum_{k\in[L]}(w_{k}-c_{\varepsilon}(w))^{2}_{+}=\varepsilon^{2},$

where the last two steps hold by definition of $t_{\varepsilon}(w)$ and
$c_{\varepsilon}(w)$, while the strict inequality holds since
$t_{\varepsilon}(w)\leq
c_{\varepsilon}(w)<c_{\varepsilon}(w)+\varepsilon/\sqrt{2}$ as established
above, while $w_{j}>t_{\varepsilon}(w)$ by assumption. Therefore
$j\in\operatorname{argmax}^{\varepsilon}(w)$.

Now we check the converse. For this last step, we will need a lemma. The
following result characterizes the projection of $w$ to the set
$R_{j}^{\varepsilon}$.

######  Lemma 15.

Fix any $w\in\mathbb{R}^{L}$ and any $j\in[L]$. Then there is a unique
$a\in\mathbb{R}$ satisfying

$w_{j}=a+\varepsilon/\sqrt{2}-\sum_{k\neq j}(w_{k}-a)_{+}.$ (9)

Moreover, defining the vector $v\in\mathbb{R}^{L}$ as

$v_{k}=\begin{cases}a+\varepsilon/\sqrt{2},&k=j,\\\ a\wedge w_{k},&k\neq
j,\end{cases}$ (10)

it holds that $v=\operatorname{argmin}_{u\in R_{j}^{\varepsilon}}\|w-u\|$,
i.e., $v$ is the projection of $w$ to the set $R_{j}^{\varepsilon}$.

Next suppose $w_{j}\leq t_{\varepsilon}(w)$. We need to verify that
$j\not\in\operatorname{argmax}^{\varepsilon}(w)$. Let $a$ and $v$ be defined
as in Lemma 15 above. We can compare the equation (9) to the definition of
$t_{\varepsilon}(w)$,

$t_{\varepsilon}(w)=c_{\varepsilon}(w)+\varepsilon/\sqrt{2}-\sum_{k\in[L]}(w_{k}-c_{\varepsilon}(w))_{+}=c_{\varepsilon}(w)+\varepsilon/\sqrt{2}-\sum_{k\neq
j}(w_{k}-c_{\varepsilon}(w))_{+},$

where the last step holds since $w_{j}\leq t_{\varepsilon}(w)$ by assumption,
and $t_{\varepsilon}(w)\leq c_{\varepsilon}(w)$ as proved above. Since
$c\mapsto c+\varepsilon/\sqrt{2}-\sum_{k\neq j}(w_{k}-c)_{+}$ is an increasing
function, and $w_{j}\leq t_{\varepsilon}(w)$, this implies $a\leq
c_{\varepsilon}(w)$. We then calculate

$\displaystyle\textnormal{dist}(w,R_{j}^{\varepsilon})^{2}$
$\displaystyle=\|w-v\|^{2}$
$\displaystyle=\left(w_{j}-\left(a+\varepsilon/\sqrt{2}\right)\right)^{2}+\sum_{k\neq
j}\left(w_{k}-a\wedge w_{k}\right)^{2}\textnormal{ by~{}\eqref{eqn:v_a}}$
$\displaystyle=\left(w_{j}-\left(a+\varepsilon/\sqrt{2}\right)\right)^{2}+\sum_{k\neq
j}(w_{k}-a)^{2}_{+}$ $\displaystyle=\left(\sum_{k\neq
j}(w_{k}-a)_{+}\right)^{2}+\sum_{k\neq j}(w_{k}-a)^{2}_{+}\textnormal{
by~{}\eqref{eqn:find_a}}$ $\displaystyle\geq\left(\sum_{k\neq
j}(w_{k}-c_{\varepsilon}(w))_{+}\right)^{2}+\sum_{k\neq
j}(w_{k}-c_{\varepsilon}(w))^{2}_{+}\textnormal{ since $a\leq
c_{\varepsilon}(w)$}$
$\displaystyle=\left(\sum_{k\in[L]}(w_{k}-c_{\varepsilon}(w))_{+}\right)^{2}+\sum_{k\in[L]}(w_{k}-c_{\varepsilon}(w))^{2}_{+}\textnormal{
since $w_{j}\leq t_{\varepsilon}(w)\leq c_{\varepsilon}(w)$}$
$\displaystyle=\varepsilon^{2}\textnormal{ by definition of
$c_{\varepsilon}(w)$.}$

Therefore, $\textnormal{dist}(w,R_{j}^{\varepsilon})\geq\varepsilon$, and so
$j\not\in\operatorname{argmax}^{\varepsilon}(w)$.

###  A.6 Proofs of technical lemmas

####  A.6.1 Proof of Lemma 13

First, since $y$ is constrained to have nonnegative entries,

$\|x+y+r\mathbf{1}_{k}\|^{2}\geq\|(x+r\mathbf{1}_{k})_{+}\|^{2},$

where for a vector $v=(v_{1},\dots,v_{L})\in\mathbb{R}^{L}$, we write
$(v)_{+}$ to denote the vector with $j$th entry given by
$(v_{j})_{+}=\max\\{v_{j},0\\}$ for each $j$. The analogous bound holds for
$\|x^{\prime}+y^{\prime}+r^{\prime}\mathbf{1}_{k^{\prime}}\|^{2}$.

$\|x+y+r\mathbf{1}_{k}\|^{2}+\|x^{\prime}+y^{\prime}+r^{\prime}\mathbf{1}_{k^{\prime}}\|^{2}\geq\|(x+r\mathbf{1}_{k})_{+}\|^{2}+\|(x^{\prime}+r^{\prime}\mathbf{1}_{k^{\prime}})_{+}\|^{2}.$

We now need to show that

$\|(x+r\mathbf{1}_{k})_{+}\|^{2}+\|(x^{\prime}+r^{\prime}\mathbf{1}_{k^{\prime}})_{+}\|^{2}\geq\varepsilon^{2}.$

Next let

$\bar{r}=\frac{r+r^{\prime}}{2}=\frac{x^{\top}\mathbf{1}_{k}+x^{\prime}{}^{\top}\mathbf{1}_{k^{\prime}}-\varepsilon\sqrt{2}}{2},\quad\Delta=\frac{-r+r^{\prime}}{2},$

so that we can write

$r=\bar{r}-\Delta,\ r^{\prime}=\bar{r}+\Delta$

for some $\Delta\in\mathbb{R}$. We we therefore need to show
$\inf_{\Delta\in\mathbb{R}}f(\Delta)\geq\varepsilon^{2}$, where

$f(\Delta)=\|(x+(\bar{r}-\Delta)\mathbf{1}_{k})_{+}\|^{2}+\|(x^{\prime}+(\bar{r}+\Delta)\mathbf{1}_{k^{\prime}})_{+}\|^{2}.$

First, we observe that we can restrict the range of $\Delta$. Specifically,
for any $\Delta\geq\max_{j\in[L]}x_{j}+\bar{r}$, we have
$f(\Delta)=\|(x^{\prime}+(\bar{r}+\Delta)\mathbf{1}_{k^{\prime}})_{+}\|^{2}$,
which is a nondecreasing function; therefore,

$\inf_{\Delta\geq\max_{j\in[L]}x_{j}+\bar{r}}f(\Delta)=f\left(\max_{j\in[L]}x_{j}+\bar{r}\right),$

meaning that we do not need to consider values of $\Delta$ beyond this upper
bound. Applying a similar argument for a lower bound, we see that from this
point on we only need to verify that

$f(\Delta)\geq\varepsilon^{2}\textnormal{ \ for
$-\max_{j\in[k^{\prime}]}x^{\prime}_{j}-\bar{r}\leq\Delta\leq\max_{j\in[L]}x_{j}+\bar{r}$}.$

Moreover, for any $\Delta$, we have

$(x_{j}+(\bar{r}-\Delta))^{2}=(x_{j}+(\bar{r}-\Delta))^{2}_{+}+(x_{j}+(\bar{r}-\Delta))^{2}_{-}\leq(x_{j}+(\bar{r}-\Delta))^{2}_{+}+(\bar{r}-\Delta)^{2}$

for all $j$, by nonnegativity of $x$. Furthermore we must have

$(x_{j}+(\bar{r}-\Delta))^{2}=(x_{j}+(\bar{r}-\Delta))^{2}_{+}$

for at least one $j\in[L]$ when $\Delta\leq\max_{j\in[L]}x_{j}+\bar{r}$ as
specified above (i.e., because for $j$ maximizing the entry $x_{j}$, the value
$(x+(\bar{r}-\Delta)\mathbf{1}_{k})_{j}$ is nonnegative). Therefore,

$\displaystyle\|(x+(\bar{r}-\Delta)\mathbf{1}_{k})_{+}\|^{2}$
$\displaystyle\geq\|x+(\bar{r}-\Delta)\mathbf{1}_{k}\|^{2}-(k-1)(\bar{r}-\Delta)^{2}.$

Similarly, we can show that

$\|(x^{\prime}+(\bar{r}+\Delta)\mathbf{1}_{k^{\prime}})_{+}\|^{2}\geq\|x^{\prime}+(\bar{r}+\Delta)\mathbf{1}_{k^{\prime}}\|^{2}-(k^{\prime}-1)(\bar{r}+\Delta)^{2}.$

We then calculate

$\displaystyle f(\Delta)$
$\displaystyle=\|(x+(\bar{r}-\Delta)\mathbf{1}_{k})_{+}\|^{2}+\|(x^{\prime}+(\bar{r}+\Delta)\mathbf{1}_{k^{\prime}})_{+}\|^{2}$
$\displaystyle\geq\|x+(\bar{r}-\Delta)\mathbf{1}_{k}\|^{2}-(k-1)(\bar{r}-\Delta)^{2}+\|x^{\prime}+(\bar{r}+\Delta)\mathbf{1}_{k^{\prime}}\|^{2}-(k^{\prime}-1)(\bar{r}+\Delta)^{2}$
$\displaystyle=\|x\|^{2}+2(\bar{r}-\Delta)x^{\top}\mathbf{1}_{k}+(\bar{r}-\Delta)^{2}+\|x^{\prime}\|^{2}+2(\bar{r}+\Delta)x^{\prime}{}^{\top}\mathbf{1}_{k^{\prime}}+(\bar{r}+\Delta)^{2}$
$\displaystyle=2\varepsilon^{2}-(x^{\top}\mathbf{1}_{k})^{2}-(x^{\prime}{}^{\top}\mathbf{1}_{k^{\prime}})^{2}+2(\bar{r}-\Delta)x^{\top}\mathbf{1}_{k}+2(\bar{r}+\Delta)x^{\prime}{}^{\top}\mathbf{1}_{k^{\prime}}+2\bar{r}^{2}+2\Delta^{2},$

where the last step holds by (6).

Writing $z=x^{\top}\mathbf{1}_{k}$ and
$z^{\prime}=x^{\prime}{}^{\top}\mathbf{1}_{k^{\prime}}$ for convenience, we
have

$\displaystyle f(\Delta)$ $\displaystyle\geq
2\varepsilon^{2}-z^{2}-z^{\prime}{}^{2}+2(\bar{r}-\Delta)z+2(\bar{r}+\Delta)z^{\prime}+2\bar{r}^{2}+2\Delta^{2}$
$\displaystyle=2\varepsilon^{2}-z^{2}-z^{\prime}{}^{2}+2\left(\frac{z+z^{\prime}-\varepsilon\sqrt{2}}{2}-\Delta\right)z+2\left(\frac{z+z^{\prime}-\varepsilon\sqrt{2}}{2}+\Delta\right)z^{\prime}$
$\displaystyle\hskip
72.26999pt+2\left(\frac{z+z^{\prime}-\varepsilon\sqrt{2}}{2}\right)^{2}+2\Delta^{2}$
$\displaystyle=\varepsilon^{2}+4\left(z-\varepsilon/\sqrt{2}\right)\left(z^{\prime}-\varepsilon/\sqrt{2}\right)+2\left(\Delta-\frac{z-z^{\prime}}{2}\right)^{2}$
$\displaystyle\geq\varepsilon^{2}+4\left(z-\varepsilon/\sqrt{2}\right)\left(z^{\prime}-\varepsilon/\sqrt{2}\right).$

where the second step holds because
$\bar{r}=\frac{r+r^{\prime}}{2}=\frac{z+z^{\prime}-\varepsilon\sqrt{2}}{2}$ by
(7).

Finally, we also have $x^{\top}\mathbf{1}_{k}=\|x\|_{1}\geq\|x\|$, where the
first step holds since $x$ is nonnegative, and so by assumption (6) we have
$z=x^{\top}\mathbf{1}_{k}\geq\varepsilon/\sqrt{2}$. Similarly,
$z^{\prime}\geq\varepsilon/\sqrt{2}$. Therefore,
$4\left(z-\varepsilon/\sqrt{2}\right)\left(z^{\prime}-\varepsilon/\sqrt{2}\right)\geq
0$, and so we have shown that $f(\Delta)\geq\varepsilon^{2}$, as desired.

####  A.6.2 Proof of Lemma 14

First, suppose $w\in R_{j}^{\varepsilon}$. Then
$\operatorname{argmax}(w)=\\{j\\}$, and so
$j\in\operatorname{argmax}^{\varepsilon}(w)$ by 9. Next, let
$c=w_{j}-\varepsilon/\sqrt{2}$. Then $w_{\ell}\leq c$ for all $\ell\neq j$,
and so defining $f_{w}$ as in (8), we have

$f_{w}(c)=(w_{j}-c)^{2}+(w_{j}-c)^{2}=\varepsilon^{2}\ \Longrightarrow\
c_{\varepsilon}(w)=c=w_{j}-\varepsilon/\sqrt{2}.$

Then

$t_{\varepsilon}(w)=c_{\varepsilon}(w)+\varepsilon/\sqrt{2}-\sum_{\ell\in[L]}(w_{\ell}-c_{\varepsilon}(w))_{+}=\left(w_{j}-\varepsilon/\sqrt{2}\right)+\varepsilon/\sqrt{2}-\varepsilon/\sqrt{2}=w_{j}-\varepsilon/\sqrt{2},$

and so for all $\ell\neq j$, we have $w_{\ell}\leq
w_{j}-\varepsilon/\sqrt{2}=t_{\varepsilon}(w)$ and thus
$\ell\not\in\operatorname{argmax}^{\varepsilon}(w)$. This verifies that
$\operatorname{argmax}^{\varepsilon}(w)=\\{j\\}$.

To prove the converse, suppose
$\operatorname{argmax}^{\varepsilon}(w)=\\{j\\}$. By 9, we then have
$\operatorname{argmax}(w)=\\{j\\}$. Let $k\in\operatorname{argmax}_{\ell\neq
j}w_{\ell}$, then to show $w\in R_{j}^{\varepsilon}$ it suffices to show that
$w_{j}\geq w_{k}+\varepsilon/\sqrt{2}$. Define $v\in\mathbb{R}^{L}$ as

$v_{\ell}=\begin{cases}w_{k},&\ell=j,\\\
w_{k}+\varepsilon/\sqrt{2},&\ell=k,\\\ w_{\ell},&\ell\neq j,k.\end{cases}.$

Then $v\in R_{k}^{\varepsilon}$, and so we must have $\|w-v\|\geq\varepsilon$
since $k\not\in\operatorname{argmax}^{\varepsilon}(w)$. We calculate

$\|w-v\|^{2}=(w_{j}-v_{j})^{2}+(w_{k}-v_{k})^{2}=(w_{j}-w_{k})^{2}+(\varepsilon/\sqrt{2})^{2}.$

This implies that we must have $(w_{j}-w_{k})^{2}\geq\varepsilon^{2}/2$, which
completes the proof.

####  A.6.3 Proof of Lemma 15

First we verify existence and uniqueness of $a$. Let

$f(a)=w_{j}-a+\sum_{k\neq j}(w_{k}-a)_{+}.$

Then $f:\mathbb{R}\rightarrow\mathbb{R}$ is a continuous and strictly
decreasing bijection, so $f(a)=\varepsilon/\sqrt{2}$ must have a unique
solution.

Next let $a$ be this unique solution and let $v$ be defined as in (10). Now we
verify that $v$ is the unique solution to the optimization problem

$v=\operatorname{argmin}_{u\in\mathbb{R}^{L}}\left\\{\frac{1}{2}\|u-w\|^{2}:u_{j}\geq
u_{k}+\varepsilon/\sqrt{2}\textnormal{ for all $k\neq j$}\right\\},$

which defines the projection of $w$ to $R_{j}^{\varepsilon}$. By the
Karush–Kuhn–Tucker (KKT) conditions for first-order optimality, it is
sufficient to verify that

$v-w-\sum_{k\neq j}\lambda_{k}(\mathbf{e}_{j}-\mathbf{e}_{k})=0,$ (11)

where $\lambda_{k}\geq 0$ for all $k$, and $\lambda_{k}=0$ for any inactive
constraints, i.e., any $k$ for which $v_{j}>v_{k}+\varepsilon/\sqrt{2}$. Let
$\lambda_{k}=(w_{k}-a)_{+}$ for each $k\neq j$. Note that if
$v_{j}>v_{k}+\varepsilon/\sqrt{2}$ (i.e., an inactive constraint), then we
must have $w_{k}<a$ and so $\lambda_{k}=0$, by construction of $v$. Then, for
any $\ell\neq j$,

$\left(v-w-\sum_{k\neq
j}\lambda_{k}(\mathbf{e}_{j}-\mathbf{e}_{k})\right)_{\ell}=v_{\ell}-w_{\ell}+\lambda_{\ell}=a\wedge
w_{\ell}-w_{\ell}+(w_{\ell}-a)_{+}=0,$

where $\mathbf{e}_{k}$ denotes the $k$th canonical basis vector in
$\mathbb{R}^{L}$. This proves that the $\ell$th coordinate of the system of
equations (11) holds for each $\ell\neq j$. Now we verify the $j$th
coordinate. We have

$\left(v-w-\sum_{k\neq
j}\lambda_{k}(\mathbf{e}_{j}-\mathbf{e}_{k})\right)_{j}=v_{j}-w_{j}-\sum_{k\neq
j}\lambda_{k}\\\ =a+\varepsilon/\sqrt{2}-w_{j}-\sum_{k\neq
j}(w_{k}-a)_{+}=\varepsilon/\sqrt{2}-f(a)=0,$

by our choice of $a$ as the solution to $f(a)=\varepsilon/\sqrt{2}$. This
verifies the KKT conditions, and thus, $v$ is the projection of $w$ to
$R_{j}^{\varepsilon}$ as claimed.

##  Appendix B Extension to randomized algorithms

As mentioned in Section 2, in many applications it is common to use
randomization in the construction of a classification procedure, in the
learning algorithm $\mathcal{A}$ and/or in the selection rule $s$. In this
section we formalize this more general framework, and will see how our results
apply.

First, in the non-random setting, a learning algorithm $\mathcal{A}$ maps a
data set $\mathcal{D}\in(\mathcal{X}\times[L])^{n}$ to a fitted probability
estimate function $\hat{p}:\mathcal{X}\to\Delta_{L-1}$—we write the regression
procedure as $\hat{p}=\mathcal{A}(\mathcal{D})$. In the randomized setting, we
also allow for external randomness, $\hat{p}=\mathcal{A}(\mathcal{D};\xi)$,
where $\xi\sim\textnormal{Uniform}[0,1]$ is a random seed (e.g., we might use
$\xi$ to randomly shuffle the training data when running stochastic gradient
descent).

Next, in the non-random setting, a selection rule $s$ is a map from
$\Delta_{L-1}$ to subsets of $[L]$, resulting in a candidate set of labels
$\hat{S}=s(\hat{p}(x))$. In the randomized setting (for instance, if $s$ is
the argmax but with ties broken at random), we instead include a random seed
$\zeta\sim\textnormal{Uniform}[0,1]$, and write $\hat{S}=s(\hat{p}(x);\zeta)$.
Formally, then, we have
$s:\Delta_{L-1}\times[0,1]\to\raisebox{1.79993pt}{\Large$\wp$}([L])$. (Note
that the selection rule proposed in this work—the inflated argmax—is itself
_not_ random.)

Combining the two stages of the procedure, then, the classification algorithm
is given by $\mathcal{C}=s\circ\mathcal{A}$, where given a training set
$\mathcal{D}$ and a test point $x$, along with i.i.d. random seeds
$\xi,\zeta\sim\textnormal{Uniform}[0,1]$, we return the candidate set of
labels given by

$\mathcal{C}(\mathcal{D},x;\xi,\zeta)=s(\hat{p}(x);\zeta)\textnormal{ where
}\hat{p}=\mathcal{A}(\mathcal{D};\xi).$

Of course, we can derive the non-random setting as a special case, simply by
designing $\mathcal{A}$ and $s$ to not depend on the random seeds $\xi$ and
$\zeta$.

Now, how do the results of this paper extend to this randomized setting?
First, we need to reconsider our definition of selection stability (Definition
1). We will now define it as

$\frac{1}{n}\sum_{i=1}^{n}\Pr\left\\{\hat{S}\cap\hat{S}^{\setminus
i}=\varnothing\right\\}\leq\delta,$

where $\hat{S}=\mathcal{C}(\mathcal{D},x;\xi,\zeta)$ and $\hat{S}^{\setminus
i}=\mathcal{C}(\mathcal{D}^{\setminus i},x;\xi,\zeta)$, and where the
probabilities are taken with respect to the distribution of the random seeds.
We can also consider a notion of algorithmic stability (Definition 2) that
allows for randomization: we say that $\mathcal{A}$ has tail stability
$(\varepsilon,\delta)$ if

$\frac{1}{n}\sum_{i=1}^{n}\Pr\left\\{\|\hat{p}(x)-\hat{p}^{\setminus
i}(x)\|\geq\varepsilon\right\\}\leq\delta,$

where $\hat{p}=\mathcal{A}(\mathcal{D};\xi)$ and $\hat{p}^{\setminus
i}=\mathcal{A}(\mathcal{D}^{\setminus i};\xi)$, and where again probability is
taken with respect to the distribution of $\xi$.

With these definitions in place, we observe that the result of 4 still holds
in this more general setting for the selection rule
$\operatorname{argmax}^{\varepsilon}$: if $\mathcal{A}$ is a randomized
algorithm with tail stability $(\varepsilon,\delta)$ (under the new definition
given above), then the randomized classification algorithm given by
$\mathcal{C}=\operatorname{argmax}^{\varepsilon}\circ\mathcal{A}$ has
selection stability $\delta$ (again, under the new definition given above),
simply due to the $\varepsilon$-compatibility property of
$\operatorname{argmax}^{\varepsilon}$ (Theorem 8). Soloff et al., 2024b
proved Theorem 7 in this more general setting of randomized learning
algorithms $\mathcal{A}$, so Theorem 5 holds for randomized $\mathcal{A}$
under this more general definition of selection stability.

##  Appendix C Extension to finitely many bags

In Section 3.1, we discussed how in practice, the bagged algorithm
$\widetilde{\mathcal{A}}_{m}$ would be constructed by taking an empirical
average over some large number $B$ of sampled bags, rather than computing the
exact expectation $\mathbb{E}_{r}[\mathcal{A}(\mathcal{D}^{r})(x)]$. Now that
we have allowed for randomized learning algorithms (as in Appendix B), we can
formalize this setting: for a random seed $\xi\sim\textnormal{Uniform}[0,1]$,
we define

$\widetilde{\mathcal{A}}_{m}(\mathcal{D},\xi)(x)=\frac{1}{B}\sum_{b=1}^{B}\mathcal{A}(\mathcal{D}^{r_{b}})(x),$

where the random draw of $B$ many bags, $r_{1},\dots,r_{B}$, is generated
using the random seed $\xi$. Recall that Soloff et al., 2024b  showed tail
stability for the bagged version of any algorithm (as restated in Theorem
7)—their work also gives a result for the finite-$B$ case, as follows:

######  Theorem 16 ((Soloff et al., 2024b, )).

For any base learning algorithm $\mathcal{A}$ returning outputs in
$\Delta_{L-1}$, the bagged algorithm $\widetilde{\mathcal{A}}_{m}$ (computed
with a finite number of bags, $B$) has tail stability $(\varepsilon,\delta)$
for any $\varepsilon,\delta>0$ satisfying

$\delta=\varepsilon^{-2}\cdot(1-1/L)\left(\frac{1}{n-1}\cdot\frac{p_{n,m}}{1-p_{n,m}}+\frac{16e^{2}}{B}\right).$

In particular, combined with our result on the $\varepsilon$-compatibility of
the inflated argmax, we have the following generalization of our main result
(Theorem 5):

######  Theorem 17.

Fix any sample size $n$, any bag size $m$, any inflation parameter
$\varepsilon>0$, and any number of bags $B\geq 1$. For any base learning
algorithm $\mathcal{A}$, the classification algorithm
$\mathcal{C}=\operatorname{argmax}^{\varepsilon}\circ\widetilde{A}_{m}$,
obtained by combining the bagged version of $\mathcal{A}$ (with $B$ many bags)
together with the inflated argmax, satisfies selection stability $\delta$
where

$\delta=\varepsilon^{-2}\cdot(1-1/L)\left(\frac{1}{n-1}\cdot\frac{p_{n,m}}{1-p_{n,m}}+\frac{16e^{2}}{B}\right).$
(12)

Comparing to our main result for the derandomized case (Theorem 5, which can
be interpreted as taking $B\to\infty$), we see that this finite-$B$ result is
essentially the same as Theorem 5 once $B\gg n$.

##  Appendix D Compare to a simpler selection rule: the fixed-margin rule

In this section, we compare to the simpler fixed-margin selection rule, which
was defined earlier in (3)—for convenience, we define it again here:

$s_{\textnormal{margin}}^{\varepsilon}(w)=\left\\{j\in[L]:w_{j}>\max_{\ell\in[L]}w_{\ell}-\varepsilon/\sqrt{2}\right\\},$

the set of all indices $j$ that are within some margin of being the maximum.

In this section, we will see that this fixed-margin selection rule is
$\varepsilon$-compatible. Since this rule clearly has the advantage of being
much simpler than the inflated argmax (both in terms of its definition and
interpretability, and in terms of its computation), we might ask whether this
rule is perhaps preferable to the more complex inflated argmax. However, we
will also see that the fixed-margin selection rule can be very inefficient
compared to the inflated argmax: the inflated argmax can never return a larger
set, and will often return a substantially smaller one.

###  D.1 Theoretical results

First, we verify that this simple rule satisfies $\varepsilon$-compatibility
(Definition 3).

######  Proposition 18.

For any $\varepsilon>0$, the selection rule
$s_{\textnormal{margin}}^{\varepsilon}$ is $\varepsilon$-compatible.

In particular, since $s_{\textnormal{margin}}^{\varepsilon}$ clearly contains
the argmax (i.e., $\operatorname{argmax}(w)\subseteq
s_{\textnormal{margin}}^{\varepsilon}(w)$) and satisfies permutation
invariance (in the sense of 9), by 11 this immediately implies that

$s_{\textnormal{margin}}^{\varepsilon}(w)=\\{j\\}\ \Longrightarrow\
\operatorname{argmax}^{\varepsilon}(w)=\\{j\\}$

for all $w\in\mathbb{R}^{L}$, i.e., $\operatorname{argmax}^{\varepsilon}$ is
at least as good as $s_{\textnormal{margin}}^{\varepsilon}$ at returning a
singleton set. However, for this particular selection rule, we can state an
even stronger result:

######  Proposition 19.

For any $w\in\mathbb{R}^{L}$ and any $\varepsilon>0$,

$s_{\textnormal{margin}}^{\varepsilon}(w)\supseteq\operatorname{argmax}^{\varepsilon}(w).$

Figure 3: The fixed-margin selection rule (3), left, and the inflated argmax
(2), right, when $L=3$.

In particular, this theoretical result ensures that the set of candidate
labels $\hat{S}$ returned by the inflated argmax can never be larger than for
the fixed-margin selection rule. In low dimensions, however, the improvement
is small. Indeed, for $L=2$, we actually have
$s_{\textnormal{margin}}^{\varepsilon}(w)=\operatorname{argmax}^{\varepsilon}(w)$
for all $w$. For $L=3$, Figure 3 shows that the inflated argmax is strictly
better (i.e., the set inclusion result of 19 can, for certain values of $w$,
be a strict inclusion), but the difference in this low-dimensional setting
appears minor. Next, however, we will see empirically that as the dimension
$L$ grows, the benefit of the inflated argmax can be substantial.

###  D.2 Simulation to compare the selection rules

In this section, we compare the fixed-margin rule to the inflated argmax using
randomly generated probability weights. For our simulation, we consider the
size of the two sets for various values of $L$. To sample the weights, we draw
a standard Gaussian random vector $Z\sim\mathcal{N}(0,\mathbf{I}_{L})$ and
define $w$ as the softmax of $Z$, i.e.,

$w_{j}=\frac{e^{Z_{j}}}{\sum_{\ell\in[L]}e^{Z_{\ell}}}$

for $j\in[L]$. We then compute $\operatorname{argmax}^{\varepsilon}(w)$ and
$s_{\textnormal{margin}}^{\varepsilon}(w)$, where $\varepsilon=0.1$. In Figure
4, we present results comparing the average size of the sets returned by each
of the two methods. In this setting, we see that the inflated argmax has
substantial gains over the fixed-margin selection rule when the number of
classes $L$ is large. Even when $L=25$, the ratio of the expected sizes is
about $78\%$, so the inflated argmax has a nontrivial advantage in that
setting. When $L=100$, the ratio of expected sizes is $48\%$, meaning the
inflated argmax is on average more than twice as small.

Figure 4: Simulation to compare the selection rules
$\operatorname{argmax}^{\varepsilon}$ and
$s_{\textnormal{margin}}^{\varepsilon}$. The figure shows the average set
size, $|\operatorname{argmax}^{\varepsilon}(w)|$ and
$|s_{\textnormal{margin}}^{\varepsilon}(w)|$, avaeraged over $1,000$ random
draws of $w$ (with standard error bars shown). See Section D.2 for details.

###  D.3 Proofs for the fixed-margin selection rule

####  D.3.1 Proof of 18

Let $w,v\in\mathbb{R}^{L}$ with $s_{\textnormal{margin}}^{\varepsilon}(w)\cap
s_{\textnormal{margin}}^{\varepsilon}(v)=\varnothing$. Let
$j\in\operatorname{argmax}_{\ell}w_{\ell}$ and
$k\in\operatorname{argmax}_{\ell}v_{\ell}$. Then clearly $j\in
s_{\textnormal{margin}}^{\varepsilon}(w)$, so since
$s_{\textnormal{margin}}^{\varepsilon}(w)\cap
s_{\textnormal{margin}}^{\varepsilon}(v)=\varnothing$, this implies $j\not\in
s_{\textnormal{margin}}^{\varepsilon}(v)$. By definition, then,
$v_{j}\leq\max_{\ell}v_{\ell}-\varepsilon/\sqrt{2}=v_{k}-\varepsilon/\sqrt{2}$.
By an identical argument, we have $w_{k}\leq w_{j}-\varepsilon/\sqrt{2}$. Then

$\displaystyle\|w-v\|^{2}$
$\displaystyle\geq(w_{j}-v_{j})^{2}+(w_{k}-v_{k})^{2}$
$\displaystyle\geq(w_{j}-v_{j})^{2}_{+}+(v_{k}-w_{k})^{2}_{+}$
$\displaystyle\geq\left(w_{j}-(v_{k}-\varepsilon/\sqrt{2})\right)_{+}^{2}+\left(v_{k}-(w_{j}-\varepsilon/\sqrt{2})\right)_{+}^{2}$
$\displaystyle\geq\inf_{t\in\mathbb{R}}\left\\{(\varepsilon/\sqrt{2}-t)_{+}^{2}+(\varepsilon/\sqrt{2}+t)_{+}^{2}\right\\}$
$\displaystyle=2(\varepsilon/\sqrt{2})^{2}=\varepsilon^{2}.$

This proves $\|w-v\|\geq\varepsilon$, and therefore, we have shown that
$\varepsilon$-compatibility is satisfied.

####  D.3.2 Proof of 19

Suppose $j\in\operatorname{argmax}^{\varepsilon}(w)$. Then
$\textnormal{dist}(w,R_{j}^{\varepsilon})<\varepsilon$, so we can find some
$v\in R_{j}^{\varepsilon}$ with $\|w-v\|<\varepsilon$. Then, for any $k\neq
j$,

$\displaystyle w_{j}-w_{k}$
$\displaystyle=v_{j}-v_{k}+(w_{j}-v_{j})+(v_{k}-w_{k})$
$\displaystyle\geq\varepsilon/\sqrt{2}+(w_{j}-v_{j})+(v_{k}-w_{k})\textnormal{
since $v\in R_{j}^{\varepsilon}$}$
$\displaystyle\geq\varepsilon/\sqrt{2}-\big{(}|w_{j}-v_{j}|+|w_{k}-v_{k}|\big{)}$
$\displaystyle\geq\varepsilon/\sqrt{2}-\sqrt{2}\sqrt{(w_{j}-v_{j})^{2}+(w_{k}-v_{k})^{2}}$
$\displaystyle\geq\varepsilon/\sqrt{2}-\sqrt{2}\|w-v\|$
$\displaystyle>\varepsilon/\sqrt{2}-\sqrt{2}\cdot\varepsilon$
$\displaystyle=-\varepsilon/\sqrt{2}.$

Since this holds for all $k\neq j$, then,

$w_{j}>\max_{k\in[L]}w_{k}-\varepsilon/\sqrt{2},$

which proves that $j\in s_{\textnormal{margin}}^{\varepsilon}(w)$.

## References

  * Andonova et al., (2002) Andonova, S., Elisseeff, A., Evgeniou, T., and Pontil, M. (2002).  A simple algorithm for learning stable machines.  In European Conference on Artificial Intelligence (ECAI), pages 513–517.
  * Angelopoulos et al., (2023) Angelopoulos, A. N., Bates, S., et al. (2023).  Conformal prediction: A gentle introduction.  Foundations and Trends® in Machine Learning, 16(4):494–591.
  * Blondel et al., (2020) Blondel, M., Martins, A. F. T., and Niculae, V. (2020).  Learning with Fenchel-Young losses.  J. Mach. Learn. Res., 21:Paper No. 35, 69.
  * Bousquet and Elisseeff, (2002) Bousquet, O. and Elisseeff, A. (2002).  Stability and generalization.  The Journal of Machine Learning Research, 2:499–526.
  * (5) Breiman, L. (1996a).  Bagging predictors.  Machine learning, 24(2):123–140.
  * (6) Breiman, L. (1996b).  Heuristics of instability and stabilization in model selection.  The Annals of Statistics, 24(6):2350–2383.
  * Bridle, (1990) Bridle, J. S. (1990).  Probabilistic interpretation of feedforward classification network outputs, with relationships to statistical pattern recognition.  In Neurocomputing (Les Arcs, 1989), volume 68 of NATO Adv. Sci. Inst. Ser. F: Comput. Systems Sci., pages 227–236. Springer, Berlin.
  * Caramanis et al., (2011) Caramanis, C., Mannor, S., and Xu, H. (2011).  Robust optimization in machine learning.
  * Chzhen et al., (2021) Chzhen, E., Denis, C., Hebiri, M., and Lorieul, T. (2021).  Set-valued classification–overview via a unified framework.  arXiv preprint arXiv:2102.12318.
  * de Brébisson and Vincent, (2016) de Brébisson, A. and Vincent, P. (2016).  An exploration of softmax alternatives belonging to the spherical loss family.  In International Conference on Learning Representations (ICLR).
  * Del Coz et al., (2009) Del Coz, J. J., Díez, J., and Bahamonde, A. (2009).  Learning nondeterministic classifiers.  Journal of Machine Learning Research, 10(10).
  * Friedler et al., (2019) Friedler, S. A., Scheidegger, C., Venkatasubramanian, S., Choudhary, S., Hamilton, E. P., and Roth, D. (2019).  A comparative study of fairness-enhancing interventions in machine learning.  In Proceedings of the conference on fairness, accountability, and transparency, pages 329–338.
  * Grycko, (1993) Grycko, E. (1993).  Classification with set-valued decision functions.  In Information and Classification: Concepts, Methods and Applications Proceedings of the 16th Annual Conference of the “Gesellschaft für Klassifikation eV” University of Dortmund, April 1–3, 1992, pages 218–224. Springer.
  * Hechtlinger et al., (2018) Hechtlinger, Y., Póczos, B., and Wasserman, L. (2018).  Cautious deep learning.  arXiv preprint arXiv:1805.09460.
  * Huang and Vishnoi, (2019) Huang, L. and Vishnoi, N. (2019).  Stable and fair classification.  In International Conference on Machine Learning, pages 2879–2890. PMLR.
  * Lei, (2014) Lei, J. (2014).  Classification with confidence.  Biometrika, 101(4):755–769.
  * Martins and Astudillo, (2016) Martins, A. and Astudillo, R. (2016).  From softmax to sparsemax: A sparse model of attention and multi-label classification.  In International conference on machine learning, pages 1614–1623. PMLR.
  * Murdoch et al., (2019) Murdoch, W. J., Singh, C., Kumbier, K., Abbasi-Asl, R., and Yu, B. (2019).  Definitions, methods, and applications in interpretable machine learning.  Proceedings of the National Academy of Sciences, 116(44):22071–22080.
  * Papadopoulos et al., (2002) Papadopoulos, H., Proedrou, K., Vovk, V., and Gammerman, A. (2002).  Inductive confidence machines for regression.  In Machine learning: ECML 2002: 13th European conference on machine learning Helsinki, Finland, August 19–23, 2002 proceedings 13, pages 345–356. Springer.
  * Paszke et al., (2019) Paszke, A., Gross, S., Massa, F., Lerer, A., Bradbury, J., Chanan, G., Killeen, T., Lin, Z., Gimelshein, N., Antiga, L., et al. (2019).  PyTorch: An imperative style, high-performance deep learning library.  Advances in neural information processing systems, 32.
  * Quach et al., (2023) Quach, V., Fisch, A., Schuster, T., Yala, A., Sohn, J. H., Jaakkola, T. S., and Barzilay, R. (2023).  Conformal language modeling.  arXiv preprint arXiv:2306.10193.
  * Reuleaux, (1876) Reuleaux, F. (1876).  Kinematics of machinery: outlines of a theory of machines.  Dover, New York.  Reprinted as “The Kinematics of Machinery”.
  * Sadinle et al., (2019) Sadinle, M., Lei, J., and Wasserman, L. (2019).  Least ambiguous set-valued classifiers with bounded error levels.  J. Amer. Statist. Assoc., 114(525):223–234.
  * (24) Soloff, J. A., Barber, R. F., and Willett, R. (2024a).  Bagging provides assumption-free stability.  Journal of Machine Learning Research, 25(131):1–35.
  * (25) Soloff, J. A., Barber, R. F., and Willett, R. (2024b).  Stability via resampling: statistical problems beyond the real line.  arXiv preprint arXiv:2405.09511.
  * Vovk et al., (2005) Vovk, V., Gammerman, A., and Shafer, G. (2005).  Algorithmic learning in a random world.  Springer, New York.
  * Weisstein, (2005) Weisstein, E. W. (2005).  Reuleaux triangle.  MathWorld–A Wolfram Web Resource.  https://mathworld.wolfram.com/ReuleauxTriangle.html.
  * Xiao et al., (2017) Xiao, H., Rasul, K., and Vollgraf, R. (2017).  Fashion-MNIST: a novel image dataset for benchmarking machine learning algorithms.
  * Yu and Barter, (2024) Yu, B. and Barter, R. (2024).  Veridical Data Science: The Practice of Responsible Data Analysis and Decision Making.  MIT Press.
  * Yu and Kumbier, (2020) Yu, B. and Kumbier, K. (2020).  Veridical data science.  Proceedings of the National Academy of Sciences, 117(8):3920–3929.
  * Zhao et al., (2023) Zhao, S., Li, Q., Yang, Y., Wen, J., and Luo, W. (2023).  From softmax to nucleusmax: A novel sparse language model for chinese radiology report summarization.  ACM Transactions on Asian and Low-Resource Language Information Processing, 22(6):1–21.