# Implicit Regularization of Deep Residual Networks towards Neural ODEs

Pierre Marion, Yu-Han Wu††footnotemark:
Sorbonne Université, CNRS,
Laboratoire de Probabilités, Statistique et Modélisation,
F-75005 Paris, France
<EMAIL_ADDRESS>
<EMAIL_ADDRESS>Michael E. Sander
Ecole Normale Supérieure, CNRS,
Département de Mathématiques et Applications,
F-75005 Paris, France
<EMAIL_ADDRESS>Gérard Biau
Sorbonne Université, CNRS,
Laboratoire de Probabilités, Statistique et Modélisation,
F-75005 Paris, France
<EMAIL_ADDRESS>Equal contribution.

###### Abstract

Residual neural networks are state-of-the-art deep learning models. Their
continuous-depth analog, neural ordinary differential equations (ODEs), are
also widely used. Despite their success, the link between the discrete and
continuous models still lacks a solid mathematical foundation. In this
article, we take a step in this direction by establishing an implicit
regularization of deep residual networks towards neural ODEs, for nonlinear
networks trained with gradient flow. We prove that if the network is
initialized as a discretization of a neural ODE, then such a discretization
holds throughout training. Our results are valid for a finite training time,
and also as the training time tends to infinity provided that the network
satisfies a Polyak-Łojasiewicz condition. Importantly, this condition holds
for a family of residual networks where the residuals are two-layer
perceptrons with an overparameterization in width that is only linear, and
implies the convergence of gradient flow to a global minimum. Numerical
experiments illustrate our results.

##  1 Introduction

Residual networks are a successful family of deep learning models popularized
by breakthrough results in computer vision (He et al., 2016b). The key idea
behind residual networks, namely the presence of skip connections, is now
ubiquitous in deep learning, and can be found, for example, in Transformer
models (Vaswani et al., 2017). The main advantage of skip connections is to
allow successful training with depth of the order of a thousand layers, in
contrast to vanilla neural networks, leading to significant performance
improvement (e.g., Wang et al., 2022). This has motivated research on the
properties of residual networks in the limit where the depth tends to
infinity. One of the main explored directions is the neural ordinary
differential equation (ODE) limit (Chen et al., 2018).

Before presenting neural ODEs, we first introduce the mathematical formalism
of deep residual networks. We consider a single model throughout the paper to
simplify the exposition, but most of our results apply to more general models,
as will be discussed later. We consider the formulation

$h_{k+1}=h_{k}+\frac{1}{L\sqrt{m}}V_{k+1}\sigma\Big{(}\frac{1}{\sqrt{q}}W_{k+1}h_{k}\Big{)},\quad
k\in\\{0,\dots,L-1\\},$ (1)

where $L$ is the depth of the network, $h_{k}\in{\mathbb{R}}^{q}$ is the
output of the $k$-th hidden layer, $V_{k}\in{\mathbb{R}}^{q\times m}$,
$W_{k}\in{\mathbb{R}}^{m\times q}$ are the weights of the $k$-th layer, and
$\sigma$ is an activation function applied element-wise. Scaling with the
square root of the width is classical, although it often appears as an
equivalent condition on the variance at initialization (Glorot & Bengio, 2010;
LeCun et al., 2012; He et al., 2015). We make the scaling factors explicit to
have weights of magnitude $\mathcal{O}(1)$ independently of the width and the
depth. The $1/L$ scaling factor is less common, but it is necessary for the
correspondence with neural ODEs to hold. More precisely, if there exist
Lipschitz continuous functions $\mathcal{V}$ and $\mathcal{W}$ such that
$V_{k}=\mathcal{V}(k/L)$ and $W_{k}=\mathcal{W}(k/L)$, then the residual
network (1) converges, as $L\to\infty$, to the ODE

$\frac{dH}{ds}(s)=\frac{1}{\sqrt{m}}\mathcal{V}(s)\sigma\Big{(}\frac{1}{\sqrt{q}}\mathcal{W}(s)H(s)\Big{)},\quad
s\in[0,1],$ (2)

where $s$ is the continuous-depth version of the layer index. It is important
to note that this correspondence holds for fixed limiting functions
$\mathcal{V}$ and $\mathcal{W}$. This is especially true at initialization,
for example by setting the $V_{k}$ to zero and the $W_{k}$ to Gaussian
matrices weight-tied across the depth. The initial residual network is then
trivially equal to the neural ODE $\frac{dH}{ds}(s)=0$. Of course, more
sophisticated initializations are possible, as shown, e.g., in Marion et al.
(2022); Sander et al. (2022b). However, regardless of an ODE structure at
initialization, a more challenging question is that of the structure of the
network during and after training. Since the weights are updated during
training, there is a priori no guarantee that an ODE limit still holds after
training, even if it does at initialization.

The question of a potential ODE structure for the trained network is not a
mere technical one. In fact, it is important for at least three reasons.
First, it gives a precise answer to the question of the connection between
(trained) residual networks and neural ODEs, providing more solid ground to a
common statement in the community that both can coincide in the large-depth
limit (see, e.g., Haber & Ruthotto, 2017; E et al., 2019; Dong et al., 2020;
Massaroli et al., 2020; Kidger, 2022). Second, it opens exciting perspectives
for understanding residual networks. Indeed, if trained residual networks are
discretizations of neural ODEs, then it is possible to apply results from
neural ODEs to the large family of residual networks. In particular, from a
theoretical point of view, the approximation capabilities of neural ODEs are
well understood (Teshima et al., 2020; Zhang et al., 2020) and it is
relatively easy to obtain generalization bounds for these models (Hanson &
Raginsky, 2022; Marion, 2023). From a practical standpoint, advantages of
neural ODEs include memory-efficient training (Chen et al., 2018; Sander et
al., 2022b) and weight compression (Queiruga et al., 2021). This is important
because in practice memory is a bottleneck for training residual networks
(Gomez et al., 2017). Finally, our analysis is a first step towards
understanding the implicit regularization (Neyshabur et al., 2014; Vardi,
2023) of gradient descent for deep residual networks, that is, characterizing
the properties of the trained network among all minimizers of the empirical
risk.

Throughout the document, it is assumed that the network is trained with
gradient flow, which is a continuous analog of gradient descent. The
parameters $V_{k}$ are updated according to an ODE of the form
$\frac{dV_{k}}{dt}(t)=-L\frac{\partial\ell}{\partial V_{k}}(t)$ for
$t\geqslant 0$, where $\ell$ is an empirical risk (the exact mathematical
context and assumptions are detailed in Section 3), and similarly for $W_{k}$.
The scaling factor $L$ is the counterpart of the factor $1/L$ in (1), and
prevents vanishing gradients as $L$ tends to infinity. Note that the gradient
flow is defined with respect to a time index $t$ different from the layer
index $s$.

##### Contributions.

Our first main contribution (Section 4.1) is to show that a neural ODE limit
holds after training up to time $t$, i.e., there exists a function
$\mathcal{V}(s,t)$ such that the residual network converges, as $L$ tends to
infinity, to the ODE

$\frac{dH}{ds}(s)=\frac{1}{\sqrt{m}}\mathcal{V}(s,t)\sigma\Big{(}\frac{1}{\sqrt{q}}\mathcal{W}(s,t)H(s)\Big{)},\quad
s\in[0,1].$

This large-depth limit holds for any finite training time $t\geqslant 0$.
However, the convergence of the optimization algorithm as $t$ tends to
infinity, which we refer to as the long-time limit to distinguish it from the
large-depth limit $L\to\infty$, is not guaranteed without further assumptions,
due to the non-convexity of the optimization problem. We attack the question
(Section 4.2) when the width is large enough by proving a Polyak-Łojasiewicz
(PL) condition, which is now state of the art in analyzing the properties of
optimization algorithms for deep neural networks (Liu et al., 2022). The main
assumption for our PL condition to hold is that the width $m$ of the hidden
layers should be greater than some constant times the number of data $n$. As a
second main contribution, we show that the PL condition yields the long-time
convergence of the gradient flow for residual networks with linear
overparameterization. Finally, we prove the convergence with high probability
in the long-time limit, namely the existence of functions
$\mathcal{V}_{\infty}$ and $\mathcal{W}_{\infty}$ such that the discrete
trajectory defined by the trained residual network (1) converges as both $L$
and $t$ tend to infinity to the solution of the neural ODE (2) with
$\mathcal{V}=\mathcal{V}_{\infty}$ and $\mathcal{W}=\mathcal{W}_{\infty}$. In
addition, our approach points out that this limiting ODE interpolates the
training data. Finally, our results are illustrated by numerical experiments
(Section 5).

##  2 Related work

##### Deep residual networks and neural ODEs.

Several works study the large-depth convergence of residual networks to
differential equations, but without considering the training dynamics (Thorpe
& van Gennip, 2023; Cohen et al., 2021; Marion et al., 2022; Hayou, 2023).
Closer to our setting, Cont et al. (2022) and Sander et al. (2022b) analyze
the dynamics of gradient descent for deep residual networks, as we do, but
with significant differences. Cont et al. (2022) consider a
$\nicefrac{{1}}{{\sqrt{L}}}$ scaling factor in front of the residual branch,
resulting in a limit that is not a neural ODE. In addition, only $W$ is
trained. Furthermore, to obtain convergence in the long-time limit, it is
assumed that the data points are nearly orthogonal. Sander et al. (2022b)
prove the existence of an ODE limit for trained residual networks, but in the
simplified case of a linear activation and under a more restricted setting.

##### Long-time convergence of wide residual networks.

Polyak-Łojasiewicz conditions are a modern tool to prove long-time convergence
of overparameterized neural networks (Liu et al., 2022). These conditions are
a relaxation of convexity, and mean that the gradients of the loss with
respect to the parameters cannot be small when the loss is large. They have
been applied to residual networks with both linear (Bartlett et al., 2018; Wu
et al., 2019; Zou et al., 2020) and nonlinear activations (Allen-Zhu et al.,
2019; Frei et al., 2019; Barboni et al., 2022; Cont et al., 2022; MacDonald et
al., 2022). Building on the proof technique of Nguyen & Mondelli (2020) for
non-residual networks, we need only a linear overparameterization to prove our
PL condition, i.e., we require $m=\Omega(n)$. This compares favorably with
results requiring polynomial overparameterization (Allen-Zhu et al., 2019;
Barboni et al., 2022) or assumptions on the data, either a margin condition
(Frei et al., 2019) or a sample size smaller than the dimension of the data
space (Cont et al., 2022; MacDonald et al., 2022).

##### Implicit regularization.

Our paper can be related to a line of work on the implicit regularization of
gradient-based algorithms for residual networks (Neyshabur et al., 2014). We
show that the optimization algorithm does not just converge to any residual
network that minimizes the empirical risk, but rather to the discretization of
a neural ODE. Note that most implicit regularization results state that the
optimization algorithm converges to an interpolator that minimizes some
complexity measure, which can be a margin (Lyu & Li, 2020), a norm (Boursier
et al., 2022), or a matrix rank (Li et al., 2021). Thus, an interesting next
step is to understand if the neural ODE found by gradient flow actually
minimizes some complexity measure, and to characterize its generalization
properties.

##  3 Definitions and notation

This section is devoted to specifying the setup outlined in Section 1. Proofs
are given in the appendix.

##### Residual network.

A (scaled) residual network of depth $L\in\mathbb{N}^{*}$ is defined by

$\displaystyle h_{0}^{L}$ $\displaystyle=A^{L}x$ (3) $\displaystyle
h_{k+1}^{L}$
$\displaystyle=h_{k}^{L}+\frac{1}{L\sqrt{m}}V_{k+1}^{L}\sigma\Big{(}\frac{1}{\sqrt{q}}W_{k+1}^{L}h_{k}^{L}\Big{)},\quad
k\in\\{0,\dots,L-1\\},$ $\displaystyle F^{L}(x)$
$\displaystyle=B^{L}h_{L}^{L}.$

To allow the hidden layers $h_{k}^{L}\in{\mathbb{R}}^{q}$ to have a different
dimension than the input $x\in{\mathbb{R}}^{d}$, we first map $x$ to
$h_{0}^{L}$ with a weight matrix $A^{L}\in{\mathbb{R}}^{q\times d}$. We assume
that the hidden layers belong to a higher dimensional space than the input and
output, i.e., $q\geqslant\max(d,d^{\prime})$. The residual transformations are
two-layer perceptrons parameterized by the weight matrices
$V_{k}^{L}\in{\mathbb{R}}^{q\times m}$ and $W_{k}^{L}\in{\mathbb{R}}^{m\times
q}$. This is standard in the literature (e.g., He et al., 2016a; Dupont et
al., 2019; Barboni et al., 2022). The last weight matrix
$B^{L}\in{\mathbb{R}}^{d^{\prime}\times q}$ maps the last hidden layer to the
output $F^{L}(x)$ in ${\mathbb{R}}^{d^{\prime}}$. Also,
$\sigma:{\mathbb{R}}\to{\mathbb{R}}$ is an element-wise activation function
assumed to be $\mathcal{C}^{2}$, non-constant, Lipschitz continuous, bounded,
and such that $\sigma(0)=0$. The convenient shorthand
$Z_{k}^{L}=(V_{k}^{L},W_{k}^{L})$ is occasionally used, and we denote
$\|Z_{k}^{L}\|_{F}$ the sum of the Frobenius norms
$\|V_{k}^{L}\|_{F}+\|W_{k}^{L}\|_{F}$.

##### Data and loss.

The data is a sample of $n$ pairs $(x_{i},y_{i})_{1\leqslant i\leqslant
n}\in(\mathcal{X}\times\mathcal{Y})^{n}$ where $\mathcal{X}\times\mathcal{Y}$
is a compact set of $\mathbb{R}^{d}\times\mathbb{R}^{d^{\prime}}$. The
empirical risk is the mean squared error
$\ell^{L}=\frac{1}{n}\sum_{i=1}^{n}\|F^{L}(x_{i})-y_{i}\|^{2}$.

##### Initialization.

We initialize $A^{L}=(I_{{\mathbb{R}}^{d\times
d}},0_{{\mathbb{R}}^{(q-d)\times d}})$ as the identity matrix in
${\mathbb{R}}^{d\times d}$ concatenated row-wise with the zero matrix in
${\mathbb{R}}^{(q-d)\times d}$, to act as a simple projection of the input
onto the higher dimensional space ${\mathbb{R}}^{q}$, and similarly
$B^{L}=(0_{{\mathbb{R}}^{d^{\prime}\times(q-d^{\prime})}},I_{{\mathbb{R}}^{d^{\prime}\times
d^{\prime}}})$. The weights $V_{k}^{L}$ are initialized to zero and the
$W_{k}^{L}$ as weight-tied standard Gaussian matrices, i.e., for all
$k\in\\{1,\dots,L\\}$, $W_{k}^{L}=W\sim\mathcal{N}(0,1)^{\otimes(m\times q)}$.
Initializing outer matrices to zero is standard practice (Zhang et al., 2019),
while taking weight-tied matrices instead of i.i.d. ones is less common. We
show in Section 5 that it is still possible to learn with this initialization
scheme on real world data. As explained in Section 4.3, other initialization
choices are possible, provided they correspond to the discretization of a
Lipschitz continuous function, but we focus on this one in the main text for
simplicity.

##### Training algorithm.

Gradient flow is the limit of gradient descent as the learning rate tends to
zero. The parameters are set at time $t=0$ by the initialization, and then
evolve according to the ODE

$\frac{dA^{L}}{dt}(t)=-\frac{\partial\ell^{L}}{\partial
A^{L}}(t),\quad\frac{dZ_{k}^{L}}{dt}(t)=-L\frac{\partial\ell^{L}}{\partial
Z_{k}^{L}}(t),\quad\frac{dB^{L}}{dt}(t)=-\frac{\partial\ell^{L}}{\partial
B^{L}}(t),\quad t\geqslant 0,$ (4)

for $k\in\\{1,\dots,L\\}$. In the following, the dependence in $t$ is made
explicit when necessary, e.g., we write $h_{k}^{L}(t)$ instead of $h_{k}^{L}$,
and $F^{L}(x;t)$ instead of $F^{L}(x)$.

It turns out that, without further assumptions, the gradient flow can diverge
in finite time. This is because the dynamics are not (globally) Lipschitz
continuous, breaking the conditions of the Picard-Lindelöf theorem (see Lemma
16) for existence and uniqueness of ODE solutions. A common practice
(Goodfellow et al., 2016, Section 10.11.1) is to consider instead a clipped
gradient flow

$\frac{dA^{L}}{dt}(t)=\pi\Big{(}-\frac{\partial\ell^{L}}{\partial
A^{L}}(t)\Big{)},\quad\frac{dZ_{k}^{L}}{dt}(t)=\pi\Big{(}-L\frac{\partial\ell^{L}}{\partial
Z_{k}^{L}}(t)\Big{)},\quad\frac{dB^{L}}{dt}(t)=\pi\Big{(}-\frac{\partial\ell^{L}}{\partial
B^{L}}(t)\Big{)},$ (5)

where $\pi$ is a generic notation for a bounded Lipschitz continuous operator.
For example, clipping each coordinate of the gradient at some $C>0$ amounts to
taking $\pi$ as the projection on the ball centered at $0$ of radius $C$ for
the $\ell_{\infty}$ norm. Clipping ensures that the dynamics are well defined,
as shown in the next proposition that is a consequence of the Picard-Lindelöf
theorem.

######  Proposition 1.

The (clipped) gradient flow (5) has a unique solution for all $t\geqslant 0$.

In Section 4.2, we make additional assumptions to prove the long-time
convergence of the gradient flow. We then prove that these assumptions ensure
that the dynamics of the gradient flow (4) are well defined, eliminating the
need for clipping.

##### Neural ODE.

The neural ODE corresponding to the residual network (3) is defined by

$\displaystyle\begin{split}H(0)&=Ax\\\
\frac{dH}{ds}(s)&=\frac{1}{\sqrt{m}}\mathcal{V}(s)\sigma\Big{(}\frac{1}{\sqrt{q}}\mathcal{W}(s)H(s)\Big{)},\quad
s\in[0,1],\\\ F(x)&=BH(1),\end{split}$ (6)

where $x\in{\mathbb{R}}^{d}$ is the input, $H\in{\mathbb{R}}^{q}$ is the
variable of the ODE, $\mathcal{V}:[0,1]\to{\mathbb{R}}^{q\times m}$ and
$\mathcal{W}:[0,1]\to{\mathbb{R}}^{m\times q}$ are Lipschitz continuous
functions, $A\in{\mathbb{R}}^{q\times d}$ and
$B\in{\mathbb{R}}^{d^{\prime}\times q}$ are matrices, and the output is
$F(x)\in{\mathbb{R}}^{d^{\prime}}$. The following proposition shows that the
neural ODE is well defined. In addition, its output is close to the residual
network (3) provided the weights are discretizations of $\mathcal{V}$ and
$\mathcal{W}$.

######  Proposition 2.

The neural ODE (6) has a unique solution $H:[0,1]\to{\mathbb{R}}^{q}$.
Consider, moreover, the residual network (3) with $A^{L}=A$,
$V_{k}^{L}=\mathcal{V}(k/L)$ and $W_{k}^{L}=\mathcal{W}(k/L)$ for
$k\in\\{1,\dots,L\\}$, and $B^{L}=B$. Then there exists $C>0$ such that, for
all $L\in\mathbb{N}^{*}$,
$\sup_{x\in\mathcal{X}}\|F(x)-F^{L}(x)\|\leqslant\frac{C}{L}$.

Clearly, our choices of $V_{k}^{L}$ and $W_{k}^{L}$ at initialization are
discretizations of the Lipschitz continuous (in fact, constant) functions
$\mathcal{V}(s)\equiv 0$ and $\mathcal{W}(s)\equiv
W\sim\mathcal{N}(0,1)^{\otimes(m\times q)}$. Thus, Proposition 2 holds at
initialization, and the residual network is equivalent to the trivial ODE
$\frac{dH}{ds}(s)=0$. The next section shows that after training we obtain
non-trivial dynamics, which still discretize neural ODEs.

##  4 Large-depth limit of residual networks

We study the large-depth limit of trained residual networks in two settings.
In Section 4.1, we consider the case of a finite training time. We move in
Section 4.2 to the case where the training time tends to infinity, which is
tractable under a Polyak-Łojasiewicz condition. Proofs are given in the
appendix.

###  4.1 Clipped gradient flow and finite training time

We first consider the case where the neural network is trained with clipped
gradient flow (5) on some training time interval $[0,T]$, $T>0$. This allows
us to prove large-depth convergence to a neural ODE without further
assumptions. We emphasize that stopping training after a finite training time
is a common technique in practice, referred to as early stopping (Goodfellow
et al., 2016, Section 7.8). It is considered as a form of implicit
regularization, and our result sheds light on this intuition by showing that
the complexity of the trained networks increases with $T$.

The following proposition is a key step in proving the main theorem of this
section.

######  Proposition 3.

There exist $M,K>0$ such that, for any $t\in[0,T]$, $L\in\mathbb{N}^{*}$, and
$k\in\\{1,\dots,L\\}$,

$\max\big{(}\left\|A^{L}(t)\right\|_{F},\left\|V_{k}^{L}(t)\right\|_{F},\left\|W_{k}^{L}(t)\right\|_{F},\left\|B^{L}(t)\right\|_{F}\big{)}\leqslant
M,$

and, for $k\in\\{1,\dots,L-1\\}$,

$\max\big{(}\left\|V_{k+1}^{L}(t)-V_{k}^{L}(t)\right\|_{F},\left\|W_{k+1}^{L}(t)-W_{k}^{L}(t)\right\|_{F}\big{)}\leqslant\frac{K}{L}.$

Moreover, with probability at least $1-\exp\big{(}-\frac{3qm}{16}\big{)}$, the
following expressions hold for $M$ and $K$:

$M=TM_{\pi}+2\sqrt{qm},\quad K=\beta Te^{\alpha T},$ (7)

where $M_{\pi}$ is the supremum of $\pi$ in Frobenius norm, and $\alpha$ and
$\beta$ depend on $\mathcal{X}$, $\mathcal{Y}$, $M$, and $\sigma$.

This proposition ensures that the size of the weights and the difference
between successive weights remain bounded throughout training. We can now
state the main result, which states the convergence, for any training time in
$[0,T]$, of the neural network to a neural ODE as $L\to\infty$. Recall that a
sequence of functions $f^{L}$ of some variable $u$ is said to converge
uniformly over $u\in U$ to $f$ if $\sup_{u\in U}\|f^{L}(u)-f(u)\|\to 0$.

######  Theorem 4.

Consider the residual network (3) with the training dynamics (5). Then the
following statements hold as $L$ tends to infinity:

  1. (i)

There exist functions $A:[0,T]\to{\mathbb{R}}^{q\times d}$ and
$B:[0,T]\to{\mathbb{R}}^{d^{\prime}\times q}$ such that $A^{L}(t)$ and
$B^{L}(t)$ converge uniformly over $t\in[0,T]$ to $A(t)$ and $B(t)$.

  2. (ii)

There exists a Lipschitz continuous function
$\mathcal{Z}:[0,1]\times[0,T]\to{\mathbb{R}}^{q\times
m}\times\mathbb{R}^{m\times q}$ such that

$\mathcal{Z}^{L}:[0,1]\times[0,T]\to\mathbb{R}^{q\times
m}\times\mathbb{R}^{m\times q},\
(s,t)\mapsto\mathcal{Z}^{L}(s,t)=Z_{\left\lfloor(L-1)s\right\rfloor+1}^{L}(t)$
(8)

converges uniformly over $s\in[0,1]$ and $t\in[0,T]$ to
$\mathcal{Z}:=(\mathcal{V},\mathcal{W})$.

  3. (iii)

Uniformly over $s\in[0,1]$, $t\in[0,T]$, and $x\in\mathcal{X}$, the hidden
layer $h_{\left\lfloor Ls\right\rfloor}^{L}(t)$ converges to the solution at
time $s$ of the neural ODE

$\displaystyle\begin{split}H(0,t)&=A(t)x\\\ \frac{\partial H}{\partial
s}(s,t)&=\frac{1}{\sqrt{m}}\mathcal{V}(s,t)\sigma\Big{(}\frac{1}{\sqrt{q}}\mathcal{W}(s,t)H(s,t)\Big{)},\quad
s\in[0,1].\end{split}$ (9)

  4. (iv)

Uniformly over $t\in[0,T]$ and $x\in\mathcal{X}$, the output $F^{L}(x;t)$
converges to $B(t)H(1,t)$.

Let us sketch the proof of statement $(ii)$, which is the cornerstone of the
theorem. A first key idea is to introduce in (8) the piecewise-constant
continuous-depth interpolation $\mathcal{Z}^{L}$ of the weights, whose ambient
space does not depend on $L$, in contrast to the discrete weight sequence
$Z_{k}^{L}$. Since the weights remain bounded during training by Proposition
3, the Arzelà-Ascoli theorem guarantees the existence of an accumulation point
for $\mathcal{Z}^{L}$. We show that the accumulation point is unique because
it is the solution of an ODE satisfying the conditions of the Picard-Lindelöf
theorem. The uniqueness of the accumulation point then implies the existence
of a limit for the weights.

There are two notable byproducts of our proof. The first one is an explicit
description of the training dynamics of the limiting weights $A$, $B$, and
$\mathcal{Z}$, as the solution of an ODE system, as presented in Appendix A.5.
The second one, which we now describe, consists of norm bounds on the weights.
Proposition 3 bounds the discrete weights and the difference between two
consecutive weights respectively by some $M,K>0$. The proof of Theorem 4 shows
that this bound carries over to the continuous weights, in the sense that
$A(t)$, $\mathcal{V}(s,t)$, $\mathcal{W}(s,t)$, and $B(t)$ are uniformly
bounded by $M$, and $\mathcal{V}(\cdot,t)$ and $\mathcal{W}(\cdot,t)$ are
uniformly Lipschitz continuous with Lipschitz constant $K$. Formally, this
last property means that, for any $s,s^{\prime}\in[0,1]$ and $t\in[0,T]$,

$\|\mathcal{V}(s^{\prime},t)-\mathcal{V}(s,t)\|_{F}\leqslant
K|s^{\prime}-s|\quad\text{and}\quad\|\mathcal{W}(s^{\prime},t)-\mathcal{W}(s,t)\|_{F}\leqslant
K|s^{\prime}-s|.$

The boundedness and Lipschitz continuity of the weights are important features
because they limit the statistical complexity of neural ODEs (Marion, 2023).
More generally, norm-based bounds are a common approach in the statistical
theory of deep learning (see, e.g., Bartlett et al., 2017, and references
therein). Looking at the formula (7) for $M$ and $K$, one can see in
particular that the bounds diverge exponentially with $T$, providing an
argument in favor of early stopping.

Our approach so far characterizes the large-depth limit of the neural network
for a finite training time $T$, but two questions remain open. A first
challenge is to characterize the value of the loss after training. A second
one is to provide insight into the convergence of the optimization algorithm
in the long-time limit, i.e., as $T$ tends to infinity. To answer these
questions, we move to the setting where the width of the network is large
enough, which allows us to prove a Polyak-Łojasiewicz condition and thereby
the long-time convergence of the training loss to zero.

###  4.2 Convergence in the long-time limit for wide networks

Proving convergence of gradient-based optimization algorithms for neural
networks is a major difficulty in deep learning theory. One direction recently
explored considers sufficiently wide neural networks, and leverages the
Polyak-Łojasiewicz (PL) condition. In our setting, they are written as follows
(with the notation $Z^{L}=(V_{k}^{L},W_{k}^{L})_{k\in\\{1,\dots,L\\}}$):

######  Definition 1.

For $M,\mu>0$, the residual network (3) is said to satisfy the $(M,\mu)$-local
PL condition around a set of parameters
$(\bar{A}^{L},\bar{Z}^{L},\bar{B}^{L})$ if, for every set of parameters
$(A^{L},Z^{L},B^{L})$ such that

$\|A^{L}-\bar{A}^{L}\|_{F}\leqslant
M,\quad\sup_{k\in\\{1,\dots,L\\}}\|Z_{k}^{L}-\bar{Z}_{k}^{L}\|_{F}\leqslant
M,\quad\|B^{L}-\bar{B}^{L}\|_{F}\leqslant M,$

one has

$\Big{\|}\frac{\partial\ell^{L}}{\partial
A^{L}}\Big{\|}_{F}^{2}+L\sum_{k=1}^{L}\Big{\|}\frac{\partial\ell^{L}}{\partial
Z_{k}^{L}}\Big{\|}_{F}^{2}+\Big{\|}\frac{\partial\ell^{L}}{\partial
B^{L}}\Big{\|}_{F}^{2}\geqslant\mu\ell^{L},$

where the loss $\ell^{L}$ is evaluated at the set of parameters
$(A^{L},Z^{L},B^{L})$.

The next important point is to observe that, under the setup of Section 3 and
some additional assumptions, the residual network satisfies the local PL
condition of Definition 1.

######  Proposition 5.

Assume that the sample points $(x_{i},y_{i})$ are i.i.d. such that
$\|x_{i}\|_{2}=\sqrt{q}$. Then there exist $c_{1},\ldots,c_{4}>0$ (depending
only on $\sigma$) and $\delta>0$ such that, if

$q\geqslant d+d^{\prime},\quad m\geqslant c_{1}n,\quad L\geqslant
c_{2}\sqrt{nq},$

then, with probability at least $1-\delta$, the residual network (3) satisfies
the $(M,\mu)$-local PL condition around its initialization, with
$\displaystyle M=\frac{c_{3}}{\sqrt{nq}}$ and
$\displaystyle\mu=\frac{c_{4}}{n\sqrt{n}q}$.

We emphasize that Proposition 5 requires the width $m$ to scale only linearly
with the sample size $n$, which improves on the literature (see Section 2).
The other assumptions are mild. Note that our proof shows that the parameter
$\delta$ is small if $n$ grows at most polynomially with $d$ (see Appendix
B.5).

We are now ready to state convergence in the long-time and large-depth limits
to a global minimum of the empirical risk, when the local PL condition holds
and the norm of the targets $y_{i}$ is small enough.

######  Theorem 6.

Consider the residual network (3) with the training dynamics (4), and assume
that the assumptions of Proposition 5 hold. Then there exist $C,\delta>0$ such
that, if $\frac{1}{n}\sum_{i=1}^{n}\|y_{i}\|^{2}\leqslant C$, then, with
probability as least $1-\delta$, the gradient flow is well defined on
${\mathbb{R}}_{+}$, and, for $t\in{\mathbb{R}}_{+}$ and $L\in\mathbb{N}^{*}$,

$\ell^{L}(t)\leqslant\exp\big{(}-\frac{C^{\prime}t}{n\sqrt{n}q}\big{)}\ell^{L}(0),$
(10)

for some $C^{\prime}>0$ depending on $\sigma$. Moreover, the following
statements hold as $t$ and $L$ tend to infinity:

  1. (i)

There exist matrices $A_{\infty}\in{\mathbb{R}}^{q\times d}$ and
$B_{\infty}\in{\mathbb{R}}^{d^{\prime}\times q}$ such that $A^{L}(t)$ and
$B^{L}(t)$ converge to $A_{\infty}$ and $B_{\infty}$.

  2. (ii)

There exists a Lipschitz continuous function
$\mathcal{Z}_{\infty}:[0,1]\to{\mathbb{R}}^{q\times
m}\times\mathbb{R}^{m\times q}$ such that

$\mathcal{Z}^{L}:[0,1]\times{\mathbb{R}}_{+}\to\mathbb{R}^{q\times
m}\times\mathbb{R}^{m\times q},\
(s,t)\mapsto\mathcal{Z}^{L}(s,t)=Z_{\left\lfloor(L-1)s\right\rfloor+1}^{L}(t)$

converges uniformly over $s\in[0,1]$ to
$\mathcal{Z}_{\infty}:=(\mathcal{V}_{\infty},\mathcal{W}_{\infty})$.

  3. (iii)

Uniformly over $s\in[0,1]$ and $x\in\mathcal{X}$, the hidden layer
$h_{\left\lfloor Ls\right\rfloor}^{L}(t)$ converges to the solution at time
$s$ of the neural ODE

$\displaystyle\begin{split}H(0)&=A_{\infty}x\\\
\frac{dH}{ds}(s)&=\frac{1}{\sqrt{m}}\mathcal{V}_{\infty}(s)\sigma\Big{(}\frac{1}{\sqrt{q}}\mathcal{W}_{\infty}(s)H(s)\Big{)},\quad
s\in[0,1].\end{split}$

  4. (iv)

Uniformly over $x\in\mathcal{X}$, the output $F^{L}(x;t)$ converges to
$F_{\infty}(x)=B_{\infty}H(1)$. Furthermore, $F_{\infty}(x_{i})=y_{i}$ for all
$i\in\\{1,\dots,n\\}$.

This theorem proves two important results of separate interest. On the one
hand, equation (10) shows the long-time convergence of the gradient flow for
deep residual networks under the linear overparameterization assumption
$m\geqslant c_{1}n$ of Proposition 5. On the other hand, when both $t$ and $L$
tend to infinity, the network converges to a neural ODE that further
interpolates the training data. Note that the order in which $t$ and $L$ tend
to infinity does not matter by uniform convergence properties.

###  4.3 Generalizations to other architectures and initialization

To simplify the exposition, we have so far considered a particular residual
architecture defined in (3). However, most of our results hold for a more
general residual network of the form

$h_{k+1}^{L}=h_{k}^{L}+\frac{1}{L}f(h_{k}^{L},Z_{k+1}^{L}),\quad
k\in\\{0,\dots,L-1\\},$ (11)

where $f:{\mathbb{R}}^{q}\times{\mathbb{R}}^{p}\to{\mathbb{R}}^{q}$ is a
$\mathcal{C}^{2}$ function such that $f(0,\cdot)\equiv 0$ and $f(\cdot,z)$ is
uniformly Lipschitz for $z$ in any compact. All our results are shown in the
appendix for this general model, except the PL condition of Proposition 5,
which we prove only for the specific setup of Section 3. In particular, the
conclusions of Theorem 4 hold for the general model (11), as well as those of
Theorem 6 if the network satisfies a $(M,\mu)$-local PL condition with $\mu$
sufficiently large (see Appendix B for details).

It is easy to see that our residual network of interest (3) is a special case
of the general model (11) if $\sigma$ satisfies the assumptions of Section 3.
However, other choices are possible, such as convolutional layers or a
Lipschitz continuous version of Transformer (Kim et al., 2021). This latter
application is particularly interesting in the light of the literature
analyzing the Transformer architecture from a neural ODE point of view (Lu et
al., 2019; Sander et al., 2022a; Geshkovski et al., 2023).

Moreover, the initialization assumption made in Section 3 can also be relaxed
to include any so-called smooth initialization of the weights (Marion et al.,
2022). A smooth initialization corresponds to taking $V_{k}^{L}(0)$ and
$W_{k}^{L}(0)$ as discretizations of some Lipschitz continuous functions
$\mathcal{V}_{0}:[0,1]\to{\mathbb{R}}^{q\times m}$ and
$\mathcal{W}_{0}:[0,1]\to{\mathbb{R}}^{m\times q}$, that is, for
$k\in\\{1,\dots,L\\}$, $V_{k}^{L}(0)=\mathcal{V}_{0}(\frac{k}{L})$ and
$W_{k}^{L}(0)=\mathcal{W}_{0}(\frac{k}{L})$. A typical concrete example is to
let the entries of $\mathcal{V}_{0}$ and $\mathcal{W}_{0}$ be independent
Gaussian processes with expectation zero and squared exponential covariance
$K(x,x^{\prime})=\exp(-\frac{(x-x^{\prime})^{2}}{2\ell^{2}})$, for some
$\ell>0$. As shown by Proposition 2, a smooth initialization means that the
network discretizes a neural ODE.

##  5 Numerical experiments

We now present numerical experiments to validate our theoretical findings,
using both synthetic and real-world data. Experimental details are given in
Appendix E.

###  5.1 Synthetic data

We consider the residual network (3) with the initialization scheme of Section
3. So, the $V_{k}^{L}$ are initialized to zero and the $W_{k}^{L}$ to weight-
tied standard Gaussian matrices. To ease the presentation, we consider the
case where $q=d=d^{\prime}$, and we do not train the weights $A^{L}$ and
$B^{L}$, which therefore stay equal to the identity. The activation function
is GELU (Hendrycks & Gimpel, 2016), which is a smooth approximation of ReLU:
$x\mapsto\max(x,0)$. The sample points $(x_{i},y_{i})_{1\leqslant i\leqslant
n}$ follow independent standard Gaussian distributions. Note that it does not
hurt to take $x$ and $y$ independent since, in this subsection, our focus is
on optimization results only and not on statistical aspects. The mean-squared
error is minimized using full-batch gradient descent. The following
experiments exemplify the large-depth ($t\in[0,T]$, $L\to\infty$) and long-
time ($t\to\infty$, $L$ finite) limits.

Figure 1: Left: $1/L$ convergence of the maximum distance between two
successive weight matrices $\mathrm{max}_{1\leqslant k\leqslant
L,t\in[0,T]}(\|Z^{L}_{k}(t)-Z^{L}_{k+1}(t)\|_{F})$. Right: uniform convergence
of $\mathcal{Z}^{L}$ to its large-depth limit $\mathcal{Z}$. Here, for a
matrix-valued function $f$, $\|f\|$ denotes
$(\int_{0}^{1}\|f(s)\|^{2}_{F}ds)^{1/2}$.

##### Large-depth limit.

We illustrate key insights of Proposition 3 and Theorem 4, with $T=500$. In
Figure 1 (left), we plot the maximum distance between two successive weight
matrices, i.e., $\mathrm{max}_{1\leqslant k\leqslant
L,t\in[0,T]}(\|Z^{L}_{k}(t)-Z^{L}_{k+1}(t)\|_{F})$, for different values of
$L$ and training time $t\in[0,T]$. We observe a $1/L$ convergence rate, as
predicted by Proposition 3. Moreover, for a fixed $L$, the distance between
two successive weight matrices increases with the training time, however at a
much slower pace than the exponential upper bound on $K$ given in identity
(7). Figure 1 (right) depicts the uniform convergence of $\mathcal{Z}^{L}$ to
its large-depth limit $\mathcal{Z}$, illustrating statement $(ii)$ of Theorem
4. The function $\mathcal{Z}$ is computed using $\mathcal{Z}^{L}$ for
$L=2^{14}$. Note that the convergence is slower for larger training times.

##### Long-time limit.

We now turn to the long-time training setup, training for $80{,}000$
iterations with $L=64$. In Figure 2, we plot a specific (randomly-chosen)
coefficient of matrices $V_{k}^{L}$ and $W_{k}^{L}$ across layers, for
different training times. This illustrates Theorem 6 in a practical setting
since, visually, the weights behave as a Lipschitz continuous function for any
training time and converge to a Lipschitz continuous function as $t\to\infty$.
We also display the loss as a function of the training time, corroborating the
convergence of the loss to zero in Theorem 6.

Figure 2: Left: Randomly-chosen coefficient of the weight matrices across
layers ($x$-axis) for various training times $t$ (lighter color indicates
higher training time). Right: Loss against training time.

###  5.2 Real-world data

We now investigate the properties of deep residual networks on the CIFAR 10
dataset (Krizhevsky, 2009). We deviate from the mathematical model (3) by
using convolutions instead of fully connected layers. More precisely, $A^{L}$
is replaced by a trainable convolutional layer, and the residual layers write
$h_{k+1}^{L}=h_{k}^{L}+\frac{1}{L}\mathrm{bn}^{L}_{2,k}(\mathrm{conv}^{L}_{2,k}(\sigma(\mathrm{bn}^{L}_{1,k}(\mathrm{conv}^{L}_{1,k}(h_{k}^{L}))))),\quad
k\in\\{0,\dots,L-1\\},$ where $\mathrm{conv}^{L}_{i,k}$ are convolutions and
$\mathrm{bn}^{L}_{i,k}$ are batch normalizations. The output of the residual
layers is mapped to logits through a linear layer $B^{L}$. We initialize
$\mathrm{bn^{L}_{2,k}}$ to $0$, and $\mathrm{bn^{L}_{1,k}}$ and
$\mathrm{conv^{L}_{i,k}}$ either to weight-tied or to i.i.d. Gaussian. Table 1
reports the accuracy of the trained network, and whether it has Lipschitz
continuous (or smooth) weights after training, depending on the activation
function $\sigma$ and on the initialization scheme. To assess the smoothness
of the weights, we simply resort to visual inspection. For example, Figure 3
(left) shows two random coefficients of the convolutions across layers with
GELU and a weight-tied initialization: the smoothness is preserved after
training. Smooth weights indicate that the residual network discretizes a
neural ODE (see, e.g., Proposition 2). On the contrary, if an i.i.d.
initialization is used, smoothness is not preserved after training, as shown
in Figure 3 (right), and the residual network does not discretize a neural
ODE.

Figure 3: Random coefficients of the convolutions across layers ($x$-axis) after training. Left: Weight-tied initialization leads to smooth weights. Right: i.i.d. initialization leads to non-smooth weights. Act. function | Init. scheme | Train Acc. | Test Acc. | Smooth trained weights
---|---|---|---|---
Identity | Weight-tied | $56.5\pm 0.1$ | $59.8\pm 0.7$ | ✓
i.i.d. | $56.1\pm 0.3$ | $59.6\pm 0.7$ | $\bm{\times}$
GELU | Weight-tied | $80.5\pm 0.7$ | $79.9\pm 0.2$ | ✓
i.i.d. | $89.8\pm 0.5$ | $85.7\pm 0.1$ | $\bm{\times}$
ReLU | Weight-tied | $97.4\pm 0.6$ | $88.1\pm 0.1$ | $\bm{\times}$
i.i.d. | $98.4\pm 0.1$ | $88.4\pm 0.5$ | $\bm{\times}$
Table 1: Accuracy and smoothness of the trained weights depending on the
choice of activation function $\sigma$ and initialization scheme. We display
the median over $5$ runs and the interquartile range between the first and
third quantile. Smooth weights correspond to a neural ODE structure.

Table 1 conveys several important messages. First, in accordance with our
theory (Theorem 4), we obtain a neural ODE structure when using a smooth
activation function and weight-tied initialization (lines $1$ and $3$ of Table
1). This is not the case when using the non-smooth ReLU activation and/or
i.i.d. initialization. In fact, we prove in Appendix D that the smoothness of
the weights is lost when training with ReLU in a simple setting, confirming
this experimental observation. Furthermore, the third line of Table 1 shows
that it is possible to obtain a reasonable accuracy with a neural ODE
structure, which, as emphasized in Section 1, also comes with theoretical and
practical advantages. Nevertheless, we obtain an improvement in accuracy in
the cases corresponding to non-smooth weights, i.e., to a residual network
that does not discretize an ODE. Extending our theory to such cases is left
for future work.

##  6 Conclusion

We study the convergence of deep residual networks to neural ODEs. When
properly scaled and initialized, residual networks trained with fixed-horizon
gradient flow converge to neural ODEs as the depth tends to infinity. This
result holds for very general architectures. In the case where both training
time and depth tend to infinity, convergence holds under a local Polyak-
Łojasiewicz condition. We prove such a condition for a family of deep residual
networks with linear overparameterization.

The setting of neural ODE-like networks comes with strong guarantees, at the
cost of some performance gap when compared with i.i.d. initialization as
highlighted by the experimental section. Extending the mathematical large-
depth study to the i.i.d. case is an interesting problem for future research.
Previous work suggests that the correct limit object is then a stochastic
differential equation instead of an ODE (Cohen et al., 2021; Cont et al.,
2022; Marion et al., 2022).

#### Acknowledgments

P.M. is supported by a grant from Région Île-de-France and by MINES Paris -
PSL. P.M. and Y.-H.W. are funded by a Google PhD Fellowship. M.S. is supported
by the “Investissements d’avenir” program, reference ANR19-P3IA-0001, and by
the European Research Council (ERC project NORIA). This work was granted
access to the HPC resources of IDRIS under the allocation 2020-[AD011012073]
made by GENCI.

## References

  * Allen-Zhu et al. (2019) Zeyuan Allen-Zhu, Yuanzhi Li, and Zhao Song.  A convergence theory for deep learning via over-parameterization.  In K. Chaudhuri and R. Salakhutdinov (eds.), _Proceedings of the 36th International Conference on Machine Learning_ , volume 97, pp. 242–252. PMLR, 2019.
  * Arnold (1992) Vladimir I. Arnold.  _Ordinary Differential Equations_.  Springer, Berlin, 1992.
  * Barboni et al. (2022) Raphaël Barboni, Gabriel Peyré, and François-Xavier Vialard.  On global convergence of ResNets: From finite to infinite width using linear parameterization.  In S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh (eds.), _Advances in Neural Information Processing Systems_ , volume 35, pp. 16385–16397. Curran Associates, Inc., 2022.
  * Bartlett et al. (2017) Peter L. Bartlett, Dylan J. Foster, and Matus J. Telgarsky.  Spectrally-normalized margin bounds for neural networks.  In I. Guyon, U. Von Luxburg, S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett (eds.), _Advances in Neural Information Processing Systems_ , volume 30, pp. 6240–6249. Curran Associates, Inc., 2017\.
  * Bartlett et al. (2018) Peter L. Bartlett, Dave P. Helmbold, and Philip M. Long.  Gradient descent with identity initialization efficiently learns positive definite linear transformations by deep residual networks.  In J. Dy and A. Krause (eds.), _Proceedings of the 35th International Conference on Machine Learning_ , volume 80, pp. 521–530. PMLR, 2018.
  * Boursier et al. (2022) Etienne Boursier, Loucas Pillaud-Vivien, and Nicolas Flammarion.  Gradient flow dynamics of shallow ReLU networks for square loss and orthogonal inputs.  In S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh (eds.), _Advances in Neural Information Processing Systems_ , volume 35, pp. 20105–20118. Curran Associates, Inc., 2022.
  * Chen et al. (2018) Ricky T.Q. Chen, Yulia Rubanova, Jesse Bettencourt, and David K. Duvenaud.  Neural ordinary differential equations.  In S. Bengio, H. Wallach, H. Larochelle, K. Grauman, N. Cesa-Bianchi, and R. Garnett (eds.), _Advances in Neural Information Processing Systems_ , volume 31, pp. 6572–6583. Curran Associates, Inc., 2018.
  * Cohen et al. (2021) Alain-Sam Cohen, Rama Cont, Alain Rossier, and Renyuan Xu.  Scaling properties of deep residual networks.  In M. Meila and T. Zhang (eds.), _Proceedings of the 38th International Conference on Machine Learning_ , volume 139, pp. 2039–2048. PMLR, 2021.
  * Cont et al. (2022) Rama Cont, Alain Rossier, and RenYuan Xu.  Convergence and implicit regularization properties of gradient descent for deep residual networks.  _arXiv:2204.07261_ , 2022.
  * Debnath & Bhatta (2014) Lokenath Debnath and Dambaru Bhatta.  _Integral Transforms and Their Applications_.  CRC press, Boca Raton, third edition, 2014.
  * Dong et al. (2020) Chengyu Dong, Liyuan Liu, Zichao Li, and Jingbo Shang.  Towards adaptive residual network training: A neural-ODE perspective.  In H. Daumé III and A. Singh (eds.), _Proceedings of the 37th International Conference on Machine Learning_ , volume 119, pp. 2616–2626. PMLR, 2020.
  * Dragomir (2003) Sever S. Dragomir.  _Some Gronwall Type Inequalities and Applications_.  Nova Science Publishers, 2003.
  * Dupont et al. (2019) Emilien Dupont, Arnaud Doucet, and Yee Whye Teh.  Augmented neural ODEs.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett (eds.), _Advances in Neural Information Processing Systems_ , volume 32, pp. 3140–3150. Curran Associates, Inc., 2019.
  * E et al. (2019) Weinan E, Jiequn Han, and Qianxiao Li.  A mean-field optimal control formulation of deep learning.  _Research in the Mathematical Sciences_ , 6:10, 2019.
  * Filippov (1988) Aleksei F. Filippov.  _Differential Equations with Discontinuous Righthand Sides_.  Springer, Dordrecht, 1988.
  * Frei et al. (2019) Spencer Frei, Yuan Cao, and Quanquan Gu.  Algorithm-dependent generalization bounds for overparameterized deep residual networks.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett (eds.), _Advances in Neural Information Processing Systems_ , volume 32, pp. 14769–14779. Curran Associates, Inc., 2019.
  * Geshkovski et al. (2023) Borjan Geshkovski, Cyril Letrouit, Yury Polyanskiy, and Philippe Rigollet.  The emergence of clusters in self-attention dynamics.  _arXiv:2305.05465_ , 2023.
  * Glorot & Bengio (2010) Xavier Glorot and Yoshua Bengio.  Understanding the difficulty of training deep feedforward neural networks.  In Y.W. Teh and M. Titterington (eds.), _Proceedings of the Thirteenth International Conference on Artificial Intelligence and Statistics_ , volume 9, pp. 249–256. PMLR, 2010.
  * Goel et al. (2020) Surbhi Goel, Aravind Gollakota, Zhihan Jin, Sushrut Karmalkar, and Adam Klivans.  Superpolynomial lower bounds for learning one-layer neural networks using gradient descent.  In H. Daumé III and A. Singh (eds.), _Proceedings of the 37th International Conference on Machine Learning_ , volume 119, pp. 3587–3596. PMLR, 2020.
  * Gomez et al. (2017) Aidan N. Gomez, Mengye Ren, Raquel Urtasun, and Roger B. Grosse.  The reversible residual network: Backpropagation without storing activations.  In I. Guyon, U. Von Luxburg, S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett (eds.), _Advances in Neural Information Processing Systems_ , volume 30, pp. 2214–2224. Curran Associates, Inc., 2017\.
  * Goodfellow et al. (2016) Ian Goodfellow, Yoshua Bengio, and Aaron Courville.  _Deep Learning_.  MIT Press, 2016.
  * Haber & Ruthotto (2017) Eldad Haber and Lars Ruthotto.  Stable architectures for deep neural networks.  _Inverse Problems_ , 34:014004, 2017.
  * Hanson & Raginsky (2022) Joshua Hanson and Maxim Raginsky.  Fitting an immersed submanifold to data via Sussmann’s orbit theorem.  In _2022 IEEE 61st Conference on Decision and Control (CDC)_ , pp. 5323–5328, 2022.
  * Hayou (2023) Soufiane Hayou.  On the infinite-depth limit of finite-width neural networks.  _Transactions on Machine Learning Research_ , 2023.
  * He et al. (2015) Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Delving deep into rectifiers: Surpassing human-level performance on ImageNet classification.  In _Proceedings of 2015 IEEE Conference on Computer Vision and Pattern Recognition_ , pp. 1026–1034. IEEE, 2015.
  * He et al. (2016a) Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Identity mappings in deep residual networks.  In B. Leibe, J. Matas, N. Sebe, and M. Welling (eds.), _Computer Vision – ECCV 2016_ , pp. 630–645. Springer International Publishing, 2016a.
  * He et al. (2016b) Kaiming He, Xiangyu Zhang, Shaoqing Ren, and Jian Sun.  Deep residual learning for image recognition.  In _Proceedings of 2016 IEEE Conference on Computer Vision and Pattern Recognition_ , pp. 770–778. IEEE, 2016b.
  * Hendrycks & Gimpel (2016) Dan Hendrycks and Kevin Gimpel.  Gaussian Error Linear Units (GELUs).  _arXiv:1606.08415_ , 2016.
  * Kidger (2022) Patrick Kidger.  _On Neural Ordinary Differential Equations_.  PhD thesis, University of Oxford, 2022.
  * Kim et al. (2021) Hyunjik Kim, George Papamakarios, and Andriy Mnih.  The Lipschitz constant of self-attention.  In M. Meila and T. Zhang (eds.), _Proceedings of the 38th International Conference on Machine Learning_ , volume 139, pp. 5562–5571. PMLR, 2021.
  * Krizhevsky (2009) A. Krizhevsky.  Learning multiple layers of features from tiny images.  Technical report, University of Toronto, 2009.
  * Laurent & Massart (2000) Béatrice Laurent and Pascal Massart.  Adaptive estimation of a quadratic functional by model selection.  _The Annals of Statistics_ , 28:1302–1338, 2000.
  * LeCun et al. (2012) Yann A. LeCun, Léon Bottou, Genevieve B. Orr, and Klaus-Robert Müller.  Efficient backprop.  In G. Montavon, G.B Orr, and K.-R. Müller (eds.), _Neural Networks: Tricks of the Trade: Second Edition_ , pp. 9–48, Berlin, 2012. Springer.
  * Li et al. (2021) Zhiyuan Li, Yuping Luo, and Kaifeng Lyu.  Towards resolving the implicit bias of gradient descent for matrix factorization: Greedy low-rank learning.  In _International Conference on Learning Representations_ , 2021.
  * Liu et al. (2022) Chaoyue Liu, Libin Zhu, and Mikhail Belkin.  Loss landscapes and optimization in over-parameterized non-linear systems and neural networks.  _Applied and Computational Harmonic Analysis_ , 59:85–116, 2022.
  * Loyka (2015) Sergey Loyka.  On singular value inequalities for the sum of two matrices.  _arXiv:1507.06630_ , 2015.
  * Lu et al. (2019) Yiping Lu, Zhuohan Li, Di He, Zhiqing Sun, Bin Dong, Tao Qin, Liwei Wang, and Tie-Yan Liu.  Understanding and improving transformer from a multi-particle dynamic system point of view.  _arXiv:1906.02762_ , 2019.
  * Luk (2017) Jonathan Luk.  Notes on existence and uniqueness theorems for ODEs, 2017.  URL: http://web.stanford.edu/~jluk/math63CMspring17/Existence.170408.pdf.
  * Lyu & Li (2020) Kaifeng Lyu and Jian Li.  Gradient descent maximizes the margin of homogeneous neural networks.  In _International Conference on Learning Representations_ , 2020.
  * MacDonald et al. (2022) Lachlan E. MacDonald, Hemanth Saratchandran, Jack Valmadre, and Simon Lucey.  A global analysis of global optimisation.  _arXiv:2210.05371_ , 2022.
  * Marion (2023) Pierre Marion.  Generalization bounds for neural ordinary differential equations and deep residual networks.  _arXiv:2305.06648_ , 2023.
  * Marion et al. (2022) Pierre Marion, Adeline Fermanian, Gérard Biau, and Jean-Philippe Vert.  Scaling ResNets in the large-depth regime.  _arXiv:2206.06929_ , 2022.
  * Massaroli et al. (2020) Stefano Massaroli, Michael Poli, Jinkyoo Park, Atsushi Yamashita, and Hajime Asama.  Dissecting neural ODEs.  In H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin (eds.), _Advances in Neural Information Processing Systems_ , volume 33, pp. 3952–3963. Curran Associates, Inc., 2020.
  * Neyshabur et al. (2014) Behnam Neyshabur, Ryota Tomioka, and Nathan Srebro.  In search of the real inductive bias: On the role of implicit regularization in deep learning.  _arXiv:1412.6614_ , 2014.
  * Nguyen & Mondelli (2020) Quynh N. Nguyen and Marco Mondelli.  Global convergence of deep networks with one wide layer followed by pyramidal topology.  In H. Larochelle, M. Ranzato, R. Hadsell, M.F. Balcan, and H. Lin (eds.), _Advances in Neural Information Processing Systems_ , volume 33, pp. 11961–11972. Curran Associates, Inc., 2020.
  * Paszke et al. (2019) Adam Paszke, Sam Gross, Francisco Massa, Adam Lerer, James Bradbury, Gregory Chanan, Trevor Killeen, Zeming Lin, Natalia Gimelshein, Luca Antiga, Alban Desmaison, Andreas Kopf, Edward Yang, Zachary DeVito, Martin Raison, Alykhan Tejani, Sasank Chilamkurthy, Benoit Steiner, Lu Fang, Junjie Bai, and Soumith Chintala.  PyTorch: An imperative style, high-performance deep learning library.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett (eds.), _Advances in Neural Information Processing Systems_ , volume 32, pp. 8024–8035. Curran Associates, Inc., 2019.
  * Queiruga et al. (2021) Alejandro Queiruga, N. Benjamin Erichson, Liam Hodgkinson, and Michael W. Mahoney.  Stateful ODE-nets using basis function expansions.  In M. Ranzato, A. Beygelzimer, Y. Dauphin, P.S. Liang, and J. Wortman Vaughan (eds.), _Advances in Neural Information Processing Systems_ , volume 34, pp. 21770–21781. Curran Associates, Inc., 2021.
  * Sander et al. (2022a) Michael E. Sander, Pierre Ablin, Mathieu Blondel, and Gabriel Peyré.  Sinkformers: Transformers with doubly stochastic attention.  In G. Camps-Valls, F.J.R. Ruiz, and I. Valera (eds.), _Proceedings of The 25th International Conference on Artificial Intelligence and Statistics_ , volume 151, pp. 3515–3530. PMLR, 2022a.
  * Sander et al. (2022b) Michael E. Sander, Pierre Ablin, and Gabriel Peyré.  Do residual neural networks discretize neural ordinary differential equations?  In S. Koyejo, S. Mohamed, A. Agarwal, D. Belgrave, K. Cho, and A. Oh (eds.), _Advances in Neural Information Processing Systems_ , volume 35, pp. 36520–36532. Curran Associates, Inc., 2022b.
  * Teshima et al. (2020) Takeshi Teshima, Koichi Tojo, Masahiro Ikeda, Isao Ishikawa, and Kenta Oono.  Universal approximation property of neural ordinary differential equations.  _arXiv:2012.02414_ , 2020.
  * Thorpe & van Gennip (2023) Matthew Thorpe and Yves van Gennip.  Deep limits of residual neural networks.  _Research in the Mathematical Sciences_ , 10:6, 2023.
  * Tropp (2012) Joel A. Tropp.  User-friendly tail bounds for sums of random matrices.  _Foundations of Computational Mathematics_ , 12:389–434, 2012.
  * Vardi (2023) Gal Vardi.  On the implicit bias in deep-learning algorithms.  _Communications of the ACM_ , 66:86–93, 2023.
  * Vaswani et al. (2017) Ashish Vaswani, Noam Shazeer, Niki Parmar, Jakob Uszkoreit, Llion Jones, Aidan N. Gomez, Łukasz Kaiser, and Illia Polosukhin.  Attention is all you need.  In I. Guyon, U. Von Luxburg, S. Bengio, H. Wallach, R. Fergus, S. Vishwanathan, and R. Garnett (eds.), _Advances in Neural Information Processing Systems_ , volume 30, pp. 5998–6008. Curran Associates, Inc., 2017\.
  * Vershynin (2018) Roman Vershynin.  _High-Dimensional Probability: An Introduction with Applications in Data Science_.  Cambridge Series in Statistical and Probabilistic Mathematics. Cambridge University Press, Cambridge, 2018.
  * Wang et al. (2022) Hongyu Wang, Shuming Ma, Li Dong, Shaohan Huang, Dongdong Zhang, and Furu Wei.  DeepNet: Scaling transformers to 1,000 layers.  _arXiv:2203.00555_ , 2022.
  * Wu et al. (2019) Lei Wu, Qingcan Wang, and Chao Ma.  Global convergence of gradient descent for deep linear residual networks.  In H. Wallach, H. Larochelle, A. Beygelzimer, F. d'Alché-Buc, E. Fox, and R. Garnett (eds.), _Advances in Neural Information Processing Systems_ , volume 32, pp. 13389–13398. Curran Associates, Inc., 2019.
  * Zhang et al. (2020) Han Zhang, Xi Gao, Jacob Unterman, and Tom Arodz.  Approximation capabilities of neural ODEs and invertible residual networks.  In H. Daumé III and A. Singh (eds.), _Proceedings of the 37th International Conference on Machine Learning_ , volume 119, pp. 11086–11095. PMLR, 2020.
  * Zhang et al. (2019) Hongyi Zhang, Yann N. Dauphin, and Tengyu Ma.  Fixup initialization: Residual learning without normalization.  In _International Conference on Learning Representations_ , 2019.
  * Zou et al. (2020) Difan Zou, Philip M. Long, and Quanquan Gu.  On the global convergence of training deep linear ResNets.  In _International Conference on Learning Representations_ , 2020.

Appendix

##### Organization of the Appendix.

In Section A, we give some results on the general residual network (11). In
Section B, these results are instantiated in the specific case of the residual
network (3), thus proving the results of the paper. Section C contains some
lemmas that are useful for the proofs. We present in Section D a counter-
example showing that a residual network with the ReLU activation can move away
from the neural ODE structure during training. Finally, Section E presents
some experimental details.

##  Appendix A Some results for general residual networks

##### Lipschitz continuity.

Let $(\mathscr{U},\|\cdot\|)$, $(\mathscr{V},\|\cdot\|)$, and
$(\mathscr{W},\|\cdot\|)$ be generic normed spaces. Then a function of two
variables $g:\mathscr{U}\times\mathscr{V}\to\mathscr{W}$ is:

  1. $(i)$

(Globally) Lipschitz continuous if there exists $K\geqslant 0$ such that, for
$(u,v),(u^{\prime},v^{\prime})\in\mathscr{U}\times\mathscr{V}$,

$\|g(u,v)-g(u^{\prime},v)\|\leqslant K\|u-u^{\prime}\|+K\|v-v^{\prime}\|.$

  2. $(ii)$

Locally Lipschitz continuous in its first variable if, for any compacts
$E\subset\mathscr{U},E^{\prime}\subset\mathscr{V}$, there exists $K\geqslant
0$ such that, for $(u,v),(u^{\prime},v)\in E\times E^{\prime}$,

$\|g(u,v)-g(u^{\prime},v)\|\leqslant K\|u-u^{\prime}\|.$

Equivalent definitions hold for a function of one variable. Moreover,
$g(\cdot,v)$ is said to be uniformly Lipschitz continuous for $v$ in
$\mathscr{V}$ if there exists $K\geqslant 0$ such that, for
$(u,v),(u^{\prime},v)\in\mathscr{U}\times\mathscr{V}$,

$\|g(u,v)-g(u^{\prime},v)\|\leqslant K\|u-u^{\prime}\|,$

and uniformly Lipschitz continuous for $v$ in any compact if, for any compact
$E^{\prime}\subset\mathscr{V}$, there exists $K\geqslant 0$ such that, for
$(u,v),(u^{\prime},v)\in\mathcal{U}\times E^{\prime}$,

$\|g(u,v)-g(u^{\prime},v)\|\leqslant K\|u-u^{\prime}\|.$

Throughout, we refer to a Lipschitz continuous function with Lipschitz
constant $K\geqslant 0$ as $K$-Lipschitz.

##### Model.

As explained in Section 4.3, most of our results are proven for the general
residual network

$\displaystyle h_{0}^{L}(t)$ $\displaystyle=A^{L}(t)x$ $\displaystyle
h_{k+1}^{L}(t)$
$\displaystyle=h_{k}^{L}(t)+\frac{1}{L}f(h_{k}^{L}(t),Z_{k+1}^{L}(t)),\quad
k\in\\{0,\dots,L-1\\},$ (12) $\displaystyle F^{L}(x;t)$
$\displaystyle=B^{L}(t)h_{L}^{L}(t),$

where $Z^{L}(t)=(Z_{1}^{L}(t),\dots,Z_{L}^{L}(t))\in({\mathbb{R}}^{p})^{L}$
and $f:{\mathbb{R}}^{q}\times{\mathbb{R}}^{p}\to{\mathbb{R}}^{q}$ is a
$\mathcal{C}^{2}$ function such that $f(0,\cdot)\equiv 0$ and $f(\cdot,z)$ is
uniformly Lipschitz for $z$ in any compact. Let us introduce the
backpropagation equations, which are instrumental in the study of the gradient
flow dynamics. These equations define the backward state
$p_{k}^{L}(t)\in{\mathbb{R}}^{q}$ through the backward recurrence

$\displaystyle p_{L}^{L}(t)$ $\displaystyle=2B^{L}(t)^{\top}(F^{L}(x;t)-y)$
$\displaystyle p_{k}^{L}(t)$
$\displaystyle=p_{k+1}^{L}(t)+\frac{1}{L}\partial_{1}f(h_{k}^{L}(t),Z_{k+1}^{L}(t))p_{k+1}^{L}(t),\quad
k\in\\{0,\dots,L-1\\},$ (13)

where $\partial_{1}f\in{\mathbb{R}}^{q\times q}$ stands for the Jacobian
matrix of $f$ with respect to its first argument. Similarly, we let
$\partial_{2}f\in{\mathbb{R}}^{q\times p}$ be the Jacobian matrix of $f$ with
respect to its second argument. For a sample $(x_{i},y_{i})_{1\leqslant
i\leqslant n}\in(\mathcal{X}\times\mathcal{Y})^{n}$, we let $h_{k,i}^{L}(t)$
and $p_{k,i}^{L}(t)$ be, respectively, the hidden layer $h_{k}^{L}(t)$ and the
backward state $p_{k}^{L}(t)$ associated with the $i$-th input $x_{i}$.
Denoting the mean squared error associated with the sample by $\ell^{L}$, we
have, by the chain rule,

$\displaystyle\frac{\partial\ell^{L}}{\partial A^{L}}(t)$
$\displaystyle=\frac{1}{n}\sum_{i=1}^{n}p_{0,i}^{L}(t)x_{i}^{\top}$ (14)
$\displaystyle\frac{\partial\ell^{L}}{\partial Z_{k}^{L}}(t)$
$\displaystyle=\frac{1}{nL}\sum_{i=1}^{n}\partial_{2}f(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))^{\top}p_{k-1,i}^{L}(t),\quad
k\in\\{1,\dots,L\\},$ (15) $\displaystyle\frac{\partial\ell^{L}}{\partial
B^{L}}(t)$
$\displaystyle=\frac{2}{n}\sum_{i=1}^{n}(F^{L}(x_{i};t)-y_{i})h_{L,i}^{L}(t)^{\top}.$
(16)

##### Initialization.

The parameters $(Z_{k}^{L}(t))_{1\leqslant k\leqslant L}$ are initialized to
$Z_{k}^{L}(0)=Z^{\textnormal{init}}\big{(}\frac{k}{L}\big{)}$, where
$Z^{\textnormal{init}}:[0,1]\to{\mathbb{R}}^{p}$ is a Lipschitz continuous
function. Furthermore, we initialize $A^{L}(0)$ to some matrix
$A^{\textnormal{init}}\in{\mathbb{R}}^{q\times d}$ and
$B^{L}(0)=B^{\textnormal{init}}\in{\mathbb{R}}^{d^{\prime}\times q}$. Note
that this initialization scheme is a generalization of the one presented in
Section 3.

##### Additional notation.

For a vector $x$, $\left\|x\right\|$ denotes the Euclidean norm. For a matrix
$A$, the operator norm induced by the Euclidean norm is denoted by
$\left\|A\right\|_{2}$, and the Frobenius norm is denoted by
$\left\|A\right\|_{F}$. Finally, we use the notation $A^{L}$ (resp.
$Z_{k}^{L}$, $B^{L}$) to denote the function $t\mapsto A^{L}(t)$ (resp.
$t\mapsto Z_{k}^{L}(t)$, $t\mapsto B^{L}(t)$), since the parameters are
considered as functions of the training time throughout this appendix.

##### Overview of Appendix A.

First, in Section A.1, we study the case of the (clipped) gradient flow (5).
We show that the weights and the difference between successive weights are
bounded during the entire training. Section A.2 shows a similar result for the
standard gradient flow (4) under a PL condition. In Section A.3, we show a
generalized version of the Arzelà-Ascoli theorem, which allows us to prove the
existence of a converging subsequence of the weights in the large-depth limit.
Section A.4 is devoted to the convergence of the Euler scheme for
parameterized ODEs. We then proceed to prove in Section A.5 our main result,
i.e., the large-depth convergence of the gradient flow. The key step is to
establish the uniqueness of the adherence point of the weights. Finally, in
Section A.6, we prove the existence of a double limit for the weights and the
hidden states when both the depth and the training time tend to infinity.

###  A.1 The trained weights are bounded in the finite training-time setup

Before stating the result, let us introduce the notation
$\partial_{22}f(h,z)\in{\mathbb{R}}^{q\times p\times p}$, which is the third-
order tensor of second partial derivatives of $f$ with respect to $z$. We
endow the space ${\mathbb{R}}^{q\times p\times p}$ with the operator norm
$\|\cdot\|_{2}$ induced by the Euclidean norm in ${\mathbb{R}}^{p}$ and the
$\|\cdot\|_{2}$ norm in ${\mathbb{R}}^{q\times p}$. In other words,

$\|\partial_{22}f(h,z)\|_{2}=\sup_{u\in{\mathbb{R}}^{p},\|u\|=1}\|\partial_{22}f(h,z)u\|_{2},$

where $\partial_{22}f(h,z)u\in{\mathbb{R}}^{q\times p}$ is the tensor product
of $\partial_{22}f(h,z)$ against $u$. Similarly,
$\partial_{21}f(h,z)\in{\mathbb{R}}^{q\times p\times q}$ denotes the third-
order tensor of cross second partial derivatives of $f$, and the space
${\mathbb{R}}^{q\times p\times q}$ is endowed with the operator norm
$\|\cdot\|_{2}$ induced by the Euclidean norm in ${\mathbb{R}}^{q}$ and the
$\|\cdot\|_{2}$ norm in ${\mathbb{R}}^{q\times p}$.

######  Proposition 7.

Consider the residual network (12) initialized as explained in Appendix A and
trained with the gradient flow (5) on $[0,T]$, for some $T\in(0,\infty)$. Let

$\displaystyle M_{\pi}$
$\displaystyle=\max\Big{(}\max_{A\in{\mathbb{R}}^{q\times
d}}\|\pi(A)\|_{F},\,\max_{Z\in{\mathbb{R}}^{p}}\|\pi(Z)\|,\,\max_{B\in{\mathbb{R}}^{d^{\prime}\times
q}}\|\pi(B)\|_{F}\Big{)},$ $\displaystyle M_{0}$
$\displaystyle=\max\Big{(}\|A^{\textnormal{init}}\|_{F},\,\sup_{s\in[0,1]}\|Z^{\textnormal{init}}(s)\|,\,\|B^{\textnormal{init}}\|_{F}\Big{)}\quad\text{and}\quad
M=M_{0}+TM_{\pi}.$

Then the gradient flow is well defined on $[0,T]$, and, for $t\in[0,T]$,
$L\in\mathbb{N}^{*}$, and $k\in\\{1,\dots,L\\}$,

$\|A^{L}(t)\|_{F}\leqslant M,\quad\|Z_{k}^{L}(t)\|\leqslant
M,\quad\text{and}\quad\|B^{L}(t)\|_{F}\leqslant M.$ (17)

Moreover, there exist $\alpha,\beta>0$ such that, for $t\in[0,T]$ and
$k\in\\{1,\dots,L-1\\}$,

$\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|\leqslant\Big{(}\|Z_{k+1}^{L}(0)-Z_{k}^{L}(0)\|+\frac{\beta
T}{L}\Big{)}e^{\alpha T}.$

The following expressions for $\alpha$ and $\beta$ hold:

$\alpha=2e^{K}K^{\prime}M(e^{K}M^{2}M_{X}+M_{Y})\quad\text{and}\quad\beta=2Ke^{K}M(K+e^{K}K^{\prime}MM_{X})(e^{K}M^{2}M_{X}+M_{Y}),$

where

$\displaystyle M_{X}$ $\displaystyle=\sup_{x\in\mathcal{X}}\|x\|,\quad
M_{Y}=\sup_{y\in\mathcal{Y}}\|y\|,\quad K_{1}=\sup_{\|z\|\leqslant
M}\big{\|}\partial_{1}f(h,z)\big{\|}_{2}$ (18) $\displaystyle E$
$\displaystyle=\\{(h,z)\in\mathbb{R}^{d}\times\mathbb{R}^{p},\|h\|\leqslant
e^{K_{1}}MM_{X},\,\|z\|\leqslant M\\}$ (19) $\displaystyle K_{2}$
$\displaystyle=\sup_{(h,z)\in E}\big{\|}\partial_{2}f(h,z)\big{\|}_{2},\quad
K=\max(K_{1},K_{2})$ $\displaystyle K^{\prime}$ $\displaystyle=\sup_{(h,z)\in
E}\Big{(}\max\big{(}\big{\|}\partial_{22}f(h,z)\big{\|}_{2},\big{\|}\partial_{21}f(h,z)\big{\|}_{2}\big{)}\Big{)}.$

###### Proof.

The time-independent dynamics

$(A^{L},Z_{k}^{L},B^{L})\mapsto\Big{(}\pi\Big{(}-\frac{\partial\ell^{L}}{\partial
A^{L}}\Big{)},\pi\Big{(}-L\frac{\partial\ell^{L}}{\partial
Z_{k}^{L}}\Big{)},\pi\Big{(}-\frac{\partial\ell^{L}}{\partial
B^{L}}\Big{)}\Big{)}$

defining the gradient flow (5) are locally Lipschitz continuous, hence the
gradient flow is defined on a maximal interval $[0,T_{\max})$ by the Picard-
Lindelöf theorem (see Lemma 16). Let us show by contradiction that
$T_{\max}=T$. Assume that $T_{\max}<T$. If this is true, again by the Picard-
Lindelöf theorem, we know that the parameters diverge to infinity at
$T_{\max}$. However, for any $t\in[0,T_{\max})$, we have

$\displaystyle\|A^{L}(t)\|_{F}\leqslant\|A^{L}(0)\|_{F}+\int_{0}^{t}\Big{\|}\frac{dA^{L}}{dt}(\tau)\Big{\|}_{F}d\tau\leqslant
M_{0}+\int_{0}^{t}M_{\pi}d\tau\leqslant M_{0}+TM_{\pi}=M.$

Bounds on $B^{L}$ and $Z_{k}^{L}$ by $M$ can be shown similarly. This
contradicts the divergence of the parameters at $t=T_{\max}$. We conclude that
the gradient flow is well defined on $[0,T]$ and that the bounds (17) hold.

It remains to bound the difference $\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|$. We have,
for $t\in[0,T]$ and $k\in\\{1,\dots,L-1\\}$,

$\displaystyle\Big{\|}\frac{dZ_{k+1}^{L}}{dt}(t)$
$\displaystyle-\frac{dZ_{k}^{L}}{dt}(t)\Big{\|}=L\Big{\|}\frac{\partial\ell^{L}}{\partial
Z_{k+1}^{L}}(t)-\frac{\partial\ell^{L}}{\partial Z_{k}^{L}}(t)\Big{\|}$
$\displaystyle\leqslant\sum_{i=1}^{n}\frac{1}{n}\big{\|}\partial_{2}f(h_{k,i}^{L}(t),Z_{k+1}^{L}(t))^{\top}p_{k,i}^{L}(t)-\partial_{2}f(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))^{\top}p_{k-1,i}^{L}(t)\big{\|}$
$\displaystyle\leqslant\frac{1}{n}\sum_{i=1}^{n}\left\|\partial_{2}f(h_{k,i}^{L}(t),Z_{k+1}^{L}(t))\right\|_{2}\left\|p_{k,i}^{L}(t)-p_{k-1,i}^{L}(t)\right\|$
$\displaystyle\quad+\left\|p_{k-1,i}^{L}(t)\right\|\left\|\partial_{2}f(h_{k,i}^{L}(t),Z_{k+1}^{L}(t))-\partial_{2}f(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))\right\|_{2}$
(20)

Furthermore, for $t\in[0,T]$, $k\in\\{0,\dots,L-1\\}$, and
$i\in\\{1,\dots,n\\}$,

$\displaystyle\|h_{k+1,i}^{L}(t)\|=\|h_{k,i}^{L}(t)+\frac{1}{L}f(h_{k,i}^{L}(t),Z_{k+1}^{L}(t))\|\leqslant(1+\frac{K_{1}}{L})\|h_{k,i}^{L}(t)\|,$

since $f(\cdot,Z_{k+1}^{L}(t))$ is $K_{1}$-Lipschitz, where $K_{1}$ is defined
by (18), and $f(0,Z_{k+1}^{L}(t))=0$. Therefore, for any
$k\in\\{1,\dots,L\\}$,

$\|h_{k,i}^{L}(t)\|\leqslant
e^{K_{1}}\|h_{0,i}^{L}(t)\|=e^{K_{1}}\|A^{L}(t)x_{i}\|\leqslant
e^{K_{1}}MM_{X}.$ (21)

This bound shows that the pair $(h_{k,i}^{L}(t),Z_{k+1}^{L}(t))$ belongs to
the compact $E$ defined in (19) for every $t\in[0,T]$, $k\in\\{1,\dots,L\\}$,
and $i\in\\{1,\dots,n\\}$. In particular,
$\|\partial_{2}f(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))\|_{2}\leqslant K$, and

$\displaystyle\big{\|}\partial_{2}f(h_{k,i}^{L}(t),Z_{k+1}^{L}(t))$
$\displaystyle-\partial_{2}f(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))\big{\|}_{2}$
$\displaystyle\leqslant
K^{\prime}\|h_{k,i}^{L}(t)-h_{k-1,i}^{L}(t)\|+K^{\prime}\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|.$

Returning to (20), we obtain

$\displaystyle\Big{\|}\frac{dZ_{k+1}^{L}}{dt}(t)-\frac{dZ_{k}^{L}}{dt}(t)\Big{\|}$
$\displaystyle\leqslant\frac{1}{n}\sum_{i=1}^{n}K\|p_{k,i}^{L}(t)-p_{k-1,i}^{L}(t)\|$
$\displaystyle\quad+K^{\prime}\|p_{k-1,i}^{L}(t)\|\big{(}\|h_{k,i}^{L}(t)-h_{k-1,i}^{L}(t)\|+\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|\big{)}.$

For $k\in\\{1,\dots,L\\}$ and $i\in\\{1,\dots,n\\}$,

$\left\|p_{k,i}^{L}(t)-p_{k-1,i}^{L}(t)\right\|=\frac{1}{L}\big{\|}\partial_{1}f(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))p_{k,i}^{L}(t)\big{\|}\leqslant\frac{K}{L}\left\|p_{k,i}^{L}(t)\right\|,$

and, similarly,

$\left\|h_{k,i}^{L}(t)-h_{k-1,i}^{L}(t)\right\|=\frac{1}{L}\|f(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))\|\leqslant\frac{K}{L}\left\|h_{k-1,i}^{L}(t)\right\|\leqslant\frac{Ke^{K}MM_{X}}{L}.$

Thus,

$\displaystyle\Big{\|}\frac{dZ_{k+1}^{L}}{dt}(t)-\frac{dZ_{k}^{L}}{dt}(t)\Big{\|}\leqslant\frac{1}{n}\sum_{i=1}^{n}{\|p_{k,i}^{L}(t)\|}\Big{(}\frac{K^{2}}{L}+\frac{K^{\prime}K}{L}e^{K}MM_{X}+K^{\prime}\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|\Big{)}.$

Moreover, for $k\in\\{0,\dots,L\\}$ and $i\in\\{1,\dots,n\\}$,

$\displaystyle\|p_{k,i}^{L}(t)\|$
$\displaystyle\leqslant\|p_{k+1,i}^{L}(t)\|+\frac{1}{L}\big{\|}\partial_{1}f(h_{k,i}^{L}(t),Z_{k+1}^{L}(t))p_{k+1,i}^{L}(t)\big{\|}\leqslant\|p_{k+1,i}^{L}(t)\|+\frac{K}{L}\|p_{k+1,i}^{L}(t)\|.$

Hence

$\displaystyle\|p_{k,i}^{L}(t)\|$ $\displaystyle\leqslant
e^{K}\|p_{L,i}^{L}(t)\|=2e^{K}\|B^{L}(t)^{\top}(F^{L}(x_{i};t)-y_{i})\|$
$\displaystyle\leqslant
2e^{K}M\big{(}\|B^{L}(t)h_{L,i}^{L}(t)\|+\|y_{i}\|\big{)}\leqslant
2e^{K}M(e^{K}M^{2}M_{X}+M_{Y}),$

where we use (17) and (21) for the last inequality. Putting all the pieces
together, we obtain

$\Big{\|}\frac{dZ_{k}^{L}}{dt}(t)-\frac{dZ_{k+1}^{L}}{dt}(t)\Big{\|}\leqslant\alpha\|Z_{k}^{L}(t)-Z_{k+1}^{L}(t)\|+\frac{\beta}{L}.$

Integrating between $0$ and $t$, we see that

$\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|\leqslant\|Z_{k+1}^{L}(0)-Z_{k}^{L}(0)\|+\frac{\beta
t}{L}+\int_{0}^{t}\alpha\|Z_{k}^{L}(\tau)-Z_{k+1}^{L}(\tau)\|d\tau.$

Applying Grönwall’s inequality (see, e.g., Dragomir, 2003), we conclude that
$\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|\leqslant(\|Z_{k+1}^{L}(0)-Z_{k}^{L}(0)\|+\frac{\beta
T}{L})e^{\alpha T}$, as desired. ∎

###  A.2 The trained weights are bounded under the local PL condition

######  Proposition 8.

Consider the residual network (12) initialized as explained in Appendix A and
trained with the gradient flow (4) on $[0,\infty]$. Then, for $M>0$, there
exists $\mu>0$ such that, if the residual network satisfies the
$(M,\mu)$-local PL condition (1) around its initialization for any
$L\in\mathbb{N}^{*}$, then:

  1. (i)

The gradient flow is well defined on $\mathbb{R}_{+}$, and, for
$t\in\mathbb{R}_{+}$, $L\in\mathbb{N}^{*}$, and $k\in\\{1,\dots,L\\}$,

$\|A^{L}(t)\|_{F}\leqslant M_{A},\quad\|Z_{k}^{L}(t)\|\leqslant
M_{Z},\quad\text{and}\quad\|B^{L}(t)\|_{F}\leqslant M_{B},$

where

$M_{A}=\|A^{\textnormal{init}}\|_{2}+M,\quad
M_{Z}=\sup_{s\in[0,1]}\|Z^{\textnormal{init}}(s)\|+M,\quad\text{and}\quad
M_{B}=\|B^{\textnormal{init}}\|_{2}+M.$

  2. (ii)

There exists $\tilde{K}>0$ such that, for $t\in{\mathbb{R}}_{+}$,
$L\in\mathbb{N}^{*}$, and $k\in\\{1,\dots,L\\}$,

$\|Z_{k}^{L}(t)-Z_{k+1}^{L}(t)\|\leqslant\frac{\tilde{K}}{L}.$

  3. (iii)

There exists a bounded integrable function $b:{\mathbb{R}}_{+}\to{\mathbb{R}}$
such that, for $t\in{\mathbb{R}}_{+}$, $L\in\mathbb{N}^{*}$, and
$k\in\\{1,\dots,L\\}$,

$\max\Big{(}\Big{\|}\frac{dA^{L}}{dt}(t)\Big{\|},\Big{\|}\frac{dZ_{k}^{L}}{dt}(t)\Big{\|},\Big{\|}\frac{dB^{L}}{dt}(t)\Big{\|}\Big{)}\leqslant
b(t)$

  4. (iv)

$A^{L}(t)$, $B^{L}(t)$, and $Z_{k}^{L}(t)$ admit a limit uniformly over
$L\in\mathbb{N}^{*}$ and $k\in\\{1,\dots,L\\}$ as $t\to\infty$.

  5. (v)

For $t\in\mathbb{R}_{+}$ and $L\in\mathbb{N}^{*}$, $\ell^{L}(t)\leqslant
e^{-\mu t}\ell^{L}(0)$.

Moreover, the following expression for $\mu$ hold:

$\mu=\max(M_{B}K,M_{B}M_{X},M_{A}M_{X})\frac{8e^{K}}{M}\sup_{L\in\mathbb{N}^{*}}\sqrt{\ell^{L}(0)},$
(22)

where

$\displaystyle M_{X}$ $\displaystyle=\sup_{x\in\mathcal{X}}\|x\|,\quad
K_{1}=\sup_{\|z\|\leqslant M_{Z}}\big{\|}\partial_{1}f(h,z)\big{\|}$
$\displaystyle E$
$\displaystyle=\\{(h,z)\in\mathbb{R}^{d}\times\mathbb{R}^{p},\|h\|\leqslant
e^{K_{1}}M_{A}M_{X},\,\|z\|\leqslant M_{Z}\\}$ $\displaystyle K_{2}$
$\displaystyle=\sup_{(h,z)\in E}\big{\|}\partial_{2}f(h,z)\big{\|},\quad
K=\max(K_{1},K_{2}).$

###### Proof.

Let $M>0$, $\mu$ defined by (22), and assume that the residual network
satisfies the $(M,\mu)$-local PL condition (1) around its initialization for
any $L\in\mathbb{N}^{*}$.

The time-independent dynamics

$(A^{L},Z_{k}^{L},B^{L})\mapsto\Big{(}-\frac{\partial\ell^{L}}{\partial
A^{L}},-L\frac{\partial\ell^{L}}{\partial
Z_{k}^{L}},-\frac{\partial\ell^{L}}{\partial B^{L}}\Big{)}$

defining the gradient flow (5) are locally Lipschitz continuous, hence the
gradient flow is defined on a maximal interval $[0,T_{\max})$ by the Picard-
Lindelöf theorem (see Lemma 16). Let us show by contradiction that
$T_{\max}=\infty$. Assume that $T_{\max}<\infty$. If this is true, again by
the Picard-Lindelöf theorem, we know that the parameters diverge to infinity
at $T_{\max}$. In particular, there exist $t\in(0,T_{\max})$ and
$k\in\\{1,\dots,L\\}$ such that

$\|A^{L}(t)-A^{L}(0)\|_{F}>M\text{ or }\|Z_{k}^{L}(t)-Z_{k}^{L}(0)\|>M\text{
or }\|B^{L}(t)-B^{L}(0)\|_{F}>M.$

Let $t^{*}\in(0,T_{\max})$ be the infimum of such times $t$. Then, for
$t<t^{*}$ and $k\in\\{1,\dots,L\\}$,

$\|A^{L}(t)-A^{L}(0)\|_{F}\leqslant M\text{ and
}\|Z_{k}^{L}(t)-Z_{k}^{L}(0)\|\leqslant M\text{ and
}\|B^{L}(t)-B^{L}(0)\|_{F}\leqslant M,$ (23)

and, by continuity of $A^{L}$, $B^{L}$, and $Z_{k}^{L}$, these inequalities
also hold for $t=t^{*}$. By definition, this means that the $(M,\mu)$-local PL
condition is satisfied for $t\leqslant t^{*}$, and ensures that

$\Big{\|}\frac{\partial\ell^{L}}{\partial
A^{L}}(t)\Big{\|}_{F}^{2}+L\sum_{k=1}^{L}\Big{\|}\frac{\partial\ell^{L}}{\partial
Z_{k}^{L}}(t)\Big{\|}^{2}+\Big{\|}\frac{\partial\ell^{L}}{\partial
B^{L}}(t)\Big{\|}_{F}^{2}\geqslant\mu\ell^{L}(t).$

Therefore, by definition of the gradient flow (4),

$\displaystyle\frac{d\ell^{L}}{dt}(t)$
$\displaystyle=\Big{\langle}\frac{\partial\ell^{L}}{\partial
A^{L}}(t),\frac{dA^{L}}{dt}(t)\Big{\rangle}+\sum_{k=1}^{L}\Big{\langle}\frac{\partial\ell^{L}}{\partial
Z_{k}^{L}}(t),\frac{dZ_{k}^{L}}{dt}(t)\Big{\rangle}+\Big{\langle}\frac{\partial\ell^{L}}{\partial
B^{L}}(t),\frac{dB^{L}}{dt}(t)\Big{\rangle}$
$\displaystyle=-\Big{\|}\frac{\partial\ell^{L}}{\partial
A^{L}}(t)\Big{\|}_{F}^{2}-L\sum_{k=1}^{L}\Big{\|}\frac{\partial\ell^{L}}{\partial
Z_{k}^{L}}(t)\Big{\|}^{2}-\Big{\|}\frac{\partial\ell^{L}}{\partial
B^{L}}(t)\Big{\|}_{F}^{2}$ $\displaystyle\leqslant-\mu\ell^{L}(t).$

Thus, by Grönwall’s inequality, for $t\leqslant t^{*}$,

$\ell^{L}(t)\leqslant e^{-\mu t}\ell^{L}(0).$ (24)

Furthermore, by (23) and the definition of $M_{A}$, $M_{B}$, $M_{Z}$, we have,
for $t\leqslant t^{*}$ and $k\in\\{1,\dots,L\\}$,

$\|A^{L}(t)\|_{F}\leqslant M_{A},\quad\|Z_{k}^{L}(t)\|\leqslant
M_{Z},\quad\text{and}\quad\|B^{L}(t)\|_{F}\leqslant M_{B}.$

A quick scan through the proof of Proposition 7 reveals that by similar
arguments, we have, for $t\leqslant t^{*}$, $k\in\\{1,\dots,L\\}$, and
$i\in\\{1,\dots,n\\}$,

$(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))\in
E\quad\text{and}\quad\|p_{k-1,i}^{L}(t)\|\leqslant
2e^{K}\|p_{L,i}^{L}(t)\|\leqslant 2e^{K}M_{B}\|F^{L}(x_{i};t)-y_{i}\|.$

Thus, for $k\in\\{0,\dots,L\\}$,

$\frac{1}{n}\sum_{i=1}^{n}\|p_{k,i}^{L}(t)\|\leqslant\frac{2e^{K}M_{B}}{n}\sum_{i=1}^{n}\|F^{L}(x_{i};t)-y_{i}\|\leqslant
2e^{K}M_{B}\sqrt{\ell^{L}(t)}\leqslant 2e^{K}M_{B}e^{-\frac{\mu
t}{2}}\sqrt{\ell^{L}(0)},$ (25)

where the second inequality is a consequence of the Cauchy-Schwartz
inequality. Let us now bound $\|Z_{k}^{L}(t^{*})-Z_{k}^{L}(0)\|$. We have, for
$k\in\\{1,\dots,L\\}$,

$\displaystyle\|Z_{k}^{L}(t^{*})-Z_{k}^{L}(0)\|$
$\displaystyle\leqslant\int_{0}^{t^{*}}\Big{\|}\frac{dZ_{k}^{L}}{dt}(t)\Big{\|}dt$
$\displaystyle\leqslant\frac{1}{n}\sum_{i=1}^{n}\int_{0}^{t^{*}}\big{\|}\partial_{2}f(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))^{\top}p_{k-1,i}^{L}(t)\big{\|}dt$
(by (15)).
$\displaystyle\leqslant\frac{K}{n}\sum_{i=1}^{n}\int_{0}^{t^{*}}\|p_{k-1,i}^{L}(t)\|dt,$

since $(h_{k-1,i}^{L}(t),Z_{k}^{L}(t))\in E$ and
$\|\partial_{2}f(h,z)\|\leqslant K$ for $(h,z)\in E$. Therefore, by (25),

$\displaystyle\|Z_{k}^{L}(t^{*})-Z_{k}^{L}(0)\|\leqslant
2Ke^{K}M_{B}\int_{0}^{t^{*}}e^{-\frac{\mu
t}{2}}\sqrt{\ell^{L}(0)}dt\leqslant\frac{4Ke^{K}M_{B}}{\mu}\sqrt{\ell^{L}(0)}\leqslant\frac{M}{2},$

where the last inequality is a consequence of the definition of $\mu$.
Similarly, by (14) and (25),

$\displaystyle\|A^{L}(t^{*})-A^{L}(0)\|_{F}$
$\displaystyle\leqslant\int_{0}^{t^{*}}\Big{\|}\frac{dA^{L}}{dt}(t)\Big{\|}_{F}dt$
$\displaystyle\leqslant\int_{0}^{t^{*}}\frac{1}{n}\sum_{i=1}^{n}\big{\|}p_{0,i}^{L}(t)x_{i}^{\top}\big{\|}_{F}dt$
$\displaystyle\leqslant
2e^{K}M_{B}M_{X}\sqrt{\ell^{L}(0)}\int_{0}^{t^{*}}e^{-\frac{\mu t}{2}}dt$
$\displaystyle\leqslant\frac{4e^{K}M_{B}M_{X}}{\mu}\sqrt{\ell^{L}(0)}$
$\displaystyle\leqslant\frac{M}{2}.$

Finally, by (16),

$\displaystyle\|B^{L}(t^{*})-B(0)\|_{F}$
$\displaystyle\leqslant\int_{0}^{t^{*}}\Big{\|}\frac{dB^{L}}{dt}(t)\Big{\|}_{F}dt$
$\displaystyle\leqslant\int_{0}^{t^{*}}\frac{2}{n}\sum_{i=1}^{n}\|(F^{L}(x_{i};t)-y_{i})h_{L,i}^{L}(t)^{\top}\|_{F}dt$
$\displaystyle\leqslant
2e^{K}M_{A}M_{X}\sqrt{\ell^{L}(0)}\int_{0}^{t^{*}}e^{-\frac{\mu t}{2}}dt$
$\displaystyle\leqslant\frac{4e^{K}M_{A}M_{X}}{\mu}\sqrt{\ell^{L}(0)}$
$\displaystyle\leqslant\frac{M}{2},$

where the third inequality is a consequence of the Cauchy-Schwartz inequality
and of the fact that $\|h_{L,i}^{L}(t)\|\leqslant e^{K}M_{A}M_{X}$. By
continuity of $A^{L}$, $Z_{k}^{L}$, and $B^{L}$, these three bounds contradict
the definition of $t^{*}$. We conclude that $T_{\max}=\infty$ and that the
parameters stay within a ball of radius $M$ of their initialization, yielding
the inequalities, for $t\in\mathbb{R}_{+}$, $L\in\mathbb{N}^{*}$, and
$k\in\\{1,\dots,L\\}$,

$\|A^{L}(t)\|_{F}\leqslant M_{A},\quad\|B^{L}(t)\|_{F}\leqslant
M_{B},\quad\|Z_{k}^{L}(t)\|\leqslant M_{Z}.$

This proves statement $(i)$ of the proposition. Moreover, the analysis above
show that the derivatives of $A^{L}$, $Z_{k}^{L}$, and $B^{L}$ are bounded by
a bounded integrable function independent of $L$ and $k$. This shows $(iii)$,
together with the fact that the functions $A^{L}(t)$, $Z_{k}^{L}(t)$, and
$B^{L}(t)$ admit limits as $t\to\infty$. Furthermore, the convergence towards
their limit is uniform over $L$ and $k$, as we show for example for
$A^{L}(t)$. If we denote by $A_{\infty}^{L}$ its limit, and apply the same
steps as for bounding $\|A^{L}(t^{*})-A^{L}(0)\|_{F}$, we obtain, for any
$t\geqslant 0$,

$\displaystyle\|A_{\infty}^{L}-A^{L}(t)\|_{F}$
$\displaystyle\leqslant\int_{t}^{\infty}\Big{\|}\frac{dA^{L}}{d\tau}(\tau)\Big{\|}_{F}d\tau$
$\displaystyle\leqslant
2e^{K}M_{B}M_{X}\sqrt{\ell^{L}(0)}\int_{t}^{\infty}e^{\frac{-\mu\tau}{2}}d\tau$
$\displaystyle=\frac{4e^{K}M_{B}M_{X}}{\mu}e^{\frac{-\mu
t}{2}}\sqrt{\ell^{L}(0)}$ $\displaystyle\leqslant\frac{M}{2}e^{\frac{-\mu
t}{2}},$

where the last inequality comes from the definition of $\mu$. The bound is
independent of $L$, proving statement $(iv)$. Statement $(v)$ readily follows
from (24).

To complete the proof, it remains to prove statement $(ii)$ by bounding the
differences $\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|$. Now that we know that the
weights are bounded, we can follow the same steps as in the proof of
Proposition 7 and show the existence of $C_{1}$, $C_{2}>0$ such that

$\displaystyle\Big{\|}\frac{dZ_{k+1}^{L}}{dt}(t)$
$\displaystyle-\frac{dZ_{k}^{L}}{dt}(t)\Big{\|}\leqslant\frac{1}{n}\sum_{i=1}^{n}{\|p_{k,i}^{L}(t)\|}\Big{(}\frac{C_{1}}{L}+C_{2}\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|\Big{)}.$

Using (25), we obtain

$\displaystyle\Big{\|}\frac{dZ_{k+1}^{L}}{dt}(t)-\frac{dZ_{k}^{L}}{dt}(t)\Big{\|}\leqslant
2e^{K}M_{B}e^{-\frac{\mu
t}{2}}\sqrt{\ell^{L}(0)}\Big{(}\frac{C_{1}}{L}+C_{2}\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|\Big{)}.$

Integrating between $0$ and $t$, we obtain

$\displaystyle\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|$
$\displaystyle\leqslant\|Z_{k+1}^{L}(0)-Z_{k}^{L}(0)\|+\int_{0}^{t}2e^{K}M_{B}e^{-\frac{\mu\tau}{2}}\sqrt{\ell^{L}(0)}\frac{C_{1}}{L}d\tau$
$\displaystyle\quad+\int_{0}^{t}2e^{K}M_{B}e^{-\frac{\mu\tau}{2}}\sqrt{\ell^{L}(0)}C_{2}\|Z_{k+1}^{L}(\tau)-Z_{k}^{L}(\tau)\|d\tau$
$\displaystyle\leqslant\|Z_{k+1}^{L}(0)-Z_{k}^{L}(0)\|+\frac{C_{1}M}{2M_{X}L}$
$\displaystyle\quad+\int_{0}^{t}2e^{K}M_{B}e^{-\frac{\mu\tau}{2}}\sqrt{\ell^{L}(0)}C_{2}\|Z_{k+1}^{L}(\tau)-Z_{k}^{L}(\tau)\|d\tau,$

where the second inequality uses the definition of $\mu$. By Grönwall’s
inequality,

$\displaystyle\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|$
$\displaystyle\leqslant\Big{(}\|Z_{k+1}^{L}(0)-Z_{k}^{L}(0)\|+\frac{C_{1}M}{2M_{X}L}\Big{)}\exp\Big{(}\int_{0}^{t}2e^{K}M_{B}e^{-\frac{\mu\tau}{2}}\sqrt{\ell^{L}(0)}C_{2}d\tau\Big{)}$
$\displaystyle\leqslant\Big{(}\|Z_{k+1}^{L}(0)-Z_{k}^{L}(0)\|+\frac{C_{1}M}{2M_{X}L}\Big{)}\exp\Big{(}\frac{C_{2}M}{2M_{X}}\Big{)},$

again by definition of $\mu$. Finally, since
$Z_{k}^{L}(0)=Z^{\text{init}}(\frac{k}{L})$ and $Z^{\text{init}}$ is Lipschitz
continuous, this proves the existence of $\tilde{K}>0$ (independent of $L$,
$t$ and $k$) such that
$\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|\leqslant\frac{\tilde{K}}{L}$, which yields
statement $(ii)$. ∎

###  A.3 Generalized Arzelà–Ascoli theorem

######  Proposition 9 (Generalized Arzelà–Ascoli theorem).

Let $I\subseteq\mathbb{R}_{+}$ be an interval, and
$(Z_{k}^{L})_{L\in\mathbb{N}^{*},1\leqslant k\leqslant L}$ be a family of
$\mathcal{C}^{1}$ functions from $I$ to $\mathbb{R}^{p}$. Define

$\mathcal{Z}^{L}:[0,1]\times I\to\mathbb{R}^{p},\
(s,t)\mapsto\mathcal{Z}^{L}(s,t)=Z_{\lfloor(L-1)s\rfloor+1}^{L}(t).$

Assume that there exist a constant $C>0$ and a bounded integrable function
$b:I\to\mathbb{R}$ such that the following statements hold for any $t\in I$
and $L\in\mathbb{N}^{*}$:

  1. (i)

For $k\in\\{1,\dots,L-1\\}$,
$\|Z_{k+1}^{L}(t)-Z_{k}^{L}(t)\|\leqslant\frac{C}{L}$,

  2. (ii)

For $k\in\\{1,\dots,L\\}$, $\|Z_{k}^{L}(t)\|\leqslant C$ and
$\|\frac{dZ_{k}^{L}}{dt}(t)\|\leqslant b(t)$.

Then there exist a subsequence $(\mathcal{Z}^{\phi(L)})_{L\in\mathbb{N}^{*}}$
of $(\mathcal{Z}^{L})_{L\in\mathbb{N}^{*}}$ and a Lipschitz continuous
function $\mathcal{Z}^{\phi}:[0,1]\times I\to\mathbb{R}^{p}$ such that
$\mathcal{Z}^{\phi(L)}(s,t)$ tends to $\mathcal{Z}^{\phi}(s,t)$ uniformly over
$s$ and $t$.

Note that if $I$ is a compact interval, then the existence of a (uniformly)
convergent subsequence is guaranteed by the standard Arzelà–Ascoli theorem.
Indeed, the uniform equicontinuity is a consequence of assumptions $(i)$ and
$(ii)$, while $(ii)$ provides a uniform bound. However, if $I$ is not compact,
more involved arguments are needed.

###### Proof.

Assume, without loss of generality, that $b$ is also bounded by $C$. According
to assumption $(i)$, for $t\in I$ and $i,j\in\\{1,\dots,L\\}$,

$\displaystyle\|Z_{i}^{L}(t)-Z_{j}^{L}(t)\|\leqslant\frac{C|i-j|}{L}.$

Also, according to $(ii)$, for $t,t^{\prime}\in I$ and $k\in\\{1,\dots,L\\}$,

$\displaystyle\|Z_{k}^{L}(t)-Z_{k}^{L}(t^{\prime})\|=\Big{\|}\int_{t^{\prime}}^{t}\frac{dZ_{k}^{L}}{d\tau}(\tau)d\tau\Big{\|}\leqslant
C|t-t^{\prime}|.$

It follows that, for $s,s^{\prime}\in[0,1]$ and $t,t^{\prime}\in I$,

$\displaystyle\|\mathcal{Z}^{L}(s,t)-\mathcal{Z}^{L}(s^{\prime},t^{\prime})\|$
$\displaystyle\leqslant\|\mathcal{Z}^{L}(s,t)-\mathcal{Z}^{L}(s,t^{\prime})\|+\|\mathcal{Z}^{L}(s,t^{\prime})-\mathcal{Z}^{L}(s^{\prime},t^{\prime})\|$
$\displaystyle\leqslant
C|t-t^{\prime}|+\frac{C|\lfloor(L-1)s\rfloor-\lfloor(L-1)s^{\prime}\rfloor|}{L}.$

Therefore, with some simple algebra, we obtain

$\displaystyle\|\mathcal{Z}^{L}(s,t)-\mathcal{Z}^{L}(s^{\prime},t^{\prime})\|\leqslant
C|t-t^{\prime}|+C|s-s^{\prime}|+\frac{C}{L}.$ (26)

The statement of the proposition is then a consequence of the next three
steps.

##### There exists a convergent subsequence of
$(\mathcal{Z}^{L}(s,t))_{L\in\mathbb{N}^{*}}$.

First, let
$((s_{i},t_{i}))_{i\in\mathbb{N}}=(\mathbb{Q}\cap[0,1])\times(\mathbb{Q}\cap
I)$. By $(ii)$, the sequence
$(\mathcal{Z}^{L}(s_{i},t_{i}))_{L\in\mathbb{N}^{*},i\in\mathbb{N}}$ is
bounded. It is therefore possible to construct by a diagonal procedure a
subsequence $(\mathcal{Z}^{\phi(L)})_{L\in\mathbb{N}^{*}}$ such that, for each
$i\in\mathbb{N}$, $(\mathcal{Z}^{\phi(L)}(s_{i},t_{i}))_{L\in\mathbb{N}^{*}}$
is a convergent sequence.

Let us now show that $(\mathcal{Z}^{\phi(L)}(s,t))_{L\in\mathbb{N}^{*}}$
converges for any $s\in[0,1]$ and $t\in I$, by proving that it is a Cauchy
sequence in the complete metric space $\mathbb{R}^{p}$. Let $\varepsilon>0$,
$s\in[0,1]$, and $t\in I$. Since $((s_{i},t_{i}))_{i\in\mathbb{N}}$ is dense
in $[0,1]\times I$, there exists some $j\in\mathbb{N}$ such that
$|s_{j}-s|\leqslant\varepsilon$ and $|t_{j}-t|\leqslant\varepsilon$. Then, for
$L,M\in\mathbb{N}^{*}$, we have

$\displaystyle\|\mathcal{Z}^{\phi(L)}(s,t)$
$\displaystyle-\mathcal{Z}^{\phi(M)}(s,t)\|$
$\displaystyle\leqslant\|\mathcal{Z}^{\phi(L)}(s,t)-\mathcal{Z}^{\phi(L)}(s_{j},t_{j})\|+\|\mathcal{Z}^{\phi(L)}(s_{j},t_{j})-\mathcal{Z}^{\phi(M)}(s_{j},t_{j})\|$
$\displaystyle\quad+\|\mathcal{Z}^{\phi(M)}(s_{j},t_{j})-\mathcal{Z}^{\phi(M)}(s,t)\|$
$\displaystyle\leqslant
2C\varepsilon+\frac{C}{\phi(L)}+\|\mathcal{Z}^{\phi(L)}(s_{j},t_{j})-\mathcal{Z}^{\phi(M)}(s_{j},t_{j})\|+2C\varepsilon+\frac{C}{\phi(M)},$

where we used inequality (26) twice. Since
$(\mathcal{Z}^{\phi(L)}(s_{j},t_{j}))_{L\in\mathbb{N}^{*}}$ is a convergent
sequence, it is a Cauchy sequence. Thus, the bound can be made arbitrarily
small for $L,M$ large enough. This shows that
$(\mathcal{Z}^{\phi(L)}(s,t))_{L\in\mathbb{N}^{*}}$ is also a Cauchy sequence.
It is therefore convergent, and we denote by $\mathcal{Z}^{\phi}(s,t)$ its
limit.

##### The function $\mathcal{Z}^{\phi}$ is Lipschitz continuous.

By considering (26) for the subsequence $\phi(L)$ and letting $L\to\infty$, we
have that, for any $s,s^{\prime}\in[0,1]$ and $t,t^{\prime}\in I$,

$\|\mathcal{Z}^{\phi}(s,t)-\mathcal{Z}^{\phi}(s^{\prime},t^{\prime})\|\leqslant
C(|s-s^{\prime}|+|t-t^{\prime}|).$ (27)

##### The convergence of $(\mathcal{Z}^{\phi(L)}(s,t))_{L\in\mathbb{N}^{*}}$
to $\mathcal{Z}^{\phi}(s,t)$ is uniform over $s$ and $t$.

Let $\varepsilon>0$, $s\in[0,1]$, and $t\in I$. Then, by (26) and (27), it is
possible to find $\delta>0$ such that, for any
$s^{\prime},s^{\prime\prime}\in[0,1]$ and $t^{\prime},t^{\prime\prime}\in I$
satisfying $|s^{\prime}-s^{\prime\prime}|\leqslant\delta$ and
$|t^{\prime}-t^{\prime\prime}|\leqslant\delta$,

$\displaystyle\begin{split}\|\mathcal{Z}^{\phi(L)}(s^{\prime},t^{\prime})-\mathcal{Z}^{\phi(L)}(s^{\prime\prime},t^{\prime})\|&\leqslant\varepsilon+\frac{C}{\phi(L)}\quad\text{and}\quad\|\mathcal{Z}^{\phi}(s^{\prime},t^{\prime})-\mathcal{Z}^{\phi}(s^{\prime\prime},t^{\prime})\|\leqslant\varepsilon,\end{split}$
(28)

and

$\displaystyle\begin{split}\|\mathcal{Z}^{\phi(L)}(s^{\prime},t^{\prime})-\mathcal{Z}^{\phi(L)}(s^{\prime},t^{\prime\prime})\|\leqslant\varepsilon+\frac{C}{\phi(L)}\quad\text{and}\quad\|\mathcal{Z}^{\phi}(s^{\prime},t^{\prime})-\mathcal{Z}^{\phi}(s^{\prime},t^{\prime\prime})\|\leqslant\varepsilon.\end{split}$
(29)

Furthermore, there exists a finite set $\\{s_{1},\dots,s_{S}\\}\subset[0,1]$
such that

$\displaystyle[0,1]\subset\bigcup_{i=1}^{S}(s_{i}-\delta,s_{i}+\delta).$

In the sequel, we denote by $s^{*}$ an element of $\\{s_{1},\dots,s_{S}\\}$
that is at distance at most $\delta$ from $s$.

If $I$ is unbounded, then, by assumption $(ii)$ and since $b$ is integrable,
there exists some $t_{0}>0$ such that, for $t\geqslant t_{0}$,

$\displaystyle\|\mathcal{Z}^{\phi(L)}(s,t)-\mathcal{Z}^{\phi(L)}(s,t_{0})\|$
$\displaystyle\leqslant\int_{t_{0}}^{t}\Big{\|}\frac{d}{dt}Z_{\lfloor(\phi(L)s-1)\rfloor+1}^{\phi(L)}(\tau)\Big{\|}d\tau\leqslant\int_{t_{0}}^{t}b(\tau)d\tau\leqslant\varepsilon.$
(30)

The same inequality holds for $\mathcal{Z}^{\phi}$ by letting $L$ tend to
infinity. If $I$ is bounded, we simply let $t_{0}=\sup I$.

We may then pick a finite set $\\{t_{1},\dots,t_{T}\\}\subset[0,t_{0}]$ such
that

$\displaystyle[0,t_{0}]\subset\bigcup_{i=1}^{T}(t_{i}-\delta,t_{i}+\delta).$

Two cases may arise depending on the value of $t$. If $t\in[0,t_{0}]$, then
there exists an element of the set $\\{t_{1},\dots,t_{T}\\}$ at distance at
most $\delta$ from $t$, and we denote it by $t^{*}$. If $t>t_{0}$, we let
$t^{*}=t_{0}$. According to (29) and (30), we then have in both cases that

$\displaystyle\begin{split}\|\mathcal{Z}^{\phi(L)}(s,t)-\mathcal{Z}^{\phi(L)}(s,t^{*})\|\leqslant\varepsilon+\frac{C}{\phi(L)}\quad\text{and}\quad\|\mathcal{Z}^{\phi}(s,t)-\mathcal{Z}^{\phi}(s,t^{*})\|\leqslant\varepsilon.\end{split}$
(31)

To conclude, we have to bound the term
$\|\mathcal{Z}^{\phi(L)}(s,t)-\mathcal{Z}^{\phi}(s,t)\|$ uniformly over $s$
and $t$. We first have

$\displaystyle\|\mathcal{Z}^{\phi(L)}(s,t)$
$\displaystyle-\mathcal{Z}^{\phi}(s,t)\|$
$\displaystyle\leqslant\|\mathcal{Z}^{\phi(L)}(s,t)-\mathcal{Z}^{\phi(L)}(s,t^{*})\|+\|\mathcal{Z}^{\phi(L)}(s,t^{*})-\mathcal{Z}^{\phi}(s,t^{*})\|$
$\displaystyle\quad+\|\mathcal{Z}^{\phi}(s,t^{*})-\mathcal{Z}^{\phi}(s,t)\|$
$\displaystyle\leqslant
2\varepsilon+\frac{C}{\phi(L)}+\|\mathcal{Z}^{\phi(L)}(s,t^{*})-\mathcal{Z}^{\phi}(s,t^{*})\|,$

where the last inequality is a consequence of (31). The last term can be
bounded as follows:

$\displaystyle\|\mathcal{Z}^{\phi(L)}(s,t^{*})$
$\displaystyle-\mathcal{Z}^{\phi}(s,t^{*})\|$
$\displaystyle\leqslant\|\mathcal{Z}^{\phi(L)}(s,t^{*})-\mathcal{Z}^{\phi(L)}(s^{*},t^{*})\|+\|\mathcal{Z}^{\phi(L)}(s^{*},t^{*})-\mathcal{Z}^{\phi}(s^{*},t^{*})\|$
$\displaystyle\quad+\|\mathcal{Z}^{\phi}(s^{*},t^{*})-\mathcal{Z}^{\phi}(s,t^{*})\|$
$\displaystyle\leqslant
2\varepsilon+\frac{C}{\phi(L)}+\max_{i\in\\{1,\dots,S\\}}\|\mathcal{Z}^{\phi(L)}(s_{i},t^{*})-\mathcal{Z}^{\phi}(s_{i},t^{*})\|,$

by using (28) and the fact that $s^{*}\in\\{s_{1},\dots,s_{S}\\}$. Putting all
the pieces together, we finally obtain

$\displaystyle\|\mathcal{Z}^{\phi(L)}(s,t)$
$\displaystyle-\mathcal{Z}^{\phi}(s,t)\|\leqslant
4\varepsilon+\frac{2C}{\phi(L)}+\max_{i\in\\{1,\dots,S\\},j\in\\{1,\dots,T\\}}\|\mathcal{Z}^{\phi(L)}(s_{i},t_{j})-\mathcal{Z}^{\phi}(s_{i},t_{j})\|.$

By taking $L$ large enough, independent of $s$ and $t$, the sum of the last
two terms can be made less than $\varepsilon$. Since $\varepsilon$ is
arbitrary, this concludes the proof. ∎

A consequence of this result is a simplified version for sequences of
functions only indexed by $L$ and not $k$, as follows.

######  Corollary 10.

Let $I\subseteq\mathbb{R}_{+}$ be an interval, and
$(Z^{L})_{L\in\mathbb{N}^{*}}$ be a family of $\mathcal{C}^{1}$ functions from
$I$ to $\mathbb{R}^{p}$. Assume that there exist a constant $C>0$ and a
bounded integrable function $b:I\to\mathbb{R}$ such that, for any $t\in I$ and
$L\in\mathbb{N}^{*}$, $\|Z^{L}(t)\|\leqslant C$ and
$\|\frac{dZ^{L}}{dt}(t)\|\leqslant b(t)$. Then there exist a subsequence
$(Z^{\phi(L)})_{L\in\mathbb{N}^{*}}$ of $(Z^{L})_{L\in\mathbb{N}^{*}}$ and a
function $Z^{\phi}:I\to\mathbb{R}^{p}$ such that $Z^{\phi(L)}(t)$ tends to
$Z^{\phi}(t)$ uniformly over $t$.

###  A.4 Consistency of the Euler scheme for parameterized ODEs

######  Proposition 11 (Consistency of the Euler scheme for parameterized
ODEs.).

Let $(\theta_{k}^{L})_{L\in\mathbb{N}^{*},1\leqslant k\leqslant L}$ be a
bounded family of vectors of ${\mathbb{R}}^{p}$, and let

$\Theta^{L}:[0,1]\to\mathbb{R}^{p},\
s\mapsto\theta_{\left\lfloor(L-1)s\right\rfloor+1}^{L}.$

Assume that there exists $\Theta:[0,1]\to\mathbb{R}^{p}$ a Lipschitz
continuous function such that $\Theta^{L}(s)$ tends to $\Theta(s)$ uniformly
over $s$. Let $(a^{L})_{L\in\mathbb{N}^{*}}$ be a sequence of vectors in some
compact $E\subset{\mathbb{R}}^{d}$ converging to $a\in E$. Let
$g:\mathbb{R}^{d}\times\mathbb{R}^{p}\to\mathbb{R}^{d}$ be a $\mathcal{C}^{1}$
function such that $g(0,\cdot)\equiv 0$ and $g(\cdot,\theta)$ is uniformly
Lipschitz continuous for $\theta$ in any compact of ${\mathbb{R}}^{p}$.
Consider the discrete scheme

$\displaystyle\begin{split}u_{0}^{L}&=a^{L}\\\
u_{k+1}^{L}&=u_{k}^{L}+\frac{1}{L}g(u_{k}^{L},\theta_{k+1}^{L}),\quad
k\in\\{0,\dots,L-1\\}.\end{split}$ (32)

Then $u_{\left\lfloor Ls\right\rfloor}^{L}$ tends to $U(s)$ uniformly over
$s\in[0,1]$, where $U$ is the unique solution of the ODE

$\displaystyle\begin{split}U(0)&=a\\\
\frac{dU}{ds}(s)&=g(U(s),\Theta(s)),\quad s\in[0,1].\end{split}$ (33)

Moreover, the convergence only depends on the sequence
$(a^{L})_{L\in\mathbb{N}^{*}}$ and on its limit $a\in E$ through
$(\|a^{L}-a\|)_{L\in\mathbb{N}^{*}}$.

###### Proof.

Let $M$ be a bound of the sequence
$(\theta_{k}^{L})_{L\in\mathbb{N}^{*},1\leqslant k\leqslant L}$. By definition
of $\Theta^{L}$, the sequence $(\Theta^{L})_{L\in\mathbb{N}^{*}}$ is also
uniformly bounded by $M$, and the same is true for $\Theta$. Then the function
$g(\cdot,\Theta(s))$ is uniformly Lipschitz for $s\in[0,1]$. Furthermore,
$(U,s)\mapsto g(U,\Theta(s))$ is continuous in $s$ because $g$ and $\Theta$
are continuous. Thus the ODE (33) has a unique solution on $[0,1]$ by the
Picard-Lindelöf theorem (see Lemma 16).

Denote by $C$ the uniform Lipschitz constant of $g(\cdot,\theta)$ for
$\|\theta\|\leqslant M$. Since $g(0,\cdot)\equiv 0$ and $g(\cdot,\Theta(s))$
is $C$-Lipschitz, one has

$\Big{\|}\frac{dU}{ds}(s)\Big{\|}=\|g(U(s),\Theta(s))\|\leqslant C\|U(s)\|.$

Therefore, by Grönwall’s inequality,

$\|U(s)\|\leqslant\|U(0)\|\exp(C)=\|a\|\exp(C)\leqslant D_{E}\exp(C),$

where $D_{E}=\sup_{x\in E}\|x\|<\infty$. A similar reasoning applies to the
discrete scheme (32), using the discrete version of Grönwall’s inequality.
More precisely, for any $k\in\\{0,\dots,L-1\\}$,

$\|u_{k+1}^{L}\|\leqslant\|u_{k}^{L}\|+\frac{1}{L}\|g(u_{k}^{L},\theta_{k+1}^{L})\|\leqslant\Big{(}1+\frac{C}{L}\Big{)}\|u_{k}^{L}\|.$

Thus,

$\|u_{k}^{L}\|\leqslant\|u_{0}^{L}\|\exp(C)=\|a^{L}\|\exp(C)\leqslant
D_{E}\exp(C).$

Overall, we can consider a restriction of $g$ to a compact set depending only
on $M$, $C$, and $E$, which we will still denote by $g$ with a slight abuse of
notation.Since $g$ is $\mathcal{C}^{1}$, it is therefore bounded and Lipschitz
continuous, and we still let $C$ be its Lipschitz constant.

For $L\in\mathbb{N}^{*}$ and $k\in\\{0,\dots,L\\}$, we denote by
$\Delta_{k}^{L}$ the gap between the continuous and the discrete schemes,
i.e.,

$\Delta_{k}^{L}=\Big{\|}U\Big{(}\frac{k}{L}\Big{)}-u_{k}^{L}\Big{\|}.$

The next step is to recursively bound the size of this gap, first observing
that $\Delta_{0}^{L}=\|a^{L}-a\|$. We have that

$s\mapsto\frac{dU}{ds}(s)=g(U(s),\Theta(s))$ (34)

is a Lipschitz continuous function with some Lipschitz constant $\tilde{C}$.
To see this, just note that $U$ itself is Lipschitz continuous in $s$, since
$g$ is bounded, and therefore the function (34) is a composition of Lipschitz
continuous functions. In particular, $\frac{dU}{ds}$ is almost everywhere
differentiable, and its derivative $\frac{d^{2}U}{ds^{2}}(s)$ is bounded in
the supremum norm by $\tilde{C}$. As a consequence, for
$k\in\\{0,\dots,L-1\\}$, the Taylor expansion of $U$ on
$[\frac{k}{L},\frac{k+1}{L}]$ takes the form

$\displaystyle U\Big{(}\frac{k+1}{L}\Big{)}$
$\displaystyle=U\Big{(}\frac{k}{L}\Big{)}+\frac{1}{L}\frac{dU}{ds}\Big{(}\frac{k}{L}\Big{)}+\int_{k/L}^{(k+1)/L}\Big{(}\frac{k+1}{L}-s\Big{)}\frac{d^{2}U}{ds^{2}}(s)ds,$

where the norm of the remainder term is less than $\tilde{C}/L^{2}$.
Therefore,

$\displaystyle\Delta_{k+1}^{L}$
$\displaystyle=\Big{\|}U\Big{(}\frac{k+1}{L}\Big{)}-u_{k+1}^{L}\Big{\|}$
$\displaystyle=\Big{\|}U\Big{(}\frac{k}{L}\Big{)}+\frac{1}{L}g\Big{(}U\Big{(}\frac{k}{L}\Big{)},\Theta\Big{(}\frac{k}{L}\Big{)}\Big{)}+\int_{k/L}^{(k+1)/L}\Big{(}\frac{k+1}{L}-s\Big{)}\frac{d^{2}U}{ds^{2}}(s)ds$
$\displaystyle\quad-
u_{k}^{L}-\frac{1}{L}g(u_{k}^{L},\theta_{k+1}^{L})\Big{\|}$
$\displaystyle\leqslant\Big{\|}U\Big{(}\frac{k}{L}\Big{)}-u_{k}^{L}\Big{\|}+\Big{\|}\frac{1}{L}g\Big{(}U\Big{(}\frac{k}{L}\Big{)},\Theta\Big{(}\frac{k}{L}\Big{)}\Big{)}-\frac{1}{L}g(u_{k}^{L},\theta_{k+1}^{L})\Big{\|}$
$\displaystyle\quad+\int_{k/L}^{(k+1)/L}\Big{(}\frac{k+1}{L}-s\Big{)}\Big{\|}\frac{d^{2}U}{ds^{2}}(s)\Big{\|}ds$
$\displaystyle\leqslant\Delta_{k}^{L}+\frac{C}{L}\Delta_{k}^{L}+\frac{C}{L}\Big{\|}\Theta\Big{(}\frac{k}{L}\Big{)}-\theta_{k+1}^{L}\Big{\|}+\frac{\tilde{C}}{L^{2}}.$

In the last inequality, we used the fact that $g$ is $C$-Lipschitz. Since,by
definition, $\theta_{k+1}^{L}$ = $\Theta^{L}(\frac{k}{L-1})$, we obtain, for
$k\in\\{0,\dots,L-1\\}$,

$\displaystyle\Delta_{k+1}^{L}$
$\displaystyle\leqslant\Big{(}1+\frac{C}{L}\Big{)}\Delta_{k}^{L}+\frac{C}{L}\Big{\|}\Theta\Big{(}\frac{k}{L}\Big{)}-\Theta^{L}\Big{(}\frac{k}{L-1}\Big{)}\Big{\|}+\frac{\tilde{C}}{L^{2}}$
$\displaystyle\leqslant\Big{(}1+\frac{C}{L}\Big{)}\Delta_{k}^{L}+\frac{C}{L}\sup_{s\in[0,1]}\|\Theta(s)-\Theta^{L}(s)\|+\frac{C}{L}\Big{\|}\Theta\Big{(}\frac{k}{L}\Big{)}-\Theta\Big{(}\frac{k}{L-1}\Big{)}\Big{\|}+\frac{\tilde{C}}{L^{2}}$
$\displaystyle\leqslant\Big{(}1+\frac{C}{L}\Big{)}\Delta_{k}^{L}+\frac{C}{L}\sup_{s\in[0,1]}\|\Theta(s)-\Theta^{L}(s)\|+\frac{CC_{\Theta}}{L^{2}}+\frac{\tilde{C}}{L^{2}},$

where $C_{\Theta}$ is the Lipschitz constant of $\Theta$. By the discrete
Grönwall’s inequality, we deduce that, for $k\in\\{0,\dots,L-1\\}$,

$\displaystyle\Delta_{k+1}^{L}$
$\displaystyle\leqslant\Big{(}\Delta_{0}^{L}+\sup_{s\in[0,1]}\|\Theta(s)-\Theta^{L}(s)\|+\frac{C_{\Theta}}{L}+\frac{\tilde{C}}{LC}\Big{)}e^{C}$
$\displaystyle=\Big{(}\|a^{L}-a\|+\sup_{s\in[0,1]}\|\Theta(s)-\Theta^{L}(s)\|+\frac{C_{\Theta}}{L}+\frac{\tilde{C}}{LC}\Big{)}e^{C}.$
(35)

This shows that the gaps $\Delta_{k}^{L}$ converge to zero uniformly over
$k\in\\{0,\dots,L\\}$ as $L$ tends to infinity.

We conclude by observing that, for any $s\in[0,1]$,

$\|U(s)-u_{\left\lfloor
Ls\right\rfloor}^{L}\|\leqslant\Big{\|}U(s)-U\Big{(}\frac{\lfloor
Ls\rfloor}{L}\Big{)}\Big{\|}+\Big{\|}U\Big{(}\frac{\lfloor
Ls\rfloor}{L}\Big{)}-u_{\left\lfloor
Ls\right\rfloor}^{L}\Big{\|}\leqslant\frac{C_{U}}{L}+\Delta_{\lfloor
Ls\rfloor}^{L},$ (36)

where $C_{U}$ is the Lipschitz constant of $U$. Both terms converge to zero
uniformly over $s$ as $L$ tends to infinity. Finally, an inspection of our
bounds shows that the convergence only depends on
$(a^{L})_{L\in\mathbb{N}^{*}}\in E^{\mathbb{N}^{*}}$ through $\|a^{L}-a\|$. ∎

The results of Proposition 11 can be extended without much effort to two other
related cases. First, the parameters $\theta_{k}^{L}$ may depend on some other
variable $t$, as long as all assumptions are verified uniformly over $t$.
Second, these parameters may converge to some limit parameters as both $L$ and
$t$ go to infinity. This is encapsulated in the following two corollaries.

######  Corollary 12.

Let $I\subseteq\mathbb{R}_{+}$ be an interval. Let
$(\theta_{k}^{L})_{L\in\mathbb{N}^{*},1\leqslant k\leqslant L}$ be a uniformly
bounded family of functions from $I$ to $\mathbb{R}^{p}$, and let

$\Theta^{L}:[0,1]\times I\to\mathbb{R}^{p},\
(s,t)\mapsto\theta_{\left\lfloor(L-1)s\right\rfloor+1}^{L}(t).$

Assume that there exists a function $\Theta:[0,1]\times I\to\mathbb{R}^{p}$
such that $\Theta^{L}(s,t)$ tends to $\Theta(s,t)$ uniformly over $s$ and $t$,
and $\Theta(\cdot,t)$ is uniformly Lipschitz continuous for $t\in I$. Let
$(a^{L})_{L\in\mathbb{N}^{*}}$ be a family of functions from $I$ to some
compact $E\subset{\mathbb{R}}^{d}$, uniformly converging to $a:I\to E$. Let
$g:\mathbb{R}^{d}\times\mathbb{R}^{p}\to\mathbb{R}^{d}$ be a $\mathcal{C}^{1}$
function such that $g(0,\cdot)\equiv 0$ and $g(\cdot,\theta)$ is uniformly
Lipschitz continuous for $\theta$ in any compact of ${\mathbb{R}}^{p}$.
Consider the discrete scheme, for $t\in I$,

$\displaystyle u_{0}^{L}(t)$ $\displaystyle=a^{L}(t)$ $\displaystyle
u_{k+1}^{L}(t)$
$\displaystyle=u_{k}^{L}(t)+\frac{1}{L}g(u_{k}^{L}(t),\theta_{k+1}^{L}(t)),\quad
k\in\\{0,\dots,L-1\\}.$

Then $u_{\left\lfloor Ls\right\rfloor}^{L}(t)$ tends to $U(s,t)$ uniformly
over $s\in[0,1]$ and $t\in I$, where $U(\cdot,t)$ is the unique solution of
the ODE

$\displaystyle U(0,t)$ $\displaystyle=a(t)$ $\displaystyle\frac{\partial
U}{\partial s}(s,t)$ $\displaystyle=g(U(s,t),\Theta(s,t)),\quad s\in[0,1].$

Moreover, the convergence only depends on the sequence
$(a^{L})_{L\in\mathbb{N}^{*}}$ and on its limit $a\in E^{I}$ through
$(\sup_{t\in I}\|a^{L}(t)-a(t)\|)_{L\in\mathbb{N}^{*}}$.

######  Corollary 13.

Let $I\subseteq\mathbb{R}_{+}$ be an interval. Let
$(\theta_{k}^{L})_{L\in\mathbb{N}^{*},1\leqslant k\leqslant L}$ be a uniformly
bounded family of functions from $I$ to $\mathbb{R}^{p}$, and let

$\Theta^{L}:[0,1]\times\mathbb{R}_{+}\to\mathbb{R}^{p},\
(s,t)\mapsto\theta_{\left\lfloor(L-1)s\right\rfloor+1}^{L}(t).$

Assume that there exists a function $\Theta_{\infty}:[0,1]\to\mathbb{R}^{p}$
such that $\Theta^{L}(s,t)$ tends to $\Theta_{\infty}(s)$ uniformly over $s$
as $L,t\to\infty$, and $\Theta_{\infty}$ is Lipschitz continuous. Let
$(a^{L})_{L\in\mathbb{N}^{*}}$ be a family of functions from $I$ to some
compact $E\subset{\mathbb{R}}^{d}$, and converging to $a_{\infty}\in E$ as
$L,t\to\infty$. Let $g:\mathbb{R}^{d}\times\mathbb{R}^{p}\to\mathbb{R}^{d}$ be
a $\mathcal{C}^{1}$ function such that $g(0,\cdot)\equiv 0$ and
$g(\cdot,\theta)$ is uniformly Lipschitz continuous for $\theta$ in any
compact of ${\mathbb{R}}^{p}$. Consider the discrete scheme, for $t\in I$,

$\displaystyle u_{0}^{L}(t)$ $\displaystyle=a^{L}(t)$ $\displaystyle
u_{k+1}^{L}(t)$
$\displaystyle=u_{k}^{L}(t)+\frac{1}{L}g(u_{k}^{L}(t),\theta_{k+1}^{L}(t)),\quad
k\in\\{0,\dots,L-1\\}.$

Then $u_{\left\lfloor Ls\right\rfloor}^{L}(t)$ tends to $U(s)$ uniformly over
$s\in[0,1]$ as $L,t\to\infty$, where $U$ is the unique solution of the ODE

$\displaystyle U(0)$ $\displaystyle=a_{\infty}$
$\displaystyle\frac{dU}{ds}(s)$
$\displaystyle=g(U(s),\Theta_{\infty}(s)),\quad s\in[0,1].$

Moreover, the convergence only depends on the sequence
$(a^{L})_{L\in\mathbb{N}^{*}}$ and on its limit $a\in E^{I}$ through
$(\sup_{t\in I}\|a^{L}(t)-a(t)\|)_{L\in\mathbb{N}^{*}}$.

###  A.5 Large-depth convergence of the gradient flow

This section is devoted to proving the main result of Appendix A, namely the
large-depth convergence of the gradient flow. The setting we consider
encompasses both Section 4.1 (finite training time and clipped gradient flow)
and Section 4.2 (arbitrary training time and standard gradient flow). To this
end, we consider a training interval $I=[0,T]\subseteq{\mathbb{R}}_{+}$, for
$T\leqslant\infty$, and the gradient flow formulation (5), which is equivalent
to the standard gradient flow (4) if $\pi$ equals the identity. Note that we
do not need to assume in the following proof that $\pi$ is bounded (but only
Lipschitz continuous). Therefore, the proof also holds in the case where $\pi$
equals the identity.

######  Theorem 14.

Consider the residual network (12) initialized as explained in Appendix A and
trained with the gradient flow (5) on $I=[0,T]\subseteq{\mathbb{R}}_{+}$, for
some $T\in(0,\infty]$. Assume that there exists a unique solution to the
gradient flow, such that $(A^{L})_{L\in\mathbb{N}^{*}}$ and
$(B^{L})_{L\in\mathbb{N}^{*}}$ each satisfies the assumptions of Corollary 10,
and $(Z_{k}^{L})_{L\in\mathbb{N}^{*},1\leqslant k\leqslant L}$ satisfies the
assumptions of Proposition 9. Then the following four statements hold as $L$
tends to infinity:

  1. (i)

There exist functions $A:I\to{\mathbb{R}}^{q\times d}$ and
$B:I\to{\mathbb{R}}^{d^{\prime}\times q}$ such that $A^{L}(t)$ and $B^{L}(t)$
converge uniformly over $t\in I$ to $A(t)$ and $B(t)$.

  2. (ii)

There exists a Lipschitz continuous function $\mathcal{Z}:[0,1]\times
I\to\mathbb{R}^{p}$ such that

$\mathcal{Z}^{L}:[0,1]\times I\to\mathbb{R}^{p},\
(s,t)\mapsto\mathcal{Z}^{L}(s,t)=Z_{\left\lfloor(L-1)s\right\rfloor+1}^{L}(t)$

converges uniformly over $s\in[0,1]$ and $t\in I$ to $\mathcal{Z}(s,t)$.

  3. (iii)

Uniformly over $s\in[0,1]$, $t\in I$, and $x\in\mathcal{X}$, the hidden layer
$h_{\left\lfloor Ls\right\rfloor}^{L}(t)$ converges to the solution at time
$s$ of the neural ODE

$\displaystyle\begin{split}H(0,t)&=A(t)x\\\ \frac{\partial H}{\partial
s}(s,t)&=f(H(s,t),\mathcal{Z}(s,t)),\quad s\in[0,1].\end{split}$

  4. (iv)

Uniformly over $t\in I$ and $x\in\mathcal{X}$, the output $F^{L}(x;t)$
converges to $B(t)H(1,t)$.

###### Proof.

According to Proposition 9, there exists a subsequence
$(\mathcal{Z}^{\phi(L)})_{L\in\mathbb{N}^{*}}$ of
$(\mathcal{Z}^{L})_{L\in\mathbb{N}^{*}}$ and a Lipschitz continuous function
$\mathcal{Z}^{\phi}:[0,1]\times I\to\mathbb{R}^{p}$ such that
$\mathcal{Z}^{\phi(L)}(s,t)$ tends to $\mathcal{Z}^{\phi}(s,t)$ uniformly over
$s$ and $t$. Similarly, by Corollary 10, there exists subsequences of
$(A^{L})_{L\in\mathbb{N}^{*}}$ and $(B^{L})_{L\in\mathbb{N}^{*}}$ that
converge uniformly. With a slight abuse of notation, we still denote these
subsequences by $\phi$, and the corresponding limits by $A^{\phi}$ and
$B^{\phi}$.

In the remainder, we prove the uniqueness of the accumulation point
$(Z^{\phi},A^{\phi},B^{\phi})$ by showing that it is the solution of an ODE
that satisfies the assumptions of the Picard-Lindelöf theorem. The statements
$(i)$ to $(iv)$ then follow easily.

Consider a general input $(x,y)\in\mathcal{X}\times\mathcal{Y}$, and let
$H^{L}(s,t)=h_{\left\lfloor Ls\right\rfloor}^{L}(t)$ (recall that
$h_{k}^{L}(t)$ is defined by the forward propagation (12)). Corollary 12, with
$\theta_{k}^{L}=Z_{k}^{\phi(L)}$, $\Theta=\mathcal{Z}^{\phi}$,
$a^{L}=A^{\phi(L)}x$, $g=f$, ensures that $H^{\phi(L)}(s,t)$ converges
uniformly (over $s$ and $t$) to $H^{\phi}(s,t)$ that is the solution at time
$s$ of the ODE

$\displaystyle H^{\phi}(0,t)$ $\displaystyle=A^{\phi}(t)x$
$\displaystyle\frac{\partial H^{\phi}}{\partial s}(s,t)$
$\displaystyle=f(H^{\phi}(s,t),\mathcal{Z}^{\phi}(s,t)),\quad s\in[0,1].$

By inspecting the proof of the corollary, we also have that
$(h_{k}^{\phi(L)})_{L\in\mathbb{N}^{*},1\leqslant k\leqslant\phi(L)}$ and
$(H^{\phi(L)})_{L\in\mathbb{N}^{*}}$ are uniformly bounded and that
$H^{\phi}(\cdot,t)$ is uniformly Lipschitz continuous for $t\in I$.

We now turn our attention to the backpropagation recurrence (13), which
defines the backward state $p_{k}^{L}(t)$. First observe that the convergence
of $H^{\phi(L)}$ implies that

$p_{\phi(L)}^{\phi(L)}(t)=2B^{\phi(L)}(t)^{\top}(B^{\phi(L)}(t)h_{\phi(L)}^{\phi(L)}(t)-y)=2B^{\phi(L)}(t)^{\top}(B^{\phi(L)}(t)H^{\phi(L)}(1,t)-y)$

converges uniformly to
$2B^{\phi}(t)^{\top}(B^{\phi}(t)H^{\phi}(1,t)-y)\in{\mathbb{R}}^{d}$. Now, let
$P^{L}(s,t)=p_{\left\lfloor Ls\right\rfloor}^{L}(t)$. We apply again Corollary
12, this time to the backpropagation recurrence (13), with
$\theta_{k}^{L}=(h_{k}^{\phi(L)},Z_{k}^{\phi(L)})$,
$\Theta=(H^{\phi},\mathcal{Z}^{\phi})$,
$g:(p,(h,Z))\mapsto\partial_{1}f(h,Z)p$, and
$a^{L}=2(B^{\phi(L)})^{\top}(B^{\phi(L)}H^{\phi(L)}(1,\cdot)-y)$. Let us
quickly check that the conditions of the corollary are met:

  * •

The sequence $(h_{k}^{\phi(L)})_{L\in\mathbb{N}^{*},1\leqslant
k\leqslant\phi(L)}$ is bounded, as noted previously, and the same holds for
$(Z_{k}^{\phi(L)})_{L\in\mathbb{N}^{*},1\leqslant k\leqslant\phi(L)}$ by the
assumptions of Theorem 14.

  * •

The function $H^{\phi}(\cdot,t)$ is uniformly Lipschitz continuous for $t\in
I$, as noted previously, and the same is true for $Z^{\phi}(\cdot,t)$ since
$Z^{\phi}$ is Lipschitz continuous.

  * •

The function $h_{\left\lfloor(\phi(L)-1)s\right\rfloor+1}^{\phi(L)}(t)$ tends
to $H^{\phi}(s,t)$ uniformly over $s$ and $t$, as seen in the beginning of the
proof. More precisely, we know that
$H^{\phi(L)}(s,t)=h_{\left\lfloor\phi(L)s\right\rfloor}^{\phi(L)}(t)$ tends to
$H^{\phi}(s,t)$. Simple algebra and the fact that two successive iterates of
(12) are separated by a distance proportional to $1/L$ show that both
statements are equivalent. Furthermore, $\mathcal{Z}^{\phi(L)}(s,t)$ tends to
$\mathcal{Z}^{\phi}(s,t)$ uniformly over $s$ and $t$ as noted above.

  * •

The sequence $(a^{L})_{L\in\mathbb{N}^{*}}$ is uniformly bounded, since
$B^{\phi(L)}$ and $H^{\phi(L)}(1,\cdot)$ are. It also converges uniformly to
$a:t\mapsto 2B^{\phi}(t)^{\top}(B^{\phi}(t)H^{\phi}(1,t)-y)$.

  * •

The function $g$ is $\mathcal{C}^{1}$ since $f$ is $\mathcal{C}^{2}$. We
clearly have $g(0,\cdot)\equiv 0$. Finally, $g(\cdot,(h,Z))$ is uniformly
Lipschitz continuous for $(h,Z)$ in any compact since $\partial_{1}f$ is
continuous.

Overall, we obtain that $P^{\phi(L)}(s,t)$ converges uniformly (over $s$ and
$t$) to $P^{\phi}(s,t)$, the solution at time $s$ of the backward ODE

$\displaystyle P^{\phi}(1,t)$
$\displaystyle=2B^{\phi}(t)^{\top}(B^{\phi}(t)H^{\phi}(1,t)-y)$
$\displaystyle\frac{\partial P^{\phi}}{\partial s}(s,t)$
$\displaystyle=\partial_{1}f(H^{\phi}(s,t),\mathcal{Z}^{\phi}(s,t))P^{\phi}(s,t),\quad
s\in[0,1].$

Furthermore, the proof of the corollary shows that
$(P^{\phi(L)})_{L\in\mathbb{N}^{*}}$ is uniformly bounded. Now, recall that
the gradient flow for $Z_{k}^{\phi(L)}(t)$, given by (5) and (15), takes the
following form, for $t\in I$ and $k\in\\{1,\dots,\phi(L)\\}$,

$\frac{\partial Z_{k}^{\phi(L)}(t)}{\partial
t}=\pi\Big{(}-\frac{1}{n}\sum_{i=1}^{n}\partial_{2}f(h_{k-1,i}^{\phi(L)}(t),Z_{k}^{\phi(L)}(t))^{\top}p_{k-1,i}^{\phi(L)}(t)\Big{)},$

where the $i$ subscript corresponds to the $i$-th input $x_{i}$. By
definition, for $s\in[0,1]$,
$\mathcal{Z}^{\phi(L)}(s,t)=Z_{\left\lfloor(\phi(L)-1)s\right\rfloor+1}^{\phi(L)}(t)$.
Thus, the equation above can be rewritten, for $s\in[0,1]$ and $t\in I$,

$\frac{\partial\mathcal{Z}^{\phi(L)}(s,t)}{\partial
t}=\pi\Big{(}-\frac{1}{n}\sum_{i=1}^{n}\partial_{2}f(h_{\lfloor(\phi(L)-1)s\rfloor,i}^{\phi(L)}(t),Z_{\lfloor(\phi(L)-1)s\rfloor+1}^{\phi(L)}(t))^{\top}p_{\lfloor(\phi(L)-1)s\rfloor,i}^{\phi(L)}(t)\Big{)}.$
(37)

The term inside $\pi$ can be rewritten as

$-\frac{1}{n}\sum_{i=1}^{n}\partial_{2}f\Big{(}H_{i}^{\phi(L)}\Big{(}\frac{\lfloor(\phi(L)-1)s\rfloor}{\phi(L)},t\Big{)},\mathcal{Z}^{\phi(L)}(s,t)\Big{)}^{\top}P_{i}^{\phi(L)}\Big{(}\frac{\lfloor(\phi(L)-1)s\rfloor}{\phi(L)},t\Big{)}.$

Since $f$ is $\mathcal{C}^{2}$, $\partial_{2}f$ is locally Lipschitz
continuous. Applying the first part of the proof to the specific case of
$x_{i}$, we know that $H_{i}^{\phi(L)}$ and $P_{i}^{\phi(L)}$ uniformly
bounded, and that $H_{i}^{\phi(L)}(s,t)$ and $P_{i}^{\phi(L)}(s,t)$ converge
uniformly to $H_{i}^{\phi}(s,t)$ and $P_{i}^{\phi}(s,t)$. Therefore, the
right-hand side of (37) converges uniformly over $s$ and $t$ to

$\pi\Big{(}-\frac{1}{n}\sum_{i=1}^{n}\partial_{2}f(H_{i}^{\phi}(s,t),\mathcal{Z}^{\phi}(s,t))^{\top}P_{i}^{\phi}(s,t)\Big{)}.$

We have just shown the uniform convergence of the derivative in $t$ of
$\mathcal{Z}^{\phi(L)}(s,t)$. Furthermore, we know that, for $s\in[0,1]$, the
sequence $(t\mapsto\mathcal{Z}^{\phi(L)}(s,t))_{L\in\mathbb{N}^{*}}$ converges
to $\mathcal{Z}^{\phi}(s,\cdot)$. These two statements imply that
$\mathcal{Z}^{\phi}$ is differentiable with respect to $t$ and that, for
$s\in[0,1]$, its derivative satisfies the ordinary differential equation

$\frac{\partial\mathcal{Z}^{\phi}(s,t)}{\partial
t}=\pi\Big{(}-\frac{1}{n}\sum_{i=1}^{n}\partial_{2}f(H_{i}^{\phi}(s,t),\mathcal{Z}^{\phi}(s,t))^{\top}P_{i}^{\phi}(s,t)\Big{)}.$
(38)

Moreover, by our initialization scheme,

$\mathcal{Z}^{\phi}(s,0)=Z^{\textnormal{init}}(s).$ (39)

A similar approach reveals that $A^{\phi}(t)$ and $B^{\phi}(t)$ are
differentiable and that they verify the equations

$\displaystyle\frac{dA^{\phi}}{dt}(t)$
$\displaystyle=\pi\Big{(}-\frac{1}{n}\sum_{i=1}^{n}P_{i}^{\phi}(0,t)x_{i}^{\top}\Big{)},\quad$
$\displaystyle A^{\phi}(0)=A^{\textnormal{init}},$ (40)
$\displaystyle\frac{dB^{\phi}}{dt}(t)$
$\displaystyle=\pi\Big{(}-\frac{2}{n}\sum_{i=1}^{n}(B^{\phi}(t)H_{i}^{\phi}(1,t)-y_{i})H_{i}^{\phi}(1,t)^{\top}\Big{)},\quad$
$\displaystyle B^{\phi}(0)=B^{\textnormal{init}}.$ (41)

The equations (38) to (41) can be seen as an initial value problem whose
variables are the function
$\mathcal{Z}^{\phi}(\cdot,t):[0,1]\to{\mathbb{R}}^{p}$ and the matrices
$A^{\phi}(t)\in{\mathbb{R}}^{q\times
d},B^{\phi}(t)\in{\mathbb{R}}^{d^{\prime}\times q}$. To complete the proof, it
remains to show, using the Picard-Lindelöf theorem (see Lemma 16), that there
exists a unique solution to this problem. First, note that the space
$\mathcal{B}([0,1],{\mathbb{R}}^{p})$ of bounded functions from $[0,1]$ to
${\mathbb{R}}^{p}$ endowed with the supremum norm is a Banach space, which is
the proper space in which to apply the Picard-Lindelöf theorem. We therefore
endow the space of parameters
$\mathcal{B}([0,1],{\mathbb{R}}^{p})\times{\mathbb{R}}^{q\times
d}\times{\mathbb{R}}^{d^{\prime}\times q}$ with the norm

$\|(\mathcal{Z},A,B)\|:=\sup_{s\in[0,1]}\|\mathcal{Z}(s)\|+\|A\|_{2}+\|B\|_{2},$

which makes it a Banach space. We have to show that the mapping

$\displaystyle\begin{split}(\mathcal{Z},A,B)\mapsto\bigg{(}&s\mapsto\pi\Big{(}-\frac{1}{n}\sum_{i=1}^{n}\partial_{2}f(H_{i}(s),\mathcal{Z}(s))^{\top}P_{i}(s)\Big{)},\\\
&\pi\Big{(}-\frac{1}{n}\sum_{i=1}^{n}P_{i}(0)x_{i}^{\top}\Big{)},\;\pi\Big{(}-\frac{2}{n}\sum_{i=1}^{n}(BH_{i}(1)-y_{i})H_{i}(1)^{\top}\Big{)}\bigg{)}\end{split}$
(42)

is locally Lipschitz continuous with respect to this norm, where we recall
that $H_{i}(s)$ in (42) is the solution at time $s$ of the initial value
problem

$\displaystyle\begin{split}H_{i}(0)&=Ax_{i}\\\
\frac{dH_{i}}{ds}(s)&=f(H_{i}(s),\mathcal{Z}(s)),\quad s\in[0,1],\end{split}$
(43)

and $P_{i}(s)$ is the solution at time $s$ of the initial value problem

$\displaystyle\begin{split}P_{i}(1)&=2B^{\top}(BH_{i}(1)-y_{i})\\\
\frac{dP_{i}}{ds}(s)&=\partial_{1}f(H_{i}(s),\mathcal{Z}(s))P_{i}(s),\quad
s\in[0,1].\end{split}$ (44)

To prove that the mapping (42) is locally Lipschitz continuous, we first check
that it is well defined. Since $\mathcal{Z}$ is assumed to be only bounded
(and not continuous), the solutions of the initial value problems (43) and
(44) are well defined in the sense of the Caratheodory conditions, which are
given in Lemma 17.

Next, we can show that $(\mathcal{Z},A,B)\mapsto H_{i}$ is locally Lipschitz
continuous for $i\in\\{1,\dots,n\\}$. To do this, consider two sets of
parameters $(\mathcal{Z},A,B)$ and $(\tilde{\mathcal{Z}},\tilde{A},\tilde{B})$
belonging to a compact set $D$. Let $H_{i}$ and $\tilde{H}_{i}$ denote the
corresponding hidden states. As in the proof of Proposition 11, it holds that
$H_{i}$ and $\tilde{H}_{i}$ belong to some compact set $E$ that depends only
on $D$ and $f$. Let $K_{f}$ be the Lipschitz constant of the $\mathcal{C}^{1}$
function $f$ on $E\times D$. Then,

$\displaystyle\|\tilde{H}_{i}(s)-H_{i}(s)\|$
$\displaystyle\leqslant\|\tilde{H}_{i}(0)-H_{i}(0)\|+\int_{0}^{s}\Big{\|}\frac{d\tilde{H}_{i}}{dr}(r)-\frac{dH_{i}}{dr}(r)\Big{\|}dr$
$\displaystyle\leqslant\|\tilde{H}_{i}(0)-H_{i}(0)\|+\int_{0}^{s}\|f(\tilde{H}_{i}(r),\tilde{\mathcal{Z}}(r))-f(H_{i}(r),\mathcal{Z}(r))\|dr.$

The norm inside the integral can be bounded by

$\displaystyle\|f(\tilde{H}_{i}(r),\tilde{\mathcal{Z}}(r))$
$\displaystyle-f(\tilde{H}_{i}(r),\mathcal{Z}(r))\|+\|f(\tilde{H}_{i}(r),\mathcal{Z}(r))-f(H_{i}(r),\mathcal{Z}(r))\|$
$\displaystyle\leqslant
K_{f}\sup_{r\in[0,1]}\|\tilde{\mathcal{Z}}(r)-\mathcal{Z}(r)\|+K_{f}\|\tilde{H}_{i}(r)-H_{i}(r)\|.$

Therefore,

$\displaystyle\|\tilde{H}_{i}(s)-H_{i}(s)\|\leqslant\|\tilde{A}-A\|_{2}\|x_{i}\|+K_{f}\sup_{r\in[0,1]}\|\tilde{\mathcal{Z}}(r)-\mathcal{Z}(r)\|+\int_{0}^{s}K_{f}\|\tilde{H}_{i}(r)-H_{i}(r)\|dr.$

Using Grönwall’s inequality, we obtain, for any $s\in[0,1]$,

$\displaystyle\|\tilde{H}_{i}(s)-H_{i}(s)\|\leqslant\Big{(}\|\tilde{A}-A\|_{2}\|x_{i}\|+K_{f}\sup_{r\in[0,1]}\|\tilde{\mathcal{Z}}(r)-\mathcal{Z}(r)\|\Big{)}\exp(K_{f}).$

This shows that the function $(\mathcal{Z},A,B)\mapsto H_{i}$ is locally
Lipschitz continuous. One proves by similar arguments that the function
$(\mathcal{Z},A,B)\mapsto P_{i}$ is locally Lipschitz continuous. Thus,
overall, the mapping (42) is locally Lipschitz continuous as a composition of
locally Lipschitz continuous functions.

The Picard-Lindelöf theorem guarantees the uniqueness of the maximal solution
of the initial value problem (38)–(41) in the space
$\mathcal{B}([0,1],{\mathbb{R}}^{p})\times{\mathbb{R}}^{d\times
q}\times{\mathbb{R}}^{d^{\prime}\times q}$. Since any accumulation point
$(\mathcal{Z}^{\phi},A^{\phi},B^{\phi})$ is a solution belonging to this
space, this proves the uniqueness of the accumulation point, which we
therefore denote as $(\mathcal{Z},A,B)$.

The uniform convergence of $(\mathcal{Z}^{L},A^{L},B^{L})$ to
$(\mathcal{Z},A,B)$ is then easily shown by contradiction. Suppose that
uniform convergence does not hold. If this is true, then there exists a
subsequence that stays at distance $\varepsilon>0$ from $(\mathcal{Z},A,B)$
(in the sense of the uniform norm). Then arguments similar to the beginning of
the proof show the existence of a second accumulation point, which is a
contradiction. This shows the uniform convergence, yielding statements $(i)$
and $(ii)$ of the theorem.

Finally, reapplying Corollary 12 with $\theta_{k}^{L}=Z_{k}^{L}$,
$\Theta=\mathcal{Z}$, $a^{L}=A^{L}x$, $g=f$, completes the proof by proving
statements $(iii)$ and $(iv)$. ∎

##### Training dynamics of the limiting weights.

Interestingly, the proof of Theorem 14 provides us with an explicit
description of the evolution of the continuous-depth limiting weights during
training. With the notation of the proof, the continuous weights satisfy the
training dynamics:

$\displaystyle\frac{dA}{dt}(t)$
$\displaystyle=\pi\Big{(}-\frac{1}{n}\sum_{i=1}^{n}P_{i}(0,t)x_{i}^{\top}\Big{)}$
$\displaystyle\frac{\partial\mathcal{Z}}{\partial t}(s,t)$
$\displaystyle=\pi\Big{(}-\frac{1}{n}\sum_{i=1}^{n}\partial_{2}f(H_{i}(s,t),\mathcal{Z}(s,t))^{\top}P_{i}(s,t)\Big{)}$
$\displaystyle\frac{dB}{dt}(t)$
$\displaystyle=\pi\Big{(}-\frac{2}{n}\sum_{i=1}^{n}(B(t)H_{i}(1,t)-y_{i})H_{i}(1,t)^{\top}\Big{)},$

where we recall that $H_{i}(s,t)$ is the solution at time $s$ of the initial
value problem

$\displaystyle H_{i}(0,t)$ $\displaystyle=A(t)x_{i}$
$\displaystyle\frac{\partial H_{i}}{\partial s}(s,t)$
$\displaystyle=f(H_{i}(s,t),\mathcal{Z}(s,t)),\quad s\in[0,1],$

and $P_{i}(s,t)$ is the solution at time $s$ of the problem

$\displaystyle P_{i}(1,t)$ $\displaystyle=2B(t)^{\top}(B(t)H_{i}(1,t)-y_{i})$
$\displaystyle\frac{\partial P_{i}}{\partial s}(s,t)$
$\displaystyle=\partial_{1}f(H_{i}(s,t),\mathcal{Z}(s,t))P_{i}(s,t),\quad
s\in[0,1].$

These equations can be thought of as the continuous-depth equivalent of the
backpropagation equations.

###  A.6 Existence of the double limit when $L,t$ tend to infinity

######  Proposition 15.

Consider the residual network (12), and assume that:

  1. (i)

$A^{L}(t)$, $Z_{\left\lfloor Ls\right\rfloor}^{L}(t)$, and $B^{L}(t)$ converge
uniformly over $L\in\mathbb{N}^{*}$ and $s\in[0,1]$ as $t\to\infty$.

  2. (ii)

$A^{L}(t)$, $Z_{\left\lfloor Ls\right\rfloor}^{L}(t)$, and $B^{L}(t)$ converge
uniformly over $t\in\mathbb{R}_{+}$ and $s\in[0,1]$ as $L\to\infty$.

  3. (iii)

The loss $\ell^{L}(t)$ converges to $0$ uniformly over $L\in\mathbb{N}^{*}$ as
$t\to\infty$.

Then the following four statements hold as $t$ and $L$ tend to infinity:

  1. (i)

There exist matrices $A_{\infty}\in{\mathbb{R}}^{q\times d}$ and
$B_{\infty}\in{\mathbb{R}}^{d^{\prime}\times q}$ such that $A^{L}(t)$ and
$B^{L}(t)$ converge to $A_{\infty}$ and $B_{\infty}$.

  2. (ii)

There exists a Lipschitz continuous function
$\mathcal{Z}_{\infty}:[0,1]\to{\mathbb{R}}^{p}$ such that $Z_{\left\lfloor
Ls\right\rfloor}^{L}(t)$ converges to $\mathcal{Z}_{\infty}(t)$ uniformly over
$s\in[0,1]$.

  3. (iii)

Uniformly over $s\in[0,1]$ and $x\in\mathcal{X}$, the hidden layer
$h_{\left\lfloor Ls\right\rfloor}^{L}(t)$ converges to the solution at time
$s$ of the ODE

$\displaystyle\begin{split}H(0)&=A_{\infty}x\\\
\frac{dH}{ds}(s)&=f(H(s),\mathcal{Z}_{\infty}(s)),\quad s\in[0,1].\end{split}$
(45)

  4. (iv)

Uniformly over $x\in\mathcal{X}$, the output $F^{L}(x;t)$ converges to
$F_{\infty}(x)=B_{\infty}H(1)$. Furthermore, $F_{\infty}(x_{i})=y_{i}$ for
$i\in\\{1,\dots,n\\}$.

###### Proof.

The existence of limits $A_{\infty}$ and $B_{\infty}$ to $A^{L}(t)$ and
$B^{L}(t)$ as $L$ and $t$ tend to infinity is given by Lemma 19. The same
argument applies to $Z_{\left\lfloor sL\right\rfloor}^{L}(t)$, which provides
a limit $\mathcal{Z}_{\infty}(s)$ to the sequence. Furthermore, following the
proof of the lemma, we see that the convergence of $Z_{\left\lfloor
sL\right\rfloor}^{L}(t)$ to $\mathcal{Z}_{\infty}(s)$ is uniform over
$s\in[0,1]$. Corollary 13, applied with $\theta_{k}^{L}=Z_{k}^{L}$,
$\Theta_{\infty}=\mathcal{Z}_{\infty}$, $a^{L}=A^{L}x$, $g=f$, then ensures
that $h_{\left\lfloor Ls\right\rfloor}^{L}(t)$ converges uniformly (over
$s\in[0,1]$ and $x\in\mathcal{X}$) to $H(s)$ that is the solution at time $s$
of (45), as $L$ and $t$ tend to infinity. As a consequence, $F^{L}(x;t)$
converges uniformly over $x$ to $F_{\infty}(x)$ as $L,t\to\infty$.
Furthermore, recall that

$\ell^{L}(t)=\frac{1}{n}\sum_{i=1}^{n}\|F^{L}(x_{i};t)-y_{i}\|_{2}^{2}.$

The left-hand side converges as $L,t\to\infty$ to $0$ by assumption of the
proposition, while the right-hand side converges to

$\frac{1}{n}\sum_{i=1}^{n}\|F_{\infty}(x_{i})-y_{i}\|_{2}^{2}.$

Therefore, $F_{\infty}(x_{i})=y_{i}$ for $i\in\\{1,\dots,n\\}$, and the proof
is complete. ∎

##  Appendix B Proofs of the results of the main paper

In this section, we prove the results of the main paper. Most of these results
follow from those presented in Section A. The only substantial proof is that
of Proposition 5, which shows the local PL condition. It uses a result of
Nguyen & Mondelli (2020) involving the Hermite transform and the sub-Gaussian
variance proxy, which we define briefly. We refer to Debnath & Bhatta (2014,
Chapter 17) and Vershynin (2018, Sections 2.5.2 and 3.4.1), respectively, for
more detailed explanations.

##### Hermite transform.

The $r$-th normalized probabilist’s Hermite polynomial is given by

$h_{r}(x)=\frac{1}{\sqrt{r!}}(-1)^{r}e^{x^{2}/2}\frac{d^{r}}{dx^{r}}e^{-x^{2}/2},\quad
r\geqslant 0.$

This family of polynomials forms an orthonormal basis of square-integrable
functions for the inner product

$\langle
f_{1},f_{2}\rangle=\frac{1}{\sqrt{2\pi}}\int_{-\infty}^{\infty}f_{1}(x)f_{2}(x)e^{-x^{2}/2}dx.$

Therefore, any function $\sigma$ such that
$\frac{1}{\sqrt{2\pi}}\int_{-\infty}^{\infty}\sigma^{2}(x)e^{-x^{2}/2}dx<\infty$
can be decomposed on this basis. The $r$-th coefficient of this decomposition
is denoted by $\eta_{r}(\sigma)$.

##### Sub-Gaussian random vector.

A random vector $x\in{\mathbb{R}}^{d}$ is sub-Gaussian with variance proxy
$v_{x}>0$ if, for every $y\in{\mathbb{R}}^{d}$ of unit norm,

$\mathbb{P}(|\langle x,y\rangle|\geqslant t)\leqslant
2\exp\Big{(}-\frac{t^{2}}{2v_{x}^{2}}\Big{)}.$

##### Additional notation.

For a matrix $A$, we let $s_{\min}$ and $s_{\max}$ its minimum and maximum
singular values, and similarly, $\lambda_{\min}$ and $\lambda_{\max}$ its
minimum and maximum eigenvalues (whenever they exist).

Before delving into the proofs, we briefly describe the parts of this section
that make use of the specific model (3). The most important one is the proof
of Proposition 5, i.e., the proof that the residual network satisfies the
$(M,\mu)$-local PL condition. Additionally, in the proof of Proposition 3, the
expressions for $M$ and $K$ are valid only for the specific model (3).
Finally, in the proof of Theorem 6, the beginning of the proof reveals that
condition (22) of Proposition 8 on $\mu$ can be expressed as a condition on
the norm of the labels $y_{i}$. This applies only to the specific model (3).
Observe that, if one assumes that the general residual network of Section A
satisfies the $(M,\mu)$-local PL condition with $\mu$ given by (22), then the
rest of the proof of Theorem 6 unfolds, and the conclusions of the theorem
hold for the general model.

###  B.1 Proof of Proposition 1

Proposition 1 is a consequence of Proposition 7 with
$f(h,(V,W))=\frac{1}{\sqrt{m}}V\sigma(\frac{1}{\sqrt{q}}Wh)$.

###  B.2 Proof of Proposition 2

Proposition 11, with $\theta_{k}^{L}=(V_{k}^{L},W_{k}^{L})$,
$\Theta=(\mathcal{V},\mathcal{W})$, $a^{L}=Ax$,
$g(h,(V,W))=\frac{1}{\sqrt{m}}V\sigma(\frac{1}{\sqrt{q}}Wh)$, gives the
existence and uniqueness of the solution of the neural ODE (6). Moreover,
inspecting the proof of Proposition 11, equations (35) gives that, for any
input $x\in\mathcal{X}$, the difference between the last hidden layer
$h_{L}^{L}$ of the discrete residual network (3) and its continuous
counterpart $H(1)$ in the neural ODE (6) is bounded by

$C^{\prime}\Big{(}\frac{1}{L}+\sup_{s\in[0,1]}\|\Theta(s)-\Theta^{L}(s)\|\Big{)},$

where $C^{\prime}>0$ is independent of $L$ and $x\in\mathcal{X}$, and
$\Theta^{L}(s)=\theta_{\left\lfloor(L-1)s\right\rfloor+1}^{L}$. The function
$\Theta^{L}$ is a piecewise-constant interpolation of $\Theta$ with pieces of
length $\frac{1}{L-1}$. Since $\Theta$ is Lipschitz continuous, the distance
between $\Theta$ and $\Theta^{L}$ decreases as
$\nicefrac{{C^{\prime\prime}}}{{L}}$ for some $C^{\prime\prime}>0$ depending
on $\Theta$ but not on $L$. This yields
$\|h_{L}^{L}-H(1)\|\leqslant\frac{C^{\prime}(1++C^{\prime\prime})}{L}$, where
$C^{\prime}$ and $C^{\prime\prime}$ are independent of $L$ and
$x\in\mathcal{X}$. Since $F^{L}(x)=Bh_{L}^{L}$ and $F(x)=BH(1)$, the result is
proven.

###  B.3 Proof of Proposition 3

We apply Proposition 7 with
$f(h,(V,W))=\frac{1}{\sqrt{m}}V\sigma(\frac{1}{\sqrt{q}}Wh)$. Recall that the
parameters $Z=(V,W)$ are considered in Proposition 7 as a vector. In
particular, $\|Z\|=\|V\|_{F}+\|W\|_{F}$. Therefore, Proposition 7 shows that,
for $t\in[0,T]$, $L\in\mathbb{N}^{*}$, and $k\in\\{1,\dots,L\\}$,

$\|A^{L}(t)\|_{F}\leqslant
M,\quad\|V_{k}^{L}(t)\|_{F}+\|W_{k}^{L}(t)\|_{F}\leqslant
M,\quad\text{and}\quad\|B^{L}(t)\|_{F}\leqslant M,$

where

$\displaystyle M$ $\displaystyle=M_{0}+TM_{\pi}$ $\displaystyle M_{0}$
$\displaystyle=\max\Big{(}\|A^{L}(0)\|_{F},\,\|V_{0}^{L}(0)\|_{F}+\|W_{0}^{L}(0)\|_{F},\,\|B^{L}(0)\|_{F}\Big{)}$
$\displaystyle M_{\pi}$
$\displaystyle=\max\Big{(}\max_{A\in{\mathbb{R}}^{q\times
d}}\|\pi(A)\|_{F},\,\max_{Z\in{\mathbb{R}}^{q\times
m}\times{\mathbb{R}}^{m\times
q}}\|\pi(Z)\|,\,\max_{B\in{\mathbb{R}}^{d^{\prime}\times
q}}\|\pi(B)\|_{F}\Big{)}.$

Furthermore, due to our initialization scheme described in Section 3,

$\|A^{L}(0)\|_{F}=\sqrt{d},\quad\|V_{0}^{L}(0)\|_{F}=0,\quad\|W_{0}^{L}(0)\|_{F}\leqslant
2\sqrt{qm},\quad\|B^{L}(0)\|_{F}=\sqrt{d^{\prime}},$

where the third inequality holds with probability at least
$1-\exp(-\frac{3qm}{16})$ by Lemma 20. Since we take
$q\geqslant\max(d,d^{\prime})$, this implies that, with high probability,
$M_{0}\leqslant 2\sqrt{qm}$, yielding the formula for $M$ in Proposition 3.
Finally, the existence of $K=\beta Te^{\alpha T}$ such that the difference
between two successive weight matrices is bounded by $\nicefrac{{K}}{{L}}$, as
well as the dependence of $\alpha$ and $\beta$ on $\mathcal{X}$,
$\mathcal{Y}$, $M$, and $\sigma$, follows easily from Proposition 7, given
that our initialization scheme ensures that $Z_{k}^{L}(0)=Z_{k+1}^{L}(0)$ for
all $L\in\mathbb{N}^{*}$ and $k\in\\{1,\dots,L\\}$.

###  B.4 Proof of Theorem 4

By Proposition 3 and the fact that $\pi$ is bounded, the sequences
$(A^{L})_{L\in\mathbb{N}^{*}}$ and $(B^{L})_{L\in\mathbb{N}^{*}}$ each satisfy
the assumptions of Corollary 10, and
$(Z_{k}^{L})_{L\in\mathbb{N}^{*},1\leqslant k\leqslant L}$ satisfies the
assumptions of Proposition 9. Theorem 4 then follows directly from Theorem 14,
by taking, as previously,
$f(h,(V,W))=\frac{1}{\sqrt{m}}V\sigma(\frac{1}{\sqrt{q}}Wh)$.

###  B.5 Proof of Proposition 5

We drop the $L$ superscripts for this proof, since $L$ is fixed. Denote by
$\bar{A},\bar{B},\bar{V}_{k},\bar{W}_{k}$ parameters sampled according to the
initialization scheme of Section 3, which means in particular that
$\bar{V}_{k}=0$ and $\bar{W}_{k}=\bar{W}\sim\mathcal{N}^{\otimes(m\times q)}$.
Since, by assumption, the activation function $\sigma$ is bounded and not
constant, it cannot be a polynomial function. As a consequence, there are