Fractional Stochastic Volatility Pricing of European Option Based on Self-Adaptive Differential Evolution

Yue Hu, Hongling Dong^*, Le Fu, Jiayang Zhai

School of Science, Zhejiang University of Science and Technology, Hangzhou, China.

DOI: 10.4236/jmf.2022.123029                

Abstract

The option pricing model can predict the future trend of the financial market. In order to more accurately describe the changing process of the financial market, the Hurst index which can describe the characteristics of long-term memory is introduced into the traditional Heston model. Under the assumption that the underlying asset price follows fractional Brownian motion, the fractional stochastic volatility pricing of European option pricing model (Hurst-Heston model) is constructed, and the closed solution of the model is obtained according to the partial differential equation satisfied by the model. By analyzing the relationship between Hurst index and asset price, it is found that the movement process of asset price under the hypothesis of this model is more consistent with the real market change law, which verifies the rationality of the model. In the process of empirical analysis of SSE 50ETF put option data, Self-adaptive Differential Evolution algorithm is used to estimate parameters. The results showed that the error of Hurst-Heston model is smaller than other models, and the prediction error for consecutive trading days is similar. It showed that the pricing results of Hurst-Heston model are more accurate and stable.

Keywords

Fractional Brownian Motion, Stochastic Volatility, European Option Pricing, Self-adaptive Differential Evolution Algorithm

1. Introduction

The option derived in the late 18th century is a financial tool used by financial institutions and enterprises to avoid risks [1]. It first appeared in the European and American financial markets. As the contract, it gives holders the right to purchase or sell an asset at the agreed price at a specific time or at any time before that time. On February 9, 2015, Shanghai Stock Exchange pilot listed the first stock option “SSE 50ETF option” in China. After that, China’s option market developed rapidly, and many different types of options were also listed and traded.

A reasonable option pricing model can accurately predict the future trend of the market, which is of great significance to the steady development of the financial market. In 1973 [2], the first complete option pricing model “B-S model” created by Black and Scholes was used publicly, and since then this model has been widely used in European options pricing [3] [4]. In its basic assumptions, such as the stock price follows a log-normal distribution, the risk-free interest rate is known, the stock price volatility is constant, there is no transaction cost in the hedging portfolio, etc. These conditions cannot be met in the real financial market transactions, so phenomena such as “volatility smile” often occur in the actual application. In order to make up for the deficiency of B-S model, Merton first proposed the lognormal jump diffusion model in 1976 [5], it added a jump process into the B-S model, which was more in line with the real market with sudden risks. Subsequently, [6] proposed a call option replication strategy based on proportional transaction costs, and gave the option pricing formula with transaction costs. In addition, because stock price volatility is random, some scholars introduced stochastic volatility model to improve the option pricing model: [7] constructed Heston model on the basis of traditional B-S model and assumed that volatility of underlying asset followed O-U process; [8] proposed fast mean regression volatility model, which reduced model parameters to facilitate parameter estimation, and at the same time, it better explained the characteristics of volatility smile and return on assets peak thick tail; [9] [10] used non-affine stochastic volatility model for option pricing, which is more accurate than other models; [11] generalized the nonlinear partial differential equation when the underlying asset follows the stochastic single-factor interest rate model, regarded the nonlinear term in the equation as the transaction cost, and proved the existence of the classical solution of the model; [12] proposed the volatility decomposition model, , it divided the stochastic volatility process of the underlying asset into two risk sources based on the Heston model, which could describe the relevant characteristics of the underlying asset more accurately.

The Hurst index H was discovered by British hydrologist H.E. Hurst, when he studied the relationship between water flow and storage capacity [13], and took it as an indicator to judge whether timing data followed the process of biased random walk. $H=0.5$, means that the time series is a standard Brownian motion process; $0.5<H<1$, indicating positive correlation of time series and has long memory; $0\le H<0.5$ is the negative correlation of time series, and is the mean recovery process; In other words, if $H\ne 0.5$, time series data can be described by fractional Brownian motion. As for the characteristics of China’s stock market, a large number of research results show that there is a significant long memory in China’s stock market [14] [15] [16] [17] [18], so many scholars have introduced the Hurst index into the financial market as an indicator to describe the long memory. For example, [19] improved the European option pricing model by combining the switch mode and the fractal Black-Scholes market hypothesis; [20] obtained the explicit integral representation of the early exercise premium and the critical exercise price of the American look back option under the assumption that the stock price follows the mixed jump diffusion fractional Brownian motion; [21] constructed a fractional fuzzy option pricing model based on the B-S model using fractal theory and fuzzy set theory, which improved the accuracy of pricing results.

In order to obtain more accurate prediction results on option pricing, we should not only consider the randomness of stock price fluctuations, but also not ignore the long memory of financial markets. However, the existing option pricing models cannot meet these two properties at the same time, which is the original intention of this paper.

This paper introduced some basic knowledge in Chapter 2. Chapter 3 introduced the newly constructed option pricing model and the method of parameter optimization in the model. Chapter 4 analyzed the relationship between Hurst index and asset prices. Chapter 5 used SSE50ETF put option data to do empirical analysis. Chapter 6 summarized the results of this paper. The main contents are as follows:

1) Based on the traditional Heston model, the Hurst index His introduced, so the fractional Brownian motion is used to describe the process of stock price change, a long-memory fractional stochastic volatility European option pricing model is constructed, and the corresponding option pricing formula is obtained. This model makes up for the deficiency that the existing stochastic volatility option pricing models cannot reflect the characteristics of the fractal market.

2) By analyzing the Hurst index in the model, it is found that the introduction of H into the option pricing model can more accurately reflect the real change process of the financial market, which verifies the rationality of the model.

3) Using the SSE 50ETF put option data for empirical analysis, due to the large number of parameters in this model, the traditional linear least squares parameter estimation method is invalid. Referring to the estimation methods of multi-parameter models in other literatures, this paper will choose the Self-adaptive Differential Evolution algorithm to estimate the parameters of the model. The final analysis results show that, compared with the other pricing model, the fractional stochastic volatility model has more accurate pricing results and less error, which further illustrates the feasibility and accuracy of the pricing model proposed in this paper.

2. Preliminary Knowledge

When constructing the option pricing model, in order to describe the long-term memory of the financial market more intuitively, it is necessary to introduce Fractional Brownian Motion (FBM), the concept of which is given by definition 1, and the Fractional Ito Lemma under the FBM can be seen as lemma 1.

Definition 1 [22]: If the continuous Gaussian process $\left\{B_H\left(t\right),t\in R\right\}$ with parameter H satisfies the following conditions:

1) $B_H\left(0\right)=0$ and for any t and $\Delta t>0$, the incremental expectation of FBM is 0;

2) Fort ands at different moments, their covariance function is:

$E\left[B_H\left(t\right)B_H\left(s\right)\right]=\frac{1}{2}\left({\left|t\right|}^{2H}+{\left|s\right|}^{2H}-{\left|t-s\right|}^{2H}\right)$. (1)

Then the Gaussian process $\left\{B_H\left(t\right),t\in R\right\}$ is called fractional Brownian motion. In formula (1), $H\left(0<H<1\right)$ is the Hurst exponent that describes the relationship between the motion increments. If $H=0.5$, $B_H\left(t\right)$ is the standard Brownian motion.

The core property of FBM is that the increment has stationarity, autocorrelation and self-similarity, and the stock price fluctuation process that conforms to fractal characteristics also has similar properties. So FBM can be used to describe the stock price fluctuation process.

Lemma 1 [23]: Assuming that the price of the derivative is $f=f\left(S_t,t\right)$, and the underlying asset price S_t follows a biased FBM, that is $d{S}_t=\mu S_{t}dt+\sigma S_{t}d{B}_H\left(t\right)$, then, at any time t, there is the following relationship:

$df=\left(\mu S_t\frac{\partial f}{\partial S_t}+\frac{\partial f}{\partial t}+H{\sigma }^2{t}^{2H-1}{S}_t^2\frac{{\partial }^{2}f}{\partial S_t^2}\right)dt+\sigma S_t\frac{\partial f}{\partial S_t}d{B}_H\left(t\right)$. (2)

Equation (2) is called fractional Ito Lemma, where $\mu$ is the drift rate and $\sigma$ is the volatility.

3. Pricing Model Construction and Parameter Optimization

The option pricing model constructed in this paper is an improvement on the basis of the Heston model. It not only considers the dynamic process of the financial market as a stochastic fluctuation, but also combines its characteristics of long memory, which makes up for the defects of previous models that only consider a single influence.

3.1. Hurst-Heston Option Pricing Model

In the financial market, stock market volatility is a dynamic process with long memory characteristics. In order to reflect the nature of the financial market more comprehensively in the option pricing model, this paper improves Heston model and introduces FBM on the basis of this model to obtain the Hurst-Heston model as shown in Definition 2.

Definition 2: It is assumed that under the neutral probability measure Q, the underlying asset S_t follows fractional Brownian motion and the volatility V_t follows the O-U process, that is, the price and volatility of the underlying asset satisfy the following differential equations respectively:

$d{S}_t=\mu S_{t}dt+\sqrt{V_t}{S}_{t}d{B}_H^1\left(t\right)$, (3)

$d{V}_t=\kappa \left(\theta -V_t\right)+\sigma \sqrt{V_t}d{B}_H^2\left(t\right)$. (4)

Equation (3) and Equation (4) used to solve option price is called Hurst-Heston option pricing model. Among them $\mu$ is the drift rate, V_t is the volatility, $\theta$ is the long-term mean value of volatility, $\kappa$ is the mean reversion rate, $\sigma$ is the coefficient of variance variation, that is, the volatility of volatility. Both $B_H^1\left(t\right)$ and $B_H^2\left(t\right)$ are fractional Brownian motion, and $d{B}_H^1\left(t\right)d{B}_H^2\left(t\right)=\rho dt$.

In this model, it is assumed that there are risk-free assetsS_t which can be freely traded in the financial market and meet the following conditions:

1) The price of stock index options in the fractal market only depends on the stock price and the risk-free interest rate;

2) There are no taxes and transaction costs in the market, and no dividends are paid during the transaction;

3) The underlying asset priceS_t satisfies the fractional Brownian motion;

4) The underlying asset fluctuates randomly, and its volatility V_t follows the O-U process.

3.2. Analytical Pricing of European Options under the Hurst-Heston Model

Theorem 1: The partial differential equation satisfied by European option price $f=f\left(t,S,V\right)$ under Hurst-Heston model is:

$\begin{array}{l}\frac{\partial f}{\partial t}+V{S}^{2}H{t}^{2H-1}\frac{{\partial }^{2}f}{\partial S^2}+{\sigma }^{2}VH{t}^{2H-1}\frac{{\partial }^{2}f}{\partial V^2}+\frac{1}{2}SV\sigma \rho \frac{{\partial }^{2}f}{\partial S\partial V}+r\frac{\partial f}{\partial S}S\\ +\left[\kappa \left(\theta -V\right)+\lambda v\right]\frac{\partial f}{\partial V}-rf=0\end{array}$ (5)

In Equation (5), t means any moment, V is the volatility of stock price, S is the stock value, H is the Hurst index, $\sigma$ is the volatility of volatility, $\rho$ is the correlation coefficient of fractional Brownian motion $d{B}_H^1\left(t\right)$ and $d{B}_H^2\left(t\right)$, r is the risk-free interest rate, $\kappa$ is the mean regression rate, $\theta$ is the long-term mean of volatility, $\lambda$ is the ratio of the market price of volatility to volatility.

Prove: According to the standard arbitrage theory, assuming that there is no arbitrage opportunity, there is a risk-free investment portfolio I, which contains an option with a value of $f\left(S,V,t\right)$, $\Delta$ stocks with a value of S, and $\phi$ options with another value of $Z\left(S,V,t\right)$, then the value of the portfolio I is:

$I=f+\Delta S+\phi Z$. (6)

According to differential Equation (3), Equation (4) and Lemma 1, due to $d{B}_H\left(t\right)=\epsilon \sqrt{d{t}^{2H}}$, then in the time interval $dt$, the change of I is:

$\begin{array}{l}dI=df+\Delta dS+\phi dZ\\ =\left[\frac{\partial f}{\partial S}\mu S+\frac{\partial f}{\partial t}+\frac{{\partial }^{2}f}{\partial S^2}V{S}^{2}H{t}^{2H-1}+\frac{{\partial }^{2}f}{\partial V^2}V{\sigma }^{2}H{t}^{2H-1}+\frac{1}{2}\frac{{\partial }^{2}f}{\partial S\partial V}VS\sigma \rho +\Delta \mu S\right]dt\\ +\phi \left[\frac{\partial Z}{\partial S}\mu S+\frac{\partial Z}{\partial t}+\frac{{\partial }^{2}Z}{\partial S^2}V{S}^{2}H{t}^{2H-1}+\frac{{\partial }^{2}Z}{\partial V^2}V{\sigma }^{2}H{t}^{2H-1}+\frac{1}{2}\frac{{\partial }^{2}Z}{\partial S\partial V}VS\sigma \rho \right]dt\end{array}$

$+\left[\frac{\partial f}{\partial V}+\phi \frac{\partial Z}{\partial V}\right]dV+\left[\frac{\partial f}{\partial S}\sqrt{V}S+\Delta \sqrt{V}S+\phi \frac{\partial Z}{\partial S}\sqrt{V}S\right]d{B}_H^1\left(t\right)$ (7)

Since the assumption is under the risk-neutral measure, the $dV$ and $d{B}_H^1\left(t\right)$ terms in Equation (7) are eliminated, so that:

$\{\begin{array}{l}\frac{\partial f}{\partial V}+\phi \frac{\partial Z}{\partial V}=0;\\ \frac{\partial f}{\partial S}\sqrt{V}S+\Delta \sqrt{V}S+\phi \frac{\partial Z}{\partial S}\sqrt{V}S=0.\end{array}$ (8)

Then we can get:

$\{\begin{array}{l}\phi =-\frac{\partial f}{\partial V}/\frac{\partial Z}{\partial V};\\ \Delta =-\phi \frac{\partial Z}{\partial S}-\frac{\partial f}{\partial S}.\end{array}$ (9)

And because $dI=rIdt$, where r is the risk-free interest rate, put it together with Equation (9) into Equation (7), then Equation (7) can be written as follows:

$\begin{array}{l}\frac{\frac{\partial f}{\partial t}+V{S}^{2}H{t}^{2H-1}\frac{{\partial }^{2}f}{\partial S^2}+{\sigma }^{2}VH{t}^{2H-1}\frac{{\partial }^{2}f}{\partial V^2}+\frac{1}{2}SV\sigma \rho \frac{{\partial }^{2}f}{\partial S\partial V}-rf+r\frac{\partial f}{\partial S}S}{\frac{\partial f}{\partial V}}\\ =\frac{\frac{\partial Z}{\partial t}+V{S}^{2}H{t}^{2H-1}\frac{{\partial }^{2}Z}{\partial S^2}+{\sigma }^{2}VH{t}^{2H-1}\frac{{\partial }^{2}Z}{\partial V^2}+\frac{1}{2}SV\sigma \rho \frac{{\partial }^{2}Z}{\partial S\partial V}-rf+r\frac{\partial Z}{\partial S}S}{\frac{\partial Z}{\partial V}}.\end{array}$ (10)

According to the research of Heston [6], it is known that the function $Z\left(S,V,t\right)=-\kappa \left(\theta -V\right)+\lambda \left(S,V,t\right)$, where $\lambda \left(S,V,t\right)$ is the market price of volatility risk. Meanwhile, according to Brecden’s assumption that the market price of volatility is its linear function, that is $\lambda \left(S,V,t\right)=\lambda v$, then Equation (5) can be obtained from Equation (10), which is proved.

According to Theorem 1, analogy with the solution process of Heston model, the option price is simulated by numerical method under the condition of $f\left(t,S,V\right)=\max \left\{S-K,0\right\}$. Using Gaussian quadrature and fast Fourier transform [24], the following Theorem 2 can be obtained by making $x=\ln \left(S_t\right)$, $\varphi =\ln \left(r/r_0\right)$.

Theorem 2: The Pricing formula of European call option in the Hurst-Heston model is:

$\begin{array}{c}C=e^{-r\left(T-t\right)}E\left[\max \left(S_{T-t}-K,0\right)\right]\\ =S_{T-t}{P}_1-K{e}^{-r\left(T-t\right)}{P}_2.\end{array}$ (11)

Among them,

$P_j=\frac{1}{2}+\frac{1}{\pi }\mathrm{Re}\left[\frac{e^{-i\varphi \ln K}{f}_j\left(x,V,t,T;\varphi \right)}{i\varphi }\right]d\varphi ,j=1,2;$

$f_j\left(x,V,t;\varphi \right)=e^{A\left(t,T;\varphi \right)+B\left(t,T;\varphi \right)V+i\varphi x};$

$\begin{array}{c}A\left(t,T;\varphi \right)=ir\varphi \left(T-t\right)+\frac{\kappa \theta }{{\sigma }^2}[\left(b_j-\rho \sigma H{t}^{2H-1}\varphi i+h_j\right)\left(T-t\right)\begin{array}{c}{}_{}{}^{}\\ \stackrel{}{}\end{array}\\ -2\ln \left(\frac{1-g_j{e}^{h_j\left(T-t\right)}}{1-g_j}\right)];\end{array}$

$B\left(t,T;\varphi \right)=\frac{b_j-\rho \sigma H{t}^{2H-1}\varphi i+h_j}{{\sigma }^2}\left(\frac{1-e^{h_j\left(T-t\right)}}{1-g_j{e}^{h_j\left(T-t\right)}}\right);$

$g_j=\frac{b_j-\rho \sigma H{t}^{2H-1}\varphi i+h_j}{b_j-\rho \sigma H{t}^{2H-1}\varphi i-h_j};$

$h_j=\sqrt{{\left(\rho \sigma H{t}^{2H-1}\varphi i\right)}^2-{\sigma }^2\left(2u_j\varphi i-{\varphi }^2\right)};$

$u_1=\frac{1}{2},u_2=-\frac{1}{2};$

$b_1=\kappa -\rho \sigma H{t}^{2H-1},b_2=\kappa ,j=1,2.$

The put option pricing formula is:

$P=S_{T-t}\left(P_1-1\right)+K{e}^{-r\left(T-t\right)}\left(1-P_2\right)$. (12)

Note: The parity formula between the European call option price and the put option price is:

$C-P=S_{T-t}-K{e}^{-r\left(T-t\right)}$. (13)

3.3. Parameter Estimation for Pricing Model

Because the option pricing model proposed in this paper is a multi-parameter model, the common parameter estimation methods such as the least squares method and the maximum likelihood estimation method are not applicable, so the Self-adaptive Differential Evolution (SaDE) algorithm will be used to estimate the parameters of this model. The evolution process of the SaDE algorithm is equivalent to the Differential Evolution algorithm [25] [26] [27], which is improved and optimized on the basis of DE algorithm and solves the problems of premature convergence and search stagnation of DE algorithm. The DE algorithm is an optimization algorithm based on the theory of swarm intelligence. It continuously evolves through cooperation and competition among individuals in the group, retains good individuals, eliminates inferior individuals, and guides the search to approach the optimal solution. Its essence is a multi-objective optimization algorithm used to solve the overall optimal solution in multi-dimensional space. The evolution process of the DE algorithm mainly includes population initialization, mutation, crossover and selection and its evolution process is shown in Figure 1, the population is initialized at first, NP D-dimensional random variables are generated, and the fitness function value of each vector is calculated. If the conditions are met, the optimal parameter value is obtained. Otherwise, the mutation, crossover, and selection operations are continued. For the newly

Figure 1. Flow chart of adaptive differential evolution algorithm.

generated vector Calculate the value of the fitness function again, repeat these steps until the fitness function value satisfies the conditions, then obtain the optimal parameter value. After improving the mutation and crossover operations in the evolution process, the operation steps of the SaDE algorithm are as follows:

1) Population initialization

In the D-dimensional space, NP real-valued parameter vectors are randomly generated as the population of each generation, and the general initial population fits the uniform probability distribution. Let $X_{i,j}^G\left(i=1,2,\cdots ,NP;j=1,2,\cdots ,D\right)$ denote the value of the jth dimension of the ith vector.

2) Mutation operation

The basic mutation vector is generated according to formula (14):

$H_i^{G+1}=X_1^G+F\times \left(X_2^G-X_3^G\right)$. (14)

In Equation (14), $H_i^{G+1}$ is the ith vector in the next generation, $F\in \left[0,2\right]$ is the scaling factor, also known as the mutation operator, which is used to control the amplification of the deviation vector, and this factor can reflect the global optimization ability of the algorithm. The smaller the F value, the better the local searching ability of the algorithm. The larger the F value, the more the fitness function value can jump out of the local minimum point, and the slower the convergence speed. $X_1^G,X_2^G,X_3^G$ are the three different individuals in the population whose fitness function value is optimal.

In the basic differential evolution algorithm, the mutation operator often takes a constant, but its value is difficult to determine accurately. If the mutation rate is too large, the global optimal solution will be low; If the value of mutation rate is too small, the diversity of the population will decline, and it is easy to appear the phenomenon of “premature”. Therefore, the improved adaptive mutation operator $\lambda$ is adopted in this paper:

$\lambda ={\text{e}}^{1-\frac{G_m}{G_m+1-G}},F=F_0\times 2^{\lambda }$. (15)

In Equation (15), F₀ is the initial mutation operator, G_m is the current evolution algebra, and G is the maximum evolution algebra. The mutation operator of this mutation form is 2F₀ at the beginning, which can maintain the diversity of the population and prevent premature maturity. With the development of evolution, the mutation operator is gradually reduced toF₀, which can effectively avoid the destruction of the optimal solution.

3) Crossover operation

In order to increase the diversity of vectors in the population, the following crossover operation is introduced:

$\begin{array}{l}{V}_{i,j}^{G+1}=\{\begin{array}{l}{H}_{i,j}^{G+1},rand\left(0,1\right)\le CR\text{or}j=rand\left(1,D\right)\\ X_{i,j}^{G+1},\text{else}\end{array}\\ \left(i=1,2,\cdots ,NP;j=1,2,\cdots ,D\right)\end{array}$ (16)

In Equation (16), $V_{i,j}^{G+1}$ is the value of the jth dimension of the ith vector in the next generation, $rand\left(0,1\right)$ means that random numbers are generated between [0, 1], $CR\in \left[0,1\right]$ is a crossover operator. The larger the value of CR, the faster the convergence speed of the algorithm. This paper adopts the crossover operator of random range as follows:

$CR=0.5\left[1+rand\left(0,1\right)\right]$. (17)

This method can keep the mean value of the crossover operator at about 0.75, which ensures the diversity of the population.

4) Selection operation

In order to determine whether the vector in the population can become a member of the next generation, we need to compare the test vector with the current target vector, and calculate the fitness function $F\left(k\right)$ value of each vector, then the following selection operation is performed:

$X_i^{G+1}=\{\begin{array}{l}{V}_i^{G+1},\text{if}F\left(V_i^{G+1}\right)<F\left(X_i^G\right)\\ X_i^G,\text{otherwise}\end{array}$ (18)

The fitness function selected in this paper is:

$F\left(k\right)=\frac{1}{n}{\displaystyle \underset{i=1}{\overset{n}{\sum }}{\left(C_i^y-C_i^M\right)}^2}$. (19)

In Equation (19), n is the sample size, $C_i^y,C_i^M$ are the option price predicted by the ith data and the actual option price respectively, $F\left(k\right)$ is the fitness value of the kth individual, and the vector with the minimum fitness function value will appear in the next generation emerges.

5) Boundary condition processing

The vectors beyond the bounds are replaced by randomly generated parameter vectors in the feasible region.

4. The Relationship between Hurst Index and Asset Price

In this paper, the Hurst index H, which measures the degree of long memory, is introduced into the stochastic volatility option pricing model, and then studies the distribution characteristics of the asset returns subscripted by the Hurst-Heston model withH. By observing the changes of the probability density curve of stock price returns under different values of H, we can more intuitively understand the relationship between Hurst index and stock price changes in the model proposed in this paper.

Theorem 3: Under the assumption of Definition 2, the underlying asset price approximately obeys the following distribution:

$\ln S_t~N\left[\ln S_0+rt-H{t}^{2H}\theta +\left(V_0-\theta \right)H{t}^{2H}{e}^{-\kappa t},{\sigma }^2{V}_0\frac{1-e^{-2\kappa t}}{2\kappa }{t}^{2H}\right]$. (20)

Prove: According to Equation (2), Equation (3) and theorem 1, we can get the following formula:

$\begin{array}{c}\text{d}\ln S_t=\frac{1}{S_t}d{S}_t+\frac{1}{2}\left(-\frac{1}{S_t^2}\right){\left(d{S}_t\right)}^2\\ =\frac{1}{S_t}{S}_t\left[\mu dt+\sqrt{V_t}d{B}_H^1\left(t\right)\right]+\frac{1}{2}\left(-\frac{1}{S_t^2}\right)S_t^2{\left[\mu dt+\sqrt{V_t}d{B}_H^1\left(t\right)\right]}^2\\ =\left(\mu -H{t}^{2H-1}{V}_t\right)dt+\sqrt{V_t}d{B}_H^1\left(t\right).\end{array}$

So $\ln S_t=\ln S_0+{\displaystyle {\int }_0^{t}d\ln S_u}=\ln S_0+\left(\mu -H{t}^{2H-1}{V}_t\right)t+\sqrt{V_t{t}^{2H}}$. (21)

According to Equation (3), let $f\left(V_t\right)=V_t{e}^{\kappa t}$, then from Theorem 1 we can obtain:

$\begin{array}{c}df\left(V_t,t\right)=\left(\frac{\partial f}{\partial t}+\frac{\partial f}{\partial V_t}\kappa \left(\theta -V_t\right)+\frac{{\partial }^{2}f}{\partial V_t^2}{\sigma }^2{V}_{t}H{t}^{2H-1}\right)dt+\frac{\partial f}{\partial V_t}\sigma \sqrt{V_t}d{B}_H^2\left(t\right)\\ =\kappa \theta e^{\kappa t}dt+\sigma \sqrt{V_t}{e}^{\kappa t}d{B}_H^2\left(t\right).\end{array}$

So:

$V_t=\left(1-e^{-\kappa t}\right)\theta +V_0+{\displaystyle {\int }_0^t\sigma \sqrt{V_t}}{e}^{-\kappa \left(t-s\right)}d{B}_H^2\left(t\right)$. (22)

According to Equation (21) and Equation (22), it can be deduced that the underlying asset price S_t approximately obeys the following distribution (Theorem 3):

$\ln S_t~N\left[\ln S_0+rt-H{t}^{2H}\theta +\left(V_0-\theta \right)H{t}^{2H}{e}^{-\kappa t},{\sigma }^2{V}_0\frac{1-e^{-2\kappa t}}{2\kappa }{t}^{2H}\right]$.

By using Theorem 3, basic parameters $S_0=2$, $r=0.02$, $\theta =2$, $V_0=0.3$, $\kappa =1.5$, $\sigma =2.5$ are selected to obtain density curves of asset prices at different maturity times, as shown in Figure 2. It can be seen from Figure 2(a) that when $t=0.5$, the probability density curve gradually shows a sharp peak as the value of H increases; from Figure 2(b), it can be seen that when $t=0.25$, the probability density curve shifts to the right as the value of H increases, and the phenomenon of sharp peak becomes more and more obvious. Mathematically, the distribution of this characteristic can be described by the Levy distribution, and two important parameters $\alpha$ (depicting kurtosis) and $\beta$ (descriptive skewness) are used to describe the distribution. When the distribution of a random variable satisfies power-law attenuation and $\alpha <2$, the distribution of the random variable shows a sharp spike. For the return series $H>0.5$ with long memory, and the relationship between $\alpha$ and Hurst index H is $\alpha =1/H$, then the parameter $\alpha <2$ in the return distribution. Thus, it can be seen from this

(a)(b)

Figure 2. Probability density curve of underlying yield. (a) $t=0.5$ (b) $t=0.25$.

that the yield series with long memory will show a peak distribution, and as the value of H gradually increases, the peak phenomenon observed in the yield distribution will become more and more obvious. In the real financial market, the yield of the underlying assets has the characteristics of peak and thick tail. In the “peak” period, the market is stable as a whole, and the price trend conforms to the technical analysis and statistical laws; but in the “thick tail” period, the market changed dramatically due to the impact of special events. It can be seen that the larger the value of the Hurst index, which measures long-term memory, the more it can indicate that the market is stable at this stage, and the “peak” state of the stock return density curve is more obvious. In conclusion, it can be verified that the change process of the underlying assets under the assumption of the Hurst-Heston model is closer to the real situation, which can better describe the changing trend of stock prices in the real financial market.

5. Empirical Analysis

Before pricing options using the pricing model, it is necessary to estimate the parameters in the model, and the implicit parameter estimation method will be adopted in this paper [28]. This method uses the actual price of the market to infer the model parameter value, and the parameter value obtained is the parameter value under the risk neutral probability measure. This method requires a relatively small amount of data, and the estimated value is more effective.

In order to estimate the parameters of the Hurst-Heston model, this paper selects the daily closing data of the SSE 50ETF put option on August 26, 2021 and August 27, 2021, and uses the SaDE algorithm to estimate the model parameters, all data in this paper comes from https://option.eastmoney.com/. A total of 98 contracts were traded over two days, with expiry dates in September, October, December and March. The closing prices of the SSE 50ETF in the past two days were 3.128 yuan and 3.159 yuan respectively, and the execution prices ranged from 2.85 yuan to 4.4 yuan. In this paper, the risk-free rate selects the daily yield of 1-year Treasury bonds. In the algorithm optimization process, the initial population number is set as 60, the maximum evolutionary generation is 100, the initial mutation operator $F_0=0.4$, and the parameter value ranges are $v_0\in \left(0,1\right)$, $\kappa \in \left(0,1\right)$, $\theta \in \left(0,2\right)$, $\sigma \in \left(0,1\right)$, $\rho \in \left(-1,1\right)$, $H\in \left(0.5,1\right)$ respectively. In the process of parameter optimization of Hurst-Heston model, the optimal value of fitness function is 1.7133 × 10⁻⁵, and the estimated results of model parameters $v_0,\kappa ,\theta ,\sigma ,\rho ,H$ are [0.04608806, 0.36575, 0.32713226, 0.99363029, 0.69351875, 0.60894724] respectively.

The parameter optimization result of the model by SaDE algorithm can be obtained with the value of H is 0.60894724, $H>0.5$ indicating that the fluctuation of SSE Composite index is not a random walk, but a partial random walk, that is, it has long memory. Using the above parameters, the Hurst-Heston model is used to predict the option price of the next 10 trading days from August 30 to September 10. The expiration date of the trading contract on each trading day includes September, October, December 2021 and March 2022, and the strike price includes a variety of different prices ranging from 2.85 yuan to 4.4 yuan. Then calculated the Mean Relative Error (MRE), Mean Absolute Error (MAE), Mean Square Error (MSE), and Root Mean Square Error (RMSE) between the daily valuation and the actual market price as the model evaluation indicators, the results were shown in Table 1.

As can be seen from the results in Table 1 that the value of MRE fluctuates around 0.04 in the first 5 days, indicating that the daily error of the parameters obtained by the implicit parameter estimation method when used to estimate the option price in the next five trading days is relatively close, and the fluctuation range is small, indicating that the model is relatively stable. From the sixth day, the MRE value of the predicted value becomes larger. If this group of parameters is used to predict the option price of the next 10 trading days, the longer the interval is, the greater the error will be. This result indicates that the parameters of the pricing model are dynamic, and the related parameters are updated daily in the financial market for this reason, so the same group of parameters cannot be used for long-term prediction.

In order to evaluate the effect of fractional stochastic volatility model in option pricing more intuitively, we used the closing data of SSE 50ETF put options on September 10, 2021 to make pricing predictions. There were 50 trading contracts on that day, the closing price of SSE 50ETF was 3.292, and the 1-year Treasury bond yield was $r=0.022861$ on that day. The pricing results of this model are compared with B-S model, Heston model and fractional B-S model with long-memory features. The daily closing data of SSE 50ETF put options on September 8, 2021 and September 9, 2021 are selected to optimize the parameters of each model using SaDE algorithm. Each iteration of parameter optimization needs to solve the problem of fitness function minimization. The iterative curves of the Hurst-Heston model and Heston model in the optimization process are shown in Figure 3, and the final parameter estimation results are shown in

Table 1. Daily error comparison of HH model.

(a)(b)

Figure 3. Parametric optimization iteration curve. (a) HH model (b) Heston model.

Table 2. The prediction results of each model and the errors with the actual market price are shown in Table 3, and the absolute error comparison diagram of the three models is shown in Figure 4.

In Table 2, HH represents the Hurst-Heston model proposed in this paper, and BSH represents the fractional B-S model. In Table 3, SE_HH, SE_BS, SE_HS and SE_BSH respectively represent the square error between the pricing results of the Hurst-Heston model, B-S model, Heston model, fractional B-S model and the actual market price. It can be seen from Table 3 that the MSE of the pricing result of the Hurst-Heston model and the actual market price is 0.00002, the MSE of B-S model, Heston model and fractional B-S model are 0.0003, 0.00004 and 0.00008 respectively. In addition, the Mean Relative Error of fractional random

Table 2. Model parameters.

Table 3. Comparison of model pricing results.

Figure 4. Absolute error of model pricing results.

volatility model and market price is 0.1267, the Mean Relative Error of B-S model, Heston model and fractional B-S model in pricing are 0.2095, 0.1776 and 0.2264 respectively. By comparison, it is obvious that the error of the Hurst-Heston model is smaller than that of the other three basic models. Due to the limited space, daily forecast data cannot be listed one by one. But after many experiments, the results showed that the pricing error of the Hurst-Heston model proposed in this paper is smaller than that of other comparison models, indicating that this model has better pricing effect and higher accuracy.

It can be observed from Figure 4 that the absolute error of Hurst-Heston model proposed in this paper for the prediction of different option contracts are all smaller than that of the other three models, and the fluctuation is relatively stable, among which the error fluctuation range of B-S model is the largest. In addition, the MAE between the option price predicted by Hurst-Heston model optimized based on SaDE algorithm and the actual market price is 0.00298, the pricing result errors of B-S model, fractional B-S model and Heston model are 0.01262, 0.005133 and 0.00748 respectively. It is obvious that the MAE of this model is smaller than that of B-S model, fractional B-S model and Heston model, it shows that the option pricing results of the Hurst-Heston model optimized by SaDE algorithm is more accurate.

6. Conclusion

In this paper, the Hurst index H is introduced on the basis of the Heston option pricing model, and the Hurst-Heston European option pricing model with long memory characteristics is constructed, then obtained Hurst-Heston option pricing formula by solving the partial differential equation under the model by Fourier transform method. In addition, through the distribution function of the underlying asset price, the influence of H on the asset price is analyzed. When the value of H is between 0.5 and 1, as the value of H increases, the probability density curve of S_t gradually presents the phenomenon of peak and thick tail, which is more in line with the reality. In the empirical analysis, the sample data is used to estimate the model parameters, and then the option prices in the next few trading days are predicted. The prediction results of Hurst-Heston model are compared with those of Heston model, B-S model and fractional B-S model. By comparing the absolute error, Mean Square Error and other evaluation indicators between each model and the market price, it can be seen that the error of Hurst-Heston model is relatively smaller. The empirical results showed that the Hurst-Heston model is reasonable and feasible in European option pricing, and the pricing results are more accurate. The Hurst-Heston pricing model can not only describe the peak and thick tail nature of the underlying asset rate of return, but also reflect the long-memory characteristics of the financial market. The change process of the underlying asset depicted by this model is closer to the real situation and can better describe the stock price change trend in the real financial market. In this paper, we assumed that the stock price fluctuates randomly with a single factor in the process of research, but in practice, the stock price fluctuation can be affected by a variety of factors, so we can consider the situation that the stock price change is a multi-factor random fluctuation in the follow-up research. In the empirical analysis, stock option data are used to verify the model, and the model has not been applied to other financial markets, so we can try to extend the model to real option pricing in the future.

Conflicts of Interest

The authors declare no conflicts of interest.

References

[1]	Qi, Y., Sun, J. and Li, F. (2019) Research Progress of Digital Option Theory. Economic Dynamics, No. 5, 119-134.
[2]	Black, F. and Myron, S. (1973) The Pricing of Options and Corporate Liabilities. Journal of Political Economy, 81, 639-654. https://doi.org/10.1086/260062
[3]	Zhang, J., Wang, Y. and Zhang, S. (2022) A New Homotopy Transformation Method for Solving the Fuzzy Fractional Black-Scholes European Option Pricing Equations under the Concept of Granular Differentiability. Fractal and Fractional, 6, Article No. 286. https://doi.org/10.3390/fractalfract6060286
[4]	He, X.-J. and Lin, S. (2021) A Fractional Black-Scholes Model with Stochastic Volatility and European Option Pricing. Expert Systems with Applications, 178, Artilce No. 114983. https://doi.org/10.1016/j.eswa.2021.114983
[5]	Merton, R.C. (1976) Option Pricing When Underlying Process of Stock Returns Are Discontinuous. Journal of Financial Economics, 3, 124-144. https://doi.org/10.1016/0304-405X(76)90022-2
[6]	Leland, H. (1985) Option Pricing and Replication with Transactions Costs. The Journal of Finance, 40, 1283-1301. https://doi.org/10.1111/j.1540-6261.1985.tb02383.x
[7]	Heston, S. (1993) A Closed-Form Solution for Options with Stochastic Volatility with Applications to Bond and Currency Options. The Review of Financial Studies, 6, 327-343. https://doi.org/10.1093/rfs/6.2.327
[8]	Li, P., Yang, J. and Lin, Y. (2020) Parameter Estimation and Application of Fast Mean Reversion Stochastic Volatility Model. Operations Research and Management Science, 29, 137-143.
[9]	He, J. and Wei, Z. (2019) The Pricing of Carry Options Based on Jump Diffusion Model of Non-Affine Stochastic Volatility. Mathematics in Practice and Theory, 49, 132-139.
[10]	Wu, X., Li, X. and Mam C, (2019) Research on SSE 50ETF Option Pricing Based on Stochastic Volatility Model. Journal of Applied Statistics and Management, 38, 115-131.
[11]	SenGupta, I. (2014) Option Pricing with Transaction Costs and Stochastic Interest Rate. Applied Mathematical Finance, 21, 399-416. https://doi.org/10.1080/1350486X.2014.881263
[12]	Zhou, R. (2018) Option Pricing Based on Volatility Decomposition. Systems Engineering Theory and Practice, 38, 1919-1929.
[13]	Li, D. (2017) Research on Long Memory of Financial Market Based on Fractal Method. Ph.D. Thesis, Foreign Economic and Trade University, China.
[14]	Yuan, Y., Zhuang, X. and Jin, X. (2016) Intraday Effect, Long Memory and Multifractality of Chinese Stock Market: Based on the Perspective of Price-Quantity Cross-Correlation. Journal of Systems Management, 25, 28-35.
[15]	Geng, K. and Zhang, S. (2008) Research on Long Memory of UHF Duration Sequence in Chinese Stock Market. Chinese Journal of Management Science, No. 2, 7-13.
[16]	Xue, C. and Li, X. (2008) Long-Memory Characteristics of China’s Stock Market Returns and Volatility. Statistics and Decision, No. 4, 117-119.
[17]	Yuan, Y. and Zhuang, X. (2008) The Long Memory and Market Development Status of China’s Stock Market. Journal of Applied Statistics and Management, No. 1, 156-163.
[18]	Li, Y., He, Y. and Zhu, H. (2003) An Empirical Study on China’s Stock Market Return and Volatility Long Memory. Systems Engineering Theory and Practice, No. 1, 9-15.
[19]	Shi, Y., Tao, X. and Zhang, S. (2011) On-Off Hurst Index Fractal Black-Scholes Market Pricing of European Options. Economic Mathematics, 28, 40-44.
[20]	Yang, Z. and Caselli, F. (2017) Optimal Exercise Boundary of American Fractional Lookback Option in a Mixed Jump-Diffusion Fractional Brownian Motion Environment. Mathematical Problems in Engineering, 2017, Article ID: 5904125. https://doi.org/10.1155/2017/5904125
[21]	Qin, X., Lin, X. and Wang, W. (2019) Research on Fuzzy Pricing of European Options Based on Long Memory Characteristics. Systems Engineering Theory and Practice, 39, 3073-3083.
[22]	Xue, D., Zhu, G., Zhu, Y. and Xiong, Y. (1995) Stationarity Analysis and Description of Fractional Brownian Motion. Journal of Huazhong University of Science and Technology, 10, 78-81.
[23]	Huang, Y. (2009) Time Axis Transformation Method for Option Pricing under Geometric Fractional Brownian Motion. Statistics and Decision, No. 5, 10-12.
[24]	Fadugba, S.E. and Nwozo, C.R. (2016) Valuation of European Call Options via the Fast Fourier Transform and the Improved Mellin Transform. Journal of Mathematical Finance, 6, 338-359. https://doi.org/10.4236/jmf.2016.62028
[25]	Ding, Q. and Yin, X. (2017) Overview of Differential Evolution Algorithms. Journal of Intelligent Systems, 12, 431-442.
[26]	Yang, Q., Cai, L. and Xue, Y. (2008) Overview of Differential Evolution Algorithms. Pattern Recognition and Artificial Intelligence, 21, 506-513.
[27]	Hu, F., Dong, Q. and Lu, L. (2021) Improved Differential Evolution Algorithm for Adaptive Quadratic Mutation and Its Application. Computer Applications and Software, 38, 271-280.
[28]	Bakshi, G., Cao, C. and Chen, Z. (1997) Empirical Performance of Alternative Option Pricing Models. The Journal of Finance, 52, 2003-2049. https://doi.org/10.1111/j.1540-6261.1997.tb02749.x