# Large deviations rates for stochastic gradient descent with strongly convex
functions

Dragana Bajovic
Faculty of Technical Sciences, University of Novi Sad, Novi Sad, Serbia
<EMAIL_ADDRESS>Dusan Jakovetic
Faculty of Sciences, University of Novi Sad, Novi Sad, Serbia
<EMAIL_ADDRESS>Soummya Kar
Carnegie Mellon University, Pittsburgh, PA, USA
<EMAIL_ADDRESS>This work is supported by the European Union’s Horizon
2020 Research and Innovation program under grant agreement No 957337. The
paper reflects only the view of the authors and the Commission is not
responsible for any use that may be made of the information it contains.

###### Abstract

Recent works have shown that high probability metrics with stochastic gradient
descent (SGD) exhibit informativeness and in some cases advantage over the
commonly adopted mean-square error-based ones. In this work we provide a
formal framework for the study of general high probability bounds with SGD,
based on the theory of large deviations. The framework allows for a generic
(not-necessarily bounded) gradient noise satisfying mild technical
assumptions, allowing for the dependence of the noise distribution on the
current iterate. Under the preceding assumptions, we find an upper large
deviations bound for SGD with strongly convex functions. The corresponding
rate function captures analytical dependence on the noise distribution and
other problem parameters. This is in contrast with conventional mean-square
error analysis that captures only the noise dependence through the variance
and does not capture the effect of higher order moments nor interplay between
the noise geometry and the shape of the cost function. We also derive exact
large deviation rates for the case when the objective function is quadratic
and show that the obtained function matches the one from the general upper
bound hence showing the tightness of the general upper bound. Numerical
examples illustrate and corroborate theoretical findings.

##  1 Introduction

The large deviations theory represents a well-established principled approach
for studying _rare events_ that occur with stochastic processes, e.g., [7].
Typically, we are concerned with a sequence of rare events $E_{k}$ related
with the stochastic process of interest, indexed by, e.g., time $k$. In this
setting, the probability of event $E_{k}$, $k=1,2,...$ typically decays
exponentially in $k$; the large deviations theory then enables to quantify
this exponential rate. Such an approach has found many applications in
statistics [4], mechanics [24], communications [22], and information theory
[23].

To be more concrete, consider an example of a sequence of random vectors
$X_{k}$ taking values in ${\mathbb{R}}^{d}$ that converge, e.g., almost
surely, to a (deterministic) limit point $x^{\star}\in{\mathbb{R}}^{d}$. The
rare event of interest $E_{k}$ can then be, for example,
$E_{k}=\\{\|X_{k}-x^{\star}\|\geq\delta\\}$, for some positive quantity
$\delta$, with $\|\cdot\|$ denoting the Euclidean norm. Equivalently, $E_{k}$
can be represented as $\\{X_{k}\in C_{\delta}\\}$, where $C_{\delta}$ is the
complement of the $l_{2}$ ball of radius $\delta$ centered at $x^{\star}$.
Large deviations analysis then aims at discovering the corresponding rate of
decay, i.e., the inaccuracy rate $\mathbf{I}(C_{\delta})$:

$\mathbb{P}\left(X_{k}\in
C_{\delta}\right)=e^{-k\,\mathbf{I}(C_{\delta})+o(k)},$ (1)

where $o(k)$ denotes terms growing slower than linearly with $k$. The
inaccuracy rate $\mathbf{I}(C_{\delta})$ can usually be expressed via the so
called _rate function_ $I:\mathbb{R}^{d}\mapsto\mathbb{R}$, according to the
following formula [1]:

$\mathbf{I}(C_{\delta})=\inf_{x\in C_{\delta}}I(x).$ (2)

Differently from the set function $\mathbf{I}$, the rate function $I$ does not
depend on the region $C_{\delta}$; that is, when $C_{\delta}$ changes, only
the region over which we minimize in (2) changes, while the function remains
unchanged. Furthermore, this is true for arbitrary set $C_{\delta}$. This
means that, once the rate function is computed, the corresponding inaccuracy
rate can be obtained via (2) for a new given region of interest.

In this paper, we are interested in applying the large deviations theory to
analyzing the stochastic gradient descent (SGD) method. SGD is a simple but
widely used optimization method that finds numerous practical applications,
such as training machine learning and deep learning models, e.g., [21, 11,
17]. More precisely, we consider unconstrained optimization problems where the
goal is to minimize a smooth, strongly convex function
$f:\,{\mathbb{R}}^{d}\rightarrow\mathbb{R}$, via the SGD method of the form:

$X_{k+1}=X_{k}-\alpha_{k}\,(\nabla f(X_{k})-Z_{k}).$ (3)

Here, $k=1,2,...$ is the iteration counter, $\alpha_{k}=a/k$, $a>0$ is the
step-size, and $Z_{k}$ is a zero-mean gradient noise that may depend on
$X_{k}$. In this context, we are interested in solving for (1) and (2) for the
SGD method (3), where now $x^{\star}$ is interpreted as the (deterministic)
global minimizer of $f$. In other words, we are interested in finding (or
approximating) the rate function $I(x)$ that quantifies the “tails” or “rare
events” of how the SGD sequence iterates $X_{k}$ deviate from the solution
$x^{\star}$.

Clearly, evaluating (2) for SGD is of significant interest. It readily
provides insights into the high-probability bounds for SGD that have been
subject of much research effort recently, [8, 9, 16, 10, 6]. However, unlike
the typical high probability bound studies, the large deviations approach here
is fully flexible with respect to the choice of set $C_{\delta}$; e.g., the
$l_{2}$-ball complement may be replaced with an arbitrary open set, such as
$l_{p}$ norm complement of an arbitrary $l_{p}$-norm. While large deviations
theory is a well-established field, there has been a limited body of work that
applies large deviations to the analysis of SGD. Reference [26] is concerned
with large deviations analysis for a scalar stochastic process equivalent to
SGD in one dimension. The authors of [14] study large deviations of SGD when
the step-size converges to zero; however, they are not concerned with large
deviations when the iteration counter $k$ increases – the case of our interest
here.

Contributions. In this paper, we are interested in evaluating the large
deviations rates in (1) and (2) for the SGD method, when the objective
function $f$ is smooth and strongly convex. Our main contributions are as
follows. When $f$ is a (strongly convex) quadratic function, we establish the
so-called full large deviations principle for the sequence $X_{k}$. This means
that we evaluate rate function $I(x)$ exactly, i.e., the corresponding rare
event probability is computed exactly, with upper and lower bounds matched, up
to exponentially decaying factors. We further explicitly quantify the rate
function $I(x)$ as a function of the distribution of the gradient noise. This
reveals a significant influence of higher order moments on the performance (in
the sense of rare event probabilities) of SGD. This is in contrast with
conventional SGD analyses, that typically capture only the dependence on the
gradient noise variance. The large deviations principle for quadratic
functions is established under a very general class of gradient noise
distributions that are essentially only required to have a finite moment
generating function. Next, for generic smooth and strongly convex costs $f$,
we establish a large deviations upper bound (a lower bound on function $I(x)$)
that certifies an exponential decay of the rare event probabilities in (1)
with SGD. This is achieved when the distribution of the gradient noise is sub-
Gaussian. We further show that the obtained large deviations upper bound is
tight, as the corresponding rate function actually matches, up to higher order
factors, the exact rate function that we formerly establish for the quadratic
costs.

Our results are related with high probability bounds-type studies of SGD and
related stochastic methods [12, 8, 9, 16, 10]. Therein, for a given $\delta>0$
and a confidence level $1-\beta$, $\beta\in(0,1)$, the goal is to find
$K(\delta,\beta)$ such that $f(X_{k})-f(x^{\star})\leq\delta$ with probability
at least $1-\beta$, for all $k\geq K(\delta,\beta)$. The works [8, 9, 16, 10]
provide estimates of $K(\delta,\beta)$ that depend _logarithmically_ on
$\beta$. In more detail, [8, 9] establish high probability bounds for the
stochastic gradient methods therein assuming sub-Gaussian gradient noises. The
work [16] calculates the corresponding bounds for the basic SGD and the mirror
descent that utilize a gradient truncation technique, while relaxing the noise
sub-Gaussianity. The work [10] establishes high probability bounds for an
accelerated SGD that also utilizes a clipping nonlinearity. The large
deviations rates in (1) and (2) \- give estimates of $K(\delta,\beta)$ that
also depend logarithmically on $\beta$, when $\beta$ is small (goes to
zero).111It is easy to see this by noting that, for $\mu$-strongly convex
costs, we have $f(x)-f(x^{\star})\geq\frac{\mu}{2}\|x-x^{\star}\|^{2}$, for
all $x\in{\mathbb{R}}^{d}$, requiring that the the right hand side of (1) be
less than $\beta$, and reverse-engineering the smallest iterate $k$ for which
the latter holds.

Compared with existing high probability bound works, our results give the
_exact_ (tight) exponential decay rate in (2), and for an _arbitrary set_ that
does not contain $x^{\star}$, not only the Euclidean ball complements. To be
concrete, the closest results to ours are obtained in [12]. While they are not
directly concerned with obtaining large deviations rates, their results (with
some additional work) lead to an exponential decay rates for Euclidean ball
complements. In contrast, our results work for arbitrary open sets.
Furthermore, focusing only on Euclidean ball complements, our results provide
much tighter exponential rate bounds. Specifically, as we show in the paper,
the exponential rate that we provide captures the interplay between the noise
geometry and the cost function curvature, see Section 4.2 for details. From
the technical perspective, this is achieved by working directly with the SGD
iterates, as opposed to working with the distance of the iterates from the
solution. To do so, we derive a novel set of techniques that build upon the
large deviations theory rather than on martingale concentration inequalities.

The current paper is also related with large deviations analyses of stochastic
processes that arise with distributed inference, such as estimation and
detection. Distributed detection has been studied in [5], for Gaussian
observations, and in [3], for generic observations. The work [19] evaluates
large deviations of the local states with a distributed detection method, when
the step size parameter decreases. Reference [20] further analyzes the non-
exponential terms and consider directed networks for a similar problem. The
paper [18] considers distributed detection with 1-bit messages. [13] consider
social learning problems. Reference [2] analyzes large deviations for
distributed estimation and social learning. Unlike these works on distributed
inference, we are not directly concerned with distributed systems; also, the
cost functions that we consider are more general and, unlike the works above,
do not result in linear (distributed averaging) dynamics; hence, novel tools
for large deviations analysis are required here.

The rest of the paper is organized as follows. Section 2 explains the problem
that we consider and gives the required preliminaries. Section 3 provides the
main results of the paper – a large deviations upper bound for generic costs,
and the full (exact) large deviations rates for quadratic costs. Specializing
to the Gaussian noise, Section 4 provides analytical, closed-form expressions
for the large deviations rate function. Section 5 gives the proof of the main
lemma underlying the upper bound for the general functions. Finally, we
conclude in Section 6. Appendix contains additional insights and examples,
numerical results, and missing proofs.

##  2 Setup and preliminaries

We consider unconstrained optimization problem of the form

$\min_{x\in\mathbb{R}^{d}}f(x).$ (4)

We assume that $f$ is $L$-smooth and $\mu$-strongly convex, and that the
stepsize in algorithm (3) is of the form $\alpha_{k}=a/(k+b),$ where
$a,\,b>0.$

######  Assumption 1.

We assume that $f$ is twice differentiable, $L$-smooth and $\mu$-strongly
convex, where $0<\mu\leq L$.

Strong convexity implies uniqueness of the solution of (4), which is denoted
by $x^{\star}.$ We make the following assumption regarding the stepsize
parameter $a.$

######  Assumption 2.

The stepsize parameter $a$ satisfies $a\mu>1$.

Assumption 1 is standard in the analysis of optimization methods, i.e., it
corresponds to a standard class of functions over which an optimization method
analysis is carried out. Assumption 2 is required for some asymptotic
arguments ahead, as $k\rightarrow\infty$. In practice, it may be restrictive
that the constant $a$ is too large in the step-size choice $a/k$, as at the
initial iterations (small $k$’s), we would have very large step-sizes. This is
alleviated by having an appropriately chosen constant $b>1$.

We denote by $\tilde{g}(X_{k})$ the stochastic gradient of $f$ returned by the
gradient oracle at the current iterate $X_{k},$ and by $g(X_{k})$ the (exact)
gradient of $f$ at the current iterate $X_{k}.$ The difference between
$\tilde{g}(X_{k})$ and $g(X_{k})$ (the gradient “noise”) is denoted by
$Z_{k}=g(X_{k})-\tilde{g}(X_{k})$. We make the following assumptions on
$Z_{k}.$

######  Assumption 3.

  1. 1.

For each $k$, $Z_{k}$ depends on the past iterates only through $X_{k}$.

  2. 2.

For each $k$, the distribution of $Z_{k}$ given $X_{k}$ depends on $X_{k}$
only through its realization and does not depend on the current iterate index,
$k$.

  3. 3.

For any given $x,$ $\mathbb{E}[Z_{k}|X_{k}=x]=0,$ i.e., conditioned on the
current iterate, the noise is zero-mean.

Assumption 3 allows for a general gradient noise that may actually depend on
the current iterate $X_{k}$. This is a more general setting than the
frequently studied case when $Z_{k}$ is i.i.d. and independent of $X_{k}$.
Item 3. of Assumption 3 says that, conditioned on the current iterate, the
noise is zero-mean on average. This is also a standard bias-free noise
assumption. Finally, note that items 1. and 2. in Assumption 3 typically hold
in machine learning settings. Therein, the goal is typically to minimize a
population loss $f(x)=\mathbb{E}[\phi(x,v)]$ where the expectation is taken
over the distribution of the data $v$, and $\phi$ is an instantaneous loss
function. Given that, at some iteration $k$, $X_{k}$ takes a value $x$, the
gradient noise equals $\nabla\phi(x,v_{k})-\mathbb{E}[\nabla\phi(x,v)]$, where
$v_{k}$ is the data point sampled at iteration $k$. Then, items 1. and 2. are
clearly satisfied, provided that the data sampling process is independent of
the evolution of $X_{k}$.

For $x\in\mathbb{R}^{d},$ we denote by $H(x)$ the Hessian matrix of $f$
computed at $x.$ For short, we denote $H^{\star}=H(x^{\star})$, i.e.,
$H^{\star}$ is the Hessian matrix of $f$ computed at $x^{\star}.$ For any
$x\in\mathbb{R}^{n}$, define $h:\mathbb{R}^{n}\mapsto\mathbb{R}^{d}$ as the
residual of the first order Taylor’s approximation of the gradient $g$ at
$x^{\star}$,

$h(x)=g(x)-H^{\star}(x-x^{\star}),$ (5)

for $x\in\mathbb{R}^{n}$. For each $\delta>0,$ define also

$\overline{h}(\delta)=\sup_{x\in\mathbb{B}_{x^{\star}}(\delta)}\|h(x)\|,$ (6)

where $B_{x}(\delta)$ denotes the Euclidean ball in $\mathbb{R}^{d}$ of radius
$\delta\geq 0,$ centered at $x.$ The following result holds by a well-known
corollary of Taylor’s remainder theorem.

######  Lemma 1.

There holds $\overline{h}(\delta)=o(\delta),$ i.e., $\lim_{\delta\rightarrow
0}\frac{\overline{h}(\delta)}{\delta}=0.$

######  Remark 1.

Clearly, when $f$ is quadratic, $H(x)$ is constant for all
$x\in\mathbb{R}^{d}$ and equal to $H^{\star},$ implying $h(x)\equiv 0$ and
also $\overline{h}(\delta)\equiv 0.$

######  Remark 2.

Quantity $h(x)$ can be explicitly characterized if, in addition, it is assumed
that the Hessian of function $f$ is Lipschitz continuous, i.e., if
$\left\|H(x)-H(y)\right\|\leq L_{H}\,\|x-y\|$, for all
$x,y\in{\mathbb{R}}^{d}$, for some nonnegative constant $L_{H}$. It is easy to
show that, in this case, we have $\|h(x)\|\leq L_{H}\,\|x-x^{\star}\|^{2}$,
for any $x\in{\mathbb{R}}^{d}.$ The latter implies a quadratic upper bound in
$\delta$ on $\overline{h}(\delta)$, i.e., $\overline{h}(\delta)\leq
L_{H}\delta^{2},$ for each $\delta\geq 0.$

###  2.1 Distance to solution recursion

For analytical purposes, it is of interest to study the squared distance to
solution of the current iterates $\|X_{k}-x^{\star}\|^{2}.$ To characterize
the evolution of this quantity, we use standard arguments that follow from
strong convexity and Lipschitz smoothness:

$\displaystyle\|X_{k+1}-x^{\star}\|^{2}$
$\displaystyle\leq\left(1-2\alpha_{k}\mu+2\alpha_{k}^{2}L^{2}\right)\|X_{k}-x^{\star}\|^{2}+2\alpha_{k}(X_{k}-x^{\star})^{\top}Z_{k}$
$\displaystyle\,\,\,+2\alpha_{k}^{2}\|Z_{k}\|^{2};$ (7)

details of the derivations can be found in Appendix A.

We introduce function $\beta_{k}:\mathbb{R}^{2}\mapsto\mathbb{R},$ defined by
$\beta_{k}(u,v)=1-\alpha_{k}u+\alpha_{k}^{2}v.$ Similarly, for any two
iteration indices $l\leq k,$ we define
$\beta_{k,l}:\mathbb{R}^{2}\mapsto\mathbb{R}$ by
$\beta_{k,l}(u,v)=\beta_{k}(u,v)\cdots\beta_{l}(u,v).$ The following technical
lemma providing bounds on the product functions $\beta_{k,l}$ will be useful
for the study of recursion (2.1) as well as other similar recursions that will
emerge from the analysis.

######  Lemma 2.

Let $l$ and $k$ be two iteration indices such that $l<k$. For any nonnegative
$u,\,v\in\mathbb{R},$ and $\alpha_{k}=a/(k+b),$ where $b\geq 1,$ there holds:

  1. 1.

$\beta_{k,l}(u,v)\leq\left(\frac{l+b}{k+b+1}\right)^{au}e^{\frac{a^{2}v}{l+b-1}}$;

  2. 2.

for each $l$ such that $l+b\geq\frac{5au}{2}$, there holds
$\beta_{k,l}(u,v)\geq\left(\frac{l+b-1}{k+b}\right)^{au}e^{-\frac{a^{2}u^{2}}{l+b-1}};$

The proof of Lemma 2 is given in Appendix A.

Finally, for each iteration index $k$, we denote by $\mu_{k}$ the Borel
measure on $\mathbb{R}^{d}$ induced by $X_{k}$. Similarly, we denote by
$\nu_{k}$ the Borel measure induced by $\|X_{k}-x^{\star}\|.$

###  2.2 Large deviations preliminaries

We next give a definition of the rate function and the large deviations
principle.

Rate function $I$ and the large deviations principle.

######  Definition 1 (Rate function $I$ [7]).

Function $I:\mathbb{R}^{d}\mapsto[0,+\infty]$ is called a _rate function_ if
it is lower semicontinuous, or, equivalently, if its level sets are closed.
If, in addition, the level sets of $I$ are compact (i.e., closed and bounded),
then $I$ is called a good rate function.

######  Definition 2 (The large deviations principle [7]).

Suppose that $I:\mathbb{R}^{d}\mapsto[0,+\infty]$ is lower semicontinuous. A
sequence of measures $\mu_{k}$ on
$\left(\mathbb{R}^{d},\mathcal{B}\left(\mathbb{R}^{d}\right)\right)$, $k\geq
1$, is said to satisfy the large deviations principle (LDP) with rate function
$I$ if, for any measurable set $D\subseteq\mathbb{R}^{d}$, the following two
conditions hold:

  1. 1.

$\displaystyle\limsup_{k\rightarrow+\infty}\,\frac{1}{k}\,\log\mu_{k}(D)\leq-\,\inf_{x\in\overline{D}}I(x);$

  2. 2.

$\displaystyle\liminf_{k\rightarrow+\infty}\,\frac{1}{k}\,\log\mu_{k}(D)\geq-\,\inf_{x\in
D^{\mathrm{o}}}I(x).$

Log-moment generating functions of the noise $Z_{k}$ and the iterates $X_{k}.$
Following Assumption 3, we define the conditional LMGF of $Z_{k}$ given the
last iterate $X_{k}.$

######  Definition 3 (Conditional LMGF of $Z_{k}$).

We denote by $\Lambda(\cdot;x)$ the log-moment generating function (LMGF) of
$Z_{k}$ given $X_{k}=x$,

$\Lambda(\lambda;x):=\log\mathbb{E}\left[\left.e^{\lambda^{\top}Z_{k}}\right|X_{k}=x\right],\>\>\mathrm{for\>\>}\lambda,x\in\mathbb{R}^{d}.$
(8)

It will also be useful to define the conditional moment-generating function of
$\|Z_{k}\|^{2}$, which we denote by $M(\cdot;x)$:

$M(\nu;x):=\mathbb{E}\left[\left.e^{\nu\|Z_{k}\|^{2}}\right|X_{k}=x\right],$
(9)

for $\nu\in\mathbb{R},$ $x\in\mathbb{R}^{d}$. By the inequality $e^{x}\leq
x+e^{x^{2}},$ which holds for all $x\in\mathbb{R},$ we have
$\mathbb{E}\left[\left.e^{\lambda^{\top}Z_{k}}\right|X_{k}\right]\leq\mathbb{E}[\lambda^{\top}Z_{k}|X_{k}]+\mathbb{E}\left[e^{(\lambda^{\top}Z_{k})2}|X_{k}\right]\leq\mathbb{E}\left[\left.e^{\|\lambda^{2}\|\|Z_{k}\|^{2}}\right|X_{k}\right],$
where we used the Cauchy-Schwartz inequality, for the second term, and the
fact that $Z_{k}$ is zero-mean, for the first term. Thus,

$\Lambda(\lambda;x)\leq\log M(\|\lambda^{2}\|;x)$ (10)

for any realization $x$ of $X_{k}.$

Lemma 3 lists properties of $\Lambda$ that will be used in the paper.

######  Lemma 3 (Properties of $\Lambda$).

For any given $x\in\mathbb{R}^{d}$ the following properties hold:

  1. 1.

$\Lambda(\cdot;x)$ is convex and differentiable in the interior of its domain;

  2. 2.

$\Lambda(0;x)=0$ and $\nabla\Lambda(0;x)=\mathbb{E}[Z_{k}|X_{k}=x]=0$;

  3. 3.

$\Lambda(\lambda;x)\geq 0,$ for each $\lambda.$

###### Proof.

Convexity and differentiability are general properties of log-moment
generating functions [7], as well as the zero value at the origin property and
also that the gradient at the origin equals the mean vector;
$\nabla\Lambda(0;x)=0$ follows by the assumption that the noise is zero-mean,
Assumption 3. The non-negativity from Part 3 follows by invoking convexity and
exploiting the two properties from part 2, i.e., for any $x\in\mathbb{R}^{d}$:
$\Lambda(\lambda;x)\geq\Lambda(0;x)+\nabla\Lambda(0;x)^{\top}\lambda=0.$ ∎

######  Example 1.

To illustrate the LMGF function $\Lambda,$ we consider the case when,
conditioned on an arbitrary realization $X_{k}=x,$ the gradient noise $Z_{k}$
is Gaussian, with mean vector equal to zero vector and covariance matrix
$\Sigma(x).$ Using standard formula for the LMGF of a Gaussian multivariate,
we have

$\Lambda(\lambda;x)=\frac{1}{2}\lambda^{\top}S(x)\lambda,$ (11)

for $\lambda\in\mathbb{R}^{d}.$ We note that when the gradient noise $Z_{k}$
is independent of the current iterate $X_{k}$, the indices $X_{k}$ in the
preceding formula can be omitted, i.e., the expression for $\Lambda$
simplifies to $\Lambda(\lambda;X_{k})=\frac{1}{2}\lambda^{\top}S\lambda,$ for
all realizations $X_{k}.$

It will also be of interest to define the (unconditional) log-moment
generating function of the iterates $X_{k}$.

######  Definition 4 (LMGF of $X_{k}-x^{\star}$).

We let $\Gamma_{k}$ denote the (unconditional) moment generating function of
$X_{k}$,

$\Gamma_{k}(\lambda):=\mathbb{E}\left[e^{\lambda^{\top}(X_{k}-x^{\star})}\right],$
(12)

for $\lambda\in\mathbb{R}^{d}$. The (unconditional) log-moment generating
function of $X_{k}$ is then given by $\log\Gamma_{k}$.

We assume that the initial iterate $X_{1}$ is deterministic222We note that
this assumption can be relaxed to allow for random initial iterate; see
Appendix D for details.. Hence, $\Gamma_{1}$ is finite for all
$\lambda\in\mathbb{R}^{d}$.

We assume that the family of functions $\Lambda(\cdot;x)$ satisfy the
following regularity conditions.

######  Assumption 4 (Lipschitz continuity in $x$).

There exists a constant $L_{\Lambda}$ such that for every $\lambda$, $x$,
$y\in\mathbb{R}^{d}$, there holds:

$|\Lambda(\lambda;x)-\Lambda(\lambda;y)|\leq
L_{\Lambda}\|\lambda\|^{2}\|x-y\|.$ (13)

######  Remark 3.

We note that Assumption 4 is trivially satisfied when the noise distribution
does not depend on the current iterate. For another illustration, consider
Gaussian random noise distribution from Example 1, for which we have:

$\displaystyle\Lambda(\lambda;x)-\Lambda(\lambda;y)$
$\displaystyle=\frac{1}{2}\lambda^{\top}(S(x)-S(y))\lambda$ (14)
$\displaystyle\leq\frac{1}{2}\|\lambda\|^{2}\|S(x)-S(y)\|.$ (15)

Comparing with the condition in (13), we see that (13) is satisfied when
entries of the covariance matrix $S$, as functions of $x$, are Lipschitz
continuous.

The assumption below will be used for the proof of the main result of the
paper, when the case of general convex functions is considered.

######  Assumption 5 (Sub-Gaussian noise).

There exists a constant $C_{1}>0$ such that, for each
$\lambda,\,x\in\mathbb{R}^{d}$

$\Lambda(\lambda;x)\leq C_{1}\frac{\|\lambda\|^{2}}{2}.$ (16)

######  Remark 4.

Assumption 5 means that the gradient noise has “light tails,” i.e., there
exist positive constants $c_{1},c_{2}$, such that the probability that the
magnitude of the norm of the noise vector is above $\epsilon$ is upper bounded
by $c_{1}\,e^{-c_{2}\epsilon^{2}}$, for any $\epsilon>0$. Clearly, a Gaussian
zero-mean multivariate distribution satisfies this property, and also any
noise distribution with compact support.

This assumption also ensures that, for each given $\lambda,$ the value of the
variance “proxy” $C_{1}$ cannot grow without bound as the domain of iterates
$x$ enlarges. For a Gaussian distribution, this means that the variance, as a
function of the current iterate should be uniformly bounded over the domain of
the iterates, which is a typical assumption in related works.

We also use the following implications of Assumption 5.

######  Proposition 1.

  1. 1.

There exists $C_{2}>0$ such that

$\mathbb{E}\left[\left.\exp\left(\frac{\|Z_{k}^{2}\|}{C_{2}}\right)\right|X_{k}\right]\leq
e.$ (17)

  2. 2.

For any $\nu\in[0,1/C_{2}]$ there holds

$M(\nu;X_{k})\leq\exp(\nu C_{2}).$ (18)

###### Proof.

The proof of part 1 can be derived by applying properties of sub-Gaussian
random variables to $\|Z_{k}\|;$ see, e.g., Proposition 2.5.2 in [25] and also
[15] for a treatment of sub-Gaussian random vectors.

To show part 2, fix $\nu\in[0,1/C_{2}].$ By Hölder’s inequality (applied for
$\text{``}p\text{''}=1/(\nu C_{2})\geq 1$)

$\displaystyle M(\nu;X_{k})$
$\displaystyle\leq\left(\mathbb{E}\left[\left.\exp\left(1/C_{2}\|Z_{k}\|^{2}\right)\right|X_{k}\right]\right)^{\nu
C_{2}}$ (19) $\displaystyle\leq\exp(\nu C_{2})$ (20)

where in the second inequality we used part 1. ∎

######  Remark 5.

When the distribution of $Z_{k}$ is Gaussian, zero mean and with covariance
matrix $\Sigma,$ and independent of the current iterate, we have

$\Lambda(\lambda)=\frac{1}{2}\lambda^{\top}\Sigma\lambda\leq\frac{1}{2}\sigma_{\max}^{2}\|\lambda\|^{2},$
(21)

where $\sigma_{\max}^{2}$ is the maximal eigenvalue of $\Sigma.$ Comparing
with Assumption 5, we see that condition (16) holds with
$C_{1}=\sigma_{\max}^{2}.$ It can also be shown that part 1. of Proposition 1
holds for $C_{2}\geq 2\sigma_{\max}^{2}.$

###  2.3 Key technical lemma

######  Definition 5.

The Fenchel-Legendre transform, or the conjugate, of a given function
$\Psi:\mathbb{R}^{d}\mapsto\mathbb{R}$ is defined by

$I(x)=\sup_{\lambda\in\mathbb{R}^{d}}x^{\top}\lambda-\Psi(\lambda),\>\>\mathrm{for\>\>}x\in\mathbb{R}^{d}.$
(22)

######  Lemma 4.

Let $\Psi_{k}$ be a sequence of log-moment generating functions associated to
a given sequence of measures
$\mu_{k}:\mathcal{B}(\mathbb{R}^{d})\mapsto[0,1].$ Suppose that, for each
$\lambda\in\mathbb{R}^{d},$ the following limit exists:

$\limsup_{k\rightarrow+\infty}\frac{1}{k}\Psi_{k}(k\lambda)\leq\Psi(\lambda).$
(23)

If $\Psi(\lambda)<\infty$ for each $\lambda\in\mathbb{R}^{d},$ then the
sequence $\mu_{k}$ satisfies the LDP upper bound with the rate function $I$
equal to the Fenchel-Legendre transform of $\Psi.$ If, in addition, (23) holds
as a limit and with equality, then the sequence of measures satisfies the LDP
with rate function $I.$

The second part of the lemma follows by the Gärtner-Ellis theorem. The first
part can be proven by similar arguments as in the proof of the upper bound of
the Gärtner-Ellis theorem; for details, see also the proof of Lemma 35 in [2].

##  3 Large deviations rates for SGD iterates $X_{k}$

###  3.1 Large deviations rates for $\|X_{k}-x^{\star}\|$

To derive the main result – the large deviations rate function for the SGD
sequence $X_{k},$ we first study large deviations properties of the sequence
$\|X_{k}-x^{\star}\|,$ $k=1,2,...$ For the latter, we first exploit the idea
from [12] to obtain a high probability bound for the (scaled) quantity
$\|X_{k}-x^{\star}\|^{2},$ via its moment generating function. We then use
this bound to derive a rate function (bound) for $\|X_{k}-x^{\star}\|.$ Since
our assumptions are distinct than those in [12] (e.g., the recursive form that
we work with here contains factors that require special treatment than the one
in [12], also we do not assume bounded noisy gradient, as is the case with the
proof available in [12]), we provide full proof details, see appendix.

######  Lemma 5.

For any $k,$ there holds

$\mathbb{P}\left(\|X_{k}-x^{\star}\|\geq\delta\right)\leq
ee^{-(k+k_{0})B\delta^{2}},$ (24)

where
$B=\min\\{\frac{1}{k_{0}\|X_{1}-x^{\star}\|},\frac{2a\mu-1}{4\max\\{C_{1},2C_{2}\\}a^{2}}\\}$
and $k_{0}=4a^{2}L^{2}/(2a\mu-1).$

######  Remark 6.

The preceding theorem establishes a large deviations upper bound for the
sequence of squared distance to solution iterates $\xi_{k},$ by exploiting
noise sub-Gaussianity. By its nature, this result is a rough characterization
of the large deviations rate function for the sequence $X_{k}.$ In addition to
being a result of independent interest, the utility consists in bounding the
tails of distribution $\mu_{k},$ as an enabling step towards deriving a fine,
close to exact rate function for the SGD iterates $X_{k},$ as the main
contribution of this paper. The latter is the subject of the next section.

###  3.2 Main result: Large deviations rates for $X_{k}$

We now present our result for general convex functions satisfying assumptions
from Section 2. The pillar of the analysis is the limit of the sequence of
log-moment generating functions $\log\Gamma_{k}$ of the SGD iterates.

######  Lemma 6.

Suppose that Assumptions 1-5 hold and that the stepsize is given by
$\alpha_{k}=a/(k+k_{0}).$ For any $\lambda\in\mathbb{R}^{d},$

$\limsup_{k\rightarrow+\infty}\frac{1}{k}\log\Gamma_{k}(k\lambda)\leq\overline{\Psi}(\lambda):=\Psi^{\star}(\lambda)+r(\lambda),$
(25)

where $\Psi^{\star}$ is defined by

$\Psi^{\star}(\lambda)=\int_{0}^{1}\Lambda(aQD(\theta)Q^{\top}\lambda;x^{\star})d\theta,$
(26)

where $H^{\star}=QDQ^{\top}$, $QQ^{\top}=I$,
$D=\mathrm{diag}\\{\rho_{1},...,d_{n}\\}$,
$D(\theta)=\mathrm{diag}\\{\theta^{a\rho_{1}-1},...,\theta^{ad_{n}-1}\\}$,
$r(\lambda)=\frac{4a^{2}\overline{\gamma}^{2}L_{\Lambda}}{B^{2}}\|\lambda\|^{4}+a\|\lambda\|\overline{h}\left(\frac{2\overline{\gamma}\|\lambda\|}{B}\right),$
and $\overline{\gamma}=\max\\{1,\sqrt{(1-a\mu)^{2}+a^{2}(L^{2}-\mu^{2})}\\}.$

The proof of Lemma 6 is given in section 5. Having the limit in (25), LDP
upper bound follows by Lemma 4.

######  Theorem 1.

Suppose that Assumptions 1-5 hold and that the stepsize is given by
$\alpha_{k}=a/(k+k_{0}).$ Then, the sequence of iterates $X_{k}$ satisfies the
LDP upper bound with rate function $\overline{I}$ given as the Fenchel-
Legendre transform of $\overline{\Psi}$ from Lemma 6, i.e., for any closed set
$F$:

$\limsup_{k\rightarrow+\infty}\frac{1}{k}\log\mathbb{P}\left(X_{k}\in
F\right)\leq-\inf_{x+x^{\star}\in F}\overline{I}(x).$ (27)

######  Remark 7.

The rate function $\overline{I}$ depends on the Hessian matrix at the
solution, $H(x^{\star})$. However, coarser exponential rate bounds can be
obtained by uniformly bounding the eigenvalues of $H(x^{\star}),$ as by our
assumptions they are all confined in the interval $[\mu,L].$ See Appendix D
for details.

###  3.3 Discussions and interpretations

####  3.3.1 Positivity of $\overline{I}$ and exponential decay

From the fact that $\Psi^{\star},\,r\geq 0,$ and that both $\Psi^{\star}$ and
$r$ are finite on $\mathbb{R}^{d},$ it can be shown that $\overline{I}\geq 0$
and that $\overline{I}$ is a good rate function. Specifically,
$\overline{I}(0)=0$ and $\overline{I}(x)>0$ for any $x\neq 0.$ Therefore, for
any closed set $F$ such that $x^{\star}\notin F,$ we have

$\inf_{x+x^{\star}\in F}\overline{I}(x)>0,$ (28)

that is, the exponent in (27) is strictly positive ensuring the exponential
decay of the probabilities $\mathbb{P}\left(X_{k}\in F\right)$. To illustrate
this in intuitive terms, we take as a special case the set
$F=B_{x^{\star}}^{\mathrm{c}}(\delta),$ for some $\delta>0.$ Then, the event
of interest becomes $\\{X_{k}\in F\\}=\\{\|X_{k}-x^{\star}\|\geq\delta\\}.$
Thus, for any $\delta>0,$ Theorem 1 implies that

$\limsup_{k\rightarrow+\infty}\frac{1}{k}\log\mathbb{P}\left(\|X_{k}-x^{\star}\|\geq\delta\right)\leq-R(\delta),$
(29)

where $R(\delta)=\inf_{\|x\|\geq\delta}I(x)>0.$

####  3.3.2 Remainder term $r$

Recalling Lemma 1, it is easy to see that $r(\lambda)=o(\|\lambda\|^{2}),$
i.e., $\lim_{\|\lambda\|\rightarrow 0}\frac{r(\lambda)}{\|\lambda\|^{2}}=0.$
Also, for a function $f$ that has Lipschitz Hessian, see Remark 2, the
residual function $r$ behaves roughly as $\sim\|\lambda\|^{3}$.

Further, for the special case when $f$ is quadratic, $\overline{h}(\delta)=0,$
and hence $r$ contains only the first term, and thus
$r(\lambda)\sim\|\lambda\|^{4}$. Similarly, when the noise distribution does
not depend on the current iterate, we have that $L_{\Lambda}=0,$ and hence
$r(\lambda)=o(\|\lambda\|^{2}).$ Finally, for the case when both of the
preceding conditions hold, the residual term is zero at all points: $r\equiv
0,$ and hence the rate function $\overline{I}=I^{\star},$ where $I^{\star}$ is
the Fenchel-Legendre transform of $\Psi^{\star}.$

####  3.3.3 Small deviations regime

When high precision estimates are sought, or equivalently, for small $\delta$
in (29), the candidate values of $\overline{I}$ in the minimization are very
close to $0.$ By the fact that the remainder term
$r(\lambda)=o(\|\lambda^{2}\|),$ it can be shown that, in the small deviations
regime, $\overline{I}$ is determined by $\Psi^{\star}$ only, i.e.,
$\overline{I}\approx I^{\star},$ and, also, its behaviour is dominantly
characterized by the noise variance.

###  3.4 LDP for quadratic functions

In this section we provide the full LDP for the case when $f$ is a quadratic
function. The proof of Theorem 2 is given in Appendix E.

######  Theorem 2.

Suppose that the objective function $f$ is quadratic, that Assumptions 2-3
hold, with the step size given by $\alpha_{k}=a/k$. Suppose also that the
noise distribution does not depend on the current iterate and that it has a
finite log-moment generating function $\Lambda$. Then, the sequence $X_{k}$
satisfies the large deviations principle with the rate function $I^{\star}$
given as the conjugate of $\Psi^{\star}$ defined in (26), with
$\Lambda(\cdot;x^{\star})$ replaced by $\Lambda.$

The rate function $I^{\star}$ depends on the distribution of $Z_{k}$ and fully
captures all moments of this distribution. In particular, for non-Gaussian
distributions, it captures exactly the dependence not only on the variance,
but also on higher order moments.

######  Remark 8.

We note that, in contrast with Theorem 1, for Theorem 2 the conditional
distribution of $Z_{k}$ can be arbitrary, as long as $\Lambda$ is finite. In
particular, it allows for distributions that are not light-tailed, such as
Laplacian.

######  Remark 9.

Recalling the discussion from subsection 3.3.2, we see that the upper bound
rate function from Theorem 1 and the rate function from Theorem 2 match, hence
showing that the bound in Theorem 1 is tight.

##  4 Gaussian noise: analytical characterization of the rate function

If the noise $Z_{k}$ has a Gaussian distribution with mean value zero and
covariance matrix $\Sigma$, then $\Psi^{\star}$ is computed by

$\Psi^{\star}(\lambda)=\frac{a^{2}}{2}\int_{0}^{1}\lambda^{\top}QD(\theta)Q^{\top}\Sigma
QD(\theta)Q^{\top}\lambda d\theta.$ (30)

To simplify the notation, let $S=Q^{\top}\Sigma Q,$ and
$M(\theta)=D(\theta)SD(\theta).$ It is easy to verify that
$M_{ij}(\theta)=S_{ij}\theta^{a(\rho_{i}+\rho_{j})-2},$ for any $i,j=1,...,d,$
and thus $\int_{0}^{1}M_{ij}(\theta)d\theta=S_{ij}/(a(\rho_{i}+\rho_{j})-1).$
Hence, we obtain the following closed-form expression for $\Psi^{\star}:$

$\Psi^{\star}(\lambda)=\frac{a^{2}}{2}\lambda^{\top}QS^{\star}Q\lambda,$ (31)

where $S_{ij}^{\star}=S_{ij}/(a(\rho_{i}+\rho_{j})-1),$ for $i,j=1,...,d.$

Recalling the Definition 5, it can be shown that the Fenchel-Legendre
transform $I^{\star}$ of $\Psi^{\star}$ is given by

$I^{\star}(z)=\frac{1}{2a^{2}}z^{\top}Q^{\top}{S^{\star}}^{-1}Qz.$ (32)

To obtain further intuition about the rate function $I^{\star}$, we consider
the special case when the Hessian matrix $H^{\star}$ and the covariance matrix
$\Sigma$ share the same eigenspace (given by the columns of the matrix $Q$).
Intuitively, the latter means that the orientation of the quadratic
approximation of $f$ at the origin is aligned with the gradient noise
distribution in each of the axes. In this case, it follows that
$S=Q^{\top}\Sigma Q$ is diagonal with $S_{ii}=\sigma_{ii}^{2},$ where
$\sigma_{ii}^{2}$ is the $i$-th eigenvalue of $\Sigma$ (i.e., the eigenvalue
of $\Sigma$ corresponding to its eigenvector given by the $i$-th column of
matrix $Q$). It follows that $S^{\star}$ is also diagonal with
$S_{ii}^{\star}=\sigma_{ii}^{2}/(2a\rho_{i}-1).$ Thus, the following neat
expression for the rate function $I^{\star}$ emerges:

$I(z)=\frac{1}{2a^{2}}z^{\top}Q^{\top}\mathrm{diag}\left(\frac{2a\rho_{1}-1}{\sigma_{11}^{2}},...,\frac{2a\rho_{d}-1}{\sigma_{dd}^{2}}\right)Qz.$
(33)

###  4.1 Decay rates with $l_{2}$ balls

We consider the case when in the large deviations event of interest
$\\{X_{k}\in F\\}$ the set $F$ is given as the complement of an $l_{2}$ ball
around the solution $x^{\star}:$ $F=B^{\mathrm{c}}_{x^{\star}}(\delta),$ i.e.,
$\\{X_{k}\in F\\}=\\{\|X_{k}-x^{\star}\|\geq\delta\|.$ Assuming that the
residual is zero (see the result for quadratic functions in Section 3.4), by
Theorem 1, we have

$\limsup_{k\rightarrow+\infty}\frac{1}{k}\log\mathbb{P}\left(\left\|X_{k}-x^{\star}\right\|\geq\delta\right)\leq\inf_{\|z\|\geq\delta}I(z)=:\mathbf{I}(B^{\mathrm{c}}_{x^{\star}}(\delta)).$
(34)

For the Gaussian noise assumed in this section, we have:

$\displaystyle\mathbf{I}(B^{\mathrm{c}}_{x^{\star}}(\delta))$
$\displaystyle=\inf_{\|z\|\geq\delta}\frac{1}{2a^{2}}z^{\top}Q^{\top}{S^{\star}}^{-1}Qz$
$\displaystyle=\frac{\delta^{2}}{2a^{2}}\inf_{\|w\|\geq
1}w^{\top}Q^{\top}{S^{\star}}^{-1}Qw$
$\displaystyle=\frac{\delta^{2}}{2a^{2}}\frac{1}{\lambda_{\max}(S^{\star})},$
(35)

where $\lambda_{\max}(S^{\star})$ is the largest eigenvalue of the matrix
$S^{\star}.$ Hence, to find the value of the exponent $\mathbf{I}$ for any
given ball-shaped set, it suffices to find (once) the maximal eigenvalue of
$S^{\star},$ and the exponent $\mathbf{I}$ would be easily computed by the
quadratic function (4.1).

We close the analysis with a particularly elegant solution for the special
case when $H^{\star}$ and $\Sigma$ are axes-aligned. As detailed at the
beginning of the section, in the latter case, $S^{\star}$ is diagonal, with
$S_{ii}^{\star}=\sigma_{ii}^{2}/(2a\rho_{i}-1),$ and the rate function is
given by (33). Thus, to find the maximal eigenvalue of $S^{\star}$ reduces to
finding the index $i$ for which $\frac{\sigma^{2}_{ii}}{2a\rho_{i}-1}$ is
highest, or, equivalently, $\frac{2a\rho_{i}-1}{\sigma^{2}_{ii}}$ the lowest,
which then yields:

$\displaystyle\mathbf{I}(B^{\mathrm{c}}_{x^{\star}}(\delta))=\frac{\delta^{2}}{2a^{2}}\min\\{\frac{2a\rho_{i}-1}{\sigma^{2}_{ii}}:i=1,...d\\},$
(36)

where, we recall, $\rho_{i}$ is the $i$-th eigenvalue of $H^{\star}.$ What the
expression above is saying is that, in order to find the exponential decay
rate for an $l_{2}$ ball, we should search for the direction $i$ in which the
value $\frac{\sigma^{2}_{ii}}{2\rho_{i}-1}$ is highest. In a sense, the latter
quantity can be thought of as the effective noise variance, capturing the
interplay between the noise distribution and the shape of the function at the
solution. Specifically, if along the direction where the noise variance is
highest, say $i^{\star},$ the function has a high curvature (i.e., large
$\rho_{i^{\star}}$), this will effectively alleviate the effects of noise and
increase the rate function, in comparison to the case when the curvature along
$i$ is lower, and therefore result in faster convergence.

Finally, when the noise is isotropic, i.e., such that
$\sigma_{ii}^{2}=\sigma^{2},$ for all $i,$ exploiting the fact that the
spectrum of $H^{\star}$ lies inside the interval $[\mu,L],$ the rate function
is found by:

$\displaystyle\mathbf{I}(B^{\mathrm{c}}_{x^{\star}}(\delta))=\frac{\delta^{2}}{2a^{2}}\frac{2a\mu-1}{\sigma^{2}}.$
(37)

###  4.2 Comparison with the rate from Lemma 5

We now compare the rate function bounds obtained from Lemma 5 and Theorem 1.
To gain deeper insights, we will assume that the residual term $r$ equals zero
(compare with Section 3.4). We also assume that the noise is Gaussian and
axes-aligned with the matrix $H^{\star}$ (see the preceding subsection). The
exponent $B$ from 5 can be upper bounded by333The dependence in $B$ on $X_{1}$
in Lemma 5 seems to be an artifact of the conducted proof method, rather than
an essential property of the exponential rate that Lemma 5 pursues. Hence, for
unbiased comparison, we omit this factor in the analysis of the rate $B$.

$B\leq\frac{2a\mu-1}{4{\sigma_{\max}}^{2}a^{2}},$

where we exploited the fact that, for Gaussian noise,
$C_{1}=\sigma_{\max}^{2},$ see Remark 5. Hence, for an $l_{2}$ ball of radius
$\delta,$ the exponent that Lemma 5 provides is bounded by

$B\delta^{2}\leq\frac{\delta^{2}}{4a^{2}}\frac{2a\mu-1}{\sigma_{\max}^{2}}.$
(38)

The counterpart obtained from Theorem 1 is given by expression (36). To show
direct comparison with (38), we further upper bound this value by decoupling
the minimization over $i:$

$\displaystyle\mathbf{I}(B^{\mathrm{c}}_{x^{\star}}(\delta))$
$\displaystyle=\frac{\delta^{2}}{2a^{2}}\min\\{\frac{2a\rho_{i}-1}{\sigma^{2}_{ii}}:i=1,...d\\}$
$\displaystyle\geq\frac{\delta^{2}}{2a^{2}}\frac{\min\\{2a\rho_{i}-1:i=1,..,d\\}}{\max\\{\sigma_{ii}^{2}:i=1,..,d\\}}$
$\displaystyle=\frac{\delta^{2}}{2a^{2}}\frac{2a\mu-1}{\sigma_{\max}^{2}}.$
(39)

Comparing with (38) (and ignoring the scaling constant $2$), the following
important point can be noted: on intuitive level, the derivation of the rate
$B$ is equivalent to that of decoupling the effects of the noise distribution
and the shape of the function $f$ at the origin with the rate $I^{\star}$.
Hence, in contrast with $I^{\star},$ the rate $B$ is oblivious to the
interplay between these two quantities – from a purely technical perspective,
this distinction is a consequence of relying on recursions on the iterates’
distance to solution, $\|X_{k}-x^{\star}\|,$ as opposed to working directly
with the iterates $X_{k},$ as is the case in the proof of Theorem 1.

##  5 Proof of Lemma 6

This section provides the main elements of the proof of the limit in (25); the
proofs of omitted results can be found in Appendix C. Fix
$\lambda\in\mathbb{R}^{d}.$ Fix $k\geq 1.$ Define $\eta_{l}=B_{k,l}\eta_{k},$
$B_{k,l}=(I-\alpha_{l}H^{\star})\cdots(I-\alpha_{k}H^{\star}),$
$\eta_{k}=k\lambda.$ By Lemma 2,

$\|\eta_{l}\|\leq
k\left(\frac{l+k_{0}}{k+k_{0}+1}\right)^{a\mu}\|\lambda\|\leq(l+k_{0})\|\lambda\|.$
(40)

For an arbitrary $l\leq k,$ there holds

$\displaystyle\Gamma_{l+1}(\eta_{l})$
$\displaystyle=\mathbb{E}\left[\exp\left(\eta_{l}^{\top}(X_{l+1}-x^{\star})\right)\right]$
$\displaystyle=\mathbb{E}\left[\mathbb{E}\left[\exp\left(\eta_{l}^{\top}(X_{l}-\alpha_{l}g(X_{l})+\alpha_{l}Z_{l}-x^{\star})\right)|X_{l}\right]\right]$
$\displaystyle=\mathbb{E}\left[\exp\left(\Lambda(\alpha_{l}\eta_{l};X_{l})+\eta_{l}^{\top}(X_{l}-\alpha_{l}g(X_{l})-x^{\star})\right)\right]$
$\displaystyle=\int_{x\in\mathbb{R}^{d}}\Gamma_{l+1|l}(\eta_{l};x)\mu_{l}(dx),$
(41)

where $\Gamma_{l+1|l}(\cdot;x)$ denotes the conditional moment generating
function of $X_{l+1},$ given $X_{l}=x.$ We now fix $\delta>0$ (the exact value
to be chosen later) and split the analysis in two cases: 1)
$\mathcal{A}_{l,\delta}=\left\\{X_{l}\in B_{x^{\star}}(\delta)\right\\};$ and
2) $\mathcal{A}^{\mathrm{c}}_{l,\delta}=\left\\{X_{l}\in
B^{\mathrm{c}}_{x^{\star}}(\delta)\right\\}.$

Introduce

$\displaystyle\Gamma_{l+1|\mathcal{A}_{l,\delta}}(\eta_{l})$
$\displaystyle:=\mathbb{E}\left[1_{\|X_{l}-x^{\star}\|\leq\delta}\Gamma_{l+1|l}(\eta_{l};X_{l})\right]$
$\displaystyle=\int_{\|x-x^{\star}\|\leq\delta}\Gamma_{l+1|l}(\eta_{l};x)\mu_{l}(dx)$
(42) $\displaystyle\Gamma_{l+1|\mathcal{A}^{\mathrm{c}}_{l,\delta}}(\eta_{l})$
$\displaystyle:=\mathbb{E}\left[1_{\|X_{l}-x^{\star}\|>\delta}\Gamma_{l+1|l}(\eta_{l};X_{l})\right]$
$\displaystyle=\int_{\|x-x^{\star}\|>\delta}\Gamma_{l+1|l}(\eta_{l};x)\mu_{l}(dx);$
(43)

note that

$\Gamma_{l+1}(\eta_{l})=\Gamma_{l+1|A_{l,\delta}}(\eta_{l})+\Gamma_{l+1|A^{\mathrm{c}}_{l,\delta}}(\eta_{l}).$
(44)

_Case 1: $x\in\mathcal{A}_{l,\delta}$._ Fix $x\in\mathbb{R}^{d}$ such that
$\|x\|\leq\delta.$ We have:

$\displaystyle\Gamma_{l+1|l}(\eta_{l};x)$
$\displaystyle=\exp\left(\Lambda(\alpha_{l}\eta_{l};x)+\eta_{l}^{\top}(x-\alpha_{l}g(x)-x^{\star})\right)$
$\displaystyle=\exp\left(\Lambda(\alpha_{l}\eta_{l};x)+\eta_{l}^{\top}((I-\alpha_{l}H^{\star})(x-x^{\star})-\alpha_{l}h(x)\right)$
$\displaystyle\leq\exp\left(\Lambda(\alpha_{l}\eta_{l};x^{\star})+L_{\Lambda}\|\eta_{l}^{2}\|\|x-x^{\star}\|+\alpha_{l}\|\eta_{l}\|\|h(x)\|\right)$
$\displaystyle\,\,\,\,\,\times\exp\left(\eta_{l}^{\top}((I-\alpha_{l}H^{\star})(x-x^{\star})\right)$
(45)
$\displaystyle\leq\exp\left(\Lambda(\alpha_{l}\eta_{l};x^{\star})+L_{\Lambda}\alpha_{l}^{2}\|\eta_{l}^{2}\|\delta^{2}+\alpha_{l}\|\eta_{l}\|\overline{h}(\delta)+\eta_{l-1}^{\top}(x-x^{\star})\right),$
(46)

where in (45) we used Lipschitz continuity of $\Lambda$ in $x$, Assumption 4,
and in (46) we used the fact that $\|x-x^{\star}\|\leq\delta$. It follows that

$\Gamma_{l+1|\mathcal{A}_{\delta}}(\eta_{l})\leq\exp\left(\Lambda(\alpha_{l}\eta_{l};x^{\star})+r_{0}(\lambda,\delta)\right)\Gamma_{l}(\eta_{l-1}),$
(47)

where
$r_{0}(\lambda,\delta)=L_{\Lambda}a^{2}\|\lambda\|^{2}\delta^{2}+a\|\lambda\|\overline{h}(\delta).$

_Case 2: $x\in\mathcal{A}^{\mathrm{c}}_{l,\delta}$_. By strong convexity and
Lipschitz smoothness of $f$ in Assumption 1, for each $l\geq 1,$ the following
holds:

$\displaystyle\|X_{l}-g(X_{l})-x^{\star}\|$
$\displaystyle\leq\gamma_{l}\|X_{l}-x^{\star}\|,$ (48)
$\displaystyle\leq\overline{\gamma}\|X_{l}-x^{\star}\|$ (49)

where $\gamma_{l}=(1-2\alpha_{l}\mu+\alpha_{l}^{2}L^{2})^{1/2}$, see Appendix
A for the proof, and $\overline{\gamma}=\sup\\{\gamma_{l}:l=1,2,...\\};$ it is
easy to verify that
$\overline{\gamma}=\max\\{1,\sqrt{(1-a\mu)^{2}+a^{2}(L^{2}-a^{2})}\\}.$

For an arbitrary $x\in\mathbb{R}^{d},$ we have:

$\displaystyle\Gamma_{l+1|l}(\eta_{l};x)$
$\displaystyle=\exp\left(\Lambda(\alpha_{l}\eta_{l};x)+\eta_{l}^{\top}(x-\alpha_{l}g(x)-x^{\star})\right)$
$\displaystyle\leq\exp\left(\frac{C_{1}\alpha_{l}^{2}\|\eta_{l}\|^{2}}{2}\right)\exp\left(\overline{\gamma}\|\eta_{l}\|\|x-x^{\star}\|\right),$
(50)
$\displaystyle\leq\exp\left(\frac{C_{1}a^{2}\|\lambda\|^{2}}{2}\right)\exp\left(\overline{\gamma}(l+k_{0})\|\lambda\|\|x-x^{\star}\|\right),$
(51)

where in (50) we used the assumption that $Z_{k}$ is sub-Gaussian, Assumption
5, for the first term, together with (48) and Cauchy-Schwartz, for the second
term, while in (51) we exploited (40). Recalling the induced measure
$\nu_{l},$ we now have

$\displaystyle\Gamma_{l+1|\mathcal{A}^{\mathrm{c}}_{l,\delta}}(\eta_{l})$
$\displaystyle\leq\exp\left(\frac{C_{1}a^{2}\|\lambda\|^{2}}{2}\right)\int_{z\geq\delta}e^{(l+k_{0})\overline{\gamma}\|\lambda\|z}\nu_{l}(dz).$
(52)

The idea of analysing the “tail” term
$\Gamma_{l+1|\mathcal{A}^{\mathrm{c}}_{l,\delta}}(\eta_{l})$ is the following:
by Theorem 5, we know that the probability density $\nu_{l}$ at a given point
$z$ behaves roughly as $e^{-(l+k_{0})Bz^{2}}.$ If $\delta$ is sufficiently
large, then, for all $z\geq\delta,$ the negative exponential rate of the
measure $\nu_{l}(z)$ is in absolute terms higher than the exponent
$(l+k_{0})\overline{\gamma}\|\lambda\|z.$ Integrating by parts, we obtain that
for $\delta=\frac{2\overline{\gamma}\|\lambda\|}{B},$ the integral on the
right hand-side of (52) is upper bounded by a constant $K.$ Thus:

$\Gamma_{l+1|\mathcal{A}^{\mathrm{c}}_{l,\delta}}(\eta_{l})\leq
K\exp\left(\frac{C_{1}a^{2}\|\lambda\|^{2}}{2}\right).$ (53)

Combining with (47) and recalling (44),

$\Gamma_{l+1}(\eta_{l})\leq\exp\left(\Lambda(\alpha_{l}\eta_{l};x^{\star})+r(\lambda)\right)\Gamma_{l}(\eta_{l-1})\\\
+K\exp\left(\frac{C_{1}a^{2}\|\lambda\|^{2}}{2}\right),$

where
$r(\lambda)=r_{0}\left(\lambda,\frac{2\overline{\gamma}\|\lambda\|}{B}\right).$
Iterating the preceding recursion, where we exploit the nonnegativity of
$\Lambda,$ property 3. from Lemma 3, we obtain:

$\displaystyle\Gamma_{k+1}(k\lambda)\leq\exp\left(\sum_{l=1}^{k}\left(\Lambda(\alpha_{l}\eta_{l};x^{\star})+r(\lambda)\right)\right)\Gamma_{1}(\alpha_{1}\eta_{1})$
$\displaystyle+K\exp\left(\frac{C_{1}a^{2}\|\lambda\|^{2}}{2}\right)\sum_{l=1}^{k}e^{\sum_{j=l}^{k}\left(\Lambda(\alpha_{j}\eta_{j};x^{\star})+r(\lambda)\right)}$
$\displaystyle\leq(k+1)K\exp\left(\frac{C_{1}a^{2}\|\lambda\|^{2}}{2}+a\|\lambda\|\|X_{1}-x^{\star}\|\right)$
$\displaystyle\;\;\;\;\;\;\;\times
e^{\sum_{l=1}^{k}\left(\Lambda(\alpha_{l}\eta_{l};x^{\star})+r(\lambda)\right)}.$
(54)

Taking the limit, dividing by $k,$ and taking the $\limsup$

$\displaystyle\limsup_{k\rightarrow+\infty}\frac{1}{k}\log\Gamma_{k+1}(k\lambda)\leq$
$\displaystyle\;\;\;\;\;\;\;\;\;\;\;\;\;\;\;r(\lambda)+\limsup_{k\rightarrow+\infty}\frac{1}{k}\sum_{l=1}^{k}\Lambda(\alpha_{l}\eta_{l};x^{\star}).$
(55)

Finally, it can be shown that

$\lim_{k\rightarrow+\infty}\frac{1}{k}\sum_{l=1}^{k}\Lambda(\alpha_{l}\eta_{l};x^{\star})=\int_{0}^{1}\Lambda(aQD(\theta)Q^{\top}\lambda;x^{\star})d\theta.$
(56)

The proof of (56) is provided in Appendix C. This completes the proof of Lemma
6.

##  6 Conclusions

We developed large deviations analysis for the stochastic gradient descent
(SGD) method, when the objective function is smooth and strongly convex. For
(strongly convex) quadratic costs, we establish the full large deviations
principle. That is, we derive the exact exponential rate of decay of the
probability that the iterate sequence generated by SGD stays within an
arbitrary set that is away from the problem solution. This is achieved for a
very general class of gradient noises, that may be iteration-dependent and are
required to have a finite log-moment generating function. For generic costs,
we derive a tight large deviations upper bound that, up to higher order terms,
matches the exact rate derived for the quadratics.

## References

  * [1] R.. Bahadur  “On the Asymptotic Efficiency of Tests and Estimates”  In _Sankhya: The Indian Journal of Statistics, 1933-1960_ 22.3/4 Springer, 1960, pp. 229–252  URL: http://www.jstor.org/stable/25048458
  * [2] Dragana Bajovic  “Inaccuracy rates for distributed inference over random networks with applications to social learning”  arXiv, 2022  DOI: 10.48550/ARXIV.2208.05236
  * [3] D. Bajović, D. Jakovetić, J… Moura, J. Xavier and B. Sinopoli  “Large Deviations Performance of Consensus+Innovations Distributed Detection with Non-Gaussian Observations”  In _IEEE Transactions on Signal Processing_ 60.11, 2012, pp. 5987–6002
  * [4] J.. Bucklew  “Large Deviations Techniques in Decision, Simulation and Estimation”  New York: Wiley, 1990
  * [5] D., D., J., B. and J… Moura  “Distributed Detection via Gaussian Running Consensus: Large Deviations Asymptotic Analysis”  In _IEEE Transactions on Signal Processing_ 59.9, 2011, pp. 4381–4396
  * [6] Damek Davis, Dmitriy Drusvyatskiy, Lin Xiao and Junyu Zhang  “From Low Probability to High Confidence in Stochastic Convex Optimization.”  In _J. Mach. Learn. Res._ 22, 2021, pp. 49–1
  * [7] A. Dembo and O. Zeitouni  “Large Deviations Techniques and Applications”  Boston, MA: JonesBarlett, 1993
  * [8] S. Ghadimi and G. Lan  “Optimal stochastic approximation algorithms for strongly convex stochastic composite optimization I: A generic algorithmic framework”  In _SIAM J. Optim._ 22.4, 2012, pp. 1469–1492
  * [9] S. Ghadimi and G. Lan  “Optimal stochastic approximation algorithms for strongly convex stochastic composite optimization, II: Shrinking procedures and optimal algorithms”  In _SIAM J. Optim._ 23.4, 2013, pp. 2061–2089
  * [10] Eduard Gorbunov, Marina Danilova and Alexander Gasnikov  “Stochastic optimization with heavy-tailed noise via accelerated gradient clipping”  In _arXiv preprint arXiv:2005.10785_ , 2020
  * [11] Eduard Gorbunov, Filip Hanzely and Peter Richtárik  “A unified theory of SGD: Variance reduction, sampling, quantization and coordinate descent”  In _International Conference on Artificial Intelligence and Statistics_ , 2020, pp. 680–690  PMLR
  * [12] Nicholas J.. Harvey, Christopher Liaw, Yaniv Plan and Sikander Randhawa  “Tight analyses for non-smooth stochastic gradient descent”  In _32nd Annual Conference on Learning Theory_ 99, 2019, pp. 1–35  PMLR
  * [13] Ping Hu, Virginia Bordignon, Stefan Vlaski and Ali H. Sayed  “Optimal Aggregation Strategies for Social Learning over Graphs”  arXiv, 2022  DOI: 10.48550/ARXIV.2203.07065
  * [14] Wenqing Hu, Chris Junchi Li, Lei Li and Jian-Guo Liu  “On the diffusion approximation of nonconvex stochastic gradient descent”  In _Annals of Mathematical Sciences and Applications_ 4.1, 2019
  * [15] Chi Jin, Praneeth Netrapalli, Rong Ge, Sham M. Kakade and Michael I. Jordan  “A Short Note on Concentration Inequalities for Random Vectors with SubGaussian Norm”  arXiv, 2019  DOI: 10.48550/ARXIV.1902.03736
  * [16] A. Juditsky, A. Nazin, A. Nemirovsky and A. Tsybakov  “Algorithms of robust stochastic optimization based on mirror descent method”  In _arXiv:1907.02707_ , 2019
  * [17] Lihua Lei and Michael I Jordan  “On the adaptivity of stochastic gradient-based optimization”  In _SIAM Journal on Optimization_ 30.2 SIAM, 2020, pp. 1473–1500
  * [18] S. Marano and A.. Sayed  “Detection under one-bit messaging over adaptive networks”  In _IEEE Trans. Information Theory_ 65.10, 2019, pp. 6519–6538
  * [19] V. Matta, P. Braca, S. Marano and A.. Sayed  “Diffusion-based adaptive distributed detection: Steady-state performance in the slow adaptation regime”  In _IEEE Trans. Information Theory_ 62.8, 2016, pp. 4710–4732
  * [20] V. Matta, P. Braca, S. Marano and A.. Sayed  “Distributed detection over adaptive networks: Refined asymptotics and the role of connectivity”  In _IEEE Trans. Signal and Information Processing over Networks_ 2.4, 2016, pp. 442–460
  * [21] Feng Niu, Benjamin Recht, Christopher Ré and Stephen J Wright  “Hogwild!: A lock-free approach to parallelizing stochastic gradient descent”  In _arXiv preprint arXiv:1106.5730_ , 2011
  * [22] Adam Shwartz and Alan Weiss  “Large Deviations for Performance Analysis: Queues, Communications, and Computing”  New York: ChapmanHall, 1995
  * [23] T.. and J..  “Elements of Information Theory”  New York: John WileySons, 1991
  * [24] Hugo Touchette  “The large deviation approach to statistical mechanics”  In _Physics Reports_ 478.1, 2009, pp. 1–69
  * [25] Roman Vershynin  “High-Dimensional Probability: An Introduction with Applications in Data Science”, Cambridge Series in Statistical and Probabilistic Mathematics  Cambridge University Press, 2018  DOI: 10.1017/9781108231596
  * [26] Michael Woodroofe  “Normal Approximation and Large Deviations for the Robbins-Monro Process”  In _Z. Wahrscheinlichkeitstheorie verw. Geb._ 21, 1972, pp. 329–338

## Appendix A.

###### Proof of recursion 2.1.

For any $k\geq 1,$ we have:

$\displaystyle\|X_{k+1}-x^{\star}\|$
$\displaystyle=\|X_{k}-\alpha_{k}g(X_{k})+\alpha_{k}Z_{k}-x^{\star}\|$
$\displaystyle=\|X_{k}-x^{\star}\|^{2}-2\alpha_{k}(X_{k}-x^{\star})^{\top}(g(X_{k})-Z_{k})$
$\displaystyle+\alpha_{k}^{2}\|g(X_{k})-Z_{k}\|^{2}$
$\displaystyle\leq(1-2\alpha_{k}\mu)\xi_{k}+2\alpha_{k}(X_{k}-x^{\star})^{\top}Z_{k}$
$\displaystyle+2\alpha_{k}^{2}\|g(X_{k})\|+2\alpha_{k}^{2}\|Z_{k}\|^{2}$
$\displaystyle\leq\left(1-2\alpha_{k}\mu+2\alpha_{k}^{2}L^{2}\right)\xi_{k}+2\alpha_{k}(X_{k}-x^{\star})^{\top}Z_{k}$
$\displaystyle+2\alpha_{k}^{2}\|Z_{k}\|^{2},$ (57)

where the first inequality follows from the strong convexity of $f,$
Assumption 1, and the fact that, for $a,b\in\mathbb{R}^{d},$ $\|a-b\|^{2}\leq
2\|a\|^{2}+2\|b\|^{2},$ and the second inequality follows from the Lipschitz
smoothness of $f$, Assumption 1. ∎

###### Proof of Lemma 2.

Fix $l$ and $k$ where $1\leq l\leq k.$ Fix $u,\,v\geq 0.$ From the upper and
the lower Darboux sum for the logarithmic function applied to the interval
$[l,k]$, we obtain:

$\log\frac{k+1}{l}\leq\frac{1}{l}+\ldots+\frac{1}{k}\leq\log\frac{k}{l-1}.$
(58)

For the $2$-sum we use the following simple bound $1/l^{2}\leq
1/(l(l-1))=1/(l-1)-1/l$ to obtain:

$\frac{1}{l^{2}}+\ldots+\frac{1}{k^{2}}\leq\frac{1}{l-1}-\frac{1}{l}+\ldots+\frac{1}{k-1}-\frac{1}{k}\leq\frac{1}{l-1}.$
(59)

To prove part 1, we use that $1+x\leq e^{x}$ applied to each of the terms in
the product $\beta_{k,l},$ together with the left hand-side inequality of (58)
and the right hand-side inequality of (59):

$\displaystyle\beta_{k,l}(u,v)$ $\displaystyle\leq
e^{-au\sum_{j=l}^{k}\frac{1}{j+b}+a^{2}v\sum_{j=l}^{k}\frac{1}{(j+b)^{2}}}$
$\displaystyle\leq
e^{-au\log\left(\frac{k+b+1}{l+b}\right)+\frac{a^{2}v}{l+b-1}}$
$\displaystyle=\left(\frac{l+b}{k+b+1}\right)^{au}e^{\frac{a^{2}v}{l+b-1}}.$
(60)

To prove part 2, we first note that, since $v\geq 0,$ there holds
$\beta_{k,l}(u,v)\geq\beta_{k,l}(u,0),$ i.e.,
$\beta_{k,l}(u,v)\geq(1-\alpha_{k}u)\cdots(1-\alpha_{l}u).$ We now use that,
for $x\leq\frac{2}{5},$ $1-x\geq e^{-x-x^{2}}:$

$\displaystyle\beta_{k,l}(u,v)$ $\displaystyle\geq
e^{-au\sum_{j=l}^{k}\frac{1}{j+b}-a^{2}u^{2}\sum_{j=l}^{k}\frac{1}{(j+b)^{2}}}$
$\displaystyle\geq\left(\frac{l+b-1}{k+b}\right)^{au}e^{-\frac{a^{2}u^{2}}{l+b-1}}.$
(61)

This completes the proof of the lemma. ∎

###### Proof of (48).

Here we prove an alternative recursion on $\|X_{k}-x^{\star}\|$, used within
the proof of Lemma 6. Specifically, we show that, for any $k,$

$\|X_{k+1}-x^{\star}\|\leq\gamma_{k}\|X_{k}-x^{\star}\|+\alpha_{k}\|Z_{k}\|,$
(62)

where, we recall,
$\gamma_{k}=\left(1-2\alpha_{k}\mu+\alpha_{k}^{2}L^{2}\right)^{1/2}.$

From the triangle inequality applied to the Euclidean norm,

$\displaystyle\|X_{k+1}-x^{\star}\|$
$\displaystyle=\|X_{k}-\alpha_{k}g(X_{k})+\alpha_{k}Z_{k}-x^{\star}\|$
$\displaystyle\leq\|X_{k}-\alpha_{k}g(X_{k})-x^{\star}\|+\alpha_{k}\|Z_{k}\|$
(63)

Exploiting $L$-smoothness and $\mu$\- convexity of $f$ for the second term:

$\displaystyle\|X_{k}-\alpha_{k}g(X_{k})-x^{\star}\|^{2}\leq$
$\displaystyle\|X_{k}-x^{\star}\|^{2}-2\alpha_{k}(X_{k}-x^{\star})^{\top}g(X_{k})+\alpha_{k}^{2}\|g(X_{k})\|$
$\displaystyle\leq\|X_{k}-x^{\star}\|^{2}-2\alpha_{k}\mu\|X_{k}-x^{\star}\|^{2}+\alpha_{k}^{2}L^{2}\|X_{k}-x^{\star}\|^{2}$
$\displaystyle=\gamma_{k}^{2}\|X_{k}-x^{\star}\|^{2}.$ (64)

Taking the square root and replacing in (Proof of (48).) yields (62). ∎

## Appendix B.

###### Proof of Lemma 5..

First, we transform the recursion in (2.1) by defining
$Y_{k+1}=(k+k_{0})\|X_{k+1}-x^{\star}\|^{2},$ to obtain:

$Y_{k+1}\leq
a_{k}Y_{k}-b_{k}\sqrt{k+k_{0}-1}(X_{k}-x^{\star})^{\top}Z_{k}+c_{k}\|Z_{k}\|^{2},$
(65)

where

$\displaystyle a_{k}$
$\displaystyle=\frac{k+k_{0}}{k+k_{0}-1}(1-2\alpha_{k}\mu+2\alpha_{k}^{2}L^{2})$
(66) $\displaystyle b_{k}$ $\displaystyle=\frac{a}{\sqrt{k+k_{0}-1}}$ (67)
$\displaystyle c_{k}$ $\displaystyle=\frac{a^{2}}{k+k_{0}}.$ (68)

The key technical result behind Lemma 5 is the following upper bound on the
tail probability of the $Y_{k}$ iterates:

$\mathbb{P}\left(Y_{k}\geq\epsilon\right)\leq ee^{-B\epsilon},$ (69)

which holds for each $k\geq 1,$ and $\epsilon\geq 0.$ The result of Lemma 5
directly follows from (69) by taking $\epsilon_{k}=k\delta^{2},$ for each $k.$

Thus, in the remainder of the proof we focus on proving (24). It can be easily
verified that, for each $k,$

$a_{k}=1-\frac{2a\mu-1}{k+k_{0}-1}\left(1-\frac{2a^{2}L^{2}}{(2a\mu-1)(k+k_{0}-1)}\right).$
(70)

Recalling Assumption 2 and the value of $k_{0},$ we see that the above
quantity is smaller than $1$ for each $k.$

Denote by $\Phi_{k}$ the moment generating function of $Y_{k},$ and by
$\Phi_{k+1|k}(\cdot;X_{k})$ the moment generating function of $Y_{k}$
conditioned on $X_{k}:$

$\displaystyle\Phi_{k}(\nu)$ $\displaystyle:=\mathbb{E}\left[\exp(\nu
Y_{k})\right]$ (71) $\displaystyle\Phi_{k+1|k}(\nu;X_{k})$
$\displaystyle:=\mathbb{E}\left[\left.\exp(\nu Y_{k})\right|X_{k}\right],$
(72)

for $\nu\in\mathbb{R};$ note that
$\Phi_{k+1}(\nu)=\mathbb{E}\left[\Phi_{k+1|k}(\nu;X_{k})\right],$ for each
$\nu\in\mathbb{R}.$ From the recursion (2.1), we have:

$\displaystyle\Phi_{k+1|k}(\nu;X_{k})=$ $\displaystyle\exp(a_{k}\nu
Y_{k})\mathbb{E}\left[\left.\exp(-b_{k}\sqrt{k+k_{0}-1}(X_{k}-x^{\star})^{\top}Z_{k}+c_{k}\|Z_{k}\|^{2})\right|X_{k}\right]$
$\displaystyle\leq\exp(a_{k}\nu
Y_{k})\,\left(\mathbb{E}\left[\left.\exp(-2b_{k}\nu\sqrt{k+k_{0}-1}(X_{k}-x^{\star})^{\top}Z_{k})\right|X_{k}\right]\right)^{1/2}\times$
$\displaystyle\,\,\,\,\,\,\,\left(\mathbb{E}\left[\left.\exp(2c_{k}\nu\|Z_{k}\|^{2})\right|X_{k}\right]\right)^{1/2}$
$\displaystyle\leq\exp(a_{k}\nu
Y_{k})\,\exp(2b_{k}^{2}\nu^{2}Y_{k})\,\left(\mathbb{E}\left[\left.\exp(2c_{k}\nu\|Z_{k}\|^{2})\right|X_{k}\right]\right)^{1/2}$
(73)

Recalling (2), the last term is finite for $\nu\leq 1/(2a^{2}C_{2})=:B_{0},$
and for such $\nu$, the corresponding value is equal to $\exp(C_{2}c_{k}\nu)$.
Thus, for each $\nu\leq B_{0},$

$\Phi_{k+1|k}(\nu;X_{k})\leq\exp(\nu(a_{k}+2b_{k}^{2}\nu)Y_{k}+C_{2}c_{k}\nu).$
(74)

It is easy to see that $B\leq B_{0}.$ Consider $\nu\leq B$. Taking the
expectation on both sides of (74), the following recursive inequality on
$\Phi_{k}$ is obtained for any $\nu\leq B$ and any $k\geq 1:$

$\Phi_{k+1}(\nu)\leq\Phi((a_{k}+2b_{k}^{2}B)\nu)\exp(C_{2}c_{k}\nu).$ (75)

From this point, the proof proceeds similarly as in [12], i.e., by induction,
and using $k=1$ as the base, it can be shown that, for each $\nu\leq B,$

$\Phi_{k}(\nu)\leq e^{\frac{\nu}{B}}.$ (76)

By exponential Markov, from (76), for each $\nu\leq B,$

$\mathbb{P}\left(Y_{k}\geq\epsilon\right)\leq\mathbb{E}\left[\exp{\nu
Y_{k}}e^{-\nu\epsilon}\right].$ (77)

Taking $\nu=B$ yields the desired result. ∎

## Appendix C.

###### Proof of (56) .

Introduce step-wise constant function $s_{k}:[0,1]\mapsto\mathbb{R}$, defined
by

$s_{k}(\theta)=\left\\{\begin{array}[]{cc}\Lambda(\alpha_{l}\eta_{l};x^{\star}),&\mathrm{for\,}\frac{l-1}{k}<\theta\leq\frac{l}{k}\\\
0,&\mathrm{for\,}\theta=0\end{array}\right..$ (78)

It is easy to verify that the integral of $s_{k}$ over $[0,1]$ equals the
desired sum in the right hand-side of (5), i.e.,

$\int_{0}^{1}s_{k}(\theta)=\frac{1}{k}\sum_{l=1}^{k}\Lambda(\alpha_{l}\eta_{l};x^{\star}).$
(79)

We next show that

$\lim_{k\rightarrow+\infty}s_{k}(\theta)=\Lambda(aQD(\theta)Q^{\top}\lambda;x^{\star}),$
(80)

where $D(\theta)$ is as defined in the claim of the theorem. To show the
preceding limit, note that, for each $\theta\in(0,1]$,

$s_{k}(\theta)=\Lambda(k\alpha_{l_{k}}QD_{k,l_{k}}Q^{\top}\lambda;x^{\star})$
(81)

where $[D_{k,l_{k}}]_{ii}=\beta_{k,l_{k}}(\rho_{i},0),$ $l_{k}$ is the index
of the interval in the definition of $s_{k}$ to which $\theta$ belongs,
$l_{k}=\min\\{l=1,..,k:\theta\leq\frac{l}{k}\\},$ and $\rho_{i}$ is, we
recall, the $i$-th eigenvalue of $H^{\star}.$

Using the bounds from Lemma 2, it is easy to establish the by sandwiching
argument that

$\lim_{k\rightarrow\infty}k\alpha_{l_{k}}\beta_{k,l_{k}}(\rho_{i},0)=a\theta^{a\rho_{i}-1}.$
(82)

The limit in (80) now follows by the continuity of $\Lambda(\cdot;x^{\star}),$
which follows by convexity of $\Lambda(\cdot;x^{\star}),$ Lemma 3.

Using the fact that $s_{k}$ can be uniformly bounded for all $k$ and
$\theta\in[0,1]$, we can exchange the order of the limit and the integral, to
obtain:

$\lim_{k\rightarrow+\infty}\int_{0}^{1}s_{k}(\theta)d\theta=\int_{0}^{1}\lim_{k\rightarrow+\infty}s_{k}(\theta)d\theta=\int_{0}^{1}\Lambda(aQD(\theta)Q^{\top}\lambda;x^{\star}),$
(83)

establishing the claim of the lemma. ∎

## Appendix D.

###### Derivations with Remark 6.

Consider function $\overline{\Psi}(\lambda)=\Psi^{\star}(\lambda)+r(\lambda)$
in Lemma 6. We derive here a lower bound on rate function $\overline{I}$ in
Theorem 1 that does not explicitly depend on $H(x^{\star})$. In view of the
fact that $\overline{I}$ is the Fenchel-Legendre transform of
$\overline{\Psi}$, a lower bound on $\overline{I}$ is readily obtained by
deriving an upper bound on $\overline{\Psi}(\lambda).$ Note that $r(\lambda)$
does not explicitly depend on $H(x^{\star})$, hence we only need to derive an
upper bound on $\Psi^{\star}$. By Assumption 5, we have, for any
$\theta\in[0,1]$, that
$\Lambda(aQD(\theta)Q^{\top}\lambda;x^{\star})\leq\frac{C_{1}\,a^{2}}{2}\lambda^{\top}(QDQ^{\top})^{2}\lambda$
$\leq\frac{C_{1}\,a^{2}}{2}\|\lambda\|^{2}\,\|D(\theta)\|^{2}$, where we
recall that $\|\cdot\|$ denotes the 2-norm of its vector or matrix argument.
Next, note that $\|D(\theta)\|\leq\theta^{a\,\mu-1}$, for all
$\theta\in[0,1]$, because all eigenvalues $\rho_{i}$’s of $H(x^{\star})$
belong to the interval $[\mu,L]$. Therefore, we obtain:

$\displaystyle\Psi^{\star}(\lambda)\leq\frac{C_{1}\,a^{2}}{2}\|\lambda\|^{2}\int_{0}^{1}\theta^{2a\,\mu-2}d\theta=\frac{C_{1}\,a^{2}}{2(2\,a\mu-1)}\|\lambda\|^{2}.$

∎

###### The case of random initial iterate $X_{1}$.

Recall that, by definition,
$\Gamma_{1}(\lambda)=\mathbb{E}\left[e^{\lambda^{\top}(X_{1}-x^{\star})}\right].$
When $X_{1}$ is random, $\Gamma_{1}$, as a function of $\lambda,$ is therefore
the log-moment generating function of $X_{1}-x^{\star}.$ Provided its domain
is $\mathbb{R}^{d},$ all arguments in the proof of Theorem 1 remain the same.
In particular, in eq. (5), the factor $e^{\|\lambda\|\|X_{1}-x^{\star}\|}$
would be replaced by a (finite-valued) function (of $\lambda$), and the
subsequent results would be unaltered; a similar comment applies for the
statement and the proof of Lemma 5. ∎

## Appendix E.

###### Proof of Theorem 2.

It is easy to show that for the assumed quadratic form, the iterates $X_{k}$
have the following representation:

$X_{k+1}=A_{k0}X_{1}+\sum_{l=1}^{k}\alpha_{l}A_{k,l+1}Z_{l},$ (84)

where $A_{k,l}=\prod_{j=l}^{k}(I-\alpha_{j}H)$. By the assumption that the
noise realizations at different times are independent and with a constant
distribution, we obtain:

$\Gamma_{k+1}(\lambda)=e^{\lambda^{\top}X_{1}}e^{\sum_{l=1}^{k}\Lambda(\alpha_{l}A_{k,l+1}\lambda)}.$
(85)

The proof now follows from Lemma 4 and the limit established in 56. ∎

## Appendix F. Numerical results

We now illustrate the achieved results through a numerical simulation. We
consider a strongly convex quadratic cost function
$f:{\mathbb{R}}^{d}\rightarrow\mathbb{R}$, defined by
$f(x)=\frac{1}{2}x^{\top}Ax+bx$, $d=10$, where the symmetric $d\times d$
matrix $A$ and the $d\times 1$ vector $b$ are generated randomly.
Specifically, we generate the entries of $b$ mutually independently, according
to the standard normal distribution. The matrix $A$ is generated as follows.
We let $A=Q\Lambda Q^{\top}$, where $Q$ is the matrix whose columns are the
orthonormal eigenvectors of matrix $(B+B^{\top})/2$, and the entries of $B$
are drawn mutually independently from the standard normal distribution; the
matrix $\Lambda$ is the diagonal matrix whose diagonal entries are drawn from
the uniform distribution on the interval $[1,2]$. Clearly, the optimal
solution for the problem equals $x^{\star}=A^{-1}b$.

We consider the gradient noise that is generated in an i.i.d. manner over
iterations and over the gradient noise vector elements, independently from the
solution iterate sequence. Two different noise distributions per gradient
noise entry are considered, such that the per-entry noise variance is kept
equal for the two distributions, equal to $\sigma^{2}$. In this way, we
evaluate the effects of higher order moments on the performance of SGD. The
first distribution is zero-mean Gaussian with variance $\sigma^{2}$. The
second distribution is the zero-mean Laplacian with the same variance. We set
$\sigma^{2}=0.04.$

We numerically estimate, via Monte Carlo simulations, the probability
$P\left(\|X_{k}-x^{\star}\|>\delta\right)$ along iterations $k=1,2,...$ We
denote the corresponding numerical estimate by ${p}_{k}$. Two different values
of $\delta$ are considered, $\delta=0.3$, and $\delta=0.03$. For each Monte
Carlo run, $X_{1}$ is set to the zero vector. For the numerical example here,
$\|x^{\star}\|=2.342$, and hence $\delta=0.3$ corresponds to the relative
error level $\delta/\|x^{\star}\|\approx 0.13$, while $\delta=0.03$
corresponds to $\delta/\|x^{\star}\|\approx 0.013$. Figure 1 plots ${p}_{k}$
versus iteration counter $k$ (in linear scale for the horizontal axis, and
$\mathrm{log}_{10}$-scale for the vertical axis) for the Gaussian noise case
(blue line) and the Laplacian noise case (red line). The top Figure is for
$\delta=0.3$, and the bottom Figure is for $\delta=0.03$. We can see that, for
a large value of $\delta$, the two curves are very different: the Laplacian
gradient noise case leads to a worse performance. This is because, for large
$\delta$, the argument $\lambda$ of the LMGF $\Lambda$ that corresponds to the
minimizer in the rate function value $I^{\star}$ is large (see Theorem 4), and
hence higher order polynomial coefficients ($\sim\lambda^{3}$ and higher) play
a significant role. As the higher order moments of the Gaussian and Laplace
distributions are very different (equal to zero for the Gaussian and strictly
positive for the Laplacian), the result is the different large deviations
performance (worse for the Laplacian case) as seen in Figure 1, top. On the
other hand, for a small value of $\delta$ (bottom Figure), the argument
$\lambda$ of the LMGF that corresponds to the minimizer in the rate function
expression $I^{\star}$ is small, and hence only the first two order polynomial
coefficients of $\Lambda$ play a significant role. As the two distributions
here are both zero mean and have equal variance (hence having equal first and
second order moments), the large deviation performance for the two noises
matches, as seen in Figure 1, bottom. This behavior is in accordance with the
theory derived.

Figure 1: Monte Carlo estimate of $P\left(\|X_{k}-x^{\star}\|>\delta\right)$
along iterations $k=1,2,...$ for SGD with Gaussian (blue line) and Laplacian
(red line) gradient noise with equal per-entry variance $\sigma^{2}=0.04$. Top
Figure: $\delta=0.3$; Bottom Figure: $\delta=0.03$.