1 Introduction↩︎

It is classically known that by using certain deterministic stepsize schedules \(\{\alpha_t\}\), the Gradient Descent algorithm (GD) \[\begin{align} x_{t+1} = x_t - \alpha_t \nabla f(x_t), \qquad t=0,1,2,\dots \label{eq-intro:gd} \end{align}\tag{1}\] solves convex optimization problems \(\min_{x \in \mathbb{R}^d} f(x)\) to arbitrary precision from an arbitrary initialization \(x_0\). This result can be found in any textbook on convex optimization, e.g., [1]–[8]. The present paper investigates the following question: \[\begin{align} \text{Is it possible to accelerate the convergence of GD by using \emph{random} stepsizes?} \label{eq-intro:motivation} \end{align}\tag{2}\] This is different from the traditional approach for obtaining accelerated convergence, which modifies the GD algorithm beyond just stepsizes (e.g., using momentum [9]); see the survey [10].

Our main result is that for separable convex optimization, the answer is yes. Specifically, we provide an i.i.d.stepsize schedule that makes GD converge at the optimal accelerated rate. See §1.1 for a precise statement; here we first contextualize with the relevant literature on random stepsizes.

The special case of quadratics.

This question is understood (only) for optimizing convex quadratic functions. In particular, the line of work [11]–[13] shows that for quadratic functions \(f\) that are \(m\)-strongly convex and \(M\)-smooth (i.e., \(m I \preceq \nabla^2 f \preceq MI\)), there is a unique optimal distribution for i.i.d.stepsizes that makes GD converge at a rate of \[\begin{align} \left( \frac{\|x_n - x^*\|}{\|x_0 - x^*\|} \right)^{1/n} \longrightarrow\; \frac{\sqrt{\kappa} - 1}{\sqrt{\kappa} + 1}\,, \label{eq-intro:rate} \end{align}\tag{3}\] where \(x^*\) denotes the minimizer of \(f\), \(\kappa := \tfrac{M}{m}\) denotes the condition number, and this convergence is in an appropriate probabilistic sense (e.g., almost-sure or \(L^1\)). This stepsize distribution is most naturally stated in terms of the inverse stepsizes \(\beta_t := 1/\alpha_t\) having the Arcsine distribution on the interval \((m,M)\), i.e., density \[\begin{align} \frac{d \mu}{d\beta} = \frac{\mathbb{1}_{\beta \in (m,M)}}{\pi \sqrt{(M-\beta)(\beta-m)}} \,. \label{eq-intro:arcsine} \end{align}\tag{4}\] Briefly, this Arcsine distribution arises naturally for the special case of quadratics because: (1) the convergence rate is invariant w.r.t.the stepsize order and thus depends only on the empirical measure of the stepsizes; and (2) the Arcsine distribution is the limit as \(n \to \infty\) of the empirical measure of the optimal \(n\) (deterministic) stepsizes for GD [14]. See §2.1 for details.

The resulting convergence rate 3 requires roughly \(\Theta(\sqrt{\kappa} \log \tfrac{1}{\varepsilon})\) iterations to produce an \(\varepsilon\)-accurate solution. Such a rate is called “accelerated” because its dependence on \(\kappa\) improves over the standard \(\Theta(\kappa \log \tfrac{1}{\varepsilon})\) convergence rate required by GD with constant stepsizes. Moreover, this convergence rate \(\Theta(\sqrt{\kappa} \log \tfrac{1}{\varepsilon})\) is information-theoretically optimal for any first-order algorithm [3], [15], and matches the worst-case convergence rate of accelerated algorithms that change GD beyond just its stepsizes, e.g., the Conjugate Gradient method [16], Polyak’s Heavy Ball method [17], Nesterov-style momentum methods [9], and more [10].

Beyond quadratics?

For general convex optimization, momentum-based algorithms still achieve these optimally accelerated rates \(\Theta(\sqrt{\kappa} \log \tfrac{1}{\varepsilon})\) [9]. However, existing acceleration results for random stepsizes do not extend beyond the special case of quadratic functions \(f\)—otherwise \(\nabla f\) is non-linear, making GD a non-linear map, which completely breaks all previous analyses. This failure of the analyses is due to fundamental differences between quadratic optimization and convex optimization:

• The optimal deterministic stepsizes for quadratic optimization [14] (which the random stepsizes are built upon) can lead to divergent behavior for convex optimization [18].

• The convergence rate is invariant w.r.t.the stepsize order for quadratic optimization (which motivates the opportunity for random stepsizes), but this is false for convex optimization [18].

• Using time-varying stepsizes is equivalent to using momentum for quadratic optimization (which enables GD to achieve accelerated rates), but this is false for convex optimization [10].

Due to these challenges, it remained unclear whether random stepsizes could lead to the accelerated rate \(\Theta(\sqrt{\kappa} \log \tfrac{1}{\varepsilon})\) for convex optimization. In fact, it was unknown even if random stepsizes could provide an arbitrary small improvement over the \(\Theta(\kappa \log \tfrac{1}{\varepsilon})\) rate of constant-stepsize GD in any setting beyond quadratic optimization.

1.1 Contribution↩︎

This paper shows that random stepsizes provably accelerate GD beyond the setting of quadratics. Specifically, our main result (Theorem 2) shows that GD converges at the optimally accelerated rate for separable convex functions. This separability property is defined below and, for \(C^2\) functions \(f\), is equivalent to commutativity of the Hessians \(\{\nabla^2 f(x)\}_{x \in \mathbb{R}^d}\), see Appendix 8.2.

1.2 Discussion↩︎

Optimality.

The rate in Theorem 2 is optimal since there is an exactly matching lower bound: there is an explicit \(\kappa\)-conditioned quadratic function for which any (Krylov) first-order algorithm cannot contract to the optimum faster than this rate \(\tfrac{\sqrt{\kappa}-1}{\sqrt{\kappa}+1}\) [3].

Relation with other stepsize schedules.

Here we provide a brief historical perspective to contextualize the random stepsizes studied in this paper. As mentioned earlier, the classic result of [14] showed that GD can be accelerated for the special case of quadratics by using deterministic stepsizes related to Chebyshev polynomials. Since then, many alternative stepsize schedules have been investigated, some adaptive, e.g., exact line search [4], Armijo-Goldstein schedules [3], Polyak-type schedules [4], and Barzilai-Borwein-type schedules [20]. While some of these schedules can be quite effective in practice, none has led to an analysis for any setting beyond quadratics that outperforms the “unaccelerated” contraction factor \((\tfrac{\|x_n - x^*\|}{\|x_0 - x^*\|})^{1/n} \leqslant\tfrac{\kappa - 1}{\kappa + 1}\) of constant stepsize GD, which corresponds to a convergence rate of \(\Theta(\kappa \log \tfrac{1}{\varepsilon})\). Starting with [18], a recent line of work has shown partially accelerated convergence rates for general convex optimization using deterministic stepsize schedules suggested by semidefinite-programming-based estimates of convergence rates [18], [19], [21]–[32], see the recent survey [33]. In particular, [19] showed the silver convergence rate \(O(\kappa^{\log_{1 + \sqrt{2}} 2} \log \tfrac{1}{\varepsilon}) \approx O(\kappa^{0.78} \log \tfrac{1}{\varepsilon})\) by using a certain fractal-like deterministic stepsize schedule called the silver stepsizes. This rate is conjecturally optimal among all deterministic stepsize schedules [19]. Theorem 2 is incomparable: it shows a fully accelerated rate \(O(\kappa^{1/2} \log \tfrac{1}{\varepsilon})\), but under the additional assumption of separability. See Table 1. This incomparability raises interesting questions about the possible benefit of randomization for GD, see the discussion of future work in §6.

Generality.

The phenomenon of acceleration via random stepsizes in Theorem 2 extends to other structural assumptions on \(f\) and other notions of convergence (see §4), and to settings where GD is performed inexactly (see §5.3). Note also that since GD does not depend on the choice of the basis \(U\) in the separable decomposition 5 , Theorem 2 only requires the existence of some such \(U\) (which is the definition of separability). In particular, the algorithm need not know \(U\).

Variability.

Randomized algorithms have random trajectories. The best run can be significantly better than the typical run (whose rate converges almost surely to the accelerated rate, see §5.2) and the worst run (whose rate can diverge but only with exponentially small probability, see §5.1).

Game-theoretic connections and lower bounds.

This paper considers random stepsizes in order to accelerate GD. This randomized perspective is also helpful for the dual problem of constructing lower bounds, i.e., exhibiting functions for which GD converges slowly. In particular, lifting to distributions over functions leads to a game between an algorithm (who chooses a distribution over stepsizes) and an adversary (who chooses a distribution over functions) which is symmetric in that the Arcsine distribution is the optimal answer for both. This leads to new perspectives on classical lower bounds for GD; details in §5.4.

The equalization property and logarithmic potential theory.

A distinguishing feature of our analysis is a connection to potential theory: the variational characterization for the optimal stepsize distribution mirrors the variational characterization for a certain equilibrium distribution which configures a unit of charged particles so as to minimize the logarithmic potential energy. The optimal distribution satisfies a remarkable “equalization property” which in the physical context amounts to a constant potential over space, and in the optimization context amounts to an identical convergence rate over all quadratic functions. A key technical insight is that martingale arguments extend this phenomenon to all separable convex functions. Links between Krylov-subspace algorithms and potential theory are classical [34]; we point out that these connections extend beyond quadratic objectives and to random stepsizes. See §2 for a conceptual overview of our analysis techniques and how we exploit this equalization property to prove Theorem 2.

1.3 Organization↩︎

§2 provides a conceptual overview of our analysis approach. §3 proves our main result, Theorem 2. §4 discusses the generality of our result. §5 provides auxiliary results about notions of convergence, stability, and game-theoretic interpretations. §6 describes several directions for future work that are motivated by our results. For brevity, several proofs are deferred to the appendix.

2 Conceptual overview: acceleration by random stepsizes↩︎

Here we provide informal derivations of our main result in order to overview the key conceptual ideas. We begin by briefly recalling in §2.1 how random stepsizes can accelerate GD for the special case of quadratic optimization. Then in §2.2, we describe the challenges for acceleration beyond quadratics, and our new techniques for circumventing these challenges.

2.1 Quadratic optimization↩︎

Consider running GD on the class \(\mathcal{Q}_{m,M}\) of quadratic functions \(f\) that are \(m\)-strongly convex and \(M\)-smooth. What is the optimal distribution from which to draw i.i.d.stepsizes? In order to understand this, we first recall the optimal deterministic stepsizes.

Optimal deterministic stepsizes.

Without loss of generality after translating the solution \(x^* = 0\), we have \[\begin{align} f(x) = \frac{1}{2} x^T H x\, \qquad \text{ where } \qquad mI \preceq H \preceq MI\,. \label{eq:hi-quad:obj} \end{align}\tag{6}\] Since the function is quadratic, its gradient is linear \(\nabla f(x) = Hx\), hence GD is a linear map \[\begin{align} x_{t+1} = x_t - \alpha_t \nabla f(x_t) = (I - \alpha_t H) x_t\,, \label{eq:hi-quad:alg} \end{align}\tag{7}\] and therefore the final iterate is \[\begin{align} x_n = \prod_{t=0}^{n-1} (I - \alpha_t H) x_0\,. \label{eq:hi-quad:poly} \end{align}\tag{8}\] The convergence rate for a given stepsize schedule—in the worst case over objective functions in \(\mathcal{Q}_{m,M}\) and initializations \(x_0\)—is therefore \[\begin{align} \sup_{f \in \mathcal{Q}_{m,M}, \, x_0 \neq x^*} \left( \frac{\|x_n - x^*\|}{\|x_0 - x^*\|} \right)^{1/n} &= \sup_{mI \preceq H \preceq MI, \, x_0 \neq x^*} \left( \frac{\|\prod_{t=0}^{n-1} (I - \alpha_t H) x_0\|}{\|x_0\|} \right)^{1/n} \nonumber \\ &= \sup_{mI \preceq H \preceq MI} \| \prod_{t=0}^{n-1} (I - \alpha_t H) \|^{1/n} \tag{9} \\ &= \sup_{\lambda \in [m,M]} \prod_{t=0}^{n-1} |1 - \alpha_t \lambda|^{1/n}\,. \tag{10} \end{align}\] Above, the penultimate step is by definition of the spectral norm, and the final step is by diagonalizing. It is a simple and celebrated fact, dating back to Young in 1953 [14], that minimizing this convergence rate 10 over stepsizes is equivalent to a certain extremal polynomial problem, and that as a consequence, the optimal stepsizes minimizing 10 are the inverses of the roots of a (translated and scaled) Chebyshev polynomial, in any order. Said more precisely: for any horizon length \(n\), the optimal \(n\) stepsizes are any permutation of \[\begin{align} \left\{ \left( \frac{M+m}{2} + \frac{M-m}{2} \cos \left( \frac{2t+1}{2n}\pi \right) \right)^{-1} \right\}_{t=0}^{n-1} \label{eq:quad-cheb} \end{align}\tag{11}\] and the corresponding \(n\)-step convergence rate in 10 is \[\begin{align} \sup_{f \in \mathcal{Q}_{m,M}, \, x_0 \neq x^*} \left( \frac{\|x_n - x^*\|}{\|x_0 - x^*\|} \right)^{1/n} = \left[ 2\frac{(\sqrt{\kappa} + 1)^n(\sqrt{\kappa} - 1)^n}{(\sqrt{\kappa} + 1)^{2n} + (\sqrt{\kappa} - 1)^{2n}} \right]^{1/n} \overset{n \to \infty}{\longrightarrow} \frac{\sqrt{\kappa} - 1}{\sqrt{\kappa} + 1}\,. \label{eq:quad-rate} \end{align}\tag{12}\]

Optimal random stepsizes.

Much less well-known is the fact that this accelerated rate can also be asymptotically achieved by simply using i.i.d.stepsizes from an appropriate distribution [11]–[13]. The intuition is that the convergence rate 10 is invariant in the order of the stepsizes, instead depending only on their empirical measure, and as \(n \to \infty\) one can achieve the same empirical measure (and thus the same convergence rate) by replacing the carefully chosen deterministic Chebyshev stepsize schedule 11 with i.i.d.stepsizes from their limiting empirical measure. See [18] for a discussion of this through the lens of exchangeability and De Finetti’s Theorem.

Concretely, if the inverse stepsizes \(\beta_t := 1/\alpha_t\) are i.i.d.from some distribution \(\mu\), then (forgoing measure-theoretic technicalities) the Law of Large Numbers ensures that the \(n\)-step convergence rate in 10 converges almost surely and in expectation to \[\begin{align} \sup_{f \in \mathcal{Q}_{m,M}, \, x_0 \neq x^*} \left( \frac{\|x_n - x^*\|}{\|x_0 - x^*\|} \right)^{1/n} &= \sup_{\lambda \in [m,M]} \prod_{t=0}^{n-1} |1 - \alpha_t \lambda|^{1/n} \nonumber \\ &= \exp\left( \sup_{\lambda \in [m,M]} \frac{1}{n} \sum_{t=0}^{n-1} \log |1 - \beta_t^{-1} \lambda| \right) \nonumber \\ &\overset{n \to \infty}{\longrightarrow} \exp\left( \sup_{\lambda \in [m,M]} \mathbb{E}_{\beta \sim \mu}\log |1 - \beta^{-1} \lambda| \right) \, \label{eq:quad-rand-rate} \end{align}\tag{13}\] This informal derivation suggests that the distribution \(\mu\) minimizing the asymptotic convergence rate 13 should be the appropriate limit as \(n \to \infty\) of the empirical measure for the inverse Chebyshev stepsize schedule (this limit is the \(\mathop{\mathrm{{\mathrm{Arcsine}}(m,M)}}\) distribution, since \(a + b\cos (\Theta)\) has the Arcsine distribution on \([a-b,a+b]\) if \(\Theta\) is uniformly distributed on \([0,\pi],\) see Figure 2), and moreover the optimal value should be the accelerated rate \(\frac{\sqrt{\kappa}-1}{\sqrt{\kappa}+1}\). This intuition is formalized by the following lemma.

This lemma was proven in the optimization literature via a long, tour-de-force integral calculation [13]. In Appendix 7, we present an alternative, short proof from Altschuler’s thesis [18] that argues by connecting this extremal problem for stepsize distributions to a classical extremal problem in potential theory about electrostatic capacitors. The key feature of this connection is the following equalization property of the Arcsine distribution.

This equalization property uniquely defines the Arcsine distribution in that there is no other distribution satisfying ?? ; see Appendix 7. In words, this property shows that the Arcsine distribution “equalizes” the convergence rate over all quadratic functions (parameterized by their curvatures \(\lambda\)). We emphasize that this is not just an algebraic identity: this equalization property also has a natural interpretation in terms of a game between an optimization algorithm (who chooses a stepsize to minimize the convergence rate) and an adversary (who chooses a quadratic function to maximize the convergence rate), see §5.4. In Appendix 7 we elaborate on the connections to potential theory, and leverage that connection to provide short proofs of these lemmas. For pedagogical reasons and to provide further intuition, in Appendix 8.1 we also provide a short semi-rigorous proof using only elementary single-variable calculus.

2.2 Beyond quadratic optimization↩︎

We now turn to accelerating GD beyond the special case of quadratic functions \(f\). The central challenge is that the curvature, as measured by the Hessian \(\nabla^2 f\), is no longer constant. This breaks the above quadratic analysis in several ways since it is no longer possible to summarize or re-parameterize the Hessians along the GD trajectory in terms of a single scalar \(\lambda \in [m,M]\). To explain these issues, let us attempt to adapt this analysis to the general convex setting, assuming for simplicity here that \(f\) is twice differentiable (our actual analysis in §3 avoids any such assumption).

Time-varying curvature.

Clearly, if \(f\) is not quadratic, then the gradient \(\nabla f\) is not linear, and as a consequence the GD map is not linear. Assuming as before that \(x^* = 0\) without loss of generality, one can generalize the linear GD map 8 to the non-quadratic setting as \[\begin{align} x_n = \prod_{t=0}^{n-1} (I - \alpha_t H_t) x_0 \qquad \text{where} \qquad H_t := \int_0^1 \nabla^2 f( ux_t) du\,, \label{eq:hi-convex:poly} \end{align}\tag{14}\] since by the Fundamental Theorem of Calculus, \(\nabla f(x_t) = H_tx_t\) and thus \(x_{t+1} = x_t - \alpha_t \nabla f(x_t) = (I - \alpha_t H_t) x_t\). The interpretation of \(H_t\) is that it is the average Hessian between \(x^* = 0\) and \(x_t\), i.e., it measures the relevant local curvature in the \(t\)-th iteration of GD. In the quadratic setting, \(H_t = H\) for all \(t\) since the Hessian \(\nabla^2f \equiv H\) is constant, which recovers 8 from 14 . However, for the general non-quadratic setting, the curvature \(H_t\) varies in \(t\). This is the overarching challenge mentioned above. For the purpose of accelerating GD, this challenge manifests in two interrelated ways: the lack of commutativity and predictability of the curvature \(H_t\) along the trajectory of GD. Below we describe both challenges and how we circumvent them.

2.2.1 Issue 1: commutativity (deriving the extremal problem)↩︎

The first issue is that for convex functions, the Hessians at different points do not necessarily commute. Thus, in general, \(\{H_t\}\) do not commute and therefore are not simultaneously diagonalizable. This precludes the key step 10 in the quadratic analysis above, which decomposes the convergence rate across eigendirections, thereby reducing the (multivariate) problem of convergence analysis to an analytically tractable (univariate) extremal problem. This is the only reason we assume separability: it is equivalent to commuting Hessians (Lemma 12). Concretely, separability allows us to show an analog of 10 in which the convergence rate of a stepsize schedule \(\{\alpha_t\}\) is bounded by \[\begin{align} \prod_{t=0}^{n-1} |1 - \alpha_t \lambda_t|^{1/n}\,, \label{eq:hi-convex:rate} \end{align}\tag{15}\] the difference to the quadratic case being that there all \(\lambda_t \in [m,M]\) are the same, whereas now \(\lambda_t \in [m,M]\) can be time-varying and can depend on all previous stepsizes \(\{\alpha_s\}_{s < t}\).

By taking an appropriate¹ limit of the \(n\)-step convergence rate 15 , it can be shown that the (logarithm of the) asymptotic rate for i.i.d.inverse stepsizes \(\beta_t := \alpha_t^{-1}\) from a distribution \(\mu\) is governed by the extremal problem: \[\begin{align} \inf_\mu \sup_{\{\mathcal{A}_t\}-\text{adapted process }\{\lambda_t \in [m,M] \}} \limsup_{n \to \infty} \frac{1}{n} \sum_{t=0}^{n-1} \mathbb{E}_{\beta_t \sim \mu} \log \left\lvert 1 - \beta_t^{-1} \lambda_t\right\rvert \label{eq:extremal-convex} \end{align}\tag{16}\] where \(\{\mathcal{A}_t\}\) denotes the natural filtration corresponding to the inverse stepsizes \(\{\beta_t\}\). Observe that 16 recovers the extremal problem 13 for the quadratic case if \(\lambda_t\) are constant, but in general is different because \(\lambda_t\) can be time-varying. This brings us to the second issue.

2.2.2 Issue 2: unpredictability (solving the extremal problem)↩︎

The second issue is that the curvature is time-varying (as measured by \(\lambda_t\)), thus less predictable, which makes it difficult to choose stepsizes that lead to fast convergence (as measured by the extremal problem 16 ). The key upshot of our random stepsizes is that the randomness of \(\beta_t\) prevents an adversarial choice of \(\lambda_t \in [m,M]\) since \(\lambda_t\) must be “chosen” before \(\beta_t\) is drawn. More precisely, we exploit the equalization property of the Arcsine distribution (Lemma 2) to argue that for any choice of \(\lambda_t\), the distribution of \(\beta_t\) ensures that \[\begin{align} \mathbb{E}_{\beta_t \sim \mu} \log |1 - \beta_t^{-1} \lambda_t| = \log \frac{\sqrt{\kappa} - 1}{\sqrt{\kappa} + 1}\,. \label{eq:equalization-unpredictability} \end{align}\tag{17}\] Intuitively, this can be interpreted as using random stepsizes to hedge between the possible curvatures \(\lambda_t\): the expected “progress” from an iteration is the same regardless of the curvature. Quantitatively, this equalization property 17 lets us argue an accelerated convergence rate using a martingale version of the Law of Large Numbers. In particular, this implies that the Arcsine distribution—which was optimal for the extremal problem ?? in the quadratic setting—achieves the same value for the extremal problem 16 in the convex setting, and therefore is also optimal for this setting (as it is more general).

3 Proof of main result↩︎

Here we prove Theorem 2. We build up to the full result by first proving it for univariate objectives, and then extending to the general (separable) multivariate setting. See §2 for a high-level overview of the analysis techniques and the two key conceptual challenges regarding unpredictability of the Hessian (dealt with in the univariate proof via the equalization property of the Arcsine distribution) and commutativity of the Hessian (dealt with in the multivariate proof via separability). For shorthand, we denote the accelerated rate by \[\begin{align} \mathop{\mathrm{\mathrm{R_{acc}}}}:= \frac{\sqrt{\kappa} - 1}{\sqrt{\kappa} + 1 }\,. \end{align}\]

3.1 Univariate setting↩︎

We begin by remarking that unforeseen complications can occur even in the univariate setting. For example, Polyak’s Heavy Ball method fails to converge for such functions [37].

We first state a simple helper lemma that lets us avoid assumptions of twice-differentiability when expanding the convergence rate of GD in the form 15 described in the overview §2.2.

For intuition, if \(f\) is twice-differentiable, then this lemma follows immediately from the Intermediate Value Theorem since \(f'(x) = f''(z) (x - x^*)\) for some \(z\) between \(x\) and \(x^*\), and the Hessian \(f''(z) \in [m,M]\) by \(m\)-strong convexity and \(M\)-smoothness. The proof is not much more difficult without the assumption of twice-differentiability.

3.2 Multivariate setting↩︎

Here we prove Theorem 2 in the general multivariate setting. The proof follows quickly from the univariate setting proved above and the definition of separability.

4 Generalizations↩︎

The purpose of this paper is to show that random stepsizes can accelerate GD on classes of strongly convex and smooth functions that are more general than quadratics. For concreteness, Theorem 2 focuses on functions that are separable, strongly convex, and smooth. But this result extends to other classes of functions. We describe two extensions below, thematically splitting modifications of the separability assumption (§4.1) from the strong convexity and smoothness assumptions (§4.2). These two types of modifications are complementary in the sense that they individually deal with the two key issues discussed in §2.2: commutativity and unpredictability of the Hessian, respectively.

4.1 Beyond separability↩︎

Here we consider one example of a structural assumption that is different from separability yet still allows for the acceleration phenomenon via random stepsizes. To introduce this, recall from Definition 1 that separability amounts to decomposability into a certain type of building block: univariate functions. Here we consider another building block: radial functions.

We begin by showing acceleration for the case of a single radial function, before moving to functions built from multiple radial functions. In words, the key idea behind the proof is that although radial functions do not have commutative Hessians—which we exploited in our main result to reduce the (multivariate) problem of convergence analysis to an analytically tractable (univariate) extremal problem, see the overview §2.2—radial functions are univariate in nature (in particular all gradients in a GD trajectory are colinear), and this enables a similar reduction.

It follows that acceleration via random stepsizes generalizes to functions which can decompose into radial functions and univariate convex functions.

4.2 Beyond strong convexity and smoothness↩︎

Recall that for twice-differentiable functions \(f\), the assumptions of \(m\)-strong convexity and \(M\)-smoothness are equivalent to the assumption that the spectrum of \(\nabla^2 f(x)\) lies in the interval \([m,M]\) for all points \(x\). The algorithmic techniques developed in this paper extend to functions \(f\) whose Hessians have more structured spectra, lying within some set \(E \subset \mathbb{R}_{> 0}\) that is not necessarily an interval. Note that strong convexity of \(f\) ensures a strictly positive real spectrum. Such extensions to more structured Hessians have led to improved results for other algorithms in quadratic settings and beyond, see e.g., [39]–[41].

This extension is rooted in connections to logarithmic potential theory, and proceeds as follows. Recall from the overview §2 that the asymptotic convergence rate of GD with inverse stepsizes i.i.d.from a distribution \(\mu\), is related to the maximum logarithmic potential of \(\mu\) over \([m,M]\), i.e., \[\begin{align} \sup_{\lambda \in [m,M]} \mathbb{E}_{\beta \sim \mu} \log |1 - \beta^{-1} \lambda|\,, \end{align}\] and that the unique distribution minimizing this quantity is the Arcsine distribution (Lemma 1). This extends to spectral sets \(E \subset \mathbb{R}_{> 0}\) beyond intervals, in particular the finite union of non-singleton compact intervals. Indeed, by a similar derivation, the convergence rate in this case is related to the analogous extremal problem \[\begin{align} \sup_{\lambda \in E} \mathbb{E}_{\beta \sim \mu} \log |1 - \beta^{-1} \lambda|\,, \end{align}\] and the unique distribution minimizing this quantity is the unique optimal stepsize distribution. What is this optimal distribution? This is a classical topic in potential theory and there is a well-developed machinery (typically based on computing conformal mappings, Green’s functions, or Fekete sets) for determining this distribution explicitly in closed-form when possible (typically this only occurs for very simple sets) or otherwise approximating it; we refer the interested reader to the surveys and textbooks [34], [42]–[44].

5 Discussion and related results↩︎

5.1 Divergent runs and notion of convergence↩︎

When considering randomized optimization algorithms, there are several notions for measuring the convergence rate. Two natural ones are: \[\begin{align} \mathbb{E}\Bigg[ \Bigg( \frac{\|x_{n} - x^*\|}{\|x_0 - x^*\|} \Bigg)^{1/n} \Bigg] \qquad \text{ vs } \qquad \mathbb{E}\Bigg[\frac{\|x_{n} - x^*\|}{\|x_0 - x^*\|} \Bigg]^{1/n} \,. \label{eq:notions} \end{align}\tag{24}\] These are different. Indeed, our main result shows that by running GD with random stepsizes, the former quantity is asymptotically equal to the accelerated rate \(\mathop{\mathrm{\mathrm{R_{acc}}}}= \tfrac{\sqrt{\kappa} - 1}{\sqrt{\kappa} + 1} < 1\). Whereas the latter quantity can be larger than \(1\) (details in Appendix 8.4). Why the discrepancy between these two notions of convergence rate? Which is relevant for optimization?

The discrepancy arises because the error \(\|x_{n} - x^*\|\) is exponentially small with overwhelming probability (leading to our convergence result for the first notion), but can be exponentially large with exponentially small probability (leading to the instability for the second notion). The key point is that for optimization algorithms, it is very desirable to succeed with overwhelming probability, even if performance is poor in the unlikely complement event. Simply re-run the algorithm if that happens. This is why we study the quantity on the left hand side of 24 . For intuition, we conclude this section with an illustrative toy example that explains this discrepancy in its simplest form.

This distinction is analogous to the classical Kelly criterion in information theory and gambling: for repeated multiplicative bets, the objective governing long-run typical performance is not the arithmetic expectation of terminal wealth, but rather the expected logarithmic growth rate, or equivalently the asymptotic geometric mean [45]–[47]. Here the sign is reversed—we seek fast multiplicative decay of error rather than fast growth of wealth—but the principle is the same. The normalized logarithm of the error captures the almost-sure exponential rate, whereas the expectation of the unnormalized error can be controlled by exponentially rare trajectories.

5.2 Typical runs and high probability bounds↩︎

Theorem 2 establishes that random stepsizes lead to an accelerated rate in the limit as the number of iterations \(n \to \infty\). The same analysis arguments can also establish non-asymptotic convergence rates that apply for any finite \(n\) and hold with high probability. Indeed, recall from 18 that our argument in §3 relates the convergence rate (for any finite \(n\)) to an empirical sum, namely \[\begin{align} \Bigg|\frac{x_{n} - x^*}{x_0 - x^*} \Bigg|^{1/n} = \exp\left( \frac{1}{n} \sum_{t=0}^{n-1} \log |1 - \beta_t^{-1} \lambda_t| \right)\,. \label{eq:rate-hp} \end{align}\tag{25}\] Now, instead of using the Law of Large Numbers to argue that this empirical sum converges to its expectation (giving the accelerated rate), one can use concentration-of-measure to argue about the deviations from its expectation.

Below we state such a result. This is established by bounding the variance of the aforementioned empirical sum via a martingale-type argument and applying Chebyshev’s inequality; and then extending this convergence rate to general multivariate objectives \(f\) by a union bound. Full proof details are deferred to Appendix 8.5. We note that one can show more refined bounds by bounding higher-order moments or the moment generating function (see, e.g., the textbooks [48]–[50]), but that requires computing involved logarithmic integrals and thus for simplicity of exposition, here we state only this simple illustrative result.

This implies that for a typical run, the convergence rate is close to the accelerated rate; see Figure 3 for an illustration.

5.3 Stability↩︎

Throughout, we have assumed for simplicity that all calculations are performed with exact gradients and in exact arithmetic. It is of practical importance to understand how inexactness affects convergence, and there is a long line of work on such questions for first-order algorithms, see, e.g., [51]–[57] and the references within. There are many models of inexact computation. Here, as an illustrative example, we study one common model: the multiplicative error model for gradients [51], [52]. Specifically, suppose that all calculations are still performed in exact arithmetic, but GD performs updates using an estimate \(\tilde{g}_t\) of the true gradient \(\nabla f(x_t)\), where \[\begin{align} \|\tilde{g}_t - \nabla f(x_t)\| \leqslant\varepsilon\|\nabla f(x_t)\|\,, \label{eq:inexact} \end{align}\tag{26}\] for some inexactness parameter \(\varepsilon\in (0,1)\). This model of inexact gradients does not align with separability, hence we state the following results for univariate objectives for simplicity. It is an interesting question to extend beyond.

The proof is deferred to Appendix 8.6, since it is similar to that of Theorem 2. Conceptually, the only difference is that the average Hessians \(\lambda_t \in [m,M]\) are now perturbed multiplicatively, and thus the effective curvatures \(\tilde{\lambda}_t\) lie in the enlarged interval \([(1-\varepsilon)m,(1+\varepsilon)M]\). The convergence rate is therefore governed by \(\exp( \mathbb{E}_{\beta \sim \mathop{\mathrm{{\mathrm{Arcsine}}(m,M)}}} \log | 1- \beta^{-1} \tilde{\lambda}|)\) where \(\tilde{\lambda}\) can now be slightly outside the interval \([m,M]\). When the inverse stepsizes are still drawn from \(\mathop{\mathrm{{\mathrm{Arcsine}}(m,M)}}\), the worst perturbation occurs at the upper endpoint \((1+\varepsilon)M\), which gives the bound in Theorem 9.

We remark that for small error \(\varepsilon\), the slowdown factor in Theorem 9 is of the order \[\begin{align} \frac{1}{1 + \tilde{\varepsilon} - \sqrt{(1 + \tilde{\varepsilon})^2 - 1}} = 1 + \sqrt{2 \tilde{\varepsilon}} + \Theta(\tilde{\varepsilon}) = 1 + \sqrt{\frac{4M\varepsilon}{M-m}} + \Theta(\varepsilon)\,. \label{eq:slowdown-1} \end{align}\tag{27}\] In other words, for large condition number \(\kappa \gg 1\), the number of iterations required by GD increases by a multiplicative factor of roughly \(1 + 2\sqrt{\varepsilon}\) when using inexact gradients. Note that this slowdown is more than the \((1 + \Theta(\varepsilon))\) slowdown for constant-stepsize schedules [51]. This is a manifestation of the general phenomenon that accelerated algorithms are often less stable to model mis-specification, see e.g., [53].

The above discussion analyzes the effect of inexact gradients for the i.i.d.stepsize schedule that is optimized for exact gradients; one can also ask for the i.i.d.stepsize schedule that is optimized for this setting of inexact gradients. This follows from a minor adaptation of our techniques.

The key idea is that the inexact gradients effectively worsen the parameters for strong convexity (from \(m\) to \(\tilde{m}\)) and smoothness (from \(M\) to \(\tilde{M}\)). Theorem 10 then follows by appealing to our main result, Theorem 2, for this enlarged range. Proof details are deferred to Appendix 8.6. By a Taylor expansion in \(\varepsilon\), the resulting convergence rate \(\tilde{R}_{\mathrm{acc}} := \tfrac{\sqrt{\tilde{\kappa}}-1}{\sqrt{\tilde{\kappa}}+1}\) has a slowdown factor of \[\begin{align} \frac{\tilde{R}_{\mathrm{acc}}}{\mathop{\mathrm{\mathrm{R_{acc}}}}} = 1 + \frac{2\sqrt{\kappa}}{\kappa-1}\varepsilon+ O(\varepsilon^2) = 1 + \frac{2\sqrt{Mm}}{M-m}\varepsilon+ O(\varepsilon^2)\,. \label{eq:slowdown-2} \end{align}\tag{28}\] This \(1+\Theta(\varepsilon)\) slowdown achieved by re-optimizing the stepsizes improves over the \(1+\Theta(\sqrt{\varepsilon})\) slowdown for using the original \(\mathop{\mathrm{{\mathrm{Arcsine}}(m,M)}}\) stepsizes, c.f., 27 .

5.4 Game-theoretic connections and lower bounds↩︎

This paper considers random stepsizes for the algorithmic purpose of accelerating GD. It turns out that this randomized perspective is also helpful for the dual problem of constructing algorithmic lower bounds. We describe this duality here.

Recall from §2 that the optimal distribution \(\mu\) over i.i.d.inverse stepsizes \(\beta = \alpha^{-1}\) has the following extremal characterization: \[\begin{align} \inf_{\mu \in \mathcal{M}(E)} \; \sup_{\lambda \in E}\; \mathbb{E}_{\beta \sim \mu} \log \left\lvert 1 - \beta^{-1}\lambda\right\rvert\,. \label{eq:games:1} \end{align}\tag{29}\] Above, \(\mu\) is constrained to \(E := [m,M]\) without loss of generality (since this holds at optimality), \(\lambda \in E\) should be interpreted as the curvature of a “worst-case” univariate quadratic function \(f(x) = \tfrac{\lambda}{2}x^2\), and \(\mathbb{E}_{\beta \sim \nu} \log |1 - \beta^{-1} \lambda|\) is the logarithm of the asymptotic convergence rate.

To bridge this problem of accelerating GD with the problem of constructing lower bounds, we consider a standard two-player zero-sum game where players choose their actions simultaneously and independently. In our case, one player chooses inverse stepsizes \(\beta^{-1}\) and the other player chooses a quadratic function (parameterized by its curvature \(\lambda)\). By lifting to mixed (i.e., randomized) strategies, one player chooses a distribution \(\mu\) over stepsizes, the other chooses a distribution \(\nu\) over quadratic functions, with payoff given by \({\cal R}(\mu,\nu) := \mathbb{E}_{\beta \sim \mu, \, \lambda \sim \nu} \log \left\lvert 1 - \beta^{-1}\lambda\right\rvert\).

The key observation is that due to the specific payoff function, the game is essentially symmetric between the two players (modulo inversion), and as a consequence, so are the optimal mixed strategies for both players. Below, let \(\mathop{\mathrm{\mathrm{FlippedArcsine}(m,M)}}\) denote the law of \(\lambda\) when \(\lambda^{-1} \sim \mathrm{Arcsine}(M^{-1},m^{-1})\).

As a corollary, Lemma 6 gives a new proof that \(\tfrac{\sqrt{\kappa} - 1}{\sqrt{\kappa} + 1}\) is the fastest possible convergence rate for GD. This result is classically due to [15], who proved this by constructing a single, infinite-dimensional “worst-case function” for which any Krylov-subspace algorithm converges slowly. Here we construct a “worst-case distribution” over univariate functions for which any GD stepsize schedule converges slowly in expectation. By the probabilistic method, for any specific GD stepsize schedule, we can then conclude the existence of a single univariate function for which convergence is slow. For simplicity, we state this lower bound here for i.i.d.stepsizes; this can be extended to arbitrary GD stepsize schedules and also to non-asymptotic lower bounds, see [18].

6 Conclusion and future work↩︎

This paper shows that randomizing stepsizes provides an alternative mechanism for accelerated convex optimization beyond the special case of quadratics. This raises several natural questions about the possible benefit of randomness for GD—or said equivalently, whether the Arcsine Stepsize Schedule can be derandomized. The current results differ on two fronts: Theorem 2 shows that random stepsizes achieve the fully accelerated rate \(O(\kappa^{1/2} \log \tfrac{1}{\varepsilon})\) on separable functions, whereas [19] shows that deterministic stepsize achieves the partially accelerated “Silver” rate \(O(\kappa^{\log_{1+\sqrt{2}} 2} \log \tfrac{1}{\varepsilon}) \approx O(\kappa^{0.78} \log \tfrac{1}{\varepsilon})\) on general convex functions. These two sources of differences motivate the questions:

• Is it possible to de-randomize Theorem 2, i.e., construct a deterministic stepsize schedule which achieves the fully accelerated rate for separable functions? If so, a natural conjecture would be some ordering of the Chebyshev stepsize schedule.

• Does the Arcsine Stepsize Schedule fully accelerate on general (non-separable) convex functions? If not, can any randomized stepsize schedule improve over the Silver Convergence Rate for general convex functions?

There are also many questions about extensions to other optimization settings of interest. For other settings, what is the optimal stepsize distribution if not Arcsine? In any setting, is there a provable benefit of using non-i.i.d.stepsizes? Can random stepsizes be alternatively understood through discretization analyses of continuous-time gradient flows, or via the continuized framework of [59]? While other mechanisms for acceleration (e.g., momentum) are now well understood [10], little is known about acceleration via randomness.

Acknowledgements.

We are grateful to Abraham Wyner for suggesting connections to repeated gambling, discussed in §5.1.

7 Connections to logarithmic potential theory↩︎

Recall from the overview section §2 that the optimal distribution \(\mu\) for inverse stepsizes \(\beta\) is given by the extremal problem \[\begin{align} \inf_{\mu \in \mathcal{M}(E)} \;\sup_{\lambda \in E} \;\mathbb{E}_{\beta \sim \mu} \log \left\lvert 1 - \beta^{-1} \lambda \right\rvert\,, \label{eq:extremal-app} \end{align}\tag{30}\] where \(E := [m,M]\) and \(\mathcal{M}(E)\) denotes the probability distributions supported on \(E\).² Indeed, for quadratic optimization, this extremal problem exactly characterizes the asymptotic convergence rate (see §2.1); and for the more general setting of separable convex optimization, even though the relevant extremal problem is somewhat different, the solution is identical (see §2.2).

In this appendix section, we connect this extremal problem over stepsize distributions with the extremal problem for electrostatic capacitors that is classically studied in potential theory. As we explain below, although these two problems seek probability distributions minimizing different objectives, the two problems are intimately linked because they have the same optimality condition: the optimal distribution \(\mu\) must satisfy the equalization property from Lemma 2 (i.e., the objective \(\mathbb{E}_{\beta \sim \mu} \log |1 - \beta^{-1}\lambda|\) has equal value for all \(\lambda \in [m,M]\)). Because this equalization property uniquely characterizes the distribution (a classical fact from potential theory, see Lemma 11), this means that the two problems have the same unique optimal solution. In this way, the equalization property plays the central role connecting these two problems.

This appendix section is organized as follows. In §7.1 we recall the electrostatics problem and highlight the syntactic similarities and differences to the stepsize problem. In §7.2 we explain the equalization property in the context of both problems and how it arises as the optimality condition. This is the key link between the two problems. In §7.3 we exploit this connection concretely by using the machinery of potential theory to solve the stepsize problem. This provides alternative proofs of Lemmas 1 and 2 that are significantly shorter than the integral calculations of [11], [13].

7.1 Extremal problem for electrostatic capacitors↩︎

The extremal problem central to logarithmic potential theory is \[\begin{align} \inf_{\mu \in \mathcal{M}(E)} \iint \log \frac{1}{\left\lvert\lambda - \omega\right\rvert} d\mu(\omega) d\mu(\lambda)\,, \label{eq:extremal-pot} \end{align}\tag{31}\] given a set \(E \subset \mathbb{C}\). This problem is sometimes called the electrostatics problem for a capacitor due to its physical interpretation: how should you place a distribution of repelling particles on \(E\) so as to minimize the total repulsive energy? Here the repulsive energy is measured via \[I(\mu) := \iint \log \frac{1}{\left\lvert\lambda - \omega\right\rvert} d\mu(\omega) d\mu(\lambda) = \mathbb{E}_{\lambda \sim \mu} \phi_{\mu}(\lambda)\,,\] which is the \(\mu\)-average of the logarithmic potential energy at \(\lambda\), denoted \[\phi_{\mu}(\lambda) := \mathbb{E}_{\omega \sim \mu} \log \frac{1}{|\lambda - \omega|}\,.\]

Standard arguments based on lower semi-continuity and strict convexity of the objective ensure that an optimal distribution exists and is unique. For certain sets \(E\), this optimal distribution can be computed in closed form. For example it is well-known that for the interval \(E = [m,M]\), the unique optimal distribution is the \(\mathop{\mathrm{{\mathrm{Arcsine}}(m,M)}}\) distribution, see e.g., [44].

Syntactic similarities and differences.

Let \(E = [m,M]\). Then the stepsize extremal problem 30 clearly bears similarities to the electrostatics extremal problem 31 , namely \[\begin{align} \inf_{\mu \in \mathcal{M}(E)} \sup_{\lambda \in E} \phi_{\mu}(0) - \phi_{\mu}(\lambda) \qquad \text{versus} \qquad \inf_{\mu \in \mathcal{M}(E)} \mathbb{E}_{\lambda\sim \mu} \phi_{\mu}(\lambda)\,. \end{align}\] Yet there are two differences. First, the stepsize problem evaluates the logarithmic potential \(\phi_{\mu}(\lambda)\) at a specific point \(\lambda \in [m,M]\), whereas the electrostatics problem averages this potential over \(\lambda \sim \mu\). Second, the stepsize problem has an additional term of \(\phi_{\mu}(0)\).

In general, changing the objective or constraints in an extremal problem can significantly alter the optimal solution. Yet, remarkably, the optimal solution \(\mu\) is identical for both problems. This is because both problems have the same optimality condition, described next.

7.2 Equalization property↩︎

Here we describe how the equalization property arises as the optimality condition in both problems.

7.2.1 Equalized potential for the electrostatics problem↩︎

A celebrated fact—sometimes called the Fundamental Theorem of Potential Theory [42]—is that the optimal distribution \(\mu\) to the electrostatics problem 31 satisfies an equalization property: the potential \(\lambda \mapsto \phi_{\mu}(\lambda)\) is constant on \(E\). (Note that this is equivalent to the constancy on \(E\) of the objective for the stepsize problem, namely \(\lambda \mapsto \mathbb{E}_{\beta} \log |1 - \beta^{-1}\lambda| = \phi_{\mu}(0) - \phi_{\mu}(\lambda)\).)

This equalization property has an intuitive physical interpretation: if the potential is not constant on \(E\), then particles will flow to points of \(E\) with lower potential, so the current state is not at equilibrium and therefore cannot be optimal. This is why the optimal solution \(\mu\) is often called the equilibrium distribution. This result is attributed to Frostman; here we state a simple version of it that suffices for our setting, namely a compact interval \(E = [m,M]\), although it extends to much more general sets \(E \subset \mathbb{C}\) [42].

A proof can be found in any textbook on logarithmic potential theory, e.g., [42], [43], [60]. Below we use calculus of variations to give an informal back-of-the-envelope calculation that provides intuition for how this equalization property arises as the optimality condition. Let \[\begin{align} F(\mu) := \iint \log \frac{1}{|x-y|} d\mu(x) d\mu(y) + \gamma \left( \int d\mu(x) - 1\right) \end{align}\] denote the original objective plus a Lagrangian term for the normalization constraint on \(\mu\). Then a quick calculation shows that the first variation is \[\begin{align} \delta F(\mu,h) &= \frac{d}{d\varepsilon} F(\mu+ \varepsilon h) \big|_{\varepsilon= 0} = \int \left[ \gamma + 2 \int \log \frac{1}{|x-y|} d\mu(x) \right] dh(y)\,. \end{align}\] In order for this first variation to vanish along any direction \(h\), it must hold that \[\begin{align} \gamma + 2 \int \log \frac{1}{|x-y|} d\mu(x) = 0, \qquad \forall y \in E\,. \end{align}\] This is the equalization property.

7.2.2 Equalized convergence rate for the stepsize problem↩︎

As discussed in §2.1 and formally stated in Lemma 1, the stepsize problem has the same equalization property at optimality. This is rigorously proven in §7.3.2 using the machinery of potential theory; to provide intuition for how this arises as the optimality condition, here we give an informal back-of-the-envelope argument based on the game-theoretic perspective developed in §5.4.

Suppose \(\mu\) is a solution to the stepsize problem 30 . By the calculation in §5.4, there is a distribution \(\nu\) (the Flipped Arcsine distribution) for which \[\begin{align} \mathbb{E}_{\beta \sim \mu,\, \lambda \sim \nu} \log \left\lvert 1 - \beta^{-1}\lambda\right\rvert = \log \mathop{\mathrm{\mathrm{R_{acc}}}}\,. \end{align}\] Observe that this is the \(\nu\)-average of \(P_{\mu}(\lambda) := \mathbb{E}_{\beta \sim \mu} \log |1 - \beta^{-1}\lambda| = \phi_{\mu}(0) - \phi_{\mu}(\lambda)\). Thus, if \(\mu\) does not equalize, then (forgoing measure-theoretic technicalities) there exists \(\lambda \in \operatorname{supp}(\nu) = [m,M]\) for which \(P_{\mu}(\lambda)\) is strictly larger than the \(\nu\)-average. That is, \(\max_{\lambda \in [m,M]} P_{\mu}(\lambda) > \log \mathop{\mathrm{\mathrm{R_{acc}}}}\). But this contradicts the optimality of \(\mu\) since the value of the stepsize problem 30 is \(\log \mathop{\mathrm{\mathrm{R_{acc}}}}\). Therefore the equalization property must be satisfied by the optimal distribution \(\mu\).

7.3 Solution to the extremal problem over stepsizes↩︎

So far we have focused on developing conceptual connections between the stepsize problem and electrostatics problem. We now exploit these connections concretely by using the machinery of potential theory to solve the stepsize problem. This provides short proofs of Lemmas 1 and 2.

7.3.1 Proof of Lemma 2↩︎

Here we prove the equalization property. In the optimization literature, this was proven by computing the integral directly using trigonometric identities [11], or by decomposing the logarithmic integrand into the Chebyshev polynomial basis [13]. Here we mention an alternative approach by appealing to the following well-known identity from potential theory, which gives the potential for the equilibrium distribution on \([-1,1]\) at any point \(z\); see e.g., [43].

7.3.2 Proof of Lemma 1↩︎

Here we prove that the Arcsine distribution is the unique optimal solution for the stepsize problem 30 by mirroring the standard argument from potential theory that the Arcsine distribution is the unique optimal solution to the electrostatic problem 31 , see e.g., [44]. Specifically, our approach is to show that any optimal distribution \(\mu\) satisfies the equalization property (see §7.2.2 for a heuristic argument of this). Indeed, it then follows that the optimal distribution \(\mu\) is unique and is the Arcsine distribution by the converse of Frostman’s Theorem (stated below), which ensures that the Arcsine distribution is the unique distribution satisfying equalization.

We begin with a simple observation.

We now show that any optimal distribution, constrained to this interval, must satisfy the equalization property. For this it is helpful to recall several classical facts from potential theory. These can be found, e.g., in page 171 and Theorems 3.2, 2.8, 2.2, and 3.1, respectively, of [44].

8 Deferred proofs↩︎

8.1 Alternative intuitive proof of Lemma 2 without potential theory↩︎

For pedagogical reasons and to provide further intuition for the equalizing property (Lemma 2), here we present a short semi-rigorous proof that uses only elementary calculus (c.f., the rigorous proof using potential theory in Appendix 7). Begin by re-parameterizing: \(\beta \sim \mathop{\mathrm{{\mathrm{Arcsine}}(m,M)}}\) has distribution \(\beta = c + d \cos \theta\) where \(\theta \sim \mathrm{Unif}(0,\pi)\), \(c := \tfrac{M+m}{2}\), \(d := \tfrac{M-m}{2}\). The desired integral is then \[F(\lambda) := \frac{1}{\pi} \int_0^{\pi} \log \left|1 - \frac{\lambda}{c + d \cos \theta}\right| d\theta\,.\]

Step 1: equalization.

To prove \(F\) is constant on \([m,M]\), we argue that \(F'\) is zero on \((m,M)\). Forgoing technical details about swapping the derivative and integral, \(F'\) is calculated to be \[\begin{align} F'(\lambda) = \frac{1}{\pi} \int_0^{\pi} \frac{d}{d\lambda} \log \left|1 - \frac{\lambda}{c + d \cos \theta}\right| d\theta = - \frac{1}{\pi} \mathrm{PV} \int_0^{\pi} \frac{d\theta}{c + d\cos \theta - \lambda} \,. \end{align}\] Here the principal value arises from the simple pole at the value of \(\theta\) where \(c+ d \cos \theta = \lambda\). This vanishes since it is well-known³ that, more generally, \[\begin{align} \mathrm{PV} \int_0^{\pi} \frac{d\theta}{a + b\cos \theta} = 0 \end{align}\] for any \(|a| < |b|\) and \(b> 0\). (For us, \(a = c - \lambda\) and \(b = d\).)

Step 2: value.

Since \(F\) is constant, we can evaluate at the midpoint \(\lambda = c = (M+m)/2\). Compute \[\begin{align} F(\lambda) &= \frac{1}{\pi} \int_0^{\pi} \log \Big| 1 - \frac{c}{c + d \cos \theta} \Big| d \theta = \frac{1}{\pi} \int_0^{\pi} \Big[ \log d + \log |\cos \theta| - \log (c + d\cos \theta) \Big] d\theta \\ &= \log d - \log 2 - \log \frac{c + \sqrt{c^2-d^2}}{2} = \log \frac{\sqrt{\kappa}-1}{\sqrt{\kappa}+1}\,. \end{align}\] Above, the first step is by definition of \(F\), the second step is by expanding the logarithm, and the final step is by collecting the logarithms and simplifying. It remains to justify the three integral calculations in the penultimate step. The first integral is trivial as the integrand \(\log d\) is constant in \(\theta\); the second and third are due to the standard integral formulas [61].

8.2 Equivalence of separability and commutative Hessians↩︎

Here we prove an alternative characterization of separability. We remark that this provides a simple way to certify that a function is not separable: it suffices to exhibit points \(x\) and \(y\) for which the Hessians \(\nabla^2 f(x)\) and \(\nabla^2 f(y)\) do not commute.

Proof of (i) \(\Longrightarrow\) (ii). The function \(h(y) := f(U^Ty) = \sum_{i=1}^d f_i([y]_i)\) is separable, hence its Hessian \(\nabla^2 h\) is diagonal. The chain rule implies \(\nabla^2 f = U^T \nabla^2 h U\), thus the Hessians of \(f\) are simultaneously diagonalizable in the basis \(U\), hence they commute.

Proof of (ii) \(\Longrightarrow\) (i). Since the Hessians commute, they are simultaneously diagonalizable, i.e., there exists an orthogonal matrix \(U\) such that \(\nabla^2 f(x) = U^T D(x) U\) for all \(x\), where \(D(x)\) is diagonal. Defining \(h(y) := f(U^Ty)\), it follows by the chain rule that \(\nabla^2 h(y) = D(U^Ty)\) is diagonal. Thus all mixed partial derivatives vanish, hence \(\tfrac{\partial}{\partial [y]_i}h(y)\) is a function only of \([y]_i\). Integrating implies \(h(y) = \sum_{i=1}^d f_i([y]_i)\) for some functions \(f_i\). Thus \(f(x) = h(Ux) = \sum_{i=1}^d f_i([Ux]_i)\) as desired.

8.3 Proof of Lemma 4↩︎

We first isolate a useful identity: the random variables \(Z_t = \log |1 - \beta_t^{-1} \lambda_t|\) defined in §3.1 are uncorrelated—despite not being independent. Intuitively, this is due to the equalization property of the Arcsine distribution (from which \(\beta_t^{-1}\) are independently drawn) since it ensures that each \(Z_t\) has the same (conditional) expectation regardless of the dependencies injected through the non-stationary process \(\lambda_t\).

8.4 Details for §5.1↩︎

Here we observe that random stepsizes can make GD unstable in alternative notions of convergence. Instability in the following lemma occurs for \(\kappa > 4\) since then \(\sqrt{\kappa} - 1 > 1\). See §5.1 for a discussion of this phenomenon and why the notion of convergence studied in the rest of this paper is more relevant to optimization practice.

To prove this, we make use of the following integral representation of the geometric mean of two positive scalars \(a, b > 0\) in terms of the weighted harmonic mean \(H_{t}(a,b) := (t a^{-1} + (1-t )b^{-1})^{-1}\). This identity is well-known (e.g., it follows from [62], or directly from [61] by taking \(q = 1/b\) and \(p = 1/b-1/a\)).

This lets us compute the expectation of the proposed random stepsizes.

We are now ready to prove Lemma 14.

8.5 Proof of Theorem 8↩︎

For shorthand, let \(V := \sup_{\lambda \in [m,M]} \mathop{\mathrm{\mathrm{Var}}}_{\beta \sim \mathop{\mathrm{{\mathrm{Arcsine}}(m,M)}}} \left( \log \left\lvert 1 - \beta^{-1} \lambda\right\rvert \right)\) denote a uniform upper bound on the variance of any \(Z_t\). By the integral calculations in [13], \[\begin{align} V \leqslant\pi^2 + o_{\kappa}(1)\,. \end{align}\] Thus because \(\mathop{\mathrm{\mathrm{Cov}}}(Z_s, Z_t) = 0\) for \(s \neq t\) by Lemma 13, \[\begin{align} \mathop{\mathrm{\mathrm{Var}}}\left( \frac{1}{n} \sum_{t=0}^{n-1} Z_t \right) = \frac{1}{n^2}\sum_{t=0}^{n-1} \mathop{\mathrm{\mathrm{Var}}}\left( Z_t \right) \leqslant \frac{V}{n} \leqslant\frac{\pi^2 + o_{\kappa}(1)}{n}\,. \end{align}\] Thus by Chebyshev’s inequality and the equalization property of the Arcsine distribution (Lemma 2), \[\begin{align} \mathbb{P}\left( \frac{1}{n} \sum_{t=0}^{n-1} \log|1 - \beta_t^{-1} \lambda_t| \geqslant\log \mathop{\mathrm{\mathrm{R_{acc}}}}+ \varepsilon\right) \leqslant\frac{\pi^2 + o_{\kappa}(1)}{n \varepsilon^2}\,, \end{align}\] where \(\beta_t\) are i.i.d.\(\mathop{\mathrm{{\mathrm{Arcsine}}(m,M)}}\) and \(\lambda_t \in [m,M]\) is any process adapted to the filtration \(\sigma(\{\beta_s\}_{s < t})\). Thus by 18 , the average convergence rate \(R_{n,i} := (\tfrac{| [x_{n}]_i - [x^*]_i | }{ | [x_0]_i - [x^*]_i|})^{1/n}\) for the function \(f_i\) satisfies \[\begin{align} \mathbb{P}\left( R_{n,i} \geqslant e^{\varepsilon} \, \mathop{\mathrm{\mathrm{R_{acc}}}}\right) \leqslant\frac{\pi^2 + o_{\kappa}(1)}{n \varepsilon^2} \end{align}\] for any coordinate \(i \in S := \{i \in [d] : [x_0]_i \neq [x^*]_i\}\), i.e., for any \(i\) for which the initialization is not already optimal. Thus the overall convergence rate \(R_n := (\tfrac{\|x_{n} - x^*\| }{ \|x_0 - x^*\|})^{1/n}\) for \(f\) satisfies \[\begin{align} \mathbb{P}\left(R_n \geqslant e^{\varepsilon} \, \mathop{\mathrm{\mathrm{R_{acc}}}}\right) \leqslant \mathbb{P}\left(\max_{i \in S} R_{n,i} \geqslant e^{\varepsilon} \, \mathop{\mathrm{\mathrm{R_{acc}}}}\right) \leqslant \frac{d(\pi^2 + o_{\kappa}(1))}{n \varepsilon^2}\,. \end{align}\] Above the first step is because \(R_n \leqslant\max_{i \in S} R_{n,i}\) (since \(R_n^2 = \tfrac{\sum_{i \in S} |[x_{n}]_i - [x^*]_i|^2}{\sum_{i \in S} |[x_{0}]_i - [x^*]_i|^2} \leqslant\max_{i \in S} \tfrac{|[x_{n}]_i - [x^*]_i|^2}{ |[x_{0}]_i - [x^*]_i|^2} = \max_{i \in S} R_{n,i}^2\)). The second step is by a union bound and the previous display. The proof is then complete by setting \(\delta\) equal to the right hand side of the above display.

8.6 Details for §5.3↩︎

References↩︎

---

1. A subtle but important technical issue is that one cannot simply use the Law of Large Numbers to derive the asymptotic convergence rate from the \(n\)-step rate (c.f., 13 for how that worked in the quadratic case). Indeed, in the convex case here, the relevant empirical sum \(\frac{1}{n} \sum_{t=0}^{n-1} \log |1 - \beta_t^{-1} \lambda_t|\) has non-i.i.d.summands: they are neither identically distributed (since \(\lambda_t\) is time-varying) nor independent (since \(\lambda_t\) can depend on previous stepsizes). Briefly, we overcome the identically distributed issue by using the equalization property of the Arcsine distribution (Lemma 2) to show that these summands still have the desired expectation, and we overcome the independence issue by using a martingale-based argument to show that the fresh randomness in \(\alpha_t\) makes these summands “independent enough” to prove a specially tailored version of Kolmogorov’s Strong Law of Large Numbers for this setting. Details in §3.1.↩︎

2. Although the stepsize problem was not originally written with this constraint in §2, it is easily seen that this constraint is satisfied at optimality and therefore does not change the extremal problem. See Observation 17.↩︎

3. This can be proved with the substitution \(t = \tan(\theta/2)\). Indeed, after some algebra, this substitution gives \(\mathrm{PV} \int_0^{\pi} \frac{d\theta}{a + b\cos \theta} = \mathrm{PV}\int_0^{\infty} \frac{2dt}{P-Qt^2}\) where \(P := a+b\) and \(Q := b-a\) and the simple pole is now at \(t_0 = \sqrt{\tfrac{P}{Q}}\). The antiderivative \(\int \frac{dt}{P-Qt^2} = \frac{1}{2\sqrt{PQ}} \log |\tfrac{\sqrt{P} + \sqrt{Q}t}{\sqrt{P}-\sqrt{Q}t}| := g(t)\), hence the PV integral is \(\lim_{\varepsilon\searrow 0} (g(t_0 - \varepsilon) - g(0) + g(\infty) - g(t_0 + \varepsilon)) = 0\).↩︎