Black-Scholes Theory and Diffusion Processes on the Cotangent Bundle of the Affine Group.

The Black-Scholes partial differential equation (PDE) from mathematical finance has been analysed extensively and it is well known that the equation can be reduced to a heat equation on Euclidean space by a logarithmic transformation of variables. However, an alternative interpretation is proposed in this paper by reframing the PDE as evolving on a Lie group. This equation can be transformed into a diffusion process and solved using mean and covariance propagation techniques developed previously in the context of solving Fokker-Planck equations on Lie groups. An extension of the Black-Scholes theory with coupled asset dynamics produces a diffusion equation on the affine group, which is not a unimodular group. In this paper, we show that the cotangent bundle of a Lie group endowed with a semidirect product group operation, constructed in this paper for the case of groups with trivial centers, is always unimodular and considering PDEs as diffusion processes on the unimodular cotangent bundle group allows a direct application of previously developed mean and covariance propagation techniques, thereby offering an alternative means of solution of the PDEs. Ultimately these results, provided here in the context of PDEs in mathematical finance may be applied to PDEs arising in a variety of different fields and inform new methods of solution.

1. Introduction

The Nobel Prize-winning Black-Scholes Equation [1] is arguably the most well-known partial differential equation in mathematical finance. The equation rests on a parsimonious option pricing model that has been fairly successful in informing banks and portfolio managers of the construction of risk-free hedges. Additionally, the model has since provided the framework for the import of a variety of tools, such as the theory of stochastic processes, from physics and mathematics. In this paper, we offer a new Lie group-theoretic interpretation of the Black-Scholes equation by reformulating the original equation and extensions of it as diffusion processes on Lie groups.

Group-theoretic approaches have been used extensively in the analysis of symmetries of partial differential equations (PDEs) in mathematical finance [2,3,4]. One of the central questions there is to identify the group of transformations of variables that can be applied to the equations that preserve the structure of the equations while reducing it to a simpler form for analysis and solution. For instance, the one-asset [2] and in general, the multi-asset Black Scholes equation [5] can be reduced to a heat equation through a logarithmic transformation of variables.

In this paper, instead of analysing the symmetry properties of PDEs, we reformulate the PDE as a diffusion equation on the Lie group. Functions of real parameters are upgraded to functions that take group elements as their arguments. That is, we start with a linear (parabolic) PDE of the form, where L is the linear partial differential operator and is a vector of coordinates that may be chosen based on how the differential equation is derived. In the case of the Black-Scholes equation and related asset models considered in this paper, L and are originally defined in an N-dimensional Euclidean space. By matching the derivatives in L to Lie derivatives of correctly chosen groups, one can rewrite (1) as the following linear PDE over the Lie group G, where parameterizes the group element and . The operator is now a differential operator consisting of Lie directional derivatives. Such differential equations arise in bevel-tip needle steering [6,7], error propagation [8,9], DNA statistical mechanics [10], multi-robot localization [11], stochastic kinematic cart models in SLAM [12] and in image contour completion and enhancement [13,14,15].

The benefit of reframing differential equations this way would be that variable-coefficient PDEs of N variables in may reduce to constant-coefficient PDEs on the Lie group G of dimension of at least N, which can be analysed using techniques developed in [6,7,8,9,10,11,12]. These applications motivated the development of mean and covariance propagation techniques to approximate the solution of diffusion equations in the special Euclidean group and more generally, in unimodular groups [16,17]. Other numerical techniques to solve differential equations on Lie groups, through a generalisation of Euler and Runge-Kutta schemes have been developed in [18,19,20].

In this paper, we extend the regime of applicability of mean and covariance propagation techniques by first considering the groups and that arise in the one-asset and two-asset Black-Scholes equations, respectively, and then the affine group, that arises in a coupled-asset dynamics extension of the Black-Scholes theory. The application of mean and covariance propagation to becomes the main subject of the paper. There exist other methods to approximately solve PDEs in mathematical finance, such as finite difference methods [21,22], finite element methods [23] and the Adomian decomposition method [24,25]. In this paper, we also make a comparison between the mean and covariance propagation technique and a standard finite difference scheme used to solve the governing equations.

The Lie group is not unimodular, thereby precluding the direct application of the theory of mean and covariance propagation used in the solution of diffusion equations on the group. However, we show that the cotangent bundle group of a Lie group (with a trivial center) is unimodular. Thus, we analyse diffusion processes on the affine group by matching the diffusion on the affine group with a degenerate diffusion on the cotangent bundle to which mean and covariance propagation techniques can be applied.

This paper considers three types of asset dynamics models: one-asset, two-asset and a coupled asset model. Each of these models give rise to a PDE describing the evolution of the option price as a function of the asset prices, and are introduced in Section 2. These PDEs in Euclidean space are reframed as diffusion processes on Lie groups in Section 3. The theory of mean and covariance propagation in the solution of diffusion processes on Lie groups is reviewed in Section 4. The one-asset and two-asset Black-Scholes models, which can be trivially solved using the logarithmic transformation of variables, are used as examples to illustrate the techniques, and thereby sets the ground for the non-trivial case of a coupled asset model. Since the coupled asset model leads to a diffusion equation over a non-unimodular group, Section 5 proves that a related structure, the cotangent bundle group of a Lie group with trivial center, is always unimodular. Section 6 describes the mean and covariance propagation over the unimodular cotangent bundle of the affine group, thereby solving the option price dynamics for the coupled asset model. Section 7 provides some numerical results for the mean and covariance propagated solution in comparison with finite difference methods. Finally, Section 8 seeks to demonstrate backward compatibility and completes the analysis by applying mean and covariance propagation on the cotangent bundle group to the solution of the one-asset Black-Scholes equation.

2. Asset Dynamics Models

In this section, we review the derivation of the well-known one-asset Black-Scholes equation, the two-asset Black-Scholes equation that evolve with correlated Wiener process, and follow a similar derivation to develop a coupled asset model.

2.1. One-Asset Black-Scholes Equation

The one-asset Black-Scholes equation is derived from the following Itô stochastic differential equation [26] governing the dynamics of an asset value a, where is the increment of a Wiener process, corresponding to random draws from a Gaussian probability distribution [27,28]. The increment of the Wiener process satisfies the following relations: and in general for an N-dimensional Wiener increment we have, where in both cases denotes the expected value with respect to the underlying probability distribution function and is the Kronecker-delta ( if and only if , otherwise ). The and in Equation (3) respectively describe a drift and volatility for the asset price evolution. The volatility is defined such that is the variance of the random price fluctuation. Equation (3) is a stochastic differential equation whose solution is a geometric Brownian motion. We emphasise that the equation in (3) is to be interpreted as an Itô stochastic differential equation. Stratonovich stochastic differential equations also appear in this paper (see Equation (110) for instance) and will be distinguished from Itô equations with a Ⓢ, i.e., for a stochastically varying quantity x, a general one-dimensional Stratonovich stochastic differential equation with Wiener noise would be, for a drift and diffusion coefficient . Nevertheless the properties of in (4) and of in (5) hold for both Itô and Stratonovich equations.

The option price is a function of time t and the underlying asset price a. Using Itô’s Lemma [29], we have for correct to order ,

In the Black-Scholes model, there also exists another asset, the bond price, which evolves deterministically (for a non-stochastic interest rate) as where r is the risk-free interest rate [26]. A portfolio is constructed as a linear combination of the two assets and can be written as,

This form of is due to the fact that the portfolio is assumed to be self-financing with no external money flows [30]. Banks choose in order to obtain a zero-risk portfolio [31], i.e., so that would have no Wiener noise terms. Assets without stochasticity would evolve at the risk-free interest rate r and therefore we obtain, and using , a hedge of the form would cancel out the term in (6) and therefore make the dynamics of risk-free. Substituting this hedge into (8) and using the form of in (6) we obtain the one-asset Black-Scholes equation as,

The parameters and do not feature in this equation.

2.2. Two-Asset Black-Scholes Equation

In the two-asset model, one has asset values a and b that evolve with correlated Wiener processes—this is a specific case of the multi-asset model in [5]. That is, where and represents the expectation. In the two-asset problem, the option price is . Applying Itô’s Lemma to now gives:

The portfolio is now where and ; and are both chosen such that experiences risk-free dynamics—very similar to the procedure in the derivation of the one-asset Black-Scholes equation in the previous section. The two-asset Black-Scholes equation [32] will then be, where r is the risk-free interest rate. Like the one-asset Black-Scholes equation, and do not feature in the two-asset equation.

2.3. Option Price Evolution with Coupled Assets

The evolution equations of the two assets a and b in the two-asset Black-Scholes model (10) featured correlated Wiener processes. In this section, we imagine a different form of coupling of the form:

The evolution of the value of asset a is independent of the evolution of b; a is therefore the ‘leading’ asset value. On the other hand, the drift and variance in the stochastic differential equation for are dependent on the current values of the leading asset. Hence b is the ‘trailing’ asset value. One may imagine this form of dependency when a represents a raw material and b represents a finished good that makes use of this raw material.

Both and are forced by the same Wiener process ; this implies that for the risk-free interest rate r. This can be derived by using the risk-neutral measure with the knowledge that all assets evolve with a risk-free interest rate of r in this measure [31]. The option price is and applying Itô’s Lemma to now gives,

We now assume a general hedge of the form . Since there is only one Wiener process , removing uncertainty would not be a sufficient condition to determine both and . Thus, for simplicity, we let which leads to the relationship . The governing partial differential equation for with two coupled assets would be, where the leading asset/trailing asset coupling between the variables breaks the symmetry between the assets a and b that existed in (12).

3. Reframing Partial Differential Equations as Diffusion Processes on Lie Groups

The governing PDEs for the three types of asset dynamics models in (9), (12) and (15) can be re-expressed in terms of Lie derivatives of , and , respectively, by matching the group parameters to the asset variables. A review of the concept of Lie derivatives is presented in Section 3.1.

3.1. Preliminary Definitions

Let G be an N-dimensional matrix Lie group with Lie algebra . Then, let an element be parameterized as where , using the notation in [29,33]. The ‘right’ Jacobian of the group is defined [34] as the following matrix, and the ‘left’ Jacobian is the matrix, where the square brackets reinforce that we are dealing with a matrix. This is not to be confused with the ‘right’ and ‘left’ Jacobian determinants that arise in the volume forms of the group, which would be the determinants of the matrices in (16) and (17), respectively. Note that the ‘right’ Jacobian has the appearing on the right whereas the ‘left’ Jacobian has the term on the left. However, the ‘right’ Jacobian is left invariant, i.e., and the ‘left’ Jacobian is right invariant, i.e., , assuming that parameterizes the whole group G and that these shifts are permitted in the function domain. Finally, the ∨ operator is defined as a bijection mapping to , and vectorizes the matrix element in . The inverse of the ∨ operator is a ∧ that maps to .

We make an additional remark regarding the parameterization of the group with . This is to say that the whole group (except for a set of measure zero) is parameterized by one coordinate chart. For instance, for groups such as or where one may use Euler angles to parameterize the rotations, there exists a set of measure zero corresponding to the set of Euler angles where the Jacobian matrices for the parameterization becomes singular. For denoting the ZXZ or ZYZ Euler angles, singularities occur at and . Additionally, since a rotation by describes the same rotation as that by , the open coordinate chart will be , which has a one-to-one correspondence with a subset of rotations that excludes the rotations at and , and the rotations at and . The closure of this coordinate chart will however establish a many-to-one map with the group and parameterize all group elements. A similar issue exists in the case of using the Iwasawa decomposition to parameterize with where parameterizes the 2D rotation, and parameterize the upper triangular matrix. Here, the coordinate chart would be ; yet again, the closure of this chart will establish a many-to-one map with and parameterize all group elements. For other cases, such as using a vector drawn from to parameterize as well as for the affine group and the cotangent bundle group of considered in this paper, the coordinate chart used parameterizes the entire group.

The operator of the adjoint representation of group G at is given by where for any . When expressed as a matrix, we denote the operator as . If for forms an N-dimensional orthonormal basis for , is given as,

Orthonormality of the Lie algebra basis is defined here with respect to an inner product of the form where , effectively fixing a metric for G. Using to define the Lie bracket, we can also represent the “little ad” operator, where for as, where and .

For a differentiable function on the group , one can construct the right and left Lie directional derivatives as,

In parametric form, where , and is the gradient operator for . Here, ‘’ represents the inverse transpose operation, i.e., the inverse of the transpose of the matrix. The Lie directional derivative operators are and . In the sequel, we drop the tildes, as it will be clear from the arguments whether or is considered, and likewise for the Jacobians.

3.2. One-Asset Black-Scholes as a Diffusion on

The set of positive real numbers equipped with the multiplication operation forms a commutative (Abelian) group. This is a subgroup of the general linear group of one dimension and is represented by . The basis of the Lie algebra is the number 1 and for a group element , the corresponding element in the Lie algebra is . Using (21) for a differentiable function , we obtain, where both left and right Lie derivatives yield the same result since the group is Abelian.

Using this relationship, we can rewrite the one-asset Black-Scholes Equation (9) as, using the shorthand . The variable-coefficient Black-Scholes equation has thus been transformed to a constant-coefficient equation on . Additionally, using time reversal to convert the backward parabolic equation to a forward parabolic equation and setting we have, which is a diffusion equation with drift and diffusivity . Here, is interpreted as a function over . Henceforth the term initial condition will be with respect to , which due to time reversal, refers to a final condition with respect to t.

We note that (9) by itself is not a diffusion process. Instead, we obtain a diffusion equation in (24) only after the time reversal and making the exponential transformation . Therefore, solving for indirectly solves for and solves the Black-Scholes equations. In subsequent sections, we apply similar transformations, while noting that solutions of such diffusion equations for u indirectly provides a solution for V.

3.3. Two-Asset Black-Scholes as a Diffusion on

An element parameterized as can be represented by a diagonal matrix as,

The Lie algebra of can be represented by the orthonormal basis,

The group is Abelian and much like , the left and right Lie derivatives coincide for a differentiable function as, where . This allows us to write the two-asset Black-Scholes Equation (12) in terms of Lie derivatives as, where and . Hence, the solution of (28) allows one to construct the solution to the two-asset Black-Scholes Equation (12).

3.4. Option Price Evolution with Coupled Assets as a Diffusion on

The affine group of the positive real line consists of all that transforms the scalar to that is, . An element can be expressed as,

The group action is a matrix multiplication of g with (where x is expressed in homogeneous coordinates as ). The Lie algebra of is two dimensional and spanned by the orthonormal basis elements and ,

Since , for in the Lie algebra of , . Or equivalently for , we would have . Using (18), the adjoint representation of would be,

An important feature of is that the determinant of is a, which is generally not equal to 1. This implies that the group is not unimodular. For non-unimodular groups, there exist distinct left-invariant and right-invariant Haar measures [35,36]. Using (16) and (17), the left and right Jacobians and of expressed as matrices are, and since the Lie group is not unimodular.

The left and right Lie directional derivatives can be evaluated using (21) for a function as, where . Therefore, we can now rewrite (15) in terms of Lie derivatives of as where and . Hence, the solution of (34) allows one to construct the solution to the coupled asset model in Equation (15).

4. Mean and Covariance Propagation for Unimodular Lie Groups

Diffusion equations on unimodular Lie groups, such as (24) and (28) can be solved approximately using mean and covariance propagation techniques developed in [6,7,8,9,10,11,12,16]. This technique would not be applicable directly on the affine group in the solution of Equation (34) since the group is not unimodular; instead, the later sections will show how the technique can be modified by converting (34) to a diffusion over the cotangent bundle of the affine group, which is unimodular. This section will review the theory of covariance propagation for a general N-dimensional unimodular matrix Lie group G with Lie algebra and apply the technique to solve (24) and (28).

Consider a general diffusion equation for on G in terms of right Lie derivatives as, where the drift vector and diffusivity matrix are independent of g but can be time-dependent in general. Since the equation is a Fokker–Planck equation over a Lie group [17,33], one can interpret as a probability density (see Appendix B). This analogy also assumes that is square-integrable in G. Additionally, , which also holds true for all values of time. Here, is the Haar measure [35,36].

Given an initial condition , we aim to propagate the solution to the next time-step () to obtain . In general, knowing we seek to obtain the solution at . The solution at can be obtained by a group convolution [33] from the solution at t as, where is the fundamental solution/Green’s function that propagates the solution over . We recognise that over a sufficiently small the fundamental solution evolves in . Since can be bijectively mapped to the Euclidean space , an approximate expression for would be a multivariate Gaussian that is a function of the N-dimensional used to parameterize the group element . This Gaussian would evolve in with a drift of . This drift can be injected on the group G by an evolution of the form where ∧ maps an element of to [12,17]. The covariance of this Gaussian would be , which would equal the group-theoretic covariance when restricted to small . The definitions of the group-theoretic mean, , and covariance, , would be provided later in (39) and (40), respectively, but here it suffices to note that the group-theoretic mean and covariance match the mean and covariance defined in Euclidean space (see Appendix C) for small , thereby allowing us to construct in terms of the group parameters [17,33] as, where the log maps the group element to and ∨ maps the element of the Lie algebra to . This form naturally extends to a Gaussian over the group in (47). If , which is a group Dirac delta distribution centred about , the distribution at remains tightly focused (the rigorous definition of such a distribution is provided in [9]) and from (36). The propagated solution would be, which is a specialized case of propagation of probability density using transition probabilities in a continuous-time Markovian process. We can define the group-theoretic mean and covariance using the following equations [17,33]: where defines a distance from the mean. Note that is usually only approximately equal to [8]. Similar definitions of mean and covariance are constructed in [37,38,39] and alternative definitions of an algebraic covariance have also been discussed in [17,40,41].

It is possible to obtain expressions for and without explicitly solving for using the concepts of mean and covariance propagation developed in [6,7,8,9,10,11,12,16] by using the relation in (38). For , we have the Baker-Campbell-Hausdorff (BCH) formula: where [8,11],

Here, , is the Lie bracket in . Since and making use of the fact that we have, where the expansion in (41) is truncated to second-order with respect to , i.e., higher order terms involving and , and above are neglected since once substituted in (39) and (40), they give rise to products of covariances or higher moments, which are assumed to be negligible. Using the definitions of mean and covariance in (39) and (40) for a convolved probability density , we obtain the following expressions for the mean and covariance of the convolved distribution correct to second order: where F captures the second order propagation in the covariance. Its exact form is derived for in [8] and more generally for unimodular groups in [11,33]. In this section, we concentrate on first order propagation since it admits a closed-form solution.

Successive compositions of using (43) provides the following closed form solution for the mean, which is obtained by integrating the exponent over time [12] (assuming that the initial condition is the identity element of the group): where is defined in (35). Setting the initial condition to be the identity element does not lead to a loss of generality since the solution at any other initial condition can be obtained by convolving the fundamental solution (37) with a group Dirac delta function at the initial condition. Using Equation (41) we see that a successive composition becomes an integration in the exponent when the Lie bracket [7]; this holds in all examples considered in this paper since is time-independent. Approximating Equation (44) to first order and recursively applying Equation (44) by discretising a domain from 0 to t into n segments with step-size and taking the limits and , we obtain the following integral [12]:

Then, an approximate solution to (35) can be constructed as, where is the initial condition. While initial conditions such as that for European put options [22] can also be considered, we focus on a Dirac delta initial condition centred at the group identity, since the solution for any other initial condition can be constructed by a convolution using (47). Moreover, the initial condition should be such that is square-integrable in G.

4.1. Mean and Covariance Propagation for Diffusion Processes on

In this subsection, we apply the propagation technique to (23). An initial condition of is set (here is the Euclidean Dirac delta distribution defined for ). Using (45) and (46), we have and . The expressions are in terms of reversed time . Substituting these expressions into (47) and noting that this is the fundamental solution , we have, and the corresponding solution for the one-asset Black-Scholes Equation (9) would be . The solution obtained from covariance propagation exactly matches the standard analytical solution obtained after applying a logarithmic transformation to (9), as described in Appendix A.

4.2. Mean and Covariance Propagation for Diffusion Processes on

We now apply the propagation technique to (28) with an initial condition of . Comparing (35) with (28), we have,

Then, the mean propagation equation would be, and . Substituting these expressions in (47), we obtain the fundamental solution as, for and . Propagation yields a solution that is equivalent to the standard solution obtained after applying the logarithmic transformations and to (12) as described in Appendix A. In summary, mean and covariance propagation applied to the one-asset and two-asset Black-Scholes models yields results that exactly match analytical solutions.

4.3. Mean and Covariance Propagation Requires Unimodularity of the Lie Group

Both and are unimodular Lie groups. That is, the Haar measure for a group G is bi-invariant and is the only ‘natural’ measure (upto an arbitrary scaling by a constant) for which the following holds, for a fixed . However is not a unimodular group, and one can define a right-invariant as well as a left-invariant Haar measure (but no bi-invariant measure exists). The theory of mean and covariance propagation in this section implicitly relies on the group being unimodular; this ensures that there is a unique bi-invariant Haar measure with respect to which a probability density function can be defined. If the group was not unimodular, a probability density defined with respect to one measure may not be a density with respect to the other. Therefore, to apply mean and covariance propagation to obtain the solution for (34), one approach would be to map the diffusion on the group to a diffusion on a related but unimodular group. In the next section, we show that that cotangent bundle of a group, when equipped with semidirect product group operation, is unimodular. This property will be exploited thereon by reinterpreting (34) as a diffusion on the cotangent bundle of the affine group.

5. Unimodularity of the Cotangent Bundle Group

Let G be an N dimensional matrix Lie group and the Lie algebra is the tangent space of the group at identity with N basis vectors, represented by for . The tangent bundle can be constructed and equipped with a group operation □ as,

Similarly, the cotangent bundle, which is the dual space to the tangent bundle is , and can be equipped with the group operation ▪ as, where there exists a bijection between and and similarly between and . Note that the two expressions indicate that the tangent and cotangent bundles have been endowed with a semidirect product. For an element and , the corresponding element in the tangent bundle would be and the group operation in the tangent bundle will be, where ∘ is the group operation in G and is the adjoint representation of where is defined such that for . If we were to represent as matrix, we have and the coadjoint representation would be . This is because for and , we can construct a bijection such that and where . Note that the ∨ operator for elements from and elements from are different but its usage will be clear from the object to which it is applied to. The inner product is defined to be , . If we use the adjoint representation of the group, a typical element of the tangent space would be and . Then, the image of the cotangent space on would transform as since . This motivates the use of as the dual to the adjoint representation in (55). Then, for , we have, where for and . Tangent and cotangent bundles, endowed with these operations have been used before in [42].

A matrix representation can be constructed for both the tangent and cotangent bundles incorporating the semidirect product as, and,

The dimension of the cotangent and tangent bundle group would be provided that the adjoint representation is N-dimensional. This may not hold in certain cases, such as for , where is zero dimensional. More generally, these constructions will result in a -dimensional semidirect product group if the group G has a trivial center. For the affine group , whose cotangent bundle group will be considered in the paper, this is true. However for , the center is one-dimensional (consisting of scalar multiples of the identity matrix ) and would be dimensional. We postpone the discussion of constructing the cotangent and tangent bundle groups for groups with non-trivial centers to a future paper.

5.1. Properties of the Adjoint Operator

The proofs of unimodularity will depend on the properties of the adjoint operator in the tangent and cotangent bundle group. We present a few properties that will be useful in constructing the proofs in the later subsections.

5.1.1.

Since and , we can convert the matrix representation of the operator into its coordinate free form through the ∧ operator. For where , we have,

The right-hand side can be simplified as,

Hence,

This proves that since it must hold for all .

5.1.2.

The relationship would follow if and . The ∨ operator defined on maps the Lie algebra of the tangent bundle to whereas the ∨ operator acting on maps to . The overloading of the ∨ symbol should not lead to confusion since the object on which it acts will determine its interpretation.

5.2. Lie Algebra of and

The basis vectors spanning the Lie algebra of the tangent bundle can be obtained by differentiating an element representation (57) with respect to each of the N variables parameterizing it, at identity. The Lie algebra basis can be obtained for the parameters parameterizing the group G as well as for the parameters parameterizing as,

Note that where is a basis of .

For defined to be a basis element of the Lie algebra of , the ∨ operator is defined such that where this time and the ∨ operator maps from this Lie algebra to . A similar definition can also be created for the Lie algebra of the cotangent bundle, . In this case, the basis of the Lie algebra can be deduced as,

Here, but , the latter of which is related to the cotangent space at identity. In general for an element of the Lie algebra of , the vee operator is defined such that, where and . A similar definition exists for the cotangent bundle, where and . The specific ∨ used will be clear from the context in terms of which argument it is applied to. The Lie algebra for the (co-)tangent bundle group defined this way is only -dimensional for groups with trivial centers.

The tangent bundle group

This can be proven by considering the adjoint representation of an element within the tangent bundle. The adjoint, where , is given through matrix notation as,

Solving this by substituting the two different forms of obtained from (61) and using the relationship that from Section 5.1.1 and Section 5.1.2 we obtain, where * denotes a matrix of the same size as . The determinant of the adjoint is . A group is unimodular if and only if the absolute value of the determinant of the adjoint is equal to 1. In this case, it can be seen that would be unimodular if and only if and if and only if where e is the identity element of G.

If G has a non-trivial center , then for a , one can decompose such that h and k commute. This is true since one can partition G such that any can be decomposed into where , , and is the ‘fundamental domain’ [33] associated with the quotient group . The dimension of will then be . Therefore, the current construction of using to represent the group G will no longer be faithful when implying that since the theorem and proof are restricted to the special construction , they only hold in their current form for the case where G has a trivial center and is dimensional.  □

The cotangent bundle group

We construct a proof in a very similar way as that in Theorem 1. Here, we obtain the adjoint representation of the cotangent bundle in the following form for ,

Since , this ensures that the cotangent bundle group is unimodular, independent of the unimodularity of G. The restriction to trivial centers follows from the fact that is only N-dimensional for groups with trivial centers and therefore will only have the same number of dimensions as the cotangent bundle with such a restriction.  □

6. Option Price Evolution with Coupled Assets as a Diffusion Process on the Cotangent Bundle of the Affine Group

Now we apply the results from the previous section to construct the cotangent bundle of the affine group. The motivation is to re-express Equation (34) in terms of the Lie derivatives of the cotangent bundle of . An element h from the cotangent bundle of the affine group, , can be expressed in the following form using (31) and (58), where . Using (62), the orthonormal basis of the Lie algebra of this group, expressed in matrices, is

For , we have , and hence we can construct a bijection from the Lie algebra to such that for . We argue that although the number of dimensions is doubled when we consider the cotangent bundle, this has the benefit of making the group unimodular. To see this, consider the left and right Jacobians of the cotangent bundle group:

Then, we observe that , implying unimodularity. The adjoint of the affine cotangent bundle is, where , as expected. For completeness, for is,

The left and right Lie derivatives are, where . We now observe from (33) and (73) that and , which allows us to rewrite (34) in terms of the cotangent bundle group Lie derivatives as, which is a degenerate diffusion over the affine cotangent bundle group. We denote the ‘master’ function as, , and as the marginalised version of this ‘master’ function marginalised over the variables x and y.

6.1. Mean and Covariance Propagation for Diffusion Processes in the Cotangent Bundle of

Since is unimodular, the theory of mean and covariance propagation can be applied directly from Section 4 to solve (74). The initial condition is where would be the identity matrix and the mean propagation equation would be, for and .

Comparing (74) with (35), we have, which is a singular matrix, implying that this is a degenerate diffusion over the cotangent bundle group. The covariance propagation equation is then, such that, where , and . Noticeably is singular, but we treat this issue by writing, where is the 2 × 2 identity matrix. Then, can be written as:

For the affine cotangent bundle, it is possible to obtain a closed-form expression for . Any can be expressed in the following form, using the expression of h in (68) and where and are functions of . Then, one has, where,

The derivations of these expressions are provided in Appendix D. For convenience of notation, we express the first two components of as and the third and fourth components with so that we have,

Substituting this form into (80) and taking the limit , we obtain, where is the Dirac delta function. Additionally, we re-express in terms of as, which can be derived by calculating the Jacobian determinant of the system of Equations in (83) and (84). We note that the variables x and y are extraneous to the original problem and we aim to marginalise the distribution from (80) to as,

Using (87) and (86), we have, which is the fundamental solution to (74).

6.2. Normalisation of Probability Distribution Functions on the Affine Cotangent Bundle Group

For the affine cotangent bundle group, we can define a general Gaussian as, for , where a general normalisation factor , which is a function of the covariance , is used. In the previous section, we used a factor of the form, for in the affine cotangent bundle group (prior to marginalisation; see Equation (80)). However, this is only correct for small covariances and assuming that the determinant of the group Jacobian (70) is sufficiently close to 1. The goal of this section is to derive a higher order correction to this normalisation factor for the affine cotangent bundle group. To do so, we first convert to exponential coordinates so that for an arbitrary element , where is the exponential coordinate parameterization of the group. Then, the right Jacobian of the group defined in (16) for this coordinate system would be where the use of makes contact with the earlier parameterization in (68); expressing and as functions of we have,

Then, we have, where since the cotangent bundle group is unimodular, we have dropped the R subscript noting that independent of the parameterization. Expanding to small we have,

We can also write this relation as, for

Since to , we have

Now if we consider the integral of (90) over the cotangent bundle group, and make the substitution , we can rewrite (100) as, by making use of the unimodularity of the cotangent bundle and noting that for a parameterization in terms of the exponential coordinates . Thus, where we use the relationship for the determinant of the Jacobian in (99). Equating (102) to the number 1 and noting that the integral evaluates to , we have,

In the context of the degenerate diffusions on the affine cotangent bundle group that arise in the coupled asset model, we have, for

We also note that the sign of the non-zero element of is negative. Whereas, all terms of are positive. This suggests that at a sufficiently large covariance, it is possible that as a consequence of the current approximation. While this was not observed for the small covariances used in this paper, this suggests an opportunity to use a different method of approximating the integral, by using the following general relation for [29,43]: where is the trace operator. Nevertheless, in this paper we do not pursue this approximation and instead use the result in (104) assuming that the eigenvalues of are sufficiently large (which is the case for the range of parameters used in the numerical simulations in the subsequent sections).

7. Numerical Results for Option Price Evolution with Coupled Assets

We solve the PDE in (74) by four methods: (1) Finite difference method (implicit and explicit), (2) first order propagation, (3) second order propagation and (4) Euler–Maruyama integration of the underlying stochastic differential equations. The theory of second order covariance propagation has not yet been introduced for , and will also be described in this section. All simulations were performed using MATLAB R2019a on a 2.7 GHz Dual-Core Intel Core i5 processor. The CPU time to run 10 time steps of second order propagation, explicit finite difference method and implicit finite difference method was 18.94 s, 0.99 s and 11.55 s (assuming that the finite difference matrices are constructed before-hand), respectively; however, if the matrix logarithms are evaluated analytically rather than numerically—which is possible for the affine cotangent bundle group (82)—the CPU time for first order propagation reduces to 0.30 s for 10 time steps.

7.1. Finite Difference Method

The 2D finite difference scheme was implemented on a rectangular a–b grid with second-order accuracy. Explicit and implicit schemes were constructed, using the forward and backward Euler scheme, respectively. A time-step of units was used for the explicit scheme and a time-step of units for the implicit scheme. A smaller time-step was used for the explicit method to ensure stability. A grid spacing of approximately units along a and units along b was used. The simulation domain was chosen to best capture the distribution while ensuring that the boundaries were sufficiently far from the mode of the distribution. This was to ensure that the Dirichlet boundary condition of could be set. A Dirac delta initial condition was approximated using a circular Gaussian with a standard deviation equal to twice the grid spacing along a.

7.2. Euler–Maruyama Integration of Underlying Stochastic Differential Equations

It is possible to convert a Fokker–Planck equation over an N-dimensional Lie group to a stochastic differential equation in the Lie algebra. A Fokker–Planck equation for a probability density of the form in the group G with a drift of and diffusivity is, where is a right directional Lie derivative in G. From [17,29], we recognise that the stochastic differential equation described by this Fokker-Planck equation on the Lie algebra would be of the form, where is a vector of increments of N uncorrelated Wiener processes, , corresponding to random draws from a Gaussian with zero mean and variance , and . Equation (108) can be interpreted either as a Stratonovich or Itô equation since the diffusion term B is independent of . Since , this process occurs in the Lie algebra of G. This process can then be injected on to the group [17] as, thereby defining a stochastic process on G. If we were to parameterize the group elements with a vector such that where is the right Jacobian matrix, we obtain a Stratonovich stochastic differential equation in the parameter space as, where Ⓢ emphasises that this is a Stratonovich equation. For the cotangent bundle group of the affine group, and the form of is provided in (70). Matching (107) with (74), we obtain, and,

Using the form of from (70) and the expressions for and B from (111) and (112), we have the following Stratonovich stochastic differential equations for : where and highlight the degenerate nature of the stochastic differential equation. However, to implement an Euler-Maruyama integration in parameter space, we require the stochastic differential equations to be expressed as Itô equations; the distinction between the Itô and Stratonovich form of a stochastic differential equation is especially important here since the diffusivity is now a function of a. The Itô version of the equations is thus, where the terms and correct for the drift. Equivalently, instead of solving (114), it is also possible to obtain an evolution in parameter space by projecting the stochastic evolution of in the cotangent bundle (109) on to the space of parameters . This is because for any we have, which determines a unique point in a stochastic trajectory in parameter space . Due to the degenerate nature of the diffusion process, there is no evolution in the parameters . This method of obtaining the stochastic process in parameter space, although equivalent in principle to that obtained by integrating (114), is henceforth referred to as an Itô-Gangolli method since it makes use of the McKean-Gangolli injection [17] and numerically solves the stochastic differential equations in the Lie algebra of the group (108) by an Euler-Maruyama integration.

The group-theoretic covariance and mean of the probability density function were deduced from the ensemble generated by (109) on the group using the methods in [8,33]. That is, the discrete version of (39) can be obtained by setting the sampled probability density to be for a total of samples. Here, is the group Dirac delta function. A similar substitution in (40) gives the sample covariance. Thus, the sampled mean and covariance are,

In the context of the Euler–Maruyama integration, is the total number of sample paths and the averaging is performed at each time slice t. The value of is obtained from the evolution process in G described in (109). One would expect that for large , and will approximate the corresponding mean and covariance of the solution of (74) to the extent that the exact solution remains a Gaussian on the cotangent bundle group. Hence, the results from the Euler-Maruyama integration of a large number of sample paths can serve as a baseline truth to compare results against. Note that (116) needs to be solved iteratively by beginning with the following guess solution, and the mean at the iteration can be computed using, where we observe that quantifies an error that goes to zero once the mean converges. This converged value is substituted as in (117) to obtain the sampled covariance. In the simulations, = 30,000 and a time-step of units was used. After generating the ensemble of sample paths, a continuous probability distribution was created in space by kernel density estimation using a Gaussian kernel. This kernel density estimated probability distribution was used as a baseline to compare against finite difference and propagation solutions (shown later in the contour plots of Figure 3 and Figure 6).

7.3. Second Order Propagation for

The mean propagates by (43) where for using the definitions in (75). Using (72) and (42) and substituting it in the definition of covariance in (40), we have the following expression for second-order covariance propagation [11,33]: where, for

Much like the case for first order propagation, only the first subspace of is non-zero, corresponding to . A time-step of units was used in second order propagation. Figure 1 depicts the convergence of the error in with respect to time-step. The baseline truth in this case was sampled from the ensemble of paths generated in the Euler-Maruyama integration (117). The parameters used in this study were: , , , and . To compare the covariances, the metric in (122) was used.

Figure 1

Convergence of second order propagated with reducing time-step, relative to the sample standard deviation . The horizontal axis shows the reversed time from 0.1 to 0.7 units.

7.4. Results

For the numerical simulations, we consider the PDE in (74) with two sets of parameters, and ; the set of parameters describes a scenario where diffusion of asset price is low and describes a scenario where the diffusion is relatively large. Since the values of , and r are same in both cases, the values of and will be used to distinguish these two sets. The constraint was used to fix the ratio in both cases. Additionally, the initial condition is a Dirac delta distribution at the identity element of the cotangent bundle group of .

Firstly, we show the differences between first and second order group-theoretic covariance propagation for the large diffusion scenario () in Figure 2. The error in covariance is plotted as a relative deviation using the following error metric, where denotes the Frobenius matrix norm, is the sampled covariance from (117) and is either the first or second order propagated covariance.

Figure 2

Relative error comparing at a given time-step with the sampled covariance for parameter values and , representing a scenario with large diffusion.

The relative error in mean was evaluated in a similar fashion: where is the sampled mean from (116) and is either the first or second order propagated mean. For the set of parameters describing ‘large diffusion’ (parameter set ) the relative error in the mean was approximately 0.1 to 0.2% and first order and second order mean propagation were indistinguishable. Furthermore, in the case of ‘small diffusion’ (parameter set ), the error in covariance was approximately 0.4 to 1.2% and with minimal difference between first order and second order propagated results; finally, the relative error in mean for this set of parameters was in the order of 0.02% and again with no difference between first and second order propagation. Only approximate values are given since this range of error is within the variability of the sampled mean and covariance itself, which are used as the baseline.

We now proceed to compare the results from first order and second order propagation with those from finite-difference methods, relative to the probability density function obtained from the Itô-Gangolli method used to indirectly solve the stochastic differential equations in (114). The probability distribution corresponding to the ensemble of points in parameter space was considered as the ground truth for the following numerical studies. One can then estimate the mean and covariance of this ground truth at time by, for sample paths and . Note that these are the expressions for the mean and covariance defined in Euclidean space (see Appendix C). It is important to emphasise that here we are not comparing on the basis of the group-theoretic mean and covariance but rather on the basis of a mean and covariance defined in . Since the solution from propagation (followed by marginalisation) or finite difference methods would yield a probability density over the affine cotangent bundle group, it is important to convert these results to an equivalent probability density function on parameter space. That is, if is a probability distribution on the affine cotangent bundle, and would be a probability distribution in parameter space . In the special case of degenerate diffusion for the coupled asset model, is a probability distribution in the Euclidean half-space . The mean and covariance of such a distribution can then be evaluated as,

For ; we also have from (70). The mean and covariance defined this way is evaluated for the finite difference and propagation solution (based on (89) but using the higher-order normalisation factor from (104)) and compared against the sampled mean and covariance obtained from the ground truth in (124,125).

7.4.1. Small Diffusion: ,

Contour plots were generated for first order propagation, second order propagation, explicit and implicit finite difference methods. The baseline probability density was smoothed by a kernel density estimation procedure using a Gaussian kernel, which was used as the ground truth for the contour plots to qualitatively assess the shape of the distribution. These are shown at 300 time steps into the simulation in Figure 3.

Figure 3

Contour plots at 300 time steps into the simulation () for the small diffusion scenario, showing the close match between the first order and second order propagation against the ground truth (kernel estimated probability density) but a worse match for the finite difference solutions.

The poor performance of the finite difference solution relative to the propagation solution is also evident in Figure 4. In this figure, the relative error is measured with respect to the sampled covariance (125) and calculated using the same formula in (122), but where is the covariance in terms of parameters from (128).

Figure 4

Relative error in covariance: Comparison between first order propagation, second order propagation and explicit and implicit finite difference (inset) for the small diffusion scenario. The propagation results nearly coincide and therefore cannot be distinguished in the plot.

The relatively poor performance of the finite difference solution is because the covariance is very small (but not zero) such that one observes spurious oscillations in the finite difference solution (see Figure 5). Negative values of are artefacts of the discretisation. In simulating the evolution of distributions with low covariance, a finite difference solution requires a very fine mesh near the mode of the distribution to avoid such artefacts whereas a propagated solution requires no such discretisation in space and automatically avoids these issues. Nevertheless, the relative error in mean, measured through a Euclidean norm, was approximately 0.2 to 0.4% for both propagation and finite difference simulations.

Figure 5

Finite difference solution (explicit) at 450 time steps () showing spurious oscillations due to the small covariances in the small diffusion scenario.

7.4.2. Large Diffusion: ,

Contour plots were generated for first order propagation, second order propagation, explicit and implicit finite difference methods. These are shown at 300 time steps into the simulation in Figure 6.

Figure 6

Contour plots at 300 time steps into the simulation () for the large diffusion scenario, with the kernel density estimated probability distribution used as the ground truth.

Figure 7 and Figure 8 show the mean and covariance evaluated using (127) and (128) compared relative to (124) and (125). The error in covariance was measured by a relative Frobenius norm and the relative error in mean was measured using a Euclidean norm defined as, for and is the sampled mean (124).

Figure 7

Relative error in mean: Comparison between first order, second order propagation and finite difference methods (explicit and implicit) for the large diffusion scenario. The propagation results nearly coincide and therefore cannot be distinguished in the plot.

Figure 8

Relative error in covariance: Comparison between first order, second order propagation and finite difference methods (explicit and implicit) for the large diffusion scenario.

Mean and covariance propagation when applied to the coupled asset model tend to yield results with lower error when the covariances are small, in line with the original assumptions made to derive the technique. Additionally, covariance propagation is more suitable in dealing with Dirac delta initial conditions than a finite difference method. Another advantage of the propagation technique as opposed to a finite difference solution is that there is no grid involved: one effectively has a solution that is not discretized in the domain and only requires a temporal discretization (for the second order propagation) or no discretization at all (for the first order propagation where the solution can be written in closed-form). Finally, a synergy of the two methods would be useful in spanning a broader range of covariances than either method can handle on its own.

8. Backward Compatibility of Propagation on the Cotangent Bundle with the One-Asset Black-Scholes Equation

In this final section, we consider reframing the one-asset Black-Scholes equation as a diffusion process on the affine cotangent bundle group. This is possible because the Lie derivative of (22) can also be represented by from (73). Hence, we can write (24) in terms of the affine cotangent bundle group Lie derivatives as,

By first order propagation, we know that the mean is (assuming an identity initial condition),

The diffusion matrix is now, and, since and D are both diagonal. Similar results can also be obtained numerically using second order propagation for both mean and covariance. Following the discussion in Section 6.1, we construct a Gaussian on the cotangent bundle of the form, where using the form of in (81) we have, so that, and . Then we can show that, by constructing a Jacobian matrix from the system of equations in (136)–(138) and using its determinant to normalise the Dirac delta function. Marginalising over the variables , we obtain the solution for the one-asset Black-Scholes equation given by propagation on the affine cotangent bundle as,

It is important to note that the solution in (140) differs from (48) by the Jacobian factor and therefore is not an exact solution. However we see that for small values of , to . In this case, it is also possible to compare the propagated result with the analytical solution for the Black-Scholes equation. We make this comparison for representing a small diffusion scenario and representing a large diffusion scenario, at and show the plots in Figure 9 and Figure 10. Moreover, the higher-order normalisation factor in (104) is used to normalise the propagated result in (140) instead of . We see a closer match with the analytical solution for the small diffusion scenario.

Figure 9

Plot of the analytical and propagated solution to the converted 1D Black-Scholes equation in (130), , for the small diffusion case (, and ). Both first order and second order propagation results coincide.

Figure 10

Plot of the analytical and propagated solution to the converted 1D Black-Scholes equation in (130), , for the large diffusion case (, and ). Both first order and second order propagation results coincide.

9. Conclusions

Reframing PDEs in as diffusion processes on Lie groups offers an alternative method to solve PDEs by using mean and covariance propagation techniques developed previously in the context of Fokker–Planck equations on Lie groups [6,7,8,9,10,12]. In the case of asset dynamics from mathematical finance, the method yields the exact solution for the one-asset and two-asset problems by matching the Lie derivatives of the one-asset and two-asset Black-Scholes equations with the Lie derivatives of and , respectively; this trivially reduces to the logarithmic coordinate transformation that converts these equations to heat equations.

While using the apparatus of mean/covariance propagation on Lie groups is undue for the one-asset and two-asset Black-Scholes equations, the matching is especially useful for the model of option price evolution under coupled asset dynamics introduced in the paper where the logarithmic coordinate transformation characteristic of the one-asset and two-asset Black-Scholes PDE can no longer be applied. Instead, we solve the equation by matching the derivatives with the Lie derivatives of . We provide proofs of the unimodularity of the cotangent bundle of a Lie group, and exploit this property to perform a mean/covariance propagation on the cotangent bundle of .

The mathematical apparatus developed can be applied to different PDEs in mathematical finance, such as those arising from stochastic volatility models or other forms of asset coupling in multi-asset models, and to linear convection-diffusion equations in transport theory, to name a few. Due to the unimodularity of the cotangent bundle group, the Lie group to which the derivatives are matched with need not be unimodular. Additionally, subsequent research can also be directed towards extending the cotangent bundle group construction presented here to groups with non-trivial centers, and to a general stability analysis of propagation schemes.