Sayak Ray Chowdhury<EMAIL_ADDRESS>
Microsoft Research, India and Patrick Saux<EMAIL_ADDRESS>
Univ. Lille, Inria, CNRS, UMR 9189 - CRIStAL, F-59000, Lille, France and
Odalric-Ambrym Maillard<EMAIL_ADDRESS>
Univ. Lille, Inria, CNRS, UMR 9189 - CRIStAL, F-59000, Lille, France and
Aditya Gopalan<EMAIL_ADDRESS>
Indian Institute of Science, Bangalore, India

# Bregman Deviations of Generic Exponential Families

###### Abstract

We revisit the method of mixtures, or Laplace method, to study the
concentration phenomenon in generic (possibly multidimensional) exponential
families. Using duality properties of the Bregman divergence to construct
nonnegative martingales, we establish a generic inequality controlling the
deviation between the parameter of the family and a finite sample estimate.
This bound is expressed in the local geometry induced by the Bregman pseudo-
metric. Moreover, it is time-uniform and involves a quantity extending the
classical information gain to exponential families, which we call the Bregman
information gain. For the practitioner, we instantiate this novel bound to
several classical families, e.g., Gaussian (including with unknown variance or
multivariate), Bernoulli, Exponential, Weibull, Pareto, Poisson and Chi-
square, yielding explicit forms of the confidence sets and the Bregman
information gain. We further compare the resulting confidence bounds to state-
of-the-art time-uniform alternatives and show this novel method yields
competitive results. Finally, we apply our result to the design of generalized
likelihood ratio tests for change detection, capturing new settings such as
variance change in Gaussian families.

###### keywords:

Exponential families, Bregman divergence, concentration bounds.

##  1 Introduction

Concentration inequalities are a powerful set of methods in statistical theory
with key applications in machine learning. Often in machine learning
applications, a learner estimates some quantity solely based on samples from
an unknown distribution and would like to know the magnitude of the estimation
error. The typical example is that of the mean $\mu$ of some real-valued
random variable $X$, estimated by its empirical mean built from a sample of
$n$ independent and identically distributed (i.i.d.) observations. We refer
the interested reader to the monographs of [Boucheron et al. (2013), Raginsky
and Sason (2018), Zeitouni and Dembo (1998)] for standard results and related
topics.

In many situations, one may want to further estimate some vector parameter, as
in, e.g., linear bandits (Abbasi-Yadkori et al., 2011) and logistic bandits
(Faury et al., 2020). A closely related problem is to estimate the parameter
of a distribution coming from a parametric family (Chowdhury et al., 2021).
_Exponential families_ are a flexible way to formalize such distributions over
a set $\mathcal{X}$ by describing densities of the form $p_{\theta}(x)\propto
h(x)\exp(\langle\theta,F(x)\rangle)$, for some given feature function
$F:\mathcal{X}\mapsto\mathbb{R}^{d}$ and base function $h$ (Section 2). Most
classical distributions fall into this category, e.g., Gaussian (possibly
multivariate, with or without known variance), Exponential, Gamma, Chi-square,
Weibull, Pareto, Poisson, Bernoulli and Multinomial distributions (Amari,
2016).

A number of problems in current-day machine learning involve sequential,
active data-sampling strategies (Cesa-Bianchi and Lugosi, 2006). This includes
multi-armed bandits, reinforcement learning, active learning and federated
learning, to mention a few application domains. Since the decision to sample a
novel observation results from the interaction between the learning algorithm
and the environment, and depends on past observations, one needs to design
concentration inequalities working with a _random_ number of observations
typically at a random stopping time (Durrett, 2019). A natural way to handle
this difficulty is to derive time-uniform concentration inequalities,
producing sequences of confidence sets valid uniformly over all number of
observations with high probability, as opposed to being valid for a single
number of observation.

A popular method in bandit theory is to combine supermartingale techniques
with union bound arguments over a geometric time grid, a technique known as
time peeling (or stitching) – see Bubeck (2010); Cappé et al. (2013) for early
uses in bandits, as well as Garivier (2013), or more recently Maillard
(2019b); see also Howard et al. (2020, 2021) for a recent, complementary
survey of the history of this field, and Kuchibhotla and Zheng (2021) for an
extension of Bentkus’ concentration bounds (Bentkus, 2004) for bounded
distributions using time peeling. The method of mixtures, initiated by Robbins
and Pitman (1949); Robbins (1970) and popularized in Peña et al. (2008) is a
powerful alternative to peeling for developing time-uniform confidence sets.
It has been applied to sub-Gaussian families in Abbasi-Yadkori et al. (2011),
leading to a variety of applications (Chowdhury and Gopalan, 2017; Durand et
al., 2018; Kirschner and Krause, 2018). In Kaufmann and Koolen (2021), a
generalization to handle one-dimensional exponential families is considered,
with applications to Gaussian (known variance) and Gamma distributions (known
shape). A fairly different ‘capital process’ construction technique has been
recently developed for bounded distributions in Shafer and Vovk (2019), and
popularized further in Waudby-Smith and Ramdas (2023).

In this work, we revisit the method of mixtures for parametric exponential
families of arbitrary dimension, expressing deviations in the natural
(Bregman) divergence of the family. The setting of exponential families is
convenient for integration in a Bayesian setup, thanks to the notion of
conjugate prior that enables us to reduce computation of tedious integrals to
simple parameter updates. Here, we exploit this property to obtain explicit
mixtures of martingales. Exponential families are largely used in modern
machine learning, yet concentration tools available to the practitioner are
comparatively scarce beyond the Gaussian case. To help close this gap, we
obtain both sharp and computationally tractable confidence sets, especially in
the small sample regime.

#### Outline and contributions.

In Section 2, we first recall some background material on exponential families
and their associated Bregman divergences. Section 3.1 states our main result
(Theorem 3): a _time-uniform_ concentration inequality for exponential
families. Specifically, we control the Bregman deviations associated with the
log-partition function of the family using a novel information-theoretic
quantity, the _Bregman information gain_. On a high level, this quantifies the
information gain about the parameter of the family after observing i.i.d.
samples from it, which is measured in terms of the natural Bregman divergence
of the family. To illustrate the utility of this general result, we detail in
Section 3.2 how Bregman information gain and deviation inequalities specialize
for well-known exponential families, resulting in fully explicit confidence
sets (see Table 1). To the best of our knowledge, we are the first to derive
an explicit _time-uniform_ deviation inequality for two-parameter Gaussian
(i.e. both mean and variance are unknown), Chi-square, Weibull, Pareto, and
Poisson distributions. Our result is an adaptation of the method of mixtures
technique, and a proof sketch is outlined in Section 3.3. In Section 4, we
numerically evaluate the high-probability confidence sets built from our
method for classical families, and achieve state-of-the-art time-uniform
bounds. Finally, in Section 5, we generalize Theorem 3 to obtain a _doubly
time-uniform_ concentration inequality for generic exponential families, which
could be of independent interest (Theorem 5). We present an application of
both results in controlling the false alarm probability of the Generalized
Likelihood Ratio (GLR) test, used for change detection in the exponential
family model.

##  2 Exponential Families and Bregman Divergence

In this section, we introduce exponential families and the link between their
Kullback-Leibler (KL) and Bregman divergences, as well as useful properties of
these divergences.

#### Exponential families.

We consider an exponential family of distributions
$\\{p_{\theta}\\}_{\theta\in\Theta}$ over some set $\mathcal{X}$,
parameterized in some open set $\Theta\subset\mathbb{R}^{d}$, whose density or
mass function has the form
$p_{\theta}(x)\\!=\\!h(x)\exp(\langle\theta,F(x)\rangle\\!-\\!\mathcal{L}(\theta))$.
Here, $F\\!:\\!\mathcal{X}\\!\to\\!\mathbb{R}^{d}$ is the feature function,
$h\\!:\\!\mathcal{X}\\!\to\\!\mathbb{R}_{+}$ is the base function, and
$\mathcal{L}$ represents the normalization term (a.k.a. log-partition
function, convex w.r.t. $\theta$) given by
$\mathcal{L}(\theta)=\log\int_{\mathcal{X}}h(x)\exp(\langle\theta,F(x)\rangle)dx$.
We denote by
$\Theta_{\mathcal{D}}=\big{\\{}\theta\in\mathbb{R}^{d}:\mathcal{L}(\theta)<\infty\big{\\}}$
the domain of $\mathcal{L}$ and by
$\Theta_{I}=\big{\\{}\theta\in\Theta_{\mathcal{D}}:\det\nabla^{2}\mathcal{L}(\theta)>0\big{\\}}$
the set on which its Hessian is invertible. We assume that
$\Theta\subset\Theta_{I}$, which is tantamount to assuming that the family is
minimal, and ensures we only consider non-degenerate distributions
($\nabla\mathcal{L}$ is one-to-one on its domain). We use notations
$\mathbb{P}_{\theta},\mathbb{E}_{\theta}$ to explicitly refer to the
probability and expectation associated to the distribution $p_{\theta}$.

#### Bregman divergence of an exponential family.

A fundamental property of exponential families is the following form for the
KL divergence between two distributions with parameters
$\theta,\theta^{\prime}\in\Theta_{\mathcal{D}}$:

$\displaystyle\texttt{KL}(p_{\theta},p_{\theta^{\prime}})=\langle\theta-\theta^{\prime},\mathbb{E}_{\theta}(F(X))\rangle-\mathcal{L}(\theta)+\mathcal{L}(\theta^{\prime})\,.$

Here, $\mathbb{E}_{\theta}(F(X))$ is called the vector of expectation
parameters (a.k.a. dual parameters), and is equal to
$\nabla\mathcal{L}(\theta)$. Hence, it holds that
$\texttt{KL}(p_{\theta},p_{\theta^{\prime}})=\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)$,
where $\mathcal{B}_{\mathcal{L}}$ is known as the Bregman divergence (Bregman,
1967) with potential function $\mathcal{L}$, defined by

$\displaystyle\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\stackrel{{\scriptstyle\rm
def}}{{=}}\mathcal{L}(\theta^{\prime})-\mathcal{L}(\theta)-\langle\theta^{\prime}-\theta,\nabla\mathcal{L}(\theta)\rangle\,.$

#### Tail and duality properties.

The canonical Bregman divergence of an exponential family enjoys two
fundamental properties. The first one links it to the log-moment generating
function of the random variable $F(X)$, which makes Bregman divergences well
suited to control the tail behavior of random variables appearing in
concentration inequalities. The second one highlights duality properties that
enables convenient algebraic manipulations. To this end, for any
$\lambda\in\mathbb{R}^{d}$, we define the function
$\mathcal{B}_{{\mathcal{L}},\theta}(\lambda)\\!=\\!\mathcal{B}_{\mathcal{L}}(\theta\\!+\\!\lambda,\theta)\\!=\\!{\mathcal{L}}(\theta\\!+\\!\lambda)\\!-\\!{\mathcal{L}}(\theta)\\!-\\!\langle\lambda,\nabla{\mathcal{L}}(\theta)\rangle$.
Also, we introduce the Legendre-Fenchel dual operator $\star$ associating a
function $G$ to its dual
$G^{\star}(x)=\sup_{\lambda}\langle\lambda,x\rangle-G(\lambda)$.

######  Lemma 1 (Properties of Bregman divergences)

For all $\theta\\!\in\\!\Theta_{\mathcal{D}}$ and
$\lambda\\!\in\\!\mathbb{R}^{d}$ such that
$\theta\\!+\\!\lambda\in\Theta_{\mathcal{D}}$,

$\displaystyle\log\mathbb{E}_{\theta}\left[\exp\left(\langle\lambda,F(X)-\mathbb{E}_{\theta}[F(X)]\rangle\right)\right]=\mathcal{B}_{\mathcal{L},\theta}(\lambda)\,.$

Furthermore, if $\nabla\mathcal{L}$ is one-to-one, the following Bregman
duality relations hold for any
$\theta,\theta^{\prime}\in\Theta_{\mathcal{D}}:$

$\displaystyle\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\mathcal{B}_{\mathcal{L},\theta^{\prime}}^{\star}(\nabla{\mathcal{L}}(\theta)-\nabla{\mathcal{L}}(\theta^{\prime}))=\mathcal{B}_{\mathcal{L}^{\star}}(\nabla{\mathcal{L}}(\theta),\nabla{\mathcal{L}}(\theta^{\prime})).$

More generally, the following holds for any $\alpha\in[0,1]:$

$\displaystyle\mathcal{B}_{{\mathcal{L}},\theta^{\prime}}^{\star}\left(\alpha(\nabla{\mathcal{L}}(\theta)-\nabla{\mathcal{L}}(\theta^{\prime}))\right)=\mathcal{B}_{\mathcal{L}}\big{(}\theta^{\prime},\theta_{\alpha}\big{)},\,\,\text{where}\,\,\theta_{\alpha}=\nabla{\mathcal{L}}^{-1}\big{(}\alpha\nabla{\mathcal{L}}(\theta)+(1-\alpha)\nabla{L}(\theta^{\prime})\big{)}.$

The second half of this technical lemma is essentially a change of variable
formula to move back and forth between two representations of an exponential
family: in natural parameters (measured by the Bregman divergence between
$\theta^{\prime}$ and $\theta$) and in expectation parametrization (measured
by the dual Bregman divergence between
$\nabla\mathcal{L}(\theta)=\mathbb{E}_{\theta}[F(X)]$ and
$\nabla\mathcal{L}(\theta^{\prime})=\mathbb{E}_{\theta^{\prime}}[F(X)]$ ). For
more background on this, which forms the basis of the information geometry
field, we refer to Amari (2016, Sections 2.1 and 2.7). This result is at the
root of the martingale construction behind Theorem 3 in the next section. For
completeness, the proof of this classical lemma is given in Appendix B.

##  3 Time-uniform Bregman Concentration

In this section, we are interested in controlling the deviation between a
parameter $\theta\in\Theta$ and its estimate $\theta_{n}$ built from $n$
observations from distribution $p_{\theta}$. We naturally measure this
deviation in terms of the canonical Bregman divergence of the family. Further,
we would like to control this deviation not only for a single sample number
$n$, but simultaneously for all $n\in\mathbb{N}$. Namely, we would like to
upper bound quantities of the form $\mathbb{P}\left[\exists
n\\!\in\\!\mathbb{N}:\mathcal{B}_{\mathcal{L}}(\theta,\theta_{n})\geqslant\dots\right]$.
Such a control is very useful in contexts when observations are gathered
sequentially (either actively or otherwise), and especially when the number of
observations $n$ is unknown beforehand. Classical examples include multi-armed
bandits (Auer et al., 2002), or model-based reinforcement learning (Jaksch et
al., 2010).

###  3.1 A generic deviation inequality

In this section, we consider the problem of controlling the Bregman-deviations
of a parameter estimate of $\theta$. To this end, we adapt the method of
mixtures (a.k.a. Laplace method) from (Peña et al., 2008). The method is
originally designed in the context of Gaussian distributions, where it yields
simple closed-form expressions, even though it can be applied more generally.
We state below a generic extension of the method to parametric exponential
families and introduce a quantity that measures a form of _information gain_
about $\theta$ after observing $n$ samples, but expressed in terms of the
natural Bregman divergence. For this reason, we call this quantity the Bregman
information gain.

######  Definition 2 (Bregman information gain)

Let $X_{1},\ldots,X_{n}\\!\sim\\!p_{\theta}$ be i.i.d. samples generated from
$p_{\theta}$, where $\theta\in\Theta\subset\Theta_{I}$, and let
$\theta_{0}\\!\in\\!\Theta$ be a reference parameter. For any constant
$c\\!>\\!0$, the Bregman information gain about $\theta$ after observing
$X_{1},\ldots,X_{n}$ from $\theta_{0}$ is defined as

$\displaystyle\gamma_{n,c}(\theta_{0})\\!=\\!\log\\!\left(\\!\frac{\int_{\Theta}\exp\big{(}\\!-\\!c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{0})\big{)}d\theta^{\prime}}{\int_{\Theta}\\!\exp\\!\Big{(}\\!-\\!(n\\!+\\!c)\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta_{0}))\Big{)}d\theta^{\prime}}\\!\right)\\!,$

where
$\theta_{n,c}(\theta_{0})\\!=\\!(\nabla\mathcal{L})^{-1}\\!\bigg{(}\\!\frac{\sum_{t=1}^{n}F(X_{t})+c\nabla\mathcal{L}(\theta_{0})}{n+c}\\!\bigg{)}$
denotes a parameter estimate of $\theta$.111$\theta_{n,c}(\theta_{0}\\!)$ is
actually a maximum a posteriori estimate under a conjugate prior on $\theta$,
depending on the reference point $\theta_{0}$.

#### Dependence on $\theta_{0}$ and example.

The acute reader can note that the considered parameter estimate
$\theta_{n,c}(\theta_{0})$ and Bregman information gain
$\gamma_{n,c}(\theta_{0})$ involve a reference parameter $\theta_{0}$. It
makes sense to have such a local reference point since the Bregman divergence
is typically linked to metrics with _local_ (non-constant) curvature. Hence,
the (information) geometry seen from the perspective of different points
$\theta_{0}$ may be different, unlike in the Gaussian case. Specifically, for
a Gaussian $\mathcal{N}(\mu,\sigma^{2})$ with known variance $\sigma^{2}$, the
Bregman information gain w.r.t. a reference point $\mu_{0}$ reads

$\displaystyle\gamma_{n,c}^{\mathcal{N}}(\mu_{0})\\!=\\!\frac{1}{2}\\!\log\frac{2\pi\sigma^{2}}{c}\\!-\\!\frac{1}{2}\\!\log\frac{2\pi\sigma^{2}}{n\\!+\\!c}\\!=\\!\frac{1}{2}\\!\log\frac{n\\!+\\!c}{c}\,.$

It is independent of the reference parameter $\mu_{0}$, which is a consequence
of the fact that the Bregman divergence in this case is proportional to the
squared Euclidean distance; in other words, Gaussian distributions with known
variance exhibit invariant geometry. However, for other exponential families,
the geometries are inherently different, and hence Bregman information gain
depend explicitly on local reference points (see Section 3.2 for details). For
reference, the classical Gaussian information gain, i.e., the mutual
information between a prior $\mu\sim\mathcal{N}(\mu_{0},\sigma^{2})$ and the
average of an i.i.d. sample $X_{1},\dots,X_{n}$ drawn from
$\mathcal{N}(\mu,c\sigma^{2})$ is $\frac{1}{2}\log\frac{n\\!+\\!c}{c}$, which
matches the Bregman information gain.

We now present the main result of this paper – a time uniform confidence bound
for $\theta$ connecting the Bregman divergence geometry of the exponential
family with its Bregman information gain.

######  Theorem 3 (Main result: Laplace method for generic exponential
families)

Fix any $\delta\in(0,1]$ and $n\\!\in\\!\mathbb{N}$. Under the hypothesis of
Definition 2, consider the confidence set

$\displaystyle\Theta_{n,c}(\delta)=\left\\{\theta_{0}\in\Theta:(n+c)\mathcal{B}_{\mathcal{L}}\left(\theta_{0},\theta_{n,c}(\theta_{0})\right)\leqslant\log\frac{1}{\delta}+\gamma_{n,c}(\theta_{0})\right\\}.$

The following time-uniform control holds whenever the Bregman information gain
is well-defined:

$\displaystyle\mathbb{P}_{\theta}\left[\exists
n\in\mathbb{N}:\theta\notin\Theta_{n,c}(\delta)\right]\leqslant\delta\,.$

Note the implicit definition of the confidence set $\Theta_{n,c}(\delta)$.
where the parameter of interest $\theta_{0}$ appears in both arguments of the
Bregman divergence, as well as in the Bregman information gain
$\gamma_{n,c}(\theta_{0})$. Because of this, computing this confidence set
from the equation in Theorem 3 may seen non-trivial at first glance. However,
we show in Section 3.2 how these sets simplify for many classical families,
revealing how the computation can be made efficiently. Moreover, we observe
that these confidence sets are actually tighter than those of prior work, and
as such are especially well suited to be used when $n$ is small (for large
$n$, most methods produce essentially equivalent sets). Note that our
concentration bound holds uniformly over all $n$. Equivalently (Howard et al.,
2020, Lemma 3), it also holds that
$\mathbb{P}_{\theta}\left[\theta\notin\Theta_{\tau,c}(\delta)\right]\leqslant\delta$
for any random stopping time $\tau$.

#### Comparison with prior work.

Similar to this work, Kaufmann and Koolen (2021) extend the method of mixtures
technique to derive time-uniform concentration bounds for exponential
families. However, their proof technique is fairly different, relying instead
on discrete mixtures and stitching, which only works for single parameter
families, and involves case-specific calculations that are difficult to
generalize beyond Gaussian (with known variance) and Gamma (with known shape).
In contrast, our method applies to generic exponential families, including
distributions with more than one parameter such as Gaussian when both mean and
variance are unknown. Moreover, their discrete prior construction leads to
technical constants seemingly unrelated to the exponential family model. Our
prior is naturally induced by the exponential family, leveraging key
properties of Bregman divergences, yields more intrinsic quantities (Bregman
information gain) and perhaps a more elegant and shorter proof. Hence, our
results are not only more general, but also of fundamental interest.

#### Asymptotic behavior.

The asymptotic width of $\Theta_{n,c}(\delta)$ depends on the behavior of
$\gamma_{n,c}(\theta_{0})$ as $n\rightarrow+\infty$. Standard arguments show
that $\theta_{n,c}(\theta_{0})\rightarrow\theta$ and the Taylor expansion of
$\theta^{\prime}\mapsto\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)$
around $\theta^{\prime}=\theta$ is
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\frac{1}{2}(\theta^{\prime}-\theta)^{\top}\nabla^{2}\mathcal{L}(\theta)(\theta^{\prime}-\theta)+o\left(\lVert\theta^{\prime}-\theta\rVert^{2}\right)$
(note that $\nabla^{2}\mathcal{L}(\theta)$ is positive definite for
$\theta\in\Theta_{I}$). Laplace’s method for integrals (chapter 20 in
Lattimore and Szepesvári (2019), Shun and McCullagh (1995)) then gives the
following estimate:

$\displaystyle\int_{\Theta}\\!\exp\\!\Big{(}\\!-\\!(n\\!+\\!c)\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta_{0}))\Big{)}d\theta^{\prime}\approx\int_{\mathbb{R}^{d}}\\!\exp\\!\Big{(}\\!-\\!\frac{n\\!+\\!c}{2}(\theta-\theta^{\prime})^{\top}\nabla^{2}\mathcal{L}(\theta)(\theta-\theta^{\prime})\Big{)}d\theta^{\prime}\,,$

which is a simple Gaussian integral, and thus
$\gamma_{n,c}(\theta_{0})=\frac{d}{2}\log(1\\!+\\!n/c)+\mathcal{O}(1)$. This
asymptotic scaling is worse than the $\log\log n$ rate of the law of iterated
logarithm. However, this is a standard feature of the method of mixtures
compared to stitching (Maillard, 2019b; Howard et al., 2020), which is
compensated by its improved nonasymptotic sharpness, as evidenced in Section
4.

#### Dependence on the parameter $c$.

The parameter $c$ is usually chosen to be $1$ in the Gaussian case. However,
it is useful to study its influence over the bounds. We provide in Section 4 a
detailed study of this parameter, revealing first that the confidence bounds
are not significantly altered over a large range of values, and then
explaining how to pick an optimized value of $c$, e.g., for a specific horizon
$n=n_{0}$. Note that choosing a variable $c$ changing with $n$ is not allowed
by the theory, as this would break the martingale property used in our proof
(see Section 3.3). Moreover, it would be incompatible with the time-uniform
lower bound $\lvert\Theta_{n,c}\left(\delta\right)\rvert=\Omega(\sqrt{\log\log
n/n})$ provided by the law of iterated logarithm. Indeed, the Bregman
information gain with $c_{n}\\!\propto\\!n$ would be asymptotically
$\gamma_{n,c_{n}}(\theta_{0})=\frac{d}{2}\log(1\\!+\\!\mathcal{O}(1))+\mathcal{O}(1)=\mathcal{O}(1)$,
leading to
$\lvert\Theta_{n,c}\left(\delta\right)\rvert=\mathcal{O}(1/\sqrt{n})$ when
$n\rightarrow+\infty$.

#### Application to bandits.

One can apply the technique developed in proving Theorem 3 to build confidence
sets in standard $K$-armed bandit problems. In such problems, we are given $K$
distributions $\\{p_{\theta_{i}}\\}_{i=1}^{K}$ from a generic exponential
family (e.g., Gaussian) with parameters $\\{\theta_{i}\\}_{i=1}^{K}$, from
which we can draw samples, interpreted as rewards we want to maximize. This
setting is standard to analyze regret-optimal bandit algorithms (Cappé et al.,
2013; Korda et al., 2013; Baudry et al., 2020). Specifically, at each time
$n$, we choose an arm $i_{n}\in[K]$ based on past observations, and draw a
sample $X_{n}$ from its distribution $p_{\theta_{i_{n}}}$. The samples are
then used to update knowledge about the parameters
$\\{\theta_{i}\\}_{i=1}^{K}$ by building confidence sets for each of them. Let
$N_{i}(n)$ denote the number of times that we have chosen action $i$ up to
time $n$. Also, let $\theta_{i,N_{i}(n),c}(\theta_{0})$ and
$\gamma_{i,N_{i}(n),c}(\theta_{0})$ denote the parameter estimate and Bregman
information gain for arm $i$, respectively (similar to Definition 2 with $n$
replaced by $N_{i}(n)$ ). We construct the confidence set for arm $i$ at time
$n$ as

$\displaystyle\Theta_{i,n,c}(\delta)=\bigg{\\{}\theta_{0}\in\Theta:(N_{i}(n)+c)\mathcal{B}_{\mathcal{L}}\left(\theta_{0},\\!\theta_{i,N_{i}(n),c}(\theta_{0})\right)\leqslant\log\frac{1}{\delta}+\gamma_{i,N_{i}(n),c}(\theta_{0})\bigg{\\}}.$

Then, similar to Theorem 3, it holds that the true parameter $\theta_{i}$ lies
in the set $\Theta_{i,n,c}(\delta)$ for all time-steps $n\in\mathbb{N}$ with
probability at least $1-\delta$. Finally, we take a union bound over $i\in[K]$
to obtain confidence sets for all arms (with widths inflated by an additive
$\log K$ factor). Such construction is standard and is used in Abbasi-Yadkori
et al. (2011) for (sub)-Gaussian families. Possible applications include UCB
algorithms for regret minimization and designing GLR stopping rules for
tracking algorithms in pure exploration (Garivier and Kaufmann, 2016) in the
context of generic exponential families (see also Section 5 for another
application of GLR tests using the parameter estimate
$\theta_{n,c}(\theta_{0})$).

###### Remark 4

We provide in Appendix A.2 a complementary result (Corollary 8) using a
Legendre function $\mathcal{L}_{0}$ instead of the regularizing parameter $c$,
which eschews the use of a local reference $\theta_{0}$. However, we argue
that such a global regularization is actually less convenient to use except
for univariate and multivariate Gaussian distributions that anyway exhibit
invariant geometry. Furthermore, in Theorem 7, we prove a more general result
that handles the case of a sequence $(X_{t})_{t\in\mathbb{N}}$ of random
variables that are not independent, and having possibly different
distributions from each other, which we apply to build confidence sets in
linear bandits (see Appendix F).

###  3.2 Specification to classical families

In this section, we specify the result of Theorem 3 to some classical
exponential families. Interestingly, the literature on time-uniform
concentration bounds outside of random variables that are bounded, gamma with
fixed shape or Gaussian with known variance is significantly scarce, even
though many more distributions are commonly used in machine learning models.
We derive below explicit confidence sets for a range of distributions, which
we believe will be of interest for the wider machine learning and statistics
community. For instance, consider active learning in bandit problems
(Carpentier et al., 2011), where one targets upper confidence bounds on the
variance (rather than the mean); in the Gaussian case, this can be achieved
with Chi-square concentration. Hao et al. (2019) studies the classical UCB
algorithm for bandits under a weaker assumption that sub-Gaussianity,
involving the Weibull concentration. In differential privacy, concentration of
Laplace distribution (symmetrized exponential) is often used to study the
utility of differentially private mechanisms (Dwork et al., 2014). Finally,
heavy-tailed distributions such as Pareto have recently been of interest to
study risk-averse or corruption in bandit problems (Holland and Haress, 2021;
Basu et al., 2022).

We now make explicit the Bregman information gains and confidence sets for
some illustrative families (more examples and full derivations are provided in
Appendix C).

#### Gaussian (unknown mean and variance).

Let $X\sim\mathcal{N}(\mu,\sigma^{2})$. Given samples $X_{1},\ldots,X_{n}$, we
define $S_{n}=\sum_{t=1}^{n}X_{t}$ and $\widehat{\mu}_{n}=S_{n}/n$. Further,
for $\mu,\sigma\in\mathbb{R}\times\mathbb{R}^{+}$, we define the normalized
sum of squares
$Z_{n}(\mu,\sigma)=\frac{1}{\sigma^{2}}\sum_{t=1}^{n}(X_{t}-\mu)^{2}$. Then,
for reference parameters
$\mu_{0}\in\mathbb{R},\sigma_{0}\in\mathbb{R}^{+},c>0$, the Bregman
information gain reads

$\displaystyle\gamma_{n,c}^{\mathcal{N}}(\mu_{0},\sigma_{0})=\\!\frac{3}{2}\log\left(\frac{n}{n\\!+\\!c}Z_{n}(\widehat{\mu}_{n},\sigma_{0})+\frac{c}{n\\!+\\!c}Z_{n}(\mu_{0},\sigma_{0})+c\right)\\!+\\!f_{n,c}\,,$

where
$f_{n,c}\\!=\\!\left(\frac{n+c+1}{2}\right)\log\\!\left(n\\!+\\!c\right)\\!-\\!\left(\frac{c}{2}\\!+\\!2\right)\log
c\\!-\\!\frac{n}{2}(1\\!+\\!\log
2)\\!+\\!\log\Gamma\\!\left(\frac{c+3}{2}\right)\\!-\\!\log\Gamma\\!\left(\frac{n+c+3}{2}\right)$.

#### Bernoulli.

Let $X\sim\text{Bernoulli}(\mu)$, with unknown mean $\mu\in[0,1]$. Define, for
reference parameter $\mu_{0}\in\mathbb{R}$, $c>0$, the estimate
$\mu_{n,c}(\mu_{0})\\!=\\!\frac{S_{n}+c\mu_{0}}{n+c}$. Then, the Bregman
information gain is given by

$\displaystyle\gamma_{n,c}^{\text{Bernoulli}}(\mu_{0})=c\,\mathbb{H}(\mu_{0})\\!-\\!(n\\!+\\!c)\mathbb{H}(\mu_{n,c}(\mu_{0}))+\log\frac{\mathbf{B}(c\mu_{0},c(1-\mu_{0}))}{\mathbf{B}((n+c)\mu_{n,c}(\mu_{0}),(n+c)(1-\mu_{n,c}(\mu_{0}))},$

where
$\mathbf{B}(\alpha,\beta)\\!=\\!\int_{0}^{1}u^{\alpha-1}(1\\!-\\!u)^{\beta-1}du$
is the Beta function and $\mathbb{H}(\cdot)$ the Bernoulli entropy function.

#### Exponential.

Let $X\sim\text{Exp}(1/\mu)$, with unknown mean $\mu>0$. For reference
parameters $\mu_{0},c>0$, the Bregman information gain reads

$\displaystyle\gamma_{n,c}^{\text{Exp}}(\mu_{0})=\log\left(S_{n}/\mu_{0}+c\right)+\log\frac{\Gamma(c)}{\Gamma(n+c)}+(n+c-1)\log(n+c)-c\log
c-n~{}.$

#### Pareto.

Let $X\sim\text{Pareto}\left(\alpha\right)$, with unknown shape $\alpha>0$.
Define $L_{n}\\!=\\!\sum_{t=1}^{n}\log X_{t}$. Then, for reference parameters
$\alpha_{0},c>0$, the Bregman information gain is given by

$\displaystyle\gamma_{n,c}^{\text{Pareto}}(\alpha_{0})\\!=\\!\log\left(\alpha_{0}L_{n}\\!+\\!c\right)\\!+\\!\log\frac{\Gamma(c)}{\Gamma(n+c)}\\!+\\!(n\\!+\\!c\\!-\\!1)\log(n+c)\\!-\\!c\log
c\\!-\\!n.$

#### Chi-square.

Let $X\\!\sim\\!\chi^{2}(k)$, where $k\\!\in\\!\mathbb{N}$ is unknown. Define
$K_{n}\\!=\\!\sum_{t=1}^{n}\log\frac{X_{t}}{2}$. For reference points
$k_{0}\in\mathbb{N},c>0$, introduce an estimate $k_{n,c}(k_{0})$ satisfying
$\psi_{0}\\!\left(\\!\frac{k_{n,c}(k_{0})}{2}\\!\right)\\!=\\!\frac{K_{n}+c\;\psi_{0}\left(\frac{k_{0}}{2}\right)}{n+c}$.
In this case, the Bregman information gain is given by

$\displaystyle\gamma_{n,c}^{\chi^{2}}(k_{0})$
$\displaystyle=\frac{k_{n,c}(k_{0})}{2}\Big{(}K_{n}\\!+\\!c\psi_{0}\Big{(}\frac{k_{0}}{2}\Big{)}\Big{)}\\!-\\!(n\\!+\\!c)\log\Gamma\Big{(}\frac{k_{n,c}(k_{0})}{2}\Big{)}\\!+\\!c\log\Gamma\Big{(}\frac{k_{0}}{2}\Big{)}\\!-\\!c\;\frac{k_{0}}{2}\psi_{0}\Big{(}\frac{k_{0}}{2}\Big{)}$
$\displaystyle\\!+\\!\log\frac{J(c,c\psi_{0}(\frac{k_{0}}{2}))}{J(n\\!+\\!c,K_{n}\\!+\\!c\psi_{0}(\frac{k_{0}}{2}))}\,,\,\text{where}\,\,J(a,b)\\!:=\\!\sum\limits_{k^{\prime}=1}^{\infty}\exp\left(\\!-\\!a\log\Gamma\left(\\!\frac{k^{\prime}}{2}\\!\right)\\!+\\!b\frac{k^{\prime}}{2}\\!\right)\\!.$

The function $J(a,b)$ can be estimated using numerical methods (see Section
4). This paves a way to build high-probability confidence sets for Chi-square
distribution which has not been adequately captured in prior work. The sums
over $k^{\prime}$ above derive from the martingale construction of Theorem 3
with discrete mixture (i.e., w.r.t the counting measure). A continuous (i.e.,
w.r.t the Lebesgue measure) mixture similar to the other families is also
possible, as detailed in the appendix in Remark 9.

#### Explicit confidence sets.

We now turn to illustrate the confidence sets in the same exponential families
as above. They are obtained by specifying the generic form and simplifying the
resulting expression. We provide the confidence sets for two-parameter
Gaussian (i.e., unknown mean and variance), Bernoulli, Exponential, Pareto and
Chi-square distributions, respectively, in Table 1. The technical details of
the derivation of the specific forms for each illustrative family is postponed
to Appendix C, along with other distributions (Gamma, Poisson, Weibull) in
Table 2.

Table 1: Bregman confidence sets given by Theorem 3 for representative families Distribution | Parameters | Confidence Set
---|---|---
Gaussian |  $\mu\in\mathbb{R}$ $\sigma\in\mathbb{R}_{+}$ |  $\frac{1}{2}Z_{n}(\mu,\sigma)\\!-\\!\frac{n+c+3}{2}\log\left(\frac{n}{n\\!+\\!c}Z_{n}(\widehat{\mu}_{n},\sigma)\\!+\\!\frac{c}{n\\!+\\!c}Z_{n}(\mu,\sigma)\\!+\\!c\right)$ $\quad\leqslant\log\frac{1}{\delta}\\!-\\!\frac{n}{2}\log 2\\!-\\!\left(\frac{c}{2}\\!+\\!2\right)\log c\\!+\\!\frac{1}{2}\log\left(n\\!+\\!c\right)\\!+\\!\log\frac{\Gamma\left(\frac{c+3}{2}\right)}{\Gamma\left(\frac{n+c+3}{2}\right)}$
Bernoulli | $\mu\in[0,1]$ |  $S_{n}\log\\!\frac{1}{\mu}\\!+\\!(n\\!-\\!S_{n})\log\\!\frac{1}{1-\mu}\\!+\\!\log\\!\frac{\Gamma(S_{n}+c\mu)\Gamma(n-S_{n}+c(1-\mu))}{\Gamma(c\mu)\Gamma(c(1-\mu))}$ $\quad\leqslant\log\\!\frac{1}{\delta}\\!+\\!\log\\!\frac{\Gamma(n+c)}{\Gamma(c)}$
Exponential | $\mu\in\mathbb{R}_{+}$ |  $\frac{S_{n}}{\mu}\\!-\\!(n\\!+\\!c\\!+\\!1)\log\left(\frac{S_{n}}{\mu}+c\right)$ $\quad\leqslant\log\frac{1}{\delta}\\!+\\!\log\frac{\Gamma(c)}{\Gamma(n+c)}\\!-\\!\log(n\\!+\\!c)\\!-\\!c\log c$
Pareto | $\alpha\in\mathbb{R}$ |  $\alpha L_{n}-(n\\!+\\!c\\!+\\!1)\log\left(\alpha L_{n}+c\right)$ $\quad\leqslant\log\frac{1}{\delta}+\log\frac{\Gamma(c)}{\Gamma(n+c)}-\log(n+c)-c\log c$
Chi-square | $k\in\mathbb{N}$ |  $n\log\Gamma\left(\frac{k}{2}\right)\\!-\\!\frac{k}{2}K_{n}\\!-\\!\log J\left(c,c\psi_{0}\left(\frac{k}{2}\right)\right)$ $\quad+\log J\left(n\\!+\\!c,K_{n}\\!+\\!c\psi_{0}\left(\frac{k}{2}\right)\right)\leqslant\log\frac{1}{\delta}$

We now provide a set of illustrative numerical experiments to display the
confidence envelopes resulting from Theorem 3 in the case of classical
exponential families. We plot for a given value of the confidence level
$\delta$ and regularization parameter $c$ the convex sets
$n\mapsto\bigcap\nolimits_{n^{\prime}\leqslant
n}\Theta_{n^{\prime},c}(\delta)$ (taking the running intersection is standard
in sequential testing and was pioneered in Darling and Robbins (1967); it
ensures that the upper (resp. lower) confidence envelopes are nonincreasing
(resp. nondecreasing) as the sample size grows, thus providing tighter
bounds). In dimension $d=1$, we report the extremal points of these intervals,
which we call the upper and lower envelopes respectively. We refer to Figure 1
for Pareto, Chi-square and Gaussian (with unknown $\mu,\sigma$) and Appendix D
for many other families. Because we exploit the Bregman geometry, our bounds
capture a larger setting than typical mean estimation; for instance, we are
able to concentrate around the exponent $\alpha$ of a Pareto distribution even
when the distribution is not integrable ($\alpha\\!<\\!1$). Furthermore, for
two-parameter Gaussian, apart form being _anytime_ , our confidence sets are
convex and bounded in contrast to the one based on Chi-square quantiles with a
crude union bound (see Appendix D.1).

(a) Gaussian (mean and variance)

(b) Bernoulli

(c) Pareto

(d) Chi-square

Figure 1: Example of Gaussian time-uniform joint confidence sets for
$(\mu,\sigma)=(1,1)$ with sample sizes $n\in\left\\{10,25,50,100\right\\}$
observations (smaller confidence sets correspond to larger sample sizes), and
examples of time-uniform confidence envelopes for $\text{Bernoulli}(0.8)$,
$\text{Pareto}(0.5)$ and $\chi^{2}(5)$ on several realizations as a function
of the number of observations $n$. Thick lines indicate the median curve over
1000 replicates.

###  3.3 Proof Sketch: Theorem 3

We now sketch the proof of Theorem 3, and refer the interested reader to the
appendix for details.

#### Step 1: Martingale construction.

For any $\lambda\in\mathbb{R}^{d}$, we introduce the quantity

$\displaystyle M_{n}^{\lambda}$ $\displaystyle=$
$\displaystyle\exp\left(\langle\lambda,n(\mu_{n}\\!-\\!\mu)\rangle\\!-\\!n\mathcal{B}_{\mathcal{L},\theta}(\lambda)\right)\,,$

where $\mu_{n}\\!=\\!\frac{1}{n}\sum_{t=1}^{n}\\!F(X_{t})$ and
$\mu\\!=\\!\mathbb{E}_{\theta}[F(X)]$. By Lemma 1, $M_{n}^{\lambda}$ is a
martingale such that $\mathbb{E}[M_{n}^{\lambda}]\\!=\\!1$. We now introduce
the distribution
$q(\theta|\alpha,\beta)=\exp(\langle\theta,\alpha\rangle-\beta\mathcal{L}(\theta))H(\alpha,\beta)$
where $H(\cdot,\cdot)$ is the normalization function, and define the mixture
martingale $M_{n}=\int M_{n}^{\lambda}q(\theta+\lambda|\alpha,\beta)d\lambda$.
Note that $M_{n}$ also satisfies $\mathbb{E}[M_{n}]=1$.

####  Step 2. Choice of parameters and duality properties.

Choosing $\alpha=c\nabla\mathcal{L}(\theta)$ and $\beta=c$ yields

$\displaystyle
M_{n}\\!=\\!G(\theta,c)\\!\int\\!\\!\exp\\!\left(\\!\langle\lambda,n(\mu_{n}\\!-\\!\mu)\rangle\\!-\\!(n\\!+\\!c)\mathcal{B}_{\mathcal{L},\theta}(\lambda)\\!\right)d\lambda\,,$

where
$G(\theta,c)=1/\int\exp\big{(}-c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\big{)}d\theta^{\prime}$.
Now, we consider the function
$\mathcal{B}_{\mathcal{L},\theta}^{\star}(x)=\max_{\lambda}\langle\lambda,x\rangle-\mathcal{B}_{\mathcal{L},\theta}(\lambda)$.
Note that its maximal point $\lambda^{\star}_{x}$ satisfies
$x+\nabla\mathcal{L}(\theta)=\nabla\mathcal{L}(\theta+\lambda^{\star}_{x})$
for every $x$. In particular, for the choice $x=\frac{n}{n+c}(\mu_{n}-\mu)$,
we have
$\frac{n}{n+c}\mu_{n}+\frac{c}{n+c}\nabla\mathcal{L}(\theta)=\nabla\mathcal{L}(\theta+\lambda^{\star}_{x})$,
which justifies introducing the following regularized estimate

$\displaystyle\theta_{n,c}(\theta)=(\nabla\mathcal{L})^{-1}\left(\frac{n}{n+c}\mu_{n}+\frac{c}{n+c}\nabla\mathcal{L}(\theta)\right),$

provided that $\nabla\mathcal{L}$ is one-to-one in $\Theta$.

####  Step 3. Martingale rewriting and conclusion.

Note that the duality property also yields
$\langle\lambda,x\rangle-\mathcal{B}_{\mathcal{L},\theta}(\lambda)-\langle\lambda^{\star}_{x},x\rangle+\mathcal{B}_{\mathcal{L},\theta}(\lambda^{\star}_{x})=\langle\lambda-\lambda^{\star}_{x},\nabla\mathcal{L}(\theta+\lambda^{\star}_{x})\rangle+\mathcal{L}(\theta+\lambda^{\star}_{x})-\mathcal{L}(\theta+\lambda)$
for each $x$. Using this property, we obtain

$\displaystyle
M_{n}=\exp\big{(}(n\\!+\\!c)\mathcal{B}_{\mathcal{L},\theta}^{\star}(x)\big{)}\frac{G(\theta,c)}{G(\theta_{n,c}(\theta),n+c)}\,.$

After applying Markov’s inequality, we obtain the following inequality:

$\displaystyle\mathbb{P}_{\theta}\bigg{[}\mathcal{B}_{\mathcal{L},\theta}^{\star}\bigg{(}\frac{n}{n\\!+\\!c}(\mu_{n}-\mu)\bigg{)}\geqslant\frac{1}{n\\!+\\!c}\log\frac{G(\theta_{n,c}(\theta),n\\!+\\!c)}{\delta\,G(\theta,c)}\bigg{]}\leqslant\delta\,.$

Now, we use duality property of Bregman divergence to obtain
$\mathcal{B}_{\mathcal{L},\theta}^{\star}\\!\left(\\!\frac{n}{n+c}(\mu_{n}\\!-\\!\mu)\right)\\!=\\!\mathcal{B}_{\mathcal{L}}\big{(}\theta,\theta_{n,c}(\theta)\big{)}$,
which yields the desired form. Finally, by properties of the martingale
$M_{n}$, the deviation bound also holds for any random stopping time $N$, and
hence we obtain the time uniform bound by employing a stopping time
construction similar to that of Peña et al. (2008); Abbasi-Yadkori et al.
(2011).

##  4 Numerical Experiments

In this section, we provide illustrative numerical results to show the
resulting confidence bands built from Theorem 3 are promising alternatives to
existing competitors (detailed in Appendix D). All the confidence sets
presented here are implemented in the open source concentration-lib Python
package (https://pypi.org/project/concentration-lib/).

#### Comparison.

We compare our confidence envelopes to state-of-the-art time-uniform bounds in
Figure 2 and Appendix D. Interestingly, our upper and lower bounds are not
necessarily symmetrical, and in the Bernoulli case, they fit within the
distribution support without clipping, contrary to most other methods, thus
adapting to the local geometry of the family. Our bounds are comparable to
those of Kaufmann and Koolen (2021), if slightly tighter; however, we
emphasize again that our scope is wider and captures more distributions.

#### Numerical complexity.

(a) Gaussian

(b) Exponential

(c) Bernoulli

(d) Bernoulli (tuning $c$)

Figure 2: Comparison of median confidence envelopes around the mean for
$\mathcal{N}(0,1)$, Bernoulli$(0.8)$ and $\text{Exp}(1)$, as a function of the
sample size $n$, over $1000$ independent replicates. Grey lines are
trajectories of empirical means $\widehat{\mu}_{n}$. Bottom right: confidence
envelopes for varying regularization $c\\!\in\\![0.1,30]$ for Bernoulli$(0.5)$
distribution. Dotted black line: heuristic $c_{n}\\!\approx\\!0.12n$.

Some of our formulas (Poisson, Chi-square) require the evaluation of integrals
for which no closed-form expression exist. However, these can be estimated up
to arbitrary precision by numerical routines, see Remarks 10 and 11 for
details. Similarly, the use of special functions (digamma, ratio of Gamma) may
lead to numerical instabilities or overflows for large $n$; we recommend
instead to use the log-Gamma function and an efficient implementation of the
log-sum-exp operator. We report results for $n\\!\leqslant\\!200$ as we think
the small sample regime is where our bounds shine (most reasonable methods
produce similar confidence sets for large $n$). Finally, most confidence sets
reported in Figure 2 are implicitly defined (after reorganizing terms) as
level sets $\\{F(\theta)\leqslant 0\\}$ for some function $F$. For one-
dimensional families, we use a root search routine for a fast estimation of
the boundaries of such sets; in higher dimension, we evaluate $F$ on a uniform
grid (e.g, in Figure 1(a), we use a $1024\\!\times\\!1024$ grid over
$[\\!-\\!2,4]\\!\times\\![0.1,4]$).

#### Influence of the regularization parameter $c$.

Here, we fix $\delta\\!=\\!0.05$ (note that all bounds exhibit the typical
$\log\frac{1}{\delta}$ dependency). The parameter $c$ is fixed to $c\\!=\\!1$
and can be interpreted as a number of virtual prior samples, thus introducing
a bias. The choice $c\\!=\\!1$ is a typical value considered in literature on
Gaussian concentration. To see the influence of the parameter $c$ on the
bounds, we provide a second set of experiments, in which we compute confidence
envelopes for $c$ varying between $0.1$ and $30$ (Figure 2(d)). We find that
$c$ has limited influence on the width of the confidence interval. To get more
insight on the tuning of $c$, we perform additional experiments in Appendix D,
where we optimize the confidence width w.r.t. to $c$ for varying sample sizes
$n$. The resulting choice
$c^{*}_{n}\\!=\\!\mathop{\mathrm{argmin}}_{c}\lvert\Theta_{n,c}\left(\delta\right)\rvert$
exhibits a linear trend $c^{*}_{n}\\!\approx\\!0.12n$, which seems consistent
across the tested distributions. Of note, the constant $0.12$ is reminiscent
of the optimal tuning of the sub-Gaussian Laplace method (Howard et al.,
2020). This is consistent with the observation in Section 3.1 that the
confidence width $\lvert\Theta_{n,c}\left(\delta\right)\rvert$ is equivalent
to the Gaussian confidence width when $n\rightarrow+\infty$. We report in
Figure 2(d) this heuristic for $c$. However, we restate that sample a size-
dependent regularization parameter $c$ is not supported by the theory (as it
would violate the law of iterated logarithms). Now, this does not prevent from
choosing a specific horizon of interest $n_{0}$ and set $c\\!=\\!0.12n_{0}$ to
promote sharpness around $n=zn_{0}$.

##  5 Application: Generalized Likelihood Ratio Test in Exponential Families

In this section, we apply our result to control the false alarm rate of
Generalized Likelihood Ratio (GLR) tests when detecting a change of measure
from a sequence of observations. The GLR test in the exponential family model
$\mathcal{E}$ is defined, for a threshold $\alpha>0$, as

$\displaystyle\tau(\alpha;\mathcal{E})=\min\\{t\in\mathbb{N}:\max_{s\in[0,t)}G^{\mathcal{E}}_{1:s:t}\geqslant\alpha\\},$

where the GLR is defined to be
$G^{\mathcal{E}}_{1:s:t}\\!=\\!\inf_{\theta^{\prime}}\sup_{\theta_{1},\theta_{2}}\log\bigg{(}\\!\frac{\prod_{t^{\prime}=1}^{s}p_{\theta_{1}}(X_{t^{\prime}})\\!\prod_{t^{\prime}=s+1}^{t}p_{\theta_{2}}(X_{t^{\prime}})}{\prod_{t^{\prime}=1}^{t}p_{\theta^{\prime}}(X_{t^{\prime}})}\\!\bigg{)}$.

Hereafter, we consider that all the observations
$(X_{t^{\prime}})_{t^{\prime}\in\mathbb{N}}$ come from a distribution with
same parameter $\theta\in\Theta$. Then, the false alarm rate of the GLR test
can be bounded using

$\displaystyle\mathbb{P}_{\theta}(\tau(\alpha;\mathcal{E})<\infty)\\!=\\!\mathbb{P}_{\theta}\big{(}\exists(t,s)\in\mathbb{N}^{2},s<t:G^{\mathcal{E}}_{1:s:t}\geqslant\alpha\big{)}.$

Observe that solving for the inner supremums in the expression of GLR, we
obtain

$\displaystyle G^{\mathcal{E}}_{1:s:t}$ $\displaystyle=$
$\displaystyle\inf_{\theta^{\prime}}s\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{1:s})+(t-s)\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{s+1:t}),$

where we introduce parameter estimates $\theta_{1:s}$ and $\theta_{s+1:t}$
such that
$\nabla\mathcal{L}(\theta_{1:s})=\frac{1}{s}\sum_{t^{\prime}=1}^{s}F(X_{t^{\prime}})$
and
$\nabla\mathcal{L}(\theta_{s+1:t})=\frac{1}{t-s}\sum_{t^{\prime}=s+1}^{t}F(X_{t^{\prime}})$,
respectively. Hence, the false alarm rate can be controlled as

$\displaystyle\mathbb{P}_{\theta}(\tau(\alpha;\mathcal{E})\\!<\\!\infty)\\!\leqslant\\!\mathbb{P}_{\theta}\Big{(}\\!\exists
s\\!\in\\!\mathbb{N}\\!:\\!s\mathcal{B}_{\mathcal{L}}(\theta,\theta_{1:s})\\!\geqslant\\!\alpha_{1}\\!\Big{)}\\!+\\!\mathbb{P}_{\theta}\Big{(}\\!\exists(t,s)\\!\in\\!\mathbb{N}^{2},\\!s\\!<\\!t\\!:\\!(t\\!-\\!s)\mathcal{B}_{\mathcal{L}}(\theta,\theta_{s+1:t})\\!\geqslant\\!\alpha_{2}\\!\Big{)}$

for appropriate terms $\alpha_{1},\alpha_{2}$ such that
$\alpha=\alpha_{1}+\alpha_{2}$. The control of the first term comes naturally
from our time-uniform deviation result (Theorem 3) using a regularized version
of the estimate $\theta_{1:s}$. For the second term, we need to study the
concentration of $\mathcal{B}_{\mathcal{L}}(\theta,\theta_{s+1:t})$ uniformly
over $s$ and $t$.

#### Doubly time-uniform concentration.

For any $s\\!<\\!t$, $c>0$ and reference parameter $\theta_{0}\in\Theta$, we
define the regularized estimate
$\theta_{s+1:t,c}(\theta_{0})=(\nabla\mathcal{L})^{-1}\left(\frac{1}{t-s+c}\left(\sum_{t^{\prime}=s+1}^{t}F(X_{t^{\prime}})+c\nabla\mathcal{L}(\theta_{0})\right)\right)$
built from $(t\\!-\\!s)$ observations $X_{s+1},\ldots,X_{t}$. Similarly, we
introduce the corresponding Bregman information gain
$\gamma_{s+1:t,c}(\theta_{0})$ as in Definition 2. The following result
extends the Laplace method for exponential families (Theorem 3) to control the
Bregman deviation of $\theta_{s+1:t,c}(\theta_{0})$ around $\theta_{0}$.

######  Theorem 5 (Doubly time-uniform concentration)

Let $\delta\\!\in\\!(0,1]$ and
$g\colon\mathbb{N}\\!\rightarrow\\!\mathbb{R}_{+}^{*}$ such that
$\sum_{t=1}^{\infty}1/g(t)\\!\leqslant\\!1$ (e.g.,
$g(t)\\!=\\!\kappa(1\\!+\\!t)\log^{2}(1\\!+\\!t)$ where $\kappa=2.10975$).
Then, under the assumptions of Definition 2, it holds that

$\displaystyle\mathbb{P}_{\theta}\\!\left[\exists t\in\mathbb{N},\exists
s\\!<\\!t:\\!(t\\!-\\!s\\!+\\!c)\mathcal{B}_{\mathcal{L}}(\theta,\theta_{s+1:t,c}(\theta))\\!\geqslant\\!\log\frac{g(t)}{\delta}+\gamma_{s+1:t,c}(\theta)\right]\leqslant\delta\,.$

###### Remark 6

Maillard (2019a) proves a doubly uniform concentration inequality for means of
(sub)-Gaussian random variables and conjectures that it could be extended to
other types of changes, such as changes of variance in a Gaussian family.
Theorem 5 can be seen as a generalization of this result to generic
exponential families. This naturally comes with additional challenges due to
the intrinsic local geometry of exponential families. The proof of this result
is deferred to Appendix E.

#### A regularized GLR test.

For a given false alarm (i.e., false positive) probability $\delta\in(0,1]$
and parameter $c>0$, we define the regularized GLR test
$\tau^{\delta}_{c}(\mathcal{E})$ as follows:

$\displaystyle\Theta_{s,c}(\delta)=\left\\{\theta\in\Theta:(s+c)\mathcal{B}_{\mathcal{L}}(\theta,\theta_{1:s,c}(\theta))\leqslant\log(1/\delta)+\gamma_{1:s,c}(\theta)\right\\}\,,$
$\displaystyle\Theta_{s+1:t,c}(\delta)=\left\\{\theta\in\Theta:(t-s+c)\mathcal{B}_{\mathcal{L}}(\theta,\theta_{s+1:t,c}(\theta))\leqslant\log(g(t)/\delta)+\gamma_{s+1:t,c}(\theta)\right\\}\,,$
$\displaystyle\tau^{\delta}_{c}(\mathcal{E})=\min\left\\{t\in\mathbb{N}:\exists
s<t,\Theta_{s,c}(\delta/2)\cap\Theta_{s+1:t,c}(\delta/2)=\emptyset\right\\}\,.$

Note that, under the measure $\mathbb{P}_{\theta}$,
$\theta\in\Theta_{s,c}(\delta/2)$ holds with probability at least $1-\delta/2$
by time-uniform concentration over $s$ (Theorem 3), and
$\theta\in\Theta_{s+1:t,c}(\delta/2)$ holds with probability at least
$1\\!-\\!\delta/2$ by doubly time-uniform concentration over both $s$ and $t$
(Theorem 5). Then, by a union-bound, this test is guaranteed to have a false
alarm probability controlled by $\delta$. The factor $g(t)$ is essentially a
function growing slightly faster than linearly; other common choices include
$g(t)\propto t^{\eta}$ for some $\eta>1$, though it leads to a looser bound
due to the faster growth. Furthermore, if one is interested in detecting
changes up to a known horizon $T$, one can replace the function $g(t)$ by the
slightly tighter factor $g(t)=\kappa(1\\!+\\!t)\log^{2}(1\\!+\\!t)$ with
$\kappa\\!=\\!\sum_{t=1}^{T}1/g(t)$ and still guarantee a false alarm rate
under $\delta$. In Appendix E.3, we report simulations of this GLR test to
detect changes of variance in the family of centered Gaussian (Figure 14),
which answers an open question from Maillard (2019a).

##  6 Conclusion

We apply the method of mixture to derive a time-uniform deviation inequality
for generic parametric exponential families expressed in terms of their
Bregman divergences, highlighting the role of a quantity, the _Bregman
information gain_ , that is related to the geometry of the family. We
specialize this general result to build confidence sets for classical
examples. Our method compares favorably to the state-of-the art for Gaussian,
Bernoulli and exponential distributions, and extends to other families like
Chi-square, Pareto and two-parameter Gaussian, where no known time-uniform
bound exists. Our method is also general enough to design GLR tests for change
detection. An interesting direction for future work would be to consider the
case where the log-partition function (and thus the Bregman divergence) is
misspecified or unknown and has to be learned alongside the parameter.

The authors acknowledge the funding of the French National Research Agency,
the French Ministry of Higher Education and Research, Inria, the MEL and the
I-Site ULNE regarding project R-PILOTE-19-004-APPRENF and Bandits For Health
(B4H). We thank the anonymous reviewers for their careful reading of the paper
and their suggestions for improvements. Experiments presented in this paper
were carried out using the Grid’5000 testbed, supported by a scientific
interest group hosted by Inria and including CNRS, RENATER and several
universities as well as other organizations (see https://www.grid5000.fr).

## References

  * Abbasi-Yadkori et al. (2011) Yasin Abbasi-Yadkori, Dávid Pál, and Csaba Szepesvári.  Improved algorithms for linear stochastic bandits.  In _Advances in Neural Information Processing Systems_ , pages 2312–2320, 2011.
  * Amari (2016) Shun-ichi Amari.  _Information geometry and its applications_ , volume 194.  Springer, 2016.
  * Auer et al. (2002) P. Auer, N. Cesa-Bianchi, and P. Fischer.  Finite-time Analysis of the Multi-armed Bandit Problem.  _Machine Learning_ , 47(2):235–256, 2002.
  * Basu et al. (2022) Debabrota Basu, Odalric-Ambrym Maillard, and Timothée Mathieu.  Bandits corrupted by nature: Lower bounds on regret and robust optimistic algorithm.  _arXiv preprint arXiv:2203.03186_ , 2022.
  * Baudry et al. (2020) Dorian Baudry, Emilie Kaufmann, and Odalric-Ambrym Maillard.  Sub-sampling for efficient non-parametric bandit exploration.  _Advances in Neural Information Processing Systems_ , 33:5468–5478, 2020.
  * Bentkus (2004) Vidmantas Bentkus.  On hoeffding’s inequalities.  _The Annals of Probability_ , 32(2):1650–1673, 2004.
  * Berend and Kontorovich (2013) Daniel Berend and Aryeh Kontorovich.  On the concentration of the missing mass.  _Electronic Communications in Probability_ , 18:1–7, 2013\.
  * Boucheron et al. (2013) Stéphane Boucheron, Gábor Lugosi, and Pascal Massart.  _Concentration inequalities: A nonasymptotic theory of independence_.  Oxford university press, 2013.
  * Bregman (1967) L.M. Bregman.  The relaxation method of finding the common point of convex sets and its application to the solution of problems in convex programming.  _USSR Computational Mathematics and Mathematical Physics_ , 7(3):200–217, 1967.  ISSN 0041-5553.
  * Bubeck (2010) Sébastien Bubeck.  _Bandits games and clustering foundations_.  PhD thesis, 2010.
  * Cappé et al. (2013) Olivier Cappé, Aurélien Garivier, Odalric-Ambrym Maillard, Rémi Munos, and Gilles Stoltz.  Kullback-leibler upper confidence bounds for optimal sequential allocation.  _The Annals of Statistics_ , pages 1516–1541, 2013.
  * Carpentier et al. (2011) Alexandra Carpentier, Alessandro Lazaric, Mohammad Ghavamzadeh, Rémi Munos, and Peter Auer.  Upper-confidence-bound algorithms for active learning in multi-armed bandits.  In _International Conference on Algorithmic Learning Theory_ , pages 189–203. Springer, 2011.
  * Cesa-Bianchi and Lugosi (2006) Nicolo Cesa-Bianchi and Gábor Lugosi.  _Prediction, learning, and games_.  Cambridge university press, 2006.
  * Chowdhury and Gopalan (2017) Sayak Ray Chowdhury and Aditya Gopalan.  On kernelized multi-armed bandits.  In _International Conference on Machine Learning_ , pages 844–853. PMLR, 2017.
  * Chowdhury et al. (2021) Sayak Ray Chowdhury, Aditya Gopalan, and Odalric-Ambrym Maillard.  Reinforcement learning in parametric mdps with exponential families.  In _International Conference on Artificial Intelligence and Statistics_ , pages 1855–1863. PMLR, 2021.
  * Darling and Robbins (1967) DA Darling and Herbert Robbins.  Iterated logarithm inequalities.  _Proceedings of the National Academy of Sciences_ , 57(5):1188–1192, 1967.
  * Durand et al. (2018) Audrey Durand, Odalric-Ambrym Maillard, and Joelle Pineau.  Streaming kernel regression with provably adaptive mean, variance, and regularization.  _Journal of Machine Learning Research_ , 19(1):650–683, 2018.
  * Durrett (2019) Rick Durrett.  _Probability: theory and examples_ , volume 49.  Cambridge university press, 2019.
  * Dwork et al. (2014) Cynthia Dwork, Aaron Roth, et al.  The algorithmic foundations of differential privacy.  _Found. Trends Theor. Comput. Sci._ , 9(3-4):211–407, 2014.
  * Faury et al. (2020) Louis Faury, Marc Abeille, Clément Calauzènes, and Olivier Fercoq.  Improved optimistic algorithms for logistic bandits.  In _International Conference on Machine Learning_ , pages 3052–3060. PMLR, 2020.
  * Garivier and Kaufmann (2016) Aurélien Garivier and Emilie Kaufmann.  Optimal best arm identification with fixed confidence.  In _Conference on Learning Theory_ , pages 998–1027. PMLR, 2016.
  * Garivier (2013) Aurélien Garivier.  Informational confidence bounds for self-normalized averages and applications.  In _2013 IEEE Information Theory Workshop (ITW)_ , pages 1–5, 2013\.  10.1109/ITW.2013.6691311.
  * Hao et al. (2019) Botao Hao, Yasin Abbasi Yadkori, Zheng Wen, and Guang Cheng.  Bootstrapping upper confidence bound.  _Advances in Neural Information Processing Systems_ , 32:12123–12133, 2019.
  * Hoeffding (1963) Wassily Hoeffding.  Probability inequalities for sums of bounded random variables.  _Journal of the American Statistical Association_ , 58(301):13–30, 1963.  ISSN 01621459.
  * Holland and Haress (2021) Matthew Holland and El Mehdi Haress.  Learning with risk-averse feedback under potentially heavy tails.  In _International Conference on Artificial Intelligence and Statistics_ , pages 892–900. PMLR, 2021.
  * Howard et al. (2020) Steven R Howard, Aaditya Ramdas, Jon McAuliffe, and Jasjeet Sekhon.  Time-uniform chernoff bounds via nonnegative supermartingales.  _Probability Surveys_ , 17:257–317, 2020.
  * Howard et al. (2021) Steven R Howard, Aaditya Ramdas, Jon McAuliffe, and Jasjeet Sekhon.  Time-uniform, nonparametric, nonasymptotic confidence sequences.  _The Annals of Statistics_ , 2021.
  * Jaksch et al. (2010) Thomas Jaksch, Ronald Ortner, and Peter Auer.  Near-optimal regret bounds for reinforcement learning.  _Journal of Machine Learning Research_ , 11(4), 2010.
  * Kaufmann and Koolen (2021) Emilie Kaufmann and Wouter M Koolen.  Mixture martingales revisited with applications to sequential tests and confidence intervals.  _Journal of Machine Learning Research_ , 22:246–1, 2021\.
  * Kearns and Saul (1998) Michael Kearns and Lawrence Saul.  Large deviation methods for approximate probabilistic inference.  In _Proceedings of the Fourteenth conference on Uncertainty in artificial intelligence_ , pages 311–319, 1998.
  * Kirschner and Krause (2018) Johannes Kirschner and Andreas Krause.  Information directed sampling and bandits with heteroscedastic noise.  In _Conference On Learning Theory_ , pages 358–384. PMLR, 2018.
  * Korda et al. (2013) Nathaniel Korda, Emilie Kaufmann, and Remi Munos.  Thompson sampling for 1-dimensional exponential family bandits.  _Advances in neural information processing systems_ , 26, 2013.
  * Kuchibhotla and Zheng (2021) Arun K Kuchibhotla and Qinqing Zheng.  Near-optimal confidence sequences for bounded random variables.  In Marina Meila and Tong Zhang, editors, _Proceedings of the 38th International Conference on Machine Learning_ , volume 139 of _Proceedings of Machine Learning Research_ , pages 5827–5837. PMLR, 18–24 Jul 2021.
  * Lattimore and Szepesvári (2019) T. Lattimore and C. Szepesvári.  _Bandit Algorithms_.  Cambridge University Press, 2019.
  * Laurent and Massart (2000) Beatrice Laurent and Pascal Massart.  Adaptive estimation of a quadratic functional by model selection.  _Annals of Statistics_ , pages 1302–1338, 2000.
  * Maillard (2019a) O.-A. Maillard.  Sequential change-point detection: Laplace concentration of scan statistics and non-asymptotic delay bounds.  In _Algorithmic Learning Theory (ALT)_ , 2019a.
  * Maillard (2019b) Odalric-Ambrym Maillard.  _Mathematics of Statistical Sequential Decision Making_.  Habilitation à diriger des recherches, Université de Lille Nord de France, February 2019b.
  * Peña et al. (2008) Victor H Peña, Tze Leung Lai, and Qi-Man Shao.  _Self-normalized processes: Limit theory and Statistical Applications_.  Springer Science & Business Media, 2008.
  * Raginsky and Sason (2018) Maxim Raginsky and Igal Sason.  _Concentration of Measure Inequalities in Information Theory, Communications, and Coding_.  Now Publishers, 2018.
  * Robbins (1970) Herbert Robbins.  Statistical methods related to the law of the iterated logarithm.  _The Annals of Mathematical Statistics_ , 41(5):1397–1409, 1970.
  * Robbins and Pitman (1949) Herbert Robbins and EJG Pitman.  Application of the method of mixtures to quadratic forms in normal variates.  _The annals of mathematical statistics_ , pages 552–560, 1949.
  * Shafer and Vovk (2019) Glenn Shafer and Vladimir Vovk.  _Game-Theoretic Foundations for Probability and Finance_ , volume 455\.  John Wiley & Sons, 2019.
  * Shun and McCullagh (1995) Zhenming Shun and Peter McCullagh.  Laplace approximation of high dimensional integrals.  _Journal of the Royal Statistical Society: Series B (Methodological)_ , 57(4):749–760, 1995.
  * Ville (1939) Jean Ville.  Etude critique de la notion de collectif.  _Bull. Amer. Math. Soc_ , 45(11):824, 1939.
  * Waudby-Smith and Ramdas (2023) Ian Waudby-Smith and Aaditya Ramdas.  Estimating means of bounded random variables by betting.  _Journal of the Royal Statistical Society: Series B (Statistical Methodology)_ , 2023.
  * Zeitouni and Dembo (1998) A Dembo O Zeitouni and A Dembo.  Large deviations techniques and applications.  _Applications of Mathematics_ , 38, 1998.

##  Appendix A Proof of the Main Result About Bregman Deviations

In this section, we detail the proof of the main results regarding Bregman
deviation inequalities, stated with its two variants Theorem 3 and 7, We start
in Section A.1 with the proof of Theorem 3 when considering a numerical
constant $c$ to build the regularizer. In Section A.2, we then provide the
proof for Theorem 7 when considering instead a Legendre function
$\mathcal{L}_{0}$ as a regularizer, and also extend to the non i.i.d. case.

###  A.1 Proof of Bregman concentration using mixture parameter $c$

The proof of this first result follows an adaptation of the method of mixture,
combined with properties of Bregman divergences recalled in Lemma 1.

See 3

Proof of Theorem 3:

Step 1. Martingale and mixture martingale construction.

Let us note that $\mathbb{E}_{\theta}[F(X)]=\nabla\mathcal{L}(\theta)$ and
$\log\mathbb{E}_{\theta}\left[\exp\langle\lambda,F(X)\rangle\right]=\mathcal{L}(\theta+\lambda)-\mathcal{L}(\theta)$.
Hence, we deduce that

$\displaystyle\log\mathbb{E}_{\theta}\left[\exp\langle\lambda,F(X)-\mathbb{E}_{\theta}[F(X)]\rangle\right]=\mathcal{B}_{\mathcal{L},\theta}(\lambda)\stackrel{{\scriptstyle\rm
def}}{{=}}\mathcal{L}(\theta+\lambda)-\mathcal{L}(\theta)-\langle\lambda,\nabla\mathcal{L}(\theta)\rangle\,.$

Now, let $\mu_{n}=\frac{1}{n}\sum_{t=1}^{n}F(X_{t})$ and
$\mu=\mathbb{E}_{\theta}[F(X)]$. For any $\lambda\in\mathbb{R}^{d}$, the
following quantity

$\displaystyle M_{n}^{\lambda}$ $\displaystyle=$
$\displaystyle\exp\bigg{(}\langle\lambda,n(\mu_{n}-\mu)\rangle-n\mathcal{B}_{\mathcal{L},\theta}(\lambda)\bigg{)}\,,$

is thus a martingale such that $\mathbb{E}[M_{n}^{\lambda}]=1$.

We now introduce the distribution
$q(\theta|\alpha,\beta)=H(\alpha,\beta)\cdot\exp(\langle\theta,\alpha\rangle-\beta\mathcal{L}(\theta))$
where $H$ is the normalization term. We further introduce the following
quantity

$\displaystyle
M_{n}=\int_{\Lambda_{\theta}}M_{n}^{\lambda}q(\theta+\lambda|\alpha,\beta)d\lambda\,,$

where $\Lambda_{\theta}=\\{\lambda:\theta+\lambda\in\Theta\\}$. This also
satisfies $\mathbb{E}[M_{n}]=1$. Further, we have the rewriting

$\displaystyle
M_{n}=H(\alpha,\beta)\int_{\Lambda_{\theta}}\exp\bigg{(}\langle\lambda,n(\mu_{n}-\mu)\rangle-n\mathcal{B}_{\mathcal{L},\theta}(\lambda)+\langle\lambda+\theta,\alpha\rangle-\beta\mathcal{L}(\theta+\lambda)\bigg{)}d\lambda\,.$

Step 2. Choice of parameters and duality properties. Considering in particular
the choice $\alpha=c\nabla\mathcal{L}(\theta)$ and $\beta=c$, we get

$\displaystyle M_{n}$ $\displaystyle=$ $\displaystyle
H(c\nabla\mathcal{L}(\theta),c)\int_{\Lambda_{\theta}}\exp\bigg{(}\langle\lambda,n(\mu_{n}-\mu)\rangle-n\mathcal{B}_{\mathcal{L},\theta}(\lambda)+c\langle\lambda+\theta,\nabla\mathcal{L}(\theta)\rangle-c\mathcal{L}(\theta+\lambda)\bigg{)}d\lambda$
$\displaystyle=$ $\displaystyle
G(\theta,c)\int_{\Lambda_{\theta}}\exp\bigg{(}\langle\lambda,n(\mu_{n}-\mu)\rangle-(n+c)\mathcal{B}_{\mathcal{L},\theta}(\lambda)\bigg{)}d\lambda~{},$

where we introduce the following quantity

$\displaystyle
G(\theta,c)=H(c\nabla\mathcal{L}(\theta),c)\exp\bigg{(}c\langle\theta,\nabla\mathcal{L}(\theta)\rangle-c\mathcal{L}(\theta)\bigg{)}$
$\displaystyle=$
$\displaystyle\frac{\exp\bigg{(}c\langle\theta,\nabla\mathcal{L}(\theta)\rangle-c\mathcal{L}(\theta)\bigg{)}}{\int_{\Theta}\exp\bigg{(}c\langle\theta^{\prime},\nabla\mathcal{L}(\theta)\rangle-c\mathcal{L}(\theta^{\prime})\bigg{)}d\theta^{\prime}}$
$\displaystyle=$
$\displaystyle\bigg{[}\int_{\Theta}\exp\big{(}-c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\big{)}d\theta^{\prime}\bigg{]}^{-1}~{}.$

At this point, let us introduce $x=\frac{n}{n+c}(\mu_{n}-\mu)$. We also
consider the Legendre-Fenchel dual function
$\mathcal{B}_{\mathcal{L},\theta}^{\star}(x)=\max_{\lambda}\langle\lambda,x\rangle-\mathcal{B}_{\mathcal{L},\theta}(\lambda)$
and denote $\lambda^{\star}_{x}$ its maximal point. Now, we note that
$x+\nabla\mathcal{L}(\theta)=\nabla\mathcal{L}(\theta+\lambda^{\star}_{x})$
and thus

$\displaystyle\langle\lambda,x\rangle-\mathcal{B}_{\mathcal{L},\theta}(\lambda)-\langle\lambda^{\star}_{x},x\rangle+\mathcal{B}_{\mathcal{L},\theta}(\lambda^{\star}_{x})$
$\displaystyle=$
$\displaystyle\langle\lambda-\lambda^{\star}_{x},x+\nabla\mathcal{L}(\theta)\rangle+\mathcal{L}(\theta+\lambda^{\star}_{x})-\mathcal{L}(\theta+\lambda)$
$\displaystyle=$
$\displaystyle\langle\lambda-\lambda^{\star}_{x},\nabla\mathcal{L}(\theta+\lambda^{\star}_{x})\rangle+\mathcal{L}(\theta+\lambda^{\star}_{x})-\mathcal{L}(\theta+\lambda)\,.$

Also, since it holds that
$x=\frac{n}{n+c}(\mu_{n}-\nabla\mathcal{L}(\theta))$, then the quantity
$x+\nabla\mathcal{L}(\theta)=\nabla\mathcal{L}(\theta+\lambda^{\star}_{x})$
rewrites
$\frac{n}{n+c}\mu_{n}+\frac{c}{n+c}\nabla\mathcal{L}(\theta)=\nabla\mathcal{L}(\theta+\lambda^{\star}_{x})$,
which justifies to introduce the following regularized parameter estimate

$\displaystyle\theta_{n,c}(\theta)=(\nabla\mathcal{L})^{-1}\bigg{(}\frac{n}{n+c}\mu_{n}+\frac{c}{n+c}\nabla\mathcal{L}(\theta)\bigg{)}\,.$

Step 3. Martingale rewriting and conclusion. Using this property, we note that

$\displaystyle M_{n}$ $\displaystyle=$
$\displaystyle\exp\big{(}(n\\!+\\!c)\mathcal{B}_{\mathcal{L},\theta}^{\star}(x)\big{)}\cdot
G(\theta,c)\\!\\!\int_{\Lambda_{\theta}}\\!\\!\exp\\!\bigg{(}\\!(n\\!+\\!c)\Big{[}\langle\lambda\\!-\\!\lambda^{\star}_{x},\nabla\mathcal{L}(\theta\\!+\\!\lambda^{\star}_{x})\rangle\\!-\\!\mathcal{L}(\theta\\!+\\!\lambda)\\!+\\!\mathcal{L}(\theta\\!+\\!\lambda^{\star}_{x})\Big{]}\\!\bigg{)}d\lambda$
$\displaystyle=$
$\displaystyle\exp\big{(}(n\\!+\\!c)\mathcal{B}_{\mathcal{L},\theta}^{\star}(x)\big{)}\cdot
G(\theta,c)\\!\int_{\Lambda_{\theta}}\\!\exp\\!\bigg{(}\langle\lambda\\!+\\!\theta,(n\\!+\\!c)\nabla\mathcal{L}(\theta\\!+\\!\lambda^{\star}_{x})\rangle\\!-\\!(n\\!+\\!c)\mathcal{L}(\theta\\!+\\!\lambda)\bigg{)}d\lambda\times$
$\displaystyle\exp\bigg{(}\\!-(n\\!+\\!c)\langle\lambda^{\star}_{x}\\!+\\!\theta,\nabla\mathcal{L}(\theta\\!+\\!\lambda^{\star}_{x})+(n\\!+\\!c)\mathcal{L}(\theta\\!+\\!\lambda^{\star}_{x})\rangle\bigg{)}$
$\displaystyle=$
$\displaystyle\exp\big{(}(n\\!+\\!c)\mathcal{B}_{\mathcal{L},\theta}^{\star}(x)\big{)}\frac{G(\theta,c)}{G(\theta+\lambda^{\star}_{x},n+c)}=\exp\big{(}(n\\!+\\!c)\mathcal{B}_{\mathcal{L},\theta}^{\star}(x)\big{)}\frac{G(\theta,c)}{G(\theta_{n,c}(\theta),n+c)}\,.$

We then apply simple Markov inequality, which yields for all constant $c$,

$\displaystyle\forall\delta\in(0,1],\quad\mathbb{P}\bigg{(}\mathcal{B}_{\mathcal{L},\theta}^{\star}\bigg{(}\frac{n}{n\\!+\\!c}(\mu_{n}\\!-\\!\mu)\bigg{)}\geqslant\frac{1}{n+c}\log\bigg{(}\frac{G(\theta_{n,c}(\theta),n+c)}{G(\theta,c)}\frac{1}{\delta}\bigg{)}\bigg{)}\leqslant\delta.$

To conclude, we use the duality property of the Bregman divergence (Lemma 1),
considering that $\nabla\mathcal{L}$ is invertible. Indeed, denoting
$\alpha=n/(n+c)$, and $\theta_{n}=(\nabla\mathcal{L})^{-1}(\mu_{n})$, it comes

$\displaystyle\mathcal{B}_{\mathcal{L},\theta}^{\star}\bigg{(}\alpha(\mu_{n}-\mu)\bigg{)}=\mathcal{B}_{\mathcal{L}}\big{(}\theta,\theta_{\alpha}\big{)}\,\qquad\text{where}\,\,\theta_{\alpha}=\nabla\mathcal{L}^{-1}\big{(}\alpha\nabla\mathcal{L}(\theta_{n})+(1-\alpha)\nabla\mathcal{L}(\theta)\big{)}=\theta_{n,c}(\theta)\,.$

Last, we note that this extends from a single $n$ to being time-uniform thanks
to Doob’s maximal inequality for nonnegative supermartingale applied to
$M_{n}$ (also known as Vile’s inequality Ville (1939)). $\hfill\square$

###  A.2 Bregman concentration using Legendre function $\mathcal{L}_{0}$

In this section, we state a more general result that can handle the case of
sequence $(X_{t})_{t\geqslant 1}$ of random variables that are not
independent, each having possibly different distribution from others.

######  Theorem 7 (Laplace method for non i.i.d. samples)

Suppose $\\{\mathcal{H}_{t}\\}_{t=0}^{\infty}$ is a filtration such that for
each $t\geqslant 1$, (i) $X_{t}$ is $\mathcal{H}_{t}$-measurable, (ii) $F_{t}$
and $h_{t}$ are $\mathcal{H}_{t-1}$-measurable and (iii) given
$\mathcal{H}_{t-1}$, $X_{t}\sim p_{\theta_{0},t}$ where $p_{\theta_{0},t}$
belongs to an exponential family with parameter $\theta_{0}$, feature function
$F_{t}$ and base function $h_{t}$. Let $\mathcal{L}_{t}$ be the log-partition
function corresponding to $F_{t}$. For any Legendre (i.e., strictly convex and
continuously differentiable) function $\mathcal{L}_{0}:\Theta\to\mathbb{R}$
such that $\exp(-\mathcal{L}_{0})$ is integrable and any $n\in\mathbb{N}$, we
introduce the parameter estimate and Bregman information gain

$\displaystyle\theta_{n,\mathcal{L}_{0}}$
$\displaystyle=\Big{(}\sum_{t=0}^{n}\nabla\mathcal{L}_{t}\Big{)}^{-1}\Big{(}\sum_{t=1}^{n}F_{t}(X_{t})\Big{)}\,,$
$\displaystyle\gamma_{n,\mathcal{L}_{0}}$
$\displaystyle=\log\bigg{(}\frac{\int_{\Theta}\exp\Big{(}-\mathcal{L}_{0}(\theta^{\prime})\Big{)}d\theta^{\prime}}{\int_{\Theta}\exp\Big{(}-\sum_{t=0}^{n}\mathcal{B}_{\mathcal{L}_{t}}(\theta^{\prime},\theta_{n,\mathcal{L}_{0}})\Big{)}d\theta^{\prime}}\bigg{)}\,,$

and then, for any $\delta\in(0,1]$, the set

$\displaystyle\Theta_{n,\mathcal{L}_{0}}(\delta)=\left\\{\theta:\sum_{t=0}^{n}\mathcal{B}_{\mathcal{L}_{t}}\left(\theta,\theta_{n,\mathcal{L}_{0}}\right)-\mathcal{L}_{0}(\theta)\leqslant\log\frac{1}{\delta}+\gamma_{n,\mathcal{L}_{0}}\right\\}\,.$

Then, the following time-uniform concentration holds:

$\displaystyle\mathbb{P}\big{[}\exists
n\in\mathbb{N}:\theta_{0}\notin\Theta_{n,\mathcal{L}_{0}}(\delta)\big{]}\leqslant\delta\,.$

In contrast to Theorem 3 that involves a local regularization using the true
parameter $\theta_{0}\in\Theta$ of the family and a constant $c\\!>\\!0$, here
we make use of a global regularization, in the form of the Legendre function
$\mathcal{L}_{0}$. Concretely, the regularized parameter estimate
$\theta_{n,\mathcal{L}_{0}}$ and the Bregman information gain
$\gamma_{n,\mathcal{L}_{0}}$ do not depend on the true parameter $\theta_{0}$,
thus making for a more explicit confidence set $\Theta_{n,\mathcal{L}_{0}}$.
The trade-off here is that the choice of regularizer is limited by the
integrability assumption on $\exp(-\mathcal{L}_{0})$, which is critical to
build an appropriate prior for the method of mixtures.

The proof of this result follows a similar line of proof as before. However
the use of $\mathcal{L}_{0}$ instead of $c$ induces a few changes that we
detail below. In particular, we use a different prior to build the mixture of
martingales.

Proof of Theorem 7:

Step 1. Martingale and mixture martingale construction. For any
$\lambda\in\mathbb{R}^{d}$, we define

$\displaystyle
M_{n}^{\lambda}=\exp\left(\sum_{t=1}^{n}\left(\lambda^{\top}\left(F_{t}(X_{t})-\mathbb{E}_{\theta}\left[F_{t}(X_{t})\right]\right)-B_{\mathcal{L}_{t},\theta}(\lambda)\right)\right)~{},$

where we introduced
$B_{\mathcal{L}_{t},\theta}(\lambda)=\mathcal{L}_{t}(\theta+\lambda)-\mathcal{L}_{t}(\theta)-\langle\lambda,\nabla\mathcal{L}_{t}(\theta)\rangle$
for convenience.

Note that $M_{n}^{\lambda}>0$ and

$\displaystyle\mathbb{E}_{\theta}\left[\exp\left(\lambda^{\top}F_{t}(X_{t})\right)\right]$
$\displaystyle=\exp\left(\mathcal{L}_{t}(\theta+\lambda)-\mathcal{L}_{t}(\theta)\right)~{},$
$\displaystyle\mathbb{E}_{\theta}\left[\exp\left(\lambda^{\top}\mathbb{E}_{\theta}[F_{t}(X_{t})\right)\right]$
$\displaystyle=\exp\left(\lambda^{\top}\mathcal{L}_{t}(\theta)\right)~{}.$

Note that $M_{n}^{\lambda}$ is $\mathcal{H}_{n}$-measurable and in fact
$\mathbb{E}[M_{n}^{\lambda}|\mathcal{H}_{n-1}]=M_{n-1}^{\lambda}$. Therefore
$\\{M_{n}^{\lambda}\\}_{n=0}^{\infty}$ is a nonnegative martingale adapted to
the filtration $\\{\mathcal{H}_{n}\\}_{n=0}^{\infty}$ and actually satisfies
$\mathbb{E}\left[M_{n}^{\lambda}\right]=1~{}$. For any prior density
$q(\theta)$ for $\theta$, we now define a mixture of martingales

$\displaystyle
M_{n}=\int_{\Lambda_{\theta}}M_{n}^{\lambda}q(\theta+\lambda)d\lambda~{}.$ (1)

where $\Lambda_{\theta}=\\{\lambda:\theta+\lambda\in\Theta\\}$. Then $M_{n}$
is also a martingale and actually satisfies
$\mathbb{E}\left[M_{n}\right]=1~{}$. Now consider the prior density

$\displaystyle
q(\theta+\lambda)=c_{0}\cdot\exp\left(-\mathcal{L}_{0}(\theta+\lambda)\right)~{},$
(2)

where
$c_{0}=\frac{1}{\int_{\Theta}\exp\left(-\mathcal{L}_{0}(\theta^{\prime})\right)d\theta^{\prime}}$
(which is well-defined since $\exp(-\mathcal{L}_{0})$ is integrable over
$\Theta$). We then have

$\displaystyle
M_{n}=c_{0}\int_{\Lambda_{\theta}}\exp\left(\lambda^{\top}S_{n}-\sum_{t=1}^{n}\mathcal{B}_{\mathcal{L}_{t},\theta}(\lambda)-\mathcal{L}_{0}(\theta+\lambda)\right)d\lambda~{},$

where we denote
$S_{n}=\sum_{t=1}^{n}\left(F_{t}(X_{t})-\mathbb{E}_{\theta}\left[F_{t}(X_{t})\right]\right)$.
Now, from the formula of parameter estimate, we have

$\displaystyle\quad\sum_{t=1}^{n}\nabla\mathcal{L}_{t}(\theta_{n,\mathcal{L}_{0}})+\nabla\mathcal{L}_{0}(\theta_{n,\mathcal{L}_{0}})=\sum_{t=1}^{n}F_{t}(X_{t})~{}.$
(3)

This yields
$S_{n}=\sum_{t=1}^{n}\left(\nabla\mathcal{L}_{t}(\theta_{n,\mathcal{L}_{0}})-\nabla\mathcal{L}_{t}(\theta)\right)+\nabla\mathcal{L}_{0}(\theta_{n,\mathcal{L}_{0}})$.

Step 2. Legendre function and Bregman properties. We now introduce the
function $\mathcal{L}_{1:n}(\theta)=\sum_{t=1}^{n}\mathcal{L}_{t}(\theta)$.
Note that $\mathcal{L}_{1:n}$ is a also Legendre function and its associated
Bregman divergence satisfies

$\displaystyle\mathcal{B}_{\mathcal{L}_{1:n}}(\theta^{\prime},\theta)=\sum_{t=1}^{n}\left(\mathcal{L}_{t}(\theta^{\prime})-\mathcal{L}_{t}(\theta)-(\theta^{\prime}-\theta)^{\top}\nabla\mathcal{L}_{t}(\theta)\right)=\sum_{t=1}^{n}\mathcal{B}_{\mathcal{L}_{t}}(\theta^{\prime},\theta).$

In this notation, we can rewrite
$\sum_{t=1}^{n}\mathcal{B}_{\mathcal{L}_{t},\theta}(\lambda)=\mathcal{B}_{\mathcal{L}_{1:n},\theta}(\lambda)$
and
$S_{n}=\nabla\mathcal{L}_{n}(\theta_{n,\mathcal{L}_{0}})-\nabla\mathcal{L}_{n}(\theta)+\nabla\mathcal{L}_{0}(\theta_{n,\mathcal{L}_{0}})$.
We then obtain

$\displaystyle
M_{n}=c_{0}\cdot\exp\left(-\mathcal{L}_{0}(\theta)\right)\int_{\Lambda_{\theta}}\exp\left(\lambda^{\top}x_{n}-\mathcal{B}_{\mathcal{L}_{1:n},\theta}(\lambda)+\lambda^{\top}x_{0}-\mathcal{B}_{\mathcal{L}_{0},\theta}(\lambda)\right)d\lambda~{},$
(4)

where we have introduced
$x_{0}=\nabla\mathcal{L}_{0}(\theta_{n,\mathcal{L}_{0}})-\nabla\mathcal{L}_{0}(\theta)$
and
$x_{n}=\nabla\mathcal{L}_{1:n}(\theta_{n,\mathcal{L}_{0}})-\nabla\mathcal{L}_{1:n}(\theta)$.

We now have from the Bregman-duality property recalled in Lemma 1 that

$\displaystyle\sup_{\lambda\in\mathbb{R}^{d}}\left(\lambda^{\top}x_{n}-\mathcal{B}_{\mathcal{L}_{1:n},\theta}(\lambda)\right)$
$\displaystyle=\mathcal{B}_{\mathcal{L}_{1:n},\theta}^{\star}(x_{n})$
$\displaystyle=\mathcal{B}_{\mathcal{L}_{1:n},\theta}^{\star}(\nabla\mathcal{L}_{1:n}(\theta_{n,\mathcal{L}_{0}})-\nabla\mathcal{L}_{1:n}(\theta))=\mathcal{B}_{\mathcal{L}_{1:n}}(\theta,\theta_{n,\mathcal{L}_{0}})~{}.$

Further, any optimal $\lambda$ must satisfy

$\displaystyle\nabla\mathcal{L}_{1:n}(\theta+\lambda)-\nabla\mathcal{L}_{1:n}(\theta)=x_{n}\implies\nabla\mathcal{L}_{1:n}(\theta+\lambda)=\nabla\mathcal{L}_{1:n}(\theta_{n,\mathcal{L}_{0}})~{}.$

One possible solution is $\lambda^{\star}=\theta_{n,\mathcal{L}_{0}}-\theta$.
We then have

$\displaystyle\lambda^{\top}x_{n}-\mathcal{B}_{\mathcal{L}_{1:n},\theta}(\lambda)$
$\displaystyle=\lambda^{\top}x_{n}-\mathcal{B}_{\mathcal{L}_{1:n},\theta}(\lambda)+\mathcal{B}_{\mathcal{L}_{1:n}}(\theta,\theta_{n,\mathcal{L}_{0}})-\left({\lambda^{\star}}^{\top}x_{n}-\mathcal{B}_{\mathcal{L}_{1:n},\theta}(\lambda^{\star})\right)$
$\displaystyle=\mathcal{B}_{\mathcal{L}_{1:n}}(\theta,\theta_{n,\mathcal{L}_{0}})+(\lambda-\lambda^{\star})^{\top}\nabla\mathcal{L}_{n}(\theta+\lambda^{\star})+\mathcal{B}_{\mathcal{L}_{1:n},\theta}(\lambda^{\star})-\mathcal{B}_{\mathcal{L}_{1:n},\theta}(\lambda)-(\lambda-\lambda^{\star})^{\top}\nabla\mathcal{L}_{1:n}(\theta)$
$\displaystyle=\mathcal{B}_{\mathcal{L}_{1:n}}(\theta,\theta_{n,\mathcal{L}_{0}})+(\lambda-\lambda^{\star})^{\top}\nabla\mathcal{L}_{1:n}(\theta+\lambda^{\star})+\mathcal{L}_{1:n}(\theta+\lambda^{\star})-\mathcal{L}_{1:n}(\theta+\lambda)~{}.$
(5)

Similarly, we have

$\displaystyle\lambda^{\top}x_{0}-\mathcal{B}_{\mathcal{L}_{0},\theta}(\lambda)=\mathcal{B}_{\mathcal{L}_{0}}(\theta,\theta_{n,\mathcal{L}_{0}})+(\lambda-\lambda^{\star})^{\top}\nabla\mathcal{L}_{0}(\theta+\lambda^{\star})+\mathcal{L}_{0}(\theta+\lambda^{\star})-\mathcal{L}_{0}(\theta+\lambda)\,.$
(6)

Step 3. Martingale rewriting and conclusion. Plugin-in (6) and (5) in (4), we
now obtain

$\displaystyle M_{n}$
$\displaystyle=c_{0}\exp\Big{(}\sum_{j\in\\{0,1:n\\}}\mathcal{B}_{\mathcal{L}_{j}}(\theta,\theta_{n,\mathcal{L}_{0}})-\mathcal{L}_{0}(\theta)\Big{)}$
$\displaystyle\quad\quad\times\int_{\Lambda_{\theta}}\exp\Big{(}\sum_{j\in\\{0,1:n\\}}\left((\lambda-\lambda^{\star})^{\top}\nabla\mathcal{L}_{j}(\theta+\lambda^{\star})+\mathcal{L}_{j}(\theta+\lambda^{\star})-\mathcal{L}_{j}(\theta+\lambda)\right)\Big{)}d\lambda$
$\displaystyle=c_{0}\exp\Big{(}\sum_{j\in\\{0,1:n\\}}\mathcal{B}_{\mathcal{L}_{j}}(\theta,\theta_{n,\mathcal{L}_{0}})-\mathcal{L}_{0}(\theta)\Big{)}\cdot\exp\Big{(}-\sum_{j\in\\{0,1:n\\}}\left((\theta+\lambda^{\star})^{\top}\nabla\mathcal{L}_{j}(\theta+\lambda^{\star})-\mathcal{L}_{j}(\theta+\lambda^{\star})\right)\Big{)}$
$\displaystyle\quad\quad\times\int_{\Lambda_{\theta}}\exp\Big{(}\sum_{j\in\\{0,1:n\\}}\left((\theta+\lambda)^{\top}\nabla\mathcal{L}_{j}(\theta+\lambda^{\star})-\mathcal{L}_{j}(\theta+\lambda)\right)\Big{)}d\lambda$
$\displaystyle=\frac{c_{0}}{c_{n}}\exp\Big{(}\sum_{j\in\\{0,1:n\\}}\mathcal{B}_{\mathcal{L}_{j}}(\theta,\theta_{n,\mathcal{L}_{0}})-\mathcal{L}_{0}(\theta)\Big{)}$
$\displaystyle\quad\quad\times\frac{\int_{\Lambda_{\theta}}\exp\left(\sum_{j\in\\{0,1:n\\}}\left((\theta+\lambda)^{\top}\nabla\mathcal{L}_{j}(\theta+\lambda^{\star})-\mathcal{L}_{j}(\theta+\lambda)\right)\right)d\lambda}{\int_{\Theta}\exp\left(\sum_{j\in\\{0,1:n\\}}\left((\theta^{\prime})^{\top}\nabla\mathcal{L}_{j}(\theta+\lambda^{\star})-\mathcal{L}_{j}(\theta^{\prime})\right)\right)d\theta^{\prime}}$
$\displaystyle=\frac{c_{0}}{c_{n}}\exp\Big{(}\sum_{t=1}^{n}\mathcal{B}_{\mathcal{L}_{t}}(\theta,\theta_{n,\mathcal{L}_{0}})+\mathcal{B}_{\mathcal{L}_{0}}(\theta,\theta_{n,\mathcal{L}_{0}})-\mathcal{L}_{0}(\theta)\Big{)}~{},$

where we have introduced
$c_{n}=\frac{\exp\left(\sum_{j\in\\{0,1:n\\}}\left((\theta+\lambda^{\star})^{\top}\nabla\mathcal{L}_{j}(\theta+\lambda^{\star})-\mathcal{L}_{j}(\theta+\lambda^{\star})\right)\right)}{\int_{\Theta}\exp\left(\sum_{j\in\\{0,1:n\\}}\left((\theta^{\prime})^{\top}\nabla\mathcal{L}_{j}(\theta+\lambda^{\star})-\mathcal{L}_{j}(\theta^{\prime})\right)\right)d\theta^{\prime}}~{}.$

We then obtain the following result by a simple application of Markov’s
inequality, which yields

$\displaystyle\mathbb{P}\left[\sum_{t=1}^{n}\mathcal{B}_{\mathcal{L}_{t}}(\theta,\theta_{n,\mathcal{L}_{0}})+\mathcal{B}_{\mathcal{L}_{0}}(\theta,\theta_{n,\mathcal{L}_{0}})-\mathcal{L}_{0}(\theta)\geqslant\log\left(\frac{c_{n}}{c_{0}}\frac{1}{\delta}\right)\right]$
$\displaystyle=\mathbb{P}\left[M_{n}\geqslant\frac{1}{\delta}\right]\leqslant\delta\cdot\mathbb{E}\left[M_{n}\right]=\delta~{}.$
(7)

Finally, since $\lambda^{\star}=\theta_{n,\mathcal{L}_{0}}-\theta$, we have
the more explicit form

$\displaystyle\frac{c_{n}}{c_{0}}$
$\displaystyle=\frac{\int_{\Theta}\exp\Big{(}-\mathcal{L}_{0}(\theta^{\prime})\Big{)}d\theta^{\prime}}{\int_{\Theta}\exp\Big{(}-\sum\limits_{j\in\\{0,1:n\\}}\mathcal{B}_{\mathcal{L}_{j}}(\theta^{\prime},\theta+\lambda^{\star})\Big{)}d\theta^{\prime}}$
$\displaystyle=\frac{\int_{\Theta}\exp\Big{(}-\mathcal{L}_{0}(\theta^{\prime})\Big{)}d\theta^{\prime}}{\int_{\Theta}\exp\Big{(}-\sum_{t=1}^{n}\mathcal{B}_{\mathcal{L}_{t}}(\theta^{\prime},\theta_{n,\mathcal{L}_{0}})-\mathcal{B}_{\mathcal{L}_{0}}(\theta^{\prime},\theta_{n,\mathcal{L}_{0}})\Big{)}d\theta^{\prime}}~{}.$

We obtain the time-uniform bound by applying Doob’s maximal inequality for
nonnegative supermartingale. $\hfill\square$

The general result of Theorem 7 specifies straightforwardly to the case when
all observations have same law, yielding the following corollary that we state
now for completeness.

###### Corollary 8

Let $\mu_{n}=\frac{1}{n}\sum_{t=1}^{n}F(X_{t})$. For all Legendre function
$\mathcal{L}_{0}$ such that $\exp(-\mathcal{L}_{0})$ is integrable over
$\Theta$, let

$\displaystyle\theta_{n,\mathcal{L}_{0}}\\!=\big{(}\nabla\mathcal{L}\\!+\\!\frac{1}{n}\nabla\mathcal{L}_{0}\big{)}^{-1}(\mu_{n})\,\text{
and
}\,\gamma_{n,\mathcal{L}_{0}}\\!=\log\left(\frac{\int_{\Theta}\exp\Big{(}-\mathcal{L}_{0}(\theta^{\prime})\Big{)}d\theta^{\prime}}{\int_{\Theta}\exp\Big{(}-\\!n\mathcal{B}_{\mathcal{L}}\left(\theta^{\prime},\theta_{n,\mathcal{L}_{0}}\right)\\!-\\!\mathcal{B}_{\mathcal{L}_{0}}\left(\theta^{\prime},\theta_{n,\mathcal{L}_{0}}\right)\Big{)}d\theta^{\prime}}\right)\,.$

Then for all $\delta\in(0,1],$

$\displaystyle\mathbb{P}\bigg{(}\exists
n\in\mathbb{N},\,\,n\mathcal{B}_{\mathcal{L}}\left(\theta,\theta_{n,\mathcal{L}_{0}}\right)+\mathcal{B}_{\mathcal{L}_{0}}\left(\theta,\theta_{n,\mathcal{L}_{0}}\right)\geqslant\log\frac{1}{\delta}+\mathcal{L}_{0}(\theta)+\gamma_{n,\mathcal{L}_{0}}\bigg{)}\leqslant\delta\,.$

#### Choice of Legendre function $\mathcal{L}_{0}$

A natural choice for the Legendre regularizer $\mathcal{L}_{0}$ is to use the
log-partition function of the exponential family at hand. If
$\mathcal{L}_{0}=c\mathcal{L}$, where $c>0$ is a scaling coefficient, the
formula above simplify to

$\displaystyle\theta_{n,\mathcal{L}_{0}}\\!=\frac{\theta_{n}}{1+\frac{c}{n}}\,\text{
and
}\,\gamma_{n,\mathcal{L}_{0}}\\!=\log\left(\frac{\int_{\Theta}\exp\Big{(}-c\mathcal{L}(\theta^{\prime})\Big{)}d\theta^{\prime}}{\int_{\Theta}\exp\Big{(}-\\!(n+c)\mathcal{B}_{\mathcal{L}}\left(\theta^{\prime},\theta_{n,\mathcal{L}_{0}}\right)\Big{)}d\theta^{\prime}}\right)\,,$

where $\theta_{n}=\left(\nabla\mathcal{L}\right)^{-1}(\mu_{n})$ is the
standard maximum likelihood estimate of the parameter $\theta_{0}$.

In the special case of the one-dimensional Gaussian family
$\mathcal{N}(\theta_{0},\sigma^{2})$ where the variance $\sigma^{2}$ is known,
the log-partition function regularizer defined by
$\mathcal{L}_{0}(\theta)=\frac{\theta^{2}}{2\sigma^{2}}$ satisfies the
integrability assumption
$\int_{\mathbb{R}}e^{-c\mathcal{L}_{0}(\theta)}d\theta<\infty$ with $c>0$.
Straightforward calculations show that the resulting confidence set is the
same as the one derived from Theorem 3 with local regularization $c$ (see
Appendix C.1). For many other standard families however, this integrability
assumption may fail, as in the case of the exponential distributions
$\mathcal{E}(-\theta)$, for which $\mathcal{L}(\theta)=-\log(-\theta)$ (see
Appendix C.5).

For other choices of Legendre function, computing the regularized parameter
estimate $\theta_{n,\mathcal{L}_{0}}$ and the information gain
$\gamma_{n,\mathcal{L}_{0}}$ requires inverting the function
$\nabla\mathcal{L}+\frac{1}{n}\nabla\mathcal{L}_{0}$ and computing a tedious
integral, both of which seldom result in closed-form expressions. For these
reasons, we recommend as a first step the use of local regularization (Theorem
3) over Legendre regularizers. In Appendix F, we detail an application of
Theorem 7 to Gaussian contextual bandits using a quadratic regularizer.

##  Appendix B Properties of Bregman Divergences

In this section, we detail the derivation of some useful properties of the
Bregman divergence that are used in the proof of the Bregman deviation
inequalities.

See 1

Proof of Lemma 1:

The first equality is immediate, since
$\mathbb{E}_{\theta}[F(X)]=\nabla{\mathcal{L}}(\theta)$ and

$\displaystyle\log\mathbb{E}_{\theta}\left[\exp(\langle\lambda,F(X)\rangle)\right]$
$\displaystyle=$
$\displaystyle\log\int\exp(\langle\lambda,F(x)\rangle+\langle{\theta},F(x)\rangle-{\mathcal{L}}(\theta))dx$
$\displaystyle=$
$\displaystyle{\mathcal{L}}(\lambda+\theta)-{\mathcal{L}}(\theta)\,.$

We now turn to the duality formula. Using the definition of each terms, we
have

$\displaystyle\mathcal{B}_{{\mathcal{L}},\theta^{\prime}}^{\star}(\nabla{\mathcal{L}}(\theta)-\nabla{\mathcal{L}}(\theta^{\prime}))$
$\displaystyle=$
$\displaystyle\sup_{\lambda}\langle\lambda,\nabla{\mathcal{L}}(\theta)-\nabla{\mathcal{L}}(\theta^{\prime})\rangle-\Big{[}{\mathcal{L}}(\theta^{\prime}+\lambda)-{\mathcal{L}}(\theta^{\prime})-\langle\lambda,\nabla{\mathcal{L}}(\theta^{\prime})\rangle\Big{]}$
$\displaystyle=$
$\displaystyle\sup_{\lambda}\langle\lambda,\nabla{\mathcal{L}}(\theta)\rangle-{\mathcal{L}}(\theta^{\prime}+\lambda)+{\mathcal{L}}(\theta^{\prime})\,.$

An optimal $\lambda$ must satisfy
$\nabla{\mathcal{L}}(\theta)=\nabla{\mathcal{L}}(\theta^{\prime}+\lambda)$.
Hence, provided that $\nabla{\mathcal{L}}$ is invertible, this means
$\lambda=\theta-\theta^{\prime}$. Plugin-in this value, we obtain

$\displaystyle\mathcal{B}_{{\mathcal{L}},\theta^{\prime}}^{\star}(\nabla{\mathcal{L}}(\theta)-\nabla{\mathcal{L}}(\theta^{\prime}))=\langle\theta-\theta^{\prime},\nabla{\mathcal{L}}(\theta)\rangle-{\mathcal{L}}(\theta)+{\mathcal{L}}(\theta^{\prime})=\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\,.$

The remaining equality
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\mathcal{B}_{{\mathcal{L}}^{\star}}(\nabla{\mathcal{L}}(\theta),\nabla{\mathcal{L}}(\theta^{\prime}))$
is a standard result. Finally, regarding the generalization, we note that

$\displaystyle\mathcal{B}_{{\mathcal{L}},\theta^{\prime}}^{\star}\Big{(}\alpha(\nabla{\mathcal{L}}(\theta)-\nabla{\mathcal{L}}(\theta^{\prime}))\Big{)}$
$\displaystyle=$
$\displaystyle\sup_{\lambda}\langle\lambda,\alpha\nabla{\mathcal{L}}(\theta)+(1-\alpha)\nabla\mathcal{L}(\theta^{\prime})\rangle-{\mathcal{L}}(\theta^{\prime}+\lambda)+{\mathcal{L}}(\theta^{\prime})\,.$

Hence, an optimal $\lambda$ must now satisfy
$\alpha\nabla{\mathcal{L}}(\theta)+(1-\alpha)\nabla\mathcal{L}(\theta^{\prime})=\nabla{\mathcal{L}}(\theta^{\prime}+\lambda)$,
that is $\lambda=\theta_{\alpha}-\theta^{\prime}$. This further yields

$\displaystyle\mathcal{B}_{{\mathcal{L}},\theta^{\prime}}^{\star}\Big{(}\alpha(\nabla{\mathcal{L}}(\theta)-\nabla{\mathcal{L}}(\theta^{\prime}))\Big{)}$
$\displaystyle=$
$\displaystyle\langle\theta_{\alpha}-\theta^{\prime},\nabla{\mathcal{L}}(\theta_{\alpha})\rangle-{\mathcal{L}}(\theta_{\alpha})+{\mathcal{L}}(\theta^{\prime})=\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{\alpha})\,.$

$\hfill\square$

##  Appendix C Specification to Illustrative Exponential Families

In this section, we provide the technical derivations to specify our generic
concentration result to a few classical distributions. The results are
summarized in Table 2.

###  C.1 Gaussian with unknown mean, known variance

Consider $X\sim\mathcal{N}(\mu,\sigma^{2})$, where $\mu$ is unknown and
$\sigma$ is known. This corresponds to an exponential family of distributions,
with parameter $\theta=\mu$, feature function $F(x)=\frac{1}{\sigma^{2}}x$ and
log-partition function $\mathcal{L}(\theta)=\frac{\theta^{2}}{2\sigma^{2}}$.
The Bregman divergence between two parameters $\theta^{\prime}$ and $\theta$
associated with $\mathcal{L}$ is given by
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\texttt{KL}(P_{\mu},P_{\mu^{\prime}})=\frac{(\mu-\mu^{\prime})^{2}}{2\sigma^{2}}$.
We further have
$\mathcal{L}^{\prime}(\theta)=\frac{\theta}{\sigma^{2}}=\frac{\mu}{\sigma^{2}}$
and $\mathcal{L}^{\prime\prime}(\theta)=\frac{1}{\sigma^{2}}$, implying that
$\mathcal{L}^{\prime}$ is invertible. Denoting $S_{n}=\sum_{t=1}^{n}X_{t}$, we
obtain that

$\displaystyle\mu_{n,c}(\mu)=\theta_{n,c}(\theta)=\frac{S_{n}+c\mu}{n+c}~{}.$

The Bregman deviations simplify as follows

$\displaystyle(n+c)\cdot\mathcal{B}_{\mathcal{L}}(\theta,\theta_{n,c}(\theta))$
$\displaystyle=(n+c)\frac{(\theta_{n,c}(\theta)-\theta)^{2}}{2\sigma^{2}}=\frac{1}{n+c}\frac{\left(S_{n}-n\mu\right)^{2}}{2\sigma^{2}}~{}.$

Now, on the other hand, let us see that the information gain is explicitly
given by

$\displaystyle\gamma_{n,c}(\mu)$
$\displaystyle=\frac{1}{2}\log\frac{n+c}{c}~{}.$

We obtain from Theorem 3 w.p. $\geqslant 1-\delta$, $\forall n\in\mathbb{N}$,

$\displaystyle\frac{1}{n+c}\frac{(S_{n}-n\mu)^{2}}{2\sigma^{2}}\leqslant\log(1/\delta)+\frac{1}{2}\log\frac{n+c}{c}~{}.$

###  C.2 Gaussian with known mean, unknown variance

We consider $X\sim\mathcal{N}\left(\mu,\sigma^{2}\right)$, where $\mu$ is
known and $\sigma$ is unknown. This corresponds to a one-dimensional
exponential family model with parameter
$\theta=-\frac{1}{2\sigma^{2}}\in\mathbb{R}^{*}_{-}$, feature function
$F(x)=\left(x-\mu\right)^{2}$ and log-partition function
$\mathcal{L}(\theta)=-\frac{1}{2}\log\left(-2\theta\right)$. The first two
derivatives of $\mathcal{L}$ are given by
$\mathcal{L}^{\prime}(\theta)=-\frac{1}{2\theta}$ and
$\mathcal{L}^{\prime\prime}(\theta)=\frac{1}{2\theta^{2}}$, which shows that
$\mathcal{L}^{\prime}$ is invertible on the domain
$\Theta=\mathbb{R}^{*}_{-}$. The Bregman divergence between two parameters
$\theta^{\prime}$ and $\theta$ is therefore
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\frac{1}{2}\log\left(\theta/\theta^{\prime}\right)+\frac{1}{2}\left(\theta^{\prime}/\theta-1\right)$.

Let $Q_{n}=\sum_{t=1}^{n}\left(X_{t}-\mu\right)^{2}$. A short calculation
shows that the following expression holds:

$\displaystyle\theta_{n,c}(\theta)$
$\displaystyle=-\frac{n+c}{c/\theta-2Q_{n}}\,.$

To compute the Bregman information gain, note that the above expression and a
change of variable implies

$\displaystyle\int_{-\infty}^{0}\exp\left(-c\mathcal{B}_{\mathcal{L}}\left(\theta^{\prime},\theta\right)\right)d\theta^{\prime}$
$\displaystyle=-\theta\left(2/c\right)^{c/2+1}e^{c/2}\Gamma\left(c/2+2\right)\,.$

Combining the last two lines yields

$\displaystyle\gamma_{n,c}(\theta)$
$\displaystyle=\log\frac{\theta}{\theta_{n,c}(\theta)}-(1\\!+\\!\log
2)\frac{n}{2}-(\frac{c}{2}\\!+\\!1)\log
c+(\frac{n\\!+\\!c}{2}\\!+\\!1)\log(n\\!+\\!c)-\frac{\log\Gamma\left(\frac{n\\!+\\!c}{2}\\!+\\!2\right)}{\log\Gamma\left(\frac{c}{2}\\!+\\!2\right)}$
$\displaystyle=\log\frac{c\\!-\\!2\theta Q_{n}}{n\\!+\\!c}-(1\\!+\\!\log
2)\frac{n}{2}-(\frac{c}{2}\\!+\\!1)\log
c+(\frac{n\\!+\\!c}{2}\\!+\\!1)\log(n\\!+\\!c)-\frac{\log\Gamma\left(\frac{n\\!+\\!c}{2}\\!+\\!2\right)}{\log\Gamma\left(\frac{c}{2}\\!+\\!2\right)}$
$\displaystyle=\log\frac{Q_{n}/\sigma^{2}+c}{n\\!+\\!c}-(1\\!+\\!\log
2)\frac{n}{2}-(\frac{c}{2}\\!+\\!1)\log
c+(\frac{n\\!+\\!c}{2}\\!+\\!1)\log(n\\!+\\!c)-\frac{\log\Gamma\left(\frac{n\\!+\\!c}{2}\\!+\\!2\right)}{\log\Gamma\left(\frac{c}{2}\\!+\\!2\right)}\,.$

Moreover, we deduce from the expression of $\theta_{n,c}(\theta)$ and the
Bregman divergence that:

$\displaystyle\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta))$
$\displaystyle=-\frac{1}{2}\log\left(\frac{c\theta^{\prime}/\theta-2\theta^{\prime}Q_{n}}{n\\!+\\!c}\right)+\frac{1}{2}\left(\frac{c\theta^{\prime}/\theta-2\theta^{\prime}Q_{n}}{n\\!+\\!c}-1\right)\,.$

So, we have

$\displaystyle(n+c)\mathcal{B}_{\mathcal{L}}(\theta,\theta_{n,c}(\theta))$
$\displaystyle=-\frac{n+c}{2}\log\left(\frac{Q_{n}/\sigma^{2}+c}{n\\!+\\!c}\right)+\frac{Q_{n}/\sigma^{2}+c}{2}-\frac{n}{2}\,.$

After some simple algebra and using the natural parametrization
$\theta=-\frac{1}{2\sigma^{2}}$, we derive from Theorem 3 that w.p. $\geqslant
1-\delta$, $\forall n\in\mathbb{N}$,

$\displaystyle-(\frac{n\\!+\\!c}{2}\\!+\\!1)\log\frac{Q_{n}/\sigma^{2}+c}{n\\!+\\!c}+\frac{Q_{n}}{2\sigma^{2}}$
$\displaystyle\leqslant\log\frac{1}{\delta}-\frac{n}{2}\log
2-(\frac{c}{2}\\!+\\!1)\log
c+(\frac{n\\!+\\!c}{2}\\!+\\!1)\log(n\\!+\\!c)+\log\Gamma\left(\frac{c}{2}\\!+\\!2\right)-\log\Gamma\left(\frac{n\\!+\\!c}{2}\\!+\\!2\right)\,.$
(8)

###  C.3 Gaussian with unknown mean and variance

We consider $X\sim\mathcal{N}\left(\mu,\sigma^{2}\right)$, where both
$\mu\in\mathbb{R}$ and $\sigma>0$ are unknown. These distributions form a two-
dimensional exponential family with parameter
$\theta=\begin{pmatrix}\theta_{1}\\\
\theta_{2}\end{pmatrix}=\begin{pmatrix}\frac{\mu}{\sigma^{2}}\\\
-\frac{1}{2\sigma^{2}}\end{pmatrix}$ belonging to the domain
$\Theta=\mathbb{R}\times\mathbb{R}^{*}_{-}$. The corresponding feature and
log-partition functions writes $F(x)=\begin{pmatrix}x\\\ x^{2}\end{pmatrix}$
and
$\mathcal{L}(\theta)=-\frac{1}{2}\log(-\theta_{2})-\frac{\theta_{1}^{2}}{4\theta_{2}}$
respectively. The gradient and Hessian of $\mathcal{L}$ follows from
straightforward calculations and write:

$\displaystyle\nabla\mathcal{L}(\theta)=\begin{pmatrix}-\frac{\theta_{1}}{2\theta_{2}}\\\
-\frac{1}{2\theta_{2}}+\frac{\theta_{1}^{2}}{4\theta_{2}^{2}}\end{pmatrix}\,\,,\nabla^{2}\mathcal{L}(\theta)=\begin{pmatrix}-\frac{1}{2\theta_{2}}&\frac{\theta_{1}}{2\theta_{2}^{2}}\\\
\frac{\theta_{1}}{2\theta_{2}^{2}}&\frac{1}{2\theta_{2}^{2}}-\frac{\theta_{1}^{2}}{2\theta_{2}^{3}}\end{pmatrix}\,.$

In particular, $\text{Tr
}\nabla^{2}\mathcal{L}(\theta)=-\frac{1}{2\theta_{2}}+\frac{1}{2\theta_{2}^{2}}-\frac{\theta_{1}^{2}}{2\theta_{2}^{3}}>0$
and $\det\nabla^{2}\mathcal{L}(\theta)=-\frac{1}{4\theta_{2}^{3}}>0$,
therefore the symmetric matrix $\nabla^{2}\mathcal{L}(\theta)$ is positive
definite (both its eigenvalues are positive). The inverse of the gradient
mapping $\nabla\mathcal{L}\colon\Theta\rightarrow\mathcal{P}$, where the image
set is $\mathcal{P}=\left\\{(x,y)\in\mathbb{R}^{2},y>x^{2}\right\\}$, can be
explicitly computed as:

$\displaystyle\left(\nabla\right)^{-1}\mathcal{L}\begin{pmatrix}x\\\
y\end{pmatrix}=\frac{1}{y-x^{2}}\begin{pmatrix}x\\\
-\frac{1}{2}\end{pmatrix}\,.$

The expression of the Bregman divergence between two parameters $\theta$ and
$\theta^{\prime}$ can be calculated from the above results and reads:

$\displaystyle\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)$
$\displaystyle=\frac{1}{2}\log\frac{\theta_{2}}{\theta^{\prime}_{2}}+\frac{\theta^{\prime}_{2}}{2\theta_{2}}-\theta^{\prime}_{2}\left(\frac{\theta^{\prime}_{1}}{2\theta^{\prime}_{2}}-\frac{\theta_{1}}{2\theta_{2}}\right)^{2}-\frac{1}{2}.$

Let $S_{n}=\sum_{t=1}^{n}X_{t}$ and $Q_{n}=\sum_{t=1}^{n}X_{t}^{2}$ and
$V_{n,c}=\frac{1}{n+c}Q_{n}-\left(\frac{1}{n+c}S_{n}\right)^{2}$. We have that

$\displaystyle\theta_{n,c}(\theta)$
$\displaystyle=\frac{1}{V_{n,c}-\frac{c}{n\\!+\\!c}\frac{1}{2\theta_{2}}+\frac{nc}{(n\\!+\\!c)^{2}}\frac{\theta_{1}^{2}}{4\theta_{2}^{2}}+\frac{c}{(n\\!+\\!c)^{2}}\frac{\theta_{1}}{\theta_{2}}S_{n}}\begin{pmatrix}\frac{1}{n\\!+\\!c}S_{n}-\frac{c}{2(n\\!+\\!c)}\frac{\theta_{1}}{\theta_{2}}\\\
-\frac{1}{2}\end{pmatrix}\,.$

To compute the Bregman information gain, we will need to evaluate the integral
$\int_{\Theta}\exp\left(-c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\right)d\theta^{\prime}$.
As the integrand is nonnegative, we can integrate first along
$\theta_{1}^{\prime}$ (Fubini’s theorem), which writes:

$\displaystyle\int_{\Theta}\exp\left(-c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\right)d\theta^{\prime}$
$\displaystyle=\int_{-\infty}^{0}\int_{-\infty}^{\infty}\exp\left(-\frac{c}{2}\log\left(\frac{\theta_{2}}{\theta^{\prime}_{2}}\right)-\frac{c}{2}\frac{\theta^{\prime}_{2}}{\theta_{2}}+c\theta^{\prime}_{2}\left(\frac{\theta^{\prime}_{1}}{2\theta^{\prime}_{2}}-\frac{\theta_{1}}{2\theta_{2}}\right)^{2}+\frac{c}{2}\right)d\theta_{1}^{\prime}d\theta_{2}^{\prime}$
$\displaystyle=e^{\frac{c}{2}}\int_{-\infty}^{0}\left(\frac{\theta^{\prime}_{2}}{\theta_{2}}\right)^{\frac{c}{2}}e^{-\frac{c}{2}\frac{\theta^{\prime}_{2}}{\theta_{2}}}\int_{-\infty}^{\infty}\exp\left(c\theta^{\prime}_{2}\left(\frac{\theta^{\prime}_{1}}{2\theta^{\prime}_{2}}-\frac{\theta_{1}}{2\theta_{2}}\right)^{2}\right)d\theta_{1}^{\prime}d\theta_{2}^{\prime}\,.$

Let $m=\frac{\theta_{1}}{2\theta_{2}}$ and
$s=\frac{1}{\sqrt{-2c\theta^{\prime}_{2}}}$. The integral w.r.t
$\theta^{\prime}_{1}$ can be rewritten using the change of variable
$y=\frac{\theta^{\prime}_{1}}{2\theta^{\prime}_{2}}$ as:

$\displaystyle\int_{-\infty}^{\infty}\exp\left(c\theta^{\prime}_{2}\left(\frac{\theta^{\prime}_{1}}{2\theta^{\prime}_{2}}-\frac{\theta_{1}}{2\theta_{2}}\right)^{2}\right)d\theta_{1}^{\prime}$
$\displaystyle=-2\theta_{2}^{\prime}\sqrt{2\pi}s\underbrace{\int_{-\infty}^{\infty}\frac{1}{\sqrt{2\pi}s}e^{-\frac{\left(y-m\right)^{2}}{2s^{2}}}dy}_{=1}\,.$

Therefore, the above calculation simplifies to:

$\displaystyle\int_{\Theta}\exp\left(-c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\right)d\theta^{\prime}$
$\displaystyle=2\sqrt{\frac{\pi}{c}}e^{\frac{c}{2}}\int_{-\infty}^{0}\left(\frac{\theta^{\prime}_{2}}{\theta_{2}}\right)^{\frac{c}{2}}e^{-\frac{c}{2}\frac{\theta^{\prime}_{2}}{\theta_{2}}}\sqrt{-\theta_{2}^{\prime}}d\theta_{2}^{\prime}\,,$

which after a linear change of variable on $\theta^{\prime}_{2}$ can be
related to the Gamma function as follows:

$\displaystyle\int_{\Theta}\exp\left(-c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\right)d\theta^{\prime}$
$\displaystyle=2\sqrt{\frac{\pi}{c}}\left(\frac{2}{c}\right)^{\frac{c+1}{2}+1}\Gamma\left(\frac{c+1}{2}+1\right)\left(-\theta_{2}\right)^{\frac{3}{2}}\,.$

The Bregman information gain is thus:

$\displaystyle\gamma_{n,c}(\theta)$ $\displaystyle=-\left(\frac{1\\!+\\!\log
2}{2}\right)n-\left(\frac{c}{2}+2\right)\log
c+\left(\frac{n\\!+\\!c}{2}+2\right)\log(n\\!+\\!c)$
$\displaystyle\quad+\log\Gamma\left(\frac{c\\!+\\!3}{2}\right)-\log\Gamma\left(\frac{n\\!+\\!c\\!+\\!3}{2}\right)+\frac{3}{2}\log\frac{\theta_{2}}{\theta_{n,c}(\theta)_{2}}\,,$

with
$\frac{\theta_{2}}{\theta_{n,c}(\theta)_{2}}=-2\theta_{2}V_{n,c}+\frac{c}{n\\!+\\!c}-\frac{nc}{(n\\!+\\!c)^{2}}\frac{\theta_{1}^{2}}{2\theta_{2}}-\frac{2c\theta_{1}}{(n\\!+\\!c)^{2}}S_{n}$.
Now, applying the result of Theorem 3 shows that w.p. $\geqslant 1-\delta$,
$\forall n\in\mathbb{N}$,

$\displaystyle-\frac{n\\!+\\!c\\!+\\!3}{2}\log\left(\frac{\theta_{2}}{\theta_{n,c}(\theta)_{2}}\right)-\theta_{2}(n\\!+\\!c)V_{n,c}-\frac{n\theta_{1}^{2}}{4\theta_{2}}-\frac{\theta_{2}}{n\\!+\\!c}S_{n}^{2}-\theta_{1}S_{n}$
$\displaystyle\quad\leqslant\log\frac{1}{\delta}-\frac{n}{2}\log
2-\left(\frac{c}{2}\\!+\\!2\right)\log
c+\left(\frac{n\\!+\\!c}{2}\\!+\\!2\right)\log\left(n\\!+\\!c\right)+\log\Gamma\left(\frac{c\\!+\\!3}{2}\right)-\log\Gamma\left(\frac{n\\!+\\!c\\!+\\!3}{2}\right)\,.$

After expanding the ratio $\frac{\theta_{2}}{\theta_{n,c}(\theta)_{2}}$ and
substituting the natural parametrization in terms of $\mu$ and $\sigma$, we
finally obtain that w.p. $\geqslant 1-\delta$, $\forall n\in\mathbb{N}$,

$\displaystyle-\frac{n\\!+\\!c\\!+\\!3}{2}\log\left(\frac{1}{\sigma^{2}}V_{n,c}+\frac{c}{n\\!+\\!c}+\frac{nc}{(n\\!+\\!c)^{2}}\frac{\mu^{2}}{\sigma^{2}}-\frac{2c}{(n\\!+\\!c)^{2}}\frac{\mu}{\sigma^{2}}S_{n}\right)+\frac{n\\!+\\!c}{2\sigma^{2}}V_{n,c}+\frac{n\mu^{2}}{2\sigma^{2}}+\frac{1}{2(n\\!+\\!c)\sigma^{2}}S_{n}^{2}-\frac{\mu}{\sigma^{2}}S_{n}$
$\displaystyle\qquad\leqslant\log\frac{1}{\delta}-\frac{n}{2}\log
2-\left(\frac{c}{2}\\!+\\!2\right)\log
c+\left(\frac{n\\!+\\!c}{2}\\!+\\!2\right)\log\left(n\\!+\\!c\right)+\log\Gamma\left(\frac{c\\!+\\!3}{2}\right)-\log\Gamma\left(\frac{n\\!+\\!c\\!+\\!3}{2}\right)\,.$

To simplify this formula, we introduce the standardized sum of squares
$Z_{n}(m,s)=\sum\limits_{t=1}^{n}\left(\frac{X_{t}-m}{s}\right)^{2}$ for
$(m,s)\in\mathbb{R}\times\mathbb{R}_{+}^{*}$. After rearranging terms and
denoting $\widehat{\mu}_{n}=S_{n}/n$, the above formula reads:

$\displaystyle\frac{1}{2}Z_{n}(\mu,\sigma)-\frac{n+c+3}{2}\log\left(\frac{n}{n\\!+\\!c}Z_{n}(\widehat{\mu}_{n},\sigma)+\frac{c}{n\\!+\\!c}Z_{n}(\mu,\sigma)+c\right)$
$\displaystyle\leqslant\log\frac{1}{\delta}-\frac{n}{2}\log
2-\left(\frac{c}{2}\\!+\\!2\right)\log
c+\frac{1}{2}\log\left(n\\!+\\!c\right)+\log\Gamma\left(\frac{c\\!+\\!3}{2}\right)-\log\Gamma\left(\frac{n\\!+\\!c\\!+\\!3}{2}\right)\,.$

###  C.4 Bernoulli

We consider $X\sim\text{Bernoulli}(\mu)$. This corresponds to an exponential
family model, with parameter $\theta=\log\frac{\mu}{1-\mu}$, feature function
$F(x)=x$ and log-partition function
$\mathcal{L}(\theta)=\log(1+\exp(\theta))$. The Bregman divergence between two
parameters $\theta^{\prime}$ and $\theta$ associated with $\mathcal{L}$ is
given by
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\texttt{kl}(\mu,\mu^{\prime})=\mu\log\frac{\mu}{\mu^{\prime}}+(1-\mu)\log\frac{1-\mu}{1-\mu^{\prime}}$.
We further have
$\mathcal{L}^{\prime}(\theta)=\frac{\exp(\theta)}{1+\exp(\theta)}=\mu$ and
$\mathcal{L}^{\prime\prime}(\theta)=\frac{\exp(\theta)}{(1+\exp(\theta))^{2}}=\mu(1-\mu)$.
Therefore $\mathcal{L}^{\prime}$ is invertible and we have the expression
$\theta_{n,c}(\theta)=(\mathcal{L}^{\prime})^{-1}\left(\frac{1}{n+c}\sum_{t=1}^{n}X_{t}+\frac{c}{n+c}\frac{\exp(\theta)}{1+\exp(\theta)}\right)$.
Then, denoting $S_{n}=\sum_{t=1}^{n}X_{t}$, we get

$\displaystyle\mu_{n,c}(\mu):=\mathcal{L}^{\prime}(\theta_{n,c}(\theta))=\frac{\sum_{t=1}^{n}X_{t}+c\mu}{n+c}=\frac{S_{n}+c\mu}{n+c}~{}.$

Therefore, the Bregman deviation specifies to the following closed-form
formula

$\displaystyle(n+c)\cdot\mathcal{B}_{\mathcal{L}}(\theta,\theta_{n,c}(\theta))$
$\displaystyle=(n+c)\cdot\texttt{kl}(\mu_{n,c}(\mu),\mu)$
$\displaystyle=(S_{n}+c\mu)\log\frac{\mu_{n,c}(\mu)}{\mu}+\left(n-S_{n}+c(1-\mu)\right)\log\frac{1-\mu_{n,c}(\mu)}{1-\mu}~{}.$

Now, we observe that

$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}-\\!c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\big{)}d\theta^{\prime}$
$\displaystyle=\frac{\text{B}(c\mu,c(1-\mu))}{\mu^{c\mu}(1-\mu)^{c(1-\mu)}}~{},$
$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}-\\!(n+c)\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta))\big{)}d\theta^{\prime}$
$\displaystyle=\frac{\text{B}((n+c)\mu_{n,c}(\mu),(n+c)(1-\mu_{n,c}(\mu)))}{\mu_{n,c}(\mu)^{(n+c)\mu_{n,c}(\mu)}(1-\mu_{n,c}(\mu))^{(n+c)(1-\mu_{n,c}(\mu))}}~{}.$

Therefore, we deduce that the Bregman information gain rewrites

$\displaystyle\gamma_{n,c}(\mu)=$
$\displaystyle(S_{n}+c\mu)\log\mu_{n,c}(\mu)+\left(n-S_{n}+c(1-\mu)\right)\log(1-\mu_{n,c}(\mu))-c\mu\log\mu-c(1-\mu)\log(1-\mu)$
$\displaystyle+\log\frac{\Gamma(c\mu)\Gamma(c(1-\mu))}{\Gamma(S_{n}+c\mu)\Gamma(n-S_{n}+c(1-\mu))}+\log\frac{\Gamma(n+c)}{\Gamma(c)}~{}.$

Combining the above and using Theorem 3, we obtain that w.p. $\geqslant
1-\delta$, $\forall n\in\mathbb{N}$,

$\displaystyle
S_{n}\log\frac{1}{\mu}+(n-S_{n})\log\frac{1}{1-\mu}+\log\frac{\Gamma(S_{n}+c\mu)\Gamma(n-S_{n}+c(1-\mu))}{\Gamma(c\mu)\Gamma(c(1-\mu))}\leqslant\log\frac{1}{\delta}+\log\frac{\Gamma(n+c)}{\Gamma(c)}~{}.$

###  C.5 Exponential

We consider $X\sim\text{Exp}(1/\mu)$ with unknown mean $\mu$. The distribution
of $X$ is supported on $[0,+\infty)$ with density
$p_{\mu}(x)=\frac{1}{\mu}e^{-\frac{x}{\mu}}$. This corresponds to an
exponential family model with parameter $\theta=-\frac{1}{\mu}$, feature
function $F(x)=x$ and log-partition function
$\mathcal{L}(\theta)=\log(-\frac{1}{\theta})$. The Bregman divergence between
two parameters $\theta^{\prime}$ and $\theta$ associated with $\mathcal{L}$ is
given by
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\texttt{KL}(P_{\mu},P_{\mu^{\prime}})=\frac{\mu}{\mu^{\prime}}-1-\log\frac{\mu}{\mu^{\prime}}$.

We have $\mathcal{L}^{\prime}(\theta)=-\frac{1}{\theta}=\mu$ and
$\mathcal{L}^{\prime\prime}(\theta)=\frac{1}{\theta^{2}}=\mu^{2}$. Therefore
$\mathcal{L}^{\prime}$ is invertible and we have

$\displaystyle\mu_{n,c}(\mu):=\mathcal{L}^{\prime}(\theta_{n,c}(\theta))=\frac{S_{n}+c\mu}{n+c}~{}.$

Therefore, we deduce that the Bregman divergence takes the following form

$\displaystyle(n+c)\cdot\mathcal{B}_{\mathcal{L}}(\theta,\theta_{n,c}(\theta))$
$\displaystyle=(n+c)\cdot\texttt{KL}(P_{\mu_{n,c}(\mu)},P_{\mu})=\frac{S_{n}}{\mu}-(n+c)\log\frac{\mu_{n,c}(\mu)}{\mu}-n~{}.$

Now, we observe that

$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}-\\!c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\big{)}d\theta^{\prime}$
$\displaystyle=\frac{\Gamma(c+1)}{c\mu}\left(\frac{e}{c}\right)^{c}=\frac{\Gamma(c)}{\mu}\left(\frac{e}{c}\right)^{c},$
$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}-\\!(n+c)\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta))\big{)}d\theta^{\prime}$
$\displaystyle=\frac{\Gamma(n+c+1)}{(n+c)\mu_{n,c}(\mu)}\left(\frac{e}{n+c}\right)^{n+c}=\frac{\Gamma(n+c)}{\mu_{n,c}(\mu)}\left(\frac{e}{n+c}\right)^{n+c}~{}.$

Therefore, the Bregman information gain writes explicitly as follows

$\displaystyle\gamma_{n,c}(\theta)=\log\frac{\mu_{n,c}(\mu)}{\mu}+\log\left(\Gamma(c)\left(\frac{e}{c}\right)^{c}\right)-\log\left(\Gamma(n+c)\left(\frac{e}{n+c}\right)^{n+c}\right).$

We can now specify the inequality
$(n+c)\cdot\mathcal{B}_{\mathcal{L}}(\theta,\theta_{n,c}(\theta))\leqslant\log(1/\delta)+\gamma_{n,c}(\theta)$.
Combining the above, we obtain w.p. $\geqslant 1-\delta$, $\forall
n\in\mathbb{N}$,

$\displaystyle\frac{S_{n}}{\mu}-\left(n+c+1\right)\log\left(\frac{S_{n}+c\mu}{(n+c)\mu}\right)\leqslant(n+c)\log(n+c)+\log\frac{\Gamma(c)}{\Gamma(n+c)}+\log\frac{1}{\delta}-c\log
c~{}.$

###  C.6 Gamma with fixed shape

We consider $X\sim\text{Gamma}(\lambda,k)$ with fixed shape $k>0$ and unknown
scale $\lambda>0$.222Note that $\text{Gamma}(\lambda,k)$ with $k=1$ is
$\text{Exp}(1/\lambda)$. The distribution of $X$ is supported on $[0,\infty)$
with density
$p_{\lambda}(x)=\frac{1}{\Gamma(k)\lambda^{k}}x^{k-1}e^{-\frac{x}{\lambda}}$.
This corresponds to an exponential family model with parameter
$\theta=-\frac{1}{\lambda}$, feature function $F(x)=x$ and log-partition
function $\mathcal{L}(\theta)=\log(\lambda^{k})=k\log(-\frac{1}{\theta})$. The
Bregman divergence between two parameters $\theta^{\prime}$ and $\theta$
associated with $\mathcal{L}$ is given by
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\texttt{KL}(p_{\lambda},p_{\lambda^{\prime}})=k\left(\frac{\lambda}{\lambda^{\prime}}-1-\log\left(\frac{\lambda}{\lambda^{\prime}}\right)\right)$.

Note that $\mathcal{L}^{\prime}(\theta)=-\frac{k}{\theta}=k\lambda$ and
$\mathcal{L}^{\prime\prime}(\theta)=\frac{k}{\theta^{2}}=k\lambda$. Therefore
$\mathcal{L}^{\prime}$ is invertible, and we get

$\displaystyle
k\lambda_{n,c}(\lambda):=\mathcal{L}^{\prime}(\theta_{n,c}(\theta))=\frac{\sum_{t=1}^{n}X_{t}+ck\lambda}{n+c}$

Therefore, the Bregman divergence takes the following form

$\displaystyle(n+c)\cdot\mathcal{B}_{\mathcal{L}}(\theta,\theta_{n,c}(\theta))$
$\displaystyle=(n+c)\cdot\texttt{KL}(p_{\lambda_{n,c}(\lambda)},p_{\lambda})$
$\displaystyle=k(n+c)\left(\frac{\sum_{t=1}^{n}X_{t}+ck\lambda}{(n+c)k\lambda}\right)-k(n+c)\log\left(\frac{\sum_{t=1}^{n}X_{t}+ck\lambda}{(n+c)k\lambda}\right)-k(n+c)~{}.$

Now, we observe that

$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}-\\!c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\big{)}d\theta^{\prime}$
$\displaystyle=\frac{\Gamma(ck)}{\lambda}\left(\frac{e}{ck}\right)^{ck},$
$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}-\\!(n+c)\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta))\big{)}d\theta^{\prime}$
$\displaystyle=\frac{\Gamma((n+c)k)}{\lambda_{n,c}(\lambda)}\left(\frac{e}{(n+c)k}\right)^{(n+c)k}~{}.$

Therefore, the Bregman information gain writes explicitly as follows

$\displaystyle\gamma_{n,c}(\theta)=\log\left(\frac{\sum_{t=1}^{n}X_{t}+ck\lambda}{(n+c)k\lambda}\right)+\log\left(\Gamma(ck)\left(\frac{e}{ck}\right)^{ck}\right)-\log\left(\Gamma((n+c)k)\left(\frac{e}{(n+c)k}\right)^{(n+c)k}\right).$

Using the above with Theorem 3, we obtain w.p. $\geqslant 1-\delta$, $\forall
n\in\mathbb{N}$,

$\displaystyle
k(n+c)\left(\frac{\sum_{t=1}^{n}X_{t}+ck\lambda}{(n+c)k\lambda}\right)-(k(n+c)+1)\log\left(\frac{\sum_{t=1}^{n}X_{t}+ck\lambda}{(n+c)k\lambda}\right)$
$\displaystyle\leqslant\log\frac{1}{\delta}+\log\frac{\Gamma(ck)}{\Gamma((n+c)k)}+(n+c)k\log((n+c)k)+ck-
ck\log ck.$ (9)

###  C.7 Weibull with fixed shape

We consider $X\sim\text{Weibull}(\lambda,k)$ with fixed shape $k>0$ and
unknown scale $\lambda>0$.333Note that $\text{Weibull}(\lambda,k)$ with $k=1$
is $\text{Exp}(1/\lambda)$. The distribution of $X$ is supported on
$[0,\infty)$ with density
$p_{\lambda}(x)=\frac{k}{\lambda}\left(\frac{x}{\lambda}\right)^{k-1}e^{-\left(\frac{x}{\lambda}\right)^{k}}$.
This corresponds to an exponential family model with parameter
$\theta=-\frac{1}{\lambda^{k}}$, feature function $F(x)=x^{k}$ and log-
partition function
$\mathcal{L}(\theta)=\log(\lambda^{k})=\log(-\frac{1}{\theta})$. The Bregman
divergence between two parameters $\theta^{\prime}$ and $\theta$ associated
with $\mathcal{L}$ is given by
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\texttt{KL}(P_{\lambda},P_{\lambda^{\prime}})=\left(\frac{\lambda}{\lambda^{\prime}}\right)^{k}-1-\log\left(\frac{\lambda}{\lambda^{\prime}}\right)^{k}$.

Note that $\mathcal{L}^{\prime}(\theta)=-\frac{1}{\theta}=\lambda^{k}$ and
$\mathcal{L}^{\prime\prime}(\theta)=\frac{1}{\theta^{2}}=\lambda^{2k}$.
Therefore $\mathcal{L}^{\prime}$ is invertible, and we get

$\displaystyle(\lambda_{n,c}(\lambda))^{k}:=\mathcal{L}^{\prime}(\theta_{n,c}(\theta))=\frac{\sum_{t=1}^{n}X_{t}^{k}+c\lambda^{k}}{n+c}.$

Therefore, the Bregman divergence takes the following form

$\displaystyle(n+c)\cdot\mathcal{B}_{\mathcal{L}}(\theta,\theta_{n,c}(\theta))$
$\displaystyle=(n+c)\cdot\texttt{KL}(p_{\lambda_{n,c}(\lambda)},p_{\lambda})$
$\displaystyle=(n+c)\left(\frac{\sum_{t=1}^{n}X_{t}^{k}+c\lambda^{k}}{(n+c)\lambda^{k}}\right)-(n+c)\log\left(\frac{\sum_{t=1}^{n}X_{t}^{k}+c\lambda^{k}}{(n+c)\lambda^{k}}\right)-(n+c)~{}.$

Now, we observe that

$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}-\\!c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\big{)}d\theta^{\prime}$
$\displaystyle=\frac{\Gamma(c)}{\lambda^{k}}\left(\frac{e}{c}\right)^{c},$
$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}-\\!(n+c)\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta))\big{)}d\theta^{\prime}$
$\displaystyle=\frac{\Gamma(n+c)}{(\lambda_{n,c}(\lambda))^{k}}\left(\frac{e}{n+c}\right)^{n+c}~{}.$

Therefore, the Bregman information gain writes explicitly as follows

$\displaystyle\gamma_{n,c}(\theta)=\log\left(\frac{\sum_{t=1}^{n}X_{t}^{k}+c\lambda^{k}}{(n+c)\lambda^{k}}\right)+\log\left(\Gamma(c)\left(\frac{e}{c}\right)^{c}\right)-\log\left(\Gamma(n+c)\left(\frac{e}{n+c}\right)^{n+c}\right).$

Using the above with Theorem 3, we obtain w.p. $\geqslant 1-\delta$, $\forall
n\in\mathbb{N}$,

$\displaystyle(n+c)\left(\frac{\sum_{t=1}^{n}X_{t}^{k}+c\lambda^{k}}{(n+c)\lambda^{k}}\right)-(n+c+1)\log\left(\frac{\sum_{t=1}^{n}X_{t}^{k}+c\lambda^{k}}{(n+c)\lambda^{k}}\right)$
$\displaystyle\leqslant\log\frac{1}{\delta}+\log\frac{\Gamma(c)}{\Gamma(n+c)}+(n+c)\log(n+c)+c-c\log
c.$ (10)

###  C.8 Pareto with fixed scale

We consider $X\sim\text{Pareto}\left(\alpha\right)$, where $\alpha>0$ is
unknown.444We assume scale is fixed to the value $1$. The distribution of $X$
is supported in $[1,+\infty)$, with density
$p_{\alpha}(x)=\frac{\alpha}{x^{\alpha\\!+\\!1}}$, which corresponds to a one-
dimensional exponential family model with parameter
$\theta=-\alpha-1\in(-\infty,-1)$, feature function $F(x)=\log x$ and log-
partition function $\mathcal{L}(\theta)=-\log\left(-1\\!-\\!\theta\right)$.
The first two derivatives of $\mathcal{L}$ are given by
$\mathcal{L}^{\prime}(\theta)=-\frac{1}{1+\theta}$ and
$\mathcal{L}^{\prime\prime}(\theta)=\frac{1}{(1+\theta)^{2}}$, therefore
$\mathcal{L}^{\prime}$ is invertible on the domain $(-\infty,-1)$. Using these
expressions, the Bregman divergence between two parameters $\theta^{\prime}$
and $\theta$ writes
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=-\log\left(-1-\theta^{\prime}\right)+\log\left(-1-\theta\right)+\frac{\theta^{\prime}\\!-\\!\theta}{1\\!+\\!\theta}$.
Using the shorthand $L_{n}=\sum_{t=1}^{n}\log X_{t}$, it follows from the
definition that:

$\displaystyle\theta_{n,c}(\theta)$
$\displaystyle=-1+\frac{n\\!+\\!c}{\frac{c}{1\\!+\\!\theta}-L_{n}}\,.$

To compute the Bregman information gain, we rewrite the following integral
thanks to an affine change of variable in order to relate it to the Gamma
function:

$\displaystyle\int_{-\infty}^{-1}\exp\left(-c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\right)d\theta^{\prime}$
$\displaystyle=\left(-1\\!-\\!\theta\right)^{-c}e^{\frac{c\theta}{1\\!+\\!\theta}}\int_{-\infty}^{-1}\left(-1\\!-\\!\theta^{\prime}\right)^{c}e^{-\frac{c\theta^{\prime}}{1\\!+\\!\theta}}d\theta^{\prime}$
$\displaystyle=\left(-1\\!-\\!\theta\right)\Gamma(c)\left(\frac{e}{c}\right)^{c}\,.$

The expression of the Bregman information then follows immediately:

$\displaystyle\gamma_{n,c}(\theta)$
$\displaystyle=-\log\left(\frac{n\\!+\\!c}{(-1-\theta)L_{n}+c}\right)-n-c\log
c+(n\\!+\\!c)\log(n\\!+\\!c)+\log\Gamma(c)-\log\Gamma(n\\!+\\!c)\,.$

Moreover, we deduce from the expression of $\theta_{n,c}(\theta)$ and the
Bregman divergence that:

$\displaystyle\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta))$
$\displaystyle=-\log(-1-\theta^{\prime})+\log\left(\frac{n\\!+\\!c}{L_{n}-\frac{c}{1\\!+\\!\theta}}\right)-\frac{\left(L_{n}\\!-\\!\frac{c}{1\\!+\\!\theta}\right)\left(\theta^{\prime}\\!+\\!1\right)}{n\\!+\\!c}-1\,.$

Therefore, Theorem 3 combined with the natural parameter $\theta=-1-\alpha$
yields that w.p. $\geqslant 1-\delta$, $\forall n\in\mathbb{N}$,

$\displaystyle\alpha L_{n}-(n\\!+\\!c\\!+\\!1)\log\left(\alpha L_{n}+c\right)$
$\displaystyle\leqslant\log\frac{1}{\delta}-c\log
c-\log(n\\!+\\!c)+\log\Gamma(c)-\log\Gamma(n\\!+\\!c)\,.$ (11)

###  C.9 Chi-square

We finally consider $X\sim\chi^{2}(k)$ or, equivalently,
$X\sim\texttt{Gamma}\left(\frac{k}{2},\frac{1}{2}\right)$, i.e.,
$p_{k}(x)=\frac{(\frac{1}{2})^{\frac{k}{2}}}{\Gamma\left(\frac{k}{2}\right)}x^{\frac{k}{2}-1}e^{-\frac{x}{2}}$,
$x\geqslant 0$, $k\in\mathbb{N}$ (or $k\in\mathbb{R}_{+}\setminus\\{0\\}$ if
one considers Gamma distributions). This corresponds to an exponential family
model with parameter $\theta=\frac{k}{2}-1$, feature function $F(x)=\log x$
and log-partition function $\mathcal{L}(\theta)=(\theta+1)\log
2+\log\Gamma(\theta+1)$. The Bregman divergence between two parameters
$\theta^{\prime}$ and $\theta$ associated with $\mathcal{L}$ is given by

$\displaystyle\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=\texttt{KL}(P_{k},P_{k^{\prime}})=\frac{1}{2}(k-k^{\prime})\psi_{0}(k/2)-\log\frac{\Gamma(k/2)}{\Gamma(k^{\prime}/2)}~{},$

where $\psi_{0}(x)=\frac{d}{dx}\log\Gamma(x)$ denotes the digamma function.

We further have $\mathcal{L}^{\prime}(\theta)=\log 2+\psi_{0}(\theta+1)$ and
$\mathcal{L}^{\prime\prime}(\theta)=\psi_{1}(\theta+1)$, where
$\psi_{1}(x)=\frac{d^{2}}{dx^{2}}\log\Gamma(x)$ denotes the trigamma function.
Therefore $\mathcal{L}^{\prime}$ is invertible, and the parameter estimate is
given by

$\displaystyle\log
2+\psi_{0}(\theta_{n,c}(\theta)+1)=\frac{1}{n+c}\sum_{t=1}^{n}\log
X_{t}+\frac{c}{n+c}(\log 2+\psi_{0}(\theta+1))~{},$

yielding

$\displaystyle k_{n,c}(k):=2(1+\theta_{n,c}(\theta))$
$\displaystyle=2\psi_{0}^{-1}\left(\frac{1}{n+c}\sum_{t=1}^{n}\log
X_{t}+\frac{c}{n+c}\psi_{0}(k/2)-\frac{n}{n+c}\log 2\right)$
$\displaystyle=2\psi_{0}^{-1}\left(\frac{1}{n+c}K_{n}+\frac{c}{n+c}\psi_{0}(k/2)\right)\,,$

where $K_{n}=\sum\limits_{t=1}^{n}\log\frac{X_{t}}{2}$. Therefore, the Bregman
divergence rewrites as follows

$\displaystyle(n+c)\cdot\mathcal{B}_{\mathcal{L}}(\theta,\theta_{n,c}(\theta))$
$\displaystyle=(n+c)\cdot\texttt{KL}(P_{k_{n,c}(k)},P_{k})$
$\displaystyle=\frac{1}{2}(k_{n,c}(k)-k)\left(K_{n}+c\;\psi_{0}\left(\frac{k}{2}\right)\right)-(n+c)\log\frac{\Gamma\left(\frac{k_{n,c}(k)}{2}\right)}{\Gamma\left(\frac{k}{2}\right)}~{}.$

Now, we see that

$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}\\!-\\!c\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)\big{)}d\theta^{\prime}\\!$
$\displaystyle=\\!\int\frac{1}{2}\exp\left(\\!-\frac{c}{2}(k-k^{\prime})\psi_{0}(k/2)\right)\left(\frac{\Gamma(k/2)}{\Gamma(k^{\prime}/2)}\right)^{c}dk^{\prime}~{},$
$\displaystyle\int_{\mathbb{R}^{d}}\exp\big{(}\\!-\\!(n+c)\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta))\big{)}d\theta^{\prime}\\!$
$\displaystyle=\\!\int\frac{1}{2}\exp\left(\\!-\frac{n\\!+\\!c}{2}(k_{n,c}(k)-k^{\prime})\psi_{0}(k_{n,c}(k)/2)\right)\left(\frac{\Gamma(k_{n,c}(k)/2)}{\Gamma(k^{\prime}/2)}\right)^{n\\!+\\!c}\\!dk^{\prime}~{}.$

Therefore, the Bregman information gain writes

$\displaystyle\gamma_{n,c}(\theta)$
$\displaystyle=c\log\Gamma\left(\frac{k}{2}\right)-(n+c)\log\Gamma\left(\frac{k_{n,c}(k)}{2}\right)-c\;\frac{k}{2}\psi_{0}\left(\frac{k}{2}\right)+(n+c)\frac{k_{n,c}(k)}{2}\psi_{0}\left(\frac{k_{n,c}(k)}{2}\right)$
$\displaystyle\quad+\log\frac{\int\left(\Gamma\left(k^{\prime}/2\right)\right)^{-c}\exp\left(c\;\frac{k^{\prime}}{2}\psi_{0}\left(\frac{k}{2}\right)\right)dk^{\prime}}{\int\left(\Gamma\left(k^{\prime}/2\right)\right)^{-(n+c)}\exp\left((n+c)\frac{k^{\prime}}{2}\psi_{0}\left(\frac{k_{n,c}(k)}{2}\right)\right)dk^{\prime}}$
$\displaystyle=c\log\Gamma\left(\frac{k}{2}\right)-(n+c)\log\Gamma\left(\frac{k_{n,c}(k)}{2}\right)-c\;\frac{k}{2}\psi_{0}\left(\frac{k}{2}\right)+\frac{k_{n,c}(k)}{2}\left(K_{n}+c\psi_{0}\left(\frac{k}{2}\right)\right)$
$\displaystyle\quad+\log\frac{\int\left(\Gamma\left(k^{\prime}/2\right)\right)^{-c}\exp\left(c\;\frac{k^{\prime}}{2}\psi_{0}\left(\frac{k}{2}\right)\right)dk^{\prime}}{\int\left(\Gamma\left(k^{\prime}/2\right)\right)^{-(n+c)}\exp\left(\frac{k^{\prime}}{2}\left(K_{n}+c\psi_{0}\left(\frac{k}{2}\right)\right)\right)dk^{\prime}}$
$\displaystyle=c\log\Gamma\left(\frac{k}{2}\right)-(n+c)\log\Gamma\left(\frac{k_{n,c}(k)}{2}\right)-c\;\frac{k}{2}\psi_{0}\left(\frac{k}{2}\right)+\frac{k_{n,c}(k)}{2}\left(K_{n}+c\psi_{0}\left(\frac{k}{2}\right)\right)$
$\displaystyle\quad+\log
J\left(c,c\psi_{0}\left(\frac{k}{2}\right)\right)-\log
J\left(n\\!+\\!c,K_{n}\\!+\\!c\psi_{0}\left(\frac{k}{2}\right)\right)\,,$

where we define the auxiliary function
$J(a,b)=\int\exp\left(-a\log\Gamma\left(\frac{k^{\prime}}{2}\right)+b\frac{k^{\prime}}{2}\right)dk^{\prime}$.
Combining the above with Theorem 3, and after some simple algebra, we obtain
that w.p. $\geqslant 1-\delta$, $\forall n\in\mathbb{N}$,

$\displaystyle n\log\Gamma\left(\frac{k}{2}\right)-\frac{k}{2}K_{n}-\log
J\left(c,c\psi_{0}\left(\frac{k}{2}\right)\right)+\log
J\left(n\\!+\\!c,K_{n}\\!+\\!c\psi_{0}\left(\frac{k}{2}\right)\right)\leqslant\log\frac{1}{\delta}\,.$
(12)

###### Remark 9

The integral terms $\int dk^{\prime}$ in the above derive from the martingale
construction in A.1 and the mixture distribution $q(\theta|\alpha,\beta)$ over
the parameter $\theta\in\Theta$ of the exponential family. In the case of
$\texttt{Gamma}\left(\frac{k}{2},\frac{1}{2}\right)$ with unknown shape
$\frac{k}{2}>0$, we have $\theta=\frac{k}{2}$ and
$\Theta=\left(0,+\infty\right)$, therefore $dk^{\prime}$ corresponds to the
Lebesgue measure over $\left(0,+\infty\right)$. When restricted to the Chi-
square family, $\Theta=\mathbb{N}$, and $dk^{\prime}$ is instead the counting
measure, effectively turning integrals into discrete sums. In Appendix D, we
report figures using both versions, see Figure 6 and Figure 11.

###### Remark 10

The ratio of integrals (or infinite sums) in (12) can be efficiently
implemented using a simple integration scheme (or by truncation). Indeed, for
a given $k_{\max}\in\mathbb{N}$, $b\ll 1$ and $B\gg 1$, we define
$x_{k^{\prime}}=b+\frac{k^{\prime}}{k_{\max}}\left(B-b\right)$ for
$k^{\prime}=0,\dots,k_{\max}$, so that

$\displaystyle\log\frac{\int_{0}^{+\infty}\left(\Gamma\left(\frac{k^{\prime}}{2}\right)\right)^{-c}\exp\left(\frac{k^{\prime}}{2}c\psi_{0}\left(\frac{k}{2}\right)\right)dk^{\prime}}{\int_{0}^{+\infty}\left(\Gamma\left(\frac{k^{\prime}}{2}\right)\right)^{-(n+c)}\exp\left(\frac{k^{\prime}}{2}\left(c\psi_{0}\left(\frac{k}{2}\right)+\sum_{t=1}^{n}\log\frac{X_{t}}{2}\right)\right)dk^{\prime}}$
$\displaystyle\approx\underset{k^{\prime}=1,\dots,k_{\max}}{\emph{logsumexp}}\left(-c\log\Gamma\left(\frac{x_{k^{\prime}}}{2}\right)+\frac{cx_{k^{\prime}}}{2}\psi_{0}\left(\frac{k}{2}\right)+\log\left(x_{k^{\prime}}-x_{k^{\prime}-1}\right)\right)$
$\displaystyle\quad-\underset{k^{\prime}=1,\dots,k_{\max}}{\emph{logsumexp}}\left(-\left(n+c\right)\Gamma\left(\frac{x_{k^{\prime}}}{2}\right)+\frac{x_{k^{\prime}}}{2}\left(c\psi_{0}\left(\frac{k}{2}+\sum_{t=1}^{n}\log\frac{X_{t}}{2}\right)\right)+\log\left(x_{k^{\prime}}-x_{k^{\prime}-1}\right)\right).$

Similarly, we have

$\displaystyle\log\frac{\sum\limits_{k^{\prime}=1}^{+\infty}\left(\Gamma\left(\frac{k^{\prime}}{2}\right)\right)^{-c}\exp\left(\frac{k^{\prime}}{2}c\psi_{0}\left(\frac{k}{2}\right)\right)}{\sum\limits_{k^{\prime}=1}^{+\infty}\left(\Gamma\left(\frac{k^{\prime}}{2}\right)\right)^{-(n+c)}\exp\left(\frac{k^{\prime}}{2}\left(c\psi_{0}\left(\frac{k}{2}\right)+\sum_{t=1}^{n}\log\frac{X_{t}}{2}\right)\right)}$
$\displaystyle\approx\underset{k^{\prime}=1,\dots,k_{\max}}{\emph{logsumexp}}\left(-c\log\Gamma\left(\frac{k^{\prime}}{2}\right)+\frac{k^{\prime}}{2}c\psi_{0}\left(\frac{k}{2}\right)\right)$
$\displaystyle\quad-\underset{k^{\prime}=1,\dots,k_{\max}}{\emph{logsumexp}}\left(-\left(n+c\right)\Gamma\left(\frac{k^{\prime}}{2}\right)+\frac{k^{\prime}}{2}\left(c\psi_{0}\left(\frac{k}{2}+\sum_{t=1}^{n}\log\frac{X_{t}}{2}\right)\right)\right),$

The final steps correspond to the right-rectangular scheme over $(b,B)$ with
$k_{\max}$ steps and the truncation to the first $k_{\max}$ terms
respectively. Note the use of
$\underset{k^{\prime}=1,\dots,k_{\max}}{\emph{logsumexp}}(z)=\log\sum\limits_{k^{\prime}=1}^{k_{\max}}\exp\left(z_{k^{\prime}}\right)$
for $z\in\mathbb{R}^{k_{\max}}$, which is efficiently implemented in many
libraries for scientific computing and better handles summation of large
numbers. Empirically, we found that $\log b=-10,\log B=10$ and $k_{\max}=2000$
provided sufficient accuracy and that using finer approximation schemes did
not significantly impact the numerical results.

###  C.10 Poisson

We consider $X\sim\text{Poisson}\left(\lambda\right)$, where $\lambda>0$ is
unknown. We recall that the distribution of $X$ is supported on $\mathbb{N}$
with probability mass function
$\mathbb{P}\left(X=k\right)=\frac{\lambda^{k}e^{-\lambda}}{k!}$. This
corresponds to a one-dimensional exponential family model with parameter
$\theta=\log\lambda\in\mathbb{R}$, feature function $F(x)=x$ and log-partition
function $\mathcal{L}(\theta)=e^{\theta}$ (which is invertible on
$\mathbb{R}$). The Bregman divergence between two parameters $\theta^{\prime}$
and $\theta$ is therefore
$\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta)=e^{\theta^{\prime}}-e^{\theta}-(\theta^{\prime}-\theta)e^{\theta}$.

Using the shorthand $S_{n}=\sum_{t=1}^{n}X_{t}$, it follows from the
definition that:

$\displaystyle\theta_{n,c}(\theta)$
$\displaystyle=\log\left(\frac{S_{n}\\!+\\!ce^{\theta}}{n\\!+\\!c}\right)\,.$

The Bregman information gain is expressed using the auxiliary function
$I(a,b)=\int_{-\infty}^{+\infty}e^{-ae^{\theta}+b\theta}d\theta$ as:

$\displaystyle\gamma_{n,c}(\theta)$
$\displaystyle=c(1-\theta)e^{\theta}-(n\\!+\\!c)\left(1-\theta_{n,c}(\theta)\right)e^{\theta_{n,c}(\theta)}+\log
I\left(c,ce^{\theta}\right)-\log
I\left(n\\!+\\!c,(n\\!+\\!c)e^{\theta_{n,c}(\theta)}\right)$
$\displaystyle=c(1-\theta)e^{\theta}-\left(1-\log\left(\frac{S_{n}\\!+\\!ce^{\theta}}{n\\!+\\!c}\right)\right)\left(S_{n}\\!+\\!ce^{\theta}\right)+\log
I\left(c,ce^{\theta}\right)-\log
I\left(n\\!+\\!c,S_{n}\\!+\\!ce^{\theta}\right)\,.$

Moreover, we deduce from the expression of $\theta_{n,c}(\theta)$ and the
Bregman divergence that:

$\displaystyle\mathcal{B}_{\mathcal{L}}(\theta^{\prime},\theta_{n,c}(\theta))$
$\displaystyle=e^{\theta^{\prime}}-\left(\frac{S_{n}\\!+\\!ce^{\theta}}{n\\!+\\!c}\right)-\left(\theta^{\prime}-\log\left(\frac{S_{n}\\!+\\!ce^{\theta}}{n\\!+\\!c}\right)\right)\left(\frac{S_{n}\\!+\\!ce^{\theta}}{n\\!+\\!c}\right)\,.$

Therefore, Theorem 3 combined with the natural parametrization
$\theta=\log\lambda$ yields that w.p. $\geqslant 1-\delta$, $\forall
n\in\mathbb{N}$,

$\displaystyle n\lambda-S_{n}\log\lambda$
$\displaystyle\leqslant\log\frac{1}{\delta}+\log I\left(c,c\lambda\right)-\log
I\left(n\\!+\\!c,S_{n}\\!+\\!c\lambda\right)\,.$

###### Remark 11

Although, to the best of our knowledge, the integral $I(a,b)$ does not have a
closed-form expression, it can be numerically estimated up to arbitrary
precision. We recommend the same implementation as discussed in Remark 10,
using the _logsumexp_ operator for stability, and refer to the code for
further details.

###  C.11 Summary table

We summarize the high probability confidence sets of Theorem 3 derived in this
Table 2.

Table 2: Summary of Bregman confidence sets given by Theorem 3 for
representative families. Throughout, the following notations are used:
$S_{n}=\sum\limits_{t=1}^{n}X_{t}$, $\widehat{\mu}_{n}=\frac{S_{n}}{n}$,
$Q_{n}(\mu)=\sum\limits_{t=1}^{n}\left(X_{t}-\mu\right)^{2}$,
$Z_{n}(\mu,\sigma)=\frac{Q_{n}(\mu)}{\sigma^{2}}$,
$S_{n}^{(k)}=\sum\limits_{t=1}^{n}X_{t}^{k}$, $L_{n}=\sum\limits_{t=1}^{n}\log
X_{t}$, $K_{n}=\sum\limits_{t=1}^{n}\log\frac{X_{t}}{2}$,
$I(a,b)=\int_{-\infty}^{+\infty}e^{-ae^{\theta}+b\theta}d\theta$,
$J(a,b)=\int_{0}^{\infty}\exp\left(-a\log\Gamma\left(\frac{k}{2}\right)\\!+\\!\frac{bk}{2}\right)dk$
($dk$ is the Lebesgue measure if $k\in\mathbb{R}_{+}$ and the counting measure
if $k\in\mathbb{N}$).
Name | Parameters | Formula
---|---|---
Gaussian | $\mu\in\mathbb{R}$ | $\frac{1}{n+c}\frac{\left(S_{n}-n\mu\right)^{2}}{2\sigma^{2}}\leqslant\log\frac{1}{\delta}+\frac{1}{2}\log\frac{n+c}{c}$
Gaussian | $\sigma\in\mathbb{R}_{+}$ |  $\frac{Q_{n}(\mu)}{2\sigma^{2}}-\left(\frac{n\\!+\\!c}{2}\\!+\\!1\right)\log\left(\frac{Q_{n}(\mu)+c}{2\sigma^{2}}\right)$ $\quad\leqslant\log\frac{1}{\delta}\\!+\\!\log\frac{\Gamma\left(\frac{c}{2}\\!+\\!2\right)}{\Gamma\left(\frac{n\\!+\\!c}{2}\\!+\\!2\right)}\\!-\\!\frac{n}{2}\log 2\\!-\\!(\frac{c}{2}\\!+\\!1)\log c\\!+\\!(\frac{n\\!+\\!c}{2}\\!+\\!1)\log(n\\!+\\!c)$
Gaussian |  $\mu\in\mathbb{R}$ $\sigma\in\mathbb{R}_{+}$ |  $\frac{1}{2}Z_{n}(\mu,\sigma)-\frac{n+c+3}{2}\log\left(\frac{n}{n\\!+\\!c}Z_{n}(\widehat{\mu}_{n},\sigma)+\frac{c}{n\\!+\\!c}Z_{n}(\mu,\sigma)+c\right)$ $\quad\leqslant\log\frac{1}{\delta}-\frac{n}{2}\log 2-\left(\frac{c}{2}\\!+\\!2\right)\log c+\frac{1}{2}\log\left(n\\!+\\!c\right)$ $\qquad+\log\Gamma\left(\frac{c\\!+\\!3}{2}\right)-\log\Gamma\left(\frac{n\\!+\\!c\\!+\\!3}{2}\right)$
Bernoulli | $\mu\in[0,1]$ |  $S_{n}\log\\!\frac{1}{\mu}\\!+\\!(n\\!-\\!S_{n})\log\\!\frac{1}{1\\!-\\!\mu}\\!+\\!\log\\!\frac{\Gamma(S_{n}\\!+\\!c\mu)\Gamma(n\\!-\\!S_{n}\\!+\\!c(1\\!-\\!\mu))}{\Gamma(c\mu)\Gamma(c(1\\!-\\!\mu))}$ $\quad\leqslant\log\\!\frac{1}{\delta}\\!+\\!\log\\!\frac{\Gamma(n\\!+\\!c)}{\Gamma(c)}$
Exponential | $\mu\in\mathbb{R}_{+}$ |  $\frac{S_{n}}{\mu}\\!-\\!(n\\!+\\!c\\!+\\!1)\log\left(\frac{S_{n}}{\mu}+c\right)$ $\quad\leqslant\log\frac{1}{\delta}\\!+\\!\log\frac{\Gamma(c)}{\Gamma(n\\!+\\!c)}\\!-\\!\log(n\\!+\\!c)\\!-\\!c\log c$