Surprising Instabilities in Training Deep Networks and a Theoretical Analysis 
Authors: Yuxin Sun
Dong Lao
Ganesh Sundaramoorthi
Anthony Yezzi
Venue: NeurIPS 2022 Accept (2022)
Keywords: Deep Learning Theory， Stability Analysis

Abstract: We empirically demonstrate numerical instabilities in training standard deep networks with SGD. Specifically, we show numerical error (on the order of the smallest floating point bit) induced from floating point arithmetic in training deep nets can be amplified significantly and result in significant test accuracy variance, comparable to the test accuracy variance due to stochasticity in SGD. We show how this is likely traced to instabilities of the optimization dynamics that are localized over iterations and regions of the weight tensor space. We do this by presenting a theoretical framework using numerical analysis of partial differential equations (PDE), and analyzing the gradient descent PDE of a one-layer convolutional neural network, which is sufficient to illustrate these instabilities. We show that it is stable only under certain conditions on the learning rate and weight decay. We reproduce the localized instabilities in the PDE for the one-layer network, which arise when the conditions are violated.

Surprising Instabilities in Training Deep Networks and a Theoretical Analysis
Yuxin Sun! Dong Lao² Ganesh Sundaramoorth Anthony Yezzi!
GGeorgia Institute of Technology, RUCLA, *Raytheon Technologies
{syuxin3, syezzi)@gatech. edu, lao8cs.ucla ucla. edu, ganesh. sundaramoorthi@rtt..
Abstract
We discover restrained numerical instabilities in current training practices of deep networks with stochastic gradient descent (SGD). We show numerical eπror (on the order of the smallest floating point bit) induced from floating point arithmetic in training deep nets can be amplified significantly and result in significant test accuracy variance, comparable to the test accuracy variance due stochasticity in SGD. We show how this is likely traced instabilities of the optimization dynamics that are restrained, i.e., localized over iterations and regions of the weight tensor space. We do this by presenting theoretical framework using numerical analysis of partial differential equations (PDE), and analyzing the gradient descent PDE of a simplified convolutional neural network (CNN). We show that it is stable only under certain conditions on the learning rate and weight decay. We reproduce the localized instabilities in the PDE for the simplified network, which arise when the conditions are violated.
Introduction
Deep learning is now standard practice in number of application areas including computer vision and natural language processing. However, the practice has out-paced the theory. Theoretical advances in deep learning could improve practice, and could also speed the adoption the technology in e.g., a number of safety critical application areas.
In this paper, we discover phenomena of restrained instabilities in current deep learning training optimization. Numerical instabilities are present but are restrained so that there is not full global divergence, these restrained instabilities still cause tiny numerical finite precision errors to amplify, which leads to significant variance in the test accuracy. We introduce theoretical framework using tools from numerical partial diffcrential equations (PDE) to provide insight into the phenomena and analytically show how such restrained instabilities can arise decay for simplified CNN model. This study is a step toward a theory for principally choosing learning rates in deep network training from the perspective of numerical conditioning.
function of learning rate and weight
In deep learning practice, is well known that large learning rates can cause divergence and small learning rates can be slow and cause the optimization to be stuck at high training error (T), (2). Learning rates and schedules are often chosen heuristically to address this trade-off. On the other hand, in numerical PDEs, there is mature theory for choosing the learning rate (3] or step size for discrete schemes for solving PDEs. There are known conditions on the step size (e.g., the Courant- Friedrichs-Lent (CFL) condition (4) ensure that the resulting scheme is stable and converges. These methods would thus be natural to study learning rate schemes for training deep networks with SGD, which in the continuum limit is PDE. However, to the best of our knowledge. numerical PDEs have not been applied to this domain. This might be for two reasons. First, the complicated structure of deep networks with many layers and non-linearities makes it difficult to apply these analytic techniques. Second, the stochastic element of SGD complicates applying these techniques,
36th Conference on Neural Information Processing Systems (NcurIPS 2022).

---

which primarily do not consider stochastic elements. While we do not yet claim to overcome these issues, we do show the relevance of PDE analysis in analyzing training of deep networks with SGD, which gives promise for further developing this area.
Our partieular contributions 1. We found convincing evidence of and thus discovered the presence of restrained numerical instabilities in standard training practices for deep neural networks. 2. We show that this results in inherent floating points arithmetic errors having divergence effects on the training optimization that are even comparable to the divergence due to stochastic variability in SGD. 3. We present and study a simplified model, which allows us to exploit numerical PDE analysis, in order to offer an explanation on how these instabilities might arise.
These insights suggest us that the learning rate should be chosen adaptively over localized space- time regions of the weight tensor space in order to ensure stability while not compromising speed. While we understand this is complex task, without a simple solution, our study motivates the need for investigating such challenging strategies. Note there are a number of heuristically driven methods for adaptive learning rates (e.g.. (5)-(9), but they do not address numerical stability.
2 Related Work
As stated earlier, we are not aware of related work in applying numerical PDE analysis to study stability and choice of learning rates of deep network training There has been work on studying convergence of SGD and learning rates (e-g.. [TO]- (T2), but they are often based on convexity and Lipschitz assumptions that could be difficult to apply to deep networks. For deep networks, there is
work on global stability training in relation to curvature of the loss, e.g., (13). [3 is more focused
on generalization rather than numerical stability and does not address restrained instabilites.
More broadly, although not directly related to our work, there has been number of recent works in connecting PDEs with neural networks for various end goals (see (14) (16) for detailed surveys).
For example, applying neural networks to solve PDEs - 21 J and interpreting the forward
pass of neural networks as approximations to PDEs [22 (24) Since SGD can be interpreted as
a discretization of a PDE, PDEs also play a role developing optimization schemes for neural networks 25 31.
proposed to reduce variance of final accuracy induced by the stochastic gradient, we find quite surprisingly that even if all random seeds in optimization are fixed, variance induced by epsilon-small numerical noise (perturbing only the last floating point bit) and hence floating point arithmetic is nevertheless greater or similar to variance across independent trials (see Table ID. This observation potentially indicates that numerical stability is a critical but long-ignored factor affecting neural network optimization.
Despite various methods [34-334] 34
There has been work training deep networks using limited numerical precision (e.g.. half precision) (35), (36). The purpose of these works is to introduce schemes to retain accuracy, which is well known to be reduced significantly when lowering precision these levels. Different from this literature, show an unknown fact: errors just in the last bit of floating point representations result in significant divergence. This is to provide evidence of numerical instability resulting from learning rate choices.
3 Background: PDEs, Discretization, and Numerical Stability
We start by providing background on analyzing stability of linear PDE discretizations. This method- ology will be applied to neural networks in later sections as a PDE is the continuum limit of SGD. We illustrate concepts using basic PDE, i.e., the heat equation (used to model a diffusion process):
Ozu(t, x) = dra u(t, w(0, x) = m0(x)),
(1)
where u T) x R + RR, O.u denotes the partial with respect to time I, Onn' is the second derivative with respect to the spatial dimension x, and R - RR is the initial condition. To solve
'While this first step to study deep networks using more matured tools for the numerical analysis of discretized PDE's is focused on lllayer CNN, we show in the appendix how this starting point still yields useful insight into the more complicated dynamics of multi-layer networks, especially when highly regularized, thereby offering strong promise for further developing this theory.
"In deep learning. "stability" often refers to vanishing/exploding gradients, network dynamics, or model robustness. These topics are different from numerical stability of the training optimization studied in this paper.

---

244
2.4
1.5
1.5
200
1.0
800
1.0
0.5
30
10 15 20
Figure 1: Illustration of instability in discretizing the heat equation. The initial condition ug is triangle (blue), with boundary conditions u(0) u(30) 0. [Left]] When the CFL condition met, i.e., Az/ 0.4 < 5. the method is stable and approximates the solution of the PDE. Note the true steady state is O, which matches the plot. [Right]: When the CFL condition is not met, i.e., 1 small numerical errors are amplified and the scheme diverges.
At/(Az)2 = 0.8 2
this numerically, discretizes the equation, using finite difference approximations, A standard discretization is through ut, x) ult.a+
Ar] - (((Ar)) and 0,u(t, =
(Aa)²
Att+ + O(A),, where Ax is the spatial increment and At is the step size, which gives At the following update scheme:
At
u1 nー (x) = (x)
+
[u" (x + Ar) - 2u" v= Ar)] + (x),
(2)
Arr)
where denotes the iteration number, u² is the approximation of u at nAt, and denotes the eтror of the approximation (i.e., due to discretization and finite precision arithmetic).
A key question is whether u" remains bounded as n co, i.e., numerically stable. It is typically easier to understand stability in the frequency domain, which is referred to as Von Neumann analysis 3 [37]. Computing the spatial Discrete Fourier Transform (DFT) yields:
n-1
(w) Alw=-" (w) + ("(W,, or (w(w) = (w)(ow) + Σ A'(w)e"
(3)
where the hat denotes DFT, denotes frequency, and A(w) = 1
Ar [1m cos(war)) is the
amplifier function. Note that u" is stable if and only if [A(w))
Notice that if A(w)) 1, u²
converges and errors e" are attenuated over iterations. [A(w]] 2 1, u" diverges and the errors e are amplified. The case, [A(w)] This is known as the CFL condition. Notice the restriction on the time step (learning rate) for numerical stability. If the discretization error of the operator on the right hand side of the PDE is less than 0(Ar), stability also implies convergence to the PDE solution as At + 3], 37
for all w, implies the conditions that At 0 and AV (Ax)²²
See Figure I for simulation of the discretization scheme for the heat PDE. The instability causes blowup the solution globally. For standard optimization of CNNs, we will show rather that the instability is restrained (localized in time and space), which nevertheless causes divergence.
4 Empirical Evidence of Instability in Training Deep Networks
We discover and provide empirical evidence of restrained instabilities in current deep lcarning training practice. To this end, show that optimization paths in current practice of training convolutional neural networks (CNNs) can diverge significantly due to the smallest errors from finite precision arithmetic. We show that the divergence can be eliminated with learning rate choice.
4.1 A Perturbed SGD by Introducing the Smallest Floating Point Perturbations
We introduce a modified version of SGD that introduces perturbations the original SGD update on the order of the smallest possible machine representable perturbation, modeling noise on the order of

---

error due to finite precision floating point arithmetic:
D.+1 02 - 179t
(SGD)
D.t1 = (w/kk) x (kgg))
(PERTURBED SGD)
where v) >> 0 the learning rate, 9t is the stochastic gradient or the momentum vector if momentum is used, and k > 1 is an odd integer. Notice that mathematically, both updates are exactly equivalent (division by cancels the multiplication by k), however considering floating point arithmetic, these may not be equivalent. A floating point number is represented
significand X base expoeen (s bits significand, e bits exponent)
(4)
where typically base 2 (binary) is used. Thus, the modified update introduces potential perturbation
in the last significant bit of 1/9r. hence a relative difference of the right hand sides on the order of 2- according IEEE standards on floating point computation. In other words, Relative y) := |ar /l/|=|| < 2 2- where x f((()))) y f((n//k) /k) fl(kgt), and fi(x)) is the floating point representation of æ. Note that this is the same order of relative error between 1/9; and fl(1g.)) ) the inherent machine error due to floating point representation. The perturbation induces the smallest possible perturbation of 17: representable by the machine. Thus, perturbed SGD is way to show the effects of machine error. Notice that when is a power of two, as floating point numbers are typically stored in binary form, the division just subtracts the exponent and the multiplication adds to the exponent so both updates are equivalent and no perturbation is introduced. However, choosing odd can create changes at the last significant bit, with each choice of k yielding a different change.
4.2 Perturbed SGD on Common CNNs and Demonstration of Restrained Instabilities
In SGD, there are several sources of randomness, i.e., random initialization, random shuffled batches/selection random data augmentation, and there are also sources of non-determinism in implementations of the deep learning frameworks. To ensure the variance merely originated from the floating point perturbation, in the following experiments, we eliminate all randomness by fixing the initialization, the random seed for batch selection, and the non non-determinism flag. Appendix empircally validates that all seeds are fixed. We use Pytorch using default 32 bit floating point precision (similar conclusions for 64 bit hold).
Perturbed SGD with ReLU: Our first experiment uses the ResNet-56 architecture [38] which we train on CIFAR-10 (39) using perturbed SGD. We use the standard parameters for training this network 38 : 1r=0.1 1, batch size 128, weight decay
5e-4, momentum=0.9 9. The standard step decay learning rate schedule is used: the learning rate is divided by 10 every 40 epochs for total of 200 epochs. We report the final test accuracies for 6 different values of kt including k = (standard SGD), and compute the standard deviation of the accuracies (STD). Table [ (left) shows the results. This results test accuracy variability over different k. We repeat this experiment over different seeds for the stochastic batch selection (shown in the next columns) fixing the initialization. The standard deviation for each seed over k on average 0.16%. To illustrate the significance of this value, we compare to the variability due to batch selection on the right of Table We empirically estimate the relative fluctuation of the gradient estimate from batch selection, which is 26. 72.72. This large compared to 2 23 for the floating point perturbation, yet this smallest fluctuation that is machine representable (and models finite precision error) results in the about the same accuracy variance as batch selection. Thus, finite precision error, inherent in the system, surprisingly results in about the same test accuracy variance as stochastic batch selection.
We first conjectured that the variance due to floating point errors could be caused by the ReLU activation, the activation standard in ResNets. This because the gradient of the ReLU is discontinuous at zero. Hence small positive values (c.g., 1c-6) would give a gradient of I and small negative values (e.g., -1e-6) would give gradient of zero. Hence the ReLU gradient is sensitive to numerical fluctuation around zero. Note common deep learning implementations choose a value zero, among O. 0.5 or 1: but this has the same sensitivity. Therefore, we adopted the Swish activation function (proposed in (40). Swish(x) := x [1 + exp (-Bz)] (3 > 0. as 3 - oc Swish approaches ReLUU,, which is an approximation of the ReLU that has a continuous gradient at zero. Recent work has shown the Swish function improves accuracy [40] as well. Thus, we repeated the previous experiment (under the same settings) replacing the ReL. with the Swish activation. Results are shown in Appendix, and show similar results as the previous experiment with ReLU: the
Perturbed SGD with Swish:

---

Saup
2
5.
STD
93.36
93.40
93.10
93.14
93.34
93.33
0.11
3
93.49
93.37
93.08
93.68
93.16
93.12
0.22
==
5
93.64
93.22
93.39
93.17
93.26
93.42
0.16
METiON
ERaOR
Sm
=== 7
93.36
93.31
93:14
93.28
0.09
FLOATING PT
2
0.16
93.87
93.55
93.42
BATCH SELECTION
26.73
0.17
03.31
93.14
0.21
Table 1: [Left]] Final test accuracy different seeds (batch selections) and different floating point perturbations (rows) for Resnet56 trained on CIFAR-10. STD is standard deviation. [Right]: Comparison of errors induced in the gradient and resulting perturbation vs batch selection. Floating point error perturbs the gradient by an amount 8 orders of magnitude smaller than batch selection yet surprisingly yields nearly the same test accuracy variation!
accuracy STD for floating point
standard deviation test accuracy due to floating point perturbation (93.72 0.15) is significant and is comparable (even exceeds) the variance due to stochastic batch selection (0.09).
Other Architecturee/Datassts To verify that this phenomenon not specific the choice of network architecture or dataset, we repeated the same experiment for a different network (VGGI6
41 and different dataset (Pashion-MNISTT (42). We summarize the results in the Appendix, which confirms that the phenomena is still present across different networks and datasets.
Divergence in Optimization Paths: We plot the difference in network weights (using the measure described next) between the original SGD weights, 0; (i is the index for location in the vector), and the perturbed SGD weights, 0f over epochs and show that the errors amplified, suggesting instability. We use the average relative I.1 absolute difference between weights: RelL. (0, 0') == N E. 0, ( + (0]1). where N is the number of network parameters. This measure used as the numerical eтror is multiplicative and relative to the weight. Note that RelL is bounded by 2. Results are shown for multiple networks in Figure (left). The difference rapidly grows, indicating errors are amplified, and then the growth is restramed at larger epochs. Although the errors are restrained, the initial error growth is enough to result different test accuracies.
Average Relative L1 Difference
Average Relstive L1 Difference
1.5
S
2.5e-2
1.25e
6 6.25a-3
0.08
11125-3
125e-4
0.04
Keenethil Citar10 Resnetté FMNIST
Vec1 Cilar10
0.02
0.5
Vag IMNIST
50
100 Epoch
150
200
0
20
40
60
80
100
Scaled lteration
Figure 2: [Left]]: Relative LI difference in weights across epochs for SGD (k = 1) and modified
SGD (k = 3). A decayed learning rate schedule used. The errors quickly build up, but then are contained at higher epochs. The initial build up of errors is enough to result in different test accuracy. [Right]: Relative L1 difference in weights over initial iterations for SGD
and modified SGD (k = 3) with various fixed learning rates. Lower than some minimum learning rate, the floating point perturbation gets attenuated and higher than that value errors are amplified, suggesting an instability.
Evidence Divergence is Due to Instability: We provide evidence that this divergence phenomenon due an instability by showing that decreasing the learning rate small enough value can eliminate the divergence, even initially. We experiment with various fixed learning rates. We use Resnet56 under the same settings as before. For smaller learning rates, more iterations would have to be run to move equivalent amount for larger rate. We thus choose the iterations for each run such that the learning rate times the number of iterations constant. The result shown Figure 2 (right), where the axis is scaled mentioned. When the learning rate is chosen small enough, the
5

---

floating point errors are attenuated over iterations, indicating the process is stable (note it is not stuck local minima see Appendix C). Above this level, the errors are amplified.
5 PDE Stability Analysis of a Simplified CNN
In this section, we use numerical PDE tools introduced in Section 3 to analyze one-layer CNN to show analytically how the instability phenomenon described the previous section can arise. To do this, we start with a continuum model of the network, derive the analytical formula for the gradient descent, which a PDE, discretize it, which serves as a model for the optimization algorithm, and analyze the stability of the discretization through Von-Neumann analysis.
5.1 Network, Loss and Gradient Flow PDE
Let 7 IR2 + R be an input image the network, K R R be convolution kernel (we assume and K to be of finite support), r - R be the Swish activation, and 8 - IR is the sigmoid function. We consider the following CNN, which is the simplest version of VGG:
f()) = J r(K (x)
(5)
This network consists of convolution layer, activation layer, a global pooling layer and sigmoid layer. As in deep learning packages, will denote the cross-correlation (which we call convolution). The bias is not included in the network it does not impact our analysis nor conclusions.
We assume dataset that consists of just single image I and its label y E {0, 1}. This is enough to demonstrate the instability phenomena. We consider the following regularized cross-entropy loss:
1 L(K) = e(v. + K
where = f())
(6)
É denotes cross-entropy, the second term is regularization (IL2 norm squared of the kernel), and >> 0 a parameter (weight decay). Conclusions are not are restricted to cross-entropy, but it makes some expressions simpler. The gradient descent PDE is given by (see the Appendix for proof):
Theorem 5.1 (Gradient Descent PDE). The gradient descent PDE with respect to the loss (6 is
0,K - VKL(K) - - y)r (K 1) aK,
where r' denotes the derivative of the activation, and parametrizes the evolution. If we maintain the constraint that K is finite support that K is zero outside 1- /2/ then the constrained gradient descent is
0,K = (- - y)r * I aK) | W,
(8)
where W a windowing function (1 inside [-w /2, w/222 and zero outside).
5.2 Stability Analysis of Optimization
We now analyze the stability of the standard discretization the gradient descent PDEs. Note that the PDEs non-linear, and in order perform the stability analysis, study linearizations of the PDE around the switching regime of the activation, i.e., K 0. As we will see in the next sub-section, away from the transitioning regime, a similar analysis applies.
We first study the linearization (7), and then we adapt the analysis to 8 The linearization the PDE results in the following (see Appendix for the derivation):
Theorem 5.2 (Linearized PDE). The linearization of the PDE around K i given by
a (i B[K 1) 1]) 0,k = -
-ak] W, where a == s'Ogl = - y,
is considered constant with respect to K, denotes the integral (sum) of the values of (assumed finite support), and B > 0 is the parameter of the Swish activation, r(x) Ite pr.

---

Note the sign of a y can be either positive or negative; in the case a >> 0, the second condition automatically met (since 3 > 0). However in the case that O, the weight decay must be chosen large enough to be stable. The first condition shows that the weight decay and the step size satisfy relationship the optimization is to be stable. Note that the conditions also depend on the structure of the network, the current state of the network (i.e., through dependence of a), and the frequency content of the input. This is different than current deep learning practice, where the weight decay is In Appendix D, we also consider a non-constant linear model of K for a. This leads to the similar conclusions as the constant model, and so we present the constant model for simplicity.
We discretize the linearization, which results a discrete optimization algorithm, equivalent to the algorithm employed in standard deep learning packages. Using forward Fuler discretization gives
K K"
A + B((K* T) 7]) - Atak")
W
(10)
where n denotes the iteration number, and At denotes the step size (learning rate SGD). To derive the stability conditions, we compute the spatial DFT of (10) (see the Appendix D for details):
Theorem 5.3 (DFT of Discretization). The DFT of the discretization (T0) of the linearized PDE is
km++1 R = a AI( + (R sinc
w 2
Atak"
(11)
2
where K" denotes the DFT of K and sine denotes the sinc function. In the case that w + oo (the window support becomes large), the DFT of K can be written in terms of the amplifier A as
Alw)=1-Ar (w) = A(W) A(W)R - "A,,, where A(w) = Att Ar(a (a++as((()P))) (a)=)i)a((ta
(12)
In order for the updates to be a stable process, the magnitude of the amplifier should be less than one, i... |A| 1. This results in the following stability criteria (see the Appendix D for derivation): Theorem 5.4 (Stability Conditions). The discretization of the linearized PDE 7 whose DFT is given by [IT stable (in the case w + oo, i.e., the window gets large) if and only if the following conditions on the weight decay a and the step size (learning rate) At met:
a Omma
2
2 Ja min mimmmnn/eee [[(w)]2 and a O'mnn == aB max ((w))2
(13)
Al
typically chosen constant through the evolution and does not adapt with the state of the network, and the step size is typically chosen empirically without dependence on the weight decay.
5.3 Empirical Validation of Stability Bounds for Linear PDE
Figure 3 empirically validates the existence of upper and lower bounds on a. We chose the learning rate At 1e-8, a -0.5, 8 1, initialize the weights K (32 x 32) according to a normal distribution (mean 0, variance 1), and is 256 x 256 checkerboard pattern filled with 1 and -1.
The above stability conditions support. Thus, the evolution could cause the initial kernel's (finite) support to grow over iterations, which is not representative of fixed, finite support kernels used in deep network practice. In the case that W is finite and fixed support, it is not possible to analytically find simple multiplicative amplifier factor as in (2) since the convolution with the sinc function in (TI) does not separate from K" In the case that K is a constant (i.e., K is delta function or support of size 1), (f² (î(?) sinc (% Ren ()7 * sine (% -)), which case one can write the amplifier function as
under the condition that the windowing function, W, had infinite
1 A(w) = At a ( {
sinc u.
(ww))
(14)
which has the effect that the upper bound for increases and the lower bound decreases in (13) but maintains existence of upper and lower bounds.
We conjecture this is also true when K has support larger than 1, and the bounds on converge to (13) as the support of the kernel is increased. We provide empirical evidence. We find the lower and upper bounds for by running the PDE and noting the values of a when the scheme first becomes
7

---

a 1.1c7 dess
10a
a= (=11.6 1.2e7 < Genn
10ah
Gea < 0 1.9e0
Kernel
Ormin
Comin
Ommo
10
a= 0e6 dewe
16x16
3.8c6
1.07c9
2c8
2c8
10
2.1c8
a 2.2
D Gave
32x32
1.3e7
1.07e9
2e8
2e8
10tc
64x64
3.8c7
1.07c9
2c8
2c8
10a
102
109
10
10
Iteration
Figure 3: Empirical validation of stability bounds for the linearized PDE (9) When the weight decay is chosen such that a E (Aiii amaa). ). the PDE is stable and the loss converges, otherwise the PDE can be unstable and diverge. The pink and brown lines clipped the loss exceeded the largest float.
Table 2:
Comparison of the bounds on weight decay, a. between the non- windowed linear PDE, i.c. amin and amea in ((3) and the windowed linear PDE in ()) found empirically, i.e., amin and o mar Verifies bounds exist in windowed case.
unstable. We compare it to the bounds in (13 We use the same settings as the previous experiment. See results Table which verifies the bounds approach the infinite support case.
In this section, we analyzed the linearized gradient descent PDE. Momentum is often used. See Appendix E for analysis of the momentum case, which leads similar conclusions this section.
5.4 Restrained Instabilities in the Non-linear PDE
We now explain how the linear analysis around K 0 in the previous sections applies to analyzing full non-linear PDE. In particular, we show that the non-linear PDE is also not stable unless similar conditions on the weight decay and step size are met. However, the non-linear PDE can transition to regimes that have more generous stability conditions where it could be stable, and the PDE can move back and forth between stable and unstable regimes, resulting in restrained instabilities.
The non-line non-linear PDE could move between three different regimes: when the activation is "activated" (r' = 1), "not activated" (x) = 0), and "transitioning" (r'(z) is a non-constant linear function). In the transitioning regime, K x 0, in which case the linear analysis of the previous section applies, and in this regime, there are conditions on the weight decay and step size for stability. In the not-activated regime, K is negative and away from zero, the non-linear PDE driven only by the regularization term, which is stable when a < 37 This is typically true in current deep learning practice where weight decay is chosen on the order of 1e-4 and the learning rate is on order of <le-1. the activated regime, K * I is positive and away from zero and the non- lonelinear PDE can be approximated by a linear PDE (see Appendix D.5), which gives the stability condition: condition in the transitioning regime.
Approximation of Non-Linear PDE by Linear PDEs and Stability Analysis:
2
I' similar to the stability
At
Formation of Restrained Instabilities: Note that different regions in the kemne domain could each be different regimes and therefore governed by different linear PDEs discussed earlier. See Table Such regions could then switch between different regimes throughout the evolution of the PDE, and thus the cvolution could go between stable and unstable regimes. For example, in the transitioning regime, the instability quickly drives the magnitude of the kernel up (if the stability conditions not met) and thus the region to possibly stable (c.g., activated) regime. Data-driven terms could then kick the kernel back into the (unstable) transitioning regime, and this switching between regimes could happen indefinitely. As a result instabilities occur, but can be "restrained" by the activation. The optimization nevertheless amplifies and accumulates small noise in the transitioning regime.
Empirical Validation:
We now empirically demonstrate that restrained instabilites are present and error amplification occurs in the (non-linear) gradient descent PDE of the one layer CNN even when the dataset consists of a single image. We choose a 20, 3 1, y and the same input image and kernel K initialization specified in Section 5.3 Figure 4 (left) shows the loss function over
8

---

Activated
Transitioning
Not Activated
Region
{x: K x T(x) > 0}
{a K * [(x) R 0)
{ K * [(x) < 0}
0,K -s((0) VV
0.K -917 M((K Kx 1])
Linear PDE
0,K -ak
1 oK
-ak
Stability
f
<
2 - 20 min |[]]] Ai
<. A7
28 max [ww]
Table 3: The non-linear PDE is approximated by three linear PDEs in different regions of the tensor (kernel) space ("activated", "not activated" and "transitioning"). The non-linear PDE can transition between these regions in which different linear PDEs approximate the behavior.
iterations for various learning rate (At) choices, and the right plot shows the evolution of the L1 relative difference between weights of SGD (k = 1) and perturbed SGD (k 3) as discussed in Section From the previous discussion and (T3), we know that is necessary that A/ 0.1 in order for all regimes to be stable. This is confirmed by the plot on the left: when A/ > 0.1, the loss blows up. When At < 0.03, which is the empirically determined threshold for the transitioning region, all regimes are stable. and the kernel converges as shown in the plot. When 0.03 < At S. 0.1. the PDE transitions between the (unstable) transitioning regime and the stable one (activated), which introduces oscillations. This is consistent with our theory outlined previously. On the plot on the right of Figure we show divergence caused by the floating perturbations introduced in Section for for various At. When 0.03 < At S 0.1 (in the restrained instability regime), the error builds up quickly and then levels off. When At < 0.03, errors are not amplified. This is consistent with the behavior of multi-layer networks in Section validating our theory.
Average Retative Difference
0.01
10
0.05
10
0.07
0. 1itaaallee thresholvi
10a
0r00
0.5
0-0.0
10
0.05
0.07
0009
10P
10a
10a
10
1000
Iteration
Scaled Iteration
Figure 4: [Left]]: Loss iterations of the non-linear PDE for various choices of learning rates (At). The loss for various learning rates are consistent with the theory predicted for fully stable, fully unstable and restrained instabilities. [Right]: LI error accumulation in the non-linear PDE. The plot is consistent with expected error accumulation for restrained instabilities, and fully stable regimes.
Remarks Deeper Networks: This section formulated theory how the divergence phenomena for practical deep network optimization reported in Section Section/cooul coull arise. We found that even the simplest (one-layer) CNN exhibits this phenomena through analytical techniques. In multi-layer networks, there are more possibilities for such restrained instabilities to arise in multiple layers and multiple kernels of the network. Such instabilities may only arise in some kernels/layers. The analysis in this section is relevant to any one layer of such networks by treating the part of the network before the layer input; within small time period, one treat this nearly constant input the layer.
6 Conclusion
We discovered restrained numerical instabilities in standard training procedures of deep networks. In particular, we showed that epsilon-small emrrr inherent floating point arithmetic can amplify and lead to divergence in final test accuracies, comparable to stochastic batch selection. Given such variations are comparable accuracy gains reported in many papers, such instabilities could inflate or even suppress such gains, and further investigation is needed. Using simplified CNN, which is

---

amenable to linear analysis, we employed PDEs to explain how these instabilites might arise. We derived CFL conditions, which imposed conditions on the learning rate and weight decay.
This study provided a step in principally choosing learning rates for stability. We showed that numerical PDEs have promise to shed insight. The analysis showed how restrained instabilities can arise even in simplified CNN, which suggests the same phenomena exists in larger networks. Although empirically
showed how restrained instabilites could be eliminated with globally small learning rates, such instabilities are localized to space-time regions of the tensor space, suggesting in practice the need to develop adaptive learning rate schemes tailored have small rates in these local regions to ensure numerical stability, while maintaining speed.
Acknowledgements
This research was supported in part by Army Research Labs (ARL) W911NF-22-1-0261 and Raytheon Technologies Research Center.
References
[1] Y. Bengio, "Practical recommendations for gradient-based training of deep architectures," in Neural networks: Tricks of the trade, Springer, 2012, pp. 437-478.
[2] I. Goodfellow, Y. Bengio, and A. Courville, Deep learning. MIT press, 2016.
[3] L. N. Trefethen, "Finite difference and spectral methods for ordinary and partial differential equations," 1996.
[4] Courant, K. Friedrichs, and H. Lewy, "On the partial difference equations of mathematical physics," IBM journal of Research and Development, vol. 11, no. 2, pp. 215-234, 1967.
[5] J. Duchi, E. Hazan, and Y. Singer. "Adaptive subgradient methods for online learning and stochastic optimization," Journal of machine learning research, vol. 12, no. Jul, pp. 2121-2159,
2011.
[6] M. D. Zeiler, "Adadelta: An adaptive learning rate method," arXiv preprint 1212::2112.5011
2012.
[7] D. P. Kingma and J. Ba, "Adam: A method for stochastic optimization." arXiv preprint
arXiv:1412.660 2014.
[8] Y. Bengio, "Rmsprop and equilibrated adaptive learning rates for nonconvex optimization," corr htl//502.044390 2015.
[9] L. Luo, Y. Xiong, and Y. Liu, "Adaptive gradient methods with dynamic bound of learning rate." in International Conference on Learning Representations, 2019. [Online]. Available:
https://openrevier net/ wwwoormddduggggg/
[10] J. Latz, "Analysii of stochastic gradient descent in continuous time," " Statistics and Computing, vol. 31, 4. pp. 1-25, 2021.
[11] Rakhlin, O. Shamir, and K. Sridharan, "Making gradient descent optimal for strongly convex stochastic optimization," arXiv preprint arXiv: 1rXiv:110.55644 56477647, 2011.
[12] P. Toulis, D. Tran, and E. Airoldi, "Towards stability and optimality in stochastic gradient
descent," in Artificial Intelligence and Statistics, PMLR, 2016, pp. 1290-1298.
[13] J. M. Cohen, S. Kaur, Y. Li, J. 7. Kolter, and A. Talwalkar, "Gradient descent on neural networks typically occurs the edge of stability." CoRR, vol. 0tt/22133.0065. 2021. arXiv:
2103. 00065 [Online]. Available: https: : //t///aaaaa //axxviivv.. wtrggabs/21003 200e//200000065
[14] E. Weinan, C. Ma, S. Wojtowytsch, and L. Wu, "Towards a mathematical understanding of neural network- based machine learning: What we know and what we don't." arXiv preprint arXiv:2009. 10713, 2020.
[15] G. E. Karniadakis, G. Kevrekidis, L. Lu, P. Perdikaris, S. Wang, and L. Yang. "Physics- informed machine learning," Nature Reviews Physics, vol. 3, no. 6, pp. 422-440, 2021. [16] M. Burger, L. Ruthotto, S. Osher, et al., "Connections between deep learning and partial differential equations," European Journal of Applied Mathematics, vol. 32, no. 3, pp. 395-396. 2021.
[17] J. Sirignano and K. Spiliopoulos, "Dgm: A deep learning algorithm for solving partial differ- ential equations," Journal of computational physics, vol. 375, pp. 1339-1364, 2018.
10

---

[18] R. Т. Chen, Y. Rubanova, J. Bettencourt, and D. Duvenaud, "Neural ordinary differential equations," arXiv preprint 1806: 0866607666 2018.
[19] Z. Long, Y. Lu, X. Ma, and B. Dong. "Pde-net: Learning pdes from data," in International Conference on Machine Learning, PMLR, 2018, pp. 3208-3216.
[20] Y. Khoo, J. Lu, and L. Ying, "Solving parametric pde problems with artificial neural networks," European Journal of Applied Mathematics, vol. 32, no. 3, pp. 421-435, 2021.
[21] J. Han, A. Jentzen, and E. Weinan, "Solving high-dimensional partial differential equations using deep learning." Proceedings of the National Academy of Sciences, vol. 115, no. 34, pp. 8505-85 1000,, 2018.
[22] I. Ruthotto and F. Haber, "Deep neural networks motivated by partial differential equations,"
Journal of Mathematical Imaging and Vision, vol. 62, no. 3. pp. 352-364, 2020.
[23] E. Haber and L. Ruthotto, "Stable architectures for deep neural networks," Inverse problems, vol. 34, no. 1, P. 014 004, 2017.
[24] E. Weinan, "A proposal on machine learning via dynamical systems," Communications in Mathematics and Statistics, vol. 5, no. 1, pp. 1-11, 2017.
[25] P. Chaudhari and S. Soatto, "Stochastic gradient descent performs variational inference, converges to limit cycles for deep networks," in 2018 Information Theory and Applications
Workshop (ITA), IEEE, 2018, pp. 1-10.
[26] P. Chaudhari, A. Oberman, S. Osher, S. Soatto, and G. Carlier, "Deep relaxation: Partial differential equations for optimizing deep neural networks," Research in the Mathematical
Sciences, vol. 5, no. 3. pp. 1-30, 2018.
[27] A. C. Wilson, B. Recht, and M. I. Jordan, "A Iyapunov analysis of momentum methods in optimization," arXiv preprint 1.000633,, 2016.
[28] S. Osher, B. Wang, P. Yin, et al., "Laplacian smoothing gradient descent," arXiv preprint 2018.
arXiv: 10Xiv:8006.06311
[29] D. Lao, P. Zhu, P. Wonka, and G. Sundaramoorthi, "Channel-directdd gradients for optimization of convolutional neural networks," arXiv preprint arXiv:2008. 10766, 2020.
[30] G. Sundaramoorthi and A. Yezzi, "Variational pdes for acceleration on manifolds and appli- cation to diffeomorphisms." Advances in Neural Information Processing Systems, vol. 31, 2018.
[31] Y. Sun, D. Lao, G. Sundaramoorthi, and A. Yezzi, "Accelerated pdes for construction and theoretical analysis of sgd extension," in The Symbiosis of Deep Learning and Differential Equations, 2021.
[32] A. Defazio, F. Bach, and S. Lacoste-Julien, "Saga: A fast incremental gradient method with support for non-strongly convex composite objectives," in Advances in neural information processing systems, 2014, pp. 1646-1654.
[33] R. Johnson and T. Zhang, "Accelerating stochastic gradient descent using predictive variance reduction," Advances in neural information processing systems, vol. 26, pp. 315-323, 2013. [34] L. Lei, C. Chen, and M. I. Jordan, "Non-convex finite-sum optimization via sesg methods." arXiv preprint arXiv: 001.00966,,9 2017.
[35] M. Courbariaux, Y. Bengio, and J.-P. David, "Training deep neural networks with low precision multiplications." arXiv preprint arXiv: 1412.7024 2014.
[36] S. Gupta, A. Agrawal, K. Gopalakrishnan, and P. Narayanan, "Deep learning with limited numerical precision," in International conference machine learning, PMLR, 2015, pp. 1737- 1746.
[37] R. D. Richtmyer and K. W. Morton, "Difference methods for initial-value problems," Malabar, 1994.
[38] K. He, X. Zhang, S. Ren, and Sun, "Deep residual learning for image recognition," in Proceedings of the IEEE Conference on Computer Vision and Pattern Recognition (CVPR), 2016.
[39] A. Krizhevsky, G. Hinton, al., "Learning multiple layers of features from tiny images." 2009.
[40] P. Ramachandran, B. Zoph, and Q V. Le, "Searching for activation functions," arXiv preprint arXiv:/710.05941, 2017.
[41] K. Simonyan and A. Zisserman, "Very deep convolutional networks for large-scale image recognition." arXiv preprint arXiv:1409.1556 2014.

---

[42]
H. Xiao, K. Rasul, and R. Vollgraf. "Fashion minst: A novel image dataset for benchmarking
machine learning algorithms," arXiv preprint 1rXvv:1708.077477 2017.
12


Reviews:
[{"review_id": "tsYwKqrCNS", "reviewer": "Reviewer_rvJB", "summary": "This paper models optimization algorithms as PDEs and studies their stability to obtain bounds on the hyperparameters of the optimization algorithms (such as the learning rate). The authors first empirically demonstrate instabilities in training of deep networks induced by floating point approximations, and use that to make the case for choosing the learning rate carefully.", "questions": "The authors take a discrete time system (gradient descent optimization), and analyze a continuous analog of it (the gradient flow PDE), and then study ways to discretize the PDE. This is a curious approach - is there an advantage to this approach over studying the stability of the discrete time system in the first place?", "limitations": "Yes", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "The paper is well supported by theory and experiments, and is easy to follow. Section 3 provides a nice introductory example to the technique that the authors are using. While the paper is technically sound in terms of the technique presented, there is no comparison to other bounds on the learning rate obtained from convex optimization, or BIBO stability of control systems. A more compelling case could be made by either showing that the PDE technique obtains bounds that are in line with previous observations, or that the technique obtains better bounds on the learning rate.\n\nThe paper would also be stronger if the authors could extend their approach to other optimization algorithms like Adam/Adagrad, etc.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "hysmqGssCz6", "reviewer": "Reviewer_ux6y", "summary": "In this paper, the authors show that there are numerical instabilities in the training of neural networks that inherently result from floating point errors. They empirically observe the phenomenon using a simple perturbation to the learning rate of SGD, which leads to some variations in the test performance. The authors further analyze one-layer CNN in which they can analytically derive an PDE of the gradient flow and use linear analysis to derive CFL conditions for stabilities. ", "questions": "1) What does “restrained” mean exactly?\n\n2) Theorem 5.2: What can we say about the PDE with RELU activation?\n\n3) In Section 5, you consider the regularized cross-entropy loss. How about unregularized loss?\n", "limitations": "N/A", "rating": 6, "confidence": 3, "soundness": 3, "presentation": 3, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "This paper discovers some interesting empirical observations on the numerical instabilities of training neural networks with SGD. The author presents a theoretical analysis in a simplified setting and implications in choosing principle to choosing learning rates. \n\nHowever, while there are variations in test accuracy shown in Table 1, I do not think they are that surprising. We do not know whether there is any other source of non-determinism in the implementation, which might cause the variations. I wonder if it is worth checking with a simple network on toy data using batch gradient descent?\n\nOn the theoretical side, the analysis, while revealing some conditions on the learning rate and weight decay, is rather limited to a simple architecture and the optimization requires regularization albeit common.\n", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "gS6WUWJYn4", "reviewer": "Reviewer_4tXB", "summary": "This paper observes that the training dynamics of deep neural networks is unstable in the sense that a tiny perturbation of the updates lead to significant changes in the dynamics (at the level of the trajectory), and also to a change in the final test accuracy (although rather benign). The authors propose to study this phenomenon from the point of view of numerical analysis and carry this program for a simplified model of a one-hidden-layer convnet with a single input and an infinite dimensional conv filter. For this toy model, a linearized analysis shows that there is a range of stable and unstable stepsizes and weight decay regularization. This observation is confirmed numerically on a deeper model.\n", "questions": "- In Section 5, I would rather call the system an (infinite dimensional) ODE rather than PDE because only a derivative in time is involved (while PDE is usually reserved for equations which involve partial derivatives in different variables).", "limitations": "N/A", "rating": 6, "confidence": 3, "soundness": 4, "presentation": 4, "contribution": 3, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "I enjoyed reading this paper: the presentation of the results is fair and with an appreciated scientific humbleness, the experiments are carried out with care and with a proper control of the various effects that could influence the outcome and the simple mathematical model is insightful.\n\nThe main weakness is that it is not clear whether this observation/analysis is an isolated remark or whether will lead to a fruitful line of research which will potentially lead to useful practical recommendations for deep learning (admittedly, this weakness applies to many papers that propose a new point of view). So the main question is whether the conclusions and recommandations remain relevant for more complex architectures (with more that one input, more than one filter, etc). In particular, it is not clear whether it is worth the computational cost (longer training time) to attempt to avoid this instability in practice by taking smaller step-sizes (because the variance on the test error is not that big).", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}, {"review_id": "IRw9oKPg_ID", "reviewer": "Reviewer_ojkU", "summary": "The submission shows the numerical instabilities in the current training of deep networks (DNN) with SGD and claims that the numerical error is due to floating point arithmetic. \nThen, it shows the condition of the learning rate and the weight decay to perform stable training from the perspective of numerical PDE.\n", "questions": "See the first point in Weaknesses.", "limitations": "Not related to the potential negative societal impact.", "rating": 5, "confidence": 3, "soundness": 3, "presentation": 2, "contribution": 2, "strengths": "", "weaknesses": "", "quality": null, "clarity": null, "significance": null, "originality": null, "strengths_and_weaknesses": "Strengths \n- The submission shows some evidence of numerical instabilities in reproducing the results of DNN training.\n- The submission studies a simplified CNN by the numerical PDE analysis and offers some principles for choosing learning rate and weight decay.\n\nWeaknesses\n- The current experiment lacks some details, so the non-reproducibility of the training trajectory cannot be strictly guaranteed to come from floating point arithmetic. As we all know, if we perform GPU-based training, too many seeds need to be fixed (involving cuda, cudnn, etc.). Small negligence will lead to inconsistent trajectory and output in the training process. Therefore if it is claimed that the instability comes from the floating-point calculation, it should be ensured that the dynamic trajectories of weights can strictly coincide for several tries when k=1.\n- The experimental findings and theoretical results seem to be separated. No matter what the unstable factors are, as long as the conditions given by the theorem are not satisfied, it appears that the optimal trajectory will diverge.\n- The authors provide too many linear PDE-related validation experiments. Since the title and intro try to present the contribution from the perspective of DNN training, it is better to show the application of stability conditions in neural network training rather than the numerical simulation of PDE.\n- The traditional convergence analysis for the optimizer can also get the conditions related to LR and weight. It's NOT that there is ``no theoretical result related to learning rate'', but it's not as good as the heuristic methods in practical use. Secondly, in practice, it seems that the large LR (violating the stability condition) can still be used to train DNN without the divergence described by the theorem. From this perspective, the conclusion of this submission is not very practical.", "was_revised": false, "final_justification": "", "rebuttal": "", "claims_and_evidence": "", "methods_and_evaluation_criteria": "", "theoretical_claims": "", "experimental_designs_or_analyses": "", "relation_to_broader_scientific_literature": "", "essential_references_not_discussed": "", "other_strengths_and_weaknesses": "", "other_comments_or_suggestions": "", "questions_for_authors": "", "summary_of_paper": "", "summary_of_recommendation": "", "technical_quality": "", "clarity_of_presentation": "", "potential_impact": "", "robotics_focus": "", "quality_of_limitations": "", "questions_for_rebuttal": "", "reasons_to_accept": "", "reasons_to_reject": "", "questions_to_authors": "", "main_review": "", "summary_of_the_paper": "", "summary_of_the_review": "", "strength_and_weaknesses": "", "limitations_and_societal_impact": "", "review_text": "", "extra_scores": {"ethics_review_area": ""}, "extra_text": {}}]

Author rebuttal:

Decision: Accept