A path-dependent PDE solver based on signature kernels

Introduction

Apart from the realm of rough volatility, stochastic Volterra processes have appeared in the modeling of electricity prices and turbulence , in principal-agent problems , climate modeling and in the study of certain asymptotic regimes . More broadly, fractional Brownian motions and related processes are present in a variety of situations where time-correlation is observed.

While the PPDE representation finds applications via studying its properties (e.g. weak error rates in ), this article is concerned with their numerical computation. The solution to a PPDE is a function of a path (and other finite dimensional inputs). Contrary to the wide literature on numerical schemes for high-dimensional PDEs, we are here interested in purely infinite dimensional inputs that we will treat as such instead of discretising them. Despite the central place of such PPDEs in mathematical finance, numerical methods pertaining to them are scarse, see Section 1.6 for more details. Our motivation is then twofold: design a numerical solver for such PDEs and overcome the technical challenges of learning a function on paths.

2. Signature kernels

Kernel methods are a well-established class of algorithms that are the key components of several machine learning models such as support vector machines (SVMs). The central idea of kernel methods is to lift the (possibly unstructured) input data to a (possibly infinite dimensional) Hilbert space by means of a nonlinear feature map defined in terms of a kernel. The advantage of doing so lies in the fact that a non-linear optimisation tasks in the original space becomes linear once lifted to the feature space. Thus, the solution of the original optimisation problem can often be expressed entirely in terms of evaluations of the kernel on training points by means of celebrated representer theorems. The kernel can often be evaluated efficiently with no reference to the feature map, a property known as kernel trick. The selection of an effective kernel will usually be a task-dependent problem, and this challenge is particularly difficult when the data is sequential. Signature kernels are a class of universal kernels on sequential data which have received attention in recent years thanks to their efficiency in handling path-dependent problems .

3. Kernels vs neural networks for PDEs

In recent years, PDE solvers based on neural networks and kernel methods have been introduced as alternatives to classical finite difference and finite elements schemes. The most famous examples of neural PDE solvers are the Deep Galerkin model (DGM) and Physics informed neural networks (PINNs) . These techniques essentially parameterise the solution of the PDE as a neural network which is then trained using the PDE and its boundary conditions as loss function. More recently, DGM and PINNs have been extended to path-dependent PDEs by . PDE solvers based on kernel methods have been investigated in various works , but the development of kernel methods for path-dependent PDEs has remained open. PDE solvers based on kernel methods hold potential for considerable advantages over neural PDE solvers, both in terms of theoretical analysis and numerical implementation. In effect, the theoretical analysis of neural PDE solvers is often limited to density/universal approximation results aimed at showing the existence of a network of a requisite size achieving a certain error rate, without guaranteeing whether this network is computable in practice via stochastic gradient descent. On the contrary, as we will demonstrate in this work, kernel-based PDE solvers are provably convergent and amenable to rigorous numerical analysis based on a rich arsenal of functional analytic results from the kernel literature.

4. Monte Carlo methods for pricing under rough volatility

The notion of rough volatility has recently disrupted stochastic modelling in quantitative finance , leaving no one indifferent. In this setting, the instantaneous volatility of assets is modelled as a fractional Brownian motion with small Hurst parameter as in (1.1) or as a stochastic Volterra process. This feature has been shown empirically to be consistent with historical market trajectories and to capture important stylised facts of equity options in a parsimonous way. However, this new paradigm comes with an increased computational cost compared to classical techniques. With the notable exception of the rough Heston model, the absence of Markovianity of volatility prevents any pricing tools other than Monte Carlo simulations. However, the new promises–for estimation and calibration–of these rough volatility models have encouraged deep and fast innovations in numerical methods for pricing, in particular the now standard Hybrid scheme . The complexity of a Monte Carlo scheme scales roughly as $\mathcal{O}(h^{-2}Nd^{2})$, where $N$ is the number of sample paths, $h$ is the discretisation step and $d$ is the dimension of the process. This complexity can be further reduced to $\mathcal{O}((h\log(h))^{-1}Nd^{2})$ with hybrid schemes . In this paper, we propose a simulation-free alternative to Monte Carlo pricing constructing a kernel-based solver for the corresponding PPDE. The complexity of the resulting solver scales as $\mathcal{O}(h^{-2}\hat{N}d^{2})$ where this time $\hat{N}$ is the number of collocation paths. This cost can be further reduced to $\mathcal{O}(h^{-1}\hat{N}d^{2})$ if the operations are carried out on a GPU, see for details. As our numerical experiments demonstrate, a transparent comparison of complexities between our approach and Monte Carlo for pricing under rough volatility can only be achieved by obtaining error rates for the convergence of the approximations, which we intend to explore as future work.

5. Contributions

In this paper we design a novel numerical solver for PPDEs based on signature kernel and accompanied by convergence guarantees. Our numerical scheme does not use any classical probabilistic representation or any other pre-existing solver; it constitutes a valid alternative to the ubiquitous Monte Carlo methods, in particular in the context of option pricing under rough volatility. The hurdles in the theoretical analysis are many-fold and their resolution contributes to the development of the theory of signature kernels. Challenges include the translation of the PPDE framework into one amenable for the use of signature kernels, the proof of existence and uniqueness of a global minimum in the optimisation problem as well as an analytic form of the solution by means of a representer theorem (Theorem 4.2 and Corollary 4.4), and a consistency result based on the adaptation of convergence results from to this infinite-dimensional setting via a compactness argument for functions on path spaces (Theorem 4.6 and Lemma 4.8). Finally, several numerical examples are presented in Section 5 and the code is open-sourced at https://github.com/crispitagorico/sigppde.

6. Related literature

Several deep learning and signature-based approaches have been proposed in recent years. In , the authors use a probabilistic approach based on the martingale representation theorem to design an objective function for solving time-discretised PPDEs via a parametric model combining RNNs and signatures. The martingale representation theorem is also central in , where the authors make use of a continuous-time RNNs developed in . In both these papers, the signature has to be truncated for the implementation. The authors of extend the deep Galerkin framework proposed by using an LSTM architecture at each time step. The only paper applicable to Volterra processes develops a backward discretisation algorithm using neural networks, drawing inspiration from the deep learning methodology introduced in . We note that none of the aforementioned papers present theoretical results on the proposed numerical schemes and they all rely on the training of neural networks to approximate the solution. Aloof from the learning paradigm, convergence guarantees of monotone schemes can be found in for the viscosity solutions to first generation PPDEs. The present paper is, to our knowledge, the first to propose a numerical scheme for second generation PPDEs with convergence guarantees.

7. Outline

Section 2 describes the rigorous framework of PPDEs and how we intent to approximate them. This is followed in Section 3 by an introduction to signature kernels and the definitions we will need in the remainder of the paper. The main theoretical analysis and results are presented in Section 4 while the numerical experiments can be found in Section 5. We give an outlook on future research in Section 6. Finally, the Appendix A gathers some technical proofs.

Path-dependent PDEs

As the example from the introduction suggests, in situations where the underlying process is fractional, such as fractional Brownian motion, rough volatility models or stochastic Volterra processes, conditional expectations are functions of paths. In this section, we will present several instances where these functions are solutions to well-posed path-dependent PDEs. Beforehand, we must set the mathematical scene in a rigorous way. We introduce the notations and definitions needed to formalise equations such as (1.2). This framework enables to present central examples in subsection 2.3, then in subsection 2.4 we adapt this setup to the signature realm and finally outline the numerical method is subsection 2.5.

The solutions to the path-dependent PDEs of interest are functionals of time, space and paths, hence we define suitable space and distance, borrowed from :

2. Functions on continuous paths

and similarly for the higher derivatives. The generalisation of the space $\mathcal{C}^{\mathcal{N}}$ to this type of derivative is more involved and requires regularity estimates tailored to the speed of explosion of the kernel. For precise definitions, we will point towards the appropriate references in the examples below. We are now ready to state the equations we wish to solve.

3. The target PDEs

More precisely, this paper will focus on three examples of parabolic second order PDEs which arise from the Feynman-Kac formula. The solution to these PPDEs thus have a probabilistic representation, a crucial advantage to assess the accuracy of our method. However we emphasise that the numerical scheme does not use this a priori knowledge as in supervised learning tasks. Much like deep Galerkin or PINN algorithms , our method can be deployed to solve PPDEs for which no numerical method is known.

In the three subsequent subsections, we present (i) the model (ii) the value function $u$ and (iii) the PPDE, in that order. We point to the suitable references for well-posedness and for more details. We emphasise that all three examples are particular cases of [10, Proposition 2.14]; this result ensures, under general assumptions, existence and uniqueness of a classical solution to the PPDE associated to a Volterra process.

In this first toy example we consider test functions of $(\widehat{W}_{s})_{s\in[t,T]}$ which corresponds to the state space $\widetilde{\Omega}$ with $d=0$ and $e=1$. Let us define the function

Based on [52, Theorem 4.1], under some standard regularity assumptions on $\phi,\psi$ and $K$, $u$ is a solution to the path-dependent heat equation

It can also be proved that this solution is unique in the space where the functional Itô formula [52, Theorem] is applicable.

3.2. VIX options under rough Bergomi model

3.3. Options under rough Bergomi model

In this last example, which was already presented in the introduction, the underlying is the (log-)asset price in a rough Bergomi-type model. The value function in the rough Bergomi model is indexed on the state space $\widetilde{\Omega}$ with $d=1$. As we have already seen,

for all $(t,x,\gamma)\in\widetilde{\Omega}$ and with boundary condition $u(T,x,\gamma)=\phi(x)$. Notice the boundary condition only applies to $x$ and is not path-dependent.

In the rough Bergomi model, the paths $\Theta$ have a financial interpretation as they are linked to the forward variance curves $\xi$, a quantity (more or less) observed on the markets. We refer to [40, Remark 4.7 (5)] for more details.

3.4. An important remark on the singularity of the kernel

We emphasise that the PPDEs presented in these examples fit in (2.3) albeit with a non-continuous direction $K^{t}$. The pathwise derivative is thus defined as the limit of the regularised derivative with the approximation $K^{t,\delta}$ as in (2.2). One can show thanks to the probabilistic representations that the solution to the regularised PPDE converges to the solution of the singular one as $\delta$ goes to zero. Our numerical scheme approximates the PPDE (2.3) with a regularised kernel, for a fixed $\delta>0$.

4. Functions on bounded variation paths

Making the framework of subsection 2.2 amenable to the application of signature kernels requires some carefulness. Firstly, the latter are defined for bounded variation paths, a subset of the continuous paths considered in the previous subsection. In this subsection, functions are thus defined over $\Omega$ instead of $\widetilde{\Omega}$.

where the sum is taken over multi-indices $n\in\mathcal{N}$.

5. A kernel method for PPDEs

We saw that we cannot solve a PPDE of the type $\widetilde{\mathcal{L}}_{t}u=f$ on $\widetilde{\Omega}$ directly with signature kernel methods, hence why we defined an analogous weaker PPDE $\mathcal{L}_{t}u=f$ on $\Omega\subset\widetilde{\Omega}$ which solution can be suitably represented with signature kernels. Functions admitting Fréchet derivatives also admit Gateaux derivatives and they are equal, therefore solutions to the former PPDE are also solutions to the latter. This will allow us in Theorem 4.6 to bridge the gap between the equation we solve and the equation we target. Now that we have this picture in mind, we can present the gist of our method, inspired by and extended to the path-dependent case.

In other words, we search for the function with minimal RKHS norm that satisfies the PDE constraints at all the collocation points. This problem formulation raises a number of questions, of both theoretical and practical importance:

Does the domain of $\mathcal{L}_{t}$ contain $\mathcal{H}$ ?

Does the minimisation problem (2.13) have a unique solution $u_{m,n}$?

How does one solve this optimisation over an infinite-dimensional space ?

Does the PPDE (2.3) have a unique classical solution $u^{\star}$ ?

Does $u_{m,n}$ converge to $u^{\star}$ as $m,n\to\infty$ ?

After designing our kernel and RKHS in Section 3, we answer positively to all these questions under appropriate conditions. In order to successively execute limit arguments, we exploit that $\mathfrak{K}$ is a compact set of $\Omega$, which we define explicitly under the topology induced by $\bm{d}$. Exploiting the robust signatures introduced in may allow to lift this assumption; we will investigate this idea in the future. The interested reader may find a summary of the answers as follows.

Yes, Proposition 4.1 proves this assertion.

Yes, as this is a type of optimal recovery problem, see Theorem 4.2.

As detailed in Corollary 4.4, the optimal recovery problem reduces to a finite dimensional quadratic optimisation problem with linear constraints. The latter can be solved via gradient descent, aloof from sophisticated learning algorithms.

As pointed out in Section 2.3, existence and uniqueness of this PPDE holds under some conditions, which are thoroughly checked in the given references for the examples of interest.

Theorem 4.6 ensures convergence towards the true solution under two different sets of assumptions, both requiring that $u^{\star}$ is the unique solution of the PPDE (2.3) and that the collocation points form a countable dense set of $\mathcal{H}$.

Let us summarise the several levels of approximation. We numerically solve the optimisation problem (2.13), which is a discretisation of the PPDE (2.12) over $\Omega$. Under the right set of assumptions, the solution to the optimisation problem actually approximates the solution to the PPDE (2.3) over $\widetilde{\Omega}$. Finally, the latter is an approximation of the PPDE with a singular direction $\bm{\eta}$. While our setup is designed to match PPDEs of second generation (arising from Volterra processes), it could be adapted to the first generation of Dupire (corresponding to path-dependent payoffs).

Signature kernels

Recall that a Hilbert space $\mathcal{H}$ of functions defined on a set $\Omega$ is a reproducing kernel Hilbert space (RKHS) over $\Omega$ if, for each $\omega\in\Omega$, the point evaluation functional at $\omega$, $f\mapsto f(\omega)$, is a continuous linear functional, i.e. there exists a constant $C_{\omega}\geq 0$ such that

The classical Moore-Aronszajn Theorem ensures that if $\kappa$ is a positive semidefinite kernel, then there exists a unique RKHS $\mathcal{H}$ such that $\kappa$ has the following reproducing property

2. Signature kernels

where $a=(a_{1},...,a_{n})$ and $b=(b_{1},...,b_{n})$ are tensors in $V^{\otimes n}$. Let $T((V))=\prod_{n=1}^{\infty}V^{\otimes n}$ the vector space of formal tensor series and $T(V)\subset T((V))$ the space of formal tensor polynomials defined as

There are many ways of extending by linearity the inner product in (3.3) to an inner product on $T(V)$. The simplest way to achieve this is by defining the inner product as follows

for any $a=(a_{0},a_{1},...),b=(b_{0},b_{1}...)$ in $T(V)$.

We denote by $\overline{T(V)}$ be the Hilbert space obtained by completing $T(V)$ with respect to $\langle\cdot,\cdot\rangle_{T((V))}$.

The principal ingredient needed to define the signature kernel is a classical transform in stochastic analysis known as the signature. Next we recall its definition.

When it is clear from the context, we will suppress the dependence on the interval and instead denote the signature of $\gamma$ by $S(\gamma)$. We refer the interested reader to for an account and examples on the use of signature methods in machine learning.

The signature of a BV path $\gamma$ can be equivalently defined as the solution of a $\overline{T(V)}$-valued differential equations controlled by $\gamma$, as stated in the following lemma.

admits $Y_{s}=A\cdot S(\gamma)_{[t,s]}$ as unique solution, where the product $\cdot$ on $\overline{T(V)}$ is defined for any $v=\left(v_{0},v_{1},...\right)$ and $w=\left(w_{0},w_{1},...\right)$ in $\overline{T(V)}$ as the element $v\cdot w=\left(z_{0},z_{1},...\right)$ in $\overline{T(V)}$ such that

and where $d\gamma_{t}$ is identified with $(0,d\gamma_{t},0,...)\in\overline{T(V)}$. In particular, the signature $S(\gamma)_{t,t^{\prime}}$ is the unique solution to (3.7) with initial condition $A=\mathbf{1}:=(1,0,0,...)\in\overline{T(V)}$.

By construction, the iterated integrals in (3.6) admits the following recursive structure

Thus the coordinates of the signature when taken in isolation, do not satisfy a differential equation.

We can now define the signature kernel as the inner product of two signatures.

It was proved in that signature kernel can be realised as the solution of path-dependent PDE as stated in the next lemma; this effectively provides a kernel trick for the signature kernel, i.e. a way of computing it without reference to the signature map.

and with boundary conditions $\kappa_{\text{sig}}(\gamma,\tau)_{t,q}=\kappa_{\text{sig}}(\gamma,\tau)_{r,s}=1$ for all $q\in[s,s^{\prime}],r\in[t,t^{\prime}]$.

Several finite difference schemes are available for numerically evaluating solutions to equation (3.8), see [44, Section 3.1] for details.

3. Signature kernels from static kernels

where the inner product is taken in the closure of $T(\mathcal{H}_{g})$.

Lemma 3.6 directly yields a kernel trick for the $g$-lifted signature kernel $\kappa^{g}_{\text{sig}}$:

with the boundary conditions $\kappa^{g}_{\text{sig}}(\gamma,\tau)_{t,q}=\kappa^{g}_{\text{sig}}(\gamma,\tau)_{r,s}=1$ for all $q\in[s,s^{\prime}],r\in[t,t^{\prime}]$.

The integration with respect to $\mathcal{H}_{g}$ paths is understood in the sense of

4. Derivatives of static kernels

To approximate solutions to PPDEs, it will be necessary to differentiate a kernel with respect to its input variables. The RBF kernel $g_{\sigma}$, defined in (3.2), is clearly $C^{\infty}$ and its $n^{th}$ order derivative can be obtained using elementary calculus. For example, the first two derivatives, which are the only derivatives we will need to consider in the experimental section, admits the following expressions

By the results in , for all $x\in V$, the derivatives of $g(x,\cdot)$ are elements of $\mathcal{H}_{g}$. Furthermore, [55, Theorem 1] states that, for any $h\in\mathcal{H}_{g}$ and any multi-index $I$,

We make the following standing assumption on the static kernel.

This assumption is satisfied for the Gaussian kernel (3.2).

5. Derivatives of signature kernels

It was shown in that the directional derivative in equation (3.14) is such that

The proofs of all the results of this section are postponed to Appendix A.2 to ease the flow of the paper.

Expliciting the partial derivatives in (3.16) yields the equivalent formulation:

Equation (3.16) is more conducive for numerical computations while Equation (A.2) will be used for the theoretical analysis. One can perform the same expansion for the second derivative, as can be seen in the proof (Equation (A.2)).

This result is analogous to [25, Theorem 4.4] albeit in $\mathcal{H}_{g}$ valued paths.

6. Computations of derivatives

satisfies the following linear system of hyperbolic PDEs

7. Universal product kernels

[9, Lemma 5.2] Let $\kappa_{1},\kappa_{2}$ be two cc-universal kernels indexed on $\Omega_{1},\Omega_{2}$ respectively. Then the product kernel defined for any $\omega_{1},\omega_{1}^{\prime}\in\Omega_{1}$ and $\omega_{2},\omega_{2}^{\prime}\in\Omega_{2}$ as

is cc-universal on $\Omega_{1}\times\Omega_{2}$.

However, signature kernels are only universal when restricted to a class of paths where the signature is an injective map, as we explain next.

Moreover, as a consequence of Lemmas 3.16 and 3.17, the product kernel $\kappa$ is cc-univeral on $\mathfrak{K}\times\mathfrak{K}$. In other words, $\mathcal{H}$ is dense in $\mathcal{C}(\mathfrak{K})$ with respect to the topology of uniform convergence.

Main theoretical results

This section aims at answering the questions (1)-(3) of Section 2.5 regarding the well-posedness of the optimisation problem (2.13). Our first result verifies that elements of $\mathcal{H}$ are indeed twice differentiable, and therefore belong to the domain of the PPDE operator $\mathcal{L}_{t}$. This answers Question (1) in 2.5. The proof builds on the estimates derived through blood and sweat in Lemma A.2. The technical hurdle consists in interchanging the infinite series and the derivative operator. As advertised earlier, this operation requires the restriction of $\mathcal{H}$ to a compact subset of $\Omega$ called $\mathfrak{K}$. Even though $\mathfrak{K}$ is closed, the time derivative of $h\in\mathcal{H}$ is well-defined since it only perturbs to the right and we set $\partial_{t}u(T,x,\gamma)=0$. Moreover the spatial and pathwise derivatives are also well-defined over $\mathfrak{K}$ since the kernel $\kappa$ is defined over $\Omega\times\Omega$.

In particular the domain of $\mathcal{L}_{t}$ is included in $\mathcal{H}$ for all $t\in[0,T)$.

This highlights the link between the two types of pathwise derivatives $\widehat{\partial}^{\eta}_{\gamma}$ and $\partial^{\eta}_{\gamma}$.

which is uniformly bounded over all $\omega\in\ \mathfrak{K}$ by Lemma A.1. We used dominated convergence in the last line which follows from (A.2). This entails the convergence of $\partial^{\eta}_{\gamma}h_{n}$ towards $h^{\prime}$ is uniform in $\mathfrak{K}$; therefore $h$ is differentiable and $\partial^{\eta}_{\gamma}h=h^{\prime}=\sum_{i=1}^{\infty}\alpha_{i}\partial^{\eta}_{\gamma}\kappa(\omega_{i},\cdot)$.

and we will prove that $\lim_{n\to\infty}\partial^{\eta{\bar{\eta}}}_{\gamma}h_{n}=h^{\prime\prime}$. First note that

Applying the same approach as for the first derivative, we obtain by Cauchy-Schwarz inequality

where we used dominated convergence as for the first derivative. Lemma A.1 thus ensures that $\partial^{\eta{\bar{\eta}}}_{\gamma}h_{n}$ converges towards $h^{\prime\prime}$ uniformly in $\mathfrak{K}$ and thus $h$ is twice differentiable with $\partial^{\eta{\bar{\eta}}}_{\gamma}h=h^{\prime\prime}=\sum_{i=1}^{\infty}\alpha_{i}\partial^{\eta{\bar{\eta}}}_{\gamma}\kappa(\omega_{i},\cdot)$. ∎

Now that we have checked that the constraints of the optimal recovery problem (2.13) make sense, we arrive at the following representer theorem, which gives an answer to Question (2) of our list 2.5.

associated to the RKHS $\mathcal{H}$ is well-posed and has a unique minimiser $u_{m,n}\in\mathcal{H}$ of the form

The proof is similar to the proof of classical representer theorems with the addition of pointwise observation of higher order derivatives. The existence and uniqueness of the minimiser of the optimal recovery problem (4.2)) is guaranteed by the fact that the Lagrangian functional

where the last equality follows from the fact that $u^{\perp}\in\mathcal{A}^{\perp}$. Hence, the following equality holds

Since for any $u\in\mathcal{H}$ the functional $u\mapsto\frac{1}{2}\left\lVert u\right\rVert_{\mathcal{H}}^{2}$ is strictly monotonically increasing, one has $\frac{1}{2}\left\lVert u\right\rVert_{\mathcal{H}}^{2}\geq\frac{1}{2}\left\lVert\hat{u}\right\rVert_{\mathcal{H}}^{2}$, which leads to the required inequality

Note that $(\alpha_{i})_{i=1}^{m+n}$ crucially depend on the collocation points and on $m$ and $n$; in particular they have no reason to remain stable for different values of $m,n$. This makes comparison among $(u_{m,n}$ much trickier. As another avenue for future research, we note that orthogonalising the collocation points in a way that the $(\alpha_{i})_{i}$ remain constant for all $m,n$ could lead to better control of the norm.

Note that because the functional $u_{m,n}$ defined in (4.3) is an element of the RKHS $\mathcal{H}$, by the reproducing property we have the following expression of the squared norm

This lifts the main computational obstacle and opens a clear path to solving the optimisation problem, thereby resolving Question (3) of 2.5.

The matrix $\widetilde{\mathcal{K}}$ is well-defined by the differentiability of the product kernel.

2. Consistency

We associate to it the range of derivatives $\mathcal{D}:=\{{\rm id},\partial_{t},\partial_{x},\partial^{2}_{x},\partial^{\eta}_{\gamma}\partial_{x},\partial_{\gamma}^{\eta},\partial^{\eta{\bar{\eta}}}_{\gamma}\}$ such that the norm also reads

Our main result, inspired by [12, Theorem 1.2], follows. Its proof is given at the end of the section.

There exists a unique solution $u^{\star}$ to the PPDE (2.3) in $\mathcal{H}$.

The assumptions of the theorem all deserve separate discussion, studied in inverse order.

On the collocation points. It is natural to ask the collocation points to be dense in $\mathfrak{K}$, however this leaves a lot of flexibility for choosing them in an optimal way, especially in the space of paths. This issue is directly related to the rate of convergence of the method which lies beyond the scope of this paper. In practice though, one is only interested in covering the support of the evaluation points. During the course of our experiments, we found that sampling the collocation points randomly (as Brownian motion trajectories) led to similar accuracy as sampling them from the same distribution as the evaluation points.

On the well-posedness of (2.3) Existence and uniqueness of (2.3) is discussed in several examples of interest in Section 2.3, and in the more general Volterra case in . Uniqueness allows us to show that all the converging subsequences of $(u_{m,n})_{m,n}$ have the same limit, thus implying the sequence does converge. Unfortunately, known uniqueness results hold on a subset of $\mathcal{C}(\widetilde{\Omega})$ whereas we require uniqueness in (a subset of) $\mathcal{C}(\mathfrak{K})$ where $\mathfrak{K}$ is strictly included in $\widetilde{\Omega}$.

On the assumption $u^{\star}\in\mathcal{H}$.

Proving convergence of $(u_{m,n})_{m,n}$ at the very least requires to find a norm in which this family is uniformly bounded. The reason for supposing $u^{\star}\in\mathcal{H}$ is that it implies both the RKHS and $\mathcal{C}^{\mathcal{N}}$ norms of $(u_{m,n})_{m,n}$ are uniformly bounded, which in turn yields the relative compactness of this family in $\mathcal{C}^{\mathcal{N}}$, by Lemma 4.8. It is still slightly different from condition (i) because uniqueness holds in a different space.

For $u^{\star}$ to be in $\mathcal{H}$ would require at minima to be in $\mathcal{C}^{\infty}$. In this direction, Theorem 3.10 in shows how to expand a path (semimartingale) functional with respect to terms of the signature. It is conceivable that such a representation also holds in the fractional setting for sufficiently smooth functionals, which could then be identified to an element of the signature RKHS. When the underlying process is a semimartingale, smoothness of the payoff function (and the coefficients) essentially ensures smoothness of the conditional expectation (i.e. the solution to the Kolmogorov equation). However when the direction of the pathwise derivative is singular (only square integrable) one can only prove twice differentiabiliy of the solution (see [10, Proposition 2.23]). This remains a glass ceiling until one unveils how to exploit the regularisation properties of the expectation.

Several other works in the literature assume that $u^{\star}\in\mathcal{H}$. make in addition the classical assumption that $\mathcal{H}$ is embedded in a Sobolev space, a condition we drop because we are able to prove sufficient regularity only with the help of the bounded RKHS norm. In finite dimensions, this condition allows to prove convergence rates, see [54, Proposition 11.30]. See also who provide error bounds under additional assumptions.

The proof strategy for Theorem 4.6 using condition (ii) consists in showing that $(u_{m,n})$ is relatively compact and then extracting converging subsequences. Hence before presenting it we need to identify compact subsets of $\mathcal{H}$ with respect to the $\mathcal{C}^{\mathcal{N}}$-topology. The RKHS norm is linked to the regularity of the function, therefore a family of $\mathcal{H}$ bounded in RKHS norm is equicontinuous. If this family is also bounded in $\mathcal{C}^{\mathcal{N}}$ norm we can conclude by Arzelà-Ascoli theorem that it is relatively compact.

Let Assumption 3.10 hold for $g$. Let $\mathcal{Y}\subset\mathcal{H}$ be bounded under both the RKHS and $\mathcal{C}^{\mathcal{N}}$ norms. Then $\mathcal{Y}$ is relatively compact with respect to the $\mathcal{C}^{\mathcal{N}}$ topology.

We take evaluation points over a compact and $\partial\kappa(\omega,\cdot)$ is continuous hence $\sup_{\omega\in\mathfrak{K}}\left\lVert\partial\kappa(\omega,\cdot)\right\rVert_{\mathcal{H}}^{2}<\infty$. Thus there exists $C>0$ such that for all such $\eta,{\bar{\eta}}$,

Under condition (i), the same conclusion holds except that $u_{\infty}$ is an element of the closure of $\mathcal{H}$ with respect to the $\mathcal{C}^{\mathcal{N}}$ norm; the universality property of the RKHS entails that $\overline{\mathcal{H}}=\mathcal{C}^{\mathcal{N}}(\mathfrak{K})$. We now need to show that $u_{\infty}$ solves the PPDE (2.3).

In addition, $v_{m,n}$ converges uniformly to $v$ because $u_{m,n}$ converges to $u_{\infty}$ in $\mathcal{C}^{\mathcal{N}}$. Since $\varepsilon>0$ was arbitrary, Equation (4.7) thus entails that $v(\omega)=f(\omega)$. Following a similar argument it can be shown that $u_{\infty}(T,x,\gamma)=\phi(T,x,\gamma)$ for any $x,\gamma$.

Numerical experiments

In this section, we present numerical experiments benchmarking the signature kernel PPDE solver presented in the previous section against either analytic solutions (if they are available) or classical Monte-Carlo solvers. We consider the two examples mentioned in the introduction, namely the path-dependent heat equation (2.7)) and the rough Bergomi PPDE (2.9). For both experiments, we consider two types of errors between the predicted prices and true prices, namely the mean squared error (MSE), and the mean absolute error (MAE). The optimal kernel hyperparameters were determined by cross-validation. We will conclude the section with a discussion to compare our kernel approach with recent neural networks techniques for solving PPDEs. Our code is available at https://github.com/crispitagorico/sigppde.

In this first example, we consider a one dimensional fractional Brownian motion $X$ with Hurst exponent $H=0.1$. We recall that by [52, Theorem 4.1], under appropriate regularity assumotions on assumptions on the functions $\phi,f$, the conditional expectation

is realised as the solution of the path-dependent heat equation

In our experiments we choose $f\equiv 0$ and consider the following three instances of function $\phi$, for which analytic expressions of the conditional expectation are available:

Because analytic prices are available, we limit ourselves to asses the performance of our kernel algorithm to recover these true prices. Collocation points for the kernel method were chosen uniformly on $$for the time variable and sampled from the process$\Theta^{t}$ for the path variable.

As it can be observed from Figure 1, our kernel approach is capable to recovering with a high level of accuracy the prices for all considered payoff profiles.

2. Rough Bergomi

In this section, we make use of the following formulation of the rough Bergomi model introduced by and already written down in the introduction

To evaluate the conditional expectation in the Monte Carlo benchmark we remark that, setting

if $s>t$ and otherwise, the conditional expectation reduces to the following unconditional expectation as we showed in the introduction

since $\Theta^{t}$ is an $\mathcal{F}_{t}$ measurable process and $I^{t}$ is independent from $\mathcal{F}_{t}$. Collocation points for the kernel method were sampled uniformly at random on $$for the time variable, uniformly at random on$[x_{min},x_{max}]$for the price variables, where$x_{min},x_{max}$were chosen using ground truth prices, and sampled from the process$\Theta^{t}$ for the path variable.

We consider two values of $H=0.1$ and $H=0.3$ and two strike values $K=0.1$ and $K=1.0$. We take the same number of sample paths and collocation points for Monte Carlo and the signature kernel solver respectively, and consider three different increasing values, namely $20$, $100$ and $500$. As it can be observed from the tables, the numerical performance in terms of MSE and MAE of the two algorithms are comparable across the different settings.

The time complexity of a Monte Carlo scheme scales roughly as $\mathcal{O}(h^{-2}Nd^{2})$, where $N$ is the number of sample paths, $h$ is the discretisation step and $d$ is the dimension of the process. This complexity can further reduced to $\mathcal{O}((h\log(h))^{-1}Nd^{2})$ using schemes such as the one proposed by . The complexity for evaluating our solver (once it has been trained offline) scales as $\mathcal{O}(h^{-2}\hat{N}d^{2})$ where $\hat{N}$ is the number of collocation points. This cost can be further reduced $\mathcal{O}(h^{-1}\hat{N}d^{2})$ if the operations are carried out on a GPU, see for additional details. It is clear that a transparent comparison of complexities between our approach and Monte Carlo for pricing under rough volatility can only be achieved by obtaining error rates for the convergence of the approximations governing the values of $N$ and $\hat{N}$ respectively in the above complexities, which we intend to explore as future work.

Nevertheless, the PDE approach provides by design advantages over Monte Carlo methods. Indeedn, finite-dimension greeks are prone to instability, let alone pathwise ones, while all derivatives of the PPDE solution appear from (the linear combination of) standard PDEs.

3. Comparison with neural PDE solvers

As antiticipated, we will conclude this section with some remarks aimed at comparing our kernel approach to the numerous neural network techniques to solve PPDE which have been proposed in the literature in recent years. Neural PDE solvers such as DGM and PINNs essentially parameterise the solution of the PDE as a neural network which is then trained using the PDE and its boundary conditions as loss function. One drawback of these approaches in the path-dependent setting is that they require a time-discretisation of the solution. The time grids over which solutions and derivatives are evaluated during training and testing need to agree, usually followed by an arbitrary interpolation. On the contrary, our kernel approach operates directly on continuous paths, and it is therefore mesh-free. Furthermore, derivatives can be accessed directly by differentiating kernels, with no need to prematurely discretise in time. We note that concomitant to this paper is , which uses a neural rough differential equation (Neural RDE) model to parameterise the solution of a PPDE and is therefore also mesh-free. Another drawback of neural PDE solvers is that their theoretical analysis is often limited to density/universal approximation results; showing the existence of a network of a requisite size achieving a certain error rate, without guarantees whether this network is computable in practice. In fact the optimisation, done with variants of stochastic gradient descent, is never convex so there is no guarantee to attain a global minimum. Our kernel approach boils down to a convex finite dimensional optimisation, guaranteed to attain the unique global minimum, and that is consistent with the original problem when the number of collocation points is sent to infinity.

Outlook

The numerical solver presented in this paper competes with Monte Carlo methods in terms of the balance between complexity and accuracy. Our next objective will be to derive error estimates for the proposed method with respect to the number of collocation points in order to carry out a fair comparison. This step is essential for choosing the collocation points in an optimal way and better assess the efficiency of the algorithm. For this purpose we may follow the approach of or decide to regularise our optimisation problem as in . Further extensions include PPDEs arising from general Volterra processes, non-linear PPDEs, and path-dependent payoffs as the VIX example presented in Subsection 2.3.2. Finally, the inverse problem is particularly relevant for financial applications. We expect that the calibration of model parameters can be considerably sped up if the kernel is enhanced with this additional data.

Appendix A Proofs of signature and signature kernel estimates

2) Furthermore, [25, Proposition 2.9] gives the $1$-variation bound

Thus, applying the inequality $\left\lVert\left\lVert\gamma^{g}-\bar{\gamma}^{g}\right\rVert_{\mathcal{H}_{g}}\right\rVert_{0}\leq\overline{g_{11}}\left\lVert\gamma-\bar{\gamma}\right\rVert_{0}$ and Grönwall’s lemma we obtain

3) We start with the following observation. For all $u\in[t,T]$ we have

Further, we note that, for all $\varepsilon>0$,

Hence, by dominated convergence and (A.8), we get

Note that $\left\lVert\partial_{x^{ij}}^{2}g(\gamma_{u},\cdot)\right\rVert_{\mathcal{H}_{g}}=\partial_{x^{ij}y^{ij}}^{4}g(\gamma_{u},\gamma_{u})$, and hence

4) Noting that, for all $t\leq s^{\prime}\leq t^{\prime}\leq T$,

Then, from the integration by parts formula we observe that

This allows us to study the following regularity:

5) The latest estimate allows to use dominated convergence one more time and compute the second derivative in the direction of $\eta,{\bar{\eta}}$:

Using the same arguments as in (A.8) one deduces that

Hence the difference of the source term for two paths $\gamma,\bar{\gamma}$ can be bounded as follows

A.2. Differentiability and continuity estimates for the signature kernel

Then for any $\varepsilon>0$, Cauchy-Schwarz inequality and (A.2) yield

Therefore, dominated convergence, Cauchy-Schwarz and (A.3) entail

The same ideas coupled with the bounds (A.5) and (A.6) lead to the last estimate. ∎

1) In virtue of (A.2), dominated convergence allows to push the limit inside the integrals:

For all $u\in[t,T]$ and $v\in[s,T]$ we have

The expressions (A.15) and (A.14) yield (3.16) and (A.2) respectively.

2) We take the terms of (A.2) one by one to compute

Thanks to Lemma A.1, we can use dominated convergence to swap limits and integrals. After rearranging terms we obtain

Eventually this yields the more compact formulation

References