Efficient Simulation for Pricing Barrier Options with Two-Factor Stochastic Volatility and Stochastic Interest Rate

This paper presents an extension of the double Heston stochastic volatility model by combiningHull-White stochastic interest rates. By the change of numeraire and quadratic exponential scheme, this paper develops a new simulation scheme for the extended model. By combining control variates and antithetic variates, this paper provides an efficient Monte Carlo simulation algorithm for pricing barrier options. Based on the differential evolution algorithm the extended model is calibrated to S&P 500 index options to obtain themodel parameter values. Numerical results show that the proposed simulation scheme outperforms the Euler scheme, the proposed simulation algorithm is efficient for pricing barrier options, and the extended model is flexible to fit the implied volatility surface.


Introduction
A barrier option is a path-dependent option which is exterminated (knocked out) or initiated (knocked in) if the underlying spot price hits the specified barrier level during the life of the option.Because of this supplementary risk, barrier options are cheaper than plain vanilla options and thus are widely traded in exchanges worldwide.One-factor stochastic volatility models can generate "smile, " leverage effects, and term structure effects which cannot be explained by the Black-Scholes model [1,2].Consequently many papers [3][4][5][6][7][8][9] evaluate barrier options under one-factor stochastic volatility models.
However, one-factor stochastic volatility models including the Heston model [10] present the poor performance when fitting the stiff volatility skews [11].One extension by using multiple stochastic volatility factors has been presented in some literatures [2,12,13].Christoffersen et al. [2] confirm that the double Heston model significantly improves the flexibility of the one-factor stochastic volatility model in capturing the volatility term structure.In addition, stochastic interest rate is crucial for option pricing because it ensures proper discounting of future payoffs.In recent literatures [14][15][16][17], Hull-White stochastic interest rate which is analytically tractable has been incorporated into one-factor stochastic volatility model for pricing path-dependent options.Therefore, the model which incorporates multifactor stochastic volatility and stochastic interest rate may be more reasonable for pricing barrier options.
Barrier options with less stochastic factors can be efficiently evaluated by partial differential equation (PDE) methods [3,4,18,19].However, the evaluation of barrier options with multiple stochastic factors is to solve a highdimensional PDE which makes PDE methods quite complex and potentially prone to accuracy and stability problems [1].A more efficient method for pricing barrier options in this case is the Monte Carlo method; see [20,21].The main purpose of this paper is to provide a Monte Carlo method for pricing barrier options under a two-factor stochastic volatility and stochastic interest rate model.
The main contributions of this paper are threefold.Firstly, this paper extends the double Heston model to stochastic interest rate.Secondly, this paper provides a new simulation scheme for the extended model.Thirdly, the paper develops an efficient Monte Carlo algorithm for pricing barrier options.The rest of the paper is organized as follows.Section 2 presents the extended model.Section 3 details the simulation scheme for the proposed model.Section 4 develops the

Simulation Scheme for the DHHW Model
To reduce the dimension of the Monte Carlo simulation, we change from the measure  to the -forward measure  by using (0, ) as numeraire.Set By the Itô formula, we rewrite (6) as follows: 3.1.Variance Simulation.Based on the fact that a noncentral  2 distribution with high noncentrality parameter can be well approximated by a normal distribution, we use quadratic exponential scheme [24] for discrete variance process   ().By (2) and simple calculation, we have Set   =  2  / 2  .Provided that   ≤ 2 and   → ∞, we approximate   ( + Δ) (the segment of high value) by where and   are independent standard normal random variables.
Provided that   ≥ 1 and   → 0, we approximate   ( + Δ) (the segment of low value) by where   are independent uniform random numbers.

Asset Price Simulation.
Under the -forward measure , by integrating (11) and applying the Cholesky decomposition, we rewrite asset price process (6) as follows: where   () are Brownian motions independent of    () and   ().
By the drift interpolation method [25], we approximate the integral of the variance process   () ( = 1, 2) by By ( 7) and ( 8), we have The integral ∫ where   are independent standard normal random variables.

Simulation Algorithm for Pricing Barrier Options
Under the -forward measure , we evaluate barrier options by the following formula: where Λ((  ),   ) is a payoff function,  is barrier level, and   is the first time when barrier is hit.For a down and out call (DOC) option   is defined as follows: Λ((  ),   ) is given by where  is exercise price at maturity time .

The Basic Algorithm.
Based on the simulation scheme for the DHHW model, we evaluate a DOC option by Algorithm 1.

The Hybrid Algorithm Based on Variance Reduction Techniques.
We combine antithetic variate and control variate techniques to improve the efficiency of the Algorithm 1.We use a European option price as control variate and modify step (19) of Algorithm 1 as follows: where   (0) is the exact price of the European call option which can be efficiently calculated by the fast Fourier transform method [26,27].Given the characteristic function,   (0) can be approximated by the following formula: where  = log ,  is damping factor,  is an imaginary unit,  = ℎ/2,  = 2/ℎ (ℎ is a regular spacing),   is the Kronecker delta function that is unity for  = 0 and zero otherwise, and where Based on the results of Grzelak et al. [14] and simple calculation, we have the discounted -forward characteristic function of the DHHW model as follows: where For computing Φ −1 (  ) in step (7) of Algorithm 1, we use the approximation algorithm of Wichura [28]  (29)

Calibration of the Model.
Before applying the DHHW model to option pricing, we need to estimate the model parameters.For this purpose we use financial market data to calibrate the DHHW model.We define the following relative mean squares error (RMSE): where  Θ  and   are the model and market prices, respectively.We calibrate the DHHW model by solving the following nonlinear least-squares optimization problem: We choose the S&P 500 index call options on March 28, 2017, available online at http://www.cboe.com/.The options have maturities between 366 days and 997 days and strike prices ranging from 2225 to 2500.The closing price of underlying is 2358.57.For the sake of simplicity, we set the dividend yield to zero.We calculate the model prices by (26).We use differential evolution algorithm (global optimum) to seek the optimal parameter set Θ * .Table 1 lists the calibrated parameter values and associated RMSE.

Implementation of the Simulation Algorithm.
With the obtained parameter values we first test the performance of the simulation scheme.We apply the simulation scheme to European options and take the fast Fourier transform solutions as benchmark.In (26) we use  = 4096, ℎ = /300 = 0.01047, and  = 1.25 and set  = 100,  = 100, and  = 2.For our simulation scheme, we specify   = 1.5 throughout the paper.We use 1000, 10000, and 100000 Table 2 shows that although there is almost no difference in the sample standard error between the two schemes, significant relative error remains in Euler scheme compared to our simulation scheme.Moreover, it is observed that standard error of our simulation scheme decreases with the increase of the number of simulations and is almost unaffected by time steps.Table 2 implies that our simulation scheme outperforms the Euler scheme.
With the same parameter setting we evaluate the DOC options with barrier  = 90 by the hybrid algorithm.We still use 1000, 10000, and 100000 simulation trials.Since our simulation scheme is almost unaffected by time steps we only use 100 time steps.Table 3 compares the hybrid algorithm with basic algorithm, control variates, and antithetic variates techniques.The criterion used is sample standard errors produced by the above methods.Table 3 shows that both control variates and antithetic variates can reduce standard errors.The hybrid algorithm which combines the two techniques significantly reduces the standard errors.Table 3 implies that the hybrid algorithm is efficient for pricing barrier options.
Furthermore, based on the hybrid algorithm, we plot the implied volatility surface of the DHHW model.Figure 1 examines the effects of interest rate parameters   ,   and the second variance parameters  2 ,  2 ,  2 ,  2 on the implied volatility surface.To examine the effects of the above parameters we change each parameter by setting three different values while fixing all other parameter values.
Figure 1 summarizes our findings.Both   and   only affect the long end regarding the maturity of the volatility surface.However, decreasing   leads to a significant rise of the slope of the implied volatility smile, while the contrary occurs by decreasing   .Increasing  2 results in a fall of volatility, and this effect on the short-term volatility is more significant, which leads to the slight decrease of the curvature of the smile.Increasing  2 results in the rise of overall volatility, but the slope of the smile hardly changed.The curvature of the smile slightly increases by increasing  2 .By changing  2 , the slope of the smile significantly changes and a negative slope is displayed in the implied volatility.Figure 1 implies that interest rates and the second variance process have important effects on the implied volatility surface of the DHHW model.

Conclusion
We propose the DHHW model by combining the double Heston stochastic volatility and Hull-White stochastic interest rate.Under the -forward measure, this paper develops a simulation scheme for the DHHW model.proposed scheme outperforms Euler scheme and the hybrid algorithm is efficient and easy to implement in pricing barrier options.Extensive implied volatility experiments show that implied volatility surfaces present different shapes by varying the parameter values of interest rate and the second variance process, which verify that the DHHW model is flexible to fit the implied volatility smile.

[
not by the Matlab function norminv.m.This can greatly decrease the simulation time.With the obtained   , we consider antithetic variable −  and compute   () by   (  −   ) 2 and thus obtain the second path of the DHHW model.By combining control variate technique we modify the last step of Algorithm   (0) −  (0, ) max ( () − , 0)] .
Combining the control variates and antithetic variates, this paper provides a hybrid Monte Carlo algorithm for pricing barrier options under the DHHW model.Numerical results show that the

Table 1 :
Estimated parameter values for the DHHW model calibrated to S&P 500 index call options onMarch 28, 2017.

Table 2 :
Comparison of the accuracy of our simulation scheme and Euler scheme for pricing European options.The exact price  = 13.209.

Table 3 :
Comparison of the standard errors among the different Monte Carlo techniques for pricing DOC options, where time steps  = 100 and barrier  = 90.