Option Pricing under Double Stochastic Volatility Model with Stochastic Interest Rates and Double Exponential Jumps with Stochastic Intensity

We present option pricing under the double stochastic volatility model with stochastic interest rates and double exponential jumps with stochastic intensity in this article. We make two contributions based on the existing literature. First, we add double stochastic volatility to the option pricingmodel combining stochastic interest rates and jumps with stochastic intensity, and we are the first to fill this gap. Second, the stochastic interest rate process is presented in the Hull–White model. Some authors have concentrated on hybrid models based on various asset classes in recent years. )erefore, we build a multifactor model with the term structure of stochastic interest rates. We also approximated the pricing formula for European call options by applying the COS method and fast Fourier transform (FFT). Numerical results display that FFTand the COS method are much faster than the numerical integration approach used for obtaining the semi-closed form prices. )e COS method shows higher accuracy, efficiency, and stability than FFT. )erefore, we use the COS method to investigate the impact of the parameters in the stochastic jump intensity process and the existence of the process on the call option prices. We also use it to examine the impact of the parameters in the interest rate process on the call option prices.


Introduction
An abundance of empirical studies show the existence of the asymmetric leptokurtic features and the volatility smile after Black and Scholes [1] did some experimental and pioneering work in European option pricing. Allowing the volatility and the interest rate to change, allowing for the existence of jumps, and the change of the jump intensity over time represent reasonable dynamics of the asset returns.
Stochastic volatility models have been playing a significant part in European option modelling since volatility should be a random variable based on extensive empirical studies. Some authors proposed several representative stochastic volatility models [2][3][4][5]. Heston [6] specified the variance (the square of volatility) with a Cox-Ingersoll-Ross (CIR) process which is more proper for application than other models. He contributed to the existing literature mainly by modelling the variance with the CIR process which displays mean-reverting and nonnegative properties, deriving the formulae for the characteristic functions with PDE approach and applying the Fourier transform for obtaining the closed-form valuation formula for European options since the density function is expressed with the characteristic function using inverse Fourier transform. Schöbel and Zhu [7] developed a model with stochastic volatility which is also proper for application. e volatility follows an Ornstein-Uhlenbeck process in their model with the asset returns, and its volatility being correlated with each other. ey contributed to the existing literature mainly by using the expectation approach for deriving the formulae for the characteristic functions instead of PDE approach. Lewis [8] developed a stochastic variance model that the variance follows a 3/2 nonaffine stochastic process because the option prices under this model are local martingales instead of martingales, and 3/2 process is not stationary and it is improper for application. Grasselli [9] proposed that the variance follows a 4/2 process which is a mix of the 1/2 and the 3/2 terms. Since they quoted that the 4/2 model shares the same properties as the 3/2 model, so it is not proper for practical application as well. Christoffersen et al. [10] modelled the variance with a two-factor mean-reverting square root process that provides more flexibility than the Heston model. eir empirical study shows that their model works better than the one-factor model. Based on the existing forms of processes used to describe the dynamics of the volatility, we decide to study option pricing under twofactor stochastic volatility in this article since it is more applicable for practical application. e interest rate is also time varying in the real economy. Meanwhile, stochastic interest rate models have a longer history than stochastic volatility models. ey are initially used to study the zero-coupon bond and interest rate options and derive the formulae for them for application. ere are four typical stochastic interest rate models. Vasicek [11] modelled the interest rate with an Ornstein-Uhlenbeck process for describing the change of it. Cox et al. [12] modelled the interest rate with mean-reverting square root process; henceforth, it was applied for reference to develop the stochastic volatility model. ey contributed to the existing literature by expanding Vasicek's model that they added a term to the diffusion coefficient, and it maintains the mean-reverting and nonnegative properties that make it more proper for application than Vasicek's model. Since then, this model has been a benchmark to specify the dynamic change of the variance and the interest rate for decades. Longstaff [13] proposed a mean-reverting double square root model, and compared to the one square root model, it has some special features that it requires the parameters in the model to satisfy a specific condition. Since the dynamic change of the variance and the interest rate shares some common features, Zhu [14] expanded Longstaff's model to the stochastic variance model that they modelled the variance with double square root process. Hull and White [15] proposed a special process to specify the change of the interest rate with all the parameters in the model being time varying. Since this model cannot capture the market shapes very well in reality, they noted that the calibration of this model needs to be carefully dealt with. To make it perform better for practical application, Hull and White [16] improved their model to be a more reliable and applicable one that only one parameter in the process is time varying. It can be transformed to another form which is generally called the Hull-White interest rate process [17].
Because of the contribution of the authors who developed the stochastic interest rate models, some authors began to introduce the stochastic interest rate processes into European option pricing models to make them more reasonable for practical application. eir empirical work supports the significant improvement of stochastic interest rates [18,19]. Meanwhile, several authors focused on option pricing combining stochastic volatility and interest rates to build hybrid models. CIR and Hull-White models are generally used to display the dynamics of the interest rate. Grzelak and Oosterlee [20] proposed an option pricing model with stochastic volatility and stochastic interest rates. e interest rate follows Hull-White and CIR processes in their model. Grzelak et al. [21] proposed option pricing combining stochastic volatility with Schöbel-Zhu and CIR processes and stochastic interest rate with Hull-White process, and they used some techniques for obtaining the formula for discounted characteristic function.
Jumps are used to describe the discontinuous behavior of the asset returns, and adding jumps to the option pricing model is also an extension. Merton [22] added the lognormal jumps to the option pricing model since he mentioned that the changes in stock prices appear to be jumps. However, the volatility is still a constant parameter in his model. Kou [23] proposed a double exponential jump-diffusion model for capturing discontinuous nature that the asset returns have; jump size is double exponentially distributed in his paper, and it can capture the leptokurtic feature and the volatility smile. Empirical studies also support the improvement of the option pricing model with jumps. Bates [24] established a model combining stochastic volatility and lognormal jumps, and his empirical work shows that his model can improve on fitting option prices. Bakshi et al. [25] also found that adding jumps to the option pricing model with stochastic volatility can improve the performance on pricing options, especially for short time to maturity. Some authors also make some expansions that they focused on the option pricing modelling combining stochastic volatility, stochastic interest rates, and jumps, and several authors think that this kind of model is more reasonable and appropriate for application. Scott [26] supported better performance of this kind of the option pricing model. Jiang [27] also indicated that though their empirical study demonstrates that the dynamics of the interest rate has little impact on option pricing, option pricing modelling combining the three factors is still robust over time. Besides assuming the interest rate follows a CIR process, the Hull-White interest rate process also deserves studying [28].
In addition, empirical studies also support the existence of the change in the jump intensity over time. Santa-Clara and Yan [29] proposed that the volatility and the jump intensity change over time in their model. eir empirical results show that the volatility varies over time, and the jump intensity varies much wider than the former. Chang et al. [30] used ten years of stock returns data to affirm the existence of the change in the jump intensity over time and capture the switching of the jump intensity. Huang et al. [31] proposed a model combing stochastic volatility and jumps with stochastic intensity. ey derived the characteristic function and did some numerical study based on it by applying FFT. Several authors presented models combing stochastic volatility, stochastic interest rates, and jumps with stochastic intensity [32,33].
We present option pricing combining double stochastic volatility, stochastic interest rates, and double exponential jumps with stochastic intensity in this article. We derive the semi-closed form pricing formula and use it as the benchmark to examine some properties of two numerical approaches generally used to approximate the pricing formula for European options. Fast Fourier transform (FFT) and the COS method are two accurate and efficient numerical approaches generally used to approximate the formula for European option prices. Carr and Madan [34] developed a straightforward and efficient expression for the Fourier transform and used FFT to approximate the pricing formula for the options numerically. It makes computation simple and efficient and has a significant reduction in computation time. Fang and Oosterlee [35] developed the COS method for approximating the pricing formula for European options. Since the density function is expressed with the characteristic function via inverse Fourier transform, they used the Fourier cosine series to replace it for approximating the pricing formula. eir numerical results show that it is a highly efficient approach. In the numerical analysis part, we use FFT and the COS method to approximate the formula for European call option prices and compare the computation speed, the accuracy, the efficiency, and the stability between the two approaches. We examine the impact of the parameters in the jump intensity process and the existence of the process on call option prices. We also examine the impact of the parameters in the interest rate process on the call option prices.
We make two contributions based on the existing literature. First, we add double stochastic volatility to the option pricing model combining stochastic interest rates and jumps with stochastic intensity since the double stochastic volatility model is more applicable for practical application than the one-factor stochastic volatility model [10], and we are the first to fill this gap. Second, the stochastic interest rate process is presented in the Hull-White model. In recent years, some authors have concentrated on hybrid models based on various asset classes [20,21], and the option pricing model with all these features have the possibility to develop more proper option prices [21]. erefore, we build a multifactor model with the term structure of stochastic interest rates. It is an extension of the work by Grzelak et al. [21]. Since our model is in the class of the affine jumpdiffusion (AJD) process which was introduced by Duffie et al. [36], we use their results to obtain the formula for the discounted characteristic function. Duffie et al. [36] did the pioneering research on the AJD process, the securities were modelled under an equivalent martingale measure, they extended the Heston model to be a multidimensional model, and the Fourier transform of the security price is in closed form which means that the coefficients in the expression of the Fourier form need to satisfy some ordinary differential equations by applying the Feynman-Kac theorem [36].
We form the structure of this article with the following sections. We develop option pricing modelling combining double stochastic volatility, stochastic interest rates, and double exponential jumps with stochastic intensity and derive the semiclosed form pricing formula in Section 2. In Section 3, we discuss and investigate some stochastic differential equations relevant to the Hull-White interest rate process and derive the discounted characteristic function. In Section 4, we do some numerical work. Section 5 provides the conclusion. T] , Q)be a complete probability space with a filtration and Q presents a risk-neutral measure. e stock price S t is expressed by the following dynamic system:

The Model and Semi-Closed Form Formula
, v jt is the volatility, V jt is its square which is called the variance, and λ t is the jump intensity. θ j and θ λ are their mean-reversion levels, κ j and κ λ are their mean-reversion rates, and σ j and σ λ are their volatilities, respectively. r t is the instantaneous spot interest rate, δ is its mean-reversion speed, η is its volatility, and ϑ t is a timevarying drift term, and it is used to match the initial term structure of the interest rates.
N t represents Poisson process with intensity λ t and J represents the jump size, and we assume that lnJ has an asymmetric double exponential distribution with density function p df u (z): where η 1 ≻ 1, η 2 ≻ 0, p, q ≻ 0, and p + q � 1, where p and q represent the probabilities for positive and negative jumps, respectively; therefore, we can obtain that We set X t � ln S t , τ � T − t, Y � ln J, and k � ln K, where T is the maturity date, τ is the time to maturity, and K is the strike price. Under the risk-neutral measure Q, the price of a call option C(S, V 1 , V 2 , r, λ, t) at time t ∈ [0, T] with strike price K and maturity date T is given by we can rewrite it as Mathematical Problems in Engineering 3 e derivation of the semi-closed form formula is presented by applying Radon-Nikodym derivatives. We consider switching Q to the measure Q S and the T forward measure Q T for the first and second expectation parts respectively in (4). We give the following Radon-Nikodym derivatives: is the price at time t of a zerocoupon bond which matures at time T. en, (4) can be rewritten as Since the density function f(x) and the characteristic function f(u) form a Fourier pair, and we define where φ S (u) denotes the characteristic function under Q S , φ T (u) denotes the characteristic function under Q T , and φ(u) denotes the discounted characteristic function under Q.
We can obtain the following equations using Radon-Nikodym derivatives: en, (4) can be rewritten as erefore, we can get the semi-closed form formula once we obtain the formula for the discounted characteristic function φ(u). e formula for φ(u) is derived in the next section.

The Discounted Characteristic Function
e derivation of the formula for the discounted characteristic function is presented in this section. To be specific, first, we discuss and investigate some stochastic differential equations relevant to the Hull-White interest rate process, present the Hull-White decomposition, and enumerate some relevant formulae including the pricing formula for a zero-coupon bond. Second, we use the results given by Duffie et al. [36] to obtain the formula for the discounted characteristic function.
Applying Itô's lemma to the Hull-White model we obtain that d e δt r t � δe δt ϑ t dt + ηe δt dW r t .
We integrate (15) to obtain that erefore, r T is a normally distributed conditional on F t with 4 Mathematical Problems in Engineering e interest rate process in (1) can be decomposed into r t � r t + ψ t , and this is well known as the Hull-White decomposition. ψ t and r t are given by We can obtain the following stochastic differential equation using Itô's lemma: Theorem 1. If the dynamics of r t is given by the stochastic process (21), we define f(t, T, r): it takes the following form: where Proof. Since (21) is an Ornstein-Uhlenbeck process and it possesses an affine term structure, we conjecture that where f(t, T, r) satisfies a partial differential equation by applying the Feynman-Kac theorem [36]: with boundary condition f(T, T, r) � 1.
We can obtain two ordinary differential equations by rearranging (25) in terms of (24) with boundary conditions C(T, T) � 0 and D(T, T) � 0.
We can obtain the solutions for C(t, T) and D(t, T) by solving the above two ordinary differential equations, thus the proof is complete. erefore, we can obtain the pricing formula for P(0, T): us, ψ T can be given by where f(0, T) � − (zlnP(0, T)/zT), and we denote f(0, T) as the instantaneous forward interest rate. According to (19), ϑ t can be expressed as ϑ t � (1/δ)(zψ t /zt) + ψ t , and substituting (28) into it yields Setting is the market instantaneous forward rates, the superscript M represents that the value is calculated according to a set yield curve [37]. □ Lemma 1. If the dynamics of r t is governed by the stochastic process (21), we can have the following equation: Proof. e proof of Lemma 1 is in Appendix A.

Theorem 2. Since r T is a normally distributed conditional on
Proof. e Hull-White decomposition leads to We can obtain the formula for the first integral part in (33) by using Lemma 1, and we only need to derive the formula for the second integral part to obtain mean μ R and variance σ 2 R . According to (27), we can express P M (0, T) as erefore, we can obtain that which leads to erefore, we can obtain that Mathematical Problems in Engineering 5 (0, t)).

(37)
We can get equations (31) and (32) where Proof. eorem 2 can lead us to the following equation: where ψ � ψ t . We can obtain (38) by rearranging the above equation. e proof is complete.
□ Theorem 4. If the asset price is governed by the dynamic system (1), the discounted characteristic function φ(u; X, V 1 , V 2 , r, λ, τ) takes the following form: where Proof. Although our model is in the class of the AJD process, we still need to separate X t into two parts and use the Hull-White decomposition before using the results given by Duffie et al. [36]. According to Grzelak et al. [21], we define X t : � X t + Ψ t , with Ψ t � t 0 ψ s ds, and then we can obtain the following system by using the Hull-White decomposition r t � r t + ψ t : en, we use the results given by Duffie et al. [36] to derive the discounted characteristic function.
We can get a system of six ordinary differential equations by rearranging (45) in terms of (46) and (47): We can obtain the formulae for C A (u, τ), D X (u, τ), D V j (u, τ), D r (u, τ), and D λ (u, τ) by solving the above ordinary differential equations, and thus we can obtain the formula for φ(u; X, V 1 , V 2 , r, λ, τ). e discounted characteristic function φ(u; X, V 1 , V 2 , r, λ, τ) can be expressed as the following form: According to (36), we can obtain that (0,t)) .

Numerical Discussion
e approximated pricing formula for the call options using the COS method is derived in this section. We use it to do some numerical analyses to compare the computation speed, the accuracy, the efficiency, and the stability between FFT and the COS method. We also investigate the impact of the parameters in the jump intensity process and the interest rate process and the existence of the jump intensity process on the call option prices.

Theorem 5.
e pricing formula for the call options with the COS method V(t, x) at time t ∈ [0, T] is approximated on a bounded interval [a, b]: where the apostrophe on the right side of the summation symbol means that the first item is weighted by 1/2, R(·) is the real part, x � ln(S/K), and ϕ(u; x) is the characteristic function of x, and it satisfies the following equation [35,38]: , a and b are specific constants that satisfy a≺b, and W k are the cosine series coefficients of the call option payoff.
Proof. e pricing formula under the T forward measure Q T is given by wheref(y|x) is the density function of y given x with respect to Q T , y � ln(S T /K), and V(T, y) is the call option payoff. According to Fang and Oosterlee [35], f(y | x) decays fast to zero as y ⟶ ± ∞, so we can truncate the integration to make the difference between the true value and approximation negligible by choosing a specific interval [a, b], thus its approximation can be given by Since f(y | x) decays fast, we can approximate the characteristic function ϕ(u; x) on the interval [a, b]: To find an analytical approximation formula for f(y | x), we can express it with the following Fourier cosine series expansion: then A k can be given by thus the approximated pricing formula can be rewritten by where the cosine coefficient W k is given by Since the coefficients decay fast, the summation can be truncated to obtain that e call option payoff is given by then W k can be rewritten by Using the basic integration rules straightforward, we can obtain that where Mathematical Problems in Engineering To make the difference between the true value and approximation negligible, we need to determine the interval [a, b] appropriately, and the range of interval [a, b] is determined by [35] where c n (n � 1, . . . , 4) are the n-th cumulant of ln(S T /K) and L is the truncation parameter. c n can be given by c n � z n (ln(ϕ(u))) zu n u�0 . (66) Since we cannot obtain the cumulants directly under the condition that the interest rate follows the Hull-White process, we use some specific and suitable range for approximation according to Grzelak et al. [21]: We use the formulae derived above and set the values of the parameters to do some numerical analyses.
e parameters are given in Table 1. Table 2 shows the numerical results. We use the integration approach to calculate the semi-closed form prices obtained by (14), and it takes a large amount of time for calculation. e COS method and FFT are much faster than the integration approach which means that the two approaches make big improvement in computation speed. e numerical results demonstrate that the price differences between the COS method and semi-closed form prices are negligible compared to the price differences between FFT and semi-closed form prices which means that the COS method shows higher accuracy than FFT.
We compare the error convergence between the COS method and FFT with different grid points. We set the semiclosed form price as the benchmark, the values of the grid points with N � 2 n (6, . . . 10), the strike price K � 80, 100, and 120, T � 1, and Table 3 shows the result. e differences between the semi-closed form price and the prices computed by applying the COS method are negligible, and the COS method converges much faster than FFT which means that the COS method shows higher efficiency than FFT.
We examine the relative differences of the call option prices with different grid points to compare the COS method and FFT in terms of the stability. We set the value of grid points N � 2 n (n � 15) as the benchmark, grid points with N � 2 n (6, . . . 10), the strike prices K � 80, 100, and 120, T � 1, and Table 4 shows the result. e relative differences computed using the COS method are lower than those computed using FFT for all the chosen values of grid points, respectively. It demonstrates that the COS method is more stable than FFT which means that the COS method shows higher stability than FFT.
Since the COS method is more accurate, efficient, and stable than FFT to approximate the option prices, we use it for approximation to investigate the impact of the parameters in the jump intensity process and the interest rate process and the existence of the jump intensity process on the call option prices. Figure 1 illustrates that the change of the mean-reversion rate κ λ has little impact on call option prices, and the change of the mean-reversion level θ λ has important impact on call option prices. As θ λ increases, call option price also increases, and the change of the call option price is an increasing function of θ λ . e reason that this phenomenon happens is possibly when the mean level of the jump intensity is high and the cost of investing stocks becomes higher; therefore, the investors tend to buy call options to lower the cost of their investment which makes the call option price become higher.
We investigate the impact of the existence of the jump intensity process on the call option prices, and Figure 2 shows the result. It illustrates that the call option prices with stochastic jump intensity are higher than the call option prices with constant jump intensity. e reason that this phenomenon happens is possibly when the jump intensity changes over time, and there is a great opportunity that the cost of investing stocks becomes higher; therefore, the investors tend to buy call options to lower the cost of their investment which makes the call option price become higher.
We investigate the impact of the parameters in the interest rate process on the call option prices, and Figure 3 shows the result. It illustrates that the change of the meanreversion rate δ has little impact on call option prices and the      change of the volatility η has important impact on call option prices. As η increases, call option price also increases, and the change of call option price is an increasing function of η. e possible reason that this phenomenon happens is that greater volatility of the interest rate means a greater opportunity of the increase of it. e interest rate is the opportunity cost of the investment in stocks, options, and other financial products. When the interest rates increase, the cost of investing stocks becomes higher; therefore, the investors tend to buy call options to lower the cost of their investment which makes the call option price become higher. e investors can obtain the same profits by investing in the call options instead of investing in the stocks.

Conclusion
We addressed European option pricing under a double stochastic volatility model with stochastic interest rates and double exponential jumps with stochastic intensity in this article.
In theoretical part, we used Radon-Nikodym derivatives to derive the semi-closed form valuation formula with the expression of the discounted characteristic function which means we only need to derive the formula for the discounted characteristic function for obtaining the semi-closed form valuation formula. We used the results given by Duffie et al. [36] to derive the discounted characteristic function.
In the numerical analysis part, we derived the approximated pricing formula by applying the COS method and FFT and compared the calculation speed, the accuracy, the efficiency, and the stability between the COS method and FFT. e numerical results demonstrate that it takes a large amount of time to calculate the semi-closed form prices using the integration approach. Both the COS method and FFT takes less time, and they improve in computation speed. e price differences between the COS method and semi-closed form prices are negligible compared to the price differences between FFT and semiclosed form prices which means that the COS method shows higher accuracy than FFT. We compare the COS method and FFT in terms of the convergence, and it demonstrates that the COS method shows higher efficiency than FFT. We examined the relative differences of call option prices with different grid points to compare the COS method and FFT in terms of the stability. e result demonstrates that the COS method shows higher stability than FFT. Because of the higher accuracy, efficiency, and stability of the COS method, we use it to investigate the impact of the parameters in the jump intensity process on call option prices. e numerical results illustrate that the change of the mean-reversion rate has little impact on call option prices and the change of the mean-reversion level has important impact on call option prices. We examine the impact of the existence of the jump intensity on the call option prices with the COS method, and it illustrates that the call option prices with stochastic jump intensity are higher than the call option prices with constant jump intensity. We also investigate the impact of the parameters in the interest rate process on the call option prices with the COS method. It illustrates that the change of the mean-reversion rate has little impact on call option prices, and the change of the volatility has important impact on call option prices.

A. Proof of Lemma 1
In this appendix, we derive Lemma 1, using the method initially presented by Brigo and Mercurio [39]. We can get the following equation by applying stochastic integration by parts: