Survey of Optimization Algorithms in Modern Neural Networks

Abstract

The main goal of machine learning is the creation of self-learning algorithms in many areas of human activity. It allows a replacement of a person with artificial intelligence in seeking to expand production. The theory of artificial neural networks, which have already replaced humans in many problems, remains the most well-utilized branch of machine learning. Thus, one must select appropriate neural network architectures, data processing, and advanced applied mathematics tools. A common challenge for these networks is achieving the highest accuracy in a short time. This problem is solved by modifying networks and improving data pre-processing, where accuracy increases along with training time. Bt using optimization methods, one can improve the accuracy without increasing the time. In this review, we consider all existing optimization algorithms that meet in neural networks. We present modifications of optimization algorithms of the first, second, and information-geometric order, which are related to information geometry for Fisher–Rao and Bregman metrics. These optimizers have significantly influenced the development of neural networks through geometric and probabilistic tools. We present applications of all the given optimization algorithms, considering the types of neural networks. After that, we show ways to develop optimization algorithms in further research using modern neural networks. Fractional order, bilevel, and gradient-free optimizers can replace classical gradient-based optimizers. Such approaches are induced in graph, spiking, complex-valued, quantum, and wavelet neural networks. Besides pattern recognition, time series prediction, and object detection, there are many other applications in machine learning: quantum computations, partial differential, and integrodifferential equations, and stochastic processes.

1. Introduction

The question of increasing the accuracy of neural networks remains relevant. There exist many approaches to solve this problem, including the data augmentation [1], improving the mathematical model of neurons [2], adding the complement neural network [3], and so on. Indeed, all these approaches solve the problem of increasing the accuracy in neural networks, but these approaches are not generalized for all models. Thus, there exist universal methods that improve neural network training. One such is the optimization of the loss function. The main problem in neural networks is a descent to the global minimum. The first attempts to achieve the minimal value were realized using stochastic gradient descent (SGD) [4]. Such a first-order optimizer is still widely used in modern neural networks. Later, SGD was equipped with momentum and Nesterov condition (SGDM Nesterov) [5]. This modification increased the convergence rate with high accuracy in fewer iterations.

According to the structure of SGD, there has been proposals of AdaGrad in [6], RMSProp [7], Adadelta [8], and Adam [9], which allow attaining the minimum of functions that contain local extremes. In the case of the Rastrigin, Rosenbrouck, and Ackley functions from [10], these approaches do not reach the global minimum. This happens by achieving the minimum, considering only gradient directions. Exponential moving averages do not suffice for the minimization of such test functions. The same technique is used in modifications of Adam and SGD. Updating the exponential moving averages, Adam-type algorithms are built to improve the minimization for multi-extreme functions. The advantage of these approaches is that the initialization of biases can be counteracted, resulting in bias-corrected estimates. However, they can minimize the loss function but not always at the global minimum. Therefore, second-order optimization should be applied for higher accuracy attaining.

The main goal of second-order optimization algorithms [11] is to achieve the global minimum in a short time, because they are slower that first-order optimization algorithms. Second-order algorithms consider the convexity (curvature) of objective functions by the Hessian matrix. This approach is called the Newton optimization method [12]. Computations of an inverted Hessian matrix make second-order optimization more complex. In machine learning, Newton optimization is an ineffective tool due to the number of neurons, which can exceed one hundred. However, the approximation of inverted Hessian matrix allows the loss function to be minimized within the required time consumption. Such a technique is called the quasi-Newton optimization method [13]. Quasi-Newton optimization algorithms analyze the loss function from a functional point of view, which allows the local extremes to be avoided and rapid convergence in the global minimum. For the improvement of minimization, we observe the objective function from geometric and probabilistic points of view.

The first attempt to geometrically minimize smooth functions was applied in [14], where authors engaged the properties of Riemannian geometry. Such an optimization used the gradient flow for receiving the global minimum. Later, there was an idea to engage information geometry in [15], which is an intersection of Riemannian geometry and probability theory. Such an approach is divided into two branches, which utilize Fisher–Rao [16] or Bregman [17] metrics. The optimization with Fisher–Rao metrics (Fisher information matrix) is called a natural gradient descent. In the case of utilizing Bregman metrics, one receives the mirror descent. These optimization algorithms have demonstrated abilities in pattern recognition [18] and time series prediction [19]. They mainly find application in physics-informed neural networks [20], which are devoted to solving partial differential and integrodifferential equations.

The main goal of this article is to classify optimization algorithms, identify their properties, and provide minimization techniques for further research. After gradient-based optimization methods, we demonstrate alternative approaches to loss function minimization: gradient-free and bilevel optimization.

All provided optimization algorithms are presented in Figure 1, which indicates the contents covered in this paper.

In

Table 1. Summary of review papers on optimization algorithms in the theory of neural networks.

Reference	Summary of Work	Limitations
[21]	A survey demonstrating first-order optimization algorithms in convolutional neural networks.	This survey presents only gradient-based optimization algorithms. They are well known and do not give any novel ideas for neural networks of different architecture.
[25]	A study reviewing the Grasshopper optimization algorithm, which is used for various problems in machine learning, image processing, wireless networking, engineering design, and control systems.	Regardless of the good presentation of the local, evolutionary, and swarm optimization algorithms, the author does not present global optimizers, which can improve the work of gradient-free neural networks, or other models in machine learning.
[26]	A survey studying various models, datasets, and gradient-based optimization techniques in deep meta-learning.	There are only stochastic gradient and adaptive moment estimation approaches, customized under conditions of meta-learning systems.
[27]	This survey describes the information-geometric optimization approach from probabilistic and geometrical points of view. There are new type manifolds with described Riemannian metrics and connections.	This review does not give exact application domains and lacks the mirror descent, which is a duality of the natural gradient descent.
[23]	This survey describes fractional optimization. There are are implied mathematical formulations of fractional error backpropagation.	This review mentions only Riemann–Liouville, Caputo, and Grunwald–Letnikov derivatives.
[28]	This review studies fractional, gradient-free, and information-geometric optimization algorithms. The author shows new types of manifolds with implied Riemannian metrics and connections.	Such a review does not mention other types of fractional derivatives. It considers only particle swarm optimization algorithms from gradient-free approaches and briefly describes the natural gradient descent with Fisher matrix approximation.

, reviews have provided an incomplete classification of all existing optimization methods, such as first-order (SGD and Adam type), second-order (Newtonian and quasi-Newtonian), information-geometric, fractional, and gradient-free methods. In [21], authors present only SGD, NAG, AdaGrad, AdaDelta, RMSProp, Adam, NAdam, and AMSGrad. These optimizers are well known and do not provide new ideas for further research. In many other articles, fractional optimizers use the well-known Riemann–Liouville, Caputo, and Grunwald–Letnikov derivatives. In the review [22], the author classifies all existing fractional derivatives with the corresponding generalizations: Tarasov, Hadamard, Marchaux, Liouville–Weil, and others. In [23], the researchers introduced first-order, info-geometric, fractional, and gradient-free optimizers. The fractional approach is well described for activation functions in recurrent, convolutional, and complex-valued neural networks. Particle swarm optimization and information-geometric Levenberg–Marquardt algorithm are briefly presented. In this review, we have shown all existing types of optimizers. We have described the most promising approaches in machine learning, namely the natural gradient and mirror descent algorithms, which are well described in [24].

The main contribution of this review is the classification of optimization algorithms in the modern theory of neural networks. Figure 1 shows the evolution of first-order, second-order, and information-geometric optimization algorithms. AdaPNM and Adan algorithms use the extended versions of exponential moving averages, which give better convergence than Adam-type approaches. In addition to the basic Newton optimizer, we considered the conjugate gradient with all known update parameters. In addition to the BFGS and SR-1 algorithms, there are modern quasi-Newtonian Apollo and AdaHessian algorithms. In this review, we describe promising information-geometric optimizers. The new natural gradient descent is geometrically defined using the Fisher information matrices of Gaussian, polynomial, and Dirichlet distributions. The mirror descent optimizer is a duality of the natural gradient descent with an equivalent convergence rate. These methods consider not only gradient directions but also the curvature and duality of the loss function. For further research, we mentioned gradient-free optimization algorithms: local, global, population, sequential model, evolutionary, and genetic algorithms. Afterward, we present fractional optimizers provided with fractional derivatives containing kernel differences, logarithms, and exponentials.

Conclusions about provided optimizers are valid, considering the corresponding mathematical properties, experiments, exploitation, and comparative analysis in various neural network types. All parameters in optimizers are controllable at the request of a network. These parameters are learning rate, weight decay, dampening, Hessian approximation, Fisher information matrix, and other variables. These variables are applied according to existing problems of global minimum reaching, avoiding local extremes, and rectifying the updated value in case of exploding and vanishing gradient. This research is built on implying the convergence rate (regret bound), testing on Rosenbrock, Rastrigin, Ackley, Lévi, and other multimodal functions, and their exploitation on networks with any accessible datasets. Conclusions about provided optimizers can be generalized by extending exponential into positive–negative moving averages in Adam-type algorithms, and Fisher matrix modification replaces beta-distribution with Dirichlet distribution. The main feature of provided optimizers is their universality. They can be used in any gradient-based neural network with arbitrary available datasets on any CPU and GPU platform.

The rest of the paper is structured as follows. Section 2 is devoted to first-order optimization algorithms, which are described in Section 2.1–Section 2.3: SGD-type, Adam-type, and PNM-type optimization algorithms. Section 3 presents second-order optimization techniques, such as Newton and Quasi-newton algorithms. Section 4 consists of the representation of information-geometric methods, which are divided into Fisher–Rao and Bregman metrics. Part of this paper contains the most progressive approaches in the modern theory of optimization methods. Section 5 contains applications of all optimization algorithms introduced in Section 2, Section 3 and Section 4. In Section 6, we present conclusions on the provided optimization algorithms and proposals for their modification and further implementation in neural networks.

2. First-Order Optimization Algorithms

2.1. SGD-Type Algorithms

As we can see from Figure 1, the earliest first-order optimization algorithm is SGD [4], which can be described by the iterative formula

$${\theta }_{t+1}={\theta }_t-{\alpha }_t\nabla f\left({\theta }_t\right),$$

(1)

where ${\theta }_t$ is a weight, $f\left({\theta }_t\right)$ is a loss function with its gradient $\nabla f\left({\theta }_t\right)$, and ${\alpha }_t$ is a learning rate. Later, SGD was modified to SGDM with Nestreov condition [5] presented by

$${\theta }_{t+1}={\theta }_t-{\alpha }_t(\nabla f\left({\theta }_t\right)+\mu b_{t+1}),$$

(2)

where ${\theta }_0$ is an initial point, $b^{(k+1)}=\mu b^{\left(k\right)}+(1-\tau )(\nabla f\left({\theta }^{\left(k\right)}\right)+\lambda {\theta }^{\left(k\right)})$, where $\tau$ is a damping parameter and $\mu$ is a momentum. The algorithm of SGDM with Nesterov condition is constructed in [29]. This optimization method is used in approximation theory and machine learning. Various backpropagation methods ([30,31,32,33]) are based on calculating local partial derivatives, which rectify the value of weights of neural networks using (2). Such an approach can be modified to other more-effective versions, which more rapidly converge to a minimum. Therefore, because of the relatively low convergence rate of SGDM with the Nesterov condition and its step-size updating, the authors in [6] introduced AdaGrad.

The AdaGrad differs from SGDM with the Nesterov condition according to the adaptive step size. This advantage allows for an increase in learning rate and reduces time consumption. AdaGrad is described using the following formula

$${\theta }_{t+1}={\theta }_t-\frac{{\alpha }_t}{\sqrt{G_{t+1}+ϵ}}\nabla f\left({\theta }_t\right),$$

where the sum of gradients is

$$G_{t+1}=G_t+\nabla f\left({\theta }_{t-1}\right).$$

As was noted above, AdaGrad updates the step size ${\alpha }_t$ on the given information of all previous gradients observed along the way. However, its disadvantage is similar, like in SGD, where minimization is based on gradient directions and step-size regulation, which does not guarantee convergence in the global minimum neighborhood. Afterward, Adagrad was equipped with the gradient norm information, which improves the convergence rate. Such modification is called AdaGrad-Norm [6], presented as follows:

$${\theta }_{t+1}={\theta }_t-\frac{\alpha }{1+(t-1)\eta }\frac{g_t}{\sqrt{G_t}+ϵ}.$$

(3)

For the assessment of the optimizer, researchers use the regret bound $R\left(T\right)$, which measures the effectiveness of reinforcement learning systems. In other words, it verifies the validity of the optimizer convergence. The regret bound formula is described as follows:

$$R\left(T\right)=\sum _{t=1}^T[f_t\left({\theta }_t\right)-f_t\left({\theta }^{\ast }\right)],$$

where ${\theta }^{\ast }=argmi{n}_{\theta }{\sum }_{t=1}^T{f}_t\left(\theta \right)$ and $f_t\left({\theta }_t\right)$ is a loss function. First-order optimizers satisfy $O\left(\sqrt{T}\right)$ regret bound. Such a measure is comparable to convex online learning methods. We consider the following notations: $g_{1:t,i}=(g_{1,i},\cdots .,g_{T,i})\in {\mathbb{R}}^T$ is the gradient vector in the i-th dimension on $1,\cdots ,T$ iterations. The regret bound of AdaGrad is

$$R\left(T\right)\le \left(\frac{D^2}{2\alpha }+\alpha \right)\sum _{i=1}^d{\parallel g_{1:T,i}\parallel }_2,$$

where gradient $g_{t,\theta }$ satisfies the following conditions for any positive numbers $D,D_{\infty },G,G_{\infty }$:

$$\parallel g_{t,\theta }{\parallel }_2={\left({\int }_{\mathbb{R}}{g}_{t,\theta }^{2}d{\theta }_t\right)}^{\frac{1}{2}}\le G,{\parallel g_{t,\theta }\parallel }_{\infty }=\underset{{\theta }_t\in \mathbb{R}}{sup}{g}_{t,\theta }\le G_{\infty },$$

$$\parallel {\theta }_n-{\theta }_m{\parallel }_2\le D,{\parallel {\theta }_n-{\theta }_m\parallel }_{\infty }\le D_{\infty },$$

(4)

for all $m,n\in \{1,\cdots ,T\}$ and $\theta \in {\mathbb{R}}^d$. Algorithm estimation by the regret bound is used in the three cases where information about the sequence is not known in advance. The first case describes a flat domain, where the optimizer expects a large update but the $g_t$ gradient is very small. The second case depicts a large gradient domain, where the optimizer expects a large update supported by a large $g_t$. The third case represents a steep and narrow valley domain that mimics the minimum function. In such a domain, the optimizer expects a small update supported by a small $g_t$, which leads to descent into a local minimum. AdaGrad can solve only the third case because the previous gradient accumulations increase the update.

AdaGrad and AdaGrad-Norm adapt the step size more accurately than SGD and SGDM with Nesterov condition. It increases the learning rate to raise the probability of reaching the global minimum. The accumulation of previous gradients does not allow convergence in the neighborhood of the global minimum. This approach, like SGD, descends toward negative directions of gradients. Such a disadvantage caused researchers to come up with the idea of adapting the step size using mean moments, which implies the root mean square propagation (RMSProp) optimization algorithm.

The RMSprop optimizer [7] partially uses the same technique as SGDM and limits the oscillations in the upright direction, which increases the learning rate. It takes substantial step sizes in the horizontal direction, with convergence of $v_t$. Such an approach is

$$v_t=\gamma v_{t-1}+(1-\gamma )g_t^2,$$

$${\theta }_{t+1}={\theta }_t-\alpha \frac{g_t}{\sqrt{v_t}+ϵ},$$

(5)

where $\gamma \in (0,1)$ is a moment. The regret bound is

$$R\left(T\right)\le \left(\frac{D^2}{2\alpha }+\frac{\alpha (2-\gamma )}{\gamma }\right)\sqrt{T}\sum _{i=1}^d\left(\sqrt{v_{T,i}}+ϵ\right).$$

The parameter $v_t$ lets us solve the second possible case in the optimization process. The RMSProp is still an actual algorithm in many neural networks, such as SGDM with Nesterov condition, because it contains prerequisites of exponential moving averages. According to techniques in AdaGrad and RMSProp, AdaDelta was introduced in [8].

AdaDelta, based on techniques in RMSProp and AdaGrad, which separate dynamic learning rate per dimension, requires no manual setting of a learning rate and takes minimal computation over gradient descent. Additionally, it is insensitive to hyperparameters and robust to blow-up gradients, noise, and network choices. Unlike AdaGrad, such a method reduces aggressive, monotonically decreasing learning rates. AdaDelta limits the accumulated window of past gradients to a fixed size without accumulating all squares of past gradients. The running average $E{\left[g^2\right]}_t$ depends on the previous average and current gradient. Initial accumulation variables $E{\left[g^2\right]}_0$, $E{\left[\mathrm{\Delta }{\theta }^2\right]}_0$ are equal to 0. The AdaDelta algorithm is described as follows. Accumulate gradient:

$$E{\left[g^2\right]}_t=\rho E{\left[g^2\right]}_{t-1}+(1-\rho )g_t^2,$$

where $\rho$ is a decay rate parameter. Compute update:

$$\mathrm{\Delta }{\theta }_t=-\frac{RMS{\left[\mathrm{\Delta }\theta \right]}_{t-1}}{RMS{\left[g\right]}_t}{g}_t,$$

where

$$RMS{\left[g\right]}_t=\sqrt{E{\left[g^2\right]}_t+ϵ}.$$

Accumulate updates:

$$E{\left[\mathrm{\Delta }{\theta }^2\right]}_t=\rho E{\left[\mathrm{\Delta }{\theta }^2\right]}_{t-1}+(1-\rho )\mathrm{\Delta }{\theta }_t^2,$$

$${\theta }_{t+1}={\theta }_t+\mathrm{\Delta }{\theta }_t.$$

(6)

The SGD was modified in other variations, which, like AdaGrad, AdaDelta, and RMSProp, increase the test accuracy and accelerate the algorithm implementation. Such modifications improve recognition, prediction, generation, and decision-making processes, which develops the theory of neural networks. Modification (2)–(5) does not achieve higher accuracy compared with SGDM with Nesterov condition in deep convolutional neural networks, in [34]. Therefore, one needs algorithms that consider gradient directions and weight values. The valuable modification of SGDM is Nesterov accelerated gradient (NAG), presented in [35].

NAG is widely applied in deep networks and other supervised learning models. It often provides improvements over SGDM Nesterov. Strictly speaking, fast gradient methods only have provable corrections over regular gradient descent for the deterministic case with exact gradients. In the stochastic case, fast gradients partially mimic their exact gradient counterparts, resulting in some practical gain. Such an approach is described as

$$v_t=\mu v_{t-1}+\alpha f(\theta -\mu v_{t-1}),$$

$${\theta }_t={\theta }_{t-1}+v_t.$$

(7)

As said before, there is an efficient technique of utilizing the Lebesgue norm of the weights for improving the resulting test accuracy. Such an approach is called $L^2$ regularization, which is described in [36].

The $L^2$ regularization proposed in [37] acts similarly to the usual weight decay that is used in SGD. Indeed, both approaches evaluate weights closer to zero at the same rate. The $L^2$ regularized SGD (SGDW) is a modification of the usual SGD, which differs by recovering original weight decay using decoupling weight decay from the optimization steps for the loss function. For adaptive gradient algorithms, the main difference is the presence of adapted sums of gradients of the loss function. With $L^2$ regularization, both types of gradients are normalized by their summed magnitudes, and weights with large gradient magnitudes are regularized by a smaller relative amount than other weights. The introduced SGDW is described as follows:

$$g_t=\nabla f\left({\theta }_t\right)+\lambda {\theta }_{t-1},$$

where $\lambda$ is a weight decay parameter ($L_2$ regularization factor),

$$m_t={\beta }_1{m}_{t-1}+{\eta }_t\alpha g_t,$$

$${\theta }_t={\theta }_{t-1}-m_t-{\eta }_t\lambda {\theta }_{t-1},$$

(8)

where ${\eta }_t$ is a schedule multiplier. Unfortunately, the SGDW does not significantly improve the usual SGD, especially in deep learning. The convolutional neural networks confuse the $L^2$ regularization and make minimization too difficult for SGDW. However, the two techniques can improve minimization of the loss function using projection [38] and hyperparameter methods [39].

Before introducing the projection technique for SGD, it is necessary to recall batch normalization [40]. In practice, normalization techniques play a main role in modern deep learning. They allow weights to converge more rapidly with better generalization performances. The normalization-induced scale invariance among the weights gives advantages to SGD, such as the effective and automatic regularization of step size and stabilization of the training procedure. In practice, one can notice that including momentum in SGD-type optimizers reduces the step sizes for weight scaling. Such a phenomenon has not yet been studied and causes side effects in the minimization process. This is a critical issue because modern deep learning widely uses SGD- and Adam-type optimizers, which contain momentum and scale-invariant parameters. Therefore, SGD with projection (SGDP) was proposed, which removes the radial component, or the norm-increasing direction, at every iteration. Due to scale invariance, such a modification only alters the effective step without changing the update directions, thus satisfying the original convergence properties of the gradient descent optimizer.

The SGDP can be presented as the following iterative formulas:

$$\begin{array}{c}\hfill p_t=\mu p_{t-1}+\nabla f\left({\theta }_{t-1}\right),\\ \hfill {\theta }_t={\theta }_{t-1}-\alpha p_t.\end{array}$$

(9)

This approach, compared with SGDW, increases the convergence rate, but it is not sufficient for significantly advancing global loss function minimization. It can be seen applied minimizing the Rastrigin function in [41]. The hyperparameter method is called the quasi-hyperbolic momentum algorithm (QHM).

Momentum-based acceleration of SGD is used in deep learning. In [42], the author introduces QHM as a modification of SGDM Nesterom. This approach has an immediate $\nu$ discount factor that encapsulates SGD ($\nu =0$) and SGDM ($\nu =1$). A self-explanatory interpretation of QHM is the $\nu$-weighted average of the momentum update step and the simple SGD update step. The expressive power of QHM comes intuitively from separating the pulse buffer discount factor $\beta$ from the contribution of the current gradient to the $1-\nu \beta$ update rule. On the contrary, momentum is closely related to the discount factor $\beta$ and the contribution of the current gradient $1-\beta$.

Let us present QHM as the following iterative process:

$$g_t=\beta g_{t-1}+(1-\beta ){\widehat{g}}_t,$$

$${\theta }_t={\theta }_{t-1}-\alpha \left[(1-\nu )·{\widehat{g}}_t+\nu g_t\right].$$

(10)

Such a method can surmount the local minimums and has advantages over the abovedemonstrated algorithms. Like other optimizers, QHM does not suffice for achieving the highest accuracy in deep neural networks because it still takes into account the gradient directions, parameter amendments, and gradient normalization.

SGD-type algorithms are used in convolutional, recurrent, autoencoder, and graph neural networks. Their significant disadvantage is insufficient information about the properties of the loss function to reach the global minimum. The SGD is led by the gradient’s direction, which is not enough for the test accuracy to increase. In modern neural networks, more preferred approaches are Adam and its modifications (Adam-type algorithms). They contain exponential moving averages of gradient and squared gradient, which significantly improves the neural network training. These methods enforce the optimization by moment estimation, which gives more information about the global minimum. Moreover, such an improvement enhances pattern recognition, time series prediction, and object classification accuracy.

2.2. Adam-Type Algorithms

The adaptive moment estimate (Adam) [43] is a modified version of the gradient descent, which uses exponential moving averages of gradient and its square. Taking into account gradient directions and means of moments, loss function minimization with higher frequency converges in the global minimum. The iterative formula is presented as

$$m_t={\beta }_1{m}_{t-1}+(1-{\beta }_1)g_t,$$

$$v_t={\beta }_2{v}_{t-1}+(1-{\beta }_2)g_t^2,$$

$${\widehat{m}}_t=m_t/(1-{\beta }_1^t),{\widehat{v}}_t=v_t/(1-{\beta }_2^t),$$

$${\theta }_t={\theta }_{t-1}-\alpha {\widehat{m}}_t/(\sqrt{{\widehat{v}}_t}+ϵ).$$

(11)

Note that ${\beta }_1,{\beta }_2$ are called moments, and $m_t$ and $v_t$ are exponential moving averages of $g_t$ and $g_t^2$, respectively. The regret bound of Adam is

$$R\left(T\right)\le \frac{D^2\sqrt{T}}{2\alpha (1-{\beta }_1)}\sum _{i=1}^d\sqrt{{\widehat{v}}_{T,i}}+\frac{\alpha (1+{\beta }_1)G_{\infty }^3{G}^{-2}}{(1-{\beta }_1)\sqrt{1-{\beta }_2}{(1-\gamma )}^2}\sum _{i=1}^T{\parallel g_{1:T,i}\parallel }_2+\frac{D_{\infty }^2{G}_{\infty }\sqrt{1-{\beta }_2}}{2\alpha (1-{\beta }_1){(1-\lambda )}^2},$$

$\gamma =\frac{{\beta }_1^2}{\sqrt{{\beta }_2}}$ and $\lambda \in (0,1)$. The exponential moving averages let Adam solve the second and third possible cases in the optimization process, which happens in Rosenbrock function global minimization.

According to the adaptive moment estimate algorithm (11), there exist modifications, distinct by step size adaptation and manipulation with exponential moving averages. Like the SGDW, there is a modification in [44], called the Adam with $L^2$ regularization. In the case of the usual Adam, weights are not regularized as much as they would with decoupled weight decay since the gradient of the regularizer is scaled. The AdamW with the same exponential moving averages is described as the following iterative formula:

$${\theta }_t={\theta }_{t-1}-\eta +t\left[\alpha ·{\widehat{m}}_t/(\sqrt{{\widehat{v}}_t}+ϵ)+\lambda {\theta }_t\right],$$

(12)

where ${\eta }_t$ is a schedule multipliers and $\lambda$ is a $L^2$ regularization factor. This split weight reduction modification gives significantly better generalization performance than Adam optimizer with $L^2$ regularization in [45]. Another modification of Adam is projection Adam, which is called AdamP [46].

AdamP, such as SGDP in (8), is based on a sum projection of gradient and momentum vectors onto the tangent space of weights. Such an approach allows accelerating the effective step sizes for scale-invariant weights. This technique, together with its applications, is described in [47] and has the following representation:

$$g_t=\nabla f\left({\theta }_{t-1}\right)+\lambda {\theta }_{t-1},$$

$$m_t={\beta }_1{m}_{t-1}+(1-{\beta }_1)g_t,$$

$$v_t={\beta }_2{v}_{t-1}+(1-{\beta }_2)g_t^2,$$

$$p_t=\frac{m_t}{\sqrt{v_t}+ϵ}.$$

Then, according to the cosine condition, we obtain

$$ifcos\left(g_i\right)<\delta /\sqrt{dim\left(\theta \right)}:q_i=p_i-(\frac{\theta }{{\parallel \theta \parallel }_2}·p_i)\frac{\theta }{{\parallel \theta \parallel }_2}$$

$$else:q_i=p_i$$

Afterward, one obtains

$${\theta }_t={\theta }_{t-1}-{\eta }_t\left(\frac{\alpha m_t}{\sqrt{v_t}+ϵ}+\lambda {\theta }_{t-1}\right).$$

(13)

According to the notes in [48], this implies that momentum-based optimizers induce the excessive growth of scale-invariant weight norms, which causes premature decay of the effective optimization steps, leading to sub-optimal performances. The resulting AdamP, like SGDP in Section 2.1, successfully suppresses the weight norm growth and trains a model at an unobstructed speed. Another approach to raise the convergence rate is applying the quasi-hyperbolic momentum, similarly as for QHM in (9).

As for the QHM, there is the proposed QHAdam [49], in which Adam’s moment estimations are replaced with quasi-hyperbolic terms. This approach is described as

$$g_{t+1}={\beta }_1{g}_t+(1-{\beta }_1)\nabla f\left({\theta }_t\right),g_{t+1}^{\prime }={(1-{\beta }_1^{t+1})}^{-1}{g}_{t+1},$$

$$s_{t+1}={\beta }_2{s}_t+(1-{\beta }_2){(\nabla f\left({\theta }_t\right))}^2,s_{t+1}^{\prime }={(1-{\beta }_2^{t+1})}^{-1}{s}_{t+1},$$

$${\theta }_{t+1}={\theta }_t-\alpha \left[\frac{(1-{\nu }_1)\nabla f\left({\theta }_t\right)+{\nu }_1{g}_{t+1}^{\prime }}{\sqrt{(1-{\nu }_2{(\nabla f\left({\theta }_t\right))}^2)+{\nu }_2{s}_{t+1}^{\prime }}+ϵ}\right].$$

(14)

As was noted in the previous subsection, there is NAG, which accelerates the usual SGD in [50]. The same modification propagates to the Adam and transforms it to the NAdam [51]. This technique converges to the neighborhood of global minimum with higher frequency and a lower number of iterations, compared with Adam and its previous modifications, and it is constructed more simply than AdamW, AdamP, and QHAdam. For cases of Rastrigin and Rosenbrock functions, however, (11)–(13) do not converge in the global minimum. $L^2$ regularization, projection, and quasi-hyperbolic parameter techniques influence the step size, taking into account gradient directions. Such modifications can accelerate the convergence process but unnecessarily at the global minimum. However, there exist two techniques that make the minimization process “smoother”. Such approaches are called Nesterov-accelerated (NAdam) and rectified (RAdam) adaptive moment estimation.

Let us present the iterative formula of NAdam:

$${\theta }_{t+1}={\theta }_t-\frac{\eta }{\sqrt{{\widehat{v}}_t}+ϵ}\left({\beta }_1{\widehat{m}}_t+\frac{(1-{\beta }_1)g_t}{1-{\beta }_1^t}\right)$$

(15)

This method is the continuation of NAG due to the addition of exponential moving averages. Compared with Adam and its previous modifications, the NAdam increases the accuracy of converging in deep convolutional neural networks and makes the minimization process faster by “smoothing” the descent. Such a technique is ineffective in physics-informed neural networks because of the smoothing descent trajectory, which gives more deviations for partial differential equation solutions. This disadvantage is resolved by rectifying Adam (RAdam).

Proposed in [52], RAdam differs from other optimization methods by introducing a term to rectify the variance of the adaptive language modeling and learning rate. This modification proved its ability to receive higher test accuracy. Such optimization method has the following iterative formula:

$${\rho }_t={\rho }_{\infty }-2t{\beta }_2^t/(1-{\beta }_2^t)$$

$${\rho }_{\infty }=\frac{2}{1-{\beta }_2}-1$$

If the variance is tractable (${\rho }_t>4$) then the adaptive learning rate is computed as

$$l_t=\sqrt{(1-{\beta }_2^t)/v_t},$$

the variance rectification term is calculated as

$$r_t=\sqrt{\frac{({\rho }_t-4)({\rho }_t-2){\rho }_{\infty }}{({\rho }_{\infty }-4)({\rho }_{\infty }-2){\rho }_t}}$$

and we update parameters with adaptive momentum:

$${\theta }_t={\theta }_{t-1}-{\alpha }_t{r}_t{\widehat{m}}_t{l}_t.$$

(16)

If the variance is not tractable we update instead with

$${\theta }_t={\theta }_{t-1}-{\alpha }_t{\widehat{m}}_t.$$

(17)

Such a method overtakes the NAdam and other algorithms, especially in deep neural networks [53], such as AlexNet [54], ResNet [55], InceptionNet [56], GoogLeNet [57], and Res-Next [58]. However, RAdam adapts an overly complex learning rate and, like previous analogs, can not converge in the global minimum of the Rastrigin function. Moreover, there exist approaches that make the minimization process faster and more accurate. One such is the difference gradient approach, which is called DiffGrad [59].

The difference gradient approach is based on the moment estimate technique and calculates only a DiffGrad friction coefficient (DFC). The main distinction of DiffGrad is that such an approach is based on the change in short-term gradients and adjusts the learning rate according to the need for dynamic learning rate adjustment. Thus, DiffGrad makes the smaller parameter update in regions of low gradient change. The DiffGrad friction coefficient (DFC) controls the learning rate using information about the short-term gradient behavior. The DFC is represented by ${\xi }_t$ and defined as

$${\xi }_t=\frac{1}{1+exp-|\mathrm{\Delta }{g}_t|},$$

(18)

where $\mathrm{\Delta }{g}_t$ is the change between previous and current gradients, given as

$$\mathrm{\Delta }{g}_t=g_{t-1}-g_t.$$

In DiffGrad, the calculation of the first-order moment $m_t$ with offset correction and the second-order moment $v_t$ with offset correction is conducted in the same way as in Adam [60]. This optimization algorithm updates ${\theta }_{t+1}$, using the following update rule:

$${\theta }_{t+1}={\theta }_t-\frac{{\alpha }_t{\xi }_t{\widehat{m}}_t}{\sqrt{{\widehat{v}}_t}+ϵ}.$$

(19)

The regret bound is

$$R\left(T\right)\le \frac{D^2\sqrt{T}}{2\alpha (1-{\beta }_1)}\sum _{i=1}^d\frac{\sqrt{{\widehat{v}}_{T,i}}}{{\xi }_{T,i}}+\frac{\alpha (1+{\beta }_1)G_{\infty }^3{G}^{-2}}{(1-{\beta }_1)\sqrt{1-{\beta }_2}{(1-\gamma )}^2}\sum _{i=1}^T{\parallel g_{1:T,i}\parallel }_2+\frac{D_{\infty }^2{G}_{\infty }\sqrt{1-{\beta }_2}}{2\alpha (1-{\beta }_1){(1-\lambda )}^2}.$$

This optimizer can solve the second and third possible cases in the optimization process. For solving the first case, there is induced DFC, which lets us avoid the local minimums. Therefore, (19) can solve all three possible scenarios in loss function optimization, which Rastrigin and Rosebrock’s functions contain. The DiffGrad generates a high learning rate if the optimization is far from the optimum solution and a low learning rate if the optimization is near the optimum solution. Moreover, such a technique allows avoiding some local minimums, which occur in the minimization of Rastrigin and Rosenbrock functions [61]. This approach is suitable for deep convolutional neural networks due to moment estimation and analyzing previous and current gradients. Moreover, DiffGrad adjusts the learning rate very accurately, avoiding overshooting the global minimum and reducing oscillation around it. Such algorithms are tested on ResNet50 in [62], solving the pattern recognition problem on images from CIFAR10 and CIFAR100. According to the results of pattern recognition, it became clear that DiffGrad functions better than SGDM (2), AdaDelta (5), and Adam (11). However, this approach does not guarantee high accuracy in other neural networks, especially in quantum, spiking, complex-valued, and physics-informed neural networks. This disadvantage is explained by the lack of analysis of the curvature of the loss function during minimization. For that reason, the optimization algorithm in [63] was introduced, which is called Yogi.

The Yogi algorithm, like Adam, relies on gradient scaling by the square root of exponential moving averages of past squared gradients and controls the update in the effective learning rate, leading to even better performance with similar theoretical guarantees on convergence in [64]. This lets us solve the problem of convergence failure in simple convex optimization settings, which Adam-type algorithms, such as AdamW, AdamP, QHAdam, NAdam, and RAdam, cannot handle. The difference between $v_t$ and $v_{t-1}$ and its magnitude depend on $v_{t-1}$ and $g_t^2$. When $v_{t-1}$ is much larger than $g_t^2$, Yogi, like Adam, increases the effective learning rate, but such a procedure is more controllable. This approach has the following description:

$$m_t={\beta }_1{m}_{t-1}+(1-{\beta }_1)g_t,$$

$$v_t=v_{t-1}-(1-{\beta }_2)sign(v_{t-1}-g_t^2)g_t^2,$$

$${\widehat{m}}_t=m_t/(1-{\beta }_1^t),{\widehat{v}}_t=v_t/(1-{\beta }_2^t),$$

$${\theta }_t={\theta }_{t-1}-\gamma {\widehat{m}}_t/(\sqrt{{\widehat{v}}_t}+ϵ).$$

(20)

Compared with DiffGrad and other previous Adam-type algorithms, this method shows better results in deep convolutional neural networks. The regret bound of Yogi is the same as that of Adam. The performed ${\widehat{v}}_t$ solves three possible cases in loss function optimization. The authors of [65,66] defined new optimization methods for deep learning, such as AdaBelief and AdaBound.

The main feature of AdaBelief is adapting the learning rate according to the “belief” in the current gradient direction. There is a difference between AdaBelief and Adam in parameters $v_t$ and $s_t$, which are defined as exponential moving averages of $g_t^2$ and ${(g_t-m_t)}^2$, respectively. According to $s_t$ as the prediction of the gradient at the next time step, if the observed gradient greatly deviates from the prediction, then one distrusts the current observation and takes a small step; therefore, if the observed gradient is close to the prediction, one trusts it and take a large step. This allows us to achieve high test accuracy in convolutional neural networks in [66] on ImageNet and CIFAR10.

$$m_t={\beta }_1{m}_{t-1}+(1-{\beta }_1)g_t,$$

$$s_t={\beta }_2{s}_{t-1}+(1-{\beta }_2){(g_t-m_t)}^2+ϵ,$$

$$\widehat{m_t}=\frac{m_t}{1-{\beta }_1^t},\widehat{s_t}=\frac{s_t}{1-{\beta }_2^t},$$

$${\theta }_t={\theta }_{t-1}-\gamma {\widehat{m}}_i/(\sqrt{{\widehat{s}}_i}+ϵ).$$

(21)

The regret bound is

$$R\left(T\right)\le \frac{D^2\sqrt{T}}{2\alpha (1-{\beta }_1)}\sum _{i=1}^d\sqrt{{\widehat{s}}_{T,i}}+\frac{\alpha (1+{\beta }_1)(1+logT)}{2\sqrt{c}{(1-{\beta }_1)}^3}\sum _{i=1}^T{\parallel g_{1:T,i}\parallel }_2+\frac{D_{\infty }^2{G}_{\infty }{\beta }_1}{2\alpha (1-{\beta }_1){(1-\lambda )}^2},$$

where $s_{t,i}\ge c>0$, for all $t\in [1,T]$. Such an approach solves all three possible cases in loss function optimization. AdaBelief has a modified version with the fast gradient sign method (FGSM) presented in [67,68].

The AdaBound method allows restriction of the learning rate between the upper and lower continuous functions, called clips. Such a technique reduces the probability of vanishing and blow-up gradient. This approach is defined as

$$m_t={\beta }_1{m}_{t-1}+(1-{\beta }_1)g_t,$$

$$v_t={\beta }_2{v}_{t-1}+(1-{\beta }_2)g_t^{2}and{V}_t=diag\left(v_t\right),$$

$${\widehat{\eta }}_t=Clip(\alpha /\sqrt{V_t},{\eta }_l\left(t\right),{\eta }_u\left(t\right))and{\eta }_t={\widehat{\eta }}_t/\sqrt{t},$$

$${\theta }_{t+1}={\theta }_t-{\eta }_t{\widehat{m}}_t/(\sqrt{{\widehat{v}}_i}+ϵ).$$

(22)

The regret bound is

$$R\left(T\right)\le \frac{D^2\sqrt{T}}{2\alpha (1-{\beta }_1)}\sum _{i=1}^d\frac{1}{{\widehat{\eta }}_{T,i}}+\frac{{\beta }_1{D}_{\infty }^{{}^2}}{2(1-{\beta }_1){(1-\lambda )}^2{L}_{\infty }}+(2\sqrt{T}-1)\frac{G^2{R}_{\infty }}{(1-{\beta }_1)},$$

where $L_{\infty }={\eta }_l\left(1\right)$ and $R_{\infty }={\eta }_u\left(1\right)$. Such an approach solves all three possible cases in loss function optimization. This approach is more complex compared with DiffGrad, Yogi, and AdaBelief but can converge faster in the global minimum. There were experiments on minimization of Rastrigin and Rosenbrock functions [41], where AdaBound descent is the global minimum. However, such an approach is too complex for optimization, and there exists a more simple method, which is called AdamInject [69], which reduces time consumption while preserving the convergence rate.

AdamInject is one of the most recent approaches in first-order optimization algorithms, which, unlike the AdaBelief algorithm, modifies $m_t$, which is the exponential moving average of $g_t$, into $s_t$. Such a parameter is equipped with the difference between the previous parameters ${\theta }_{t-1}$ and ${\theta }_{t-2}$. AdamInject has the following description.

If $t=1$:

$$s_t={\beta }_1{s}_{t-1}+(1-{\beta }_1)g_t.$$

Else:

$$\mathrm{\Delta }\theta ={\theta }_{t-2}-{\theta }_{t-1}$$

$$s_t={\beta }_1{s}_{t-1}+(1-{\beta }_1)(g_t+\mathrm{\Delta }{g}_t^2)/k.$$

Afterward, one makes the usual calculations

$$v_t={\beta }_2{v}_{t-1}+(1-{\beta }_2)g_t^2,$$

$$\widehat{s_t}=\frac{s_t}{1-{\beta }_1^t},\widehat{v_t}=\frac{v_t}{1-{\beta }_2^t},$$

$${\theta }_t={\theta }_{t-1}-\gamma \widehat{s_t}/(\sqrt{\widehat{v_t}}+ϵ).$$

The regret bound is

$$R\left(T\right)\le \frac{D^2\sqrt{T}}{2\alpha (1-{\beta }_1)}\sum _{i=1}^d\frac{\sqrt{{\widehat{v}}_{T,i}}}{{\xi }_{T,i}}+\frac{\alpha (1+{\beta }_1)G_{\infty }^3{G}^{-2}}{(1-{\beta }_1)\sqrt{1-{\beta }_2}{(1-\gamma )}^2}\sum _{i=1}^T{\parallel g_{1:T,i}\parallel }_2$$

$$+\frac{D_{\infty }^2{G}_{\infty }\sqrt{1-{\beta }_2}}{2\alpha (1-{\beta }_1){(1-\lambda )}^2}+4\sum _{i=1}^d{G}_{\infty }^2{D}_{\infty }{\parallel g_{1:T,i}\parallel }_2^2.$$

Such an approach solves all three possible cases in loss function optimization. This algorithm is tested on the Rastrigin and Rosenbrock functions. Afterward, there AdamInject is trained on CIFAR10 VGG 16, ResNet, ResNext, SENet, and DenseNet and presented the best results compared with known analogs.

The introduced Adam-type optimization algorithms are used in deep convolutional neural networks, such as Res-Net, Res-Next, InceptionNet, GoogLeNet, and others. They also find application in recurrent and spiking neural networks due to their described advantages, which SGD-type algorithms do not have. In quantum, complex-valued, and quaternion-valued neural networks, however, such approaches suffer a loss in accuracy due to the usual SGDM with Nesterov condition. This problem caused researchers in [70] to devise extending the number of moments from 2 to 3. This approach is called positive–negative momentum (PNM) and Nesterov’s adaptation.

2.3. Positive–Negative Momentum

The development of SGD- and Adam-type optimization algorithms cannot be infinite. There have to be other approaches and techniques for extending the class of first-order methods. This issue has made researchers consider approaches that allow more than two exponential moving averages $m_t$ and $v_t$, because step-size regularization, according to moment estimation and modified Section 2.1 and Section 2.2, has its limits, which can be seen in convolutional neural networks, such as ResNet18, GoogLeNet, and DenseNet. In these cases, there has to be an additional exponential moving average that allows descending to the global minimum neighborhood.

In their paper [71], the authors introduce the conventional momentum method, also called heavy ball (HB) in [72]. Later, positive–negative momentum (PNM) optimization algorithms were proposed. In this approach, the main feature is positive–negative averaging, an exponential moving average analog. This averaging is

$$m_t=\sum _{k=0}^t{\beta }_3{\beta }_1^{t-k}{g}_k.$$

Inspired by this simple idea, there was a proposal to combine positive–negative averaging with the conventional momentum method in [71]. The positive–negative average is

$$m_t=(1+{\beta }_0)m_t^{\left(odd\right)}+{\beta }_0{m}_t^{even}=(1+{\beta }_0)\sum _{k=1,3,\cdots ,t}{\beta }_3{\beta }_0^{t-k}{g}_k+{\beta }_0\sum _{k=0,2,\cdots ,t}{\beta }_3{\beta }_0^{t-k}{g}_k.$$

The stochastic gradient descent equipped with the conventional momentum estimates, which adjust the learning rate and value of the gradient, is written in the following formula

$$m_t={\beta }_1^2{m}_{t-2}+(1-{\beta }_1^2)g_t,$$

$${\theta }_{t+1}={\theta }_t-\frac{\eta }{\sqrt{{(1+{\beta }_0)}^2+{\beta }_0^2}}[(1+{\beta }_0)m_t-{\beta }_0{m}_{t-1}].$$

(23)

This approach is an analog of the SGD-type algorithm, which differs by the presence of positive–negative average $m_t$ and step-size adaptation. The Adam-type analog of the PNM algorithm is also proposed, which is called AdaPNM, described as

$$m_t={\beta }_1^2{m}_{t-2}+(1-{\beta }_1^2)g_t,$$

$${\widehat{m}}_t=\frac{(1+{\beta }_0)m_t-{\beta }_0{m}_{t-1}}{1-{\beta }_1^t},$$

$$v_t={\beta }_2{v}_{t-1}+(1-{\beta }_2)g_t^2,$$

$$v_{max}=max(v_t,v_{max}),{\widehat{v}}_t=\frac{v_{max}}{1-{\beta }_2^t},$$

$${\theta }_{t+1}={\theta }_t-\frac{\eta }{\sqrt{{(1+{\beta }_0)}^2+{\beta }_0^2}}{\widehat{m}}_t.$$

(24)

The advantage of PNM and AdaPNM can be seen in [73], where these approaches are deep neural networks, such as ResNet, VGG, DenseNet, and GoogLeNet tested on image bases CIFAR10 and CIFAR100. Such techniques give higher test accuracy than advanced Adam-type optimization algorithms, such as Yogi and AdaBound. If one equips AdaPNM and PNM with techniques contained in known analogs, then it can improve the optimization process. Moreover, there exists the adaptive Nesterov momentum algorithm (Adan).

In [74], the authors propose the adaptive Nesterov momentum algorithm devoted to effectively accelerating neural network training. Adan revises vanilla Nesterov acceleration to develop a new method for estimating Nesterov momentum, reducing the amount of computation and memory cost of computing the gradient at the extrapolation point. Adan uses first- and second-order moments in adaptive gradient algorithms to speed up convergence. This approach has the following form:

$$m_t=(1-{\beta }_0)m_{t-1}+{\beta }_0{g}_t,$$

$$v_t=(1-{\beta }_1)v_{t-1}+{\beta }_1(g_t-g_{t-1}),$$

$$n_t=(1-{\beta }_2)n_{t-1}+{\beta }_2{[g_t+(1-{\beta }_1)(g_t-g_{t-1})]}^2,$$

$${\eta }_t=\eta /\left(\sqrt{n_t+ϵ}\right),$$

$${\theta }_{t+1}={(1+\lambda \eta )}^{-1}[{\theta }_t-{\eta }_t(m_t+(1-{\beta }_1))v_t].$$

(25)

This method is the generalization of SGD and Adam as the AdaPNM. There are demonstrated advantages of this approach over AdaBelief, which show the third highest test accuracy after the usual SGD. Regardless of the modifications of Adam-type algorithms, the usual SGD can give even better results, which claims to search other methods that approve their modifications.

First-order optimization methods are suitable for pattern recognition, time series prediction, and object classification. They do not consume much execution time and power, which makes them a realistic option in modern neural networks. However, first-order optimization algorithms, except PNM, AdaPNM, and Adan, cannot significantly increase the accuracy of neural networks with complex architecture, such as graph, complex-valued, and quantum neural networks. For solving differential equations, they cannot overtake the results of Adam because physics-informed neural networks contain automatic differentiation, which works after multilayer perceptron. Therefore, one needs second-order optimization algorithms, which can significantly improve loss function minimization.

3. Second-Order Optimization Algorithms

First-order optimization algorithms are ineffective for receiving the global minimum of objective smooth function. Such minimization approaches generally take into account the gradient directions in every iteration. The customization and performance of first-order optimization algorithms can only increase the accuracy and avoid some local minimums. Thus, there are provided second-order optimization algorithms. They reach the minimum taking into account not only the directions of gradients but the convexity (curvature) of the objective function, which is measured with the Hessian. Such an approach is called the Newton method.

3.1. Newton Algorithms

The main idea of the Newton optimization method [75] is based on gradient descent, containing the calculation of inverse Hessian of the smooth objective function. Such an approach increases the minimization accuracy of functions with multiple numbers of local minimums. Newton’s optimization algorithm is described as

$${\theta }_{t+1}={\theta }_t-Hess{\left({\theta }_t\right)}^{-1}\nabla f\left({\theta }_t\right),$$

(26)

where $t=0,\dots ,M>0$ and

$$Hess\left({\theta }_t\right)={\nabla }^{2}f\left({\theta }_t\right).$$

(27)

This iterative formula is received from second-order Taylor expansion.

In [76], the authors introduce the Newton minimum residual (Newton-MR) optimization method, which calculates the Hessian matrix by Moore–Penrose inverse [77] operator ${[·]}^{†}$ such as in

$${\theta }_{t+1}={\theta }_t+p_t={\theta }_t+{\left[{\nabla }^{2}f\left({\theta }_t\right)\right]}^{†}\nabla f\left({\theta }_t\right).$$

(28)

According to Newton and Newton-MR methods, OverSketch Newton Fast convex optimization was suggested. Such an approach is described in [78] with appropriate algorithms, which compute the update direction considering the case of strong convexity. However, there are other modifications of Newton’s optimization methods that differ in their simplicity in realization and implementation in neural networks.

There are Krylov methods, such as conjugate gradients (CG) [79], the minimal residual method (MINRES) [80], and the generalized minimal residual method (GMRES) [81]. GMRES applies to indefinite matrices, and MINRES applies to symmetric indefinite matrices. Besides GD, MINRES, and GMRES, there exists a generic stochastic inexact Newton-Krylov method, described in [82]. Such an approach can be implemented and applied in physics-informed neural networks because of suitable approximations for time-dependent equations with divergence operator $(\nabla ·)$. For the majority of neural networks, CG is the most preferred.

The most preferred Newton’s optimization method is the conjugate gradient. This approach includes an unrestricted class of optimizers with low memory requirements and high convergence rates. Such optimization algorithms are described as the following iterative formula:

$${\theta }_{t+1}={\theta }_t+{\eta }_t{d}_t,$$

$$d_{t+1}=-g_{t+1}+{\beta }_t{d}_t,d_0=-g_0.$$

The main part of this formula is CG update parameter ${\beta }_t$, which is presented in

Table 2. List of conjugate gradient update parameters.

CG Update Parameter	Authors	Year
${\beta }_t^{HS}=\frac{g_{t+1}^T{y}_t}{d_t^T{y}_t}$	Hestenes and Stiefel [83]	1952
${\beta }_t^{FR}=\frac{\parallel g_{t+1}{\parallel }_2^2}{\parallel g_t{\parallel }_2^2}$	Fletcher and Reeves [84]	1964
${\beta }_t^D=\frac{g_{t+1}^T{\nabla }^{2}f\left({\theta }_t\right)d_t}{d_t^T{\nabla }^{2}f\left({\theta }_t\right)d_t}$	Daniel [85]	1967
${\beta }_t^{PRP}=\frac{g_{t+1}^T{y}_t}{\parallel g_t{\parallel }_2^2}$	Polak and Ribière [86] and Polyak [87]	1969
${\beta }_t^{CD}=\frac{\parallel g_{t+1}{\parallel }_2^2}{-d_t^T{g}_t}$	Fletcher [88], CD stands for “Conjugate Descent”	1987
${\beta }_t^{LS}=\frac{g_{t+1}^T{y}_t}{-d_t^T{g}_t}$	Liu and Storey [89]	1991
${\beta }_t^{DY}=\frac{\parallel g_{t+1}{\parallel }_2^2}{d_t^T{g}_t}$	Dai and Yuan [90]	1999
${\beta }_t^N={\left(y_{t}2d_t\frac{\parallel y_{t+}{\parallel }_2^2}{d_t^T{g}_t}\right)}^T\frac{g_{t+1}}{d_t^T{y}_t}$	Hager and Zhang [91]	2005

.

Despite the increasing minimization accuracy, Newton’s method is slower than the first-order optimizers. This disadvantage significantly impacts time consumption, which slows the training of deep neural networks. For accelerating the minimization process, many approximations of the Hessian matrix have been developed. Second-order optimization algorithms, which contain the approximation of the Hessian matrix, are called quasi-Newton.

3.2. Quasi-Newton Algorithms

The class of quasi-Newton optimization algorithms shows approximately the same accuracy as Newton optimization algorithms, but the ability to converge faster makes them useful in machine learning. The first quasi-Newton optimization algorithms are BFGS [92] and L-BFGS [93].

The BFGS method and its limited memory version, at the t-th iteration, are presented as the following iterative formula:

$${\theta }_{t+1}={\theta }_t-{\alpha }_t{H}_t\nabla f\left({\theta }_t\right),$$

(29)

where ${\alpha }_t$ is the step length, $\nabla f\left({\theta }_t\right)$ is the gradient, and $H_t$ is the inverse BFGS Hessian approximation, which is updated at every iteration using the formula

$$H_{t+1}=V_k^T{H}_t{V}_t+{\rho }_t{s}_t{s}_t^T,$$

(30)

$${\rho }_t=\frac{1}{y_t^T{s}_t},V_t=I-{\rho }_t{y}_t{s}_t^T,$$

(31)

where the curvature pairs $(s_t,y_t)$ are defined as

$$s_t={\theta }_t-{\theta }_{t-1},y_t=\nabla f\left({\theta }_t\right)-\nabla f\left({\theta }_{t-1}\right).$$

(32)

The curvature pairs $(s_t,y_t)$ are constructed sequentially at every iteration, and inverse Hessian approximation at the t-th iteration $H_t$ depends on iterate and gradient information from past iterations.

The inverse BFGS Hessian approximations have to satisfy secant and curvature conditions:

$$H_{t+1}{y}_t=s_t,s_t^T{y}_t>0.$$

(33)

As a result, as long as the initial inverse Hessian approximation is positive definite, then all subsequent inverse BFGS Hessian approximations are also positive definite. Note the inverse Hessian approximation $H_{k+1}$ distinct from the approximation $H_k$ by a rank-2 matrix.

In the limited memory version, the matrix $H_k$ is defined at each iteration as the result of applying m BFGS updates to a multiple of the identity matrix using the set of m most recent curvature pairs $\{s_t,y_t\}$ kept in storage. Thus, one does not need to store the approximations of the four dense inverse Hessian matrices; rather, one can store two $m\times d$ matrices and compute the matrix-vector product via two-loop recursion [94]. After the step has been computed, the oldest pair $(s_t,y_t)$ is discarded, and the new curvature pair is stored.

Analogically, the symmetric rank-one (SR-1) [95] formula and its limited memory version were proposed in [96]. At the t-th iteration, the SR1 method computes a new iterate by the formula

$${\theta }_{k+1}={\theta }_k+p_k,$$

(34)

where $p_k$ is the minimizer of the following equation

$$\underset{p}{min}{m}_k\left(p\right)=f\left({\theta }_k\right)+\nabla f{\left({\theta }_k\right)}^{T}p+\frac{1}{2}{p}^T{B}_{k}p,$$

(35)

$${\parallel p\parallel }_2\le {\mathrm{\Delta }}_k,$$

where ${\mathrm{\Delta }}_k$ is the trust region and $B_k$ is the SR1 Hessian approximation computed as

$$B_{k+1}=B_k+\frac{(y_k-B_k{s}_k){(y_k-B_k{s}_k)}^T}{{(y_k-B_k{s}_k)}^T{s}_k}.$$

(36)

Similar to L-BFGS, the limited memory version of SR1 defines the matrix $B_k$ as the result of applying m SR1 updates to the identity matrix using a set of m correction pairs $\{s_i,y_i\}$ kept in storage.

The provided BFGS and SR-1 with low-memory versions preserve the convergence rate and reduce the computational complexity. However, these methods remain ineffective in deep neural networks, which make calculations of matrices of high dimensions. Therefore, one needs to come up with simplified versions of quasi-Newton algorithms. One of the most used quasi-Newton optimization methods in neural networks is Apollo [97]. Apollo is a non-convex stochastic optimization algorithm. It includes the loss function curvature by approximating the Hessian with a diagonal matrix. This algorithm updates and stores the diagonal approximation of the Hessian as efficiently as Adam-type optimizers with linear complexity in both time and memory. Dealing with non-convexity, the Hessian is replaced by its positive definite rectified absolute value. Experiments with deep neural networks, which recognize images from CIFAR10 and ImageNet, show that Apollo achieves a global minimum with time and memory consumption similar to first-order optimizers. Such an approach is described as

$$g_{t+1}=\nabla f\left({\theta }_t\right),$$

$$m_{t+1}=\frac{\beta (1-{\beta }^t)}{1-{\beta }^{t+1}}{m}_t+\frac{1-\beta }{1-{\beta }^{t+1}}{g}_{t+1},$$

$$\alpha =\frac{d_t^T(m_{t+1}-m_t)+d_t^T{B}_t{d}_t}{(\parallel d_t{{\parallel }_4+ϵ)}^4},$$

$$B_{t+1}=B_t-\alpha ·diag\left(d_t^2\right),$$

$$D_{t+1}=rectify(B_{t+1},1),$$

$$d_{t+1}=D_{t+1}^{-1}{m}_{t+1},$$

$${\theta }_{t+1}={\theta }_t-{\eta }_{t+1}{d}_{t+1},$$

(37)

The regret bound is

$$R\left(T\right)\le \frac{D^2{\beta }_2}{2\eta (1-\beta ){(1-{\lambda }^2)}^2}+\frac{3\eta G^2}{2(1-\beta )}(2\sqrt{T}-1).$$

Such an approach solves all three possible cases in loss function optimization. Apollo is a method devoted to accelerating the minimization process. There is, however, a problem of accuracy loss in achieving lower time consumption. Because of this disadvantage, the authors in [98] proposed another quasi-Newton optimization method in neural networks, which is called AdaHessian.

AdaHessian is a second-order stochastic optimization algorithm that incorporates the loss function curvature with adaptive estimates of the Hessian. Second-order optimizers are the most advanced approaches with convergence properties superior to those of the above-presented first-order methods. The main disadvantage of traditional second-order methods is their heavier computation for each iteration and low accuracy. However, AdaHessian includes a fast Hutchinson method for curvature matrix approximation with low computational cost, root mean square exponential moving average for smoothing variations of the Hessian diagonal across different iterations, and block diagonal averaging for reducing the variance of Hessian diagonal elements. This approximation of the Hessian matrix makes the learning process faster while preserving the convergence rate.

Before presenting the AdaHessian algorithm, the gradient $g_t=\nabla f\left({\theta }_t\right)$ and Hessian $H_t={\nabla }^{2}f\left({\theta }_t\right)$ need to be noted. Let us present the following matrix diagonalization of the Hessian matrix as

$$D_t=diag\left(H\right)=\mathbb{E}[{\theta }_t\odot \left(H_t{\theta }_t\right)],$$

where ⊙ is a componentwise multiplication of vector. Perform a simple spatial averaging on the Hessian diagonal as follows:

$$D^{\left(s\right)}[ib+j]=\frac{1}{b}\sum _{k=1}^{b}D[ib+k],$$

for $1\le j\le b$, $0\le i\le \frac{d}{b}-1$. The Hessian diagonal with momentum is

$${\overline{D}}_t=\sqrt{\frac{(1-{\beta }_2){\sum }_{i=1}^t{\beta }_2^{t-i}{D}_i^{\left(s\right)}{D}_i^{\left(s\right)}}{1-{\beta }_2^t}},$$

$$m_t=\frac{(1-{\beta }_1){\sum }_{i=1}^t{\beta }_1^{t-i}{g}_i}{1-{\beta }_1^t},$$

$$v_t={\left({\overline{D}}_t\right)}^k,$$

$${\theta }_t={\theta }_{t-1}-\eta m_t/v_t.$$

(38)

Second-order optimization algorithms have a higher convergence rate. Unfortunately, they are not suitable in deep convolutional, recurrent, and spiking neural networks because of their high complexity. Therefore, one needs methods that satisfy the requirements for high convergence rate and low time consumption. In 2012, there was a proposition to develop second-order optimization algorithms using the smooth manifolds in [99]. In [100], the authors came up with the idea of using probability distribution manifolds instead of smooth. Such a technique appears at the intersection of geometry, probability, statistics, and optimization.

4. Information-Geometric Optimization Methods

Second-order optimization methods converge faster to the neighborhood of the global minimum compared with first-order methods. Unfortunately, their time consumption is not sufficiently fast for deep neural networks. Therefore, it is necessary to make quasi-Newton optimization algorithms faster. The most recent way involves applying information geometry. In this section, we present a concise and modern view of the basic structures of information geometry and report some applications in machine learning (statistical mixture clustering).

Initially, the idea of using geometric means in optimization problems appeared in 2010 in [101]. The Hessian matrix, presented in second-order optimization algorithms, does not contain full information about the surface of the loss function. Moreover, it makes the descent too complex from a computational point of view. Taking into account the geometric structure of the surface and corresponding spaces (Euclidean, hyperbolic, parabolic), however, allows us to define the shortest way for descent and can reduce unnecessary computations. In [102], for defining the shortest way, the authors applied gradient flow, which improved the quality of searching global minimum. This method applied smooth manifolds and a generalized definition of a gradient. Later, however, researchers came up with the idea to use probability distribution manifolds, in which their advantage over smooth manifolds was demonstrated in the optimization problem. The intersection of probability theory and statistics with Riemannian geometry produces information geometry.

Analogically to information theory, considering the communication of messages over noisy transmission channels, one defines information sciences as the fields of studying the connection between data and families of models, i.e., information sciences create methods to transform information from data to models. Such transformation is made by probabilistic, statistical, and geometric means, which allows us to include it in machine learning.

There is another definition of information geometry, given by F. Nilsson in [103], as decision-making geometry. It involves model fitting (inference), which can be interpreted as a decision problem, namely the decision of which model parameter to choose from a family of parametric models. In [104,105], Abraham Wald proposed this scope, considering all statistical problems as statistical decision problems. Dissimilarities (also commonly referred to as distances among others) play a crucial role not only in measuring the degree of fit of model data (probability in statistics, classifier loss functions in machine learning, objective functions in mathematical programming or operations research, etc.) but also when measuring the discrepancy (or deviation) between models.

In this section, we show two distinct optimization algorithms based on information-geometric approaches—natural gradient descent [106] and mirror descent [107]. These methods differ from other optimization algorithms by the ability to measure distances in non-Euclidean domains using Kullback–Leibler and Bregman divergences, respectively.

4.1. Natural Gradient Descent

Initially, gradient descent with gradient flow was proposed in [108]. This approach is a generalization of the second-order optimization method on Riemannian manifolds, which can be Euclidean or non-Euclidean. Let $({\mathcal{M}}^n,g)$ be a Riemannian manifold, where ${\mathcal{M}}^n$ is a topological space, which can be expressed in the local coordinate system of an atlas $\mathcal{A}={\left\{(U_i,x_i)\right\}}_i$ of charts $(U_i,x_i)$, and for tangent bundle $T{\mathcal{M}}^n$ Riemannian metric $g:T{\mathcal{M}}^n\otimes T{\mathcal{M}}^n\to \mathbb{R}$.

The dynamics of the Riemannian gradient flow $\theta \left(t\right)$ for the optimization problem ${min}_{\theta \in {\mathcal{M}}^n}f\left(\theta \right)$ is obtained by searching for an infinitesimal change in $\theta \left(t\right)$, which would lead to better improvement in the objective value while controlling the length of the change in terms of the geometry of the manifold, i.e.,

$$\theta (t+dt)=argmi{n}_{\theta }f\left(\theta \right)dt+\frac{1}{2}d{\theta }^{T}g(\theta ,\theta +d\theta )d\theta .$$

(39)

For $dt$, utilizing $d\theta \left(t\right)=\theta (t+dt)-\theta \left(t\right)$, one can substitute $f\left(\theta \right(t\left)\right)+\langle d\theta ,\nabla f\left(\theta \right(t\left)\right)\rangle$ instead of $f\left(\theta \right)$ and receive

$$d\theta \left(t\right)=argmi{n}_{d\theta }\langle d\theta ,\nabla f\left(\theta \left(t\right)\right)\rangle dt+\frac{1}{2}d{\theta }^{T}g(\theta ,\theta +d\theta )d\theta .$$

Solving for $d\theta$, one obtains

$$\frac{d\theta \left(t\right)}{dt}=-g{(\theta ,\theta +d\theta )}^{-1}\nabla f\left(\theta \left(t\right)\right).$$

(40)

For better understanding, the examples provided of the Riemannian metric is the standard Euclidean manifold that has the corresponding metric $g=I$, which reduces the gradient flow to the usual gradient descent. There is, however, another important example—probability distribution manifold with K-L divergence metric [109,110], which is the beginning of quantum neural networks.

There is a noted advantage of using probability distribution manifolds in [111], which can be equipped not necessarily with KL-divergence but with Bergman, Jensen, and others. Extensions of Rimannian manifolds with Levi-Civita connection $({\mathcal{M}}^n,g,{\nabla }_{LC})$ to conjugate connection manifolds $({\mathcal{M}}^n,g,\nabla ,{\nabla }^{\ast })$ have been presented, where ${\nabla }_{LC}=\frac{\nabla +{\nabla }^{\ast }}{2}$. Conjugate connection manifolds are the particular case of divergence manifolds, denoted as $({\mathcal{M}}^n,D^g,{\nabla }^D),{\nabla }^{\ast D}$. In such manifolds, D can be Kullback–Leibler or Bregman divergences. There are two ways to imply natural gradient descent using direct K-L divergence, such as in [100], and a more fundamental way, presented in [103]. Let $f\in C^{\infty }\left({\mathcal{M}}^n\right)$ be a smooth loss function and ${exp}_{\theta }:T_{\theta }{\mathcal{M}}^n\to {\mathcal{M}}^n$ be the Riemannian exponential map for updating the sequence of points ${\theta }_t$ on the manifold as follows:

$${\theta }_{t+1}={exp}_{{\theta }_t}(-{\alpha }_t{\nabla }_{M}f\left({\theta }_t\right)),$$

(41)

where

$${\nabla }_{M}f\left(\theta \right)={\nabla }_v\left(f\left({exp}_{\theta }\left(v\right)\right)\right){|}_{v=0}=\underset{h\to 0}{lim}\frac{f(\theta +hv)-f\left(\theta \right)}{h}.$$

(42)

Returning to the formula of gradient flow, we obtain

$${\theta }_{t+1}=R_{{\theta }_t}(-{\alpha }_t{\nabla }_{\theta }f\left({\theta }_t\right)).$$

Utilizing the retraction $R_{\theta }\left(v\right)=\theta +v$, which is the first-order Taylor approximation of the exponential map, follows the natural gradient descent:

$${\theta }_{t+1}={\theta }_t-{\alpha }_t{g}^{-1}(\theta ,\theta +d\theta ){\nabla }_{\theta }f\left({\theta }_t\right).$$

(43)

Here $g(\theta ,\theta +d\theta )=F\left({\theta }_t\right)$ is called Fisher information matrix, and the natural gradient is defined as follows:

$${\nabla }^{NG}f\left(\theta \right)=F^{-1}\left(\theta \right){\nabla }_{\theta }f\left(\theta \right).$$

(44)

Afterward, we come to the natural gradient descent:

$${\theta }_{t+1}={\theta }_t-{\alpha }_t{\nabla }_{\theta }^{NG}f\left({\theta }_t\right).$$

(45)

The corresponding regret bound is

$$R\left(T\right)\le \frac{D^2}{2\alpha }Tr\left(F\right)+\frac{\alpha }{2}\sum _{t=1}^T{\parallel g_t\parallel }_F^2,$$

where $Tr\left(F\right)=F_{1,1}+\cdots +F_{n,n}$ and $\parallel g_t{\parallel }_F=\sqrt{{\sum }_{i=1}^n{\sum }_{j=1}^n\left|g_{t,ij}\right|}$ for the Riemannian metric $g_t$. Such an approach solves all three possible cases in loss function optimization. Natural gradient descent differs from the first- and second-order optimization algorithms presented in Section 2 and Section 3, respectively, by the ability to converge in global minimum for time consumption, suitable for deep learning. As was said, such an approach creates a new branch in the theory of artificial intelligence—quantum machine learning. The application of neural networks for training quantum circuits necessitates a new structure for the mathematical model of neurons appropriate for quantum computations. Afterward, vanilla natural gradient descent was replaced with quantum in [112]. Such networks are part of developing quantum information theory and quantum computers, which makes the training significantly faster. In actual research, however, vanilla natural gradient descent with Dirichlet distributions has already been tested on Rastrigin and Rosebrock functions in [113] and exploited in convolutional, recurrent neural networks in [114]. The important thing in this method is selecting the appropriate probability distribution, which can increase the convergence rate and even accelerate the learning process. The probability distributions are presented in

Table 3. List of Fisher information matrices for different probability distributions.

Probability Density Function	Fisher Information Matrix	Probability Distribution
$p(x;\mu ,\sigma )=\frac{2{\pi }^{-\frac{n}{2}}}{\sqrt{{\sigma }_1^2\cdots {\sigma }_n^2}}{e}^{-\frac{{(x-\mu )}^{T}diag{({\sigma }_1^2,\cdots ,{\sigma }_n^2)}^{-1}(x-\mu )}{2}}$	$F_G=\left(\begin{array}{ccccc}\frac{1}{{\sigma }_1^2}& 0& \cdots & 0& 0\\ 0& \frac{2}{{\sigma }_1^2}& \cdots & 0& 0\\ \cdots & \cdots & \cdots & \cdots & \cdots \\ 0& 0& \cdots & \frac{1}{{\sigma }_n^2}& 0\\ 0& 0& \cdots & 0& \frac{2}{{\sigma }_n^2}\end{array}\right)$	Gauss [115]
$p(x;t)=\frac{n!}{x_1!\cdots x_n!}{t}_1^{x_1}\cdots t_n^{x_n}$	$F_M=\left(\begin{array}{cccc}\frac{t_1+t_n}{t_1{t}_n}& \frac{1}{t_n}& \cdots & \frac{1}{t_n}\\ \frac{1}{t_n}& \frac{t_2+t_n}{t_2{t}_n}& \cdots & \frac{1}{t_n}\\ \cdots & \cdots & \cdots \\ \frac{1}{t_n}& \frac{1}{t_n}& \cdots & \frac{t_1+t_n}{t_{n-1}{t}_n}\end{array}\right)$	Multinomial [116]
	$F_D^{ii}={\psi }^{\prime }\left({\alpha }_i\right)-{\psi }^{\prime }\left({\sum }_i{\alpha }_i\right),$
$p(x;\alpha )=\frac{1}{B\left(\alpha \right)}{\prod }_{i=1}^n{x}_i^{{\alpha }_i-1},B\left(\alpha \right)=\frac{{\prod }_i\Gamma \left({\alpha }_i\right)}{\Gamma \left({\sum }_i{\alpha }_i\right)}$	$F_D^{ij}=-{\psi }^{\prime }\left({\sum }_i{\alpha }_i\right)$	Dirichlet [113,117]
	$i\ne j,i=1,..,n$
$p(x;\alpha ,\beta )={\prod }_{i=1}^n\frac{x_i^{{\alpha }_i-1}}{B({\alpha }_i,{\beta }_i)}{\left(1-{\sum }_{j=1}^i{x}_j\right)}^{{\gamma }_i}$	$F_{GD}=diag\left(F_D({\alpha }_1,{\beta }_1),\cdots ,F_D({\alpha }_1,{\beta }_1)\right)$
${\gamma }_i={\beta }_i-{\alpha }_{i+1}-{\beta }_{i+1}$ for $i=1,\cdots ,n-1$		Generalized Dirichlet [113]
and ${\gamma }_n={\beta }_{n-1}$	$\mathcal{O}$-zero matrix

, which can improve the quality of minimization of the loss functions.

As can be seen, this is a second-order optimization algorithm. However, by selecting appropriate probability distribution, such as Gauss and Dirichlet, we can reduce the variable $\theta$ in the Fisher information matrix, which makes it possible to avoid its calculation in every iteration. Such approach is realized in [113,114,115,116]. The natural gradient descent can replace second-order optimization algorithms due to convergence rate and time consumption. However, there is another approach that uses Bregman metrics, called mirror descent.

4.2. Mirror Descent

The Bregman metric is an alternative to Fisher–Rao metrics based on dual Hessian manifolds. This approach considers not only the gradient directions and curvature of the loss function but also the duality. The dual spaces are defined as the set of maps from the topological spaces to their underlying field, which is usually presented as $\mathbb{R}$. This feature has impacted the name of the method. This technique is a reason for continuation of second-order optimization algorithms toward the information geometry. Moreover, it can be used in physics-informed neural networks, where dual averaging procedures can increase the accuracy of final solutions.

Recall that the gradient descent can be extended by proximity function $\Phi (·,·)$ as follows:

$${\theta }_{t+1}=argmi{n}_{\theta }\left\{\langle \theta ,\nabla f\left({\theta }_t\right)\rangle +\frac{1}{{\alpha }_t}\Phi (\theta ,{\theta }_t)\right\}.$$

(46)

If $\Phi (\theta ,{\theta }_t)=\frac{1}{2}{\parallel \theta -{\theta }_t\parallel }_2^2$, then we obtain the usual gradient descent. For Bregman divergence, the proximity function is $\Phi (\theta ,{\theta }_t)=D_{\psi }(\theta ,{\theta }_t)$, and we get mirror descent

$${\theta }_{t+1}=argmi{n}_{\theta }\left\{\langle \theta ,\nabla f\left({\theta }_t\right)\rangle +\frac{1}{{\alpha }_t}{D}_{\psi }(\theta ,{\theta }_t)\right\},$$

(47)

where

$$D_{\psi }(\theta ,{\theta }_t)=\psi \left(\theta \right)-\psi \left({\theta }_t\right)-\langle \theta -{\theta }_t,\nabla \psi \left({\theta }_t\right)\rangle .$$

(48)

According to the above divergence manifolds, it is possible to imply equivalence between natural gradient descent and mirror descent. In addition, this fact was presented in [103], but without geometric tools. Mirror descent is formulated as a duality of natural gradient descent. Natural gradient descent on the dual Hessian manifold $(M,g^{\ast }={\nabla }^{2}F\left(\eta \right))$ is equivalent to Bregman mirror descent on the Hessian manifold $(M,g={\nabla }^{2}F\left(\theta \right))$, where F is the Bregman generator, $\eta =\nabla F\left(\theta \right)$, and $\theta =\nabla F^{\ast }\left(\eta \right)$.

The mirror descent gives the following natural gradient update rule:

$$N{G}^{\ast }:{\xi }_{t+1}={\xi }_t-\eta {\left(g_{\xi }^{\ast }\right)}^{-1}\left({\xi }_t\right){\nabla }_{\xi }{f}_{\theta }\left(\theta \left({\xi }_t\right)\right),$$

(49)

$$N{G}^{\ast }:{\xi }_{t+1}={\xi }_t-\eta {\left(g_{\xi }^{\ast }\right)}^{-1}\left({\xi }_t\right){\nabla }_{\xi }{f}_{\xi }\left({\xi }_t\right)$$

where $g_{\xi }^{\ast }\left(\xi \right)={\nabla }^2{F}^{\ast }\left(\xi \right)={\left({\nabla }_{\theta }^{2}F\left(\theta \right)\right)}^{-1}$ and $\theta \left(\xi \right)=\nabla F^{\ast }$.

The method is called mirror descent because of performing the gradient step in the dual space, which plays the role of “mirror”. This means that mirror descent seeks the global minimum according to the duality of the probability distribution manifold, which is equivalent to the Fisher information matrix for natural gradient descent.

Besides the usual mirror descent, there is the stochastic mirror descent (SMD) proposed by Nemirovski and Yudin in [118]. This method presented high accuracy on ResNet18 in recognizing images from Cifar10. For a strictly convex differentiable function $\psi (·)$, called the potential (proximity) function, stochastic mirror descent is presented by the following iterative formula

$$\nabla \psi \left({\theta }_{t+1}\right)=\nabla \psi \left({\theta }_t\right)-\eta \nabla f\left({\theta }_t\right),$$

(50)

which is equivalent to the following expression

$${\theta }_{t+1}=arg\underset{\theta }{min}\left\{\eta {\theta }_t^T\nabla f\left({\theta }_t\right)+D_{\psi }(\theta ,{\theta }_t)\right\},$$

(51)

where Bregman divergence can be presented as

$$D_{\psi }(\theta ,{\theta }_t)=\psi \left(\theta \right)-\psi \left({\theta }_t\right)-\nabla \psi {\left({\theta }_t\right)}^T(\theta -{\theta }_t),$$

(52)

which is the Bregman divergence to the potential function $\psi$. Note that $D_{\psi }$ is nonnegative, convex in its first argument, and that due to strict convexity, $D_{\psi }(\theta ,{\theta }^{\prime })=0$ if and only if $\theta ={\theta }^{\prime }$. The regret bound is

$$R\left(T\right)\le \frac{D^2}{2\eta }+\frac{\eta }{2}\sum _{t=1}^T{\parallel g_t\parallel }_{\ast }^2,$$

where $\parallel g_t{\parallel }_{2,\ast }$ is a dual norm of $\parallel g_t{\parallel }_2$ for the gradient $g_t$. Such an approach solves all three possible cases in loss function optimization.

Different choices of the potential function $\psi$ yield different optimization algorithms, which will potentially have different implicit biases. A few examples follow in

Table 4. Potential function with corresponding Bregman divergences.

Potential Function $\mathit{\psi }\left(\mathit{\theta }\right)$	Bregman Divergence ${\mathit{D}}_{\mathit{\psi }}(\mathit{\theta },{\mathit{\theta }}^{\prime })$	Algorithm
$\frac{1}{2}{\parallel \theta \parallel }_2^2$	$\frac{1}{2}{∥\theta -{\theta }^{\prime }∥}_2^2$	Gradient Descent
${\sum }_j{\theta }_{j}log{\theta }_j$	${\sum }_j{\theta }_{j}log\frac{{\theta }_j}{{\theta }_j^{\prime }}-{\sum }_j{\theta }_j+{\sum }_j{\theta }_j^{\prime }$	Exponentiated Gradient Descent

.

Potentially, this method can improve the loss function minimization in the convolutional, graph, and recurrent neural networks of huge architectures. Moreover, mirror descent can be equipped with adaptive moment estimation and other modifications proper for first-order optimization methods.

5. Application of Optimization Methods in Modern Neural Networks

The introduced optimization methods find application in various artificial neural networks. All presented optimization algorithms show good results in some determined problems, which can be solved using machine learning.

First-order optimization methods are reliable in problems of pattern recognition, time series prediction, voice detection, and text analysis. In such cases, one needs deep convolutional and recurrent neural networks, equipped with meta-data [119]. The provided architectures take too long in training neural networks, and then one needs first-order optimizers, which are classified in SGD- and Adam-type algorithms. Such approaches do not consume much time and power, simplifying the training. For convolutional neural networks, such as AlexNet, GoogLeNet, ResNet, SqueezeNet [120], and VGG [121], it suffices to use SGD-type algorithms for achieving high test accuracy. In the case of DenseNet [122], Xception [123], ShuffleNet [124], and GhostNet [125], it is necessary to use advanced Adam-type algorithms, such as DiffGrad, Yogi, AdaBelief, AdaBound, AdamInject, and AdaPNM. The more complex the deep neural network that is applied in recognition problems, the more advanced first-order optimization methods it requires. The same proposition holds for recurrent neural networks, where for SVR [126], XGBoost [127], LSTM [128], and GRU [129], it is enough to apply SGD-type and Adam-type algorithms, such as Adam, RAdam, NAdam, and QHAdam. Recurrent neural networks, containing CNN-LSTM [130], CNN-GRU [131], and TCN-LSTM [132] layers, require advanced Adam-type algorithms.

Second-order optimization methods have a higher rate of convergence in the neighborhood of global minimum, but their time consumption is higher. Even quasi-Newton methods, which are dedicated to reducing time consumption and preserving the convergence rate, still consume more time resources compared with first-order optimization algorithms. They can be used in some convolutional neural networks, such as AlexNet, Res-Net, VGG, and SqueezeNet. For these architectures, time and power consumption is not critical, and the most appropriate algorithms are Apollo and AdaHessian. In recurrent neural networks, second-order optimization algorithms show better results than first-order, but they increase the training time. The second-order methods are, however, good in physics-informed neural networks (PINN). For finding the solution of some partial differential equations with initial and boundary conditions, one needs to analyze the loss function. In this case, L-BFGS, SR1, Apollo, and AdaHessian can reach a higher accuracy than first-order optimization methods. Simple PINN, DeepONet [133], and MFNN [134] are not large like CNN, which allows using of second-order optimization methods, which gives solutions with high accuracy. For Riemannian neural networks, which raise the accuracy by determining the geodesic, one prefers to apply Apollo and AdaHessian, because with gradient directions, they analyze the curvature of the loss function. For cases of Riemannian convolutional neural networks [135], however, they cannot reduce time and power consumption, and one loses the test accuracy by applying a first-order optimizer. Therefore, it is necessary to engage algorithms based on information geometry.

Optimization methods based on information geometry have advantages in speed and accuracy. They can be used in convolutional, recurrent, physics-informed, and Riemannian neural networks. According to their principle of work, there are provided quantum neural networks, which correlate with complex-valued neural networks. Natural gradient descent and mirror descent can achieve the same accuracy and consume low time resources as second- and first-order optimization methods, respectively. Such possibilities allow the inclusion of such algorithms in convolutional, recurrent, graph, and auto-encoder neural networks. Natural gradient and mirror descent extend the application of neural networks in quantum computations, where vanilla natural gradient descent is modified into its quantum analog. For physics-informed neural networks opportunity to converge in the neighborhood of global minimum, consuming time and power similar to first-order methods, can improve the process of solving partial differential and integral equations. This allows it to compete with traditional finite difference, element, and volume methods.

A summary of the above areas of applications is in

Table 5. Applications of optimization algorithm, divided into four types.

Type of Optimization Algorithm	Optimizer	Application	Advantages	Disadvantages
SGD-type	SGD	PINN [136], SNN [137], CVNN [138], AE [139]	These optimizers are fast and can easily be customized. They still meet in many modern neural networks. The convergence rate is from $O\left(\sqrt{T}\right)$ to $O(logT).$	These optimizers cannot reach the global minimum of the loss function. As the consequence, the training accuracy is decreasing. The majority do not have regret bound estimation.
	AdaGrad	CNN [140]
	AdaDelta	CNN [141], RNN [141], SNN [142]
	RMSProp	CNN [143], RNN [144], SNN [145]
	SGDW	CNN [37]
	SGDP	CNN [38]
	QHM	CNN [39]
	NAG	CNN [146], RNN [147]
Adam-type	Adam	CNN [148], RNN [149], SNN [150], PINN [151], GNN [152] CVNN [153]	ADue to exponential moving averages and their modification, the optimization process is more accurate and rapid. The convergence rate can be improved from $O\left(\sqrt{T}\right)$ to $O(logT)$. The number of parame ters is extended, which makes the optimization more controllable.	AOnly DiffGrad, Yogi, AdaBelief, AdaBound, and AdamInject can reach the global minimum of the loss function. They are appropriate for CNN and RNN.
	AdamW	CNN [44]
	AdamP	CNN [46]
	QHAdam	CNN [49]
	Nadam	CNN [51]
	Radam	CNN [52]
	DiffGrad	CNN [154], RNN [60], GNN [155]
	Yogi	CNN [156], RNN [157]
	AdaBelief	CNN [158], RNN [158], GNN [159]
	AdaBound	CNN [160], RNN [161]
	AdamInject	CNN [69]
PNM-type	PNM	CNN [73]	These optimizers are improved by positive–negative moment estimations, which help to reach the global extreme.	They are appropriate only for CNN.
	AdaPNM	CNN [73]
	Adan	CNN [74]
Newton	Newton approach	CNN [162]	These optimizers can reach the global minimum using less iterations. They can be extended on non-Euclidean domains.	They are appropriate only for deep CNN, GNN, and PINN. The optimization process is too long.
	CG	CNN [163], GNN [164]
Quasi-Newton	(L-)BFGS	PINN [165]	These optimzers are faster than Newton approaches. Their main ability is to achieve the global minimum in a short time.	These algorithms are not fast as first-order approaches. They are not useful for deep CNN. Only Apollo has a regret bound.
	SR1	CNN [166]
	Apollo	CNN [97]
	AdaHessian	CNN [98]
Information geometry	NGD	CNN [167], RNN [114], GNN [168], PINN [169], QNN [169]	These optimizers are novel and can be useful in neural networks of any type. The set of hyperparameters is controllable and wider than that in first-order approaches.	The mathematical model of these optimization methods is too complex for customization. Not many probability distributions and potential functions have been investigated for NGD and MD, respectively.
	MD	CNN, RNN [170]

.

Note that CNN, RNN, GNN, PINN, SNN, CVNN, and QNN are convolutional, recurrent, graph, physics-informed, spiking, complex-valued, and quantum neural networks, respectively. These neural networks have a proven ability to solve various problems related to recognition, prediction, generation, processing, detection, and so on. All of them belong to the set of neural networks with gradient-based architectures. Recall that machine learning is the theory that studies self-learning algorithms, which also can be classified. Gradient-based neural networks present one of the classes in machine learning. Meanwhile, there are many advanced gradient and gradient-free learning methods.

6. Challenges and Potential Research

6.1. Promising Approaches in Optimization

Despite the advancements in optimization methods, there exist problems, which concern the fundamental theory of machine learning. All algorithms presented above are utilizable for neural networks with gradient backpropagation, which is not a unique way to rectify weights. Potential methods are presented in [171], which allow reducing gradient calculation, i.e., making the gradient-free error backpropagation process. Therefore, for loss function minimization, one needs to engage alternative optimization methods. For example, the alternating direction method of multipliers in [172] or ensemble neural network in [173], where the gradient is reduced. For such models, one needs to use the following gradient-free optimization methods in

Table 6. Types of gradient-free optimization algorithm.

Type of Optimization Algorithm	Optimizer
Local optimization	Hill Climbing [174],
	Stochastic Hill Climbing [175],
	Simulated Annealing [176],
	Downhill Simplex Optimization [177]
Global optimization	Random Search [178],
	Grid Search [179],
	Random Restart Hill Climbing [180],
	Random Annealing [181],
	Pattern Search [182],
	Powell’s Method [183]
Population-based optimization	Parallel Tempering [184],
	Particle Swarm Optimization [185],
	Spiral Optimization [186],
	Evolution Strategy [187]
Sequential model-based optimization	Bayesian Optimization [188],
	Lipschitz Optimization [189],
	Tree of Parzen Estimators [190]

.

Hybrid optimization algorithms in recent research attract the attention of data scientists. Such approaches are available not only for single neural networks but also for extended ensemble models. Ensemble learning can include neural networks, decision trees, evolutionary algorithms, and other machine learning models. Hybrid optimizers are combinations of different types of algorithms. Hybrid algorithm traversing particle swarm and gradient descent excited by [191]. A combination of gradient descent and genetic optimization is found in [192]. Bio-inspired techniques, traversing population-based and gradient-based optimizers, are presented in [193]. The proposed methods have an application in modern neural networks and show a higher convergence rate than non-hybrid optimizers. The provided gradient-based and gradient-free optimizers can be used in the parallel optimization problems. This question remains in the multi-disciplinary optimization of a typical transport aircraft wing as an example. In the distributed optimization of multi-agent systems, agents cooperate for the global function minimization, which is the sum of local objective functions. Depending on the applications, including energy systems, smart buildings, smart manufacturing, and sensor networks, many distributed optimization algorithms have been developed. In these algorithms, gradient-based and gradient-free optimization algorithms can play an important role.

Gradient-free optimization methods have attracted the interest of many researchers because of their computational simplicity, absence of unwanted differentiation operators, and convergence rate, which is not lower than that of the usual gradient descent. Such an approach has allowed mathematicians to develop alternative mathematical models of neural networks, which show their effectiveness over gradient-based neural networks in [194]. However, this does not mean that one needs to refuse gradient-based optimization algorithms. They have worthy continuation, which is an engaging fractional derivative for computing the gradient.

Fractional calculus, a reasonable continuation of classical calculus, influences theories of partial differential and integral equations, approximations, signal processing, and optimizations. Attempts to generalize differential operators have implicated various properties and helpful propositions concerning machine learning. This can be seen in extending gradient-based optimization methods by fractional derivatives.

Table 7. Types of fractional derivatives on finite interval $[a,b]$ for gradient descent.

Type of Fractional Derivatives	Formulas
Riemann–Liouville	${(}_{RL}{D}_{a+}^{\alpha }f)\left(t\right)=\frac{1}{\Gamma (n-\alpha )}{\left(\frac{d}{dt}\right)}^{\kappa }{\int }_a^t\frac{f\left(\tau \right)}{{(\tau -t)}^{\alpha -\kappa +1}}d\tau$,
	${(}_{RL}{D}_{b-}^{\alpha }f)\left(t\right)=\frac{1}{\Gamma (n-\alpha )}{\left(\frac{d}{dt}\right)}^{\kappa }{\int }_t^b\frac{f\left(\tau \right)}{{(\tau -t)}^{\alpha -\kappa +1}}d\tau$,
	where $Re\left(\alpha \right)>0$
Liouville–Sonine–Caputo	${(}_{LSC}{D}_{a+}^{\alpha }f)\left(t\right)=\frac{1}{\Gamma (\kappa -\alpha )}{\int }_a^t\frac{f^{\left(\kappa \right)}\left(\tau \right)}{{(\tau -t)}^{\alpha -\kappa +1}}d\tau$,
	${(}_{LSC}{D}_{b-}^{\alpha }f)\left(t\right)=\frac{{(-1)}^{\kappa }}{\Gamma (\kappa -\alpha )}{\int }_t^b\frac{f^{\left(\kappa \right)}\left(\tau \right)}{{(\tau -t)}^{\alpha -\kappa +1}}d\tau$,
	where $Re\left(\alpha \right)>0$ and $\kappa =\left[Re\right(\alpha \left)\right]+1$
Tarasov	${(}_T{D}_{a+}^{\alpha }f)\left(t\right)=\frac{\alpha }{\Gamma (1-\alpha )}{\int }_0^{\infty }\frac{f\left(\tau \right)-f(t-\tau )}{{\left(\tau \right)}^{\alpha +1}}d\tau$,
	${(}_T{D}_{b-}^{\alpha }f)\left(t\right)=\frac{\alpha }{\Gamma (1-\alpha )}{\int }_0^{\infty }\frac{f\left(\tau \right)-f(t+\tau )}{{\left(\tau \right)}^{\alpha +1}}d\tau$,
	where $0<\alpha <1$ and $a=0$, $b=\infty$
Hadamard	${(}_H{D}_{a+}^{\alpha }f)\left(t\right)=\frac{1}{\Gamma \left(\alpha \right)}{\left(t\frac{d}{dt}\right)}^{\kappa }{\int }_a^t{\left(log\frac{t}{\tau }\right)}^{\sigma }f\left(\tau \right)\frac{d\tau }{tau}$,
	${(}_H{D}_{b-}^{\alpha }f)\left(t\right)=\frac{1}{\Gamma \left(\alpha \right)}{\left(t\frac{d}{dt}\right)}^{\kappa }{\int }_a^t{\left(log\frac{\tau }{t}\right)}^{\sigma }f\left(\tau \right)\frac{d\tau }{\tau }$,
	where $Re\left(\alpha \right)>0$, $\sigma \in \mathbb{C}$ and $\kappa =\left[Re\right(\alpha \left)\right]+1$
Marchaud	${(}_M{D}_{+}^{\alpha }f)\left(t\right)=\frac{\alpha }{\Gamma (1-\alpha )}{\int }_{-\infty }^t\frac{f\left(t\right)-f\left(\tau \right)}{{(t-\tau )}^{\alpha +1}}d\tau$,
	${(}_M{D}_{-}^{\alpha }f)\left(t\right)=\frac{\alpha }{\Gamma (1-\alpha )}{\int }_t^{\infty }\frac{f\left(t\right)-f\left(\tau \right)}{{(t-\tau )}^{\alpha +1}}d\tau$,
	where $0<Re\left(\alpha \right)<1$
Liouville–Weyl	${(}_{LW}{D}_{+}^{\alpha }f)\left(t\right)=\frac{1}{\Gamma (\kappa -\alpha )}\frac{d^{\kappa }}{d{t}^{\kappa }}{\int }_{-\infty }^t\frac{f\left(\tau \right)}{{(t-\tau )}^{\alpha -\kappa +1}}d\tau$,
	${(}_{LW}{D}_{+}^{\alpha }f)\left(t\right)=\frac{1}{\Gamma (\kappa -\alpha )}\frac{d^{\kappa }}{d{t}^{\kappa }}{\int }_t^{\infty }\frac{f\left(\tau \right)}{{(\tau -t)}^{\alpha -\kappa +1}}d\tau$,
	where $0<Re\left(\alpha \right)$ and $-\infty <x<b<\infty$
Sabzikar–Meerschaert–Chen	${(}_{SMC}{D}_{+}^{\alpha ,\lambda }f)\left(t\right)=\frac{\alpha }{\Gamma (1-\alpha )}{\int }_0^{\infty }\frac{e^{-\lambda \tau }(f\left(t\right)-f(t-\tau ))}{{\tau }^{\alpha +1}}d\tau$,
	${(}_{SMC}{D}_{+}^{\alpha ,\lambda }f)\left(t\right)=\frac{\alpha }{\Gamma (1-\alpha )}{\int }_0^{\infty }\frac{e^{-\lambda \tau }(f\left(t\right)+f(t-\tau ))}{{\tau }^{\alpha +1}}d\tau$,
	where $0<Re\left(\alpha \right)$ and $\lambda \in \mathbb{C}$
Katugampola	${(}_K{D}_{a+,\sigma }^{\alpha ,\lambda }f)\left(t\right)=\frac{{\sigma }^{\alpha -\kappa +1}{e}^{\lambda t}}{\Gamma (\kappa -\alpha )}{\left(t^{1-\sigma }\frac{d}{dt}\right)}^{\kappa }{\int }_a^t\frac{{\tau }^{\sigma -1}{e}^{-\lambda \tau }f\left(\tau \right)}{{(t^{\sigma }-{\tau }^{\sigma })}^{\alpha -\kappa +1}}$,
	${(}_K{D}_{b-,\sigma }^{\alpha ,\lambda }f)\left(t\right)=\frac{{\sigma }^{\alpha -\kappa +1}{e}^{\lambda t}}{\Gamma (\kappa -\alpha )}{\left(t^{1-\sigma }\frac{d}{dt}\right)}^{\kappa }{\int }_t^b\frac{{\tau }^{\sigma -1}{e}^{-\lambda \tau }f\left(\tau \right)}{{({\tau }^{\sigma }-t^{\sigma })}^{\alpha -\kappa +1}}$,
	where $0<Re\left(\alpha \right)$ and $\lambda \in \mathbb{C}$

summarizes all known fractional derivatives [195], which can improve the minimization of loss functions in neural networks.

The efficiency of fractional optimization methods has been demonstrated in [196], but there is the following question: does the chain rule for fractional derivatives make the error propagation modeling more complex? No, because there are implied generalized chain and Leibniz rules in [197], which do not differ much from the usual versions. Therefore, it is possible to generalize first-order optimization methods from SGD- to Adam-type algorithms.

Simultaneously, there is the problem of developing neural networks containing bilevel optimization [198]. For neural networks equipped with meta-learning, utilization of the provided optimization methods does not suffice. There have to be used bilevel optimization algorithms, such as BSA [199], TTSA [200], HOAG [201], AID-FP [202], AID-CG [203], and stocBiO [204]. Such approaches provide the theoretical guarantee in hyperparameter optimization, meta, and ensemble learning. Potentially, these methods can be improved using information geometry or techniques from advanced first-order optimization algorithms.

Besides meta and ensemble learning, there exist soft computing and federated analogs. Soft computing is a set of techniques aimed at modeling and solving real-world problems that are difficult to solve mathematically. These approaches, approximately equivalent to human decision-making, are designed to allow for partial truth, approximate reasoning, uncertainty, and imprecision. This is different from the usual hard computing model, which relies entirely on numerical calculations and logical reasoning. SC is an integrative field in which there is a new phase of artificial intelligence called computational intelligence. Federated learning is an effective learning strategy for disparate datasets [205]. It prevents leakage of sensitive information when training a model on data from multiple devices. FL has received a lot of attention, which has served as motivation for several useful initiatives to build learning apps on a wide variety of decentralized devices. FL provides decentralized learning that does not require the transfer of raw data between nodes, thereby protecting user information. Moreover, FL guarantees a reduction in communication costs between the server and the client. The communication between the client and the server has been shortened, as training-related client data is not sent to the server. The advantages of FL include increased privacy and lower communication costs. FL is used in situations where respect for confidentiality and privacy is paramount [206].

6.2. Open Problems in the Modern Theory of Neural Network

Another problem in physics-informed neural networks is their extension on delay differential equations. The problem of solving such equations with various delays has not yet been solved. There is a model for solving fractional differential equations in [207]. In this case, one raises the question about using different activation functions, which, consequently, impacts the loss function optimization. Therefore, one can use optimizers based on information geometry. Afterward, they can be induced in delay physics-informed models. In particular, information-geometric methods are relevant for complex-valued neural networks.

Complex-valued neural networks demonstrate their efficiency in engineering areas such as antenna design, radar imaging, optical/lightwave information processing, and quantum devices such as superconductive devices. All these areas of applications contain mathematical models with rotational points, wave functions, and integral transforms such as Fourier, Laplace, and Hilbert. This model suggests the future realization of intrinsically non-von Neumann computers, including pattern-information-representing devices. Conventional quantum computation is strictly limited in its treatable problems. In contrast, CVNN-based quantum computing can solve more general tasks, leading to an application of quantum informatics. Therefore, it is necessary to imply optimization methods based on quantum and tensor computing. For example, Shampoo [208] can perform optimization at the tensor level.

Quantum neural networks have attracted the attention of many researchers who study and develop advanced machine learning methods. With the appearance of quantum computers, it was necessary to expand the theory of neural networks on quantum devices, which allowed analyzing quantum processes and tensor-network circuits. For problems of image recognition, time series prediction, and moving object detection, one needs to use convolutional neural networks [209], which engage a quantum natural gradient. This optimization algorithm presents a quantum analog of vanilla natural gradient descent, but it uses the Fubini–Study tensor instead of the Fisher information matrix. Unlike the usual Fisher and Hessian matrices, such a tensor considers quantum computations with wave functions, bra- and ket-vectors. Such a method can be equipped with exponential moving and positive–negative averages. There is the proposition of combining quantum neural networks with spiking neurons in recent research. This will lead to other tasks for developing optimization methods with memory.

Recall that neural networks of the third generation are based on weights with memory. Such a model is close to biological neural networks and makes the learning process more accurate. For such models, one can use corresponding memory-based optimization algorithms, such as mixed stochastic gradient descent, denoted MEGA I and MEGA II in [210]. The ability to remember past experiences while being trained on a continuum of tasks makes it possible to raise the test accuracy of spiking neural networks. This kind of neural network has experienced extension from multilayer perceptrons to various advanced networks, such as feedback SNN, lattice map SNN, deep convolutional SNN, spike-timing-dependent plasticity, and many other networks. Howeover, one of the most recent architectures is graph networks, which can successfully train models without spiking neurons.

Graph neural networks, whose structures are reminiscent of simple graphs, are utilized in many tasks related to medicine and biology. As can be seen in neuroscience resources, biological neural networks do not present a linear, sequential, and organized determined-order model but construct other more progressive connections, which significantly impact mathematical models of neurons. If it is possible to create a model with the corresponding structure, then one receives the network without unwanted layers and computations. However, such a model yields the error propagation problem, which differs from classical error backpropagation. Subsequently, alternative error rectification methods were proposed, called neighbor aggregation and information update in [211]. For these methods, it is suitable to apply the first-order information-geometric optimizers. One of these is mirror descent. Such an algorithm can be modified with new proximal functions and extended versions of Bregman divergence in further research. Nevertheless, there are divergence formulas that can produce new information-geometric optimization methods. For modifying the structure of neural networks, however, one may use wavelet decompositions, which can potentially process input data more accurately.

Neural networks, based on a wavelet, are dispersed in signal processing. Moreover, such a tool is used in numerical solving of partial differential equations, which seems an important advance of physics-informed neural networks. Wavelet decomposition has been used in convolution, LSTM, and GRU layers, which makes data processing more accurate, especially in cases of scaling input information. Moreover, a model of graph-wavelet neural networks was present that gave the best results of node classification in [212], compared with spectral CNN, ChebyNet, GCN, and MoNet. For wavelet neural networks, one can use algorithms presented in Table 5. However, there are whale and butterfly optimization methods, which are presented in [213]. In the case of binary neural networks [214], the provided optimizers can be useful. For such a network, researchers often utilize the particle swarm method, which is gradient-free. Therefore, there is a possibility of comparing gradient-free optimization methods between each other and developing new approaches.

In summary, one needs to say that developing optimization algorithms is correlated with the evolution of neural network architectures and challenges posed before research. Moreover, such a relationship works backward, which happens in the case of quantum natural gradient descent and quantum neural networks. Note that neural networks can be improved only in terms of the pattern recognition problem. Therefore, one needs to consider the moving detection problem. In the case of time series predictions, the stochastic process and Brownian motion should be considered, which come from statistical physics and thermodynamics. Such challenges have an influence on the theory of not only neural networks but also of artificial intelligence.

7. Conclusions

In the conclusion of this review, we can say that fundamental development of the theory of neural networks allows us to simplify the work of humans on additional scales. Optimization methods have allowed us to achieve not only higher test accuracy in a short time but to imply the necessary features and disadvantages of existing models, which through time has allowed research to provide advanced architectures. First-order optimization algorithms, presented in Section 2, have improved and are changing simultaneously with the growth of neural networks in size and quality. SGD-, Adam- and PNM-type algorithms are well described, as are their properties, improvements to the techniques, and evolution. We presented second-order optimization algorithms, Newton and quasi-Newton methods, in Section 3. We also presented their applications in neural networks, simplifying their computation complexity using approximation theory, which produced other modifications and variations. Such an approximation led to changing the Hessian by gradient flow tensor and Fisher information matrix. This modification brought us to information geometry, which allowed us to provide natural gradient and mirror descents, whose equivalence is introduced in Section 4. Afterward, we summarized in Table 5 all types of neural networks suitable for the introduced optimization methods. Modern networks are provided in which the introduced optimization algorithms have already been used and those that could be used for increasing test accuracy. In Section 6, we demonstrated further ways of developing and applications of optimization methods. Also provided are gradient-free, fractional order, and bilevel optimization algorithms, their existing versions, and potentials for increasing test accuracy for various types of neural networks. In the end, this survey has presented to readers all types of existing optimization methods, their modifications, and applications, which can help in studies to comprehend the state of the modern theory of optimization and machine learning. Such information allows us to create advanced fundamental modifications and extend the area of applications for modern neural networks of all types. Optimization algorithms are constructed according to the practice and neural network architectures. However, their approach to minimizing the loss functions impacts the neural network architectures. The evolution of convolutional neural networks has passed along with the development of Adam-type optimizers, and conversely. Newton and quasi-Newton optimization methods have extended the neural networks to their geometrical analogs, such as Riemannian and Kahlerian, for real- and complex-valued models, respectively. Information-geometric optimization algorithms influence the development of quantum neural networks, where the loss function minimization is based on quantum natural gradient descent, which calculates the Fisher–Rao metrics instead of the Fisher matrix. Gradient-based and gradient-free optimizers impact the development of the modern theory of ensemble learning.