Bismut–Elworthy–Li formula for subordinated Brownian motion applied to hedging financial derivatives

Abstract

The objective of the paper is to extend the results in Fournié, Lasry, Lions, Lebuchoux, and Touzi (1999), Cass and Fritz (2007) for continuous processes to jump processes based on the Bismut–Elworthy–Li (BEL) formula in Elworthy and Li (1994). We construct a jump process using a subordinated Brownian motion where the subordinator is an inverse $\alpha$-stable process ${(L_t)}_{t\ge 0}$ with $(0,1]$. The results are derived using Malliavin integration by parts formula. We derive representation formulas for computing financial Greeks and show that in the event when $L_t\equiv t$, we retrieve the results in Fournié et al. (1999). The purpose is to by-pass the derivative of an (irregular) pay-off function in a jump-type market by introducing a weight term in form of an integral with respect to subordinated Brownian motion. Using MonteCarlo techniques, we estimate financial Greeks for a digital option and show that the BEL formula still performs better for a discontinuous pay-off in a jump asset model setting and that the finite-difference methods are better for continuous pay-offs in a similar setting. In summary, the motivation and contribution of this paper demonstrates that the Malliavin integration by parts representation formula holds for subordinated Brownian motion and, this representation is useful in developing simple and tractable hedging strategies (the Greeks) in jump-type derivatives market as opposed to more complex jump models.

Keywords:

• stable distributions
• Malliavin calculus
• integration by parts
• subordinated Brownian motion
• Bismut–Elworthy–Li formula

AMS Subject Classifications:

• C61
• C63
• G11
• G12
• G13

View correction statement:

Erratum

Public Interest Statement

This paper entitled Bismut–Elworthy–Li Formula for Subordinated Brownian Motion Applied to Hedging Financial Derivatives provides pricing and risk management methods usable by quantitative risk analysts. It provides knowledge on how one can design and determine prices of new derivative products in a jump-type financial market. In addition, market risk can be managed using Greek representations to determine for instance, proportions in one’s investment portfolio, make informed decision based on the behaviour of volatility and/or potential extreme shocks in the market. It can be argued that simpler methods could be used instead to achieve the above objectives. Agreed, however, care must be taken to ensure that simple model of choice represents the market as closely as possible. This current paper is geared towards achieving exactly that. Stochastic volatility models with similar modelled jumps could be further innovations.

1. Introduction

The problem of computating the Greeks of derivatives with smooth pay-off functions has been extensively studied. The problem where the pay-off function is not necessarily regular poses a different level of difficult and requires a different approach. See for instance Fournié, Lasry, Lions, Lebuchoux, and Touzi (1999) and Cass and Fritz (2007). Notice that existing techniques avoid differentiating the pay-off function by introducing a weight function. The Bismut–Elworthy–Li (BEL) representation formula (Elworthy & Li, 1994) is one scenario of such innovations. In this paper we show that the known relationship between the Malliavin derivative and the first variation process still holds for an alpha-stable subordinated Brownian motion and results in an explicit martingale weight factor. This allows for an extension of the BEL formula to subordinated Brownian motion and as a result, Greeks can be easily computed in jump-type and/or emerging markets. The subordinator belongs to the Lévy family of four parameter alpha-stable distributions. Price dynamics of almost all financial instruments in the different financial markets are observed to deviate from the Gaussian distribution. Various models in literature have been developed to closely estimate the dynamics of these markets. The rich and robust family of alpha-stable distributions has proven successful over most traditional models in capturing skewed and heavy tailed distributions. An application to estimate the densities of solutions to stochastic differential equation with subordinated Brownian motion under Malliavin framework is discussed in Kusuoka (2010). A rather different approach is discussed in Wyłomaǹska (2012) where the authors, in addition to investigating the densities of subordinated Brownian motion, also discuss some properties related to transforms and averaged mean squared displacements of the process. They consider and compare both cases of stable processes and their inverses. In addition they provide some parameter estimation techniques. An intuitive study related to the work in this current paper is Elworthy and Li (1994). However, this is limited to the Delta. The authors derived derivatives of solutions of diffusion equations and demonstrate that these derivatives exhibit and allow for the estimation of a diffusion equations’ smoothing properties. In addition, they use their results to study the logarithmic gradient of heat kernels. Their results can be extended to deriving representation formulas for more Greeks. A less related but still interesting work is by Song and Vondraček (2003) where properties of a killed subordinated Brownian motion by an $\alpha /2$-stable are compared with those of an $\alpha$-stable subordinated Brownian motion. They show possible comparability in their killing measures and propose bounds on the Green function and a jumping kernel of a subordinated ($\alpha /2$) process. The work by Zhang (2012) on the derivative and gradient estimate for SDEs driven by stable processes is another intuitive literature vital to our study. Other related works that are recent include Takeuchi (2010), Khedher (2012) and Kawai and Kohatsu-Higa (2010). There is more existing literature on stable distributions that we cannot exhaustively discuss. Interested readers can refer to Sun and Xie (2014) for a comprehensive literature on stable distributions.

In summary, the motivation and contribution of this paper demonstrates that Malliavin integration by parts representation formula holds for subordinated Brownian motion and, this representation is useful in developing simple and tractable hedging strategies (the Greeks) in jump-type derivatives market as opposed to more complex jump models.

The rest of the paper is organised as follows: In Section 2 we review some traditional techniques for computing the Greeks. In Section 3 we derive some useful results about stable processes and Malliavin Calculus with regards to subordinated Brownian motion. In Section 4 we analyse the differential calculus of our choice of process and show that the integration by parts formula exists. In Section 5 we present the BEL formula for subordinated Brownian motion and use the result to derive the main results of the paper. In Section 6 we discuss some applications. Section 7 concludes.

2. Sensitivity analysis

The theory of risk-neutral valuation asserts that given a complete filtered probability space $(\mathrm{\Omega },\mathcal{F},{({\mathcal{F}}_t)}_{t\ge 0},\mathtt{P})$ with a subjective probability measure $\mathtt{P}$ and a filtration ${({\mathcal{F}}_t)}_{t\ge 0}$ satisfying the usual conditions, we can construct a pay-off function $\mathrm{\Phi }$ of an option under an equivalent martingale measure whose price is given by(1) $$\begin{array}{c}\hfill V_t=\mathtt{E}[\tilde{\mathrm{\Phi }}(X_T)|{\mathcal{F}}_t],\end{array}$$(1)

where $\tilde{\mathrm{\Phi }}=e^{-rT}\mathrm{\Phi }$ is a discounted pay-off function by a risk free interest rate, r and $X_T$ is the value of the underlying at maturity time T. A Greek is a derivative of (1) with respect to a certain model parameter p which could be the initial value $X_0$ of the underlying, the volatility parameter $\sigma$, time to maturity $\tau =T-t$, the option strike, E or the interest rate, r i.e.(2) $$\begin{array}{c}\hfill \text{Greek}=\frac{\partial V_t}{\partial p}.\end{array}$$(2)

Clearly, (2) poses a problem if $\tilde{\mathrm{\Phi }}(X_T)$ is not differentiable. Different techniques in literature have been explored to resolve this challenge. For instance the Likelihood method (Fournié et al., 1999) suitable for known distribution of the underlying. It takes the form(3) $$\begin{array}{c}\hfill \frac{\partial V_t}{\partial p}=\mathtt{E}[\tilde{\mathrm{\Phi }}(p){\partial }_{p}ln(\rho (p))|{\mathcal{F}}_t],\end{array}$$(3)

where $\rho$ denotes the density of the underlying. Another widely known approach is the method of Malliavin calculus (see Bavouzet & Messaoud, 2006; Bayazit & Nolder, 2013; Bermin, 2000; Di Nunno, Øksendal, & Proske, 2008; Nualart, 1995). It eliminates differentiation of the pay-off function by introducing a weight term composed of Malliavin derivative and the Ornstein–Uhlenbeck operator terms, i.e.(4) $$\begin{array}{c}\hfill \frac{\partial V_t}{\partial p}=\mathtt{E}[\tilde{\mathrm{\Phi }}(X(p))\mathcal{E}(X(p),G(p))|{\mathcal{F}}_t],\end{array}$$(4)

where the weight factor $\mathcal{E}$ consists of Malliavin derivatives of random variables X and G belonging to some space with nice properties e.g. $L^2$. This approach is more flexible than the previous one in the sense that the distribution of the underlying is irrelevant in computing the Greeks. However, it is computationally expensive. Another approach is the Bismut–Elworthy–Li representation formula Elworthy and Li (1994), Cass and Fritz (2007) and Baños, Duedahl, Meyer-Brandis, and Proseke (2015),(5) $$\begin{array}{c}\hfill \frac{\partial V_t}{\partial x}=\mathtt{E}\left[\mathrm{\Phi }(X_T)\underset{0}{\overset{T}{\int }}{a}_s\frac{\partial X_s}{\partial x}d{B}_s\left|\right.{\mathcal{F}}_t\right];X_0=x,\end{array}$$(5)

where $p=x$, the underlying initial spot price and $a_t$ is some bounded function satisfying(6) $$\begin{array}{c}\hfill \underset{0}{\overset{T}{\int }}{a}_{s}ds=1\text{e.g.}a=\frac{1}{T}.\end{array}$$(6)

The Bismut–Elworthy–Li Formula (5) applies to continuous diffusion processes but can be adapted to finite (see Cass & Fritz, 2007) and infinite (see Chen, Song, & Zhang, 2015) jump processes. The usefulness of this formula is its allowance for an explicit representation of the Delta of a financial derivative.

For instance by employing the theory of Malliavin calculus, the weight $\mathcal{E}$ can be obtained explicitly for different Brownian motion functionals. Malliavin calculus for both continuous and jump diffusion processes has been extensively discussed in literature and there exist enormous applications on the subject (see Baños et al., 2015; Bavouzet & Messaoud, 2006, Bayazit & Nolder, 2013; Cass & Fritz, 2007; Di Nunno et al., 2008; Elworthy & Li, 1994; Kusuoka, 2010). The focus of this paper is to compute the Greeks for a wide range of pay-offs irrespective of their structure by employing the BEL formula in the framework of subordinated Brownian motion by $\alpha /2$-stable distributions, with $\alpha \in (0,2)$.

3. Stable distributions

In this section we shall adapt some of the definitions and notations from Kateregga, Mataramvura, and Taylor (2017).

Definition 3.1

A stable distribution is a four-parameter family of distributions denoted by $S(\alpha ,\beta ,\nu ,\mu )$ where

1. $\alpha \in (0,2]$ is the characteristic exponent responsible for the shape of the distribution.

2. $\beta \in [-1,1]$ is responsible for skewness of the distribution.

3. $\nu >0$ is the scale parameter (it narrows/extends the distribution around $\mu$).

4. $\mu \in \mathbb{R}$ is the location parameter (it shifts the distribution to the left or the right).

In this paper we are only interested in positive, non-decreasing, stable cádlág processes denoted by ${(S_t)}_{t\ge 0}$ whose density functions are defined on ${\mathbb{R}}_{+}$. For simplicity, we shall denote by $\mathcal{K}$ the set of all such processes. More so, we shall denote the inverse of $S_t$ by $L_t,t\ge 0$. The purpose is to introduce jumps in asset returns through subordinated Brownian motion and in turn, provide pricing and hedging (the Greeks) representation formulas.

However, the density of $S_t$ does not have a closed-form representation formula except for Gaussian ($\alpha =2$), Cauchy ($\alpha =1,\beta =0$) and Inverse Gaussian or Pearson ($\alpha =0.5,\beta =\pm 1$) distributions.

Essentially, densities of stable processes are estimated using characteristic functions through Fourier transforms.¹ The characteristic function of $S_t$ is obtained using its domain of attraction and the Lévy-Khinchine representation formula (see Applebaum, 2004) and it’s given by(7) $$\begin{array}{c}\hfill \mathrm{\Phi }(\theta )=\mathtt{E}[exp(i\theta S)]=\left\{\begin{array}{cc}exp\left(-{\nu }^{\alpha }{|\theta |}^{\alpha }\left(1-i\beta \text{sign}(\theta )tan\left(\frac{\pi \alpha }{2}\right)\right)+i\mu \theta \right);\hfill & \text{for}\alpha \ne 1.\hfill \\ exp\left(-\nu |\theta |\left(1+i\beta \text{sign}(\theta )\frac{2}{\pi }log|\theta |\right)+i\mu \theta \right);\hfill & \text{for}\alpha =1.\hfill \end{array}\right.\end{array}$$(7)

Consequently, the density of $S_t$ is given by(8) $$\begin{array}{c}\hfill h_{S_t}(u)=\frac{1}{\pi }\underset{0}{\overset{\infty }{\int }}{e}^{-iut}\mathrm{\Phi }(t)dt.\end{array}$$(8)

Definition 3.2

Let $S_t\in \mathcal{K}$. Its inverse $L_t$, $t\in [0,T]$ is defined by (see Meerschaert & Straka, 2013):(9) $$\begin{array}{c}\hfill L_s:=\left\{\begin{array}{cc}inf\left\{t:S_t>s\right\}\hfill & \text{if}s\in [0,S_t).\hfill \\ T\hfill & \text{if}s=S_T.\hfill \end{array}\right.\end{array}$$(9)

For $l\in [0,T]$, it is easy to see that the following equivalence relation holds:(10) $$\begin{array}{c}\hfill S_l<t⟺L_t\ge l.\end{array}$$(10)

The process $L_t$ is also interpreted as the first passage time of S. Moreover(11) $$\begin{array}{c}\hfill L_{S_t}=t\text{and}{S}_{L_{t-}}\le t\le S_{L_t}.\end{array}$$(11)

Let $h_L(l;t)$ denote the density function of $L_t$. Using relation (10) we deduce $F(t;l):=\mathtt{P}(S_l<t)=\mathtt{P}(L_t\ge l)=\underset{0}{\overset{\infty }{\int }}{h}_L(\tau ;t)d\tau$ which implies(12) $$\begin{array}{c}\hfill h_L(l;t)=-\frac{\partial F(t;l)}{\partial l}=-\frac{\partial }{\partial l}\underset{-\infty }{\overset{t}{\int }}{h}_S(u,l)du.\end{array}$$(12)

We can therefore approximate $h_L(l;t)$ by estimating the integral in (12) using the density $h_S(t;l)$ in (8) and its characteristic function (7).

3.1. Subordination

Let ${\left(B_t\right)}_{t\ge 0}$, denote standard Brownian motion and $(L_t,\alpha \in (0,1])$ be an inverse stable process independent of B. Then the process $B_{L_t}$ is referred to as stable subordinated Brownian motion.

Definition 3.3

(Joint Probability Space)    Suppose ${\left\{B_t\right\}}_{t\ge 0}$ is standard Brownian motion defined on $(W,\mathcal{B}(W),{\mu }^W)$ and let $(U,\mathcal{B}(U),{\mu }^U)$ be the space of non-decreasing càdlàg processes ${\left(L_t\right)}_{t\ge 0}$, $t\in [0,T]$ starting at zero. We define a joint product probability space(13) $$\begin{array}{c}\hfill (\mathrm{\Omega },\mathcal{F},\mathtt{P}):=(W\times U,\mathcal{B}(W)\times \mathcal{B}(U),{\mu }^W\times {\mu }^U),\end{array}$$(13)

where $\mathrm{\Omega }=\mathcal{C}([0,L_T])$ endowed with the natural filtrations:(14) $$\begin{array}{c}\hfill {\mathcal{F}}_t:=\sigma \left\{B_{L_{\tau }}:\tau \le L_t\right\},{\mathcal{F}}_t^W:=\sigma \left\{B_{\tau }:\tau \le L_t\right\}.\end{array}$$(14)

We further introduce a separable Hilbert space defined by(15) $$\begin{array}{c}\hfill \mathcal{H}:=\left\{h\in \mathcal{C}(\mathrm{\Omega };\mathbb{R}),h\text{is absolutely continuous and}\dot{h}\in {\mathcal{L}}^2(\mathrm{\Omega };\mathbb{R})\right\},\end{array}$$(15)

to obtain a complete abstract probability space $(\mathcal{H},\mathrm{\Omega },\mathtt{P},{({\mathcal{F}}_t)}_{t\ge 0})$.

The first and second moments of subordinated Brownian motion follow by applying semigroup properties. Literature on semigroups and infinitesimal generators is detailed in Song and Vondraček (2003) and Applebaum (2004).

Definition 3.4

Suppose $B:=(B_t,{\mathtt{P}}_x)=\sqrt{2}{B}_t^{\ast }$ where ${(B_t^{\ast })}_{t\ge 0}$ is standard Brownian motion in $\mathbb{R}$. The transition density $p(x,y;t)$ of B is given by(16) $$\begin{array}{c}\hfill p(x,0;t)=p(x;t)=\frac{1}{2\sqrt{\pi t}}{exp(-|x|}^2/4t),t>0,x,y\in \mathbb{R},\end{array}$$(16)

and its semigroup ${(P_t)}_{t\ge 0}$ given by(17) $$\begin{array}{c}\hfill P_{t}f(x)={\mathtt{E}}_x[f(B_t)]=\underset{\mathbb{R}}{\int }p(x,y;t)f(y)dy,\end{array}$$(17)

where f is a nonnegative Borel function on $\mathbb{R}$.

Definition 3.5

Suppose $X_t:=B_{L_t},t\ge 0$ is a subordinated Brownian motion. Then we define the operator $(Q_t:t\ge 0)$ satisfying(18) $$\begin{array}{c}\hfill Q_{t}f(x;t)={\mathtt{E}}_x[f(X_t;t)]={\mathtt{E}}_x[f(B_{L_t})]=\underset{0}{\overset{\infty }{\int }}{P}_{l}f(x;l)h_L(l;t)dl.\end{array}$$(18)

The operator ${(Q_t)}_{t\ge 0}$ has a transition density $q(x,y;t)$ defined by(19) $$\begin{array}{c}\hfill q(x;t)=\underset{0}{\overset{\infty }{\int }}p(x;l)h_L(l;t)dl.\end{array}$$(19)

Lemma 3.6

The mean and variance of $B_{L_t}$ are computed by(20) $$\begin{array}{cc}& {\mathtt{E}}_x[B_{L_t}]=0.\hfill \end{array}$$(20) (21) $$\begin{array}{cc}& {\mathtt{E}}_x[B_{L_t}^2]={\mathtt{E}}_x[L_t].\hfill \end{array}$$(21)

Proof

Suppose f in Definitions 3.4 and 3.5 is such that $f(z)=z$. Using (17) and (18) we immediately observe that$$\begin{array}{c}\hfill {\mathtt{E}}_x[B_{L_t}]=\underset{0}{\overset{\infty }{\int }}{\mathtt{E}}_x[B_l]h_L(l;t)dl=0.\end{array}$$

Equation (21) follows immediately by definition of variance of $B_{L_l}$. That is$$\begin{array}{c}\hfill {\mathtt{E}}_x[B_{L_t}^2]=\underset{0}{\overset{\infty }{\int }}\mathtt{E}[B_l^2]h_L(l;t)dl=\underset{0}{\overset{\infty }{\int }}l{h}_L(l;t)dl={\mathtt{E}}_x[L_t].\end{array}$$$\square$

Note that for $\alpha \in (0,1]$, $\mathtt{E}[L_t]$ exists and can be computed. If $L_t\equiv t$ in (20) and (21) we recover standard Brownian motion with mean 0 and variance t. In general, the k-th moment of $L_t$ is given by(22) $$\begin{array}{c}\hfill ⟨L_t^k⟩=\frac{\mathrm{\Gamma }(k+1)t^{k\alpha }}{\mathrm{\Gamma }(k\alpha +1)},\end{array}$$(22)

where $k\ge 1,(k\in \mathbb{R})$ and $\mathrm{\Gamma }(a)=\underset{0}{\overset{\infty }{\int }}{x}^{a-1}{e}^{-x}dx$ is the Gamma function.

Lemma 3.7

The covariance of $B_{L_t}$ is given by$$\begin{array}{c}\hfill Cov[B_{L_t},B_{L_s}]=min(\mathtt{E}[L_t],\mathtt{E}[L_s]).\end{array}$$

Proof

Let $s\le t$ then $L_s\le L_t$ and we have$$\begin{array}{cc}\hfill B_{L_t}& =B_{L_s}+(B_{L_t}-B_{L_s}).\hfill \\ \hfill B_{L_t}{B}_{L_s}& =B_{L_s}^2+B_{L_s}(B_{L_t}-B_{L_s}).\hfill \end{array}$$

Since for all $t\in {\mathbb{R}}_{+}$, $B_{L_t}$ has independent increments with zero mean, we have$$\begin{array}{cc}\hfill Cov[B_{L_t},B_{L_s}]& =\mathtt{E}[B_{L_t}{B}_{L_s}]\hfill \\ \hfill & =\mathtt{E}[B_{L_s}^2]+\mathtt{E}[B_{L_s}(B_{L_t}-B_{L_s})]\hfill \\ \hfill & =\mathtt{E}[B_{L_s}^2]+\mathtt{E}[B_{L_s}]\mathtt{E}[(B_{L_t}-B_{L_s})]\hfill \\ \hfill & =\mathtt{E}[B_{L_s}^2]\hfill \\ \hfill & =\mathtt{E}[L_s].\hfill \end{array}$$

Similarly for $L_t\le L_s$ we have the covariance as $\mathtt{E}[L_t]$. Then, we write$$\begin{array}{c}\hfill Cov[B_{L_t},B_{L_s}]=min(\mathtt{E}[L_t],\mathtt{E}[L_s]).\end{array}$$$\square$

4. Malliavin derivative in the direction of jump processes

In this section we explore the differential calculus of $B_{L_t}$ using directional derivatives. Little has been done on this, a few references include Kusuoka (2010) and Zhang (2012).

Lemma 4.1

Following Definition 3.3, we let f be an $({\mathcal{F}}_t)$-adapted right-continuous process with left limits satisfying(23) $$\begin{array}{c}\hfill \mathtt{E}\left[\underset{0}{\overset{T}{\int }}{|f(\tau -)|}^{2}d{L}_{\tau }\right]<\infty .\end{array}$$(23)

Then its $({\mathcal{F}}_t)$-martingale stochastic integral exists, it’s well defined and can be expressed as an $({\mathcal{F}}_t^W)$-martingale stochastic integral i.e.(24) $$\begin{array}{c}\hfill \underset{0}{\overset{T}{\int }}f(t-)d{B}_{L_t}=\underset{0}{\overset{L_T}{\int }}f(S_{\tau }-)d{B}_{\tau },\end{array}$$(24)

where $S_t$ is the inverse stable process of $L_t$.

Proof

This follows from the standard change of time. $\square$

Following Definition 3.3, we denote by D the Malliavin derivative operator defined on $\mathcal{H}$ such that $\dot{h}$ represents differentiation of h with respect to $L_t$. In addition, $D_h$ will denote Malliavin differentiation in direction h.

Lemma 4.2

Let ${\left(B_{L_t}\right)}_{t\ge 0}$ be subordinated Brownian motion associated with $(\mathcal{H},\mathrm{\Omega },{({\mathcal{F}}_t)}_{t\ge 0},\mathtt{P})$. Then $D_h{B}_{L_t}=h(L_t)$, for all $h\in H$ and $t\in [0,T]$. Moreover, for a Càdlàg process f we have(25) $$\begin{array}{c}\hfill D_h\underset{0}{\overset{T}{\int }}f(t-)d{B}_{L_t}=\underset{0}{\overset{L_T}{\int }}f(S_t-)dh(t),\end{array}$$(25)

where $S_t$ is the inverse stable process of $L_t$.

Proof

It is easy to see that $D_h{B}_{L_t}=h(L_t)$ since by definition(26) $$\begin{array}{c}\hfill (D_{h}F)(L_t)=\underset{\epsilon \to 0}{lim}\frac{F(L_t+\epsilon h)-F(L_t)}{\epsilon },h\in H,\end{array}$$(26)

for $F\in (\mathcal{H},\mathrm{\Omega },\mathcal{F},{({\mathcal{F}}_t)}_{t\ge 0},\mathtt{P})$, provided the limit exists.

For the second part of the lemma we notice that since $L_t$ is of bounded variation, the contribution of its small jumps is almost negligible and the number of jumps is finite. We partition $[0,T]$ as $0={\tau }_0<{\tau }_1<\cdots <{\tau }_{n-1}<{\tau }_n=T$ where ${\tau }_i,i=1,\dots ,n-1$ are the jump times of $B_{L_t}$ and let $\left\{t_{i,j};j=0,1,2,\dots ,N_i\right\}$ be a partion of $[{\tau }_{i-1},{\tau }_i)$. Suppose $\mathrm{\Delta }:={max}_{i,j}(t_{i,j}-t_{i,j-1})$, we have(27) $$\begin{array}{cc}\hfill D_h\underset{0}{\overset{T}{\int }}f(t)d{B}_{L_t}& =\underset{\mathrm{\Delta }\to 0}{lim}\sum _{i=1}^n\left[\sum _{j=1}^{n-1}f(t_{i,j-1}-)\left[h(L_{t_{i,j}})-h(L_{t_{i,j-1}})\right]\right.\hfill \\ \hfill & +f(t_{i,n_i-1}-)[h(L_{t_{i,n_i}}-)-h(L_{t_{i,n_i-1}})]\hfill \\ \hfill & \left.+f(t_{i,n_i}-)[h(L_{t_{i,n_i}})-h(L_{t_{i,n_i}}-)]\right].\hfill \end{array}$$(27)

Now if we let $u_{i,j}:=L_{t_{i,j}}$ then $S_{u_{i,j}}=t_{i,j}$ for $j=0,1,2,\dots ,N_i-1$ and as a consequence we write$$\begin{array}{cc}\hfill D_h\underset{0}{\overset{T}{\int }}f(t)d{B}_{L_t}& =\underset{\mathrm{\Delta }\to 0}{lim}\sum _{i=1}^n\left[\sum _{j=1}^{n-1}f(S_{u_{i,j-1}}-)[h(u_{i,j})-h(u_{i,j-1})]\right.\hfill \\ \hfill & +f(S_{u_{i,n_i-1}}-)[h(u_{i,n_i}-)-h(u_{i,n_i-1})]\hfill \\ \hfill & \left.+f(S_{u_{i,n_i}}-)[h(u_{i,n_i})-h(u_{i,n_i}-)]\right].\hfill \\ \hfill & =\sum _{i=1}^n\underset{L_{{\tau }_{i-1}}}{\overset{L_{{\tau }_i-}}{\int }}f(S_u-)dh(u)+\sum _{i=1}^{n}f(S_{u_i,n_i}-)[h({\tau }_i)-h({\tau }_i-)].\hfill \end{array}$$

Note that $S_u$ is constant on $u\in [L_{{\tau }_i-},L_{{\tau }_i}]$. Therefore$$\begin{array}{c}\hfill f(S_{u_i,n_i}-)[h(L_{{\tau }_i})-h(L_{\tau }-)]=\underset{L_{{\tau }_{i-1}}}{\overset{L_{{\tau }_i-}}{\int }}f(S_u-)dh(u),i=1,2,\dots ,n.\end{array}$$

The result follows immediately, i.e.(28) $$\begin{array}{c}\hfill D_h\underset{0}{\overset{T}{\int }}f(t-)d{B}_{L_t}=\underset{0}{\overset{L_T}{\int }}f(S_t-)dh(t),\end{array}$$(28) $\square$

Lemma 4.3

Suppose f is a right continuous function with left limits, then(29) $$\begin{array}{c}\hfill \underset{0}{\overset{T}{\int }}f(\tau -)d{L}_{\tau }=\underset{0}{\overset{L_T}{\int }}f(S_{\tau }-)d\tau ,\end{array}$$(29)

where $S_t$ is the inverse stable process of $L_t$.

Proof

This follows from standard change of time computations.

$\square$

Lemma 4.4

Let g be an $({\mathcal{F}}_t)$-adapted right-continuous process with left limits satisfying,(30) $$\begin{array}{c}\hfill \mathtt{E}\left[\underset{0}{\overset{T}{\int }}{|g(\tau -)|}^{2}d{L}_{\tau }\right]<\infty ,\end{array}$$(30)

then(31) $$\begin{array}{c}\hfill \underset{0}{\overset{T}{\int }}g(\tau -)d{B}_{L_{\tau }}=\underset{0}{\overset{L_T}{\int }}g(S_{\tau }-)d{B}_{\tau }a.s.\end{array}$$(31)

This holds for all $t\in [0,T]$.

Proof

This follows a similar argument as in the previous Lemma. $\square$

4.1. Discrete multiple stochastic integral

In this section we derive the integration by parts formula associated with our choice of random process $B_{L_t}$. We denote by ${\zeta }_k$ the time between the k-th and $(k+1)$-th jumps of $B_{L_t}$.

We start with the triplet $(\mathcal{H},\mathrm{\Omega },\mathtt{P})$, a joint probability space introduced in Definition 3.3 with a real seperable Hilbert space $\mathcal{H}$, with a scalar product ${⟨·,·⟩}_{\mathcal{H}}$ and the norm for an element $g\in \mathcal{H}$, denoted by ${\Vert g\Vert }_{\mathcal{H}}$, $\mathrm{\Omega }$ is the completion of $\mathcal{H}$ and $\mathtt{P}$ is the extension to the Borel $\sigma -$algebra of $\mathrm{\Omega }$ of a cylindrical measure. We define independent stable random variables ${\zeta }_k:=B_{L_k}-B_{L_{k-1}}$ which are canonical projections from $\mathrm{\Omega }$ to $\mathbb{R}$. We assume $L_k$ is a form of subordinator introduced in Janicki and Weron (1993, p. 33) which is an $\alpha /2$-stable totally skewed Lévy motion with increasing sample paths $(\alpha \in (0,2),\beta =1)$. This is a symmetric alpha-stable process ($S\alpha S$) with positive Poisson jumps. Therefore, $B_{L_k}$ belongs to a class of $S\alpha S$ Lévy motion processes with its jumps only at the jump times of $L_k$. As a consequence, we use Charlier polynomials to define the multiple stochastic integrals with respect to our process.

Definition 4.5

The Charlier polynomials are defined as(32) $$\begin{array}{c}\hfill C_n^{(\lambda )}(x)={(-1)}^n{\lambda }^{-x}{e}^{\lambda }\frac{d^n}{d{\lambda }^n}(e^{-\lambda }{\lambda }^x).\end{array}$$(32)

They can also be expressed explicitly as$$\begin{array}{c}\hfill C_n^{(\lambda )}(x)=\sum _{k=0}^n\left(\begin{array}{c}n\\ k\end{array}\right){(x)}_k{(-\lambda )}^{n-k},x\in \mathrm{\Omega },\end{array}$$

where ${(x)}_k:=x(x-1)\cdots (x-k+1)$.

The Charlier polynomials form an orthogonal basis of ${\mathcal{L}}^2(\mathrm{\Omega },\mathcal{F},\mathtt{P})$ with respect to the Poisson measure $\mu (dx)=\frac{{\lambda }^x}{x!}{e}^{\lambda }dx$. Moreover, we have:$$\begin{array}{c}\hfill \sum _{x=0}^{\infty }\frac{{\lambda }^x}{x!}{C}_n^{(\lambda )}(x)C_m^{(\lambda )}(x)=\frac{n!e^{-\lambda }}{{\lambda }^n}{\delta }_{\mathit{nm}},\lambda >0,\end{array}$$

where ${\delta }_{\mathit{nm}}=0$ when $n\ne m$ and ${\delta }_{\mathit{nm}}=1$ for $n=m$. Therefore any function $F\in {\mathcal{L}}^2(\mathrm{\Omega },\mathcal{F},\mathtt{P})$ can be uniquely represented as$$\begin{array}{c}\hfill F(x)=\sum _{n\ge 0}{f}_n{C}_n^{(\lambda )}(x),f_n\in {\mathbb{R}}_{+}.\end{array}$$

with the corresponding norm given by ${\Vert f\Vert }^2=\sum _{n\ge 0}{|f_n|}^2{\lambda }^{n}n!$.

Next, we construct the discrete multiple stochastic integral using the Wick product following similar arguments in Privault (1990).

Suppose ${\mathcal{P}}^{\ast }$ is a set of all functionals of the form $Q({\zeta }_0,\dots ,{\zeta }_{n-1})$ where Q is a real polynomial and $n\in \mathbb{N}$, we regard ${\mathcal{P}}^{\ast }$ as an algebra generated by $\left\{C_n^{(\lambda )}({\zeta }_k):k,n\in \mathbb{N}\right\}$ and define the Wick product in the following manner.

Definition 4.6

The Wick product of two elements $F,G\in {\mathcal{P}}^{\ast }$ denoted by $F\diamond G$ is defined (relaxing $\lambda$ for simplicity) as:$$\begin{array}{cc}\hfill (C_{n_1}({\zeta }_{k_1})\cdots C_{n_d}({\zeta }_{k_d}))& \diamond (C_{m_1}({\zeta }_{k_1})\cdots C_{m_d}({\zeta }_{k_d}))\hfill \\ \hfill & =\frac{(n+m)!}{n!m!}{C}_{n_1+m_1}({\zeta }_{k_1})\cdots C_{n_d+m_d}({\zeta }_{k_d}),\hfill \end{array}$$

where for $a\in {\mathbb{N}}^d$, $a!=a_1!\dots a_d!$ and $n=(n_1,\dots ,n_d)$, $m=(m_1,\dots ,m_d)$ and $k_1\ne \cdots \ne k_d$.

Let $\mathcal{H}=l^2(\mathbb{N})$ be a space of square-summable sequences. There exists a discrete chaotic decomposition of ${\mathcal{L}}^2(\mathrm{\Omega },\mathtt{P})$ whose elements F can each be represented as a sum of multiple stochastic integrals of kernels of ${\mathcal{H}}^{\circ n}=l^2{(\mathbb{N})}^{\circ n}$. That is,(33) $$\begin{array}{c}\hfill F=\sum _{n=0}^{\infty }{I}_n(f_n),\end{array}$$(33)

where $f_n\in {\mathcal{H}}^{\circ n},n\in \mathbb{N}$ and $I_n(f_n)$ is the discrete multiple stochastic integral of symmetric functions of discrete variable.

The stochastic integral of $f\in l^2(\mathbb{N})$ is an isometry from $\mathcal{H}=l^2(\mathbb{N})$ to ${\mathcal{L}}^2(\mathrm{\Omega },\mathcal{F},\mathtt{P})$ as follows:

Proposition 4.7

Let $f\in l^2(\mathbb{N})$ and define $I_1(f):=\underset{0}{\overset{+\infty }{\int }}f(B_{L_t})d{B}_{L_t}$. Then(34) $$\begin{array}{c}\hfill \mathtt{E}[I_1{(f)}^2]=\sum _{k,l\in \mathbb{N}}{f}_k{f}_l\mathtt{E}[{\zeta }_k{\zeta }_l]=\sum _0^{\infty }{f}_k^2={\Vert f^2\Vert }_{{\mathcal{L}}^2([0,L_T],\mathrm{\Omega })}.\end{array}$$(34)

Proof

Consider a partition ${\tau }_1<\cdots <{\tau }_{m-1}$ where ${\tau }_i,i=1,\dots ,m-1$ are the jump times of $B_{L_t}$ and let $\left\{t_{i,j};j=0,1,2,\dots ,n_i\right\}$ be a partion for each $[{\tau }_{i-1},{\tau }_i)$. Suppose $\mathrm{\Delta }:={max}_{i,j}(t_{i,j}-t_{i,j-1})$, we have$$\begin{array}{c}\hfill I_1(f)=\underset{\mathrm{\Delta }\to 0}{lim}\sum _{\genfrac{}{}{0pt}{}{i=1}{0\le j\le n_i}}^{i=m-1}f(B_{L_{t_{i,j}}})(B_{L_{t_{i,j+1}}}-B_{L_{t_{i,j}}})+f(B_{L_{t_{n_i}}})(B_{L_1}-B_{L_{t_{n_i}}}),\end{array}$$

where $n_i=max\left\{j:t_j\le 1\right\}$. Consequently$$\begin{array}{cc}\hfill \mathtt{E}[I_1{(f)}^2]& =\underset{\mathrm{\Delta }\to 0}{lim}\sum _{\genfrac{}{}{0pt}{}{i=1}{0\le j_1,j_2\le n_i}}^M\mathtt{E}\left[f(B_{L_{t_{i,j_1}}})f(B_{L_{t_{i,j_2}}})\right.\hfill \\ \hfill & \left.\times (B_{L_{t_{i,j_1+1}}}-B_{L_{t_{i,j_1}}})(B_{L_{t_{i,j_2+1}}}-B_{L_{t_{i,j_2}}})\right].\hfill \end{array}$$

We notice that for $j_1<j_2$ on each $[{\tau }_i,{\tau }_{i+1})$, we have$$\begin{array}{c}\hfill \mathtt{E}[f(B_{L_{t_{i,j_1}}})f(B_{L_{t_{i,j_2}}})(B_{L_{t_{i,j_1+1}}}-B_{L_{t_{i,j_1}}})(B_{L_{t_{i,j_2+1}}}-B_{L_{t_{i,j_2}}})]=0,\end{array}$$

which is arrived at by conditioning with respect to ${\mathcal{F}}_{t_{i,j_2}}$ and applying the tower property. Meanwhile for $j_1=j_2=j$ we have$$\begin{array}{cc}& \mathtt{E}[f^2(B_{L_{t_{i,j}}}){(B_{L_{t_{i,j+1}}}-B_{L_{t_{i,j}}})}^2].\hfill \\ \hfill & =\mathtt{E}[f(B_{L_{t_{i,j}}}^2)\mathtt{E}[{(B_{L_{t_{i,j+1}}}-B_{L_{t_{i,j}}})}^2|{\mathcal{F}}_{t_{i,j}}]].\hfill \\ \hfill & =\mathtt{E}[f^2(B_{L_{t_{i,j}}})\mathtt{E}[(L_{t_{i,j+1}}-L_{t_{i,j}})|{\mathcal{F}}_{t_{i,j}}].\hfill \\ \hfill & =\mathtt{E}[f^2(B_{L_{t_{i,j}}})(L_{t_{i,j+1}}-L_{t_{i,j}})].\hfill \end{array}$$

The last equation follows from the law of total expectation. The result follows immediately by combining both cases. $\square$

The discrete multiple stochastic integral $I_n(f_n)$, $f_n$ symmetric in $l^2({\mathbb{N}}^n)$ with finite support can be defined directly using the Wick product.

Definition 4.8

The symmetric tensor product $f_1\circ \cdots \circ f_n$ is defined as(35) $$\begin{array}{c}\hfill f_1\circ \cdots \circ f_n=\frac{1}{n!}\sum _{\sigma \in {\mathrm{\Sigma }}_n}{f}_{\sigma (1)}\otimes \cdots \otimes f_{\sigma (n)}\end{array}$$(35)

where $f_1,\dots ,f_n\in \mathcal{H}$ and ${\mathrm{\Sigma }}_n$ is the set of all permutations of $\left\{1,\dots ,n\right\}$. Moreover, suppose $g_1,\dots ,g_n\in l^2(\mathbb{N})$ with finite supports, we have$$\begin{array}{c}\hfill I_n(g_n\circ \cdots \circ g_n)=I_1(g_1)\diamond \cdots \diamond I_1(g_n),\end{array}$$

where $I_1(g_i)=\sum _{k=0}^{\infty }{g}_i(k)C_1(k)$, $1\le i\le n$.

The definitions above suggest the results (Privault, 1990, Prop. 2 and 3) for the Poisson process, also hold for our choice of process $B_{L_t}$ with similar proofs.

Lemma 4.9

 

(1)

Suppose ${(e_n)}_{n\in \mathbb{N}}$ is a canonical basis in $l^2{(\mathbb{N})}^{\circ n}$. Then (a) For $k_1\ne \cdots \ne k_d$ and $n_1+\cdots +n_d=n$, we have $$\begin{array}{c}\hfill I_n(e_{k_1}^{\circ n_1}\circ \cdots \circ e_{k_d}^{\circ n_d})=n_1!\cdots n_d!C_{n_1}({\zeta }_{k_1})\cdots C_{n_d}({\zeta }_{k_d}).\end{array}$$ b Suppose $f=\sum _{k=0}^{\infty }{f}_k{e}_k\in l^2(\mathbb{N})$ has finite support, then $$\begin{array}{c}\hfill I_n(f^{\circ n})=n!\sum _{\genfrac{}{}{0pt}{}{\genfrac{}{}{0pt}{}{k_1\ne \cdots \ne k_d}{n_1+\cdots +n_d=n}}{n_1,\dots ,n_d>0}}{f}_{k_1}^{n_1}\cdots f_{k_d}^{n_d}{C}_{n_1}({\zeta }_1)\cdots C_{n_d}({\zeta }_d).\end{array}$$ (c) If $f_n\in l^2({\mathbb{N}}^n)$, $g_m\in l^2({\mathbb{N}}^m)$ are symmetric with finite supports, then $$\begin{array}{c}\hfill I_n(f_n)\diamond I_m(g_m)=I_{n+m}(f_n\circ g_m).\end{array}$$

(2)

Let ${\mathcal{D}}_n=\left\{(k_1,\dots ,k_n)\in {\mathbb{N}}^n:\exists i\ne j\text{such that}{k}_i=k_j\right\}$ which represents the diagonals in ${\mathbb{N}}^n$ and let ${\mathcal{X}}_n={\mathbb{N}}^n\backslash {\mathcal{D}}_n$. Suppose $f_n\in l^2({\mathbb{N}}^n)$ and $g_m\in l^2({\mathbb{N}}^m)$ are symmetric with finite supports, then $$\begin{array}{c}\hfill {⟨I_n(f_n),I_m(g_m)⟩}_{{\mathcal{L}}^2(\mathrm{\Omega })}=\left\{\begin{array}{cc}n!{⟨f_n,g_m⟩}_{l^2({\mathcal{X}}_n)}+{(n!)}^2{⟨f_n,g_m⟩}_{l^2({\mathcal{D}}_n)}\hfill & \text{if}n=m\hfill \\ {⟨I_n(f_n),I_m(g_m)⟩}_{{\mathcal{L}}^2(\mathrm{\Omega })}=0\hfill & \text{if}n\ne m.\hfill \end{array}\right.\end{array}$$

Proof

See Privault (1990). $\square$

Lemma 4.10

Suppose ${\mathcal{C}}_n$ defines chaos of order $n\in \mathbb{N}$ in ${\mathcal{L}}^2(\mathrm{\Omega })$ given by ${\mathcal{C}}_n=\left\{I_n(f_n):f_n\in l^2({\mathbb{N}}^{\circ n})\right\}$. Then ${\mathcal{L}}^2(\mathrm{\Omega })$ has a chaotic decomposition:$$\begin{array}{c}\hfill \mathcal{L}(\mathrm{\Omega },\mathcal{F},\mathtt{P})=\underset{n=0}{\overset{\infty }{⨁}}{\mathcal{C}}_n.\end{array}$$

Moreover, if $K_n$ is the tensor product ${\mathcal{H}}^{\circ n}$, endowed with the norm$$\begin{array}{c}\hfill |\Vert f_n{\Vert |}_n^2=n!{⟨f_n,f_n⟩}_{l^2({\mathcal{X}}_n)}+{(n!)}^2{⟨f_n,f_n⟩}_{l^2({\mathcal{D}}_n)},\end{array}$$

equivalent to ${\Vert ·\Vert }_{l^2(\mathbb{N})}$ then the Fock space $\mathcal{F}(\mathcal{H}):={⨁}_{n=0}^{\infty }{K}_n$ is isometrically isomorphic to ${\mathcal{L}}^2(\mathrm{\Omega })$.

Proof

The ${\mathcal{C}}_n$’s are orthogonal according to Lemma 4.9. Secondly, Q is dense in ${\mathcal{L}}^2(\mathrm{\Omega })$ since the polynomials of ${\mathcal{P}}^{\ast }$ are dense in ${\mathcal{L}}^2({\mathbb{R}}_{+},\frac{{\lambda }^x}{x!}{e}^{\lambda }dx)$: Suppose $F\in {\mathcal{L}}^2(\mathrm{\Omega })$ and $\mathtt{E}[F{Q}_n({\zeta }_0)\cdots Q_n({\zeta }_n)]=0,$, for any $Q_0,\dots ,Q_n,n\in \mathbb{N}$, then $\mathtt{E}[F|{\zeta }_0,\dots ,{\zeta }_n]=0,n\in \mathbb{N}$. This implies$$\begin{array}{c}\hfill \underset{n\to \infty }{lim}\mathtt{E}[F|{\zeta }_0,\dots ,{\zeta }_n]=F\mathtt{P}-a.s.,\end{array}$$

since $\mathtt{E}[F|{\zeta }_0,\dots ,{\zeta }_n]$ is a discrete-time martingale. Therefore $F=0,\mathtt{P}-a.s.$ and ${\mathcal{P}}^{\ast }$ is dense in ${\mathcal{L}}^2(\mathrm{\Omega })$. $\square$

The annihilation operator defined in (26) has an equivalent in the discrete chaotic decomposition of ${\mathcal{L}}^2(\mathrm{\Omega },\mathcal{F},\mathtt{P})$ given by(36) $$\begin{array}{c}\hfill D_k(I_n(f_n))=\sum _{p=0}^{n-1}\frac{n!}{p!}{I}_p\left(f_n(\ast ,\underset{n-p\text{times}}{\underbrace{k,\dots ,k}})\right);k\in \mathbb{N}.\end{array}$$(36)

Moreover, the following lemma holds.

Lemma 4.11

Suppose $\mathcal{U}$ denotes a dense set of elements of $u\in {\mathcal{L}}^2(\mathrm{\Omega })\otimes l^2(\mathbb{N})$ such that $u_k={\zeta }_k{h}_k,k\in \mathbb{N}$ where $h:\mathbb{N}\to Q$ has finite support in $\mathbb{N}$ and define an operator $\delta :\mathcal{U}\to {\mathcal{L}}^2(\mathrm{\Omega })$ by$$\begin{array}{c}\hfill \delta (u)=-\sum _{k=0}^{\infty }(u_k+D_k{u}_k),\end{array}$$

Then for any $F\in Dom(D)$ and $u\in Dom(\delta )$, we have(37) $$\begin{array}{c}\hfill \mathtt{E}[{(DF,u)}_{l^2(\mathbb{N})}]=\mathtt{E}[F\delta (u)],\end{array}$$(37)

where Dom represents domain and the operators D and $\delta$ are also closable and adjoint to each other.

Proof

See Privault (1990). $\square$

The above results can be extended ${\mathcal{L}}^2(\mathrm{\Omega })\otimes l^2(\mathbb{N})$ to ${\mathcal{L}}^2(\mathrm{\Omega })\otimes {\mathcal{L}}^2({\mathbb{R}}_{+})$ yielding:(38) $$\begin{array}{c}\hfill \underset{0}{\overset{L_T}{\int }}Fu(B_{\tau })d{B}_{\tau }=F\underset{0}{\overset{L_T}{\int }}u(\tau )d{B}_{\tau }-\underset{0}{\overset{L_T}{\int }}(D_{t}F)u(\tau )d\tau ,t\le \tau .\end{array}$$(38)

This in turn leads to the following duality relation:(39) $$\begin{array}{c}\hfill \mathtt{E}\left[\underset{0}{\overset{L_T}{\int }}(D_{t}F)u(\tau )d\tau \right]=\mathtt{E}\left[F\underset{0}{\overset{L_T}{\int }}u(\tau )d{B}_{\tau }\right],t\le \tau .\end{array}$$(39)

The duality formula is an important tool for finding alternative representations of derivatives of expectations of irregular functions. We discuss this later.

Note that for the rest of the article we shall use abbreviations SDE and SSDE to represent stochastic differential equation and subordinated stochastic differential equation respectively.

4.2. Malliavin derivative of solutions to subordinated SDEs

In the following, K is defined as a Hilbert space, $D^{1,2}(K)$ as a Sobolev space of K-valued functions associated with the H-derivative, ${\mathcal{L}}_2(H;K)$ as a total set of a K-valued linear operator of Hilbert-Schmidt class on H, ${\mathcal{L}}^{1,2}(d{B}_{L_t};K)$ as the total set of $({\mathcal{F}}_t)$-predictable $(\mathbb{R}\times K)$-valued functions $\sigma$ such that for $\sigma (t,X)\in D^{1,2}(\mathbb{R}\times K)$, $t\in [0,T]$ with(40) $$\begin{array}{c}\hfill {\Vert \sigma \Vert }_{\mathcal{L}(d{B}_{L_t};K)}:=\mathtt{E}{\left[|D_h{\sigma (t-,X)|}_{{\mathcal{L}}_2(H;\mathbb{R}\times K)}^{2}d{L}_t\right]}^{1/2}<\infty ,\end{array}$$(40)

where h is given in (15). We denote $\mathcal{L}(dt;K)$ as the total set of $({\mathcal{F}}_t)$-predictable $(\mathbb{R}\times K)$-valued functions b satisfying $b(t,X)\in D^{1,2}(\mathbb{R}\times K)$, $t\in [0,T]$ with(41) $$\begin{array}{c}\hfill {\Vert b\Vert }_{\mathcal{L}(dt;K)}:=\underset{0}{\overset{T}{\int }}{\left|D_{h}b(t-,X)\right|}_{{\mathcal{L}}_2(H;K)}^{2}dt<\infty .\end{array}$$(41)

Let $(\mathrm{\Omega },\mathcal{F},\mathtt{P})$ represent a joint probability space introduced in the previous sections and consider the following stochastic differential equation:(42) $$\begin{array}{c}\hfill d{X}_t=b(t,X_t)dt+\sigma (t,X_t)d{B}_{L_t};X_0=x,\end{array}$$(42)

A scenario of (14) with Lipschitz coefficients and standard Brownian motion has been discussed in Fournié et al. (1999) and Bayazit and Nolder (2013) to compute financial Greeks. The case of non-Lipschitz coefficients with subordinated Brownian motion is handled in Bavouzet and Messaoud (2006), Sun and Xie (2014) and much less related in Cass and Fritz (2007) and Di Nunno et al. (2008). For non-Lipschitz with standard Brownian motion, (see for instance Baños et al., 2015). In the current paper we are interested in a model with Lipschitz coefficients and subordinated Brownian motion.

Proposition 4.12

Suppose in the stochastic differential equation (42) has Lipschitz coefficients and assume $\sigma (t,X)\in D^{1,2}(\mathbb{R}\times K)$ and $b(t,X)\in D^{1,2}(\mathbb{R}\times K)$. Then the solution $X_t$ to (42) exists, is unique and belongs to $D^{1,2}$ for all $t\in [0,T]$.

Proof

The proof is based on Picard’s successive approximation and it follows similar steps of (Kusuoka, 2010, Thm 3.1 and Prop 2.1). $\square$

Let $[·,·]$ denote the dot product endowed on $\mathcal{H}$ in (15). A representation of the derivative of $X_t$ follows in the following proposition.

Proposition 4.13

Denote the directional derivative of $X_t$ by $D_r{X}_t[h]$, $r\le t$ where $D_r$ is the Malliavin derivative operator. Then from (42) we have(43) $$\begin{array}{c}\hfill D_r{X}_t[h]=\underset{r}{\overset{t}{\int }}{D}_{r}b(s)[h]ds+\underset{r}{\overset{t}{\int }}{D}_r\sigma (s-)[h]d{B}_{L_s}+\underset{r}{\overset{t}{\int }}\sigma (s-)dh(L_s).\end{array}$$(43)

If we assume $\left\{h\right\}$ is a complete orthonormal basis in H, then(44) $$\begin{array}{cc}\hfill D_r{X}_t& =\underset{r}{\overset{t}{\int }}{D}_{r}b(s)ds+\underset{r}{\overset{t}{\int }}{D}_r\sigma (s-)d{B}_{L_s}\hfill \\ \hfill & +\sum _{i=1}^{\infty }{h}^i(s)\otimes \underset{L_r}{\overset{L_t}{\int }}\sigma (S_s-){\dot{h}}^i(s)ds.\hfill \end{array}$$(44)

Moreover, if $\sigma \equiv 1$, then(45) $$\begin{array}{c}\hfill D_r{X}_t=\underset{r}{\overset{t}{\int }}{b}^{\prime }(X(s))D_{r}X(s)ds+\sum _{i=1}^{\infty }{h}^i(t)\otimes \underset{L_r}{\overset{L_t}{\int }}{\dot{h}}^i(s)ds.\end{array}$$(45)

Proof

First, represent (42) in its integral form:(46) $$\begin{array}{c}\hfill X_t=x+\underset{0}{\overset{t}{\int }}b(\tau ,X_{\tau })d\tau +\underset{0}{\overset{t}{\int }}\sigma (\tau ,X_{\tau })d{B}_{L_{\tau }}.\end{array}$$(46)

Then apply the product rule and (26) to the second term on the RHS of (46) to obtain (43).

If $\left\{h\right\}$ is an orthogonal basis, we can express $\sigma$ as(47) $$\begin{array}{c}\hfill \sigma =\sum _{i=1}^{\infty }\frac{[\sigma ,h^i]}{\Vert h^i\Vert }{h}^i,\end{array}$$(47)

where $h=\Vert h\Vert \widehat{h}$ and $\widehat{h}$ is a unit vector of h. For h orthonormal gives (44). $\square$

If $\sigma \equiv 1$ the last term of (45) is zero and by the Grönwall’s inequality, we have:(48) $$\begin{array}{c}\hfill D_r{X}_t=exp\left(\underset{r}{\overset{t}{\int }}{b}^{\prime }(X_s)ds\right).\end{array}$$(48)

Also note from (46) that the first variation process can be deduced as:(49) $$\begin{array}{c}\hfill \frac{\partial X_t}{\partial x}=exp\left(\underset{0}{\overset{t}{\int }}{b}^{\prime }(X_s)ds\right).\end{array}$$(49)

Combining (48) and (49) results into the following useful relation(50) $$\begin{array}{c}\hfill D_r{X}_t=\frac{\partial X_t}{\partial x}exp\left(-\underset{0}{\overset{r}{\int }}{b}^{\prime }(X_s)ds\right).\end{array}$$(50)

Alternatively,(51) $$\begin{array}{c}\hfill \frac{\partial X_t}{\partial x}=D_r{X}_t\frac{\partial X_r}{\partial x},r\le t.\end{array}$$(51)

5. BEL formula for subordinated stochastic differential equations

Bismut–Elworthy–Li formula for general Lévy processes is studied in Cass and Fritz (2007). We derive representations for subordinated Brownian motion based on (Elworthy & Li, 1994).

Proposition 5.1

Let $X_t$ be the solution to (42) on the horizon $[0,T]$ where $L_t\equiv t$ (see for instance Baños et al., 2015; Cass & Fritz, 2007) and let $\mathrm{\Phi }:\mathbb{R}\to \mathbb{R}$ denote some bounded function. Suppose we can define a functional $V_t(X_T)$ of $X_T$ by(52) $$\begin{array}{c}\hfill V_t(X_T)=\mathtt{E}[\mathrm{\Phi }(X_T];X_0=x.\end{array}$$(52)

Then the derivative of V with respect to x is given by(53) $$\begin{array}{c}\hfill \frac{\partial V_t}{\partial x}=\mathtt{E}\left[\mathrm{\Phi }(X_T)\underset{0}{\overset{T}{\int }}{a}_s\frac{\partial X_s}{\partial x}d{B}_s\right],\end{array}$$(53)

where $a_t$ is some bounded function satisfying(54) $$\begin{array}{c}\hfill \underset{0}{\overset{T}{\int }}{a}_{s}ds=1.\end{array}$$(54)

Proof

Apply the classical chain rule on $\mathtt{E}[{\partial }_x\mathrm{\Phi }(X_T)]$ and use the relation (51) followed by the chain rule in the Malliavin sense. Finally apply the duality relation (39), in that order. A similar proof is provided in Sturm (2004) using the identity $D_t{X}_T=J_T{J}_t^{-1}\sigma (t,X_{t-})1_{\left\{t\le T\right\}}$, where $J_t:=\frac{\partial X_t}{\partial X_0}$, $X_0=x$. $\square$

Equation (53) is Bismut–Elworthy–Li formula for Geometric Brownian motion.

Proposition 5.2

Suppose (46) has Lipschitz coefficients. Let $R(t,X_t)$ denote the right inverse of $\sigma (t,X_t)$ where $\sigma (t,X_t)$ is elliptic. For any function $\mathrm{\Phi }\in C_b^1(\mathbb{R})$ and $h\in \mathbb{R}$, we have (the x argument is relaxed for simplicity)(55) $$\begin{array}{c}\hfill D_h\mathtt{E}[\mathrm{\Phi }(X_t)]=\frac{1}{L_t}\mathtt{E}\left[\mathrm{\Phi }(X_t)\underset{0}{\overset{t}{\int }}R(\tau )·D_h{X}_{\tau }d{B}_{L_{\tau }}\right],\end{array}$$(55)

where(56) $$\begin{array}{c}\hfill D_h{X}_t=h(t)+\underset{0}{\overset{t}{\int }}{b}^{\prime }(s,X_s)·D_h{X}_{s}ds.\end{array}$$(56)

Proof

Working backwards and using the results obtained above we have$$\begin{array}{c}\hfill D_h\mathtt{E}[\mathrm{\Phi }(X_T)]=\mathtt{E}[D_r\mathrm{\Phi }(X_t)[h]]=\mathtt{E}[{\mathrm{\Phi }}^{\prime }(X_t)D_r{X}_t[h]].\end{array}$$

Next we apply the duality relation (39), Equation (50), Grönwall’s inequality and (45) for some arbitrary h not necessarily an orthonormal basis and $\sigma \ge 1$ in that order, where we have chosen $a_t=T\text{for all}t\in [0,T].$$\square$

We provide a detailed analysis on the above result in the following section.

5.1. Main results

This section presents the main results of the paper. The idea is to extend the results by Fournié et al. (1999) to a subordinated stochastic differential equation model by deriving the first- and second-order derivative representation for the expectation of a function that is not necessarily regular. Specifically, the idea is to by-pass the derivative of the expected (irregular) function by introducing a weight term in form of an integral with respect to subordinated Brownian motion. The results in this section are employed in the following section to estimate the Greeks using MonteCarlo simulations.

In this section, operators $\delta$ and D will be used interchangeably to represent weak derivatives, and $J_t$ shall denote the first variation process given by(57) $$\begin{array}{c}\hfill J_t={\left.\frac{\partial X_t}{\partial x},\text{and}{J}_0=\frac{\partial X_t}{\partial x}\right|}_{t=0}.\end{array}$$(57) $\mathcal{L}$ shall denote a space of bounded integrable functions, $(Q_t:t\ge 0)$ shall denote the semigroup of the solution $X_t$ to (46). Lastly, ${\left\{S_t\right\}}_{t\ge 0}$ shall denote a non-decreasing cádlág $\alpha$-stable process and ${\left\{L_t\right\}}_{t\ge 0}$ its inverse with $\alpha \in (0,1]$. We require $X_t$ to be complete to enforce some integrability conditions on $D{X}_t$.

Corollary 5.3

Let $U:\mathbb{R}\to \mathcal{L}(\mathbb{R})$ with bounded first derivative such that $\delta Q_t(U):\mathbb{R}\to \mathcal{L}(\mathbb{R})$ and define $U_X:={(df)}_X=Df(X_t)$. Then the weak derivative of $Q_t$ with respect to the initial state x is given by(58) $$\begin{array}{c}\hfill \delta Q_t{(U)}_x[J_0]=\mathtt{E}[U_X[J_t]]=\mathtt{E}[Df(X_t)[J_t]],\end{array}$$(58)

provided the last term exists. Moreover(59) $$\begin{array}{c}\hfill {(\delta Q_t(df))}_x[J_0]=d{(Q_{t}f)}_x[J_0].\end{array}$$(59)

Since $X_t$ is non-explosive, the weak derivative $\delta Q_t$ is well defined.

Proof

This follows directly from the application of a weak chain rule. $\square$

Corollary 5.4

Suppose $X_t\in \mathbb{R}$ is non-degenerate and elliptic, there exist an inverse $R(t,X_t)$ of $\sigma (t,X_t)$ smooth in $X_t$ such that $\left|R(X)[Y]\right|\le \epsilon {|Y|}^2$ for all $X,Y\in \mathbb{R}$ for some $\epsilon >0$. Moreover, if(60) $$\begin{array}{c}\hfill \underset{0}{\overset{L_t}{\int }}\mathtt{E}[|J_{S_{\tau }}{|}^2]d\tau <\infty .\end{array}$$(60)

then,(61) $$\begin{array}{c}\hfill \underset{0}{\overset{L_t}{\int }}\mathtt{E}[|R(S_{\tau },X_{S_{\tau }})J_{S_{\tau }}|]d\tau <\infty .\end{array}$$(61)

Proof

Recall that the result holds for the case of continuous processes (Cass and Fritz, 2007; Elworthy & Li, 1994). We can therefore employ similar arguments of partitioning as in the second part of the Proof of Lemma 4.2 and apply similar steps as in the continuous case but piece-wise, to arrive at the required result. $\square$

Theorem 5.5

Let $\mathrm{\Phi }:\mathbb{R}\to \mathbb{R}$ with its first derivative bounded and continuous:(62) $$\begin{array}{c}\hfill \delta Q_t(d\mathrm{\Phi })=d(Q_t\mathrm{\Phi })a.s.t\ge 0,\end{array}$$(62)

Moreover, for $x\in \mathbb{R}$, $T>0$, the derivative with respect to x is given by(63) $$\begin{array}{c}\hfill {\partial }_x\mathtt{E}[\mathrm{\Phi }(X_T)]=\frac{1}{L_T}\mathtt{E}\left[\mathrm{\Phi }(X_T)\underset{0}{\overset{L_T}{\int }}R(S_{\tau },X_{S_{\tau }})J_{S_{\tau }}d{B}_{\tau }\right],\end{array}$$(63)

where $\underset{0}{\overset{L_T}{\int }}R(S_{\tau },X_{S_{\tau }})J_{S_{\tau }}d{B}_{\tau },T\ge 0$ is a martingale.

Proof

From Corollary 5.4, let $T>0$, applying Itô’s formula to(64) $$\begin{array}{c}\hfill (t,X_t)\mapsto Q_{T-t}\mathrm{\Phi }(X_t),0\le t<T,\end{array}$$(64)

yields(65) $$\begin{array}{c}\hfill Q_{T-t}\mathrm{\Phi }(X_t)=Q_T\mathrm{\Phi }(x)+\underset{0}{\overset{t}{\int }}d{(Q_{T-\tau }\mathrm{\Phi })}_{X_{\tau }}\sigma (\tau ,X_{\tau })d{B}_{L_{\tau }},\text{for}t\in [0,T].\end{array}$$(65)

As $t\to T$, and applying the knowledge from Lemmas 4.3 and 4.4 yields(66) $$\begin{array}{c}\hfill \mathrm{\Phi }(X_T)=Q_T\mathrm{\Phi }(x)+\underset{0}{\overset{L_T}{\int }}d{(Q_{S_T-S_{\tau }}\mathrm{\Phi })}_{X_{S_{\tau }}}(\sigma (S_{\tau },X_{S_{\tau }})d{B}_{\tau }.\end{array}$$(66)

Multiplying (66) by a martingale $\underset{0}{\overset{L_T}{\int }}R(S_{\tau },X_{S_{\tau }})(J_{S_{\tau }})d{B}_{\tau }$ yields(67) $$\begin{array}{cc}\hfill \mathtt{E}\left[\mathrm{\Phi }(X_T)\underset{0}{\overset{L_T}{\int }}R(S_{\tau },X_{S_{\tau }})(J_{S_{\tau }})d{B}_{\tau }\right]& =\mathtt{E}\left[\underset{0}{\overset{L_T}{\int }}d{(Q_{S_T-S_{\tau }}\mathrm{\Phi })}_{X_{S_{\tau }}}{J}_{S_{\tau }}d\tau \right].\hfill \\ \hfill & =\mathtt{E}\left[\underset{0}{\overset{L_T}{\int }}{((\delta Q_{S_T-S_{\tau }})(d\mathrm{\Phi }))}_{X_{S_{\tau }}}(J_{S_{\tau }})d\tau \right].\hfill \\ \hfill & =\underset{0}{\overset{L_T}{\int }}{((\delta Q_{S_{\tau }})((\delta Q_{S_T-S_{\tau }})(d\mathrm{\Phi })))}_x(J_0)d\tau .\hfill \\ \hfill & =\underset{0}{\overset{L_T}{\int }}{(\delta Q_{S_T}(d\mathrm{\Phi }))}_x(J_0)d\tau .\hfill \\ \hfill & =L_T\delta Q_T{(d\mathrm{\Phi })}_x(J_0).\hfill \end{array}$$(67)

Since $J_0=1$ and $\delta Q_T{(d\mathrm{\Phi })}_x=\mathtt{E}[{(d\mathrm{\Phi })}_x]={\partial }_x\mathtt{E}[\mathrm{\Phi }]$. The required result follows. $\square$

Corollary 5.6

Suppose $X_t\in {\mathbb{R}}^2,t\ge 0$ and indexes $0\le j,k\le m$, then(68) $$\begin{array}{c}\hfill \underset{0}{\overset{L_T}{\int }}\mathtt{E}[|D{X}_{S_{\tau },x}(J_0^j){|}^2]d\tau \le \epsilon |J_0^j{|}^2,\epsilon >0.\end{array}$$(68)

and(69) $$\begin{array}{c}\hfill \underset{0\le S_{\tau }\le t}{sup}\underset{x\in {\mathbb{R}}^2}{sup}\mathtt{E}[|D^2{X}_{S_{\tau },x}(J_0^j,J_0^k)|]\le \epsilon |J_0^j||J_0^k|,\end{array}$$(69)

and(70) $$\begin{array}{c}\hfill \underset{0\le S_{\tau }\le t}{sup}\underset{x\in {\mathbb{R}}^2}{sup}\mathtt{E}[|tD{X}_{S_{\tau },x}|]\le \epsilon .\end{array}$$(70)

where $|·|$ denotes the Euclidean norm.

Proof

The proof follows directly by applying (5) where $h\equiv J_0$. $\square$

Theorem 5.7

From Corollaries 5.4 and 5.6, let $\mathrm{\Phi }:{\mathbb{R}}^2\to \mathbb{R}$ be such that its first and second derivatives are bounded and continuous and,(71) $$\begin{array}{c}\hfill d(Q_t\mathrm{\Phi })=\delta Q_t(d\mathrm{\Phi })a.s.t\ge 0\end{array}$$(71)

such that for almost all $x_j,x_k\in x$ for $0\le j,k\le 2$,(72) $$\begin{array}{cc}\hfill D^2{Q}_t\mathrm{\Phi }(x)(J_0^j)(J_0^k)& =\mathtt{E}\left[D^2\mathrm{\Phi }(X_t)(D{X}_{S_t,x}{J}_0^j,D{X}_{S_t,x}{J}_0^k)\right]\hfill \\ \hfill & +\mathtt{E}\left[D\mathrm{\Phi }(X_t)\sigma (S_t,X_{S_t})(J_0^j,J_0^k)\right],\hfill \end{array}$$(72)

where(73) $$\begin{array}{c}\hfill X_0=x,J_0^j={\left.\frac{\partial X_{S_t}}{\partial x_j}\right|}_{t=0}\text{and}{J}_0^k={\left.\frac{\partial X_{S_t}}{\partial x_k}\right|}_{t=0}.\end{array}$$(73)

Then$$\begin{array}{cc}\hfill \frac{{\partial }^2}{\partial x_j{x}_k}\mathtt{E}[\mathrm{\Phi }(X_T)]& =\frac{4}{L_T^2}\mathtt{E}\left[\mathrm{\Phi }(X_T)\underset{L_T/2}{\overset{L_T}{\int }}\left(R(S_t,X_{S_t-})\frac{\partial X_{S_t-}}{\partial x_j}\right)d{B}_t\right.\hfill \\ \hfill & \left.\times \underset{0}{\overset{L_T/2}{\int }}\left(R(S_t,X_{S_t-})\frac{\partial X_{S_t-}}{\partial x_k}\right)d{B}_t\right]\hfill \\ \hfill & +\frac{2}{L_T}\mathtt{E}\left[\mathrm{\Phi }(X_T)\underset{0}{\overset{L_T/2}{\int }}\left(DR(S_t,X_{S_t-})\frac{\partial X_{S_t-}}{\partial x_j}\frac{\partial X_{S_t-}}{\partial x_k}\right)d{B}_t\right]\hfill \\ \hfill & +\frac{2}{L_T}\mathtt{E}\left[\mathrm{\Phi }(X_T)\underset{0}{\overset{L_T/2}{\int }}\left(R(S_t,X_{S_t-})\frac{{\partial }^2{X}_{S_t-}}{\partial x_j{x}_k}{|}_{t=0}\right)d{B}_t\right].\hfill \end{array}$$

Proof

From Equation (72) we deduce$$\begin{array}{cc}\hfill L_T\left(D^2{Q}_T\mathrm{\Phi }(x)(J_0^j,J_0^k)\right)& =\mathtt{E}\left[D\mathrm{\Phi }(X_T)(J_T)\underset{0}{\overset{L_T}{\int }}R(S_{\tau },X_{S_{\tau }-})(J_{S_{\tau }})d{B}_{\tau }\right]\hfill \\ \hfill & -\mathtt{E}\left[\underset{0}{\overset{L_T}{\int }}D(Q_{S_T-S_{\tau }}\mathrm{\Phi })(X_{S_{\tau }-})\right.\hfill \\ \hfill & \left.\times (D\sigma (S_{\tau },X_{S_{\tau }-})(J_{S_{\tau }}^k)R(S_{\tau },X_{S_{\tau }-})(J_{S_{\tau }}^j)d\tau \right]\hfill \\ \hfill & +\mathtt{E}\left[\underset{0}{\overset{L_T}{\int }}(Q_{S_T-S_{\tau }}\mathrm{\Phi })(X_{S_t-})(D^2{X}_{S_t-,x})(J_0^j,J_0^k)d\tau \right].\hfill \end{array}$$

Suppose $L_T=L_t/2$ and $0\le \tau \le t/2$ then$$\begin{array}{cc}\hfill D^2{Q}_{t/2}\mathrm{\Phi }(x)(J_0^j,J_0^k)& =\frac{4}{L_t^2}\mathtt{E}\left[\mathrm{\Phi }(X_t)\underset{L_t/2}{\overset{L_t}{\int }}R(S_{\tau },X_{S_{\tau }})J_{S_{\tau }}^{j}d{B}_{\tau }\right.\hfill \\ \hfill & \left.\times \underset{0}{\overset{L_t/2}{\int }}R(S_{\tau },X_{S_{\tau }-})J_{S_{\tau }}^{k}d{B}_{\tau }\right]\hfill \\ \hfill & -\frac{2}{L_t}\mathtt{E}\left[\underset{0}{\overset{L_t/2}{\int }}D(Q_{S_t/2-S_{\tau }}\mathrm{\Phi })(X_{S_{\tau }-})\right.\hfill \\ \hfill & \left.\times (D\sigma (S_{\tau },X_{S_{\tau }-})(J_{S_{\tau }}^j)(R(S_{\tau },X_{S_{\tau }-})(J_{S_{\tau }}^k)d\tau \right]\hfill \\ \hfill & +\frac{2}{L_t}\mathtt{E}\left[\underset{0}{\overset{L_t/2}{\int }}D(Q_{S_t/2-S_{\tau }}\mathrm{\Phi })(X_{S_{\tau }})(D^2{X}_{S_{\tau },\tau })(J_0^j,J_0^k)d\tau \right].\hfill \end{array}$$

Applying Itô formula to $\left\{Q_{t-\tau }\mathrm{\Phi }(X_{\tau }):0\le L_{\tau }<L_t\right\}$ at $L_{\tau }=L_t/2$, yields(74) $$\begin{array}{c}\hfill Q_{t/2}\mathrm{\Phi }(X_t)=Q_t\mathrm{\Phi }(x)+\underset{0}{\overset{L_t/2}{\int }}D(Q_{S_t/2-S_{\tau }}\mathrm{\Phi })(X_{S_{\tau }})(\sigma (S_{\tau },X_{S_{\tau }-})d{B}_{\tau }.\end{array}$$(74)

Multiply (74) by ${\int }_0^{L_t/2}DR(S_{\tau },X_{S_{\tau }})J_{S_{\tau }}^j{J}_{S_{\tau }}^{k}dB\tau$ and ${\int }_0^{L_t/2}R(S_{\tau },x)(J_0^j,J_0^k)d{B}_{\tau }$. Next, taking expectations and applying the identity (see Elworthy & Li, 1994)(75) $$\begin{array}{c}\hfill D\sigma J^{j}R{J}^k+\sigma DR{J}^j{J}^k=0,\end{array}$$(75)

yields the required result.

$\square$

6. Applications

This section is dedicated to estimating the Delta, Gamma and Vega from two stochastic models of the asset price namely; the subordinated stochastic differential equation (SSDE) and the Geometric Brownian motion (GBm).

6.1. Hedging discontinuous-pay-off type options in Black–Scholes framework

We focus only on the digital option but the analysis could be extended to other discontinuous-pay-off type or irregular pay-off options (see Fournié et al., 1999 for instance). The Delta and Gamma follow directly from Theorems 5.5 and 5.7 respectively.

Let $X_t$ satisfy the stochastic differential equation(76) $$\begin{array}{c}\hfill d{X}_t=r{X}_{t}d{L}_t+\sigma X_{t}d{B}_{L_t};X_0=x,\end{array}$$(76)

with the solution given by(77) $$\begin{array}{c}\hfill X_t=xexp(r{L}_t+\sigma B_{L_t}).\end{array}$$(77)

Figure 1 shows the dynamics of solutions to SSDE and GBm.

Figure 1. Price evolution from SSDE and GBm models.

The Greeks are computed by conditioning on $L_T=\tau$, as follows :

6.1.1. Delta

According to Theorem 5.5(78) $$\begin{array}{c}\hfill {\partial }_x\mathtt{E}[\mathrm{\Phi }(X_T)]=\frac{1}{\tau }\mathtt{E}\left[\mathrm{\Phi }(X_T)\underset{0}{\overset{T}{\int }}\frac{1}{\sigma X_t}\frac{\partial X_t}{\partial x}d{B}_{L_t}\right].\end{array}$$(78)

Assume a discounted pay-off of a digital option i.e. $\mathrm{\Phi }(X_T)=e^{-rT}{\mathbf{1}}_{E<X_T}$ where E is the strike price, then we can express the Delta by(79) $$\begin{array}{c}\hfill {\partial }_x\mathtt{E}[\mathrm{\Phi }(X_T)]=\mathtt{E}\left[e^{-rT}{\mathbf{1}}_{\left\{E<X_T\right\}}\frac{B_{L_T}}{x\sigma \tau }\right].\end{array}$$(79)

Figure 2. Digital option delta: $\alpha =0.8,r=0.1,\sigma =0.2,S_0=110,E=100,T=1$.

Figure 2 shows the digital option delta from both the SSDE and GBm models. Observe that the delta from SSDE is slightly higher due to existence of jumps and its convergence is slightly slower.

6.1.2. Gamma

Consider $n=2$ and let $X_t^1=X_t^2$ in Theorem 5.7 then we deduce(80) $$\begin{array}{cc}\hfill \frac{{\partial }^2}{\partial x^2}\mathtt{E}[\mathrm{\Phi }(X_T]& =\frac{4e^{-rT}}{x^2{\sigma }^2{\tau }^2}\mathtt{E}[\mathrm{\Phi }(X_T)B_{L_T/2}(B_{L_T}-B_{L_T/2})]\hfill \\ \hfill & -\frac{2e^{-rT}}{x^2\sigma \tau }\mathtt{E}[\mathrm{\Phi }(X_T)B_{L_T/2}]\hfill \\ \hfill & +\frac{2e^{-rT}}{x^2\sigma \tau }\mathtt{E}[\mathrm{\Phi }(X_T)B_{L_T/2}].\hfill \end{array}$$(80)

Suppose $L_T/2\equiv L_T$. We obtain a similar expression of the Gamma as in the case of continuous Brownian motion by applying the identity $\mathtt{E}[B_{2L_T}{B}_{L_T}]=\mathtt{E}[B_{L_T}(B_{2L_T}-B_{L_T})+B_{L_T}^2]$ and simplying. That is(81) $$\begin{array}{c}\hfill \frac{{\partial }^2}{\partial x^2}\mathtt{E}[\mathrm{\Phi }(X_T]=e^{-rT}\mathtt{E}\left[\mathrm{\Phi }(X_T)\frac{1}{x^2\sigma \tau }\left(\frac{B_{L_T}^2}{\sigma \tau }-B_{L_T}-\frac{1}{\sigma }\right)\right].\end{array}$$(81)

Figure 3. Digital option gamma: $\alpha =0.8,r=0.1,\sigma =0.2,S_0=110,E=100,T=1$.

Figure 3 shows digital option gamma from the SSDE and GBm models. Again we observe that the gamma from SSDE is slightly higher than that from GBm.

6.2. Vega

The Vega can be deduced similarly using integration by parts. That is(82) $$\begin{array}{c}\hfill \frac{\partial }{\partial \sigma }\mathtt{E}[\mathrm{\Phi }(X_T]=e^{-rT}\mathtt{E}\left[\mathrm{\Phi }(X_T)\left(\frac{B_{L_T}^2}{\sigma \tau }-B_{L_T}-\frac{1}{\sigma }\right)\right].\end{array}$$(82)

Figure 4. BEL formula with $\alpha =0.8,r=0.1,\sigma =0.2,S_0=110,E=100$.

Figure 4 shows the Vega from the SSDE and GBm models. Note that the Vega from SSDE is slightly higher. Observe that despite fact the SSDE model has jumps, the convergence rate for the estimation of the Greeks from the model is as good as in the GBm model.

As a matter of interest, we apply the finite-difference method on the subordinated Brownian motion model to estimate the Greeks for a call option from the SSDE model. Recall from Fournié et al. (1999) that the finite-difference method is recommended for computing the Greeks from European options compared to the BEL formula, it performs better in this case. Figure 5 shows the estimation of the Greeks of a European call option using the SSDE model.(83) $$\begin{array}{c}\hfill {\partial }_x\mathtt{E}[\mathrm{\Phi }(X_T)]=\mathtt{E}\left[e^{-rT}max(X_T-E,0)\frac{B_{L_T}}{x\sigma L_T}\right].\end{array}$$(83)

Figure 5. Finite-difference method for Call Greeks from SSDE with $\alpha =0.7,r=0.1,\sigma =0.2,S_0=110,E=100$.

7. Conclusion

We extended the integration by parts formula approach to computing the Greeks of options with discontinuous pay-offs presented in Fournié et al. (1999) to markets with jumps. As an application, we estimated Greeks from the SSDE model and observed that BEL formula still performs well for SSDE as in the continuous diffusion models. As a concrete practical application, our model can be applied by investors in emerging/illiquid markets to construct hedge portfolios.

8. Cover image

Source: Authors.

Additional information

Funding

This work was supported by funds from the National Research Foundation of South Africa (NRF), the African Institute for Mathematical Sciences (AIMS) and the African Collaboration for Quantitative Finance and Risk Research (ACQuFRR) which is the research section of the African Institute of Financial Markets and Risk Management (AIFMRM), which delivers postgraduate education and training in financial markets, risk management and quantitative finance at the University of Cape Town in South Africa.

Notes on contributors

M. Kateregga

M. Kateregga is a PhD student at the University of Cape Town in South Africa. His research is in the field of mathematical finance and his PhD thesis is entitled Stable Distributions with Applications in Finance. The current paper is a chapter in his thesis. He is also a researcher at the African Collaboration for Quantitative Finance and Risk Research (ACQuFRR) which is the research section of the African Institute of Financial Markets and Risk Management (AIFMRM), which delivers postgraduate education and training in financial markets, risk management and quantitative finance. Mr Kateregga also works with Mira Networks in South Africa as a Software Engineer.

Notes

1 Note that the characteristic function will always exist.