# Local Risk Bounds for Statistical Aggregation

Jaouad Mourtada Tomas Vaškevičius Nikita Zhivotovskiy CREST, ENSAE, Institut
Polytechnique de Paris, France<EMAIL_ADDRESS>of
Mathematics, EPFL, Switzerland<EMAIL_ADDRESS>of
Statistics, University of California Berkeley, United States,
<EMAIL_ADDRESS>
###### Abstract

In the problem of aggregation, the aim is to combine a given class of base
predictors to achieve predictions nearly as accurate as the best one. In this
flexible framework, no assumption is made on the structure of the class or the
nature of the target. Aggregation has been studied in both sequential and
statistical contexts. Despite some important differences between the two
problems, the classical results in both cases feature the same global
complexity measure. In this paper, we revisit and tighten classical results in
the theory of aggregation in the statistical setting by replacing the global
complexity with a smaller, local one. Some of our proofs build on the PAC-
Bayes localization technique introduced by Catoni. Among other results, we
prove localized versions of the classical bound for the exponential weights
estimator due to Leung and Barron [Trans. Inf. Theory., 52(8):3396–3410, 2006]
and deviation-optimal bounds for the $Q$-aggregation estimator. These bounds
improve over the results of Dai, Rigollet and Zhang [Ann. Statist.
40(3):1878-1905, 2012] for fixed design regression and the results of Lecué
and Rigollet [Ann. Statist. 42(1):211-224, 2014] for random design regression.

##  1 Introduction

Aggregation††footnotetext: Accepted for presentation at the Conference on
Learning Theory (COLT) 2023. is a fundamental task in statistical learning
theory and sequential prediction. It can be informally described as follows:
given a family of predictors, combine their predictions to produce forecasts
almost as accurate as those of the best predictor in the provided class. This
basic problem can be instantiated in a variety of situations, and provides a
flexible building block to construct accurate and theoretically sound
procedures. The aggregation problem has been explored both in the statistical
setting [56, 72, 16, 63, 44, 22, 7, 41, 43, 13, 53] and in the sequential
setting (see, e.g., [33, 28, 68, 36, 11, 69] and the classic survey [19]). It
is also related to the classical _model selection_ problem in statistics; see,
e.g., the lecture notes [48] for an introduction to the modern non-asymptotic
treatment of the subject.

In this work, we focus on statistical variants of the aggregation problem,
primarily within the frameworks of _statistical learning/random-design
regression_ and _denoising/fixed-design regression_. We aim to relate and
contrast these statistical problems with their sequential counterparts,
emphasizing suitable complexity measures associated with each problem setting.
Before highlighting the differences between the two variants of aggregation,
we start by discussing some important links and commonalities between them.

##### Sequential and statistical aggregation.

For concreteness and to emphasize the complexity aspect of the problem, we
will concentrate on regression with the squared loss. Henceforth, let
$\mathcal{X}$ be an input space and $\\{f_{\theta}:\theta\in\Theta\\}$ be a
class of functions $\mathcal{X}\to[-1,1]$, called base predictors or experts.
The goal is to approximate the performance of the best function in this class.
An important special case occurs when $\Theta$ is a finite but otherwise
arbitrary set. In the statistical context, such a setting is called _model
aggregation_ [56, 63], also known as _expert aggregation_ in the sequential
setting.

In the sequential setting, we are given a sequence of observations
$(x_{1},y_{1}),\dots,(x_{n},y_{n})\in\mathcal{X}\times[-1,1]$ revealed one at
a time. The aim is to construct a sequence of predictors
$\widehat{f}_{0},\dots,\widehat{f}_{n-1}:\mathcal{X}\to[-1,1]$ such that
$\widehat{f}_{i}$ depends only on the first $i$ observations and for any data
sequence, the _regret_

$\frac{1}{n}\sum_{i=1}^{n}\big{(}\widehat{f}_{i-1}(x_{i})-y_{i}\big{)}^{2}-\min_{\theta\in\Theta}\frac{1}{n}\sum_{i=1}^{n}\big{(}f_{\theta}(x_{i})-y_{i}\big{)}^{2}$
(1)

is upper bounded by an appropriate function of $M$ and $n$.

In contrast, in the statistical learning framework, the input-output pairs
$(X_{i},Y_{i})_{1\leqslant i\leqslant n}$ are assumed to be $n$ i.i.d. random
variables with the same distribution $P$ on $\mathcal{X}\times[-1,1]$ as the
generic pair $(X,Y)$. The quality of a regression function
$f:\mathcal{X}\to\mathbf{R}$ is then measured by its _risk_
$R(f)={\mathbf{E}}[(f(X)-Y)^{2}]$. The distribution $P$ is unknown, but the
sample $(X_{i},Y_{i})_{1\leqslant i\leqslant n}$ is available. The aim is to
produce a predictor $\widehat{f}_{n}:\mathcal{X}\to\mathbf{R}$ depending on
the random sample that achieves a small _excess risk_

$R(\widehat{f}_{n})-\min_{\theta\in\Theta}R(f_{\theta}),$

either in expectation or with high probability.

The fact that the sequential problem deals with arbitrary (deterministic) data
suggests that guarantees in this setting can be transported to the statistical
setting. This is indeed partly the case, as shown by a standard conversion
procedure (see [12, 46, 72, 16, 18]): for any sequence of predictors
$\widehat{f}_{0},\dots,\widehat{f}_{n}$ (with $\widehat{f}_{i}$ based on the
observations $(X_{j},Y_{j})_{1\leqslant j\leqslant i}$) and any reference
function $f:\mathcal{X}\to\mathbf{R}$, we have

$\frac{1}{n+1}\sum_{i=1}^{n+1}{\mathbf{E}}[R(\widehat{f}_{i-1})]-R(f)={\mathbf{E}}\bigg{[}\frac{1}{n+1}\sum_{i=1}^{n+1}(\widehat{f}_{i-1}(X_{i})-Y_{i})^{2}-\frac{1}{n+1}\sum_{i=1}^{n+1}(f(X_{i})-Y_{i})^{2}\bigg{]},$
(2)

which holds because $\widehat{f}_{i-1}$ is independent of $(X_{i},Y_{i})$.
Hence, the empirical regret with respect to $f_{0}$ controls the average
excess risk of the predictors $\widehat{f}_{i}$ for $1\leqslant i\leqslant
n+1$. This conversion is, however, only partial for two reasons. First, the
obtained guarantee is only on the _expected_ risk rather than its deviations;
as shown by Audibert [7], there may be a real difference between the expected
and high-probability behavior of the risk. Second, the excess risk bound
obtained this way features the complexity measure of the harder, sequential
version of the problem. This leads to suboptimal bounds, as discussed below.

##### Limitation of selection rules and the need for aggregation.

An important aspect of aggregation theory is that it allows us to treat rather
generic problems without imposing strong structural assumptions. For instance,
the class of base predictors may have unfavorable geometry due to being non-
convex, and likewise, the target to be predicted is not assumed to lie in some
favorable position relative to the class of predictors. An interesting
consequence of this lack of structure is that _selection/proper learning_
rules, namely, rules that select an element of the class of base predictors,
cannot achieve optimal prediction accuracy. In particular, this eliminates the
most natural procedure – that of empirical risk minimization – as a candidate
for optimal aggregation.

The necessity of improper learning is a fundamental facet of the theory of
aggregation, shared by both the sequential and statistical versions of the
problem. More precisely, no proper/selection rule (deterministic or
randomized) can improve upon the suboptimal worst-case regret or excess risk
rate of order $\sqrt{\log M/n}$. This lower bound holds due to the possible
non-convexity of the class and the arbitrary position of the target $Y$ with
respect to it. Roughly speaking, non-convexity can induce an erratic behavior
of selection rules and prevents one from exploiting the curvature of the
squared loss. In contrast, when the base class of predictors is convex,
curvature of the loss function can be exploited to obtain a faster rate
scaling as $\log M/n$ for empirical risk minimization. For relevant lower
bounds in the statistical setting, we refer to [16, §3.5], [3, §20.2] and [37,
7]. A corresponding lower bound in the sequential setup, which can be deduced
from the statistical lower bound by (2), can be found in [19, p. 72].

##### Aggregation by exponential weights.

Perhaps the most important aggregation procedure, both in the sequential and
statistical settings, is the method of _exponential weights_. This estimator
is useful and nontrivial even in the case of model aggregation, but it is most
naturally defined and analyzed in a more general context described below.

Let $\\{f_{\theta}:\theta\in\Theta\\}$ denote the class of base predictors,
where the parameter set $\Theta$ is now general and not necessarily finite. As
before, we make no special assumptions on the structure of the class, which
may be non-convex. A natural idea to artificially introduce some structure to
the class is to “convexify” it. Specifically, let $\mathcal{P}(\Theta)$ denote
the set of probability measures on $\Theta$. This set possesses a convex
structure and $\Theta$ naturally embeds into it by associating each
$\theta\in\Theta$ to the Dirac mass $\delta_{\theta}$. In addition, each
function $h:\Theta\to\mathbf{R}$ defines a linear function on
$\mathcal{P}(\Theta)$ by integration:
$\rho\mapsto\langle\rho,h\rangle=\int_{\Theta}h(\theta)\rho(\mathrm{d}\theta)$.
Now, let $\pi$ be a probability measure on $\Theta$, called “prior” in analogy
with Bayesian statistics; this distribution encodes the “inductive bias” of
the procedure, and should generally put mass on “simpler” regions of the
class. The prior $\pi$ endows $\mathcal{P}(\Theta)$ with a convex regularizer
given by the _relative entropy_ or _Kullback-Leibler divergence_
$\rho\mapsto\mathrm{KL}({\rho},{\pi})=\int_{\Theta}\log(\frac{\mathrm{d}\rho}{\mathrm{d}\pi})\mathrm{d}\rho$,
with the convention that $\mathrm{KL}({\rho},{\pi})=+\infty$ when $\rho$ is
not absolutely continuous with respect to $\pi$. The _exponential weights_
rule is obtained by optimizing the linearized objective on
$\mathcal{P}(\Theta)$ with an entropic regularization. Specifically, the
exponential weights posterior
$\widehat{\rho}^{\mathsf{ew}}_{n}\in\mathcal{P}(\Theta)$ with prior $\pi$ and
inverse temperature (learning rate) $\beta>0$ is given by

$\widehat{\rho}^{\mathsf{ew}}_{n}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\rho\in\mathcal{P}(\Theta)}\big{\\{}\langle\rho,R_{n}\rangle+\beta^{-1}\cdot\mathrm{KL}({\rho},{\pi})\big{\\}}$
(3)

where $R_{n}(\theta)=n^{-1}\sum_{i=1}^{n}(f_{\theta}(X_{i})-Y_{i})^{2}$
denotes the empirical risk of the predictor $f_{\theta}$. The solution to this
optimization problem can be expressed explicitly using the Donsker-Varadhan
variational formula:

$\frac{\mathrm{d}\widehat{\rho}^{\mathsf{ew}}_{n}}{\mathrm{d}\pi}(\theta)=\frac{\exp\\{-\beta
R_{n}(\theta)\\}}{\int_{\Theta}\exp\\{-\beta
R_{n}(\vartheta)\\}\pi(\mathrm{d}\vartheta)}.$ (4)

Finally, each distribution $\rho$ gives rise to the “averaged” or “mixture”
predictor $f_{\rho}(x)=\int_{\Theta}f_{\theta}(x)\rho(\mathrm{d}\theta)$; in
particular, the exponential weights predictor is
$\widehat{f}^{\mathsf{ew}}_{n}=f_{\widehat{\rho}^{\mathsf{ew}}_{n}}$.

In the model aggregation context with $|\Theta|=M$, the exponential weights
algorithm with uniform prior on $\Theta$ and inverse temperature
$\beta_{i}=i/8$ for $\widehat{f}^{\mathsf{ew}}_{i}$ has a regret (1) bounded
by $8\log(M)/n$ [19, pp. 45-46]. By identity (2) and convexity of the risk,
this implies that the following averaged estimator, proposed in [12, 16, 72]
and known as _progressive mixture_ rule:

$\widehat{f}^{\mathsf{pm}}_{n}=\frac{1}{n+1}\sum_{i=0}^{n}\widehat{f}^{\mathsf{ew}}_{i}$
(5)

achieves the excess risk bound
${\mathbf{E}}[R(\widehat{f}^{\mathsf{pm}}_{n})]-\min_{\theta\in\Theta}R(f_{\theta})\leqslant
8\log(M)/(n+1)$, which is known to be the optimal rate of aggregation [63].
Remarkably, while working in the space of measures (the convex hull of the
class) improves the geometry of the class by ensuring convexity, the above
bound features the complexity of the original class rather than that of its
convex hull. Indeed, the best achievable excess risk over the convex hull of
$M$ functions is substantially larger than the $\log M/n$ rate for model
aggregation [63].

In fact, the progressive mixture was the only known rate-optimal model
aggregation procedure until the works of [7, 41], which introduced procedures
with optimal high-probability bounds (that the progressive mixture does not
achieve).

##### Global and local entropic complexities.

While the exponential weights/progressive mixture estimators achieve the
optimal rates of aggregation for finite classes, entropic penalization is
well-defined in a broader context. The class of “entropy-regularized”
procedures leads to a natural “entropic” complexity measure, which is
different from the more traditional measures of complexity such as covering
numbers, Rademacher complexity, or VC dimension [65, 39, 31, 70].

Entropic complexities depend on the prior $\pi$. For any $\pi$, it follows
from the regret bound for exponential weights in the sequential setting
together with the identity (2), that if $\beta_{i}=i/c$ with $c\geqslant 8$
then

$\frac{1}{n+1}\sum_{i=0}^{n}{\mathbf{E}}[R(\widehat{f}^{\mathsf{ew}}_{i})]\leqslant\inf_{\gamma\in\mathcal{P}(\Theta)}\bigg{\\{}\langle\gamma,R\rangle+c\cdot\frac{\mathrm{KL}({\gamma},{\pi})}{n+1}\bigg{\\}}.$
(6)

The quantity on the right-hand side plays an important role in the study of
exponential weights and entropy-regularized aggregation. We call this quantity
the _global entropic complexity_ [69, 73, 8]. It is defined as

$\mathsf{C}_{\mathsf{glob}}(\beta,\pi)=\inf_{\gamma\in\mathcal{P}(\Theta)}\Big{\\{}\langle\gamma,R\rangle+\beta^{-1}\cdot\mathrm{KL}({\gamma},{\pi})\Big{\\}}=-\frac{1}{\beta}\log\bigg{(}\int_{\Theta}e^{-\beta
R(\theta)}\pi(\mathrm{d}\theta)\bigg{)}.$ (7)

It is worth noting that the complexity $\mathsf{C}_{\mathsf{glob}}(\beta,\pi)$
only depends on $\beta$ and on the distribution of $R(\theta)$ under
$\theta\sim\pi$, and not explicitly on the “size” or dimension of the class.
In addition, the more mass the prior puts on parameters $\theta$ with small
risk, the smaller the complexity; this quantity, therefore, measures the
quality of the prior $\pi$ for the actual risk $R$.

The global entropic complexity appears in classical results in aggregation,
both in the sequential and statistical settings (and within the statistical
framework, both for learning and denoising). For instance, in the statistical
learning setting, the bound (6) together with convexity of the risk $R$
implies that the progressive mixture (5) satisfies the following bound:
${\mathbf{E}}[R(\widehat{f}^{\mathsf{pm}}_{n})]\leqslant\mathsf{C}_{\mathsf{glob}}((n+1)/c,\pi)$.

At the same time, sequential prediction seems genuinely harder than
statistical estimation, for the following reason. In sequential prediction,
one must issue predictions at each time step $i=1,\dots,n$, having seen only
the first $i-1$ observations, and the error is averaged over sample sizes.
Clearly, less information is available at early rounds, which should result in
a lower precision. In contrast, in statistical estimation, one has at hand a
full sample of $n$ observations. That the global complexity
$\mathsf{C}_{\mathsf{glob}}(\beta,\pi)$ bounds the hardness of sequential
prediction suggests that it provides an overly pessimistic account of the
hardness of the easier, statistical problem. This suboptimality manifests
itself in different ways, such as additional logarithmic factors, inability to
handle unbounded parameter domains, or it can provide misleading indication on
the qualitative effect of some problem parameters.

The refinement of the global entropic complexity that we consider is the
following _local entropic complexity_ : letting
$\pi_{h}=e^{h}\pi/\langle\pi,e^{h}\rangle\in\mathcal{P}(\Theta)$ for a
function $h:\Theta\to\mathbf{R}$,

$\mathsf{C}_{\mathsf{loc}}(\beta,\pi)=\frac{\int_{\mathbf{R}^{d}}R(\theta)e^{-\beta
R(\theta)}\pi(\mathrm{d}\theta)}{\int_{\mathbf{R}^{d}}e^{-\beta
R(\theta)}\pi(\mathrm{d}\theta)}=\langle\pi_{-\beta R},R\rangle.$ (8)

The local complexity is the average risk under the exponential weights
posterior with respect to the _true_ (rather than empirical) risk. In analogy
with thermodynamics, the global complexity corresponds to the “free energy”,
while the local complexity corresponds to the “average energy”. The link
between both complexities is given by the relations

$\frac{\mathrm{d}}{\mathrm{d}\beta}\\{\beta\cdot\mathsf{C}_{\mathsf{glob}}(\beta,\pi)\\}=\mathsf{C}_{\mathsf{loc}}(\beta,\pi)\quad\text{{i.e.,}\@}\quad\mathsf{C}_{\mathsf{glob}}(\beta,\pi)=\frac{1}{\beta}\int_{0}^{\beta}\mathsf{C}_{\mathsf{loc}}(\eta,\pi)\mathrm{d}\eta.$
(9)

In addition, by a standard property of log-Laplace transforms,
$-\beta\cdot\mathsf{C}_{\mathsf{glob}}(\beta,\pi)$ is convex in $\beta$, so
that $\mathsf{C}_{\mathsf{loc}}(\beta,\pi)$ is decreasing in $\beta$.
Combining this with the second identity in (9) shows that
$\mathsf{C}_{\mathsf{loc}}(\beta,\pi)\leqslant\mathsf{C}_{\mathsf{glob}}(\beta,\pi)$.
Another essentially equivalent quantity is

$\inf_{\gamma\in\mathcal{P}(\Theta)}\Big{\\{}\langle\gamma,R\rangle+\beta^{-1}\cdot\mathrm{KL}({\gamma},{\pi_{-\beta
R}})\Big{\\}}=\mathsf{C}_{\mathsf{glob}}(\beta,\pi_{-\beta
R})=\frac{1}{\beta}\int_{\beta}^{2\beta}\mathsf{C}_{\mathsf{loc}}(\eta,\pi)\mathrm{d}\eta,$

(as since $\pi_{-\beta R})_{-\eta R}=\pi_{-(\beta+\eta)R}$. In particular, by
monotonicity of $\mathsf{C}_{\mathsf{loc}}(\eta,\pi)$ we have

$\mathsf{C}_{\mathsf{loc}}(2\beta,\pi)\leqslant\mathsf{C}_{\mathsf{glob}}(\beta,\pi_{-\beta
R})\leqslant\mathsf{C}_{\mathsf{loc}}(\beta,\pi).$ (10)

Combining inequality (6) with the relation (9) leads to the bound

$\frac{1}{n+1}\sum_{i=0}^{n}{\mathbf{E}}[R(\widehat{f}_{i})]\leqslant\frac{1}{n+1}\int_{0}^{n+1}\mathsf{C}_{\mathsf{loc}}(c\,\eta,\pi)\,\mathrm{d}\eta\leqslant\frac{1}{n+1}\sum_{i=0}^{n}\mathsf{C}_{\mathsf{loc}}(c\,i,\pi).$
(11)

In words, the averaged (over sample sizes $i=0,\dots,n$) risk of the
exponential weights estimator $\widehat{f}_{i}$ (based on the first $i$
observations and with inverse temperature $\beta_{i}=c\,i$) is bounded by the
averaged local complexity. Given this, it is natural to aim for the following
strengthening of (11): control the risk of exponential weights (or another
aggregation procedure) with $n$ samples by the local complexity
$\mathsf{C}_{\mathsf{loc}}(c\,n,\pi)$. However, this inequality does not
follow from (11); instead, the classical approach is to deduce (using
convexity of $R$) the global bound
${\mathbf{E}}[R(\widehat{f}^{\mathsf{pm}}_{n})]\leqslant\mathsf{C}_{\mathsf{glob}}(c(n+1),\pi)$
for the averaged estimator $\widehat{f}^{\mathsf{pm}}_{n}$. One of the main
contributions of this work is to show that the risk of suitable entropy-
regularized aggregation methods can indeed be controlled by the local
complexity $\mathsf{C}_{\mathsf{loc}}(c\,n,\pi)$, at the same level of
generality as the global bounds.

##### Related work.

Local priors and entropic complexities were first introduced by Catoni [16] in
the context of PAC-Bayesian bounds for classification, with a similar approach
proposed by Zhang [74]. Alquier [1] also obtained local PAC-Bayesian bounds
for regression. However, previous works [16, 74, 1, 2] either produced
generalization error bounds at a slow rate or assumed the target was in a
favorable position relative to the class (Bernstein assumption) rather than in
an arbitrary position as in the general aggregation problem. This
distinguishes our results from those of prior works.

Another relevant reference is [13], where the problem of aggregation is
considered in the specific case of tuning the regularization parameter for
ridge regression. In this case, the authors obtain guarantees with an excess
risk term that does not depend on the number of possible parameters in the
discrete grid, or on the range of possible parameters. These improvements are
related to those obtained from local bounds. However, the approach from [13]
seems to rely on the specific structure of the considered class, in particular
on the fact that it is comprised of ordered (affine) spectral smoothers.

Guarantees featuring the relative entropy as a complexity measure are
sometimes referred to as “PAC-Bayesian” in the literature [49, 17, 2].
Specifically, PAC-Bayesian inequalities, first used by McAllester [49, 50],
are obtained by combining the Chernoff method for establishing deviation
inequalities with duality for the relative entropy. We note that this
methodology has also been applied to study the generalization error
(difference between empirical and population risks) of learning algorithms,
which is distinct from the excess risk studied in this paper. We refer to a
sample of works in this direction [49, 50, 17, 25, 61, 47] and to the recent
surveys [2, 57] for further references.

The present work belongs to the literature on aggregation in the statistical
setting, which has been considered for both fixed and random design.
Guarantees in expectation for the exponential weights estimator were obtained
in [12, 72, 16, 37] (with additional averaging) for random design, and in [29,
44, 22] for fixed design. Sharp bounds in deviation for model aggregation were
first obtained in [7, 41]; closer to our work, the method of $Q$-aggregation
was shown in [59, 21, 43, 13] to achieve the worst-case optimal rate of
aggregation in deviation. We also mention that additional aggregation methods
have been proposed in [71, 64], using a sequential procedure featuring a
second-order correction related to that used in $Q$-aggregation; a related
approach has been considered in the context of online learning in [8, 40, 26].
Finally, extensions to heavy-tailed situations have been considered in [52,
53, 55].

###  1.1 Summary of main results

The rest of the paper is split into two parts. First, in Section 2, we state
our localized PAC-Bayesian bounds in three different statistical learning
settings. Then, in Section 3, we instantiate our abstract results to concrete
problems. Our main contributions can be summarized as follows:

  * •

In the fixed design regression setting, we prove two localized bounds. Theorem
1 provides a localized version of the classical risk bound for the exponential
weights estimator due to Leung and Barron [44]. Theorem 2 proves an analogous
bound for the $Q$-aggregation estimator, this time optimal in deviation.

  * •

Our fixed design local risk bound for the $Q$-aggregation estimator is
extended to the random design setting in Theorem 3. A special case of this
bound resolves a question raised by Lecué and Mendelson [42, Question 1.2] on
model aggregation guarantees that adapt to the hardness of the problem instead
of always paying for the worst-case minimax rate.

  * •

Section 3.2 deals with Gaussian process priors, making explicit some
distinctions between local and global complexities and connecting the
$Q$-aggregation and ridge regression estimators.

  * •

Section 3.3 is devoted to the analysis of improper ridge-type estimators
within the context of random design regression. Although we utilize a leave-
one-out type of analysis in this section, we are able to recover a complexity
measure similar to those that appear in Section 3.2.

  * •

Theorem 5 in Section 3.4 provides a deviation-optimal localized bound for the
$Q$-aggregation estimator in the transductive learning setup with almost
exchangeable priors. This result improves over the results from [58] for
Vapnik-Chervonenkis classes by accounting for their local structure in the
resulting bounds.

The proofs and further technical details are deferred to the appendix.

###  1.2 Notation

This section summarizes the notation used throughout the paper. We let
$\Theta$ and $\mathcal{X}$ be measurable spaces, and
$f:\Theta\times\mathcal{X}\to\mathbf{R}$ be a measurable function. For any
$\theta\in\Theta$, let $f_{\theta}=f(\theta,\cdot):\mathcal{X}\to\mathbf{R}$.
We denote by $\mathcal{P}(\Theta)$ the space of all probability distributions
on $\Theta$. Henceforth, we leave measurability assumptions implicit to
lighten the presentation, and adopt the convention that all introduced
functions between measurable spaces are measurable. In order to distinguish
between expectations over $\Theta$ (which come from the fact that we work on
the probability simplex $\mathcal{P}(\Theta)$) and randomness that comes from
random samples $(X,Y)$ in $\mathcal{X}\times\mathbf{R}$, we denote
expectations over $\Theta$ using integral and measure-theoretic notation,
while we use probabilistic notation for random variables such as random
samples. Specifically, for $\mu$ a measure on $\Theta$ and
$g:\Theta\to\mathbf{R}$, we denote by
$\langle\mu,g\rangle=\int_{\Theta}g(\theta)\mu(\mathrm{d}\theta)$ the integral
of $f$ with respect to the measure $\mu$, assuming it is well-defined. For any
distribution $\rho\in\mathcal{P}(\Theta)$, and any function
$f:\Theta\times\mathcal{X}\to\mathbf{R}$, the notation $\langle\rho,f\rangle$
averages over $\Theta$ and results in a function
$f_{\rho}:\mathcal{X}\to\mathbf{R}$ defined for any $x\in\mathcal{X}$ by
$f_{\rho}(x)=\langle\rho,f(\cdot,x)\rangle=\int_{\Theta}f_{\theta}(x)\rho(\mathrm{d}\theta).$
We introduce the shorthand notation $R(\rho)=R(f_{\rho})$ and
$\langle\rho,R_{n}\rangle=\int_{\Theta}R(f_{\theta})\rho(\mathrm{d}\theta).$
The empirical norm $\|\cdot\|_{n}^{2}$ is defined by
$\|f\|_{n}^{2}=\frac{1}{n}\sum_{i=1}^{n}f(X_{i})^{2}$. For any function
$g:\mathcal{X}\to\mathbf{R}$, and any probability measure
$\rho\in\mathcal{P}(\Theta)$, and any function
$f:\Theta\times\mathcal{X}\to\mathbf{R}$ we also introduce the following
shorthand notation
$\langle\rho,\|f-g\|_{n}^{2}\rangle=\int_{\Theta}\|f(\theta,\cdot)-g(\cdot)\|_{n}^{2}\rho(\mathrm{d}\theta)$.
For any functions $f,g$, we write $f\lesssim g$ if there is an absolute
constant $c>0$ such that $f\leqslant cg$.

##  2 Main results

We now state the main results of this paper: localized PAC-Bayesian bounds for
various statistical aggregation problems. In Section 2.1, we treat the fixed
design (or denoising) setting, where the covariates are assumed to be
deterministic. The case of random design regression (statistical learning) is
then considered in Section 2.2. We also mention that a third setting, namely
transductive learning, will be considered in Section 3.4.

While each of the learning formulations listed above admits a different
population risk functional $R$, its empirical counterpart $R_{n}$ remains the
same across the three problems. Given a sample $(X_{i},Y_{i})_{i=1}^{n}$ of
$n$ observations, the empirical risk $R_{n}(f)$ of a function
$f:\mathcal{X}\to\mathbf{R}$ is defined by

$R_{n}(f)=\frac{1}{n}\sum_{i=1}^{n}(f(X_{i})-Y_{i})^{2}.$ (12)

Our results apply to two entropy-regularized aggregation procedures. The first
one is the _exponential weights_ estimator. Given a prior distribution
$\pi\in\mathcal{P}(\Theta)$ and an inverse temperature parameter $\beta>0$, it
is defined by

$\widehat{f}^{\mathsf{ew}}_{n}=\langle\widehat{\rho}^{\mathsf{ew}}_{n},f\rangle,\text{
where
}\widehat{\rho}^{\mathsf{ew}}_{n}=\widehat{\rho}^{\mathsf{ew}}_{n}(\pi,\beta)=\mathop{\mathrm{arg}\,\mathrm{min}}_{\rho\in\mathcal{P}(\Theta)}\left\\{\langle\rho,R_{n}\rangle+\frac{1}{\beta}\mathrm{KL}(\rho,\pi)\right\\}.$
(13)

The estimator $\widehat{f}^{\mathsf{ew}}_{n}$ is analyzed in Section 2.1 in
the classical Gaussian fixed design setting. The purpose of this analysis is
to demonstrate some challenges associated with PAC-Bayesian localization in
the simplest possible setting and also, to provide an improvement over the
classical non-localized bound of Leung and Barron [44].

Despite optimal in-expectation performance in the denoising setting, the
exponential weights estimator $\widehat{f}^{\mathsf{ew}}_{n}$ is known to
suffer from sub-optimal tail behavior [21]. In the fixed design setting, a
known remedy is to replace the linear term $\langle\rho,R_{n}\rangle$ in (13)
by the convex combination
$\alpha\langle\rho,R_{n}\rangle+(1-\alpha)R_{n}(\rho)$ for any
$\alpha\in(0,1)$. The intuition for this manipulation is as follows. The
choice $\alpha=1$ yields the exponential weights estimator, which as mentioned
above is deviation sub-optimal. While minimizing the objective with the term
$\langle\rho,R_{n}\rangle$ introduces inductive bias towards “sparse”
probability measures, yielding good complexity regularization for the model
selection aggregation problem, this regularizing effect is completely lost at
the other extreme $\alpha=0$. On the other hand, any $\alpha\in(0,1)$ is a
good choice: the resulting objective still benefits from the inductive bias of
the linear term $\alpha\langle\rho,R_{n}\rangle$, while the term
$(1-\alpha)R_{n}(\alpha)$ allows us to benefit from the curvature of the
quadratic loss that is lost with the choice $\alpha=1$; this curvature of the
quadratic loss is crucial to regularize variance terms that render the
estimator $\widehat{f}^{\mathsf{ew}}$ deviation suboptimal despite its in-
expectation optimality for the denoising problem. While any $\alpha\in(0,1)$
would work, in order to simplify our results we define the _$Q$ -aggregation_
estimator with the particular choice $\alpha=1/2$:

$\widehat{f}^{\mathsf{q}}_{n}=\langle\widehat{\rho}^{\mathsf{q}}_{n},f\rangle,\text{
where
}\widehat{\rho}^{\mathsf{q}}_{n}=\widehat{\rho}^{\mathsf{q}}_{n}(\pi,\beta)=\mathop{\mathrm{arg}\,\mathrm{min}}_{\rho\in\mathcal{P}(\Theta)}\left\\{\frac{1}{2}\langle\rho,R_{n}\rangle+\frac{1}{2}R_{n}(\rho)+\frac{1}{\beta}\mathrm{KL}(\rho,\pi)\right\\}.$
(14)

This estimator was first proposed by Audibert [5] in the context of convex
aggregation (where one aims to compare to convex combinations of functions in
the class). It was first recognized as an optimal model aggregation procedure
by Rigollet [59], and variants thereof were further analyzed in [21, 43, 13].
In particular, Dai, Rigollet and Zhang [21] obtain a global PAC-Bayes bound in
the fixed design setup.

###  2.1 Fixed design

In this section, let $X_{1},\dots,X_{n}\in\mathcal{X}$ be arbitrary
_deterministic_ covariate vectors and consider the model:

$Y_{i}=f^{\star}(X_{i})+Z_{i},$ (15)

where $Y=(Y_{1},\dots,Y_{n})\in\mathbf{R}^{d}$ is the vector of response
variables, $Z=(Z_{i},\dots,Z_{n})$ is a random vector of zero-mean noise, and
$f^{\star}$ is the regression function. Observe that the only randomness in
the model (15) comes from the noise vector $Z$. Define the population risk
functional $R$ by

$R(f)=\mathbf{E}_{Y}\left[\frac{1}{n}\sum_{i=1}^{n}(f(X_{i})-Y_{i})^{2}\right]=\|f-f^{\star}\|_{n}^{2}+\frac{1}{n}\sum_{i=1}^{n}\mathrm{Var}(Z_{i}),$
(16)

while its empirical counterpart $R_{n}$ is defined in (12). With a slight
abuse of notation, in what follows we identify $f^{\star}$ with the vector
$(f^{\star}(X_{i}))_{1\leqslant i\leqslant n}\in\mathbf{R}^{n}$.

This section contains two results. The following result, proved in Section
A.3, provides the first localized in-expectation bound for the exponential
weights estimator under the assumption that the noise vector $Z$ is standard
normal. In particular, this yields an improvement over the classical _global_
PAC-Bayesian bound obtained in the same setup by Leung and Barron [44]. The
assumption that $Z$ is Gaussian, in particular, allows us to leverage Stein’s
Unbiased Risk Estimator (SURE) in the analysis of $\widehat{f}^{\mathsf{ew}}$.
Such an approach in the context of aggregation via exponential weights was
first considered by George [29], and later systematically applied by Leung and
Barron [44], Dalalyan and Tsybakov [22] and Dalalyan and Salmon [24], among
others. For a treatment of more general noise distributions, we refer to [23]
and references therein.

######  Theorem 1.

Consider the fixed design setting described above and assume that
$(Z_{i})_{i=1}^{n}$ are i.i.d. centered Gaussian random variables with
variance $\sigma^{2}$. Let $\\{f_{\theta}:\theta\in\Theta\\}$ be a family of
vectors in $\mathbf{R}^{n}$ and $\pi\in\mathcal{P}(\Theta)$ a fixed
probability distribution. Then, the exponential weights estimator
$\widehat{f}^{\mathsf{ew}}_{n}=\widehat{f}^{\mathsf{ew}}_{n}(\pi,\beta)$
defined in (13), with prior $\pi$ and inverse temperature
$\beta=n/(8\sigma^{2})$, satisfies

${\mathbf{E}}\Big{[}\|\widehat{f}^{\mathsf{ew}}_{n}-f^{\star}\|_{n}^{2}\Big{]}\leqslant\int_{\Theta}\lVert{f_{\theta}-f^{\star}}\rVert_{n}^{2}\,\pi_{-nR/(16\sigma^{2})}(\mathrm{d}\theta).$
(17)

Equivalently and with the notation (8) from the introduction, the bound (17)
writes

${\mathbf{E}}\big{[}R(\widehat{f}^{\mathsf{ew}}_{n})\big{]}\leqslant\mathsf{C}_{\mathsf{loc}}(n/(16\sigma^{2}),\pi).$

As already mentioned at the beginning of Section 2, the exponential weights
estimator is deviation suboptimal. On the other hand, it was shown in [21]
showed that the $Q$-aggregation estimator satisfies a _global_ exponential-
tail bound in the fixed design setting with sub-Gaussian noise; we say that a
centered real random variable $Z$ is $\sigma^{2}$-sub-Gaussian if for any
$\lambda\in\mathbf{R}$ we have ${\mathbf{E}}\exp(\lambda
Z)\leqslant\exp(\lambda^{2}\sigma^{2}/2)$. The following result shows that in
fact, the $Q$-aggregation estimator satisfied a deviation optimal _local_
bound. See Appendix A.4 for the proof.

######  Theorem 2.

Consider the fixed design setting and assume that the noise random variables
$Z_{i}$ are i.i.d., zero-mean and $\sigma^{2}$-sub-Gaussian. Let
$\\{f_{\theta}:\theta\in\Theta\\}$ be a family of vectors in $\mathbf{R}^{n}$
and $\pi\in\mathcal{P}(\Theta)$ a fixed probability distribution. Let
$\widehat{f}^{\mathsf{q}}_{n}=\widehat{f}^{\mathsf{q}}_{n}(\pi,\beta)$ denote
the $Q$-aggregation estimator defined in (14), with prior $\pi$ and inverse
temperature $\beta=n/(12\sigma^{2})$. Then, for any $\delta\in(0,1)$, with
probability at least $1-\delta$ it holds that

$\displaystyle\|\widehat{f}^{\mathsf{q}}_{n}-f^{\star}\|_{n}^{2}$
$\displaystyle\leqslant\int_{\Theta}\lVert{f_{\theta}-f^{\star}}\rVert_{n}^{2}\,\pi_{-nR/(36\sigma^{2})}(\mathrm{d}\theta)+\frac{18\sigma^{2}\log(1/\delta)}{n}.$
(18)

###  2.2 Random design

We now to to the _random design_ setting, where we observe $n$ input-output
pairs $(X_{i},Y_{i})$, $i=1,\dots,n$, each drawn independently from some
distribution $P$. Unlike in the fixed design setting considered in Section
2.1, we work under boundedness assumptions on both the response variable $Y$
and the functions in the class. Specifically, we assume in this section that,
for some constant $b>0$, one has $|Y|\leqslant b$ almost surely, and for every
$\theta\in\Theta$,

$\lVert{f_{\theta}}\rVert_{\infty}=\sup_{x\in\mathcal{X}}|f_{\theta}(x)|\leqslant
b.$ (19)

When working with arbitrary classes of functions, such boundedness conditions
on the response variable $Y$ and the reference functions are natural and have
been considered by many authors; see, e.g., [63, 7, 43]. Finally, the
population risk functional $R$ is defined by $R(f)={\mathbf{E}}_{(X,Y)\sim
P}\left[(f(X)-Y)^{2}\right],$ while its empirical version $R_{n}$ is defined
in (12). The following result in proved in Appendix A.5.

######  Theorem 3.

Consider the random design setup described above. Let
$\pi\in\mathcal{P}(\Theta)$ be any fixed probability distribution and let
$\beta=c_{1}n/b^{2}$ be the inverse temperature parameter. Let
$\widehat{f}^{\mathsf{q}}_{n}=\widehat{f}^{\mathsf{q}}_{n}(\pi,\beta)$ denote
the $Q$-aggregation estimator defined in (14). Then, for any $\delta\in(0,1)$,
with probability at least $1-\delta$ it holds that

$R(\widehat{f}^{\mathsf{q}}_{n})\leqslant\int_{\Theta}R(f_{\theta})\,\pi_{-c_{1}nR/(3b^{2})}(\mathrm{d}\theta)+\frac{c_{2}b^{2}\log(3/\delta)}{n},$
(20)

where $c_{1},c_{2}>0$ are universal constants.

Finally, it is important to note that the estimator
$\widehat{f}^{\mathsf{q}}_{n}$, within the framework of model aggregation
involving random design, has not been previously analyzed. The question of
analyzing the excess risk associated with this particular estimator was posed
in [51].

##  3 Applications

In this section, we instantiate the local bounds in several different
contexts.

###  3.1 Model aggregation

As discussed in the introduction, a problem where
$\Theta=\\{\theta_{1},\dots,\theta_{M}\\}$ is a finite set is called _model
aggregation_. While the results we are about to discuss are applicable in all
three learning setups considered in Section 2, for the sake of concreteness,
let us consider the random design setup described in Section 2.2. In
particular, the response variable $Y$ is almost surely contained in $[-b,b]$
and the reference class of functions satisfies the boundedness condition (19).

Let $\widehat{f}=\widehat{f}(S_{n})$ denote the output of any statistical
estimator, where $S_{n}=(X_{i},Y_{i})_{i=1}^{n}$ is the observed data sample.
In the model aggregation problem, we wish to obtain bounds on

$R(\widehat{f})-\min_{\theta\in\Theta}R(f_{\theta}).$

The following minimax lower bound is well-known [63]:

$\inf_{\widehat{f}}\sup_{P,|\Theta|=M}{\mathbf{E}}R(\widehat{f})-\min_{\theta\in\Theta}R(f_{\theta})\geqslant
c\frac{b^{2}\log M}{n},$ (21)

where $c>0$ is a universal constant, the expectation is taken with respect to
$S_{n}$, the infimum is taken over all statistical estimators $\widehat{f}$,
and the supremum is taken over all data-generating distributions and reference
classes of functions of size $M$ that respect the above-described constraints.

Lecué and Mendelson [42] asked whether it is possible to prove an upper bound
for the model aggregation problem that exhibits some form of “instance
adaptivity”, that is, adaptivity to the hardness of the specific problem at
hand. This adaptivity should be quantitatively reflected in a risk bound that
does not depend on the full size of the reference set $\Theta$, but only on
the number of functions that are close to the best function in the class. The
proposition stated below provides such a guarantee. We remark that partial
results were obtained in [42] under the so-called Bernstein assumption.

######  Proposition 1.

Consider the random design model aggregation setup described above and let
$\theta^{\star}\in\mathop{\mathrm{arg}\,\mathrm{min}}_{j=1,\dots,M}R(\theta_{j})$.
Let $\pi\in\mathcal{P}(\Theta)$ be the uniform distribution over the finite
set $\Theta$ and set the inverse temperature parameter $\beta=n/(576b^{2})$.
Let $\widehat{f}^{\mathsf{q}}_{n}=\widehat{f}^{\mathsf{q}}_{n}(\pi,\beta)$
denote $Q$-aggregation estimator defined in (14). Then, for any
$\delta\in(0,1)$, with probability at least $1-\delta$ it holds that

$R(\widehat{f}^{\mathsf{q}}_{n})-R(f_{\theta^{\star}})\leqslant
c_{1}\frac{b^{2}}{n}\bigg{[}\log\bigg{(}\sum_{j=1}^{M}\exp\left(-c_{2}\frac{n}{b^{2}}\big{(}R(\theta_{j})-R(\theta^{\star})\big{)}\right)\bigg{)}+\log(1/\delta)\bigg{]},$
(22)

where $c_{1},c_{2}>0$ are universal constants.

###### Proof.

The result follows from Theorem 3 by upper bounding the local complexity using
the inequality (10) and by plugging in $\gamma=\delta_{\theta^{\star}}$, a
Dirac mass at $\theta^{\star}$, into the resulting bound. ∎

###  3.2 Gaussian priors

In this section, we consider the case where the prior $\pi$ on functions
$f:\mathcal{X}\to\mathbf{R}$ is a centered Gaussian process. In order to
simplify notation, we assume that $\pi$ is supported on a finite-dimensional
space of functions; this is not critical to our analysis, which would apply to
the infinite-dimensional case with minor changes in terminology. Furthermore,
up to a change the of variable $x$, we may assume that
$\mathcal{X}=\mathbf{R}^{d}$ and $\mathcal{F}$ is the class of linear
functions $f_{\theta}(x)=\langle\theta,x\rangle$ for
$\theta\in\mathbf{R}^{d}$; up to a linear transformation, we may assume that
$\pi=\mathcal{N}(0,\gamma^{-1}I_{d})$ for some $\gamma>0$.

We first introduce some notation for both fixed and random design. The sample
covariance matrix is
$\widehat{\Sigma}_{n}=n^{-1}\sum_{i=1}^{n}X_{i}X_{i}^{\mathsf{T}}$, while for
$\lambda>0$ the _ridge estimator_ $\widehat{\theta}_{\lambda}$ is:

$\widehat{\theta}_{\lambda}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\theta\in\mathbf{R}^{d}}\big{\\{}R_{n}(f_{\theta})+\lambda\lVert{\theta}\rVert^{2}\big{\\}}=(\widehat{\Sigma}_{n}+\lambda
I_{d})^{-1}\cdot\frac{1}{n}\sum_{i=1}^{n}Y_{i}X_{i}.$

We now determine the $Q$-aggregation estimator under Gaussian prior for the
linear class.

######  Proposition 2.

Under the prior $\pi=\mathcal{N}(0,\gamma^{-1}I_{d})$ and with inverse
temperature $\beta/2$, the $Q$-aggregation posterior is
$\mathcal{N}(\widehat{\theta}_{\lambda},\beta^{-1}(\frac{1}{2}\widehat{\Sigma}_{n}+\lambda
I_{d})^{-1})$ while the exponential weights posterior is
$\mathcal{N}(\widehat{\theta}_{\lambda},\beta^{-1}(\widehat{\Sigma}_{n}+\lambda
I_{d})^{-1})$, where $\lambda=\gamma/\beta$. The associated prediction
functions both coincide with that of the ridge regression estimator
$f_{\widehat{\theta}_{\lambda}}$.

Exponential weights is well-known in this context, and can be computed
directly thanks to its explicit expression. In contrast, the solution to the
$Q$-aggregation optimization problem does not have a closed form in general.
In this case, the solution is obtained by first showing that it is Gaussian
(using Proposition 5), and then solving the optimization problem among
Gaussian posteriors.

Next, we express and compare the global and local complexities.

######  Proposition 3.

Assume that $\mathcal{X}=\mathbf{R}^{d}$ and let

$\mathcal{F}=\\{x\mapsto\langle\theta,x\rangle:\theta\in\mathbf{R}^{d}\\}$

be the class of linear functions. Set
$\theta_{\lambda}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\theta\in\mathbf{R}^{d}}\\{R(\theta)+\lambda\lVert{\theta}\rVert^{2}\\}$,
and $\Sigma={\mathbf{E}}[XX^{\mathsf{T}}]$ (assumed finite) for random design,
and $\Sigma=\widehat{\Sigma}_{n}=n^{-1}\sum_{i=1}^{n}X_{i}X_{i}^{\mathsf{T}}$
for fixed design. Set $\pi=\mathcal{N}(0,\gamma^{-1}I_{d})$. Then, for any
$\beta>0$, letting $\lambda=\gamma/\beta$,

$\displaystyle\mathsf{C}_{\mathsf{glob}}(\beta/2,\pi)$
$\displaystyle=R(\theta_{\lambda})+\lambda\lVert{\theta_{\lambda}}\rVert^{2}+\frac{1}{\beta}\log\det\big{(}I_{d}+\lambda^{-1}\Sigma\big{)}$
(23) $\displaystyle\mathsf{C}_{\mathsf{loc}}(\beta/2,\pi)$
$\displaystyle=R(\theta_{\lambda})+\frac{1}{\beta}\cdot\operatorname{Tr}[(\Sigma+\lambda
I_{d})^{-1}\Sigma].$ (24)

In addition, assuming for the second inequality below that
$0<\lambda\leqslant\|{\Sigma}\|_{\mathrm{op}}$,

$\operatorname{Tr}\big{[}(\Sigma+\lambda
I_{d})^{-1}\Sigma\big{]}\leqslant\log\det\big{(}I_{d}+\lambda^{-1}\Sigma\big{)}\leqslant
2\log\bigg{(}1+\frac{\|{\Sigma}\|_{\mathrm{op}}}{\lambda}\bigg{)}\cdot\operatorname{Tr}\big{[}(\Sigma+\lambda
I_{d})^{-1}\Sigma\big{]}.$ (25)

We note that the global entropic complexity already appears in the seminal
work of [69] in online linear regression, where it bounds the regret of the
so-called Vovk-Azoury-Warmuth forecaster. Here, the local complexity improves
the global complexity in two ways. First, it replaces the log determinant by a
trace, allowing in particular to take $\lambda\to 0$; this allows to remove a
$\log n$ factor as $n\to\infty$, and also to obtain bounds for arbitrary large
norms of the comparison parameter $\theta$ (by taking $\lambda$ small enough
that $\lambda\lVert{\theta}\rVert^{2}\ll d/n$). Second, it removes an
approximation term of $\lambda\lVert{\theta_{\lambda}}\rVert^{2}$, which
allows to obtain better rates in some regimes (see, e.g., [14]). We also note
that in the context of robust estimation, versions of local entropic
complexities have been computed in [10, 9].

Theorems 1 and 2 in the fixed design setting apply to the present context, and
provide bounds for the ridge regression estimator featuring the local
complexity of Proposition 3. However, the random-design result does not apply
to the Gaussian setting, since it requires boundedness of the base functions.
In Theorem 4, we indicate a different approach giving bounds for random
design.

###  3.3 Random-design bounds for ridge-type improper estimators

In this section, we consider an alternative approach to obtaining excess risk
bounds featuring related complexities as the one from the previous section, in
the random design setting.

At this point, it is worth discussing some difficulties that are specific to
the random design case. Unlike for fixed design, we now aim to predict for a
new, unseen point $X$. This difference has the technical implication that,
when analyzing an estimator like ridge regression, which is based on empirical
risk minimization and corresponds to $Q$-aggregation as per Proposition 2, we
must relate two distinct covariance matrices: the sample covariance matrix
$\widehat{\Sigma}_{n}$ (which measures the closeness of predictions of
different parameters on the sample) and the population covariance $\Sigma$
(which controls closeness on a new random sample $X$). To obtain guarantees of
this nature, it is not hard to show that some assumptions on the distribution
of $X$ are required beyond just the existence of $\Sigma$, such as controlled
moments of order $p>2$.

Here, a key point is that precise excess risk bounds can nonetheless be
obtained without additional assumptions on the distribution of $X$, provided
that different estimators than ridge regression are considered; we refer to
[31, 27, 63, 55] for more information on this topic, in the case of linear
regression (that is, with guarantees uniform over the linear class).

An important aspect of the problem at hand is that, in order to obtain
nontrivial guarantees without additional assumptions on the distribution of
$X$, one must consider aggregation procedures that return a function outside
the class, that is, a nonlinear function [60, 55]. While this is similar to
the need for aggregation rules for finite classes discussed in the
introduction, it calls for some explanation. Indeed, for finite classes, the
need for aggregation rules arises due to nonconvexity of the class; in
contrast, the class of linear functions is convex. In this case, the need for
aggregation methods stems from the unbounded nature of linear functions,
leading to large errors when $X$ distribution can take large values with small
probabilities [55]. Hence, the guarantees of Section 2.2 do not apply as
Theorem 2 needs bounded functions in the class. To ensure boundedness, we can
truncate the original functions [31, 55], projecting their predictions to a
bounded interval. However, this operation makes the resultant truncated linear
functions class nonconvex, even if the original class was convex.

A natural approach to obtain guarantees is then to perform $Q$-aggregation on
the class of _truncated_ linear functions, under a Gaussian prior on
parameters. One could then apply Theorem 3, since its boundedness condition
does apply. However, this approach has two main limitations. First, the
resulting estimator involves a non-convex minimization problem, whose
computational complexity could depend exponentially on the dimension. Second,
it is unclear that the resulting procedure would satisfy a risk bound
featuring the local complexity from Proposition 3. Indeed, if one applies
Theorem 2, the resulting bound would depend on the local entropic complexity
(under Gaussian prior) of _truncated_ linear functions. Assuming that
$|Y|\leqslant b$, truncating linear functions at $b$ only reduces their risk,
and thus also the associated global entropic complexity. Hence, a bound based
on the global complexity of Proposition 3 holds for this method. In contrast,
replacing functions by functions with smaller risk does not necessarily
decrease the _local_ complexity, and as such the local entropic complexity of
truncated linear functions could be larger than that of linear functions.

We now describe an alternative approach, which partly circumvents these
limitations. This approach is based on leave-one-out and exchangeability
considerations, and has been used to design aggregation procedures for linear
regression in [27], and for density estimation and logistic regression in
[54]. This approach allows to obtain bounds of “local” nature, without the
additional logarithmic factors associated with “global” bounds. In addition,
it generally leads to estimators that are simple efficiently computable. On
the other hand, it only leads to guarantees in expectation rather than in
deviation [55]. In addition, the resulting bound is slightly larger than the
local entropic complexity of Proposition 3.

We now describe the considered procedures and their guarantees. Recall that
the _ridge_ estimator is given by

$\widehat{\theta}_{\lambda}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\theta\in\mathbf{R}^{d}}\big{\\{}R_{n}(\theta)+\lambda\lVert{\theta}\rVert^{2}\big{\\}}=(\widehat{\Sigma}_{n}+\lambda
I_{d})^{-1}\cdot\frac{1}{n}\sum_{i=1}^{n}Y_{i}X_{i}.$ (26)

For any $\lambda>0$ and $x\in\mathbf{R}^{d}$, define the _ridge leverage
score_ of $x\in\mathbf{R}^{d}$ by

$h_{\lambda}(x)=\bigg{\langle}\bigg{(}\sum_{i=1}^{n}X_{i}X_{i}^{\mathsf{T}}+\lambda
nI_{d}+xx^{\mathsf{T}}\bigg{)}^{-1}x,x\bigg{\rangle},$ (27)

and define (letting $\lambda^{\prime}=(1+1/n)\lambda$)

$\widehat{f}_{\lambda}(x)=(1-h_{\lambda^{\prime}}(x))^{2}\cdot\langle\widehat{\theta}_{\lambda^{\prime}},x\rangle.$
(28)

This estimator is a ridge-regularized analogue of the Forster-Warmuth
procedure introduced in [27], which corresponds to the case $\lambda=0$. We
therefore call it _ridge FW_. It can also be seen a linear version of the
ridge-regularized logistic regression procedure from [55]. Also, for $b>0$
define the truncation function $\psi_{b}:\mathbf{R}\to\mathbf{R}$ by
$\psi_{b}(x)=\max(-b,\min(b,x))$ for $x\in\mathbf{R}$. Then, let

$\widehat{g}_{\lambda,b}(x)=\psi_{b}\big{(}\langle\widehat{\theta}_{\lambda^{\prime}},x\rangle\big{)}.$
(29)

This procedure is a ridge-regularized version of the truncated least squares
estimator considered in [31, 55].

These estimators satisfy the following excess risk bounds:

######  Theorem 4.

For any distribution $P$ on $\mathbf{R}^{d}\times\mathbf{R}$ such that
$|Y|\leqslant b$ almost surely and $\Sigma={\mathbf{E}}XX^{\mathsf{T}}$
exists,

$\displaystyle{\mathbf{E}}\big{[}R(\widehat{f}_{\lambda})\big{]}$
$\displaystyle\leqslant\inf_{\theta\in\mathbf{R}^{d}}\big{\\{}R(\theta)+\lambda\lVert{\theta}\rVert^{2}\big{\\}}+2m^{2}\,\frac{\operatorname{Tr}[(\Sigma+\lambda
I_{d})^{-1}\Sigma]}{n+1}\,;$ (30)
$\displaystyle{\mathbf{E}}\big{[}R(\widehat{g}_{\lambda,b})\big{]}$
$\displaystyle\leqslant\inf_{\theta\in\mathbf{R}^{d}}\big{\\{}R(\theta)+\lambda\lVert{\theta}\rVert^{2}\big{\\}}+8m^{2}\,\frac{\operatorname{Tr}[(\Sigma+\lambda
I_{d})^{-1}\Sigma]}{n+1}.$ (31)

In addition, let $\widehat{b}=\max_{1\leqslant i\leqslant n}|Y_{i}|$ and
$\widehat{g}_{\lambda}(x)=\widehat{g}_{\lambda,\widehat{b}}(x)$. Then, we have

${\mathbf{E}}\big{[}R(\widehat{g}_{\lambda})\big{]}\leqslant\inf_{\theta\in\mathbf{R}^{d}}\big{\\{}R(\theta)+\lambda\lVert{\theta}\rVert^{2}\big{\\}}+8b^{2}\,\frac{\operatorname{Tr}[(\Sigma+\lambda
I_{d})^{-1}\Sigma]}{n+1}+\frac{b^{2}}{n+1}.$ (32)

The three estimators $\widehat{f}_{\lambda},\widehat{g}_{\lambda,b}$ and
$\widehat{g}_{\lambda}$ achieve the same risk bound up to constant factors.
The main differences are:

  * •

$\widehat{f}_{\lambda}$ is more expensive computationally than truncated ridge
at prediction time. Indeed, truncated ridge only requires to store the vector
$\widehat{\theta}_{\lambda^{\prime}}\in\mathbf{R}^{d}$, which takes $O(d)$
memory (in the case of infinite- or high-dimensional reproducing kernel
Hilbert spaces, using the so-called kernel trick one can replace $d$ by $n$ in
all this discussion); computing predictions involves a scalar product, which
takes $O(d)$ time. In contrast, $\widehat{f}_{\lambda}$ requires storing both
$\widehat{\theta}_{\lambda^{\prime}}$ and the precision matrix
$S_{\lambda}=(\sum_{i=1}^{n}X_{i}X_{i}^{\mathsf{T}}+\lambda(n+1)I_{d})^{-1}$,
which takes $O(d^{2})$ memory. Computation of the leverage
$h_{\lambda^{\prime}}(x)=\langle S_{\lambda}x,x\rangle/(1+\langle
S_{\lambda}x,x\rangle)$ (by Sherman-Morrison’s identity) involves computing
$\langle S_{\lambda}x,x\rangle$, which takes $O(d^{2})$ time.

  * •

On the other hand, from a statistical perspective, $\widehat{g}_{\lambda,b}$
requires prior knowledge of the parameter $b$, while $\widehat{f}_{\lambda}$
is adaptive to $b$. This shortcoming of $\widehat{g}_{\lambda,b}$ is addressed
by the adaptive truncated estimator $\widehat{g}_{\lambda}$.

In fact, one could replace $\operatorname{Tr}[(\Sigma+\lambda
I_{d})^{-1}\Sigma]$ by the (smaller) quantity
${\mathbf{E}}\operatorname{Tr}[(\widehat{\Sigma}_{n+1}+\lambda)^{-1}\widehat{\Sigma}_{n+1}]$,
which is always defined without any assumption on the distribution of $X$. We
stated the bound in terms of the population covariance since it is more
explicit. In any case, it is worth noting that this bound is true under only a
second moment assumption on $X$.

We also mention that a remarkable procedure for linear regression in the
sequential setting was proposed by Vovk [69] and Azoury and Warmuth [11].
Combining its regret guarantee with the conversion described in the
introduction, the averaged version $\widehat{g}_{\lambda}$ of this procedure
satisfies the following risk bound:

${\mathbf{E}}R(\widehat{g}_{\lambda})\leqslant\inf_{\theta\in\mathbf{R}^{d}}\big{\\{}R(\theta)+\lambda\lVert{\theta}\rVert^{2}\big{\\}}+\frac{b^{2}\cdot\log\det(\lambda^{-1}\Sigma+I_{d})}{n+1}.$

Compared to this bound, the guarantees from Theorem 4 allow to replace the
log-determinant by a smaller trace term.

###  3.4 Aggregation of classifiers and almost exchangeable priors

In the classic paper [36], the following question is presented: Given a finite
set $\mathcal{F}$ comprised of $\\{0,1\\}$-valued functions and assuming that
$y_{i}\in\\{0,1\\}$, how can one demonstrate the optimal sequential regret
bound (1) for a particular loss function with curvature as a function of $n$
and the cardinality of $\mathcal{F}$? An important particular case of their
main result is the squared loss function, which is also the subject of our
focus. Another interesting analysis of classes of binary-valued functions in
the context of statistical aggregation appears in [58]. To present their
result, we first need some standard definitions. Let
$\mathcal{F}=\\{f_{\theta}:\theta\in\Theta\\}$ be a class of
$\\{0,1\\}$-valued functions with the domain $\mathcal{X}$. The Vapnik-
Chervonenkis (VC) dimension of $\mathcal{F}$ is the largest integer $d$ such
that

$\max\limits_{x_{1},\ldots,x_{d}\in\mathcal{X}}\left|\\{(f(x_{1}),\ldots,f(x_{d})):f\in\mathcal{F}\\}\right|=2^{d}.$

Classes of $\\{0,1\\}$-valued functions with a finite VC dimension will be
referred to as VC classes. Recall that for the squared loss we have
$R(f)={\mathbf{E}}[(f(X)-Y)^{2}]$. The result in [58, Theorem 3 and Remark 4]
can be stated as follows. For any VC class $\mathcal{F}$ whose VC dimension if
equal to $d$ there is an aggregation procedure $\widehat{f}_{\textrm{RST}}$
such that, with probability at least $1-\delta$,

$R(\widehat{f}_{\textrm{RST}})-\inf_{f\in\mathcal{F}}R(f)\leqslant
c\left(\frac{d}{n}\log\left(\frac{en}{d}\right)+\frac{1}{n}\log\left(\frac{1}{\delta}\right)\right),$
(33)

where $c>0$ is some absolute constant provided that $n\geqslant d$. The
primary difference from binary classification is the use of squared loss
instead of the standard binary loss. This aligns with the work [36]. As a
result, one can get the fast for the improper estimator
$\widehat{f}_{\textrm{RST}}$ without making any low noise/margin assumptions.
Rakhlin, Sridharan and Tsybakov [58] also present an example of VC class
showing that this bound cannot be improved (including the
$\log\left(en/d\right)$ term) up to multiplicative constant factors. That is,
at least in some cases, the inequality (33) can be reversed (but with a
different absolute constant).

Recent efforts in the field of binary classification have focused on improving
the $\log\left(en/d\right)$ factor in the bounds. The authors of [34] and [75]
leveraged additional complexity measures of VC classes to capture the refined
local structure of the classes. A key complexity measure discussed in these
papers is the _star number_ , which is defined as the largest integer $s$ such
that there exists a subset $S=\\{x_{1},\dots,x_{s}\\}\subseteq\mathcal{X}$ and
functions $f_{0},f_{1},\dots,f_{s}\in\mathcal{F}$ such that for all
$i=1,\ldots,s$, $f_{0}$ and $f_{i}$ disagree on a unique point in $S$. That
is, for all $i=1,\ldots,s$, we have $\\{x_{i}\\}=\\{x\in S:f_{0}(x)\neq
f_{i}(x)\\}$. In binary classification, it is shown in [34, 75] that the term
$\log\left(en/d\right)$ can usually be replaced in the risk bounds with
$\log\left(\min\\{en/d,s\\}\right)$ if the star number of $\mathcal{F}$ is
finite. A natural question arises as to whether the star number, which
characterizes the local structure of VC classes, can be applied in the context
of bound (33). It appears that the methods of this paper allow us to do so. In
order to provide such a setup, we extend our analysis to handle the so-called
_transductive algorithms_. This approach was proposed by Vapnik [65] and has
been demonstrated to have useful connections with PAC-Bayesian arguments in
subsequent works [17, 6]. In the simplest form of the transductive setup, the
learner works with an i.i.d sample of size $2n$. Among these points, the first
$n$ points have their labels $Y_{i}$, while for the remaining $n$ points we
only observe the objects $X_{i}$, but not their labels. We are now aiming to
prove the high probability upper bound on the empirical excess risk computed
with respect to the second part of the sample:

$R^{\prime}_{n}(\widehat{f})-\min\limits_{f\in\mathcal{F}}R^{\prime}_{n}(f)=\frac{1}{n}\bigg{(}\sum_{i=n+1}^{2n}(\widehat{f}(X_{i})-Y_{i})^{2}-\min\limits_{f\in\mathcal{F}}\sum_{i=n+1}^{2n}(f(X_{i})-Y_{i})^{2}\bigg{)},$

where the estimator $\widehat{f}$ depends on $n$ labeled instances and all
$2n$ unlabeled points. The idea is to replace the population risk $R(f)$ by
its empirical counterpart $R^{\prime}_{n}(f)$ computed with respect to the
second half of the sample. Furthermore, for the inverse temperature parameter
$\beta>0$, we redefine the (data-dependent) localized prior as follows:

$\pi_{\beta(R_{n}+R^{\prime}_{n})}(\mathrm{d}\theta)=\frac{\exp(-\beta(R_{n}(f_{\theta})+R^{\prime}_{n}(f_{\theta})))}{\int_{\Theta}\exp(-\beta(R_{n}(f_{\vartheta})+R^{\prime}_{n}(f_{\vartheta})))\pi(\mathrm{d}\vartheta)}.$
(34)

The next theorem is a direct extension of Theorem 3 to the transductive setup.

######  Theorem 5 (A local bound with almost exchangeable priors).

Consider the transductive classification setup discussed above. Assume that
the boundedness assumption (19) holds. Let $\pi\in\mathcal{P}(\Theta)$ be any
probability distributions that may depend on the sample of size $2n$. Let
$\widehat{f}^{\mathsf{q}}_{n}$ be the output of the $Q$-aggregation estimator
defined in (14) with inverse temperature $\beta=c_{1}n/b^{2}$. Then, for any
$\delta\in(0,1)$, with probability at least $1-\delta$, it holds that

$R_{n}^{\prime}(\widehat{f}^{\mathsf{q}}_{n})\leqslant\inf_{\gamma\in\mathcal{P}(\Theta)}\left\\{\langle\gamma,R^{\prime}_{n}\rangle+\frac{3b^{2}}{2c_{1}n}\mathrm{KL}(\gamma,\pi_{-\frac{c_{1}n}{6b^{2}}(R_{n}+R^{\prime}_{n})})\right\\}+c\frac{b^{2}\log(3/\delta)}{n}.$
(35)

where $c,c_{1}>0$ are universal constants.

The fact that our localized prior distribution now depends on the entire
sample of size $2n$ makes the setup of Theorem 5 closely related to the idea
of _almost exchangeable priors_ ††The prior distribution $\pi$ depends on the
entire sample of size $2n$, while two part of the sample are exchangeable due
to independence, hence the terminology. introduced in [17, 6] and to
_conditional mutual information_ appearing in the recent works [61, 30, 32].
Our focus is on the aggregation problem with squared loss without the
realizability assumption, whereas prior work primarily addressed the standard
binary classification setup.

We now demonstrate the application of Theorem 5 in the context of binary
classification. Consider a VC class $\mathcal{F}$. Given an i.i.d. sample
$S_{2n}=\\{X_{1},\ldots,X_{2n}\\}$, we first consider the (random) projection
class
$\mathcal{F}|_{S_{2n}}=\\{(f(X_{1}),\ldots,f(X_{2n})):f\in\mathcal{F}\\}$
sample of unlabeled instances and choose the prior distribution $\pi_{S_{2n}}$
to be the uniform measure on $\mathcal{F}_{S_{2n}}$. When working with our
$Q$-Aggregation estimator, we identify $\mathcal{F}|_{S_{2n}}$ with the set
$\Theta$. Observe that $\mathcal{F}|_{S_{2n}}$ is a finite class, so that our
uniform distribution is correctly defined.

######  Proposition 4.

Consider the transductive setup of Theorem 5 and let $\mathcal{F}$ be a class
of $\\{0,1\\}$-valued functions with VC dimension $d$ and star number $s$. Let
$\pi=\pi(S_{2n})$ be a uniform measure on the projection class
$\mathcal{F}|_{S_{2n}}$, where $S_{2n}$ denotes the full sample of unlabeled
points of size $2n$. Let $\widehat{f}^{\mathsf{q}}_{n}$ be the output of the
$Q$-aggregation estimator defined in (14) with the inverse temperature
parameter $\beta=c_{1}n$ and prior $\pi$. Then,

$R^{\prime}_{n}(\widehat{f}^{\mathsf{q}}_{n})-\min\limits_{f\in\mathcal{F}}R^{\prime}_{n}(f)\leqslant
c\left\\{\frac{d}{n}\log\left[es\left(\frac{nR(f^{\star})}{d}+\log\left(\frac{en}{d}\right)\right)\right]+\frac{\log(4/\delta)}{n}\right\\},$

where $c,c_{1}>0$ are absolute constants and
$R(f^{\star})=\inf_{f\in\mathcal{F}}R(f)$.

Our result provides the first improvement over the bound (33) in [58], but
only in the transductive setting. An interesting regime is when
$R(f^{\star})\lesssim\frac{s^{\kappa}}{n}$, where $\kappa>1$ is some absolute
constant. The reason for including $R(f^{\star})$ in the logarithmic factor is
that we need to account for the possibility of multiple minimizers of the
population risk. This distinguishes our setup from all previous classification
bounds that involved the star number, where the low-noise assumption implied
the existence of a unique risk minimizer, making the local structure of the
class significantly simpler.

##### Acknowledgments.

Jaouad Mourtada and Nikita Zhivotovskiy acknowledge the support of the France-
Berkeley Fund.

## References

  * [1] P. Alquier.  PAC-Bayesian bounds for randomized empirical risk minimizers.  Mathematical Methods of Statistics, 17:279–304, 2008.
  * [2] P. Alquier.  User-friendly introduction to PAC-Bayes bounds.  arXiv preprint arXiv:2110.11216, 2021.
  * [3] M. Anthony and P. L. Bartlett.  Neural Network Learning: Theoretical Foundations.  Cambridge University Press, 2009.
  * [4] N. Aronszajn.  Theory of reproducing kernels.  Transactions of the American Mathematical Society, 68(3):337–404, 1950.
  * [5] J.-Y. Audibert.  Aggregated estimators and empirical complexity for least square regression.  Annales de l’Institut Henri Poincaré, Probabilités et Statistiques, 40(6):685–736, 2004.
  * [6] J.-Y. Audibert.  PAC-Bayesian statistical learning theory.  PhD thesis, Université Paris VI, 2004.
  * [7] J.-Y. Audibert.  Progressive mixture rules are deviation suboptimal.  In Advances in Neural Information Processing Systems 20, pages 41–48, 2008.
  * [8] J.-Y. Audibert.  Fast learning rates in statistical inference through aggregation.  The Annals of Statistics, 37(4):1591–1646, 2009.
  * [9] J.-Y. Audibert and O. Catoni.  Linear regression through PAC-Bayesian truncation.  Preprint arXiv:1010.0072, 2010.
  * [10] J.-Y. Audibert and O. Catoni.  Robust linear least squares regression.  The Annals of Statistics, 39(5):2766–2794, 2011.
  * [11] K. S. Azoury and M. K. Warmuth.  Relative loss bounds for on-line density estimation with the exponential family of distributions.  Machine Learning, 43(3):211–246, 2001.
  * [12] A. R. Barron.  Are Bayes rules consistent in information?  In Open Problems in Communication and Computation, pages 85–91. Springer, 1987.
  * [13] P. C. Bellec.  Optimal bounds for aggregation of affine estimators.  The Annals of Statistics, 46(1):30–59, 2018.
  * [14] G. Blanchard and N. Mücke.  Optimal rates for regularization of statistical inverse learning problems.  Foundations of Computational Mathematics, 18(4):971–1013, 2018\.
  * [15] E. Carlen.  Trace inequalities and quantum entropy: an introductory course.  Entropy and the quantum, 529:73–140, 2010.
  * [16] O. Catoni.  Statistical Learning Theory and Stochastic Optimization: Ecole d’Eté de Probabilités de Saint-Flour XXXI - 2001, volume 1851 of Lecture Notes in Mathematics.  Springer-Verlag, 2004.
  * [17] O. Catoni.  PAC-Bayesian Supervised Classification: The Thermodynamics of Statistical Learning, volume 56 of IMS Lecture Notes Monograph Series.  Institute of Mathematical Statistics, 2007.
  * [18] N. Cesa-Bianchi, A. Conconi, and C. Gentile.  On the generalization ability of on-line learning algorithms.  IEEE Transactions on Information Theory, 50(9):2050–2057, 2004\.
  * [19] N. Cesa-Bianchi and G. Lugosi.  Prediction, Learning, and Games.  Cambridge University Press, Cambridge, New York, USA, 2006.
  * [20] T. M. Cover and J. A. Thomas.  Elements of Information Theory.  Wiley Series in Telecommunications and Signal Processing. Wiley-Interscience, New York, USA, 2nd edition, 2006.
  * [21] D. Dai, P. Rigollet, and T. Zhang.  Deviation optimal learning using greedy $Q$-aggregation.  The Annals of Statistics, 40(3):1878–1905, 2012.
  * [22] A. Dalalyan and A. B. Tsybakov.  Aggregation by exponential weighting, sharp PAC-Bayesian bounds and sparsity.  Machine Learning, 72(1):39–61, 2008.
  * [23] A. S. Dalalyan.  Simple proof of the risk bound for denoising by exponential weights for asymmetric noise distributions.  arXiv preprint arXiv:2212.12950, 2022.
  * [24] A. S. Dalalyan and J. Salmon.  Sharp oracle inequalities for aggregation of affine estimators.  The Annals of Statistics, 40(4):2327–2355, 2012.
  * [25] G. K. Dziugaite and D. M. Roy.  Computing nonvacuous generalization bounds for deep (stochastic) neural networks with many more parameters than training data.  In Proceedings of Uncertainty in Artificial Intelligence, 2017\.
  * [26] T. V. Erven, W. M. Koolen, and D. Van der Hoeven.  Metagrad: Adaptation using multiple learning rates in online learning.  Journal of Machine Learning Research, 22(161):1–61, 2021.
  * [27] J. Forster and M. K. Warmuth.  Relative expected instantaneous loss bounds.  Journal of Computer and System Sciences, 64(1):76–102, 2002.
  * [28] D. P. Foster.  Prediction in the worst case.  The Annals of Statistics, 19:1084–1090, 1991.
  * [29] E. I. George.  Minimax multiple shrinkage estimation.  The Annals of Statistics, 14(1):188–205, 1986.
  * [30] P. Grünwald, T. Steinke, and L. Zakynthinou.  PAC-Bayes, MAC-Bayes and conditional mutual information: Fast rate bounds that handle general VC classes.  In Conference on Learning Theory, pages 2217–2247, 2021.
  * [31] L. Györfi, M. Kohler, A. Krzyzak, and H. Walk.  A distribution-free theory of nonparametric regression.  Springer, 2002.
  * [32] M. Haghifam, G. K. Dziugaite, S. Moran, and D. Roy.  Towards a unified information-theoretic framework for generalization.  Advances in Neural Information Processing Systems, 34:26370–26381, 2021.
  * [33] J. Hannan.  Approximation to Bayes risk in repeated play.  Contributions to the Theory of Games, 3:97–139, 1957.
  * [34] S. Hanneke.  Refined error bounds for several learning algorithms.  Journal of Machine Learning Research, 17(1):4667–4721, 2016.
  * [35] S. Hanneke and L. Yang.  Minimax analysis of active learning.  Journal of Machine Learning Research, 16(1):3487–3602, 2015.
  * [36] D. Haussler, J. Kivinen, and M. K. Warmuth.  Sequential prediction of individual sequences under general loss functions.  IEEE Transactions on Information Theory, 44(5):1906–1925, 1998\.
  * [37] A. Juditsky, P. Rigollet, and A. B. Tsybakov.  Learning by mirror averaging.  The Annals of Statistics, 36(5):2183–2206, 2008.
  * [38] V. Kanade, P. Rebeschini, and T. Vaskevicius.  Exponential tail local Rademacher complexity risk bounds without the Bernstein condition.  arXiv preprint arXiv:2202.11461, 2022.
  * [39] V. Koltchinskii.  Oracle Inequalities in Empirical Risk Minimization and Sparse Recovery Problems, volume 2033 of École d’Été de Probabilités de Saint-Flour.  Springer-Verlag Berlin Heidelberg, 2011.
  * [40] W. M. Koolen and T. Van Erven.  Second-order quantile methods for experts and combinatorial games.  In Conference on Learning Theory, pages 1155–1175, 2015.
  * [41] G. Lecué and S. Mendelson.  Aggregation via empirical risk minimization.  Probability Theory and Related Fields, 145(3):591, 2009.
  * [42] G. Lecué and S. Mendelson.  On the optimality of the aggregate with exponential weights for low temperatures.  Bernoulli, 19(2):646–675, 2013.
  * [43] G. Lecué and P. Rigollet.  Optimal learning with Q-aggregation.  The Annals of Statistics, 42(1):211–224, 2014.
  * [44] G. Leung and A. R. Barron.  Information theory and mixing least-squares regressions.  IEEE Transactions on Information Theory, 52(8):3396–3410, 2006\.
  * [45] T. Liang, A. Rakhlin, and K. Sridharan.  Learning with square loss: localization through offset Rademacher complexity.  In Proceedings of The 28th Conference on Learning Theory, pages 1260–1285, 2015.
  * [46] N. Littlestone.  From on-line to batch learning.  In Proceedings of the 2nd annual workshop on Computational Learning Theory, pages 269–284. Morgan Kaufmann Publishers Inc., 1989.
  * [47] G. Lugosi and G. Neu.  Online-to-PAC conversions: Generalization bounds via regret analysis.  arXiv preprint arXiv:2305.19674, 2023.
  * [48] P. Massart.  Concentration inequalities and model selection. Ecole d’Eté de Probabilités de Saint-Flour XXXIII, volume 1896 of Lecture Notes in Mathematics.  Springer, 2007.
  * [49] D. A. McAllester.  Some PAC-Bayesian theorems.  Machine Learning, 37(3):355–363, 1999.
  * [50] D. A. McAllester.  PAC-Bayesian stochastic model selection.  Machine Learning, 51(1):5–21, 2003.
  * [51] N. Mehta.  Fast rates with high probability in exp-concave statistical learning.  In Proceedings of the 20th International Conference on Artificial Intelligence and Statistics, pages 1085–1093, 2017.
  * [52] S. Mendelson.  On aggregation for heavy-tailed classes.  Probability Theory and Related Fields, 168:641–674, 2017.
  * [53] S. Mendelson.  An unrestricted learning procedure.  Journal of the ACM, 66(6):1–42, 2019.
  * [54] J. Mourtada and S. Gaïffas.  An improper estimator with optimal excess risk in misspecified density estimation and logistic regression.  Journal of Machine Learning Research, 23(31):1–49, 2022.
  * [55] J. Mourtada, T. Vaškevičius, and N. Zhivotovskiy.  Distribution-free robust linear regression.  Mathematical Statistics and Learning, 4(3):253–292, 2021.
  * [56] A. Nemirovski.  Topics in non-parametric statistics.  Lectures on Probability Theory and Statistics: Ecole d’Ete de Probabilites de Saint-Flour XXVIII-1998, 28:85–277, 2000.
  * [57] M. Raginsky, A. Rakhlin, and A. Xu.  Information-theoretic stability and generalization.  Information-Theoretic Methods in Data Science, 10:302–329, 2021\.
  * [58] A. Rakhlin, K. Sridharan, and A. B. Tsybakov.  Empirical entropy, minimax regret and minimax risk.  Bernoulli, 23(2):789–824, 2017.
  * [59] P. Rigollet.  Kullback-Leibler aggregation and misspecified generalized linear models.  The Annals of Statistics, 40(2):639–665, 2012.
  * [60] O. Shamir.  The sample complexity of learning linear predictors with the squared loss.  Journal of Machine Learning Research, 16(108):3475–3486, 2015.
  * [61] T. Steinke and L. Zakynthinou.  Reasoning about generalization via conditional mutual information.  In Conference on Learning Theory, pages 3437–3452. PMLR, 2020.
  * [62] I. Steinwart and A. Christmann.  Support Vector Machines.  Information Science and Statistics. Springer-Verlag New York, 2008.
  * [63] A. B. Tsybakov.  Optimal rates of aggregation.  In Proceedings of the 16th Conference on Computational Learning Theory, pages 303–313, 2003.
  * [64] D. Van der Hoeven, N. Zhivotovskiy, and N. Cesa-Bianchi.  A regret-variance trade-off in online learning.  Advances in Neural Information Processing Systems, 35, 2022.
  * [65] V. Vapnik.  The Nature of Statistical Learning Theory.  Springer Science & Business Media, 1999.
  * [66] V. N. Vapnik and A. Y. Chervonenkis.  Theory of Pattern Recognition: Statistical Learning Problems (in Russian).  Nauka, Moscow, 1974.
  * [67] S. Vijaykumar.  Localization, convexity, and star aggregation.  Advances in Neural Information Processing Systems, 34:4570–4581, 2021.
  * [68] V. Vovk.  A game of prediction with expert advice.  Journal of Computer and System Sciences, 56(2):153–173, 1998.
  * [69] V. Vovk.  Competitive on-line statistics.  International Statistical Review, 69(2):213–248, 2001.
  * [70] M. J. Wainwright.  High-dimensional statistics: A non-asymptotic viewpoint.  Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, 2019.
  * [71] O. Wintenberger.  Optimal learning with Bernstein online aggregation.  Machine Learning, 106(1):119–141, 2017.
  * [72] Y. Yang.  Mixing strategies for density estimation.  The Annals of Statistics, 28(1):75–87, 2000.
  * [73] T. Zhang.  From $\varepsilon$-entropy to KL-entropy: Analysis of minimum information complexity density estimation.  The Annals of Statistics, 34(5):2180–2210, 2006.
  * [74] T. Zhang.  Information-theoretic upper and lower bounds for statistical estimation.  IEEE Transactions on Information Theory, 52(4):1307–1321, 2006\.
  * [75] N. Zhivotovskiy and S. Hanneke.  Localization of VC classes: Beyond local Rademacher complexities.  Theoretical Computer Science, 742:27–49, 2018.

##  Appendix A Proofs

This appendix contains the proofs of all the results stated in the main text.

###  A.1 Basic identities and inequalities

In this section, we collect some basic facts related to the $\mathrm{KL}$
divergence, quadratic loss and the curvature of the empirical the
$Q$-aggregation objective defined in (14). These facts are repeatedly used
throughout the proofs of our main theorems stated in Section 2. We will
repeatedly use the convex-conjugate duality of the moment generating function
and the $\mathrm{KL}$ divergence written in the following form: for any
function $h:\Theta\to\mathbf{R}$ such that the quantities below are well-
defined,

$\sup\limits_{\rho\in\mathcal{P}(\Theta)}\left\\{\langle\rho,h\rangle-\mathrm{KL}(\rho,\pi)\right\\}=\log\int_{\Theta}\exp\left(h(\theta)\right)\pi(\mathrm{d}\theta).$
(36)

In order to obtain localized PAC-Bayes bounds, at some point in our proofs we
will need to replace the difference of localized $\mathrm{KL}$ divergences
$\mathrm{KL}(\rho,\pi_{-f})-\mathrm{KL}(\gamma,\pi_{-f})$ by the difference
$\mathrm{KL}(\rho,\pi)-\mathrm{KL}(\gamma,\pi)$. The following lemma shows
that the price of this interchange is equal to $\langle\rho-\gamma,f\rangle$,
for any function $f$ used to localize the prior $\pi$.

######  Lemma 1.

For any probability distributions $\rho,\gamma,\pi\in\mathcal{P}(\Theta)$ such
that $\rho\ll\pi$ and $\gamma\ll\pi$, and for any function
$f:\Theta\to\mathbf{R}$ we have

$\mathrm{KL}(\rho,\pi_{-f})-\mathrm{KL}(\gamma,\pi_{-f})=\mathrm{KL}(\rho,\pi)-\mathrm{KL}(\gamma,\pi)+\langle\rho-\gamma,f\rangle.$
(37)

###### Proof.

A direct computation yields

$\displaystyle\mathrm{KL}(\gamma,\pi_{-f})$
$\displaystyle=\int_{\Theta}\log\bigg{(}\frac{\mathrm{d}\gamma}{\mathrm{d}\pi_{-f}}(\theta)\bigg{)}\gamma(\mathrm{d}\theta)$
(38)
$\displaystyle=\int_{\Theta}\log\left(\frac{\mathrm{d}\gamma}{\mathrm{d}\pi}(\theta)\cdot\frac{\langle\pi,\exp(-f)\rangle}{\exp(-f(\theta))}\right)\gamma(\mathrm{d}\theta)$
(39)
$\displaystyle=\mathrm{KL}(\gamma,\pi)+\langle\gamma,f\rangle+\log\langle\pi,\exp(-f)\rangle.$
(40)

Similarly, it holds that

$\mathrm{KL}(\rho,\pi_{-f})=\mathrm{KL}(\rho,\pi)+\langle\rho,f\rangle+\log\langle\pi,\exp(-f)\rangle.$
(41)

Subtracting (41) from (40) concludes the proof of this lemma. ∎

Next, we recall some basic properties of the squared loss. Define

$V_{n}(\rho)=\langle\rho,R_{n}\rangle-
R_{n}(\rho)=\langle\rho,\|f-f_{\rho}\|_{n}^{2}\rangle=\int_{\Theta}\|f_{\theta}-f_{\rho}\|_{n}^{2}\,\rho(\mathrm{d}\theta),$
(42)

where we remind that $R_{n}$ is defined in (12). Then, for any
$\gamma\in\mathcal{P}(\Theta)$ we have

$V_{n}(\rho)+\|f_{\rho}-f_{\gamma}\|_{n}^{2}=\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle.$
(43)

The two identities stated above are repeatedly applied in the rest of the
proofs without explicitly referring to this section.

The following lemma highlights the main property of the $Q$-aggregation
estimator that distinguishes its analysis from the exponential weights. This
lemma bears similarity with the analysis using offset Rademacher complexities
[45, 75, 67, 38], where the curvature of the loss function allows to obtain a
certain “complexity regularizing” term. However, the term obtained in the
lemma stated below is different from the one present in the offset Rademacher
complexity analysis, which would correspond to a term proportional to
$-\|f_{\widehat{\rho}}-f_{\gamma}\|_{n}^{2}$.

######  Lemma 2.

Let $\pi\in\mathcal{P}(\Theta)$ and $\beta>0$ be arbitrary. Let
$\widehat{\rho}=\widehat{\rho}^{\mathsf{q}}_{n}(\beta,\pi)$ be the minimizer
of the $Q$-aggregation objective defined in (14). Then, for any
$\gamma\in\mathcal{P}(\Theta)$, the following deterministic inequality holds:

$R_{n}(\widehat{\rho})-R_{n}(\gamma)-\frac{1}{\beta}\mathrm{KL}(\gamma,\pi)+\frac{1}{\beta}\mathrm{KL}(\widehat{\rho},\pi)\leqslant-\frac{1}{2}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{1}{2}V_{n}(\gamma).$
(44)

###### Proof.

Let
$Q_{n}(\rho)=\frac{1}{2}R_{n}(\rho)+\frac{1}{2}\langle\rho,R_{n}\rangle+\frac{1}{\beta}\mathrm{KL}(\rho,\pi)$.
Thus, $\widehat{\rho}$ is the unique minimizer of $Q_{n}$. Observe that

$\displaystyle
R_{n}(\widehat{\rho})-R_{n}(\gamma)=\left(Q_{n}(\widehat{\rho})-Q_{n}(\gamma)\right)-\frac{1}{2}V_{n}(\widehat{\rho})+\frac{1}{2}V_{n}(\gamma)-\frac{1}{\beta}\mathrm{KL}(\widehat{\rho},\pi)+\frac{1}{\beta}\mathrm{KL}(\gamma,\pi).$
(45)

Strong convexity of the squared loss combined with the optimality of
$\widehat{\rho}$ for the objective $Q_{n}$ yields:

$\displaystyle
Q_{n}(\widehat{\rho})-Q_{n}(\gamma)\leqslant-\frac{1}{2}\|f_{\widehat{\rho}}-f_{\gamma}\|_{n}^{2}.$
(46)

Plugging the above into (45) yields

$\displaystyle R_{n}(\widehat{\rho})-R_{n}(\gamma)$ (47)
$\displaystyle\leqslant-\frac{1}{2}\left(\|f_{\widehat{\rho}}-f_{\gamma}\|_{n}^{2}+V_{n}(\widehat{\rho})\right)+\frac{1}{2}V_{n}(\gamma)-\frac{1}{\beta}\mathrm{KL}(\widehat{\rho},\pi)+\frac{1}{\beta}\mathrm{KL}(\gamma,\pi).$
(48)
$\displaystyle=-\frac{1}{2}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{1}{2}V_{n}(\gamma)-\frac{1}{\beta}\mathrm{KL}(\widehat{\rho},\pi)+\frac{1}{\beta}\mathrm{KL}(\gamma,\pi),$
(49)

where the final identity follows from the variance identity (43). ∎

The distribution $\widehat{\rho}^{\mathsf{ew}}_{n}$ defined in (13) admits a
closed form analytic expression. This considerably simplifies the analysis of
the exponential weights estimator. In contrast, there is no simple closed form
expression for the distribution $\widehat{\rho}^{\mathsf{q}}$ defined in (14).
Thus, in general, the analysis of the $Q$-aggregation estimator poses more
challenges than the exponential weights estimator. A priori,
$\widehat{\rho}^{\mathsf{q}}_{n}$ is an arbitrary probability distribution on
$\Theta$, and therefore it takes values in the infinite-dimensional space
$\mathcal{P}(\Theta)$ whenever $\Theta$ is infinite-dimensional. Drawing
similarities with the _representer theorem_ and the _kernel trick_ commonly
employed in the context of reproducing kernel Hilbert spaces [4, 62], we
obtain the following result proved in Section A.2. It states that
$\widehat{\rho}^{\mathsf{q}}_{n}$ lies in a data-dependent finite-dimensional
exponential family of distributions. This characterization is particularly
useful when the prior distribution $\pi$ is Gaussian as we shall later explain
in Section 3.2.

######  Proposition 5.

There exist (data-dependent) constants $(a_{i})_{1\leqslant i\leqslant
n}\in\mathbf{R}^{n}$, $(b_{i})_{1\leqslant i\leqslant n}\in\mathbf{R}^{n}$ and
$c_{0}\in\mathbf{R}$ such that, for $\pi$-almost all $\theta\in\Theta$,

$\frac{\mathrm{d}\widehat{\rho}^{\mathsf{q}}_{n}}{\mathrm{d}\pi}(\theta)=\exp\bigg{\\{}\sum_{i=1}^{n}a_{i}f_{\theta}(X_{i})^{2}+\sum_{i=1}^{n}b_{i}f_{\theta}(X_{i})+c_{0}\bigg{\\}}.$

###  A.2 Proof of Proposition 5

By definition of the $Q$-aggregation estimator, we have
$\widehat{\rho}_{n}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\rho}\\{\beta
Q_{n}(\rho)+\mathrm{KL}({\rho},{\pi})\\}$, with

$\displaystyle Q_{n}(\rho)$
$\displaystyle=(1-\alpha)\langle\rho,R_{n}\rangle+\alpha R_{n}(\rho)$
$\displaystyle=\frac{1-\alpha}{n}\sum_{i=1}^{n}\int_{\Theta}(f_{\theta}(X_{i})-Y_{i})^{2}\rho(\mathrm{d}\theta)+\frac{\alpha}{n}\sum_{i=1}^{n}\bigg{(}\int_{\Theta}f_{\theta}(X_{i})\rho(\mathrm{d}\theta)-Y_{i}\bigg{)}^{2}$
$\displaystyle=\Phi\bigg{\\{}\bigg{(}\int_{\Theta}f_{\theta}(X_{i})^{2}\rho(\mathrm{d}\theta)\bigg{)}_{1\leqslant
i\leqslant
n},\bigg{(}\int_{\Theta}f_{\theta}(X_{i})\rho(\mathrm{d}\theta)\bigg{)}_{1\leqslant
i\leqslant n}\bigg{\\}}\,,$

where $\Phi:\mathbf{R}^{n}\times\mathbf{R}^{n}\to\mathbf{R}$ is defined by:

$\Phi\big{\\{}(u_{i})_{1\leqslant i\leqslant n},(v_{i})_{1\leqslant i\leqslant
n}\big{\\}}=\frac{1-\alpha}{n}\sum_{i=1}^{n}\big{(}u_{i}-2v_{i}Y_{i}+Y_{i}^{2}\big{)}+\frac{\alpha}{n}\sum_{i=1}^{n}(v_{i}-Y_{i})^{2}.$

Now, let
$u_{i}=\int_{\Theta}f_{\theta}(X_{i})^{2}\widehat{\rho}_{n}(\mathrm{d}\theta)$
and $v_{i}=\int_{\Theta}f_{\theta}(X_{i})\widehat{\rho}_{n}(\mathrm{d}\theta)$
for $i=1,\dots,n$. Define the set:

$\Delta=\bigg{\\{}\rho\in\mathcal{P}(\Theta):\int_{\Theta}f_{\theta}(X_{i})^{2}\rho(\mathrm{d}\theta)=u_{i},\
\int_{\Theta}f_{\theta}(X_{i})\rho(\mathrm{d}\theta)=v_{i},\ 1\leqslant
i\leqslant n\bigg{\\}}.$

Since $\widehat{\rho}_{n}\in\Delta$, we have
$\widehat{\rho}_{n}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\rho\in\Delta}\\{\beta
Q_{n}(\rho)+\mathrm{KL}({\rho},{\pi})\\}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\rho\in\Delta}\mathrm{KL}({\rho},{\pi})$
since $Q(\rho)=Q(\widehat{\rho}_{n})$ for $\rho\in\Delta$ by the above. The
conclusion then follows from the maximum entropy theorem; see, e.g., [20,
Theorem 12.1.1, p. 410]. ∎

###  A.3 Proof of Theorem 1

Under the Gaussian noise assumption, we have

$Y\sim\mathcal{N}(f^{\star},\sigma^{2}I_{n}),$

where with a slight abuse of notation, $f^{\star}$ denotes the vector
$(f^{\star}(X_{1}),\dots,f^{\star}(X_{n}))\in\mathbf{R}^{n}$.

We now recall some basic facts about Stein’s Unbiased Risk Estimator (SURE).
Let $\widehat{f}:\mathbf{R}^{n}\to\mathbf{R}^{n}$ be any estimator of the mean
vector $f^{\star}$ as a function of the observations vector $Y$. Then, SURE
for the estimator $\widehat{f}$ is given by the identity

$\mathrm{SURE}(\widehat{f},Y)=\|\widehat{f}(Y)-Y\|_{2}^{2}+2\sigma^{2}\sum_{i=1}^{n}\frac{\partial}{\partial
Y_{i}}\widehat{f}(Y)_{i}-n\sigma^{2},$ (50)

where $\widehat{f}(Y)_{i}$ denotes the $i$-th component of the $n$-dimensional
vector $\widehat{f}(Y)$. It suffices to assume the existence, continuity, and
integrability of the above partial derivatives so that a simple Gaussian
integration by parts argument yields

${\mathbf{E}}_{Y}\left[\|\widehat{f}-f^{\star}\|_{2}^{2}\right]={\mathbf{E}}_{Y}\left[\mathrm{SURE}(\widehat{f},Y)\right].$
(51)

For the exponential weights estimator
$\widehat{f}^{\mathsf{ew}}_{n}=\langle\widehat{\rho}^{\mathsf{ew}}_{n}(\beta,\pi),f\rangle$,
SURE admits a closed-form expression. To simplify the notation, in what
follows let $\widehat{\rho}=\widehat{\rho}^{\mathsf{ew}}_{n}(\beta,\pi)$ and
$\widehat{f}=\widehat{f}^{\mathsf{ew}}_{n}$. Then, we have

$\widehat{\rho}(\mathrm{d}\theta)=\frac{\exp(-\beta
R_{n}(f_{\theta}))\pi(d\theta)}{\int_{\Theta}\exp(-\beta
R_{n}(f_{\theta}))\pi(d\theta)}.$ (52)

A direct computation of the partial derivatives yields

$\frac{\partial}{\partial
Y_{i}}(\widehat{f}(Y))_{i}=\frac{2\beta}{n}\int_{\Theta}\left((f_{\theta})_{i}-(\widehat{f})_{i}\right)^{2}\widehat{\rho}(\mathrm{d}\theta)$
(53)

and hence,

$\frac{1}{n}\mathrm{SURE}(\widehat{f},Y)=\langle\widehat{\rho},R_{n}\rangle+\bigg{(}\frac{4\beta\sigma^{2}}{n}-1\bigg{)}V_{n}(\widehat{\rho})-\sigma^{2}.$
(54)

In previous analyses of the exponential weights estimator [44, 22], the second
term in the right hand side of the above identity is dropped by taking
$\beta\leqslant\frac{n}{4\sigma^{2}}$ . This is the part where our analysis
departs from the proofs that yield global bounds: taking a smaller $\beta$ and
keeping this term—that now becomes strictly negative—is crucial in order to
prove the localized bound stated in Theorem 1. We now proceed to finish the
proof of this theorem.

The constraints put on $\beta$ in the theorem statement guarantee that
$\frac{4\beta\sigma^{2}}{n}-1\leqslant-\frac{1}{2}$. Hence, we may bound the
SURE estimator (54) by

$\frac{1}{n}\mathrm{SURE}(\widehat{f},Y)\leqslant\langle\widehat{\rho},R_{n}\rangle-\frac{1}{2}V_{n}(\widehat{\rho})-\sigma^{2}.$
(55)

Adding $V_{n}(\widehat{\rho})+\sigma^{2}$ to both sides of the above
inequality and recalling that $\widehat{\rho}$ minimizes the
$\mathrm{KL}$-divergence penalized empirical risk minimization objective (13),
for any $\gamma\in\mathcal{P}(\Theta)$ it holds that

$\displaystyle\frac{1}{n}\mathrm{SURE}(\widehat{f},Y)+\frac{1}{2}V_{n}(\widehat{\rho})+\sigma^{2}$
(56)
$\displaystyle\leqslant\langle\widehat{\rho},R_{n}\rangle+\frac{1}{\beta}\mathrm{KL}(\rho,\pi)-\frac{1}{\beta}\mathrm{KL}(\rho,\pi)$
(57)
$\displaystyle\leqslant\langle\gamma,R_{n}\rangle+\frac{1}{\beta}\mathrm{KL}(\gamma,\pi)-\frac{1}{\beta}\mathrm{KL}(\rho,\pi)$
(58)
$\displaystyle=\langle\gamma,R_{n}\rangle+\frac{1}{\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{2}R})-\frac{1}{\beta}\mathrm{KL}(\rho,\pi_{-\frac{\beta}{2}R})+\frac{1}{2}\langle\widehat{\rho}-\gamma,R\rangle$
(59)
$\displaystyle\leqslant\langle\gamma,R_{n}\rangle+\frac{1}{\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{2}R})+\frac{1}{2}\langle\widehat{\rho}-\gamma,R\rangle,$
(60)

where the equality step follows from Lemma 1. Now observe that

$\langle\widehat{\rho},R\rangle=\langle\widehat{\rho},\|f-f^{\star}\|_{n}^{2}\rangle+\sigma^{2}=\|\widehat{f}-f^{\star}\|_{n}^{2}+V_{n}(\widehat{\rho})+\sigma^{2}.$
(61)

Plugging the above identity into (60) and rearranging yields

$\displaystyle\frac{1}{n}\mathrm{SURE}(\widehat{f},Y)-\frac{1}{2}\|\widehat{f}-f^{\star}\|_{n}^{2}\leqslant\langle\gamma,R_{n}-\frac{1}{2}R\rangle+\frac{1}{\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{2}R})-\frac{1}{2}\sigma^{2}.$
(62)

Taking expectations on both sides and using the fact that $\gamma$ does not
depend on the sample $Y=(Y_{1},\dots,Y_{n})^{\mathsf{T}}$, we have

$\displaystyle{\mathbf{E}}_{Y}\left[\frac{1}{n}\mathrm{SURE}(\widehat{f},Y)\right]-\frac{1}{2}{\mathbf{E}}_{Y}\left[\|\widehat{f}-f^{\star}\|_{n}^{2}\right]$
$\displaystyle\leqslant\frac{1}{2}\langle\gamma,R\rangle+\frac{1}{\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{2}R})-\frac{1}{2}\sigma^{2}$
(63)
$\displaystyle=\frac{1}{2}{\mathbf{E}}_{\theta\sim\gamma}\left[\|f_{\theta}-f^{\star}\|_{n}^{2}\right]+\frac{1}{\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{2}R}).$
(64)

Finally, because SURE is unbiased (51), it follows that

${\mathbf{E}}_{Y}\left[\|\widehat{f}-f^{\star}\|_{n}^{2}\right]\leqslant{\mathbf{E}}_{\theta\sim\gamma}\left[\|f_{\theta}-f^{\star}\|_{n}^{2}\right]+\frac{2}{\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{2}R}).$

Recalling that the choice of $\gamma\in\mathcal{P}(\Theta)$ is arbitrary and
using the inequality (10), we finish the proof. ∎

###  A.4 Proof of Theorem 2

To simplify the notation, let $\widehat{f}=\widehat{f}^{\mathsf{q}}_{n}$ and
$\widehat{\rho}=\widehat{\rho}^{\mathsf{q}}_{n}(\beta,\pi)$. Moreover, we
identify all the functions $f:\mathcal{X}\to\mathbf{R}$ as $n$-dimensional
vectors $(f(X_{1}),\dots,f(X_{n}))^{\mathsf{T}}\in\mathbf{R}^{n}$. Let
$\gamma$ be the minimizer of
$\gamma^{\prime}\mapsto\langle\gamma^{\prime},\|f-f^{\star}\|_{n}^{2}\rangle+\frac{3}{2\beta}\mathrm{KL}(\gamma^{\prime},\pi_{-\frac{\beta}{3}R}).$
Then, we have

$\displaystyle\|f_{\widehat{\rho}}-f^{\star}\|_{n}^{2}-\inf_{\gamma^{\prime}\in\mathcal{P}(\Theta)}\left\\{\langle\gamma^{\prime},\|f-f^{\star}\|_{n}^{2}\rangle+\frac{3}{2\beta}\mathrm{KL}(\gamma^{\prime},\pi_{-\frac{\beta}{3}R})\right\\}$
(65)
$\displaystyle=\|f_{\widehat{\rho}}-f^{\star}\|_{n}^{2}-\langle\gamma,\|f-f^{\star}\|_{n}^{2}\rangle-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{3}R})$
(66)
$\displaystyle=\|f_{\widehat{\rho}}-f^{\star}\|_{n}^{2}-\|f_{\gamma}-f^{\star}\|_{n}^{2}-V_{n}(\gamma)-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{3}R})$
(67)
$\displaystyle=R_{n}(\widehat{\rho})-R_{n}(\widehat{\gamma})+\frac{2}{n}\langle
f_{\widehat{\rho}}-f_{\gamma},Z\rangle-
V_{n}(\gamma)-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{3}R})$
(68)
$\displaystyle=R_{n}(\widehat{\rho})+\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})-R_{n}(\widehat{\gamma})-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{3}R})$
(69) $\displaystyle\quad\quad+\frac{2}{n}\langle
f_{\widehat{\rho}}-f_{\gamma},Z\rangle-
V_{n}(\gamma)-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R}).$
(70)

Next, applying Lemma 1 we have

$\displaystyle\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{3}R})$
(71)
$\displaystyle=\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi)-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi)+\frac{3}{2\beta}\langle\widehat{\rho}-\gamma,\frac{\beta}{3}R\rangle$
(72)
$\displaystyle=\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi)-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi)+\frac{1}{2}R(\widehat{\rho})+\frac{1}{2}V_{n}(\widehat{\rho})-\frac{1}{2}R(\gamma)-\frac{1}{2}V_{n}(\gamma)$
(73)
$\displaystyle=\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi)-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi)+\frac{1}{2}R_{n}(\widehat{\rho})+\frac{1}{2}V_{n}(\widehat{\rho})-\frac{1}{2}R_{n}(\gamma)-\frac{1}{2}V_{n}(\gamma)+\frac{1}{n}\langle
f_{\widehat{\rho}}-f_{\gamma},Z\rangle.$ (74)

Plugging the above identity into (70) results in the following equality

$\displaystyle\|f_{\widehat{\rho}}-f^{\star}\|_{n}^{2}-\inf_{\gamma^{\prime}\in\mathcal{P}(\Theta)}\left\\{\langle\gamma^{\prime},\|f-f^{\star}\|_{n}^{2}\rangle+\frac{3}{2\beta}\mathrm{KL}(\gamma^{\prime},\pi_{-\frac{\beta}{3}R})\right\\}$
(75)
$\displaystyle=\frac{3}{2}\left(R_{n}(\widehat{\rho})+\frac{1}{\beta}\mathrm{KL}(\widehat{\rho},\pi)-R_{n}(\widehat{\gamma})-\frac{1}{\beta}\mathrm{KL}(\gamma,\pi)\right)$
(76) $\displaystyle\quad\quad+\frac{3}{n}\langle
f_{\widehat{\rho}}-f_{\gamma},Z\rangle-\frac{3}{2}V_{n}(\gamma)+\frac{1}{2}V_{n}(\widehat{\rho})-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R}).$
(77)

Bounding the term (76) via the deterministic inequality proved in Lemma 2, the
right-hand side of the above identity is upper bounded by

$\displaystyle-\frac{3}{4}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{3}{4}V_{n}(\gamma)+\frac{3}{n}\langle
f_{\widehat{\rho}}-f_{\gamma},Z\rangle-\frac{3}{2}V_{n}(\gamma)+\frac{1}{2}V_{n}(\widehat{\rho})-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})$
(78)
$\displaystyle\leqslant-\frac{3}{4}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{3}{n}\langle
f_{\widehat{\rho}}-f_{\gamma},Z\rangle+\frac{1}{2}V_{n}(\widehat{\rho})-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})$
(79)
$\displaystyle\leqslant-\frac{1}{4}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{3}{n}\langle
f_{\widehat{\rho}}-f_{\gamma},Z\rangle-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R}),$
(80)

where in the final inequality we used the fact that
$\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle\geqslant
V_{n}(\widehat{\rho})$. Putting everything together results in the following
deterministic inequality:

$\displaystyle\frac{2\beta}{3}\|f_{\widehat{\rho}}-f^{\star}\|_{n}^{2}-\frac{2\beta}{3}\inf_{\gamma^{\prime}\in\mathcal{P}(\Theta)}\left\\{\langle\gamma^{\prime},\|f-f^{\star}\|_{n}^{2}\rangle+\frac{3}{2\beta}\mathrm{KL}(\gamma^{\prime},\pi_{-\frac{\beta}{3}R})\right\\}$
(81)
$\displaystyle\leqslant-\frac{\beta}{6}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{2\beta}{n}\langle
f_{\widehat{\rho}}-f_{\gamma},Z\rangle-\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})$
(82)
$\displaystyle\leqslant\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{-\frac{\beta}{6}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{2\beta}{n}\langle
f_{\rho}-f_{\gamma},Z\rangle-\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}$
(83)
$\displaystyle=\log\left({\mathbf{E}}_{\theta\sim\pi_{-\frac{\beta}{3}R}}\left[\exp\left(-\frac{\beta}{6}\|f_{\theta}-f_{\gamma}\|_{n}^{2}+\frac{2\beta}{n}\langle
f_{\theta}-f_{\gamma},Z\rangle\right)\right]\right)=U,$ (84)

where the last line follows from (36). We will now show that the right-hand
side of the above inequality is small with high probability by bounding its
moment-generating function. Indeed, we have

$\displaystyle{\mathbf{E}}_{Z}\left[\exp(U)\right]$
$\displaystyle={\mathbf{E}}_{Z}{\mathbf{E}}_{\theta\sim\pi_{-\frac{\beta}{3}R}}\left[\exp\left(-\frac{\beta}{6}\|f_{\theta}-f_{\gamma}\|_{n}^{2}+\frac{2\beta}{n}\langle
f_{\theta}-f_{\gamma},Z\rangle\right)\right]$ (85)
$\displaystyle={\mathbf{E}}_{\theta\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{Z}\left[\exp\left(-\frac{\beta}{6}\|f_{\theta}-f_{\gamma}\|_{n}^{2}+\frac{2\beta}{n}\langle
f_{\theta}-f_{\gamma},Z\rangle\right)\right]$ (86)
$\displaystyle={\mathbf{E}}_{\theta\sim\pi_{-\frac{\beta}{3}R}}\left[\exp\left(-\frac{\beta}{6}\|f_{\theta}-f_{\gamma}\|_{n}^{2}\right){\mathbf{E}}_{Z}\left[\exp\left(\frac{2\beta}{n}\langle
f_{\theta}-f_{\gamma},Z\rangle\right)\right]\right]$ (87)
$\displaystyle={\mathbf{E}}_{\theta\sim\pi_{-\frac{\beta}{3}R}}\left[\exp\left(-\frac{\beta}{6}\|f_{\theta}-f_{\gamma}\|_{n}^{2}\right)\prod_{i=1}^{n}{\mathbf{E}}_{Z_{i}}\left[\exp\left(\frac{2\beta}{n}(f_{\theta}(X_{i})-f_{\gamma}(X_{i}))Z_{i}\right)\right]\right],$
(88)

where the last line follows by the independence of the noise variables
$Z_{1},\dots,Z_{n}$. Further, by the sub-Gaussianity assumption on the noise
we have

$\displaystyle{\mathbf{E}}_{Z}\left[\exp(U)\right]$
$\displaystyle\leqslant{\mathbf{E}}_{\theta\sim\pi_{-\frac{\beta}{3}R}}\left[\exp\left(-\frac{\beta}{6}\|f_{\theta}-f_{\gamma}\|_{n}^{2}\right)\prod_{i=1}^{n}\exp\left(\frac{2\beta^{2}(f_{\theta}(X_{i})-f_{\gamma}(X_{i}))^{2}\sigma^{2}}{n^{2}}\right)\right]$
(89)
$\displaystyle={\mathbf{E}}_{\theta\sim\pi_{-\frac{\beta}{3}R}}\left[\exp\left(\left[\frac{2\beta^{2}\sigma^{2}}{n}-\frac{\beta}{6}\right]\|f_{\theta}-f_{\gamma}\|_{n}^{2}\right)\right].$
(90)

In particular, for any $\beta\in(0,\frac{n}{12\sigma^{2}}]$ it holds that
${\mathbf{E}}_{Z}[\exp(U)]\leqslant 1$. Hence, by Markov’s inequality and the
bound (84), for any $\delta\in(0,1)$, with probability at least $1-\delta$, it
holds that

$\frac{2\beta}{3}\|f_{\widehat{\rho}}-f^{\star}\|_{n}^{2}-\frac{2\beta}{3}\inf_{\gamma^{\prime}\in\mathcal{P}(\Theta)}\left\\{\langle\gamma^{\prime},\|f-f^{\star}\|_{n}^{2}\rangle+\frac{3}{2\beta}\mathrm{KL}(\gamma^{\prime},\pi_{-\frac{\beta}{3}R})\right\\}\leqslant
U\leqslant\log(1/\delta).$ (91)

By re-scaling the above inequality by $3/(2\beta)$ and using the inequality
(10), we finish the proof. ∎

###  A.5 Proof of Theorem 3

Let $P_{X}$ be the marginal distribution of the covariates $X$ and for any
$f\in L_{2}(P_{X})$ introduce the notation

$\|f\|_{2}^{2}=\|f\|_{L_{2}(P_{X})}^{2}={\mathbf{E}}_{X\sim P_{X}}[f(X)^{2}].$

Then, for any $\gamma\in\mathcal{P}(\Theta)$ let
$V(\gamma)=\langle\gamma,\|f-f_{\gamma}\|_{2}^{2}\rangle$; note that this is a
population counterpart to the empirical variance term $V_{n}$ defined in (42).

Because our proof involves classical empirical process arguments based on
symmetrization and contraction, it is convenient to introduce notation
commonly used in the empirical processes theory. Recall that $P$ is the data-
generating probability measure and let $P_{n}$ be its empirical counterpart;
that is, $P_{n}=n^{-1}\sum_{i=1}^{n}\delta_{(X_{i},Y_{i})}$, where
$\delta_{(X_{i},Y_{i})}$ denotes the Dirac measure concentrated on the $i$-th
data point. Then, for any function $h:\mathcal{X}\times[-b,b]\to\mathbf{R}$,
we set

$Ph={\mathbf{E}}_{(X,Y)\sim P}[h(X,Y)]\quad\text{and}\quad
P_{n}h=\frac{1}{n}\sum_{i=1}^{n}h(X_{i},Y_{i}).$

In particular, for any probability measure $\gamma\in\mathcal{P}(\Theta)$ we
introduce the function $\ell_{\gamma}(X,Y)=(f_{\gamma}(X)-Y)^{2}$. Then, we
have

$P\ell_{\gamma}=R(\gamma)\quad\text{and}\quad
P_{n}\ell_{\gamma}=R_{n}(\gamma).$

Finally, to simplify the notation let
$\widehat{\rho}=\widehat{\rho}^{\mathsf{q}}_{n}(\beta,\pi)$. We are now ready
to proceed with the proof of the theorem. Let $\gamma$ be the minimizer of
$R(\cdot)+\frac{3}{4}V(\cdot)+\frac{3}{2\beta}\mathrm{KL}(\cdot,\pi_{-\frac{\beta}{3}R})$
over $\mathcal{P}(\Theta)$. Then, we have

$\displaystyle
R(\widehat{\rho})-\inf_{\gamma^{\prime}\in\mathcal{P}(\Theta)}\left\\{\langle\gamma^{\prime},R\rangle+\frac{3}{2\beta}\mathrm{KL}(\gamma^{\prime},\pi_{-\frac{\beta}{3}R})\right\\}$
(92) $\displaystyle\leqslant
R(\widehat{\rho})-\inf_{\gamma^{\prime}\in\mathcal{P}(\Theta)}\left\\{R(\gamma^{\prime})+\frac{3}{4}V(\gamma^{\prime})+\frac{3}{2\beta}\mathrm{KL}(\gamma^{\prime},\pi_{-\frac{\beta}{3}R})\right\\}$
(93)
$\displaystyle=R(\widehat{\rho})-R(\gamma)-\frac{3}{4}V(\gamma)-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{3}R})$
(94)
$\displaystyle=P(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{3}{4}V(\gamma)+\frac{3}{2\beta}\Big{\\{}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})-\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{3}R})\Big{\\}}-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})$
(95)
$\displaystyle=P(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{3}{4}V(\gamma)+\frac{3}{2\beta}\Big{\\{}\mathrm{KL}(\widehat{\rho},\pi)-\mathrm{KL}(\gamma,\pi)\Big{\\}}+\frac{1}{2}\langle\widehat{\rho}-\gamma,R\rangle-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{2}R})$
(96)
$\displaystyle=\frac{3}{2}P(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{5}{4}V(\gamma)+\frac{1}{2}V(\widehat{\rho})-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi)+\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi)-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{2}R}),$
(97)

where the penultimate equality follows from Lemma 1. Adding and subtracting
$\frac{3}{2}P_{n}(\ell_{\widehat{\rho}}-\ell_{\gamma})$ we obtain the
following deterministic identity:

$\displaystyle
R(\widehat{\rho})-\inf_{\gamma\in\mathcal{P}(\Theta)}\left\\{R(\gamma)+\frac{3}{4}V(\gamma)+\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{3}R})\right\\}$
(98)
$\displaystyle=\frac{3}{2}(P-P_{n})(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})-\frac{5}{4}V(\gamma)+\frac{1}{2}V(\widehat{\rho})$
(99)
$\displaystyle\quad\quad+\frac{3}{2}\left(P_{n}(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{1}{\beta}\mathrm{KL}(\gamma,\pi)+\frac{1}{\beta}\mathrm{KL}(\widehat{\rho},\pi)\right).$
(100)

Plugging in the result of Lemma 2 into (100), and multiplying both sides by
$\frac{2}{3}$, we obtain

$\displaystyle\frac{2}{3}\left[R(\widehat{\rho})-\inf_{\gamma\in\mathcal{P}(\Theta)}\left\\{R(\gamma)+\frac{3}{4}V(\gamma)+\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{3}R})\right\\}\right]$
(101)
$\displaystyle\leqslant(P-P_{n})(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{1}{\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})-\frac{5}{6}V(\gamma)+\frac{1}{3}V(\widehat{\rho})-\frac{1}{2}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{1}{2}V_{n}(\gamma)$
(102)
$\displaystyle=(P-P_{n})(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{1}{18}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{1}{18}\langle\widehat{\rho},\|f-f_{\gamma}\|_{2}^{2}\rangle-\frac{1}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})$
(103)
$\displaystyle\quad\quad+\frac{1}{3}V(\widehat{\rho})+\frac{1}{18}\langle\widehat{\rho},\|f-f_{\gamma}\|_{2}^{2}\rangle-\frac{8}{18}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{1}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})$
(104) $\displaystyle\quad\quad+\frac{1}{2}V_{n}(\gamma)-\frac{5}{6}V(\gamma).$
(105)

It remains to obtain high-probability upper bounds on each of the three terms
(103), (104) and (105) and use the union bound. We show how to bound each term
separately in the three lemmas stated below and proved at the end of this
section.

To bound the term (103), we combine Talagrand’s contraction lemma with the
usual PAC-Bayes analysis based on the convex-conjugate duality of the moment
generating function and the $\mathrm{KL}$ divergence. This part of the
analysis closely follows that of Lecué and Rigollet [43].

######  Lemma 3.

Fix any $\gamma\in\mathcal{P}(\Theta)$. Then, for any
$\beta\in(0,n/(576b^{2})]$ and for any $\delta\in(0,1)$, with probability at
least $1-\delta/3$, it holds that

$\displaystyle\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{(P-P_{n})(\ell_{\rho}-\ell_{\gamma})-\frac{1}{18}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{1}{18}\langle\rho,\|f-f_{\gamma}\|_{2}^{2}\rangle-\frac{1}{2\beta}\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}$
(106) $\displaystyle\lesssim\frac{\log(3/\delta)}{\beta}.$ (107)

We bound the term (104) using the relation (36) combined with Bernstein’s
concentration inequality.

######  Lemma 4.

Fix any $\gamma\in\mathcal{P}(\Theta)$. Then, for any
$\beta\in(0,n/(50b^{2})]$ and for any $\delta\in(0,1)$, with probability at
least $1-\delta/3$, it holds that

$\displaystyle\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{\frac{1}{3}V(\rho)+\frac{1}{18}\langle\rho,\|f-f_{\gamma}\|_{2}^{2}\rangle-\frac{8}{18}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{1}{2\beta}\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}\lesssim\frac{\log(3/\delta)}{\beta}.$
(108)

Finally, since $\gamma\in\mathcal{P}(\Theta)$ is a fixed vector that does not
depend on the data, we bound (105) via a single application of Bernstein’s
inequality.

######  Lemma 5.

Fix any $\gamma\in\mathcal{P}(\Theta)$. Then, for any $\delta\in(0,1)$, with
probability at least $1-\delta/3$ it holds that

$\displaystyle\frac{1}{2}V_{n}(\gamma)-\frac{5}{6}V(\gamma)\lesssim\frac{b^{2}\log(3/\delta)}{n}.$
(109)

The proof of Theorem 3 follows by bounding the terms (103), (104) and (105)
via Lemma 3, 4 and 5, respectively, and applying the union bound. ∎

####  A.5.1 Proof of Lemma 3

Define

$\displaystyle U$
$\displaystyle=2\beta\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{(P-P_{n})(\ell_{\rho}-\ell_{\gamma})-\frac{1}{18}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{1}{18}\langle\rho,\|f-f_{\gamma}\|_{2}^{2}\rangle-\frac{1}{2\beta}\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}$
(110)
$\displaystyle=\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{2\beta(P-P_{n})(\ell_{\rho}-\ell_{\gamma})-\frac{\beta}{9}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{\beta}{9}\langle\rho,\|f-f_{\gamma}\|_{2}^{2}\rangle-\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}.$
(111)

Let $\sigma_{1},\dots,\sigma_{n}$ denote i.i.d. Rademacher signs. And let
$X_{i}^{\prime},Y_{i}^{\prime}$ be an independent copy of $X_{i},Y_{i}$. Let

$P_{n}^{\sigma}h=\frac{1}{n}\sum_{i=1}^{n}\sigma_{i}h(X_{i}).$

Throughout the remained of the proof we use $X$ and $Y$ as shorthand notations
for the sequences $X_{1},\dots,X_{n}$ and $Y_{1},\dots,Y_{n}$, respectively.
Then, symmetrization and contraction gives us

$\displaystyle{\mathbf{E}}_{(X,Y)}\left[\exp\left(U\right)\right]$ (112)
$\displaystyle\leqslant{\mathbf{E}}_{(X,Y)}{\mathbf{E}}_{\sigma}\left[\exp\left(2\cdot\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{2\beta
P_{n}^{\sigma}(\ell_{\rho}-\ell_{\gamma})-\frac{\beta}{9}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{1}{2}\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}\right)\right]$
(113)
$\displaystyle\leqslant{\mathbf{E}}_{(X,Y)}{\mathbf{E}}_{\sigma}\left[\exp\left(\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{4\beta
P_{n}^{\sigma}(\ell_{\rho}-\ell_{\gamma})-\frac{2\beta}{9}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle-\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}\right)\right]$
(114)
$\displaystyle\leqslant{\mathbf{E}}_{(X,Y)}{\mathbf{E}}_{\sigma}\left[\exp\left(\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{16b\beta
P_{n}^{\sigma}(f_{\rho}-f_{\gamma})-\frac{2\beta}{9}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle-\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}\right)\right],$
(115)

where the last inequality follow via contraction principle as in [43, Eq.
(3.11)]. In particular, the quadratic loss is Lipschitz over bounded domain.
That is, for any $x,y$:

$\displaystyle\left|\left(\ell_{\rho}(x,y)-\ell_{\gamma}(x,y)\right)-\left(\ell_{\rho^{\prime}}(x,y)-\ell_{\gamma}(x,y)\right)\right|$
(116)
$\displaystyle=\left|(f_{\rho}(x)-f_{\rho^{\prime}}(x))(f_{\rho}(x)+f_{\rho^{\prime}}(x)-2y)\right|$
(117) $\displaystyle\leqslant 4b\left|f_{\rho}(x)-f_{\rho^{\prime}}(x)\right|$
(118)
$\displaystyle=4b\left|(f_{\rho}(x)-f_{\gamma}(x))-(f_{\rho^{\prime}}(x)-f_{\gamma}(x))\right|.$
(119)

We are ready to apply the PAC-Bayesian machinery to (115). Using (36), we have

$\displaystyle{\mathbf{E}}_{(X,Y)}\left[\exp\left(U\right)\right]$ (120)
$\displaystyle\leqslant{\mathbf{E}}_{(X,Y)}{\mathbf{E}}_{\sigma}{\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}\left[\exp\left(\frac{16b\beta}{n}\sum_{i=1}^{n}\sigma_{i}(f_{J}(X_{i})-f_{\gamma}(X_{i}))-\frac{2\beta}{9}\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\right]$
(121)
$\displaystyle={\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{(X,Y)}\left[\exp\left(-\frac{2\beta}{9}\|f_{J}-f_{\gamma}\|_{n}^{2}\right){\mathbf{E}}_{\sigma}\left[\exp\left(\frac{16b\beta}{n}\sum_{i=1}^{n}\sigma_{i}(f_{J}(X_{i})-f_{\gamma}(X_{i}))\right)\right]\right].$
(122)

Now we will upper bound the expectation with respect to the Rademacher signs
via Hoeffding’s lemma. We have

$\displaystyle{\mathbf{E}}_{\sigma}\left[\exp\left(\frac{16b\beta}{n}\sum_{i=1}^{n}\sigma_{i}(f_{J}(X_{i})-f_{\gamma}(X_{i}))\right)\right]$
$\displaystyle=\prod_{i=1}^{n}{\mathbf{E}}_{\sigma_{i}}\left[\exp\left(\frac{16b\beta}{n}\sigma_{i}(f_{J}(X_{i})-f_{\gamma}(X_{i}))\right)\right]$
(123)
$\displaystyle\leqslant\prod_{i=1}^{n}\exp\left(\frac{128b^{2}\beta^{2}}{n^{2}}(f_{J}(X_{i})-f_{\gamma}(X_{i})^{2})\right)$
(124)
$\displaystyle=\exp\left(\frac{128b^{2}\beta^{2}}{n}\|f_{J}-f_{\gamma}\|_{n}^{2}\right).$
(125)

Plugging the above into (122) yields

$\displaystyle{\mathbf{E}}_{(X,Y)}\left[\exp\left(U\right)\right]\leqslant{\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{(X,Y)}\left[\exp\left(-\frac{2\beta}{9}\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\exp\left(\frac{128b^{2}\beta^{2}}{n}\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\right].$
(126)

For any $\beta\in(0,n/(576b^{2})]$, the above is at most equal to one. It
follows that for any $\beta\in(0,n/(576b^{2})]$ we have,

${\mathbf{E}}\exp(U)\leqslant 1$

and therefore, via Markov’s inequality, with probability at least
$1-\delta/3$, it holds that

$U\leqslant\log(3/\delta),$

which concludes the proof of this lemma. ∎

####  A.5.2 Proof of Lemma 4

Define

$\displaystyle U$
$\displaystyle=2\beta\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{\frac{1}{3}V(\rho)+\frac{1}{18}\langle\rho,\|f-f_{\gamma}\|_{2}^{2}\rangle-\frac{8}{18}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{1}{2\beta}\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}$
(127)
$\displaystyle=\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{\frac{2\beta}{3}V(\rho)+\frac{\beta}{9}\langle\rho,\|f-f_{\gamma}\|_{2}^{2}\rangle-\frac{8\beta}{9}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle-\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{3}R})\right\\}.$
(128)

Next, observe that

$V(\rho)=\langle\rho,\|f-f_{\gamma}\|_{2}^{2}\rangle-\|f_{\rho}-f_{\gamma}\|_{2}^{2}\leqslant\langle\rho,\|f-f_{\gamma}\|_{2}^{2}\rangle.$

Applying the above inequality, it follows that

$U\leqslant\sup_{\rho\in\mathcal{P}(\Theta)}\left\\{\frac{7\beta}{9}\langle\rho,\|f-f_{\gamma}\|_{2}^{2}\rangle-\frac{8\beta}{9}\langle\rho,\|f-f_{\gamma}\|_{n}^{2}\rangle-\mathrm{KL}(\rho,\pi_{-\frac{\beta}{3}R})\right\\}.$

Let $X$ be a shorthand for $X_{1},\dots,X_{n}$ and let
$\sigma=(\sigma_{1},\dots,\sigma_{n})$ be a sequence of independent Rademacher
random signs. Then, computing the above supremum explicitly and applying the
symmetrization inequality (cf. the proof of Lemma 3) yields

$\displaystyle{\mathbf{E}}_{X}\left[\exp\left(U\right)\right]$ (129)
$\displaystyle={\mathbf{E}}_{X}{\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}\left[\exp\left(\frac{\beta}{18}\left(14\|f_{J}-f_{\gamma}\|_{2}^{2}-16\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\right)\right]$
(130)
$\displaystyle={\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{X}\left[\exp\left(\frac{\beta}{18}\left(14\|f_{J}-f_{\gamma}\|_{2}^{2}-16\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\right)\right]$
(131)
$\displaystyle={\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{X}\left[\exp\left(\frac{\beta}{18}\left(15\|f_{J}-f_{\gamma}\|_{2}^{2}-15\|f_{J}-f_{\gamma}\|_{n}^{2}-\|f_{J}-f_{\gamma}\|_{2}^{2}-\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\right)\right]$
(132)
$\displaystyle\leqslant{\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{X}{\mathbf{E}}_{\sigma}\left[\exp\left(\frac{\beta}{9}\left(\frac{15}{n}\sum_{i=1}^{n}\sigma_{i}(f_{J}(X_{i})-f_{\gamma}(X_{i}))^{2}-\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\right)\right]$
(133)
$\displaystyle={\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{X}\exp\left(-\frac{\beta}{9}\|f_{J}-f_{\gamma}\|_{n}^{2}\right){\mathbf{E}}_{\sigma}\left[\exp\left(\frac{15\beta}{9n}\sum_{i=1}^{n}\sigma_{i}(f_{J}(X_{i})-f_{\gamma}(X_{i}))^{2}\right)\right].$
(134)

We can now bound the expectation with respect to the Rademacher signs via
Hoeffding’s inequality. In particular, we have

$\displaystyle{\mathbf{E}}_{X}\left[\exp(U)\right]$ (135)
$\displaystyle\leqslant{\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{X}\exp\left(-\frac{\beta}{9}\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\prod_{i=1}^{n}{\mathbf{E}}_{\sigma_{i}}\left[\exp\left(\frac{15\beta}{9n}\sigma_{i}(f_{J}(X_{i})-f_{\gamma}(X_{i}))^{2}\right)\right]$
(136)
$\displaystyle\leqslant{\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{X}\exp\left(-\frac{\beta}{9}\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\prod_{i=1}^{n}\exp\left(\frac{225\beta^{2}}{162n^{2}}(f_{J}(X_{i})-f_{\gamma}(X_{i}))^{4}\right)$
(137)
$\displaystyle\leqslant{\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{X}\exp\left(-\frac{\beta}{9}\|f_{J}-f_{\gamma}\|_{n}^{2}\right)\prod_{i=1}^{n}\exp\left(\frac{225\beta^{2}}{162n^{2}}\cdot
4b^{2}(f_{J}(X_{i})-f_{\gamma}(X_{i}))^{2}\right)$ (138)
$\displaystyle={\mathbf{E}}_{J\sim\pi_{-\frac{\beta}{3}R}}{\mathbf{E}}_{X}\exp\left(\left[\frac{900b^{2}\beta^{2}}{162n}-\frac{\beta}{9}\right]\|f_{J}-f_{\gamma}\|_{n}^{2}\right).$
(139)

For any $\beta\in(0,n/(50b^{2})]$ the above is bounded by one. Hence, by
Markov’s inequality, for any $\delta\in(0,1)$ and for any
$\beta\in(0,n/(50b^{2})]$ we have

$\mathbf{P}(U\geqslant\log(3/\delta))\leqslant\frac{\delta}{3},$

which is what we wanted to show. ∎

####  A.5.3 Proof of Lemma 5

Define the function $g$ by

$g(X)=\langle\gamma,(f(X)-f_{\gamma}(X))^{2}\rangle={\mathbf{E}}_{J\sim\gamma}[(f_{J}(X)-f_{\gamma}(X))^{2}].$

Then,

$V_{n}(\gamma)=\frac{1}{n}\sum_{i=1}^{n}g(X_{i})\quad\text{and}\quad
V(\gamma)={\mathbf{E}}[g(X)].$

Note that $|g(X)|\leqslant 4b^{2}$, and

$\displaystyle\operatorname{Var}(g(X))$
$\displaystyle\leqslant{\mathbf{E}}g(X)^{2}={\mathbf{E}}_{X}\left({\mathbf{E}}_{J\sim\gamma}[(f_{J}(X)-f_{\gamma}(X))^{2}]\right)^{2}$
(140)
$\displaystyle\leqslant{\mathbf{E}}_{X}{\mathbf{E}}_{J\sim\gamma}[(f_{J}(X)-f_{\gamma}(X))^{4}]\leqslant
4b^{2}{\mathbf{E}}_{X}{\mathbf{E}}_{J\sim\gamma}[(f_{J}(X)-f_{\gamma}(X))^{2}]$
(141) $\displaystyle=4b^{2}V(\gamma).$ (142)

Hence, by Bernstein’s inequality, with probability at least $1-\delta/3$ we
have

$\displaystyle V_{n}(\gamma)-V(\gamma)\leqslant\sqrt{\frac{2\cdot
4b^{2}V(\gamma)\cdot\log(3/\delta)}{n}}+\frac{4b^{2}\log(3/\delta)}{3n}.$
(143)

It follows that with probability at least $1-\delta/3$ we have

$\displaystyle\frac{1}{2}V_{n}(\gamma)-\frac{5}{6}V(\gamma)$
$\displaystyle=\frac{1}{2}\left(V_{n}(\gamma)-V(\gamma)\right)-\frac{1}{3}V(\gamma)$
(144) $\displaystyle\leqslant\frac{1}{2}\sqrt{\frac{2\cdot
4b^{2}V(\gamma)\cdot\log(3/\delta)}{n}}-\frac{1}{3}V(\gamma)+\frac{1}{2}\frac{4b^{2}\log(3/\delta)}{3n}$
(145)
$\displaystyle\leqslant\sup_{x\in\mathbf{R}}\bigg{\\{}\frac{1}{2}x\sqrt{\frac{2\cdot
4b^{2}\cdot\log(3/\delta)}{n}}-\frac{1}{3}x^{2}\bigg{\\}}+\frac{1}{2}\frac{4b^{2}\log(3/\delta)}{3n}$
(146) $\displaystyle\lesssim\frac{b^{2}\log(3/\delta)}{n},$ (147)

which is what we wanted to show. ∎

###  A.6 Proof of Theorem 5

The proof follows the proof of Theorem 3 almost verbatim. Let us discuss the
necessary changes. First, replacing $R(\cdot)$ with $R_{n}^{\prime}(\cdot)$,
$V(\cdot)$ with $V_{n}^{\prime}(\cdot)$, $P$ with $P_{n}^{\prime}$, and
$\|\cdot\|_{2}$ with $\|\cdot\|^{\prime}_{n}$ (the versions with primes
correspond to empirical means with respect to the second half of the sample)
in the first part of the proof of Theorem 3, we obtain for any
$\gamma\in\mathcal{P}(\Theta)$ that possibly depends on the entire sample of
size $2n$, We have

$\displaystyle
R^{\prime}_{n}(\widehat{\rho})-\left(R_{n}^{\prime}(\gamma)+\frac{3}{4}V_{n}^{\prime}(\gamma)+\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})\right)$
(148)
$\displaystyle=P_{n}^{\prime}(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{3}{4}V_{n}^{\prime}(\gamma)+\frac{3}{2\beta}\Big{\\{}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})-\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})\Big{\\}}$
(149)
$\displaystyle\qquad-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})$
(150)
$\displaystyle=P_{n}^{\prime}(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{3}{4}V_{n}^{\prime}(\gamma)+\frac{3}{2\beta}\Big{\\{}\mathrm{KL}(\widehat{\rho},\pi)-\mathrm{KL}(\gamma,\pi)\Big{\\}}$
(151)
$\displaystyle\qquad+\frac{1}{4}\langle\widehat{\rho}-\gamma,R_{n}^{\prime}\rangle+\frac{1}{4}\langle\widehat{\rho}-\gamma,R_{n}\rangle-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})$
(152)
$\displaystyle=\frac{5}{4}(P_{n}^{\prime}-P_{n})(\ell_{\widehat{\rho}}-\ell_{\gamma})-V_{n}^{\prime}(\gamma)+\frac{1}{4}V_{n}^{\prime}(\widehat{\rho})-\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi)+\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi)$
(153)
$\displaystyle\qquad+\frac{3}{2}P_{n}(\ell_{\widehat{\rho}}-\ell_{\gamma})+\frac{1}{4}V_{n}(\widehat{\rho})-\frac{1}{4}V_{n}(\gamma)-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})}).$
(154)

where the penultimate equality follows from Lemma 1. By Lemma 2 we have

$R_{n}(\widehat{\rho})-R_{n}(\gamma)-\frac{1}{\beta}\mathrm{KL}(\gamma,\pi)+\frac{1}{\beta}\mathrm{KL}(\widehat{\rho},\pi)\leqslant-\frac{1}{2}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{1}{2}V_{n}(\gamma).$

Plugging this inequality into (154), we have

$\displaystyle
R^{\prime}_{n}(\widehat{\rho})-\left(R_{n}^{\prime}(\gamma)+\frac{3}{4}V_{n}^{\prime}(\gamma)+\frac{3}{2\beta}\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})\right)$
$\displaystyle\leqslant\frac{5}{4}(P_{n}^{\prime}-P_{n})(\ell_{\widehat{\rho}}-\ell_{\gamma})-V_{n}^{\prime}(\gamma)+\frac{1}{4}V_{n}^{\prime}(\widehat{\rho})-\frac{3}{2\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})$
$\displaystyle\quad+\frac{1}{2}V_{n}(\gamma)-\frac{3}{4}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle+\frac{1}{4}V_{n}(\widehat{\rho}).$
$\displaystyle=\frac{5}{4}(P_{n}^{\prime}-P_{n})(\ell_{\widehat{\rho}}-\ell_{\gamma})-\frac{3}{16}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{3}{16}\langle\widehat{\rho},(\|f-f_{\gamma}\|_{n}^{\prime})^{2}\rangle-\frac{3}{4\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})$
$\displaystyle\quad+\frac{1}{4}V_{n}^{\prime}(\widehat{\rho})+\frac{1}{4}V_{n}(\widehat{\rho})+\frac{3}{16}\langle\widehat{\rho},(\|f-f_{\gamma}\|_{n}^{\prime})^{2}\rangle-\frac{9}{16}\langle\widehat{\rho},\|f-f_{\gamma}\|_{n}^{2}\rangle-\frac{3}{4\beta}\mathrm{KL}(\widehat{\rho},\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})$
$\displaystyle\quad+\frac{1}{2}V_{n}(\gamma)-V_{n}^{\prime}(\gamma).$

As in the proof of Theorem 3, we individually bound each term appearing in the
last three lines. To do so, we repeat the lines of the proof of Lemma 3, Lemma
4 and Lemma 5. There are two observations that allow us to work with the prior
distribution $\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})}$ that depends on
the entire sample of size $2n$. In fact, one can still use the symmetrization
since the random measure $\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})}$ does
not change if we swap any two observations between the the first and the
second half of the sample of size $2n$. Furthermore, after we apply the
symmetrization step, all the randomness in our analysis comes from the
Rademacher random signs. The remaining technical steps are exactly the same as
in the proofs of Lemma 3, and Lemma 4. It is only left to modify the argument
of Lemma 5. The only difference is that $\gamma$ depends on the entire sample
of size $2n$ and we cannot use the Bernstein inequality as before.
Nevertheless, we can still use the exchangeability of our sample. For any
$h>0$ we have for sequence of independence Rademacher random variables (using
the notation of the proof of Lemma 5)

$\displaystyle{\mathbf{E}}\exp\left(hn\left(\frac{1}{2}V_{n}(\gamma)-V_{n}^{\prime}(\gamma)\right)\right)$
$\displaystyle={\mathbf{E}}{\mathbf{E}}_{\sigma}\exp\left(\frac{3h}{4}\sum\limits_{i=1}^{n}\left(\sigma_{i}(g(X_{i})-g(X_{n+i}))-\frac{1}{3}g(X_{i})-\frac{1}{3}g(X_{n+i})\right)\right)$
$\displaystyle\leqslant{\mathbf{E}}{\mathbf{E}}_{\sigma}\exp\left(\frac{3h}{2}\sum\limits_{i=1}^{n}\left(\sigma_{i}g(X_{i})-\frac{1}{3}g(X_{i})\right)\right)$
$\displaystyle\leqslant{\mathbf{E}}\exp\left(\frac{36h^{2}}{8}\sum\limits_{i=1}^{n}b^{2}g(X_{i})-\frac{1}{2}\sum\limits_{i=1}^{n}g(X_{i})\right)=1,$

for $h=\frac{1}{3b}$. Here we used that $g(X)\geqslant 0$. By Markov’s
inequality, with probability at least $1-\delta/3$,

$\frac{1}{2}V_{n}(\gamma)-V_{n}^{\prime}(\gamma)\leqslant\frac{3b\log(3/\delta)}{n}.$

The claim follows as before. ∎

###  A.7 Proof of Theorem 4

Before presenting the proof, we make some remarks. It appears that the
existing analysis of aggregation procedures [41, 52, 53] goes through the
complexity of the class formed by pairwise differences between functions in
the reference class. That is, the existing upper bounds are written in term of
the localized complexity of the class
$\mathcal{F}-\mathcal{F}=\\{f-g:f,g\in\mathcal{F}\\}$. While for convex
classes $\mathcal{F}-\mathcal{F}$ has essentially the same local complexity as
$\mathcal{F}$ itself, the gap can be significant for non-convex classes as
shown by the following example.

######  Example 1 (The local complexity of the class $\mathcal{F}-\mathcal{F}$
can be larger than that of $\mathcal{F}$).

Consider the class of thresholds in $\mathbb{Z}$. That is, define
$\mathcal{F}=\\{x\mapsto\bm{1}(x\geqslant t):t\in\mathbb{Z}\\}$. It is
immediate to verify that the star number $s$ of $\mathcal{F}$ is equal to $2$.
We further observe that the class $\mathcal{F}-\mathcal{F}$ contains the class
of singletons $\mathcal{G}=\\{x\mapsto\bm{1}(x=t):t\in\mathbb{Z}\\}$. However,
the star number of $\mathcal{G}$ is equal to infinity. Informally, while each
function in $\mathcal{F}$ has at most two neighbors when restricted to a
finite sample, the class $\mathcal{F}-\mathcal{F}$ is more expressive and each
function has a rich neighborhood.

Thus, if one wants to incorporate the star number in the bound such as (33),
we need to avoid the existing techniques that focus on the class
$\mathcal{F}-\mathcal{F}$. Our techniques are particularly suitable for this
purpose.

###### Proof.

The key point is to upper bound
$\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})$ in Theorem
5. Fix $\beta=c_{1}n$. Observe that our data dependent projection class
$\mathcal{F}|_{S_{2n}}=\\{(f(X_{1}),\ldots,f(X_{2n})):f\in\mathcal{F}\\}$ is
finite and can be seen as a subset of $\\{0,1\\}^{2n}$. Moreover, its VC
dimension is bounded by $d$ and its star number is bounded by $s$. Observe
that $\beta(R_{n}(f)+R^{\prime}_{n}(f))$ is proportional to the number of
mistakes that $f$ does on the sample of size $2n$. Assume that $\gamma$ puts
its whole mass on any $f^{\prime}\in\mathcal{F}|_{S_{2n}}$ that minimizes
$R_{n}(f)+R_{n}^{\prime}(f)$. We pick one of these minimizers and denote it by
$f^{\prime}$. We have

$\displaystyle\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})$
(155)
$\displaystyle=\log\left(\sum\limits_{f\in\mathcal{F}|_{S_{2n}}}\exp(-c_{1}n(R_{n}(f)+R^{\prime}_{n}(f)-R_{n}(f^{\prime})-R^{\prime}_{n}(f^{\prime}))/6)\right)$
(156)
$\displaystyle=\log\left(\sum\limits_{k=0}^{2n}\big{|}\\{f\in\mathcal{F}|_{S_{2n}}:n(R_{n}(f)+R^{\prime}_{n}(f)-R_{n}(f^{\prime})-R^{\prime}_{n}(f^{\prime}))=k\\}\big{|}\cdot\exp(-c_{1}k/6)\right)$
(157)

We now want to bound

$|\\{f\in\mathcal{F}|_{S_{2n}}:n(R_{n}(f)+R^{\prime}_{n}(f)-R_{n}(f^{\prime})-R_{n}^{\prime}(f^{\prime}))=k\\}|=(I)_{k},$

for all $k=0,\ldots,2n$. Using the bound in [35, Theorem 10], one can verify
that

$\displaystyle\bigl{|}\\{x\in S_{2n}:\textrm{there exist}\
f,g\in\mathcal{F}|_{S_{2n}}\ \textrm{such that}\ f(x)\neq g(x),$
$\displaystyle\qquad\qquad\qquad\textrm{and}\
R_{n}(f)+R^{\prime}_{n}(f)=R_{n}(g)+R^{\prime}_{n}(g)=k/n+R_{n}(f^{\prime})+R^{\prime}_{n}(f^{\prime})\\}\bigr{|}$
$\displaystyle\qquad\leqslant(s+1)(n(R_{n}(f^{\prime})+R^{\prime}_{n}(f^{\prime}))+k).$

Thus, by the standard Sauer-Vapnik-Chervonenkis lemma we have

$(I)_{k}\leqslant\left(\frac{e(s+1)(n(R_{n}(f^{\prime})+R^{\prime}_{n}(f^{\prime}))+k)}{d}\right)^{d}.$

Plugging this into (157) we have for some $c_{2},c_{3},c_{4},c_{5}>0$,

$\displaystyle\mathrm{KL}(\gamma,\pi_{-\frac{\beta}{6}(R_{n}+R^{\prime}_{n})})$
$\displaystyle\leqslant\log\left(\sum\limits_{k=0}^{2n}\left(\frac{e(s+1)(n(R_{n}(f^{\prime})+R^{\prime}_{n}(f^{\prime}))+k)}{d}\right)^{d}\exp(-c_{1}k/6)\right)$
$\displaystyle\leqslant\log\Bigg{(}c_{2}\left(e(s+1)\right)^{d}\Bigg{(}\left(\frac{c_{3}n(R_{n}(f^{\prime})+R^{\prime}_{n}(f^{\prime})+d)}{d}\right)^{d}$
$\displaystyle\qquad+\sum\limits_{k=c_{4}(R_{n}(f^{\prime})+R^{\prime}_{n}(f^{\prime}))+d}^{\infty}\exp(c_{5}d\log(k/d)-c_{1}k/6)\Bigg{)}\Bigg{)}.$

Using $\log(k/d)\leqslant k/d$ for large enough $k$, it is easy to see that
the last summand does not play any role and is absorbed by absolute constants.
Finally, we use the standard relative Vapnik-Chervonenkis bound [66] to prove
that, with probability at least $1-\delta$,

$R_{n}(f^{\prime})+R^{\prime}_{n}(f^{\prime})\leqslant
4\inf\limits_{f\in\mathcal{F}}R(f)+\frac{c(d\log(en/d)+\log(1/\delta))}{n}.$

where $c>0$ is some absolute constant. Combining the bounds and using the
union bound, we obtain the bound

$\mathbf{P}\left(R^{\prime}_{n}(\widehat{\rho})-\min\limits_{f\in\mathcal{F}}R^{\prime}_{n}(f)\geqslant
c_{6}\left(\frac{d}{n}\log\left(\frac{es\left(nR(f^{\star})+d\log(\frac{n}{d})+\log(\frac{4}{\delta})\right)}{d}\right)+\frac{\log(\frac{4}{\delta})}{n}\right)\right)\leqslant\delta,$

where $c_{6}>0$ is some constant. Using [34, Lemma 20], we simplify the
expression under the logarithm and conclude the proof. ∎

###  A.8 Proofs from Section 3.2

We introduce additional notation. Let $\lVert{x}\rVert_{A}=\langle
Ax,x\rangle^{1/2}$ for $A$ a positive semi-definite matrix and
$x\in\mathbf{R}^{d}$. In the random-design setting, letting
$\theta^{\star}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\theta\in\mathbf{R}^{d}}R(f_{\theta})$,
we have
$R(f_{\theta})-R(f_{\theta^{\star}})={\mathbf{E}}[(f_{\theta}(X)-f_{\theta^{\star}}(X))^{2}]=\lVert{\theta-\theta^{\star}}\rVert_{\Sigma}^{2}$
for every $\theta\in\mathbf{R}^{d}$. Likewise, in the fixed-design setting,
$R(f_{\theta})-R(f_{\theta^{\star}})=\lVert{\theta-\theta^{\star}}\rVert_{\widehat{\Sigma}_{n}}^{2}$
Recall, for $\lambda>0$, the optimal regularized parameter is given by

$\theta_{\lambda}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\theta\in\mathbf{R}^{d}}\big{\\{}R(\theta)+\lambda\lVert{\theta}\rVert^{2}\big{\\}}=\mathop{\mathrm{arg}\,\mathrm{min}}_{\theta\in\mathbf{R}^{d}}\big{\\{}\lVert{\theta-\theta^{\star}}\rVert_{\Sigma}^{2}+\lambda\lVert{\theta}\rVert^{2}\big{\\}}=(\Sigma+\lambda
I_{d})^{-1}\Sigma\theta^{\star},$

such that

$\inf_{\theta\in\mathbf{R}^{d}}\big{\\{}R(\theta)+\lambda\lVert{\theta}\rVert^{2}\big{\\}}-R(\theta^{\star})=R(\theta_{\lambda})-R(\theta^{\star})+\lambda\lVert{\theta_{\lambda}}\rVert^{2}=\lambda\langle(\Sigma+\lambda
I_{d})^{-1}\Sigma\theta^{\star},\theta^{\star}\rangle,$

and

$R(\theta_{\lambda})-R(\theta^{\star})=\lambda^{2}\langle(\Sigma+\lambda
I_{d})^{-2}\Sigma\theta^{\star},\theta^{\star}\rangle.$

We use these expressions in what follows.

####  A.8.1 Proof of Proposition 3

Letting $\lambda=\gamma/\beta$, we have for any $\theta\in\mathbf{R}^{d}$,

$\displaystyle\beta
R(f_{\theta})+\gamma\lVert{\theta}\rVert^{2}=\beta\big{\\{}R(\theta)+\lambda\lVert{\theta}\rVert^{2}\big{\\}}=\beta\big{\\{}\langle(\Sigma+\lambda
I_{d})(\theta-\theta_{\lambda}),\theta-\theta_{\lambda}\rangle+R(\theta_{\lambda})+\lambda\lVert{\theta_{\lambda}}\rVert^{2}\big{\\}}$

where the last equality stems from the fact that the two functions of $\theta$
are quadratic, with the same Hessian, minimizer and minimum. As a result,
denoting by $\pi_{\beta}$ the distribution
$\mathcal{N}(\theta_{\lambda},\beta^{-1}(\Sigma+\lambda
I_{d})^{-1})=\mathcal{N}(\theta_{\lambda},(\beta\Sigma+\gamma I_{d})^{-1})$,
we have

$\displaystyle e^{-\beta R(\theta)/2}\pi(\mathrm{d}\theta)$
$\displaystyle=\frac{e^{-(\beta
R(f_{\theta})+\gamma\lVert{\theta}\rVert^{2})/2}\mathrm{d}\theta}{({2\pi}\gamma^{-1})^{d/2}}$
$\displaystyle=e^{-\beta\\{R(\theta_{\lambda})+\lambda\lVert{\theta_{\lambda}}\rVert^{2}\\}/2}\frac{(2\pi)^{d/2}\det[(\beta\Sigma+\gamma
I_{d})^{-1}]^{1/2}}{({2\pi}\gamma^{-1})^{d/2}}\pi_{\beta}(\mathrm{d}\theta)$
$\displaystyle=e^{-\beta\\{R(\theta_{\lambda})+\lambda\lVert{\theta_{\lambda}}\rVert^{2}\\}/2}\det\big{(}I_{d}+{\lambda^{-1}}\Sigma\big{)}^{-1/2}\pi_{\beta}(\mathrm{d}\theta).$

This implies that the global complexity writes

$\displaystyle-\frac{2}{\beta}\log\bigg{(}\int_{\mathbf{R}^{d}}e^{-\beta
R(\theta)/2}\pi(\mathrm{d}\theta)\bigg{)}$
$\displaystyle=-\frac{2}{\beta}\log\bigg{(}e^{-\beta\\{R(\theta_{\lambda})+\lambda\lVert{\theta_{\lambda}}\rVert^{2}\\}/2}\det\big{(}I_{d}+{\lambda^{-1}}\Sigma\big{)}^{-1/2}\bigg{)}$
$\displaystyle=R(\theta_{\lambda})+\lambda\lVert{\theta_{\lambda}}\rVert^{2}+\beta^{-1}\log\det(I_{d}+\lambda^{-1}\Sigma).$