On the RND under Heston's stochastic volatility model
aa r X i v : . [ q -f i n . P R ] J a n On the RND under Heston’s stochasticvolatility model
Ben BoukaiDepartment of Mathematical Sciences, IUPUIIndianapolis, IN 46202 , USAJanuary 12, 2021
Abstract
We consider Heston’s (1993) stochastic volatility model for valuationof European options to which (semi) closed form solutions are avail-able and are given in terms of characteristic functions. We provethat the class of scale-parameter distributions with mean being theforward spot price satisfies Heston’s solution. Thus, we show thatany member of this class could be used for the direct risk-neutralvaluation of the option price under Heston’s SV model. In fact, wealso show that any RND with mean being the forward spot pricethat satisfies Hestons’ option valuation solution, must be a memberof a scale-family of distributions in that mean. As particular exam-ples, we show that one-parameter versions of the
Log-Normal, Inverse-Gaussian, Gamma, Weibull and the
Inverse-Weibull distributions areall members of this class and thus provide explicit risk-neutral den-sities (RND) for Heston’s pricing model. We demonstrate, via exactcalculations and Monte-Carlo simulations, the applicability and suit-ability of these explicit RNDs using already published Index data witha calibrated Heston model (
S&P500 , Bakshi, Cao and Chen (1997), and
ODAX,
Mr´azek and Posp´ıˇsil (2017)), as well as current option marketdata (
AMD ). Keywords : Heston model, option pricing, risk-neutral valuation,calibration. Introduction
The stochastic volatility model for option valuation of Heston (1993) is widelyaccepted nowadays by both, academics and practitioners. It prescribes, undera risk-neutral probability measure Q , say, the dynamics of the spot’s (stock,index) price process S = { S t , t ≥ } , in relation to a corresponding, thoughunobservable (untradable ) volatility process V = { V t , t ≥ } via a systemof stochastic deferential equations. This system is given by dS t = rS t dt + p V t S t dW ,t dV t = κ ( θ − V t ) + η p V t dW ,t , (1)where r is the risk-free interest rate, κ, θ and η are some constants (to bediscussed below) and where W = { W ,t , t ≥ } and W = { W ,t , t ≥ } aretwo Brownian motion processes under Q with d ( W W ) = ρdt .The quest to incorporate a non-constant volatility in the option valuationmodel, has risen in the literature (e.g., Wiggins (1987) or Stein and Stein(1991)) ever since the seminal work of Black and Scholes (1973) and of Merton(1973), (abbreviated here as the BSM) in modeling the price of a Europeancall option when the spot’s price was assumed to evolve, with a constantvolatility of the spot’s returns, σ , as a geometric Brownian motion, dS t = rS t dt + σS t dW ,t . (2)Coupled with an ingenious argument of instantaneous portfolio hedging (alongwith other assumptions such as self-financing, no-cost trading/carry, etc.)and an application of Ito’s Lemma to the underlying PDE, the BSM modelprovides an exact solution for the price of an European call option C ( · ).Specifically, given the current spot price S τ = S and the risk-free interestrate r , the price of the corresponding call option with price-strike K andduration T , C S ( K ) = S Φ( d ) − K e − rt Φ( d ) , (3)where t = T − τ is the remaining time to expiry. Here, using the conventionalnotation, Φ( · ) and φ ( · ) denote the standard Normal cdf and pdf , respectively,and d := − log( KS ) + ( r + σ ) tσ √ t and d := d ( k ) − σ √ t. (4)2n similarity to the form of the BSM solution in (3), Heston (1993) ob-tained that the solution to the system of PDE resulting from the stochasticvolatility model, (1), is given by C S ( K ) = S P − K e − rt P , (5)where P j j = 1 ,
2, are two related (under a risk-neutral probability measure Q ) conditional probabilities that the option will expire in-the-money, con-ditional on the given current stock price S τ = S and the current volatility, V τ = V . However, unlike the explicit BSM solution in (3) which is givenin terms of the normal (or log-normal) distribution, Heston (1993) provided(semi) closed-form solutions to these two probabilities, P and P in termsof their characteristic functions (for more details, see the Appendix). Hence, C S ( K ) in (5) is readily computable, via complex integration, for any choiceof the parameters ϑ = ( κ, θ, η, ρ ) in (1), all in addition to S, V and r . Theseparameters have particular meaning in context of the SV model (1): r is theprevailing risk-free interest rate; ρ is the correlation between the two Brown-ian motions comprising it; θ is the long-run average, κ is the mean-reversionspeed and η is the variance of the volatility V (see also Section 4 below). Itshould be noted that different choices of ϑ will lead to different values C S ( K )in (5) and hence, the value ϑ = ( κ, θ, η, ρ ) must be appropriately ‘calibrated’ first for C S ( K ) to actually match the option market data.The role of the risk-neutral probability measure Q in option valuation ingeneral and in determining the specific solution in (5) (or in (3)) in particular,cannot be overstated (in the ‘risk-neutral’ world). As was established by Coxand Ross (1976), the risk-neutral equilibrium requires that for T > τ (with t = T − τ ), E ( S T | S τ ) = Z S T d Q ( S T )= Z S T · q ( S T ) dS T = S τ e rt (6)and that (in the case of a European call option), C S τ ( K ) must alo satisfy C S τ ( K ) = e − rt E (max( S T − K, | S τ )= e − rt Z ( S T − K ) + d Q ( S T )= e − rt Z ∞ K ( S T − K ) q ( S T ) dS T , (7)3here, for any x ∈ R , x + := max ( x, q ( · ) is the risk-neutral density(RND) under Q , reflective of the conditional distribution of the spot price S T at time T , given the spot price, S τ at time τ < T . The risk-neutralprobability Q links together the option evaluation and the distribution ofspot price S T and its stochastic dynamics governing the model (as in (1),say). As was mentioned earlier, in the case of the BSM in (2) the RNDis unique and is given by the log-normal distribution. However, since theHeston (1993) model involves the dynamics of two stochastic processes, oneof which (the volatility, V) is untradable and hence not directly observable,there are innumerable many choices of RNDs, q ( · ), that would satisfy (6)-(7)and hence, the general solutions of P an P in (5) by means of characteristicfunctions (per each choice of ϑ = ( κ, θ, η, ρ )).Needless to say, there is an extensive body of literature dealing with the ‘extraction’, ‘recovery’, ‘estimation’ or ‘approximation’ , in parametric or non-parametric frameworks, of the RND, q ( · ) from the available (market) optionprices; see Jackwerth (2004), Figlewski (2010), Grith and Kr¨atschmer (2012)and Figlewski (2018), for comprehensive reviews of the subject. With theparametric approach in particular, one strives to estimate by various means(maximum likelihood, method of moments, least squares, etc.) the parame-ters of some assumed distribution so as to approximate available option dataor implied volatilities (c.f. Jackwerth and Rubinstein (1996)). This type of assumed multi-parameter distributions includes some mixtures of log-normaldistribution (Mizrach (2010), Grith and Kr¨atschmer (2012)), generalizedgamma (Grith and Kr¨atschmer (2012)), generalized extreme value (Figlewski(2010)), the gamma and the Weibull distributions (Savickas (2005)), amongothers. While empirical considerations have often led to suggesting theseparametric distributions as possible RNDs, the motivation to these consid-erations did not include direct link to the governing pricing model and itdynamics, as was the case in the BSM model, linking directly the log-normaldistribution and the price dynamics of model (2).In this paper, we present in our Theorem 2 a more direct approach (andhence, a link) to the RND quest in the case of Heston’s (1993) SV model (1).By expanding the last term in (7), with S τ = S , we obtain, C S ( K ) = e − rt Z ∞ K S T · q ( S T ) dS T − Ke − rt Q ( S T > K ) . (8)Clearly, by comparing (5) to (8), it follows that P = Q ( S T > K ) ≡ − Q ( K )(the risk-neutral probability of the option expiring in the money), whereas4y, (6) and (5), P ≡ Z ∞ K S T Se rt · q ( S T ) dS T = Z ∞ K S T E ( S T | S τ ) · q ( S T ) dS T . (9)We note that since by (6), we have, Z ∞ S T E ( S T | S τ ) · q ( S T ) dS T = 1 , the probability P is also being interpreted (see for example xxxx) as theprobability of the option expiring in the money, but under the so-calledphysical probability measure that is being dominated by Q . However here, inthe case of of the Heston’s (1993) model, we consider a different interpretationof this term P , which enables us to characterize a class of RND candidatesthat satisfy (5).It is a standard notation to denote by ∆( K ) the so-call delta function (orhedging fraction) in the option valuation, as defined by∆( K ) = ∂C S ( K ) ∂S . (10)In the Appendix, we show that for Heston’s call option price C S ( K ) as givenin (5), one has (see also Bakshi, Cao and Chen (1997)), P ≡ ∆( K ) . (11)Hence, under model (1), Heston’s solution for the option price in (5) can bewritten in an equivalent form as: C S ( K ) ≡ S · ∆( K ) − K e − rt · (1 − Q ( K )) . (12)We point out in passing that this presentation (12) also trivially applies tothe BSM option price in (3) since in that case, Φ( d ) ≡ ∆( K ) and Φ( d ) isjust the probability that option will expire in the money (calculated underthe log-normal distribution).In Section 2, we identify the class of distributions (and hence of RNDs)that admit the presentation in (12) of Heston’s (1993) option price as as isgiven in (5). Specifically, we show that any risk-neutral probability distribu-tion Q that satisfies (6)-(7) with a scale parameter µ = S · e rt would admit thepresentation in (12) and hence would satisfy Hestons’ (1993) option pricing5odel in (5) . In fact, we also show in the Appendix that the RNDs that maybe calculated directly from Heston’s characteristic functions (correspondingto P and P ) are members of this class of distributions as well. In Theorem2 below we establish the direct link, through Heston’s (1993) solution in (5)(or (12)) between this class of RNDs and the assumed stochastic volatilitymodel in (1) governing the spot price dynamics.In Section 3 we provide some specific examples of well known distribu-tions that satisfy (12). These include the Gamma, Inverse Gaussian, Log-Normal , the
Weibull and the
Inverse Weibull distributions (all with a partic-ular parametrization) as possible RND solutions for option valuation underHeston’s SV model. The extent agreement between each of these five partic-ular distributions as a possible RND for the Heston’s model, and the actualHeston’s RND (calculated from P in the Appendix) and the simulated dis-tribution of the spot prices obtained under (a discretized version of) model(1) is illustrated numerically in Section 4. We demonstrated the applicabil-ity and suitability of these explicit RNDs using already published Index datawith a calibrated Heston model ( S&P500 , Bakshi, Cao and Chen (1997), and
ODAX,
Mr´azek and Posp´ıˇsil (2017)), as well as current option market data(
AMD ). In the Appendix, we present the expressions for for Heston’s charac-teristic functions and discuss some of the immediate properties leading toour main results as stated in Theorem 2.
In this section we identify the class of distributions, and therefore of possibleRNDs that admit the presentation in (12) for the price of a European calloption. Specifically, we show that any RND candidate that satisfies (6)-(7)with a scale parameter µ = S · e rt would admit the presentation in (12) andhence in light of the result in (11) (see Claim 7 in the Appendix), wouldequivalently satisfy Hestons’ (1993) option pricing model in (5).To that end and to simplify the presentation, we consider a continuouspositive random variable X with mean µ > Q ). We denote by Q µ ( · ) and q µ ( · ) the cdf and pdf of X , respectively, to emphasize their dependency on µ , as a parameter.Similarly, for a given µ >
0, we write E µ ( · ) for the expectation of X (or6unctions thereof) under Q µ so that, µ := E µ ( X ) = Z ∞ xq µ ( x ) dx ≡ Z ∞ (1 − Q µ ( x )) dx. In similarity to (7), we define, for each s ≥ c µ ( s ) := E µ [( X − s ) + ] . Clearly, c µ (0) = µ . Note that c µ ( · ) is merely the undiscounted version of C S ( · ) in (7), so that with µ = S · e rt as in (6), we have, c µ ( K ) ≡ e rt · C S ( K ).It is straightforward to see that, as in (8), c µ ( s ) = Z ∞ s ( x − s ) q µ ( x ) dx = Z ∞ s xq µ ( x ) d x − s (1 − Q µ ( s )) , (13)or equivalently, c µ ( s ) ≡ Z ∞ s (1 − Q µ ( x )) dx. (14)Hence, it follows immediately from (14) that for each s > c ′ µ ( s ) := ∂c µ ( s ) ∂s = − (1 − Q µ ( s )) . (15)As we proceed to explore more of the basic properties of the function c µ ( · ),we add the simple assumption that µ is a scale-parameter of the underlyingdistribution of X . Assumption A:
We assume that Q := { Q µ , µ > } is a scale-family ofdistributions (under Q ), so that for any given µ > Q µ ( x ) ≡ Q ( x/µ ) and q µ ( x ) ≡ µ q ( x/µ ) , ∀ x > , for some cdf Q ( · ) with a pdf q ( · ) satisfying R ∞ xq ( x ) dx = 1 and R ∞ x q ( x ) dx < ∞ . In Lemma 1 below we establish the linear homogeneity of c µ ( s ) under Assumption A and provide in Lemma 2 the implied re-scaling property ofthis function and the consequential specific derived form of c µ ( s ) as presentedin Theorem 1 below. For the linear homogeneity property of the Europeanoptions, in general, see Theorems 6 & 9 of Merton (1973) or Theorem 2.3 inJiang (2005). 7 emma 1 Suppose that Q µ ∈ Q and thus it satisfies the conditions of As-sumption A, then the function c µ ( s ) as defined in (14) (or (13) is homoge-neous of degree one in s and in µ . That is, for s ′ = α s and µ ′ = α µ with α > , we have c µ ′ ( s ′ ) ≡ α c µ ( s ) . Proof.
With a simple change of variable, it follows immediately from
As-sumption A and (14) that with s ′ = α s , µ ′ = α µ , α >
0, we have c µ ′ ( s ′ ) = Z ∞ s ′ (1 − Q µ ′ ( x )) d x ≡ α Z ∞ s ′ /α (1 − Q µ ( u )) d u = α c µ ( s ) . Now, by applying the results of Lemma 1 with α = 1 /µ , we immediatelyobtain the following useful result. Lemma 2
Suppose that Q µ ∈ Q and thus it satisfies the conditions of As-sumption A, then it holds that c µ ( s ) = µ c ( s/µ ) , (16) for any s > , where c is as defined in (14), but with respect to Q , c ( s ) = Z ∞ s (1 − Q ( u )) du with c ′ ( s ) = − (1 − Q ( s )) . (17)It should be clear from the above results that this function, c µ ( s ), can bere-scaled or ”standardized” so that c µ ( s ) /µ is independent of µ . In particular,if s = b µ for some b >
0, then again by (16), c µ ( b µ ) = µ c ( b ).Next, as in (10), we define the ’Delta’ -function corresponding to the func-tion c µ ( s ) in (13) or (14), as ∆ µ ( s ) := ∂c µ ( s ) /∂µ . In the next theorem weshow that under Assumption A , ∆ µ ( · ) may be expressed in terms of thetruncated mean of X and the consequential representation of c µ ( · ). Theorem 1
Suppose that Q µ ∈ Q and thus it satisfies the conditions ofAssumption A, then for each s > , ∆ µ ( s ) = 1 µ Z ∞ s xq µ ( x ) dx. (18) Further, ∆ µ ( s ) ≡ ∆ ( s/µ ) , where ∆ ( s ) := R ∞ s uq ( u ) du ≤ for any s > .Hence, c µ ( s ) in (13) may be written as as c µ ( s ) = µ ∆ µ ( s ) − s (1 − Q µ ( s )) (19)8 roof. To prove (18), note that by Lemma 2, (17) and (13),∆ µ ( s ) = ∂∂µ [ µc ( s/µ )] = c ( s/µ ) − sµ c ′ ( s/µ ) == Z ∞ s/µ uq ( u ) du ≡ µ Z ∞ s xq µ ( x ) dx. (20)The second part follows immediately from the first part and AssumptionA and noting that ∆ ( s ) ≤ R ∞ uq ( u ) du = 1. Finally, since by (18), R ∞ s xq µ ( x ) d x = µ · ∆ µ ( s ), the main result in (19) follows directly from (13).An immediate conclusion of Theorem 1 is that if Q µ is a member of thescale-family of distribution Q (by Assumption A ) then the functions c µ ( s )can easily be evaluated by calculating first the values of c ( · ) for the ratio b = s/µ . Specifically, c µ ( s ) = µ c ( s/µ ) ≡ µ ∆ ( s/µ ) − s (1 − Q ( s/µ ) (21)The results of the Theorem, either as given in (19) or in (21), can be useddirectly for the risk-neutral valuation of European call option with a strike K and a current spot price S τ ≡ S , providing the expression for C S ( K ) asis given in (12). That is, if Q ∈ Q , then with µ = S e rt applied to (21), wehave C S ( K ) := e − rt E (( X − K ) + | S ) = e − rt c µ ( K )= S ∆ ( K/µ ) − K e − rt (1 − Q ( K/µ )) . (22)We summarize these findings in the following Corollary. Corollary 1
For any risk-neutral distribution Q µ that satisfies, in additionto (6)-(7), also the conditions of Assumption A, with µ = S e rt , so that Q µ ∈ Q (for some Q ( · ) ), we have as in (12) that C S ( K ) = S ∆ µ ( K ) − K e − rt (1 − Q µ ( K )) , where ∆ µ ( K ) := ∆ ( K/µ ) and Q µ ( K ) := Q ( K/µ ) . Hence Q µ also satisfiesHeston’s (1993) option pricing model (and solution) as given in (5). Remark 1
In the case in which the risk-neutral evaluation of the optionincludes a dividend with a rate q , then E ( S T | S ) = S e ( r − q ) t in (6), in whichcase, by applying µ = S e ( r − q ) t to (21) we obtain C S ( K ) = e − rt c µ ( K ) = S e − qt ∆ µ ( K ) − K e − rt (1 − Q µ ( K )) .
9t should be clear from (22) that since the probability distribution Q µ isassumed here to be a member of a scale family Q , its values depend on K and S only through the ratio K/S . In the Appendix, we assert in Proposition1 that any risk neutral probability distribution Q µ that satisfies the solution(5) for Heston’s option pricing model, must also be a member of a scalefamily of distributions, with a scale parameter µ = S e rt (or µ = S e ( r − q ) t , inthe case of a dividend yielding spot). This assertion follows directly from thespecific form of Heston’s RND established in the Appendix which is givenin terms of characteristic function corresponding to the term P (see (31)and the subsequent comments there). Hence combined, the statements ofCorollary 1 and Proposition 1, can be summarized in the following theorem. Theorem 2
Let Q µ ( · ) be any risk-neutral distribution that satisfies (6)-(7)with a corresponding RND q µ ( · ) and with µ = S e rt . Then Q µ ( · ) satisfiesHeston’s option pricing solution in (5) (and equivalently in (12)) if and onlyif Q µ ( · ) is member of a scale-family of distributions with a scale parameter µ . In view of Theorem 2, the quest for finding appropriate RND for Heston’ SVmodel for a particular parametrization of ϑ = ( κ, θ, η, ρ ) must be focused onlyon those members of a scale-family of distributions with a scale parameter µ = S e rt . Accordingly, we provide in this section, five particular exam-ples of well-known distributions, that satisfy the conditions of AssumptionA and hence admit, per Corollary 1, the presentation (12) for the Heston’soption pricing model in (5). These well-known distributions, namely, the
Log-Normal , the
Gamma , the
Inverse Gaussian , the
Weibull and the
InverseWeibull distributions, are re-parametrized under
Assumption A to a stan-dardized, one-parameter version having mean 1 and a second moment thatdepends on a single parameter ν > ν ≡ σ √ t , for some σ > inexpensive RNDs, easy to obtain, to calculateand calibrate as compared to the alternatives approaches available in theliterature. We note that while the gamma and the Weibull distribution wereconsidered by Savickas (2002, 2005) for deriving ‘alternative’ option pricing10ormulas, the motivation for the parametrization there was made without re-gard to the spot price dynamics (but rather for fitting kurtosis and skewness)and therefore are different.With these standardized distributions in hand and the corresponding ex-plicit expressions for Q ( · ) as obtained under Assumption A , we utilize (21)(with µ ≡
1) to first calculate in each case the expression for c ( s ) = ∆ ( s ) − s (1 − Q ( s )) , which is then used to obtain, with µ >
0, the expression for the undiscounted option price as, c µ ( s ) = µ × (cid:20) ∆ ( s/µ ) − sµ × (1 − Q ( s/µ )) (cid:21) . Finally, as in (22), the corresponding expression for the call option price isobtained as C S ( K ) = e − rt c µ ( K ) (with µ = Se rt , s = K and ν = σ √ t ). Wepoint out again, that each of these five distributions would satisfy as RND,Heston’s (1993) general solution for the valuation of a European call optionas is given in (5). We begin with the log-normal distribution which resultswith the classical Back-Scholes option pricing model (as given in (3)-(4)) . Suppose that the random variable U has the ’standard’ (one-parameter) log-normal distribution having mean E ( U ) = 1 and variance V ar ( U ) = e ν − ν >
0, so that W = log( U ) ∼ N ( − ν / , ν ). Accordingly, the pdf of U is given by; q ( u ) = 1 uν × φ (cid:18) log( u ) + ν / ν (cid:19) , u > , and its cdfQ ( u ) = P r ( U ≤ u ) = Z u q ( s ) ds = Φ (cid:18) log( u ) + ν / ν (cid:19) , ∀ u > . It is straightforward to verify that if X ≡ µU for some µ >
0, then the pdf of X is the ’scaled’ version of q , namely, q µ ( x ) = µ q ( x/µ ), so this distributionsatisfies Assumption A . 11ext, we calculate the expression of ∆ ( s ) which upon using the relation U ≡ e W , becomes∆ ( s ) := Z ∞ s uq ( u ) dx = Z ∞ log( s ) e w φ (cid:18) w + ν / ν (cid:19) dwν = Z ∞ log( s ) φ (cid:18) w − ν / ν (cid:19) dwν = 1 − Φ (cid:18) log( s ) − ν / ν (cid:19) . Hence, for the ’standardized’ model we have that c ( s ) =∆ ( s ) − s (1 − Q ( s )) ≡ (cid:20) − Φ (cid:18) log( s ) − ν / ν (cid:19)(cid:21) − s (cid:20) − Φ (cid:18) log( s ) + ν / ν (cid:19)(cid:21) Accordingly, by Lemma 2 and (21), c µ ( s ) ≡ µ × c ( s/µ ) and we thereforeimmediately obtain the following expression for c µ ( s ) as, c µ ( s ) = µ × (cid:20) ∆ ( s/µ ) − sµ × (1 − Q ( s/µ )) (cid:21) = µ × (cid:20) − Φ (cid:18) log( s/µ ) − ν / ν (cid:19)(cid:21) − s × (cid:20) − Φ (cid:18) log( s/µ ) + ν / ν (cid:19)(cid:21) ≡ µ × Φ (cid:18) log( µ/s ) + ν / ν (cid:19) − s × Φ (cid:18) log( µ/s ) − ν / ν (cid:19) , where the last equality utilized the symmetry of the normal distribution.Finally, to calculate under the log-normal RND the price of a call optionat a strike K when the current price of the spot is S , we utilize the aboveexpression, c µ ( s ), with µ ≡ S e rt , s ≡ K and ν ≡ σ √ t to obtain, C S ( K ) = e − rt c µ ( K ), which matches exactly the Black-Scholes call option price as isgiven in (3)-(4). We begin with some standard notations. We write W ∼ G ( α, λ ) to indi-cate that the random variable W has the gamma distribution with a scaleparameter λ > α >
0, in which case we write g ( · ; α, λ ) and G ( · ; α, λ ) for the corresponding pdf and cdf of W , respectively.Recall that E ( W ) = α/λ and V ar ( W ) = α/λ . Additionally, we denote12y Γ( α ) := R ∞ y α − e − y dy the gamma function whose incomplete version isΓ( ξ, α ) := R ξ y α − e − y dy , is defined for any ξ > U has the ’standard’ (one-parameter)Gamma distribution having mean E ( U ) = 1 and variance V ar ( U ) = ν , forsome ν >
0, so that U ∼ G ( a, a ) where we substituted a ≡ /ν . Accord-ingly, the pdf of U is given by q ( u ) := g ( u ; a, a ) = a ( au ) a − e − au Γ( a ) , u > . and its cdf , by Q ( u ) = P r ( U ≤ u ) := G ( u ; a, a ) = Γ( au, a )Γ( a ) , for any u >
0. It is straightforward to verify that if X ≡ µU for some µ > pdf , q µ ( · ) of X is the ’scaled’ version of q ( · ), and that AssumptionA holds in this case too. Next, we calculate the expression for ∆ ( s ),∆ ( s ) = Z ∞ s uq ( u ) du = Z ∞ s u g ( u ; a, a ) du = Z ∞ s a ( au ) a e − au a Γ( a ) du = 1 − Γ( as, a + 1)Γ( a + 1) ≡ − G ( s ; a + 1 , a ) . Accordingly, we obtain for the ’standardized’ Gamma model that c ( s ) =∆ ( s ) − s (1 − Q ( s )) ≡ = (cid:20) − Γ( as, a + 1)Γ( a + 1) (cid:21) − s (cid:20) − Γ( as, a )Γ( a ) (cid:21) = [1 − G ( s ; a + 1 , a )] − s [1 − G ( s ; a, a )]Again, by Lemma 2 and (21), c µ ( s ) ≡ µ × c ( s/µ ) and we therefore imme-diately obtain the following expression for c µ ( s ) in this case of the Gammamodel as, c µ ( s ) = µ × (cid:20) ∆ ( s/µ ) − sµ × (1 − Q ( s/µ )) (cid:21) = µ × [1 − G ( s/µ ; a + 1 , a )] − s × [1 − G ( s/µ ; a, a )] (23)Finally, to calculate under this Gamma RND the price of a call option at astrike K when the current price of the spot is S , we will utilize (23) with µ ≡ S e rt , s ≡ K and ν ≡ σ √ t (so that a ≡ /σ t ) to obtain, C S ( K ) = e − rt c µ ( K ). 13 .3 The Inverse Gaussian RND Using standard notation we write W ∼ IN ( α, λ ) to indicate that the ran-dom variable W has the Inverse Gaussian distribution with mean E ( W ) = α and V ar ( W ) = α /λ . Now suppose that a random variable U has the ’stan-dard’ (one-parameter) Inverse Gaussian distribution having mean E ( U ) = 1and variance V ar ( U ) = ν , for some ν >
0, so that U ∼ IN (1 , /ν ).Accordingly, the pdf and cdf of U are given by; q ( u ) = 1 νu / × φ (cid:18) u − ν √ u (cid:19) , u > , and Q ( u ) = Φ (cid:18) u − ν √ u (cid:19) + e /ν × Φ (cid:18) − u + 1 ν √ u (cid:19) ∀ u > . Again, one can verify that if X ≡ µU for some µ >
0, then the pdf , q µ ( · ) of X is the ’scaled’ version of q ( · ) above, so that Assumption A holds in thiscase too.In the case of this distribution, the values of of ∆ ( s ) = R ∞ s uq ( u ) du must be evaluated numerically which, when combined with the expression of Q ( s ) given above, provide the values of c µ ( s ) = µ × (cid:20) ∆ ( s/µ ) − sµ × (1 − Q ( s/µ )) (cid:21) , for any µ >
0. Here again, the corresponding values of the call option C S ( K )may be obtained,exactly along the same lines as in the previous examples,with µ ≡ S e rt , s ≡ K and ν = σ √ t . Using standard notation we write W ∼ W ( ξ, λ ) to indicate that the randomvariable W has the Weibull distribution with cdf and pdf which are of theform, F W ( w ) = 1 − e − ( w/λ ) ξ , and f W ( w ) = ξλ ( wλ ) ξ e − ( w/λ ) ξ , w > , respectively, where λ > ξ > W are given by, E ( W ) = λh ( ξ ) and V ar ( W ) = λ ( h ( ξ ) − h ( ξ ) ) , h j ( ξ ) := Γ(1 + j/ξ ) , j = 1 , . . . , . Now suppose that a random vari-able U has the ’standard’ (one-parameter) Weibull distribution having mean E ( U ) = 1 and variance V ar ( U ) = ν , for some ν >
0. That is, for a given ν >
0, we let ξ ∗ ≡ ξ ( ν ) be the (unique) solution of the equation h ( ξ ) h ( ξ ) = 1 + ν , (24)in which case, h ∗ j ≡ h j ( ξ ∗ ) , j = 1 , λ ∗ ≡ /h ∗ and U ∼ W ( ξ ∗ , λ ∗ ).Accordingly, the pdf and cdf of U are given by, Q ( u ) = 1 − e − ( u/λ ∗ ) ξ ∗ , and q ( u ) = ξ ∗ λ ∗ ( uλ ∗ ) ξ ∗ e − ( u/λ ∗ ) ξ ∗ , u > , Again, it can be easily verified that if X ≡ µU for some µ >
0, then the pdf , q µ ( · ) of X is the ’scaled’ version of q ( · ) above, so that Assumption A holdsin this case too. For this RND, the values of of ∆ ( s ) can be obtained in aclosed form as∆ ( s ) = Z ∞ s uq ( u ) du = 1 − Γ(( s/λ ∗ ) ξ ∗ ; 1 + 1 /ξ ∗ )Γ(1 + 1 /ξ ∗ ) , which, together with the expression of Q ( · ) given above, provide the valuesof c µ ( s ) = µ × (cid:20) ∆ ( s/µ ) − sµ × (1 − Q ( s/µ )) (cid:21) , for any µ >
0. Here again, the corresponding values of the call option C S ( K )may be obtained,exactly along the same lines as in the previous examples,with µ ≡ S e rt , s ≡ K and ν = σ √ t . In similarilty to the above example, we write W ∼ IW ( ξ, α ) to indicate thatthe random variable W has the Inverse Weibull distribution (see for example,de Gusm˜ao at. el. (2009) ) with cdf and pdf which are of the form, F W ( w ) = e − ( α/w ) ξ , and f W ( w ) = ξα ( αw ) ξ +1 e − ( α/w ) ξ , w > , (25)respectively, where α > ξ > W are given by, E ( W ) = α ˜ h ( ξ ) and V ar ( W ) = α (˜ h ( ξ ) − ˜ h ( ξ )) , h j ( ξ ) ≡ h j ( − ξ ) = Γ(1 − j/ξ ) , j = 1 , . . . , . Here too, we let U have the’standard’ (one-parameter) Inverse Weibull distribution with mean E ( U ) = 1and variance V ar ( U ) = ν , for some ν >
0. That is, for a given ν >
0, welet ξ ∗ ≡ ξ ( ν ) be the (unique) solution of the equation˜ h ( ξ )˜ h ( ξ ) = 1 + ν , (26)in which case, U ∼ IW ( ξ ∗ , α ∗ ) with α ∗ = 1 / ˜ h ( ξ ∗ ). Accordingly, the pdf, q ( u ), and cdf, Q ( u ), of U are as given in (25), but with ξ ∗ and α ∗ .Hence, we may proceed exactly along the lines of the previous example tocalculate c ( s ), and c µ ( s ) and hence, C S ( K ) = e − rt c µ ( K ) (with µ = Se rt , s = K and ν = σ √ t ). As can be seen, the distribution in each of these five examples satisfies theconditions of
Assumption A and hence by Corollary 1 could potential serve asRND for Heston’s SV model (1). These distributions are defined by a singleparameter, namely ν ≡ σ √ t , that affects their features, such as kurtosis and skewness , and hence their suitability as RND for particular scenariosof the SV model (1)– as are determined by the structural model parameter ϑ = ( κ, θ, η, ρ ) (more on this point in the next section). However, for sake ofcompleteness and for future reference we provide in Table 1 the expressionfor the kurtosis and skewness for these five distributions.It is interesting to note that, in the relevant parametric domain, all butthe Weibull example, have positive Skewness measure. It can be numericallyverified that in the Weibull case, γ s ( ξ ( ν )) changes it’s sign and is negativeonce ν < . γ k ( ξ ( ν )) < . < ν < . ν > . Having introduced in the previous section several examples of distributionsthat serve as possible RND for Heston’s (1993) option valuation (5) under16able 1: The skewness and excess kurtosis measures of the RND Examples3.1-3.5 as functions of the single parameter ν ≡ σ √ t .Distribution E ( U ) V ar ( U ) Skew Exc.Kurtosis G (1 /ν , /ν ) 1 ν ν ν IN (1 , /ν ) 1 ν ν ν W ( ξ, /g ( ξ )) ∗ ν γ s ( ξ ) γ k ( ξ ) − IW ( ξ, /h ( ξ )) ∗∗ ν γ s ( − ξ ) γ k ( − ξ ) − LN ( − ν / , ν ) 1 e ν − e ν + 2) √ e ν − e ν + 2 e ν + 3 e ν − ∗ Here ξ ≡ ξ ( ν ) solves equation (24); ∗∗ Here ξ ≡ ξ ( ν ) solves equation (26) and it is assumed that ν is such that ξ ( ν ) > h j ( ξ ) = Γ(1 + j/ξ ) , j =1 , . . . , γ s ( ξ ) = h ( ξ ) − h ( ξ ) h ( ξ ) + 2 h ( ξ ) (cid:2) h ( ξ ) − h ( ξ ) (cid:3) / and γ k ( ξ ) = h ( ξ ) − h ( ξ ) h ( ξ ) + 6 h ( ξ ) h ( ξ ) − h ( ξ ) (cid:2) h ( ξ ) − h ( ξ ) (cid:3) . the stochastic volatility model (1), we dedicate this section to illustrationof their applicability and their relative comparison. In the Appendix, weprovide the closed-form expressions for Heston’s P and P as are given interms of their characteristic functions (see Heston (1993)). These termsenable us to compute, for given S τ = S , V τ = V and r , and for each choiceof ϑ = ( κ, θ, η, ρ ), Heston’s call price C S ( K ) as in (5) as well as Heston’s RNDas derived from the characteristic functions of P (see Appendix for details).Features of this distribution, such as Kurtosis and
Skewness as are largelydetermined by η and ρ , respectively (see Bakshi, Cao and Chen (1997) ),would serve as guide for matching a particular proposed RND from amongour five examples (see also Table 1). For instance, in cases which admit17n RND with a distinct negative skew, the Weibull distribution could beconsidered, whereas, in those cases with a distinct positive skew, the
InverseWeibull or the other distributions discussed in Section 3 could be considered.Additionally, we may simulate observations on ( S T , V T ) from a discretizedversion of Heston’s stochastic volatility process (1) to obtain the simulatedrendition of the marginal distribution of S T . In light of the scaling propertyof the RND, we present, whenever convenient, the results in terms of therescaled spot priced, S ∗ = S T /µ , where µ = Se rt (see Corollary 4). In thesimulations, we employed either the (reflective version of) Milstein’s (1975)discretization scheme or Alfonsi’s (2010) implicit discretization scheme alldepending on whether the so-called Feller condition, ζ := κθ/η >
1, holdsor not (see Gatheral (2006) for a discussion). We note that ζ is intimatelyrelated to the conditional distribution of V T implied by the SV model (1), (seeProposition 2 of Andersen (2008) for details). In all cases we also includeda comparison of the Monte-Carlo distribution of S ∗ to the actual Heston’sRND as was numerically calculated using (31) under the ‘calibrated’ valuesof ϑ = ( κ, θ, η, ρ ).In the first two examples, we use values of the structural parameters, ϑ = ( κ, θ, η, ρ ), already calibrated to market data on as can be found fromBakshi, Cao and Chen (1997) (on the S&P 500) and Mr´azek and Posp´ıˇsil(2017) (on the ODAX). These two examples which involve market data ontraded Indexes are are used to illustrate the applicability of the Weibulldistribution to situations in which the RND is negatively skewed and largelyleptokurtic one. Other similar illustrations using calibrated parameter values.such as from Lemaire , Montes and Pag`es (2020) ( on the E URO S TOXX
Example: S&P 500:
Bakshi, Cao and Chen (1997) presented an exten-sive market data study for comparing several competing stochastic volatilitymodels, including that of Heston’s (1993). The data used covered optionsand spot prices for the S&P 500 Index starting from June 1, 1988 throughMay 31, 1991. From Table III there we find that in addition to r = 0 . .7 0.8 0.9 1.0 1.1 1.2 . . . . . Joint distribution of (S*,V)
S*−price V − v o l a t ili t y Marginal RND of S*
S*−Price D en s i t y Heston’s RNDWeibull RND
Figure 1:
Simulated joint ( S ∗ , V ) -distribution, and the Heston’s and WeibullRNDs for the S&P 500 data based on the calibrated parameter ˆ ϑ =(1 . , . , . , − . as provided by Bakshi, Cao and Chen (1997). Heston’s SV model, isˆ ϑ = (1 . , (0 . / . , . , − . . In this case, ˆ ζ = 0 . <
1, hence we used Milstein (reflective) scheme toobtain, for a short contract duration with t = 56 /
365 = 0 .
153 year, a Monte-Carlo sample of M = 10 ,
000 simulated pairs ( S ∗ , V ) with (standardized) spotprice and volatility, according to the SV model (1). Their joint distribution ispresented in Figure 1a, where we have superimposed the matching 16%, 50%,68% , 95% and 99.5% contour lines. In Figure 1b we present the histogramof the simulated marginal distribution of the spot price S ∗ . The mean and19tandard deviation of these M simulated spot price values are ¯ S ∗ = 0 . σ √ t = 0 . ϑ ) Heston’s RND as was computed directly by using (31). As isexpected in the case of (risk-neutrality) modeling the spot prices of an Index,the implied RND is negatively skewed ( sk = − . ν to the ‘observed’ value of ˆ σ √ t and used it to obtain the numericalsolution of equation (24) as ˆ ξ = 17 . h ( ˆ ξ ) = 0 . W ( ˆ ξ, /h ( ˆ ξ ))RND. As can be seen, the two RND curves are almost indistinguishable.To further illustrate the applicability of the proposed RND to situationswith distinct positive Skewness , we considered again Bakshi, Cao and Chen(1997) calibrated parametrization used for Figure 1, but now with a ( hypo-thetically ) positive correlation, so that ˆ ϑ = (1 . , (0 . / . , . , . . The simulated Monte-Carlo distributions (joint and marginal) are presentedin Figure 2, exhibiting the distinct positive
Skewness of Heston’s RND (cal-culated from (31)). This suggested a comparison to the Inverse Weibull dis-tribution discussed in Section 3.5. The mean and standard deviation of these M simulated spot price values are ¯ S ∗ = 0 . σ √ t = 0 . ν was match to ˆ σ √ t to obtain the solution of equation(26) as ˆ ξ = 18 . h ( ˆ ξ ) = 1 . IW ( ˆ ξ, / ˜ h ( ˆ ξ )) RND to Figure 2b, il-lustrating the extent of the agreement between Heston’s (implied) RND andthe Inverse Weibull distribution in this case (with a distinct positives skew). Example: ODAX:
Mr´azek and Posp´ıˇsil (2017) studies various optimizationtechniques for calibrating and simulating the the Heston model. For demon-strate their results, they used the ODAX option Index with 5 blended matu-rities of three and six months and with 107 strikes, as were recored on March19, 2013. The calibrated results of the structural parameter ϑ = ( κ, θ, η, ρ )are provided in Table 4 there,ˆ ϑ = (1 . , . , . , − . . with r = 0 . S = 7962 .
31 and with V = 0 . t = 64 / M =10 ,
000 pairs of ( S ∗ , V ) to obtain from the discretized process, the renditions20 .6 0.8 1.0 1.2 1.4 1.6 . . . . Simulated Joint distribution of (S*, V)
S*−price V − v o l a t ili t y Marginal RND of S*
S*−Price D en s i t y Heston’s RNDInv Weibull RND
Figure 2:
Illustrating the Inverse Weibull pdf for the
S&P 500 data of Bakshi,Cao and Chen (1997) calibrated parameter but with hypothetically positivecorrelation ( ρ = 0 . ) resulting with a positively skewed RND. of their joint distribution as well as the marginal distribution of S ∗ . Theseare presented in Figure 3. The mean and standard deviation of these M sim-ulated spot price values are ¯ S ∗ = 0 . σ √ t = 0 . negatively skewed ( sk = − . ν was match to ˆ σ √ t to obtain the solution of equation(26) as ˆ ξ = 19 . h ( ˆ ξ ) = 0 . W ( ˆ ξ, /h ( ˆ ξ )), is also displayed in the figure, indicating21 .6 0.8 1.0 1.2 1.4 1.6 . . . . . . . Simulated Joint distribution of (S*, V)
S*−price V − v o l a t ili t y Marginal RND of S*
S*−Price D en s i t y Heston’s RNDWeibull RND
Figure 3:
Simulated joint ( S ∗ , V ) -distribution, and the (conditional) Heston’sand Weibull RND for the ODAX data based on the calibrated parameter ˆ ϑ =(1 . , . , . , − . as provided by Mr´azek and Posp´ıˇsil (2017). the excellent agreement, in this case two, between the RNDs. Example: AMD:
This example is based on real and current option datathat we retrieved from
Yahoo Finance as of the closing of trading on Decem-ber 31, 2020. The closing price of his stock, on that day, was 91 .
71 and itpays no dividend, so that q = 0 to add to the prevailing (risk-free) interestrate of r = 0 . t = 47 /
365 and some N = 39 strikes, K , . . . , K with corresponding call option (market) prices C , . . . , C areavailable (we actually recorded the option prices as the average between thebid and ask). As standard measure of the goodness-of-fit between the model-calculated option price C Model ( K i ) and the option market price C i , we usedthe Mean Squared Error , MSE,
M SE ( M odel ) = 1 N N X i =1 ( C Model ( K i ) − C i ) To calibrate the Heston SV model, we used the optim( · ) function of R, tominimize M SE ( Heston ) over the model’s parameter , ϑ = ( κ, θ, η, ρ ) withthe initial values of (2 , . , . ,
0) and with V = 0 .
25. The results of thecalibrated values areˆ ϑ = (1 . , . , . , . . This calibrated parameter, ˆ ϑ , was then used to calculate, using Heston’scharacteristic function, the option prices according to Heston’s SV model (5).These values are displayed in Table 2, along with the actual market prices.Next, we obtained, as in the previous examples, a Monte-Carlo sample of( S ∗ , V ) whose results are displayed in Figure 4. The mean and standarddeviation of these simulated stock prices are ¯ S ∗ = 1 . sd ( S ∗ ) =0 . positively skewed ( sk = 1 . ν = σ √ t , defining these distribution directly,using the ‘standard’ Black-Scholes implied volatility, namely ˆ ν = IV BS √ t .This entails using the the optim( · ) function again to minimize the M SE ( BS )with respect to the single parameter σ . This standard estimation procedureyielded IV BS = 0 . ν = 0 . Gamma , Inverse Gaussian and the
Log-Normal
RND, as in Table 1. The extent of their agreement with the Heston’simplied RND is self-evident. To further demonstrate that very point, wecalculated the option prices under each one of these modeled RND, andcalculated the corresponding
M SE ( M odel ). The results of this comparisonare presented, side-by-side in Table 2. In Figure 5 we display the option pricecurve for each of these pricing models– they are ‘virtually’ almost identicalin this example. 23 .5 1.0 1.5 2.0 2.5 3.0
Simulated Joint distribution of (S*, V)
S*−price V − v o l a t ili t y Marginal RND of S*
S*−Price D en s i t y . . . . . . . Heston’s RNDGamma RNDInvGauss RNDLogNorm RND
Figure 4:
Simulated joint ( S ∗ , V ) -distribution, and Heston’s, Gamma, Log-Normaland Inv. Gaussian RNDs calculated based on the AMD data of Table 2.
Heston (1993) provided (semi) closed form expressions to the probabilities P and P that comprise the solution C S ( K ) in (5) to the option valuation underthe stochastic volatility model (1). Starting from a ‘guess’ of the Black-Sholesstyle solution, C = SP − Ke − rt P , (27)he has shown that with x := log ( S ), this solution must satisfy the SDEresulting from the SV model in (1), ∂P j ∂t = 12 v ∂ P j ∂x + ρηv ∂ P j ∂x∂v + 12 η v ∂ P j ∂v + ( r + u j v ) ∂P j ∂x + ( a − b j v ) ∂P j ∂v , Call Option Prices for AMD with 47 Days to Maturity
DataStrike C a ll O p t i on P r i c e MarketHestonGammaInvGauss RNDBlack−Scholes
Figure 5:
Comparing the the option prices obtained by the Heston, Gamma, In-verse Gaussian and the Black-Scholes Models for the February 19, 2021 OptionSeries of
AMD (with 47 DTE) as was quoted on the closing December 31, 2020business day for j = 1 ,
2, where u = 1 / , u = − / , b = κ − ρη, b = κ and a = κθ .These closed form expressions are given by P j = 12 + 1 π Z ∞ R e (cid:20) e − iωk ψ j ( ω, t, v, x ) iω (cid:21) dω, (28)where k := log ( K ) and ψ j ( · ) is the characteristics function ψ j ( ω, t, v, x ) := Z ∞−∞ e iωs p j ( s ) ds ≡ e B j ( ω,t )+ D j ( ω,t ) v + iωx + iω rt , (29)where p j ( · ) is the pdf of s T = log( S T ) corresponding to the probability P j , j = 1 , B j ( ω, t ) = κθη { ( b j + d j − iωρη ) t − − g j e d j t − g j ) } j ( ω, t ) = b j + d j − iωρηη ( 1 − e d j t − g j e d j t ) g j = b j − iω ρη + d j b j − iωρη − d j d j = q ( iωρη − b j ) − η (2 iωu j − ω )We point out that d j above is taken to be the positive root of the Riccatiequation involving D j . However using instead the negative root, namely d ′ j = − d j , was shown to provide an equivalent, but yet more stable solutionfor ψ j above -see Albrecher, Mayer, Schoutens, and Tistaer, (2007) for moredetails on this so-called “Heston Trap”. In either case, efficient numericalroutines such as the cfHeston and callHestoncf functions of the NMOFpackage of R, are readily available to accurately compute the values of ψ j and hence of P j and the call option values, for given t, s and v and any choiceof ϑ = ( κ, θ, η, ρ ).Now, having established (29), the standard application of the Fouriertransform provides (see for example Schmelzle (2010)) that the pdf p j ( · ) of s T = log( S T ), can be obtained, for s ∈ R , as p j ( s ) = 1 π Z ∞ R e (cid:2) e − iωs ψ j ( ω, t, v, x ) (cid:3) dω. (30)Hence, it follows immediately that the pdf ˜ p j ( · ) of S T is given by˜ p j ( u ) = 1 u × p j (log( u )) ≡ π Z ∞ R e (cid:20) e − iω log( u ) ψ j ( ω, t, v, x ) u (cid:21) dω, u > . Further, since the characteristic functions ψ j in (29) are affine in x + rt =log( S ) + rt ≡ log( µ ), where as in Corollary 1, µ = Se rt , we may rewrite ˜ p j ( u )above as ˜ p j ( u ) = 1 µπ Z ∞ R e " e − iω log( u/µ ) ˜ ψ j ( ω, t, v ) u/µ dω, (31)where log( ˜ ψ j ( ω, t, v )) := log( ψ j ( ω, t, v, x ) − iωx − iω rt.
26e point out that in light of (8) that ˜ p ( · ) in (31) is the RND (under Q ) forthe Heston’s (1993) model and can similarly be easily evaluated numericallyalong-side of evaluating P . Indeed we have, P ≡ Z ∞ K ˜ p ( u ) du = Q ( S T > K ) . It should be noticed from expression (31) that any RND, ˜ p ( · ) of the HestonModel, and the corresponding risk neutral distribution Q µ ( · ) of S T , consti-tutes a scale-family of distributions in µ = Se rt , so that it satisfies the termsof Assumption A . This assertion is summarized in Proposition 1 below.
Proposition 1
Let q µ ( · ) be any RND with a corresponding risk neutral dis-tribution Q µ ( · ) that satisfies Heston’s solution in (5), with µ = S e rt then q µ ( · ) is of the form given in (31) and therefore Q µ ( · ) must be a member of ascale-family of distributions in µ . The result stated in the next claim is known, but the details are instructiveto proving (11).
Claim 1
Let ∆( K ) = ∂C/∂S as in (10), then for the Heston solution (5)(or (27)) with P j , j = 1 , as are given in (28), we have ∆( K ) = P . Proof.
By (27) we have ∂C∂S = P + S ∂P ∂S − Ke − rt ∂P ∂S . (32)Since we have x = log( S ), it follows from (28) that, ∂P j ∂S = 1 π Z ∞ R e (cid:20) e − iωk ψ j ( ω, t, v, x ) S (cid:21) dω. Hence, by (30) we have with µ = Se rt , S ∂P j ∂S ≡ p ( k ) , and Ke − rt ∂P ∂S ≡ Kµ p ( k ) , which by (9) implies that S∂P /∂S − Ke − rt ∂P /∂S = 0 in (32).27 eferences [1] Albrecher, H., Mayer, P., Schoutens, W. and Tistaer, J., (2007). TheLittle Heston Trap. The Wilmott Magazine , 83–92.[2] Alfonsi A., (2010) High order discretization schemes for the CIR pro-cess: Application to affine term structure and Heston models.
Math.Comp.
J. Comput. Finance , 11(3), 2008, 1–42.[4] Black F., and Scholes M., (1973). The pricing of options and corporateliabilities.
The Journal of Political Economy , 637-654.[5] Bakshi G., Cao C. and Chen Z., (1997) Empirical Performance ofAlternative Option Pricing Models.
The Journal of Finance , Vol. LII,No. 5, 2003-2049.[6] Cox J.C. and Ross S., (1976) The valuation of options for alternativestochastic processes.
J. Fin. Econ.
Ann. Math.
Volatility and Time Series Economet-rics: Essays in Honor of Robert F. Engle, (eds. Tim Bollerslev, JeffreyRussell and Mark Watson)
Oxford University Press, Oxford, U.K.[9] Figlewski S., (2018). Risk Neutral Densities: A Review. Available at
SSRN http://ssrn.com/abstract=3120028 .[10] Gatheral, J., (2006).
The Volatility Surface,
John Wiley and Sons, NJ.[11] Grith M. and Kr¨atschmer V., (2012) “Parametric Estimation of RiskNeutral Density Functions”, in:
Ed: Duan JC., H¨ardle W., GentleJ. (eds) Handbook of Computational Finance
Springer Handbooks ofComputational Statistics. Springer, Berlin, Heidelberg.2812] Heston S.L., (1993) A closed-form solution for options with stochasticvolatility with applications to bond and currency options.
Rev. Financ.Stud.
Option-Implied Risk-Neutral Distributionsand Risk Aversion
Research Foundation of AIMR, Charlotteville, NC[14] Jackwerth, J. C. and Rubinstein, M., (1996). Recovering ProbabilityDistributions from Option Prices.
The Journal of Finance,
51, no. 5:1611-631.[15] Jiang, L. (2005).
Mathematical Modeling and Methods of Option Pric-ing,
Translated from Chinese by Li. C, World Scientific, Singapore.[16] Lemaire, V., Montes, T. and Pag`es, G., (2020). Stationary Hestonmodel: Calibration and Pricing of exotics using Product RecursiveQuantization. Available at arXiv [q-fin.MF]: arXiv:2001.03101v2 .[17] Merton, R., (1973). Theory of rational option pricing.
The Bell Journalof Economics and Management Science , 141-183.[18] Mil’shtein, G. N. (1975). Approximate Integration of Stochastic Dif-ferential Equations.
Theory of Probability & Its Applications , 19 (3),557–562.[19] Mizrach, B. (2010). “Estimating Implied Probabilities from OptionPrices and the Underlying” in
Handbook of Quantitative Finance andRisk Management (C.-F. Lee A. Lee and J. Lee (eds.)) . Springer Sci-ence Business Media.[20] Mr´azek, M. and Posp´ıˇsil, J. (2017). Calibration and sim-ulation of Heston model.
Open Mathematics— Vol. https://doi.org/10.1515/math-2017-0058 .[21] R Core Team, (2017).
R: A Language and Environment for StatisticalComputing . Vienna, Austria, [22] Stein J. and Stein E., (1991). Stock price distributions with stochasticvolatility: An analytic approach.
Rev. Financ. Stud.
Technical Report , Available on line at https://pfadintegral.com/articles/option-pricing-formulae-using-fourier-transform/ [24] Savickas, R., (2002). A simple option formula.
The Financial Review ,Vol 37, 207-226.[25] Savickas, R., (2005). Evidence on delta hedging and implied volatilitiesfor the Black-Scholes, gamma, and Weibull option pricing models.
TheJournal of Financial Research , Vol 18:2, 299-317.[26] Wiggins, B. J., (1987). Option values under stochastic volatility: The-ory and empirical estimates.
Journal of Financial Economics , Vol.19(2), 351-372 30able 2:
Comparing the the option prices obtained by the Heston, Gamma, InverseGaussian and the Black-Scholes Models for the February 19, 2021 Option Series of
AMD (with 47 DTE) as was quoted on the closing December 31, 2020 business day.
MSE 0.004410 0.032725 0.018126 0.016748
Strike Market Price Heston Gamma InvGaussan Black.Scholes40.0 51.775 51.720 51.738 51.737 51.73742.5 49.275 49.222 49.239 49.238 49.23845.0 46.775 46.726 46.741 46.739 46.73947.5 44.200 44.231 44.244 44.240 44.24050.0 41.825 41.741 41.751 41.742 41.74355.0 36.875 36.779 36.783 36.758 36.76260.0 31.950 31.869 31.871 31.816 31.82465.0 27.150 27.058 27.073 26.977 26.99370.0 22.450 22.423 22.476 22.339 22.36572.5 20.200 20.203 20.285 20.132 20.16475.0 17.975 18.070 18.184 18.022 18.05977.5 16.025 16.039 16.186 16.022 16.06480.0 14.050 14.127 14.302 14.145 14.19082.5 12.250 12.349 12.544 12.399 12.44985.0 10.800 10.719 10.917 10.793 10.84687.5 9.275 9.242 9.428 9.330 9.38590.0 7.925 7.923 8.077 8.010 8.06692.5 6.850 6.760 6.866 6.830 6.88895.0 5.800 5.746 5.790 5.786 5.84397.5 4.925 4.871 4.844 4.870 4.927100.0 4.100 4.120 4.021 4.073 4.129105.0 2.835 2.939 2.706 2.799 2.852110.0 2.065 2.096 1.766 1.880 1.928115.0 1.410 1.498 1.119 1.237 1.278120.0 1.025 1.075 0.688 0.798 0.833125.0 0.765 0.776 0.412 0.506 0.535130.0 0.605 0.563 0.240 0.315 0.338135.0 0.550 0.411 0.136 0.194 0.211140.0 0.205 0.302 0.076 0.117 0.131145.0 0.265 0.224 0.041 0.070 0.080150.0 0.215 0.167 0.022 0.042 0.049155.0 0.185 0.125 0.011 0.024 0.029160.0 0.110 0.094 0.006 0.014 0.017165.0 0.135 0.072 0.003 0.008 0.010170.0 0.120 0.055 0.001 0.005 0.006175.0 0.135 0.042 0.001 0.003 0.004180.0 0.095 0.032 0.000 0.001 0.002185.0 0.070 0.025 0.000 0.001 0.001190.0 0.040 0.020 0.000 0.000 0.001
Source: Yahoo Financial: