Neural ODEs as the Deep Limit of ResNets with constant weights

TL;DR

This paper provides a theoretical foundation showing that ResNets with shared weights converge to Neural ODEs in the deep limit, connecting discrete deep networks with continuous models.

Contribution

It proves that stochastic gradient descent on shared-weight ResNets converges to that on Neural ODEs, establishing a formal link between the two.

Findings

ResNets with constant weights approximate Neural ODEs in the deep limit

The convergence of loss functions is established

Fokker-Planck equations are used for the proof

Abstract

In this paper we prove that, in the deep limit, the stochastic gradient descent on a ResNet type deep neural network, where each layer shares the same weight matrix, converges to the stochastic gradient descent for a Neural ODE and that the corresponding value/loss functions converge. Our result gives, in the context of minimization by stochastic gradient descent, a theoretical foundation for considering Neural ODEs as the deep limit of ResNets. Our proof is based on certain decay estimates for associated Fokker-Planck equations.

Equations

$$h_{t+1}=\sigma (h_t+K_t^{(1)}\sigma (K_t^{(2)}{h}_t+b_t^{(2)})+b_t^{(1)}),t=0,\dots ,N-1$$

$$h_{t+1}=h_t+K_t^{(1)}\sigma (K_t^{(2)}{h}_t+b_t^{(2)})+b_t^{(1)},t=0,\dots ,N-1$$

$${\dot{h}}_t=f(h_t,{\theta }_t),t=[0,N],$$

$$f(\cdot ,\theta ):{\mathbb{R}}^d\to {\mathbb{R}}^d,f(x,\theta ):=K^{(1)}\sigma (K^{(2)}x+b^{(2)})+b^{(1)},$$

$${\dot{h}}_t=f(h_t,\theta ),t=[0,1].$$

$$h_{t_{i+1}}=h_{t_i}+\frac{1}{N}f(h_{t_i},\theta ),i=0,...,N-1,t_i=i/N.$$

$$\underset{\theta }{\min }\mathcal{R}(\theta ),\mathcal{R}(\theta ):={\mathbb{E}}_{(x,y)\sim \mu }\left[L(h(x,\theta ),y)\right],$$

$$d{\theta }_t=-{\nabla }_{\theta }\mathcal{R}dt+\sqrt{\Sigma }d{W}_t,$$

$${\dot{h}}_t=f(h_t,{\theta }_t),h_{t+1}=h_t+\frac{1}{N}f(h_t,{\theta }_t).$$

$$f_{\theta }:{\mathbb{R}}^d\to {\mathbb{R}}^d$$

$$x_{i+1}^{(N)}(x,\theta )$$

$$x_0^{(N)}(x,\theta )$$

$$f_{\theta }(x)=K^{(1)}\sigma (K^{(2)}x+b^{(2)})+b^{(1)}$$

$$\Vert f_{\theta }(x)-f_{{\theta }^{\prime }}(x)\Vert$$

$$\Vert f_{\theta }(x)-f_{\theta }(x^{\prime })\Vert$$

$$\Vert {\nabla }_{\theta }{f}_{\theta }(x)-{\nabla }_{\theta }{f}_{\theta }(x^{\prime })\Vert$$

$$\Vert {\nabla }_x{f}_{\theta }(x)-{\nabla }_x{f}_{\theta }(x^{\prime })\Vert$$

$$\Vert {\nabla }_{\theta }^2{f}_{\theta }(x)\Vert +\Vert {\nabla }_{\theta }{\nabla }_x{f}_{\theta }(x)\Vert$$

$$\Vert {\nabla }_x^{2}f\Vert$$

$${\mathcal{R}}^{(N)}(\theta ):={\mathbb{E}}_{(x,y)\sim \mu }\left[\Vert y-x^{(N)}(1,\theta ){\Vert }^2\right]+\gamma H(\theta )$$

$$\dot{x}(t)$$

$$\mathcal{R}(\theta ):={\mathbb{E}}_{(x,y)\sim \mu }\left[\Vert y-x(1,\theta ){\Vert }^2\right]+\gamma H(\theta ).$$

$$d{\theta }_t^{(N)}$$

$$d{\theta }_t$$

$$-\nabla \mathcal{R}(\theta )\cdot \theta \le c(1+\Vert \theta {\Vert }^2),-\nabla {\mathcal{R}}^{(N)}(\theta )\cdot \theta \le c(1+\Vert \theta {\Vert }^2)$$

$$\underset{t\in [0,T]}{\sup }\Vert {\theta }_t-{\theta }_t^{(N)}\Vert \to 0\text{in probability as N→∞}$$

$$\mathbb{E}[\mathcal{R}({\theta }_T)]<\infty ,\mathbb{E}[{\mathcal{R}}^{(N)}({\theta }_T^{(N)})]<\infty .$$

$${\tilde{\mathcal{R}}}^{(N)}(\theta )$$

$$\tilde{\mathcal{R}}(\theta )$$

$$d{\theta }_t^{(N)}$$

$$d{\theta }_t$$

$${\lambda }^{-1}\Vert \theta {\Vert }^2\le H(\theta )\le \lambda \Vert \theta {\Vert }^2$$

$$\displaystyle\sup_{t\in[0,T]}\bigl{\|}\mathbb{E}[\theta_{t}-\theta^{(N)}_{t}]\bigr{\|}$$

$$\displaystyle\sup_{t\in[0,T]}\bigl{|}\mathbb{E}[\widetilde{\mathcal{R}}(\theta_{t})-\widetilde{\mathcal{R}}^{(N)}(\theta^{(N)}_{t})]\bigr{|}$$

$$\displaystyle\sup_{t\in[0,T]}\bigl{\|}\mathbb{E}[\theta_{t}|\theta_{t}\in B_{\widetilde{R}_{0}}]-\mathbb{E}[\theta^{(N)}_{t}|\theta^{(N)}_{t}\in B_{\widetilde{R}_{0}}]\bigr{\|}$$

Full text

2010 Primary: 68T05, 65L20, Secondary: 34A45, 35Q84, 62F10, 60H10

Keywords and phrases: Machine Learning, Deep Neural Network, ResNet, Neural ODE, Ordinary Differential Equation, Stochastic Gradient Descent, Partial Differential Equations, Fokker-Planck.

March 2, 2024

1. Introduction

ResNet and Neural ODEs

Ever since the very popular ResNet paper [18] was published several authors have made the observation that the ResNet architecture is structurally similar to the Euler discretization of an ODE ([10, 27]). However, the original ‘ResNet block’ considered in [18] is defined as

[TABLE]

where $h_{t},b_{t}^{(1)},b_{t}^{(2)}\in\mathbb{R}^{d}$ and $K_{t}^{(1)},K_{t}^{(2)}:\mathbb{R}^{d}\to\mathbb{R}^{d}$ are linear operators, and $h_{N}$ represents the output of the network. The activation function $\sigma(x)=\max\{0,x\}$ is applied component-wise on vectors in $\mathbb{R}^{d}$. $N$ is the number of layers used in the construction. For standard neural networks the operators $K_{t}^{(1)},K_{t}^{(2)}$ are matrices, but for convolutional neural networks (CNNs) the operators are discrete convolution operators111Discrete convolutional operators can be phrased as matrices, see for instance [43]..

It is not immediately obvious that 1.1 can be viewed as the Euler discretization of an ODE as the activation function $\sigma$ is applied after the addition of $h_{t}$ and $K_{t}^{(1)}\sigma(K_{t}^{(2)}y_{t}+b_{t}^{(2)})+b_{t}^{(1)}$. However, removing222This was also done by the authors of ResNet in their improved follow-up paper [19] the outermost activation $\sigma$ and instead consider the difference equation

[TABLE]

we see that this is the Euler discretization (with time-step 1) of the ODE

[TABLE]

where

[TABLE]

and where $\theta$ is the short notation for $\theta=(K^{(1)},K^{(2)},b^{(1)},b^{(2)})$. Note that in 1.3 the time-step is 1 and hence the time horizon will be the number of layers $N$. This timescale is not optimal in the sense that if the ODE is stable then the system will attract to zero as $N\to\infty$.

In this paper we consider (general) autonomous ODEs as in 1.3, with time independent $\theta$ and a time horizon of $1$, i.e.

[TABLE]

This type of model is called a Neural ODE (NODE) [5] and has recently garnered a lot of attention as it has been proven to work very well in practice. Naturally we also consider the Euler discretization of 1.4, which can be written as the difference equation

[TABLE]

The models in 1.4 and 1.5 map $\mathbb{R}^{d}$ into $\mathbb{R}^{d}$ for $t=1$. We call the models in 1.4 and 1.5, (general) ResNet type deep neural networks defined on the interval $[0,1]$. Furthermore, if the input data is in $\mathbb{R}^{d}$ and the output data in $\mathbb{R}^{l}$, we have to complement the model with a function mapping $\mathbb{R}^{d}\to\mathbb{R}^{l}$ resulting in the final output.

Empirical risk minimization

Assume that we are given data distributed according to a probability measure $\mu$ where $(x,y)\sim\mu$, and $x\in\mathbb{R}^{d}$, $y\in\mathbb{R}^{l}$. Let $L:\mathbb{R}^{l}\times\mathbb{R}^{l}\to\mathbb{R}$ be a non-negative convex function. The learning problem for a model $h=h(\cdot,\theta):\mathbb{R}^{d}\times\mathbb{R}^{m}\to\mathbb{R}^{l}$, where $\theta\in\mathbb{R}^{m}$ indicates that the function $h$ is parameterized by weights $\theta$, can be formulated as

[TABLE]

where $\mathcal{R}$ is often referred to as the risk or the risk function. In practice, the risk we have to minimize is the empirical risk, and it is a well-established fact that for neural networks the minimization problem in 1.6 is, in general, a non-convex minimization problem [35, 2, 37, 8]. As such many search algorithms may get trapped at, or converge to, local minima which are not global minima [35]. Currently, a variety of different methods are used in deep learning when training the model, i.e. when trying to find an approximate solution to the problem in 1.6, we refer to [40] for an overview of various methods. One of the most popular methods, and perhaps the standard way of approaching the problem in 1.6, is back-propagation using stochastic gradient descent, see [20] for a more recent outline of the method. While much emphasis is put on back-propagation in the deep learning community, from a theoretical perspective it does not matter if we use a forward or a backward mode of auto-differentiation.

Continuous approximation of SGD

In [28, 29] it is proved that the stochastic gradient descent can be approximated by a continuous time process

[TABLE]

where $\Sigma$ is a covariance matrix and $dW_{t}$ is a standard $m$-dimensional Wiener process defined on a probability space. The idea of approximating stochastic gradient descent with a continuous time process has been noted by several authors, see [3, 4, 6, 13, 31, 32]. A special case of what we prove in this paper, see Theorem 2.7 below, is that the stochastic gradient descent 1.7 used to minimize the risk for the ResNet model in 1.5 converges to the stochastic gradient descent used to minimize the risk for the Neural ODE model in 1.4. This convergence is proved in the sense of expectation with respect to the random initialization of the weights in the stochastic gradient descent. Furthermore, we prove that the corresponding discrepancy errors decay as $N^{-1}$ where $N$ is the number of layers or discretization steps.

Novelty and significance

It is fair to say that in general there are very few papers making more fundamental and theoretical contributions to the understanding of deep learning and more specifically ResNet like neural networks. However, in the latter case there is a strain of recent and interesting contributions. In [39] the authors allow the parameters of the model to be layer- and time-dependent resulting in non-autonomous ODEs with corresponding Euler discretization:

[TABLE]

In particular, based on more restrictive assumptions on $f$, more restrictive compared to what we use, it is proved in [39] that as the number of layers tends to infinity in 1.8, the risk associated to 1.8 and defined as in 1.6, converges in the sense of gamma convergence to the risk associated to the corresponding (continuous) ODE in 1.8: we refer to Theorem 2.1 in [39] and to [7] for an introduction to gamma convergence. The authors obtain that the minima for finite layer networks converge to minima of the continuous limit, infinite layer, network. To prove that the limit exists and has nice properties they introduce a regularization which penalizes the norm of the difference between the weights in subsequent layers. We emphasize that in [39] the authors only consider the convergence of minimizers and not the convergence of the actual optimization procedure. In [10] the authors study the limit problem directly and reformulates the problem as an optimal control problem for an ODE acting on measures. However, the complexity of such networks can be quite substantial due to the time-dependency of the weights, and it is unclear what would be the best way to construct a penalization such that the limit has nice properties.

As we mentioned before, in [5] the authors consider the autonomous ODE in 1.4, i.e. they make the assumption that all layers share the same weights, and they develop two things. Firstly, they develop an adjoint equation that allows them to approximately compute the gradient with a depth independent memory cost. Secondly, they show through numerical examples that the approach works surprisingly well for some problems.

In general, the upshot of the ODE approach is the increased regularity, since trajectories are continuous and do not intersect, for autonomous ODEs, they are reversible, see [9]. However, the increased regularity comes with a cost as Neural ODEs can have difficulties solving certain classification problems, see Section 7.

Our main contribution is that we establish, in the context of minimization by stochastic gradient descent, a theoretical foundation for considering Neural ODEs as the deep limit ($N\to\infty$) of ResNets.

Overview of the paper

The rest of the paper is organized as follows. In Section 2 we introduce the necessary formalism and notation and state our results: Theorems 2.5, 2.7 and 2.11. In Section 3 we estimate, for $\theta$ fixed, the discretization error arising as a consequence of the Euler scheme, and we prove some technical estimates. In Section 4 we collect and develop the results concerning stochastic differential equations and Fokker-Planck equations that are needed in the proofs of Theorems 2.5, 2.7 and 2.11. Section 5 is devoted to the Fokker-Planck equations for the probability densities associated to the stochastic gradient descent for the Euler scheme and the continuous ODE, respectively. We establish some slightly delicate decay estimates for these densities, of Gaussian nature, assuming that the initialization density has compact support: see Lemma 5.1 below. In Section 6 we prove Theorems 2.5, 2.7 and 2.11. In Section 7 we discuss a number of numerical experiments. These experiments indicate that in practice the rate of convergence is highly problem dependent, and that it can be considerably faster than indicated by our theoretical bounds. Finally, Section 8 is devoted to a few concluding remarks.

2. Statement of main results

Our main results concern (general) ResNet type deep neural networks defined on the interval $[0,1]$. To outline our setup we consider

[TABLE]

where $\theta\in\mathbb{R}^{m}$ is a potentially high-dimensional vector of parameters acting as a degree of freedom. Given $N\in\mathbb{Z}_{+}$ we consider $[0,1]$ as divided into $N$ intervals each of length $N^{-1}$, and we define $x_{i}^{(N)}(x,\theta)$, $i=0,\ldots,N$, recursively as

[TABLE]

We define $x^{(N)}(t,x,\theta)=x^{(N)}_{i}(x,\theta)$ whenever $t\in[i/N,(i+1)/N)$. We will not indicate the dependency on $x$, $\theta$ when it is unambiguous.

We are particularly interested in the case when $f_{\theta}(x)$ is a general vector valued (deep) neural network having parameters $\theta$ but in the context of ResNets a specification for $f_{\theta}$ is, as discussed in the introduction,

[TABLE]

where $\theta=(K^{(1)},K^{(2)},b^{(1)},b^{(2)})$ are parameters and $\sigma$ is a globally Lipschitz activation function, for instance sigmoid, tanh, ReLU, Leaky ReLU and Swish, [33]. However, our arguments rely only on certain regularity and growth properties of $f_{\theta}$. We will formulate our results using the following classes of functions.

Definition 2.1**.**

Let $g:\mathbb{R}_{+}\to\mathbb{R}_{+}$, $g(0)\geq 1$, be a non-decreasing function. We say that the function $f_{\theta}(x):\mathbb{R}^{m}\times\mathbb{R}^{d}\to\mathbb{R}^{d}$ is in regularity class $\mathcal{A}(g)$ if

[TABLE]

whenever $x,x^{\prime}\in\mathbb{R}^{d}$, $\theta$, $\theta^{\prime}\in\mathbb{R}^{m}$.

Some remarks are in order concerning Definition 2.1. Firstly, it contains the prototype neural network in 2.2, with $g(s)=s+1$. Secondly, it is essentially closed under compositions, see Lemma 2.2 below. Therefore, finite layer neural networks satisfy Definition 2.1. We defer the proof of the following lemma to Section 3.

Lemma 2.2**.**

Let $f_{\theta},g_{\theta}\in\mathcal{A}(g)$. Then $F_{\theta}=f_{\theta_{2}}\circ g_{\theta_{1}}\in\mathcal{A}(2g^{3})$ for $\theta=(\theta_{1},\theta_{2})$.

Certain results in this paper require us to control the second derivatives of the risk. We therefore also introduce the following class of functions.

Definition 2.3**.**

Let $g:\mathbb{R}_{+}\to\mathbb{R}_{+}$, $g(0)\geq 1$, be a non-decreasing function. We say that the function $f_{\theta}(x):\mathbb{R}^{m}\times\mathbb{R}^{d}\to\mathbb{R}^{d}$ is in regularity class $\mathcal{A}_{+}(g)$ if $f_{\theta}\in\mathcal{A}(g)$ and if there exists a polynomial $P:\mathbb{R}\to\mathbb{R}$ such that

[TABLE]

whenever $x\in\mathbb{R}^{d}$, $\theta\in\mathbb{R}^{m}$.

The following lemma follows from Lemma 2.2 and an elementary calculation using Definition 2.3.

Lemma 2.4**.**

Let $f_{\theta},g_{\theta}\in\mathcal{A}_{+}(g)$, then $F_{\theta}=f_{\theta_{2}}\circ g_{\theta_{1}}\in\mathcal{A}_{+}(2g^{3})$ for $\theta=(\theta_{1},\theta_{2})$.

Given a probability measure $\mu$, $(x,y)\sim\mu$, $x\in\mathbb{R}^{d}$, $y\in\mathbb{R}^{d}$, and with $x^{(N)}(\cdot,\theta)$ defined as in Section 2, we consider the penalized risk

[TABLE]

where $\gamma\in\mathbb{R}_{+}$ is a hyper-parameter and $H=H(\theta)$ is a non-negative and convex regularization. The finite layer model in Section 2 is, as described in Section 1, the forward Euler discretization of the autonomous system of ordinary differential equations

[TABLE]

where $x(t)=x(t,\theta)\in\mathbb{R}^{d}$. Given data from the distribution $\mu$ and with $x(\cdot,\theta)$ solving the system of Neural ODEs in 2.4, we consider the penalized risk

[TABLE]

Throughout this paper we will assume that moments of all orders are finite for the probability measure $\mu$. By construction the input data, $x$, as well as $x(t,\theta)$ are vectors in $\mathbb{R}^{d}$. In the case when the output data is $\mathbb{R}^{l}$, we need to modify the $\mathcal{R}^{(N)}(\theta)$ and $\mathcal{R}(\theta)$ by performing final transformations of $x^{(N)}(1,\theta)$ and $x(1,\theta)$ achieving outputs in $\mathbb{R}^{l}$. These modifications are trivial to incorporate and throughout the paper we will in our derivations therefore simply assume that $l=d$.

We will base our analysis on the following continuous in time approximations of the stochastic gradient descent, see [28, 29],

[TABLE]

for $t\in[0,T]$. Throughout the paper we will assume, for simplicity, that the constant covariance matrix $\Sigma$ has full rank something which, in reality, may not be the case, see [4]. We want to understand in what sense $\theta^{(N)}_{t}$ approximates $\theta_{t}$ as $N\to\infty$. To answer this, we first need to make sure that $\theta^{(N)}_{t}$ and $\theta_{t}$ exist. Note that $\nabla\mathcal{R}^{(N)}(\theta)$ and $\nabla\mathcal{R}(\theta)$ can, as functions of $\theta$, grow exponentially: simply consider the ODE $\dot{x}=\theta x$ which has $e^{\theta t}$ as a solution. This creates problems as the setting does not fit the standard theory of SDEs, see [14], a theory which most commonly requires that the growth of the drift term is at most linear. However, if the drift terms are confining potentials, i.e.

[TABLE]

then we have existence and uniqueness for the SDEs in Section 2, see Section 4. In particular, if we have a bound on the growth of $\nabla_{\theta}x^{(N)}(\cdot,\theta)$, $\nabla_{\theta}x(\cdot,\theta)$ then, as we will see, we can choose $H$ to have sufficiently strong convexity to ensure the existence of a large constant $c$ such that the drift terms are confining potentials in the sense stated.

If we choose $H$ so that $\mathcal{R}^{(N)}$ and $\mathcal{R}$ are convex outside some large ball then $\mathcal{R}^{(N)}$ and $\mathcal{R}$ can be seen as bounded perturbations of strictly convex potentials, see the proof of Theorem 2.5. Using this we can use the log-Sobolev inequality and hyper-contractivity properties of certain semi-groups, see Section 4, to obtain very good tail bounds for the densities of $\theta_{t}^{(N)}$ and $\theta_{t}$. Actually these tail bounds are good enough for us to prove that the expected risks are bounded, expectation is over trajectories, and that $\theta^{(N)}_{t}\to\theta_{t}$ in probability. However, these bounds do not seem to be strong enough to give us quantitative convergence estimates for the difference $|\mathbb{E}[\mathcal{R}(\theta_{T})]-\mathbb{E}[\mathcal{R}^{(N)}(\theta^{(N)}_{T})]|$. The main reason for this is that even though we have good tail bounds for the densities of $\theta^{(N)}_{t}$ and $\theta_{t}$ we do not have good estimates for $\theta^{(N)}_{t}-\theta_{t}$ in terms of $N$. The following is our first theorem.

Theorem 2.5**.**

Let $g:\mathbb{R}_{+}\to\mathbb{R}_{+}$, $g(0)\geq 1$, be a non-decreasing function and assume that $f_{\theta}\in\mathcal{A}_{+}(g)$. Given $d$, $m$, there exists a regularizing function $H=H(g)$ such that if we consider the corresponding penalized risks $\mathcal{R}^{(N)}$ and $\mathcal{R}$, defined using $H$ and with $\gamma=1$, then $\mathcal{R}^{(N)}$ and $\mathcal{R}$ are bounded perturbations of strictly convex functions. Furthermore, given $T>0$ and a compactly supported initial distribution $p_{0}$ for $\theta_{0}$, we have

[TABLE]

and

[TABLE]

Remark 2.6**.**

Theorem 2.5 remains true but with a different rate of convergence, if we replace $\gamma=1$ with $0<\gamma<1$.

There are a number of ways to introduce more restrictive assumptions on the risk in order to strengthen the convergence $\theta^{(N)}_{t}\to\theta_{t}$ and in order to obtain quantitative bounds for the difference $|\mathbb{E}[\mathcal{R}(\theta_{T})]-\mathbb{E}[\mathcal{R}^{(N)}(\theta^{(N)}_{T})]|$. Our approach is to truncate the loss function. This can be done in several ways, but a very natural choice is to simply restrict the hypothesis space by not allowing weights with too large norm. Specifically, we let $\Lambda>0$ be a large degree of freedom, and we consider

[TABLE]

instead of $\mathcal{R}^{(N)}$ and $\mathcal{R}$, where $T_{\Lambda}(\theta)$ is a smooth function such that $T_{\Lambda}(\theta)=\theta$ if $\|\theta\|\leq{\Lambda}$ and $T_{\Lambda}(\theta)=\frac{\theta}{\|\theta\|}2{\Lambda}$ if $\|\theta\|\geq 2{\Lambda}$. It should be noted that, the choice of considering only the quadratic loss function is purely for simplicity and our technique works with minor modification for any other convex loss.

Having truncated the loss functions we run continuous forms of SGDs

[TABLE]

to minimize the modified risks. Using this setup, assuming also that $H(\theta)\approx\|\theta\|^{2}$ when $\|\theta\|$ is large, the modified risks in Section 2 will satisfy quadratic growth estimates at infinity, and the modified risk will be globally Lipschitz. As a consequence all the tools from classical stochastic differential equations are at our disposal, see Section 4. This allows us to prove that $\theta^{(N)}_{t}\to\theta_{t}$ in the sense of mean square convergence. However, still the classical SDE theory does not easily seem to allow us to conclude that $\widetilde{\mathcal{R}}^{(N)}(\theta_{t}^{(N)})$ converges in any meaningful way to $\widetilde{\mathcal{R}}(\theta_{t})$. To overcome this difficulty we develop a PDE based approach to obtain further estimates for the densities of $\theta^{(N)}_{t}$ and $\theta_{t}$, and their differences. In particular, we prove the following theorem.

Theorem 2.7**.**

Let $g:\mathbb{R}_{+}\to\mathbb{R}_{+}$, $g(0)\geq 1$, be a non-decreasing function and assume that $f_{\theta}\in\mathcal{A}(g)$. Let $\Lambda\gg 1$ and $T>0$ be fixed and assume that the initial weight randomization distribution $p_{0}$ has compact support in $B(0,R_{0})$, $R_{0}\geq 1$. Assume also that

[TABLE]

on $\mathbb{R}^{m}\setminus B(0,\rho_{0})$, $\rho_{0}\geq 1$ and for some $\lambda\geq 1$. Then for given $\gamma>0$, $\Lambda>0$ and positive definite $\Sigma$, there exists a positive and finite constant $c$, depending on the function $g$ as well as $d$, $m$, $\Lambda$, $\gamma$, $\Sigma$, $T$, $R_{0}$, $\lambda$ and $\rho_{0}$, such that

[TABLE]

Furthermore, if $\widetilde{R}_{0}\geq 2R_{0}$, then

[TABLE]

Remark 2.8**.**

The requirement that the regularization term $H$ grows quadratically is crucial for our quantitative estimates, and as such our result only covers $L^{2}$ regularization (weight decay) and not sparse type regularizations like $L^{1}$.

To prove Theorem 2.7 we develop certain estimates for $p^{(N)}(\theta,t)$, $p(\theta,t)$, i.e. for the to $\theta^{(N)}(t)$, $\theta(t)$, associate probability densities, by exploring the associated Fokker-Planck equations: see Section 5. In fact, we prove several estimates which give that $p^{(N)}(\theta,T)\to p(\theta,T)$ in a very strong sense, stronger than initially is recognized from the statement of Theorem 2.7. Particular consequences of our proofs are the estimates

[TABLE]

and

[TABLE]

whenever $\widetilde{R}_{0}\geq 2R_{0}$, $k\in\mathbb{N}$ and where $c=c(g,d,m,\Lambda,\gamma,\Sigma,T,R_{0},\lambda,\rho_{0})$. In particular, these estimates indicate that $p(\theta,T)$, $p^{(N)}(\theta,T)$ and $|p(\theta,T)-p^{(N)}(\theta,T)|$ have Gaussian tails away from the (compact) support of $p_{0}$.

Remark 2.9**.**

The estimates 2.9 and 2.10 can be interpreted from a probabilistic point of view. The bound 2.9 for $p$ is equivalent to $\mathbb{E}[e^{\gamma H(\theta_{T})/4}]\leq c\|p_{0}\|_{2}$, which implies that all moments of $H(\theta_{T})$ are finite. Secondly we can interpret 2.10 as the total variation distance between $p$ and $p^{(N)}$ being of order $N^{-1}$ with respect to the measure $e^{\gamma H(\theta)/4}d\theta$.

A direct consequence of 2.10 is the following corollary which states that $\theta_{t}^{(N)}$ is a weak order 1 approximation of $\theta_{t}$.

Corollary 2.10**.**

Let $\theta_{t},\theta_{t}^{(N)}$ be as in Theorem 2.7. Let $\varphi:\mathbb{R}^{m}\to\mathbb{R}$ be a continuous function satisfying the growth condition

[TABLE]

for some polynomial $P:\mathbb{R}^{m}\to\mathbb{R}$ of order $k$. Then there exists a constant $c$ depending on $g$ as well as $d$, $m$, $\Lambda$, $\gamma$, $\Sigma$, $T$, $R_{0}$, $\lambda,\rho_{0}$ and $P$ such that

[TABLE]

Using our estimates we can also derive the following theorem from the standard theory.

Theorem 2.11**.**

Let $g:\mathbb{R}_{+}\to\mathbb{R}_{+}$, $g(0)\geq 1$, be a non-decreasing function and assume that $f_{\theta}\in\mathcal{A}(g)$. Let $\Lambda\gg 1$ and $T>0$ be fixed and assume that $\mathbb{E}[\|\theta_{0}\|^{2}]<\infty$. Then

[TABLE]

3. ResNet and the Neural ODE: error analysis

In this section we will estimate the error formed by the discretization of the ODE. This will be done on the level of the trajectory as well as on the level of the first and second derivatives with respect to $\theta$. Note that by construction, see Sections 2 and 2.4 we have

[TABLE]

for all $t\in[0,1]$. We assume consistently that the paths $x^{(N)}(t,x,\theta)$ and $x(t,x,\theta)$ start at $x$ for $t=0$, and are driven by the parameters $\theta$. We will in the following, for simplicity, use the notation $x^{(N)}(t):=x^{(N)}(t,x,\theta)$, $x_{i}^{(N)}(t)=x^{(N)}(i/N)$, $x(t):=x(t,x,\theta)$. Recall that $\widetilde{\mathcal{R}}^{(N)}$ and $\widetilde{\mathcal{R}}$ are introduced in Section 2 using the cut-off parameter $\Lambda$. Let $\hat{\mathcal{R}}^{(N)}(\theta):=\widetilde{\mathcal{R}}^{(N)}(\theta)-\gamma H(\theta)$ and $\hat{\mathcal{R}}(\theta):=\widetilde{\mathcal{R}}(\theta)-\gamma H(\theta)$. In this section we prove estimates on the discrepancy error between the trajectories of $x(t)$ and the discrete trajectories $x_{i}^{(N)}$. We begin by bounding the difference.

Lemma 3.1**.**

Let $g:\mathbb{R}_{+}\to\mathbb{R}_{+}$, $g(0)\geq 1$, be a non-decreasing function and assume that $f_{\theta}\in\mathcal{A}(g)$ as in Definition 2.1. Then there exists a function $\widetilde{g}:\mathbb{R}_{+}\to\mathbb{R}_{+}$ such that

[TABLE]

holds for all $i\in\{0,\ldots,N\}$, and such that

[TABLE]

holds for all $i\in\{0,\ldots,N\}$.

Proof.

To start the proof we first write

[TABLE]

To bound the second term in 3.3 we note that

[TABLE]

by Definition 2.1. By the triangle inequality, the integral on the right-hand side of 3.4 is bounded by

[TABLE]

By the definition of $x(t)$ and Definition 2.1 we have that the second term in 3.5 is bounded as

[TABLE]

whenever $t\in[(i-1)/N,i/N]$. Assembling 3.3, 3.4, 3.5 and 3.6 we arrive at

[TABLE]

where $A_{i}=\left\|x({i}/{N})-x_{i}^{(N)}\right\|$ for $i=0,1,\ldots,N$, and $A_{0}=0$. 3.7 can be rewritten as $A_{i}\leq C_{0}A_{i-1}+C_{1}$. Iterating this inequality we see that

[TABLE]

for any $k\in\{1,\ldots,N\}$. Elementary calculations give us

[TABLE]

To finalize the proof of 3.1 we need to establish a bound for $\|x\|_{L^{\infty}([0,1])}$ in terms of the initial value $\|x\|$. By the definition of $x(t)$ and Definition 2.1 we have for any $t_{1}\in[0,1]$ such that $t_{1}+t\leq 1$

[TABLE]

Choosing $t$ small enough i.e. $t=\frac{1}{2g(\|\theta\|)}$ we deduce

[TABLE]

Iterating the above inequality $N=2\lceil g(\|\theta\|)\rceil$ times we obtain

[TABLE]

Combining 3.9 and 3.10 proves 3.1.

It remains to prove 3.2 and to start the proof note that the first term on the left in 3.2 is already bounded by 3.10. Thus, we only need to establish that $\|x_{i}^{(N)}\|$ is bounded. We first note, using the definition of $x_{i}^{(N)}$ and Definition 2.1, that

[TABLE]

By the triangle inequality we get by rearrangement that

[TABLE]

This is again an estimate of the form $A_{i}\leq C_{0}A_{i-1}+C_{1}$, where $A_{0}=\|x\|$. Iterating this we obtain as in 3.8 that

[TABLE]

and by elementary calculations as in 3.9 we can conclude that

[TABLE]

This proves the final estimate stated in 3.2 and finishes the proof. ∎

We next upgrade the previous lemma to the level of the gradient of the trajectories with respect to $\theta$.

Lemma 3.2**.**

Let $g:\mathbb{R}_{+}\to\mathbb{R}_{+}$, $g(0)\geq 1$, be a non-decreasing function and assume that $f_{\theta}\in\mathcal{A}(g)$ as in Definition 2.1. Then there exists a function $\widetilde{g}:\mathbb{R}_{+}\to\mathbb{R}_{+}$ such that

[TABLE]

holds for all $i\in\{0,\ldots,N\}$, and such that

[TABLE]

holds for all $i\in\{0,\ldots,N\}$.

Proof.

We begin by proving 3.11 under the assumption that 3.12 holds. First note that

[TABLE]

As in the proof of Lemma 3.1 we want to build iterative inequalities, and our first goal is to bound the second term in 3.13. First note that

[TABLE]

Now using Definition 2.1 multiple times we see that

[TABLE]

We want to bound the terms on the right-hand side in Section 3 and to do this we first note that

[TABLE]

The second term appearing in 3.15 is what we want to bound in an iterative scheme. Focusing on the first term in 3.15, using Definition 2.1, we see that

[TABLE]

Again, since we assume 3.12 we can use 3.13, 3, 3, 3.15 and 3.1 to get

[TABLE]

for some non-decreasing function $g_{1}$, $A_{i}:=\|\nabla_{\theta}x(i/N)-\nabla_{\theta}x^{(N)}_{i}\|$, and $A_{0}:=0$. The above is an iterative inequality of the same type as 3.7, and we get

[TABLE]

where $g_{2}$ is another non-decreasing function. This completes the proof of 3.11 under the assumption 3.12.

We now prove 3.12. Let $t\in[0,1]$ and $\epsilon>0$, and note that from Definition 2.1 we have

[TABLE]

Fix $\epsilon=1/(2g(\|\theta\|))$. Then

[TABLE]

Iterating this inequality, using also 3.2, we get

[TABLE]

for some non-decreasing function $\widetilde{g}:\mathbb{R}_{+}\to\mathbb{R}_{+}$. Similar bounds can be deduced for $\nabla_{\theta}x^{(N)}$. This establishes 3.12 and completes the proof of the lemma. ∎

In the proof of Theorem 2.5 we need to control the second derivatives of the risk and thus we need to bound the second derivatives of the trajectories w.r.t. $\theta$. We have collected what is needed in the following lemma.

Lemma 3.3**.**

Let $g:\mathbb{R}_{+}\to\mathbb{R}_{+}$, $g(0)\geq 1$, be a non-decreasing function and assume that $f_{\theta}\in\mathcal{A}_{+}(g)$. Then there exists a non-decreasing positive function $\hat{g}:\mathbb{R}_{+}\to\mathbb{R}_{+}$ such that

[TABLE]

for some polynomial $\widetilde{P}$.

Proof.

We will proceed similarly as in the proof of Lemma 3.2. Starting with $t\in[0,1]$ and $\epsilon>0$, using Definition 2.3 we have

[TABLE]

for a polynomial $\widetilde{P}$. This gives, with $\epsilon=1/(2g^{3}(\|\theta\|))$ and after absorption, using Lemmas 3.1 and 3.2, that

[TABLE]

Iterating this inequality we deduce that

[TABLE]

for some function $\hat{g}:\mathbb{R}_{+}\to\mathbb{R}_{+}$. Similar bounds can be deduced for $\nabla^{2}_{\theta}x^{(N)}$. ∎

The next lemma of this section uses Lemmas 3.1 and 3.2 to bound the difference of the risks for the discrete and the continuous systems, as well as the difference at the gradient level. Recall that $\hat{\mathcal{R}},\hat{\mathcal{R}}^{(N)}$ are the un-penalized truncated risks as defined at the beginning of this section.

Lemma 3.4**.**

Let $g:\mathbb{R}_{+}\to\mathbb{R}_{+}$, $g(0)\geq 1$, be a non-decreasing function and assume that $f_{\theta}\in\mathcal{A}(g)$. Then

[TABLE]

Proof.

We will only prove the estimates

[TABLE]

as all the other estimates are proved similarly. To get started we first note that

[TABLE]

where now $x^{(N)}_{N}=x^{(N)}_{N}(T_{\Lambda}(\theta))$, $x(1)=x(1,T_{\Lambda}(\theta))$. As $T_{\Lambda}(\theta)$ is constant when $\|\theta\|\geq 2\Lambda$, and $\|\nabla_{\theta}T_{\Lambda}(\theta)\|$ is bounded, we see that it suffices to simply derive the bounds when $\|\theta\|\leq\Lambda$. In this case $T_{\Lambda}(\theta)=\theta$. Using Section 3 we see that to estimate the terms in Section 3 it suffices to bound $\|x(1)\|,\|x_{N}^{(N)}\|,\|\nabla_{\theta}x^{(N)}_{N}\|$, $\|\nabla_{\theta}x(1)\|$, $\|x(1)-x^{(N)}_{N}\|$, and $\|\nabla_{\theta}x(1)-\nabla_{\theta}x^{(N)}_{N}\|$ which are all provided by Lemmas 3.1 and 3.2. ∎

We end this section with the proof of Lemma 2.2.

Proof of Lemma 2.2

The proof is just a matter of verifying each part of 2.3. We begin with the Lipschitz continuity in $\theta$. From the triangle inequality, and repeatedly applying Definition 2.1, we get

[TABLE]

The Lipschitz continuity in $x$ is a simple consequence of Definition 2.1

[TABLE]

For the $\theta$ gradient of $F$ the Lipschitz continuity requires a bit more work, and we note that

[TABLE]

From Definition 2.1 we see that the first term in Section 3 is bounded as

[TABLE]

The second term in Section 3 is slightly simpler to bound using Definition 2.1, and we get

[TABLE]

The last part is to prove the Lipschitz continuity of the $x$ gradient of $F_{\theta}$. To do this we simply note that

[TABLE]

4. SDEs and Fokker-Planck equations

In this section we collect and develop the results concerning stochastic differential equations and Fokker-Planck equations that are needed in the proofs of Theorems 2.5, 2.7 and 2.11.

4.1. Existence, Uniqueness and Stability of SDEs

We here prove results about a type of SDEs including the ones in Section 2. In particular, as the risk in Theorem 2.5 does not satisfy the standard linear growth assumption, see [14], we need to verify that the processes in Theorem 2.5 exist, that they are unique and that they satisfy a stability result w.r.t. parameters.

The notation in this section deviates slightly from the rest of the paper, this is by design, as the results are general in nature. We consider

[TABLE]

which is interpreted as the stochastic integral equation

[TABLE]

and where the second integral should be interpreted in the sense of Itô. We make the following assumptions.

1. A1

$a=a(t,x)$ and $b=b(t,x)$ are jointly $L^{2}$-measurable in $(t,x)\in[t_{0},T]\times\mathbb{R}^{d}.$ 2. A2

For each compact set $\mathbf{K}\subset\mathbb{R}^{d}$ there is a Lipschitz constant $L_{\mathbf{K}}>0$ such that for $x,y\in\mathbf{K},t\in[t_{0},T]$

[TABLE] 3. A3

There exists a constant $K>0$ such that for all $(t,x)\in[t_{0},T]\times\mathbb{R}^{d}$

[TABLE] 4. A4

$x_{t_{0}}$ is independent of the Wiener process $W_{t}$, $t\geq t_{0}$, and $\mathbb{E}[x_{t_{0}}^{2}]<\infty$.

Lemma 4.1**.**

[14, Theorem 3, §6]** If A1 and A2 hold, then the solutions to the stochastic differential equations 4.2 on $[t_{0},T]$, corresponding to the same initial value and the same Wiener process, are path-wise unique.

Theorem 4.2**.**

Assuming A1, A2, A3 and A4, the stochastic differential equation 4.2 has a path-wise unique strong solution $x_{t}$ on $[t_{0},T]$ with

[TABLE]

Proof.

As we have been unable to find a reference to this result, we here outline the modifications of [14, Theorem 3,§6, p.45] needed when dropping the linear growth assumption on $a$.

To begin, let us first remark that we can take $t_{0}=0$ without loss of generality. Consider the truncation operator $\mathcal{T}^{\Lambda}(x)=\Lambda\frac{x}{\|x\|}$ for $\|x\|>\Lambda$ and $\mathcal{T}^{\Lambda}(x)=x$ otherwise. Then define the truncated version of 4.1 as follows. Let $x^{\Lambda}_{0}=\mathcal{T}^{\Lambda}(x_{0})$, $a^{\Lambda}(t,x)=a(t,\mathcal{T}^{\Lambda}(x))$, $b^{\Lambda}(t,x)=b(t,\mathcal{T}^{\Lambda}(x))$, and consider

[TABLE]

It is easily seen that the truncations $x^{\Lambda}_{0},a^{\Lambda},b^{\Lambda}$ satisfy all the requirements of [14, Theorem 1,§6, p.40]. This gives us existence and uniqueness for the solution $x^{\Lambda}_{t}$ to 4.3.

In [14, Remark 3,§6, p.48], the authors state that the existence of a solution to 4.2 can be proven with the linear growth assumption replaced with A3. They claim that it is enough to prove that

[TABLE]

where $\psi(x)=\frac{1}{1+\|x\|^{2}}$ and where $c$ is independent of the truncation parameter $\Lambda$.

The claim 4.4 can be proved with Itô’s formula together with Grönwall’s Lemma, see [14]. However, we also need to prove that

[TABLE]

To do this we first apply Itô’s formula to $\|x^{\Lambda}_{t}\|^{2}$,

[TABLE]

Next, taking the supremum, multiplying with $\psi(x_{0})$, and finally taking the expectation we get

[TABLE]

The first term on the right in 4.6 is bounded by $1$. We will now focus on the third term in 4.6. Using Hölder’s inequality together with Doob’s $L^{p}$ inequality, see for instance [14, Theorem 1,§3, p.20], and A3, we get

[TABLE]

as $\psi^{2}(x)\leq\psi(x)$. The boundedness of the right-hand side in 4.7 follows from 4.4. The second term in 4.6 can be bounded in the same way as in 4.7, due to A3, which concludes the proof of 4.5. Note that the constant $c$ in 4.5 depends only on the structural assumptions of $a$, $b$ and not on $\Lambda$. Using 4.5 one can now argue as in [14, Theorem 3, p.45] to conclude that

[TABLE]

with a rate of convergence depending only on the structural assumptions A1, A2, A3 and A4. We have now established all needed modifications and the rest of the existence proof follows as in [14, Theorem 3, p.45]. ∎

Lemma 4.3** ([14], Theorem 3, Chap 2, §7).**

Let $x_{n,t}$, $n=0,1,2,\ldots$ be solutions of

[TABLE]

where the coefficients satisfy A1 and A4, with $a$ and $b$ globally Lipschitz continuous (with a constant independent of $n$), as well as

[TABLE]

If for each $\Lambda>0$, and $s\in[0,T]$,

[TABLE]

then

[TABLE]

as $n\to\infty$.

Theorem 4.4**.**

Let $x_{n,t}$, $n=0,1,2,\ldots$ be solutions of

[TABLE]

where $a_{n}$, $b_{n}$ satisfy A1, A2, A3 and A4 with constants independent of $n$. If for each $\Lambda>0$, and $s\in[0,T]$,

[TABLE]

then

[TABLE]

in probability as $n\to\infty$.

Proof.

We will follow [14, Theorem 3, §7] and modify accordingly. Let $a_{n}^{\Lambda}(t,x)=a_{n}(t,\mathcal{T}^{\Lambda}(x))$, $b_{n}^{\Lambda}(t,x)=b_{n}(t,\mathcal{T}^{\Lambda}(x))$, $x_{0}^{\Lambda}=\mathcal{T}^{\Lambda}(x_{0})$ and consider the solution $x_{n,t}^{\Lambda}$ to

[TABLE]

Following the proof of [14, Theorem 3, §6] we first note that

[TABLE]

I.e. as long as the truncated process $x_{n,t}^{\Lambda}$ stays within a ball of radius $\Lambda$ then we have $x_{n,t}^{\Lambda}=x_{n,t}$ a.s.: for a detailed proof see [14, Theorem 3, §7]. Using 4.9 and the triangle inequality we see that the following holds for any $\epsilon>0$

[TABLE]

The first term in 4.10 goes to zero as $n\to\infty$ because our truncated equations satisfy the requirements of Lemma 4.3. The two other summands converge to 0 as $\Lambda\to\infty$ uniformly in $n$ because the rate in 4.8 only depends on the constants in A1, A2, A3 and A4 and the initial data, see the proof of 4.8. Thus, we conclude

[TABLE]

which completes the proof. ∎

4.2. Fokker-Planck equations for stochastic minimization

In general the density of a diffusion process, for us a solution to a stochastic differential equation, satisfies a partial differential equation usually called the Fokker-Planck equation or Kolmogorov’s first (forward) equation. In this section we will explore this connection for the SDEs as in Section 2 of the form

[TABLE]

with initial datum $\theta_{0}\sim p_{0}$ and for $V\in C^{2}(\mathbb{R}^{m})$. Furthermore, we have for simplicity assumed $\Sigma=\sqrt{2}I_{m}$: the generalization to general constant full rank $\Sigma$ is straightforward. If $-\nabla V$ and $\theta_{0}$ satisfy A1, A2, A3 and A4, then the associated density $p_{t}$ for $\theta_{t}$ satisfies the following Cauchy problem for the Fokker-Planck equation, see for instance [14],

[TABLE]

The equation 4.11 is sometimes called the linear Fokker-Planck equation, and for an interesting overview see for instance [41, Chapter 8]. The Cauchy problem in 4.11 has unique ‘probability’ solutions in the sense of measures, i.e. non-negative solutions that define a probability density, see [1, 36]. Using the identity

[TABLE]

we see that we can formally rewrite 4.11 as

[TABLE]

for $\widetilde{p}=e^{V}p$, and where $\langle\cdot,\cdot\rangle$ denotes the dot–product in $\mathbb{R}^{m}$. Introducing the operator $L(\cdot)=-\Delta(\cdot)+\langle\nabla V,\nabla(\cdot)\rangle$ we note that $L$ is self-adjoint in $L^{2}_{\mu}=L^{2}_{\mu}(\mathbb{R}^{m})$ for functions in $W_{\mu}^{1,2}=W_{\mu}^{1,2}(\mathbb{R}^{m})$, where $d\mu:=e^{-V}dx$. $L^{2}_{\mu}$ and $W_{\mu}^{1,2}$ are the standard $L^{2}$ and Sobolev spaces but defined with respect to $\mu$. Indeed, integrating by parts we have

[TABLE]

where $\langle\cdot,\cdot\rangle_{L^{2}_{\mu}}$ denotes the $L^{2}$ inner product w.r.t. the measure $\mu$. Assuming that $V$ satisfies a logarithmic-Sobolev inequality, see Section 4.3, it follows, see [17, Theorem 6], that $e^{tL}$ defines a strong semi-group which is hyper-contractive. Furthermore, the semi-group is positivity preserving. This establishes the existence and uniqueness of solutions to the problem in 4.11 for initial data in $L^{2}_{\mu}$.

The question is now, how do the solutions to 4.11 and 4.12 relate to each other? To answer this question, first note if $\widetilde{p}_{0}=e^{V}p_{0}\in L^{2}_{\mu}$, $p_{0}\geq 0$, then for $\widetilde{p}=C_{0}^{-1}e^{tL}\widetilde{p}_{0}\geq 0$, with $C_{0}=\int\widetilde{p}_{0}d\mu$, we have from 4.13 that

[TABLE]

Thus we have established that if the initial data $\widetilde{p}_{0}$ is a probability density w.r.t. $\mu$ then the solution $\widetilde{p}$ is also a probability density for all $t$ w.r.t. $\mu$. This implies that $\hat{p}=C_{1}^{-1}e^{-V}\widetilde{p}$ with $C_{1}=\int e^{-V}dx$ solves 4.11 and that it is a probability density w.r.t. the Lebesgue measure. The uniqueness of ‘probability’ solutions to 4.11, again see [1, 36], finally gives us that $\hat{p}=p$.

4.3. Logarithmic Sobolev inequalities and hyper-contractivity

The consequence of the previous subsection is that we can switch between solutions to 4.11 and solutions related to the semi-group $e^{tL}$ provided that the Gibbs measure $e^{-V}dx$ satisfies a logarithmic Sobolev inequality. In this section we discuss when the logarithmic Sobolev inequality holds. We begin with the following lemma which is a direct consequence of the Bakry-Émery theorem, see for instance [42, Thm 21.2].

Lemma 4.5**.**

Assume that $V(x)\in C^{2}(\mathbb{R}^{m})$ and $\nabla^{2}V\geq KI_{m}$ then $d\mu=e^{-V(x)}dx$ satisfies

[TABLE]

The inequality in Lemma 4.5 is the so called logarithmic Sobolev inequality with constant $2/K$. Often the potential $V$ is not strictly convex, but satisfies some form of ‘convexity at infinity’. For such potentials the logarithmic Sobolev inequality carries over, as is made rigorous by the Holley-Stroock perturbation lemma, see for instance [34, Prop 3.1.18]. This lemma states that if we perturb the potential with a function of bounded oscillation, then the logarithmic Sobolev inequality is preserved at the expense of a larger constant.

Lemma 4.6**.**

Assume that the probability measure $\mu$ on $\mathbb{R}^{m}$ satisfies the logarithmic Sobolev inequality 4.14 with constant $c$. Assume that $W:\mathbb{R}^{m}\to\mathbb{R}$ is a bounded and measurable function. Then the modified probability measure

[TABLE]

satisfies the log-Sobolev inequality with constant $ce^{\operatorname{osc}W}$, where $\operatorname{osc}W=\sup W-\inf W$.

Let us now consider hyper-contractivity of the semi-group $e^{tL}$. We define the $p$ to $q$ norm of the semi-group $e^{tL}$ as

[TABLE]

It then follows from a theorem of Stroock, see [12, Section 4, Theorem 4.1, Gross] or [38], that there exists, for each $p$ and $q$, where $q\geq p$, a time $t_{p\to q}$ such that if $t\geq t_{p\to q}$ then

[TABLE]

The estimate in 4.15 is the hyper-contractivity of the semi-group $e^{tL}$.

5. Estimates for Fokker-Planck equations

In the forthcoming section we prove our main results, and we will for simplicity only give the proof in the case $\Sigma=\sqrt{2}I_{m}$: the generalization to general constant full rank $\Sigma$ being straightforward. Consequently, we in this section consistently assume $\Sigma=\sqrt{2}I_{m}$.

To prove Theorem 2.7 we will first prove perturbation estimates for the linear Fokker-Planck equation 4.12 and then utilize this to prove convergence of the processes, $\theta^{(N)}(t)$, and the value, $\mathbb{E}[\widetilde{\mathcal{R}}(\theta^{(N)}(t))]$. However, we first discuss what is a priori known about solutions to 4.11 and 4.12 in the sense of regularity and integrability.

We begin with the Fokker-Planck equations, see 4.11, for the processes $\theta^{(N)}(t)$, $\theta(t)$,

[TABLE]

for $t\in(0,T]$ and where $T$ denotes the time frame used in the stochastic gradient descent. Recall that $\widetilde{\mathcal{R}}^{(N)}$ and $\widetilde{\mathcal{R}}$ are introduced in Section 2, and that $p^{(N)}$, $p$, are the probability densities to $\theta^{(N)}(t)$, $\theta(t)$, respectively. Recall from Section 4 the equations in Section 5 have unique ‘probability‘ solutions in the sense of measures. Furthermore, since $\widetilde{\mathcal{R}}^{(N)}$, $\widetilde{\mathcal{R}}$ are smooth, it follows from Schauder estimates, see [30], and a bootstrap argument that $p^{(N)},p$ are smooth functions of $(\theta,t)$. This implies that uniqueness holds in the point-wise sense for probability solutions. Let us now consider the integrability properties of solutions to Section 5. First note that Lemma 3.4 implies that the coefficients in Section 5 satisfies all the regularity and integrability assumptions of [24, Proposition 1, Remark 7], i.e.

[TABLE]

and the same holds true for $\widetilde{\mathcal{R}}^{(N)}$. Thus, there exists a unique solution

[TABLE]

to the problem for $p$ in Section 5, in the sense that

[TABLE]

for test-functions $\phi\in C^{\infty}(\mathbb{R}^{m}\times(0,T])\cap C(\mathbb{R}^{m}\times[0,T])$ with $\mbox{supp }\phi(\cdot,t)\subset K$ for all $t\in[0,T]$ and for some compact set $K\subset\mathbb{R}^{m}$. Solutions of 5.3 are called weak solutions and more info about weak solutions can be found in [11]. Now, as $p$, $\nabla p$ are both in $L^{2}$ we are allowed to use $p$ as a test-function in 5.3, resulting in the Caccioppoli estimate

[TABLE]

if $\Delta\widetilde{\mathcal{R}}\leq\rho$. By Grönwall’s lemma we get from the above Caccioppoli estimate

[TABLE]

for $0\leq t\leq T$. 5.2, 5.3, 5.4 and 5.5 also hold true with $p$, $\widetilde{\mathcal{R}}$, replaced by $p^{(N)}$, $\widetilde{\mathcal{R}}^{(N)}$. However, this type of estimate is not strong enough for our purposes.

Another approach is to switch to the semi-group $e^{tL}$ related to the linear Fokker-Planck equation 4.12 as in Section 4.2, i.e. we consider weak solutions in spaces defined with respect to the measure $d\mu:=e^{-\widetilde{\mathcal{R}}}d\theta$. As can be easily seen, the definition of $\widetilde{\mathcal{R}}$ implies that $\widetilde{\mathcal{R}}$ is, due to the truncation, a bounded perturbation of a strictly convex potential. Hence, the relations stated in Section 4.2 hold. Therefore, we have a unique $p$ for which $\widetilde{p}=e^{\widetilde{\mathcal{R}}}p$ satisfies the weak formulation of 4.12

[TABLE]

Here $\widetilde{p}\in L^{2}_{\mu}$ and the test-functions satisfy $\phi\in C^{\infty}(\mathbb{R}^{m}\times(0,T])\cap C(\mathbb{R}^{m}\times[0,T])$ and $\phi(\cdot,t)\in L^{2}_{\mu}$ for all $t\in[0,T]$. In addition, the initial data is assumed to satisfy $\widetilde{p}_{0}:=e^{\widetilde{\mathcal{R}}}p_{0}\in L^{2}_{\mu}$. In this case we deduce the a priori estimate

[TABLE]

for $0\leq t\leq T$. Again, 5.6 and 5.7 also hold true with $p$, $\widetilde{\mathcal{R}}$, replaced by $p^{(N)}$, $\widetilde{\mathcal{R}}^{(N)}$. In particular, for our problems at hand we get stronger a priori estimates compared to 5.5 by using the latter approach instead of the one in [24]. The purpose of the section is to prove the following lemma.

Lemma 5.1**.**

Let $T>0$, $t\in[0,T]$ and consider $\widetilde{\mathcal{R}}^{(N)}$, $p^{(N)}$, $\widetilde{\mathcal{R}}$, $p$. Assume $p_{0}\in L^{2}(\mathbb{R}^{m})\cap L^{\infty}(\mathbb{R}^{m})$ and let $\alpha\in\mathbb{R}$, $\psi\in C^{\infty}(\mathbb{R}^{m})$, $0\leq\psi\leq 1$. Then there exists a constant $c=c(m)$, $1\leq c<\infty$, such that the following holds

[TABLE]

Furthermore, let $E$, $F$ be compact sets of $\mathbb{R}^{m}$ such that $d(E,F)>0$, $d(E,F)$ is the distance between $E$ and $F$, and assume that $\mbox{supp }p_{0}\subset E$. Then

[TABLE]

Proof.

We will only supply the proof of the estimates for $p$ as the proofs in the case of $p^{(N)}$ are analogous. Recall that $p$ is smooth. Let $\alpha$ and $\psi$ be as in the statement of the lemma and let $\varphi(\theta)=e^{\alpha\psi(\theta)}$. We let $\widetilde{p}=e^{\widetilde{\mathcal{R}}}p$, and we introduce the measure $d\mu=e^{-\widetilde{\mathcal{R}}}d\theta$. Using this notation, and using the test-function $p\varphi^{2}\in W^{1,2}_{\mu}$ in 5.6 together with an integration by parts, we see that

[TABLE]

Using 5.10 we first deduce that

[TABLE]

and then by Young’s inequality, Cauchy-Schwartz and re-absorption we get

[TABLE]

Let $\gamma:=c\alpha^{2}\|\nabla\psi\|_{\infty}^{2}$ and note that 5.12 together with Grönwall’s lemma give

[TABLE]

where $c=c(m)>0$. In particular going back to $p$, 5.13 becomes,

[TABLE]

and this concludes the proof of Lemma 5.1 for $p$.

To prove the estimate in Lemma 5.1 for $p$, let $E$, $F$ be compact sets of $\mathbb{R}^{m}$ such that $d:=d(E,F)>0$, $d(E,F)$ is the distance between $E$ and $F$, and assume that $\mbox{supp }p_{0}\subset E$. Let $\psi=1$ on $F$ and $\psi=0$ on $E$, $0\leq\psi\leq 1$. One can construct $\psi$ so that $\|\nabla\psi\|_{\infty}\leq cd^{-1}$, $c=c(m)>0$. Using 5.14 we see that

[TABLE]

Letting $\alpha=\beta d^{2}/t$ with $\beta$ so that $(c\beta-2)=-1$ we see that

[TABLE]

for some $\widetilde{c}=\widetilde{c}(m)>0$. This concludes the proof of Lemma 5.1 and the proof of the lemma.

∎

6. Proof of the main results: Theorems 2.5, 2.7 and 2.11

In this section we prove Theorems 2.5, 2.7 and 2.11, and we will for simplicity only give the proof in the case $\Sigma=\sqrt{2}I_{m}$.

6.1. Proof of Theorem 2.5

Let us first prove

[TABLE]

We will prove 6.1 under the assumption that $f_{\theta}\in\mathcal{A}(g)$. The first step is to show that there exists a penalization $H(\theta)$ such that

[TABLE]

for some constant $c$. Let $\hat{\mathcal{R}}(\theta):=\mathcal{R}(\theta)-\gamma H(\theta)$, $\hat{\mathcal{R}}^{(N)}(\theta):={\mathcal{R}}^{(N)}(\theta)-\gamma H(\theta)$, denote the un-penalized risk. As in Section 3 we have

[TABLE]

A bound on the gradient of the risk follows from Lemma 3.1,

[TABLE]

and the same bound holds for $\nabla\hat{\mathcal{R}}^{(N)}$. Let $\phi_{R}$ be a cutoff function supported in $B(0,2R)$ such that $\phi_{R}=1$ on $B(0,R)$ and $\|\nabla\phi_{R}\|\leq\frac{c}{R}$. Let $G(\theta)$ be a convex function satisfying $\nabla G(\theta)=c\widetilde{g}(g(\|\theta\|))^{2}\theta$. With $\theta_{R}$ and $G$ as above we obtain from 6.3

[TABLE]

with a constant $c(R,\mu)\geq 1$. Letting $H(\theta):=(1-\phi_{R}(\theta))G(\theta)$ proves 6.2. Collecting what we proved in Section 4.1, and applying Theorems 4.4 and 3.4, we get 6.1.

We next prove the second part of Theorem 2.5, i.e.

[TABLE]

We will prove 6.4 under the more restrictive assumption $f_{\theta}\in\mathcal{A}_{+}(g)$. We begin by showing that there is a rotationally symmetric penalization $H$ such that $\hat{\mathcal{R}}+H$ is strictly convex on $\mathbb{R}^{m}\setminus B(0,R)$ for $R\geq 1$ large. To do this we will bound the growth of the second derivatives of $\hat{\mathcal{R}}$. We have

[TABLE]

Using Lemmas 3.2 and 3.3 we derive the bound

[TABLE]

Hence, if we penalize using a function $\widetilde{H}(\theta)$ satisfying

[TABLE]

for a large enough $C(\mu)\geq 1$, then $\hat{\mathcal{R}}(\theta)+\widetilde{H}(\theta)$ is strictly convex outside $B(0,R)$. In fact, consider a cutoff function $\phi_{R}$ as in the first part of the proof, and define $V=(1-\phi_{R})\hat{\mathcal{R}}+\widetilde{H},W=\phi_{R}\hat{\mathcal{R}}-\phi_{R}\widetilde{H}$. Then, for a possibly slightly larger constant $C=C(R,\mu)$, we see that $V$ is strictly $K$-convex and that $W$ is bounded. Thus, $\hat{\mathcal{R}}+(1-\phi_{R})\widetilde{H}$ is the sum of a strictly convex potential $V$ and a bounded perturbation $W$. Define the penalization as $H=(1-\phi_{R})\widetilde{H}$, then $\hat{\mathcal{R}}+H$ is a confining potential that satisfies A1, A2, A3 and A4. From Lemmas 4.5 and 4.6 we see that $d\mu=ce^{-(\hat{\mathcal{R}}+H)}d\theta$ satisfies a logarithmic Sobolev inequality and thus we can conclude the validity of the estimate in 4.15.

Recall from Section 4.3 that for $V=\hat{\mathcal{R}}+H$, we have $\widetilde{p}=e^{tL}c^{-1}e^{V}p_{0}=c^{-1}e^{V}p$ where $p$ is the density for the stochastic process

[TABLE]

The hyper-contractivity, i.e. 4.15 shows that $\widetilde{p}(\cdot,t)\in L^{r}_{\mu}$ for $t\geq t_{r}$, $r\geq 2$. This implies that

[TABLE]

Using Hölder’s inequality and the above we get

[TABLE]

The first term on the second line in Section 6.1 is bounded by $1$. To bound the last term in Section 6.1, note that if $V>\frac{r}{r-1}$ then the expression $e^{-V}V^{\frac{r}{r-1}}$ decreases as $V$ increases. We can thus use the estimate

[TABLE]

for $V>\frac{r}{r-1}$. The term on the right-hand side in the last display is integrable and thus by the Lebesgue dominated convergence theorem we see that $\mathbb{E}[V(\theta_{t})]<\infty$. Specifically we could use $r=2$ and thus get 6.4, completing the proof of the theorem.

Remark 6.1**.**

If $t\to\infty$ then $p_{\infty}=ce^{-V}$ and from Section 6.1 we get

[TABLE]

6.2. Proof of Theorem 2.7

To start the proof of Theorem 2.7 we let $\delta^{(N)}(x,t):=p(\theta,t)-p^{(N)}(\theta,t)$ and we note, as $\theta^{(N)}$ and $\theta$ share the same initial density, that $\delta^{(N)}$ satisfies the initial value problem

[TABLE]

Rewriting the equation Section 6.2 as in 4.12, we get

[TABLE]

Now, we would like to test the above equation with $e^{\widetilde{\mathcal{R}}}\delta^{(N)}$ but we do not a priori know that this test-function is in $L^{2}_{\mu}$, for $d\mu=e^{-\widetilde{\mathcal{R}}}d\theta$. Instead, we use the test-function $\phi^{2}(\theta)e^{\widetilde{\mathcal{R}}}\delta^{(N)}$, where $\phi\in C_{0}^{\infty}(\mathbb{R}^{m})$, $0\leq\phi\leq 1$. From the above display using integration by parts we get

[TABLE]

where $d\mu=e^{-\widetilde{\mathcal{R}}}d\theta$. Using Youngs inequality and Cauchy-Schwarz to reabsorb terms on the left-hand side we conclude that

[TABLE]

Furthermore,

[TABLE]

To get rid of the localization that $\phi$ provides we let $\phi=\phi_{R}\in C_{0}^{\infty}(B(0,2R))$ be such that $\phi_{R}=1$ on $B(0,R)$ and $\|\nabla\phi\|<\frac{C}{R}$. We will consider the limit as $R\to\infty$ in Sections 6.2 and 6.8. Specifically, using the estimate for $p$ stated in Lemma 5.1, Lemma 5.1, we see that

[TABLE]

as $R\to\infty$ since $d(\operatorname{supp}p_{0},B(0,2R)\setminus B(0,R))\approx R$ for large values of $R$. Similarly,

[TABLE]

Now, recall the definitions of $\widetilde{\mathcal{R}}^{(N)}$ and $\widetilde{\mathcal{R}}$ introduced in Section 2. Using Lemma 3.4 we see that

[TABLE]

Thus applying 6.11 in 6.10 we get

[TABLE]

as $R\to\infty$. In particular collecting Sections 6.2, 6.8, 6.9 and 6.12, and letting $R\to\infty$ for $\phi=\phi_{R}$ we can conclude that

[TABLE]

Now, to bound the right-hand side of 6.13 we use Lemma 3.4 again to see that

[TABLE]

Hence, using 6.11 and 6.14, we arrive at

[TABLE]

We can apply the estimate for $p^{(N)}$ stated in Lemma 5.1, see Lemma 5.1, and we deduce that

[TABLE]

Now collecting 6.13, 6.15 and 6.16 we get

[TABLE]

and finally

[TABLE]

Using Cauchy-Schwartz and 6.18, we get

[TABLE]

We will now use the estimate in 6.19 to conclude the proof of Theorem 2.7. Indeed, by 6.19

[TABLE]

Now by the definition of $\widetilde{\mathcal{R}}$ as the truncated risk + a quadratic penalization we see that $\widetilde{\mathcal{R}}\approx\|\theta\|^{2}$ for large $\|\theta\|$ and thus $\sup_{\|\theta\|}\|\theta\|e^{-\widetilde{\mathcal{R}}/4}\leq c(g,\Lambda,\lambda,\rho_{0})$. This together with Section 6.2 concludes the proof of 2.7 in Theorem 2.7. To prove the other estimate, 2.8, we note that

[TABLE]

Again using 6.19 we see that

[TABLE]

Furthermore, using 6.11, we see that

[TABLE]

Collecting 6.21, 6.22 and 6.23 completes the proof of 2.8 in Theorem 2.7. To prove the second set of estimates we simply have to repeat the argument and use the estimates stated in Lemma 5.1, Lemma 5.1.

6.3. Proof of Theorem 2.11

The result follows immediately from the fact that based on the truncation procedure $T_{\Lambda}$, and Lemmas 3.1 and 3.4, we can conclude that the processes $\theta$ and $\theta^{(N)}$ both satisfy the assumptions of Lemma 4.3, and this proves our theorem.

7. Numerical experiments

To investigate how the optimized networks depend on the number of residual connections, i.e. layers, we have conducted a number of numerical simulations. The idea is simply to optimize 1.5 for different values of $N$. To do so we consider a few different topologies and problems. We have conducted our experiments on the following datasets:

1. (1)

Cifar10 [23] 2. (2)

Starlike Annuli + Starlike Annuli with augmentation.

7.1. Cifar10

Cifar10 is a well-known dataset consisting of images of 10 classes, [airplanes, cars, birds, cats, deer, dogs, frogs, horses, ships, and trucks]. The task is to take a 32x32 image as input and as output the class of the image. A lot of research papers have proposed different topologies for Cifar10 and currently the accuracy of the best performing model is $96.53\%$ [16]. However, the model closest to our experiments is the Residual Network (ResNet) [18] based on which they obtained an accuracy of $\approx 93.4\%$. For comparison, the average human performance is $93.91\pm 1.52\%$ [21]. From our perspective, the interesting point about this dataset is that it is the simplest standard dataset where convolutional networks vastly outperform other models, like for instance fully connected neural networks. For this dataset we have designed a simple network that performs fairly well compared to state of the art mentioned above.

Network topology

To continue let us describe the network topology used in our experiments, see Fig. 1 for a schematic. Since networks of the type in 1.5 need to have the same input-dimension and output-dimension, channels, we first apply a convolutional layer with $3\times 3$ kernels and $256$ output channels with a stride of $2$. This is then the input to the iterated block, of type 1.5, that has $256$ input channels and $256$ output channels. The block is applied $N$ times. After this we apply a batch normalization layer, [22], as well as a $3\times 3$ kernel layer with $64$ channels with a stride of $4$ and ReLU activation. This is then followed by another batch normalization layer, a 4x4 average pooling layer and finally a softmax layer. The network in the main iterated block has the type 2.2. Specifically, it consists of a $3\times 3$ kernel layer with $256$ output channels and ReLU activation. This feeds into another $3\times 3$ kernel layer with again $256$ output channels, with no activation. Although this network topology is quite simple the layers are fairly wide, and we end up with a total number of about 1M parameters.

Experimental setup

The experimental setup is as follows. The Cifar10 dataset consists of 50k images in the training set and 10k images in the test set. In order to get a better estimate on the generalization performance, we performed 10-fold cross validation on the training set, for $N=[1,2,3,4,5,9,15,20,30,100]$. We use a weight decay of $10^{-6}$, a momentum of $0.9$, and we use the weight initialization from [15]. The models are trained with a batch size of $128$ on a single GPU (RTX 2080ti) for $N<10$ and a batch size $512$ split over $4$ GPUs (RTX 2080ti) for $N\geq 10$. We start with a learning rate of $0.1$, divide it by $10$ at $80$, $120$, $160$ and by $2$ at $180$ epochs, where an epoch consists of 45k/batchsize iterations. Training is terminated at $200$ epochs. We use the simple data augmentation depicted in [25] when training: i.e. $3$ pixels are padded on each side, and a $32\times 32$ crop is randomly sampled from the padded image or its horizontal flip. For testing, we only evaluate the single view of the original $32\times 32$ image. All images are normalized by the pixel mean.

Results

The results are depicted in Fig. 2, and we note that our best model gives us $90.5\%$ which is good for such a simple topology. It is interesting to note that the cross validation accuracy is increasing with layer count, and we saw a similar increase in the training accuracy. As the parameter count stays the same with increasing depth this means that for this problem deeper is better. Also, note that the accuracy flattens out fairly quickly. This can be interpreted as to indicate that there is limited value in using too many layers for this type of image classification problems. It may be worthwhile to implement a limit on the number of layers in the methods developed in [5]. Another observation is that the speed of convergence for the validation accuracy is close to $1/N$, albeit slightly faster. We believe that the convergence rate for the risk (or related metric) will in practice almost always be faster than our theoretical bound of $1/N$, see Theorem 2.7.

7.2. Starlike Annuli dataset

An interesting observation that was made in [9], is that certain simple functions cannot be represented as a neural ODE of the type defined in 1.4. In [9] the authors give the following example. Let $0<r_{1}<r_{2}<r_{3}$ and let $g:\mathbb{R}^{d}\to\mathbb{R}$ be such that

[TABLE]

To attempt to represent $g$ using a neural ODE $h_{t}$ as in 1.3 we compose the output of the neural ODE $h_{1}(x)$ with a linear map $\mathcal{L}:\mathbb{R}^{d}\to\mathbb{R}$. In [9] they prove that $h_{1}(x)$ is a homeomorphism. In order for $\mathcal{L}(h_{1}(x))=g(x)$ to hold, the neural ODE $h_{1}(x)$, which is a homeomporhism, have to map the ball $B_{r_{1}}$ to one side of a hyper-plane and the outer annuli to the other side, which is impossible. For the numerical simulations in [9] the authors consider $g$ as defining a binary classification problem. The data-points are sampled uniformly at random from the two classes defined by $g$ with a 50/50 class balance. When solving this binary classification problem using the neural ODE model, the authors in [9] noted that the corresponding ODEs are numerically hard to solve and require many costly function evaluations. They propose a solution to this problem by simply embedding the dataset in a bigger space, for instance $\mathbb{R}^{d}\to\mathbb{R}^{d+1}$. This solves the problem of representability, in the sense above, and the corresponding classification problem becomes easy to solve using Neural ODEs. The process of embedding the dataset in a bigger space is called augmentation.

Dataset

Our purpose is to compare the convergence rate of the risk for non-representable problems vs representable problems, in the sense above. As such, a comparison between 7.1 and its augmentation makes sense. However, the augmented version of 7.1 is very easy as it is essentially a conical solution, which is easy to represent with a neural network using ReLU activation functions. We therefore consider a related dataset that is harder than 7.1 even in the augmented case, we call this dataset the Starlike Annuli. To describe our construction, see Fig. 3, we let $r_{1}<r_{2}<r_{3}$ and define, using polar coordinates,

[TABLE]

Here the difference $r_{3}-r_{1}$ denotes the class separation. The data is generated by sampling uniformly at random from each of the two class regions. The augmented version is the dataset embedded in 3 dimensions. For our experiments we chose $r_{1}=1$, $r_{2}=1.5$, and $r_{3}=3$.

Network topology

The topology of our network is 2.2 and iterated as in 1.5, i.e.

[TABLE]

where $K^{(2)}\in\mathbb{R}^{m}\times\mathbb{R}^{2+n}$, $b^{(2)}\in\mathbb{R}^{m}$, $K^{(1)}\in\mathbb{R}^{2+n}\times\mathbb{R}^{m}$, $b^{(1)}\in\mathbb{R}^{2+n}$ and $\sigma$ is the ReLU activation function. Following [9], $n\geq 0$ is the number of augmented dimensions, $m>0$ is the size of the inner layer. In our experiments we chose $m=16$ and $n=0,1$. At the end of the network we have a softmax layer. We are using a quadratic penalization, as in the theoretical result Theorem 2.7, with constant $\gamma=0.001,\lambda=1$.

Experimental setup

The experimental setup is as follows. The Starlike Annuli dataset consists of 50k points drawn uniformly at random from 7.2 with a 50/50 class balance. As in the Cifar10 experiment, we estimate generalization by performing 10-fold cross validation on the training set itself, for $N=[1,2,5,10,20,100,200]$. We use a weight decay of [math] and a momentum of $0.9$, and we use the weight initialization in [15]. These models are trained with a batch size of $128$ on a single GPU (RTX 2080ti). We start with a learning rate of $0.1$, divide it by $10$ at $80$, $120$, $160$ and by $2$ at $180$ epochs where an epoch consists of 45k/batchsize. Training is terminated at $200$ epochs. We ran one experiment with zero augmented dimensions and one experiment with one augmented dimension.

Results

The results are depicted in Figs. 4 and 5. We chose to display the Cross Entropy instead of accuracy for visual clarity. In the non-augmented case, results displayed in Fig. 4, it is interesting to note a few things. First of all the Cross-Entropy (Risk) is increasing with the layer count, indicating that deeper is worse for this problem. This is reasonable to expect as 7.2 is of the same type as 7.1 for the limit problem (NODE). Thus, it is non-representable, in the sense stated previously, within the NODE function class, [9]. Secondly the convergence rate of the risk is indeed very close to the theoretical bound $1/N$.

On the other hand, the one augmented dimension problem, see Fig. 5, tells a different story. For $N\geq 2$ we see no noticeable difference in performance, and we seem to gain nothing from the added depth. We interpret this as to indicate that the convergence rate is very fast and that we essentially reach the NODE performance using only a small value of $N$.

7.3. Conclusion

It is interesting to observe that for some problems the performance either increase or stay the same as depth increase, see Figs. 2 and 5, while for others the performance decrease significantly, see Fig. 4.

What is the reason for this? Looking at this from the context of regularity, we note that the deeper the network is, the closer we are to the NODE, Theorem 2.7. This forces the trajectories to become smoother and smoother, and we see that depth is a regularizer for the learning problem. On one hand, the idea of regularization is supported by the observations above, especially Fig. 4. On the other hand, for the experiments regarding Cifar10 we saw a similar increase in training accuracy for increasing depth, as we saw in the test accuracy, Fig. 2. This contradicts the idea of regularization. Putting these observations together we see that regularization alone cannot explain our numerical simulations. We believe that there is likely a complicated interplay between the effect that the depth has on the ‘loss landscape’, [26], and it’s regularizing effects.

8. Concluding remarks

In this paper we have established several convergence results, Theorems 2.5, 2.7 and 2.11, which together give a theoretical foundation for considering Neural ODEs (NODEs) as the deep limit of ResNets. The rate of convergence proved is essentially sharp as it is the order of the discretization error for the ODE. Our results can be seen as a complement to the result in [39], where the authors consider a related but different problem. However, in [39] the authors only consider the convergence of minimizers and not the convergence of the actual optimization procedure. They conjecture that minimizers also converge with the rate $1/N$, which is what we prove for the stochastic gradient descent. It should also be noted that we in this paper assume that the gradient of the NODE can be computed exactly. This is in contrast to [5] where they used an adjoint method to compute an approximate gradient. On the other hand, [5] contains no theoretical justifications. The error produced when applying the adjoint method to compute the gradient is not known and it is hard to tell how it would affect the convergence considered in Theorem 2.7.

Section 4 is of independent interest, and while the results of the section may be ‘folklore’ we have provided this section as it rigorously establishes the existence and uniqueness of stochastic gradient descent as well as the connection to the corresponding Fokker-Planck and linear Fokker-Planck equation. Using what is established in Section 4 we can also make the observation that, instead of considering the rate of convergence to minima of $V$ for

[TABLE]

we can consider the convergence of the density $p_{t}$ of $\theta_{t}$

[TABLE]

Using 4.12 we see that the corresponding solution $\widetilde{p}=e^{V}p$ satisfies

[TABLE]

According to the Poincaré inequality for measures that satisfy logarithmic Sobolev inequalities, see [41],

[TABLE]

using 4.13. Now, by Grönwall’s lemma we see that

[TABLE]

The conclusion is that $p_{t}\to p_{\infty}$ exponentially fast as $t\to\infty$. So instead of minimizing the possibly non-convex function $V$ we see that stochastic gradient descent is actually minimizing a convex problem in the sense of distributions. This could explain why stochastic gradient descent performs so well in practice.

In all our proofs, and in fact even for the existence of strong solutions to the stochastic gradient descent, Section 2, we need a penalization term. In particular, this is the case for neural networks with more than one layer, in this case we cannot simply guarantee from the neural network that a condition such as A3 holds. Instead, we have to rely on bounds concerning the absolute value of the drift, which grows polynomially, the order of the polynomial being the number of layers. For the Neural ODE it is clear that for the problem in 7.1 one would need $\|\theta\|\to\infty$ in order to produce a zero risk, and hence the penalization is necessary. However, the practical implications of such a penalization is essentially non-existent. This is because we can choose a penalization that is zero on a ball of arbitrary size, for instance $B(0,1.8\cdot 10^{308})$, which is the largest value for an IEEE 754 64–bit double.

Acknowledgment The authors were partially supported by a grant from the Göran Gustafsson Foundation for Research in Natural Sciences and Medicine.