NAG-GS: semi-implicit, accelerated and robust stochastic optimizer.
Authors: Valentin Leplat
Daniil Merkulov
Aleksandr Katrutsa
Daniel Bershatsky
Ivan Oseledets
Venue: Submitted to ICLR 2023 (2023)
Keywords: Accelerated gradient methods
stochastic optimization
stochastic differential equations
semi-implicit solver
convergence analysis
deep neural networks

Abstract: Classical machine learning models such as deep neural networks are usually trained by using Stochastic Gradient Descent-based (SGD) algorithms. The classical SGD can be interpreted as a discretization of the stochastic gradient flow. In this paper we propose a novel, robust and accelerated stochastic optimizer that relies on two key elements: (1) an accelerated Nesterov-like Stochastic Differential Equation (SDE) and (2) its semi-implicit Gauss-Seidel type discretization. The convergence and stability of the obtained method, referred to as NAG-GS, are first studied extensively in the case of the minimization of a quadratic function. This analysis allows us to come up with an optimal step size (or learning rate) in terms of rate of convergence while ensuring the stability of NAG-GS. This is achieved by the careful analysis of the spectral radius of the iteration matrix and the covariance matrix at stationarity with respect to all hyperparameters of our method. We show that NAG-GS is competitive with state-of-the-art methods such as momentum SGD with weight decay and AdamW for the training of machine learning models such as the logistic regression model, the residual networks models on standard computer vision datasets, and Transformers in the frame of the GLUE benchmark.

Under review as a conference paper at ICLR 2023
NAG-GS: SEMI-IMPLICIT, ACCELERATED ROBUST STOCHASTIC OPTIMIZERS
AND
Anonymous authors Paper under double-blind review
ABSTRACT
Classical machine learning models such as deep neural networks are usually trained by using Stochastic Gradient Descent-based (SGD) algorithms. The classical SGD can be interpreted as discretization the stochastic gradient flow. In this paper we propose a novel, robust and accelerated stochastic optimizer that relies on two key elements: (1) an accelerated Nesterov-like Stochastic Differential Equation (SDE) and (2) its semi-implicit Gauss-Seidel type discretization. The convergence and stability of the obtained method, referred to as NAG-GS, are first studied extensively in the case of the minimization of quadratic function. This analysis allows us to come up with an optimal step size (or learning rate) in terms of rate of convergence while ensuring the stability of NAG-GS. This is achieved by the careful analysis of the spectral radius of the iteration matrix and the covariance matrix stationarity with respect to all hyperparameters of our method. We show that NAG-GS is competitive with state-of-the-ar methods such as momentum SGD with weight decay and AdamW for the training of machine learning models such the logistic regression model, the residual networks models on standard computer vision datasets, and Transformers in the frame of the GLUE benchmark.
INTRODUCTION
Nowadays, machine learning, and more particularly deep learning, has achieved promising results on wide spectrum of AI application domains. order to process large amounts of data, most competitive approaches rely the of deep neural networks. Such models require be trained and the process of training usually corresponds to solving a complex optimization problem. The development of fast methods urgently needed speed up the learning process and obtain efficiently trained models. In this paper, we introduce a new optimization framework for solving such problems.
Main contributions of our paper:
We propose a new accelerated gradient method of Nesterov type for convex and non-convex stochastic optimization;
We analyze the properties of the method both theoretically and experimentally;
We show that our method is robust to selection of hyperparameters. memory-efficient compared with AdamW and competitive with baseline methods in various benchmarks.
Organization of our paper:
Section gives the theoretical background for our method.
In Section 2 we propose an accelerated system of Stochastic Differential Equations (SDE) and an appropriate solver that rely on a particular discretization of the SDE's system. The obtained method, referred as NAG-GS (Nesterov Accelerated Gradient with Gauss-Seidel Splitting). is first discussed in terms of convergence in the simple but central case of quadratic functions. Moreover, we apply our method for solving a I-dimensional non-convex SDE for which we bring strong numerical evidences of the superior acceleration allowed by NAG-GS method compared to classical SDE solvers, see Appendix B
In Section 3 NAG-GS is tested to tackle stochastic optimization problems of increasing complexity and dimension, starting from the logistic regression model to the training of large machine learning models such ResNet20, ResNet50 and Transformers.

---

Under review as a conference paper at ICLR 2023
1.1 PRELIMINARIES
We start here with some general considerations in the deterministic setting for obtaining an accelerated Ordinary Differential Equations (ODE) that will be extended in the stochastic setting Section 2.1 We consider iterative methods for solving the unconstrained minimization problem:
min f(x), NVV
(1)
where V is a Hilbert space, and f - [+ is a properly closed convex extended valled function. In the following, for simplicity, consider function f smooth on the entire space. We also suppose V is equipped with the canonical inner product (ar, = Σ Tayi and the correspondingly induced norm ||<||| / (a, x)) Finally, will consider in this section the class of functions S which stands for the set of strongly convex functions of parameter μ 0 with hipschitz-continoon gradients of constant L > 0. For such class of functions, is well-known that the global minimizer exists uniquely Nesterov| (2018). One well-known approach to derive the Gradient Descent (GD) method is discretizing the so-called gradient flow:
shall consider the particular case of RR" for V and
(t) fxx(tt)
t >>0 0.
The simplest forward (explicit) Euler method with step size og leads to the GD method:
1k+1 2g ag f=xx).
In the terminology of numerical analysis, it well-known that this method is conditionally A-stable and for f E s/.1 with 0 S μ S S co, the step size a = 1/L allowed to get linear rate of convergence. Note that the highest rate of convergence is achieved for ag F+L
2 In this case:
(0/
Q
2h
|| 2²112 S
1100 - 11
with 2, L. usually referred as the condition number of function f Nesterov (2018). One can also consider the backward (implicit) Euler method:
Fk+1 f- Fy a V V (x
which is unconditionally A-stable. Here-under, we summarize the methodology proposed by Luo & Chen (2021) to come up with a general family of accelerated gradient flows by focusing on the following simple problem:
1 x= Az 2
min f(x) =
(4)
*ER"
for which the gradient flow in equation
2 reads simply as:
i(t) Ar((t) t 0,
(5)
where A is a n-by-n symmetric positive semi definief matrix ensuring that f L where μ and respectively correspond to the minimum and maximum eigenvalues matrix A, which are real and positive by hypothesis. Instead of solving directly equation 5| authors of Luo & Chen (2021) turn to a general linear ODE system:
i(t) = Gy((). t> 0.
The main idea consists in seeking such a system 6 with some asymmetric block matrix G that transforms the spectrum of A from the real line to the complex plane and reduces the condition number from x( (A) L to xGG) = (v
Afterwards, accelerated gradient methods can be constructed from A-stable methods for solving equation 6 with significant larger step size and consequently
improve the contraction rate from O ((0/1)
2k
2k
t
Furthermore, to handle
[Q
the convex case μ = authors from Luo Chen (2021 combine the transformation idea with suitable time scaling technique. In this paper we consider one transformation that relies on the embedding of into some x 2 block matrix G with a rotation built-in Luo & Chen (20II):
-1 GNAG = 1/7 A/T -P/7T
I
2

---

Under review as a conference paper at ICLR 2023
where y is positive time scaling factor which satisfies
" = 7 7(t = 0) = 70 0.
(8)
Note that, given A positive definite, have that (() < for all E (G) with a (G) denoting the spectrum of G, i.e. the set of all
easily show that for the considered transformation,
eigenvalues of G. Further, will denote by p(G) max [X] the spectral radius of matrix G. Let
now consider the NAG block Matrix and let y (z, v), the dynamical system given in equation 6 with ((0) = y/o E R2n reads:
I,
M (x v) Ar 7
with initial conditions x0)) zo and v()) vo. Before going further, let us remark that this linear ODE can be expressed as the following second-order ODE by eliminating
7it + (7 + 1a) + Ar
(10)
where Ax is therefore the gradient of f w.r.t. Thus, one could generalize this approach for any
function f E sy: by replacing A and Ax by (f(x) respectively within equation 7 equation 9 and
equation 10 Finally, some additional and useful insights are discussed in Appendix A
2 MODEL AND THEORY
2.1
ACCELERATED STOCHASTICC GRADIENT FLOW
In the previous section, we presented a family of accelerated Gradient flows obtained by an appropriate spectral transformation G' of matrix A, see equation 9 One can observe the presence of a gradient term of the smooth function f(x)) at the second differential equation. Let us recall that Ar can be replaced by V fz) for any function f E St. In the frame of this paper, function f(x) may correspond to some loss function used to train neural networks. For such setting, we assume that gradient input (x()) is contaminated by noise due to finite-sample estimate of the gradient. The study of accelerated Gradient flows now adapted to include and model the effect of the noise: to achieve this we consider the dynamics given in equation perturbed by general martingale process. This leads us to consider the following Accelcrated Stochastic Gradient (ASG) flows:
dr D.
dl
(11)
du
(x v) Ar dz dt
1
dt
?
which corresponds to an (Accelerated) system of SDE's, where Z(1) is continuous Ito martingale. We assume that (t) has the simple expression dZ odW, where W =W1,.., W.,) is standard n--imensiioal Brownian Motion. As a simple and first approach, we consider the volatility parameter constant. In the next section, present the discretizations considered for ASG flows given in equation 11
2.2 DISCRETIZATION GAUSS-SEIDER SPLITTING AND SEMI-IMPLICITNESS
In this section, we present the main strategy discretize the Accelerated SDE's system from equation II The main motivation behind the discretization method is to derive integration schemes that are, in the best case, unconditionally A-stable or conditionally A-stable with the highest possible integration step. In the classical terminology of (discrete) optimization methods, this amounts to ensure convergence of the obtained methods with the largest possible step size and consequently improve the contraction rate (or the rate of convergence). In Section [.1 we have briefly recalled that the most well-known unconditionally A-stable scheme was the backward Euler method (see equation 3), which is implicit method and hence achieve faster convergence rate. However, this requires to either solve a linear system either, in the case a general convex function, to compute the root of non-linear equation, both situations leading high computational cost. This is the
3

---

Under review as a conference paper at ICLR 2023
main reason why few implicit schemes are used in practice for solving high- dish-simensional optimization problems. But still, is expected that an explicit scheme closer to the implicit Euler method will have good stability with larger step size than the one offered by forward Euler method. Motivated by the Gauss-Seidel (GS) method for solving linear systems, we consider the matrix splitting = M + N with M being the lower triangular part of G and N G- M, we propose the following Gauss-Seidel perturbated with noise:
splitting scheme for equation 6
yk++ yk
Myk+1 + Nyk +
W W We
(12)
which for G G AA AG (see gives the following semi-implicit scheme with step size ar 0:
1k++1
1
2k
Uk Ik++1
O
(13)
Uk+1 1
Uk 1 (xx+1 Uk+1) Arx+1 + WK+1 WR Ykk
Note that due to properties of Brownian motion we can simulate its values at the selected points by: WK+1 W + AWR, k. where AW; are independent random variables with distribution (00 (r.). On
practical point of view, we will use AW; VogT where 1/s (O0 1). Furthermore, ODE corresponding to the parameter 7 is also discretized implicitly:
+ 7k
7
70
0.
(14)
0
As already mentioned earlier, heuristically, for general f $ with p 2 0, we just replace Ar in 13 with V fx) and obtain the following NAG-GS scheme:
equation
#k++1 Tk
Uk
Ik+l.
Ork:
(15)
Uk+1
P
(II+1 Uk+1)
(f((++1) + 7/k
W Ws
0kk
O%k
Finally, we come the following methods, which are addressed NAG-GS method (see Algo- rithm D. There we assume that the gradient (((++1)) is computed with some unknown noise.
Algorithm Nesterov Accelerated Gradients with
splitting (NAG-GS).
Input: Choose point Zig E RR"", some y 2 0, 7o 0.
Set Uo == Ti-
for = 1, do
Choose step size or 0.
> Update parameters and state x:
Set dソ (a 1)-1.
Set 7 = (1 T + axpL.
Set 3k+1 = an Jark + agui-
> Update state e:
Set bR == a. + 7k+1)
Set Uk+1 == (1 by) Uk + bb bb V f((++1).
end for
Moreover, the step size update can be performed with diffferntt strategies, for instance one may choose the method proposed by Nesterov (Nettrrov 2018 Method 2.2.7) which specifies to compute as E (0, 1) such that La (1 0gk + 04/44 Note that for %o μ.. hence the sequences
7k = μ and or = / for all 0. Other strategies are discussed in Luo Chen (2021 such as
La 7x [1 + a.). Finally, one could simply choose a constant step size variant Algorithm we
select and discuss this approach in Section 2.3
Let us mention that full- aiplliiti discretizations have been considered and studied by the authors, these will be briefly discussed in Appendix A.2 However, their interests are, the moment, limited for ML applications since the obtained implicit schemes are connected to a specific family of second-order methods which are intractable for real-life MI. models.
4

---

Under review as a conference paper at ICLR 2023
Table Summary on the comparison of NAG-GS the reference optimizer for different neural archi- tectures (greater is better). Target metrics are ACC @ 1 and ACC@5 5 for RESNET20 and RESNET50 respectively and average score on GLUE for ROBERTA.
MODEL
DATASET
OPTIMIZER
SCORE
RESNET20
CIFAR-10
SGD-MW
91.25
NAG-GS
91.29
RESNET50
IMAGENET
SGD-MW
92.70
NAG-GS
89.82
ROBERTA
GLUE
ADAMW NAG-GS
82.92
81.30
2.3 CONVERGENCE ANALYSIS OF QUADRATIC CASE
We mentioned in previous section that one could choose a constant step size strategy for Algorithm We propose to study how to select maximum (constant) step size that ensures optimal contraction rate while guaranteeing the convergence, or the stability of NAG-GS method once used to solve show that the choice of the optimal (constant) step size actually
SDE's system I1 Ultimately, mostly influenced by the values of L and 7,. These (hyper)parameters are central and in order to show this, we study two key quantities, namely the spectral radius of the iteration matrix and the covariance matrix associated to NAG-GS method summarized by Algorithm [] Note that this theoretical study only concerns the case f(x) 1 Az. Considering the size limitation of the paper, we present below only the main theoretical result and place its proof in Appendix A.1.4.
Theorem For GNAG 7 given 7 2 and assuming 0 < x S... A = LL < co: if
0 0<as P- //--99 (/ then the NAG-GS method summarized by Algorithm is convergent for the n-dimensional case, with 2.
Remark important to mention that theoretical result about the accelerated convergence for the stochastic setting holds exactly in the same way as for the deterministic setting.
All the steps of the convergence analysis are fully detailed in Appendix A.I and organized follows:
Appendices A. and 1.2 respectively provide the full analysis of the spectral radius of the iteration matrix associated NAG-GS method and the covariance matrix at stationarity w.r.t. all (hyper)parameters μl, L, y and for the dimensional case 2. The theoretical results obtained are summarized in Appendix 1.1 to come up with an optimal (constant) step size in terms of contraction rate.
Numerical tests performed and detailed in Appendix A.1.5 to support the theoretical results obtained for the quadratic case.
3 EXPERIMENTS
We test NAG-GS method on several neural architectures: logistic regression, transformer models for natural language processing tasks, and residual networks for computer vision tasks. Section 3.] presents the numerical results for logistic regression. Sections 3.2 and 3.3 focus on the tests carried out on transformers and residual networks. For these two neural architectures, order to perform benchmark of our method as fair as possible, replace the reference optimizers with ours and adjust only hyperparameters of our optimizer. We keep intact the model architectures and model hyperparameters such as dropout rate, schedule, batch size, number of training epoches, evaluation methodology. Moreover, we carry out ablation study on small real-world models. The results of the benchmark for transformer models and residual networks summarized in Table
Furthermore, in Section 3.4, the experiments highlight the importance of updating 7 during the training process. In Section 3.5 we present some preliminary experiments study the relations
5

---

Under review as a conference paper at ICLR 2023
between convergence and Hessian spectrum. Moreover, in Appendix C.2 bring some numerical evidences that support the theoretical constraints on the optimizer parameters that have been derived in Section 2
The implementation details are discussed in Appendix C.3
3.1 LOGISTIC REGRESSION
In this section, we benchmark NAG-GS method against state-of-the-arl optimizers on the logistic regression training problem for MNIST dataset et al. (2010). Since this problem is convex and non-quadratic, we consider this problem as the natural and next test case after the theoretical analysis and numerical tests of NAG-GS method in Section 2.3 for the quadratic convex problem. In Figure 2 we present the comparison of NAG-GS method with competitors. We confirm numerically that NAG- than SGD Momentum and Adam W optimizers, highlighting the robustness of our method w.r.t. the selection of hyperparameters. Moreover, the results indicate that the semi-implicit nature of NAG-GS method indeed ensures the acceleration effect through the use of larger learning rates while keeping good accuracy of the model, and this holds not only for the convex quadratic problems but also for non-quadratic convex ones.
GS method allows the of larger range of values for the learning
Table 2: Test accuracies for NAG- GS, SGD-Momentum and AdamW.
Test neaneyy
0.92
uecuracy
3.30
8.00
Ir
NAG-GS
SGD
Adam W
0.35
AdkaW
0.85
Adimh
10-3
0.8934
0.9190
0.9254
101
10
10
10"1
10
0.9207
0.9224
0.9069
Learning
Learning
0.1
0.9249
0.8759
0.8425
(b) μr 10-5
0.5
0.9170
0.8982
0.8638
Figure 2: Dependence of the test accuracy on the learning rates for the considered methods. NAG-GS provides the highest test accuracy for larger step size.
3.2 TRANSFORMER MODELS
this section we test NAG-GS optimizer in the frame natural language processing the tasks of fine-tuning pretrained model on GLUE benchmark datasets Wang et al. (2018). We use pretrained RoBERTa Liu et al. (2019) model from Hugging Face's TRANSFORMERS Wolf et al. (2020) library. In this benchmark, the reference optimizer is Adam Ilya al. (2019) with polynomial learning rate schedule. The training setup defined in Liu et al. (2019) is used for both NAG-GS and Adam optimizers. We search for an optimal learning rate for NAG-GS optimizer with fixed y and y get the best performance on the task at hand. Note that NAG-GS is used with constant schedule which makes it simpler tune. In terms of learning rate values, the one allowed by AdamW W around 10-5 while NAG-GS allows much bigger value of 10 2 Evaluation results on GLUE tasks are presented Table 3 Despite rather restraint search space for NAG-GS hyperparameters, it demonstrates better performance on some tasks and worse performance on others. Figure shows the behavior of loss values and target metrics on GLUE (see Appendix C).
Table Comparison of AdamW and NAG-GS optimizers in fine-tuning on GLUE benchmark. We use reported hyperparameters for Adam W. In case of NAG-GS we search hyperparameters space for the best performance metric. Search space consists of learning rate from [10 10"), factor ? from [10-2 2, 109], and momentum μ
OPTIMIZER CoLA
MNLI
MRPC
QNLI
QQP
RTE
SST2
STS-B
WNLI
ADAMW
61.60
87.56
88.24
92.62
91.69
78.34
94.95
90.68
56.34
NAG-GS
61.60
87.24
90.69
92.42
91.01
73.65
94.50
90.21
56.34

---

Under review as a conference paper at ICLR 2023
1.00
AdamW (train)
0.6
AdamW (train)
A (eval)
AdamW (eval)
Loss
0.75
NAG-GS (train)
Loss
0.5
NAG-GS (train)
NAG-GS
(eval)
NAG-GS (eval)
0.50
0.4
0.25
0.3
0
2
4
6
8
10
0
2
4
6
8
10
Epoch
Epoch
Figure 3: Cross-entropy losses on validation and train sets for CoLA (left) and MRPC (right) tasks.
Solid lines corresponds to the best trial with NAG-GS optimizer.
3.3 RESIDUAL NETWORKS
We compare NAG-GS momentum SGD with weight decay (SGD-MW) residual networks (ResNets) He et al, (016). We choose ResNet20 for versatile experimental verification properties of our optimizer and ResNet50 to test scaling of NAG-GS on large datasets and models which used on practice.
ResNet50. Due limitation of computational resources, train ResNet50 on ImageNet Deng et al. (2009) for a few learning rate values along with fixed y and μ.. We select the best overall performing parameters for optimizer while, for the reference optimizer SGD-MW, reported in the literature. Also, we found that, unlike fine-tuning on Transformers, piecewise-conste conssant schedule gives better target metrics than constant one. The best acc @ performance for NAG-GS is within three percentage points worse than for SGD-MW (see Table 1 and Figure 4b).
use the parameters
ResNet20. Since ResNet20 has lesser number of parameters than ResNet50, intensive experiments in order to evaluate more deeply the performance of NAG-GS for computer vision tasks (residual networks in particular) and to show that NAG-GS with appropriate choice of optimizer parameters is on par with SGD-MW (see Table and Figure 4a) The classification problem solved using CIFAR-10 Krizhevsky 10
carried out more
(2009). Experimental setup is the same in all experiments except optimizer and its parameters. The best test score for NAG-GS is achieved for a = 0.11,
7 17, and μ 0.01.
3.4 UPDATABLE SCALING FACTOR 7
According the theory of NAG-GS optimizer presented in Section 2 the scaling factor decays exponentially fast to μ and, in the case To μ., remains constant along iterations. So, natural question arises: the update on y necessary? Our experiments confirm that scaling factor 7 should be updated accordingly to Algorithm even in this highly non-convex setting, in order to get better metries on test sets.
We use experimental setup for ResNet20 from Section 3.3 and search for hyperparameters for NAG- GS with updatable 7 and with constant one. Common hyperoptimization library OPTUNA Akiba et al. (2019) are used with a budget of 160 iterations to sample NAG-GS parameters. Figure 5 plots the evolution of the best score value along optimization time. The final difference is about 0.5 which is significant difference in terms of final classification error.
3.5 NON-CONVEXITY AND HESSIAN SPECTRUM
Theoretical analysis of NAG-GS highlights the importance of the smallest eigenvalue of the Hessian matrix for convex and strongly convex functions. Unfortunately, the objective functions usually considered for the training of neural networks are not convex. In this section we try to address this issue. The smallest model in our experimental setup is ResNet20. However, we cannot afford to compute exactly the Hessian matrix since ResNet20 has almost 300k parameters. Instead, use Hessian-vecte vector product (HVP) Hx)) and apply matrix-free algorithms for finding the extreme eigenvalues. We estimate the extreme eigenvalues of the Hessian spectrum with power iterations
7

---

Under review as a conference paper at ICLR 2023
6
SGD-MW NAG-GS
SGD-MW
1.50
5
NAG-GS
1.25
Loss
Loss
4
1.00
3
0.75
Mmwaa
2
0.50
0.25
0.9
0.8
Acc@1
Acc
60
0.6
0.5
SGD-MW
SGD-MW
NAG-GS
40
NAG-GS
100
150 200
250
Epoch
Epoch
(a) ResNet20
(b) ResNet50
Figure 4: Evaluation of NAG-GS with SGD-MW on ResNet-20 (left) and ResNet-50 (right) on
CIFAR-10 and ImageNet respectively.
1.0
Classification Error
constant
0.6
updatable
0.4
0.2
0.0
25
50
75
100
125
150
0
25
30
75
100
125
150
Number of samples
Number of samples
Figure 5: The best @ 1 on test set for updatable and fixed scaling factor y during hyperoptimization. NAG-GS with updatable 7 gives more frequently better results than the ones obtained with constant ?. The final difference about 0.5.
(PI) along with Rayleight quotient (RQ) Golub & van Loan 2013). PI is used to get a good initial vector which used later in the optimization of RQ. In order get more useful initial vector for the estimation of the smallest eigenvalue, we apply the spectral shift H(r) Amax and use the corresponding eigenvector.
Figure O shows the extreme eigenvalues of ResNet20 Hessian at the end of each epoch for the batch
size 256 in the same setup as Section 3.3 The largest eigenvalue strictly positive while the smallest one is negative and usually oscillates around μ -1 It turns out that there is an island of hyperparameters the vicinity of that μ. We report that training ResNet20 with hyperparameters included in this island gives good target metrics. The domain of negative momenta is non-conventional and not well understood, to the best of our knowledge. Moreover, there is no theoretical guarantees for NAG-GS in the non-convex the existence of regions of convergence for SGD with negative momentum, which supports our observations. The theoretical aspects of these observations will be studied in further works.
and negative μ.. However, Velikanov et al. (2022) reports
8

---

Under review as a conference paper at ICLR 2023
200
SCD-MW NAG-GS
150
Eigenvalue
Eigenvalus
4
100
30
SGD-MW
NAG-GS
0
10
0
50
100
150
200
250
300
0
50
100
150
200
250
300
Epoch
Epoch
Figure 6: Evolution of the extreme eigenvalues (the larget and the smallest ones) during training RESNET20 on CIFAR-10 with NAG-GS.
RELATED WORKS
The approach of interpreting and analyzing optimization methods from the ODEs discretization perspective is well-known and widely used in practice Muehlebach & Jordan (019); Wilson et al. (202T): Shi al. (202T). The main advantage of this approach is to construct a direct correspondence between the properties of some classes of ODDs and their associated optimization methods. In particular, gradient descent and Nesterov accelerated methods are discussed (Su et al. (2014) as particular discretization of ODEs. In the same perspective, many other optimization methods were analyzed, we can mention the mirror descent method and its accelerated versions Krichene et al.) (2015), the proximal methods Attouch al. (2019) and ADMM Franca al. 2018). It is well known that discretization strategy is essential for transforming a particular ODE to an efficient optimization method. (Shi et (2019) ; Zhang et al. (2018) investigate the most proper discretization techniques for different classes of ODEs. A similar analysis but for stochastic first-order methods is presented
in Laborde Oberman 2020): Malladi al. (2022).
5 CONCLUSIONS AND FURTHER WORKS
We have presented a new and theoretically motivated stochastic optimizer called NAG-GS. comes from the semi-implicit Gauss-Seidel type discretization well-chosen accelerated Nesterov-like SDE. These building blocks ensure two central properties for NAG-GS: (1) the ability to accelerate the optimization process and (2) a great robustness the selection of hyperparameters, in particular for the choice of the learning rate. We demonstrate these features theoretically and provide detailed analysis of the convergence the method in the quadratic case. Moreover, we show that NAG-GS is competitive with state-of-the-arr methods for tackling wide variety of stochastic optimization problems increasing complexity and dimension, starting from logistic regression model training of large machine learning models such as ResNet20, ResNet50 and Transformers. In all tests, NAG- AAGGGS demonstrates better or competitive performance compared with standard optimizers. Further works will focus the theoretical analysis of NAG-GS in the non-convex case and the derivation efficient and tractable higher-order methods based on the full-implicit discretization of the accelerated Nesterov-like SDE.
REFERENCES
Takuya Akiba, Shotaro Sano, Toshihiko Yanase, Takeru Ohta, and Masanori Koyama. Optuna: next-generation hyperparameter optimization framework. In Proceedings of the 25rd ACM SIGKDD International Conference on Knowledge Discovery and Data Mining, 2019.
Anton Arnold, Peter Markowich, Giuseppe Toscani, and Andreas Unterreiter. On convex sobolev inequalities and the rate of convergence to equilibrium for fokker-planck type equations. 2001.
Hedy Attouch. Zaki Chbani, and Hassan Riahi. Fast proximal methods via time scaling of damped inertial dynamics. SIAM Journal Optimization, 29(3):227--2566 2019.

---

Under review as a conference paper at ICLR 2023
Igor Babuschkin, Kate Baumli, Alison Bell, Surya Bhupatiraju, Jake Bruce, Peter Buchlovsky, David Budden, Trevor Cai, Aidan Clark, Ivo Danihelka, Claudio Fantacci, Jonathan Godwin, Chris Jones, Ross Hemsley, Tom Hennigan, Matteo Hessel, Shaobo Hou, Steven Kapturowski, Thomas Keck, lurii Kemaev, Michael King, Markus Kunesch, Lena Martens, Hamza Merzic, Vladimir Mikulik, Tamara Norman, John Quan, George Papamakarios, Roman Ring, Francisco Ruiz, Alvaro Sanchez, Rosalia Schneider, Eren Sezener, Stephen Spencer, Srivatsan Srinivasan, Luyu Wang, Wojciech Stokowiec, and Fabio Viola. The DeepMind JAX Ecosystem, 2020. URL
com/ deepmind
Aleksandar Botev, Hippolyt Ritter, and David Barber. Practical Gauss-Newton optimisation for deep learning. In Proceedings of the 34th International Conference on Machine Learning, volume of Proceedings of Machine Learning Research, pp. 557-565. PMLR, 06-11 Aug 2017.
James Bradbury, Roy Frostig, Peter Hawkins, Matthew James Johnson, Chris Leary, Dougal Maclaurin, George Necula, Adam Paszke, Jake VanderPlas, Skye Wanderman and Qiao Zhang. JAX: Composable Transformations of Python+NumPy Programs, 2018. URL
http: /tpp:// com/ coom//goole/jax
Sinho Chewi, Thibaut Le Gouic, Chen Lu, Tyler Maunu, and Philippe Rigollet. Svgd as a kemel- ized wasserstein gradient flow of the chi-squared divergence. Advances in Neural Information Processing Systems, 33:2098-2109, 2020.
Jia Deng, Wei Dong, Richard Socher, Li-Jia Li, Kai Li, and Li Fei-Fei. Imagenet: A large-scale hier- archical image database. In 2009 IEEE Conference on Computer Vision and Pattern Recognition, pp. 248-255, 2009. doi: 10. N0.1109/2PP.220005520000
Guilherme Franca, Daniel Robinson, and Rene Vidal. Admm and accelerated admm as continuous dynamical systems. In International Conference on Machine Learning, pp. 1559-15000 1567. PMLR, 2018.
Gene H. Golub and Charles F. van Loan. Matrix Computations. JHU press, 2013.
Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun. Deep residual learning for image recognition, In Proceedings of the IEEE conference on computer vision and pattern recognition, pp. 770-778, 2016.
Loshchilov Ilya, Hutter Frank, et al. Decoupled weight decay regularization. Proceedings of ICLR, 2019.
Walid Krichene, Alexandre Bayen, and Peter L Bartlett. Accelerated mirror descent in continuous and discrete time. Advances in neural information processing systems, 28, 2015.
Alex Krizhevsky. Learning Multiple Layers of Features from Tiny Images. Master's thesis, University of Toronto, 2009.
Maxime Laborde and Adam Oberman. A Iyapunov analysis for accelerated gradient methods: from deterministic to stochastic case. International Conference on Artificial Intelligence and Statistics, pp. 602-612. PMLR, 2020.
Mare Lambert, Sinho Chewi, Francis Bach, Silvère Bonnabel, and Philippe Rigollet. Variational inference via wasserstein gradient flows. arXiv preprint arXiv:2205. 15902, 2022.
Jonas Latz. Analysis of stochastic gradient descent in continuous time. Statistics and Computing, 31 (4):1-25, 2021.
Yann LeCun, Corinna Cortes, and CJ Burges. Mnist handwritten digit database. ATT Labs [Online]. Available: http://yann.ll wwwplyannneecun...m www.://enncccm.eee. com/exdb/mmist, 2, 2010.
Kenneth Levenberg. A method for the solution of certain non- non-linear problems in least squares. Quart. Appl. Math., 2. 1944.
Yinhan Liu. Myle Ott, Naman Goyal, Jingfei Du, Mandar Joshi, Danqi Chen, Omer Levy, Mike Lewis, Luke Zettlemoyer, and Veselin Stoyanov. Roberta: A robustly optimized bert pretraining approach. arXiv preprint arXiv: 1907. 1rXiv 1692, 2019.
10

---

Under review as a conference paper at ICLR 2023
Hao Luo and Long Chen. From differential equation solvers to accelerated first-order methods for convex optimization. Mathematical Programming, pp. 1-47, 2021.
Sadhika Malladi, Kaifeng Lyu, Abhishek Panigrahi, and Sanjeev Arora. On the sdes and scaling rules for adaptive gradient algorithms. arXiv preprint 20l.:2220./0088) 2022.
Donald W. Marquardt. An algorithm for least-squares estimation of nonlinear parameters. Journal of the Society for Industrial and Applied Mathematics, 11(2):431-441, 1963.
Gisiro Maruyama. Continuous markov processes and stochastic equations. Rendiconti del Circolo Matematico di Palermo, 4(1):48-90, 1955.
Michael Muehlebach and Michael Jordan. A dynamical systems perspective nesterov acceleration. In International Conference on Machine Learning, pp. 4656-4662. PMLR, 2019.
Yurii Nesterov. Lectures on Convex optimization, volume 137. Springer Optimization and Its Applications, 2018.
Adam Paszke, Sam Gross, Soumith Chintala, Gregory Chanan, Edward Yang, Zachary DeVito, Zeming Lin, Alban Desmaison, Luca Antiga, and Adam Lerer. Automatic Differentiation in PyTorch, 2017.
Fernando Pérez-Cruz. Kullback-leible divergence estimation of continuous distributions. 2008 IEEE international symposium on information theory, pp. 1666-1670. IEEE, 2008
Bin Shi, Simon S Du, Weijie Su, and Michael I Jordan. Acceleration via symplectic discretization of high-resolution differential equations. Advances in Neural Information Processing Systems, 32, 2019.
Bin Shi, Simon S Du, Michael I Jordan, and Weijie J Su. Understanding the acceleration phenomenon via high-resolution differential equations. Mathematical Programming, pp. 1-70, 2021.
Weijie Su, Stephen Boyd, and Emmanuel Candes. A differential equation for modeling nesterov's accelerated gradient method: theory and insights. Advances in neural information processing systems, 27, 2014.
Maksim Velikanov, Denis Kuznedelev, and Dmitry Yarotsky. A view of mini-batch sgd via generating functions: conditions of convergence, phase transitions, benefit from negative momenta. 2022. doi: 10.48550/arxiv.22000 124.
11124
Alex Wang, Amanpreet Singh. Julian Michael, Felix Hill, Omer Levy, and Samuel R Bowman. Glue: multi-task benchmark and analysis platform for natural language understanding, arXiv preprint arXiv:/804.0740 2018.
Ashia C Wilson, Ben Recht, and Michael Jordan. A lyapunov analysis of accelerated methods in optimization. J. Mach. Learn. Res., 22:113 11-11, 2021.
Thomas Wolf, Lysandre Debut, Victor Sanh, Julien Chaumond, Clement Delangue, Anthony Moi, Pierric Cistac, Tim Rault, Rémi Louf, Morgan Funtowicz, Joe Davison, Sam Shleifer, Patrick von Platen, Clara Ma, Yacine Jernite, Julien Plu, Canwen Xu, Teven Le Scao, Sylvain Gugger, Mariama Drame, Quentin Lhoest, and Alexander M. Rush. Transformers: State-of-the ooie---ffhhhtttt art natural language processing. In Proceedings of the 2020 Conference on Empirical Methods in Natural Language Pro- cessing: System Demonstrations, pp. 38 45, Online, October 2020. Association for Computational Linguistics. htths:ts::
ac web. www.oo/ttttttoggogooooooooooottt 2020.
demos 6)
Matthew D. Zeiler. Adadelta: An adaptive learning rate method, 2012.
Jingzhao Zhang, Aryan Mokhtari, Suvrit Sra, and Ali Jadbabaie. Direct runge-kutta discretization achieves acceleration. Advances neural information processing systems, 31, 2018.

---

Under review as a conference paper at ICLR 2023
A ADDITIONAL REMARKS RELATED TO THEORETICAL BACKGROUND
An accelerated ODE has been presented in Section I.] which relied on specific spectral transforma- tion. In this brief Appendix, we add some useful insights:
Equation 10 is variant of heavy ball model with variable damping coefficients in front of # and i:
Thanks to the scaling factor y both the convex case μ = 0 and the strongly convex case / >> 0 can be handled in a unified way.
In the continuous time, u)e- 1. 1.2 O. Since "/o > 0, we have that 7(4) > 0 for all 2 0 and 7(1) converges to exponentially and monotonically as - +00. In particular, if 700 = Jr > 0, then 7(1) /
solve easily equation as follows: 7(t) μ (%o
for all t 2 O. We remark here the links between the behavior of the scaling factor 7t) and the sequence {7k 1xc introduced by Nesterov NNsterov (2018) in its analysis of optimal
first-order methods in discrete-time, see (Nesterov) 2018 Lemma 2.2.3).
Authors from Luo & Chen (2021) prove the exponential decay property C(t) S
ELoo 1> 0 for a Taylored Lyapunov function L(()= == (x((()) f(x*) 2 (e(())
201)
x²122 where x"' E argmin f is global minimizer of f. Again we note the similarity between
the Lyapunov function proposed here and the estimating sequence { { of function f
introduced by Nesterov in its optimal first-order methods analysis Nesterov (2018). In (Nes)
terov 2018 Lemma 2.2.3), this sequence that takes the form (x = (z(x) 2 114 x
where 7h+1 (1 -1-aa ayk + a and Uk+1 -
[1 a&) >x. + ax.uyk ar f
which stand for forward Buler discretization respectively of equation 8 and second ODE of equation 9
We ask the attentive reader to remember that this discussion mainly concerns the continuous time case. A second central part of our analysis was based on the methods of discretization of equation g Indeed, these discretization's 's ensure together with the spectral transformation 7 the optimal convergence rates of the methods and their particular ability to handle noisy gradients.
A.1 CONVERGENCE/STARIIING ANALYSIS OF QUADRATIC CASE: DETAILS
As briefly mentioned in Section 3 the two key elements come up with a maximum (constant) step size for Algorithm| [ are the study of the spectral radius of iteration matrix associated to NAG- GS scheme (Appendix |A.III) and the covariance matrix at stationarity (Appendix (A.1.2) w.r.t. all the significant parameters of the scheme that are the step size (integration step/time step) a, the convexity parameters 0 S μ S L oo of function f(x), the volatility of the noise a and the positive scaling parameter yy Note that this theoretical study only concerns the case f(x) = 1 Ar. Among others results, we will show in particular that for specific intervals for μu and L. higher step size for can be reached while ensuring the A-stability of NAG-GS method. Let remark finally that the two following sections consider the special case n = 2, the theoretical results obtained are summarized in Appendix 1.3 and finally generalized for ne-dimensional case in Appendix A.1.4 with global convergence results for NAG-GS method.
A.1.1 SPECTRAL RADIUS ANALYSIS
Let us assume f(x) = 1/=T Ax and since E S by hypothesis, it is diagonalizable and can be presented as A = diag(A1 (i....,, An) without loss of generality, that is to say that will consider system of coordinates composed of the eigenvectors of matrix A. Let us note that X S ... X; S X, L..
In this setting, y = (I., v) E Rd and the matrices M and N from the Gauss-Seidel Splitting of GN equation 7 Jrre:
1
0
0
Tzx2
02x2
0
-1
0
M
1/7T2x2 A/7
12/2
0
0
p /
0
u/7
L/7
0
4/
12

---

Under review as a conference paper at ICLR 2023
02x2
12x2
02 02x2
02x2
For the minimization of f(x) 1 Ar, given the property of Brownian motion AW;
Wk+1
W = Vouekk where T7s ~ |0, 1), equation 12 reads:
M641 = (14x4 aM) (a11)~((440) (14x4 + aN)ym + (14x4 - aM)
(16)
o Vana
Since matrix M is lower-triangular, matrix 14x4 aM is as well and can be factorized as follows:
1xx4 a.M DT
(1 + 22
02x2 (1 + aK 12x2
12x2 a(A (AAulaxx))
02x2
02x2
T2x2
Hence (1xx--aM)~1
where D-1 can be easily computed. It remains compute T-1;
T can be decomposed as follows: T 11x4 + Q with Q nilpotent matrix such that QQ = O4x4- For such decomposition, it well known that:
12x2
02x2]
T (11x4 Q) L1x1 Q
alulaxa-A)
12x2
(17)
((++T)))
where 7
Combining these results, equation 16 finally reads:
0
1ta
0
0
I+a 0
Mk++1
0
0
yk
+
I+- m
0
a(uuL L)
o (a-D L)
(18)
+[++1]++1]
(1+7)(1+a]
I+T
Eyk
va
1+7 7k
with E denoting the iteration matrix associated to the NAG-GS method. Hence equation 18 includes two terms, the first is the product of the iteration matrix times the current vector Vı. and the second one features the effect of the noise. For the latter, it will be studied in Appendix 1.2 on the point of view of maximum step size for the NAG-GS method through the key quantity of covariance matrix. Let us focus on the first term; it is clear that in order to get the maximum contraction rate, we should look for a that minimize the spectrum radius of E. Since the spectral radius is the maximum of absolute value of the eigenvalues of iteration matrix F,, we start by computing them. Let find the expression of E (F.. for < S that satisfy det Al4x4) as functions of the scheme's parameters. Solving
det(E M4xxa)
( a - 1i 2nd 22 a²X and adp + La - and + axPx a²)))
(19)
(a 1)2(7 or)
leads to the following eigenvalues:
A1 ?
am
A2 1
27+ opt La² a²
A
2(7 + at + ap + a²μ1)
0n 1.a2² 2L 2Lar
27 - a + 20122 + 270 + 122 + 27/14 +
2(7 + an ap a²p))
A4 2 ap La² + a²p
2(7 + 07 ap a²11)
ov 12a2 2T 21a
27 4 + 2² + 20y2 27004 + μ1² + 27/1 + 72
2(7 + or + op + a²n)
(20)
Let us first mention some general behavior or these eigenvalues: Given and μ positive, we observe that:
13

---

Under review as a conference paper at ICLR 2023
1. A1 and 22 are positive decreasing functions w.r.t a. Moreover, for bounded y and μa, have lima-+0 A1(a)] = lima-+co [A2(a)]-
2. one can show that for a E
[1t--2y/3 "+2+2/7h] functions )a(a) and (a are complex L-n values and one can easily show that both share the same absolute value. Note that the lower bound of the interval "
T is negative soon as 7 E [2T easily show that lima- (as(a) 0 and
2 T.2 pT,,, 21.
μ + 2 - S R. Moreover, lim {4(a)) با x(A) 1. The latter limit shows that eigenvalues
plays central role the convergence of the NAG-GS method since it is the one that can reach the value one and violate the convergence condition, as soon as ((A) > 2. The analysis of X also allows us come up with a good candidate for the step size a that minimize the spectral radius of matrix E, especially and obviously at critical point Qmaz μナク+2yテL which is positive since T. 2 by hypothesis. Note that the case I, + 1 gives some preliminary hints that the maximum step size can be almost "unbounded" in some particular case.
レ-
Now, let us study these eigenvalues in more details, seems that three different scenarios must be studied:
1. For any variant of Algorithm I for which
then for all k 2 0 and therefore (a) (a)). Moreover, a 24242/7E 2nt2v/ll can easily check that
L-H
|A (a)) = [2((a)] [)a(aa = (a(())). Therefore a = 2u+2v/u the step size ensuring D-E
the minimal spectral radius and hence the maximum contraction rate. Figure shows the evolution of the absolute values of the eigenvalues of iteration matrix E w.r.t. for such setting.
2. As soon as < ns, one can easily show that (a) < >²(a)) Therefore the step size a with the minimal spectral radius is such that [^л(a)) (22(a)). One can show that the /=
One can easily check that /+>+v/um=)*+*v レ-
equality holds for
u+1+2/7L (μ 7)2 0. Hence the second candidate for step size will be bigger than L-#
the first and the distance between them increasing the squared distance between and μ. Figure [ shows the evolution the absolute values of the eigenvalues of iteration
matrix E. w.r.t. for this setting.
3. For > fr: the analysis of this setting gives the same results as the previous point. According to Algorithm [ y is either constant and equal to either decreasing to μ along iterations. Hence, the case 7 > f will be considered for the theoretical analysis when 7 * μ.
As first summary, the detailed analysis of the eigenvalues iteration matrix E w.r.t. the significant parameters of the NAG-GS method leads us to come up with two candidates for the step size that minimize the spectral radius of E, hence ensuring the highest contraction rate possible. These results will be gathered with those obtained in Appendix A.1.2 dedicated to the covariance matrix analysis.
Let us now look at the behavior of the dynamics expectation; given the properties of the Brownian motion and by applying the Expectation operator E on both sides of the system of SDE's s resulting "averaged" equations identify with the "deterministic" setting studied by Luo & Chen (2021). For such setting, authors from Luo & Chen (2021) demonstrated that, if 0
11 the
Vx(A))
then Gauss-Seidel splitting-based scheme for solving equation 9 A-stable for quadratic objectives in the deterministic setting. We conclude this section by showing that our two candidates we derived
above for step size are higher than the limit (((A)
given (Luo & Chen 2021 Theorem 1). It can
be intuitively understood the case L + μ4, however we give formal proof in Lemma T
14

---

Under review as a conference paper at ICLR 2023
0.9
A1((a) A2(a)
0.7
0.6
0.5
0.4
0.3
0.2
0.1
5
10
15
20
25
30
35
40
Figure 7: Evolution of absolute values of X; w.r.t 7
<
2
(a)
A2(a))
1.6
(a(a))
[A((a)
1.4
1.2
0.8
0.6
0.4
0.2
15
Figure 8: Evolution of absolute values of X; w.r.t o;; 7 < pb.
Lemma Given y O, and assuming 0 P L. then for and > a the following
inequalities respectively hold:
2p 2/
2
L L-1 11
V (((A)
(21)
18
Vlu 7)2 44
2
L- 11
(K(AA)
where xAA) 1
15

---

Under review as a conference paper at ICLR 2023
Proof:
Let us start for the μ 7/, hence first inequality from equation 21 becomes:
2p + L
L
LLIU
( + v L L/u > (L-
pl. L L-P p
= pV -M
which holds for any positive μ, L and satisfied by hypothesis. For the case y μ4, we have:
1+17 7 + V(u 7)2 + 17L
2
L L- 11
VT//K²
2 (L 11) - -7-1
VV (p -7)3 + 471,
VL/P
=(u 7)² + 4wy (++ + 2// /411- 111--) L) (1-1))+7 772
L
Ey 2p² 1²1 2² vp/L+ pL SOLVW/L 1 wp/ Tu L.)
where second inequality hold since L 2 1 and last inequality holds since -1 1u (one can easily check this by using L > μ)). It remains to show that:
2n² + 1 + 1² Vu/L - uLLu/L
- wp/L(p L) +
which holds for any a and L positive (technical details are skipped; it mainly consists in the study of table of signs of polynomial equation in 11).
Since 7 > by hypothesis, therefore inequality
22 + 11 L + 1²2 WpIL μL H/L
-K- Vn/L(u- LL
holds for any a and positive well, conditions satisfied by hypothesis. This concludes the proof.
Furthermore, let us note that both step size candidates, that are { 24+2VAL u+++ / L-p
レーE respectively for the cases 7 = f and y > j show that NAG-GS method converges the case - l with step size that tends to oc, this behavior cannot be anticipated by the upper-bound given by Luo & Chen 2021, Theorem 1). Some simple numerical experiments are performed in Appendix A.1.5 to support this theoretical result.
Finally, based on previous discussions, let remark that a E / xo] when 7 + μ or a E 23 oo] when = μ. we have p(E(a)) = [Aa(a)) and one can show that p(E) is D-j
L-p
strictly monotonically increasing function a for all L > > 0 and y > 0 (see Appendix A. for the discussion).
A.1.2 COVARIANCE ANALYSIS
In this section, we study the contribution to the computation of maximum step size for NAG-GS method through the analysis of covariance matrix at stationarity. Let us start by computing the covariance matrix C' obtained at iteration l + from Algorithm
C++1 E(yk+19Et))
(22)
16

---

Under review as a conference paper at ICLR 2023
0
By denoting E =
let us replace Mk++1 by its expression given in equation 18 equation 22
I+7
writes:
Ck+1 = E(U++1V(+1)
= E ((Eyk (E) Ex)
(23)
=E E (EY:UR
(ExER)
which holds since expectation operator B(.) a linear operator and by assuming statistical indepen- dence between &a and Eyn. On the one hand, by using again the properties of linearity of E and since E is seen as a constant by B(.), one can show that (EynvÉ E²*) = EC&EA, One the other hand, (O, 1), then Equation equation 23
since T7
becomes:
Ck+1 = EC +Q
(24)
02x2
02x2
where Q =
02x2
2ya2 (1+Tu)= 12x2
Let us now look at limiting behavior of Equation equation 24
that is lim CR. Let be C' = limg-+oo Ck the covariance matrix reached in the asymptotic regime.
also referred to as stationary regime. Applying the limit on both sides of Equation equation 24 C
then satisfies
= ECE + +Q
(25)
Hence equation 25 is particular case of discrete Lyapunov equation. For solving such equation, the vectorization operator denoted "is applied on both sides on equation 25 this amounts to solve the following linear system:
(14x412 - E & Q
(26)
011 B
ain B
where A & B
stands for the Kronecker product. The solution is given by:
am1 B
amn B
C (Lax E E)-LQ
(27)
where l stands for the un un-vectorized operator.
Let us note that, even for the 2-dimensional case considered in this section, the dimension of matrix C' rapidly growth and cannot written plain within this paper. For the following, we will keep its symbolic expression. The stationary matrix C' quantifies the spreading of the limit of the sequence (yk), as a direct consequence of the Brownian motion effect. Now we look the directions that maximize the scattering of the points, in other words we are looking for the eigenvectors and the associated eigenvalues of C. Actually, the required information for the analysis of the step size is contained within the expression of the eigenvalues (C). The obtained eigenvalues functions w.r.t. the parameters of the schemes, while their numerator brings less interest for us (supported further), will focus on the their denominator: We obtained the following expressions:
rationale
AC(C)= = N, p(aa,,,L,,,0)
(28)
s.t.i L., = LLa03 1 "2a""=-L'aa" +La"" ALa²
0,)
L 2Lop2 + SyLap 27² La 4y +
2
A²(C) N2(a, N2(a,,,,,,7,0) L,
Dz(a, n. L.
(29)
s.1.D L. a²2² + 30 2 370 + 20p
8yap² 2 47p² 4.2p
As(C') = N3(a, μ4 o,
Dz(a, u1 H,
(30)
D, = 3²p²² + 3a 21 + 3 2 2² + 2a4
8ya4² 2 + 2 20u + 4yp² + 4-²
17

---

Under review as a conference paper at ICLR 2023
\a H, u, u, 1,
(31)
Stt LL'22 μ 7L--- 22a2222 + 4La x2= Ly + +Lapp2 + 8yLap + 2y + 4 + 2
One can observe that:
1. given 14, positive, the denominators of eigenvalues X2 and X are positive as well, unlike eigenvalues λ and X for which some vertical asymptotes may appear. The latter will be studied in more details further. Note that, even if some eigenvalues share the same denominator, is not the case for the numerator. This will be illustrated later in Figures [1 and to ease the analysis.
2. Interestingly, the volatility of the noise defined by the parameter o does not appear within the expressions of the denominators. It gives us hint that these vertical asymptotes are due the fact that spectral radius is getting close Moreover, the parameter a appears only within the numerators and based on intensive numerical tests, this parameter has pure scaling effect onto the eigenvalues AC( when studied w.r.t. without modifying the trends of the curves.
(discussed further in Appendix AA.1.3)
Let us now study in more details the denominator of λ and Xą and seek for critical step size as function 7, y and L at which vertical asymptote may appear by solving:
L²aa 1202 - 7I'a- 2L. 4² 4La²
FyLa²ut+ + 2Lau² 22 4,Lu 47 =0
(32)
=p(2p D)a + (p + (4 D)a² (2p² + 8yp + 2y²) 4y +7)
This polynomial equation has three roots:
01 -7-
1
722 -6/ + 112 4y
02 ?
(33)
L -
03 + + 72 67/ + 112 +
L 2μ
First, it obvious that first root on is negative given assumed nonnegative and therefore can be as:
disregarded. Concerning 0g and (33, those real roots
72 6yp + 12 + 444
=(7 4)² - 47p + 4 2 0
(34)
=(7 11)² 2
- 47(-1T) T.)
which is always satisfied since and 0 μ L by hypothesis. Further, it is obvious that the study must include three scenarios:
scenario I: L - 2/4 < 0, or cquivalently μ > L/2. Given μ and 7 positive by hypothesis, implics that a3 negative and hence can be disregarded. It remains to check if a2 can be positive, amounts to verify if
# + 7 - + 1²2 + 47L <
(( + 7)² <9² 6% + 112 + L
L
=p 2
which never holds by hypothesis. Therefore, for first scenario, there is no positive critical step size at which vertical asymptote for the eigenvalues may appear.
18

---

Under review as a conference paper at ICLR 2023
2. scenario 2: L 2/ O, or equivalently μl L/2. Obviously, 0g positive and hence shall considered for the analysis of maximum step size for our NAG-GS method. It remains to check if ag is positive, that is to verify if the numerator can be negative. We have at the first scenario that 02 is negative as soon as < which is verified by hypothesis. Therefore, only 0g is positive.
3. scenario 3: L 2μ = 0. For such situation, the critical step size is located at oo and can be disregarded as a potential limitation in our study.
In summary, potential critical and limiting step size only exists in the case μ L/2, equivalently if k(A) > 2. In this setting, the critical step size is positive and equal to Acrit 1- 12 67pp-p4477 Figures to 10 display the evolution of the eigenvalues A.(C' L-2g
for 1 isisa w.r.t to for the two first scenarios, that are for μ > L/2 and μ < L/2. For the first scenario, the parameters 7 and L have been respectively set to (1,3/2, 3/2/ For the second scenario, a, 7, μ and have been respectively set to 11,3/2,1, 1,3). As expected, one can observe on Vigure9tthr no vertical asymptote is present. Furthermore, one can observe (c(C") seem to convergence to some limit point when a + oo, numerically we report that this limit point zero, for all the values of y and considered.
Finally, again as expected by the results presented in this section, Figure 10 shows the presence of two vertical asymptotes for the eigenvalues a1 and X4, and none for X2 and X3. Moreover, the critical
step size approximately located at a 6, the theoretical formula Ckeri predicts 6 Finally, one can observe that, after the vertical asymptotes, all the eigenvalues converge to some limit points, again numerically report that this limit point is zero, for all the values of y and a considered.
1-2p
Evolution of A₂(C") > L/2
10
CC
CC)
As(C)
(C)
2
-8
-10
4
6
10
Figure 9: Evolution of (C')") w.r.t for scenario μ L/2, 1, 3/2, p 1, 3/2.
A.1.3 A CONCLUSION FOR THE 2-DIMENSIONAL. CASE
In Appendix I.I.I and Appendix A.1.2 several theoretical results have been derived for coming up with appropriate choices constant step size for Algorithm[1] []] Key insights and interesting values for the step size have been discussed from the study of the spectral radius of iteration matrix E and through the analysis of the covariance matrix in the asymptotic regime. Let us summarize the theoretical results obtained:
from the spectral radius analysis of iteration matrix E; two scenarios have been highlighted, that are:
19

---

Under review as a conference paper at ICLR 2023
Evolution of A,(C) μ < L/2
10
CC
(C
(C)
6.01 0.1173061
-2
-8
6
10
a
Figure 10: Evolution of (C) w.r.t for scenario μ L/2; = - 11 1, Lu = 3.
1. case 7 = μr: the step size a that minimize the spectral radius of matrix E is 24+2/aa I-
2. case 7 Jr: the step size that minimize the spectral radius of matrix E is #+7+ (1-
1-p
from the analysis of covariance matrix C' at stationarity: in the case L 2p > 0, or equivalently μ L/2, we have seen that there vertical asymptote for two eigenvalues of leading to an intractable scattering of the limit points {yx he-to generated by Algorithm I In the case μ L/2, there is positive critical step size at which vertical asymptote for the eigenvalues may appear.
C' at Olerit = M+7+ (22-6yrt++444
L-2p
Therefore, for quadratic functions such that μ L/2, we can safely choose either 2n+2vEL
L-
when 7 = fa either a #tッ+ +1u- when y >> y to get the minimal spectral radius for iteration and hence the highest contraction rate for the NAG-GS method.
For quadratic functions such that < L/2, must show that NAG-GS method is stable for both step sizes. Let us denote by aa = /2442//L
(1ng))+++1 ), two values step size for the L-p two scenarios 7 = and > p. In other to prove stability, we must show that p( (((o.)) < 1. Let us start by computing a such that p(E(a)
レ-ド
1. As proved in Appendix A.1 for E [atcr oo].
p(E(a)) = - with A given in equation 20 we then have to compute a such that:
-A 27 + az apr La² a²
2(7 + an + +ap a²2+)
a(I2a² 2La² 2Lap - Ly - 4 + a² + 2ay² + 2yap + 1² + 27/4 + 72)1/2
2( ar + az + a²n)
This leads to computing the roots of quadratic polynomial equation in the positive root is:
+ AL + 72 67p 2
(35)
L- 2μ
which not surprisingly identifies to Osm from the covariance matrix analysis Furthermore, by recalling that p( (()))) strictly monotonically increasing function over the interval [oce, oo]. showing
It explains why the critical does not include this singularity is due to the spectral radius reaching the value
20

---

Under review as a conference paper at ICLR 2023
that p(E(ae)) in Lemma
is equivalent to show that ag is strictly lower than Ocnit. The formal proof is given
Lemma Given y > O. and assuming 0 < μ < L/2. then for = μ and y > i the following
inequalities respectively hold:
222 67/p + 11² 4
2p 2v
L - 2/4
-- p
(36)
v22 -67/4 + 11² + 44
/ + 7 V/u - 7)2 4
L 24
L-
Proof:
Let us focus the case y > ,4, since 0 < μ < L/2 by hypothesis, second inequality from
equation 36 can be written as:
H)) - -L + # +47 > 0 Ey p² + - +
(L - (L + v(n
+ 74L 6p.) +2p (omu-D 11) 47
Given y. > 0, it remains to show that:
(L -= + 7(4 114L-664 6p) + (2 ())
0
(37)
In order show this, we study the conditions for 7 such that the left-hand side of equation 37 is positive. With simple manipulations, one can show that canceling the left-hand side of equation 37 boils down canceling the following quadratic polynomial:
(L p) 2² 7 6pr + (2μ L) - 44
= (-3p= + 2L) + - 4L3)7 2Lp -3p =
The two roots
212 2/ 2μ + L 4μd² + 412 + 4 4nL
2L 3/1
-22 18-227 222 + 2 2p + L 4pL3 + 442 L2 + 1 L + 4nL
72
2L 34
which are reals and distinct as soon as:
- 24 + L 4uL3 + 44 L2 + 1 0>0
=(L - 2 - -)) L μL) > 0
which holds since 0 <4 < L by hypothesis (one can easily show that -11² + T,2 - -uLL is positive
in such setting). Moreover, the denominator 2L- 3/ is strictly positive since 14 L/2. One can check that 71 is negative for all L > 0 and 0 < < L/2 (simply show that -u² 222 + 4uL is negative) and can be disregarded since 7 is positive by hypothesis. Therefore, proving that equation 37 holds is equivalent to show that:
22 - 21 + 2/ (L- - 2p)(L 1- + L - UL) + 4pL
7
(38)
2I. 3p
remains show that:
リ> -12 1-2-210 21? + +22/2/1 (T, (T.)((T [ -
11) + 12 pp.) + 41u
- -3 3p
=0> 2 VIL (I-24)((1 (1)()-1-pppp 12 ppL) 122 μIb
2 1 1L > (L (2 11)
E/ < b
which holds by hypothesis. Since y > μ by hypothesis. therefore inequality equation 38 holds for any i and I. positive as well, conditions satisfied by hypothesis.
21

---

Under review as a conference paper at ICLR 2023
L-p Finally, since w4rトt((11)4
44742/7L for any 7,44 L 0, then first inequality in equa- L-
tion 36 holds. This concludes the proof.
We conclude this section by discussing several important insights here-under:
Except for Oggiit we do not report significant information coming from the analysis of (C)) for the computation of the step size and the validity of the candidates for that are
{ 2ut22/LL 2tY+ V /+t++((u-))+++7 (( respectively for the cases 7 = f and > μ..
1-
1-
Concerning the effeet of the volatility of the noise, have mentioned earlier that the parameter a appears only within the numerators X₂(C) and based on intensive numerical tests, this parameter has a pure scaling effect onto the eigenvalues >;(C') when studied w.r.t. a without modifying the trends of the curves. For compliance purpose, Figures [1 and respectivell show the evolution of the numerators N,(o, μr, I, a) of eigenvalues expressions of C' given in Equations equation 28 to equation 31 w.r.t. a, for both scenarios μ L/2 and jr >> L/2. One can observe monotonic polynomial increasing behavior of Ni 14, L, σ) w.r.t ø for all 1 S i S 4.
The theoretical analysis summarized in this section is valid for the 2- d-mmenssiaaal case, we show in Appendix A..14 how generalize our results for the n-dimensional case. This has no impact on the our results.
Evolution of Ni(a, 1,,,) p > L/2
14000
Ni
N2
12000
N3
NA
10000
8000
6000
4000
2000
2
3
5
Figure 11: Evolution of N((a, ,1, 1, w.r.t for scenario μ L/2; 3/2, μ L- 3/2,
p+ (1r L-
A.1.4 EXTENSION N-DIMENSIONAL CASE
In this section we show that we can easily extend the results gathered for the 2-dimensional case in Appendix A. and that for NAG transformation the general SDE's system to solve for the quadratic case is:
andd Appendix A.1.2 to the n-dimensional case with 2. Let start by recalling
I Iuxa i(() 1/7(ulnxn A) -n/7inxn
Inxn
O Oax d2
v(1)
t> O.
(39)
Let recall that y = (z, v) with 2,v E R", let n be even and let consider the permuta-
tion matrix P associated
permutation indicator given here-under in two-line form:
(1 2) (3 4)
(n-1 n) (n+1
n 2) ... (2n - 2n)
(2+1-1 1) [2+ 3 2+3) ... (2n 2n - 2)
(3
4)
(2n-1 2n)
where the bottom second-half part of corresponds the complementary of the bottom first half
22

---

Under review as a conference paper at ICLR 2023
Evolution of Ni μ. 1.0,) P < L/2
3500
N1
3000
N2
N₃
2500
NA
2000
1500
1000
2
3
Figure 12: Evolution of (1 n4 L, 7, w.r.t for scenario μ L/2: 7 3/2, μ 1, ==
147+ u+ v(w-7))44711
1
w.r.t to the set 21,2.... 2n} in the increasing order. For avoiding ambiguities, the ones element of P are indices (I, (2. for S j 2n. For such convention and since permutation matrix P associated to indicator ı is orthogonal matrix, equation 39 can be equivalently written as follows:
-Inxn
Taxn 1//7Inxn
0nnx
i(t) = PPT
PPT y(t) +
1 1/>(uTax
A)
(40)
-Iaxn
Inxa -u/7hnxn
0nx1
=PT PT
PPT v(t) -PT
[1/7(1nxx A)
Since we assumed w.l.o.g. diag( AI.... An) with == A S 21 S S An L, one
can easily see that Equation equation 40 has the structure:
21
lo
21
0
0
0
1//(pd2 A.)
-4/7/
021-1
2
2x.
0
Iu
1/sluls And)
-p/y
2
2n
(41)
which boils down m -t- independent 2-dimensional SDE's systems where A, ( with 5 < m such that Àς = μ and A,. = L.
Therefore, the SDE's systems can be studied and theoretically solved independently with the schemes and the associated step sizes presented in previous sections. However, in practice, we will tackle the full SDE's system 39
Let now use the "decoupled" structure given in equation 41 to come up with a general step size that will ensure the convergence of each system and hence the convergence of the full original system given in equation 39 Let us denote by a; the maximum step size for the i-th SDE's system with 1 < i m n/2. For convenience, let consider the case 7 > μt, apply the same method as
23

---

Under review as a conference paper at ICLR 2023
detailed in Appendix A.1.1 and Appendix A.1.2 to compute the expression of a;, obtain:
a; 1 + 7 + V(u 7)2 47/2s
A2i 1
1--1 144
ensures the convergence of NAG-GS method used to solve the SDE's system 39 in the n-dimensional case for n > 2. Theorem enunciated in Section 2.3 and the proof is given here-under.
Finally, in Theorem we show that choosing a
Proof: Let us consider the SDE's system in the form given by equation 41 and let a;
++ +(u=>)++ be the step size selected for solving the i-th SDE's system with 1 S
Azi -4
m = n2/2. In order to prove the convergence of the NAG-GS method by choosing a single step size a (w- it suffices to show that:
such that 0
+ + ((p 7)2 + 4
min 0i Simm=
(42)
combined with results from Appendix A.1.3 and from Lemma sufficient to show that for any À such that O < P S
For proving that equation 42 holds,
we have:
# # + + Vln 7)² + 4 1 1 +
7)² + 4X
(43)
--
-
which is equivalent to show:
x+++ + v (n 7)² + 44
1 N+
v Vp - 7)2 +
$
A - 1
v (1s - --)²= 4TI.
v (u - 7)2 4
1
D-
L- μ
X p
(44)
Since < S < < xo by hypothesis, one can easily show that first two terms of the last inequality are negative. It remains to show that:
v (1 7)
Vlu 7)2 47d
00
1- P
(45)
2 2 44 + 22p - -1 + (4722 + 2 1+ 2 L+
72x² 1:0-222220
21²
- 2 27X22 + x²21²² 2Apz 0
Note that we can easily show that the coefficient of L is negative, hence last inequality is satisfied
as soon as T, < ==>++>==+22 1yA- 21 this concludes the proof.
- or I 2 X. The latter condition is satisfied by hypothesis,
~ 24 22 < X.
Note that one can check that
7 2
The theoretical results derived in these sections along with the key insights are validated Ap-
pendix A.1.5 A.1.5 through numerical experiments conducted for NAG-GS method in the quadratic case.
A.1.5 NUMERICAL TESTS FOR QUADRATIC CASE
In this section we report some simple numerical tests for NAG-GS method (Algorithm M used to tackle the accelerated SDE's system given in I1 where:
the objective function f() = (x - ce) ce) with $3 a all- all-ones vector of dimension 3 and c a positive scalar. For such strongly convex setting, since the feasible set is = R-, , the minimizer arg min ff uniquely exists and simply equal to ce; will be denoted further by x². The matrix A is generated as follows: A = QAQ-1 where matrix is a diagonal matrix of size 3 and Q is random orthogonal matrix. This test procedure allows specify the minimum and maximum eigenvalues of A that are respectively and L and hence it allows us to consider the two scenarios discussed in Appendix A. that are f T./2 and μ T./2.
24

---

Under review as a conference paper at ICLR 2023
The noise volatility is set to 1, we report that this corresponds to significant level of noise.
Initial parameter "/o is set to jr.
Different values for the step size a will be considered in order to empirically demonstrate the optimal choice a terms of contraction rate, but also validate the critical values for step size in the case μ < T./2 and, finally, highlight the effect of the step size in terms of scattering of the final iterates generated by NAG-GS around the minimizer of f.
On a practical point of view, we consider m 200000 points. For each of them, NAG-GS method is ran for a maximum number of iterations allowing to reach the stationarity, and the initial state zo is generated using normal Gaussian distribution. Since ((x a quadratic function, it expected that the points will converge to some Gaussian distribution around the minimizer x = ce. Furthermore, since the initial distribution is also Gaussian, then expected than the intermediate distributions (at each iteration of NAG-GS method) are Gaussian well. Therefore, in order to quantify the rate of convergence of NAG-GS method for different values of step size, we will monitor || - xx"|l ----=x1l,, that is the distance between the empirical mean of the distribution at iteration & and the minimizer x of
Figures 13 and 14 respectively show the evolution of | - x"I along iteration and the final distribution of points obtained by NAG-GS stationarity for the scenario μ L/2, for the latter the points are projected onto the three planes to have a full visualization. As expected by the theory presented in Appendix A.1.3 there no critical a, hence one may choose arbitrary large values for step size while NAG-GS method still converges. Moreover, the choice of a ae gives the highest rate of convergence. Finally, one can observe that the distribution of limit points tighten more and more around the minimizer x² of f as the chosen step increases, as expected by the analysis of Figure 9 Hence, one may choose a very large step size a so that the limit points converge to x* almost surely but at a cost of a (much) slower convergence rate. Here comes the tradeoff between the convergence rate and the limit points scattering.
Finally, Figures and 16 provide similar results for the scenario < L/2. For such scenario, the theory detailed Appendix A..13] and Appendix [A.1.4 predicts a critical a from which convergence of NAG-GS is destroyed. In order to illustrate this gradually, different values of a have been chosen within the set {a.; a./2, 2, (a. + aenit) 0. First, one can observe that the choice of a = ae gives again the highest rate of convergence, see Figure 15 Moreover, one can clearly see that for a - O the convergence starts to fail and the spreading of the limit points tends to infinity. We report that for a Oerit, NAG-GS method diverges. Again, these numerical results are fully predicted by the theory derived in previous sections.
1.6. MONOTONICITY OF SPECTRAL RADIUS OF FOR NAG-GS METHOD
Let us recall that on [ac, oo). the spectral radius p(E(a) is equal to |Ag). the expression of x4 function of parameters of interests for the convergence analysis of NAG-GS method was given in equation 20 and recalled here-under for convenience:
27 2 +y ap La? az
A
40 a
a 2La2p
-2yLa + + + 2ay2² + 2ya + 4r² + 27/44 + 72
2(7 + ar + ap + a²μ)
(46)
Let start by showing that \ is negative on [ac, oo). Firstly, one can easily observe that the denominator of X4 is positive, secondly let us compute the values for a such that:
27 + a7 + aap app- La² + a² -ー
av LVL222 - 2La2 41 2La La 47 a2p² 2ap2 + 2 + p² + 27p4 72 =
- 272 444 + 7) + a²(22 47L + 27p + ,²) - a² 22 - 4yL + 67/ + n²) =0
=(-4r))a² - 47 + 7)a 412 =
(47)
25

---

Under review as a conference paper at ICLR 2023
||
11, 1 L/2
10a
0./2
10
a. 100
10-1
10
10-3
10-4
25
50
75
100
125
150
175
200
Iteration counter k
Figure 13: Evolution of | - x²Il along iteration for the scenario L/2; 5, ====
1.9 and o for a E {ac. 2,/2, 20c, 10cc) with a L-
2u+2 5.29.
The expression above is negative as soon as a < - < since 7, 1u > by hypothesis. The latter is always satisfied since 2 > 0 by hypothesis. Therefore p(E(a)) = = for [ake 00).
To show the monotonic increasing behavior of p(E(a)) w.r.t. E [ac, o0), remains to show that: d(p(E(a)) d( d(->a) >
(48)
da da
To ease the analysis, let us decompose -a((a) +1(()) t2(a) such that:
(a) 27 + ar + op La + a²p
2(7 +0 aap a²
(2( = evLLaa2-lnnn
2Lau 27 4 a + 2ap2 + 2yau + μ1² + 2n14 - 7
2(7 + 07 + 0y + a²n)
(49)
Let us now show that
(a) da
> and
dta(a) > 0 for any > μ > We first obtain:
de
(a)
(2> + 2μ + 4ap (27 + 07 + aap La² + a²p) (2 + 20y + 20h + 20²
da
7 + |1 2T. 201z
(50)
27 + 207 + 2aj + 20
(La² +7) μ1) 20>(L
2(a + 1) a
which is strictly positive since T. >.> p 0 and > by hypothesis. Furthermore:
dtz(a)
da
(51)
( + P)LL 2(a 1)2 7)2 (12 21 (2)) 20(7 /0]] + 2y(21. (1)
30r) (or)
-1 Lu - + -
The remaining of the demonstration is significantly long and technically heavy in the case 7 > μ. Then we limit the last part of the demonstration for the case = for which we have shown In practice, with respect NAG-GS method
previously than Oc 24742 2/T
26

---

Under review as a conference paper at ICLR 2023
summarized by Algorithm 1 quickly decrease to yu and equality μ y holds for most part of the iterations of the Algorithm, hence this case is more important to detail here.
The first term the numerator of Equation 51 is positive as soon as a 2 we determine the conditions under which the second term of the the numerator of Equation 51 is positive, that is:
the case μ = 7
a ( + 2 - Lu + p(2 2p(2L 11) > 0 = aL + 2n(11² - Lu p² > μu( - 21 + 22n -1)
(52)
First one can see that:
(L(-2 + 2n((1 - Lp + 1 00 2(-2p² + 222 - p)) > 0
(53)
hold soon L > u>> 0 which is satisfied by hypothesis. Therefore, the second term of the numerator of Equation 51 is positive as soon as
/4( 1d(2pp² 2n² + 2p (4))
2/c -V V μ
(54)
(I(-2μ²) 2p( - Lp + 12))
which exists since μ 0 by hypothesis (the second root of equation 53 being negative). Finally, since a E [ate, oc] by hypothesis, dtz(a) is positive as soon as:
da
3/
OC >
L
(55)
2p
06
L-
hold with a, 20-2V/L One can easily show that both inequalities hold as soon as L > p L-
which is satisfied by hypothesis. This concludes the proof of the strict monotonicity pl (((a)) w.r.t.
for E [ac, oo] assuming L > μ > 0 and y μs.
A.2 FULLY-IMPLICI SCHEME
In this section we present an iterative method based on the NAG transformation Gy 7 along with a fully implicit discretization to tackle equation in the stochastic setting, the resulting method shall be referred to as "NAG-FT" method. We propose the following discretization for equation with noise; given step size on 0:
J perturbated
1 1
Uk41 T6+1:
(56)
uk+1
Ub
WK+1 - We
(Ix+1 Uk+1)
Ar +
ag
7
Ykk
ag
As done for the NAG-GS method, on a practical point of view, we will use W - W4 AW,
Vorile where 7/k (O, 1), by the properties of the Brownian motion.
In the quadratic case, that is f() 1/=""
T Az, solving equation i6 equivalent to solve:
Ik Uk
(1 + ax)l Dy (A - pl) (1
/ 2 y
a++1 0k+1
(57)
where Th NOO, Furthermore, the parameter equation is again discretized implicitly:
7 *
Mk+1.
%o 0.
(58)
0
As done for NAG-GS GSGGSS method, heuristically, for general f E L. with μ O, just replace Ar in equation 56 with V (**+1) and obtain the following NAG-FI scheme:
s1.1
#k++1
F1k
Uk+1 - Ik+1,
Ot
(59)
Uk+1
P (I++1 - 0k41) vf((k+1) +
Wh+1 WR
Ock
7
0%k
27

---

Under review as a conference paper at ICLR 2023
From the first equation, get Dk++1 obtain:
8x-1 -
+
+k+11 that we substitute within the second equation,
ox
+ Thakk
(f(x++) + +o Voun
Tk+1
76
(60)
1 1+T
with Tx = 1/aa + pt/p--
Computing Ik+1 is then equivalent to computing a fixed point of the operator given by the righ-hand side of equation 60 Hence, it is also equivalent finding the root of function:
g(u) = Up
ThTk
2A V f(u) + o Vaal
74
(61)
7
with g IR." -> T In order to compute the root of this function, we consider a classical Newton- In Algorithm 2 (.) denotes the Jacobean operator
Raphson procedure detailed in Algorithm 2
Algorithm 2 Newton Raathnn method
Input:
Choose the point ug E R", some Ok, Thi Th 0.
for i 0, 1....... do Compute Jg (u:4) = In + (14m (²f(u)
Compute (uu)) using equation 61
Set ui+1 = - [Jg (u)) ((u))
end for
of function, y equation 6] w.r.t. u, In denotes the identity matrix of size n and V designates the
Hessian matrix of objective function f. Note that the iterative method detailed in Algorithm corre- sponds to a particular variant second-order methods known as Bevenberg-Marquanee (1944); Marquardt (1963) applied to the unconstrained minimization problem min [(x for
Levenberg
twice-differentiable function f. Finally NAG-FI method is summarized by Algorithm 3
Algorithm 3 NAG-FI Method
Input: Choose the point zo R", Ug ro, some 0, p O, 7o 0.
for = 0, 1........ do
Sample 7 (0,
Choose 0g > 0
Set 7 24mauu
1++0x
Set T = 1/ar + 1/7
Compute the root of equation 61 using Algorithm
Set Tx+1 = u
end for
By following similar stability analysis as the performed for NAG-GS, can show that this method is unconditionnaly A-stable as expected by the theory of implicit schemes. In particular, one can show that eigenvalues of the iterations matrix are positive decreasing functions w.r.t. step size a, allowing then the choice of any positive value for a. Similarly, one can show that the eigenvalues of the covariance matrix at stationarity associated NAG-FI method are decreasing functions w.r.t. that tend to 0 as soon as lim a - 20. It implies that Algorithm is theoretically able to generate iterates that converge arg min of almost surely. even in the stochastic setting with potentially quadratic rate of converge. This theoretical result is quickly highlighted in Figure 17] shows the final distribution of points generated by NAG-FI once used in test setup detailed Appendix 1.5) in the most interesting and critical scenario L/2. As expected we report can be chosen as big as desired, we choose here a = 1000 Moreover, for increasing value of a, the final distributions of points more and more concentrated the arg min f.
Therefore, NAG-FI method constitutes a good basis for deriving efficient second-order methods for tackling stochastic optimization problems, which is hard to find in the current SOTA. Indeed, second-orde order methods and more generally some variants of preconditioned gradient methods have recently been proposed and used in the deep learning community for the training of NN for instance.
28

---

Under review as a conference paper at ICLR 2023
However, it appears that there is limited empirical success for such methods when used for training NN when compared to well-tuned Stochastic Gradient Descent schemes, see for instance Botev et al. (2017): ; Zeiler 2012). To the best of our knowledge, theoretical explanations have been brought to formally support these empirical observations. This will be part of our future research directions.
Besides these nice preliminar theoretical results and numerical observations for small dimension problems, there is a limitation of NAG-FI method that comes from the numerical feasibility for computing the root of the non-linear equation 61 that be very challenging practice. We will try to address this issue in future works.
B CONVERGENCE TO THE STATIONARY DISTRIBUTION
Another way to study convergence of the proposed algorithms to consider Fokker Planck equation for the density function p(t, z). We will consider the simple case of the scalar SDE for the stochastic gradient flow (similarly as in equation I1. Here f RR R:
dx = - Vf(d)dt + dz v f(x)dt dt + odW, xx0) ~ p(O,
It is well known, that the density function for sr() p(t, x) satisfies the corresponding Fokker Planck equation:
Op(t, x) v (p(t. x Vf(()) + o²
2 Ap(t,x)
(62)
Ōt
For the equation 62 one could write down the stationary (with - 2o) distribution
2 (-)o
exp (-221(o) dx.
lim ptt,
exp
Z
(63)
zEV
is useful to compare differrnnt optimization algorithms terms of convergence in the probability space, because it allows to study the methods in the non-convex setting. We have address two problems with this approach. Firstly, we need to specify some distance functional between current distribution Pr x(t,x) and stationary distribution p p (z). Secondly, don't ton't need to have access to the densities Pt, p" themselves.
For the first problem will consider the following distance functionals between probability distribu- tions in the scalar case:
Kullback Leibler divergence. Several studies dedicated to convergence in probability space are available Arnold et al. (2001); Chewi et al. (2020): Lambert et al. (2022). We used approach proposed in Pérez-Cruz) (2008) to estimate KL divergence between continuous distributions based their samples.
Wasserstein distance. Wasserstein distance is relatively easy to compute for scalar den- sities. Also, it was shown, that stochastic gradient process with constant learning rate is exponentially ergodic in the Wasserstein LenseLatz (202)).
Kolmogorov - Smirnov statistics. We used the two-sample Kolmogorov-Smmingg Smirnov test for goodness of fit.
To the best of our knowledge, the explicit for- mula for the stationary distribution of Fokker- Planck equations for the ASG SDE equation 11 remains unknown. That is why we've de- cided to get samples from the empirical sta- tionary distributions using Eullr-Maruyama in- tegration Maruyama (1955) with small enough stepsize of corresponding SDE with a bunch of different independent initializations.
22-
-2.5
10.0
Two pits function.
fi(x)
50 (2 log (cosh(x)) 5)2
We tested two functions, which are presented the figure T8 We initially generated 100
M
29
-15
-30
30
15
Freguently modu latedd fiunctin

---

Under review as a conference paper at ICLR 2023
points uniformly in the function domain. Then we independently solved initial value problem in equation IT for each of them with Maruyaaa (1955). Results the integration are presented on the figure 19 One see, that in relatively easy scenario 18a NAG-GS converges faster 19a than gradient flow its stationary distri- bution, while in hard scenario 18b] NAG-GS is
more robust to the large stepsize. 19b
Wasserstein Distance (0:.0")
Knimogorow Smimav statistics (p.
10P
GF Bulls NAC
0.2
0.1
0.0
(a) Results for (x)). 8 10
10
kL
Wasserstein Distonce (p..0")
Kolm ogomogo
GP Baler NAG-S5
0.4
10a1
0.1
50
100
150
200
50
100
150
200
50
100
150
200
(b) Results for fa(x). = 1.5, = 10- 10-3,
=1
Figure 19:
Convergence in probabilities of Euler integration of Gradient Flow (GF Euler) and NAG- GS for the non-convex scalar problems.
C ADDITIONAL EXPERIMENTAL DETAILS
In this section we provide additional experimental details. particular, we discuss little bit more experimental setup and give some insights about NAG-GS well.
Our computational resources are limited to a single Nvidia DGX-1 with GPUs Nvidia V100, Almost all experiments were carried out single GPU. The only exception is for the training of ResNet50 on ImageNet which used all 8 GPUs.
C.1
ROBERTA ON GLUE
After completion of experiments toy models, we applied NAG-GS for the training of large neural networks that could be used practice. We performed a grid search come up with an optimal learning rate a for RoBERTa model (see Table 4. It is worth mentioning that there is stability region for NAG-GS which visible to the naked eye, we highlighted in blue some cells from Table to support this claim. Furthermore, based these simple experiments, it appears that choosing a 10 1 a good way to go for wide variety of tasks.
In the subsequent experiments we limited ourselves to RoBERTa's training on small GLUE tasks in order to study relations between a and y in the context high-dimensionas non-convex objective functions (see Figure 20). One can that there is oblong area in which NAG-GS converges pretty well. The ridge of this area can be empirically described with linear equation log a = a log ? + in perform similar study for other types of model.
logarithmic domain. In Appendix C.2
30

---

Under review as a conference paper at ICLR 2023
Table 4: Performance in fine-tuning on GLUE benchmark for different learning rate a, fixed moment 1, and factor (1. 1.5). Largest metric values for NAG-GS optimizer is highlighted with bold typeface. Performance metric is accuracy for all tasks with exception Matthews correlation for CoLA and Pearson correlation for STS-B. Dash in cell means an absence of run. Higher is better.
7
CoLA
MNLI
MRPC
QNLI
QQP
RTE
SST2
STS-B B
WNLI
1.0
8E-01
0.00
68.38
50.54
52.71
88.42
19.86
56.34
4E-01
0.00
81.13
50.54
52.71
50.92
87.36
56.34
2E-01
51.16
70.59
88.93
52.71
92.55
89.34
56.34
1E-01
59.64
88.73
90.39
91.09
52.71
94.27
89.97
56.34
8E-02
59.41
86.27
91.49
72.56
93.92
90.48
56.34
4E-02
59.16
87.75
91.62
72.56
94.38
90.59
56.34
2E-02
57.13
88.48
91.96
67.87
94.72
90.35
56.34
1E-02
53.21
86.52
91.76
89.36
55,23
93.92
88.91
56.34
8E-03
53.94
86.27
54.51
93.58
89.01
56.34
4m-03
49.25
78. 92
53.07
93.69
86.38
56.34
2E-03
40.73
68.38
52.71
93.35
79.56
56.34
1E-03
0.00
68.38
86.88
86.01
52.71
93.35
58.97
43.66
1.5
8F-01
0.00
35.45
68.38
50.54
89.33
52.71
91.17
24.84
56.34
4E-01
0.00
80.29
68.38
50.54
90.23
52.71
91.86
88.62
56.34
2E-01
0.00
82.01
84.80
90.04
91.01
52.71
93.46
89.44
56.34
IE-01
61.48
86.94
87.75
50.54
90.91
72.56
94.38
90.07
56.34
8E-02
$7.77
87.06
86.76
90.92
90.92
73.65
94.27
90.21
56.34
4E-02
59.97
87.24
89.71
92.42
90.59
70.04
94.50
90.11
56.34
2E-02
56.50
87.25
87.75
91.65
89.89
61.73
94.15
90.03
56.34
1E-02
55.08
86.72
85.29
91.62
88.91
55.23
93.69
88.91
56.34
8E-03
50.79
86.46
83.58
90.99
88.92
56.68
93.35
87.92
56.34
4E-03
48.74
85.51
72.30
90.46
87.87
52.71
93.23
84.79
56.34
2F-03
0.00
84.64
68.38
87.94
86.80
52.71
93.46
78.15
43.66
1F-03
0.00
83.06
68.38
86.66
85.38
52.71
92.78
34.77
43.66
C.2 PHASE DIAGRAMS
In Section 3.5 we mentioned that the lowest eigenvalues of approximatd Hessian matrices evaluated during the training of ResNer20 model GS in the convex case include some conditions on the optimizer parameters a, 7, and In particular is required that p 0 and y > μ.. In order bring some insights about these remarks in the non-convex setting and inspired by Velikanov et al. (2022), we experimentally study the convergence regions of NAG-GE BAG-GS and sketch the phase diagrams of convergence for different projection planes,
negative. Furthermore, theoretical analysis of NAG-
see Figure 21
We consider the same setup in Section 3.3 and use hyperoptimization library OPTUNA Akiba et al. (2019). Our preliminary experiments RoBERTa (see Appendix C.I) shows that should be of magnitude 10-1, With the estimate Hessian spectrum of ResNet20 we define the following search space
a ~ Loguooiooonn11-- 2000), 7 LogUniform(10-1 1,10), / Uniform(-100 100).
We sample fixed number of triples and train ResNet20 model on CIFAR-10. Objective function is top-1 classification error.
We report that there is a convergence almost everywhere within the projected search space onto a-y plane (see Figure 21 The analysis of projections onto a-µ and 7-14 planes brings different conclusions: there are regions of convergence for negative f for some a Ouh and 7 >> Yeh- Also, as it was mentioned in, there is subdomain of negative comparable to domain of positive μ in sense the target metrics. Moreover, the majority of sampled points are located in the vicinity the band Amai j Amar-
31

---

Under review as a conference paper at ICLR 2023
Table 5:
The comparision of single step duration for different optimizers on RESNET20 on
CIFAR-10. ADAM-like optimizers have in twice larger state than momentum SGDs or NAG-GS.
OPTIMIZER
MEAN,
VARIANCE,
REL. MEAN
REL. VARIANCE
SGD
0.458
0.008
1.0
1.0
NAG-GS
1.648
0.045
3.6
5.5
SGD-M
3.374
0.042
7.4
5.2
SGD-MW
3.512
0.037
17.7
4.7
ADAMW
5.208
0.102
11.4
12.6
ADAM
7.919
0.169
17.3
20.8
C.3 IMPLEMENTATION DETAILS
In our work we implemented NAG-GS in PyTorch Paszke et al. (2017) and JAX Bradbury et al. (2018) Babuschkin et (2020). Both implementations online According to Algorithm the size of NAG-GS state equals to number of optimization parameters which makes NAG-GS comparable to SGD with momentum. It worth to that Adam-like optimizers has twice larger state than NAG-GS. Arithmetic complexity of NAG-GS is a comparison of computational efficiency
used our experiments and available
linear (nn) the number of parameters. 3 of common optimizers used in practice. Although forward pass and gradient computations usually give the main contribution to training step, there is a settings when efficiency of gradient updates are important (e.g. batch size or number of intermediate activations small with respect number of parameters).
32

---

Under review as a conference paper at ICLR 2023
X Initial dist.
X Initial dist.
dc
az/2
-5
-5
-5
-5
X
Initial dist.
X
Initial dist.
5
ae +2
a. 10
-5
-5
-5
5
(a) Proiection in XY Dlane.
X Initial dist.
ag/2
Initial dist.
-5
-5
-5
5
-5
X
Initial dist.
X
Initial dist.
a. *2
a. *10
-5
-5
-5
5
(b) Projection in Y plane.
X
Initial dist.
5
ad/
X
Initial dist.
-5
dc
-5
-5
-5
X Initial dist.
X
Initial dist.
a,
a. *10
-5
-5
-5
5
(c) Projection in X plane.
Figure 14: Initial (blue x) and final (red circles) distributions of points generated by NAG-GS method for the scenario μ > L/2; 5, 7 μ 1, L = 1.9 and a for a E [ac,aa/2 20,, 10a.) with ag 20+2vpL L-H 5.29.
33

---

Under review as a conference paper at ICLR 2023
1| 11, / L/2
10a
10
ac
10
ac/2
(ac + aam)/2
10-2
0.98a
10-3
10-4
25
75 100
125
150
175
200
Iteration counter k
Figure 15: Evolution of || - xI along iteration for the scenario μ < L/2; c = 5, = =
L = 3 and
for ac/2, (ae + 02 0.98 with 0e 2ut2v/L 2.73 and
1-
Oerit ++ 7/2 6
4.83.
1-22
34

---

Under review as a conference paper at ICLR 2023
X Initial dist.
X
Initial dist.
a./2
Initial dist.
Initial dist.
-5
(d. + aom
-5
0.98
5
5
-5
0
5
(a) Proiection in XY plane.
X
Initial dist.
X
Initial dist.
-5
a./2
Initial dist.
X Initial dist.
-5
(a. + aom
-5
0.98d
-5
-5
0
5
(b) Projection in Y plane.
X
Initial dist.
-5
ac
-5
a./2
5
X
Initial dist.
X Initial dist.
-5
(a. + aom)/2
-5
0.98dtore
-5
0
5
-5
0
5
(c) Projection in X plane.
Figure 16: Initial (blue x) and final (red circles) distributions of points generated by NAG-GS method for scenario μ L/2; c 5, 7 = μ = L and o = for a E {ac, 0,/2, (a. + aerit/22 2, 0.98 with a. 34t2/pL 2.73 and Cerit #tr+ /12 122-8ttt244n 1-
4.83.
T-22p
35

---

Under review as a conference paper at ICLR 2023
Initial dist.
Initial dist. Final dist.
Initial dist. Final dist.
Final dist.
5
-5 0
5
-5
0
5
-5
0
5
Figure 17: Projection onto XY, YZ and X planes (from left to right) of intial (blue x) and final (red circles) distributions of points gnerated by NAG-FI method scenario μr L/2; 5, 7 = / 1,
= and o = 1 for a = 1000a. with ag 2442yTT 2.73.
L-P
Figure 20: Landscape of accuracy metric of RoBERTa model fine-tuned CoLA (top) MRPC (bottom) with NAG-GS. Hyperparameter optimization algorithms samples learning rate a from [10-3, 10 and factor y from [10-2 10%). Both hyperparameters sampled from log-uniform distribution. Symbol * marks the best points.
102
10a
100
10
10
10
10"
10a
0
50
100
0
50
100
Y
0.0
0.2
0.4
0.6
0.8
1.0
Figure 21: Landscapes of classification error for ResNet-20 model trained on CIFAR-10 with NAG-GS after projections onto - and optimization algorithm samples learning rate a from [10- 10"). factor from [10- 10 and factor μ from | 10, 90). Hyperparameters and 7 are sampled from log-uniform distribution, hyperparameter μ are sampled from uniform distribution.
planes (from left to right). Hyperparameter
36


Reviews:
[{"review_id": "eg-sE9is3C4", "reviewer": "Reviewer_ksyw", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper starts with a variant of accelerated gradient flow considered in [Luo & Chen 2021], and added noise in a degenerate, underdamped Langevin alike fashion. The resulting SDE is then discretized semi-implicitly, resulting in a “stochastic gradient” optimization algorithm. “Convergence analysis” is conducted on quadratic loss. The method is then empirically tested, and performances comparable to existing optimizers were demonstrated on a collection of important downstream tasks.", "summary_of_the_review": "The Gauss-Seidel splitting and semi-implicit discretization is definitely an interesting idea, and I strongly encourage the authors to explore it a bit more, possibly with improved theory so that the benefits could be made quantitative. For now, I’m not very sure about either its algorithmic or theoretical impact, and thus unable to provide a very positive recommendation.", "strength_and_weaknesses": "Strength: \n\n(1) The Gauss-Seidel splitting and semi-implicit discretization is nice. It is great to see important knowledge in numerical analysis applied to optimization / machine learning problems.\n\n(2) I also like the diversity of empirical tests.\n\nWeakness:\n\n(1) Since I could have missed the point, could the authors please explain again why the new method is important, if it produces comparable performance to existing approaches?\n\n(2) Stochastic Gradient is not the same as (deterministic) batch-gradient plus constant covariance Gaussian noise. While there is still active research studying whether the noise of SG is Gaussian or heavy-tailed (e.g., [Gürbüzbalaban et al. The heavy-tail phenomenon in SGD. ICML 21]), at least its covariance should be state-dependent (e.g., [Li et al. Stochastic modified equations and adaptive stochastic gradient algorithms. ICML 17]). The assumption in this paper, namely an isotropic additive noise, is unrealistic. Normally I would suggest to soften the claim of “SGD”, but since an SGD version of the previous result in [Luo & Chen 2021] is the main point, I am a little confused about the contribution.\n\nIn addition, the derivation up to equation (15) is based on isotropic additive Gaussian noise, denoted by $\\eta_k$, but this noise disappeared in Algorithm 1. I can’t find $\\eta$ or even its amplitude $\\sigma$ and became puzzled. My guess is $\\nabla f$ in Algorithm 1 is no longer the deterministic full-batch gradient used before, but a stochastic gradient based on mini-batch. But if this is the case, the derivation is a little disconnected from SG (see my above point(2)).\n\n(3) I wish the theoretical part could be stronger. The current status of the field is well beyond convergence analysis for quadratic objectives (Thm.1). In fact, if I didn’t misunderstand, Thm.1 is for a deterministic setup. Appendix B empirically investigated the stochastic setup, based on a statement that the stationary distribution of (11) remains unknown (and thus implying theoretical analysis is out of reach). But it is a standard procedure to characterize the stationary distribution: (11) is a linear SDE and its stationary distribution is just a Gaussian, whose covariance can be solved for algebraically via Ricatti equation. Besides, I missed why the part of Appendix B in page 29 discussed overdamped Langevin instead of equation (11).\n\nBack to the convergence analysis of the full-blown nonlinear SDE, could the authors please explain why tools in, e.g., [Ma et al. Is there an analog of Nesterov acceleration for gradient-based MCMC? Bernoulli 21] and [Li et al. Hessian-free high-resolution Nesterov acceleration for sampling. ICML 22] cannot be used to quantify the speed of convergence in a more general setup?\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}, {"review_id": "Vm_MUw-xXbj", "reviewer": "Reviewer_AnLK", "summary": "", "questions": "", "limitations": "", "rating": 3, "confidence": 4, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper considers a class of accelerated algorithms and provides some analysis for noisy inputs to the algorithm. It is shown that the algorithms will converge even when presented with gradients corrupted with gaussian noise on quadratic losses. Experimental results on deep learning benchmarks show competitive performance with AdamW.\n", "summary_of_the_review": "The theoretical results here seem to be in a somewhat preliminary state and the empirical results do not seem strong enough to justify acceptance on their own.", "strength_and_weaknesses": "\nThe semi-implicit algorithm is general well-motivated and presented. However I am confused on a couple points: the equations 14 do not appear to match exactly the pseudo-code Algorithm. For example, should we not set $x_{k+1} = (1+\\alpha_k)^{-1}(x_k \\alpha_k v_k)$ rather than $x_{k+1} = (1-\\alpha_k) x_k + \\alpha_k v_k$? Similar issue may be involved with the $v_k$ update.\n\n\nIn regards to the theoretical convergence results, these seem to be asymptotic statements that the algorithm will eventually converge on quadratic losses with gaussian noise in the gradients. This seems rather weak given the current state of the art in acceleration under unbiased gradient noise (e.g. see https://proceedings.mlr.press/v119/joulani20a/joulani20a-supp.pdf), which establishes non-asymptotic convergence rates for general convex smooth function under arbitrary unbiased gradient noise with bounded variance without requiring knowledge of several problem parameters.\n\nIn regards to the experiment, it seems to be slightly worse than SGD on computer vision tasks, and comparable to AdamW on transformer fine-tuning tasks. All in all, this does not seem convincing that we should switch from AdamW. Is it possible perhaps to make an argument that this method requires less tuning, as was hinted at in the text?\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "4: All of the claims and statements are well-supported and correct.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}, {"review_id": "DPu7D7OYGu", "reviewer": "Reviewer_DiaF", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 2, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": " This paper presents the NAG-GS algorithm, a Gauss-Seidel type discretization of the Nesterov-like SDE. NAG-GS maintains the acceleration nature while improving the robustness to the section of learning rate. Numerical simulations are implemented to show that NAG-GS is competitive with state-of-the-art methods. ", "summary_of_the_review": "In summary, I think this is nice work with good intuition and theoretical background that motivates the final algorithm smoothly. I hope to see more detailed theoretical results such as convergence rate and the improvement of the learning rate, as these two properties are stated as the main advantage of the proposed NAG-GS algorithm. It would be fantastic if the theory could be extended to the more general problems such as optimizing functions in $\\mathcal{S}_{L, \\mu, \\mu}^{1,1}$.", "strength_and_weaknesses": "## Strength\n- This paper is well-organized and relatively easy to follow. The proofs seem to be theoretically sound. The theoretical background motivates the algorithm smoothly, and the proposed method is achieving great empirical performance.\n\n## Weaknesses\n- Although simulations are implemented for various settings, the theory only holds for the quadratic functions.\n- The theory only shows NAG-GS is convergent (for quadratic functions). Other properties such as the convergence rate and improvement of the range of learning rate are missing. \n\n===========================\nI have read other reviewers' responses and authors' comments. I basically agree with other reviewers. Although applying semi-implicit discretization to SDE of NAG is novel and worth more future exploration, the current theoretical results provided in this paper are preliminary. Also, the experiments cannot convince me one should switch to NAG-GS in practice in sufficiently many settings. Hence, I'll keep my score. ", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "2: The contributions are only marginally significant or novel.", "empirical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work."}, "extra_text": {}}, {"review_id": "8vJ-ktFZVt", "reviewer": "Reviewer_Dog8", "summary": "", "questions": "", "limitations": "", "rating": 5, "confidence": 3, "soundness": null, "presentation": null, "contribution": null, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "This paper proposes the NAG-GS optimizer. The algorithm is designed by considering the dynamics of an accelerated systems of SDE discretized via Gauss-Seidel (GS) splitting. NAG-GS is theoretically analyzed on quadratic problems. Numerical experiments, including logistic regression, finetuning RoBerta on downstream tasks, and image classification via ResNet, are carried out to illustrate the usefulness of the proposed method.", "summary_of_the_review": "This paper puts more weight on the theoretical side, but the performance in typical machine learning tasks confines its practical potential.", "strength_and_weaknesses": "**Strength:**\n- (S1) Applying GS for discretizing SDEs is interesting and worth careful investigation.\n- (S2) Theoretical analysis, though limited to quadratic problems, is still insightful to understand the behavior of NAG-GS.\n\n**Weakness:**\n\nThe main weakness of this work lies in the numerical experiments. Confined by computational resources, the results do not support the efficiency of proposed method well. In particular, in many test cases, NAG-GS does not catch up with AdamW or SGDm. And sometimes the gap is not small. In particular, \n- (W1) it can be beneficial to run tests on quadratic functions to support the theoretical findings.\n- (W2) The experiment setups for logistic regression in Section 3.1 is not clear. In addition, it is more clear if the test accuracy can be summarized in a table for comparison. Moreover, when the learning rate is small, SGDm seems to have better performance over NAG-GS. For the relatively simple logistic regression problem, it does not bother too much for hyperparameter tuning. In other words, a better test accuracy is more important.\n- (W3) In section 3.2, the authors write “AdamW step size is $10^{-5}$, which is much smaller than $10^{-2}$ of NAG-GS”. In practice, it is common for Adam type algorithm to use a smaller learning rate than SGD. More importantly, on more than half downstream tasks, NAG-GS fails to catch up with AdamW.\n- (W4) It is unclear for the reason of switching acc@1 to acc@5 in section 3.3. Once again, NAG-GS does not offer competitive results.\n", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"correctness": "3: Some of the paper’s claims have minor issues. A few statements are not well-supported, or require small changes to be made correct.", "technical_novelty_and_significance": "3: The contributions are significant and somewhat new. Aspects of the contributions exist in prior work.", "empirical_novelty_and_significance": "2: The contributions are only marginally significant or novel."}, "extra_text": {}}]

Author rebuttal:

Decision: Reject