Frontiers in Business, Economics and Management ISSN: 2766-824X | Vol. 7, No. 1, 2023 61 Prediction of Carbon Emission Right Price Based on XGBoost Algorithm Fuyun Zhu1, Peiyuan Liu2, Ping Hu3 1 Professor, School of Economics & Management, Shaanxi University of Science & Technology, Xi’an, China 2 Master, School of Economics & Management, Shaanxi University of Science & Technology, Xi’an, China 3Bachelor, School of Economics & Management, Shaanxi University of Science & Technology, Xi’an, China Abstract: Reasonable carbon market price prediction can facilitate the carbon market participants, such as physical producers, to achieve the goal of efficient emission reduction through the market mechanism. In this paper, we use XGBoost, an integrated algorithm in machine learning, to forecast the domestic carbon price from 2013 to 2021. Pearson coefficient is utilized to calculate the correlation of data features, and perform PCA dimensionality reduction on the features with high correlation coefficients. Before using PCA to reduce the dimensions, in order to make the feature data more suitable for later model training, the feature data is standardized. After the standardization, PCA selects super parameters with maximum likelihood estimation and outputs the dimension reduction results of features. Finally, the integrated algorithm XGBoost is used to form a prediction model for the carbon price. RMSE output from cross validation are used to evaluate the accuracy and error of the prediction results. The results show that this paper confirms that the integrated algorithm XGBoost model has a good prediction ability for carbon price, and provides a novel idea for the field of carbon emission price prediction. It is expected to provide some rational basis for the carbon market participants to make investment decisions, so as to avoid the carbon market risk caused by the change of carbon price. Keywords: Prediction of carbon price, Machine learning, Data processing and analysis, XGBoost, Bayes optimization. 1. Introduction The introduction should provide background information (iCarbon emission trading market (carbon market ) is an important policy tool for the international community to address climate change. Excessive emission of greenhouse gas, mainly carbon dioxide, has led to a gradual rise in global temperature, and an over-focus on reducing carbon dioxide emissions will certainly bring high treatment costs. Quantifying carbon emissions means that the government provides carbon emission quotas to enterprises, and enterprises with surplus carbon emissions can transfer the surplus part to enterprises with surplus carbon emissions through trading to achieve optimal allocation of resources and achieve emission reduction target at a lower cost. In June 1992, the United Nations Framework Convention on Climate Change encouraged developed countries to reduce carbon dioxide emissions and set the ultimate goal of "stabilizing the concentration of atmospheric greenhouse gases at a level that would prevent dangerous anthropogenic interference with the climate system. This level should be achieved within a time frame sufficient to make the ecosystem sustainable." In December 1997, the Third Conference of the Parties was held in Kyoto, Japan. The signing of the Kyoto Protocol forced major industrialized countries to reduce their carbon dioxide emissions from 2008 to 2012 by 5.2% from 1990. With the Kyoto Protocol and the provision of three flexible performance mechanisms, namely Joint Implementation Mechanism (JI), Clean Development Mechanism (CDM) and International Emissions Trading (IET), these three flexible mechanisms have given birth to the international carbon emission trading market. The EU established the EU Emission Trading System (EU-ETS) in April 2005, and has made remarkable achievements in environmental and economic benefits, becoming the most influential global carbon emission trading market in the world. In October 2011, the Notice of the General Office of the National Development and Reform Commission on the Pilot Project of Carbon Emission Trading identified that Beijing, Tianjin, Shanghai, Chongqing, Hubei, Guangdong and Shenzhen should carry out the pilot work of carbon emission trading. In 2019, China's carbon dioxide emissions reached 9.826 billion tons, of which 80% will come from coal combustion, and the carbon dioxide emissions of power, industry and transportation industries will account for 90% of the total emissions. At the 75th United Nations General Assembly in 2020, China pledged d to peak its carbon dioxide emissions by 2030 and strive to achieve carbon neutrality by 2060. At the National People's Congress and the Chinese People's Political Consultative Conference in 2021, Premier Li Keqiang pointed out that efforts should be made to achieve the vision of "carbon neutrality". On July 16, 2021, the national carbon market officially started online trading. On March 15, 2022, the Ministry of Ecology and Environment issued the Notice on Key Work Related to the Management of Enterprises' Greenhouse Gas Emission Reports in 2022, which established the list of key emission units in the power generation industry in the second compliance cycle of the national carbon market. According to the verification results, the annual carbon emissions in 2020 or 2021 will reach 26,000 tons of carbon dioxide equivalent, It also has key emission units in the power generation industry that meet the quota management standards for units, and is included in the list of key emission units under quota management in the national carbon emission trading market in 2022. The research on carbon trading market price prediction has focused on two main points: one is to use only the historical data of carbon trading market price to predict 62 carbon price [1-5], and the other is to predict carbon price by combining the influencing factors of carbon trading market price [6-10]. However, since there are few articles that apply machine learning to carbon price prediction in China. Based on this, this paper will use XGBoost model to predict the carbon price,which will provide new ideas for the carbon price prediction field, and provide some practical basis for the carbon market participants to make investment decisions, so as to avoid the carbon market risk caused by the carbon price change. 2. Research Method 2.1. XGBoost algorithm principle Extreme Gradient Boosting (XGboost) is an end-to-end gradient lifting tree system proposed by Dr. Tianqi Chen of Washington University. XGBoost implements a generalized Tree Boosting algorithm. One representative of this algorithm is the Gradient Boosting Decision Tree (GBDT). XGBoost adopts the idea of integration, adds regular terms on the basis of GBDT, and extends the loss function from square loss to Taylor's second-order loss. The following article will briefly describe the principle of XGBoost algorithm: 2.1.1. Boosting additive model representation (1) Where refers to the predicted value of the ith sample equal to the sum of m regression trees, and xi refers to the characteristic value of the ith sample. Forward distribution algorithm: )(ˆˆ )1()( xfyy it t i t i   (2) Where t represents the t-th tree to be optimized, t-1 represents the t-1 tree that has been optimized, and )(xf it represents the t-th tree to be optimized. objective function:    t j j N i t ii t fyyobj L 11 )( )(),( ˆ ）（ (3) Where   N i t ii yyL 1 )( ),( ˆ represents the total loss of N samples, and    t j j f 1 )( is a regular term (adding the complexity of t regression trees). The purpose of adding a regular term is to prevent over fitting. Optimization objective: objwww t T minarg )(** 2 * 1 ,, ），（  (4) Set T leaf nodes, determine their node values, and minimize the loss function. 2.1.2. Processing of regular terms    T j jt wf T 1 2 2 1 ）（ (5) Where, 、 is a super parameter and a penalty item for variables. T is equal to leaf nodes of the regression tree.   T jw 1j 2 represents the sum of the squares of each node value. When the number of leaf nodes T of the regression tree is larger, the depth of the tree is deeper, and the over fitting phenomenon is more likely to occur. Therefore, T should be punished, and the punishment intensity should be controlled by  ; When   T jw 1j 2 is larger, that is, when the node value is larger, the proportion of the node value of this regression tree in all regression trees is too large. At this time, the phenomenon of over fitting also occurs. Therefore,   T jw 1j 2 is punished, and the punishment is determined by  . 2.1.3. Treatment of    t j j f 1 )( in objective function obj t）（ )()()( 1 1 .1 fff t t j j t j j     (6) The regular term in the objective function is decomposed into the optimized first t-1 term (Constant) and the term to be optimized. In the optimization process of the objective function, that is, the minimization process, constant is not required to be considered:    T j jt t j j wff T 1 2 1 2 1 )()(  (7) Assume that the sample loss is the square error loss: )ˆ(ˆ 2 , yyyy iiL ii ）（ (8) Take the minimum loss calculation of node 1 as an example: wwyyyy wywyLLL 2 1153 2 5 2 3 22 531 )(2)( 1513 )(    ）（ 节点 (9) Where y is a real value and a constant. Minimum W can be obtained: yyw a b 531 2 min  (10) That is, when W1=y3+y5, the loss of L-node 1 is the minimum. (The same with W2, W3 and W4) Review of objective function and regular term formula: )( 1 )( ˆ i m j j m i xfy    y m i ˆ )( 63    t j j N i t ii t fyyobj L 11 )( )(),( ˆ (11)    T j jt t j j wff T 1 2 1 2 1 )()(  (12) The regular term is calculated by leaf node j. At this time, in order to ensure the consistency of the calculation method before and after the formula, the following transformation is made: Convert the loss   N i b ii yyL 1 )( ),( ˆ calculated by sample to the loss calculated by leaf node:     T j i t ii N i t ii I LL j yyyy 11 ),(),( ˆˆ ）（）（ (13) Substituted into the objective function: wyyobj j T j i t ii t I LT j 2 1 )()( 2 1 ),( ˆ      (14) where wyy j t i t i   ˆˆ )1()( , in summary: wwyyobj j T j i j t ii t I LT j 2 1 )1( 2 1 ]),([ ˆ       ）（ (15) In order to facilitate the calculation of node value W, the objective function is divided into unified leaf node calculation. The above assumption is that the sample loss is a square error loss, but the loss functions corresponding to different loss scenarios are different, and it is likely that the loss function cannot be changed into a polynomial containing W2 for solution, so Taylor second-order processing is performed on the loss function ),( ˆ )1( wyy j t ii L   : )( 0)(2 1 ))(()()( 2 0000 xxxxxx fxffxf  2.1.4. Taylor processing of ),( ˆ )1( wyy j t ii L   wyywyyyywyy j t iij t ii t iij t ii LLLL 2)1()1()1(1 ˆˆˆˆ , 2 1 ,,),( ）（）（）（ ）（   ，Where, Wj is the unique variable, and other items are constants. These constants need not be considered in the optimization of loss function, and the first step and the second step are gi and hi respectively: whwgwyy jijij t ii L 2)1( 2 1 , ˆ   ）（ (16) Substituted into the objective function: wwhwgobj jjij T j i j t I T j 22 1 )( 2 1 )] 2 1 ([       (17)    III jjj i ij i ij i jiji hwgwwhwg 22 2 1 ) 2 1 ( (18) in summary:      T j i ij i ij t II T jj hwgwobj 1 2 )]( 2 1 [  ）（ (19) In each leaf node, the sum of a step degree gi of all samples is called Gj; In each leaf node, the sum of two steps hi of all samples is called Hj. To sum up, the objective function is obtained: )]( 2 1 [ 2 1 )( HwGwobj jjj T j j t T     (20) 2.1.5. Minimize the node value of the objective function )]( 2 1 min[arg )]( 2 1 min[arg )]( 2 1 min[arg minarg,,, 2* 2 2 222 * 2 1 2 111 * 1 )(** 2 * 1 HwGww HwGww HwGww objwww jjjjj t j          ）（ (21) At this time, the optimal solution can be obtained: H G w j j j a b   2 * (22) Substitute the optimal solution into the objective function to further obtain the optimal objective function:     T j j jt H Gobjobj T t 1 2 )( 2 1 min )( *   (23) 2.2. Variable Design and Descriptive Statistics of Samples The data in this paper comes from the database of China Research Data Service Platform (CNRDS), which covers the daily data of carbon emission trading information in several provinces and cities (Beijing, Fujian, Guangdong, Hubei, Shanghai, Shenzhen, Tianjin, Chongqing) from August 5, 2013 to October 8, 2021. The data indicators include: trading date, city name, trade variety abbreviation, opening price, highest price, lowest price, average transaction price (yuan/ton) The closing price, trading volume (tons), and turnover were 10615 data samples. (Table 1) 64 Table 1. Characteristic descriptive statistics O price H price m price AT price (yuan/ton) C price TV (ton) turnover count 10615 10615 10615 10615 10615 10615 10615 mean 30.7335968 31.46141875 30.29406029 30.73066886 31.07137164 23412.55315 612664.6003 std 17.44338215 17.62580027 17.34321116 17.43692862 18.35309228 126957.9379 4539548.708 min 1 1 1 1 1 1 2 25% 18.5 19 18.255 18.405 18.52 50 949.89 50% 27.7 28.49 27.02 27.73 27.9 1648 46380 75% 39.5 40.2 39 39.5 39.995 13002.5 377307.5 max 130.9 143.99 119 122.97 555.43 8474383 354157396.3 2.3. Research Method Design In order to improve the prediction accuracy of carbon emission right price and the robustness of the model, the modeling in this paper is divided into the following three steps: ① divide the data set, take the first 10508 groups of data as the training set, and the last 107 groups of data as the test set, and predict the closing price of the last 107 groups of data. ② The Pearson correlation coefficient is utilized to analyze the correlation degree between features. In order to ensure that the data fit the model better, the data is standardized, and the features with high correlation coefficient are reduced by the MLE of PCA. ③ Building XGBoost model to make an initial prediction on data. Use Cross Validate and Bayes optimization to optimize XGBoost model and get the final prediction results. Finally, MAPE indicators are used to evaluate the robustness of the prediction model, and RMSE and MAE indicators are used to evaluate the accuracy of the prediction results. 3. Empirical Analysis 3.1. EDA Analysis In order to better fit the model training of features, this paper uses Pearson correlation coefficient index to analyze the correlation degree of feature variables (Table 2) (Figure 1). Table 2. Pearson correlation coefficient between features O price H price m price AT price (yuan/ton) C price TV (ton) turnover O price 1 0.994613748 0.996937241 0.997143961 0.958351867 -0.0402325 0.034012596 H price 0.994613748 1 0.992208845 0.995677459 0.961935042 -0.037980419 0.036354499 m price 0.996937241 0.992208845 1 0.997399303 0.958658076 -0.04180245 0.033424528 AT price (yuan/ton) 0.997143961 0.995677459 0.997399303 1 0.960931526 -0.040839603 0.034046921 C price 0.958351867 0.961935042 0.958658076 0.960931526 1 -0.038659397 0.032932777 TV (ton) -0.0402325 -0.037980419 -0.04180245 -0.040839603 -0.038659397 1 0.897452236 turnover 0.034012596 0.036354499 0.033424528 0.034046921 0.032932777 0.897452236 1 Figure 1. Heat-map of Pearson correlation coefficient among features In this paper, the P-value between features is positively correlated with the correlation coefficient, so the P-value matrix will not be described separately. It can be seen more clearly from the thermodynamic chart that the opening price, the highest price, the lowest price, and the average transaction price (yuan/ton) are highly linear (0.7 ≤| r |<1), and the transaction volume is highly linear (0.7 ≤ | r |<1) with the transaction amount. 65 Secondly, establish the relationship between the characteristics and the target closing price (Figure 2). Figure2. Relationship between characteristics and closing price As shown in Figure 2, the opening price, the highest price, the lowest price, and the average transaction price under the similar slope have a positive correlation with the closing price. There is no obvious linear relationship between the trading volume, the transaction volume, and the closing price. Based on EDA analysis, the correlation between features is relatively high. Although there are no obvious outliers, the possibility cannot be excluded . It is an indispensable application link of machine learning to solve data anomaly and data preprocessing behavior with extreme feature correlation. Therefore, feature engineering is also needed for data before building models. Table 3. Parameter Space Tuning of XGBoost Model Parameters Parameter meaning search space Tuning results learning_rate Update the contraction step in the learning process [0.01,3,0.01] 0.05 max_depth The maximum depth of the tree to avoid over fitting [2,30,1] 8 min_child_weight Minimum sample weight sum of leaf nodes to avoid over fitting [1,10,1] 1 n_estimators Control the number of weak learners [20,200,1] 171 subsample Control the proportion of random sampling for each tree [0.01,0.8,0.01] 0.39 3.2. Feature Engineering and XBGoost Model Building This experiment is performed on a computer with Intel i5- 7300HQ 2.5GHz CPU, 8 GB RAM and Windows10 operating system. Under the Anaconda environment, Python language is employed for programming. The IDE is Jupyter Notebook and GPU acceleration is not used. The optimization process of the model takes the form of seeking to minimize the loss function, and so is Bayes optimization. Generally, there are several ideas to find the minimum value of a function f (x) within the value range of x: let the first derivative of the function be 0, and determine the minimum value of f (x) by solving the value of x at this time, but this method requires that f (x) is differentiable and the differential equation has solutions; The minimum value of f (x) can also be obtained by gradient descent and other optimization methods, but this method requires that f (x) is differentiable and convex; Another method is to substitute x into f (x) to obtain all possible values of the function and then determine the minimum value. Obviously, the workload of this method is extremely large. The main idea of Bayes optimization is to randomly select points of loss function f (x), and continuously reduce the confidence interval and improve the confidence level between these points to obtain a new loss function b (x) that is infinitely close to the real loss function f (x). At this time, the minimum value of b (x) is regarded as the minimum value of f (x). Bayes optimization idea solves the problem of finding the minimum value when the loss function is smooth and uniform and complex and nondifferentiable. Before modeling, divide the data set, take the first 10,508 groups of data as the training set, and predict the closing price of the last 107 groups of data. In order to solve the problem of feature outliers and better fit the algorithm model for feature data, standardize the divided data set, that is, transform the original data into a mean value of 0, and the standard deviation is within the range of 0-1. For features with high correlation, the maximum likelihood estimation of principal component analysis (PCA) is used to reduce the feature dimension. Instantiate XGBoost Regressor, use Cross Validate to optimize the model, and KFold uses 50% cross validation to avoid disrupting the sample sequence, and use TPE (Tree Parzer Estimator) based Bayes optimization to 66 optimize the model parameters, and finally output the optimized prediction results. 3.3. Evaluating Indicator In this paper, the mean absolute error (MAE) and root mean square error (RMSE) are used to evaluate the prediction results, and the mean absolute percentage error (MAPE) is used to evaluate the stability of the prediction model.     n i yy iin RMSE 1 2 )ˆ(1 (24) (25) (26) Among them, the more novel the MAE and RMSE are, the better the prediction results are. The more novel the MAPE is, the stronger the interpretation ability of the characteristics to the target value is. The better the fitting degree of the model to the data is. 4. Experimental Results and Analysis From August 5, 2013 to October 8, 2021, there are 10615 daily time series data of the opening price, the highest price, the lowest price, the average transaction price, the closing price, the trading volume and the transaction amount of each domestic carbon emission trading market. The rise and fall of some of the closing prices are shown in Figure 4. Figure 4. Partial Up and Down Chart of Closing Price The initialized XGBoost model is adjusted by Bayes optimization based on TPE to make the model prediction results more accurate. ①Define the objective function, that is, instantiate the regression and train the model. Add a 50% cross validation that does not disrupt the order of samples to the objective function, and the output Scoring is the negative root mean square error. ② Define the parameter space. The super parameters in XGBoost algorithm are divided into three types: general parameters (for macro function control), Booster parameters (for controlling Booster details), and learning objective parameters (for determining learning strategies). This paper selects the learning that has the greatest impact on the algorithm performance_ rate、max_ depth、 min_ child_ weight、n_estimators and samples are used as the parameter space.③The optimization function is defined.The optimization agent model adopts TPE algorithm. When the loss stops for 100 consecutive times, the training is stopped. The loss is set to root mean square error (RMSE) and printed when the training is stopped. The results of optimization are shown in Table 3. Learning_ rate, max_ depth、min_ child_ weight, n_ The final tuning results of estimators and subsamples are 0.05, 8, 1, 171 and 0.39 respectively. Both the stability and prediction effect of the Bayes optimized model are significantly improved compared with the traditional XGBoost model (Table 4) and the data fitting effect is better (Figure 5). Table 4. Comparison of Algorithm Evaluation Indicators Algorithms RMSE MAE MAPE Original-XGBoost 48.42 5.15 0.02 Bayes-XGBoost 1.27 0.55 0.02 Figure 5. Final Forecast Results and Real Results 5. Conclusion a. The XGBoost model based on Bayes optimization proposed in this paper can effectively predict the carbon price, and is more accurate than the traditional XGBoost prediction method. b. The regular terms incorporated in XBGoost can reduce the complexity of the model, effectively prevent over fitting, reduce the amount of calculation, and significantly improve the prediction efficiency of the model. Therefore, it also reflects the advantages of using XGBoost combined model for prediction. To sum up, the Bayes-XGBoost model proposed in this paper provides a new idea in the field of carbon emission price prediction, which can provide some practical basis for carbon market participants to make investment decisions, so as to avoid carbon market risks caused by carbon price changes. References [1] Chevallier J. Nonparametric modeling of carbon prices [J]. Energy Economics, 2011, 33(6):1267–1282. [2] Chevallier J, Sévi B. On the realized volatility of the ECX CO2 emissions 2008 futures contract: distribution, dynamics and forecasting [J]. Ann Finance, 2011, 7:1–29. [3] Koop G, Tole L. Forecasting the European carbon market [J]. Journal of the Royal Statistical Society: Series A (Statistics in Society), 2013, 176(3):723–741. [4] Atsalakis G S. Using computational intelligence to forecast carbon prices [J]. Applied Soft Computing, 2016, 43:107-116.     n i i ii y yy n MAPE 1 ˆ1     n i i yy in MAE 1 ˆ1 67 [5] Sun W, Zhang C C. Analysis and forecasting of the carbon price using multi—resolution singular value decomposition and extreme learning machine optimized by adaptive whale optimization algorithm [J]. Applied Energy, 2018, 231:1354- 1371 [6] Wei Yu, Zhang Jiahao Research on the Prediction Method of China's Carbon Emission Trading Price Based on DMS and DMA -- Empirical Evidence from Hubei Carbon Market [J]. System Engineering,2022，40(04):3-15 [7] Zhuang Yingdong. Analysis on the influencing factors of carbon emission trading price -- taking five representative open pilot projects in China as an example [J]. Journal of Shenyang Institute of Engineering (Social Science Edition),2018，14(03):349-357. [8] Abdi A, Taghipour S. Forecasting carbon price in the Western Climate Initiative market using Bayesian networks [J]. Carbon Management, 2019, 10(3):255-268. [9] Han M, Ding L L, Zhao X, Kang W L. Forecasting carbon prices in the Shenzhen market, China: The role of mixed- frequency factors [J]. Energy, 2019, 171:69-76. [10] Zhang Yun. Research on Price Drivers of China's Carbon Trading -- A Dual Perspective Based on Market Fundamentals and Policy Information [J]. Social Science Series,2018. 01:111-120.