Option Pricing Models Driven by the Space-Time Fractional Diffusion: Series Representation and Applications

In this paper, we focus on option pricing models based on space-time fractional diffusion. We briefly revise recent results which show that the option price can be represented in the terms of rapidly converging double-series and apply these results to the data from real markets. We focus on estimation of model parameters from the market data and estimation of implied volatility within the space-time fractional option pricing models.


8
The pricing of derivatives, and notably of options, is a central subject in mathematical finance. 9 It allows the market practitioner to estimate the value of its portfolio, and to construct appropriate 10 hedging strategies. The most popular option pricing model is the one introduced by Black and Scholes 11 [4], because of its simplicity (e.g., the option price can be expressed in terms of simple mathematical The price of an option of strike K and maturity T, is a function of market parameters such as an underlying asset price S, a risk-free interest rate r and a market volatility σ. We will denote this price by V(S, K, r, σ, t). In the case of an European option, it is characterized by its payoff, that its, its value at the exercise time T; for an European call, this value is equal to V(S, K, r, σ, T) = max{S − K, 0} : For a put option, the corresponding payoff is [K − S] + . Now we recall the principles of option pricing, 38 that is, the way of determining V(S, K, r, σ, t). The risk-neutral, or risk-free approach is based on the idea that one can construct a portfolio where the (market) risk can be totally eliminated [25]. Schematically, it consists in buying an option and selling a certain quantity ∆ (to be determined) of the underlying price, so that the total value of the portfolio reads Π = V − ∆S and therefore: On the other hand, the markets are assumed to offer no arbitrage opportunity, that is, any risk-less portfolio will have the same yield as if it were capitalized at the risk-free interest rate: Equalizing (2) and (3) and making an appropriate choice for ∆ transforms the option pricing problem 41 into the resolution of a partial differential equation with terminal condition.

42
From a more theoretical point of view, within risk-neutral approach, the price can be formulated as the discounted expectations of the terminal payoff [20]: where we have introduced the time-to-maturity τ = T − t. The expectations are to be taken under the risk-neutral measure Q, which is associated to the original probability measure P via the Radon-Nikodym derivative: dQ t dP t = e S t −µt (5) The risk-neutral parameter µ can be expressed as In the BS model, the underlying asset price S is assumed to be described by a geometric Brownian motion: It follows from Itô's lemma [19] that the total differential of the option price is: Choosing ∆ = ∂V ∂S , using (7) and equalizing (2) and (3), we have shown that the call price satisfies the famous Black-Scholes equation, which is a partial differential equation (PDE) with terminal condition: It is known that, with the change of variables then the Black-Scholes PDE (9) resumes to the diffusion (or heat) equation which is a particular case of the double fractional diffusion (22) with time fractionality γ = 1 and space fractionality α = 2. It is well known that the Green function for (11) is the heat kernel and therefore, by the method of Green functions and turning back to the initial variables, we obtain the solution for the Black-Scholes PDE (in the call case): Basic manipulations on the integral (13) yield where N(.) is the normal distribution function; formula (14) is the celebrated Black-Scholes formula for the European call. The corresponding risk-neutral parameter is therefore where L α,β (t) is the Lévy process [26]. α ∈ [0, 2] and β ∈ [−1, 1] are the so-called stability and asymmetry parameters and determine the decay of the tails and the asymmetry of the probability distributions g α,β (x, t). Under the (strong) hypothesis that β = −1 (maximal negative asymmetry hypothesis) then the distribution g α,−1 := g α possesses one heavy-tail in the negative axis, and another tail in the positive axis with exponential decay as soon as α > 1, and finite exponential moments These particular Lévy distributions are sometimes called Lévy-Pareto distributions; their relevance in financial modelling has been known since the works of Mandelbrot and Fama in the 1960s [9,18]. They are known to satisfy the space-fractional equation: where D α := α−2 D α is a particular case of the Riesz-Feller operator (24) for θ = α − 2. The condition θ = α − 2 turns out to be the fractional analogue to the probabilistic condition β = −1. Note that equation (18) degenerates into the the reduced BS equation (11) when α = 2; the corresponding call option price is then where the risk-neutral parameter µ follows from (17): and reduces to µ = − σ 2 2 in the Gaussian case (α = 2). An analytic resolution of the FMLS model has been provided in [2], under the form of a quickly convergent series representation for the call price (19). The proof is based on the Mellin-Barnes representation for the solutions of the space fractional equation (18) (see [16]): if x > 0 then Let us discuss option pricing models based on space-time (double)-fractional diffusion equation, which can be expressed as * where α ∈ (0, 2], γ ∈ (0, α]. Asymmetry parameter θ is defined in the so-called Feller-Takayasu diamond |θ| ≤ min {α, 2 − α}. * 0 D γ t denotes the Caputo fractional derivative, which is defined as * and θ D α x denotes the Riesz-Feller fractional derivative, which is usually defined via its Fourier image as Let us describe the various financial models that are included in (22). The FMLS model, although far more generic than the BS one, can still be regarded as too restrictive; this is because the maximal negative asymmetry hypothesis β = −1, or equivalently θ = α − 2, does not describe well all capital markets (in particular illiquid ones, where financial assets often exhibit an almost symmetric heavy-tail). Nevertheless, it is not a priori possible to relax the maximal negative asymmetry hypothesis, because when β = −1 the expectations (6) are known to diverge [8]. The fact that the risk-neutral parameter is infinite in this case traduces the fact that the risk cannot completely be eliminated from this class of Lévy processes. Risk-minimal (instead of risk-neutral) approach has been introduced (see [6]) in this case; an interesting possibility, to generalize the the FMLS model and to remain within the risk-neutral framework, is to allow the time derivative to be also fractional (in the Caputo sense) in eq. (18): * The corresponding call option price now reads The Green functions are also known under the form of a Mellin-Barnes line integral [16]: for any x > 0. The main difference with the Lévy-stable price (19) is that the risk-neutral parameter µ γ 47 now depends on the time-fractionality γ, and is not known analytically like in the Lévy stable case 48 (20). In [3], an efficient and simple series expansion is derived for the risk neutral parameter, as well as 49 a fast converging series expansion for the call price (26). In the next section, we discuss these results in 50 detail, and test them in the real market conditions. Now, let us provide an analytic pricing formulas for the call options driven by the fractional 53 diffusion (25) (details of the proofs can be found in [3]). We assume that 1 < α ≤ 2 and 0 < γ ≤ α. The expectations in definition (6) over the probability measure P can be expressed in terms of its probability densities g α,γ (y, τ), that is: It is possible to bring the calculation back to the non time-fractional case, by writing (see details in [12]): where g γ and g θ α are solutions of single-fractional diffusion equations where M ν (z) is a function of Wright type, admitting the following Mellin-Barnes representation [17]: Inserting (29) and (33) in (28), interverting the integrals and using (17) we obtain a Mellin-Barnes representation for the risk-neutral parameter: where µ 1 = σ √ 2 α sec πα 2 is the risk-neutral parameter in the Lévy-stable case (γ = 1), cf. eq. (20). It is possible to express the integral over a vertical line (34) as a sum of residues of its analytic continuation, on the condition that the integrand decreases sufficiently fast at infinity. For Gamma function, this condition is determined by the well-known Stirling approximation [1] Γ It follows from (35) that, for a Gamma function of linear arguments of the type Γ(as + b), a, b ∈ R, its behavior at infinity depends on the sign of a, namely: (36) easily generalizes to a ratio of products of Gamma functions of linear arguments. Let us assume that a function f admits a Mellin transform f * of the form: and that this Mellin transform converges on some non-empty strip {Re(s) ∈ (c 1 , c 2 )}, so that the Mellin inversion formula holds: Introduce the characteristic quantity ∆: It follows from (36) that ∆ governs the asymptotic behavior of f * (s): Therefore, applying the residue theorem to the inversion formula (38) yields: where the choice {Re(s) > c} or {Re(s) < c} is determined by the fact that |x| −s goes to 0 at infinity in the chosen half-plane. In the case of the Mellin-Barnes representation (34), the characteristic quantity is and is negative as soon as: It follows from rule (41) that one can express (34) as the sum of the residues in the right half-plane (we choose it because µ s−1 α 1 goes to 0 at infinity in this half plane as soon as |µ 1 | < 1, which is the case in all financial applications). These poles are induced by the singularities of the Γ( 1−s α ) term, which arise at every negative integer value −n of its argument, that is at every point of the type s = 1 + αn, n ∈ N; it is well-known the residue of the Gamma function at a negative integer −n is (−1) n n! [1], and therefore we obtain: as soon as the condition (43) is fulfilled. An interesting approximation of (44) can be easily derived 60 from the Taylor approximation log(1 + u) u: which, as expected, coincides with the Lévy-stable risk-neutral parameter µ 1 when γ = 1. As a particular case, we obtain a nice approximation for the risk-neutral parameter in the fractional Black-Scholes model (α = 2): which resumes to the well-known gaussian parameter − σ 2 2 when γ = 1.

62
In Fig. 1, we plot the evolution of µ in function of the parameters γ, σ and α; thanks to the exponential 63 convergence of the series (44), it suffices to consider only the very first few terms of the series to obtain 64 an excellent level of precision.

Option price 66
In all the following we will use the notation [log] := log S K + rτ, so that the payoff in (26) can be written: [ The integral over the Green variable y becomes a particular case of a Bêta integral, which is straightforward to calculate, and one obtains the representation for the option price: The vector c := [c 1 , c 2 ] is an element of the C 2 -polyhedra P := {(t 1 , t 2 ) ∈ C 2 , 0 < Re(t 2 ) < 1, Re(t 2 − 68 t 1 ) > 0}, which generalizes the notion of convergence strip for one-dimensional Mellin transform. The double Mellin-Barnes integral (49) can also be expressed as a sum of residues. In one dimension, it is usual to sum the residues right or left to the convergence strip of the Mellin transform, like we have done for the risk-neutral parameter, where the integral (34) has been computed by right-summing the residues to obtain the series (44). In two dimensions, this procedure generalizes to a summation to a subregion of C 2 , determined by a characteristic vector associated to the integrand. The incoming residues are computed by the two-dimensional analogue to the Cauchy formula: This procedure has been introduced in [23,24]. Namely, the characteristic quantity (39) generalizes to a characteristic vector, which, in the case of the double Mellin-Barnes integral (49) reads The rule (41) generalizes to where Π ∆ is the subset of C 2 defined by in the sense of the euclidean scalar product. In the plane Re C 2 , Π ∆ is therefore the region located under the line whose slope is positive because by hypothesis γ ≤ α. In this region, poles come from functions 71 Γ(−1 − t 1 + t 2 ) and Γ(t 2 ) which are singular at every negative integer value of their argument (see 72 Fig. 2).

73
From the singular behavior of the Gamma function around a singularity [1] and the Cauchy formula (50), we obtain the series for the call price under double-fractional model: Full details of this calculation can be found in [3].

75
Let us discuss several applications of the series formula for the space-time fractional option prices. 76 We show that it can be used for estimating the market parameters of the option prices. The calculation of 77 the option price is very quick compared to the other methods (Mellin-Barnes representation, numerical 78 estimation, . . . ). We also briefly discuss the applications to implied volatility. simple series whose m-(resp. n) partial sums converges very quickly to the option price (see fig. 3,

82
where the parameters are S = 3800, K = 4000, r = 1%, σ = 20%, α = 1.7, γ = 0.9). We may observe 83 that the convergence of the m-sums are monotone, while the n-sums oscillate around the final price. In the graphs in fig. 4, we study the evolution of the option price (55) in function of different parameters.

85
In the first graph we fix S = 3800, K = 4000, r=1% σ = 20% and we plot the evolution of the price in 86 function of γ, for different stability parameters α ∈ [1.5, 2]; we choose to consider only γ > 0.33 so 87 that the condition γ > 1 − 1 α is satisfied for all stabilities, and observe that the prices are a decreasing 88 function of the time fractionality. In graph 2 we let α vary between 1 and 2 and note that when

95
The results are presented in Tab. 1, which has been taken from [13]. Table 1. Estimated values of option pricing parameters based on Black-Scholes model, FMLS model and Space-time fractional model. The estimation was done for all options and separately for call options and put options, respectively. We see that for this case is γ very close to one, which does not have to be true for illiquid markets or during the abnormal periods. AE denotes the aggregated error, which is defined as the sum of absolute differences between estimated price and market price. The process of implying the market volatility consists in finding for which volatility σ I a model-driven option price coincides with the observable price C, that is when A typical procedure is to imply a Black-Scholes volatility (by using the Black-Scholes formula for the  When the asset is "at-the-money forward", that is when then there exists an approximation for the Black-Scholes formula [5] and therefore the solution to the implied volatility equation (56) reads Such an approximation can also be derived in the double-fractional Black-Scholes model (α = 2): note that, with our notations, the ATM-forward hypothesis (57) reads [log] = 0 and therefore in this case the pricing formula (55) becomes a power series (i.e., with only positive powers of µ γ and τ): Using approximation (46) for the risk-neutral parameter in the first order term of the power series (60), we obtain the implied fractional Black-Scholes volatility (in the ATM forward case): Let us remark that the formula (62) resumes to the Black-Scholes implied volatility formula (59) when  people trained in fractional calculus. The resulting series representation can be easily grasped by any financial practitioner. We have also applied the formulas to real financial data in order to demonstrate 127 fast convergence and stability of the method. We have particularly shown numerical estimations of 128 model parameters from the real data, applications to implied volatility and presence of volatility smile.

129
Fractional models provide a fruitful field for further investigations of financial systems, including 130 portfolio management, derivative pricing, commodity pricing and many other possible applications.

131
Naturally, in these applications it is necessary to carefully define the proper fractional derivatives and 132 boundary conditions. In some cases, as e.g. in the case of fractional geometric Brownian motion, it is 133 also necessary to overcome the mathematical issues, as non-existence of moments, etc. Some of these 134 topics will be addressed in the future research.