Parameter-free accelerated gradient descent for nonconvex minimization

Introduction

This paper studies general nonconvex optimization problems:

To alleviate the problem, we propose a new first-order method. The proposed method is an AGD equipped with two restart mechanisms and enjoys the following advantages.

Our method finds an $\varepsilon$-stationary point in $O(\varepsilon^{-7/4})$ function and gradient evaluations under Section 1.

Our method estimates the Lipschitz constants $L_{f}$ and $M_{f}$ instead of requiring prior knowledge. It automatically updates the estimates, and their initial guesses are allowed to deviate from the actual Lipschitz constants, both in theory and practice.

Our method also does not need the target accuracy $\varepsilon$ as input, which implies its global convergence (Corollary 5.16).

We also conducted numerical experiments with several instances. The results show that the proposed method performs comparably to existing algorithms with similar complexity bounds, even without parameter tuning.

As shown in Figure 1, the proposed method successfully converged to a stationary point for a nonconvex test problem, regardless of the initial guesses of the Lipschitz constants, and it outperformed other state-of-the-art methods.

The main challenge in achieving the above advantages is estimating $M_{f}$ by using only first-order information, i.e., $f$ and $\nabla f$. To this end, we establish a Hessian-free analysis with two inequalities less familiar in the optimization context: a Jensen-type inequality for gradients and an error bound for the trapezoidal rule. These two inequalities do not include the Hessian itself but rather its Lipschitz constant $M_{f}$, and they enable us to estimate $M_{f}$ and simplify the analysis.

Related work

This section reviews algorithms and techniques related to this work and characterizes our method in comparison with the existing methods.

After updating the solutions, Algorithm 1 checks whether the current estimate $L$ approximates the actual value $L_{f}$ well. If

then $L$ is judged suitable, and the epoch continues; otherwise, AGD restarts from $x_{k-1}$ with a larger $L$. Since condition 9 holds in the previous iteration, starting the next epoch from $x_{k-1}$ ensures that this epoch does not increase the objective function value. Condition 9 guarantees a sufficient decrease in the objective function in the epoch and plays a crucial role in our analysis. As will be shown later, this condition is always satisfied when $L$ is sufficiently large, which justifies the restart mechanism.

This restart mechanism can be regarded as a kind of backtracking. A well-known backtracking method for estimating $L_{f}$ is based on the Armijo rule , and one might think of using the Armijo rule or its variant [8, Section 5.3] instead of condition 9. However, with the Armijo rule, the estimate $L$ changes at each iteration, complicating the complexity analysis of AGD, especially for the nonconvex case. In contrast, we fix $L$ through the epoch with the help of condition 9, preserving the simplicity of the existing analysis when $L_{f}$ is known.

Our AGD also resets the effect of acceleration by restart when $S_{k}$ becomes large; the restart condition involves an estimate $M_{k}$ of the Hessian’s Lipschitz constant $M_{f}$ (detailed later). If

the epoch continues, otherwise, AGD restarts with the computed solution $x_{k}$ and a (possibly) smaller $L$. Condition 10 is inspired by a condition in , $kM_{f}S_{k}\leq\varepsilon$; unlike the existing one, our condition does not involve $M_{f}$ or $\varepsilon$. Although $L$ does not have to be decreased on Line 11 to derive the complexity bound of $O(\varepsilon^{-7/4})$, decreasing $L$ will help the algorithm to capture the local curvature of $f$ and to converge faster empirically.

Our AGD guarantees that the gradient norm is small at a weighted average $\bar{y}_{k}$ of the solutions $y_{0},\dots,y_{k-1}$. The averaged solution $\bar{y}_{k}$ is defined by

Note that it is unnecessary to keep all $y_{0},\dots,y_{k-1}$ in memory to compute $\bar{y}_{k}$ because $Z_{k}$ and $\bar{y}_{k}$ satisfy a simple recursion:

which maintains the computational efficiency of the algorithm. Given 8, there is another form:

and this simplification of $\bar{y}_{k}$ and $Z_{k}$ is one of the aims of setting $\theta_{k}$ as in 8.

Complexity analysis

The estimate $M_{k}$ for the Hessian is updated at each iteration, unlike $L$ for the gradient. Our complexity analysis requires three inequalities on $M_{k}$: $M_{k-1}\leq M_{k}$ and

Finding the smallest $M_{k}$ that satisfies these inequalities is desirable, and fortunately, it is straightforward,

because $x_{k}$, $y_{k}$, and $\theta_{k}$ are already in hand. In other words, our $\theta_{k}$ in 8 does not depend on the estimate $M_{k}$ of $M_{f}$, which facilitates computation of $M_{k}$. Thus, we can estimate $M_{f}$ without backtracking, thereby simplifying the algorithm and improving its efficiency.

in the algorithm. For $k=1$, Eq. 18 takes the form $M_{1}=\max\{M_{0},\ast,\ast,\frac{0}{0}\}$, which should be treated as $M_{1}=\max\{M_{0},\ast,\ast\}$. This replacement ensures that

and will lead to a smaller complexity bound with lessened dependence on $M_{0}$ and $M_{f}$, as will be confirmed later in Theorem 5.13. Eq. 18 derives a theoretically better bound than 16 but requires an additional gradient evaluation at $\bar{y}_{k}$, increasing the computational cost per iteration. Therefore, we recommend using 16 in practice.

This section provides the complexity bounds of $O(\varepsilon^{-7/4})$ for Algorithm 1.

The estimates $L$ and $M_{k}$ of the Lipschitz constants should be large enough to satisfy some technical inequalities such as 9, 14, and 15, but they should not be too large. This section proves the following upper bound on $L$.

This proposition immediately follows from the following lemma.

Suppose that \theassumption(a) holds. During epochs with $L\geq L_{f}$, condition 9 always holds.

Before providing the proof, let us describe the idea behind it.

To obtain the descent condition 9, we first evaluate the decrease for one iteration. However, a common difficulty arises when analyzing AGD: the objective function value does not necessarily decrease monotonically. To deal with the problem, we introduce a potential function related to $f(x_{k})$ and show that it is nearly decreasing. The potential function $\Phi_{k}$ is defined by

inspired by 45 used in the proof. The value $\Phi_{k}$ decreases when $\left\|x_{k}-x_{k-1}\right\|$ is small enough as the following lemma shows.

Under \theassumption(a) and $L\geq L_{f}$, the following holds for all $k\geq 0$:

This lemma can be proven by putting Lemma 2.1 due to the Lipschitz gradient together with inequalities 14 and 15 due to the Lipschitz Hessian. Summing Lemma 5.3 over $k$ and doing some calculations yields Lemma 5.2. To evaluate the terms of $\left\|x_{k}-x_{k-1}\right\|^{3}$ and $\left\|x_{k}-x_{k-1}\right\|^{4}$ in Lemma 5.3 in these calculations, we use

which follows from the fact that the epoch did not end at iteration $k-1$.

Now, we prove Lemma 5.3 and use it to prove Lemma 5.2. The following proofs do not use the specific form 8 of $\theta_{k}$ but rather a more general condition on $\theta_{k}$,

which is easily verified for all $k\geq 0$ under $\theta_{0}\coloneqq\theta_{1}^{2}$ and 8 for $k\geq 1$.

to simplify the notation. From Lemmas 2.1 and 7, we have

where the last inequality uses 7 and $0\leq\theta_{k}\leq 1$. To evaluate the first term on the right-hand side, we decompose it into four terms:

Plugging the evaluations into 29 results in

Next, to bound the last term on the right-hand side, we use the following inequality from 15:

where we have used $(a+b)^{2}\leq(1+1/\theta)a^{2}+(1+\theta)b^{2}$ for $a,b,\theta>0$. Rearranging the terms yields

This inequality can be rewritten with the potential function $\Phi_{k}$ defined in 22 as

Finally, using $\theta_{k+1}^{2}-\theta_{k}\leq 0$ from 26 and

from Young’s inequality ($|\left\langle{a},{b}\right\rangle|\leq\frac{1}{2}\left\|a\right\|^{2}+\frac{1}{2}\left\|b\right\|^{2}$) yields the desired result as

Summing Lemma 5.3 over $k$ and telescoping yields

where we have used 26, 56, and $M_{i}\leq M_{k-1}$. Combining 57, 58, and 59 yields

Furthermore, we bound the last two terms as

This section proves the following upper bound on $M_{k}$.

Suppose that \theassumption(b) holds. Then, the following is true throughout Algorithm 1: $M_{k}\leq\max\{M_{0},M_{f}\}$.

Because of the update rules 16 or 18 of $M_{k}$, it suffices to show the inequalities obtained by replacing $M_{k}$ with $M_{f}$ in 14, 15, and 21. One of them has already been obtained as Lemma 3.2. Below we prove the other two; the following lemma gives the formal statement.

Under \theassumption(b), the following hold for all $k\geq 1$:

The following proofs do not use the specific form 8 of $\theta_{k}$; they only use $0\leq\theta_{k}\leq 1$.

Note that $x_{k}=\frac{1}{1+\theta_{k}}y_{k}+\frac{\theta_{k}}{1+\theta_{k}}x_{k-1}$ from 7, and Lemma 3.1 gives

Multiplying by $(1+\theta_{k})$ concludes the proof.

where we have used $x_{0}=x_{-1}$. Now, we obtain

where the last inequality follows from $p_{k,i}\leq p_{k,k-1}=1/Z_{k}$ for all $0\leq i<k$.

Next, we bound the second term on the right-hand side. Since for $0\leq i<j<k$,

Plugging this bound into 73 concludes the proof.

3 Upper bound on gradient norm

Suppose that Section 1 holds. In Algorithm 1, the following is true when $k\geq 2$:

If $M_{k}$ is computed with 18 in the algorithm, the above bound is improved to

Each term on the right-hand side can be bounded as follows:

which is the first inequality of the desired result 82. The second inequality follows from 25:

When $M_{k}$ is computed with 18, we have 21. Similarly to the proof of 82, using 21 instead of 68 yields

which implies the first inequality in 83. The second inequality is obtained similarly to the second inequality of 82.

4 Main results

Combining Lemma 5.10 with conditions 10 and 9, we obtain the following complexity bound for Algorithm 1.

In Algorithm 1, when $\bar{y}_{k}$ defined by 11 satisfies $\left\|\nabla f(\bar{y}_{k})\right\|\leq\varepsilon$ for the first time, the total iteration count $K$ is at most

If $M_{k}$ is computed with 18 in the algorithm and $\beta$ is set as $\beta=1$, the above bound is improved to

We count the number of iterations of Algorithm 1 separately for three types of epoch:

successful epoch: an epoch that does not find an $\varepsilon$-stationary point and ends at Line 11 with the descent condition 9 satisfied,

unsuccessful epoch: an epoch that does not find an $\varepsilon$-stationary point and ends at Line 9 with the descent condition 9 unsatisfied,

last epoch: the epoch that finds an $\varepsilon$-stationary point.

Let us focus on a successful epoch. Let $k$ be the iteration number of the epoch. It follows that

Plugging this bound into condition 9, we obtain

Summing this bound over all successful epochs leads to the conclusion that the total iteration number of successful epochs is at most

For later use, we will also evaluate the number of successful epochs. Combining 99 and 97 also yields

Plugging this bound into condition 9 gives

We deduce from 82 that the number of iterations for each epoch is at most

where $c_{1}$ and $c_{2}$ are defined by 94. The total iteration number of unsuccessful and last epochs is at most

Putting this bound together with 106 concludes the proof.

The proof is similar to 95 through the evaluation of the total iteration number of successful, unsuccessful, and last epochs.

Let us focus on a successful epoch. Let $k$ be the iteration number of the epoch. As a counterpart of 99, we can obtain

Plugging this bound into condition 9, we obtain

Hence, as in the proof of 95, the total iteration number of all successful epochs is at most

We deduce from 83 that the iteration number of each epoch is at most \big{(}\frac{4\bar{L}^{2}}{M_{0}\varepsilon}\big{)}^{1/4}. In addition, unsuccessful epochs occur only at most $c_{2}$ times when $\beta=1$, where $c_{2}$ is defined by 94. Thus, the total iteration number of unsuccessful and last epochs is at most

Putting this bound together with 121 concludes the proof.

The number of function and gradient evaluations is of the same order as the iteration complexity given in Theorem 5.13, because Algorithm 1 evaluates the objective function and the gradient at two or three points in each iteration.

Since Algorithm 1 does not need the target accuracy $\varepsilon$ as input, we can deduce the following global convergence property as a byproduct of the complexity analysis.

Suppose that Section 1 holds. Let $z_{K}$ denote $\bar{y}_{k}$ at the total iteration $K$ of Algorithm 1. Then, the following holds:

5 Discussion

To gain more insight into the complexity analysis provided above, let us discuss it from three perspectives.

In order to deduce from 21 that the gradient norm $\left\|\nabla f(\bar{y}_{k})\right\|$ is small, the $Z_{k}$ should be large, and hence the acceleration parameter $\theta_{k}$ should also be large, from definition 11 of $Z_{k}$. At the same time, $\theta_{k}$ should be small to obtain a significant objective decrease from 9. Our choice of $\theta_{k}=\frac{k}{k+1}$ in 8 strikes this balance, and the essence is the following two conditions: $\theta_{k}=1-\Theta(k^{-1})$ and $\theta_{k+1}^{2}\leq\theta_{k}\leq\theta_{k+1}$. As long as $\theta_{k}$ satisfies them, we can ensure a complexity bound of $O(\varepsilon^{-7/4})$ even with a different choice than in 8.

The proposed AGD restarts when $S_{k}=\sum_{i=1}^{k}\left\|x_{i}-x_{i-1}\right\|^{2}$ is large, which plays two roles in our analysis. First, as long as the epoch continues, the solutions $y_{0},\dots,y_{k-1}$ are guaranteed to be close together, and therefore Lemma 3.1 gives a good approximation $\sum_{i=0}^{k-1}p_{k,i}\nabla f\left(y_{i}\right)$ of the gradient $\nabla f(\bar{y}_{k})$. Second, the restart mechanism ensures that $\left\|x_{k}-x_{k-1}\right\|^{3}$ and $\left\|x_{k}-x_{k-1}\right\|^{4}$ are small compared to $\left\|x_{k}-x_{k-1}\right\|^{2}$ during the epoch, which helps to derive Lemma 5.2 from Lemma 5.3.

Numerical experiments

Jin et al. also use a potential function to analyze AGD under Section 1, but our proof technique differs significantly from theirs. First, our potential 22 is more complicated than the existing one. This complication is mainly caused by the dependence of our $\theta_{k}$ on $k$; $\theta_{k}$ in does not depend on $k$, while our $\theta_{k}$ does but with the benefit of not requiring knowledge of $M_{f}$ or $\varepsilon$. Second, Jin et al. prove the potential decrease only in regions where $f$ is weakly convex and use a different technique, negative curvature descent, in other regions. In contrast, our Lemma 5.3 is valid in all regions and provides a unified analysis with the potential function. This unified analysis is made possible by taking full advantage of the Lipschitz continuous Hessian through Lemmas 3.1 and 3.2.

This section compares the performance of the proposed method and six existing methods. We implemented objective functions in Python with JAX and Flax . We implemented some methods in Python and used SciPy implementation for other methods. They were executed on a computer with an Apple M1 Chip (8 cores, 3.2 GHz) and 16 GB RAM. The source code used in the experiments is available on GitHub. https://github.com/n-marumo/restarted-agd

We compared the following algorithms: Proposed, GD, JNJ2018, LL2022, OC2015, L-BFGS, and CG.

JNJ2018 [30, Algorithm 2] is an AGD method reviewed in Section 2. The parameters were set in accordance with [30, Eq. (3)]. The equation involves constants $c$ and $\chi$, whose values are difficult to determine; we set them as $c=\chi=1$.

LL2022 [32, Algorithm 2] is a state-of-the-art AGD method that makes use of the practical superiority of AGD. The parameters were set in accordance with [32, Theorem 2.2 and Section 4].

OC2015 is an AGD method with an adaptive restart scheme proposed in for convex optimization. We used the gradient scheme among the two restart schemes proposed in [37, Section 3.2]. The parameter of [37, Algorithm 1] was set as $q=0$, and the step size was determined by backtracking as in GD.

L-BFGS is the limited-memory BFGS method . We used the SciPy implementation , i.e., scipy.optimize.minimize with option method="L-BFGS-B".

CG is a nonlinear conjugate gradient algorithm, a variant of the Fletcher–Reeves method described in [36, pp.120–122]. We used the SciPy implementation , i.e., scipy.optimize.minimize with option method="CG".

The parameter setting for JNJ2018 and LL2022 requires the values of the Lipschitz constants $L_{f}$ and $M_{f}$ and the target accuracy $\varepsilon$. For these two methods, we tuned the best $L_{f}$ among $\{10^{-4},10^{-3},\dots,10^{4}\}$ and set $M_{f}=1$ and $\varepsilon=10^{-16}$ following . Note that if these values deviate from the actual values, the methods do not guarantee convergence.

2 Problem setting

We test the performance of the algorithms with three types of problem instances.

The first two instances are of training neural networks for classification and autoencoder:

The third instance is low-rank matrix completion:

3 Results

Figure 2 shows the objective function value $f(x_{k})$ and the estimates $L$ and $M_{k}$ at each iteration of Proposed. The iterations at which a restart occurred are also marked; “successful” and “unsuccessful” mean restarts at Line 11 and Line 9 of Algorithm 1, respectively.

This figure illustrates that the proposed algorithm adaptively estimates the Lipschitz constants $L_{f}$ and $M_{f}$. In particular, in Figure 1a, the values of $L$ and $M_{k}$ differ significantly between the early and final stages; the algorithm adopts a larger step size (i.e., smaller $L$) near the stationary point. This adaptability is expected to improve performance over algorithms that treat the Lipschitz constants as given values.

We can also see that the proposed method restarts frequently in the early stages but that the frequency decreases as the iterations progress. This behavior is explained by the restart condition on Line 10 of Algorithm 1; as the algorithm progresses, the iterate $x_{k}$ moves less in general, i.e., $S_{k}\coloneqq\sum_{i=1}^{k}\left\|x_{i}-x_{i-1}\right\|^{2}$ grows slower, and the restart condition holds less frequently.

The value of $M_{k}$ becomes very small at each restart because we set $M_{0}=10^{-16}$. The value is updated immediately, so changing $M_{0}$ to $10^{-12}$ or , for example, does not affect the algorithm’s behavior for the problem instances.

Figures 2a and 2d show that Proposed converges faster than the sublinear rate guaranteed by the theoretical analysis, $O(\varepsilon^{-7/4})$. This speed-up is because Proposed automatically estimates the Lipschitz constants and successfully captures the local curvature of $f$. In the setting shown in Figure 2c, LL2022 reduces the objective function value faster than Proposed, probably because the two methods converge to different stationary points; LL2022 is fortunate in this instance to have found a better stationary point in terms of the function value.

We emphasize again that we fixed all of the input parameters for Proposed, whereas we tuned $L_{f}$ for the input parameter for LL2022 and JNJ2018. Proposed achieved comparable or better performance than those of the (theoretically) state-of-the-art methods without parameter tuning.

The previous section compares the proposed algorithm with two of the theoretically superior algorithms reviewed in Section 2. This section provides further numerical comparisons with OC2015, L-BFGS, and CG, which are known to be practical but have no complexity guarantees or weaker ones than the algorithms in the previous section. We also provide the results of GD as a baseline.

The results are given in Figure 4, illustrating that OC2015, L-BFGS, and CG perform surprisingly well in practice, though not in theory. To obtain results of L-BFGS and CG, we ran the SciPy functions multiple times with the maximum number of iterations set to $2^{0},2^{1},2^{2},2^{3},\dots$ because we cannot obtain the solution at each iteration while running SciPy codes of L-BFGS and CG, but only the final result. The results are thus plotted as markers instead of lines in Figure 4. Although these methods are often ignored in the current trend of research on algorithms that achieve better complexity bounds, they are so practical that they are worth considering from theoretical perspectives. In future work, we would like to investigate whether these practical algorithms achieve the same complexity bound as the proposed algorithm or what types of instances can worsen their complexity bounds.

To simplify the notation, let $\bar{z}\coloneqq\sum_{i=1}^{n}\lambda_{i}z_{i}$. Taylor’s theorem gives

Furthermore, elementary algebra shows that $\sum_{i=1}^{n}\lambda_{i}\left\|z_{i}-\bar{z}\right\|^{2}=\sum_{1\leq i<j\leq n}\lambda_{i}\lambda_{j}\|z_{i}-z_{j}\|^{2}$ as

Appendix B Details of neural networks used in experiments

For problem 124, we used a three-layer fully connected network with bias parameters. The layers each had $784$, $32$, $16$, and $10$ nodes and had the logistic sigmoid activation. The total number of the parameters is $d=(784\times 32+32\times 16+16\times 10)+(32+16+10)=25818$.

For problem 125, we used a four-layer fully connected network with bias parameters. The layers each had $784$, $32$, $16$, $32$, and $784$ nodes and had the logistic sigmoid activation. The total number of the parameters is $d=(784\times 32+32\times 16+16\times 32+32\times 784)+(32+16+32+784)=52064$.

We initialize the parameters of the neural networks with flax.linen.Module.init method. See the documentationhttps://flax.readthedocs.io/en/latest/api_reference/flax.linen/module.html for more details of the initialization.

Acknowledgments

We are deeply grateful to the anonymous reviewers, who carefully read the manuscript and provided helpful comments. We also thank the associate editor for sharing information on relevant literature.

References