Title: ELLIPSOIDAL TRUST REGION METHODS FOR NEURAL NETWORK TRAINING

Abstract: We investigate the use of ellipsoidal trust region constraints for second-order optimization of neural networks. This approach can be seen as a higher-order counterpart of adaptive gradient methods, which we here show to be interpretable as first-order trust region methods with ellipsoidal constraints. In particular, we show that the preconditioning matrix used in RMSProp and Adam satisfies the necessary conditions for provable convergence of second-order trust region methods with standard worst-case complexities. Furthermore, we run experiments across different neural architectures and datasets to find that the ellipsoidal constraints constantly outperform their spherical counterpart both in terms of number of backpropagations and asymptotic loss value. Finally, we find comparable performance to state-ofthe-art first-order methods in terms of backpropagations, but further advances in hardware are needed to render Newton methods competitive in terms of time.

1. INTRODUCTION:

We consider finite-sum optimization problems of the form

EQUATION

which typically arise in neural network training, e.g. for empirical risk minimization over a set of data points (x i , y i ) ∈ R in × R out , i = 1, . . . , n. Here, : R out × R out → R + is a convex loss function and f : R in × R d → R out represents the neural network mapping parameterized by the concatenation of layers w ∈ R d , which is non-convex due to its multiplicative nature and potentially non-linear activation functions. We assume that L is lower bounded and twice differentiable, i.e. L ∈ C 2 (R d , R) and consider finding a first-and second-order stationary point w for which ∇L( w) ≤ g and λ min ∇ 2 L( w) ≥ -H . Non-convex optimization problems are ubiquitous in machine learning. Among the most prominent examples are present-day deep neural networks, that have achieved outstanding results on core tasks such as collaborative filtering (Wang et al., 2015) , sentence classification (Kim, 2014) and image classification (Krizhevsky et al., 2012) .

In the era of big data and deep neural networks, stochastic gradient descent (SGD) is one of the most widely used training algorithms (Bottou, 2010) . What makes SGD so attractive is its simplicity and per-iteration cost that are independent of the size of the training set (n) and scale linearly in the dimensionality (d). However, gradient descent is known to be inadequate to optimize functions that are ill-conditioned (Nesterov, 2013; Shalev-Shwartz et al., 2017) and thus adaptive gradient methods that employ dynamic, coordinate-wise learning rates based on past gradients-including Adagrad (Duchi et al., 2011) , RMSprop (Tieleman & Hinton, 2012) and Adam (Kingma & Ba, 2014)-have become a popular alternative, often providing significant speed-ups over SGD. Yet, there exist no theoretical proofs that these methods are faster than gradient descent (Li & Orabona, 2018) .

From a theoretical perspective, Newton methods provide stronger convergence guarantees by appropriately transforming the gradient in ill-conditioned regions according to second-order derivatives. It is precisely this Hessian information that allows regularized Newton methods to enjoy superlinear local convergence as well as to provably escape saddle points (Conn et al., 2000) . While second-order algorithms have a long-standing history even in the realm of neural network training (Hagan & Menhaj, 1994; Becker et al., 1988) , they were mostly considered as too computationally and memory expensive for practical applications. Yet, the seminal work of Martens (2010) renewed interest for their use in deep learning by proposing efficient Hessian-free methods that only access second-order information via matrix-vector products which can be computed at the cost of an additional backpropagation (Pearlmutter, 1994; Schraudolph, 2002) . Among the class of regularized Newton methods, trust region (Conn et al., 2000) and cubic regularization algorithms (Cartis et al., 2011) are the most principled approaches in the sense that they yield the strongest convergence guarantees. Recently, stochastic extensions have emerged (Xu et al., 2017b; Yao et al., 2018; Kohler & Lucchi, 2017; Gratton et al., 2017) , which suggest their applicability for deep learning.

We here propose a simple modification to make TR methods even more suitable for neural network training. Particularly, we build upon the following alternative view on adaptive gradient methods:

While gradient descent can be interpreted as a spherically constrained first-order TR method, preconditioned gradient methods-such as Adagrad-can be seen as first-order TR methods with ellipsoidal trust region constraint.

This observation is particularly interesting since spherical constraints are blind to the underlying geometry of the problem, but ellipsoids can adapt to local landscape characteristics, thereby allowing for more suitable steps in regions that are ill-conditioned. We will leverage this analogy and investigate the use of the Adagrad and RMSProp preconditioning matrices as ellipsoidal trust region shapes within a stochastic second-order TR algorithm (Xu et al., 2017a; Yao et al., 2018) . Since no ellipsoid fits all objective functions, our main contribution lies in the identification of adequate matrix-induced constraints that lead to provable convergence and significant practical speed-ups for the specific case of deep learning. On the whole, our contribution is threefold:

• We provide a new perspective on adaptive gradient methods that contributes to a better understanding of their inner-workings. Furthermore, we empirically find that many neural network problems exhibit diagonally dominated Hessian matrices which suggests the effectivity of diagonal preconditioning. (Section 3) • We investigate the first application of ellipsoidal TR methods for deep learning. In Theorem 1 we show that the RMSProp matrix can directly be applied as constraint inducing norm in second-order TR algorithms while preserving all convergence guarantees (Theorem 2). • Finally, we provide an experimental benchmark across different real-world datasets and architectures. We also compare against adaptive gradient methods and show results in terms of backprogations, epochs, and wall-clock time; a comparison we were not able to find in the literature. (Section 5)

Our main empirical results demonstrate that ellipsoidal constraints prove to be a very effective modification of the trust region method in the sense that they constantly outperform the spherical TR method, both in terms of number of backprogations and asymptotic loss value on a variety of tasks.

2. RELATED WORK:

First-order methods The prototypical method for optimizing Eq. ( 1) is SGD (Robbins & Monro, 1951) . While the practical success of SGD in non-convex optimization is unquestioned, the theoretical foundation of this phenomenon is still rather limited. Recent findings suggest the ability of this method to escape saddle points and reach local minima in polynomial time for general non-convex problems, but they either need to artificially add noise to the iterates (Ge et al., 2015; Lee et al., 2016) or make an assumption on the inherent noise of vanilla SGD (Daneshmand et al., 2018) . For neural network training, a recent line of research proclaims the effectiveness of SGD, but the results usually come at the cost of fairly strong assumptions such as heavy overparametrization and Gaussian inputs (Du et al., 2017; Brutzkus & Globerson, 2017; Li & Yuan, 2017; Du & Lee, 2018; Allen-Zhu et al., 2018) . Adaptive gradient methods (Duchi et al., 2011; Tieleman & Hinton, 2012; Kingma & Ba, 2014) build on the intuition that larger learning rates for smaller gradient components and smaller learning rates for larger gradient components balance their respective influences and thereby make the methods behave as if they were optimizing a more isotropic surface. Such approaches have first been suggested for neural networks by LeCun et al. (2012) . Recently, convergence guarantees for such methods are starting to appear (Ward et al., 2018; Li & Orabona, 2018) . However, these are not superior to the O( -2 g ) worst-case complexity of standard gradient descent (Cartis et al., 2012b) .

Regularized Newton methods:

The most principled class of regularized Newton methods are trust region (TR) and adaptive cubic regularization algorithms (ARC) (Conn et al., 2000; Cartis et al., 2011) , which repeatedly optimize a local Taylor model of the objective while making sure that the step does not travel too far such that the model stays accurate. While the former finds first-order stationary points within O( -2 g ), ARC only takes at most O( -3/2 g

). However, simple modifications to the TR framework allow these methods to obtain the same accelerated rate (Curtis et al., 2017) . Both methods take at most O( -3 H ) iterations to find an H approximate second-order stationary point (Cartis et al., 2012a) . These rates are optimal for second-order Lipschitz continuous functions (Carmon et al., 2017; Cartis et al., 2012a) and they can be retained even when only sub-sampled gradient and Hessian information is used (Kohler & Lucchi, 2017; Yao et al., 2018; Xu et al., 2017b; Blanchet et al., 2016; Liu et al., 2018; Cartis & Scheinberg, 2017) . Furthermore, the involved Hessian information can be computed solely based on Hessian-vector products, which are implementable efficiently for neural networks (Pearlmutter, 1994) . This makes these methods particularly attractive for deep learning, but the empirical evidence of their applicability is so far very limited. We are only aware of the works of Liu et al. (2018) and Xu et al. (2017a) , which report promising first results but these are by no means fully encompassing.

Gauss-Newton methods An interesting line of research proposes to replace the Hessian by (approximations of) the generalized-Gauss-Newton matrix (GGN) within a Levenberg-Marquardt framework1 (LeCun et al., 2012; Martens, 2010; Martens & Grosse, 2015) . These methods have been termed hessian-free since only access to GGN-vector products is required. As the GGN matrix is always positive semidefinite, they cannot leverage negative curvature to escape saddles and hence, there exist no second-order convergence guarantees. Furthermore, there are cases in neural network training where the Hessian is better conditioned than the GGN matrix (Mizutani & Dreyfus, 2008) . Nevertheless, the above works report promising preliminary results, most notably Grosse & Martens (2016) report that K-FAC can be faster than SGD on a small convnet. On the other hand, recent findings report performance at best comparable to SGD on the much larger ResNet architecture (Ma et al., 2019) . Moreover, Xu et al. (2017a) reports many cases where TR and GGN algorithms perform similarly. This line of work is to be seen as complementary to our approach since it is straight forward to replace the Hessian in the TR framework with the GGN matrix. Furthermore, the preconditioners used in Martens (2010) and Chapelle & Erhan (2011) , namely diagonal estimates of the empirical Fisher and Fisher matrix, respectively, can directly be used as matrix norms in our ellipsoidal TR framework.

3. AN ALTERNATIVE VIEW ON ADAPTIVE GRADIENT METHODS:

Adaptively preconditioned gradient methods update iterates as w t+1 = w t -η t A -1/2 t g t , where g t is a stochastic estimate of ∇L(w t ) and A t is a positive definite symmetric pre-conditioning matrix. In Adagrad, A ada,t is the un-centered second moment matrix of the past gradients computed as

EQUATION

where > 0, I is the d × d identity matrix and

G t = [g 1 , g 2 , . . . , g t ].

Building up on the intuition that past gradients might become obsolete in quickly changing non-convex landscapes, RMSprop (and Adam) introduce an exponential weight decay leading to the preconditioning matrix

EQUATION

where β ∈ (0, 1). In order to save computational efforts, the diagonal versions diag(A ada ) and diag(A rms ) are more commonly applied in practice, which in turn gives rise to coordinate-wise adaptive stepsizes that are enlarged (reduced) in coordinates that have seen past gradient components with a smaller (larger) magnitude. In that way, the optimization methods can account for gradients of potentially different scales arising from e.g. different layers of the networks.

3.1. ADAPTIVE PRECONDITIONING AS ELLIPSOIDAL TRUST REGION:

Starting from the fact that adaptive methods employ coordinate-wise stepsizes, one can take a principled view of these methods. Namely, their update steps arise from minimizing a first-order Taylor model of the function L within an ellipsoidal search space around the current iterate w t , where the diameter of the ellipsoid along a particular coordinate is implicitly given by η t and g t A -1 t . Correspondingly, vanilla (S)GD optimizes the same first-order model within a spherical constraint. Fig. 1 (top) illustrates this effect by showing not only the iterates of GD and Adagrad but also the implicit trust regions within which the local models were optimized at each step.2 Since the models are linear, the constrained minimizer is always found on the boundary. It is well known that GD struggles to progress towards the minimizer of quadratics along lowcurvature directions (see e.g., Goh (2017) ). While this effect is negligible for well-conditioned objectives (Fig. 1 , left), it leads to drastically slow-down when the problem is ill-conditioned (Fig. 1 , center). Particularly, once the method has reached the bottom of the valley, it struggles to make progress along the horizontal axis. Here is precisely where the advantage of adaptive stepsize methods comes into play. As illustrated by the dashed lines, Adagrad's search space is damped along the direction of high curvature (vertical axis) and elongated along the low curvature direction (horizontal axis). This allows the method to move further horizontally early on to enter the valley with a smaller distance to the optimizer w * along the low curvature direction which accelerates convergence.

Theorem 1 (Preconditioned gradient methods as TR). A preconditioned gradient step

w t+1 -w t = s t := -η t A -1 t g t (4)

with stepsize η t > 0, symmetric positive definite preconditioner A t ∈ R d×d and g t = 0 minimizes a first-order model around w t ∈ R d in an ellipsoid given by A t in the sense that

EQUATION

)

Corollary 1 (Rmsprop). The step s rms,t := -η t A -1/2 rms,t g t minimizes a first-order Taylor model around w t in an ellipsoid given by A 1/2 rms,t (Eq. 3) in the sense that s rms,t := arg min

EQUATION

Equivalent results can be established for Adam using g adam,t := (1 -β) t k=0 β t-k g t as well as for Adagrad by replacing the matrix A ada into the constraint in Eq. ( 6). Of course, the update procedure in Eq. ( 5) is merely a reinterpretation of the original preconditioned update, and thus the employed trust region radii are defined implicitly by the current gradient and stepsize.

3.2. DIAGONAL VERSUS FULL PRECONDITIONING:

A closer look at Fig. 1 reveals that the first two problems come with level sets that are perfectly axis-aligned, which makes these objectives particularly attractive for diagonal preconditioning. For comparison, on the right of Fig. 1 , we report another quadratic problem instance, where the Hessian is no longer zero on the off-diagonals. As can be seen, the interaction between coordinates introduces a tilt in the level sets and reduces the superiority of diagonal Adagrad over plain GD. However, using the full preconditioner A ada re-establishes the original speed up. Yet, non-diagonal preconditioning comes at the cost of taking the inverse square root of a large matrix, which is why this approach has been relatively unexplored (see Agarwal et al. (2018) for a recent exception).

Interestingly, early results by Becker et al. (1988) on the curvature structure of neural nets report a strong diagonal dominance of the Hessian matrix ∇ 2 L(w). This suggests that the loss surface is indeed somewhat axis-aligned. However, the reported numbers are only for tiny feed-forward networks of at most 256 parameters. Therefore, we generalize these findings in the following to larger networks. Furthermore, we contrast the diagonal dominance of real Hessian matrices to the expected behavior of random Wigner matrices. Of course, true Hessians do not have i.i.d. entries but the symmetry of Wigner matrices suggests that this baseline is not completely off. Furthermore, we also compare Hessians of Ordinary Least Squares problems (OLS) with random inputs. For this purpose, let δ A define the ratio of diagonal to overall mass of a matrix A, i.e. δ A := i |Ai,i| i j |Ai,j | as in (Becker et al., 1988) .

Proposition 1 (Diagonal share of Wigner matrix). For random Gaussian3 Wigner matrix W formed as

EQUATION

where ∼ stands for i.i.d. draws (Wigner, 1993) , the diagonal mass of the expected absolute matrix amounts to

EQUATION

Thus, if we suppose the Hessian at any given point w were a random Wigner matrix we would expect the share of diagonal mass to fall with O(1/d) as the network grows in size. A similar result can be derived for the large n limit in the case of OLS Hessians.

Proposition 2 (Diagonal share of OLS Hessian). Let X ∈ R d×n and assume each x i,j is generated i.i.d. with zero-mean finite second moment σ 2 > 0. Then the share of diagonal mass of the expected matrix

EQUATION

Again, we expect the diagonal mass to fall in d and interestingly, empirical simulations suggest that the result of Proposition 2 holds already in small n settings (see Fig. D .2) and it is likely that finite n results can be derived when adding assumptions such as Gaussian data. As can be seen in Fig. 2 below, even for a practical batch size of n = 32 the diagonal mass δ H of neural networks stays above both benchmarks when its input is random as well as with real-world data at random initialization, during training and after convergence. These findings are in line with Becker et al. (1988) and suggest that full matrix preconditioning is most probably not worth the additional computational cost for neural networks. Consequently, we use diagonal preconditioning for both first-and second-order methods in all of our experiments in Section 5. Further theoretical elaboration of these findings present an interesting direction of future research.

4. SECOND-ORDER TRUST REGION METHODS:

Cubic regularization (Nesterov & Polyak, 2006; Cartis et al., 2011) and trust region methods belong to the family of globalized Newton methods. Both frameworks compute parameter updates by optimizing regularized (former) or constrained (latter) second-order Taylor models of the objective L around the current iterate w t . 4 In particular, in iteration t the update step of the trust region algorithm is computed as

EQUATION

)

where ∆ t > 0 and g t and B t are either ∇L(w t ) and ∇ 2 L(w t ) or suitable approximations. The matrix A t induces the shape of the constraint set. So far, the common choice for neural networks is A t := I, ∀t which gives rise to spherical trust regions (Xu et al., 2017a; Liu et al., 2018) . By solving the constrained problem (10), TR methods overcome the problem that pure Newton steps may be ascending, attracted by saddles or not even computable. See Appendix B for more details.

4.1. CONVERGENCE OF ELLIPSOIDAL TRUST REGION METHODS:

Inspired by the success of adaptive gradient methods, we investigate the use of their preconditioning matrices as norm inducing matrices for second-order TR methods. The crucial condition for convergence is that the applied norms are not degenerate during the entire minimization process in the sense that the ellipsoids do not flatten out (or blow up) completely along any given direction. The following definition formalizes this intuition.

Definition 1 (Uniformly equivalent norms). The norms w At := (w A t w) 1/2 induced by symmetric positive definite matrices A t are called uniformly equivalent, if ∃µ ≥ 1 such that

EQUATION

We now establish a result which shows that the RMSProp ellipsoid is indeed uniformly equivalent.

Lemma 1 (Uniform equivalence). Suppose g t 2 ≤ L 2 H for all w t ∈ R d , t = 1, 2, . . . Then there always exists > 0 such that the proposed preconditioning matrices A rms,t (Eq. 3) are uniformly equivalent, i.e. Def. 1 holds. The same holds for the diagonal variant.

Consequently, the ellipsoids A rms,t can directly be applied to any convergent TR framework without losing convergence guarantees (Conn et al. (2000) , Theorem 6.6.8). 5 Interestingly, this result cannot be established for A ada,t , which reflects the widely known vanishing stepsize problem that arises since squared gradients are continuously added to the preconditioning matrix. At least partially, this effect inspired the development of RMSprop (Tieleman & Hinton, 2012) and Adadelta (Zeiler, 2012).

Why ellipsoids? There are many sources for ill-conditioning in neural networks such as un-centered and correlated inputs (LeCun et al., 2012) , saturated hidden units, and different weight scales in different layers (Van Der Smagt & Hirzinger, 1998) . While the quadratic term of model ( 10) accounts for such ill-conditioning to some extent, the spherical constraint is completely blind towards the loss surface. Thus, it is advisable to instead measure distances in norms that reflect the underlying geometry (see Chapter 7.7 in Conn et al. (2000) ). The ellipsoids we propose are such that they allow for longer steps along coordinates that have seen small gradient components in past and vice versa. Thereby the TR shape is adaptively adjusted to fit the current region of the non-convex loss landscape. This procedure is not only effective when the iterates are in an ill-conditioned neighborhood of a minimizer (Figure 1 ), but it also helps to escape elongated plateaus (see autoencoder in Section 5).

4.2. A STOCHASTIC TR FRAMEWORK FOR NEURAL NETWORK TRAINING:

Since neural network training often constitutes a large-scale learning problem in which the number of datapoints n is high, we here opt for a stochastic TR framework in order to circumvent memory issues and reduce the computational complexity. In particular, we adapt the framework of Yao et al. (2018) ; Xu et al. (2017b) to the case of iteration-dependent norm constraints (Algorithm 1). We prove that the O max -2 g -1 H , -3 H worst-case iteration complexity is retained, which is optimal if the function is second-order smooth (Cartis et al., 2012a) . Towards this end, we assume that the derivative estimates are sufficiently accurate in the following sense.

Assumption 1 (Sufficiently accurate derivatives). The approximations of the gradient and Hessian at step t satisfy

g t -∇L(w t ) ≤ δ g and B t -∇ 2 L(w t ) ≤ δ H , where δ g ≤ (1-η) g 4 and δ H ≤ min (1-η)v H 2 , 1 , for some 0 < v < 1.

For finite-sum objectives such as Eq. 1, the above condition can be met by random sub-sampling due to classical concentration results for sums of random variables (Xu et al., 2017b; Kohler & Lucchi, 2017; Tripuraneni et al., 2017) . Following these references, we assume access to the function value in each iteration for our theoretical analysis but we note that convergence can be retained even for fully stochastic trust region methods (Gratton et al., 2017; Chen et al., 2018; Blanchet et al., 2016) and indeed our experiments in Section 5 use sub-sampled function values due to memory constraints.

Algorithm 1 Stochastic Ellipsoidal Trust Region Method 1: Input: w 0 ∈ R d , γ > 1, 1 > η > 0, ∆ 0 > 0 2: for t = 0, 1, . . . , until convergence do 3:

Compute approximations g t and B t . If g t ≤ g , set g t := 0.

4::

Set A t := A rms,t or A t := diag (A rms,t ) (see Eq. ( 3)).

5::

Obtain s t by solving m t (s t ) approximately.

6::

Compute ratio of function over model decrease

EQUATION

Set 

EQUATION

5. EXPERIMENTS:

Trust region methods To validate our claim that ellipsoidal TR methods yield improved performance over spherical ones, we run a set of experiments on two image datasets and three types of network architectures. As can be seen in Figure 3 , the ellipsoidal TR methods consistently outperform their spherical counterpart in the sense that they reach full training accuracy substantially faster on all problems. Moreover, their limit points are in all cases lower than those of the uniform TR method. Interestingly, this makes an actual difference in the image reconstruction quality of autoencoders (see Figure 12 ). We thus draw the clear conclusion that the ellipsoidal trust region constraints we propose are to be preferred over their spherical counterpart when training neural networks. Both the experimental and architectural details are provided in Appendix C. Benchmark with SGD To put the previous results into context, we also benchmark several stateof-the-art gradient methods. We fix their sample size to 32 (as advocated e.g. in Masters & Luschi (2018) ) but grid search the stepsize since it is the ratio of these two quantities that effectively determines the level of stochasticity (Jastrzebski et al., 2017) . As the TR methods have a larger batch size6 of 128-512, we report results both in terms of number of backpropagations and epochs for a fair comparison. A close look at Figure 4 and 10 (Appendix) indicates that the ellipsoidal TR methods can be slightly superior in terms of backprops but at best manage to keep pace with first-order methods in terms of epochs. Furthermore, the limit points of both first-and second-order methods yield the same order of loss in most experiments. When taking gradient norms into account (plot omitted), we indeed find no spurious local minima and only the autoencoders give rise to saddle points.

6. CONCLUSION:

We investigated the use of ellipsoidal trust region constraints for neural networks. We have shown that the RMSProp matrix satisfies the necessary conditions for convergence and our experimental results demonstrate that ellipsoidal TR methods outperform their spherical counterparts significantly. We thus consider the development of further ellipsoids that can potentially adapt even better to the loss landscape such as e.g. (block-) diagonal hessian approximations (e.g. Bekas et al. (2007) ) or approximations of higher order derivatives as an interesting direction of future research.

Yet, the gradient method benchmark indicates that the value of Hessian information for neural network training is limited for mainly three reasons: 1) second-order methods rarely yield better limit points, which suggests that saddles and spurious local minima are not a major obstacle; 2) gradient methods can run on smaller batch sizes which is beneficial in terms of epoch and when memory is limited; 3) The per-iteration time complexity is noticeably lower for first-order methods (Figure 11 ). These observations suggest that advances in hardware and distributed second-order algorithms (e.g., Osawa et al. (2018) ; Dünner et al. (2018) ) will be needed before Newton-type methods can replace (stochastic) gradient methods in deep learning.

REFERENCES:

Naman Agarwal, Brian Bullins, Xinyi Chen, Elad Hazan, Karan Singh, Cyril Zhang, and Yi Zhang.

The case for full-matrix adaptive regularization. arXiv preprint arXiv:1806.02958, 2018.

B EQUIVALENCE OF PRECONDITIONED GRADIENT DESCENT AND FIRST-ORDER TRUST REGION METHODS:

Theorem 3 (Theorem 1 restated). A preconditioned gradient step

EQUATION

)

with stepsize η t > 0, symmetric positive definite preconditioner A t ∈ R d×d and g t = 0 minimizes a first-order local model around w t ∈ R in an ellipsoid given by A t in the sense that

EQUATION

)

Proof. We start the proof by noting that the optimization problem in Eq. ( 14) is convex. For η t > 0 the constraint satisfies the Slater condition since 0 is a strictly feasible point. As a result, any KKT point is a feasible minimizer and vice versa.

Let L(s, λ) denote the Lagrange dual of Eq. ( 5)

EQUATION

)

Any point s is a KKT point if and only if the following system of equations is satisfied

EQUATION

For s t as given in Eq. ( 4) we have that

EQUATION

)

and thus 17 and 18 hold with equality such that any λ ≥ 0 is feasible. Furthermore,

∇ s L(s t , λ) = ∇f (w t ) + λ s t At A t s t (4) = g t -η t λ η t g t A -1 t A t A -1 t g t = g t - λ g t A -1 t g t (21) is zero for λ = g t A -1 ≥ 0.

As a result, s t is a KKT point of the convex problem 5 which proves the assertion.

To illustrate this theoretical result we run gradient descent and Adagrad as well as the two corresponding first-order TR methods 7 on an ill-conditioned quadratic problem. While the method 1st TR optimizes a linear model within a ball in each iteration, 1st TR ada optimizes the same model over the ellipsoid given by the Adagrad matrix A ada . The results in Figure 5 show that the methods behave very similar to their constant stepsize analogues. 

C CONVERGENCE OF ELLIPSOIDAL TR METHODS:

In order to prove convergence results for ellipsoidal Trust Region methods one must ensure that the applied norms are coherent during the complete minimization process in the sense that the ellipsoids do not flatten out (or blow up) completely along any given direction. This intuition is formalized in Assumption 1 which we restate here for the sake of clarity. Definition 2 (Definition 1 restated). There exists a constant µ ≥ 1 such that

EQUATION

)

Towards this end, Conn et al. (2000) identify the following sufficient condition on the basis of which we will prove that our proposed ellipsoid A rms is indeed uniformly equivalent under some mild assumptions.

Lemma 2 (Theorem 6.7.1 in Conn et al. ( 2000)). Suppose that there exists a constant

EQUATION

then Definition 1 holds.

Having uniformly equivalent norms is sufficient to prove convergence of ellipsoidal TR methods (se AN.1 and Theorem 6.6.8 in (Conn et al., 2000) ). However, it is so far unknown how the ellipsoidal constraints influence the convergence rate itself. We here prove that the specific ellipsoidal TR method presented in Algorithm 1 preserves the rate of its spherically-constrained counterpart proposed in Yao et al. (2018) (see Theorem 2 below).

First, we show that the proposed A rms,t ellipsoid satisfies Definition 1.

Lemma 3 (Lemma 1). Suppose g t 2 ≤ L 2 H for all w t ∈ R d , t = 1, 2, . . . Then there always exists > 0 such that the proposed preconditioning matrices A rms,t (Eq. 3) are uniformly equivalent, i.e. Def. 1 holds. The same holds for the diagonal variant.

Proof. The basic building block of our ellipsoid matrix consists of the current and past stochastic gradients G t := [g 1 , g 2 , . . . , g t ].

We consider A rms which is built up as follows8 

A rms,t :=

EQUATION

)

From the construction of A rms,t it directly follows that for any unit length vector u ∈ R d \ {0}, u 2 = 1 we have

EQUATION

which proves the lower bound for ζ = 1/ . Now, let us consider the upper end of the spectrum of A rms,t . Towards this end, recall the geometric series expansion

EQUATION

and the fact that GG is a sum of exponentially weighted rank-one positive semi-definite matrices of the form g i g i . Thus

λ max (g i g i ) = Tr(g i g i ) = ∇g i 2 ≤ L 2

H , where the latter inequality holds per assumption for any sample size |S|. Combining these facts we get that

EQUATION

)

As a result we have that

≤ λ min (A rms,t ) ≤ λ max (A rms,t ) ≤ 1 -β t+1 L 2 H + (28)

Finally, to achieve uniform equivalence we need the r.h.s. of ( 28) to be bounded by 1/ . This gives rise to a quadratic equation in , namely

2 + 1 -β t+1 L 2 H -1 ≤ 0 (29)

which holds for any t and any β ∈ (0, 1) as long as

EQUATION

)

Such an always exists but one needs to choose smaller and smaller values as the upper bound on the gradient norm grows. For example, the usual value = 10 -8 is valid for all L 2 H < 9.9 • 10 7 . All of the above arguments naturally extend to the diagonal preconditioner diag(A rms ).

Second, we note that it is no necessary to compute the update step by minimizing Eq. ( 10) to global optimality. Instead, it suffices to do better than the Cauchy-and Eigenpoint simultaneously (Conn et al., 2000; Yao et al., 2018) . We here adapt this assumption for the case of iteration dependent norms (compare (Conn et al., 2000) 

EQUATION

where u t is an approximation to the corresponding negative curvature direction, i.e., for some

0 < ν < 1, u t H t u t ≤ ν ut ut t 2 λ min (B t ) and u t t = 1.

In practice, improving upon the Cauchy point is easily satisfied by any Krylov subspace method such as Conjugate Gradients, which ensures convergence to first order critical points. However, while the Steihaug-Toint CG solver can exploit negative curvature, it does not explicitly search for the most curved eigendirection and hence fails to guarantee m t (s t ) ≤ m t (s E t ). Thus more elaborate Krylov descent methods such as Lanczos method might have to be employed for second-order criticality (See also Appendix B.2 and Conn et al. (2000) Chapter 7).

We now restate two results from Conn et al. (2000) that precisely quantify the model decrease guaranteed by Assumption 2. Lemma 4 (Model decrease: Cauchy Point (Theorem 6.3.1. in Conn et al. (2000) )). Suppose that s C t is computed as in Eq. ( 32). Then

EQUATION

)

Lemma 5 (Model decrease: Eigenpoint (Theorem 6.6.1 in Conn et al. ( 2000))). Suppose that λ min (B t ) < 0 and s E t is computed as in Eq. ( 32). Then

EQUATION

)

We are now ready to prove the final convergence results. Towards this end, we closely follow the line of arguments developed in Yao et al. (2018) . First, we restate the following lemma which holds independent of the trust region constraint choice. Lemma 6 (Yao et al. (2018) ). Assume that L(w) is second-order smooth with Lipschitz constants L g and L H . Furthermore, let Assumption 1 hold. Then

EQUATION

)

Second, we show that any iterate of Algorithm 1 is eventually successful as long as either the gradient norm or the smallest eigenvalue are above (below) the critical values g and H .

Lemma 7 (Eventually successful iteration -g t ≥ g ). Assume that L(w) is second-order smooth with Lipschitz constants L g and L H . Furthermore, let Assumption 1 and 2 hold and suppose that g t ≥ g as well as

EQUATION

)

then the step s t is successful.

Therefore, recalling Eq. ( 34) as well as the fact that ut 2 ut t ≤ µ and s t 2 ≤ µ s t t ≤ µ∆ t due to Lemma 1

EQUATION

)

where the last second inequality is due to the conditions in Eq. ( 42). Therefore, ρ t ≥ η and the iteration is successful.

Together, these two results allow us to establish a lower bound on the trust region radius ∆ t . Lemma 9. Assume that L(w) is second-order smooth with Lipschitz constants L g and L H . Furthermore, let Assumption 1 and 2 hold. Suppose

δ g < 1 -η 4 g , δ H < min{ 1 -η 2 ν H , 1}.

then for Algorithm 1 we have

EQUATION

)

Proof. The proof follows directly from ∆ t ≥ ∆ t-1 /γ as well as the fact that any step is successful as soon as ∆ t falls below min g µ 1+Lg ,

(1-η) g

12L H µ 8 , (1-η) g 3µ 4 , (1-η) 2µ ν H L H

due to Lemma 7 and 8.

Lemma 10 (Number of successful iterations). Under the same setting as Lemma 9, the number of successful iterations taken by Algorithm 1 is upper bounded by

|T succ | ≤ L(x 0 ) -L(x * ) C H min{ 2 g , 2 H } , where C := η min {C 1 , C 2 }, C 1 := 1 2 min 1 1+Lg , C g , C 2 := νµ 2 2 min C 2 g , C 2 H ,C g := min g µ 1+Lg , (1-η) g 12L H µ 8 , (1-η) g 3µ 4 , C H := (1-η) 2µ ν H L H

Proof. Suppose Algorithm 1 does not terminate at iteration t. Then either g t ≥ g or λ min (B) ≤ h . If g t ≥ g , according to (33) and Lemma 1, we have

-m t (s t ) ≥ 1 2 g t min{ g t 1 + H t , ∆ t 1 µ } ≥ 1 2 g min{ g 1 + L g , C g g , C H H } ≥ C 1 g min{ g , H }.

Similarly, in the second case λ min (B t ) ≤h , from Lemma 1 and 5 we have

-m t (s t ) ≥ 1 2 ν|λ min (B t )|∆ 2 t µ 2 ≥ C 2 H min{ 2 g , 2 H }.

Let T succ denote the number of successful iterations. Since L(w) is monotonically decreasing, we have

L(w 0 ) -L(w * ) ≥ t=0 L(w t ) -L(w t+1 ) ≥ t∈Tsucc L(w t ) -L(w t+1 ) ≥ t∈Tsucc -m t (s t )η ≥ t∈Tsucc C H min{ 2 g , 2 H } ≥ |T succ |C H min{ 2 g , 2

H }, which proves the assertion.

We are now ready to prove the final result. Particularly, given the lower bound on ∆ t established in Lemma 9 we find an upper bound on the number of un-successful iterations, which combined with the result of Lemma 10 on the number of successful iterations yields the total iteration complexity of Algorithm 1.

Theorem 4 (Theorem 2 restated). Assume that L(w) is second-order smooth with Lipschitz constants L g and L H . Furthermore, let Assumption 1 and 2 hold. Then Algorithm 1 finds an O( g , H ) first-and second-order stationary point in at most O max -2 g -1

H , -3 H iterations.

Proof. The result follows by combining the lemmas 9 and 10 as in Theorem 1 of Xu et al. (2017a) .

D DIAGONAL DOMINANCE IN NEURAL NETWORKS:

In the following, we make statements about the diagonal share of random matrices. As E[ 1 x ] might not exist for a random variable x, we cannot compute the expectation of the diagonal share but rather of for computing the diagonal share of the expectation of the random matrix in absolute terms. Note that this notion is still meaningful, as the average of many non-diagonally dominated matrices with positive entries cannot become diagonally dominated.

D.1 PROOF OF PROPOSITION 1:

Proposition 3 (Proposition 1 restated). For random Gaussian Wigner matrix W formed as

EQUATION

where ∼ stands for i.i.d. draws (Wigner, 1993) , the diagonal mass of the expected absolute matrix amounts to

EQUATION

Proof.

δ E,W = d k=1 E [|W k,k |] d k=1 d l=1 E [|W k,l |] = dE [|W 1,1 |] dE [|W 1,1 |] + d(d -1)E [|W 1,2 |] = dσ 1 2/π dσ 1 2/π + d(d -1)σ 2 2/π = 1 1 + d(d-1)σ2 √ 2/π dσ1 √ 2/π = 1 1 + (d -1) σ2 σ1 (47)

which simplifies to 1 d if the diagonal and off-diagonal elements come from the same Gaussian distribution (σ 1 = σ 2 ).

For the sake of simplicity we only consider Gaussian Wigner matrices but the above argument naturally extends to any distribution with positive expected absolute values, i.e. we only exclude the Dirac delta function as probability density. 

D.2 OLS BASELINE:

When considering regression tasks, a direct competitor to neural network models is the classical Ordinary Least Squares (OLS) regression, which minimizes a quadratic loss over a linear model. In this case the Hessian simply amounts to the input-covariance matrix H ols := X X, where X ∈ R d×n . We here show that the diagonal share of the expected matrix itself also decays in d, when n grows to infinity. However, empirical simulations suggest the validity of this result even for much smaller values of n (see Figure D .2) and it is likely that finite n results can be derived when adding assumptions such as Gaussian data. Proposition 4 (Proposition 2 restated). Let X ∈ R n×d and assume each x i,k is generated i.i.d. with zero-mean finite second moment σ 2 > 0. Then the share of diagonal mass of the expected matrix

EQUATION

Where we used the fact that all x i,k are i.i.d. variables. Per assumption, we have E x 2 i,1 = σ 2 , ∀i. Furthermore, the products x i,1 x i,2 are i.i.d with expectation 0 and variance σ 4 . By the central limit theorem

Z N = 1 √ n n i=1 x i,1 x i,2 → Z in law, with Z ∼ N (0, σ 4 ), since E Z 2 n = 1 n E ( n i=1 x i,1 x i,2 2 = 1 n n i=1 E[x 2 i,1 ]E[x 2 i,2 ] = σ 4 due to the independence assumption. Then E |Z N | 1 |Z N |≥R ≤ E(|Z N | 2 /R) ≤ σ 4 /R. This implies that E (|Z N |) → E(|Z|) = 2 π σ 2

As a result, we have that in the limit of large n

EQUATION

B TRUST REGION METHODS

B.1 OUTER ITERATIONS:

Trust region methods are among the most principled approaches to overcome the above mentioned issues. These methods also construct a quadratic model m t but constrain the subproblem in such a way that the stepsize is restricted to stay within a certain radius ∆ t within which the model is trusted to be sufficiently adequate

EQUATION

Hence, contrary to line-search methods this approach finds the step s t and its length s t simultaneously by optimizing (53). Subsequently the actual decrease L(w t ) -L(w t + s t ) is compared to the predicted decrease m t (0) -m t (s t ) and the step is only accepted if the ratio ρ := L(w t ) -L(w t + s t )/(m t (0) -m t (s t )) exceeds some predefined success threshold η 1 > 0.

Furthermore, the trust region radius is decreased whenever ρ falls below η 1 and it is increased whenever ρ exceeds the "very successful" threshold η 2 0. Thereby, the algorithm adaptively measures the accuracy of the second-order Taylor model -which may change drastically over the parameter space depending on the behaviour of the higher-order derivatives9 -and adapts the effective length along which the model is trusted accordingly. See Conn et al. (2000) for more details.

As a consequence, the plain Newton step s N,t = -∇ 2 L t -1 ∇L t is only taken if it lies within the trust region radius and yields a certain amount of decrease in the objective value. Since many functions look somehow quadratic close to a minimizer the radius can be shown to grow asymptotically under mild assumptions such that eventually full Newton steps are taken in every iteration which retains the local quadratic convergence rate (Conn et al., 2000) . 

B.2 SUBPROBLEM SOLVER:

Interestingly, there is no need to optimize Eq. ( 53) to global optimality to retain the remarkable global convergence properties of TR algorithms. Instead, it suffices to do better than the Cauchy-and Eigenpoint10 simultaneously. One popular approach is to minimize m t (s) in nested Krylov subspaces. These subspaces naturally include the gradient direction as well as increasingly accurate estimates of the leading eigendirection span{g t , B t g t , B 2 t g t , . . . , B j t g t } (54) until (for example) the stopping criterion

∇m t (s j ) ≤ ∇L(w t ) min{κ K , ∇L(w t ) θ }, κ K < 1, θ ≥ 0 (55)

is met, which requires increased accuracy as the underlying trust region algorithm approaches criticality. Conjugate gradients and Lanczos method are two iterative routines that implicitly build up a conjugate and orthogonal basis for such a Krylov space respectively and they converge linearly on quadratic objectives with a square-root dependency on the condition number of the Hessian (Conn et al., 2000) . We here employ the preconditionied Steihaug-Toint CG method (Steihaug, 1983 ) in order to cope with possible boundary solutions of ( 53) but similar techniques exist for the Lanczos solver as well for which we also provide code. As preconditioning matrix for CG we use the same matrix as for the ellipsoidal constraint.

C DAMPED (GAUSS-)NEWTON METHODS An alternative approach to actively constraining the region within which the model is trusted is to instead penalize the step norm in each iteration in a Lagrangian manner. This is done by so-called damped Newton methods that add a λ > 0 multiple of the identity matrix to the second-order term in the model, which leads to the update step

EQUATION

)

This can also be solved hessian-free by conjugate gradients (or other Krylov subspace methods). The penalty parameter λ is acting inversely to the trust region radius ∆ and it is often updated accordingly. Such algorithms are commonly known as Levenberg-Marquardt algorithms and they were originally tailored towards solving non-linear least squares problems (Nocedal & Wright, 2006) but they have been proposed for neural network training already early on (Hagan & Menhaj, 1994) .

Many algorithms in the existing literature replace the use of ∇ 2 L(w t ) in ( 56) with the Generalized Gauss Newton matrix (Martens, 2010; Chapelle & Erhan, 2011) or an approximation of the latter (Martens & Grosse, 2015) . This matrix constitutes the first part of the well-known Gauss-Newton decomposition

EQUATION

where l and l are the first and second derivative of l : R out → R + assuming that out = 1 (binary classification and regression task) for simplicity here.

It is interesting to note that the GGN matrix A GGN of neural networks is equivalent to the Fisher matrix used in natural gradient descent (Amari, 1998) in many cases like linear activation function and squared error as well as sigmoid and cross-entropy or softmax and negative log-likelihood for which the extended Gauss-Newton is defined (Pascanu & Bengio, 2013) . As can be seen in ( 57) the matrix A GGN is positive semidefinite (and low rank if n < d). As a result, there exist no second-order convergence guarantees for such methods on general non-convex problems. On the other end of the spectrum, the GGN also drops possibly positive terms from the Hessian (see 57). Hence it is not guaranteed to be an upper bound on the latter in the PSD sense. Essentially, GGN approximations assume that the network is piece-wise linear and thus the GGN and Hessian matrices only coincide in the case of linear and ReLU activations or non-curved loss functions. For any other activation the GGN matrix may approximate the Hessian only asymptotically and if the (f i (•)) terms in 57 go to zero for all i ∈ {1, . . . , n}. In non-linear least squares such problems are called zero-residual problems and GN methods can be shown to have quadratic local convergence there. In any other case the convergence rate does not exceed the linear local convergence bound of gradient descent. In practice however there are cases where deep neural nets do show negative curvature in the neighborhood of a minimizer (Bottou et al., 2018 ).Finally, Dauphin et al. (2014) propose the use of the absolute Hessian instead of the GGN matrix in a framework similar to 56. This method has been termed saddle-free Newton even though its manifold of attraction to a given saddle is non-empty11 . Figure 9 : Both, the GGN method and saddle-free Newton method make a positive definite quadratic model around the current iterate and thereby overcome the abstractedness of pure Newton towards the saddle (compare Figure 8 ). However, (i) none of these methods can escape the saddle once they are in the gradient manifold of attraction and (ii) as reported in Mizutani & Dreyfus (2008) the GN matrix can be significantly less well conditioned than the absolute Hessian (here κGN = 49 487 554 and κ |H| = 1.03 so we had to add a damping factor of λ = 0.1 to make the GN step fit the plot.

C.1 COMPARISON TO TRUST REGION:

Contrary to TR methods, the Levenberg-Marquardt methods never take plain Newton steps since the regularization is always on (λ > 0). Furthermore, if a positive-definite Hessian approximation like the Generalized Gauss Newton matrix is used, this algorithm is not capable of exploiting negative curvature and there are cases in neural network training where the Hessian is much better conditioned than the Gauss-Newton matrix (Mizutani & Dreyfus, 2008 ) (also see Figure 9 ). While some scholars believe that positive-definiteness is a desirable feature (Martens, 2010; Chapelle & Erhan, 2011) , we want to point out that following negative curvature directions is necessarily needed to escape saddle points and it can also be meaningful to follow directions of negative eigenvalue λ outside a saddle since they guarantee O(|λ| 3 ) progress, whereas a gradient descent step yields at least ∇f (w) 2 progress (both under certain stepsize conditions) and one cannot conclude a-priori which one is better in general (Curtis & Robinson, 2017; Alain et al., 2018) . Despite these theoretical considerations, many methods based on GGN matrices have been applied to neural network training (see Martens (2014) and references therein) and particularly the hessian-free implementations of (Martens, 2010; Chapelle & Erhan, 2011 ) can be implemented very cheaply (Schraudolph, 2002) .

D USING HESSIAN INFORMATION IN NEURAL NETWORKS:

While many theoretical arguments suggest the superiority of regularized Newton methods over gradient based algorithms, several practical considerations cast doubt on this theoretical superiority when it comes to neural network training. Answers to the following questions are particularly unclear: Are saddles even an issue in deep learning? Is superlinear local convergence a desirable feature in machine learning applications (test error)? Are second-order methods more "vulnerable" to sub-sampling noise? Do worst-case iteration complexities even matter in real-world settings? As a result, the value of Hessian information in neural network training is somewhat unclear a-priori and so far a conclusive empirical study is still missing.

Our empirical findings indicate that the net value of Hessian information for neural network training is indeed somewhat limited for mainly three reasons: 1) second-order methods rarely yield better limit points, which suggests that saddles and spurious local minima are not a major obstacle; 2) gradient methods can indeed run on smaller batch sizes which is beneficial in terms of epoch and when memory is limited; 3) The per-iteration time complexity is noticeably lower for first-order methods. In summary, these observations suggest that advances in hardware and distributed second-order algorithms (e.g., Osawa et al. (2018) ; Dünner et al. (2018) ) will be needed before Newton-type methods can replace (stochastic) gradient methods in deep learning. The Fashion-MNIST autoencoder has the same architecture as the one used in Hinton & Salakhutdinov (2006) ; Xu et al. (2017a) ; Martens (2010) ; Martens & Grosse (2015) . The encoder structure is 784 -1000 -500 -250 -30 and the decoder is mirrored. Sigmoid activations are used in all but the central layer. The reconstructed images are fed pixelwise into a binary cross entropy loss.

The network has a total of 2 833 000 parameters. The CIFAR-10 autoencoder is taken from the implementation of https://github.com/jellycsc/PyTorch-CIFAR-10-autoencoder.

For the ResNet18, we used the implementation from torchvision for CIFAR-10 as well as a modification of it for Fashion-MNIST that adapts the first convolution to account for the single input channel.

In all of our experiments each method was run on one Tesla P100 GPU using the PyTorch 

:

This algorithm is a simplified TR method, initially tailored for non-linear least squares problems(Nocedal & Wright,

2006)

For illustrative purposes, we only plot every other trust region.

The argument naturally extends to any distribution with positive expected absolute values.

In the following we only treat TR methods, but we would like to emphasize that the use of matrix induced norms can directly be transferred to the cubic regularization framework.

Note that the assumption of bounded batch gradients, i.e. smooth objectives, is common in the analysis of stochastic algorithms(Allen-Zhu, 2017;Defazio et al., 2014;Schmidt et al., 2017;Duchi et al., 2011).

We observed weaker performance when running with smaller batch size. We hypothesize that second-order methods extract more information of each batch and are thus likely to "overfit" small batches in each step.

Essentially Algorithm 1 with mt based on a first order Taylor expansion, i.e. m 1 t (s) as in Eq. (14).

This is a generalization of the diagonal variant proposed byTieleman & Hinton (2012), which precondition the gradient step by an elementwise division with the square-root of the following estimate gt = (1 -β)gt-1 + β∇L(wt) 2 .

Note that the second-order Taylor models assume constant curvature.

which are the model minimizers along the gradient and the eigendirection associated with its smallest eigenvalue, respectively.

It is the same as that for GD, which renders the method unable to escape e.g. when initialized right on a saddle point. To be fair, the manifold of attraction for GD constitutes a measure zero set(Lee et al., 2016).

Both datasets were accessed from https://www.tensorflow.org/api_docs/python/tf/keras/datasets

Figures and Tables:

Figure 1: Top: Iterates and implicit trust regions of GD and Adagrad on three quadratic objectives with different condition number κ. Bottom: Average log suboptimality over iterations as well as 90% confidence intervals of 30 runs with random initialization

Figure 2: Diagonal mass of neural network Hessian δ H relative to δ E[|W|] and δ E[|H ols |] of corresponding dimensionality for random inputs as well as at random initialization, middle and end of training with RMSProp on CIFAR-10. Mean and 95% CI over 10 independent runs.

Figure 3: Log loss over backpropagations. Mean and 95% CI of 10 runs. Green dotted line indicates 99% acc.

Figure 5: Iterates (left) and log suboptimality (right) of GD, Adagrad and two full-featured first-order TR algorithms of which one (1st TR) is spherically constraint and the other (1st TRada) uses A ada as ellispoid.

Figure 6: Share of diagonal mass of the Hessian δ H relative to δ W of the corresponding Wigner matrix at random initialization, after 50% iterations and at the end of training with RMSprop on MNIST. Average and 95% confidence interval over 10 runs. See Figure 2 for CIFAR-10 results.

Figure 7: Validity of Proposition 2 in the small n regime for Gaussian data.

Figure8: Level sets of the non-convex, coercive objective function f (w) = 0.5w 2 0 + 0.25w 4 1 -0.5w 2 1 . Newton's Method makes a local quadratic model (blue dashed lines) and steps to its critical point. It may be thus be ascending (left) or attracted by a saddle point (right). TR methods relieve this issue by stepping to the minimizer of that model within a certain region (green dashed line).

Figure 10: Experiment comparing TR and gradient methods in terms of epochs. Average log loss as well as 95% confidence interval shown.

Figure 12: Original and reconstructed MNIST digits (left), Fashion-MNIST items (middle), and CIFAR-10 classes (right) for different optimization methods after convergence.

Chapter 6). restate this assumption here Assumption 2 (Approximate model minimization). Each update step s t yields at least as much model decrease as the Cauchy-and Eigenpoint simultaneously, i.e.

Table content: <table><tr><td colspan="4">m t (s t ) ≤ m t (s C t ) and m t (s t ) ≤ m t (s E t ),</td><td>(31)</td></tr><tr><td>where</td><td/><td/><td/></tr><tr><td>s C t := arg min 0≤α≤∆t</td><td>m t (-α</td><td>g t g t t</td><td>) and s E t := arg min</td></tr></table>

Parameters Table1reports the default parameters we consider. Only for the larger ResNet18 on CIFAR-10, we adapted the batch size to 128 due to memory constraints. Default parameters Datasets We use two real-world datasets for image classification, namely CIFAR-10 and Fashion-MNIST 12 . While Fashion-MNIST consists of greyscale 28 × 28 images, CIFAR-10 are colored images of size 32 × 32. Both datasets have a fixed training-test split consisting of 60,000 and 10,000 images, respectively.Network architecturesThe MLP architectures are simple. For MNIST and Fashion-MNIST we use a 784 -128 -10 network with tanh activations and a cross entropy loss. The networks has 101 770 parameters. For the CIFAR-10 MLP we use a 3072 -128 -128 -10 architecture also with tanh activations and cross entropy loss. This network has 410 880 parameters.

Table content: <table><tr><td>|S 0 | ∆ 0</td><td colspan="2">∆ max η 1</td><td>η 2</td><td>γ 1</td><td>γ 2</td><td>κ K (krylov tol.)</td></tr><tr><td colspan="2">TR uni 512 10 -4 10</td><td colspan="4">10 -4 0.95 1.1 1.5</td><td>0.1</td></tr><tr><td colspan="2">TR ada 512 10 -4 10</td><td colspan="4">10 -4 0.95 1.1 1.5</td><td>0.1</td></tr><tr><td colspan="2">TR rms 512 10 -4 10</td><td colspan="4">10 -4 0.95 1.1 1.75</td><td>0.1</td></tr></table>

Peer review (iclr2020):
In this paper, the authors investigate the use of ellipsoidal trust region constraints for second order optimization. The authors first show that adaptive gradient methods can be viewed as first-order trust region methods with ellipsoid constraints. The authors then show that the preconditioning matrix of RMSProp and Adam can be used as norm inducing matrices for second order trust region methods. This ellipsoidal trust region method is empirically compared with first order gradient methods and spherical second order trust region methods. Overall the paper is nicely written and very easy to read. The ideas are interesting. However, I have a number of concerns/questions about the work, that I list below. 1. Why is the preconditioning matrix of RMSProp/Adam a reasonable norm inducing matrix? One can show that the empirical Fisher is not an accurate curvature matrix in general, and so there is no reason to believe this would in fact enforce the proper ellipsoidal trust region for the method? See for example: https://arxiv.org/pdf/1905.12558.pdf. 2. I am also not convinced that Figure 2 actually shows that the curvature matrix is diagonally dominant. How do I interpret a value of 40 or 50 for this metric, and why does it imply that it is diagonally dominant? 3. The experiments also do not look very convincing to me. How sensitive is the algorithm to the hyperparameters like lambda1 and lambda2? I am also a bit confused about why different batch sizes were used for the first order gradient methods and the second order TR methods? The method overall doesn*t seem to be able to match first order gradient methods, and it is not clear whether this is because of using the RMSProp/Adam preconditioner as a curvature matrix. Given these concerns, I consider this paper to be borderline. I am happy to have a discussion with the authors and the other reviewers and change my score however. ========================================= Edit after rebuttal: I thank the reviewers for their response. While the updated paper has certainly improved, I think the paper still requires a much more thorough experimental evaluation. I am sticking to my score.