Hybrid Stochastic Gradient Descent Algorithms for Stochastic Nonconvex Optimization

TL;DR

This paper introduces a hybrid stochastic gradient estimator for nonconvex optimization, achieving improved complexity bounds over traditional methods and extending to various algorithmic variants.

Contribution

It proposes a novel hybrid SARAH-SGD algorithm with better complexity bounds for nonconvex problems and extends the approach to multiple variants.

Findings

Achieves $O( ext{max}\{rac{\sigma^3}{\varepsilon}, rac{\sigma}{\varepsilon^3} ight)$-complexity bound.

Outperforms state-of-the-art SGD in certain variance regimes.

Demonstrates effectiveness on multiple datasets with nonconvex models.

Abstract

We introduce a hybrid stochastic estimator to design stochastic gradient algorithms for solving stochastic optimization problems. Such a hybrid estimator is a convex combination of two existing biased and unbiased estimators and leads to some useful property on its variance. We limit our consideration to a hybrid SARAH-SGD for nonconvex expectation problems. However, our idea can be extended to handle a broader class of estimators in both convex and nonconvex settings. We propose a new single-loop stochastic gradient descent algorithm that can achieve $O(\max\{\sigma^3\varepsilon^{-1},\sigma\varepsilon^{-3}\})$-complexity bound to obtain an $\varepsilon$-stationary point under smoothness and $\sigma^2$-bounded variance assumptions. This complexity is better than $O(\sigma^2\varepsilon^{-4})$ often obtained in state-of-the-art SGDs when $\sigma < O(\varepsilon^{-3})$. We also consider different extensions of our method, including constant and adaptive step-size with single-loop, double-loop, and mini-batch variants. We compare our algorithms with existing methods on several datasets using two nonconvex models.

Equations

$$\min_{x\in\mathbb{R}^{p}}\Big{\{}f(x):=\mathbb{E}_{\xi}\left[f(x;\xi)\right]\Big{\}},$$

$${\mathbb{E}}_{\xi }\left[{\Vert \nabla f(x;\xi )-\nabla f(y;\xi )\Vert }^2\right]\le L^2\Vert x-y{\Vert }^2,\forall x,y\in {\mathbb{R}}^p.\text{vspace}-1.5ex$$

$${\mathbb{E}}_{\xi }\left[\Vert \nabla f(x;\xi )-\nabla f(x){\Vert }^2\right]\le {\sigma }^2,\forall x\in {\mathbb{R}}^p.$$

$$v_t:={\beta }_{t-1}{v}_{t-1}+{\beta }_{t-1}(\nabla f(x_t;{\xi }_t)-\nabla f(x_{t-1};{\xi }_t))+(1-{\beta }_{t-1})u_t,$$

$$v_t:={\beta }_{t-1}\nabla f(x_t;{\xi }_t)+(1-{\beta }_{t-1})u_t+{\beta }_{t-1}(v_{t-1}-\nabla f(x_{t-1};{\xi }_t)).$$

$${\mathbb{E}}_{({\xi }_t,{\zeta }_t)}\left[v_t\right]=\nabla f(x_t)+{\beta }_{t-1}(v_{t-1}-\nabla f(x_{t-1})).$$

$$\begin{array}[]{ll}\mathbb{E}_{(\xi_{t},\zeta_{t})}\left[\|v_{t}-\nabla{f}(x_{t})\|^{2}\right]&=\beta_{t-1}^{2}\|v_{t-1}-\nabla{f}(x_{t-1})\|^{2}-\beta_{t-1}^{2}\|\nabla{f}(x_{t-1})-\nabla{f}(x_{t})\|^{2}\vspace{1ex}\\ &+{~{}}\beta_{t-1}^{2}\mathbb{E}_{\xi_{t}}\left[\|\nabla{f}(x_{t};\xi_{t})-\nabla{f}(x_{t-1};\xi_{t})\|^{2}\right]\vspace{1ex}\\ &+{~{}}(1-\beta_{t-1})^{2}\mathbb{E}_{\zeta_{t}}\left[\|u_{t}-\nabla{f}(x_{t})\|^{2}\right].\end{array}$$

$$\mathrm{Bias}(v_t)=\Vert {\mathbb{E}}_{({\xi }_t,{\zeta }_t)}\left[v_t-\nabla f(x_t)\mid {\mathcal{F}}_t\right]\Vert ={\beta }_{t-1}\Vert v_{t-1}-\nabla f(x_{t-1})\Vert \le \Vert v_{t-1}-\nabla f(x_{t-1})\Vert .$$

$$\mathbb{E}\left[\Vert v_t-\nabla f(x_t){\Vert }^2\right]\le {\omega }_t\mathbb{E}\left[\Vert v_0-\nabla f(x^0){\Vert }^2\right]+L^2\sum _{i=0}^{t-1}{\omega }_{i,t}\mathbb{E}\left[\Vert x_{i+1}-x_i{\Vert }^2\right]+S_t,$$

$$\eta :=\frac{2}{L(\sqrt{1+4{\alpha }_m^2}+1)}\text{with}\beta :=1-\frac{c_1}{\sqrt{b(m+1)}}\text{and}{\alpha }_m^2:=\frac{{\beta }^2(1-{\beta }^{2m})}{1-{\beta }^2}.$$

$$\mathbb{E}\left[\|\nabla{f}(\widetilde{x}_{m})\|^{2}\right]\leq\frac{3b^{1/4}L\big{[}f(x^{0})-f^{\star}\big{]}}{\sqrt{c_{1}}(m+1)^{3/4}}+\left(c_{1}+\frac{1}{c_{1}}\right)\frac{\sigma^{2}}{\sqrt{b(m+1)}}.$$

$$m:=\bigg{\lfloor}\tfrac{\sigma}{\varepsilon^{3}}\Big{[}\tfrac{3Lc_{2}^{1/4}}{\sqrt{c_{1}}}\big{[}f(x^{0})-f^{\star}\big{]}+\left(c_{1}+\tfrac{1}{c_{1}}\right)\tfrac{1}{\sqrt{c_{2}}}\Big{]}^{3/2}\bigg{\rfloor}=\mathcal{O}\left(\frac{\sigma}{\varepsilon^{3}}\right).$$

$${\!\!\!\!\!\!\!\!}\begin{array}[]{ll}\mathcal{T}_{ge}&:=\dfrac{\sigma^{3}}{\varepsilon}\left[3Lc_{2}^{9/4}\big{[}f(x^{0})-f^{\star}\big{]}+2c_{2}^{3/2}\right]^{1/2}{\!\!\!}+\dfrac{3\sigma}{\varepsilon^{3}}\left[3Lc_{2}^{1/4}\big{[}f(x^{0})-f^{\star}\big{]}+\frac{2}{\sqrt{c_{2}}}\right]^{3/2}\vspace{1ex}\\ &=\mathcal{O}\left(\dfrac{\sigma^{3}}{\varepsilon}+\dfrac{\sigma}{\varepsilon^{3}}\right).\end{array}{\!\!\!\!\!\!\!\!\!}$$

$$\eta_{m}:=\frac{1}{L},~{}~{}\text{and}~{}~{}\eta_{t}:=\frac{1}{L+L^{2}\big{[}\beta^{2}\eta_{t+1}+\beta^{4}\eta_{t+2}+\cdots+\beta^{2(m-t)}\eta_{m}\big{]}}~{}~{}\text{for}~{}t=0,\cdots,m-1.$$

$$\frac{\eta }{2}\sum _{t=0}^m\mathbb{E}\left[\Vert \nabla f(x_t^{(s)}){\Vert }^2\right]\le \mathbb{E}\left[f(x_0^{(s)})\right]-\mathbb{E}\left[f(x_{m+1}^{(s)})\right]+\frac{\eta {\sigma }^2\sqrt{m+1}}{(1+\beta )\sqrt{b}}.$$

$$\frac{1}{S(m+1)}\displaystyle\sum_{s=1}^{S}\sum_{t=0}^{m}\mathbb{E}\left[\|\nabla{f}(x_{t}^{(s)})\|^{2}\right]\leq\frac{3Lb^{1/4}}{S(m+1)^{3/4}}\big{[}f(\widetilde{x}^{0})-f^{\star}\big{]}+\frac{2\sigma^{2}}{\sqrt{b(m+1)}}.$$

$$S:=\Bigg{\lfloor}\frac{3Lc_{1}^{1/4}\big{[}f(\widetilde{x}^{0})-f^{\star}\big{]}}{c_{2}^{3/4}\sigma\left(1-\frac{2}{\sqrt{c_{1}c_{2}}}\right)\varepsilon}\Bigg{\rfloor}~{}~{}~{}\text{outer iterations}.$$

$$\mathcal{T}_{ge}:=(b+3m)S=\frac{3L(c_{1}+3c_{2})c_{1}^{1/4}\big{[}f(\widetilde{x}^{0})-f^{\star}\big{]}\sigma}{c_{2}^{3/4}\left(1-\frac{2}{\sqrt{c_{1}c_{2}}}\right)\varepsilon^{3}}=\mathcal{O}\left(\frac{\sigma}{\varepsilon^{3}}\right).$$

$${\hat{v}}_t:={\beta }_{t-1}{\hat{v}}_{t-1}+\frac{{\beta }_{t-1}}{{\hat{b}}_t}\sum _{i\in {\hat{\mathcal{B}}}_t}(\nabla f(x_t;{\xi }_i)-\nabla f(x_{t-1};{\xi }_i))+(1-{\beta }_{t-1})u_t,$$

$$\eta :=\frac{2}{L\left(1+\sqrt{1+4\rho {\alpha }_m^2}\right)}\text{with}{\alpha }_m^2:=\frac{{\beta }^2(1-{\beta }^{2m})}{1-{\beta }^2}\text{and}\beta :=1-\frac{c_1}{\sqrt{\rho b(m-1)}}.$$

$$m:=\lfloor \frac{{\rho }^{1/2}\sigma }{{\epsilon }^3}{\left[\frac{3L{c}_2^{1/4}}{2c_1}\left(f(x^0)-f^{\star }\right)+\left(c_1+\frac{1}{c_1}\right)\frac{1}{2\sqrt{c_2}}\right]}^{3/2}\rfloor =\mathcal{O}\left(\frac{{\rho }^{1/2}\sigma }{{\epsilon }^3}\right).$$

$${\mathcal{T}}_{ge}:=\mathcal{O}\left(\frac{{\rho }^{1/2}{\sigma }^3}{\epsilon }+\frac{\sigma }{{\rho }^{1/2}{\epsilon }^3}\right),$$

$$\min_{x\in\mathbb{R}^{p}}\bigg{\{}f(x):=\frac{1}{n}\sum_{i=1}^{n}\Big{[}f_{i}(x):=\log(1+\exp(-a_{i}^{T}x))+\lambda\sum_{j=1}^{p}\frac{x_{i}^{2}}{1+x_{i}^{2}}\Big{]}\bigg{\}},\vspace{-0.5ex}$$

$$f^{\star}:=\min_{x\in\mathbb{R}^{p}}\Big{\{}f(x):=\frac{1}{n}\sum_{i=1}^{n}\ell(a_{i}^{\top}x,b_{i})+\frac{\lambda}{2}\|x\|^{2}\Big{\}},\vspace{-0.5ex}$$

$$\eta_{m}:=\frac{1}{L},~{}~{}~{}\text{and}~{}~{}\eta_{t}:=\frac{1}{L+L^{2}\big{[}\omega\eta_{t+1}+\omega^{2}\eta_{t+2}+\cdots+\omega^{(m-t)}\eta_{m}\big{]}},$$

$$0<{\eta }_0<{\eta }_1<\cdots <{\eta }_m=\frac{1}{L},\text{and}{\Sigma }_m:=\sum _{t=0}^m{\eta }_t\ge \frac{(m+1)\sqrt{1-\omega }}{2L}.$$

$$\begin{array}[]{ll}(m-t)\big{(}\omega\eta_{t+1}+\omega^{2}\eta_{t+2}+\cdots+\omega^{m-t}\eta_{m}\big{)}&\leq\big{(}\sum_{j=t+1}^{m}\eta_{i}\big{)}\left(\omega+\omega^{2}+\cdots+\omega^{m-t}\right)\vspace{1ex}\\ &\leq\frac{\omega}{1-\omega}\big{(}\sum_{j=t+1}^{m}\eta_{i}\big{)}.\end{array}$$

$$\left\{\begin{array}[]{ll}L^{2}\eta_{0}(\omega\eta_{1}+\omega^{2}\eta_{2}+\cdots+\omega^{m}\eta_{m})&=1-L\eta_{0}\\ L^{2}\eta_{1}(\omega\eta_{2}+\omega^{2}\eta_{3}+\cdots+\omega^{m-1}\eta_{m})&=1-L\eta_{1}\\ \cdots&\cdots\\ L^{2}\eta_{m-1}\omega\eta_{m}&=1-L\eta_{m-1}\\ 0&=1-L\eta_{m}.\end{array}\right.$$

$$\left\{\begin{array}[]{ll}\frac{\omega L^{2}}{1-\omega}\eta_{0}(\eta_{0}+\eta_{1}+\cdots+\eta_{m})&\geq m-mL\eta_{0}+\frac{\omega L^{2}}{1-\omega}\eta_{0}^{2}\\ \frac{\omega L^{2}}{1-\omega}\eta_{1}(\eta_{0}+\eta_{1}+\cdots+\eta_{m})&\geq(m-1)-(m-1)L\eta_{1}+\frac{\omega L^{2}}{1-\omega}(\eta_{1}\eta_{0}+\eta_{1}^{2})\\ \cdots&\cdots\\ \frac{\omega L^{2}}{1-\omega}\eta_{m-1}(\eta_{0}+\eta_{1}+\cdots+\eta_{m})&\geq 1-L\eta_{m-1}+\frac{\omega L^{2}}{1-\omega}(\eta_{m-1}\eta_{0}+\cdots+\eta_{m-1}^{2})\\ \frac{\omega L^{2}}{1-\omega}\eta_{m}(\eta_{0}+\eta_{1}+\cdots+\eta_{m})&\geq 1-L\eta_{m}+\frac{\omega L^{2}}{1-\omega}(\eta_{m}\eta_{0}+\cdots+\eta_{m}^{2}).\end{array}\right.$$

$$\frac{\omega L^{2}}{1-\omega}\Sigma_{m}^{2}\geq\frac{m^{2}+m+2}{2}-L(m\eta_{0}+(m-1)\eta_{1}+\cdots+\eta_{m-1}+\eta_{m})+\frac{\omega L^{2}}{2(1-\omega)}\big{(}S_{m}+\Sigma_{m}^{2}\big{)}.$$

Full text

1 Introduction

Consider the following stochastic nonconvex optimization problem of the form:

[TABLE]

where $f(\cdot;\cdot):\mathbb{R}^{p}\times\Omega\to\mathbb{R}$ is a stochastic function defined such that for each $x\in\mathbb{R}^{p}$, $f(x;\cdot)$ is a random variable in a given probability space $(\Omega,\mathbb{P})$, while for each realization $\xi\in\Omega$, $f(\cdot;\xi)$ is smooth on $\mathbb{R}^{p}$; and $\mathbb{E}_{\xi}\left[f(x;\xi)\right]$ is the expectation of $f(x;\xi)$ w.r.t. $\xi$ over $\Omega$.

Our goals and assumptions:

Since (1) is nonconvex, our goal in this paper is to develop a new class of stochastic gradient algorithms to find an $\varepsilon$-approximate stationary point $\widetilde{x}_{T}$ of (1) such that $\mathbb{E}\left[\|\nabla{f}(\widetilde{x}_{T})\|^{2}\right]\leq\varepsilon^{2}$ under mild assumptions as stated in Assumption 1.1.

Assumption 1.1**.**

The objective function $f$ of (1) satisfies the following conditions:

• $\mathrm{(a)}$

(Boundedness from below) There exists a finite lower bound $f^{\star}:=\inf_{x\in\mathbb{R}^{p}}f(x)>-\infty$.

• $\mathrm{(b)}$

($L$-average smoothness) The function $f(\cdot;\xi)$ is $L$-average smooth on $\mathbb{R}^{p}$, i.e. there exists $L\in(0,+\infty)$ such that

[TABLE]

• $\mathrm{(c)}$

(Bounded variance) There exists $\sigma\in(0,\infty)$ such that

[TABLE]

These assumptions are very standard in stochastic optimization methods [10, 14]. The $L$-average smoothness of $f$ is weaker than the smoothness of $f$ for each realization $\xi\in\Omega$. Note that our methods described in the sequel are also applicable to the finite-sum problem $\min_{x}\left\{f(x)=\frac{1}{n}\sum_{i=1}^{n}f_{i}(x)\right\}$ as long as the above assumptions hold. However, we do not specify our methods to solve this problem. In this case, $\sigma$ in (3) can be replaced by other alternatives, e.g., $\sigma_{n}^{2}:=\frac{1}{n}\sum_{i=1}^{n}\left[\|\nabla{f}_{i}(x)\|^{2}-\|\nabla{f}(x)\|^{2}\right]$.

Our key idea:

Different from existing methods, we introduce a convex combination of a biased and unbiased estimator of the gradient $\nabla{f}$ of $f$, which we call a hybrid stochastic gradient estimator. While the biased estimator exploited in this paper is SARAH in [16], the unbiased one can be any unbiased estimator. SARAH is a recursive biased and variance reduced estimator for $\nabla{f}$. Combining it with an unbiased estimator allows us to reduce the bias and variance of the hybrid estimator. In this paper, we only focus on the standard stochastic estimator as an unbiased candidate.

Related work:

Under Assumption 1.1, problem (1) covers a large number of applications in machine learning and data sciences. The stochastic gradient descent (SGD) method was first studied in [25], and becomes extremely popular in recent years. [14] seems to be the first work showing the convergence rates of robust SGD variants in the convex setting, while [15] provides an intensive complexity analysis for many optimization algorithms, including stochastic methods. Variance reduction methods have also been widely studied, see, e.g. [2, 6, 8, 11, 20, 22, 26, 27, 29]. In the nonconvex setting, [10] seems to be the first algorithm achieving $\mathcal{O}\left(\sigma^{2}\varepsilon^{-4}\right)$-complexity bound. Other researchers have also made significant progress in this direction, including [3, 4, 5, 9, 13, 17, 18, 19, 23, 28, 31]. A majority of these works, including [3, 4, 5, 13, 23], rely on SVRG estimator in order to obtain better complexity bounds. Hitherto, the complexity of SVRG-based methods remains worse than the best-known results, which is obtained in [9, 21, 28] via the SARAH estimator. However, as discussed in [21, 28], the method called SPIDER in [9, 13] does not practically perform well due to small step-size and its dependence on the reciprocal of the estimator’s norm. [28] amends this issue by using a large constant step-size, but requires large mini-batch and does not consider the single sample case and single loop variants. [21] provides a more general framework to treat composite problems where it covers (1) as special case, but it does not consider the single loop as in SGDs.

Our contribution:

To this end, our contribution can be summarized as follows:

• $\mathrm{(a)}$

We propose a hybrid stochastic estimator for a stochastic gradient of a nonconvex function $f$ in (1) by combining the SARAH estimator from [16] and any unbiased stochastic estimator such as SGD and SVRG. However, we only focus on the SGD estimator in this paper. We prove some key properties of this hybrid estimator that can be used to design new algorithms.

• $\mathrm{(b)}$

We exploit our hybrid estimator to develop a single-loop SGD algorithm that can achieve an $\varepsilon$-stationary point $\widetilde{x}_{m}$ such that $\mathbb{E}\left[\|\nabla{f}(\widetilde{x}_{m})\|^{2}\right]\leq\varepsilon^{2}$ in at most $\mathcal{O}\left(\sigma\varepsilon^{-3}+\sigma^{3}\varepsilon^{-1}\right)$ stochastic gradient evaluations. This complexity significantly improves $\mathcal{O}\left(\sigma^{2}\varepsilon^{-4}\right)$ of SGD if $\sigma<\mathcal{O}\left(\varepsilon^{-3}\right)$. We extend our algorithm to a double loop variant, which requires $\mathcal{O}\left(\max\left\{\sigma\varepsilon^{-3},\sigma^{2}\varepsilon^{-2}\right\}\right)$ stochastic gradient evaluations. This is the best-known complexity in the literature for stochastic gradient-type methods for solving (1).

• $\mathrm{(c)}$

We also investigate other variants of our method, including adaptive step-sizes, and mini-batches. In all these cases, our methods achieve the best-known complexity bounds.

Let us emphasize the following points of our contribution. Firstly, although our single-loop method requires three gradients per iteration compared to standard SGDs, it can achieves better complexity bound. Secondly, it can be cast into a variance reduction method where it starts from a “good” approximation $v_{0}$ of $\nabla{f}(x^{0})$, and aggressively reduces the variance. Thirdly, our step-size is $\eta=\mathcal{O}\left(m^{-1/3}\right)$ which is larger than $\eta=\mathcal{O}\left(m^{-1/2}\right)$ in SGDs. Fourthly, the step-size of the adaptive variant is increasing instead of diminishing as in SGDs. Finally, our method achieves the same best-known complexity as in variance reduction methods studied in [9, 21, 28]. We believe that our approach can be extended to other estimators such as SVRG [11] and SAGA [8], and can be used for Hessians to develop second-order methods as well as to solve convex and composite problems.

Paper organization:

The rest of this paper is organized as follows. Section 2 introduces our new hybrid stochastic estimator for the gradient of $f$ and investigates its properties. Section 3 proposes a single-loop hybrid SGD-SARAH algorithm and its complexity analysis. It also considers a double-loop and mini-batch variants with rigorous complexity analysis. Section 4 provides two numerical examples to illustrate our methods and compares them with state-of-the-art methods. All the proofs and additional experiments can be found in the Supplementary Document.

Notation: We work with Euclidean spaces, $\mathbb{R}^{p}$, equipped with standard inner product $\langle\cdot,\cdot\rangle$ and norm $\left\|\cdot\right\|$. For a smooth function $f$ (i.e., $f$ is continuously differentiable), $\nabla{f}$ denotes its gradient. We use $\mathbf{U}_{\mathbf{p}}(\mathcal{S})$ to denote a distribution on $\mathcal{S}$ with probability $\mathbf{p}$. If $\mathbf{p}$ is uniform, then we simply use $\mathbf{U}(\mathcal{S})$. We also use $\mathcal{O}\left(\cdot\right)$ to present big-O notion in complexity theory, and $\sigma(\cdot)$ to denote a $\sigma$-field.

2 Hybrid stochastic gradient estimators

In this section, we propose new stochastic estimators for the gradient of a smooth function $f$.

Let $u_{t}$ be an unbiased estimator of $\nabla{f}(x_{t})$ formed by a realization $\zeta_{t}$ of $\xi$, i.e. $\mathbb{E}_{\zeta_{t}}\left[u_{t}\right]=\nabla{f}(x_{t})$. We attempt to develop the following stochastic estimator for $\nabla{f}(x_{t})$ in (1):

[TABLE]

where $\xi_{t}$ and $\zeta_{t}$ are two independent realizations of $\xi$ on $\Omega$. Clearly, if $\beta_{t}=0$, then we obtain a simple unbiased stochastic estimator, and $\beta_{t}=1$, we obtain the SARAH estimator in [16]. We are interested in the case $\beta_{t}\in(0,1)$, in which we call $v_{t}$ in (4) a hybrid stochastic estimator.

Note that we can rewrite $v_{t}$ as

[TABLE]

The first two terms are two stochastic gradients estimated at $x_{t}$, while the third term is the difference $v_{t-1}-\nabla{f}(x_{t-1};\xi_{t})$ of the previous estimator and a stochastic gradient at the previous iterate. Here, since $\beta_{t-1}\in(0,1)$, the main idea is to exploit more recent information than the old ones. In fact, the hybrid estimator $v_{t}$ covers many other estimators, including SGD, SVRG, and SARAH. We can use one of the following two concrete unbiased estimators $u_{t}$ of $\nabla{f}(x_{t})$ as follows:

• •

The SGD estimator: $u_{t}:=u_{t}^{\mathrm{sgd}}=\nabla{f}(x_{t};\zeta_{t})$.

• •

The SVRG estimator: $u_{t}:=u_{t}^{\mathrm{svrg}}=\nabla{f}(\tilde{x})+\nabla{f}(x_{t};\zeta_{t})-\nabla{f}(\tilde{x};\zeta_{t})$, where $\nabla{f}(\tilde{x})$ is a full gradient evaluated at a given snapshot point $\tilde{x}$.

However, for the sake of presentation, we only focus on the SGD estimator $u_{t}:=u_{t}^{\mathrm{sgd}}$.

We first prove the following property of the estimator $v_{t}$ showing how the variance is estimated.

Lemma 2.1**.**

Let $v_{t}$ be defined by (4). Then

[TABLE]

If $\beta_{t-1}\neq 0$, then $v_{t}$ is a biased estimator. Moreover, we have

[TABLE]

Remark 2.1**.**

From (4), we can see that $v_{t}$ remains a biased estimator as long as $\beta_{t-1}\in(0,1]$. Its biased term is

[TABLE]

This shows that the bias $v_{t}$ estimator is smaller than the one in the SARAH estimator $v_{t}^{\textrm{sarah}}:=v_{t-1}^{\textrm{sarah}}+\nabla{f}(x_{t};\xi_{t})-\nabla{f}(x_{t-1};\xi_{t})$ from [16], which is $\mathrm{Bias}(v^{\mathrm{sarah}}_{t})=\|v_{t-1}^{\mathrm{sarah}}-\nabla{f}(x_{t-1})\|$.

The following lemma bounds the second moment of $v_{t}-\nabla{f}(x_{t})$ with $v_{t}$ defined in (4).

Lemma 2.2**.**

Assume that $f(\cdot,\cdot)$ is $L$-smooth and $u_{t}$ is an SGD estimator. Then, we have the following upper bound on the variance $\mathbb{E}\left[\|v_{t}-\nabla{f}(x_{t})\|^{2}\right]$ of $v_{t}$:

[TABLE]

where the expectation is taking over all the randomness $\mathcal{F}_{t}:=\sigma(v_{0},v_{1},\cdots,v_{t})$, $\omega_{t}:=\prod_{i=1}^{t}\beta_{i-1}^{2}$, $\omega_{i,t}:=\prod_{j=i+1}^{t}\beta_{j-1}^{2}$ for $i=0,\cdots,t$, and S_{t}:=\sum_{i=0}^{t-1}\big{(}\prod_{j=i+2}^{t}\beta_{j-1}^{2}\big{)}(1-\beta_{i})^{2}\sigma_{i+1}^{2} for $t\geq 0$.

Lemmas 2.1 and 2.2 provides two key properties to develop stochastic algorithm in Section 3.

3 Hybrid SARAH-SGD algorithms

In this section, we utilize our hybrid stochastic estimator $v_{t}$ in (4) to develop stochastic gradient methods for solving (1). We consider three different variants using the hybrid SARAH-SGD estimator.

3.1 The generic algorithm framework

Using $v_{t}$ defined by (4), we can develop a new algorithm for solving (1) as in Algorithm 1.

Algorithm 1 looks essentially the same as any SGD scheme with only one loop. The differences are at Step 3 with a mini-batch estimator $v_{0}$ and at Step 6, where we use our hybrid gradient estimator $v_{t}$. In addition, we will show in the sequel that it uses different step-sizes and leads to different variants. Unlike the inner loop of SARAH or SVRG, each iteration of Algorithm 1 requires three individual gradient evaluations instead of two as in these methods. The snapshot at Step 3 of Algorithm 1 relies on a mini-batch $\mathcal{B}$ of the size $b$, which is independent of $(\xi_{t},\zeta_{t})$ in the loop $t$.

3.2 Convergence analysis

We analyze two cases: constant step-size and adaptive step-size. In both cases, $\beta_{t}$ is fixed for all $t$.

3.2.a Convergence of Algorithm 1 with constant step-size $\eta$ and constant $\beta$

Assume that we run Algorithm 1 within $m$ iterations $m\geq 1$. In this case, given $0<c_{1}<\sqrt{b(m+1)}$, we choose $\eta$ and $\beta$ in Algorithm 1 as follows:

[TABLE]

The following theorem estimates the complexity of Algorithm 1 to approximate an $\varepsilon$-stationary point of (1), whose proof is given in Subsection 2.2 of the supplementary document.

Theorem 3.1**.**

Let $\{x_{t}\}$ be the sequence generated by Algorithm 1 using the step-size $\eta$ defined by (8). Let us choose $\widetilde{x}_{m}\sim\mathbf{U}(\{x_{t}\}_{t=0}^{m})$. Then

• $(\mathrm{a})$

The step-size $\eta$ satisfies \eta\geq\underline{\eta}:=\frac{2\sqrt{c_{1}}}{3{L}\big{[}b(m+1)\big{]}^{1/4}}. In addition, we have

[TABLE]

• $(\mathrm{b})$

If we choose $b:=c_{2}\sigma^{8/3}(m+1)^{1/3}$ for any $c_{2}>0$, then to guarantee $\mathbb{E}\left[\|\nabla{f}(\widetilde{x}_{m})\|^{2}\right]\leq\varepsilon^{2}$, we need to choose

[TABLE]

In particular, if we choose $c_{1}=1$, then the number of oracle calls is $\mathcal{T}_{ge}$ is

[TABLE]

Moreover, the step-size $\eta$ satisfies $\eta\geq\underline{\eta}:=\frac{2}{3{L}c_{2}^{1/4}\sigma^{2/3}(m+1)^{1/3}}=\mathcal{O}\left(m^{-1/3}\right)$.

Here, $\mathcal{T}_{ge}$ stands for the number of stochastic gradient evaluations of $f$ in (1). The complexity $\mathcal{T}_{ge}$ in (11) can be written as $\mathcal{T}_{ge}=\mathcal{O}\left(\max\left\{\sigma\varepsilon^{-3},\sigma^{3}\varepsilon^{-1}\right\}\right)$. If $\sigma<\mathcal{O}\left(\frac{1}{\varepsilon}\right)$, then our complexity is $\mathcal{T}_{ge}=\mathcal{O}\left(\sigma\varepsilon^{-3}\right)$. Even if $\sigma<\mathcal{O}\left(\frac{1}{\varepsilon^{3}}\right)$, then our complexity is still better than $\mathcal{O}\left(\sigma^{2}\varepsilon^{-4}\right)$ in SGD.

3.2.b Convergence of Algorithm 1 with adaptive step-size $\eta_{t}$ and constant $\beta$

Let $\beta:=1-\frac{c_{1}}{\sqrt{b(m+1)}}\in(0,1)$ be fixed for some $0<c_{1}<\sqrt{b(m+1)}$. Instead of fixing step-size $\eta_{t}$ as in (8), we can update it adaptively as

[TABLE]

It can be shown that $0<\eta_{0}<\eta_{1}<\cdots<\eta_{m}$. Interestingly, our step-size is updated in an increasing manner instead of diminishing as in existing SGD-type methods. Moreover, given $m$, we can pre-compute the sequence of these step-sizes $\left\{\eta_{t}\right\}_{t=0}^{m}$ in advance within $\mathcal{O}\left(m\right)$ basic operations. Therefore, it does not significantly incur the computational cost of our method.

The following theorem states the convergence of Algorithm 1 under the adaptive update (12), whose proof is given in Subsection 2.3 of the supplementary document.

Theorem 3.2**.**

Let $\{x_{t}\}$ be the sequence generated by Algorithm 1 using the step-size $\eta_{t}$ defined by (12). Let $\Sigma_{m}:=\sum_{t=0}^{m}\eta_{t}$, and $\widetilde{x}_{m}\sim\mathbf{U}_{\mathbf{p}}(\{x_{t}\}_{t=0}^{m})$ with $\mathbf{p}_{t}:=\mathbb{P}\left(\widetilde{x}_{m}=x_{t}\right)=\frac{\eta_{t}}{\Sigma_{m}}$. Then

• $(\mathrm{a})$

The sum $\Sigma_{m}$ is bounded from below as $\Sigma_{m}\geq\frac{\sqrt{c_{1}}(m+1)^{3/4}}{2Lb^{1/4}}$.

• $(\mathrm{b})$

If we choose $b:=c_{2}\sigma^{8/3}(m+1)^{1/3}$ for any $c_{2}>0$, then to guarantee $\mathbb{E}\left[\|\nabla{f}(\widetilde{x}_{m})\|^{2}\right]\leq\varepsilon^{2}$, we need to choose m:=\Big{\lfloor}\frac{\sigma}{\varepsilon^{3}}\Big{[}\frac{3Lc_{2}^{1/4}}{\sqrt{c_{1}}}\big{[}f(x^{0})-f^{\star}\big{]}+\big{(}c_{1}+\frac{1}{c_{1}}\big{)}\frac{1}{\sqrt{c_{2}}}\Big{]}^{3/2}\Big{\rfloor}=\mathcal{O}\left(\frac{\sigma}{\varepsilon^{3}}\right). Therefore, the number of stochastic gradient evaluations $\mathcal{T}_{ge}$ is at most the same as in (11).

Note that in the finite sum case, i.e. $|\Omega|=n$, we set $b:=\min\{n,c_{2}\sigma^{8/3}(m+1)^{1/3}\}$ in both Theorems 3.1 and 3.2. This complexity remains the same as in Theorem 3.1. However, the adaptive stepsize $\eta_{t}$ potentially gives a better performance in practice as we will see in Section 4.

Algorithm 1 can be considered as a single-loop variance reduction method, which is similar to SAGA [8], but Algorithm 1 aims at solving the nonconvex problem (1). It is different from standard SGD methods, where it can be initialized by a mini-batch and then update the estimator using three individual gradients. Therefore, it has the same cost as SGD with mini-batch of size $3$. As a compensation, we obtain an improvement on the complexity bound as in Theorems 3.1 and 3.2.

3.3 Convergence analysis of the double loop variant

Since the step-size $\eta_{t}$ depends on $m$, it is natural to run Algorithm 1 with multiple stages. This leads to a double-loop algorithm as SVRG, SARAH, and SPIDER, where Algorithm 1 is restarted at each outer iteration $s$. The detail of this variant is described in Algorithm 2.

To analyze Algorithm 2, we use $x^{(s)}_{t}$ to represent the iterate of Algorithm 1 at the $t$-th inner iteration within each stage $s$. From (45), we can see that each stage $s$, the following estimate holds

[TABLE]

Here, we assume that we fix the step-size $\eta_{t}=\eta>0$ for simplicity of analysis. The complexity of Algorithm 2 is given in the following theorem, whose proof is in Supplementary Document 2.4.

Theorem 3.3**.**

Let $\{x^{(s)}_{t}\}_{t=0\to m}^{s=1\to S}$ be the sequence generated by Algorithm 2 using constant step-size $\eta$ in (8). Then, the following estimate holds

[TABLE]

Let $\widetilde{x}_{T}\sim\mathbf{U}(\{x^{(s)}_{t}\}_{t=0\to m}^{s=1\to S})$. If we choose $b:=\frac{c_{1}\sigma^{2}}{\varepsilon^{2}}$ and $m+1:=\frac{c_{2}\sigma^{2}}{\varepsilon^{2}}$ for some constants $c_{1}>0$ and $c_{2}>0$ and $c_{1}c_{2}>4$, then, to guarantee $\mathbb{E}\left[\|\nabla{f}(\widetilde{x}_{T})\|^{2}\right]\leq\varepsilon^{2}$, we require at most

[TABLE]

Consequently, the total number of stochastic gradient evaluations $\mathcal{T}_{ge}$ does not exceed

[TABLE]

Note that the complexity (15) only holds if $\mathcal{O}\left(\frac{\sigma}{\varepsilon^{3}}\right)>\frac{c_{1}\sigma^{2}}{\varepsilon^{2}}$. Otherwise, the total complexity is $\mathcal{O}\left(\max\left\{\frac{\sigma}{\varepsilon^{3}},\frac{\sigma^{2}}{\varepsilon^{2}}\right\}\right)$, where other constants independent of $\sigma$ and $\varepsilon$, and are hidden. Practically, if $\beta$ is very close to $1$, one can remove the unbiased SGD term to save one stochastic gradient evaluation. In this case, our estimator reduces to SARAH but using different step-size. We observed empirically that when $\beta\approx 0.999$, the performance of our methods is not affected if we do so.

3.4 Extensions to mini-batch cases

We consider a mini-batch hybrid stochastic estimator $\hat{v}_{t}$ for the gradient $\nabla{f}(x_{t})$ defined as:

[TABLE]

where $\beta_{t-1}\in[0,1]$, and $\hat{\mathcal{B}}_{t}$ is a mini-batch of the size $\hat{b}_{t}$ and independent of the unbiased estimator $u_{t}$. Note that $u_{t}$ can also be a mini-batch unbiased estimator. For example, $u_{t}:=\frac{1}{\tilde{b}_{t}}\sum_{j\in\tilde{\mathcal{B}}_{t}}\nabla{f}(x_{t};\zeta_{j})$ is a mini-batch SGD estimator with a mini-batch $\tilde{\mathcal{B}}_{t}$ of size $\tilde{b}_{t}$, where $\tilde{\mathcal{B}}_{t}$ is independent of $\hat{\mathcal{B}}_{t}$.

Using $\hat{v}_{t}$ defined by (16), we can design a mini-batch variant of Algorithms 1 to solve (1). The following corollary is obtained as a result of Theorems 3.1 for the mini-batch variant of Algorithm 1, whose proof is in Subsection 3.2 of the supplementary document.

Corollary 3.1**.**

Let Algorithm 1 be applied to solve (1) using mini-batch update (16) for $v_{t}$ with $\hat{b}_{t}=\tilde{b}_{t}=\hat{b}\geq 1$ fixed, $0<c_{1}<\sqrt{b(m+1)}$, and the step-size

[TABLE]

If we choose $b:=c_{2}\sigma^{8/3}\left[\rho(m+1)\right]^{1/3}$ for any $c_{2}>0$, then to guarantee $\mathbb{E}\left[\|\nabla{f}(\widetilde{x}_{m})\|^{2}\right]\leq\varepsilon^{2}$, we need to choose

[TABLE]

Therefore, the number of oracle calls is $\mathcal{T}_{ge}$ is

[TABLE]

where $\rho=\rho(\hat{b}):=\frac{n-\hat{b}}{(n-1)\hat{b}}$ if $n:=|\Omega|$ is finite, and $\rho(\hat{b}):=\frac{1}{\hat{b}}$, otherwise. In particular, if we choose $\hat{b}:=\frac{\varepsilon^{2}\sigma^{2}}{c_{3}^{2}}$ for some $0<c_{3}\leq\varepsilon\sigma$, then, the overall complexity $\mathcal{T}_{ge}$ is $\mathcal{T}_{ge}:=\mathcal{O}\left(\left(c_{3}+\frac{1}{c_{3}}\right)\frac{\sigma^{2}}{\varepsilon^{2}}\right)$.

We can also develop a mini-batch variant of Algorithm 2 and estimate its complexity as in Theorem 3.3. For more details, we refer to Subsection 3.3 in the supplementary document due to space limit.

4 Numerical experiments

We verify our algorithms on two numerical examples and compare them with several existing methods: SVRG [24], SVRG+ [12], SPIDER [9], SpiderBoost [28], and SGD [10]. Due to space limit, the detailed configuration of our experiments as well as more numerical experiments can be found in Supplementary Document D. Our numerical experiments are implemented in Python and running on a MacBook Pro. Laptop with 2.7GHz Intel Core i5 and 16Gb memory.

4.1 Logistic regression with nonconvex regularizer

Our first example is the following well-known problem used in may papers including [28]:

[TABLE]

where $a_{i}\in\mathbb{R}^{p}$ are given for $i=1,\cdots,n$, and $\lambda>0$ is a regularization parameter. Clearly, problem (20) fits (1) well with $L_{f_{i}}=\frac{\|A\|^{2}}{4}+2\lambda$. In this experiment, we choose $\lambda=0.1$ and normalize the data. One can also verify Assumption 1.1 due to the bounded Hessian of $f_{i}$.

We use three datasets from LibSVM for (20): w8a ($n=49,749,p=300$), rcv1.binary ($n=20,242,p=47,236$), and real-sim $(n=72,309,p=20,958$). We run $8$ different algorithms as follows. Algorithm 1 with constant step-size (Hybrid-SGD-SL) and adaptive step-size (Hybrid-SGD-ASL) using our theoretical step-sizes (8) and (12), respectively without tuning. Hybrid-SGD-DL is Algorithm 2. SGD1 is SGD with constant step-size $\eta_{t}=\frac{0.1}{L}$, and SGD2 is SGD with adaptive step-size $\eta_{t}=\frac{0.1}{L(1+\lfloor t/n\rfloor)}$. Since the stepsize of SPIDER depends on an accuracy $\varepsilon$, we choose $\varepsilon=10^{-1}$ to get a larger step-size. Our first result in the single-sample case (i.e. when $\hat{b}=1$, not using mini-batch) is plotted in Fig. 1 after $20$ epochs.

From Fig. 1, we observe that Hybrid-SGD-SL has similar convergence behavior as SGD1, but Hybrid-SGD-ASL works better. Hybrid-SGD-DL is the best but has some oscillation. SGD2 works better than SGD1 and is comparable with Hybrid-SGD-SL/ASL in the two last datasets. SVRG performs very poorly due to its small step-size. SVRG+ works much better than SVRG, and is comparable with our methods. SPIDER is also slow even when we have increased its step-size.

Now, we run 3 single-loop algorithms with mini-batch of the size $\hat{b}:=300$. The result is shown in Fig. 2 after $20$ epochs.

Fig. 2 shows similar performance between Hybrid-SGD-SL and ASL and SGD2. Clearly, these theoretical variants of Algorithm 1 are slightly better than the adaptive SGD variant (SGD2), where a careful step-size is used.

4.2 Binary classification involving nonconvex loss and Tikhonov’s regularizer

We consider the following binary classification problem studied in [30] involving nonconvex loss:

[TABLE]

where $a_{i}\in\mathbb{R}^{p}$ and $b_{i}\in\left\{-1,1\right\}$ are given data for $i=1,\cdots,n$, $\lambda>0$ is a regularization parameter, and $\ell$ is a nonconvex loss of the forms: $\ell(\tau,s)=\left(1-\frac{1}{1+\exp(-\tau s)}\right)^{2}$ (using in two-layer neural networks). One can check that (21) satisfies Assumption 1.1 with $L\approx 0.15405\max_{i}\|a_{i}\|^{2}+\lambda$. We choose $\lambda:=0.01$, and test three variants of Algorithm 2: Hybrid-SGD-DL and compare them with SpiderBoost, SVRG, and SVRG+. Due to space limit, we only plot one experiment in Fig. 3 after $20$ epochs. Additional experiments can be found in Supplementary Document D.

As we can see from Fig. 3 that Algorithm 2 performs well and is slightly better than SpiderBoost. Note that SpiderBoost simply uses SARAH estimator with constant stepsize $\eta=\frac{1}{2L}$ but with mini-batch of the size $\lfloor\sqrt{n}\rfloor$. It is not surprise that SpiderBoost makes very good progress to decrease the gradient norms. Both SVRG and SVRG+ perform much worse than Hybrid-SGD-DL and SpiderBoost in this test, but SVRG+ is slightly better than SVRG. In our methods, due to the aid of SARAH part, they also make similar progress as SpiderBoost but using different step-sizes.

5 Conclusion

We have introduced a new hybrid SARAH-SGD estimator for the objective gradient of expectation optimization problems. Under standard assumptions, we have shown that this estimator has a better variance reduction property than SARAH. By exploiting such an estimator, we have developed a new Hybrid-SGD algorithm, Algorithm 1, that has better complexity bounds than state-of-the-art SGDs. Our algorithm works with both constant and adaptive step-sizes. We have also studied its double-loop and mini-batch variants. We believe that our approach can be extended to other choices of unbiased estimators, Hessian estimators for second-order stochastic methods, and adaptive $\beta$.

Appendix A Appendix: Properties of the hybrid stochastic estimator

This supplementary document provides the full proof of all the results in the main text. First, we need the following lemma in the sequel.

Lemma A.1**.**

Given $L>0$ and $\omega\in(0,1)$. Let $\{\eta_{t}\}_{t=0}^{m}$ be the sequence updated by

[TABLE]

for $t=0,\cdots,m-1$. Then

[TABLE]

Proof.

First, from (22) it is obvious to show that $0<\eta_{0}<\cdots<\eta_{m-1}=\frac{1}{L(1+\omega)}<\eta_{m}=\frac{1}{L}$. At the same time, since $\omega\in(0,1)$, we have $1\geq\omega\geq\omega^{2}\geq\cdots\geq\omega^{m}$. By Chebyshev’s sum inequality, we have

[TABLE]

From the update (22), we also have

[TABLE]

Using (24) into (25), we get

[TABLE]

Let $\Sigma_{m}:=\sum_{t=0}^{m}\eta_{t}$ and $S_{m}:=\sum_{t=0}^{m}\eta_{t}^{2}$. Summing up both sides of the above inequalities, we get

[TABLE]

Using again Chebyshev’s sum inequality, we have

[TABLE]

Note that $(m+1)S_{m}\geq\Sigma_{m}^{2}$ by Cauchy-Schwarz’s inequality, which shows that S_{m}+\Sigma_{m}^{2}\geq\big{(}\frac{m+2}{m+1}\big{)}\Sigma_{m}^{2}. Combining three last inequalities, we obtain the following quadratic inequation in $\Sigma_{m}$

[TABLE]

Solving this inequation with respect to $\Sigma_{m}>0$, we obtain

[TABLE]

Since $\omega\in(0,1)$, we can overestimate this as $\Sigma_{m}\geq\frac{(m+1)\sqrt{1-\omega}}{2L}$, which proves (23). ∎

1.1 The proof of Lemma 2.1: Properties of the hybrid SARAH estimator

By taking the expectation of both sides in (4) and using the fact that $\xi_{t}$ and $\zeta_{t}$ are independent, we can easily obtain (5).

To prove (6), we first write

[TABLE]

In this case, we have

[TABLE]

Let us first take expectation w.r.t. $\xi_{t}$ conditioned on $\zeta_{t}$ to obtain

[TABLE]

Now, taking the expectation over $\zeta_{t}$, and noting that $\mathbb{E}_{(\xi_{t},\zeta_{t})}\left[\cdot\right]=\mathbb{E}_{\zeta_{t}}\left[\mathbb{E}_{\xi_{t}}\left[\cdot\mid\zeta_{t}\right]\right]$ and $\mathbb{E}_{\zeta_{t}}\left[u_{t}-\nabla{f}(x_{t})\right]=0$, we get

[TABLE]

which is exactly (6). $\square$

1.2 The proof of Lemma 2.2: Bound on the variance of the hybrid estimator

We first upper bound (6) by using $\sigma_{t}^{2}:=\mathbb{E}_{\zeta_{t}}\left[\|u_{t}-\nabla{f}(x_{t})\|^{2}\right]$ and then taking the full expectation over $\mathcal{F}_{t}:=\sigma(v_{0},v_{1},\cdots,v_{t})$ as

[TABLE]

If we define $a_{t}^{2}:=\mathbb{E}\left[\|v_{t}-\nabla{f}(x_{t})\|^{2}\right]$, then the above inequality can lead to

[TABLE]

Denote $b_{t-1}^{2}:=\mathbb{E}\left[\|x_{t}-x_{t-1}\|^{2}\right]$. Then, we have from the last inequality that

[TABLE]

By induction, this inequality implies

[TABLE]

Here, we use a convention that $\prod_{i=t+1}^{t}\beta_{i}^{2}=1$. As a result, it can be rewritten in a compact form as

[TABLE]

Define $\omega_{t}:=\prod_{i=1}^{t}\beta_{i-1}^{2}$, $\omega_{i,t}:=\prod_{j=i+1}^{t}\beta_{j-1}^{2}$, and S_{t}:=\sum_{i=0}^{t-1}s_{i}=\sum_{i=0}^{t-1}\big{(}\prod_{j=i+2}^{t}\beta_{j-1}^{2}\big{)}(1-\beta_{i})^{2}\sigma^{2}_{i+1} with s_{i}:=(1-\beta_{i})^{2}\sigma^{2}_{i+1}\big{(}\prod_{j=i+2}^{t}\beta_{j-1}^{2}\big{)}. Then, we can rewrite (26) as

[TABLE]

which is exactly (7). $\square$

Appendix B Appendix: Convergence analysis of Algorithm 1 and Algorithm 2

We provide the full convergence analysis for Algorithm 1 and Algorithm 2 in the single-sample case.

2.1 The proof of Lemma B.1: One-iteration analysis

The following lemma provides a key estimate for convergence analysis of Algorithm 1.

Lemma B.1**.**

Let $\{x_{t}\}$ be the sequence generated by Algorithm 1. Then, under Assumption 1.1, we have the following estimate:

[TABLE]

where

[TABLE]

and $\omega_{t}$, $\omega_{i,t}$, and $S_{t}$ are defined in Lemma 2.2.

Proof.

First, from the $L$-smoothness of $f$, we have

[TABLE]

Taking the expectation over the randomness $(\xi_{t},\zeta_{t})$ of this estimate, we obtain

[TABLE]

Taking the full expectation over the entire history up to the $t$-th iteration, and then using (7) and noting that $x_{t}-x_{t-1}=-\eta_{t-1}v_{t-1}$, we obtain

[TABLE]

where $q_{t}^{2}:=\mathbb{E}\left[\|v_{t}\|^{2}\right]$ and $a_{t}^{2}:=\mathbb{E}\left[\|v_{t}-\nabla{f}(x_{t})\|^{2}\right]$. Here, we use $b_{t-1}^{2}:=\mathbb{E}\left[\|x_{t}-x_{t-1}\|^{2}\right]=\eta_{t-1}^{2}\mathbb{E}\left[\|v_{t-1}\|^{2}\right]=\eta_{t-1}^{2}q_{t-1}^{2}$ in the last inequality.

Summing up (29) from $t=0$ to $t=m$, we obtain

[TABLE]

Let us define $T_{m}$ as in (28), i.e.:

[TABLE]

Then, we obtain from (30) the estimate (27). ∎

2.2 The proof of Theorem 3.1: Single-loop with constant step-size

We analyze the case $\beta_{t}=\beta\in(0,1)$ fixed and the step-size $\eta_{t}=\eta>0$ fixed. From Lemma 2.2, we have $\omega_{t}=\beta^{2t}$, $\omega_{i,t}=\beta^{2(t-i)}$, and

[TABLE]

In this case, by convention that $\omega_{0}=1$, we have

[TABLE]

Now, to bound the quantity $\mathcal{T}_{m}$ defined by (28), we note that

[TABLE]

Using this expression, we can write $\mathcal{T}_{m}$ from (28) as

[TABLE]

To guarantee $\mathcal{T}_{m}\leq 0$, from (32), we need to choose

[TABLE]

Clearly, since $1-\beta^{2(m-i)}\geq 1-\beta^{2}$ for $i=0,\cdots,m-1$, if we define $\alpha_{m}^{2}:=\frac{\beta^{2}(1-\beta^{2m})}{1-\beta^{2}}$, then the condition (33) holds if $L^{2}\eta^{2}\alpha^{2}_{m}-(1-L\eta)\leq 0$. By tightening this condition, we obtain a quadratic equation $L^{2}\eta^{2}\alpha_{m}^{2}-(1-L\eta)=0$ in $\eta$, which leads to

[TABLE]

Note that since $\alpha_{m}^{2}\leq\frac{\beta^{2}}{1-\beta^{2}}$, we have $\eta\geq\underline{\eta}:=\frac{2\sqrt{1-\beta^{2}}}{L(\sqrt{1-\beta^{2}}+\sqrt{1+3\beta^{2}})}$. In that case, by using (31) and (33), (27) reduces to

[TABLE]

Note that $\mathbb{E}\left[\|v_{0}-\nabla{f}(x_{0})\|^{2}\right]\leq\frac{\sigma^{2}}{b}$ and $\mathbb{E}\left[f(x_{m+1})\right]\geq f^{\star}$, we can further bound (35) as

[TABLE]

Multiplying both sides of this inequality by $\frac{2}{\eta(m+1)}$, and then using the lower bound of $\eta$ from (34), we obtain

[TABLE]

Let us choose $\beta:=1-\frac{c_{1}}{\sqrt{b(m+1)}}$ for some $0<c_{1}<\sqrt{b(m+1)}$. In this case, the last two terms of the right-hand side of (36) become

[TABLE]

With this choice of $\beta$, (36) leads to

[TABLE]

(a) Since $\beta=1-\frac{c_{1}}{\sqrt{b(m+1)}}<1$ and $c_{1}<\sqrt{b(m+1)}$, we have

[TABLE]

and $\sqrt{1-\beta^{2}}+\sqrt{1+3\beta^{2}}\leq 1+\sqrt{1+3\beta^{2}}\leq 3$. On the other hand, from (34), we have

[TABLE]

This proves (a).

Let us define define $f^{0}:=f(x_{0})$. Then, using $\beta<1$ and (38) into (37), we get

[TABLE]

(b) Let us choose $b:=c_{2}\sigma^{8/3}(m+1)^{1/3}$ for some constant $c_{2}>0$. Then the last estimate becomes

[TABLE]

To guarantee $\frac{1}{m+1}\displaystyle\sum_{t=0}^{m}\mathbb{E}\left[\|\nabla{f}(x_{t})\|^{2}\right]\leq\varepsilon^{2}$, from (39) we need to set

[TABLE]

This leads to m+1\geq\frac{\sigma}{\varepsilon^{3}}\left[\frac{3Lc_{2}^{1/4}}{\sqrt{c_{1}}}\big{[}f^{0}-f^{\star}\big{]}+\left(c_{1}+\frac{1}{c_{1}}\right)\frac{1}{\sqrt{c_{2}}}\right]^{3/2}. Therefore, we can choose $m$ as shown in (10).

Finally, if $c_{1}=1$, then the number of stochastic gradient evaluations is $\mathcal{T}_{ge}$ is

[TABLE]

which proves (11). $\square$

2.3 The proof of Theorem 3.2: Single-loop with adaptive step-size

First, from Lemma B.1, we have

[TABLE]

where

[TABLE]

and $\omega_{t}$, $\omega_{i,t}$, and $s_{t}$ are defined in Lemma 2.2.

If we fix $\beta_{t}=\beta\in(0,1)$, then we can show that $\omega_{t}=\beta^{2t}$, $\omega_{i,t}=\beta^{2(t-i)}$, and s_{t}=(1-\beta)^{2}\Big{[}\frac{1-\beta^{2t}}{1-\beta^{2}}\Big{]}<\frac{1-\beta}{1+\beta} as in the proof of Theorem 3.1.

Now, let $u_{i}^{2}:=\mathbb{E}\left[\|v_{i}\|^{2}\right]$. To bound the quantity $\mathcal{T}_{m}$ defined by (28), we note that

[TABLE]

Using this expression, we can write $\mathcal{T}_{m}$ from (28) as

[TABLE]

To guarantee $\mathcal{T}_{m}\leq 0$, from the last expression of $\mathcal{T}_{m}$, we can impose the following condition:

[TABLE]

The condition (42) leads to the following update of $\eta_{t}$:

[TABLE]

which is exactly (12).

Next, note that \beta^{2}=\Big{(}1-\frac{c_{1}}{\sqrt{b(m+1)}}\Big{)}^{2}=1-\frac{2c_{1}}{\sqrt{b(m+1)}}+\frac{c_{1}^{2}}{b(m+1)}. Therefore, $1-\beta^{2}=\frac{2c_{1}}{\sqrt{b(m+1)}}-\frac{c_{1}^{2}}{b(m+1)}\geq\frac{c_{1}}{\sqrt{b(m+1)}}$, which implies $\sqrt{1-\beta^{2}}\geq\frac{\sqrt{c_{1}}}{(b(m+1))^{1/4}}$. Using $\sqrt{1-\omega}=\sqrt{1-\beta^{2}}\geq\frac{\sqrt{c_{1}}}{(b(m+1))^{1/4}}$ into (23) of Lemma A.1, we can show that $\Sigma_{m}\geq\frac{\sqrt{c_{1}}(m+1)^{3/4}}{2{L}b^{1/4}}$ as in the first statement (a) of Theorem 3.2.

Note that $\omega_{t}=\beta^{2t}$, by the Chebyshev sum inequality, we have

[TABLE]

Utilizing this estimate, $\mathbb{E}\left[\|v_{0}-\nabla{f}(x_{0})\|^{2}\right]\leq\frac{\sigma^{2}}{b}$, and $S_{t}\leq\frac{(1-\beta)\sigma^{2}}{1+\beta}$ into (40), and noting that $\mathcal{T}_{m}\leq 0$, we have

[TABLE]

Since $\mathbb{E}\left[f(x_{m+1})\right]\geq f^{\star}$, using this into the last estimate, and multiplying the result by $\frac{2}{\Sigma_{m}}$, we obtain

[TABLE]

Since $\left[\frac{1}{b(m+1)(1-\beta)}+(1-\beta)\right]=\left(c_{1}+\frac{1}{c_{1}}\right)\frac{1}{\sqrt{b(m+1)}}$ for $\beta=1-\frac{c_{1}}{\sqrt{b(m+1)}}$, (43) leads to

[TABLE]

The second statement (b) of Theorem 3.2 is proved similarly as in Theorem 3.1 using (44), and we omit the details. $\square$

2.4 The proof of Theorem 3.3: Double-loop with constant step-size

Similar to the proof of (37) in Theorem 3.1, we have

[TABLE]

where we use the superscript $s$ to indicate the stage $s$ in Algorithm 2. Summing up this inequality from $s=1$ to $s=S$, and then multiplying the result by $\frac{1}{(m+1)S}$ and using $\mathbb{E}\left[f(x_{m+1}^{(S)})\right]\geq f^{\star}>-\infty$, we get

[TABLE]

Here, we use the fact that \eta\geq\frac{2}{3{L}\big{[}b(m+1)\big{]}^{1/4}} from (38) in the last inequality.

If we choose $b:=\frac{c_{1}\sigma^{2}}{\varepsilon^{2}}$ and $m+1:=\frac{c_{2}\sigma^{2}}{\varepsilon^{2}}$ for some constants $c_{1}>0$ and $c_{2}>0$ and $c_{1}c_{2}>4$, then, from (46), to guarantee $\frac{1}{S(m+1)}\displaystyle\sum_{s=1}^{S}\sum_{t=0}^{m}\mathbb{E}\left[\|\nabla{f}(x_{t}^{(s)})\|^{2}\right]\leq\varepsilon^{2}$, we require

[TABLE]

Consequently, the total complexity is

[TABLE]

Since we choose $b:=\frac{c_{1}\sigma^{2}}{\varepsilon^{2}}$, the final complexity is $\mathcal{O}\left(\max\left\{\frac{\sigma}{\varepsilon^{3}},\frac{\sigma^{2}}{\varepsilon^{2}}\right\}\right)$, where other constants independent of $\sigma$ and $\varepsilon$ are hidden. $\square$

Appendix C Appendix: The convergence analysis of the mini-batch variants

In this supplementary document, we provide a full analysis of the mini-batch variants of Algorithm 1 and Algorithm 2.

3.1 Variance bound of mini-batch hybrid estimators

For $\hat{v}_{t}$ defined by (16), we have the following property.

Lemma C.1**.**

The mini-batch gradient estimator $\hat{v}_{t}$ defined by (16) satisfies

[TABLE]

where $\rho=\rho(\hat{b}):=\frac{n-\hat{b}}{(n-1)\hat{b}}$ if $n:=|\Omega|$ is finite, and $\rho(\hat{b}):=\frac{1}{\hat{b}}$, otherwise.

Proof.

Let $\hat{v}_{t}$ be defined by (16). Let $z_{t}:=\frac{1}{b_{t}}\sum_{i\in\mathcal{B}_{t}}(\nabla{f}_{\xi_{i}}(x_{t})-\nabla{f}_{\xi_{i}}(x_{t-1}))$, $\bar{z}:=\nabla{f}(x_{t})-\nabla{f}(x_{t-1})$, $\Delta_{t}:=\hat{v}_{t}-\nabla{f}(x_{t})$, and $\Delta{u}_{t}:=u_{t}-\nabla{f}(x_{t})$. Clearly, we have

[TABLE]

Moreover, we can rewrite $\hat{v}_{t}$ in (16) as

[TABLE]

Therefore, using these two expressions, we can derive

[TABLE]

For the finite-sum case, after a few elementary calculations, we can show that

[TABLE]

For the expectation case, we have

[TABLE]

In addition, under Assumption 1.1$($c$)$, we have $\mathbb{E}\left[\|\Delta{u}_{t}\|^{2}\right]\leq\rho\sigma^{2}$.

Substituting one of the two last expressions and the bound of $\mathbb{E}\left[\|\Delta{u}_{t}\|^{2}\right]$ into (48), we get (47). ∎

The following analysis is given under fixed mini-batch sizes when we choose $\hat{b}_{t}=\tilde{b}_{t}=\hat{b}$. Similar to Lemma 2.2, we can bound the variance $\mathbb{E}\left[\|\hat{v}_{t}-\nabla f(x_{t})\|^{2}\right]$ of the mini-batch hybrid estimator $\hat{v}_{t}$ from (16) in the following lemma.

Lemma C.2**.**

Assume that $f(\cdot,\cdot)$ is $L$-smooth and $u_{t}$ is an SGD estimator, $\hat{v}$ is given in (16), $\mathcal{B}_{t}$ and $\hat{\mathcal{B}}_{t}$ are mini-batches of the size $\hat{b}$. Then, we have the following upper bound on the variance $\mathbb{E}\left[\|\hat{v}_{t}-\nabla{f}(x_{t})\|^{2}\right]$:

[TABLE]

where the expectation is taking over all the randomness $\mathcal{F}_{t}:=\sigma(v_{0},v_{1},\cdots,v_{t})$, $\omega_{t}:=\prod_{i=1}^{t}\beta_{i-1}^{2}$, $\omega_{i,t}:=\prod_{j=i+1}^{t}\beta_{j-1}^{2}$ for $i=0,\cdots,t$, and S_{t}:=\sum_{i=0}^{t-1}\big{(}\prod_{j=i+2}^{t}\beta_{j-1}^{2}\big{)}(1-\beta_{i})^{2}\sigma for $t\geq 0$. $\rho=\frac{n-\tilde{b}}{\tilde{b}(n-1)}$ if $|\Omega|$ is finite and $\rho=\frac{1}{\tilde{b}}$ otherwise.

Proof.

From Lemma C.1, taking the expectation with respect to $\mathcal{F}_{t}:=\sigma(v_{0},v_{1},\cdots,v_{t})$, we have

[TABLE]

Let $a_{t}^{2}:=\mathbb{E}\left[\|\hat{v}_{t}-\nabla{f}(x_{t})\|^{2}\right]$ and $r_{t}^{2}=\mathbb{E}\left[\|x_{t+1}-x_{t}\|^{2}\right]$. By following inductive step as in the proof of Lemma 2.2, we obtain

[TABLE]

Using the definition of $\omega_{t}$, $\omega_{i,t}$, and $S_{t}$ in Lemma 2.2, the previous inequality becomes

[TABLE]

which is the same as (49). ∎

3.2 The proof of Corollary 3.1: Single loop with constant step-size and mini-batches

Using Lemma C.2 and following the same path of proof of Lemma B.1, we can show that

[TABLE]

where

[TABLE]

Clearly, we can rewrite $\widehat{\mathcal{T}}_{m}$ as

[TABLE]

where $q_{t}^{2}:=\mathbb{E}\left[\|\hat{v}_{t}\|^{2}\right]$. To guarantee $\widehat{\mathcal{T}}_{m}\leq 0$, we need to have

[TABLE]

Let $\alpha^{2}_{m}:=\frac{\beta^{2}(1-\beta^{2m})}{1-\beta^{2}}$. Since $\alpha^{2}_{1}<\alpha^{2}_{2}<\dots<\alpha^{2}_{m}$, the last condition holds if $L^{2}\eta^{2}\rho\alpha_{m}^{2}-(1-L\eta)\leq 0$. By tightening this condition, we obtain

[TABLE]

which is exactly (17). Since $\alpha_{m}^{2}\leq\frac{\beta^{2}}{1-\beta^{2}}$, we have $\eta\geq\underline{\eta}:=\frac{2\sqrt{1-\beta^{2}}}{L(\sqrt{1-\beta^{2}}+\sqrt{1+\beta^{2}(4\rho-1)})}$.

Next, we can reuse the following estimates as in the proof of Theorem 3.1:

[TABLE]

Combining these estimate into (50) and notting that $\widehat{\mathcal{T}}_{m}\leq 0$ and $\mathbb{E}\left[\|v_{0}-\nabla{f}(x_{0})\|^{2}\right]\leq\frac{\sigma^{2}}{b}$, we can show that

[TABLE]

Note that $\mathbb{E}\left[f(x_{m+1})\right]\geq f^{\star}>-\infty$, (51) can be rewritten as

[TABLE]

If we choose $\beta:=1-\frac{c_{1}}{\sqrt{\hat{\rho}b(m+1)}}$ for any $0<c_{1}<\sqrt{b(m+1)}$ such that

[TABLE]

then (52) leads to

[TABLE]

Since $\beta=1-\frac{c_{1}}{\sqrt{\rho b(m+1)}}<1$ and if we choose $b$, $\hat{b}$, and $m$ such that $\rho b(m+1)>c_{1}^{2}$, we have

[TABLE]

and $\sqrt{1-\beta^{2}}+\sqrt{1+\beta^{2}(4\rho-1)}\leq 1+\sqrt{1+3\beta^{2}}\leq 3$ since $\rho\leq 1$. Therefore, we can bound $\eta$ as

[TABLE]

Therefore, the inequality (53) can be rewritten as

[TABLE]

Let $f^{0}:=\mathbb{E}\left[f(x_{0})\right]$. We can write the bound as

[TABLE]

Let us choose $b:=c_{2}\sigma^{8/3}(\rho(m+1))^{1/3}$ for some $c_{2}>0$. Then, the last inequality leads to

[TABLE]

From (54), to guarantee $\mathbb{E}\left[\|\nabla{f}(\widetilde{x}_{m})\|^{2}\right]\leq\varepsilon^{2}$, we need to choose

[TABLE]

which leads to

[TABLE]

Hence, we can choose $m$ as in (18).

Finally, let $c_{1}=1$. Then the number of stochastic gradient evaluations $\mathcal{T}_{ge}$ is

[TABLE]

which proves (19), where $\rho\leq\frac{1}{\hat{b}}$ if $|\Omega|$ is infinite and $\rho:=\frac{n-\hat{b}}{\hat{b}(n-1)}$ if $|\Omega|$ is finite. In particular, if $|\Omega|$ is infinite and we choose $\rho:=\frac{c_{3}^{2}}{\sigma^{2}\varepsilon^{2}}$ for some $c_{3}\leq\sigma\varepsilon$, then

[TABLE]

Hence, we obtain \mathcal{T}_{ge}=\mathcal{O}\left(\big{(}c_{3}+\frac{1}{c_{3}}\big{)}\frac{\sigma^{2}}{\varepsilon^{2}}\right). $\square$

3.3 The mini-batch variant of Algorithm 2 and its complexity

Let us consider a mini-batch variant of Algorithm 2. Similar to Theorem 3.3, we can prove the following result.

Corollary C.1**.**

Let $\{x^{(s)}_{t}\}_{t=0\to m}^{s=1\to S}$ be the sequence generated by the mini-batch variant of Algorithm 2 using constant step-size $\eta$ defined in (17) with $c_{1}:=1$. Then, the following estimate holds

[TABLE]

Let $\widetilde{x}_{T}\sim\mathbf{U}(\{x^{(s)}_{t}\}_{t=0\to m}^{s=1\to S})$. If we choose $b:=\frac{c_{1}\sigma^{2}}{\varepsilon^{2}}$ and $\frac{m+1}{\hat{b}}:=\frac{c_{2}\sigma^{2}}{\hat{b}^{2}\varepsilon^{2}}$ for some constants $c_{1}>0$ and $c_{2}>0$ and $c_{1}c_{2}>4$, then, to guarantee $\mathbb{E}\left[\|\nabla{f}(\widetilde{x}_{T})\|^{2}\right]\leq\varepsilon^{2}$, we require

[TABLE]

Consequently, the total number of stochastic gradient evaluations $\mathcal{T}_{ge}$ does not exceed

[TABLE]

Proof.

First, similar to the proof of (54), we have

[TABLE]

Summing up this inequality from $s=1$ to $s=S$ and then using $\mathbb{E}\left[f(x_{m+1}^{(S)})\right]\geq f^{\star}$ and $\widetilde{x}_{0}:=x_{0}^{(1)}$, we can show that

[TABLE]

If we choose $b:=\frac{c_{1}\sigma^{2}}{\hat{b}^{2}\varepsilon^{2}}$ and $\frac{m+1}{\hat{b}}:=\frac{c_{2}\sigma^{2}}{\varepsilon^{2}}$ for some constants $c_{1}>0$ and $c_{2}>0$ and $c_{1}c_{2}>4$, then, $\rho(m+1)=\frac{c_{2}\sigma^{2}}{\varepsilon^{2}}$, $b=\frac{c_{1}\rho^{2}\sigma^{2}}{\varepsilon^{2}}$, and from (46), to guarantee $\frac{1}{S(m+1)}\displaystyle\sum_{s=1}^{S}\sum_{t=0}^{m}\mathbb{E}\left[\|\nabla{f}(x_{t}^{(s)})\|^{2}\right]\leq\varepsilon^{2}$, we require

[TABLE]

Consequently, the total complexity is

[TABLE]

Since we choose $b\hat{b}^{2}:=\frac{c_{1}\sigma^{2}}{\varepsilon^{2}}$ which shows that $b\leq\frac{c_{1}\sigma^{2}}{\varepsilon^{2}}$, the final complexity is $\mathcal{O}\left(\max\left\{\frac{\sigma}{\varepsilon^{3}},\frac{\sigma^{2}}{\varepsilon^{2}}\right\}\right)$, where other constants independent of $\sigma$ and $\varepsilon$ are hidden. ∎

Appendix D Appendix: Additional numerical experiments

In this subsection, we provide more numerical examples on two examples we tested in the main text.

4.1 Experiment setup

**Our algorithms: ** We implement the following variants of Algorithm 1 and Algorithm 2 in Python:

• •

Single-loop algorithms: We consider different variants of the single-loop algorithm, Algorithm 1. We denote them by Hybrid-SGD-SL for constant step-size variants, and Hybrid-SGD-ASL for adaptive step-size variants.

• •

Double-loop algorithms: These are variants of Algorithm 2. We denote them by Hybrid-SGD-DL[1-3] the variants corresponding to different snapshot gradient batch-sizes of $b=n^{2/3}$, $b=0.1n$, and $b=n$. We also denote Hybrid-SGD-DL as the best variants among these three choices of the batch-size for snapshot gradient.

Competitors: We also compare our methods with the most state-of-the-art candidates from the literature. We ignore other variants since their complexity bound is worse than ours and they use complicated routines for hyper-parameter selection.

• •

Stochastic gradient descent (SGD): We test two variants of SGD. The first one, called SGD1, is with constant step-size $\eta_{t}:=\frac{0.1}{L}$. The second variant, called SGD2, is with an adaptive step-size of the form $\eta_{t}:=\frac{\eta_{0}}{1+\eta^{\prime}\lfloor t/n\rfloor}$, where $\eta_{0}>0$ and $\eta^{\prime}\geq 0$ are carefully tuned to obtain the best performance. In our tests, we use $\eta_{0}:=\frac{0.1}{L}$ and $\eta^{\prime}:=1$.

• •

SVRG: This algorithmic variant is from [24], where its theoretical step-size in the single sample case is $\eta_{t}:=\frac{1}{3nL}$, and in the mini-batch case is $\eta_{t}:=\frac{1}{3L}$.

• •

SVRG+: This is a variant of SVRG studied in [12]. Its theoretical step-size in the single sample case is $\eta_{t}:=\frac{1}{6nL}$, and in the mini-batch case is $\eta_{t}:=\frac{1}{6L}$.

• •

SPIDER: SPIDER [9] is a stochastic gradient method using SARAH estimator (also called Stochastic Path-Integrated Differential EstimatoR). This method achieves the best-known complexity as Algorithm 2 but uses very different step-size $\eta_{t}:=\min\left\{\frac{\epsilon}{Ln_{0}\left\|v^{k}\right\|},\frac{1}{2Ln_{0}}\right\}$ where $n_{0}=\frac{n^{1/2}}{\hat{b}}$ with $\hat{b}$ is a given mini-batch size in the range $[1,\sqrt{n}]$.

• •

SpiderBoost: SpiderBoost [28] is a modification of SPIDER by using a large constant step-size $\eta_{t}:=\frac{1}{2L}$, but requires to set very specific mini-batch $\hat{b}=\lfloor\sqrt{n}\rfloor$ to achieve the best-known complexity as in Algorithm 2.

Problems: We consider three examples: The first one is the logistic regression with non-convex regularizer as in (20). The second example is a binary classification with non-convex loss as in (21).

Datasets: All the datasets used in this paper are downloaded from LibSVM [7] at

https://www.csie.ntu.edu.tw/ cjlin/libsvm/. We select 6 datasets: w8a ($n=49,749,p=300$), rcv1.binary ($n=20,242,p=47,236$), real-sim $(n=72,309,p=20,958$), news20.binary ($n=19,996,p=1,355,191$), url_combined ($n=2,396,130,p=3,231,961$), and epsilon ($n=400,000,p=2,000$).

4.2 Logistic regression with non-convex regularizer

In this section, we add more numerical examples to solve problem (20). Together with the convergence of the trainning loss and gradient norms in Fig. 1, the training and test accuracies are also plotted in Fig. 4 for three datasets: w8a, rcv1.binary, and real-sim.

As we can observe from Fig. 4, for w8a, all the algorithms except for SVRG achieve similar training accuracy as well as test accuracy. SVRG eventually reaches the same accuracy after around 17 epochs. For rcv1.binary, HybridSGD variants, SGD2, and SVRG+ have similar training and test accuracies, but SGD2 is more oscillated than the other methods. SGD1 performs worse than our methods in this case. Both SPIDER and SVRG still perform poorly. For real-sim, although our methods, SGD1, and SGD2 achieve lower training accuracy, they are able to reach better test accuracy than SVRG+.

In addition, the training and testing accuracies of the mini-batch case are presented in Fig. 5, where the relative residual of the train loss and the gradient norms are shown in Fig. 2. Again, our methods achieve training and test accuracies consistently with SGD2 in w8a and real-sim, while having better accuracy in rcv1.binary.

We also run SVRG, SVRG+, SpiderBoost, and our double-loop variant (Algorithm 2) on three datasets: w8a, rcv1.binary, and real-sim. The results are plotted in Fig. 6.

In this experiment, our double-loop variant and SpiderBoost outperform SVRG and SVRG+. Although the step-size of SVRG+ is $\eta=\frac{1}{6L}$ which is smaller than $\frac{1}{3L}$ in SVRG, SVRG+ still performs better than SVRG. SpiderBoost uses a large step-size $\eta=\frac{1}{2L}$ and it indeed performs slightly better than ours in the w8a dataset, but is comparable in other two. Note that our step-size $\eta$ is selected based on our theory in Theorem 3.3.

Finally, we conduct experiment on three larger datasets: epsilon, url_combined, and news20.binary. Since the the sample sizes are large, we only run mini-batch variants. The results of the single-loop variants are shown in Fig. 7.

We can observe from Fig. 7 that our single loop variants outperform SGD in all three datasets. Note that the performance of the adaptive step-size variant is similar to its fixed step-size one.

The results of the double-loop variants are also shown in Fig. 8.

Clearly, our double-loop variants achieve better performance than SVRG and SVRG+ due to better convergence rate. SpiderBoost is slightly better than ours in the dataset epsilon while they are comparable in the last two datasets since we have the same best-known convergence rate as SpiderBoost.

4.3 Binary classification involving non-convex loss and Tikhonov’s regularizer

We also conduct additional experiments to test our algorithms for solving (21). We use two different non-convex loss functions as in [30] apart from the one used in the main text, which are:

• •

Nonconvex loss in two-layer neural networks: $\ell_{1}(\tau,s)=\left(1-\frac{1}{1+\exp(-\tau s)}\right)^{2}$.

• •

Logistic difference loss: $\ell_{2}(\tau,s)=\log(1+\exp(\tau s))-\log(1+\exp(-\tau s-1))$.

These functions are smooth and satisfy Assumption 1.1.

Let us first test our algorithms and other methods on three datasets: w8a, rcv1.binary, and real-sim using single-sample setting. The results are plotted in Fig. 9 and Fig. 10. In this test, HybridSGD-DL achieves the best performance followed by HybridSGD-SL and HybridSGD-ASL. SPIDER has decent performance in the last two datasets. SGD variants also have good performance in all datasets while SGD2 is better than its fixed step-size variant. SVRG+ also has comparable performance with SGD2 whereas SVRG cannot achieve fast convergence due to its small step-size.

Next, we test mini-batch variants. On the one hand, we compare our single-loop variants HybridSGD-SL and HybridSGD-ASL with SGD. On the other hand, we compare our double-loop variants with SVRG, SVRG+, and SpiderBoost. The results for solving (21) with loss $\ell_{1}$ are shown in Fig. 11 and Fig. 12 whereas Fig. 13 and Fig. 14 present the results when using loss $\ell_{2}$.

Additionally, we repeat the experiments on three larger datasets: epsilon, news20.binary, and ulr_combined. The results are shown in Fig. 15, 16, 17, and 18.

In this experiment, although SGD2 has faster decrease during the first few epochs, our HybridSGD-SL and HybridSGD-ASL eventually achieve lower training loss and gradient norm in all datasets while reaching similar training and testing accuracies as SGD2.

Regarding the double-loop variants, our HybridSGD-DL once again has better performance than SVRG and SVRG+ while having comparable performance with SpiderBoost in terms of training loss, gradient norm, and accuracies.