In standard residual networks, representations evolve across a sequence of discrete, static layers $h_{t+1} = h_t + f(h_t, theta_t)$. In 2018, Ricky Chen, Yulia Rubanova, Jesse Bettencourt, and David Duvenaud realized that this residual difference equation is simply the forward Euler discretization of an underlying continuous process. By taking the infinitesimal step limit $Delta t to 0$, they unveiled Neural Ordinary Differential Equations (Neural ODEs)—replacing discrete layer stacks with continuous vector fields and revolutionizing memory efficiency via the continuous Adjoint Sensitivity Method.
1. From Discrete ResNets to Continuous Dynamical Flows
Recall the discrete residual update step:
$$h_{t+1} – h_t = f(h_t, theta_t)$$
Scaling step size by infinitesimal $Delta t$ and taking the continuum limit $Delta t to 0$ yields the continuous ordinary differential equation:
$$lim_{Delta t to 0} frac{h(t + Delta t) – h(t)}{Delta t} = frac{dh(t)}{dt} = f(h(t), t, theta)$$
where $f(h(t), t, theta)$ is a neural network defining a smooth vector field in state space. Computing the representation at terminal depth $T$ from initial input $h(0) = x$ is formulated as an Initial Value Problem (IVP):
$$h(T) = h(0) + int_0^T f(h(t), t, theta) dt = text{ODESolve}(h(0), f, 0, T, theta)$$
Crucially, the architecture possesses no fixed “depth”: depth is an adaptive continuum dynamically determined by modern numerical ODE solvers (such as Runge-Kutta 4/5 or Dormand-Prince with adaptive step sizes) based on user-specified error tolerances.
2. Constant Memory Backpropagation: The Adjoint State Method
Standard backpropagation through adaptive numerical ODE solvers requires storing all intermediate internal solver evaluations in memory, incurring catastrophic memory complexity scaling with the number of solver steps.
Chen et al. resurrected the Adjoint State Method—rooted in Lev Pontryagin’s 1962 Maximum Principle of optimal control theory. Define the continuous adjoint state $a(t)$ as the sensitivity of scalar loss $L$ with respect to state trajectory $h(t)$:
$$a(t) equiv frac{partial L}{partial h(t)}$$
Differentiating $a(t)$ with respect to continuous time yields the Adjoint Differential Equation:
$$frac{da(t)}{dt} = -a(t)^T frac{partial f(h(t), t, theta)}{partial h(t)}$$
To train the model, the algorithm performs a reverse-time ODE solve starting from terminal boundary condition $a(T) = frac{partial L}{partial h(T)}$. As the solver runs backward in time from $t = T$ to $t = 0$, it reconstructs the forward trajectory $h(t)$ in reverse while simultaneously accumulating the parameter gradient:
$$frac{dL}{dtheta} = -int_T^0 a(t)^T frac{partial f(h(t), t, theta)}{partial theta} dt$$
Because the forward states are recomputed in reverse time directly from $h(T)$, zero intermediate forward activations need to be stored in GPU VRAM! The training memory complexity of an arbitrarily complex Neural ODE is strictly $O(1)$ constant memory.
3. Irregularly-Sampled Time Series: Latent ODEs
Traditional sequence models (RNNs, Transformers) require observations to be aligned along uniform discrete time grids $Delta t_i = text{constant}$. Real-world telemetry—medical intensive care records, ecological sensor streams, astronomical telescope transits—is notoriously irregular, with variable gaps spanning minutes to weeks.
Neural ODEs naturally ingest irregular timestamps $(t_0, t_1, dots, t_N)$: the hidden state evolves autonomously under continuous dynamics $dh/dt$ during observation gaps, jumping discretely upon observing a new event via Bayesian filtering updates. This established Latent ODEs as the premier theoretical model for continuous physical forecasting.
4. Augmented Neural ODEs and Topological Knots
In 2019, Emilien Dupont et al. exposed a topological limitation of pure Neural ODEs: because ODE trajectories cannot intersect in phase space (by the Picard-Lindelöf Uniqueness Theorem), continuous flows cannot invert or unknot topological configurations (e.g., separating concentric spheres or learning the 1D reflection $f(x) = -x$).
Augmented Neural ODEs (ANODEs) solve this by embedding the state into an augmented dimension $h_{text{aug}}(t) = [h(t), mathbf{0}] in mathbb{R}^{d + p}$. By providing extra dimensions for trajectories to maneuver around one another without intersection, ANODEs become universal continuous function approximators while significantly smoothing the learned vector field and accelerating solver evaluation speeds.
