Neural signature kernels as infinite-width-depth-limits of controlled ResNets

Introduction

The symbiosis between differential equations and deep learning has become an active research area in recent years, notably through the introduction of hybrid models named neural differential equations (Kidger, 2022). In fact, many standard neural network architectures may be interpreted as approximations to some differential equations.

This approximation has been treated rigorously for finite-width ResNets, which in the infinite-depth limit converge in distribution to zero-drift neural stochastic differential equations (Neural SDEs) with diffusion depending on the choice of activation function (Cohen et al., 2021; Hayou, 2022; Marion et al., 2022; Cont et al., 2022).

The \saydual scenario of finite-depth and infinite-width neural networks has also been the object of many recent studies (Neal, 2012; Matthews et al., 2018; Novak et al., 2018). Notably, through the unifying algorithmic language of Tensor Programs designed by Yang (2019), many standard feedforward, convolutional and recurrent architectures of finite-depth can be shown to converge to Gaussian processes (GPs) in the infinite-width limit.

In the context of deep learning for sequential data, finite-width RNNs have been informally identified as approximations to neural controlled differential equations (Neural CDEs) introduced by Kidger et al. (2020); Morrill et al. (2021) and inspired from the homonymous class of dynamical systems studied in rough analysis, a branch of stochastic analysis providing a robust solution theory for differential equations driven by irregular signals (Lyons, 1998; Lyons et al., 2007; Friz & Hairer, 2020; Friz & Victoir, 2010).

However, contrarily to this widespread interpretation, the Euler discretization of a Neural CDE with vector fieldsTypically $f$ is taken to be a randomly initialized feedforward neural network with Gaussian weights and biases. $f$ produces a recursive relation for the hidden state $h$ in the form of (1), where the increments of the input signal $x$ enter the recursion in a multiplicative manner rather than via an additive interaction typically assumed in RNNs:

Furthermore, the addition of the previous hidden state $h_{k}$ on the right-hand side of (1), commonly referred to as a skip connection, is characteristic of ResNets and absent in classical RNNs. We will refer to architectures defined by (1) as homogeneous controlled ResNets. We will also consider their inhomogenous counterparts where the map $f=f(k,h_{k})$ depends on the iteration $k$.

Dynamical systems in the form of (1) are often called reservoirs in the paradigm of reservoir computing (Tanaka et al., 2019; Lukoševičius & Jaeger, 2009; Verstraeten et al., 2007). Contrarily to deep learning, in reservoir computing, only the final readout linear map is trained, while the function $f$ is randomly sampled but remains untrained.

It is worth noting that Neural CDEs are deep learning models that map between infinite dimensional spaces of continuous paths. Therefore, if these continuous models converge, in the infinite-width limit, to some limiting GPs, the latter should be equipped with kernel functions indexed on the same spaces of continuous paths. In the sequel we will demonstrate that controlled ResNet indeed behave like such GPs in the large width-depth regime.

Our objective here is to provide a rigorous mathematical analysis of the behavior of controlled ResNets that are randomly initialised with Gaussian weights and biases in the large width and depth regimes. More specifically:

We prove that both in the infinite-width-depth limit these architectures converge weakly to GPs with limiting kernels satisfying certain (possibly non-linear) partial differential equations varying according to the (in)homogeneity of the network and to the choice of activation function $\varphi$ (see Table 1). Moreover, we show that under some further conditions on the regularity of the driving paths the limits commute, i.e. the limiting GP is unchanged upon reversing the order of the limits. We name this new class of kernel neural signature kernels.

In the case where the system is homogeneous and $\varphi$ is the identity, we show that the equation reduces to a linear PDE and the limiting kernel is proportional to the signature kernel introduced in (Salvi et al., 2021a).

We then prove that in the infinite-depth regime, finite-width controlled ResNets converge in distribution to Neural CDEs with random vector fields. In the inhomogeneous case, these fields behave as a matrix-valued Brownian motion, while for homogeneous networks they are time-independent and Gaussian.

2 Notation

of continuous paths with a square integrable derivative.

We will consider partitions $\mathcal{D}=\{0=t_{0}<\cdots<t_{M}=1\}$ of the interval $$and write their length as$\left\lVert\mathcal{D}\right\rVert:=M$and their mesh size as$|\mathcal{D}|:=\max\limits_{i=1,\dots,\left\lVert\mathcal{D}\right\rVert}|t_{i}-t_{i-1}|$.

The explicit form of $V_{\varphi}$ changes significantly depending on the activation function. We list it for a restricted class of them in Proposition A.7 in the appendix.

The paper is organized as follows: in Section 2 we discuss some related work, in Section 3 we study the inhomogeneous version and in Section 4 we analyse the homogeneous version of controlled ResNets. We conclude in Section 5 with numerical results validating our claims. All proofs can be found in the appendix.

Related Work

Results relating infinite-width limits of neural networks to GPs have been extended from shallow networks (Neal, 2012) to richer architectures of feedforward (Lee et al., 2018; Matthews et al., 2018), convolutional (Novak et al., 2018; Garriga-Alonso et al., 2018) and recurrent (Alemohammad et al., 2020) type. This line of work culminated with the framework of Tensor Programs formulated in Yang (2019) which offers an algorithmic procedure to systematically compute the limiting GP kernels for a wide range of different architectures. One major advantage of this formalism is that it makes it possible to consider weight sharing between layers, something mostly avoided in previous literature but central in the types of systems we consider here.

The reverse scenario of infinite-depth limit of finite-width architectures has mainly been explored for ResNets. In this regime, appropriately rescaled ResNets have been shown to behave like stochastic differential equations (SDEs) (Chen et al., 2018; Cohen et al., 2021; Marion et al., 2022). In particular, Hayou (2022) considers the simpler setting where the architecture does not exhibit weight sharing (we refer to this setting as inhomogeneous) and single out limiting kernels of exponential type. In what follows, we will show that the exponential nature of the limiting kernels is somewhat kept intact even when the ResNets architectures are controlled by an external stream of information, but that the limiting kernels have more structure, particularly in the homogeneous case.

As anticipated in the introduction, ResNets that are controlled by sequential data streams are generalised forms of RNNs and correspond to Euler discretizations of Neural CDEs and variants (Kidger et al., 2020; Morrill et al., 2021; Salvi et al., 2022; Fermanian et al., 2021). These models offer a memory-efficient way to model functions of potentially irregular signals in continuous-time and have achieved state-of-the art performance on a wide range of time series tasks (Singh et al., 2022; Bellot & Van Der Schaar, 2021; Morrill et al., 2021). They stem from the well-understood mathematics of controlled differential equations, which are the central objects studied in rough analysis.

Rough path theory introduced by Lyons (1998) is a modern mathematical framework focused on making precise the interactions between highly oscillatory signals and non-linear dynamical systems. The theory provides a deterministic toolbox to recover many classical results in stochastic analysis without resorting to specific probabilistic arguments. Notably, it extends Itô’s theory of SDEs far beyond the semi-martingale setting and it has had a significant impact in the development of the theory of regularity structures by Hairer (2014), providing a mathematically rigorous description of many stochastic PDEs arising in physics.

More recently, interest has grown rapidly to develop machine learning algorithms based on rough path theoretical tools, particularly in the context of time series analysis (Kidger et al., 2019; Arribas et al., 2020; Lemercier et al., 2021b). The signature, a centrepiece of the theory, provides a top-down description of a stream; it captures crucial information such as the order of different events occurring across different channels, and filters out potentially superfluous information, such as the sampling rate of the signal.

In reservoir computing, the trajectory of a dynamical system is described through its interaction with a random dynamical system that is capable of storing information. In rough path theory the random system is replaced by a deterministic system given by the signature. Recently (Cuchiero et al., 2021a, b) have investigated empirically the idea of a continuous-time reservoir through the randomization of the signature yielding controlled residual architectures similar to the ones of interest to us.

A significant effort has been made to scale methods based on the signature to high dimensional signals. Signature kernels are defined as inner products of signatures and provide an elegant solution to this challenge thanks to the recent development of specific kernel tricks (Király & Oberhauser, 2019). Notably, Salvi et al. (2021a) establish that the signature kernel can be computed efficiently by solving a linear PDE. Algorithms based on signature kernels have been used in a wide range of applications including hypothesis testing (Salvi et al., 2021b), cybersecurity (Cochrane et al., 2021), and probabilistic forecasting (Toth & Oberhauser, 2020; Lemercier et al., 2021a) among others.

Inhomogeneous controlled ResNets

We begin by considering the case of inhomogeneous controlled ResNets. Contrarily to what one might expect, although in this setting the residual map changes at each iteration, the limiting kernels will be governed by simpler differential equations than their homogeneous counterparts, as it can be observed in Table 1. At an intuitive level, this fact can be justified by noting that sharing common random weights and biases throughout all iterations introduces a more intricate dependence structure on the dynamics of the system than if the weights and biases were independently sampled at each iteration.

with time step $\Delta t_{i}=(t_{i}-t_{i-1})>0$.

Here $\sigma_{a},\sigma_{A}>0$ and $\sigma_{b}\geq 0$ are all model hyperparameters.

The time scaling $\frac{1}{\Delta t_{i}}$ in the random weights and biases is crucial as it is exactly the scaling one needs to get an Itô diffusion in the distributional infinite-depth limit, as we will prove in Theorem 3.3 below.

2 The infinite-width-depth regime

The first problem we are interested in studying is that of characterizing the limiting behavior of these neural networks in the infinite-width-then-depth regime.

with initial condition $\kappa_{\varphi}^{x,y}(0)=\sigma_{a}^{2}$.

The statement about commutativity of limits is proved in Appendix B.3, after a characterization of the infinite-depth limit under these more stringent regularity assumptions, by proving that the distributional limit in depth is uniform in width. This generalizes the results of (Hayou & Yang, 2023) in our more complex case.

In some cases we can explicitly characterize the limiting kernels by solving analytically the differential equation (3), as stated in the following corollary We note that these characterizations expressed by means of an exponential are consistent with the results of (Hayou, 2022)..

With the same notation and assumptions as in Theorem 3.1, upon taking $\varphi=id$ the limiting kernel admits the following explicit expression

If $\varphi=ReLU$ and $x=y$, then the limiting kernel satisfies

In Lemma B.12 in the appendix we show that in the kernels governed by the dynamics (3), the parameters $\sigma_{A}$ and $\sigma_{b}$ satisfy the following path-rescaling symmetry

Next we show that infinite-depth, finite-width networks are solutions of SDEs where the vector fields are controlled by the input stream. We will then specialise to the case $\varphi=id$ and identify the limiting kernel with $\kappa_{id}$ from Corollary 3.2

3 The finite-width, infinite-depth regime

Our next result states that when their width $N$ is fixed, these networks converge to a well defined distributional limit as their depth $M$ tends to infinity. In particular, in this limit, the random weights behave like white noise, and thanks to the careful choice of time scaling we have made, the limit is in fact a zero-drift Itô diffusion with diffusion coefficient depending on the driving path.

The idea is proving that the finite difference scheme defining the inhomogeneous architecture gets closer and closer, as the mesh size of the partition becomes finer, to a Euler discretization of Equation (4). One then concludes with standard results which guarantee the convergence of Euler discretizations to the relative SDE’s solution. ∎

Equation (4) can be easily rewritten in more standard SDE form as follows (see Appendix for more details)

Passing directly to the infinite-width limit is not as easy as it could seem, Tensor Program arguments do not apply any longer since they are built for discrete layers and ”collapse” in the continuous case we have to now work with. In simpler cases the limit can be found using McKean-Vlasov arguments as in (Hayou, 2022) and we conjecture that similar results can be found in this more general setting. We leave such a study to future work.

In any case it is possible to directly prove this in the simplest case, when $\varphi=id$. The result is proved in Appendix B.2.2.

Homogeneous controlled ResNets

In this section we consider the more complex setting of networks in which the weights are shared across layers. We will see that this weight-sharing feature will yield limiting kernels governed by two-parameter, non-local partial differential differential equations. We will follow a similar structure as in the previous section, commenting on the crucial differences along the way.

As done in the homogeneous case, we now study the limiting behavior of homogeneous controlled ResNets in the infinite-width-depth limit; as in the in-homogeneous case of the previous section, we will show that the limits commute.

2 The infinite-width-depth regime

and with initial conditions for any $s,t\in$

If moreover $\varphi$ is Lipschitz the limits can be exchanged and

It is non-trivial to show not only that the problem is well-posed but even that equation (6) is well defined because the \sayinstantaneous rate of change $\partial_{s}\partial_{t}\mathcal{K}^{x,y}_{\varphi}(s,t)$ at times $s<t$ depends on the \saypast values $\mathcal{K}_{\varphi}^{x,x}(s,s)$, on the \saypresent values $\mathcal{K}_{\varphi}^{x,y}(s,t)$ and on the \sayfuture values $\mathcal{K}_{\varphi}^{y,y}(t,t)$. The nonlocal nature of these dynamics is such that it is a priori not clear that the RHS of (6) even has meaning since the matrix $\Sigma_{\varphi}^{x,y}(s,t)$ could be not positive semidefinite and $V_{\varphi}$ is only defined on PSD matrices.

Similarly to Theorem 3.1, due to the complexity of the arguments, this result is proved in several steps. In in Appendix C.1. The first step consists of showing that the infinite width limit is well defined for any choice of $\mathcal{D}_{M}$. This will be a GP defined by a kernel $\mathcal{K}_{\mathcal{D}_{M}\times\mathcal{D}_{M}}$ found as the terminal value of a finite difference scheme having the same form as that of a Euler discretization, on $\mathcal{D}_{M}\times\mathcal{D}_{M}$, of equation (6). The second step consists in proving that the kernels $\{\mathcal{K}_{\mathcal{D}_{M}\times\mathcal{D}_{M}}\}_{M}$ constitute, in a suitable metric, a Cauchy sequence as $|\mathcal{D}_{M}|\to 0$, that the limit is independent from the chosen sequence of partitions and that it does indeed uniquely solve equation (6). The final step, proved in Appendix C.3, concerns the exchange of limits. After a characterization of the infinite-depth limits, we will prove that the distributional limit as depth goes to infinity is uniform in the width, thus we will be able to use the classical Moore-Osgood theorem to justify the exchange. ∎

When the activation function $\varphi$ is the identity, equation (6) reduces to a linear hyperbolic PDE. Upon inspection, we unveil a surprising link with the signature kernel, a well-studied object in rough analysis corresponding to an inner product between two path-signatures, and that was shown by Salvi et al. (2021a) to satisfy a similar PDE. This is the content of the next corollary.

Using the same notation and assumptions as in Theorem 4.1, choosing $\varphi=id$ the limiting kernel satisfies the following identity

with initial conditions $k_{sig}^{x,y}(s,0)=k_{sig}^{x,y}(0,t)=1$.

In other words we have unveiled a novel family of kernels indexed on continuous paths which generalizes the signature kernel in (Salvi et al., 2021a). We name this new class of kernels neural signature kernels. We note that this generalization is done directly at the level of the driving PDE unlike the extensions studied in (Cass et al., 2021) which use a different inner product structure on the space where signatures live.

Analogously to the inhomogeneous case, the parameters $\sigma_{A}$ and $\sigma_{b}$ defining the neural signature kernels governed by the dynamics in equation (6) satisfy the following path-rescaling symmetry

For non-linear activation functions $\varphi$, the neural signature kernel non-linear PDE (6) and the $id$-neural signature kernel linear PDE (7) might in principle admit the same solution, which would mean essentially that linear and non-linear controlled ResNets behave in the same way in the infinite-width-depth regime. In in Figure 1 we show empirically that this is not the case in general, by comparing the limiting empirical distributions for $\varphi=id$ and $\varphi=ReLU$.

3 The finite-width, infinite-depth regime

Our next result states that when their width $N$ is fixed, homogeneous controlled ResNets converge in distribution, in $$, to a Neural CDE with random vector fields.

If we fix $a$, the $A_{k}$s and the $b_{k}$s to be the same for all $\mathcal{D}_{M}$ then we have uniform convergence by classical results. The rate of convergence can be bounded with some constants depending on the entries of $a$, $A_{k}$, $b_{k}$ and which, thanks to Gaussianity, have finite expectation. It is just a matter of applying the classical portmanteau lemma to conclude. ∎

This can be naturally extended in order to take into consideration the joint distribution for different input choices.

The solutions to equation (8) have been informally introduced in (Cuchiero et al., 2021b; Akyildirim et al., 2022) as lower dimensional approximations of path-signatures, and have been dubbed by the authors randomized signatures.

Taking directly the infinite width limit is once again far from trivial, reasoning à la Tensor Program quickly collapse and there is no clear possible future path corresponding to the McKean-Vlasov ideas for the inhomogeneous case. The problem is that the randomness is in the vector fields themselves and not in the driving paths, courtesy of the cross-layer dependencies in the homogeneous networks. This is why it’s necessary to sidestep the problem by proving the existence of uniform convergence bounds.

4 The infinite-depth-then-width regime: φ=i​d𝜑𝑖𝑑\varphi=id

As anticipated, contrary to the inhomogeneous case, in the current homogeneous setting, when $\varphi$ is the identity, we are able to prove directly that the limits in Equation 5 commute as well as explicit convergence bounds.

on $^{2}$. Moreover the convergence is of order $\mathcal{O}(\frac{1}{N})$.

Numerics

In this section, we first illustrate theoretical results established in Section 4 and then outline numerical considerations to scale the computation of signature kernels.

We start by illustrating the convergence in distribution of a homogeneous controlled ResNet to a GP endowed with neural signature kernel as per Theorem 4.1. To this aim, we consider a homogeneous ResNet $\Phi_{\varphi}^{M,N}$ with activation function $\varphi=\text{ReLU}$, and $(\sigma_{a},\sigma_{A},\sigma_{b})=(0.5,1.,1.2)$. For $R=250$ realizations of the weights and biases, we run the model on a $2$-dimensional path $x:t\mapsto(\sin(15t),\cos(30t)+3e^{t})$ observed at $100$ regularly spaced time points in $$. We then verify that, as$N$increases,$\Phi_{\varphi}^{M,N}(x)$converges to a Gaussian random variable with mean zero and variance$\mathcal{K}_{\varphi}(x,x)$. This limiting variance is computed by solving Equation 6 on a fine discretization grid. As it can be observed on Figure 3 the Gaussian fit for this one-dimensional marginal gets better as$N$ increases. Further results can be found in the appendix.

2 Scaling signature kernels

The signature kernel of two paths is typically computed by approximating the solution of the PDE in (7) on a $2$-dimensional time grid, which scales quadratically with the discretization step of the solver. Although an efficient numerical scheme leveraging GPU computations to update the solution at multiple time points on the grid in parallel has been proposed in Salvi et al. (2021a), the maximum number of threads in a GPU block imposes a hard limit on the discretization step of the solver, limiting the applicability of signature kernel methods to long time series. Theorem 4.4 offers a new way to compute the signature kernel by solving two CDEs linearly in time instead of one PDE quadratically in time; one would first run a wide and infinite-depth ResNet on the two control paths of interest, and then compute the (rescaled) dot-product between the outputs of the penultimate layer. This approach allows for more flexibility regarding the choice of path interpolations and numerical solvers, as several options are made readily available in dedicated python packages such as $\mathsf{torchcde}$ (Kidger et al., 2020). Next, we describe possible ways to increase further the scalability of this approach.

To further improve scalability of Neural CDEs for long time series Morrill et al. (2021) made use of the so-called log-ODE scheme to forward-solve the differential equation on much larger time intervals than the ones that would be expected given the sampling rate or length of the data. We leave the investigation of this numerical scheme for computing signature kernels as future work.

The forward pass of a ResNet involves several ($M\times d$ where $M$ is the number of time steps, and $d$ the dimension of the input path) matrix-vector multiplications where the entries of each $N$-by-$N$ matrix are Gaussian distributed. As remarked in (Dong et al., 2020), in the context of random RNNs, to speed-up these computations, the dense weight matrices can be replaced by structured random matrices given by the products of random (binary) diagonal matrices and Walsh-Hadamard matrices. The complexity of the matrix-vector product can be reduced to $\mathcal{O}(N\log N)$ leveraging the fast Hadamard transform algorithm (without sampling the Walsh-Hadamard matrices).

Several machine learning use cases of the signature kernel have provided empirical evidence that embedding the input paths pointwise in time in a feature space can be beneficial to increase the performance of kernel methods on sequential data. In particular, when the paths evolve in a Euclidean space, the RBF kernel often turns out to be a good choice. Although this embedding is infinite-dimensional, random Fourier features (Rahimi & Recht, 2007) make it possible to approximate it by a finite-dimensional one. One could then investigate randomly initialized ResNets, controlled by sequences of such approximate embeddings.

Conclusion and future work

In this paper we considered controlled ResNets defined as Euler-discretizations of Neural CDEs. We showed that in both the infinite-depth-then-width and in the infinite-width-then-depth limit, these converge weakly to the same GP indexed on path space endowed with neural signature kernels satisfying certain (possibly non-linear) PDEs varying according to the choice of activation function $\varphi$. In the special case where $\varphi$ is the identity, we showed that the equation reduces to a linear PDE and the limiting kernel agrees with the signature kernel. In this setting, we also provided explicit convergence rates. Finally, we showed that in the infinite-depth regime, finite-width controlled ResNets converge in distribution to Neural CDEs with random vector fields which are either time-independent and Gaussian, if the system is homogeneous, or behave like a matrix-valued Brownian motion, if the system is inhomogeneous.

We believe that a rigorous investigation of the functional analytic properties of the reproducing kernel Hilbert spaces (RKHSs) associated to the new family of neural signature kernels is also a compelling future research direction. In particular, it would allow to build an understanding of the expressivity and generalization properties of these kernels.

In the homogeneous setting, the vector fields are constant functions while in the in-homogeneous setting they are described by white noise. Investigating the intermediate regularity cases is an interesting avenue for future research; for example considering matrices and biases sampled from of Fractional Brownian Motion increments with Hurst exponent $H\in$ (the inhomogeneous case corresponds to the case $H=0.5$ while the homogeneous one to $H=1$).

Last but not least, establishing expressions and analyzing the associated Neural Tangent Kernels (NTK) (Jacot et al., 2018; Yang, 2020) would provide quantitative insights on the training mechanism of Neural CDEs by gradient descent.

All the experiments presented in this paper are reproducible following the code at https://github.com/MucaCirone/NeuralSignatureKernels

Aknowledgements

The authors would like to thank Thomas Cass, James-Micheal Lehay and David Villringer for helpful discussions.

NMC was supported by EPSRC Centre for Doctoral Training in Mathematics of Random Systems: Analysis, Modelling and Simulation (EP/S023925/1) and the Department of Mathematics, Imperial College London, through a Roth Scholarship. ML was supported by the EPSRC grant EP/S026347/1.

References

Appendix A Preliminaries

In this section we are going to state some preliminary results and considerations which we are going to refer to through the entirety of the text.

Since we are mainly interested in studying time-series, our data space will be a space of paths, more specifically we are going to consider the space

This norm is equivalent to the induced norm from \big{(}W^{1,2}()\big{)}^{d} since

A.2 Assumptions on the activation function

Here are the crucial assumptions we make on the activation:

Note how we have proved also that, under the linear boundedness assumption, using the same final bound

Let $PSD_{2}$ denote the set of $2\times 2$ positive semidefinite matrices

is $\kappa_{R}$-Lipschitz for some $\kappa_{R}>0$, i.e.

This is the content of Theorem F.4 in (Novak et al., 2018). ∎

where the first inequality follows from Jensen’s inequality, the second from Assumption A.1, the third from the triangle inequality, and the fourth from the fact that $\mathcal{N}(0,Id)$ has finite moments with $\bar{M}$ a constant incorporating $M^{2}$ and these bounds.

We end the section showing the explicit characterization of the maps $V_{\varphi}$ for some selectedby the availability in the literature. activation functions. We write $V_{\varphi}$ with the obvious meaning.

Defining $\gamma(\Sigma):=\frac{[\Sigma]_{1}^{2}}{\sqrt{[\Sigma]_{1}^{1}[\Sigma]_{2}^{2}}}$ we have

Appendix B Proofs for inhomogeneous controlled ResNets

In this section of the appendix we are going to prove all the results stated for the inhomogeneous case. The section will be subdivided in three main parts: in the first we consider the infinite-width-then-depth limit, in the second the infinite-depth-then-width one, in the final one we prove the commutativity of the integrals.

We start by recalling the defintion of the model.

with time step $\Delta t_{i}=(t_{i}-t_{i-1})>0$ and parameters $\sigma_{a},\sigma_{A}>0$ and $\sigma_{b}\geq 0$.

The main goal in this subsection is to prove the first part of Theorem 3.1, which we restate here:

with initial condition $\kappa_{\varphi}^{x,y}(0)=\sigma_{a}^{2}$.

The infinite width-convergence of the finite dimensional distributions

for a fixed depth $M$ to those of a Gaussian process $\mathcal{GP}(0,\kappa_{\mathcal{D}_{M}})$ defined by a kernel computed as the final value of a finite difference scheme on the partition $\mathcal{D}_{M}$. This will be done using Tensor Programs (Yang, 2019), and is the content of subsection B.1.1.

The infinite-depth (uniform) convergence of the discrete kernels $\kappa_{\mathcal{D}_{M}}$ to a limiting kernel $\kappa_{\varphi}$ which solves the differential equation (12). This will be done in subsection B.1.2.

In the next theorem, we show that finite width inhomogeneous ResNets converge in distribution to GPs with discrete kernels satisfying some difference equations.

with $\kappa_{\mathcal{D}_{M}}^{x,y}(0)=\sigma_{a}^{2}$ and where

We use (Yang, 2019)[Corollary 5.5] applied to the Tensor Program of Algorithm 1 where the input variables are independently sampled according to

with $\mu,\Sigma$ computed according to (Yang, 2019)[Definition 5.2] and defined on the set of all G-vars in the program i.e.

where $h=\phi((g_{k})_{k=1}^{m})$ for some function $\phi$ and $\phi(Z):=\phi((Z^{g_{k}})_{k=1}^{m})$, similarly for $g^{\prime}$.

In our setting $\mu^{in}\equiv 0$ since all Input variables are independent, from which $\mu\equiv 0$; furthermore $\Sigma^{in}(g,g^{\prime})=0$ except if $g=g^{\prime}$ when it takes values in $\{\sigma_{a}^{2},\frac{\sigma^{2}_{b}}{t_{l}-t_{l-1}},1\}$ accordingly.

Following the rules of $\Sigma$, assuming $l_{i},l_{j}\in\{1,\dots,M\}$, we obtain

In particular we see that if $l_{i}\neq l_{j}$ then $\Sigma(\mathcal{S}^{x_{i}}_{l_{i}},\mathcal{S}^{x_{j}}_{l_{j}})=\Sigma(\mathcal{S}^{x_{i}}_{l_{i\wedge j}},\mathcal{S}^{x_{j}}_{l_{i\wedge j}})$. Thus if we set, for $t_{l_{i}}\in\mathcal{D}_{M}$,

which is exactly what Equation 13 states. Then note how

There are two things to notice, done in the above proof in order to satisfy the required formalism:

In the program the output projector $v$ is sampled according to $\mathcal{N}(0,1)$ while the original $\phi\sim\mathcal{N}(0,\frac{1}{N})$. This does not pose any problems since the output of the formal programs uses ${v}/{\sqrt{N}}\sim\mathcal{N}(0,\frac{1}{N})$.

The input paths $x_{i}$ enter program 1 not as Inputs but as coefficients of LinComb, this means that for any choice of input paths we must formally consider different algorithms. In any case, for any possible choice, the result has always the same functional form; hence a posteriori it is legitimate to think about one algorithm.

Actually we have proved the even stronger statement, in the sense that the previous result holds for intermediate times too:

For all $t_{m},t_{n}\in\mathcal{D}_{M}$ one has the following distributional limit

and the matrices $\Sigma_{\mathcal{D}_{M}}^{x,y}(t_{n})$ are always in $PSD_{2}$.

B.1.2 Uniform convergence of discrete kernels

We extend in a similar way the matrix $\Sigma_{\mathcal{D}_{M}}^{x,y}$ defined in equation (14).

The main idea for proving Theorem B.6 will be to show that the sequence $\{\kappa_{\mathcal{D}_{M}}^{x,y}\}_{M}$ is uniformly bounded and uniformly equicontinuous so that Ascoli-Arzelà theorem applies. One then proves that the limit of the resulting subsequence is the unique solution of Equation 16 and that the whole sequence converges to it. Before proving Theorem B.6 we need several lemmas.

The first step is to establish a uniform lower bound.

Setting, for any $t_{m}\in\mathcal{D}$ by definition of the kernel $\kappa_{\mathcal{D}}^{x,x}$

For any $t\in$, by definition, we can write

which is $\geq 0$ by the same arguments as above. ∎

The second step is to establish a uniform upper bound.

By Gronwall inequality ((Friz & Victoir, 2010), Lemma 3.2) we have

To prove a similar inequality for $\kappa_{\mathcal{D}}^{x,x}$, consider

The following lemma shows that the kernels are in fact elements of $PSD(R)$.

We now extend the results of Lemma B.8 and Lemma B.9 to the case $x\neq y$.

Moreover for $R_{\alpha}:=C_{\alpha}\vee\sigma_{a}^{-2}$ we get

For $\kappa_{\mathcal{D}}^{x,y}(t)$ we proceed similarly to before: when $t\in[t_{m},t_{m+1})$ one has

In particular as $|\mathcal{D}|\rightarrow 0$ we have

By Lemma B.10 the sequence of functions $\{\kappa_{\mathcal{D}_{M}}^{x,y}\}_{M}$ is uniformly bounded. Assume $s<t$ then one has, with the same bounds just used, that

Here we prove that if the PDE (16) admits a solution this must be unique. Assume the existence of different solutions $K=(K_{x,x},K_{x,y},K_{y,y})$ and $G=(G_{x,x},G_{x,y},G_{y,y})$ with all the $\Sigma^{K}$ and $\Sigma^{G}$ in $PSD_{2}$ . From Eq (16) it is clear that $K_{x,x},K_{y,y},G_{x,x},G_{y,y}\geq\sigma_{a}^{2}$ and, by continuity, that they are bounded by some constant; thus all the $\Sigma_{K}$ and $\Sigma_{G}$ are in some in $PSD_{2}(\bar{R})$. Then by the Lipschitz property of $V_{\varphi}$ one sees that

thus $\left\lVert\Sigma_{K}^{x,y}(t)-\Sigma_{G}^{x,y}(t)\right\rVert_{\infty}=0$ by Gronwall for all $t\in$ i.e $K=G$.

We now need to prove that the limit $\kappa^{x,y}_{\varphi}$ of the subsequence $\{\kappa_{\mathcal{D}_{M}}^{x,y}\}_{k}$ solves the PDE, it will then follow that any sub-sequence $\{\kappa_{\mathcal{D}_{M}}^{x,y}\}$ admits a further sub-sequence converging to the same map $\kappa^{x,y}_{\varphi}$, giving us the convergence of the whole sequence. Thus without loss of generality we can assume in the sequel that the whole sequence converges.

Let us prove that the limit $\kappa_{\varphi}^{x,y}$ is, in fact, a solution of the PDE. Let $t\in[t_{m},t_{m+1})$ for some fixed $\mathcal{D}$, then

Now, considering the sequence $\mathcal{D}_{M}$, by convergence

Moreover, setting $m_{M}$ to be such that $t\in[t_{m_{M}},t_{m_{M}+1})$ in $\mathcal{D}_{M}$, we have that by Lebesgue differentiation theorem

almost surely as a function of $s$, thus using Dominated convergence we conclude that

B.1.3 Proof of Theorem 3.1: Part 1

It is finally time to prove the first part Theorem 3.1 in the main body of the paper, which we restated in the appendix as Theorem B.2.

where $\kappa_{\mathcal{D}_{M}}(x_{\alpha},x_{\beta})=\kappa_{\mathcal{D}_{M}}^{x_{\alpha},x_{\beta}}(1)$ for all $\alpha,\beta=1,\dots,n$. Thus to conclude we just have to prove that, still in distribution, it holds

This is just a matter of computation. For the case $\varphi=id$ one has

Notice then that substituting (17) for $K^{id}_{s}(x,y)$ in the integral leads to

which means, by uniqueness of solutions, that the thesis holds.

For what concerns the case $\varphi=ReLU$ notice that for one dimensional Gaussian variables centered in the origin one has

which equals $\kappa_{id}^{x,x}(t;\sigma_{a},\frac{\sigma_{A}}{\sqrt{2}},\sigma_{b})$. ∎

We conclude this section by proving Remark Remark about the path scaling symmetry mentioned in the main paper.

For all choices $(\sigma_{a},\sigma_{A},\sigma_{b})$ and for all $\varphi$ as in Theorem 3.1 we have, with abuse of notation and the obvious meaning, that

Thus the respective triplets solve the same equation and we can conclude by uniqueness. ∎

B.2 The infinite-depth-then-width regime

As mentioned in the paper, it is natural to ask what happens if the order of the width-depth limits in is reversed.

We begin by proving Theorem 3.3 which we restate for the reader’s convenience. We will follow arguments used in (Hayou, 2022), extending the results obtained therein.

We will transform Equation (19) in a usual SDE form to then use the classical result (Kloeden & Platen, 1992)[Theorem 10.2.2] to prove the convergence of Euler Discretizations to the unique solution. We first show that Equation (19) can be re-written as follows

for $i=1,\dots,N$ and $j=N(N+1)(k-1)+N(m-1)+1,\dots,N(N+1)(k-1)+Nm$. It is then just a matter of checking the required conditions to apply (Kloeden & Platen, 1992)[Theorem 10.2.2]:

thus using sublinearity of $\varphi$ for some $M>0$ one has

thus one just uses Lipschitz property of $\varphi$ with Lemma A.3.

Finally for (iii) one computes just as before

is not the Euler discretization of the above SDE. However, after classical though tedious calculations it is easy to show that

we straightforwardly extend the previous result to the case with multiple inputs considered at the same time.

B.2.2 Infinite-depth-then-width limit: φ=i​d𝜑𝑖𝑑\varphi=id

Here we directly prove that, in the case $\varphi=id$, the covariances of the infinite-depth networks converge to $\kappa_{id}$ as the width increases without bounds.

Thanks to the extension of Theorem B.13 to the multi-input case we have

where the penultimate equality follows from Itô’s isometry. Hence

meaning that $K^{N}_{t}(x_{k},x_{m})$ solves Equation (16).

Since we showed in the proof of Theorem B.6 that the unique solution to the equation is $\kappa_{id}^{x_{k},x_{m}}(1)$ we conclude.

Directly proving the convergence in distribution of the infinite-width networks to centered Gaussians with covariance function $\kappa_{id}$ is difficult since the networks are not Gaussian processes. However we believe that it is possible to employ McKean-Vlasov arguments to prove that they are at the limit, to then obtain a commmutativity of the limits; we leave the exploration of this direction to future work.

B.3 Commutativity of Limits

In this section we are going to prove the commutativity of limits in the inhomogeneous case. To do this we are going to proceed similarly to (Hayou & Yang, 2023) which proves the same result in a much restricted case.

Note that if $\varphi$ is $K$-Lip i.e. $|\varphi(x)-\varphi(y)|\leq K|x-y|$ and $\varphi(0)=0$ then

In this subsection $\varphi$ is considered Lipschitz and with $\varphi(0)=0$.

The bounding constant in equation (20) depends on the path $x$ only trough $\left\lVert\dot{x}\right\rVert_{\infty,}$. In particular the bound is uniform on bounded sets, with respect to this norm.

In particular we can interpolate $\{(t_{m},x_{t_{m}})\}_{m=0,\dots,|\mathcal{D}_{M}|}$ in such a way to have the resulting map $\bar{x}$ with

such that definitely in $M$ $\dot{x}$ is $\frac{1}{2}$-Holder continuous and

One can for example interpolate the points with a polynomial close enough to $\dot{x}$, being the polynomial defined on the compact $$it is Lipschitz on$$ hence $\frac{1}{2}$-Holder continuous. The existence of such a polynomial follows from

Then $\mathcal{S}^{M,N}(x)$ will be exactly the Euler discretisation of $\mathcal{S}^{N}(\bar{x})$ on $\mathcal{D}_{M}$, hence

Using the fact that, again (Kloeden & Platen, 1992)[10.2.2],

where $\bar{K}$ only depends on $\left\lVert\dot{\bar{x}}\right\rVert_{\infty}+\left\lVert\dot{x}\right\rVert_{\infty}$. By our choice of $\bar{x}$ we conclude. ∎

If the activation function is Lipschitz and $\varphi(0)=0$ then there is a constant $C$ depending only on $\left\lVert\dot{x}\right\rVert_{\infty,}$ in an increasing fashion such that:

where $\mu_{t}^{M,N}$ is the distribution of any coordinate of ${\mathcal{S}}^{M,N}_{t}(x)$ and $\mu_{t}^{N}(x)$ that of any coordinate of ${\mathcal{S}}^{N}_{t}(x)$ (they are identically distributed).

If the activation function is Lipschitz and $\varphi(0)=0$ then there is a constant $C$ depending only on $\left\lVert\dot{x}\right\rVert_{\infty,}$ in an increasing fashion such that:

where $\mu_{t}^{M,N}$ is the distribution of $\left\langle v^{N},{\mathcal{S}}^{M,N}_{t}(x)\right\rangle$ and $\mu_{t}^{N}(x)$ that of $\left\langle v^{N},{\mathcal{S}}^{N}_{t}(x)\right\rangle$ for some independent vector with iid entries $[v^{N}]_{\alpha}\sim\mathcal{N}(0,\frac{1}{N})$.

With the same proof we have also the convergence of the laws of the rescaled processes $\frac{1}{\sqrt{N}}{\mathcal{S}}^{M,N}_{t}(x)$.

Note that the exact arguments, being of $L^{2}$ type, can be repeated for the ”stacked” vector $(\mathcal{S}^{N,M}(x_{1}),\dots,\mathcal{S}^{N,M}(x_{N}))$ extending (qualitatively) the bounds to the multi-input case.

If the activation function is Lipschitz and $\varphi(0)=0$ then the limits in Thm. B.2 commute.

By the classical Moore-Osgood theorem we need to prove that one of the two limits is uniform in the other, for example that the limit in distribution as $M\to\infty$ is uniform in $N$ in some metric which describes convergence in distribution. But this is just the content of the previous results, extended to the multi-input case. ∎

Appendix C Proofs for homogeneous controlled ResNets

This section too will be subdivided in three main parts: in the first we consider the infinite-width-then-depth limit, in the second we reverse the order and consider the infinite-depth-then-width limit and in the third one we prove that the limits can be exchanged.

The main goal is that of proving the first part of Theorem 4.1, which we restate here:

and with initial conditions for any $s,t\in$

Once again, it is clearer to subdivide the proof of this result in two parts:

The infinite width-convergence of the finite dimensional distributions

to those of a Gaussian process $\mathcal{GP}(0,\mathcal{K}_{\mathcal{D}_{M}\times\mathcal{D}_{M}})$ defined by a Kernel computed as the final value of a difference equation on the partition $\mathcal{D}_{M}\times\mathcal{D}_{M}$ of $\times$. This will be done using Tensor Programs (Yang, 2019) in C.1.1.

The infinite-depth convergence of the discrete kernels $\mathcal{K}_{\mathcal{D}_{M}\times\mathcal{D}_{M}}$ to a limiting Kernel $\mathcal{K}_{\varphi}$ which solves the differential equation 23. This will be done in C.1.2.

We prove convergence in the finite-depth, infinite-width limit to a GP endowed with discrete kernels. using the formalism of Tensor Programs.

We use (Yang, 2019)[Corollary 5.5] applied to the Tensor Program of Algorithm 2 where the input variables are independently sampled according to

There are two things to notice, done in order to satisfy the required formalism:

In the program the output projector $v$ is sampled according to $\mathcal{N}(0,1)$ while the original $\phi\sim\mathcal{N}(0,\frac{1}{N})$. This does not pose any problems since the output of the formal programs uses ${v}/{\sqrt{N}}\sim\mathcal{N}(0,\frac{1}{N})$.

The input paths $x_{i}$ enter program 2 not as Inputs but as coefficients of LinComb, this means that for any choice of input paths we must formally consider different algorithms. In any case, for any possible choice, the result has always the same functional form; hence a posteriori it is legitimate to think about one algorithm.

Note that $S^{x}_{m}$ stands for the formal variable in program 2, not for $S^{M,N}_{t^{\mathcal{D}_{M}}_{m}}(x)$ even thogh this is the value which it ”stores” for a fixed hidden dimension $N$.

with $\mu,\Sigma$ computed according to (Yang, 2019)[Definition 5.2] and defined on the set of all G-vars in the program i.e.

where $h=\phi((g_{k})_{k=1}^{m})$ for some function $\phi$ and $\phi(Z):=\phi((Z^{g_{k}})_{k=1}^{m})$, similarly for $g^{\prime}$.

In our setting $\mu^{in}\equiv 0$ since all Input variables are independent, from which $\mu\equiv 0$; furthermore $\Sigma^{in}(g,g^{\prime})=0$ except if $g=g^{\prime}$ when it takes values in $\{\sigma^{2}_{a},\sigma^{2}_{b},\sigma^{2}_{A}\}$ accordingly.

Following the rules of $\Sigma$, assuming $m_{i},m_{j}\in\{1,\dots,\left\lVert\mathcal{D}_{M}\right\rVert\}$, we obtain

thus if we set, for $t_{m_{i}},t_{m_{j}}\in\mathcal{D}_{M}$, $\mathcal{K}_{\mathcal{D}_{M}\times\mathcal{D}_{M}}^{x_{i},x_{j}}(t_{m_{i}},t_{m_{j}}):=\Sigma(S^{x_{i}}_{m_{i}},S^{x_{j}}_{m_{j}})$ then we get

which is exactly what Theorem C.2 states. Then note how

Actually we have proved the even stronger statement that

For all $t,s\in\mathcal{D}_{M}$ one has the following distributional limit

where $[\phi^{N}]_{\alpha}\sim\mathcal{N}(0,\frac{1}{N})$ are independently sampled.

For all $t\in\mathcal{D}$, for all $s\in\mathcal{D}^{\prime}$, with a clear abouse of notation, one has the following distributional limit

where $[\phi^{N}]_{\alpha}\sim\mathcal{N}(0,\frac{1}{N})$ are independently sampled and the matrices $\Sigma_{\mathcal{D}\times\mathcal{D}^{\prime}}^{x,y}(t,s)$ are always in $PSD_{2}$.

The fact that the matrices $\Sigma_{\mathcal{D}\times\mathcal{D}^{\prime}}^{x,y}(t,s)$ are always in $PSD_{2}$ will be crucial in the following.

C.1.2 Uniform convergence to neural signature kernels

We now prove that the sequence of discrete kernels convergence uniformly to our neural signature kernels. First we need definitions for extending the discrete kernels to $^{2}$ similarly to the previous section.

as the map satisfying the following recursion

We extend in a similar way the matrix $\Sigma^{x,y}_{\mathcal{D}\times\mathcal{D}^{\prime}}$.

Next we prove that the discrete kernels converge in distribution to a limiting kernel that satisfies a two-parameters differential equation.

where $\mathcal{K}^{x,y}(s,t)$ satisfies the following differential equation

and with initial conditions for any $s,t\in$

There is a constant $C_{x}>0$ independent of $\mathcal{D}$ such that

and we can use Gronwall Inequality ((Friz & Victoir, 2010), Lemma 3.2) to obtain

In order to extend the bound to $\mathcal{K}^{x,x}$ notice that

Thus finally we can conclude, as promised, that

We are ready to extend the result to the general case.

Moreover for $R_{\alpha}:=C_{\alpha}\vee\sigma_{a}^{-2}$ we get

For $K_{x,y}(t)$ we proceed similarly to before:

The second part follows from the definition of $\operatorname{PSD}(R_{\alpha})$. ∎

In particular as $(|\mathcal{D}|\vee|\mathcal{D}^{\prime}|)\rightarrow 0$ we have

Notice that, since we can repeat the argument for $\mathcal{K}_{\mathcal{D}\times\mathcal{D}}^{x,x}$ and $\mathcal{K}^{y,y}_{\mathcal{D}^{\prime}\times\mathcal{D}^{\prime}}$, we have

Consider now another partition $\check{\mathcal{D}}\times\check{\mathcal{D}}^{\prime}$, we have

Assume to be in the case $x=y$, $\mathcal{D}=\mathcal{D}^{\prime}$, $\check{\mathcal{D}}=\check{\mathcal{D}}^{\prime}$ and define the following quantity:

One has, from the previous inequality, that

Coming back to the general case and setting

Given any sequence of partitions $\{\mathcal{D}_{n}\times\mathcal{D}^{\prime}_{n}\}$ with $|\mathcal{D}_{n}|\vee|\mathcal{D}_{n}^{\prime}|\rightarrow 0$, due to the bounds we have just proven we have

and the limit $\mathcal{K}^{x,y}_{\varphi}$ does not depend on the sequence (i.e. the limit exists and is unique).

The limit is indeed unique: assume $K^{x,y}$ and $G^{x,y}$ are limits along two different sequences of partitions $\{\mathcal{D}_{n}\times\mathcal{D}^{\prime}_{n}\}$ and $\{\mathcal{G}_{n}\times\mathcal{G}^{\prime}_{n}\}$; then the sequence $\{\mathcal{P}_{n}\times\mathcal{P}^{\prime}_{n}\}$ such that $\mathcal{P}_{2n}\times\mathcal{P}^{\prime}_{2n}:=\mathcal{D}_{n}\times\mathcal{D}^{\prime}_{n}$ and $\mathcal{P}_{2n+1}\times\mathcal{P}^{\prime}_{2n+1}:=\mathcal{G}_{n}\times\mathcal{G}^{\prime}_{n}$ is still such that $|\mathcal{P}_{n}|\vee|\mathcal{P}^{\prime}_{n}|\rightarrow 0$ hence the associated kernels have a limit which must be equal to both $K^{x,y}$ and $G^{x,y}$.

Since PSD matrices form a closed set we moreover have that the matrices $\Sigma^{x,y}_{\varphi}(s,t)$ obtained as limits using the previous result are all PSD. We can actually say more: they belong to $\operatorname{PSD}(R)$.

We can finally conclude by proving that $\mathcal{K}^{x,y}_{\varphi}$ is, in fact, a solution of the PDE:

- Part V : Uniformity in $x,y$ on bounded sets

defined on $\times$,satisfying Equation (23) and such that

We will write $|K|_{\infty}:=|K^{x,x}|\vee|K^{x,y}|\vee|K^{y,y}|$. To satisfy equation (23) the associated covariance matrices $\Sigma_{K}({\eta,\tau}),\Sigma_{G}({\eta,\tau})$ must be always $PSD_{2}$. Using the assumed bound from below the one from above given by continuity we can assume that they uniformly in time are contained in some $PSD_{2}(\bar{R})$. This true for all choices $(x,x),(x,y),(y,y)$.

Moreover this holds substituting $(x,y)$ with $(x,x)$ and $(y,y)$.

Let $\Xi_{t}:=\sup_{0\leq t_{1},s_{1}\leq\eta t}|K({t_{1},s_{1}})-G({t_{1},s_{1}})|_{\infty})$ then

By Gronwall we finally conclude that $\Xi_{t}=0$ for all $t$, which concludes the proof.z ∎

C.1.3 Proof of Theorem 4.1: Part 1

It is finally time to prove the first part of the main result of the paper, which we restate below for the reader’s convenience.

The proof is now just a matter of combining Theorem C.2 and Theorem C.7.

thus to conclude we just have to prove that, still in distribution, it holds

or, equivalently, that in $R^{n\times n}$ it holds $\lim_{M\to\infty}\mathcal{K}_{\mathcal{D}_{M}\times\mathcal{D}_{M}}(\mathcal{X},\mathcal{X})=\mathcal{K}_{\varphi}(\mathcal{X},\mathcal{X})$.

Like in the inhomogeneous case we can explicitly write the Kernel for the simplest case:

By substituting (29) for $\mathcal{K}_{id}^{x,y}(\eta,\tau)$ in the integral and using (Salvi et al., 2021a)[Theorem 2.5] we have

which means, by uniqueness of solutions (Proposition C.12), that the thesis holds. ∎

A proof similar to that given in the inhomogeneous case for $\varphi=ReLU$ does not work now since the covariance matrix is not, in general, degenerate for $x=y$ as in that case.

For all choices $(\sigma_{a},\sigma_{A},\sigma_{b})$ and for all $\varphi$ as in Theorem 3.1 we have, with abuse of notation and the obvious meaning, that

Follow the exact same steps and arguments of B.12. ∎

C.2 The infinite-depth-then-width regime

Assume the $A_{k}$ and $b_{k}$ to be fixed for all choices of $\mathcal{D}_{M}$. The system has unique solution by (Friz & Victoir, 2010)[Theorems 3.7, 3.8] with

noting that these are Lipschitz since composition of Lipschitz and linearly bounded.

By reasoning as in the proof of Theorem 4.1 one gets that

Taking the expectation over $A_{k}$,$b_{k}$ and the initial condition leads to

We call Randomized Signatures the solutions to

These are the same objects defined in (Cuchiero et al., 2021b).

C.2.2 Infinite-depth-then-width limit: φ=i​d𝜑𝑖𝑑\varphi=id

moreover the convergence is of order $\mathcal{O}(\frac{1}{N})$

C.3 Commutativity of Limits

In this section we are going to prove the commutativity of limits in the homogeneous case. The core arguments are the same as those employed for the inhomogeneous counterpart, this time however we won’t be able to take advantage of ready-made results from stochastic analysis, thus we are going to carefully obtain bounds in more direct ways.

Note how the two extensions coincide on $\mathcal{D}_{M}$.

Assume the activation function $\varphi$ is Lipschitz and linearly bounded. There is a constant $K_{x}>0$ independent of $N,M$ and increasing in $\left\lVert x\right\rVert_{1-var}$ such that

where the expectation is taken over the joint distribution of $S_{0},\{A_{k},b_{k}\}_{k=1,\dots,d}$

thus by Lemma 3.2 of (Friz & Victoir, 2010) we obtain

Assume the activation function $\varphi$ is Lipschitz and linearly bounded.

Let $t\in[t_{m},t_{m+1})$ then, using the bound in Remark Remark,

Assume the activation function $\varphi$ is Lipschitz and linearly bounded.

Taking expectations of the squares and proceeding as before we have the thesis. ∎

We can finally prove the main bound which will allow the exchange of limits:

Note first that $S^{N}(x)$ is well defined since the system has unique solution by (Friz & Victoir, 2010)[Theorems 3.7, 3.8] with

noting that these are Lipschitz since composition of Lipschitz and linearly bounded.

Proposition (C.22) then gives the sought after bound.

where $\mu_{t}^{M,N}$ is the distribution of $\left\langle v^{N},{S}^{M,N}_{t}(x)\right\rangle$ and $\mu_{t}^{N}(x)$ that of $\left\langle v^{N},{S}^{N}_{t}(x)\right\rangle$ for some independent vector with iid entries $[v^{N}]_{\alpha}\sim\mathcal{N}(0,\frac{1}{N})$.

Note that the exact arguments, being of $L^{2}$ type, can be repeated for the ”stacked” vector $(S^{N,M}(x_{1}),\dots,S^{N,M}(x_{N}))$ extending (qualitatively) the bounds to the multi-input case.

By the classical Moore-Osgood theorem we need to prove that one of the two limits is uniform in the other, for example that the limit in distribution as $M\to\infty$ is uniform in $N$ in some metric which describes convergence in distribution. But this is just the content of the previous result, extended to the multi-input case. ∎

C.4 An alternative proof for the case φ=i​d𝜑𝑖𝑑\varphi=id

In this last section we prove Theorem C.18. between Randomized Signature Kernels and the original Signature Kernel.

Our goal is that of proving the following result:

and the variance around the limit is of order $O(\frac{1}{N})$.

We know, see (Baudoin & Zhang, 2012)[Remark 2.10], that it is possible to write a closed form for $S^{N}_{t}(x)$ which decouples the effects of the vector fields and those of the driving control using the Signature:

where, with $I$ as above, $V_{I}f(x):=V_{i_{1}}(V_{i_{2}}\cdots(V_{i_{k}}f)\cdots)(x)$ with

where $y$ is another control. If we could exchange expectation with the series we would thus get

Note that if $I=()$, the empty word, then $V_{I}f(S^{N}_{0})=S^{N}_{0}$.

The most important result which will make our plan succeed is the following classical theorem:

Let $(X_{1},\dots,X_{N})$ be a zero mean multivariate normal vector, then

where the sum is over all distinct ways of partitioning $\{1,\ldots,N\}$ into pairs $\{i,j\}$, and the product is over the pairs contained in $p$.

Let us first consider the case $\sigma_{b}=0$.

Moreover, for $I=(i_{1},\dots,i_{k})$, we obtain

where $\Lambda^{N,k}_{\alpha,\beta}:=\{(\delta_{0},\dots,\delta_{k})\in\{1,\dots,N\}^{k+1}:\delta_{k}=\alpha\text{ and }\delta_{0}=\beta\}$. With this notation we can write

The time is ripe for the application of Isserlis’ Theorem.

First of all notice how the sum in Isserlis runs over the possible pairings of the index set which in our case is the set of elements of the concatenation

In particular this sum is by default if $I*J$ has an odd number of elements. This means that

Moreover, even if $|I|+|J|$ is even, the pairings must be in such a way that no factor of the product vanishes; since we are working with matrices with independent normal entries this is equivalent to requiring for each $\{a,b\}\in p\in P^{2}_{|I|+|J|}$ that $(i*j)_{a}=(i*j)_{b}$, $(\delta*\epsilon)_{a}=(\delta*\epsilon)_{b}$, and $(\delta*\epsilon)^{\prime}_{a}=(\delta*\epsilon)^{\prime}_{b}$.

in fact considering only the constraints given by $\alpha$ and $\beta$ we have $N^{|I|}$ ways to choose $\bar{\delta}\in\{1,\dots,N\}^{|I|+1}$ and $N^{|J|}$ ways to choose $\bar{\epsilon}\in\{1,\dots,N\}^{|J|+1}$, but since they must come in pairs as dictated by $p$ we actually have $N^{\frac{|I|+|J|}{2}}$ possible choices i.e. $N$ per pair.

Notice however how we have equality if and only if $I=J$, $\alpha=\beta$ and the pairings are such that $p\ni\{a,b\}=\{a,|I|+a\}$ for $a\in\{1,\dots,|I|\}$ i.e. every element in $I$ is paired to the corresponding one in $J=I$.

The if part is easy to see: the full constraints are just $\delta_{|I|}=\epsilon_{|I|}=\gamma$, $\delta_{a}=\epsilon_{a}$ for every $1<a\leq|I|$ and $\alpha=\delta_{0}=\epsilon_{0}=\beta$ thus there are

The only if is more complicated and follows from the constraint $\delta_{|I|}=\epsilon_{|J|}$. Fix $\gamma$ and assume $i_{|I|}$ is not paired with $j_{|J|}$, then the choices for $(\delta*\epsilon)_{a}=(\delta*\epsilon)_{b}$ for 2 out of the $\frac{|I|+|J|}{2}$ pairs are constrained to be $\gamma$, all in all we have $N$ choices for $\gamma$ and at most $N^{\frac{|I|+|J|}{2}-2}$ for the other entries, thus at most $N^{\frac{|I|+|J|}{2}-1}$ in total. But then we must require $i_{|I|}$ and $j_{|J|}$ to be paired. The same argument can now be repeated with $i_{|I|-1}$ and $j_{|I|-1}$ since we have established that $\{|I|,|I|+|J|\}\in p$, thus not only $\delta_{|I|}=(\delta*\epsilon)_{|I|}=(\delta*\epsilon)_{|I|+|J|}=\epsilon_{|J|}$ but also $\delta_{|I|-1}=(\delta*\epsilon)^{\prime}_{|I|}=(\delta*\epsilon)^{\prime}_{|I|+|J|}=\epsilon_{|J|-1}$. This goes on until, without loss of generality, we run out of elements in $I$. If the same happens for $J$ (i.e. $|I|=|J|$) we are done since we have proved that $\forall a\in 1,\dots,|I|$ we have $i_{a}=j_{a}$ and $\alpha=(\delta*\epsilon)^{\prime}_{1}=(\delta*\epsilon)^{\prime}_{|I|+1}=\beta$. Otherwise $|J|\geq 2+|I|$ and $(j_{|J|-|I|},\dots,j_{1})$ are paired between themselves. But since $\{1,|I|+(|J|-|I|+1)\}\in p$ we have $\alpha=(\delta*\epsilon)^{\prime}_{1}=(\delta*\epsilon)^{\prime}_{|I|+(|J|-|I|+1)}=(\delta*\epsilon)_{|I|+(|J|-|I|+1)-1}=(\delta*\epsilon)_{|I|+(|J|-|I|)}=\epsilon_{|J|-|I|}$. Which means that there is no free choice one of the remaining couples (i.e. the one containing $|I|+(|J|-|I|)$ corresponding to $j_{|J|-|I|}$) thus, reasoning just as before, we cut the number of choices of at least a factor $N$.

We have just shown that, given a pairing $p$,

except when $I=J$, $\alpha=\beta$ and the pairings are such that $p\ni\{a,b\}=\{a,|I|+a\}$ for $a\in\{1,\dots,|I|\}$, in which case

where $0<\psi(I,J,N,\alpha,\beta)\leq(|I|+|J|)!!$ i.e. $\psi(I,J,N)$ is a positive constant bounded above by the maximal number of pairings (which occur only when $I$ and $J$ are made up of the same one index). Notice how this bound depends only on $|I|$ and $|J|$ and not on $N$!

Finally, if $S^{N}_{0}$ is normally distributed then

with $0\leq\psi(I,J,N):=\frac{1}{N}\sum_{n=1}^{N}\psi(I,J,N,n,n)\leq(|I|+|J|)!!$.

Let us now look at the case with $\sigma_{b}>0$.

now we have to study, with $\hat{I}:=(i_{2},\dots,i_{|I|})$, the terms

where in the last equality we have used the independence of the terms and their 0 mean.

Concerning the second term: using independence we readily see how we must have $i_{1}=j_{1}$, then

C.4.2 Convergence to the signature kernel

Assume now that the exchange of series with limits and expectation are justified, which we will prove later, then we would like to study the variance of the expected signature kernels around their limits.

The coefficients in the expansion of the variance

are all $O(\frac{1}{N})$ when $\sigma_{b}=0$.

If $S_{0}^{N}$ is sampled from a Normal distribution as before then

is equal to $\sigma_{S_{0}}^{4}$ if $|\{\alpha,\beta,\gamma,\delta\}|=2$, to $3\sigma_{S_{0}}^{4}$ if $|\{\alpha,\beta,\gamma,\delta\}|=1$ and to otherwise.

Let us then define $T=I*J*K*L$, $\theta=(\iota*\zeta*\kappa*\lambda)$ and $\theta^{\prime}=(\iota*\zeta*\kappa*\lambda)^{\prime}$ just like before, setting $\mathcal{P}:=P^{2}_{|I|+|J|+|K|+|L|}$ we get

and once again we need to analyze these $\omega$s which, just as before, must satisfy the constraint

Since we are interested in the behavior for $N\to\infty$ and $|\mathcal{P}|$ is independent from $N$ we just need to discover when the previous inequality is an equality. This happens, as we have previously discovered, when we the pairing does not add to the possible choices of $\bar{\iota},\bar{\zeta},\bar{\kappa},\bar{\lambda}$ any more constraints than the unavoidable ones i.e. $\theta_{|I|}=\theta_{|I|+|J|}$, $\theta_{|I|+|J|+|K|}=\theta_{|I|+|J|+|K|+|L|}$, $\theta^{\prime}_{1}=\alpha$, $\theta^{\prime}_{|I|+1}=\beta$, $\theta^{\prime}_{|I|+|J|+1}=\gamma$ and $\theta^{\prime}_{|I|+|J|+|K|+1}=\delta$.

Reasoning exactly as before this can happen if and only if $I=J$, $\alpha=\beta$, $K=L$, $\gamma=\delta$, $\theta_{a}=\theta_{|I|+a}$ for $a=1,\dots,|I|$ and $\theta_{2|I|+a}=\theta_{2|I|+|K|+a}$ for $a=1,\dots,|K|$.

where $\sigma(I,J,K,L)$ is a positive constant corresponding to the maximal number of non zero pairings.

The coefficients in the expansion of the variance

are all $O(\frac{1}{N})$ with $\sigma_{b}>0$.

To ease the notation let us write $i,j,k$ and $l$ instead of, respectively, $I_{1},J_{1},K_{1}$ and $L_{1}$.

where once again we use the fact that all the $O(\frac{1}{N})$ are uniformly bounded above by some $O(\frac{1}{N})$.

Using the same arguments developed up to here all other terms in the product end up as being $N^{2}O(\frac{1}{N})$ thus dividing finally by $N^{2}$ we have the thesis. ∎

We will just do the case $(\sigma_{S_{0}},\sigma_{A},\sigma_{B})=(1,1,0)$. The arguments for the general case are the same.

First of all, to be thoroughly rigorous, we need to define a probability space over which we take all the expectations, everything is numerable thus there is no issues with this.

To justify the exchange of sum and limit we want to use Lebesgue dominated convergence, we thus need to bound the

We know, from previous considerations, that

thus, since $\omega(p,I,J,N,\alpha,\beta)\leq N^{\frac{|I|+|J|}{2}}$, we obtain

Putting everything together we have found

Finally remember how, by factorial decay,

We have now to prove, by (Tao, 2016)[8.2.1 and 8.2.2], that

As a first step assume $2|i+j$ and, writing $i\wedge j:=\min\{i,j\}$ and $i\vee j:=\max\{i,j\}$, note how

Consider then $\phi$ as the inverse of the map $(i,j)\mapsto\frac{1}{2}(i+j)(i+j+1)+j$ i.e. $\phi$ is the map enumerating pairs $(i,j)$ starting from $(0,0)$ and proceeding with diagonal motions of the form

and then $(0,m)\to(m+1,0)$. Notice that such diagonal strides have length $m+1$ and comprise all couples $(i,j)$ such that $i+j=m$.

Since there are exactly $2k+1$ choices of $(i,j)$ such that $\frac{i+j}{2}=k$, corresponding to the couples $(i,2k-i)$ for $i=0,\dots,2k$, using the above $\phi$ it suffices to prove

The exchange of sums and expectations has always been justified.

For the case with just two indices $I,J$ we need to use Fubini-Tonelli and prove that

which is proved to be $<\infty$ exactly as in the previous proof, this time taking care to consider also the case $|I|+|J|$ not even.

The case with 4 words, i.e. the variance case, goes similarly. ∎

Putting all of this together we have finally proved the theorem:

Consider randomized Signatures of the type

and the variance around the limit is of order $O(\frac{1}{N})$.

C.4.3 Convergence to a Gaussian process

$\Phi^{N}_{\mathcal{X}}$ converge in distribution to a $\mathcal{N}(0,\mathcal{K}_{id}^{\mathcal{X}})$.

By Lévy’s continuity theorem it suffices to study the limiting behavior of the characteristic functions

But if we know $S_{0}$ and the $A_{k},b_{k}$ the randomized signatures are deterministic objects, and $\phi$ is normally distributed; thus

Fortunately we are only interested in evaluating $f_{u}$ on the $G^{N}_{\mathcal{X}}$ which are all positive semidefinite matrices: