Methodology for forecasting electricity consumption by Grey and Vector autoregressive models

Highlights • The Grey and Vector autoregressive models are coupled to improve their accuracy.• Five economic and demographic parameters are included in the new hybrid model.• This new model is a reliable forecasting tool for assessing energy demand.


Method details
In this section, we have set to describe the various steps that are require for modeling and forecasting electricity consumption from each of the previous models.

The Grey model
All the steps aiming at developing the GM(1,1) Grey model are described as follows: -Step one : Exponential smoothing of raw data.
It is worth noting that forecasts based on Grey models are satisfactory when raw data show an increasing exponential trend. In case the observed trend among basic data is not exponential, Yuan C. et al. [1] suggest to use an Exponential smoothing of raw data, and to only work with actual data whose relative gap is less than 0.5. In this paper, that stage is a sinequanon condition to executing the Grey GM (1,1) model.
All other stages required for the implementation of this model are described as follows [2] : -Step two : raw data processing.
Let X (0) be the raw data set defined by We further define the X (1) series known as 1-AGO sequence as: - Step three : estimate of model parameters a and b . This is done by considering the first form of the GM(1,1) model defined by Eq. (1) : where a and b are respectively the parameters of development and control of the model, and z (1) ( k ) the background sequence of the model given by the Eq. (2) : where k = 2 , 3 , …, n . Considering the values of k , we obtain the system ( s 0 ) . .
By introducing each equation of the previous system into Eq. (1) , the following system (S) is obtained . .
The system (S) then takes the matrix form Y = Br where: . .
Thus, with ordinary least squares method, parameters a and b of the following matrix are estimated: Step four : search for the temporal response function.
For this purpose, we consider the differential equation of the Grey model defined by Eq. (5) .
Solving the latter equation yields the time response function described by Eq. (6) .
where C is a constant that depends on the initially chosen condition. With a t , and the Eq. (6) appears as shown in relation (7) .
With t = 2,3,…,n. The IAGO sequence from the relation ˆ , allows to obtain the response function restored in equation (8) .

Forecasting with Grey model
The specificity of the Grey GM(1,1) model is that it provides fairly accurate forecasts based on a small range of actual data. For this reason, the rolling mechanism makes use of the most recent data to forecast future consumption values with greater accuracy [3] .
For a series ( x (0) (1); x (0) (2); …; x (0) ( n )) of n raw data, let's assume that z data are used to construct the GM(1,1) model, and q values are estimated to test the model's forecasting capacities. The procedure for running GM (1,1) is summarized in three steps as follows: Step one : we consider sequence ( x (0) (1); x (0) (2); …; x (0) ( z )) to calibrate the model by determining values a and b of the coefficients of the model corresponding to this series. Then, values Step two : Consider a new series of raw data of z length to recalibrate the model by calculating parameters a and b . Then Step three : repeat step 2 until the desired forecast horizon is reached. The more the RMSE value is closer to zero, the better the formulated model. The NS statistic is a performance indicator which helps estimating the capacity of a model to reproduce an observed behavior. For a model to be judged satisfactory, we must have NS ≥ 0 . 7 that is, the observed models must be coherent. Those indicators are calculated based on Eqs. (9) -(13) . Table 1 gives the MAPE assessment criteria for a statistical series.

Assessment of the model accuracy
where x (0) ( k ) are the actual values, x ˆ( 0) ( k ) the forecast values, and σ 2 the variance associated to the parameter X (0) . Other statistics are used in this work for the analysis of retained parameters.
These are namely the correlation coefficient ρ, used to capture the link between two parameters, and the Average Annual Growth Rate (AAGR) which shows the evolution percentage of a parameter over consecutive years.

The Vector autoregressive model
The multilinear Vector Autoregressive (VAR(p)) model formulated in 1980 by Christopher Sims is a generalization of the ARMA models formulated by Box and Jenkins in 1978 [4] . Here, the electricity consumption known as x ( t ) at a t date is specified as a function of past values x ( t − i ) and residuals This equation is a linear one, where electricity consumption solely depends on its previous values. Christopher Sims generalizes the latter relationship in the form generated by 15, whereby the consumption x(t) can depend not only on the previous values x(t -i), but also on some exogenous parameters including GDP y p (t − i ) , GDP per capita y ph (t − i ) , population y po (t − i ) , number of households y a (t − i ) , and number of subscribers . In this case, Eq. (14) takes form 15.
In this last equation a j 1 i for j ∈ [1 , 6] are coefficients of the model, a 1 0 is the constant term and μ t is the white noise error term . Provided the optimal lag number p retained by the statistic tests is 1, the electricity consumption on the t date may depend on all variables lagged by one period. From that perspective, Eq. (15) takes form 16.
x (t) = a 0 1 + a 1 11 x (t − 1) + a 2 11 y p (t − 1) + a 3 11 y ph (t − 1) + a 4 11 y po (t − 1) + a 5 11 y m (t − 1) + a 6 11 y a (t − 1) + μ(t) This result is an equation where parameters are chronological variables which often present a non-stationary evolution. Such chronological variables in a forecasting model, imply interdependence among errors' terms, as well as the bias of forecasting results. In that case, one or several filters at differences are necessary in order to make all exogenous and endogenous variables stationary [5] . If only one filter is needed to stabilize each parameter, then the previous equation takes form 17.
Such a model must then be statistically validated before attempting any forecast.

Validation of the Vector autoregressive model
Once the coefficients a j 1 i have been estimated, we perform three statistical tests of greater significance, in order to validate the formulated model. This is the unit root test based on Dickey Fuller Augmented (ADF) and Phillips-Perron (PP) statistics, the test on the normality of residuals using the Jarque-Bera statistic as defined by Eq. (18) , and the model stability testing using the cumulative sum of the recursive residuals S ( t ) and the cumulative sum of the squares of the recursive residuals S'(t) defined by formulas (a) (19) In this equation, μ˜j is the error term, and represents the difference between the actual model and the estimated model.
In fact, Provided we have level stationary variables at the end of our tests, they are designated I(0). If on the other hand, others are stationary after one or n first-differences, they are designated I(1) or I(n). Therefore, VAR(p) modelling is only possible if all observables are stable at the same order.
The formulated model is statistically stable when the evolution of the sum of the error terms is within the confidence interval β = ±95% at the critical threshold of α = 5%. Provided the Jarque-Berra statistics are less than 5 . 99 with a P-value greater than 5%, residuals are independent and normally distributed within the model. They therefore constitute a white Gaussian noise [5] .

Forecasting with the Vector autoregressive model
In a situation where the model parameters are stationary after a first difference, the forecasting approach in difference can be adopted in order to sequentially calculate future values of electricity consumption. This is done using the relationship given by Eq. (20) , where x ˆt is the estimated value of consumption on the t date, x t −1 is the revealed value of the consumption on t-1 date, and x is the difference of consumption between years t and, t-1 and estimated by the software The error associated with each estimate is digitally estimated from Eqs. (9) and (13) .

The hybrid model
The construction of the structure of the GM(1,1)-VAR(p) hybrid model is conditioned by the characteristics of the respective residuals of each of the models [1] . If their residuals are mutually opposite in signs, and if those of the VAR(p) model are normally distributed, then it is statistically possible to construct the hybrid model as defined by Eq. (21) .
This model may be sensitive to moderate fluctuations of electricity consumption. The following flow chart in Fig. 1 shows the key stages that are useful for the execution of the GM(1,1)-VAR(p) model.
The initialization stage requires using all entering data with z length.
-Step 1 requires the designation of the endogenous variables (dependent), as well as that of all exogenous variables(independent). Within the framework of this analysis, X (0) is the endogenous parameter, while the GDP, the GDP per capita, the population, the number of households, and the number of subscribers represent all the exogeneous parameters designated by Y i . -Moreover, step 2 suggests making sure that the exponential growth trend of all variables is effective.
In fact, Yuan et al. [1] have proved in 2016, that the GM(1,1) model makes it possible to obtain better forecast-ting results than when entering energy consumption data show an exponential growth trend. This idea is kept in this work as a sine qua non condition to executing the hybrid model. All retained parameters must have an exponential growth trend.
If not, an exponential smooth of entering data will be necessary, and the after-smoothing such as variables will be considered. The authors recommend retain the smooth values whose gap with actual data is below 0.5 [1] .
Once the exponential trend of parameters is established, -Step 6 which is the last, recommends displaying estimated data X ˆ(0) in such a way that  Accuracy of the hybrid model Table 2 allows a comparison between our results and those obtained by the authors of the hybrid models that are similar to the GM(1,1)-VAR(1) model. Considering their respective MAPEs, one can conclude that the new hybrid model presented in this paper offers a better accuracy and more reliable forecast capacities. Those interesting results are due to the fact that, the new hybrid model includes additional determinants of the demand, which are all primarily characterized by an exponential growth trend. This trend improves the forecasting quality of the model. Fig. 2 presents a comparison of the provisional capacities of the GM(1,1)-VAR(1) model, with the capacities of six other models which were reproduced and tested based on our data. It can be observed that over each year, the evolution of data filtered by the GM(1,1)-VAR(1) model on top right, and which is indicated by a brown vertical line is, in most cases, identical to the evolution of actual data regarding electricity consumption of each set of lines, over the 1997-2017 period. The new formulated model surpasses hybrid models of the same nature like GM(1,1)-ARMA(2,1) (in green), and the GM(1,1)-ARIMA(2,1,1) (in dark blue). Therefore, introducing additional determinants of the endogenous variable in such type of model fairly improves their accuracy.
Additional tests and comparisons have allowed us to observe that the new hybrid model can live up with some recent artificial intelligence models. Table 3 allows us to notice that, regarding the values of the various accuracy indicators, the new model is as reliable as hybrid models: the VMD-EELM (Variational Mode Decomposition-Evolutionary Extreme Learning Machine), VMD-SVM-PSO (the VMD which is coupled to the Support Vector Machine and optimized by the PSO algorithm), and VMD-SR-SVRCBCS(the VMD which is coupled to Self-recurrent mechanism and to Support Vector Regression and improved by the CBCS algorithm), recently developed to perform water-flow rate forecasts of hydroelectrical dams [7 , 8] , and electrical loads in China [9] . Type GM(1,N)-VAR(p) hybrid models, whose genesis is given in this paper (for N = p = 1), could prove to be as efficient as complex artificial intelligence models. Therefore, they can be further developed and tested within the annual forecasting rates of hydroelectrical dams, peak loads on electrical networks, and electricity demand at various stages, depending on the determinants with exponential growth.

Declaration of Competing Interest
The authors declare that they have no known competing financial interests or personal relationships that could have appeared to influence the work reported in this paper.