# Stability via resampling: statistical problems beyond the real line

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

Model averaging techniques based on resampling methods (such as bootstrapping
or subsampling) have been utilized across many areas of statistics, often with
the explicit goal of promoting stability in the resulting output. We provide a
general, finite-sample theoretical result guaranteeing the stability of
bagging when applied to algorithms that return outputs in a general space, so
that the output is not necessarily a real-valued—for example, an algorithm
that estimates a vector of weights or a density function. We empirically
assess the stability of bagging on synthetic and real-world data for a range
of problem settings, including causal inference, nonparametric regression, and
Bayesian model selection.

##  1 Introduction

Statistical stability holds when our conclusions, estimates, fitted models,
predictions, or decisions are insensitive to small perturbations of data.
Stability is a pillar of modern data science (Yu and Kumbier, 2020).
Nonparametric tools for statistical inference—e.g. the bootstrap, subsampling,
the Jackknife—perturb data to assess variability and can thus be viewed as
aspects of stability analysis (Yu, 2013). While assessing stability is in
general statistically intractable without placing assumptions on the data (Kim
and Barber, 2023), these resampling-based approaches can also be used to
reduce the issue of instability: for instance, Breiman (1996a, b) originally
proposed bagging—which applies bootstrap resampling to construct an ensemble
of machine learning models—as an off-the-shelf stabilizer for regression and
classification. In this paper, we show how bagging can be used to stabilize
statistical methods far beyond its origins in supervised learning problems.

###  1.1 Problem settings

We begin by highlighting some specific statistical problems where provably
stable methods are desirable.

##### Problem setting 1 (background): the supervised learning setting

In the supervised learning setting, we are given a labeled data set
$(X_{1},Y_{1}),\dots,(X_{n},Y_{n})$ and are tasked with constructing a
regression function $\hat{f}$ to predict $Y$ given $X$ for new data points.
Many modern methods from the recent machine learning literature offer powerful
predictive performance, but suffer from high sensitivity to the training data,
where a prediction $\hat{f}(X)$ (at a new test point $X$) may vary immensely
if the training data is perturbed slightly. Interestingly, both high accuracy
and poor stability can often be a consequence of working with highly
overparameterized models.

Suppose that the response $Y$ is real-valued, so that our regression functions
take values in $\mathbb{R}$, i.e., $\hat{f}:\mathcal{X}\to\mathbb{R}$. Recent
work by Soloff et al. (2024) establishes that bagging the predictive
model—that is, re-fitting $\hat{f}$ on bootstrapped training sets, and then
returning the averaged prediction—yields an assumption-free guarantee on the
stability of the prediction $\hat{f}(X)\in\mathbb{R}$.

While Soloff et al. (2024)’s work is presented as a result for supervised
learning, it can be viewed more generally as a guarantee of stability for
_any_ algorithm that is trained on data and returns a real-valued output. But
even with this generalization, this result only covers a very specific
statistical scenario. For many statistical problems, we might be interested in
the stability of a procedure that returns a more complex object, such as a
function or a vector of weights, rather than a real-valued output.

##### Problem setting 2: reproducible Bayesian inference

In Bayesian data analysis, given a data set $\mathcal{D}$ containing $n$
observations drawn from some distribution $p_{\theta}$, we are interested in
the posterior distribution of the parameter $\theta$. For instance, $\theta$
might represent a particular model from a finite set of possibilities, in
which case the posterior distribution on $\theta$ is simply a vector of
weights on this finite set. Alternatively, $\theta$ might take values in a
continuous space such as $[0,1]$ (e.g., if $\theta$ represents a probability)
or $\mathbb{R}^{d}$ (e.g., coefficients for some regression model), in which
case the posterior distribution is a distribution over $[0,1]$ or over
$\mathbb{R}^{d}$—a more complex object.

Bayesian inference may sometimes lead to answers that are highly sensitive
with respect to perturbations of the data. In particular, this can occur in
settings where the model is misspecified (but nonetheless, is a reasonably
good approximation to the observed data, so that computing the posterior is
still a reasonable goal). To alleviate this instability, the BayesBag approach
(Bühlmann, 2014) averages posterior distributions over different resamplings
of the data $\mathcal{D}$. Huggins and Miller (2023) analyzes the accuracy and
stability properties of the BayesBag approach for the specific problem of
Bayesian model selection, where $\theta$ lies in a finite set of models
$\mathfrak{M}=\\{\mathfrak{m}_{1},\dots,\mathfrak{m}_{K}\\}$ (and so the
posterior on $\theta$ is given by a vector of weights of length $K$, i.e.,
with the $k$th weight indicating the posterior probability
$\mathbb{P}\\{\theta=\mathfrak{m}_{k}\mid\mathcal{D}\\}$). Their work
establishes BayesBag as an empirically effective method with asymptotic
guarantees, but has not yielded finite-sample guarantees even in the finite
model setting. Of course, we can also ask this question even more generally:
will BayesBag lead to greater stability in a general setting, where $\theta$
may lie in an arbitrary space (e.g., if the posterior is a distribution on
$[0,1]$, or on $\mathbb{R}^{d}$, etc.)?

##### Problem setting 3: synthetic controls for causal inference

The synthetic control method (SCM) is a popular approach for estimating causal
effects with panel data (Abadie and Gardeazabal, 2003; Abadie et al., 2010,
2015; Abadie, 2021). The method is particularly useful in settings with a
single treated unit and a relatively small number of control units. A
synthetic control is a weighted combination of control units, where the
weights are chosen to match the treated unit’s pre-treatment outcomes and
other covariates. SCM estimates the counterfactual—what outcomes would we
expect for the treated unit had the policy not been adopted?—using the
synthetic control’s post-treatment outcomes. For example, Abadie et al. (2015)
applies SCM to the 1990 German reunification and its impact on per capita GDP
in West Germany. Per capita GDP was increasing prior to reunification and
continued to increase in subsequent years, for this was a period of sustained
growth across industrialized countries. To assess the counterfactual, a
‘synthetic West Germany’ is constructed as follows: we find the convex
combination of other industrialized countries such that the pre-1990 economic
trend most closely resembles that of West Germany. Abadie et al. (2015) finds
that West Germany’s economy did not grow as much as it would have without
reunification.

Often SCM’s weights are of direct interest, as they say something about which
peer countries are most relevant in this counterfactual assessment. However,
in order to reliably interpret the weights, we should first ensure they are
stable (Murdoch et al., 2019). For instance, if we remove _Denmark_ from the
analysis (a country which receives zero weight in the synthetic West Germany)
then there is no longer a positive weight on the United States, which
previously received $22\%$ of the weight. This sensitivity appears to betray
our intuitive interpretation of the meaning of SCM’s weights. Can a
resampling-type approach to SCM alleviate these issues?

###  1.2 Our contributions

Algorithmic stability is desirable for statistical methods far beyond its
origins in statistical learning theory (Devroye and Wagner, 1979a, b; Bousquet
and Elisseeff, 2002). Broadly, stability is considered a ‘prerequisite for
trustworthy interpretations’ (Murdoch et al., 2019). The goal of this paper is
to establish quantitative guarantees on the stability of bagging as a general-
purpose tool across a wide range of problems where the output—the output being
estimated from the data—might not be real-valued and might even be infinite-
dimensional. Table 1 lists a range of problem settings where our framework can
be directly applied (including the settings discussed above, namely,
supervised learning, Bayesian inference, and synthetic controls). Formally,
our results will apply to any algorithm, randomized or deterministic, with
outputs in bounded subset of a separable Hilbert or Banach space—we give
details below. We highlight one particularly important finding: our main
results will show that for a Hilbert space, there is _no cost_ to high
dimensionality in terms of the stability guarantees that can be obtained (in
contrast, for the more general setting of a Banach space, this is no longer
the case).

##### Outline

The remainder of this paper is organized as follows. In Section 2, we define a
framework for algorithmic stability in any metric space, and then define
bagging for an arbitrary algorithm in Section 3. Our main results on the
stability of bagging are presented in Section 4 and Section 5. In Section 6,
we take a closer look at the implications of our results on a range of
specific settings appearing in Table 1, and discuss some consequences of our
stability guarantee. We conclude with a discussion in Section 7. All proofs
are deferred to the Appendix.

Applications | Output Space $\mathcal{W}$ | Metric $\mathsf{d}$ on $\mathcal{W}$
---|---|---
Supervised learning | Predicted value range $[a,b]\subset\mathbb{R}$ | Absolute value $|\cdot|$
| Bayesian model selection
---
Synthetic controls
Probability simplex $\Delta_{K-1}$ |  | Euclidean dist. $\|\cdot\|_{2}$ or
---
Total variation dist. $\mathsf{d}_{\mathsf{TV}}$
| Constrained Lasso ($q=1$)
---
Constrained Ridge ($q=2$)
$\ell_{q}$-ball $\mathbb{B}_{q}^{d}(R)$ |  $\ell_{q}$-norm $\|\cdot\|_{q}$
| Density estimation
---
(e.g., Bayesian posterior)
$L$-Lipschitz densities on $[0,1]$ | Total variation dist. $\mathsf{d}_{\mathsf{TV}}$
Nonparametric regression | Sobolev ball of functions $\mathbb{R}^{d}\to\mathbb{R}$ | Sobolev norm $\|\cdot\|_{s,2}$
Table 1: Example problem settings for our stability framework, including the
examples highlighted in Problem Settings 1–3 (in Section 1).

##  2 Framework: algorithmic stability in a metric space

We define a _randomized algorithm_ $\mathcal{A}$ as function that takes in a
finite sequence of data $\mathcal{D}=(Z_{i})_{i\in[n]}$, as well as some
independent randomness $\xi\sim\text{Uniform}([0,1])$, and returns an output
$\mathcal{A}(\mathcal{D};\xi)$. The term $\xi$, which we refer to as the
_random seed_ , determines all randomness used by the algorithm—for instance,
a random initialization or a random split of the data into batches. The
algorithm $\mathcal{A}$ is deterministic if it does not depend on the seed
$\xi$.

This notion of a (deterministic or randomized) algorithm covers many disparate
problems in statistical learning, such as the problem settings introduced in
Section 1.1.

  * •

Consider a supervised learning setting (Problem Setting 1) with data points
$Z_{i}=(X_{i},Y_{i})\in\mathcal{X}\times\mathbb{R}$. If we are interested in
predicting $Y$ given a new feature $X=x$, our algorithm $\mathcal{A}$ can
return $\hat{f}(x)\in\mathbb{R}$, the predicted value at $X=x$, as its output.
(Alternatively, if we are interested in building a model for $Y$ given $X$,
our algorithm $\mathcal{A}$ can return an entire fitted regression function
$\hat{f}$—a more complex output.)

  * •

In Bayesian inference (Problem Setting 2), we observe data points
$Z_{1},\dots,Z_{n}$ drawn from some distribution $p_{\theta}$, where
$\theta\in\Theta$ is the unknown parameter. To perform Bayesian inference on
$\theta$, our algorithm $\mathcal{A}$ would then return the posterior
distribution over $\theta\in\Theta$.

  * •

For the problem of synthetic controls (Problem Setting 3), our data points
$Z_{i}$ are the feature vectors for each control unit $i=1,\dots,n$, and the
algorithm $\mathcal{A}$ returns a weight vector of dimension $n$, specifying
the convex combination of control units that define the synthetic control.

In some settings it may be common to use a randomized, rather than
deterministic, algorithm. For example, in Bayesian inference, in a high-
dimensional setting we might choose to approximate the posterior through
(random) sampling rather than computing the density explicitly.

Formally, we write

$\mathcal{A}:\bigcup_{n\geq 0}\mathcal{Z}^{n}\times[0,1]\to\mathcal{W},$

where the output takes values in some space $\mathcal{W}$ equipped with a
metric $\mathsf{d}$—see Table 1 for examples. From this point on, we write
$\hat{w}=\mathcal{A}(\mathcal{D};\xi)$ to denote the output, regardless of the
nature of the output (e.g., a finite-dimensional vector of weights for the
synthetic control problem, or a posterior density function for the Bayesian
inference setting, etc.).111We assume throughout that
$t\mapsto\mathcal{A}((Z_{i})_{i\in[n]},t)$, which is a function from $[0,1]$
to $\mathcal{W}$, is measurable for any fixed
$(Z_{i})_{i\in[n]}\in\mathcal{Z}^{n}$. Measurability is defined with respect
to the Borel $\sigma$-algebra on $(\mathcal{W},\mathsf{d})$.

###  2.1 Defining stability

There are many related definitions of algorithmic stability in the literature.
We focus on stability of the algorithm’s output with respect to dropping a
data point: will the output $\hat{w}$ of the algorithm be similar if one data
point is removed at random from the training data set?

######  Definition 1.

An algorithm $\mathcal{A}$ has _mean-square stability_ $\beta^{2}$ with
respect to a metric $\mathsf{d}$ if, for all data sets
$\mathcal{D}\in\mathcal{Z}^{n}$ of size $n$,

$\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}_{\xi}\left[\mathsf{d}^{2}\left(\hat{w},\hat{w}^{\setminus
i}\right)\right]\leq\beta^{2}.$

where $\hat{w}=\mathcal{A}(\mathcal{D};\xi)$ is the output of the algorithm
run on the full data set $\mathcal{D}$, $\hat{w}^{\setminus
i}=\mathcal{A}(\mathcal{D}^{\setminus i};\xi)$ is the output of the algorithm
run on the leave-one-out data set $\mathcal{D}^{\setminus i}=(Z_{j})_{j\neq
i}$, and the expectation averages over $\xi\sim\text{Uniform}[0,1]$.

This first definition captures the average perturbation in $\hat{w}$ when one
data point is removed at random from some data set
$\mathcal{D}\in\mathcal{Z}^{n}$—that is,
$\mathsf{d}^{2}\left(\hat{w},\hat{w}^{\setminus i}\right)$, averaged over
which data point $i$ gets removed. An alternative way to define stability is
to control the tails of these perturbations
$\mathsf{d}\left(\hat{w},\hat{w}^{\setminus i}\right)$, ensuring that not too
many of them can be large.

######  Definition 2.

An algorithm $\mathcal{A}$ has _tail stability_ $(\varepsilon,\delta)$ with
respect to a metric $\mathsf{d}$ if, for all data sets
$\mathcal{D}\in\mathcal{Z}^{n}$ of size $n$,

$\frac{1}{n}\sum_{i=1}^{n}\mathbb{P}_{\xi}\left\\{\mathsf{d}(\hat{w},\hat{w}^{\setminus
i})\geq\varepsilon\right\\}\leq\delta,$

where $\hat{w}=\mathcal{A}(\mathcal{D};\xi)$ and $\hat{w}^{\setminus
i}=\mathcal{A}(\mathcal{D}^{\setminus i};\xi)$ as in Definition 1.

The terminology ‘$\mathcal{A}$ has mean-square stability $\beta^{2}$’ in
Definition 1, and ‘$\mathcal{A}$ has tail stability $(\varepsilon,\delta)$’ in
Definition 2, suppresses the dependence on the sample size $n$. Since we
perform a non-asymptotic stability analysis, we treat $n$ as a fixed positive
integer throughout.

These two notions of stability are clearly very similar in flavor, and the
following simple result shows that each implies the other.

######  Proposition 3.

If $\mathcal{A}$ has mean-square stability $\beta^{2}$ with respect to
$\mathsf{d}$, then $\mathcal{A}$ has tail stability $(\varepsilon,\delta)$
with respect to $\mathsf{d}$, for any $\varepsilon,\delta\geq 0$ satisfying
$\delta\varepsilon^{2}\geq\beta^{2}$. Conversely, if $\mathcal{A}$ has tail
stability $(\varepsilon,\delta)$ with respect to $\mathsf{d}$, then
$\mathcal{A}$ has mean-square stability
$\beta^{2}=\varepsilon^{2}+\delta\,\sup_{w,w^{\prime}\in\mathcal{W}}\mathsf{d}^{2}\left(w,w^{\prime}\right)$
with respect to $\mathsf{d}$.

###### Proof.

The first claim follows by Markov’s inequality, since

$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\mathbb{P}_{\xi}\left\\{\mathsf{d}\left(\hat{w},\hat{w}^{\setminus
i}\right)\geq\varepsilon\right\\}$
$\displaystyle=\mathbb{E}_{\xi}\left[\frac{1}{n}\sum_{i=1}^{n}\mathbf{1}\left\\{\mathsf{d}^{2}\left(\hat{w},\hat{w}^{\setminus
i}\right)\geq\varepsilon^{2}\right\\}\right]$
$\displaystyle\leq\varepsilon^{-2}\mathbb{E}_{\xi}\left[\frac{1}{n}\sum_{i=1}^{n}\mathsf{d}^{2}\left(\hat{w},\hat{w}^{\setminus
i}\right)\right]\leq\varepsilon^{-2}\beta^{2}.$

For the second claim, since $\mathsf{d}^{2}\left(\hat{w},\hat{w}^{\setminus
i}\right)\leq\varepsilon^{2}+\mathbf{1}\left\\{\mathsf{d}\left(\hat{w},\hat{w}^{\setminus
i}\right)\geq\varepsilon\right\\}\cdot\sup_{w,w^{\prime}\in\mathcal{W}}\mathsf{d}^{2}\left(w,w^{\prime}\right)$,
we have

$\displaystyle\frac{1}{n}\sum_{i=1}^{n}\mathbb{E}_{\xi}\left[\mathsf{d}^{2}\left(\hat{w},\hat{w}^{\setminus
i}\right)\right]$
$\displaystyle\leq\varepsilon^{2}+\frac{1}{n}\sum_{i=1}^{n}\mathbb{P}_{\xi}\left\\{\mathsf{d}\left(\hat{w},\hat{w}^{\setminus
i}\right)\geq\varepsilon\right\\}\cdot\sup_{w,w^{\prime}\in\mathcal{W}}\mathsf{d}^{2}\left(w,w^{\prime}\right)$
$\displaystyle\leq\varepsilon^{2}+\delta\cdot\sup_{w,w^{\prime}\in\mathcal{W}}\mathsf{d}^{2}\left(w,w^{\prime}\right),$

completing the proof. ∎

###  2.2 Related notions of stability

Our definitions of stability are stronger than many common definitions in the
literature. For instance, we require the algorithm to be stable for _any_
fixed data set $\mathcal{D}$ to certify the stability of the method without
looking at the specific data set. Definition 2 implies in particular that the
algorithm is stable when the data $(Z_{i})_{i\in[n]}$ are i.i.d. from some
distribution $P$, but our definition does not require any stochastic model.

In many cases, practitioners care specifically about the stability of the
algorithm’s output $\hat{w}$, e.g., if we want to interpret the values of some
learned weights, cluster assignments, etc. However, in other cases, we may be
more concerned with the stability of some downstream decisions or outcomes
based on an analysis. For instance, in supervised learning (Problem Setting
1), we may desire stability of a loss incurred by the output $\hat{w}$—e.g., a
predictive error $\ell(\hat{w},(x,y))=|y-\hat{w}(x)|$. These types of loss-
based stability properties can often be inferred from stability of the output
$\hat{w}$ (see Soloff et al. (2024) for further discussion). In this paper, we
only consider the stability of $\hat{w}$ itself.

##  3 Constructing a bagged algorithm via resampling

In this section, we formally define the construction of a bagged algorithm.

Stability is about (in)sensitivity to perturbed data sets. Resampling methods
deliberately perturb the data, so they are naturally a useful tool for
stabilizing algorithms: given a data set $\mathcal{D}\in\mathcal{Z}^{n}$
containing $n$ data points, we might retrain our algorithm repeatedly on
resampled data sets, and then average the resulting outputs. There are several
common schemes for constructing these resampled data sets. Here, we mention
two of the most common methods:

  * •

_Bootstrapping_ constructs each resampled data set by sampling $m$ indices
$i_{1},\dots,i_{m}$ uniformly with replacement from $[n]$ (the term
“bootstrap” commonly refers to this method with the choice $m=n$, so that the
resampled data sets are of the same size as the original $\mathcal{D}$).

  * •

Alternatively, we may instead sample data points without replacement, so that
each contains a uniformly drawn subset containing $m$ of the original $n$ data
points (e.g., with $m=n/2$). This type of approach is known as _subbagging_.

We now define some notation that unifies these two variants (as well as
allowing for other resampling schemes, as we discuss below).

Following the framework and notation of Soloff et al. (2024), we define the
set of finite sequences of indices in $[n]$:

$\textnormal{seq}_{[n]}:=\\{(i_{1},\ldots,i_{k}):k\geq
0,\,i_{1},\ldots,i_{k}\in[n]\\}.$

A sequence $r=(i_{1},\dots,i_{m})\in\textnormal{seq}_{[n]}$ is commonly
referred to as a _bag_ in the context of data set resampling. For any
$r\in\textnormal{seq}_{[n]}$, define $\mathcal{D}^{r}=(Z_{i})_{i\in
r}\in\mathcal{Z}^{m}$, the corresponding data set. 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}$.

We are now ready to define the bagged version of any base algorithm
$\mathcal{A}$, using a particular resampling distribution
$\mathcal{Q}_{n}$.222While the term “bagging” is sometimes used as a synonym
for bootstrap averaging (i.e., averaging outputs after sampling from the data
set without replacement) in the literature, in this paper we will use
“bagging” more broadly, allowing for any resampling distribution
$\mathcal{Q}_{n}$. From this point on, we assume that $\mathcal{W}$ is a
closed and compact subset of a separable Banach space—this technical condition
ensures that taking an average or an expected value of outputs in
$\mathcal{W}$ is well-defined.

######  Definition 4 (Bagging for an algorithm $\mathcal{A}$).

Given a data set $\mathcal{D}\in\mathcal{Z}^{n}$ and a distribution
$\mathcal{Q}_{n}$ on $\textnormal{seq}_{[n]}$, return the output

$\hat{w}=\frac{1}{B}\sum_{b=1}^{B}\mathcal{A}(\mathcal{D}^{r_{b}};\xi_{b}),$

where
$(r_{b},\xi_{b})\stackrel{{\scriptstyle\textnormal{iid}}}{{\sim}}\mathcal{Q}_{n}\times\textnormal{Uniform}[0,1]$.

For example, bootstrapping corresponds to choosing $\mathcal{Q}_{n}$ as the
uniform distribution over all sequences of length $m$, while subbagging takes
$\mathcal{Q}_{n}$ as the uniform distribution over sequences of $m$ _distinct_
elements.

For our theoretical analysis, it will be useful to consider taking
$B\to\infty$, so as to “derandomize” the procedure.

######  Definition 5 (Derandomized bagging for an algorithm $\mathcal{A}$).

Given a data set $\mathcal{D}\in\mathcal{Z}^{n}$ and a distribution
$\mathcal{Q}_{n}$ on $\textnormal{seq}_{[n]}$, return the output

$\hat{w}=\mathbb{E}_{r,\xi}[\mathcal{A}(\mathcal{D}^{r};\xi)],$

where the expected value is computed with respect to
$(r,\xi)\sim\mathcal{Q}_{n}\times\text{Uniform}[0,1]$.

###  3.1 A general framework for bagging

The two variants discussed above, bagging and subbagging, are special cases of
a broader framework, allowing for more flexibility in the choice of the
resampling distribution $\mathcal{Q}_{n}$. (See Soloff et al. (2024) for a
more complete discussion of this framework, and for examples of resampling
methods beyond bagging and subbagging.) Here we formalize the conditions that
$\mathcal{Q}_{n}$ needs to satisfy, in order for our main results to be
applied.

######  Assumption 1.

Fix a sample size $n\geq 1$. The following conditions must be satisfied by the
resampling distribution $\mathcal{Q}_{n}$.

__

  * •

Symmetry:_for all_ $m\geq 1$, $i_{1},\ldots,i_{m}\in[n]$, and permutations
$\sigma\in\mathcal{S}_{n}$,

$\mathcal{Q}_{n}\\{(i_{1},\ldots,i_{m})\\}=\mathcal{Q}_{n}\\{(\sigma(i_{1}),\ldots,\sigma(i_{m}))\\}.$

  * •

Nontrivial subsampling:_it holds that_ $p_{\mathcal{Q}_{n}}<1$, where we
define

$p_{\mathcal{Q}_{n}}=\mathbb{E}_{r\sim\mathcal{Q}_{n}}\left[\frac{1}{n}\sum_{i=1}^{n}\mathbf{1}\\{i\in
r\\}\right],$

_the expected fraction of data points appearing in a bag_ $r$ sampled from
$\mathcal{Q}_{n}$.

  * •

Nonpositive covariance:__

$\textnormal{Cov}_{r\sim\mathcal{Q}_{n}}(\mathbf{1}_{i\in r},\mathbf{1}_{j\in
r})\leq 0\textnormal{ for any $i\neq j\in[n]$}.$

  * •

Compatibility with the leave-one-out distribution:_the corresponding
distribution_ $\mathcal{Q}_{n-1}$ on $\textnormal{seq}_{[n-1]}$ must equal the
distribution of $r\sim\mathcal{Q}_{n}$ conditional on the event $n\not\in r$.

Before proceeding, we comment on the last part of the assumption: why does
$\mathcal{Q}_{n-1}$ appear in our assumptions on $\mathcal{Q}_{n}$? To study
the stability of the bagged algorithm (Definition 4), it is not sufficient to
consider only $\mathcal{Q}_{n}$, the distribution on
$\textnormal{seq}_{[n]}$—when we compare to the result of the algorithm on the
leave-one-out data set $\mathcal{D}^{\setminus i}$, we now have sample size
$n-1$ (instead of $n$), and consequently we need to specify a distribution
$\mathcal{Q}_{n-1}$ on $\textnormal{seq}_{[n-1]}$, as well.

The four conditions of 2 are all satisfied by the standard resampling schemes
described above, with $p_{\mathcal{Q}_{n}}=1-\left(1-\frac{1}{n}\right)^{m}$
for bootstrapping (sampling $m$ out of $n$ points uniformly _with_
replacement), and $p_{\mathcal{Q}_{n}}=\frac{m}{n}$ for subbagging (sampling
$m$ out of $n$ points uniformly _without_ replacement).

##  4 Main results for a Hilbert space

We are now ready to present our stability guarantees for bagging. Our first
main result establishes that derandomized bagging automatically provides mean-
square stability for any algorithm with outputs in a bounded subset of a
Hilbert space. In this section, we work with stability with respect to the
Hilbert norm $\|\cdot\|_{\mathcal{H}}$. For instance,

  * •

If $\mathcal{W}\subset\mathbb{R}^{d}$, then the Hilbert norm may be chosen to
be the Euclidean norm $\|\cdot\|_{2}$, or more generally, a Mahalanobis
distance $\|w\|_{\mathcal{H}}=\sqrt{w^{\top}Sw}$ for a fixed positive definite
matrix $S$.

  * •

If $\mathcal{W}\subset L_{2}(\mathbb{R}^{d})$, the space of square-integrable
functions $\mathbb{R}^{d}\to\mathbb{R}$, then the Hilbert norm is given by the
$L_{2}$ norm on functions,
$\|w\|_{\mathcal{H}}=(\int_{\mathbb{R}^{d}}w(x)^{2}\;\textnormal{d}x)^{1/2}$.
Alternatively, we can take a Sobolev space of functions
$\mathbb{R}^{d}\to\mathbb{R}$ with smoothness parameter $s$ with the Hilbert
norm being the Sobolev norm $\|\cdot\|_{s,2}$ (defined formally later on in
(3)—see Experiment 3 in Section 6.1 for details).

######  Theorem 6.

Suppose $\mathcal{A}$ takes values in a convex and closed subset
$\mathcal{W}\subset\mathcal{H}$, where $\mathcal{H}$ is a separable Hilbert
space, and define

$\mathsf{rad}_{\mathcal{H}}(\mathcal{W})=\inf_{w\in\mathcal{W}}\sup_{w^{\prime}\in\mathcal{W}}\|w-w^{\prime}\|_{\mathcal{H}}$

as the radius of the set $\mathcal{W}$. Fix a resampling distribution
$\mathcal{Q}_{n}$ on $\textnormal{seq}_{[n]}$, which satisfies 1. Then
derandomized bagging, run with base algorithm $\mathcal{A}$ and resampling
distribution $\mathcal{Q}_{n}$, has mean-square stability $\beta^{2}$ with
respect to $\|\cdot\|_{\mathcal{H}}$, with

$\displaystyle\beta^{2}:=\frac{\mathsf{rad}^{2}_{\mathcal{H}}(\mathcal{W})}{n-1}\cdot\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}.$
(1)

Combining this result with 3 immediately yields a guarantee on
$(\varepsilon,\delta)$-stability.

######  Corollary 7.

Under the notation and assumptions of Theorem 6, derandomized bagging has tail
stability $(\varepsilon,\delta)$ with respect to $\|\cdot\|_{\mathcal{H}}$ for
any $\varepsilon,\delta\geq 0$ satisfying

$\displaystyle\delta\varepsilon^{2}\geq\frac{\mathsf{rad}_{\mathcal{H}}^{2}(\mathcal{W})}{n-1}\cdot\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}.$
(2)

##### Comparing to prior bounds

The main result of (Soloff et al., 2024, Theorem 8), ensuring the stability of
algorithms that return a real-valued output, can be recovered as a special
case of Corollary 7 above, by taking
$\mathcal{W}=[0,1]\subseteq\mathcal{H}=\mathbb{R}$ (note that we then have
$\mathsf{rad}_{\mathcal{H}}(\mathcal{W})=\frac{1}{2}$). But remarkably, the
result of Theorem 6 shows that there is no price to pay for high-
dimensionality (or even for infinite dimensionality): the stability guarantee
is _exactly the same_ for any Hilbert space, regardless of dimension, as long
as the radius of the output space $\mathcal{W}$ is bounded.

###  4.1 Results for a finite number of bags

The results above apply to the derandomized bagged version of $\mathcal{A}$,
as constructed in Definition 5. In practice, of course, it is not feasible to
calculate this output $\hat{w}$ since we would need to average over all
possible bags $r$, so we would instead use the bagged estimator $\hat{w}$
constructed in Definition 4 for some finite number of bags $B$. Since we are
simply taking a Monte Carlo approximation to the derandomized bagged estimator
of Definition 5, we would expect to incur error on the order of
$\frac{1}{\sqrt{B}}$. To formalize this, we first recall a generalization of
the Azuma-Hoeffding inequality:

######  Theorem 8 ((Hayes, 2005)).

Let $\hat{w}_{1},\ldots,\hat{w}_{B}\in\mathcal{H}$ be independent random
variables in a Hilbert space $\mathcal{H}$ such that
$\|\hat{w}_{b}-\mathbb{E}[\hat{w}_{b}]\|_{\mathcal{H}}\leq C$ almost surely
for $b=1,\ldots,B$. Then for any $\delta>0$,

$\mathbb{P}\left\\{\bigg{\|}\frac{1}{B}\sum_{b=1}^{B}(\hat{w}_{b}-\mathbb{E}[\hat{w}_{b}])\bigg{\|}_{\mathcal{H}}\leq
C\sqrt{\frac{1}{2B}\log\left(\frac{2e^{2}}{\delta}\right)}\right\\}\geq
1-\delta.$

In particular,

$\mathbb{E}\left[\bigg{\|}\frac{1}{B}\sum_{b=1}^{B}(\hat{w}_{b}-\mathbb{E}[\hat{w}_{b}])\bigg{\|}_{\mathcal{H}}^{2}\right]\leq\frac{C^{2}e^{2}}{B}.$

Note that, aside from the additional factor $e^{2}$, this result is the same
as what we would obtain for the case $\mathcal{H}=\mathbb{R}$. Combining this
result with Theorem 6 yields the following stability guarantee for bagged
algorithms with a finite number of bags $B$.

######  Corollary 9.

Under the notation and assumptions of Theorem 6, let $B\geq 1$ be fixed. Then
bagging, run with base algorithm $\mathcal{A}$ and resampling distribution
$\mathcal{Q}_{n}$ and with $B$ bags, has mean-square stability $\beta^{2}$
with respect to $\|\cdot\|_{\mathcal{H}}$, with

$\beta^{2}:=\mathsf{rad}^{2}_{\mathcal{H}}(\mathcal{W})\left(\frac{1}{n-1}\cdot\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}+\frac{16e^{2}}{B}\right).$

Of course, tail stability follows as an immediate consequence of 3.

##  5 Extension: results outside the Hilbert space setting

We now turn to the question of stability in a more general setting, where the
output space is not necessarily a Hilbert space. We assume that $\mathcal{A}$
returns outputs lying in $\mathcal{W}\subset\mathcal{B}$, where $\mathcal{B}$
is a separable Banach space equipped with norm $\|\cdot\|$, and we now aim to
establish stability with respect to this norm $\|\cdot\|$ directly.

As a motivating example, consider a setting where $\mathcal{A}$ returns an
output $\hat{w}$ that specifies a distribution—for instance, this may arise in
Bayesian inference (Problem Setting 2 from Section 1.1), where $\mathcal{A}$
returns a posterior distribution over some parameter of interest. In this type
of setting, since $\hat{w}$ represents a distribution, it is natural to study
its stability with respect to $\mathsf{d}_{\mathsf{TV}}$, the total variation
distance. While $\mathsf{d}_{\mathsf{TV}}$ defines a valid norm, this is not a
Hilbert norm and so the theoretical guarantees of Section 4 cannot be directly
applied.

To allow for meaningful stability guarantees in a broader range of settings,
then, we now consider stability in this more general context, where
$\|\cdot\|$ is not necessarily a Hilbert norm.

###  5.1 Theoretical guarantee

Unlike the special case of a Hilbert space, in this more general setting we
need an additional assumption to ensure that stability will hold.

######  Assumption 2.

There exists a constant $R>0$ such that
$\sup_{w,w^{\prime}\in\mathcal{W}}\|w-w^{\prime}\|\leq R$, and moreover, there
exist some fixed elements $v_{1},\dots,v_{K}\in\mathcal{B}$ with unit norm,
$\|v_{1}\|=\dots=\|v_{K}\|=1$, and a constant $\rho>0$ such that for any
$w,w^{\prime}\in\mathcal{W}$, the difference $w-w^{\prime}$ can be
approximated by a linear combination of the $v_{j}$’s,

$\Big{\|}(w-w^{\prime})-\sum_{i=1}^{K}t_{i}v_{i}\Big{\|}\leq\rho$

for some $t=(t_{1},\dots,t_{K})\in\mathbb{R}^{K}$ with $\|t\|_{1}\leq R$.

Here $R$ is playing a role analogous to the radius,
$\mathsf{rad}_{\mathcal{H}}(\mathcal{W})$, in the Hilbert space setting. For
instance, returning to our earlier example, if $\mathcal{W}=\Delta_{K-1}$ (and
the norm is given by $\mathsf{d}_{\mathsf{TV}}$) then we can take $v_{j}$ to
be the $j$th canonical basis vector. Then 2 holds, with this same value of
$K$, and with $R=1$ and $\rho=0$.

######  Theorem 10.

Suppose $\mathcal{A}$ takes values in a convex and closed subset
$\mathcal{W}\subset\mathcal{B}$, where $\mathcal{B}$ is a separable Banach
space equipped with norm $\|\cdot\|$, and where $\mathcal{W}$ satisfies 2 for
some integer $K\geq 1$ and some $R,\rho\geq 0$. Fix a resampling distribution
$\mathcal{Q}_{n}$ on $\textnormal{seq}_{[n]}$, which satisfies 1. Then
derandomized bagging, run with base algorithm $\mathcal{A}$ and resampling
distribution $\mathcal{Q}_{n}$, has mean-square stability $\beta^{2}$ with
respect to $\|\cdot\|$, with

$\beta^{2}:=1.6\left[R\sqrt{\frac{2(1+\eta_{n})\log(2K)}{n}\left(\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}+\frac{4\eta_{n}}{(1-p_{\mathcal{Q}_{n}})^{2}}\right)}+R\sqrt{\frac{1}{n}\cdot\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}}+2\rho\cdot
p_{\mathcal{Q}_{n}}\right]^{2},$

where $\eta_{n}=\frac{H_{n}-1}{n-H_{n}}\asymp\frac{\log n}{n}$, for
$H_{n}=1+\frac{1}{2}+\dots+\frac{1}{n}$ denoting the $n$th harmonic number.

As in the Hilbert space setting, we can also extend this theorem to verify
tail stability and to prove a stability result for bagging with a finite
$B$—see Appendix B for details.

###  5.2 Examples

Now we return to Problem Setting 2, as discussed in the motivation for the
extension to Banach spaces, where $\mathcal{A}$ returns a distribution (e.g.,
a posterior), to see what the new theorem implies in this setting. To make
this more concrete, we can consider two specific scenarios: a discrete setting
and a continuous setting. We will see how, in each scenario, the new result of
Theorem 10 gives a much stronger guarantee than what can be derived from the
results of Section 4. For intuition, in this discussion, we treat
$\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}$ as a constant.

##### Problem Setting 2(a): the discrete case

First, consider a finite setting, where $\mathcal{A}$ returns a distribution
on $K$ many values—for instance, if the output of $\mathcal{A}$ specifies a
Bayesian posterior over $K$ many possible models. In this scenario, $\hat{w}$
takes values in the simplex $\Delta_{K-1}\subset\mathbb{R}^{K}$, and so the
total variation distance is given by
$\mathsf{d}_{\mathsf{TV}}(w,w^{\prime})=\frac{1}{2}\|w-w^{\prime}\|_{1}$.

Of course, since $\mathbb{R}^{K}$ is itself a Hilbert space (with norm
$\|\cdot\|_{\mathcal{H}}=\|\cdot\|_{2}$), we can apply our results from
Section 4—but this will not give a very useful bound for the stability of
derandomized bagging. Specifically, Theorem 6 tells us that mean-square
stability with respect to $\|\cdot\|_{2}$ holds, for $\beta^{2}\propto 1/n$.
Since $\|\cdot\|_{1}\leq\sqrt{K}\|\cdot\|_{2}$ (the usual inequality relating
$\ell_{1}$ norm to the $\ell_{2}$ norm in $\mathbb{R}^{K}$), this immediately
implies that

$\frac{1}{n}\sum_{i\in[n]}\mathbb{E}\left[\mathsf{d}_{\mathsf{TV}}^{2}\left(\hat{w},\hat{w}^{\setminus
i}\right)\right]\leq\frac{K}{4}\cdot\frac{1}{n}\sum_{i\in[n]}\mathbb{E}[\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{2}]\lesssim K/n,$

or in other words, mean-square stability with respect to
$\mathsf{d}_{\mathsf{TV}}$ holds, for $\beta^{2}\propto K/n$. But if the
dimension $K$ is large relative to $n$, then this is not a very useful
guarantee.

In contrast, let us now see what can be guaranteed via Theorem 10. As
mentioned above, 2 holds with $K$ equal to the dimension, and with $R=1$ and
$\rho=0$. Then Theorem 10 tells us that mean-square stability with respect to
$\|\cdot\|$ holds for $\beta^{2}\propto\frac{\log K}{n}$—a far better bound.

##### Problem Setting 2(b): the continuous case

Next, consider a continuous setting, where $\mathcal{A}$ returns a density on
$[0,1]$, given by some function $\hat{w}:[0,1]\rightarrow[0,\infty)$
satisfying $\int_{t=0}^{1}\hat{w}(t)\;\textnormal{d}t=1$. In this case we have
$\mathsf{d}_{\mathsf{TV}}(w,w^{\prime})=\frac{1}{2}\int_{t=0}^{1}|w(t)-w^{\prime}(t)|\;\textnormal{d}t$.

To formalize this setting, consider the Banach space

$\mathcal{B}=L_{1}([0,1]):=\left\\{f:[0,1]\to\mathbb{R}:\int_{t=0}^{1}|f(t)|\;\textnormal{d}t<\infty\right\\},$

equipped with the norm
$\|f\|=\frac{1}{2}\int_{t=0}^{1}|f(t)|\;\textnormal{d}t$ (note that this
differs from the usual $\ell_{1}$ norm by the factor $\frac{1}{2}$, in order
to agree with total variation distance). Suppose that our algorithm returns
densities that are further constrained to be $L$-Lipschitz, so that we have

$\mathcal{W}=\left\\{w:[0,1]\to[0,\infty):\int_{t=0}^{1}w(t)=1,\ \sup_{t\neq
t^{\prime}}\frac{|w(t)-w(t^{\prime})|}{|t-t^{\prime}|}\leq L\right\\}.$

As for Problem Setting 2(a), we again begin by asking whether the Hilbert
space result of Theorem 6 might already be sufficient in this setting. Indeed,
while $\mathcal{B}=L_{1}([0,1])$ is not a Hilbert space, the Lipschitz
constraint on $\mathcal{W}$ means that
$\int_{t=0}^{1}|w(t)|^{2}\;\textnormal{d}t<\infty$ for all $w\in\mathcal{W}$,
meaning that $\mathcal{W}$ can be viewed as a subset of the Hilbert space
$L_{2}([0,1])$, equipped with norm
$\|f\|_{\mathcal{H}}=\|f\|_{L_{2}}=(\int_{t=0}^{1}|f(t)|^{2}\;\textnormal{d}t)^{1/2}$.
A straightforward calculation verifies that
$\mathsf{rad}_{\mathcal{H}}(\mathcal{W})\propto L^{1/4}$, meaning that by
Theorem 6, derandomized bagging satisfies mean-square stability with respect
to $\|\cdot\|_{L_{2}}$ for $\beta^{2}\propto L^{1/2}/n$. Since for this
function space we have $\|\cdot\|_{L_{1}}\leq\|\cdot\|_{L_{2}}$, this
immediately implies that mean-square stability with respect to
$\mathsf{d}_{\mathsf{TV}}$ also holds with parameter $\beta^{2}\propto
L^{1/2}/n$. However, if the Lipschitz constant $L$ is large, this may not be a
satisfying result.

Now we turn to our new result, Theorem 10, to see if the stability guarantee
can be improved. First, observe that 2 holds for any $K\geq 1$ if we take
$R=1$ and $\rho\propto L/K$. This holds because we can take $v_{j}$ to be the
function $v_{j}(t)=K\cdot\mathbf{1}\\{\frac{j-1}{K}\leq t\leq\frac{j}{K}\\}$,
for each $j\in[K]$—then for any $L$-Lipschitz densities $w,w^{\prime}$, the
$2L$-Lipschitz difference $w-w^{\prime}$ can be approximated with a piecewise
constant function, i.e., a linear combination of these $v_{j}$’s. Choosing
$K\propto Ln$, we apply Theorem 10, which then leads to a guarantee on the
mean-square stability with respect to $\mathsf{d}_{\mathsf{TV}}$, with
$\beta^{2}\propto\frac{\log(Ln)}{n}$. Again, as for the discrete case, if $L$
is large (i.e., $L\gg\log^{2}n$) then this is a much stronger bound than the
one we can obtain via the Hilbert space result.

##### Summary of examples

The following table summarizes our results for the two settings, in terms of
the stability guarantee with respect to $\mathsf{d}_{\mathsf{TV}}$. The
improved scaling in the right column highlights the utility of Theorem 10 in
extending our results to a broader range of settings.

Output space $\mathcal{W}$ |  | Stability guarantee for $\mathsf{d}_{\mathsf{TV}}$
---
via Theorem 6
(Hilbert space)
| Stability guarantee for $\mathsf{d}_{\mathsf{TV}}$
---
via Theorem 10
(Banach space)
| Problem Setting 2(a):
---
distributions on $\\{1,\dots,K\\}$
$\beta^{2}\propto\frac{K}{n}$ | $\beta^{2}\propto\frac{\log K}{n}$
| Problem Setting 2(b):
---
$L$-Lipschitz densities on $[0,1]$
$\beta^{2}\propto\frac{L^{1/2}}{n}$ | $\beta^{2}\propto\frac{\log(Ln)}{n}$

######  Remark 11 (On the necessity of 2).

Our stability guarantee for general Banach spaces essentially requires that
$\mathcal{W}$ lies (approximately) in a convex hull of $2K$ many
points—namely, $Rv_{1},-Rv_{1}$, …, $Rv_{K},-Rv_{K}$. In particular, this
implies that $\mathcal{W}$ is (approximately) contained in a _finite-
dimensional_ subspace of $\mathcal{B}$. In contrast, our Hilbert space results
apply regardless of dimensionality—while we do require $\mathcal{W}$ to be
bounded via its radius, $\mathsf{rad}_{\mathcal{H}}(\mathcal{W})$, there is no
implicit cost to the dimension of $\mathcal{H}$ in the results of Theorem 6.
It turns out that some assumption of this form is actually required for the
Banach case—see 13 (in Section B.2) for an explicit counterexample verifying
stability fails without this type of assumption.

##  6 Experiments

In this section, we study the stability of subbagging in experiments across
several simulation and real data settings. Code to reproduce all experiments
is available at

https://github.com/jake-soloff/stability-resampling-experiments.

###  6.1 Data and Methods

We start by presenting an overview of our four experiments. Experiments 1, 2,
and 3 each examine stability in a Hilbert space, with respect to the
appropriate Hilbert norm $\|\cdot\|_{\mathcal{H}}$ (as studied in Section 4).
In contrast, the setting of Experiment 4 is a Banach space (as in Section 5).

##### Experiment 1 (Estimating a regression function using bagged decision
trees)

We draw data from the following data-generating process:

$\displaystyle(X_{i},\alpha_{i},\gamma_{i})$
$\displaystyle\stackrel{{\scriptstyle\text{iid}}}{{\sim}}\textnormal{Uniform}([0,1]^{d})\times\textnormal{Uniform}([-.25,.25])\times\textnormal{Uniform}([0,1]),$
$\displaystyle\tilde{Y}_{i}$
$\displaystyle=\sum_{j=1}^{d}\sin\left(\frac{X_{ij}}{j}\right)+\alpha_{i}\mathbf{1}\left\\{i=1\
(\mathrm{mod}\ 3)\right\\}+\gamma_{i}\mathbf{1}\left\\{i=1\ (\mathrm{mod}\
4)\right\\},$

with $n=500$ and $d=40$. We then set
$Y_{i}=\frac{\tilde{Y}_{i}-\min_{j}\tilde{Y}_{j}}{\max_{j}\tilde{Y}_{j}-\min_{j}\tilde{Y}_{j}}$
so that each response $Y_{i}$ is in the unit interval $[0,1]$. Note that the
algorithm has access to the observed data
$\mathcal{D}=(X_{i},Y_{i})_{i\in[n]}$, i.e., $\alpha_{i}$ and $\gamma_{i}$ are
latent variables used only to generate the data $\mathcal{D}$. We apply
`sklearn.tree.DecisionTreeRegressor` to train the regression trees, setting
`max_depth=50` and leaving all other parameters at their default values. We
evaluate the stability of the learned function in the $L_{2}$ norm
$\|\hat{w}\|_{L_{2}}=(\int_{[0,1]^{d}}\hat{w}(x)^{2}\textnormal{d}x)^{1/2}$.

##### Experiment 2 (Synthetic controls real data example)

For this experiment, we reproduce Abadie et al. (2015)’s SCM analysis on the
1990 reunification of Germany, as discussed above in Section 1.1 (see Problem
Setting 3). The data contains $n=16$ countries (aside from Germany), with $6$
predictors measured for each country (pre-treatment averages of various
macroeconomic factors such as inflation and trade openness—see Table 1 of
Abadie et al. (2015) for more details). SCM then estimates a weighted
combination of these $n$ countries that represents a ‘synthetic West Germany’,
meaning that the output of the algorithm is a vector $\hat{w}\in\Delta_{n-1}$.
When we run SCM on a subset of the $n=16$ countries, we assign a weight of $0$
to each omitted country. We run SCM using the Python package pysyncon and make
extensive use of their reproduction of the original study Abadie et al.
(2015). We evaluate the stability of these weights in the Euclidean norm
$\|\cdot\|_{2}$.

##### Experiment 3 (Least-squares spectrum analysis)

We study function estimation in a Sobolev space, corresponding to the Hilbert
space $W^{s,2}(\mathbb{R})$ with smoothness parameter $s=2$, with the
corresponding norm

$\displaystyle\|f\|_{s,2}:=\sqrt{\sum_{k=-\infty}^{\infty}(1+k^{2})^{s}\cdot|\hat{f}(k)|^{2}},$
(3)

where $\hat{f}(k)$ denotes the $k$th Fourier coefficient, indexed by
$k\in\mathbb{Z}$. We wish to recover $f$, which is equivalent to estimating
the Fourier coefficients. We restrict our estimation problem to _even
functions_ , so that we may constrain the Fourier coefficients to be real
valued.

Least squares spectral analysis (LSSA) is one method for estimating the
Fourier coefficients of a function from irregularly spaced samples Vaníček
(1969). Given a data set $\mathcal{D}=(X_{i},Y_{i})_{i\in[n]}$, we first
construct a countable set of features $\Psi_{ik}=\exp\left(jkX_{i}\right)$,
where $j=\sqrt{-1}$, $i\in[n]$ and $k\in\mathbb{Z}$. The base algorithm
$\mathcal{A}(\mathcal{D})$ returns a function $\hat{w}:[-1,1]\to\mathbb{R}$,
given by

$\hat{w}(x)=\sum_{k=-\infty}^{\infty}\hat{a}_{k}\exp(jkx),\qquad\textnormal{with}\qquad\hat{a}=\mathop{\rm
argmin}_{a\in\mathcal{C}_{n}(R)}\sum_{i=1}^{n}\left(Y_{i}-\sum_{k=-\infty}^{\infty}\Psi_{ik}a_{k}\right)^{2},$

and where the constraint set is given by

$\mathcal{C}_{n}(R):=\left\\{(\ldots,a_{-1},a_{0},a_{1},\ldots):a_{k}\in\mathbb{R},a_{k}=0\,\forall\,|k|>\lceil{n/2}\rceil,\sum_{k=-\infty}^{\infty}(1+k^{2})^{2}a_{k}^{2}\leq
R^{2}\right\\}.$

For our simulation, we sample $n=200$ random variables
$X_{1},\ldots,X_{n}\in[-1,1]$ iid from a non-uniform density
$g_{X}(x)=-\frac{1}{2}\log|x|$. Then, the corresponding responses are given by
$Y_{i}=f^{*}(X_{i})+\zeta_{i}$, where
$\zeta_{i}\stackrel{{\scriptstyle\text{iid}}}{{\sim}}\mathcal{N}(0,.01)$ and
$f^{*}(x)=\sin\left(\frac{2}{1-x}+\frac{2}{1+x}\right)$. We run the base
algorithm above using the Python package cvxpy, setting $R=10$. We evaluate
the stability of the learned function in the Sobolev norm $\|\cdot\|_{2,2}$
(see Equation 3 above).

##### Experiment 4 (Stability of the softmax function)

We simulate a binary $n\times d$ data matrix $X\in\\{0,1\\}^{n\times d}$ using
$X_{ij}\stackrel{{\scriptstyle\text{iid}}}{{\sim}}\textnormal{Bern}(0.2)$,
where $n=2000$ and $d=100$. The base algorithm takes the softmax of the column
sums
$\mathcal{A}(\mathcal{D}):=\sigma(\mathbf{1}_{n}^{\top}X)\in\Delta_{d-1}$,
i.e.

$\mathcal{A}(\mathcal{D})_{j}=\frac{\exp\left(\sum_{i=1}^{n}X_{ij}\right)}{\sum_{k=1}^{d}\exp\left(\sum_{i=1}^{n}X_{ik}\right)}.$

We evaluate stability with respect to the total variation distance
$\mathsf{d}_{\mathsf{TV}}$. (As discussed in Section 5, the
$\mathsf{d}_{\mathsf{TV}}$ norm is not a Hilbert norm—that is, unlike the
first three experiments, stability guarantees for Experiment 4 hold only due
to the extension to the Banach space setting given in Theorem 10.)

Experiment 1: Regression trees  Experiment 2: Synthetic controls, German
reunification example (Abadie et al., 2015) Experiment 3: Least-squares
spectral analysis  Experiment 4: Softmax

Figure 1: Simulation results comparing the stability of subbagging to that of
the corresponding base algorithm under the given norm. Left: Histogram of
leave-one-out perturbations. Right: for each $\varepsilon$, the smallest
$\delta$ such that the algorithm has tail stability $(\varepsilon,\delta)$ in
the sense of Definition 2. Higher curves thus represent greater instability.
In the right panels we also plot the stability guarantee of Theorem 6 for
Experiments 1–3 and Theorem 10 for Experiment 4 (the guarantees for
derandomized bagging). For all four experiments, $m=n/2$ and $B=10,000$.

###  6.2 Results

Our experiments compare the stability of the base algorithm $\mathcal{A}$ to
the stability of subbagging with $m=n/2$ and $B=10,000$. Our results are shown
in Figure 1 and Table 2. The left panels of Figure 1 show the histogram of
leave-one-out perturbations $\left\|\hat{w}-\hat{w}^{\setminus i}\right\|$ for
$i\in\left\\{1,\ldots,n\right\\}$. In the right panels of Figure 1, for a
fixed data set and algorithm, we measure stability by plotting, for each value
of $\varepsilon$, the smallest value of $\delta$ such that the algorithm is
$\left(\varepsilon,\delta\right)$-stable:

$\delta=\frac{1}{n}\sum_{i=1}^{n}1_{\|\hat{w}-\hat{w}^{\setminus
i}\|\geq\varepsilon}.$

Furthermore, in Table 2, we record the mean-square stability

$\beta^{2}:=\frac{1}{n}\sum_{i=1}^{n}\|\hat{w}-\hat{w}^{\setminus i}\|^{2}.$

For all four experiments, we observe that the subbagged algorithm has better
stability than the base algorithm, both in terms of the mean-square stability
parameter $\beta^{2}$ and in terms of the tail stability parameters (across
all, or nearly all, values of $\varepsilon$).

Setting |  | Mean-square stability $\beta^{2}$
---
for base algorithm
| Mean-square stability $\beta^{2}$
---
for subbagged algorithm
Regression trees | $0.0107$ | $0.0000047$
Synthetic controls (SCM) | $0.0532$ | $0.0125$
Spectral analysis (LSSA) | $40.033$ | $0.0793$
Softmax | $0.0186$ | $0.0006$
Table 2: Results for mean-square stability.

##  7 Conclusion

In supervised learning, algorithmic stability plays a fundamental role in
generalization theory (Bousquet and Elisseeff, 2002; Mukherjee et al., 2006;
Shalev-Shwartz et al., 2010). More broadly, various forms of stability are
fundamental across applied statistics, in part for their role in establishing
reproducibility and interpretability (Yu and Kumbier, 2020; Yu and Barter,
2024). This work provides general, finite sample theoretical guarantees
supporting the use of bagging and related techniques as general-purpose
stabilizers. Our stability guarantees represent a major generalization of
Soloff et al. (2024)—which studied algorithms with outputs on the real line—as
our results apply to infinite dimensional settings.

An important remaining open question is how to define and guarantee stability
for algorithms with discrete outputs. In the synthetic controls setting, for
example, our theory shows that bagging can stabilize the SCM _weights_.
However, our results do not imply stability in the resulting _support_ —that
is, the set of control units that receive nonzero support. In a sense, the
goal of sparsity is incompatible with stability (Xu et al., 2011), so we need
a new framework to reconcile these two desiderata. We leave these questions
for future work.

### 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.

##  Appendix A Proofs for the Hilbert space setting

###### Proof of Theorem 6.

For each $r\in\textnormal{seq}_{[n]}$, define
$\hat{w}^{r}=\mathbb{E}_{\xi}[\mathcal{A}(\mathcal{D}^{r};\xi)]$, the expected
output of the randomized algorithm for the data set given by
$\mathcal{D}^{r}$. From this point on, we write $\mathbb{E}[\cdot]$ to denote
expectation with respect to $r\sim\mathcal{Q}_{n}$. We then have
$\mathbb{E}[\hat{w}^{r}]=\hat{w}$ and

$\hat{w}^{\setminus i}=\mathbb{E}\left[\hat{w}^{r}\mid i\not\in
r\right]=\frac{1}{1-p_{\mathcal{Q}_{n}}}\mathbb{E}\left[\hat{w}^{r}\cdot\mathbf{1}\\{i\not\in
r\\}\right]$

by construction, and therefore, a straightforward calculation shows that

$\hat{w}-\hat{w}^{\setminus
i}=\frac{1}{1-p_{\mathcal{Q}_{n}}}\mathbb{E}\left[(\hat{w}^{r}-\hat{w})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right].$ (4)

Therefore,

$\displaystyle\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{\mathcal{H}}$
$\displaystyle=\sum_{i\in[n]}\left\langle\hat{w}-\hat{w}^{\setminus
i},\frac{1}{1-p_{\mathcal{Q}_{n}}}\mathbb{E}\left[(\hat{w}^{r}-\hat{w})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right]\right\rangle$
$\displaystyle=\frac{1}{1-p_{\mathcal{Q}_{n}}}\mathbb{E}\left[\left\langle\sum_{i\in[n]}(\hat{w}-\hat{w}^{\setminus
i})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right),\hat{w}^{r}-\hat{w}\right\rangle\right]$
$\displaystyle\leq\frac{1}{1-p_{\mathcal{Q}_{n}}}\mathbb{E}\left[\left\|\sum_{i\in[n]}(\hat{w}-\hat{w}^{\setminus
i})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right\|_{\mathcal{H}}^{2}\right]^{1/2}\cdot\mathbb{E}\left[\left\|\hat{w}^{r}-\hat{w}\right\|_{\mathcal{H}}^{2}\right]^{1/2}$
$\displaystyle\leq\frac{\mathsf{rad}_{\mathcal{H}}(\mathcal{W})}{1-p_{\mathcal{Q}_{n}}}\mathbb{E}\left[\left\|\sum_{i\in[n]}(\hat{w}-\hat{w}^{\setminus
i})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right\|_{\mathcal{H}}^{2}\right]^{1/2}.$

The last step holds because we have

$\mathbb{E}\left[\|\hat{w}^{r}-\hat{w}\|_{\mathcal{H}}^{2}\right]=\inf_{w\in\mathcal{W}}\mathbb{E}\left[\|\hat{w}^{r}-w\|_{\mathcal{H}}^{2}\right]\leq\inf_{w\in\mathcal{W}}\sup_{w^{\prime}\in\mathcal{W}}\|w^{\prime}-w\|_{\mathcal{H}}^{2}=\mathsf{rad}^{2}_{\mathcal{H}}(\mathcal{W}),$

where the first equality holds since $\hat{w}=\mathbb{E}[\hat{w}^{r}]$
minimizes mean-square error.

To calculate the last remaining expected value, we first note that since we
have assumed symmetry and nonpositive covariance for $\mathcal{Q}_{n}$ (see
1), we have $\textnormal{Cov}\left(\mathbf{1}\\{i\in r\\},\mathbf{1}\\{j\in
r\\}\right)=-q$ for all $i\neq j$, for some $q\geq 0$. Moreover, since the
covariance matrix of the random vector $(\mathbf{1}\\{i\in r\\})_{i\in[n]}$
must be positive semidefinite, we must have
$q\leq\frac{p_{\mathcal{Q}_{n}}(1-p_{\mathcal{Q}_{n}})}{n-1}$. We then have

$\displaystyle\mathbb{E}\left[\left\|\sum_{i\in[n]}(\hat{w}-\hat{w}^{\setminus
i})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right\|_{\mathcal{H}}^{2}\right]$
$\displaystyle=\sum_{i\in[n]}\textnormal{Var}(\mathbf{1}\\{i\in
r\\})\cdot\|\hat{w}-\hat{w}^{\setminus i}\|^{2}_{\mathcal{H}}+\sum_{i\neq
j\in[n]}\textnormal{Cov}\left(\mathbf{1}\\{i\in r\\},\mathbf{1}\\{j\in
r\\}\right)\cdot\langle\hat{w}-\hat{w}^{\setminus
i},\hat{w}-\hat{w}^{\setminus j}\rangle$
$\displaystyle=\sum_{i\in[n]}p_{\mathcal{Q}_{n}}(1-p_{\mathcal{Q}_{n}})\cdot\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{\mathcal{H}}+\sum_{i\neq
j\in[n]}(-q)\cdot\langle\hat{w}-\hat{w}^{\setminus
i},\hat{w}-\hat{w}^{\setminus j}\rangle$
$\displaystyle=\left(p_{\mathcal{Q}_{n}}(1-p_{\mathcal{Q}_{n}})+q\right)\cdot\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{\mathcal{H}}-q\cdot\left\|\sum_{i\in[n]}\hat{w}-\hat{w}^{\setminus
i}\right\|^{2}_{\mathcal{H}}$ $\displaystyle\leq\frac{n}{n-1}\cdot
p_{\mathcal{Q}_{n}}(1-p_{\mathcal{Q}_{n}})\cdot\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{\mathcal{H}}.$

Combining everything, then, we have

$\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{\mathcal{H}}\leq\frac{\mathsf{rad}_{\mathcal{H}}(\mathcal{W})}{1-p_{\mathcal{Q}_{n}}}\left(\frac{n}{n-1}\cdot
p_{\mathcal{Q}_{n}}(1-p_{\mathcal{Q}_{n}})\cdot\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{\mathcal{H}}\right)^{1/2},$

which yields

$\frac{1}{n}\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{\mathcal{H}}\leq\frac{\mathsf{rad}^{2}_{\mathcal{H}}(\mathcal{W})}{n-1}\cdot\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}$

after rearranging terms. ∎

###### Proof of Corollary 9.

Let $\hat{w}$ and $\hat{w}^{\setminus i}$ be the output of the bagged
algorithm (as in Definition 4), and let $\tilde{w}$ and $\tilde{w}^{\setminus
i}$ be the corresponding derandomized bagged estimators (as in Definition 5).
For each $i\in[n]$, we have

$\mathbb{E}\left[\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{\mathcal{H}}\right]=\mathbb{E}\left[\|(\hat{w}-\tilde{w})-(\hat{w}^{\setminus
i}-\tilde{w}^{\setminus i})+(\tilde{w}-\tilde{w}^{\setminus
i})\|^{2}_{\mathcal{H}}\right]\\\
=\mathbb{E}\left[\|(\hat{w}-\tilde{w})-(\hat{w}^{\setminus
i}-\tilde{w}^{\setminus
i})\|^{2}_{\mathcal{H}}\right]+\|\tilde{w}-\tilde{w}^{\setminus
i}\|^{2}_{\mathcal{H}},$

since $\tilde{w}-\tilde{w}^{\setminus i}$ is nonrandom, while
$(\hat{w}-\tilde{w})-(\hat{w}^{\setminus i}-\tilde{w}^{\setminus i})$ has mean
zero. Next,

$\mathbb{E}\left[\|(\hat{w}-\tilde{w})-(\hat{w}^{\setminus
i}-\tilde{w}^{\setminus i})\|^{2}_{\mathcal{H}}\right]\leq
2\mathbb{E}\left[\|\hat{w}-\tilde{w}\|^{2}_{\mathcal{H}}\right]+2\mathbb{E}\left[\|\hat{w}^{\setminus
i}-\tilde{w}^{\setminus i}\|^{2}_{\mathcal{H}}\right].$

We next have
$\hat{w}=\frac{1}{B}\sum_{b=1}^{B}\mathcal{A}(\mathcal{D}^{r_{b}};\xi_{b})$.
Note that the random variables
$\mathcal{A}(\mathcal{D}^{r_{b}};\xi_{b})\in\mathcal{H}$ are i.i.d., with mean
$\tilde{w}$, and with
$\|\mathcal{A}(\mathcal{D}^{r_{b}};\xi_{b})-\tilde{w}\|_{\mathcal{H}}\leq
2\mathsf{rad}_{\mathcal{H}}(\mathcal{W})$ almost surely. Applying Theorem 8,

$\mathbb{E}\left[\|\hat{w}-\tilde{w}\|^{2}_{\mathcal{H}}\right]\leq\frac{4e^{2}\mathsf{rad}^{2}_{\mathcal{H}}(\mathcal{W})}{B}.$

By an identical argument, the same holds for the leave-one-out outputs,
$\hat{w}^{\setminus i}$ and $\tilde{w}^{\setminus i}$. Therefore,

$\mathbb{E}\left[\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}_{\mathcal{H}}\right]\leq\|\tilde{w}-\tilde{w}^{\setminus
i}\|^{2}_{\mathcal{H}}+\frac{16e^{2}\mathsf{rad}^{2}_{\mathcal{H}}(\mathcal{W})}{B}.$

Finally, $\|\tilde{w}-\tilde{w}^{\setminus i}\|^{2}_{\mathcal{H}}$ is bounded
by Theorem 6. ∎

##  Appendix B Proofs and additional results for the Banach space setting

###  B.1 Proof of main theorem

###### Proof of Theorem 10.

Let $\|\cdot\|_{*}$ represent the dual norm to $\|\cdot\|$. First, for each
$i$, as a consequence of the Hahn–Banach theorem (see, e.g., Folland (1999),
Theorem 5.8), there is a linear operator $\hat{f}_{i}$ on $\mathcal{B}$ such
that $\|\hat{f}_{i}\|_{*}=1$, and

$\hat{f}_{i}(\hat{w}-\hat{w}^{\setminus i})=\|\hat{w}-\hat{w}^{\setminus
i}\|.$

We then have

$\displaystyle\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}=\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(\hat{w}-\hat{w}^{\setminus i})$
$\displaystyle=\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}\left(\frac{1}{1-p_{\mathcal{Q}_{n}}}\mathbb{E}\left[(\hat{w}^{r}-\hat{w})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right]\right)\textnormal{ as calculated
in~{}\eqref{eqn:w_hat_i_identity}}$
$\displaystyle=\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\frac{1}{1-p_{\mathcal{Q}_{n}}}\mathbb{E}\left[\hat{f}_{i}(\hat{w}^{r}-\hat{w})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right]\textnormal{ since $\hat{f}_{i}$ is
linear}$
$\displaystyle=\frac{1}{1-p_{\mathcal{Q}_{n}}}\cdot\mathbb{E}\left[\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(\hat{w}^{r}-\hat{w})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right].$

Next, by 2, for each $r\in\textnormal{seq}_{[n]}$ we can find some
$t^{r}=(t^{r}_{1},\dots,t^{r}_{K})\in\mathbb{R}^{K}$ with $\|t^{r}\|_{1}\leq
R$ and
$\Big{\|}(\hat{w}^{r}-\hat{w})-\sum_{j=1}^{K}t^{r}_{j}v_{j}\Big{\|}\leq\rho$.
Then, since $\hat{f}_{i}$ is linear and has unit dual norm
$\|\hat{f}_{i}\|_{*}=1$,
$\left|\hat{f}_{i}(\hat{w}^{r}-\hat{w})-\sum_{j=1}^{K}t^{r}_{j}\hat{f}_{i}(v_{j})\right|\leq\rho$
for each $i\in[n]$, and so

$\displaystyle\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus i}\|^{2}$
$\displaystyle=\frac{1}{1-p_{\mathcal{Q}_{n}}}\cdot\mathbb{E}\left[\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(\hat{w}^{r}-\hat{w})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right]$
$\displaystyle\leq\frac{1}{1-p_{\mathcal{Q}_{n}}}\cdot\mathbb{E}\left[\sum_{j\in[K]}t^{r}_{j}\cdot\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(v_{j})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right]$ $\displaystyle\hskip
72.26999pt+\frac{\rho}{1-p_{\mathcal{Q}_{n}}}\cdot\mathbb{E}\left[\sum_{i\in[n]}\left|\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right|\right]$
$\displaystyle\leq\frac{R}{1-p_{\mathcal{Q}_{n}}}\cdot\mathbb{E}\left[\max_{j\in[K]}\left|\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(v_{j})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right|\right]$ $\displaystyle\hskip
72.26999pt+\frac{\rho}{1-p_{\mathcal{Q}_{n}}}\cdot\mathbb{E}\left[\sum_{i\in[n]}\left|\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right|\right],$

where for the last step, we use the fact that $\|t^{r}\|_{1}\leq R$ for any
$r$. Since $\mathbb{E}[|\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}|]=2p_{\mathcal{Q}_{n}}(1-p_{\mathcal{Q}_{n}})$ for
each $i$, and $\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\leq\sqrt{n}\cdot\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}}$, we can rewrite this as

$\displaystyle\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}\leq\frac{R}{1-p_{\mathcal{Q}_{n}}}\cdot$
$\displaystyle\mathbb{E}\left[\max_{j\in[K]}\left|\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(v_{j})\cdot\left(\mathbf{1}\\{i\in
r\\}-p_{\mathcal{Q}_{n}}\right)\right|\right]$
$\displaystyle+2\rho\sqrt{n}\cdot
p_{\mathcal{Q}_{n}}\cdot\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}}.$

Next, define $\hat{p}^{r}=\frac{1}{n}\sum_{i\in[n]}\mathbf{1}\\{i\in r\\}$, so
that $\mathbb{E}[\hat{p}^{r}]=p_{\mathcal{Q}_{n}}$. Using 1 we can verify that
$\textnormal{Var}(\hat{p}^{r})\leq\frac{p_{\mathcal{Q}_{n}}(1-p_{\mathcal{Q}_{n}})}{n}$.
Then, since $|\hat{f}_{i}(v_{j})|\leq 1$ for all $i,j$,

$\mathbb{E}\left[\max_{j\in[K]}\left|\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(v_{j})\cdot\left(\hat{p}^{r}-p_{\mathcal{Q}_{n}}\right)\right|\right]=\mathbb{E}\left[\max_{j\in[K]}\left|\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(v_{j})\right|\cdot|\hat{p}^{r}-p_{\mathcal{Q}_{n}}|\right]\\\
\leq\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\mathbb{E}\left[|\hat{p}^{r}-p_{\mathcal{Q}_{n}}|\right]\leq\sqrt{n}\cdot\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}}\cdot\sqrt{\textnormal{Var}(\hat{p}^{r})}\\\
\leq\sqrt{p_{\mathcal{Q}_{n}}(1-p_{\mathcal{Q}_{n}})}\cdot\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}}.$

Combining with the work above, then,

$\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}\leq\frac{R}{1-p_{\mathcal{Q}_{n}}}\cdot\mathbb{E}\left[\max_{j\in[K]}\left|\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(v_{j})\cdot\left(\mathbf{1}\\{i\in
r\\}-\hat{p}^{r}\right)\right|\right]\\\
{}+R\sqrt{\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}}\cdot\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}}+2\rho\sqrt{n}\cdot
p_{\mathcal{Q}_{n}}\cdot\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}}.$

Now we bound the remaining expected value. Let $c^{(j)}\in\mathbb{R}^{n}$ have
entries $c^{(j)}_{i}=\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(v_{j})$, and observe that
$\|c^{(j)}\|^{2}_{2}\leq\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus i}\|^{2}$
and $\|c^{(j)}\|_{\infty}\leq\max_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus i}\|$,
for each $j$, because $|\hat{f}_{i}(v_{j})|\leq 1$. Moreover, for any $i$,
$\hat{w}-\hat{w}^{\setminus
i}=p_{\mathcal{Q}_{n}}(\hat{w}^{i}-\hat{w}^{\setminus i})$, where
$\hat{w}^{i}=\mathbb{E}[\hat{w}^{r}\mid i\in r]$, and so
$\|\hat{w}-\hat{w}^{\setminus i}\|\leq p_{\mathcal{Q}_{r}}R$ (by 2), and so
$\|c^{(j)}\|_{\infty}\leq p_{\mathcal{Q}_{n}}R$ for all $j$. Applying Lemma 12
(with vectors $c^{(1)},\dots,c^{(K)},-c^{(1)},\dots,-c^{(K)}$), then,

$\mathbb{E}\left[\max_{j\in[K]}\left|\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(v_{j})\cdot\left(\mathbf{1}\\{i\in
r\\}-\hat{p}^{r}\right)\right|\,\middle|\,\hat{p}^{r}\right]\\\
\leq\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}}\cdot\sqrt{2(1+\eta_{n})(\hat{p}^{r}(1-\hat{p}^{r})+4\eta_{n})\log(2K)}+\frac{2}{3}(1+\eta_{n})\log(2K)\cdot
p_{\mathcal{Q}_{n}}R.$

Finally, since $t\mapsto\sqrt{t(1-t)+c}$ is concave for $t\in[0,1]$ and for
fixed $c\geq 0$, marginalizing over $\hat{p}^{r}$, and applying Jensen’s
inequality, yields

$\displaystyle\mathbb{E}$
$\displaystyle\left[\max_{j\in[K]}\left|\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|\cdot\hat{f}_{i}(v_{j})\cdot\left(\mathbf{1}\\{i\in
r\\}-\hat{p}^{r}\right)\right|\right]$
$\displaystyle\leq\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}}\cdot\sqrt{2(1+\eta_{n})(p_{\mathcal{Q}_{n}}(1-p_{\mathcal{Q}_{n}})+4\eta_{n})\log(2K)}+\frac{2}{3}(1+\eta_{n})\log(2K)\cdot
p_{\mathcal{Q}_{n}}R.$

Combining all the above calculations proves that
$\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}\leq\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus i}\|^{2}}\cdot
A+B$, where

$A=R\left[\sqrt{2(1+\eta_{n})\left(\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}+\frac{4\eta_{n}}{(1-p_{\mathcal{Q}_{n}})^{2}}\right)\log(2K)}+\sqrt{\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}}\right]+2\rho\sqrt{n}\cdot
p_{\mathcal{Q}_{n}}$

and

$B=R^{2}\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}\cdot\frac{2}{3}(1+\eta_{n})\log(2K)\leq\frac{A^{2}}{3}.$

Then, solving the quadratic,

$\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}\leq\sqrt{\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus i}\|^{2}}\cdot
A+\frac{A^{2}}{3}\ \Longrightarrow\ \sum_{i\in[n]}\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}\leq A^{2}\cdot\left(\frac{1+\sqrt{7/3}}{2}\right)^{2}\leq 1.6A^{2}.$

This completes the proof of the theorem. ∎

######  Lemma 12.

Let $B\in\\{0,1\\}^{n}$ be a uniformly random permutation of
$(\mathbf{1}_{k},\mathbf{0}_{n-k})$, and let $p=k/n$. Let
$c^{(1)},\dots,c^{(K)}\in\mathbb{R}^{n}$ be fixed, with
$\max_{j\in[K]}\|c^{(j)}\|_{2}\leq C_{2}$ and
$\max_{j\in[K]}\|c^{(j)}\|_{\infty}\leq C_{\infty}$. Then

$\mathbb{E}\left[\max_{j\in[K]}c^{(j)\top}(B-p\mathbf{1}_{n})\right]\leq
C_{2}\sqrt{2(1+\eta_{n})(p(1-p)+4\eta_{n})\log
K}+\frac{2C_{\infty}}{3}(1+\eta_{n})\log K,$

where $\eta_{n}=\frac{H_{n}-1}{n-H_{n}}$, for
$H_{n}=1+\frac{1}{2}+\dots+\frac{1}{n}$ denoting the $n$th harmonic number.

###### Proof.

Applying (Barber, 2024, Theorem 10), for each $j$, for any
$0<\lambda<\frac{3}{2\|c^{(j)}\|_{\infty}(1+\eta_{n})}$, it holds that

$\mathbb{E}\left[\exp\left\\{\lambda
c^{(j)\top}(B-p\mathbf{1}_{n})\right\\}\right]\leq\exp\left\\{\frac{\lambda^{2}\|c^{(j)}\|^{2}_{2}(1+\eta_{n})}{2(1-\frac{2\lambda}{3}\|c^{(j)}\|_{\infty}(1+\eta_{n}))}\cdot\left(p(1-p)+4\eta_{n}\right)\right\\}.$

Therefore, from the definition of $C_{2}$ and $C_{\infty}$, for any
$0<\lambda<\frac{3}{2(1+\eta_{n})C_{\infty}}$, this implies

$\mathbb{E}\left[\exp\left\\{\lambda
c^{(j)\top}(B-p\mathbf{1}_{n})\right\\}\right]\leq\exp\left\\{\frac{\lambda^{2}C_{2}^{2}(1+\eta_{n})}{2(1-\frac{2\lambda}{3}(1+\eta_{n})C_{\infty})}\cdot\left(p(1-p)+4\eta_{n}\right)\right\\}$

for each $j$. We then have

$\displaystyle\mathbb{E}\left[\max_{j\in[K]}c^{(j)\top}(B-p\mathbf{1}_{n})\right]$
$\displaystyle\leq\lambda^{-1}\log\left(\mathbb{E}\left[\exp\left\\{\lambda\max_{j\in[K]}c^{(j)\top}(B-p\mathbf{1}_{n})\right\\}\right]\right)\textnormal{
by Jensen's inequality}$
$\displaystyle\leq\lambda^{-1}\log\left(\sum_{j\in[K]}\mathbb{E}\left[\exp\left\\{\lambda
c^{(j)\top}(B-p\mathbf{1}_{n})\right\\}\right]\right)$
$\displaystyle\leq\lambda^{-1}\log\left(K\cdot\exp\left\\{\frac{\lambda^{2}C_{2}^{2}(1+\eta_{n})}{2(1-\frac{2\lambda}{3}(1+\eta_{n})C_{\infty})}\cdot\left(p(1-p)+4\eta_{n}\right)\right\\}\right)$
$\displaystyle=\frac{\log K}{\lambda}+\frac{\lambda
C_{2}^{2}(1+\eta_{n})}{2(1-\frac{2\lambda}{3}(1+\eta_{n})C_{\infty})}\cdot\left(p(1-p)+4\eta_{n}\right).$

Taking

$\lambda=\frac{\sqrt{\frac{2(1+\eta_{n})\log
K}{C_{2}^{2}(p(1-p)+4\eta_{n})}}}{(1+\eta_{n})\left(1+\frac{2C_{\infty}}{3}\sqrt{\frac{2(1+\eta_{n})\log
K}{C_{2}^{2}(p(1-p)+4\eta_{n})}}\right)}<\frac{3}{2(1+\eta_{n})C_{\infty}},$

this yields the desired bound. ∎

###  B.2 Is 2 needed?

We return to the role of 2, as discussed in Remark 11. In particular, this
assumption requires that $\mathcal{W}$ is approximately contained in a
$K$-dimensional subspace of $\mathcal{B}$.

We will now see that no uniform stability guarantee can hold without this type
of condition. A counterexample can be found with a standard infinite-
dimensional Banach space:
$\mathcal{B}=L_{1}(\mathbb{N})=\\{(x_{1},x_{2},\dots):\sum_{j\geq
1}|x_{j}|<\infty\\}$. We will equip $\mathcal{B}$ with the norm induced by the
$\mathsf{d}_{\mathsf{TV}}$ metric, given by
$\|x\|=\frac{1}{2}\|x\|_{1}=\frac{1}{2}\sum_{j\geq 1}|x_{j}|$. This is an
infinite-dimensional and separable Banach space.

Next, define $\mathcal{W}=\Delta_{\infty}=\\{x\in\mathcal{B}:x_{i}\geq 0\
\forall i,\ \|x\|_{1}=1\\}$, the infinite-dimensional probability simplex. To
relate this to our earlier examples, this can be viewed as an infinite-
dimensional version of Problem Setting 2(b)—$\mathcal{W}$ is the space of
distributions over $\mathbb{N}$, so for instance, we might have an algorithm
$\mathcal{A}$ that returns a posterior distribution over _countably
infinitely_ many possible models. Although $\mathcal{W}$ is bounded in the
norm $\|\cdot\|$, it cannot be approximated by any finite-dimensional subspace
in $\mathcal{B}$. The following result shows stability can fail when bagging
algorithms with output in $\mathcal{W}$.

######  Proposition 13.

Let $\mathcal{Z}=\mathbb{N}$, the set of positive integers. Let
$\mathcal{B}=L_{1}(\mathbb{N})$ (equipped with the norm
$\|x\|=\frac{1}{2}\sum_{i\geq 1}|x_{i}|$), and let
$\mathcal{W}=\Delta_{\infty}$. Then there exists a (deterministic) algorithm
$\mathcal{A}:\cup_{n\geq 0}\mathcal{Z}^{n}\times[0,1]\to\mathcal{W}$, such
that for all $n\geq 1$, for any distribution $\mathcal{Q}_{n}$ satisfying 1,
derandomized bagging fails to satisfy mean-square stability with respect to
$\mathsf{d}_{\mathsf{TV}}$ for any $\beta<p_{\mathcal{Q}_{n}}$.

In other words, stability results of the type shown in Theorem 10 guarantee a
bound on $\beta^{2}$ (the stability parameter), uniformly over all algorithms
$\mathcal{A}$, that vanishes as $n\to\infty$. In contrast, the above example,
we show that—for a particular construction of the algorithm $\mathcal{A}$—the
stability parameter $\beta^{2}$ is bounded away from zero even as
$n\to\infty$.

###### Proof of 13.

Let $\mathcal{S}$ denote all finite subsets of $\mathbb{N}$, and let
$\psi:\mathcal{S}\to\mathbb{N}$ be a bijective map that assigns each finite
subset of $\mathbb{N}$ to a (unique) positive integer. For any dataset
$\mathcal{D}=(Z_{1},\dots,Z_{m})\in\mathbb{N}^{m}$, we define

$S(\mathcal{D})=\left\\{i\in\mathbb{N}:\sum_{j\in[m]}\mathbf{1}\\{Z_{j}=i\\}\geq
1\right\\},$

the set of unique values appearing in the data set. Then define
$\mathcal{A}(\mathcal{D};\xi)=\mathbf{e}_{\psi(S(\mathcal{D}))}\in\mathcal{W},$
where $\mathbf{e}_{j}$ is the $j$th canonical basis vector, with a $1$ in the
$j$th location and $0$’s elsewhere.

Now consider data set $\mathcal{D}=(1,\dots,n)$. Let $\hat{w}$ and
$\hat{w}^{\setminus i}$ denote the derandomized bagged estimators, as before,
for data sets $\mathcal{D}$ and $\mathcal{D}^{\setminus i}$, respectively. By
construction of $\mathcal{A}$, we can verify that
$\hat{w}_{j}=\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{\psi(S(\mathcal{D}^{r}))=j\\}$
and similarly $\hat{w}^{\setminus
i}_{j}=\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{\psi(S(\mathcal{D}^{r}))=j\mid
i\not\in r\\}$. Therefore,

$\displaystyle\|\hat{w}-\hat{w}^{\setminus i}\|$
$\displaystyle=\frac{1}{2}\sum_{j\in\mathbb{N}}\left|\hat{w}_{j}-\hat{w}^{\setminus
i}_{j}\right|=\frac{1}{2}\sum_{j\in\mathbb{N}}\left|\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{\psi(S(\mathcal{D}^{r}))=j\\}-\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{\psi(S(\mathcal{D}^{r}))=j\mid
i\not\in r\\}\right|$
$\displaystyle=\frac{1}{2}\sum_{j\in\mathbb{N}}\left|\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{\psi(S(\mathcal{D}^{r}))=j\\}-\frac{\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{\psi(S(\mathcal{D}^{r}))=j\\}}{1-p_{\mathcal{Q}_{n}}}\cdot\mathbf{1}\\{i\not\in\psi^{-1}(j)\\}\right|$
$\displaystyle=\frac{1}{2}\sum_{j\in\mathbb{N}}\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{\psi(S(\mathcal{D}^{r}))=j\\}\cdot\left|1-\frac{\mathbf{1}\\{i\not\in\psi^{-1}(j)\\}}{1-p_{\mathcal{Q}_{n}}}\right|$
$\displaystyle=\frac{1}{2}\sum_{j\in\mathbb{N}}\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{\psi(S(\mathcal{D}^{r}))=j\\}\cdot\frac{(1-p_{\mathcal{Q}_{n}})\cdot\mathbf{1}\\{i\in
S(\mathcal{D}^{r})\\}+p_{\mathcal{Q}_{n}}\cdot\mathbf{1}\\{i\not\in
S(\mathcal{D}^{r})\\}}{1-p_{\mathcal{Q}_{n}}}$
$\displaystyle=\frac{1}{2}\frac{1}{1-p_{\mathcal{Q}_{n}}}\left((1-p_{\mathcal{Q}_{n}})\cdot\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{i\in
S(\mathcal{D}^{r})\\}+p_{\mathcal{Q}_{n}}\cdot\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{i\not\in
S(\mathcal{D}^{r})\\}\right)$ $\displaystyle=p_{\mathcal{Q}_{n}},$

since $\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{i\in
S(\mathcal{D}^{r})\\}=\mathbb{P}_{r\sim\mathcal{Q}_{n}}\\{i\in
r\\}=p_{\mathcal{Q}_{n}}$. Since this holds for every $i$, we have
$\frac{1}{n}\sum_{i\in[n]}\|\hat{w}-\hat{w}^{\backslash
i}\|^{2}=p_{\mathcal{Q}_{n}}^{2}$, which completes the proof. ∎

####  B.2.1 Remark on dimensionality

To relate the construction in this proof back to 2 and Theorem 10, note that
in the construction of $\mathcal{A}$, at each sample size $n$ we are
essentially working in a $K$-dimensional space for $K=2^{n}$—that is, even
though $\mathcal{W}=\Delta_{\infty}$ is the space of all infinite sequences,
the algorithm $\mathcal{A}$ only ever returns sequences $\hat{w}$ with nonzero
values in the first $K=2^{n}$ entries (and in particular, $\hat{w}$ will
always lie in the convex hull of $\mathbf{e}_{1},\dots,\mathbf{e}_{K}$).

Compare this to the result of Theorem 10, which essentially suggests that
mean-square stability vanishes at the rate $\beta^{2}\propto\frac{\log K}{n}$
if $\mathcal{W}$ can be approximated as a convex hull of $\mathcal{O}(K)$ many
points. Since we have taken $\frac{\log K}{n}$ to be nonvanishing in this
example (due to $K=2^{n}$ being the “effective dimension”, as described
above), stability is no longer ensured as $n\to\infty$.

###  B.3 Additional results

First, we state a result for tail stability in the Banach space setting.

######  Corollary 14.

Under the notation and assumptions of Theorem 10, derandomized bagging has
tail stability $(\varepsilon,\delta)$ with respect to $\|\cdot\|$ for any
$\varepsilon,\delta\geq 0$ satisfying

$\delta\varepsilon^{2}\geq
1.6\left[R\sqrt{\frac{2(1+\eta_{n})\log(2K)}{n}\left(\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}+\frac{4\eta_{n}}{(1-p_{\mathcal{Q}_{n}})^{2}}\right)}+R\sqrt{\frac{1}{n}\cdot\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}}+2\rho\cdot
p_{\mathcal{Q}_{n}}\right]^{2}.$

###### Proof.

This result follows immediately from Theorem 10 together with 3. ∎

Finally, we state a finite-$B$ result for the Banach space setting.

######  Corollary 15.

Under the notation and assumptions of Theorem 10, let $B\geq 1$ be fixed. Then
bagging, run with base algorithm $\mathcal{A}$ and resampling distribution
$\mathcal{Q}_{n}$ and with $B$ bags, has mean-square stability $\beta^{2}$
with respect to $\|\cdot\|$, with

$\beta^{2}:=\inf_{\varepsilon>0}\Bigg{\\{}(1+\varepsilon)\cdot
1.6\bigg{[}R\sqrt{\frac{2(1+\eta_{n})\log(2K)}{n}\left(\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}+\frac{4\eta_{n}}{(1-p_{\mathcal{Q}_{n}})^{2}}\right)}+R\sqrt{\frac{1}{n}\cdot\frac{p_{\mathcal{Q}_{n}}}{1-p_{\mathcal{Q}_{n}}}}\\\
{}+2\rho\cdot
p_{\mathcal{Q}_{n}}\bigg{]}^{2}{}+(1+\varepsilon^{-1})\cdot\left(\frac{8KR^{2}}{B}+8\rho^{2}\right)\Bigg{\\}}$

###### Proof.

Let $\hat{w}$ and $\hat{w}^{\setminus i}$ be the output of the bagged
algorithm (as in Definition 4), and let $\tilde{w}$ and $\tilde{w}^{\setminus
i}$ be the corresponding derandomized bagged estimators (as in Definition 5).
We have

$\frac{1}{n}\sum_{i\in[n]}\mathbb{E}\left[\|\hat{w}-\hat{w}^{\setminus
i}\|^{2}\right]=\frac{1}{n}\sum_{i\in[n]}\mathbb{E}\left[\|(\hat{w}-\tilde{w})-(\hat{w}^{\setminus
i}-\tilde{w}^{\setminus i})+(\tilde{w}-\tilde{w}^{\setminus
i})\|^{2}\right]\\\
\leq\frac{1+\varepsilon}{n}\sum_{i\in[n]}\|\tilde{w}-\tilde{w}^{\setminus
i}\|^{2}+\frac{2(1+\varepsilon^{-1})}{n}\sum_{i\in[n]}\mathbb{E}\left[\|\hat{w}-\tilde{w}\|^{2}\right]+\frac{2(1+\varepsilon^{-1})}{n}\sum_{i\in[n]}\mathbb{E}\left[\|\hat{w}^{\setminus
i}-\tilde{w}^{\setminus i}\|^{2}\right],$

since
$\|a+b+c\|^{2}\leq(1+\varepsilon)\|a\|^{2}+(1+\varepsilon^{-1})\|b+c\|^{2}\leq(1+\varepsilon)\|a\|^{2}+2(1+\varepsilon^{-1})(\|b\|^{2}+\|c\|^{2})$,
for any $a,b,c\in\mathcal{B}$ (because $\|\cdot\|$ satisfies the triangle
inequality). The first term in the upper bound is bounded by Theorem 10. Now
we bound the remaining terms.

First, let $\hat{w}^{\prime}$ be an i.i.d. copy of $\hat{w}$. Under the same
notation as in the proof of Theorem 10, let $\|\cdot\|_{*}$ denote the dual
norm. Then

$\mathbb{E}[\|\hat{w}-\tilde{w}\|^{2}]=\mathbb{E}\left[\sup_{f:\|f\|_{*}\leq
1}f(\hat{w}-\tilde{w})^{2}\right]=\mathbb{E}\left[\sup_{f:\|f\|_{*}\leq
1}f\left(\hat{w}-\mathbb{E}[\hat{w}^{\prime}\mid\hat{w}]\right)^{2}\right]\\\
=\mathbb{E}\left[\sup_{f:\|f\|_{*}\leq
1}\mathbb{E}\left[f\left(\hat{w}-\hat{w}^{\prime}\right)\mid\hat{w}\right]^{2}\right]\leq\mathbb{E}\left[\mathbb{E}\left[\sup_{f:\|f\|_{*}\leq
1}f\left(\hat{w}-\hat{w}^{\prime}\right)^{2}\mid\hat{w}\right]\right]=\mathbb{E}[\|\hat{w}-\hat{w}^{\prime}\|^{2}],$

where the second step holds since $\hat{w}\perp\\!\\!\\!\perp\hat{w}^{\prime}$
and $\mathbb{E}[\hat{w}^{\prime}]=\mathbb{E}[\hat{w}]=\tilde{w}$. Now let
$\hat{w}^{(b)}=\mathcal{A}(\mathcal{D}^{r_{b}};\xi_{b})$, for each $b\in[B]$,
so that we have $\hat{w}=\frac{1}{B}\sum_{b=1}^{B}\hat{w}^{(b)}$, and
similarly write
$\hat{w}^{\prime}{}^{(b)}=\mathcal{A}(\mathcal{D}^{r_{b}^{\prime}};\xi_{b}^{\prime})$,
where the $r_{b}^{\prime}$ and $\xi^{\prime}_{b}$ samples are independently
drawn, so that
$\hat{w}^{\prime}=\frac{1}{B}\sum_{b=1}^{B}\hat{w}^{\prime}{}^{(b)}$. Then
$\hat{w}-\hat{w}^{\prime}$ is equal in distribution to
$\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}\left(\hat{w}^{(b)}-\hat{w}^{\prime}{}^{(b)}\right)$
where the $\zeta_{b}$’s are i.i.d. random signs.

Now we will bound
$\mathbb{E}[\|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}(w_{b}-w^{\prime}_{b})\|^{2}]$
for any _fixed_
$w_{1},\dots,w_{B},w^{\prime}_{1},\dots,w^{\prime}_{B}\in\mathcal{W}$. For
each $b$ we can write
$w_{b}-w^{\prime}_{b}=\sum_{j=1}^{K}t^{b}_{j}v_{j}+\Delta_{b}$ with
$\|t^{b}\|_{1}\leq R$ and $\|\Delta_{b}\|\leq\rho$, by 2. Then

$\left\|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}(w_{b}-w^{\prime}_{b})\right\|=\left\|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}\left(\sum_{j=1}^{K}t^{b}_{j}v_{j}+\Delta_{b}\right)\right\|\\\
\leq\sum_{j=1}^{K}\left\|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}t^{b}_{j}v_{j}\right\|+\left\|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}\Delta_{b}\right\|\leq\sum_{j=1}^{K}\left|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}t^{b}_{j}\right|+\rho,$

since $\|v_{j}\|=1$ for all $j$, and $\|\Delta_{b}\|\leq\rho$ for all $b$.
Therefore,

$\mathbb{E}\left[\left\|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}(w_{b}-w^{\prime}_{b})\right\|^{2}\right]\leq\mathbb{E}\left[\left(\sum_{j=1}^{K}\left|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}t^{b}_{j}\right|+\rho\right)^{2}\right]\leq
2K\sum_{j=1}^{K}\mathbb{E}\left[\left|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}t^{b}_{j}\right|^{2}\right]+2\rho^{2}\\\
\leq\frac{2K}{B^{2}}\sum_{j=1}^{K}\sum_{b=1}^{B}(t_{j}^{b})^{2}+2\rho^{2}=\frac{2K}{B^{2}}\sum_{b=1}^{B}\|t^{b}\|^{2}_{2}+2\rho^{2}\leq\frac{2K}{B^{2}}\sum_{b=1}^{B}\|t^{b}\|^{2}_{1}+2\rho^{2}\leq\frac{2KR^{2}}{B}+2\rho^{2}.$

Combining all our work we then have

$\mathbb{E}[\|\hat{w}-\tilde{w}\|^{2}]\leq\mathbb{E}[\|\hat{w}-\hat{w}^{\prime}\|^{2}]=\mathbb{E}\left[\left\|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}(\hat{w}^{b}-\hat{w}^{\prime}{}^{b})\right\|^{2}\right]\\\
=\mathbb{E}\left[\mathbb{E}\left[\left\|\frac{1}{B}\sum_{b=1}^{B}\zeta_{b}(\hat{w}^{b}-\hat{w}^{\prime}{}^{b})\right\|^{2}\
\middle|\
\hat{w}^{1},\dots,\hat{w}^{B},\hat{w}^{\prime}{}^{1},\dots,\hat{w}^{\prime}{}^{B}\right]\right]\leq\frac{2KR^{2}}{B}+2\rho^{2}.$

The same bound holds for $\mathbb{E}[\|\hat{w}^{\setminus
i}-\tilde{w}^{\setminus i}\|^{2}]$. Combining all our calculations establishes
the result. ∎

## References

  * (1)
  * Abadie (2021) Abadie, A. (2021).  Using synthetic controls: Feasibility, data requirements, and methodological aspects, Journal of Economic Literature 59(2): 391–425.
  * Abadie et al. (2010) Abadie, A., Diamond, A. and Hainmueller, J. (2010).  Synthetic control methods for comparative case studies: estimating the effect of California’s tobacco control program, J. Amer. Statist. Assoc. 105(490): 493–505.
  * Abadie et al. (2015) Abadie, A., Diamond, A. and Hainmueller, J. (2015).  Comparative politics and the synthetic control method, American Journal of Political Science 59(2): 495–510.
  * Abadie and Gardeazabal (2003) Abadie, A. and Gardeazabal, J. (2003).  The economic costs of conflict: A case study of the basque country, American economic review 93(1): 113–132.
  * Barber (2024) Barber, R. F. (2024).  Hoeffding and Bernstein inequalities for weighted sums of exchangeable random variables, arXiv preprint arXiv:2404.06457 .
  * Bousquet and Elisseeff (2002) Bousquet, O. and Elisseeff, A. (2002).  Stability and generalization, The Journal of Machine Learning Research 2: 499–526.
  * Breiman (1996a) Breiman, L. (1996a).  Bagging predictors, Machine learning 24(2): 123–140.
  * Breiman (1996b) Breiman, L. (1996b).  Heuristics of instability and stabilization in model selection, The Annals of Statistics 24(6): 2350–2383.
  * Bühlmann (2014) Bühlmann, P. (2014).  Discussion of big Bayes stories and BayesBag, Statistical science 29(1): 91–94.
  * Devroye and Wagner (1979a) Devroye, L. and Wagner, T. (1979a).  Distribution-free inequalities for the deleted and holdout error estimates, IEEE Transactions on Information Theory 25(2): 202–207.
  * Devroye and Wagner (1979b) Devroye, L. and Wagner, T. (1979b).  Distribution-free performance bounds for potential function rules, IEEE Transactions on Information Theory 25(5): 601–604.
  * Folland (1999) Folland, G. B. (1999).  Real analysis, Pure and Applied Mathematics (New York), second edn, John Wiley & Sons, Inc., New York.  Modern techniques and their applications, A Wiley-Interscience Publication.
  * Hayes (2005) Hayes, T. P. (2005).  A large-deviation inequality for vector-valued martingales, Combinatorics, Probability and Computing .
  * Huggins and Miller (2023) Huggins, J. H. and Miller, J. W. (2023).  Reproducible model selection using bagged posteriors, Bayesian Analysis 18(1): 79–104.
  * Kim and Barber (2023) Kim, B. and Barber, R. F. (2023).  Black-box tests for algorithmic stability, Inf. Inference 12(4): Paper No. iaad039, 30.
https://doi.org/10.1093/imaiai/iaad039

  * Mukherjee et al. (2006) Mukherjee, S., Niyogi, P., Poggio, T. and Rifkin, R. (2006).  Learning theory: stability is sufficient for generalization and necessary and sufficient for consistency of empirical risk minimization, Advances in Computational Mathematics 25(1): 161–193.
  * 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.
  * Shalev-Shwartz et al. (2010) Shalev-Shwartz, S., Shamir, O., Srebro, N. and Sridharan, K. (2010).  Learnability, stability and uniform convergence, The Journal of Machine Learning Research 11: 2635–2670.
  * Soloff et al. (2024) Soloff, J. A., Barber, R. F. and Willett, R. (2024).  Bagging provides assumption-free stability, Journal of Machine Learning Research 25(131): 1–35.
http://jmlr.org/papers/v25/23-0536.html

  * Vaníček (1969) Vaníček, P. (1969).  Approximate spectral analysis by least-squares fit: Successive spectral analysis, Astrophysics and Space Science 4: 387–391.
  * Xu et al. (2011) Xu, H., Caramanis, C. and Mannor, S. (2011).  Sparse algorithms are not stable: A no-free-lunch theorem, IEEE transactions on pattern analysis and machine intelligence 34(1): 187–193.
  * Yu (2013) Yu, B. (2013).  Stability, Bernoulli 19(4): 1484–1500.
  * 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, Proc. Natl. Acad. Sci 117(8): 3920–3929.