An analytical approximation formula for European option prices under a liquidity-adjusted non-a ﬃ ne stochastic volatility model

: In this paper, we investigate the pricing of European options under a liquidity-adjusted non-a ﬃ ne stochastic volatility model. An analytical European option pricing formula is successfully derived with the COS method, based on an approximation for the characteristic function of the underlying log-asset price. Numerical analysis reveals that our results are very e ﬃ cient and of reasonable accuracy, and we also present some sensitivity analysis to demonstrate the e ﬀ ects of di ﬀ erent parameters on option prices.


Introduction
The Black-Scholes (BS) model [1] has been widely adopted in the pricing of financial derivatives because of its analytical tractability. However, it is generally known that the assumptions of BS model are not appropriate with empirical results. Therefore, many researchers devote much effort to remedying these shortcomings. A great deal of scholars proposed to use the generalized Black-Scholes model to describe the dynamic process of the underlying asset and solved the generalised Black-Scholes model numerically. For example, Mohammadi [1] developed a Quintic B-spline collocation approach for solving generalized Black-Scholes partial differential equation for European option prices. Roul [3] presented a high-order compact finite difference method based on a uniform mesh to obtain a highly accurate result for the generalized Black-Scholes equation. Roul [4] constructed a new computational approach for solving the generalized Black-Scholes equation numerically by combining Crank-Nicolson scheme and sextic B-spline collocation method. Another important improvement of the BS model is to make it a fractional one since it can describe the long-range dependence displayed by financial data. It should be noted that the analytical solution to option prices is generally not available under the framework of the fractional BS model, and numerical algorithms have to be used. Roul [5] designed a high-order numerical approach for solving the fractional BS equation. Roul and Goura [6] presented a numerical technique for solving the fractional BS equation describing European options.
Another main category that is very popular in modifying the BS model is to include stochastic volatility. A breakthrough was made by Heston [7] , who proposed a model that not only possesses a range of basic properties including mean reversion and nonnegativity, but also relieve the volatility smile or smirk. Due to the popularity of the Heston model, a number of authors have worked under this framework. For example, Lee et al. [8] obtained an analytic formula for vulnerable options, while Mollapourasl et al. [9] used a radial basis function method to price American put options. He and Zhu [10] went further to introduce a stochastic interest rate into the Heston model and presented a closed-form formula for European options. In particular, the non-affine stochastic volatility model has attracted the extensive attention as an improvement of the Heston model. Christoffersen et al. [11] and Chourdakis [12] further confirmed that the non-affine volatility model represents the time series features of real option price data, which is much better than Heston model. This has prompted a number of authors to consider the non-affine volatility model. For example, Huang and Guo [13] investigated the valuation of discrete barrier options with non-affine volatility and double exponential jump, and Ma et al. [14] studied the pricing problem of vulnerable options with non-affine volatility and stochastic interest rate.
On the other hand, it needs to be pointed out that the market is assumed to be completely liquid in the above literature, whereas a lot of empirical evidence has already suggested that the market is actually not completely liquid [15,16]. Therefore, the concept of the market liquidity starts to appear in the option pricing theory. In particular, Feng et al. [17] derived an analytical solution to European option prices by incorporating a liquidity discounting factor into the option pricing model, based on which Li and Zhang [18] verified the influence of the market liquidity factor using Shanghai 50ETF options. Moreover, Li and Zhang [19] studied the pricing of European quanto options when the underlying foreign asset is in an imperfectly liquid market, while Xu et al. [20] focused on the variance and volatility swap pricing in the framework of a liquidity-adjusted underlying assets model. The pricing of vulnerable options is also considered by Wang [21] when the liquidity risk is assumed to be captured by a stochastic process. These clearly demonstrate that the effect of liquidity risk on option pricing is not negligible.
Motivated by the above research, we consider the pricing of European options under a liquidityadjusted non-affine volatility model. The advantage of this model is that it can not only describe the nonlinear characteristics of financial data, but also take the liquidity risk into account. Although it is difficult to derive a closed-form solution for the characteristic function of the underlying log-asset price, we present an analytical approximation to the characteristic function using the first-order Taylor expansion method. Based on this, an approximate price for a European option can be obtained by the COS method [22], which is recently very popular in the area of option pricing [23][24][25][26].
The rest of this article is arranged as follows. A liquidity-adjusted non-affine volatility model is briefly introduced in the Section 2. In Section 3, an approximate formula for the characteristic function of the underlying log-price is derived, based on which an analytical formula for European option prices is presented with the COS method. Some numerical examples and sensitivity analysis are provided to interpret our results in Section 4, after which we conclude in the last section.

Model specification
Let {Ω, F t , Q} denote a given complete probability space, and Q is a risk-neutral probability measure. Considering that the non-affine model can explain the nonlinear characteristics of option data in real market, and liquidity-adjusted models perform better than those without considering liquidity risk, we construct a liquidity-adjusted non-affine volatility model for the valuation of European options. Specifically, the price process S t and the process of the volatility v t under the risk-neutral world are modeled as follows: where dW 1t dW 2t = ρd t , r is the risk-free interest rate, l is the level of market liquidity, β is the sensitivity of the stock to the level of market liquidity, k donates the mean-reversion speed, θ is the long-term mean of the volatility, and σ is the instantaneous volatility of the volatility. With x t = ln(S t /K), Eqs (2.1) and (2.2) can be respectively transformed into

An analytical approximation formula
In this section, we approximate the characteristic function of the underlying log-price with the firstorder Taylor expansion and based on this, European options can be analytically valuated using the COS method.

Derivation of the characteristic function
The characteristic function of x T is defined as where T is the maturity time of the option, i is an imaginary unit, and τ = T − t. The solution to this characteristic function is presented in the following theorem. Theorem 1. If the price process and its volatility process follow the dynamics specified in Eqs (2.1) and (2.2), then the characteristic function of x T can be expressed by Proof. The Feynman-Kac theorem indicates that Φ(x, v, τ; u) satisfies the following PDE (partial differential equation) with the boundary condition given by Given that the volatility in the proposed model is of mean reversion, it is reasonable to consider the expansion of v 2 and v 3 2 at the long-term mean θ using a first-order Taylor expansion. This gives Substituting Eqs (3.2) and (3.3) into the PDE (3.1) yields If we assume that the solution to this PDE can be written in the form of with the boundary conditions A(u, 0) = B(u, 0) = 0, and substituting Eq (3.5) into Eq (3.4), we obtain As a result, the following two ordinary differential equations can be derived Since Eq (3.6) is a Riccati equation, it can be solved with some algebraic computation so that B can be formulated as Multiplying both sides of Eq (3.6) by θ/2, and substituting the result into Eq (3.7), the following expression can be derived: which further leads to This completes the proof.

European option pricing
In this subsection, an analytical approximate price formula for the European options are obtained by the COS method, the results of which are summarized in the following theorem. Theorem 2. Under the model dynamics (2.1) and (2.2), the pricing formula of European call option is expressed by: and the χ ( x 1 , x 2 ), ψ ( x 1 , x 2 ) is defined in Eqs (3.11) and (3.12).
Proof. It is well known that the price of a European option can be expressed as follows where implies the summation whose first term is divided by two, Re{·} denotes taking the real part of a complex number, and Φ(u) is the conditional characteristic function of f (y|x).
Replacing f (y|x) in (3.8) with its approximation (3.9), and interchanging summation and integration, we can derive the the price of a European option: where Clearly, the only unknown term in the pricing formula is the coefficient V k , deriving which yields a completely analytical solution. In particular, for a European call option, we have and (3.12) After some algebraic calculation, we can arrive at where ψ k and χ k are given by Eqs (3.11) and (3.12), respectively. This completes the proof. By comparing the properties of European call and put options under the same model setup, we can find that the derivation of European put options is almost similar to the proof of Theorem 2, except for the calculation of the coefficient V put (a 1 , b 1 ). Therefore, without going to details, the corresponding European put option pricing formula is presented in the following corollary. Corollary 1. Under the model dynamics (2.1) and (2.2), the pricing formula of European put option can be presented as: , 0)) .

Numerical analysis
Once an approximation formula is successfully derived, it is important to check its accuracy. Moreover, it is also interesting to investigate the influence of liquidity on option prices. Both of these will be addressed in this section.

Accuracy tests
In this subsection, numerical experiments are used to verify the pricing performance of the proposed method against the Monte Carlo simulation and the FFT method [14,27,28]. Without loss of the generality, we will only take European call options as an example. where a 0 = ln S 0 , L = 10, and c n is the n-cumulant of ln S T , whose computing formula is: We set the number of simulation paths 100, 000, the number of time steps is chosen to be 252, N = 2 10 , and the damping factor is 2.55. The values of the remaining parameters * for numerical examples are displayed in Table 1. All of the numerical examples were implemented in Matlab R2020a on a PC with an Intel Core i5 CPU. Table 2 presents our approximate prices, FFT and Monte Carlo simulation results of European call options with different strike prices and expiry times. The relative error (R.E. COS ) between the COS method and Monte Carlo simulation is calculated as: where P cos and P mc denote the prices of European options computed from Theorem 2 and Monte Carlo simulation. Similarly, the relative error (R.E. FFT ) between the FFT and Monte Carlo simulation is defined as: where P FFT denotes the price of European options computed from the FFT method. It is not difficult to find that R.E. COS and R.E. FFT are both less than 2%, while R.E. COS is no greater than (R.E. FFT in all cases. In addition, Table 2 also indicates that the CPU time to calculate the prices of European options using the different methods. We can observe that the Monte Carlo simulation generally takes about three times as long as our approach, and the FFT method is also slower than our approach. All of these clearly reveal that our approach is quite accurate and efficient.  The relative error between the COS method and Monte Carlo simulation with respect to the initial price of the underlying asset S 0 and the maturity time T − t is also shown in Figure 1 using a surface plot. One can see clearly that the relative error associated with in-the-money and at-the-money is larger than that of out-of-the money options. Moreover, the shorter the time, the larger the error. This encourages us to find more effective method to price short-tenor options in future research.

Sensitivity analysis
In this subsection, we will analyze the impact of market liquidity and stochastic volatility on European call option prices; such an impact will be studied through the influence of the following parameters: (1) the time varying level of market liquidity l t and the maturity time T −t, (2) the sensitivity of market liquidity β and the maturity time T − t, (3) long-term mean θ and risk-free interest rate r, (4) the initial price of the underlying asset S 0 and the maturity time T − t. Figure 2 shows the price of European call options with respect to the market liquidity level l t at different expiry times. Clearly, European call option price is an increasing function of the the maturity time T − t. In addition, a higher market liquidity level l t leads to higher European call option price. This is mainly because a large value of l t imply a greater level of illiquidity, which contributes to a higher underlying price and thus a larger option premium.
What are shown in Figure 3 are European call option prices with respect to the liquidity sensitivity factor β and the maturity time T − t, and the sensitivity parameter is shown to have a positive impact on option prices. This can be understood by the fact that for the same level of market liquidity, the underlying price will be affected more heavily by the market liquidity if it is more sensitive to the market liquidity, which in turn contributes to higher option prices.
Depicted in Figure 4 is how European call option price changes with respect to the long-term average of the volatility θ and risk-free interest rate r. Clearly, European call option prices are an increasing function of the risk-free interest rate, since raising the interest rate is equivalent to increasing the level of the underlying asset price. One can clearly observe that the price of European call options increases with the long-term average of volatility θ, which verifies that the price of the European call option is an increasing function of volatility. This is actually reasonable since a higher level of the long-term average of volatility tends to increase the long-term volatility level, which means that the uncertainty of the underlying asset prices increases, leading to a larger intrinsic value of the option and a higher option price.    and as we can observe, a higher level of the underlying asset price tends to increase European call option prices. The reason behind this is that a greater underlying asset price will raise the intrinsic value of the option and results in higher option prices.

Conclusions
This paper studies the pricing problem of European options under the non-affine volatility model in an illiquid market, so that the effect of market liquidity on option prices can be captured. Upon incorporating the market liquidity into the underlying dynamics, we present an approximation to the characteristic function of the underlying log-price, based on which an analytical pricing formula for European options is successfully derived using the COS method. Finally, the sensitivity analysis are used to demonstrate the validity of the model setup.