Abstract

In the paper, the pricing of the American put options under the double Heston model with Cox–Ingersoll–Ross (CIR) interest rate process is studied. The characteristic function of the log asset price is derived, and thereby Bermuda options are well evaluated by means of a state-of-the-art Shannon wavelet inverse Fourier technique (SWIFT), which is a robust and highly efficient pricing method. Based on the SWIFT method, the price of American option can be approximated by using Richardson extrapolation schemes on a series of Bermudan options. Numerical experiments show that the proposed pricing method is efficient, especially for short-term American put options.

1. Introduction

It is well known that both the academic researchers and financial engineers have been concerned with the study of option pricing theory which is one of the most important and useful research fields in financial engineering. The well-known Black–Scholes–Merton (BS) model has been widely applied in the financial industry since it was proposed in the 1970s [1]. Unfortunately, the model has many shortcomings, such as logarithm yield of the stock price are assumed to be normal distribution and volatility of stock price is assumed to be constant which is different from the real data from the financial market. This is an important reason that brings a lot of economic losses to many financial institutions and even leads to the bankruptcy of some financial institutions. And, to make matters worse, the stability and sustainability of the whole financial system are influenced. Therefore, in order to maintain the stability and sustainability of the country’s financial system, the government implemented tighter financial regulation and many scholars continue to study the option pricing theory.

Considering the abovementioned shortcomings, many researchers generalized the BS model by introducing stochastic volatility among which the popular models are the models raised by Stein and Stein [2], Heston [3], and Schöbel and Zhu [4]. However, single-factor stochastic volatility models are unable to fit the implied volatility smile very well. Many empirical research studies, such as Cont and Tankov [5], Fonseca et al. [6], Christoffersen et al. [7], and Fouque and Lorig [8] show that single-factor models can do a bad job in capturing the term structures of implicit volatilities. Christoffersen et al. [7] proposed the double Heston model by adding another uncorrelated volatility process which can provide better empirical fit to the market price and more flexibility to model the volatility term structure than the Heston model. In addition, spot interest rate is a fundamental economic variable in the financial market and it cannot be treated as a constant. Some researchers introduce stochastic interest rate into option pricing models in recent years. Grzelak and Oosterlee [9] examined the correlation issues of European product pricing with the Heston–Hull–White and Heston-CIR hybrid models [9]. Recently, Kim et al. [10] showed that incorporation of stochastic interest rates into a stochastic volatility model gave better results compared with the constant interest rate case in any maturity. It has been verified that considering stochastic interest rate into the Heston model improves the pricing performance of the model for European options in Bakshi et al. [11]. Consequently, the model which combines the two-factor stochastic volatilities and stochastic interest rate is proposed.

Because many exchange-traded derivative contracts are of American style, accurate and fast pricing of American options is important, particularly, for high frequency traders. But the pricing of American options is one of the most challenging problems. The essential difficulty is that the pricing of American options involves a difficulty-free boundary problem in mathematics. There is no exact formula for evaluating most of the American options, and some numerical pricing techniques can be developed. The existing numerical methods under stochastic volatility models can be classified into three major groups: the Monte Carlo simulation, partial integro-differential equation (PIDE) methods, and numerical integration methods. One important type of numerical integration is the so-called Fourier transform method, which has attracted more and more attention since it is efficient and often used for calibration purposes. One of the best known examples of this class is the Carr–Madan method [12] based on fast Fourier transform. Zhylyevskyy [13] also evaluated American options under stochastic volatility by the fast Fourier transform technique. However, fast Fourier transform is usually unstable affected by the dampening factor. A more efficient pricing method based on Fourier cosine series expansion is proposed by Fang and Oosterlee [14]. The main achievement of this method is that it maintains an exponential convergence rate. But, it brings forth disadvantages in the vicinity of integration boundaries because cosines exhibit periodicity behavior. And, round-off errors may accumulate near domain boundaries for long-maturity options. Gracia and Oosterlee [15] further apply the Shannon wavelet inverse Fourier technique (SWIFT) method to price European option under exponential Lévy dynamics, and SWIFT was proved to be a robust, accurate, and highly efficient financial option valuation method. In this paper, we firstly utilize the SWIFT method for pricing Bermudan options, and then the price of American option can be approximated by using Richardson extrapolation schemes on a series of Bermudan options.

An outline of the paper is as follows. In Section 2, we introduce the basic model specification. The characteristic function of the underlying log-asset price is also derived. In Section 3, we briefly introduce the SWIFT method and the procedures to price the American put option. We present a set of numerical results and analysis in Section 4 to evaluate the accuracy and efficiency of the SWIFT method benchmarked to the Monte Carlo method. We conclude the paper in Section 5.

2. Model Specification and Derivation of the Characteristic Function

Assume that denotes a complete probability space with a filtration continuous on the right, where is a risk-neutral probability. Suppose that no dividends are paid in the market and the underlying asset can be traded in continuous time and without transaction cost. Let denote the underlying asset price process. Let , and be Brownian motions. Moreover, ; any other two processes are uncorrelated with each other. The double Heston with CIR interest rate process model assumes the following dynamics:where mean the speeds of the mean reversion, long-run volatility, and instantaneous volatility of process , respectively. The stochastic volatility processes satisfy the Feller conditions [16], i.e., , to make to be positive. In addition, are the mean reversion speed, long-term interest rate, and volatility of the interest rate process , respectively. Let us suppose .

Remark 1. (1) The BS model is a special case of the above model by setting , , and . (2) The Heston model is also a special case of the above model by setting and . (3) The double Heston is another special case of the above model by setting . According to Brigo and Mercurio [17], when the short-term interest rate evolves according to equation (4), the price of a risk-free zero-coupon bond with maturity at time 0 is given bywhereThe conditional characteristic function of the log asset price process is the key component in SWIFT to price the options. We can simplify calculation and change from the numeraire measure to the T-forward measure by applying the following Radon–Nikodym derivative:Then, we derive the conditional characteristic function of the log asset price process in the following result.

Theorem 1. Suppose the price of the underlying asset follows the dynamics (1)–(4). Let . Under the T-forward measure , the conditional characteristic function of is given bywhere

Proof. By applying Feynman-Kac theorem, satisfies the following PDE:According to Duffie et al. [18], we conjecture has the following form:with initial conditions . Substituting (22) into (21) yieldsFrom the above equation, we have a system of four ordinary differential equations as follows:Then, we can obtain

3. American Option Pricing with Shannon Wavelet Expansions

3.1. Shannon Wavelet Expansions for European Options

Before we discuss the American options, we firstly introduce briefly the SWIFT method for European pricing problems as in Gracia and Oosterlee [15]. The price of a European option can be found by means of the Feynman–Kac formula as the discounted expectation of the option value at final time T. According to the risk-neutral option valuation formula,wherewhere denotes the option value, is the maturity, is the initial date, is the time to maturity, is the risk-neutral expectation, and are state variables at time and respectively, and is the probability density of given . As can be seen from the above expression, the key to price European option is to calculate the conditional density function. The SWIFT method can approximate the density function by recovering it from its characteristic function in three steps.

In the first step, is approximated by Shannon wavelet, that is,where , is the Shannon wavelet.

In the second step, the infinite summation in equation (25) has to be truncated for some , thus we obtain

The final step is to compute the density coefficients , i.e.,where . How to compute the vector of coefficients efficiently using the fast Fourier transform has been showed in by Maree [19]. To sum up, by means of the SWIFT method, the resulting option value at time t can be written as follows:

Let

It is closely related to the wavelet coefficients of the value function . For put options, the following approximate solutions can be obtained:where detailed expression of can be found in [15]. Thus, we can obtain the following result:

3.2. Shannon Wavelet Expansions for Bermudan Options

According to Chang and Chung [20], the prices of American options can be obtained by applying a Richardson extrapolation on the prices of a small number of Bermudan options. So, we need obtain the price of the Bermudan option firstly.

A Bermudan option is a financial contract that allows the holder to exercise the option for a predetermined period of times before the option expires. Without loss of generality, we can assume a Bermudan option with strike price and a set of exercise moments , and strictly ordered, , where is the option’s maturity. In the case of Bermudan option pricing, the dynamic programming principle reveals that the price of the Bermudan option at any exercise moment is the maximum of the spot payoff and the so-called continuation value. The valuation process can be regarded as a European option pricing problem between the two exercise times and can be priced by means of risk-neutral option pricing formula. For simplicity, we can divide time grid and we define and , where is the price process of the underlying asset price process. The payoff of the option is given by , and for the European option, the value of the option at maturity is given by

It is worth noting that we consider the stochastic interest rate in this article. By the dynamic programming method, the option value ahead of maturity can be expressed recursively for , byand the option value at is finally followed bywhere is the continuation value and the probability function of given is written as . The European option value is given by the integral representation in (24). According to the SWIFT pricing formula (26), we obtain

So, the option value can be known once the value coefficients are determined. According to the characterization of Bermudan options, a backward recursion will be applied to recover these coefficients . We know that the option value at maturity equals the payoff of the option, i.e., . Thus, the value coefficients for a put option is given bywhere means the payoff coefficients, and it can be efficiently computed using the fast Fourier transform as Oosterlee mentioned in [19]. Once a formula for the value coefficients at maturity is obtained, we can calculate the coefficients at any time , for , ahead of maturity recursively. When the early-exercise point is known, which is the point where the continuation value equals the payoff, i.e., , where

Once is obtained, the integral for the value coefficients can be splitted into two parts, i.e.,where denotes the continuation coefficients at time over the interval . It also can be efficiently approximated when the value coefficients at the next step are known, i.e.,wherewith

We can see that is a Hankel matrix and represents a matrix-vector product. We can compute this matrix-vector product by means of FFT since the matrix-vector product is equal to the first elements of the circular convolution , with the 2N-vectors

At the same time, a circular convolution of two vectors is equal to the inverse discrete Fourier transform of the product of the forward DFTs, i.e.,

Thus, according to the abovementioned formulas, can be obtained by the backward recursion. Between and , the valuation of the Bermudan option process can be regarded as a European option pricing problem and can be priced with the help of the risk-neutral option valuation formula (15). Thus, the price of Bermudan option at can be obtained as follows:

The prices of American options can be obtained by applying a Richardson extrapolation on the prices of a few Bermudan options [18], which is also demonstrated in [19]. Let denote the value of a Bermudan option with early exercise dates. We will use the following 4-point Richardson extrapolation scheme:where is a positive integer and denotes the approximated value of the American option. So, we can approximate the price of an American option by four Bermudan options with certain number of exercise dates.

4. Numerical Examples

In this section, some numerical examples are performed to test the performance of the SWIFT method for pricing American puts against its alternative competitor, the LSM Monte Carlo method [21], as a benchmark. In particular, we make use of the same constant parameters as specified by Zhang and Zhang [22]. These parameter values are listed in Table 1 for all our numerical examples. The computer used in the experiments equips an Intel Core i7 CPU with a 2.2 GHz processor. And, all of these numerical examples were performed in Matlab 2016a. For the SWIFT method, unless otherwise mentioned, we use L = 10, where L is the parameter used to define the truncation range [a, b] as follows:where is the nth cumulant of log-price process . It can be computed by

The relative error and CPU time information, comparing the SWIFT and the LSM method, are presented in Table 2 for the pricing of American put options. Table 2 compares relative errors between the SWIFT method and LSM for pricing American put options when in equation (27), and it shows that the SWIFT method is significantly faster than LSM. Furthermore, Table 2 compares the pricing accuracy among the two methods across a range of strike prices, and the relative percentage pricing differences of SWIFT are all less than 0.13%. Our numerical results show that the SWIFT method is fast and accurate for pricing American put options.

5. Conclusion

This article proposed an efficient method for pricing American put options under double stochastic volatility with stochastic interest rate. The conditional characteristic function of the log-asset price is derived, and an analytical approximation method for evaluating American puts is proposed. Numerical results of this article confirm that the SWIFT method is efficient and accurate for pricing the American put options. The future research is mainly about parameter calibration and hedge strategy.

Data Availability

No data were used to support this study.

Conflicts of Interest

The authors declare no conflicts of interest.

Authors’ Contributions

Shoude Huang managed the manuscript. Xunxiang Guo programmed the numerical examples and wrote the review.

Acknowledgments

This research was financially supported by the Fundamental Research Funds for the Central Universities (JBK2002003), and the authors gratefully acknowledge this support.