This paper provides an efficient pricing method for forward starting options based on Shannon wavelet expansions. Specifically, we derive the pricing results under a more realistic stock model that incorporates the double exponential jump, stochastic jump intensity, and two-factor stochastic volatilities to capture various features observed in financial markets. We obtain the characteristic function related to the payoff function; then, the options can be well evaluated by the Shannon wavelet method. Numerical experiments show that this method is fast and accurate compared to the Monte Carlo simulation. Finally, we study the influences of changing some important parameters to further illustrate the robustness and stability of the model.

1. Introduction

A forward starting option belongs to the family of path-dependent options, it is purchased and paid for now but becomes active later at a starting date with an expiration date further in the future. Forward starting options are often used by companies as employee stock options to motivate staff. In addition, insurance companies often choose forward starting options as hedging tools to reduce the risk contained in life insurance products.

A closed-form solution for the price of forward starting options was first solved by Rubinstein [1] based on the Black–Scholes (BS) framework [2]. However, many studies have shown that the BS model fails to capture some key features in financial markets, such as nonconstant volatility of stock price. Therefore, many researchers have introduced stochastic volatility in price dynamics to explain volatility smile. So far, the most common models are proposed by Heston [3], Stein and Stein [4], and Hull and White [5]. Based on this, Kruse and Nögel [6], Amerio [7], and Haastrecht and Pelsser [8] considered pricing forward starting options by adding stochastic volatility to financial markets. However, few researchers such as Da Fonseca et al. [9] have shown that single-factor stochastic volatility models do not perform well when trying to capture the long-term structure of implied volatility. Then, Christoffersen et al. [10] proposed the double Heston model, allowing the volatility dynamics to be driven by two uncorrelated stochastic processes to provide a better fit to the empirical implied volatility.

In addition to stochastic volatility, extreme price movements exist in underlying price dynamics, especially when financial crises happen. Merton [11] first incorporated log-normal jumps in asset price. Kou and Wang [12] has proposed another model assuming that the jump size follows an asymmetric double exponential distribution, it can explain extreme movements in asset price and capture the leptokurtic feature as well. Moreover, prior empirical research [13] revealed a low level of correlation between stock volatility and jump intensity; then, Huang et al. [1416] assumed in their work that the stochastic volatility and jump intensity are governed by separate processes. Therefore, in this paper, we consider combining double stochastic volatilities, asymmetric double exponential jump with stochastic intensity in the modeling framework to valuate forward starting options.

While the above stochastic structures work well in effectively capturing various market features and obtaining more accurate option prices as a result, how to find a closed-form solution under these complicated models becomes a demanding task. The numerical pricing methods for models with stochastic volatility can be classified into three main groups [14, 17]: the Monte Carlo simulation, numerical integration methods, and partial (integro-) differential equation methods. In general, the price of European-style options is determined by the discounted expectation of payoff function. Calculation of the expectation requires knowledge about the probability density function (PDF) of stock price. Although the PDF may not be available in most sophisticated price processes, the Feynman–Kac theorem connects the numerical solution of a partial (integro-) differential equation to the expectation of payoff function. Then, we can calculate the characteristic function for relevant price processes, which is the same as the Fourier transform of the PDF.

The integration methods are used when the characteristic function of the stock price is known, and a well-known example is the fast Fourier transform (FFT) method [18]. However, the convergence of FFT depends on a damping factor, so it is often unstable. Based on Fourier-cosine series expansion, Fang and Oosterlee [17] proposed a more efficient alternative called the COS method. It can achieve an exponential convergence and has been widely used in the valuation of options [1922]. One of the drawbacks of this method is that cosine expansions are formed on a global basis and exhibit periodicity behavior, and errors may accumulate near the domain boundaries.

Ortiz-Gracia and Oosterlee [23] proposed another valuation method for plain vanilla options called the Shannon Wavelet Inverse Fourier Technique (SWIFT). It performs well with more flexibility and accuracy as wavelet functions can capture local features better than cosine series. Another main improvement of the SWIFT method is that it does not require integral truncation ranges because the local wavelet basis can adaptively indicate if the error is under a defined tolerance.

The main goal of this paper is to provide an efficient method to price forward starting options. The contributions of this paper are mainly two-fold. Firstly, we utilize the SWIFT method and derive a closed-form solution. Some numerical integration methods we mentioned above have been used to price forward starting options. Kruse and Nögel [6] did the pioneering research, and they applied the FFT method under the Heston model. Zhang and Geng [20] applied the COS method to price forward starting options under the double Heston model. In this paper, we choose the SWIFT method to valuate forward starting options, due to the established accuracy and robustness of Shannon wavelets in option pricing, as demonstrated in a wealth of existing literature [14, 16, 24, 25]. Secondly, this paper proposes a more realistic model that incorporates two-factor stochastic volatilities, asymmetric double exponential jump with stochastic jump intensity to capture various features observed in financial markets. It is an extension of the work by Christoffersen et al. [10], Kou and Wang [12] and Huang et al. [16].

The rest of the paper is arranged as follows: Section 2 sets up the model of stock price dynamics. In Section 3, we derive the characteristic function of the log-asset price. Section 4 first briefly introduces the Shannon wavelet method; based on this, we derive a formula for pricing European-style forward starting options. Numerical results and sensitivity analysis are given in Section 5, and we conclude the paper in Section 6.

2. The Model

We assume that is a complete probability space where four Brownian motions , , , , and are defined, is a filtration, and is a risk-neutral measure. The stock price , two volatilities and , and jump intensity are specified as follows:where is a constant interest rate and , , and represent the mean-reverting speeds, long-term volatilities, and long-term jump intensity, respectively. The stochastic processes and satisfy the Feller conditions [26], i.e., to make the relevant process remain positive. and , and and are two pairs of Brownian motions with correlation coefficients and , respectively.

is a Poisson process with stochastic jump intensity and jump size , where is a random variable and follows an asymmetric double exponential distribution with the density function , where are up-move and down-move probabilities and represents average jump amplitude, . We further suppose that the processes , , , , and are all independent of and .

Remark 1. Model equation (1) contains some special models as follows:(1)The geometric Brownian motion by setting ;(2)The Heston model with ;(3)The double Heston model when ;(4)The double exponential jump model with .

3. Derivation of Characteristic Function

In this paper, we consider the forward starting options with starting date and expiration date . Let denote the strike price. In this section, we first obtain the characteristic function at time for log-price , then derive the characteristic function of .

Lemma 1. Assume that the underlying asset follows equation (1), then under the risk-neutral measure , the conditional characteristic function is given bywherewith

Proof. See Appendix A.

Lemma 2. Given that the process follows Cox–Ingersoll–Ross model:the expectation of is given by

Proof. See [27].

Theorem 1. Suppose that the stock price follows (1), the characteristic function of is given bywhere

Proof. See Appendix B.

4. Pricing Method

4.1. Shannon Wavelet Series Expansion

Let denote the space of square integrable functions. We start with a family of closed subspaces of with the following properties:

Definition 1. A function is said to generate a multi-resolution analysis (MRA) if it generates a sequence of closed subspaces that satisfy (i), (ii), and (iii), and the set forms an orthonormal basis of .

Definition 2. A function is called a basic scaling function or father wavelet if it generates a MRA.
We define a projection
Observe that at a chosen level of resolution , every function can be approximated:where converges to as , are called .
The scaling functions of Shannon wavelets in is computed as

4.2. Pricing Method for Forward Starting Options

Here, we consider the final payoff function of forward starting options as , with for call and for put, where is the starting time, is the maturity time, and denotes the strike price. Then, the price of forward starting options can be expressed as

Let ; then, the option price equation (12) becomeswhere represents the probability density function (PDF) of .

Based on the wavelets theory discussed above, can be expressed through Shannon wavelet expansions:the second approximately equal sign holds because if tends to zero at infinite, satisfy (see [23] for details).

Substituting equations (13) into (14) and exchanging the summation and integration,where .

4.2.1. Coefficients Computation

Proposition 1.

Proposition 2.

Proof. The results of Propositions 1 and 2 can be easily derived by two-fold duplication formula and product-to-sum formula , respectively.
According to Propositions 1 and 2, sinc(t) can be approximated asThen, we replace in (11) with (18),Consider the Fourier transform of :and , where denotes the real part of argument, the approximate computation of the scaling coefficients can be obtained:

4.2.2. Pricing Formula

For formula (15), we truncated the integration range to a finite domain without significant loss (see [23]), then

If we definewith ; through a straight forward calculation, the above integrals can be expressed as

Thereafter,for a call option, andfor a put.

Based on what we have discussed above, we can obtain the formula on the price on forward starting options:with approximated by (21) and approximated by equations (25) and (26).

5. Numerical Results

In this section, some numerical experiments are performed by utilizing the SWIFT method to price forward starting options. We first discuss the convergence of errors and find an appropriate level of resolution . Then, we display some numerical results to show the performance of the SWIFT method. Finally, we explore the price sensitivity on changing model parameters.

5.1. Error Convergence

We choose the appropriate level of resolution by two error convergence analyses. Firstly, we follow the line of [16] to compute the tail mass that the characteristic function not recovered, which can be approximated by

Figure 1 shows the loss of the tail mass for characteristic function across a range of .

In addition, Huang et al. [16] and Ortiz-Gracia and Oosterlee [23] have discussed that one advantage of utilizing Shannon wavelet is that once the scaling coefficients are derived, and the area under the approximated density function can be easily calculated as

Hence, we further investigate the lost area under probability density function (PDF) through equation (29). Figure 2 presents the error when approximating the density function with Shannon wavelet expansions.

We can observe from Figure 2 that the differences is negligible for . Meanwhile, as can be seen in Figure 1, is an appropriate level of approximation within a loss tolerance level of . Therefore, we choose for numerical experiments.

5.2. Numerical Analysis

Here, following Fang [17], the approximate truncated range can be calculated as

, where and is defined as the derivative of the underlying characteristic function, i.e., . Besides, with .

We use the pricing formula derived above to do numerical experiments and analyse the speed of computation and the accuracy and the stability of the SWIFT method.

The values of parameters are listed in Table 1 for all numerical examples. All of these numerical examples were performed in Python 3.7. Also, the computer we used equips an Intel Core i7 CPU with a 2.2 GHz processor.

Here, we use three different approaches (SWIFT, COS, and Monte Carlo simulation) to calculate the price of forward starting call options. For Monte Carlo method, we use 100,000 numbers of simulations with 200 numbers of time steps. The COS method has been proved to be highly efficient and accurate in a wealth of literature [17, 19, 22, 28]. Monte Carlo simulation is a typical numerical method in the domain of option pricing, and it can be flexibly used for various exotic options and thus becomes one of the most common approaches in practice. Therefore in this paper, we choose the results of the COS method as the benchmark and compare the performance between SWIFT and Monte Carlo methods.

Table 2 displays the resulting prices and average CPU time of the three methods. The pricing results show that the price differences between the SWIFT method and the COS method are negligible compared to the price differences between Monte Carlo simulation and COS, which means that the SWIFT method is more accurate than Monte Carlo simulation. For Monte Carlo method, we use 100,000 numbers of simulations with 200 numbers of time steps, and it takes more than 100 seconds to calculate. The SWIFT and COS make big improvement in computation speed, and the average CPU time they need is significantly less than Monte Carlo.

Then, we examine the pricing error of SWIFT and Monte Carlo with different expiration time and different strike prices . Table 3 shows the result. The error of using the SWIFT method are much lower than those computed using Monte Carlo simulation for all expiration time and strike prices. It demonstrates that the SWIFT method is more stable than Monte Carlo simulation.

5.3. Sensitivity Analysis

We will conduct sensitivity analyses to study the impact of some parameters on pricing forward starting options. These analyses can provide evidence for robustness and stability of the SWIFT method under different market conditions. We mainly study the sensitivity of option prices derived by the SWIFT method to changes in the following parameters: (1) (Figure 3) instantaneous volatility of : , (2) (Figure 4) the speeds of mean-reverting: , (3) (Figure 5) long-term jump intensity: , and (4) (Figure 6) reciprocal of the mean value for the negative jumps: .

6. Conclusion

In this paper, we use a novel valuation method to price forward starting options. We consider a model combining the asymmetric exponential jump, stochastic jump intensity, and two-factor stochastic volatilities to capture various features observed in financial markets and derive a more realistic pricing framework. We calculate the characteristic function related to the final payoff by applying Feynman–Kac formula, and then solve the differential equations. After deriving the characteristic function, the SWIFT method is applied to compute the pricing results of forward starting options. Numerical experiments are performed to show the efficiency compared to the Monte Carlo method. Finally, we investigate the impact of changing model parameters to the resulting option values, thus proving the robustness and stability of our model.


A. Proof of Lemma 1

According to Feynman–Kac theorem, satisfies the partial integrodifferential equation (PIDE) as follows:

Notice that the integral term in equation (A.1) can be written aswhere . According to Duffie et al. [29], the solution of this PIDE has the following form:with boundary conditions .

Substituting equations (A.1) and (A.2) into (A.3) yields

By matching coefficients, the problem can be simplified to four differential equations:

The authors first solve the equation (A.6). Making the substitution,

Substituting equations (A.9) into (A.6) and simplifying, the authors derive

This is a second-order homogenous linear differential equation, and the solution can be written aswhere , .

According to the boundary condition, and . The authors obtain

, then,

The equations for and are solved by analogy. Hence,with

Integrating on both sides of equation (A.5), the authors yieldwith

B. Proof of Theorem 1

Here, the authors use the same method in [20]. Notice that is the characteristic function of which can be calculated through the formula in Lemma 1. Substituting equations (2) into (A.18) yields,

Using Lemma 2 and setting , and , the authors derive

Combining equations (A.19) and the above three expression, then, Theorem 1 follows.

Data Availability

The data used to support the findings of this study are available from the corresponding author upon request.

Conflicts of Interest

The author declares no conflicts of interest.