# Forecast and Analysis of Coal Traffic in Daqin Railway Based on the SARIMA-Markov Model

Authors: Cheng Zhang; Shouchen Liu
Journal: Journal of Advanced Transportation
 (2020)

---

## Abstract

With the continuous advancement of China’s supply-side structural reform, the country’s energy consumption structure has undergone considerable changes, including an overall reduction in fossil energy use and a rapid increase in clean energy application. In the context of China’s coal overcapacity, port and rail capacities are difficult to change in the short term. This study forecasts the monthly coal traffic of Daqin Railway on the basis of the seasonal autoregressive integrated moving-average Markov model and then uses the monthly coal transport data of this railway from September 2009 to November 2019 as samples for model training and verification. Coal traffic from December 2019 to September 2020 is accurately predicted. This study also analyzes the effects of China’s industrial structure adjustment, clean energy utilization, and low-carbon usage on the coal transport volume of Daqin Railway. In addition, the characteristics of seasonal fluctuation and the development trend of Daqin Railway’s coal traffic are explored. This study provides a reference for adjusting the train operation chart of Daqin Railway’s coal transport and developing a special coal train operation plan. It can determine the time of coal transport peak warning, improve the efficiency of coal transport management, and eventually realize a reasonable allocation of resources for Daqin Railway.


---

## Body

## 1. Introduction

Accounting for 1% of the length of China’s railways, Daqin Railway comprises 20% of the country’s national railway and 13% of its national coal transport volume. Daqin Railway holds and even constantly breaks the record for the fastest train running speed, the highest running density, the largest single railway volume, and the best transportation efficiency. It is a strategic artery of China’s “west-to-east coal transport,” which continuously carries coal to the Bohai Sea at a rate of 6.3 tons per second. However, under the background of China’s coal production capacity, port and railway transport capacities are excessive and difficult to change within a short period. Moreover, the competition for China’s coal supply transport channel is becoming increasingly fierce, and consequently, Daqin Railway is expected to experience pressure from new industry competition. Historical data can be used to predict the changing trend of coal traffic in Daqin Railway. The railway transportation department can then adjust the train operation chart and formulate a coal train operation plan in accordance with the predicted coal traffic.At present, local and foreign scholars have focused on railway freight forecasts. Many research achievements have been reported in the prediction of railway freight volume. Commonly used methods include the adaptive flocking algorithm [1], support vector machines [2], the gray model [3], time series models [4], neural networks [5], combination of models [6], multiple models [7], and regression analysis models [8]. Liu and Yu [9] used the seasonal autoregressive integrated moving-average (SARIMA) model to predict and analyze railway freight volume. Zhang [10] applied a time series model and a neural network to predict the annual freight volume of the Guangzhou–Shenzhen railway line. Zhao [11] solved an autoregressive integrated moving-average (ARIMA) model by using EViews software and predicted China’s railway freight volume. Huang et al. [12] and Yuan et al. [13] analyzed the considerable error produced by the gray Verhulst model in predicting railway freight volume and used a Markov chain model to modify the prediction results of the Verhulst model, improving prediction accuracy. Zhang and Zhou [14] used the gray forecast-Markov chain-qualitative analysis method to predict railway freight volume. Tang [4] constructed an improved gray MARKOV prediction model and predicted the future freight demand of China’s railways. In Milenković et al. [15], the time series of the monthly passenger flow in Serbian railways from January 2004 to June 2014 was fitted and predicted using the SARIMA method. In Tang and Deng [16], an ARIMA model was established and R programming language was used to solve this model to make reasonable predictions of the future development trend of civil aviation passenger transport.On the basis of existing studies, local and foreign scholars have conducted considerable research on railway freight transport prediction by using different methods and from various perspectives. However, relatively few studies on railway coal transport prediction have been performed in the context of China’s coal production capacity, excessive railway transport capacity, and fierce competition among coal transport channels. In the current work, the actual coal transport volume situation of Daqin Railway is considered, and the SARIMA-Markov model is adopted to predict the monthly coal transport volume of Daqin Railway. This research analyzes the economic and market reasons for coal transportation in Daqin Railway and explores the seasonal fluctuation characteristics of coal transportation in this railway. The results of this study will provide an important reference for the managers of Daqin Railway to adjust their train operation chart and formulate a special operation plan for coal trains. The peak time of coal transport can be determined, and the efficiency of coal transport management can be improved.The remainder of this paper is organized as follows. In Section2, we provide the seasonal inspection data of monthly coal volume in Daqin Railway. In Section 3, we describe the SARIMA-Markov prediction model in detail. In Section 4, we verify the applicability of the SARIMA-Markov model by selecting the monthly coal volume of Daqin Railway from January 2009 to September 2019 as sample data for testing the model and forecasting future trends. Finally, we summarize the study in Section 5.


## 2. Seasonal Inspection of the Monthly Coal Traffic Volume of Daqin Railway

Daqin Railway provides a steady momentum for China’s sustained economic development, and it has become a landmark among China’s railways since the country’s reform and opening up 40 years ago. In the monthly forecast of the coal traffic volume of Daqin Railway, the historical data of the railway’s monthly coal traffic volume should be analyzed, and then an appropriate forecast model should be developed. In the current study, the monthly coal traffic volume of Daqin Railway from January 2009 to September 2019 is selected as the observation data, and the time chart of the coal traffic volume of Daqin Railway is drawn using EViews 10.0, as shown in Figure1.Figure 1
Time chart of the monthly coal traffic of Daqin Railway.As shown in Figure1, the monthly coal traffic of Daqin Railway exhibited a linear growth trend from 2009 to 2014. As a pilot industry of the supply-side structural reform, the State Council and the National Development and Reform Commission issued corresponding policies to cut coal capacity from 2015 to 2016, resulting in a linear decrease in the coal transport volume of Daqin Railway. After 2017, remarkable achievements were made in adjusting the industrial structure, and the overall benefit of the coal industry recovered steadily. The coal volume of Daqin Railway presented an increasing trend. In the short term, the monthly coal traffic of Daqin Railway exhibits evident seasonal fluctuation characteristics, with a cycle of 12 months. When constructing a prediction model, seasonal, trend growth, and random interference factors should be considered. In addition, the parameters of a matching prediction model should be selected to reduce prediction errors.


## 3. Modeling the SARIMA-Markov Prediction Model

### 3.1. SARIMA Prediction Model

The SARIMA model is derived from the ARIMA model. Its basic form isSARIMAp,d,qP,D,QS. The monthly coal transport volume time series of Daqin Railway is xt,t=0,±1,±2,…,. The autoregressive (AR) model is as follows:(1)xt=c+φ1xt−1+φ2xt−2+⋯+φpxt−p+εt.The moving-average (MA) model is as follows:(2)xt=c+εt−θ1εt−1−θ2εt−2−⋯−θpεt−q,where c is a constant, εt is the residual sequence, and εt∼N0,σ2. When p is the lag order of the time series, the model is referred to as the ARp model. When q is the lag order of the residual sequence, the model is referred to as the MAq model. Equations (1) and (2) can be simplified as follows:(3)φBxt=c+εt,xt=c+θBεt,where B is the backward shift operator, and(4)Bjxi=xt−j,φB=1−φ1B−φ2B2−⋯−φPBP,θB=1−θ1B−θ2B2−⋯−θpBP.ModelsARp and MAq are combined to produce ARMAp,q,(5)xt=c+φ1xt−1+φ2xt−2+⋯+φpxt−p+εt−θ1εt−1−θ2εt−2−⋯−θpεt−p.Equation (5) is simplified to the following:(6)φBxt=c+θBεt.The time series should be stationary during analysis; otherwise, false regression will be produced, leading to unreliable predicted results. Obtaining a strictly stationary time series is difficult, and thus, a time series is required to be weakly stationary. The difference method is typically used to change a weakly stationary time series into a stationary sequence. Some time series can be changed into a stationary sequence after one difference; however, other time series may require multiple differences. A stationary sequence withd-order difference is called a d-order difference sequence. A d-order difference sequence is applied to the ARMAp,q model to form the ARIMAp,d,q model, which is represented by the following:(7)φB1−Bdxt=c+θBεt.Considering the periodicity of a time series, the SARIMA model can be obtained from the seasonal difference and the seasonal parameters of the ARIMA model. The general form of the SARIMA model is as follows:(8)φBΦBS1−Bd1−BSDxt=c+θBΘBSεt,where s is the period of a time series, D is the order of seasonal difference, and BS is the seasonal shift operator. Then, ΦBS=1−Φ1B−Φ2B2−⋯−Φ2BPS, ΘBS=1−Θ1B−Θ2B2−⋯−ΦPBQS.The SARIMA model is denoted asSARIMAp,d,qP,D,Qs, where d is the difference order of each period, D is the order of seasonal difference, p is the autoregressive order, q is the moving-average order, P is the seasonal autoregressive order, and Q is the seasonal moving-average order.The steps of the SARIMA prediction model are as follows:(1)
The sample data sequence is stabilized. Historical data for the observation period are selected as the sample sequencext. The first-order difference of xt is determined to obtain dxt=Dxt,d,S. The trend component in xt is extracted, where d is the difference order and S is the number of periodic difference steps.(2)
TheSARIMAp,d,qP,D,QS model is identified. The parameters p,d,qP,D,Q of the model are determined. The correlation of the difference sequence dxt is analyzed, and the possible values of p and q are preliminarily identified in accordance with the truncated and trailing autocorrelation coefficient and partial correlation coefficient of dxt [17]. Then, the significance of the model parameters is tested. The index values of R2, the Akaike information criterion (AIC), and the Schwarz criterion (SC) are compared. The optimal model parameters are identified.(3)
Model adaptability is tested, and parameters are estimated. The correlation of the fitting residual sequenceεt is analyzed to check if it is a white noise sequence. Whether the SARIMAp,d,qP,D,QS model fully extracts the useful information contained in sequence xt is verified. If the correlation test result of εt is significant, then εt is not a white noise sequence, and the model cannot be adopted even if its evaluation index value is higher. By contrast, when the model passes the adaptability test, the least squares method is used to estimate model parameters.(4)
The model’s predictive power is evaluated.SARIMAp,d,qP,D,QS is evaluated. The evaluation index of the predictive power of the model in the sample period is MAPE=1/nx^t−xt/xt×100. If the prediction effect of MAPE≤5% is excellent, then the prediction effect of 5%≤MAPE≤10% is good, that of 10%≤MAPE≤20% is qualified, and that of 10%≤MAPE≤20% is unqualified [18].


### 3.2. Markov Prediction Model

The Markov model exerts a nonposterior effect. For a random time seriesXt,t∈T, where T is a discrete time set, the nonposterior effect is expressed as follows [19]:(9)PXn+1=in+1|X1=i1,X2=i2,…,Xn=in=PXn+1=in+1|Xn=in.The steps of the Markov prediction model are as follows:(1)
The residual sequenceei is solved as follows:(10)ei=xi−fi,i=0,1,2,…,wherexi denotes the original data, and fi is the predicted value of the SARIMA model [20].(2)
The states are divided.ei is divided into r states, with the same spacing width from large to small. The upper limit of the j state in step i is as follows:(11)Uij=minei+j−1rmaxei−minei,and lower limit is as follows:(12)Lij=minei+jrmaxei−minei.(3)
The probability of state transition ink steps is solved. Suppose the state space is Ω=i1,i2,…,ir, and the probability of state transition in step k is ii=i1,ii=i2,…,ii=ir [21].(4)
The state transition matrix is built. A one-step transfer matrix that can reflect the probability transfer of various states between systems is denoted as follows:(13)p0=p11⋅⋅⋅p1r⋮⋮pr1⋅⋅⋅prr.(5)
The predicted value of the Markov model is expressed as follows:(14)x^t+1=ft+1−∑i=1raitzi,where x^t+1 is the predicted value of the SARIMA-Markov model at time t+1 [22].aiti=1,2,…,r is the probability for a one-step state transition of the row vector of the state transition matrix [22–24]. zii=1,2,…,r is the possible predicted value of each state interval. zi=Uij+Lij/2 is selected in this study [25].


## 3.1. SARIMA Prediction Model

The SARIMA model is derived from the ARIMA model. Its basic form isSARIMAp,d,qP,D,QS. The monthly coal transport volume time series of Daqin Railway is xt,t=0,±1,±2,…,. The autoregressive (AR) model is as follows:(1)xt=c+φ1xt−1+φ2xt−2+⋯+φpxt−p+εt.The moving-average (MA) model is as follows:(2)xt=c+εt−θ1εt−1−θ2εt−2−⋯−θpεt−q,where c is a constant, εt is the residual sequence, and εt∼N0,σ2. When p is the lag order of the time series, the model is referred to as the ARp model. When q is the lag order of the residual sequence, the model is referred to as the MAq model. Equations (1) and (2) can be simplified as follows:(3)φBxt=c+εt,xt=c+θBεt,where B is the backward shift operator, and(4)Bjxi=xt−j,φB=1−φ1B−φ2B2−⋯−φPBP,θB=1−θ1B−θ2B2−⋯−θpBP.ModelsARp and MAq are combined to produce ARMAp,q,(5)xt=c+φ1xt−1+φ2xt−2+⋯+φpxt−p+εt−θ1εt−1−θ2εt−2−⋯−θpεt−p.Equation (5) is simplified to the following:(6)φBxt=c+θBεt.The time series should be stationary during analysis; otherwise, false regression will be produced, leading to unreliable predicted results. Obtaining a strictly stationary time series is difficult, and thus, a time series is required to be weakly stationary. The difference method is typically used to change a weakly stationary time series into a stationary sequence. Some time series can be changed into a stationary sequence after one difference; however, other time series may require multiple differences. A stationary sequence withd-order difference is called a d-order difference sequence. A d-order difference sequence is applied to the ARMAp,q model to form the ARIMAp,d,q model, which is represented by the following:(7)φB1−Bdxt=c+θBεt.Considering the periodicity of a time series, the SARIMA model can be obtained from the seasonal difference and the seasonal parameters of the ARIMA model. The general form of the SARIMA model is as follows:(8)φBΦBS1−Bd1−BSDxt=c+θBΘBSεt,where s is the period of a time series, D is the order of seasonal difference, and BS is the seasonal shift operator. Then, ΦBS=1−Φ1B−Φ2B2−⋯−Φ2BPS, ΘBS=1−Θ1B−Θ2B2−⋯−ΦPBQS.The SARIMA model is denoted asSARIMAp,d,qP,D,Qs, where d is the difference order of each period, D is the order of seasonal difference, p is the autoregressive order, q is the moving-average order, P is the seasonal autoregressive order, and Q is the seasonal moving-average order.The steps of the SARIMA prediction model are as follows:(1)
The sample data sequence is stabilized. Historical data for the observation period are selected as the sample sequencext. The first-order difference of xt is determined to obtain dxt=Dxt,d,S. The trend component in xt is extracted, where d is the difference order and S is the number of periodic difference steps.(2)
TheSARIMAp,d,qP,D,QS model is identified. The parameters p,d,qP,D,Q of the model are determined. The correlation of the difference sequence dxt is analyzed, and the possible values of p and q are preliminarily identified in accordance with the truncated and trailing autocorrelation coefficient and partial correlation coefficient of dxt [17]. Then, the significance of the model parameters is tested. The index values of R2, the Akaike information criterion (AIC), and the Schwarz criterion (SC) are compared. The optimal model parameters are identified.(3)
Model adaptability is tested, and parameters are estimated. The correlation of the fitting residual sequenceεt is analyzed to check if it is a white noise sequence. Whether the SARIMAp,d,qP,D,QS model fully extracts the useful information contained in sequence xt is verified. If the correlation test result of εt is significant, then εt is not a white noise sequence, and the model cannot be adopted even if its evaluation index value is higher. By contrast, when the model passes the adaptability test, the least squares method is used to estimate model parameters.(4)
The model’s predictive power is evaluated.SARIMAp,d,qP,D,QS is evaluated. The evaluation index of the predictive power of the model in the sample period is MAPE=1/nx^t−xt/xt×100. If the prediction effect of MAPE≤5% is excellent, then the prediction effect of 5%≤MAPE≤10% is good, that of 10%≤MAPE≤20% is qualified, and that of 10%≤MAPE≤20% is unqualified [18].


## 3.2. Markov Prediction Model

The Markov model exerts a nonposterior effect. For a random time seriesXt,t∈T, where T is a discrete time set, the nonposterior effect is expressed as follows [19]:(9)PXn+1=in+1|X1=i1,X2=i2,…,Xn=in=PXn+1=in+1|Xn=in.The steps of the Markov prediction model are as follows:(1)
The residual sequenceei is solved as follows:(10)ei=xi−fi,i=0,1,2,…,wherexi denotes the original data, and fi is the predicted value of the SARIMA model [20].(2)
The states are divided.ei is divided into r states, with the same spacing width from large to small. The upper limit of the j state in step i is as follows:(11)Uij=minei+j−1rmaxei−minei,and lower limit is as follows:(12)Lij=minei+jrmaxei−minei.(3)
The probability of state transition ink steps is solved. Suppose the state space is Ω=i1,i2,…,ir, and the probability of state transition in step k is ii=i1,ii=i2,…,ii=ir [21].(4)
The state transition matrix is built. A one-step transfer matrix that can reflect the probability transfer of various states between systems is denoted as follows:(13)p0=p11⋅⋅⋅p1r⋮⋮pr1⋅⋅⋅prr.(5)
The predicted value of the Markov model is expressed as follows:(14)x^t+1=ft+1−∑i=1raitzi,where x^t+1 is the predicted value of the SARIMA-Markov model at time t+1 [22].aiti=1,2,…,r is the probability for a one-step state transition of the row vector of the state transition matrix [22–24]. zii=1,2,…,r is the possible predicted value of each state interval. zi=Uij+Lij/2 is selected in this study [25].


## 4. Case Analysis

To verify the applicability of the SARIMA-Markov model, this study selected the monthly coal traffic volume of Daqin Railway from January 2009 to September 2019 as sample data to establish the prediction model. Evaluation indicators were adopted to evaluate the model. Among these, sample data from October to November 2019 were used for the model test, and trend extrapolation was performed for the period of December 2019 to September 2020.Figure 2
First-order difference sequence ofxt.


### 4.1. SARIMA Model

(1)
The data sequencext of the sample period is stabilized. EViews 10.0 software is used to find the first difference of xt. The mean value of dxt=Dxt,1,0 is approximately zero, and the trend characteristic disappears, as shown in Figure 2. The linear trend of sequence xt is fully extracted by the first-order difference [26].The autocorrelation diagram of the first-order difference sequencedxt is provided in Figure 3. The autocorrelation coefficients near the 12th and 24th orders are significantly not 0, and thus, a seasonal feature with a period of 12 exists in the first-order difference sequence dxt. This finding is consistent with the intuitive analysis results of the sequence diagram [27].After the difference operation with a period of 12 ondxt, a new sequence d′xt=Dxt,1,12 is obtained [28]. As shown in Figure 4, the 12th-order autocorrelation coefficient of the second-order difference sequence d′xt is close to 0, indicating that the periodic factor dxt is fully extracted [29–31].An augmented Dickey-Fuller (ADF) test was performed on the second-order difference sequenced′xt, and the results are provided in Table 1.The unit root statistic isADF=−7.635505, which is less than the critical value with a significance level of 1%. Thus, the null hypothesis that states that the second-order difference sequence d′xt has a unit root was rejected, indicating that the sequence d′xt is a stationary sequence. Accordingly, the parameters of the SARIMAp,d,qP,D,QS model are d=1,D=1,S=12.(2)
TheSARIMAp,d,qP,D,QS model is determined [32]. As shown in Figure 4, xt exhibits no significant correlation between different points in the same period. Thus, a simple seasonal model is established for the data sequence SARIMAp,1,q1,1,112 of the sample period. The 1st-order autocorrelation coefficient of the differential stationary d′xt is significant, and its 2nd-order partial autocorrelation coefficient is also significant. The autocorrelation and partial autocorrelation diagrams are tailed; thus, five models can be built as follows:(15)SARIMA2,1,11,1,112,SARIMA4,1,01,1,112,SARIMA1,1,31,1,112,SARIMA1,1,11,1,112,SARIMA1,1,21,1,112.The five models are tested for parameter significance, and the results are presented in Table2.R2 is the goodness of fit of the entire model; the higher the R2 value, the better the fit degree, where R2∈0,1. AIC and SC are information criteria; the lower the value, the better the fitting degree of the model [33]. The test results of five compared models are provided in Table 2. The sequence model of data xt in the sample period is SARIMA4,1,21,1,112.(3)
Model adaptability is tested, and parameters are estimated. In EViews 10.0, theSARIMA4,1,21,1,112 model is used to obtain the fitting residual sequence εt, as shown in Figure 5.A white noise test is performed on the model residualεt of SARIMA4,1,21,1,112; that is, a random test of fitting residual sequence εt. The test results are presented in Figure 6.As shown in Figure6, the statistical value of Q is 8.963, and the associated probability is 0.34>0.05 from line K=12. The null hypothesis states that no correlation exists and εt is accepted, indicating that εt is a white noise sequence. Therefore, the SARIMA4,1,21,1,112 model passes the adaptability test, and its expression is as follows:(16)1+0.6569B+0.9473B2+0.2925B3−0.0381B41−0.148B12Δ12lgxt=1−0.4657B−0.8686B21+0.87B12εt.(4)
Model predictive power is assessed. In EViews 10.0, theSARIMA4,1,21,1,112 model is used to predict the sample data of the coal traffic volume of Daqin Railway from January 2010 to September 2019. The predicted results are presented in Figure 7. An analysis model prediction level evaluation index, namely, mean absolute percentage error (MAPE) = 5.1439%, is used to determine whether the model’s predictive capability is good.TheSARIMA4,1,21,1,112 model is used to predict the monthly coal transport volume of Daqin Railway from January 2009 to September 2019, and the results are provided in Table 3. The analysis of Table 3 indicates that the error of the single prediction model is still relatively large. To avoid this situation, the Markov model is adopted to correct the residual. This model can reduce the relative residual value and find the internal regularity.Table 1
ADF test results of sequenced′xt.
t-statisticProb.∗Augmented dickey-fuller test statistic−7.6355050.0000Test critical values1% level−3.4875505% level−2.88650910% level−2.580163Table 2
Significance test of model parameters.
ModelR2AICSCSARIMA2,1,11,1,1120.397714.200714.2676SARIMA1,1,21,1,1120.407314.185114.2519SARIMA1,1,31,1,1120.406414.186614.2534SARIMA4,1,01,1,1120.008814.679514.7241SARIMA4,1,21,1,1120.676213.728113.8395Table 3
Statistical table of the SARIMA model’s predicted value and residual sequence (in million tons).
TimeACTUALSARIMA FError2010/635.538036.3685−0.83052010/736.323036.6030−0.28002010/836.211036.10510.10592010/932.576036.0014−3.42542010/1034.084032.72571.35832010/1132.880035.9553−3.07532010/1233.499034.0188−0.51982011/137.109032.58154.52752011/234.751036.8717−2.12072011/339.116036.42102.69502011/428.843035.7270−6.88402011/539.005034.04134.96372011/637.246038.0273−0.78132011/738.122038.9142−0.79222011/836.939036.76240.17662011/934.745035.8810−1.13602011/1037.969035.26462.70442011/1136.816038.5455−1.72952011/1238.331037.85820.47282012/138.174039.5214−1.34742012/234.876037.6725−2.79652012/336.843037.0296−0.18662012/431.401031.26520.13582012/539.101037.14831.95272012/634.695038.1825−3.48752012/732.731035.4198−2.68882012/831.290032.1657−0.87572012/934.953029.11895.83412012/1034.285035.6686−1.38362012/1138.157034.39113.76592012/1239.384038.94220.44182013/138.874041.9307−3.05672013/234.882036.7719−1.88992013/338.169037.81200.35702013/433.117032.65470.46232013/537.795039.3172−1.52222013/637.057035.79931.25772013/736.649337.0051−0.35582013/837.755036.32931.42572013/938.090036.75671.33332013/1035.383039.4495−4.06652013/1138.282036.56111.72092013/1239.193039.4069−0.21392014/139.741040.0323−0.29132014/234.033036.6904−2.65742014/339.399037.00872.39032014/432.122033.8711−1.74912014/540.607038.31922.28782014/638.936038.33570.60032014/738.925039.6389−0.71392014/838.154038.4263−0.27232014/937.665038.3145−0.64952014/1034.562037.6351−3.07312014/1137.535036.05871.47632014/1238.508038.17760.33042015/138.186039.0928−0.9068Figure 3
Autocorrelation of sequencedxt.Figure 4
Autocorrelation of sequenced′xt.Figure 5
Fitting of theSARIMA4,1,21,1,112 model.Figure 6
Correlation test of the residual sequenceεt in the SARIMA4,1,21,1,112 model.Figure 7
Prediction result of theSARIMA4,1,21,1,112 model in the sample period.


### 4.2. SARIMA-Markov Model

In accordance with the steps of the Markov model, the range of the residual sequenceei∈−6.884,5.834 is first determined by dividing the residual sequence into E1,E2,E3,E4 state intervals. The boundary value of state interval E1 is as follows:(17)Uij=−6.884+1−145.834−−6.884=−6.884,Lij=−6.884+145.834−−6.884=−3.705.Similarly,(18)E2∈−3.705,−0.525,E3∈−0.525,2.655,E4∈2.655,5.834.In accordance with the state distribution of the residual series, the frequency statistics of the one-step transfer fromEi to Ej are provided in Table 4.Table 4
Statistics of the one-step state transfer frequency.
yi⟶j1∑yi⟶j1j = 1j = 2j = 3j = 41003252111188383322273554157013The state transition probability matrix is as follows:(19)P0=000.60000.40000.02630.28950.47370.21050.05450.40000.49090.05450.07690.38460.53850.We obtain the following:(20)z1=−5.2945,z2=−2.115,z3=1.065,z4=4.245,from the formula zi=Uij+Lij/2. The prediction value of the SARIMA-Markov model is calculated in accordance with equations (13) and (14), and the prediction capability of the model is determined to be excellent by the analysis model’s prediction level evaluation indicator, i.e., MAPE=3.8009%. The fitting curves of the actual (ACTUAL), the predicted value of the SARIMA model (SARIMA F), and the predicted value of the SARIMA-Markov model (SARIMA-Markov F) are presented in Figure 8. As shown in the figure, the fitting effect of the SARIMA-Markov model is better than that of the SARIMA model (in million tons).Figure 8
SARIMA-Markov model fitting diagram.


### 4.3. Model Validation

To verify the prediction accuracy of the established model, the SARIMA-Markov model was used to predict the monthly coal traffic volume of Daqin Railway in October and November 2019. The predicted values for October and November 2019 were 32.5864 and 38.7377 million tons, respectively, when theSARIMA4,1,21,1,112 model was used. From Table 3, the state probability transfer vector in September 2019 is E3, and the initial state transfer vector is x0=0.05450.40.49090.0545.In accordance with the initial state probability transfer vector 1 and the state probability transfer matrix 2, the state probability transfer vectors of October and November 2019 can be calculated as follows:(21)x1=x0∗p0=0.04150.33190.49250.1328,x2=x1∗p0=0.04580.34320.49540.1133.Then, from equations (13) and (14), the predicted coal transport volume in October 2019 is 34.1619 million tons and the predicted coal transport volume in November 2019 is 32.8694 million tons. The results are provided in Table 5.Table 5
Statistics of the prediction results (in million tons).
TimeActualSARIMA FSARIMA-markov F20191036.443932.831334.161920191135.412738.601832.8694In the current study, mean absolute error (MAE) and MAPE were used to evaluate the fitting effect of the model, as indicated in Table6. Both indexes were smaller for the optimized model than those before optimization. The fitting accuracy of the SARIMA-Markov model reached 96.1991%, which was higher than the prediction accuracy of the SARIMA model (94.8561%). The prediction accuracy of the SARIMA-Markov model reached 93.2874%, which was higher than the prediction accuracy of the SARIMA model (90.5424%). Thus, the prediction accuracy of the SARIMA-Markov model is high and meets the requirements.Table 6
Model test list (in million tons).
Forecast modelModel fittingModel forecastingMAEMAPEMAEMAPE (%)SARIMA5.76115.14393.40039.4576SARIMA-markov4.25703.80092.40946.7126


### 4.4. Model Prediction

The model was built on the basis of the monthly data of the coal traffic volume of Daqin Railway from January 2009 to September 2019. The predicted sequencexf and the actual sequence x were drawn in the same diagram for comparison, as shown in Figure 9. In accordance with the state probability transition matrix P and the initial state transition vector x0, the state transition vector from December 2019 to September 2020 was determined. Simultaneously, the predicted values of SARIMA4,1,21,1,112 from December 2019 to September 2020 were obtained.Figure 9
Prediction result of theSARIMA4,1,21,1,112 model outside the sample period.The coal traffic volume of Daqin Railway from December 2019 to September 2020 was predicted using the SARIMA-Markov model. The results are provided in Table7.Table 7
Statistics of the forecast results (in million tons).
TimeCompany FSARIMA FSARIMA MARKOV F2019/12 (E)39.442940.526438.79932020/01 (E)41.583939.489440.58892020/02 (E)34.797136.617335.55342020/03 (E)41.003940.674240.73002020/04 (E)34.282734.278034.33462020/05 (E)43.410941.416141.47872020/06 (E)38.887140.640239.53332020/07 (E)40.322439.771740.71162020/08 (E)34.811339.505138.70802020/09 (E)34.597134.667834.4061


## 4.1. SARIMA Model

(1)
The data sequencext of the sample period is stabilized. EViews 10.0 software is used to find the first difference of xt. The mean value of dxt=Dxt,1,0 is approximately zero, and the trend characteristic disappears, as shown in Figure 2. The linear trend of sequence xt is fully extracted by the first-order difference [26].The autocorrelation diagram of the first-order difference sequencedxt is provided in Figure 3. The autocorrelation coefficients near the 12th and 24th orders are significantly not 0, and thus, a seasonal feature with a period of 12 exists in the first-order difference sequence dxt. This finding is consistent with the intuitive analysis results of the sequence diagram [27].After the difference operation with a period of 12 ondxt, a new sequence d′xt=Dxt,1,12 is obtained [28]. As shown in Figure 4, the 12th-order autocorrelation coefficient of the second-order difference sequence d′xt is close to 0, indicating that the periodic factor dxt is fully extracted [29–31].An augmented Dickey-Fuller (ADF) test was performed on the second-order difference sequenced′xt, and the results are provided in Table 1.The unit root statistic isADF=−7.635505, which is less than the critical value with a significance level of 1%. Thus, the null hypothesis that states that the second-order difference sequence d′xt has a unit root was rejected, indicating that the sequence d′xt is a stationary sequence. Accordingly, the parameters of the SARIMAp,d,qP,D,QS model are d=1,D=1,S=12.(2)
TheSARIMAp,d,qP,D,QS model is determined [32]. As shown in Figure 4, xt exhibits no significant correlation between different points in the same period. Thus, a simple seasonal model is established for the data sequence SARIMAp,1,q1,1,112 of the sample period. The 1st-order autocorrelation coefficient of the differential stationary d′xt is significant, and its 2nd-order partial autocorrelation coefficient is also significant. The autocorrelation and partial autocorrelation diagrams are tailed; thus, five models can be built as follows:(15)SARIMA2,1,11,1,112,SARIMA4,1,01,1,112,SARIMA1,1,31,1,112,SARIMA1,1,11,1,112,SARIMA1,1,21,1,112.The five models are tested for parameter significance, and the results are presented in Table2.R2 is the goodness of fit of the entire model; the higher the R2 value, the better the fit degree, where R2∈0,1. AIC and SC are information criteria; the lower the value, the better the fitting degree of the model [33]. The test results of five compared models are provided in Table 2. The sequence model of data xt in the sample period is SARIMA4,1,21,1,112.(3)
Model adaptability is tested, and parameters are estimated. In EViews 10.0, theSARIMA4,1,21,1,112 model is used to obtain the fitting residual sequence εt, as shown in Figure 5.A white noise test is performed on the model residualεt of SARIMA4,1,21,1,112; that is, a random test of fitting residual sequence εt. The test results are presented in Figure 6.As shown in Figure6, the statistical value of Q is 8.963, and the associated probability is 0.34>0.05 from line K=12. The null hypothesis states that no correlation exists and εt is accepted, indicating that εt is a white noise sequence. Therefore, the SARIMA4,1,21,1,112 model passes the adaptability test, and its expression is as follows:(16)1+0.6569B+0.9473B2+0.2925B3−0.0381B41−0.148B12Δ12lgxt=1−0.4657B−0.8686B21+0.87B12εt.(4)
Model predictive power is assessed. In EViews 10.0, theSARIMA4,1,21,1,112 model is used to predict the sample data of the coal traffic volume of Daqin Railway from January 2010 to September 2019. The predicted results are presented in Figure 7. An analysis model prediction level evaluation index, namely, mean absolute percentage error (MAPE) = 5.1439%, is used to determine whether the model’s predictive capability is good.TheSARIMA4,1,21,1,112 model is used to predict the monthly coal transport volume of Daqin Railway from January 2009 to September 2019, and the results are provided in Table 3. The analysis of Table 3 indicates that the error of the single prediction model is still relatively large. To avoid this situation, the Markov model is adopted to correct the residual. This model can reduce the relative residual value and find the internal regularity.Table 1
ADF test results of sequenced′xt.
t-statisticProb.∗Augmented dickey-fuller test statistic−7.6355050.0000Test critical values1% level−3.4875505% level−2.88650910% level−2.580163Table 2
Significance test of model parameters.
ModelR2AICSCSARIMA2,1,11,1,1120.397714.200714.2676SARIMA1,1,21,1,1120.407314.185114.2519SARIMA1,1,31,1,1120.406414.186614.2534SARIMA4,1,01,1,1120.008814.679514.7241SARIMA4,1,21,1,1120.676213.728113.8395Table 3
Statistical table of the SARIMA model’s predicted value and residual sequence (in million tons).
TimeACTUALSARIMA FError2010/635.538036.3685−0.83052010/736.323036.6030−0.28002010/836.211036.10510.10592010/932.576036.0014−3.42542010/1034.084032.72571.35832010/1132.880035.9553−3.07532010/1233.499034.0188−0.51982011/137.109032.58154.52752011/234.751036.8717−2.12072011/339.116036.42102.69502011/428.843035.7270−6.88402011/539.005034.04134.96372011/637.246038.0273−0.78132011/738.122038.9142−0.79222011/836.939036.76240.17662011/934.745035.8810−1.13602011/1037.969035.26462.70442011/1136.816038.5455−1.72952011/1238.331037.85820.47282012/138.174039.5214−1.34742012/234.876037.6725−2.79652012/336.843037.0296−0.18662012/431.401031.26520.13582012/539.101037.14831.95272012/634.695038.1825−3.48752012/732.731035.4198−2.68882012/831.290032.1657−0.87572012/934.953029.11895.83412012/1034.285035.6686−1.38362012/1138.157034.39113.76592012/1239.384038.94220.44182013/138.874041.9307−3.05672013/234.882036.7719−1.88992013/338.169037.81200.35702013/433.117032.65470.46232013/537.795039.3172−1.52222013/637.057035.79931.25772013/736.649337.0051−0.35582013/837.755036.32931.42572013/938.090036.75671.33332013/1035.383039.4495−4.06652013/1138.282036.56111.72092013/1239.193039.4069−0.21392014/139.741040.0323−0.29132014/234.033036.6904−2.65742014/339.399037.00872.39032014/432.122033.8711−1.74912014/540.607038.31922.28782014/638.936038.33570.60032014/738.925039.6389−0.71392014/838.154038.4263−0.27232014/937.665038.3145−0.64952014/1034.562037.6351−3.07312014/1137.535036.05871.47632014/1238.508038.17760.33042015/138.186039.0928−0.9068Figure 3
Autocorrelation of sequencedxt.Figure 4
Autocorrelation of sequenced′xt.Figure 5
Fitting of theSARIMA4,1,21,1,112 model.Figure 6
Correlation test of the residual sequenceεt in the SARIMA4,1,21,1,112 model.Figure 7
Prediction result of theSARIMA4,1,21,1,112 model in the sample period.


## 4.2. SARIMA-Markov Model

In accordance with the steps of the Markov model, the range of the residual sequenceei∈−6.884,5.834 is first determined by dividing the residual sequence into E1,E2,E3,E4 state intervals. The boundary value of state interval E1 is as follows:(17)Uij=−6.884+1−145.834−−6.884=−6.884,Lij=−6.884+145.834−−6.884=−3.705.Similarly,(18)E2∈−3.705,−0.525,E3∈−0.525,2.655,E4∈2.655,5.834.In accordance with the state distribution of the residual series, the frequency statistics of the one-step transfer fromEi to Ej are provided in Table 4.Table 4
Statistics of the one-step state transfer frequency.
yi⟶j1∑yi⟶j1j = 1j = 2j = 3j = 41003252111188383322273554157013The state transition probability matrix is as follows:(19)P0=000.60000.40000.02630.28950.47370.21050.05450.40000.49090.05450.07690.38460.53850.We obtain the following:(20)z1=−5.2945,z2=−2.115,z3=1.065,z4=4.245,from the formula zi=Uij+Lij/2. The prediction value of the SARIMA-Markov model is calculated in accordance with equations (13) and (14), and the prediction capability of the model is determined to be excellent by the analysis model’s prediction level evaluation indicator, i.e., MAPE=3.8009%. The fitting curves of the actual (ACTUAL), the predicted value of the SARIMA model (SARIMA F), and the predicted value of the SARIMA-Markov model (SARIMA-Markov F) are presented in Figure 8. As shown in the figure, the fitting effect of the SARIMA-Markov model is better than that of the SARIMA model (in million tons).Figure 8
SARIMA-Markov model fitting diagram.


## 4.3. Model Validation

To verify the prediction accuracy of the established model, the SARIMA-Markov model was used to predict the monthly coal traffic volume of Daqin Railway in October and November 2019. The predicted values for October and November 2019 were 32.5864 and 38.7377 million tons, respectively, when theSARIMA4,1,21,1,112 model was used. From Table 3, the state probability transfer vector in September 2019 is E3, and the initial state transfer vector is x0=0.05450.40.49090.0545.In accordance with the initial state probability transfer vector 1 and the state probability transfer matrix 2, the state probability transfer vectors of October and November 2019 can be calculated as follows:(21)x1=x0∗p0=0.04150.33190.49250.1328,x2=x1∗p0=0.04580.34320.49540.1133.Then, from equations (13) and (14), the predicted coal transport volume in October 2019 is 34.1619 million tons and the predicted coal transport volume in November 2019 is 32.8694 million tons. The results are provided in Table 5.Table 5
Statistics of the prediction results (in million tons).
TimeActualSARIMA FSARIMA-markov F20191036.443932.831334.161920191135.412738.601832.8694In the current study, mean absolute error (MAE) and MAPE were used to evaluate the fitting effect of the model, as indicated in Table6. Both indexes were smaller for the optimized model than those before optimization. The fitting accuracy of the SARIMA-Markov model reached 96.1991%, which was higher than the prediction accuracy of the SARIMA model (94.8561%). The prediction accuracy of the SARIMA-Markov model reached 93.2874%, which was higher than the prediction accuracy of the SARIMA model (90.5424%). Thus, the prediction accuracy of the SARIMA-Markov model is high and meets the requirements.Table 6
Model test list (in million tons).
Forecast modelModel fittingModel forecastingMAEMAPEMAEMAPE (%)SARIMA5.76115.14393.40039.4576SARIMA-markov4.25703.80092.40946.7126


## 4.4. Model Prediction

The model was built on the basis of the monthly data of the coal traffic volume of Daqin Railway from January 2009 to September 2019. The predicted sequencexf and the actual sequence x were drawn in the same diagram for comparison, as shown in Figure 9. In accordance with the state probability transition matrix P and the initial state transition vector x0, the state transition vector from December 2019 to September 2020 was determined. Simultaneously, the predicted values of SARIMA4,1,21,1,112 from December 2019 to September 2020 were obtained.Figure 9
Prediction result of theSARIMA4,1,21,1,112 model outside the sample period.The coal traffic volume of Daqin Railway from December 2019 to September 2020 was predicted using the SARIMA-Markov model. The results are provided in Table7.Table 7
Statistics of the forecast results (in million tons).
TimeCompany FSARIMA FSARIMA MARKOV F2019/12 (E)39.442940.526438.79932020/01 (E)41.583939.489440.58892020/02 (E)34.797136.617335.55342020/03 (E)41.003940.674240.73002020/04 (E)34.282734.278034.33462020/05 (E)43.410941.416141.47872020/06 (E)38.887140.640239.53332020/07 (E)40.322439.771740.71162020/08 (E)34.811339.505138.70802020/09 (E)34.597134.667834.4061


## 5. Conclusion

The SARIMA-Markov model was applied to the monthly coal traffic forecast of Daqin Railway. The SARIMA model comprehensively considered the influence of the seasonal correlation of coal traffic volume on Daqin Railway. Meanwhile, the Markov model used a residual sequence, state partition, and state transition matrix to modify the influence of the sample data on the predicted value. Compared with the simple use of the SARIMA model alone, the combination of the two models can ensure higher prediction accuracy, verifying the scientificity and feasibility of the model and providing a new method for the coal volume prediction of Daqin Railway. Moreover, the coal traffic volume of Daqin Railway from December 2019 to September 2020 was predicted via trend extrapolation, and the forecast results were analyzed as follows.(1)
China’s economy is less dependent on coal. As China’s economy moves toward high-quality development, the growth rate of the country’s gross domestic product will remain at approximately 6% in the next few years. Meanwhile, the Daqin Railway, as a strategic artery of China’s “coal transport from the west to the east,” will experience a growth rate of approximately 5%. This conclusion indicates that China’s extensive economic development model based on coal energy consumption is beginning to weaken, and the economic form is developing toward the low-carbon, efficient, and green direction. China’s coal dependence is decreasing due to the country’s effort to adjust its industrial structure, encourage innovations, and increase the intensity of scientific and technological research and development.(2)
Initial results have been achieved in optimizing the energy structure. The growth rate of the coal transport volume of Daqin Railway will gradually slow down under the influences of energy structure adjustment, industrial structure adjustment, nonfossil energy development, and other factors. This scenario shows a steady decline in the proportion of coal in energy consumption. The rapid development of low-cost nuclear, photovoltaic, and wind power will further reduce the demand for coal, which will account for a smaller share of the country’s total energy consumption, while clean, renewable energy will obtain a larger share.(3)
The process of providing heat from clean energy is proceeding in an orderly manner. China is a developing industrial country that uses coal boilers to provide heating during winter, particularly in the north of Qinling and Huai River. Under the background of the energy revolution, China has effectively promoted clean heating in the north and has replaced coal with gas and electricity for providing heat in an orderly manner, reducing the consumption of coal. In accordance with the local conditions, we expand the variety of clean heating methods to ensure the balanced development of clean heating. In our future research, the SARIMA model will be combined with the random forest model to improve the accuracy of the model prediction.



---

*Source: 1276305-2020-12-23.xml*