On Acceleration with Noise-Corrupted Gradients

Introduction

First-order methods for convex optimization play a fundamental role in the solution of modern large-scale computational problems, encompassing applications in machine learning (Bubeck, 2014), scientific computing (Spielman & Teng, 2004; Kelner et al., 2013) and combinatorial optimization (Sherman, 2017; Ene & Nguyen, 2016). A central object of study in this area is the notion of acceleration – an algorithmic technique that can be deployed when minimizing a smooth convex function $f(\cdot)$ via queries to a first-order oracle (a blackbox that on input $\mathbf{x}\in\mathcal{X}$, returns the vector $\nabla f(\mathbf{x})$ in constant time). In this setting, a function $f(\cdot)$ is $L$-smooth if it is differentiable and its gradient is $L$-Lipschitz continuous w.r.t to a pair of dual norms $\|\cdot\|,\,\|\cdot\|_{*}$, i.e.:

Notably, the idea of acceleration can be generalized beyond this notion of smoothness to various weakly smooth problems. Examples include problems in which $f(\cdot)$ has Hölder-continuous gradients, and even the problems with certain structured non-smooth objectives (Nesterov, 2005; Allen-Zhu & Orecchia, 2015; Lu et al., 2016). In this paper, we restrict our attention to the original smooth setting, to which all others can be traced back.

Acceleration is interesting because it yields faster algorithms than classical steepest-descent algorithms, often matching or closely approximating known information-theoretic lower bounds on the number of necessary queries to the oracle. In the simplest smooth setting, the optimal accelerated algorithm, Accelerated Gradient Descent (Nesterov, 1983), achieves an error that scales as $O(1/k^{2}),$ where $k$ is the number of oracle queries. This should be compared to the convergence of steepest-descent methods, which attempt to locally minimize the first-order approximation to the function and only yield $O(1/k)$-convergence (Ben-Tal & Nemirovski, 2001; Nesterov, 2013). Many of the workhorses of optimization, such as conjugate gradient and FISTA (Beck & Teboulle, 2009), are instantiations of accelerated algorithms.

Because of its generality, acceleration still proves an active topic of research. In particular, two weaknesses in the classical presentation of accelerated methods have recently attracted attention of scholars and practitioners alike: 1) the complexity and lack of underlying intuition in the convergence analysis of accelerated methods, and 2) the apparent lack of robustness to perturbations of the gradient oracle displayed by accelerated methods when compared to their non-accelerated counterparts.

Recently, some of the mystery of acceleration has faded, as different works have provided natural interpretations and alternative proofs for accelerated methods (Allen-Zhu & Orecchia, 2017; Krichene et al., 2015; Wibisono et al., 2016; Bubeck et al., 2015; Lessard et al., 2016; Hu & Lessard, 2017; Diakonikolas & Orecchia, 2017). Of particular interest to our work is the framework of (Diakonikolas & Orecchia, 2017), which completely derives accelerated algorithms from the Euler discretization of a continuous dynamics that minimizes a natural notion of duality gap.

In terms of robustness, it has long been observed empirically that a naïve application of accelerated algorithms to inexact oracles often leads to error accumulation, even in the setting of random perturbations, while standard steepest descent algorithms do not suffer from this problem (Hardt, 2014). From a theoretical point of view, a number of papers have introduced oracle models that account for inexact gradient information. For example, (d’Aspremont, 2008) proposed a restricted model of perturbations to the gradient that preserves the possibility of acceleration. More recently, (Devolder et al., 2014) proposed a more general framework that allows for larger perturbations and seems to capture the error accumulation and instability observed in practice for accelerated methods. In these works, the inexact oracle outputs an arbitrary deterministic perturbation of the true gradient oracle. In particular, (Devolder et al., 2014) shows that such perturbations can be adversarially chosen to encode non-smooth problems.

For stochastic perturbations, (Lan, 2012; Ghadimi & Lan, 2012, 2013) considered an additive-noise model, under which (Lan, 2012; Ghadimi & Lan, 2012) obtained an optimal convergence bound for the accelerated algorithm ac-sa in the smooth, non-strongly convex setting, but sub-optimal for the smooth, strongly-convex case.In particular, the deterministic term in the convergence bound in (Ghadimi & Lan, 2012) decreases as $O(1/k^{2})$ instead of the optimal $O(1-1/\sqrt{\kappa})^{k}$ convergence, where $\kappa$ is the objective function’s condition number. This bound was further improved to the optimal one in (Ghadimi & Lan, 2013) for the setting of constrained smooth and strongly convex minimization, by coupling ac-sa algorithm from (Ghadimi & Lan, 2012) with a domain-shrinking procedure. More recently, (Jain et al., 2018) completely closed this gap for the case of linear regression. Additionally, (Dvurechensky & Gasnikov, 2016) unified the deterministic model (Devolder et al., 2014), the stochastic model (Ghadimi & Lan, 2012), and the associated results. These references are the most closely related to our work.

We study the issue of robustness of accelerated methods in three steps. First, we propose a novel, simple, generic accelerated algorithm agd$+$ following the framework of (Diakonikolas & Orecchia, 2017). This algorithm has a simple interpretation and analysis, and generalizes other known accelerated algorithms.

Second, we leverage the simplicity of the analysis of agd$+$ to characterize its behavior on different models of inexact oracles. Our analysis recovers the results for the deterministic oracle models of (d’Aspremont, 2008) and (Devolder et al., 2014). More generally, we consider the more general model of noise-corrupted gradient oracle, in which the true gradient $\nabla f(\mathbf{x})$ is corrupted by additive noise $\bm{\eta}$:

where the perturbation $\bm{\eta}$ may be a random variable. Such a model captures the setting of stochastic methods, in which the gradient is only estimated from a subset of its components (Lan, 2012; Ghadimi & Lan, 2012, 2013; Atchade et al., 2014; Krichene & Bartlett, 2017; Jain et al., 2018), the setting of differentially private empirical risk minimization, in which Gaussian noise is intentionally added to the gradient to protect the privacy of the data (Bassily et al., 2014), and the setting of engineering systems in which the gradient is estimated from noisy measurements (Birand et al., 2013).

Our algorithm agd$+$ is closely related to ac-sa from (Lan, 2012) and can in fact be seen as a “lazy” (dual averaging) counterpart of ac-sa. After this paper had been submitted, Gasnikov and Nesterov independently proposed a universal method for stochastic composite optimization (Gasnikov & Nesterov, 2018). While their algorithm is defined recursively and does not explicitly account for the iterative construction of a dual solution, a simple unwinding of the recursion shows that it is identical to agd$+$. However, the fact that agd$+$ is obtained and analyzed through the use of the approximate duality gap technique (Diakonikolas & Orecchia, 2017) allows us to streamline the analysis and obtain various bounds for both deterministic and stochastic models of noise. Further, in the setting of smooth and strongly convex minimization, our analysis leads to a tighter convergence bound for a single-stage algorithm (without domain-shrinking) than previously obtained in (Ghadimi & Lan, 2012, 2013) (see Section 6 for a precise statement).

There are other models of noise that are not considered here. For example, we do not consider the model that includes both multiplicative and additive error in the gradient oracle (Hu et al., 2017). Further, stochastic methods with variance reduction (see, e.g., (Schmidt et al., 2017; Allen-Zhu, 2017) and references therein) lead to a particular structure of the gradient noise variance (e.g., Lemma 3.4 in (Allen-Zhu, 2017)) that is not explored in this work. Nevertheless, we believe that our analysis is general enough to be extended to these settings as well, which is deferred to the future version of this paper.

Our results reveal an interesting discrepancy between noise tolerance in the settings of constrained and unconstrained smooth minimization. Namely, in the setting of constrained optimization, the error due to noise does not accumulate and is proportional to the diameter of the feasible region and the expected norm of the noise. In the setting of unconstrained optimization, the bound on the error incurred due to the noise accumulates, as observed empirically by (Hardt, 2014). However, our analysis also suggests a simple restart and slow down semi-heuristic for stabilizing the noise-incurred error, which allows taking advantage of both the acceleration and the noise stability under stochastic noise.

In the case of smooth and strongly convex minimization (Section 6), the error due to noise does not accumulate even if the region is unconstrained, as long as the noise is zero-mean, independent, and has bounded variance.Obtaining similar bounds for a slightly more general model that relaxes independence (similar to (Lan, 2012; Ghadimi & Lan, 2012)) is also possible; see Section 5.2. Further, using smaller step sizes than in the standard accelerated version of the method, the error due to noise decreases at rate $1/k$ (compare this to the $1/\sqrt{k}$ rate for smooth non-strongly convex functions). This means that strong convexity of a function implies higher robustness to noise.

Finally, we verify the predictions and insights from our analysis of agd$+$ by performing numerical experiments comparing agd$+$ to other accelerated and non-accelerated methods on noise-corrupted gradient oracles. A noteworthy outcome of these experiments is the following: when a natural generic restart & slow-down semi-heuristic is applied, the accelerated algorithm axgd (Diakonikolas & Orecchia, 2017) and the algorithm agd$+$ presented in this paper seem to outperform Nesterov’s agd both in expectation and in variance in the presence of large noise. Further, we note that compared to axgd, agd$+$ reduces the oracle complexity (the number of queried gradients) by a factor of two.

Notation and Preliminaries

where $\nabla f(\cdot)$ denotes the gradient of $f(\cdot)$.

Given oracle access to (possibly noise-corrupted) gradients of $f(\cdot)$, we are interested in minimizing $f(\cdot)$. We denote by $\mathbf{x}_{*}\in\arg\min_{\mathbf{x}\in\mathcal{X}}f(\mathbf{x})$ any (fixed) minimizer of $f(\cdot)$.

We assume that there is an arbitrary (but fixed) norm $\|\cdot\|$ associated with the space, and all the statements about function properties are stated with respect to that norm. We also define the dual norm $\|\cdot\|_{*}$ in the standard way: $\|\mathbf{z}\|_{*}=\sup\{{\left\langle\mathbf{z},\mathbf{x}\right\rangle}:\|\mathbf{x}\|=1\}$. The following definitions will be useful in our analysis, and thus we state them here for completeness.

If $\psi(\cdot)$ is $\mu$-strongly convex w.r.t. a norm $\|\cdot\|$ for $\mu>0$, then $\psi^{*}(\cdot)$ is $\frac{1}{\mu}$-smooth w.r.t. the norm $\|\cdot\|_{*}$.

The Bregman divergence $D_{\psi}(\mathbf{x},\mathbf{y})$ captures the difference between $\psi(\mathbf{x})$ and its first order approximation at $\mathbf{y}.$ Notice that, for a differentiable $\psi$, we have: $\nabla_{\mathbf{x}}D_{\psi}(\mathbf{x},\mathbf{y})=\nabla\psi(\mathbf{x})-\nabla\psi(\mathbf{y}).$ The Bregman divergence $D_{\psi}(\mathbf{x},\mathbf{y})$ as a function $g_{\mathbf{y}}(\mathbf{x})$ is convex. Its Bregman divergence is itself, i.e., $D_{g_{\mathbf{y}}}(\mathbf{v},\mathbf{u})=D_{\psi}(\mathbf{v},\mathbf{u}).$

Improved Accelerated Method

In this section, we focus on the setting of smooth minimization. The case of smooth and strongly convex minimization is treated in Section 6.

To design agd$+$, we define an approximate duality gap, similar to (Diakonikolas & Orecchia, 2017, 2018), but allowing for an inexact gradient oracle according to Eq. (1.2). The construction is based on maintaining three points at each iteration $k$: $\mathbf{x}_{k}$ is the point at which the gradient is queried, while $(\mathbf{y}_{k},\mathbf{z}_{k})$ is the current primal-dual solution pair at the end of iteration $k.$ For this setup, the dual solution $\mathbf{z}_{k}$ is a conic combination of the negative gradients seen so far, taken at an initial dual point $z_{0}=\nabla\psi(\mathbf{x}_{0}),$where $\mathbf{x}_{0}$ is an arbitrary initial primal solution, i.e.,

where the sequence $a_{k}>0$, $A_{k}=\mathop{\textstyle\sum}_{i=1}^{k}a_{i}$ will be specified later. By convention, $A_{0}=0$.

The choice of sequences above immediately implies upper and lower bounds on optimum at each iteration $k.$ The upper bound is simply chosen as $U_{k}=f(\mathbf{y}_{k})$. For the lower bound, by convexity of $f(\cdot)$ (see Eq. (2.1)):

To relate the lower bound to the output of the inexact oracle, it is useful to express the gradients $\nabla f(\mathbf{x}_{i})$ as $\nabla f(\mathbf{x}_{i})=\widetilde{\nabla}f(\mathbf{x}_{i})-\bm{\eta}_{i}$. Adding and subtracting $\frac{1}{A_{k}}D_{\psi}(\mathbf{x}_{*},\mathbf{x}_{0})$ in the last equation, we have:

Finally, we can replace $\mathbf{x}_{*}$ by a minimization over $\mathcal{X}$ to obtain our final lower bound:

Applying Fact 2.4 and the definition of $\mathbf{z}_{k}$ from (3.1), we have the following characterization of the last term of $L_{k}$.

Let $\mathbf{z}_{k}$ be defined as in (agd$+$). Then:

The approximate duality gap is simply defined as $G_{k}=U_{k}-L_{k}$. Observe that, by construction of $U_{k}$ and $L_{k}$, $f(\mathbf{y}_{k})-f(\mathbf{x}_{*})\leq G_{k}$. Hence, to prove the convergence of the algorithm, it suffices to bound $G_{k}$. To do so, we will track the evolution of the quantity $A_{k}G_{k},$ i.e., we will boundFrom (Diakonikolas & Orecchia, 2017) it can be derived that $A_{k}G_{k}$ is a Lyapunov function for the continuous dynamic underlying agd$+$, i.e., $E_{k}$ is the discretization error at iteration $k$. $E_{k}=A_{k}G_{k}-A_{k-1}G_{k-1}$, so that

2 The agd++Algorithm

+Algorithm The steps of agd$+$ are defined as follows:

To seed agd$+$, we let $\mathbf{x}_{1}=\mathbf{x}_{0}$, $\mathbf{y}_{1}=\mathbf{v}_{1}=\nabla\psi^{*}(\mathbf{z}_{1})$.

3 Convergence Analysis for agd++

+ To simplify the notation, from now on we denote:

We can now bound the change $E_{k}=A_{k}G_{k}-A_{k-1}G_{k-1}$ by decomposing it into two terms: $E_{k}\leq E_{k}^{e}+E_{k}^{\eta}$, where the latter term is due to the inexact nature of the gradient oracle. The following lemma allows us to bound these terms.

Let $E_{k}^{\eta}=\left\langle\bm{\eta}_{k},\mathbf{x}_{*}-\mathbf{v}_{k}\right\rangle$ and $E_{k}^{e}=A_{k}(f(\mathbf{y}_{k})-f(\mathbf{x}_{k}))-A_{k}\left\langle\nabla f(\mathbf{x}_{k}),\mathbf{y}_{k}-\mathbf{x}_{k}\right\rangle-D_{\psi}(\mathbf{v}_{k},\mathbf{v}_{k-1}).$ Then $E_{k}\leq E_{k}^{\eta}+E_{k}^{e}$.

Let $m_{k}(\mathbf{x})=\mathop{\textstyle\sum}_{i=1}^{k}a_{i}\left\langle\widetilde{\nabla}f(\mathbf{x}_{i}),\mathbf{u}-\mathbf{x}_{i}\right\rangle+D_{\psi}(\mathbf{u},\mathbf{x}_{0})$ denote the function under the minimum in the lower bound. By Proposition 3.1, $\mathbf{v}_{k}=\nabla\psi^{*}(\mathbf{z}_{k})=\arg\min_{\mathbf{x}\in\mathcal{X}}m_{k}(\mathbf{x})$. Observe that $m_{k}(\mathbf{x})=a_{k}\left\langle\widetilde{\nabla}f(\mathbf{x}_{k}),\mathbf{x}-\mathbf{x}_{k}\right\rangle+m_{k-1}(\mathbf{x})$. By the definition of Bregman divergence:

As Bregman divergence is blind to linear and zero-order terms, we have that $D_{m_{k-1}}(\mathbf{v}_{k},\mathbf{v}_{k-1})=D_{\psi}(\mathbf{v}_{k},\mathbf{v}_{k-1})$. By Proposition 3.1, $\mathbf{v}_{k-1}=\arg\min_{\mathbf{x}\in\mathcal{X}}m_{k-1}(\mathbf{x})$, and hence $\left\langle\nabla m_{k-1}(\mathbf{v}_{k-1}),\mathbf{v}_{k}-\mathbf{v}_{k-1}\right\rangle\geq 0$. Therefore,

Using the definition of $\widetilde{\nabla}f(\mathbf{x}_{k})$, the change in the lower bound is:

For the change in the upper bound, we have:

The last piece that is needed for the analysis is the bound on the initial gap $G_{1}$, obtained in the following proposition.

$A_{1}G_{1}\leq D_{\psi}(\mathbf{x}_{*},\mathbf{x}_{0})+E_{1}^{\eta}+E_{1}^{e}$, where $E_{1}^{\eta}$ is defined as in Lemma 3.2 and $E_{1}^{e}=A_{1}(f(\mathbf{y}_{1})-f(\mathbf{x}_{1})-\left\langle\nabla f(\mathbf{x}_{1}),\mathbf{v}_{1}-\mathbf{x}_{1}\right\rangle)-D_{\psi}(\mathbf{v}_{1},\mathbf{x}_{0})$.

The proof is a straightforward application of the previously introduced definitions.

4 Convergence of agd++ with Exact Oracle

+ with Exact Oracle To prove the convergence of the method in the noiseless case, in this section we assume that $\bm{\eta}_{k}=\textbf{0}$, and, consequently, $E_{k}^{\eta}=0$. Hence, to obtain a convergence bound for agd$+$, we only need to bound $E_{k}^{e}$.

By smoothness of $f(\cdot)$, $f(\mathbf{y}_{k})-f(\mathbf{x}_{k})-\left\langle\nabla f(\mathbf{x}_{k}),\mathbf{y}_{k}-\mathbf{x}_{k}\right\rangle\leq\frac{L}{2}\|\mathbf{y}_{k}-\mathbf{x}_{k}\|^{2}$. Hence:

From (agd$+$), $\mathbf{y}_{k}-\mathbf{x}_{k}=\frac{a_{k}}{A_{k}}(\mathbf{v}_{k}-\mathbf{v}_{k-1})$. As $D_{\psi}(\mathbf{v}_{k},\mathbf{v}_{k-1})\geq\frac{\mu}{2}\|\mathbf{v}_{k}-\mathbf{v}_{k-1}\|^{2}$, it follows that:

as $\frac{{a_{k}}^{2}}{A_{k}}\leq\frac{\mu}{L}$ by the theorem assumptions. Thus: $G_{k}\leq\frac{A_{1}}{A_{k}}G_{1}$ and it remains to bound $A_{1}G_{1}$, which is just:

as $\mathbf{y}_{1}=\mathbf{v}_{1}$ and $\mathbf{x}_{1}=\mathbf{x}_{0}$. ∎

Observe that for $a_{k}=\frac{\mu}{L}\cdot\frac{k+1}{2}$ we recover the standard $1/k^{2}$ convergence rate of accelerated methods.

5 Convergence of agd++ with Inexact Oracle

+ with Inexact Oracle In this subsection, we focus on bounding the error $E_{k}^{\eta}$ that is accrued due to the additive noise $\bm{\eta}_{k}$. Additional results (including other models of noise) can be found in Section 5. From Lemma 3.2, $E_{k}^{\eta}=a_{k}\left\langle\bm{\eta}_{k},\mathbf{x}_{*}-\mathbf{v}_{k}\right\rangle$, and we have the following:

Let $\mathbf{\hat{v}}_{k}=\nabla\psi^{*}(\mathbf{z}_{k}+a_{k}\bm{\eta}_{k})=\nabla\psi^{*}(\mathbf{z}_{k-1}-a_{k}\nabla f(\mathbf{x}_{k}))$. Recall that $E_{k}^{\eta}=a_{k}\left\langle\bm{\eta}_{k},\mathbf{x}_{*}-\mathbf{v}_{k}\right\rangle$. Adding and subtracting $\mathbf{\hat{v}}_{k}$:

It is possible to relax the assumption that $\bm{\eta}_{k}$’s are independent. In fact, for Lemma 3.7 to apply, it suffices that, conditioned on the natural filtration $\mathcal{F}_{k-1}$ (all the information about the noise up to the beginning of iteration $k$), $\bm{\eta}_{k}$ is independent of $\mathbf{\hat{v}}_{k}$. More details are provided in Section 5.2.

Lemma 3.7 suggests that for unconstrained smooth minimization the sequence $a_{k}$ that leads to accelerated methods aggregates noise, as for accelerated methods $a_{k}\sim k,A_{k}\sim k^{2}$. However, if we were to resort to a slower, uniform sequence (and slower $1/k$ convergence rate), then the noise would average out, as we would have constant $a_{k}$’s and $A_{k}\sim k$. Even more, if $a_{k}\sim 1/\sqrt{k}$, then the error due to noise would decrease at rate $\log(k)/\sqrt{k}$. This is confirmed by our numerical experiments and matches the experience of practitioners, as discussed by (Hardt, 2014).

Observe that we could not get a bound on variance that is independent of $R_{\mathbf{x}_{*}}$, as the variance (unlike the expectation) of $\left\langle\bm{\eta}_{k},\mathbf{x}_{*}-\mathbf{\hat{v}}_{k}\right\rangle$ is not zero. Instead, since we upper-bound the expectation of $f(\mathbf{y}_{k})-f(\mathbf{x}_{*})$ by a non-negative quantity (and $f(\mathbf{y}_{k})-f(\mathbf{x}_{*})$ is always non-negative as $\mathbf{x}_{*}$ is the minimizer of $f(\cdot)$), we can apply Markov’s Inequality to obtain a concentration bound on $f(\mathbf{y}_{k})-f(\mathbf{x}_{*})$.

Finally, the step sizes $a_{k}$ can be chosen so as to balance the deterministic error and the error due to noise in the convergence bound. This leads to the following corollary.

Noise-Error Reduction

Based on the results of our analysis from Section 3, we now discuss how these results can be used to control the error of agd$+$ that is incurred due to the gradient oracle noise. First, we discuss how to prevent error accumulation from Lemma 3.7, which is incurred when running a vanilla version of agd$+$. The main idea is to take advantage of acceleration until the noise accumulation starts dominating the convergence, and then switch to a slower sequence $\{a_{k}\}$ for which the error averages out and the algorithm further reduces the mean. Finally, we show how, through another algorithm restart and slow down, the sequence of updates can be made convergent (i.e., the mean error is further reduced at a rate $\sim 1/\sqrt{k}$).

Observe that the result from Corollary 3.9 already gives a convergent sequence of updates. However, the choice of parameters in Corollary 3.9 is fixed and tailored to the global problem properties and worst-case effect of the additive noise. Instead, the strategy of incrementally slowing down the algorithm can take advantage of the more local, fine-grained properties of the objective function. This is confirmed by the numerical experiments provided in Section 7.

To take advantage of acceleration at the initial stage and then stabilize the mean error due to noise, we propose the following Restart+SlowDown semi-heuristics:

The only “heuristic” part of Restart+SlowDown is deciding when to switch to the slower sequence, as, due to Lemma 3.7, slower sequence is guaranteed to lead to a better bound on the approximation error due to noise in the case of unconstrained minimization. Further, switching to a slower, linearly growing sequence $A_{k}$ is guaranteed to further reduce the error mean, as discussed in the next subsection.

The intuition behind Restart+SlowDown criterion is restarting when “the signal is drowning in noise”. In particular, $\mathbf{z}_{k}$ (the weighted sum of the noisy negative gradients) is the only gradient information used in defining all steps of agd$+$ and we can interpret it as the “signal” that is used to guide algorithm updates. When the gradients are corrupted by noise, ${\mathbf{z}_{k}}=-{\sum_{i=1}^{k}a_{i}\nabla f(\mathbf{x}_{i})}-{\sum_{i=1}^{k}a_{i}\bm{\eta}_{i}}$. As the noise is assumed to be independent, the expected energy of the signal-plus-noise is equal to the sum of the energy of the signal and the expected energy of the noise:

Hence, when the criterion of Restart+SlowDown is satisfied, the energy component due to noise dominates the energy component of the signal in $\mathbf{z}_{k}$.

For constrained minimization with a small diameter, Restart+SlowDown cannot reduce the theoretical mean of the error due to noise (unless the bound from Lemma 3.7 dominates the bound from Proposition 3.5), as the noise term averages out regardless of the sequence $\{a_{i}\}$ (Proposition 3.5). Nevertheless, a slower, uniform sequence $\{a_{i}\}$ has lower variance than the accelerated sequence, and can be beneficial in the settings where the accelerated sequence produces high error variance.

2 Further Mean-Error Reduction

Quadratically-growing sequence $\{A_{i}\}$ (or linearly growing sequence $\{a_{i}\}$) is the fastest-growing sequence which guarantees that $A_{k}G_{k}$ is non-increasing in the case of smooth minimization with exact gradients. When we switch to a slower sequence $\{A_{i}\}$ by invoking Restart+SlowDown, this creates more slack in making $A_{k}G_{k}$ non-increasing in the presence of gradient noise. Hence, Restart+SlowDown reduces the mean error and keeps it bounded. However, with Restart+SlowDown alone, the mean error cannot converge to zero. To ensure that the error is converging to zero, we can perform an additional Restart+SlowDown (Restart+SlowDown-2), which uses the same criterion for restart, but slows down the sequence $a_{k}$ to $a_{k}\sim 1/\sqrt{k}$, as follows.

Thus, we have the following Corollary (of Lemma 3.7):

Finally, observe that the factor of $\log(k+1)$ in the bound from Corollary 4.1 can be removed if the number of steps $K$ is fixed in advance and $a_{k}$’s are set to $a_{k}=\frac{\mu}{L\sqrt{K}}$.

Different Models of Inexact Oracle

There are two main adversarial models of inexact gradient oracles that have been used in the convergence analysis of accelerated methods: the approximate gradient model of (d’Aspremont, 2008) and inexact first-order oracle of (Devolder et al., 2014). The approximate gradient model (d’Aspremont, 2008) defines the inexact oracle by a deterministic perturbation satisfying the following condition for all queries:

Hence, this model is only applicable to constrained optimization with bounded-diameter domain and bounded (adversarial) additive noise. Under these assumptions, (d’Aspremont, 2008) proves that it is possible to approximate $f(\mathbf{x}_{*})$ up to an error of $\delta$ achieving an accelerated rate. We can show the same asymptotic boundWe actually obtain better constants than those in Theorem 2.2 of (d’Aspremont, 2008). by applying the assumption to Equation (3.5) in Proposition 3.5. This yields:

whenever $\frac{{a_{k}}^{2}}{A_{k}}\leq\frac{\mu}{L}$ for all $k$. Setting $a_{k}=\frac{\mu}{L}\cdot\frac{k+1}{2}$ yields $A_{k}=k^{2}+O(k)$, which establishes the accelerated decrease of the first term in the error bound above.

The inexact first-order oracle (Devolder et al., 2014) is a generalization of the model from (d’Aspremont, 2008) that defines the inexact oracle by:

As stated in (Devolder et al., 2014), this model does not apply to noise-corrupted gradients per se, but rather to “non-smooth and weakly smooth convex problems”. In other words, the model was introduced to characterize the behavior of accelerated methods on objective functions that are non-smooth, but close to smooth. Our results agree with those of (Devolder et al., 2014) and lead to the same kind of error accumulation. To see this, observe that we only use the definition of smoothness when bounding $E_{k}^{e}$ in Theorem 3.4. Thus, the error from the inexact oracle would only appear as $E=E_{k}^{e}\leq A_{k}\delta$, leading to:

This is exactly the same bound as in Theorem 4 of (Devolder et al., 2014), but we obtain it through a generic algorithm with a simple analysis.

2 Generalized Stochastic Models

agd++ for Smooth and Strongly Convex Minimization

To analyze agd$+$ in this setting, we need to use a stronger lower bound $L_{k}$, which is constructed by the same arguments as before, but now using strong convexity instead of regular convexity. Such a construction gives:

While it suffices to have $\psi$ be an arbitrary function that is strongly convex w.r.t. the $\|\cdot\|_{2}$, for simplicity, we take $\psi(\mathbf{x})=\frac{\mu_{0}}{2}\|\mathbf{x}\|^{2}$, where $\mu_{0}$ will be specified later.

For $\theta_{k}=\frac{a_{k}}{A_{k}}$, the algorithm can now be stated as follows:

where, $\mathbf{x}_{1}=\mathbf{x}_{0}=\mathbf{y}_{0}=\mathbf{v}_{0}$ is an arbitrary initial point from $\mathcal{X}$.

As before, the main convergence argument is to show that $A_{k}G_{k}\leq A_{k-1}G_{k-1}$ and combine it with the bound on the initial gap $G_{1}$. We start with bounding the initial gap, as follows.

If $\psi(\mathbf{x})=\frac{\mu_{0}}{2}\|\mathbf{x}\|^{2}$, where $\mu_{0}=a_{1}(L-\mu)$, then $A_{1}G_{1}\leq\frac{A_{1}(L-\mu)}{2}\|\mathbf{x}_{*}-\mathbf{x}_{0}\|^{2}+E_{1}^{\eta}$, where $E_{1}^{\eta}=a_{1}\left\langle\bm{\eta}_{1},\mathbf{x}_{*}-\mathbf{v}_{1}\right\rangle$.

As $\mathbf{x}_{1}=\mathbf{x}_{0}$, the initial lower bound is:

As $a_{1}=A_{1}$, it follows that $\mathbf{y}_{1}=\mathbf{v}_{1}$, and hence:

where the inequality is by the smoothness of $f(\cdot)$. Combining the bounds on the initial upper and lower bounds, it follows:

as, by the initial assumption, $\mu_{0}=a_{1}(\mu-L)$ and $E_{1}^{\eta}=a_{1}\left\langle\bm{\eta}_{1},\mathbf{x}_{*}-\mathbf{v}_{1}\right\rangle$. ∎

To bound the change in the lower bound, it is useful to first bound $m_{k}(\mathbf{v}_{k})-m_{k-1}(\mathbf{v}_{k-1})$, as in the following technical proposition.

Let $\psi(\mathbf{x})=\frac{a_{1}(L-\mu)}{2}\|\mathbf{x}\|^{2}$. Then:

Observe that, by the definition of $m_{k}(\cdot)$, $m_{k}(\mathbf{v}_{k})=m_{k-1}(\mathbf{v}_{k})+a_{k}\left\langle\widetilde{\nabla}f(\mathbf{x}_{k}),\mathbf{v}_{k}-\mathbf{x}_{k}\right\rangle+a_{k}\frac{\mu}{2}\|\mathbf{v}_{k}-\mathbf{x}_{k}\|.$

The rest of the proof bounds $m_{k-1}(\mathbf{v}_{k})-m_{k-1}(\mathbf{v}_{k-1})$. Observe that, as $\mathbf{v}_{k-1}=\operatorname*{argmin}_{\mathbf{u}\in\mathcal{X}}m_{k-1}(\mathbf{u})$, it must be $\left\langle\nabla m_{k-1}(\mathbf{v}_{k-1}),\mathbf{u}-\mathbf{v}_{k-1}\right\rangle\geq 0$, $\forall\mathbf{u}\in\mathcal{X}$. As Bregman divergence is blind to linear terms:

The rest of the proof is by $\frac{a_{1}(L-\mu)}{2}\|\mathbf{v}_{k}-\mathbf{v}_{k-1}\|^{2}\geq 0$. ∎

We are now ready to move to the main part of the convergence argument, namely, to show that $A_{k}G_{k}\leq A_{k-1}G_{k-1}$ for a certain choice of $a_{k}$.

Let $\psi(\mathbf{x})=\frac{a_{1}(L-\mu)}{2}\|\mathbf{x}\|^{2}$ and $0<a_{k}\leq A_{k}\sqrt{\frac{\mu}{L}}$. Then: $A_{k}G_{k}\leq A_{k-1}G_{k-1}+E_{k}^{\eta}$, where $E_{k}^{\eta}=a_{k}\left\langle\bm{\eta}_{k},\mathbf{x}_{*}-\mathbf{v}_{k}\right\rangle$.

As $U_{k}=f(\mathbf{y}_{k})$, we have that:

Using Proposition 6.2, the change in the lower bound is:

Denote $\mathbf{w}_{k}=\frac{A_{k-1}}{A_{k}}\mathbf{v}_{k-1}+\frac{a_{k}}{A_{k}}\mathbf{x}_{k}$. By Jensen’s Inequality:

Write $a_{k}\left\langle\widetilde{\nabla}f(\mathbf{x}_{k}),\mathbf{v}_{k}-\mathbf{x}_{k}\right\rangle$ as:

As $\frac{a_{k}}{A_{k}}\leq\sqrt{\frac{\mu}{L}}$ and by smoothness of $f(\cdot):$

Combining (6.3)-(6.6), we have the following bound for the change in the lower bound:

Using the definition of $\mathbf{w}_{k}$, $\theta_{k}=\frac{a_{k}}{A_{k}},$ and ($\mu$agd$+$), it is not hard to verify that $\mathbf{y}_{k}=\mathbf{x}_{k}+\frac{a_{k}}{A_{k}}(\mathbf{v}_{k}-\mathbf{w}_{k})$ and $\frac{a_{k}}{A_{k}}(\mathbf{v}_{k-1}-\mathbf{x}_{k})=\mathbf{x}_{k}-\mathbf{y}_{k-1}$, which, using the convexity of $f(\cdot)$, gives:

Combining (6.7) and (6.2), the proof follows. ∎

Let $\psi(\mathbf{x})=\frac{L-\mu}{2}\|\mathbf{x}\|^{2}$, $a_{1}=1$, $\frac{a_{i}}{A_{i}}=\gamma_{i}\leq\sqrt{\frac{\mu}{L}}$ for $i\geq 2$, and let $\mathbf{y}_{k},\mathbf{x}_{k}$ evolve according to ($\mu$agd$+$). Then, $\forall k\geq 1$:

Applying Lemma 6.3, it follows that $G_{k}\leq\frac{A_{1}G_{1}}{A_{k}}+\frac{\sum_{i=1}^{k}E_{i}^{\eta}}{A_{k}}=\frac{A_{1}}{A_{2}}\cdot\frac{A_{2}}{A_{3}}\cdot\dots\cdot\frac{A_{k-1}}{A_{k}}G_{1}+\frac{\sum_{i=1}^{k}E_{i}^{\eta}}{A_{k}}$. As $\frac{A_{i-1}}{A_{i}}=1-\frac{a_{i}}{A_{i}}=1-\gamma_{i}$, we have $G_{k}\leq\left(\Pi_{i=1}^{k}\left(1-\gamma_{i}\right)\right)G_{1}+\frac{\sum_{i=1}^{k}E_{i}^{\eta}}{A_{k}}$. The rest of the proof is by applying Proposition 6.1 and using that $f(\mathbf{y}_{k})-f(\mathbf{x}_{*})\leq G_{k}$. ∎

Using the same arguments for bounding the noise term as in the case of smooth minimization (Section 3), we have the following corollary.

If $\bm{\eta}_{i}=\textbf{0}$ (the noiseless gradient case), setting $\gamma_{i}=\sqrt{\frac{\mu}{L}}$, we recover the standard convergence result for accelerated smooth and strongly convex minimization:

$\frac{a_{i}}{A_{i}}=\gamma_{i}=\sqrt{\frac{\mu}{L}},$

Assume that $\bm{\eta}_{i}$’s are zero-mean and independent and denote $\psi_{k}(\mathbf{x})=\sum_{i=1}^{k}a_{i}\frac{\mu}{2}\|\mathbf{x}-\mathbf{x}_{i}\|^{2}+\frac{\mu_{0}}{2}\|\mathbf{x}-\mathbf{x}_{0}\|^{2}$. Observe that the strong convexity parameter of $\psi_{k}$ is $\mu A_{k}+\mu_{0}>\mu A_{k}.$ From Fact 2.4, $\mathbf{v}_{k}=\nabla\psi^{*}_{k}(\mathbf{z}_{k})$. Similarly as for the case of smooth minimization, let $\mathbf{\hat{v}}_{k}=\nabla\psi^{*}(\mathbf{z}_{k}+a_{k}\bm{\eta}_{k})$. Then $\mathbf{\hat{v}}_{k}$ is independent of $\bm{\eta}_{k},$ and, using Fact 2.5, we have:

Note that in the setting of constrained (bounded-diameter) minimization, (Ghadimi & Lan, 2013) obtained the optimal convergence bound $O\left(\left(1-\sqrt{\frac{\mu}{L}}\right)^{k}\cdot\frac{(L-\mu)\|\mathbf{x}_{*}-\mathbf{x}_{0}\|^{2}}{2}+\cdot\frac{\sigma^{2}}{\mu k}\right)$ by coupling the algorithm from (Ghadimi & Lan, 2012) with a domain-shrinking procedure resulting in a multi-stage algorithm. We expect it is possible to obtain a similar result for $\mu$agd$+$ by coupling it with the domain-shrinking from (Ghadimi & Lan, 2013).

Numerical Experiments

In all the experiments, we used standard Python libraries to solve the considered problems to high accuracy. The resulting function value is denoted by $\hat{f}^{*}$ in the figures. In all the problems, we used $\psi(\mathbf{x})=\frac{L}{2}\|\mathbf{x}\|_{2}^{2}$ as the regularizer. For constrained problems, we implemented projected gradient descent as the “gd” algorithm.

In the graphs, to-agd$+$ denotes the “theoretically optimal” version of agd$+$; namely, it corresponds to agd$+$ with step sizes chosen according to Corollary 3.9 and Remark 3.10. In all the experiments, we compare the different accelerated algorithms (agd$+$, agd, axgd) and the non-accelerated gd under i.i.d. additive gradient noise $\bm{\eta}_{i}\sim\mathcal{N}(\textbf{0},\sigma_{\eta}I)$.

To understand the worst-case performance of agd$+$, we first compare it to Nesterov’s agd and (Diakonikolas & Orecchia, 2018)’s axgd. The instance is an unconstrained minimization problem, where $f(\mathbf{x})=\frac{1}{2}\left\langle\mathbf{A}\mathbf{x},\mathbf{x}\right\rangle-\left\langle\mathbf{b},\mathbf{x}\right\rangle$, $\mathbf{A}$ is the graph Laplacian of a cycleNamely, the difference of a tridiagonal square matrix $\mathbf{C}$ with 1’s on the main diagonal and -1’s on the remaining diagonals, and matrix $\mathbf{B}$, which is zero everywhere except for $B_{1n}=B_{n1}=1$., $b_{1}=-b_{n}=1$ and vector $\mathbf{b}$ is zero elsewhere. The initial point $\mathbf{x}_{0}$ is an all-zeros vector. The dimension of the problem is $n=100$.

The performance of agd$+$, agd, and axgd together with the performance of the slower, unaccelerated gd on the described worst-case instance is shown in Fig. 1(a)-1(d), for the exact gradient oracle (Fig. 1(a)) and noise-corrupted gradient oracle with i.i.d. $\bm{\eta}_{i}\sim\mathcal{N}(\textbf{0},\sigma_{\eta}I)$ (Fig. 1(b)-1(d)). We repeated the same experiments when the parameters for agd$+$ are chosen according to Corollary 3.9 (denoted as to-agd$+$) and when Restart+SlowDown and Restart+SlowDown-2 are employed (Fig. 1(e)-1(l)).

Without restart and slow-down, all accelerated algorithms perform similarly. In particular, as the noise standard deviation $\sigma_{\eta}$ is increased, the mean and the variance of the approximation error of all accelerated algorithms increases and the noise appears to be accumulating (see, e.g., Fig. 1(d)). On the other hand, gd generally converges to an approximation error with lower mean and variance, at the expense of converging at a slower $1/k$ rate.

When restart and slow-down are used, in the noiseless case (Fig. 1(e) and 1(i)), there is no difference compared to the vanilla case (Fig. 1(a)), which is what we want – there is no need to slow down the accelerated algorithms unless their performance is compromised by noise. In the low-noise scenario (Fig. 1(f), 1(j)), Restart+SlowDown does not change the performance of the algorithms in a noticeable way, although, in that case, the performance degradation due to noise is low. As the noise becomes higher (Fig. 1(g), 1(k), 1(f), 1(l)), restart and slow-down noticeably stabilizes all accelerated algorithms, reducing both their mean and their variance. Further, restart and slow-down generally outperforms the “theoretically optimal” agd$+$ (to-agd$+$).

“Hard Instance” over Simplex

The set of experiments in Figure 3 correspond to the minimization of the hard instance function for smooth optimization, constrained over the probability simplex. It should be compared to the unconstrained version in Figure 1. As predicted, we observe that the presence of constraints decreases the effect of error accumulation as the boundary of the feasible set limits the variance. Given the low variance due to the constraints, the effect of Restart+SlowDown is less evident for this batch of experiments.

2 Regression on Epileptic Seizure Dataset

However, the faster convergence comes at the expense of lower stability to noise as the noise becomes higher. Specifically, as the noise is increased, agd performs only marginally better and with higher variance than axgd and agd$+$ (Fig. 4(c)), and stabilizes to much higher mean and variance in the very high-noise setting (Fig. 1(d)).

Intuitively, “greedy” gradient steps that agd takes may reduce the function value significantly and lead to faster convergence in the noiseless and low-noise settings, while making the convergence very sensitive to the noise from the last iteration, as the gradient steps only depend on the last seen (noisy) gradient. In contrast, axgd and agd$+$ are more stable to noise, since both of their per-iteration steps depend on the aggregate gradient (and thus, aggregate noise) information.

As expected from the analytical results from Section 3, restart and slow-down does not noticeably improve the mean error of the algorithms (Fig. 4(e)-4(h), 4(i)-4(l)). However, in agreement with the analysis, it can reduce the error variance in the high-noise-variance setting (Fig. 4(l)). We also note that to-agd$+$ is more stable over the repeated methods’ execution, at the expense of slower initial convergence.

Logistic regression.

Finally, we evaluated the performance of the accelerated algorithms and gd for (unregularized) logistic regression on the Epileptic Seizure Recognition Dataset. The results are shown in Fig. 5.

Similar as in the case of unconstrained minimization from the beginning of this section, in the noiseless and low-noise settings (Fig. 5(a), 5(b)) all accelerated algorithms perform similarly and restart and slow-down does not lead to any noticeable improvements or degradation (Fig. 5(e), 5(i), 5(f), 5(j)). Once the noise is high enough (Fig. 5(c), 5(d)), all accelerated algorithms begin to accumulate noise, while restart and slow-down stabilize their performance to a low error mean and variance. Interestingly, in all the experiments, when Restart+SlowDown and Restart+SlowDown-2 are employed all accelerated algorithms perform at least as good as gd in terms of the error mean and variance.

Conclusion

This paper presents a new accelerated algorithm together with the analysis of its associated error bounds in the cases when the gradient oracle is corrupted by additive noise. Moreover, motivated by the analytical results, we also provide simple semi-heuristics that restart and slow down the accelerated algorithms to reduce their error mean and variance. Our numerical experiments corroborate the analytical results.

There are several interesting directions for future work that merit further investigation. For example, restart & slow-down approaches that do not require the explicit knowledge of the noise variance would be interesting for applications in engineered systems where gradients are estimated from noise-corrupted measurements.

Acknowledgements

Part of this work was done while the authors were visiting the Simons Institute for the Theory of Computing. It was partially supported by NSF grant #CCF-1718342 and by the DIMACS/Simons Collaboration on Bridging Continuous and Discrete Optimization through NSF grant #CCF-1740425. JD and LO would like to thank Guanghui Lan, Pavel Dvurechensky and the anonymous reviewers for useful comments.

The algorithm agd$+$ for the noiseless case is due to Michael B. Cohen, who termed it a “proper extension of Nesterov’s method” and shared it with JD during the Fall 2017 semester at the Simons Institute for the Theory of Computing. JD and LO dedicate this paper to the memory of Michael’s brilliance and scholarship.

References