Liquidity effects on oil volatility forecasting: From fintech perspective

Fin-tech is an emerging field, inspiring revolutionary innovations in the financial field. It may initiate the evolutionary episode of the financial research, where volatility forecasting is a crucial topic in finance. For forecasting volatility, GARCH model is a prevailing model, however, further improvement of the GARCH model is still challenging. In this paper, we demonstrate how Fintech can play a part in volatility forecasting by employing a metaheuristic procedure called Genetic Programming. On the basis, we are able to develop a new volatility forecasting model, which can beat GARCH family models (including GARCH, IGARCH and TGARCH models) in a significant way. Since genetic programming is an evolutionary algorithm based on the principles of natural selection, this innovative work will be a breakthrough point in the financial area. The innovation of this paper demonstrates how GP technology can be applied in the financial field, attempting to explore the volatility forecasting area from the combination of new technology and finance, known as fintech. More importantly, when the formula of volatility forecasting is unknown as we introduce a new factor, namely, the liquidity factor, we unveil that how GP method can be helpful in determining the specific volatility forecasting model format. We thereby exhibit the liquidity effects on volatility forecasting filed from the fintech perspective.


Motivations and aims
Fin-tech is an emerging field, which drives revolutionary innovations in the financial spectrum recently. It integrates advanced technology into financial area, provoking the evolutionary episode of the financial industry and research (Buchak et al. 2018 [1]; Chen et al. 2019 [2]). Fintech has been playing particularly important roles in two areas of financial innovations: derivative trading and high frequency trading, which demand sophistical mathematical models, computing-intensive algorithms and big data handling techniques in order to capture trading signal in microseconds. In the option markets, volatility is the dominant driving factor for option price moves. In high frequency trading, volatility forecasting is essential to manage portfolio risk. Therefore, volatility forecasting would be a key issue in the financial field and  [3]). GARCH models are prevalent in forecasting volatilities of financial markets (see Andersen et al. 2005 [4]; Tian and Hamori [5]) and oil market volatility has been intensively examined and forecasted by adopting GARCH models (see Wei et al. 2010 [6]; Efimova and Serletis 2014 [7]; Wang and Liu 2016 [8]; Xing and Wang 2019 [9]). However, existing GARCH type models are unable to capture either market micro structures, such as liquidity information or the high frequency trading signals because the volatility estimation precision from GARCH models depends on the accuracy of mean estimation models, which are in low frequency basis. Incorporating liquidity into GARCH models is necessary in order to estimate volatility more accurately, especially from high frequency trading.
In fact, the relationship between liquidity and volatility has been scrutinized in the existing literature. Fleming and Remolona (1999) [10] have investigated the relation between liquidity, volatility and public information in the US Treasury market. They declare that volatility and liquidity respond to the public information simultaneously. Feng, Hung and Wang (2014)'s [11] empirical results also support the liquidity impact on the volatility. More recently, Collin-Dufresne and Fos (2016) [12] explore the relations between liquidity and noise trading volatility. They affirm the significant impact of liquidity on noise trading volatility.
Volatility forecasting by using Artificial Intelligence (AI) and Genetic Programming (GP) can been witnessed in the literatures (see Yin et al. 2016; [13]; Ding et al. 2019 [14]; Weng et al. 2021 [15]; Mademlis and Dritsakis 2021 [16]). The goal of the paper is to adopt AI technologies to generate the best model which can comprise trading liquidity effect into volatility forecasting, and can be integrated into the existing fintech systems such as high frequency trading platforms and derivative trading platforms in the hedge fund industry.

Research contributions
The first contribution is that we adopt GP to identify the model format for forecasting oil volatility with liquidity. It is arguable that liquidity has a considerable effect on market volatility but how to integrate liquidity factor into the volatility forecasting model has no unanimous solution. Our GP system can identify the specific format of volatility forecasting with liquidity by two steps: (i) model formation search and (ii) model accuracy evolution.
In the first stage, we identify the potential form of volatility forecasting model. The GP is used to search the potential forms of the volatility forecasting model with integration of liquidity variable. GP holds the elegant characteristics that one can build the relevant performance criterion directly into the search procedure.
In the second stage, we apply Darwin's theory of evolution to model accuracy evolution. To be specific, we adopt crossover and mutation operations during our volatility forecasting model selection and development. The main advantage of crossover and mutation algorithms is that they can help models to evolve according to the historical data in order to minimize forecasting errors. As a result, the generated forecasting model can provide far better forecasting results after massive generations of model evolution.
Based on this process, our GP system can identify the most relevant terms that have effects on predicting oil volatility. The bid ask spread term has been included in the volatility forecasting model by GP system, which captures the liquidity effect as it provides the trading information of the financial market (Deuskar et al. 2011 [17]). Since the volatility forecasting model format with liquidity is unknown, our contribution relies on the specific format identification of volatility forecasting model with liquidity, which makes accurate forecasting and liquidity sensitivity analysis applicable for the oil market. Furthermore, our work is related heavily to existing literature where GP has been adopted in various financial areas. For example, Pimenta et al. (2018) [18] use GP to develop automated investing method in the stock market. More recently, Michell and Kristjanpoller (2020) [19] employ GP to develop trading rules in the US stock market. Ding et al. (2020) [20] apply GP to forecast future stock returns in different stock markets. Our paper extends the GP application to the oil market volatility forecasting, which is a crucial commodity market.
Forecasting oil price volatility is not only useful to commodity traders for speculating and hedging (Ma et al. 2019 [21]), but also important to economists for forecasting macroeconomic variables, such as inflation and industrial productions (Yin 2016 [22]; Chen et al. 2019 [23]). After the model was generated by our algorithm, we compare our model with GARCH models regarding the forecasting ability. The empirical results show that the in-sample performance of our model is more accurate than traditional GARCH models for both full-sample and most subsample predictions. Additionally, our model's out-of-sample volatility prediction is also more precise compared with traditional GARCH models for both full-sample and subsample tests. The accurate volatility prediction of our GARCH can be envisioned as a considerable improvement of the existing GARCH models and our results demonstrate that multivariable GARCH model can still deliver accurate volatility predictions.
Finally, our research also produces significant empirical contributions. We shed new insights on the concept of "Fin-tech". The term "artificial intelligence" shall be a pivotal ingredient in the fintech concept, which refers when a machine mimics cognitive functions that humans associate with other human minds, such as "learning" and "problem solving" (Russell and Norvig 2009 [24]). However, the prevailing recognition of fintech focuses on the application of technology into the financial service industry, such as online banking and mobile payment (Mackenzie 2015 [25]). This paper gives a further annotation, suggesting that fintech also involves how AI, such as genetic programming, can be integrated with financial modelling and this paper itself is a magnificent demonstration of finance with AI, which is a radical innovation in finance. This integration can provide deep insights of financial data analysis since all models are generated and evaluated purely based on financial data. Furthermore, our algorithm is also capable of dealing with big data as well as other forecasting problems, such as return forecasting. Consequently, this paper can serve as an impetus of developing other fintechs in hedge fund management including derivatives trading algorithms and high-frequency trading.
The remainder of the paper is structured as follows. In section 2, we describe the data and variables used in the paper. In section 3, we derive our GARCH model with the GP method. Section 4 presents the empirical results and the model performance in data fitting and volatility forecasting. Section 5 delivers the conclusions and implications of the paper.

The data
The GARCH(1,1) model has been widely used in the financial literatures for volatility forecasting (see Bollerslev and Wright 2001 [26]; Hansen and Lunde 2005 [27]), especially in the oil market (see Klein and Walther 2016 [28]; Lux et al. 2016 [29]). As a consequence, we also use our model to forecast the oil return volatility compared with the traditional GARCH model and GARCH family models (see Fig 1 for oil daily returns). We will use the univariate GARCH-class models, including GARCH, IGARCH and the GARCH model with asymmetric effects, namely, GJR-GARCH, as benchmark models for oil volatility predictions in our paper. The liquidity effect would be measured by the bid-ask spread (BAS), which has been demonstrated to be positively correlated with price volatility in financial markets (see Bollerslev and Melvin 1994 [30]; Haugom et al. 2014 [31]). We have not included the weekend effect in this paper, which is popular in the stock markets. Studies on the weekend and daily effect in oil futures market are scant since it has been documented that the weekend effect may not be obvious in oil futures market (Yu and Shin 2011 [32]). It might be that oil futures contract has time-to-maturity, which may have a stronger effect than the weekend effect (Geman and Kharoubi, 2008 [33]).
We obtain the WTI oil futures trading data from Thomson Datastream with daily frequency, from January 1, 2001 to December 31 2019, yielding a total of 4,752 trading days with 4,752 observations. Most literatures have adopt such database and the trading data is continuous, which would be favored for daily volatility forecasting. We use the period from January 1, 2001 to December 31, 2010 as the in-sample period and the period from January 1, 2011 to December 31, 2019 as the out-of-sample period. By "out-of-sample", we imply two things: Firstly, we use the previous year's data to estimate coefficients of the model, which forecasts this year's volatility, namely, we adopt one-year rolling window for the test.
For in-sample fitting, the estimation year and prediction year are the same. More importantly, we only use the period from January 1, 2001 to December 31, 2010 as the sample period for genetic programming to generate our model. We derive the theoretical model from genetic programming and we estimate all the model coefficients by regression. Then, we apply the generated volatility forecasting model to the period of January 1, 2011 to December 31, 2019 to testify the "out-of-sample" model performance. For the out-of-sample forecasting, we use previous year's data to forecast the volatility of this year. For example, we use the data of 2013 to estimate the model parameters and then use those parameters to forecast volatility of year 2014. Similarly, we use the data of 2014 to estimate the model parameters and then use those parameters to forecast volatility of year 2015. We thereby roll over the sample by one year.
Furthermore, we use two whole periods as the full-sample volatility forecasting test and we use one year within the whole period as the subsample forecasting tests. In addition, for both in-sample and out-of-sample tests, we use 1-day ahead prediction during the data period and the statistical test for error differences will be also employed.

Variable estimation
Because the specific computational technique we have adopted, it may not allow us to use two identical time series data even with different time lags. As a result, we adopt two different time series for modeling the variance. The first sequence is the daily simplified realized variance series. Since we adopt a relatively low frequency data, we use the low frequency daily simplified realized variance rolling average in one month to approximate the daily simplified realized variance(see Schwert 1989 [34]), which is denoted as SRV t and can be defined as: since the mean of daily realized return is near 0 with the realized return r t is defined as: (Andersen and Bollerslev 1998 [35]), where P t is the last day's settlement price and we enlarged the SRV measurement by 10,000 times and the return by 100 times to avoid extreme small values. The second sequence is the realized daily variance, which is denoted as RV t and we use the standard deviation of the oil futures daily return (Christensen and Prabhala 1998 [36]), which can be defined as: where We enlarged the measurement by 10,000 times to avoid extreme small values. For both realized daily variance and simplified realized daily variance, we use the onemonth rolling average approach and we use the updated real market data for each day prediction without iteration. We then define the simplified realized variance as the target series to be forecasted and we use other three series to fit the data, namely, Roll spread, ε t and realized variance where ε t is the residual taken from GARCH(1,1) model. As the BAS is positively related with the volatility, we use the Roll measure to approximate the daily BAS measurement. The Roll spread is an effective spread estimator for the bid-ask spread calculation and suits the low frequency data (Marshall et al. 2012 [37]; Haugom et al. 2014 [31]) and the Roll measure thereby could be severed as the market liquidity proxy. Specifically, Roll (1984) [38] assumes that the stocks have fundamental values, denoted as V t at time t and V t follows: η t is the residual term for this fundamental value random process, which differentiates from the term ε t , which is the residual taken from GARCH volatility model, and η t is the mean-zero, serially uncorrelated public information shock on day t. Next, he denotes S t as the last observed trade price on day t and presumes that S t follows: where E is the effective spread and Q t is a buy/sell indicator for the last trade that equals +1 for a buy and -1 for a sell. He further assumes that Q t is equally likely to be +1 or -1 and Q t is also serially uncorrelated, and independent of η t .
Then he takes the first difference of equation and plugs in the result from equation, which yields As a result, covðDS t ; DS tÀ 1 Þ ¼ À 1 4 E 2 or equivalently, spread ¼ 2 ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi À covðDS t ; DS tÀ 1 Þ p : Because when the auto-covariance is positive, the formula is undefined. We therefore use a modified version of the Roll estimator (Goyenko et al. 2009 [39]): ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi ffi In order to estimate the spread, we first take the difference of daily price, then we take the covariance of the price difference of t and t-1. Then, we take the square root of the negative covariance if it is greater than 0 and we take the value of 0 otherwise.
To sum up, the simplified realized variance, the realized variance, the BAS and the residual term are the four main variables concerned by this study and the residual term ε t contains the noise information that has not been captured by the model. Table 1 presents the statistical summary of the main variables used in the paper.

Preliminaries
In this section, we will develop our GARCH model based on the estimated variables in section 2. For the specific model development, we will adopt the Genetic Programming (GP) method from computer science, which is proposed by Koza (1994) [40]. GP is an evolutionary computation (EC) technique inspired by biological process Banzhaf et al. 1998 [41]; Hirsh et al. 2000 [42]; Poli et al. 2008 [43]). Since the form of volatility forecasting model with liquidity effects is uncertain, it would be beneficial to adopt GP method. One big advantage of adopting GP is that it can allow one to be agnostic about the general form of the model. In GP, a population of computer programs is evolved based on the principles of natural selection originated from Darwin's theory of evolution. After certain number of generations, GP can transform populations of programs into new and better programs. As stated in Poli et al. (2008) [43], GP has been very successful at evolving novel and unexpected ways of solving problems.
Our GP system proceeds in the following way. The system will create a number of functions as a population randomly. It firstly generates a random population of functions, and then the evaluations of every function will be performed by the system, where the forecasting accuracy will be compared with the targeted function. We define the performance of the function as the fitness of the function. Afterwards, genetic operations automatically choose one or two function(s) according the fitness evaluation. The genetic operations include two formats, namely, crossover and mutation. The crossover operation recombines two subitems from picked functions to produce new functions. The mutation operation will modify the subitems from picked functions before recombine them together to produce new functions. Then, based on crossover and mutation operations, the GP system will automatically reevaluate the newly-generated functions according to the function fitness. The probability of performing crossover and mutation operations will be pre-determined in the system. This main perquisites are the principles of evolution and the system will be discontinued after some pre-determined conditions are complied such as pre-determined threshold for function fitness. The system automatically chooses the best function as the solution according to the function fitness, which creates the new volatility forecasting model.

Genetic programming system
For our model development, we employ the following targeted functions as our forecasting task based on our GP system regarding the data sample period from January 1, 2001 to December 31, 2010: where L 2 , σ 2 , ε 2 are squared BAS, realized variance and squared residuals respectively, and r 2 denotes the simplified realized variance, known as the squared return term, which is used to approximate the volatility. Our goal is to find the most relevant terms that have effects on predicting the simplified realized variance. For the robust purpose, we also adopt our model for realized variance forecasting, which has been defined in Section 2. The empirical forecasting results for both simplified realized variance and realized variance will be presented in Section 4.
• Fitness measure: the error between the value of the individual function and the corresponding desired output (i.e. r 2 ).
• GP parameters: population = 10000, the maximum length of the program = 1000, probability of crossover operation = 0.8 and probability of mutation operation = 0.1.
• Ending conditions: when the measure of function fitness approximates 0 or the system runs up to 100 generations, the system will automatically discontinue. (For this work, the fitness measure will never reach 0, therefore the system will terminate after 100 generations.) As the number of generation increases, the average error between the value obtained by forecasting functions and the target value may decrease. In this case, GP may produce some high-order terms in order to further reduce the error. However, these high-order terms may cause the overfitting problem. Therefore, we further eliminate the forecasting functions with more than 10 terms obtained by our GP system in order to prevent the overfitting. The parameters settings used above and the choice of 10 terms are based on empirical experiences and we do not claim that these are the optimal choices. Our main purpose is to test the effectiveness of GP approach on volatility forecasting task.

Model development
With the settings stated above, we ran our GP system for 50 times in order to acquire the specific model with no ε 2 term. Since the ε 2 data series is in the program, the model solved by the program will contain the information incorporated in the term ε 2 . Therefore, the model we devise with no ε 2 term can lead to the comprehensive decomposition of the ε 2 term for the model. After simplification, the best function with no ε 2 term we obtained is the following: t Þ and L presents the BAS. We name the model as LGARCH (1,1), which is a liquidity-adjusted GARCH model. Furthermore, the GARCH (1,1) model has the following the form: where σ t is the volatility of target time series and ε t is the residual term from the return prediction equation, which is: ϕ is the conditional mean, and ε t � Nð0; s 2 t Þ. Since the conditional variance may depend on the past squared residuals of the process, namely, ε t , under the GARCH model, the ε t term plays a vital role in data fitting and conditional variance prediction. The general decomposition of ε t has been provided in the literatures, which is ε t = σ t � z t (see Lamoureux and Lastrapes 1990 [44]; Nelson 1990 [45]) and z t is the noise process. Therefore, Thus, based on this decomposition, our liquidity-adjusted GARCH model perfectly matches the decomposition of the original GARCH model. Like the GARCH-X model, the LGARCH(1,1) model uses K t to replace ε 2 t but it further elaborates the z t term embedded in ε t by adopting the BAS, which is denoted as L. It can be observed from Figs 2-4 that the BAS has the similar spike clusters with both simplified realized variance and realized variance for the periods late 2007 and early 2016. Therefore, our model implies that K tÀ 1 ¼ ε 2 tÀ 1 . Similarly, concerning the GARCH-X model, we have the general form: where x t is the exogenous variable, such as interest rate and y t denotes the residual term, and y 2 1. When L 2 t > 1, it indicates that the BAS in the market is relatively large and the liquidity of the market is insufficient. As a result, ð1 À L 2 t Þ < 0 and ðL 2 t À 2s 2 t À L 4 t Þ < 0, which generate the positive product. Since the BAS and volatility have a positive correlation, the large BAS may indicate higher volatility in the future. Therefore, the adjustment increases the volatility term by adding a positive value to the process to reflect the potential increase of volatility in the near future.
2. When L 2 t < 1, it indicates that the BAS in the market is relatively small and the market is liquid.
the product is still positive, in which case, the spread is still not small enough. Therefore, the adjustment still has positive value to unveil the potential increase of the future volatility.
(ii) When L 2 t < s 2 t þ L 4 t , the product is negative, in which case, the BAS is trivial and trading frictions in the market are negligible. Therefore, the model may foresee possible volatility decrease in the near future by adding a negative value to the previous volatility.
In order to testify our model's data fitness and volatility prediction ability, we adopt two tests for the model with both in-sample and out-of-sample tests. The first test is to predict the simplified realized variance and the second test is to predict the realized variance. We are able to show that the detailed decomposition of z t in the GARCH model substantially heightens the forecasting accuracy of the GARCH model.

Empirical results
This section gives both empirical results for regression models and model performance of volatility forecasting. In particular, we compare our data fitting results as well as prediction results with three GARCH models, namely, GARCH, IGARCH and GJR-GARCH (denoted as TGARCH in the tables).
For the model performance evaluation, we use Mean Squared Error (MSE) to measure the model performance for both in-sample fitting and out-of-sample forecasting tests since our daily data can be quite noisy (Pong et al. 2004 [46]; Bollerslev et al. 2016 [47]). More importantly, MSE is quite suitable for time-series forecasting accuracy measurement (Ismail et al. 2011 [48]). The periodic averaged MSE can be defined as: where T represents the number of observations embedded in the forecasting period, Observed t presents the observed variance from the market and Predicted t presents the variance predicted from the models. Lower MSE indicates higher forecasting accuracy as well as more powerful prediction of the model.

Empirical models
Before the model performance evaluation, we firstly run a series of regressions to certify the significant impact of BAS on both simplified realized variance and realized variance. The specification of the regression model is shown in the following equations: where SRV t is the simplified realized variance, RV t is the realized variance, L t is the bid-ask spread. We run a series of multivariable regression analysis and the results are shown in Table 2. It is observable that the impact of BAS on volatility is oftentimes significant after controlling different variables. Therefore, the empirical delivers strong support of BAS's impact on volatility and we thereby incorporate it into the GARCH model for volatility forecasting.

In-sample data fitting
For the in-sample modeling, we compare three models with our model in fitting both simplified realized variance and realized variance. Table 3 shows the in-sample fitting MSE against the simplified realized variance for all four models during the data period 2001-2010. In general, our model outperforms other three models and the improvement rate is around 25%

Out-of-sample forecasting
For the out-of-sample forecasting, we also compare three models with our model in forecasting both simplified realized variance and realized variance with the prediction errors measured by MSE. For out-of-sample forecasting, we use the normalized MSE loss function proposed by Chen and Watanabe (2019) [49]: where m is the out-of sample size. Table 5 shows out-of-sample forecasting MSE against the simplified realized variance for all four models during the data period 2011-2019. In general, our model outperforms other three models and the improvement rate is around 22% compared with other three GARCH models for the full sample test against the simplified realized variance. For the average improvement Table 4. In sample fitting with realized variance. This table presents the oil futures volatility in-sample fitting results of Liquidity-adjusted GARCH (LGARCH) with other three GARCH models against the realized variance using the Mean Squared Error (MSE). The p-values for statistical differences of the forecasting errors are also presented. The LGARCH model outweighs all other three models in both full sample and all subsample tests. Where en = � 10 n , e.g. e − 02 = � 10 −2 .

GARCH
LGARCH Improve  IGARCH  LGARCH  Improve  TGARCH LGARCH Improve rate of subsample, it is around 17% compared with the three GARCH models. Additionally, Table 6 presents the out-of-sample forecasting MSE against the realized variance for all four models during the data period 2011-2019. The results are also dominating for both full sample and subsample tests. Our model overwhelmingly surpasses all three GARCH models with the full sample improvement rate around 73% and the average improvement rate of subsample around 70%. Therefore, it is arguable that our model produces a substantially accurate forecasting results compared with all three GARCH models.
Nevertheless, for the out-of-sample forecasting, the results may not be as stable as the in sample fitting. For the data of 2013 and 2014, the oil price move trend has been changed. In 2013, the oil price had an increasing trend with the annual return of 6.90%. On the other hand, however, the oil price plummeted in 2014, with the annual return of -45.55%, which yielded Table 5. Out of sample forecasting with simplified realized variance. This table presents the oil futures volatility out-of-sample forecasting results of Liquidity-adjusted GARCH (LGARCH) with other three GARCH models against the simplified realized variance using the Mean Squared Error (MSE). The p-values for statistical differences of the forecasting errors are also presented. The LGARCH model outweighs all other three models in both full sample and all subsample tests. Where en = � 10 n , e.g. e − 02 = � 10 −2 .

GARCH
LGARCH Improve  IGARCH  LGARCH  Improve  TGARCH LGARCH Improve the lowest price of $53.45 since 2010. As a result, the return volatility movement may not follow similar pattern, which may result in the significant difference of volatility forecasting.

Robustness check
To ensure our results are robust, we also adopt the Mincer-Zarnowitz regression to verify the relative performance of our volatility forecasting model compared with other GARCH models. Following Mincer and Zarnowitz (1969) [50], we run following two regressions based on our volatility forecasting results: where SRV t and RV t are the simplified realized variance and the realized variance observed at time t respectively,v model1;t andv model2;t are the forecasted variance from model 1 and model 2 at time t respectively. For evaluating the model performance, we first run the single variable regression on the value predicted by LGARCH model as the model 1 in Eqs (12) and (13). Then, we run two-variable regression on both LGARCH as model 1 and one of three models, GARCH, IGARCH and TGRACH as model 2 in Eqs (12) and (13) sequentially. Thus, we can evaluate the model performance from two respects: Firstly, we compare the adj − R 2 from single variable regression with the adj − R 2 from twovariable regression to see whether the three GARCH models add significant explanatory power to our model. Then, for the relative performance, we investigate the significance and magnitude of the coefficients for both model 1 and model 2 in Eqs (12) and (13). Large and significant coefficients demonstrate high performance. The model performance from Mincer-Zarnowitz regression is shown in Tables 7-10. In particular, Tables 7 and 8 show in-sample model performance for both simplified realized variance and realized variance respectively. In Table 7, it can be seen that the adj − R 2 is 16.2% for our model and additional explanatory power that the three GARCH models can enhance is quite marginal. More importantly, the coefficient of our model is much larger than the coefficients for other three GARCH models. Similarly, from Table 8, it is observable that the adj − R 2 is 58.3% for our model and the explanatory power that the three GARCH models can add is quite small and the coefficients for other three GARCH models are all negative. Table 7. Regression results for simplified realized variance (in-sample). This table reports the Mincer-Zarnowitz regression result for the model performance comparison regarding the in-sample simplified realized variance. SRV0 represents the single variable regression for our model only, SRV-G represents the model 1 is our model and model 2 is GARCH model defined in Eq (12); SRV-IG represents the model 1 is our model and model 2 is IGARCH model defined in Eq (12); SRV-TG represents the model 1 is our model and model 2 is TGARCH model defined in Eq (12). � , �� , ��� indicate statistical significance at 10%, 5% and 1% levels, respectively.

SRV0 SRV-G SRV-IG SRV-TG
LGARCH Therefore, our model overwhelmingly outperforms the other three GARCH models regarding the in-sample volatility forecasting.
On the other hand, Tables 9 and 10 exhibit out-of-sample model performance for both simplified realized variance and the realized variance. Table 9 presents the adj − R 2 is 9.47% for forecasting the simplified realized variance regarding the single regression of our model. The additional adj − R 2 the other three GARCH models can add is nearly 2% and the coefficients are much smaller than our model. Similarly, Table 10 represents the adj − R 2 is 62.8% for forecasting the realized variance regarding the single regression of our model. The additional adj − R 2 the other three GARCH models can add is nearly 4%. The coefficient for our model is quite close to the IGARCH model, but is still larger than other two GARCH models. As a result, our model in general outperforms the other three GARCH models for the out-of-sample volatility forecasting.
Finally, we employ the QLIKE loss function for the robustness purpose. The QLIKE loss function could be written as: where s 2 t;j is the real value of variance andŝ 2 t;j is the predicted value of variance. From Table 11,  (12); SRV-IG represents the model 1 is our model and model 2 is IGARCH model defined in Eq (12); SRV-TG represents the model 1 is our model and model 2 is TGARCH model defined in Eq (12). � , �� , ��� indicate statistical significance at 10%, 5% and 1% levels, respectively.

SRV0 SRV-G SRV-IG SRV-TG
LGARCH where m is the out-of sample size,ŝ 2 t;j is the predicted value and s 2 t;j is the actual value. We compare our model with both ANN-GARCH and SVM-GARCH models, reporting the results in Table 12. The ANN-GARCH model is based on the Artificial Neural Networks, which is the non-parametric method. This method can be applied to uncover the nonlinear associations between the parameters of the GARCH model. Nevertheless, ANN-GARCH is less flexible compared with other models such as neural fuzzy inference system. In addition, the backpropagation algorithm of the ANN disables itself to learn from its own forecasting error. (Kristjanpoller and Minutolo, 2016 [51]).
For the ANN-GARCH method, each neural network connects a group of volatility forecasting variables with output variables and a number of hidden layers. Neurons are connected between the layers for connections that are activated by triggering a threshold. The input group of variables and output group of variables can be a combination of all of different number of neurons in each layer, where the relations between inputs and outputs are embedded (see Kristjanpoller and Hernandez, 2017 [52]; Bhattacharya and Ahmed, 2018 [53]).
The SVM-GARCH model is based on the support vector machine (SVM) stems from the statistical learning theory (see Cortes and Vapnik, 1995 [54]), For the SVM method, it develops  (13); RV-IG represents the model 1 is our model and model 2 is IGARCH model defined in Eq (13); RV-TG represents the model 1 is our model and model 2 is TGARCH model defined in Eq (13). � , �� , ��� indicate statistical significance at 10%, 5% and 1% levels, respectively.

RV-IS SRV-IS RV-OS SRV-OS
LGARCH Based on Table 12, it is thereby arguable that our GARCH outperforms other three GARCH models and we are confident that our results are robust. Further in Table 12, we adopt the QLIKE function to compare our model with the GARCH models with Artificial Neural Network (ANN) and Support Vector Machine (SVM). For SVM-GARCH, we have employed different kernel functions, including linear, polynomial and Gaussian kernels. By comparing our model with those models, we also demonstrate that our results are robust since our model has the lowest QLIKE values.

Conclusion and implications
To sum up, fin-tech is a popular financial topic more recently and with the aid from computer science, we are able to improve GARCH family models, which are widely used model in financial time series analysis. This paper proposes the improved GARCH model, which integrates the liquidity factor. Our new GARCH model, named LGARCH, takes the bid-ask spread term into the volatility forecasting, shedding insights on the understanding of the GARCH model as well as the impact of liquidity on the volatility forecasting process. Since how to integrate liquidity into the volatility forecasting model is open to discussion, we use GP method to identify the specific model format for the LGARCH model. This result also indicates that the liquidity factor plays a vital role in oil volatility forecasting since it appears in the final model format. In this vein, fintech facilitates our understanding of oil volatility forecasting as a new factor has been introduced.
Compared with GARCH, IGARCH and TGARCH, our model generally outperforms all three models for predicting both simplified realized variance and realized variance. More importantly, the improvement rate of our model compared with three GARCH models is around 20% for forecasting simplified realized variance and 70% for forecasting realized variance. We also use the Mincer-Zarnowitz regression to demonstrate our results are robust. The strong empirical results and outstanding model performance deliver significant findings on liquidity impact on volatility and yields helpful techniques for volatility forecasting. In Table 12. QLIKE loss function results for ANN and SVM models. This table reports the QLIKE loss function results for the model performance comparison regarding our model with GARCH models with ANN and SVM. SRV represents the simplified realized variance, IS represents the in-sample results (from year 2000-2010) and OS represents the out-of-sample results (from year 2011-2019). SVMl, SVMp and SVMg represent the SVM with linear, polynomial and Gaussian kernel, respectively.

RV-IS SRV-IS RV-OS SRV-OS
LGARCH addition, the more accurate forecasting model can be integrated into the existing fintech systems such as high frequency trading platforms and derivative trading platforms by providing better volatility estimation. Our sample period has 16 years oil daily data, which covered before and after financial crisis periods. Moreover, our paper has two robustness measures for both in-sample and out-of-sample periods, demonstrating the generality and validity of our results. Furthermore, our LGARCH model can also yield a couple of policy implications. Firstly, compared with traditional GARCH models, our oil volatility forecasting model can provide more reliable volatility of crude oil market. The reliable volatility prediction can help oil importing economies to determine countries' oil reserve levels in order to alleviate the negative impact on the economy. More importantly, market volatility might reflect the fragility of financial markets and economy (Celik and Ergin, 2014 [57]). Therefore, reliable forecasts of oil prices' volatility may play crucial roles for macroeconomic policy makers in setting monetary policies to stabilize the economy. Since our results have shown better volatility prediction of the crude oil market, energy economists, energy policy makers, and financial analysts may include market liquidity effect into the volatility forecasting models and utilize market liquidity as an indicator for future volatility changes.
Finally, Genetic Programming (GP) is the subdiscipline of evolutionary algorithms and it optimizes forecasting functions according to an evolutionary process with the evaluation criterion nested. We need to emphasize that, this paper is just a demonstration of using GP into volatility forecasting for the oil market, which can serve as a starting point of fintech development. The GP system in volatility forecasting can be further applied into other financial markets, such as stock market and foreign exchange market. More importantly, since GP is good at providing novel and unexpected insights, it can also be applied to other financial forecasting problems, such as return forecasting, price forecasting and correlation forecasting. As a result, it could be just a starting point of the technology integration in financial modelling for fintech.