# Using a Novel Grey System Model to Forecast Natural Gas Consumption in China

Authors: Lifeng Wu; Sifeng Liu; Haijun Chen; Na Zhang
Journal: Mathematical Problems in Engineering
 (2015)

---

## Abstract

Accurate prediction of the future energy needs is crucial for energy management. This work presents a novel grey forecasting model that integrates the principle of new information priority into accumulated generation. This grey model can better reflect the priority of the new information theoretically. The results of two practical examples demonstrate that this grey model provides very remarkable short-term predication performance compared with traditional grey forecasting model for limited data set forecasting. It is applied to Chinese gas consumption forecasting to show its superiority and applicability.


---

## Body

## 1. Introduction

Energy has an influencing role in achieving economic and social progress. Forecasting energy constitutes a vital part of energy policy of a country, especially for a developing country like China whose economy is in a stage of energy transition: from low efficiency solid fuels to oil, gas, and electric power [1]. This has motivated many researchers to focus their research on energy forecasting. For example, Jian et al. used bayesian combination model to forecast energy demand of China [2]. Ji predicted petroleum consumption in China by comparing three models [3]. Meng and Niu analysed and predicted annual electricity consumption of China [4]. Zhang et al. forecasted Chinese transport energy demand based on partial least squares regression [5]. Ünler proposed swarm intelligence to improve energy demand forecasting [6]. A large number of studies on energy consumption forecasting using grey model and improved grey model are reported, such as Meng et al. [7], Bianco et al. [8], Li et al. [9], Niu et al. [10], Pao and Tsai [11], and Wang et al. [12].Although the first-order grey model with one variable (GM(1,1)) had been widely adopted, its predicting performance still could be improved. In this paper, based on the principle of new information priority, we consider the priority of thenth component and the n-1th component, respectively, and the nth component is more prior than the n-1th component. Similarly, the 2nd component is more prior than the 1st component [13]. Two case studies show the effectiveness of the proposed modeling method, which is particularly suitable for limited sample forecasting.The rest of the paper proceeds as follows: an overview of the relevant grey models and a new GM(1,1) model with the principle of new information priority of grey system are present in Section2. The advantage of the grey model proposed in this paper over the other grey models is proved by two real cases in China in Section 3. The prediction concerning natural gas consumption in China will be conducted in Section 4. Some conclusions are discussed in final section.


## 2. Methodology

### 2.1. Grey Model Based on Autoregressive Method

AssumeX={x(1),x(2),…,x(n)}, where X stands for a data sequence. The equation
(1)xk=αxk-1+β,k=2,3,…,n,
is called the autoregressive GM(1,1) model (ARGM(1,1)). Use the ordinary least squares method to estimate the parameters as follows:
(2)α^β^=BTB-1BTY,
where
(3)Y=x(2)x(3)⋮xn-1xn,B=x11x21⋮⋮xn-21xn-11.Theorem 1.
Assume thatα^, β^ are the same as (1). Set x^(1)=x(1); then recursive function is given by
(4)x^k+1=α^kx1-β^1-α^+β^1-α^,k=1,2,…,n-1.Proof.
Substituteα^, β^ into (1); then
(5)x^k+1=α^x^k+1+β^=α^x^k+β^+β^x^k+1=⋯=α^kx^1+α^k-1+α^k-2+⋯+α^+1β^.
Letx^(1)=x(1); then
(6)x^k+1=α^kx1+1-α^k1-α^β^.
We also obtain a nonhomogenous exponential model as follows:(7)x^k=α^k-1x1-β^1-α^+β^1-α^,k=1,2,…,n-1.
Because grey autoregressive method of GM(1,1) is formulated by using the original data rather than the accumulated generation data, it does not need to inverse accumulated generating. We can see that (1) is an exponential model and (7) is a nonhomogenous exponential model. This method is called the ARGM(1,1).


### 2.2. The GM(1,1) Model with the Principle of New Information Priority via Weighted Least Squares Method

In this section, by weighted least squares method, the coefficients of (1) can be obtained as follows:
(8)Minimize:∑k=2nwkxk-αxk-1-β2,wk=k-1n-1.Then we obtain a nonhomogenous exponential model as follows:(9)x^k=α^k-1x1-β^1-α^+β^1-α^,k=1,2,…,n-1.w
k can reflect the principle of new information priority, so (9) also pays more attention to the new error and embodies the principle of new information priority. We call this method WLSGM(1,1).


### 2.3. The GM(1,1) Model with the Principle of New Information Priority

AssumeX={x(1),x(2),…,x(n)}, where X stands for a data sequence. The original data of (1) is
(10)X=x2x11x3x21⋮⋮⋮xn-1xn-21xnxn-11.The data with the principle of new information priority is(11)11⋯1101⋯11⋮⋮⋮⋮⋮00⋯1100⋯01x2x11x3x21⋮⋮⋮xn-1xn-21xnxn-11=AX=∑k=2nxk∑k=2nxk-1n-1∑k=3nxk∑k=3nxk-1n-2⋮⋮⋮∑k=n-1nxk∑k=n-1nxk-12xnxn-11,
which is called the treatment of new information priority. The equation
(12)∑i=knxi=α∑i=knxi-1+n-k+1β,k=2,3,…,n,
is referred to as grey nonhomogenous exponential model with the principle of new information priority (abbreviated as NIGM(1,1)).Theorem 2.
Assume thatX={x(1),x(2),…,x(n)} is certain. Actually, if the error is μ, then
(13)xk=αxk-1+β+μk,k=2,3,…,n.
Via the treatment of new information priority, the mean of the error is(14)EAμ=μn+n-2μn-1n-1+⋯+2μ3n-1+μ2n-1,
where A is the same as (11). Use the ordinary least squares method to estimate the parameters as follows:
(15)α^β^=BTB-1BTY,
where
(16)Y=∑k=2nxk∑k=3nxk⋮∑k=n-1nxkxn,B=∑k=2nxk-1n-1∑k=3nxk-1n-2⋮⋮∑k=n-1nxk-12xn-11.
The variable of(17)α^β^
is Var⁡α^β^=(BTATAB)-1BTATE(AμμTAT)AB(BTATAB)-1, where
(18)EAμμTAT=∑i=2nμi2μ3∑i=2nμi⋯μn[∑i=2n‍μi]μ3∑i=2nμi∑i=3nμi2⋯μn[∑i=3n‍μi]⋮⋮⋮⋮μn∑i=2nμiμn∑i=3nμi⋯μn2.Proof.
Via the treatment of new information priority, (11) becomes
(19)AY=ABαβ+Aμ.
Obviously,(20)EAμ=μn+n-2μn-1n-1+⋯+2μ3n-1+μ2n-1.
Use the ordinary least squares method to estimate the parameters as follows:(21)α^β^=BTATAB-1BTATAY.
Substituting(22)AY=ABαβ+Aμ
into (21), we obtain
(23)α^β^=BTATAB-1BTATABαβ+Aμ=αβ+BTATAB-1BTATAμ.
Then(24)Eα^β^=αβ+EBTATAB-1BTATAμ=αβ+BTATAB-1BTATEAμ.
Actually, by(25)Eα^β^
and E(Aμ), we obtain
(26)lim⁡n→∞n-2n-1-n-3n-1=0,⋯,lim⁡n→∞2n-1-1n-1=0.
The more the sample is, the less the attention we pay to the error from new information and vice versa.
We can obtain(27)Var⁡α^β^=Var⁡α^Covβ^,α^Covα^,β^Var⁡β^.
By(28)α^β^-αβ=BTATAB-1BTATAμ,
thus
(29)Var⁡α^β^=BTATAB-1BTATE×AμμTATABBTATAB-1.
In most cases, we pay attention to theVar⁡(α^) and Var⁡(β^) and overlook Cov(β^,α^); thus we denote (27) as
(30)Var⁡α^β^=BTATAB-1EAμμTAT.
So(31)EAμμTAT=∑i=2n‍μi2μ3[∑i=2n‍μi]⋯μn∑i=2n‍μiμ3∑i=2n‍μi∑i=3n‍μi2⋯μn∑i=3n‍μi⋮⋮⋮⋮μn∑i=2n‍μiμn∑i=3n‍μi⋯μn2.LetE(AμμTAT)=D; then the components from D11 to Dnn all contain μn2. The components from D11 to Dn-1n-1 all contain μn-12. The components from D11 to Dn-2n-2 all contain μn-22. The components from D11 to Dn-1n-1 all contain μn-12. Analogizing successively, the components D11 and D22 both contain μ22. Only D11 contains μ12. So
(32)Var⁡α^β^
pay more attention to the variable of the new information in the modelling process. The modeling procedure of NIGM(1,1) can be concluded as follows.Step  1. An original nonnegative series {x(1),x(2),…,x(n)} is given.Step  2. Y and B which are the same as Theorem 2 can be obtained via the treatment of new information priority.Step  3. Least squares estimation of α and β according to Theorem 2 is as follow. 								(33)α^β^=BTATAB-1BTATAY.Step  4. We obtain a nonhomogenous exponential model as follows:
(34)x^k=α^k-1x1-β^1-α^+β^1-α^,k=1,2,…,n-1.Step  5. In order to solve the effect of the initial condition, we assume
(35)xk=cα^k+b.Use the ordinary least squares method to estimate the parameters(36)c^b^=BTB-1BTY,where(37)Y=x1x2⋮xn-1xn,B=α^1α^21⋮⋮α^n-11α^n1;thus, we obtain(38)x^k=c^α^k+b^,k=1,2,…,n,which is a nonhomogenous exponential model.


## 2.1. Grey Model Based on Autoregressive Method

AssumeX={x(1),x(2),…,x(n)}, where X stands for a data sequence. The equation
(1)xk=αxk-1+β,k=2,3,…,n,
is called the autoregressive GM(1,1) model (ARGM(1,1)). Use the ordinary least squares method to estimate the parameters as follows:
(2)α^β^=BTB-1BTY,
where
(3)Y=x(2)x(3)⋮xn-1xn,B=x11x21⋮⋮xn-21xn-11.Theorem 1.
Assume thatα^, β^ are the same as (1). Set x^(1)=x(1); then recursive function is given by
(4)x^k+1=α^kx1-β^1-α^+β^1-α^,k=1,2,…,n-1.Proof.
Substituteα^, β^ into (1); then
(5)x^k+1=α^x^k+1+β^=α^x^k+β^+β^x^k+1=⋯=α^kx^1+α^k-1+α^k-2+⋯+α^+1β^.
Letx^(1)=x(1); then
(6)x^k+1=α^kx1+1-α^k1-α^β^.
We also obtain a nonhomogenous exponential model as follows:(7)x^k=α^k-1x1-β^1-α^+β^1-α^,k=1,2,…,n-1.
Because grey autoregressive method of GM(1,1) is formulated by using the original data rather than the accumulated generation data, it does not need to inverse accumulated generating. We can see that (1) is an exponential model and (7) is a nonhomogenous exponential model. This method is called the ARGM(1,1).


## 2.2. The GM(1,1) Model with the Principle of New Information Priority via Weighted Least Squares Method

In this section, by weighted least squares method, the coefficients of (1) can be obtained as follows:
(8)Minimize:∑k=2nwkxk-αxk-1-β2,wk=k-1n-1.Then we obtain a nonhomogenous exponential model as follows:(9)x^k=α^k-1x1-β^1-α^+β^1-α^,k=1,2,…,n-1.w
k can reflect the principle of new information priority, so (9) also pays more attention to the new error and embodies the principle of new information priority. We call this method WLSGM(1,1).


## 2.3. The GM(1,1) Model with the Principle of New Information Priority

AssumeX={x(1),x(2),…,x(n)}, where X stands for a data sequence. The original data of (1) is
(10)X=x2x11x3x21⋮⋮⋮xn-1xn-21xnxn-11.The data with the principle of new information priority is(11)11⋯1101⋯11⋮⋮⋮⋮⋮00⋯1100⋯01x2x11x3x21⋮⋮⋮xn-1xn-21xnxn-11=AX=∑k=2nxk∑k=2nxk-1n-1∑k=3nxk∑k=3nxk-1n-2⋮⋮⋮∑k=n-1nxk∑k=n-1nxk-12xnxn-11,
which is called the treatment of new information priority. The equation
(12)∑i=knxi=α∑i=knxi-1+n-k+1β,k=2,3,…,n,
is referred to as grey nonhomogenous exponential model with the principle of new information priority (abbreviated as NIGM(1,1)).Theorem 2.
Assume thatX={x(1),x(2),…,x(n)} is certain. Actually, if the error is μ, then
(13)xk=αxk-1+β+μk,k=2,3,…,n.
Via the treatment of new information priority, the mean of the error is(14)EAμ=μn+n-2μn-1n-1+⋯+2μ3n-1+μ2n-1,
where A is the same as (11). Use the ordinary least squares method to estimate the parameters as follows:
(15)α^β^=BTB-1BTY,
where
(16)Y=∑k=2nxk∑k=3nxk⋮∑k=n-1nxkxn,B=∑k=2nxk-1n-1∑k=3nxk-1n-2⋮⋮∑k=n-1nxk-12xn-11.
The variable of(17)α^β^
is Var⁡α^β^=(BTATAB)-1BTATE(AμμTAT)AB(BTATAB)-1, where
(18)EAμμTAT=∑i=2nμi2μ3∑i=2nμi⋯μn[∑i=2n‍μi]μ3∑i=2nμi∑i=3nμi2⋯μn[∑i=3n‍μi]⋮⋮⋮⋮μn∑i=2nμiμn∑i=3nμi⋯μn2.Proof.
Via the treatment of new information priority, (11) becomes
(19)AY=ABαβ+Aμ.
Obviously,(20)EAμ=μn+n-2μn-1n-1+⋯+2μ3n-1+μ2n-1.
Use the ordinary least squares method to estimate the parameters as follows:(21)α^β^=BTATAB-1BTATAY.
Substituting(22)AY=ABαβ+Aμ
into (21), we obtain
(23)α^β^=BTATAB-1BTATABαβ+Aμ=αβ+BTATAB-1BTATAμ.
Then(24)Eα^β^=αβ+EBTATAB-1BTATAμ=αβ+BTATAB-1BTATEAμ.
Actually, by(25)Eα^β^
and E(Aμ), we obtain
(26)lim⁡n→∞n-2n-1-n-3n-1=0,⋯,lim⁡n→∞2n-1-1n-1=0.
The more the sample is, the less the attention we pay to the error from new information and vice versa.
We can obtain(27)Var⁡α^β^=Var⁡α^Covβ^,α^Covα^,β^Var⁡β^.
By(28)α^β^-αβ=BTATAB-1BTATAμ,
thus
(29)Var⁡α^β^=BTATAB-1BTATE×AμμTATABBTATAB-1.
In most cases, we pay attention to theVar⁡(α^) and Var⁡(β^) and overlook Cov(β^,α^); thus we denote (27) as
(30)Var⁡α^β^=BTATAB-1EAμμTAT.
So(31)EAμμTAT=∑i=2n‍μi2μ3[∑i=2n‍μi]⋯μn∑i=2n‍μiμ3∑i=2n‍μi∑i=3n‍μi2⋯μn∑i=3n‍μi⋮⋮⋮⋮μn∑i=2n‍μiμn∑i=3n‍μi⋯μn2.LetE(AμμTAT)=D; then the components from D11 to Dnn all contain μn2. The components from D11 to Dn-1n-1 all contain μn-12. The components from D11 to Dn-2n-2 all contain μn-22. The components from D11 to Dn-1n-1 all contain μn-12. Analogizing successively, the components D11 and D22 both contain μ22. Only D11 contains μ12. So
(32)Var⁡α^β^
pay more attention to the variable of the new information in the modelling process. The modeling procedure of NIGM(1,1) can be concluded as follows.Step  1. An original nonnegative series {x(1),x(2),…,x(n)} is given.Step  2. Y and B which are the same as Theorem 2 can be obtained via the treatment of new information priority.Step  3. Least squares estimation of α and β according to Theorem 2 is as follow. 								(33)α^β^=BTATAB-1BTATAY.Step  4. We obtain a nonhomogenous exponential model as follows:
(34)x^k=α^k-1x1-β^1-α^+β^1-α^,k=1,2,…,n-1.Step  5. In order to solve the effect of the initial condition, we assume
(35)xk=cα^k+b.Use the ordinary least squares method to estimate the parameters(36)c^b^=BTB-1BTY,where(37)Y=x1x2⋮xn-1xn,B=α^1α^21⋮⋮α^n-11α^n1;thus, we obtain(38)x^k=c^α^k+b^,k=1,2,…,n,which is a nonhomogenous exponential model.


## 3. Results

In this section, the advantage of the NIGM(1,1) model over the other grey models is demonstrated by two real case studies of energy consumption. Mean absolute percentage error(MAPE=100%(1/n)∑k=1n(x(k)-x^(k))/x(k)) compares the real and forecasted values to evaluate the precision.Case 1 (natural gas consumption forecasting example in China [14]).
This example is from the literature [14], in order to compare the result of the literature [14] with the result of the proposed model. The sample data is the same as the literature [14]. Using the data from 1995 to 2006 to construct NIGM(1,1),
(39)xk=22.478*1.276k+144.995.
Then the value from 2007 to 2008 (out-of-sample) is predicted. Actual values and the forecasting values of five compared models are presented in Table1. As can be seen from Table 1, NIGM(1,1) yielded the lowest MAPE compared with the other models during the period from 1995 to 2008.Table 1
The fitted values (unit: 100 million cubic meters) and MAPE of five grey models.
Year
Actual value
Model (39)
ARGM(1, 1)
Reference [14]
WLSGM(1, 1)
DGM(1, 1)
1995
177.53143
173.677
177.53143
190.24323
177.53143
177.53143
1996
188.04992
181.593
187.961
183.65849
190.990
149.31
1997
176.13278
191.694
201.029
184.90875
207.706
169.38
1998
218.69985
204.583
217.403
193.99401
228.466
192.15
1999
215.3371
221.030
237.920
210.91426
254.251
217.98
2000
244.81169
242.015
263.628
235.66952
286.276
247.29
2001
274.55447
268.793
295.840
268.25976
326.051
280.53
2002
292.18069
302.961
336.202
308.68501
375.451
318.25
2003
339.45429
346.559
386.775
356.94525
436.806
361.03
2004
397.28586
402.191
450.143
413.04049
513.010
409.57
2005
486.93617
473.177
529.543
476.97073
607.654
464.63
2006
561.05120
563.756
629.031
548.73596
725.202
527.10
MAPE
3.092
8.657
4.406
18.173
6.96
2007
698.903
679.333
753.690
679.955
871.197
597.96
2008
807.857
826.811
909.888
774.429
1052.523
678.35
MAPE
2.573
10.234
10.75
27.469
15.236Case 2 (the transport energy demand forecasting in China [5]).
The example from the literature [5] which provides the sample data is considered. The data from 1990 to 2004 are used to construct NIGM(1,1):
(40)xk=11.315*1.168k+32.49.
Then the value from 2005 to 2006 (out-of-sample) is predicted. Actual values and the forecasting values of five compared models are presented in Table2.
As can be seen from Table2, NIGM(1,1) reduces the MAPE of original GM(1,1) from 4.82% to 3.7%. This indicates that NIGM(1,1) can improve the forecasting error of GM(1,1) but not the partial least squares regression method in the in-sample model. From a short-term forecasting viewpoint, NIGM(1,1) yielded the lowest MAPE from 2005 to 2006 compared with the other five models, which means that the NIGM(1,1) model reaches the objective of minimizing the forecast error and has highly accurate forecasting power. Therefore, we can conclude that NIGM(1,1) significantly enhances the precision of grey forecasting model.Table 2
Fitted values and MAPE of different models.
Year
Actual value
Model (40)
GM(1, 1)
Reference [5]
WLSGM(1, 1)
DGM(1, 1)
1990
45.41
44.71
45.41
44.07
45.41
45.41
1991
47.56
47.93
42.18
46.40
47.14
42.26
1992
50.58
50.52
46.32
49.48
49.19
46.40
1993
51.85
53.55
50.85
53.16
51.62
50.94
1994
56.25
57.09
55.83
57.49
54.50
55.93
1995
58.63
61.22
61.30
60.83
57.91
61.41
1996
59.94
66.04
67.31
60.62
61.97
67.42
1997
75.43
71.68
73.90
77.07
66.77
74.02
1998
82.45
78.27
81.14
82.15
72.48
81.27
1999
93.40
85.96
89.08
90.44
79.24
89.23
2000
100.67
94.94
97.81
98.72
87.26
97.96
2001
103.63
105.43
107.39
107.7
96.77
107.55
2002
111.71
117.68
117.91
116.93
108.05
118.09
2003
128.19
131.99
129.45
123.54
121.42
129.65
2004
151.04
148.71
142.13
146.05
137.29
142.34
MAPE
3.7
4.82
2.6
5.88
4.76
2005
166.29
168.24
156.06
164.21
156.11
156.28
2006
185.83
191.05
171.34
179.51
178.42
171.58
MAPE
1.99
3.48
2.3
5.06
3.42


## 4. Discussion and Policy Implications

Many experts have developed methods to forecast natural gas consumption [15–23] because, in recent years, the Chinese government has begun to recognize the disadvantages of the high reliance on coal. Current policies are different from previous policies. According to current policies, natural gas will be developed as a partial substitution for coal. From a short-term forecasting viewpoint, we must pay more attention to the recent data. So we choose data from 1996 to 2013 to construct NIGM(1,1):
(41)xk=122.2219*1.14736k-30.8291.We use the latest (2012-2013) value to compare the estimated value with the actual data; the result is listed in Table3.Table 3
Fitted values and MAPE of model (41).
Year
Actual value
Model (41)
2012
1463.0
1420.5
2013
1616.1
1634.3
MAPE
2.0Model (41) is used to predict the Chinese natural gas consumption based on the above comparison. The forecasted values and the trend from 2014 to 2018 are given in Table 4.Table 4
The forecasting values of model (41).
Year
Model (41)
2014
1879.7
2015
2161.2
2016
2484.3
2017
2854.9
2018
3280.2It can be seen that there is an increase in trend in the consumption data from Table4. Continued fast growth makes it important to ensure the production and availability of the amount of natural gas to meet life’s need. Further effort must be made to achieve a sustainable development of China's natural gas industry.


## 5. Conclusion

The incompleteness of information is the primary characteristic of grey system theory, so we put more emphasis on recent data. A heavier weight is assigned to the recent data with two cases which are used to examine the performance of new grey model. It may be used for other real cases for energy consumption forecasting.



---

*Source: 686501-2015-01-15.xml*