事实上，从下面的图形 3-5 中也可以看到原始数列的模拟值与原始数列之间拟合的相当好。
3.8 3.7 3.6 3.5 3.4 3.3 3.2 3.1 3 1 2 3 4 5 6 History data:real line
ˆ 图 3-5 虚线：模拟值 x
( 0)
曲线；实线：原始值 x
( 0)
曲线
在本例中，预测该市 2006 年的第三产业从业人数，则有：
ˆ ( 0 ) (7) x (1) (7) x (1) (6) 24.10 20.25 3.85 （万人） x
GM(1,1)模型参数估计、检验以及作图和预测的 Matlab 程序：
fungry1.m %GM（1,1）模型计算及检验、作图。文件名 fungry1.m function GM1=fungry1(x0) %输入原始数据 x0 T=input('T='); %从键盘输入从最后一个历史数据算起的第 T 时点 x1=zeros(1,length(x0));B=zeros(length(x0)-1,2); yn=zeros(length(x0)-1,1);Hatx0=zeros(1,length(x0)+T); Hatx00=zeros(1,length(x0));Hatx1=zeros(1,length(x0)+T); epsilon=zeros(length(x0),1);omega=zeros(length(x0),1); for i=1:length(x0) for j=1:i x1(i)=x1(i)+x0(j); end end for i=1:length(x0)-1 B(i,1)=(-1/2)*(x1(i)+x1(i+1)); B(i,2)=1; yn(i)=x0(i+1); end HatA=(inv(B'*B))*B'*yn %GM（1,1）模型参数估计 for k=1:length(x0)+T Hatx1(k)=(x0(1)-HatA(2)/HatA(1))*exp(-HatA(1)*(k-1))+HatA(2)/HatA(1); end Hatx0(1)=Hatx1(1); for k=2:length(x0)+T Hatx0(k)=Hatx1(k)-Hatx1(k-1); %累减还原得到历史数据的模拟值 end for i=1:length(x0) %开始模型检验 epsilon(i)=x0(i)-Hatx0(i); omega(i)=(epsilon(i)/x0(i))*100; end % x0；Hatx0; epsilon; omega; %必要时去掉%得到各种数据 c=std(epsilon)/std(x0); p=0; for i=1:length(x0) if abs(epsilon(i)-mean(epsilon))<0.6745*std(x0) p=p+1; end end p=p/length(x0) if p>0.95 & c<0.35 disp('The model is good,and the forecast is:'), disp(Hatx0(length(x0)+T)) elseif p>0.85 & c<0.5 disp('The model is eligibility,and the forecast is:'),
通过累加生成的数列 x
(1)
（3-16）
，计算模型参数 a 和 u 。记：
ˆ [a u ]T a
按如下公式可得模型参数 a 和 u ：
（3-17）
ˆ ( B T B ) 1 B T y n a
上式中：
（3-18）
1 (1) (1) 1 2 ( x (1) x ( 2)) 1 (1) (1) ( x ( 2 ) x ( 3 )) 1 B 2 1 ( x (1) ( n 1) x (1) ( n )) 1 2 y n [ x ( 0 ) ( 2), x ( 0 ) (3), , x ( 0 ) ( n )]T
c S1 S 0 0.0210 0.2676 0.0785
最后计算小误差概率 p ：
（3-27）
p ( 0 ) ( 0 ) 0.6745 S 0
本例中， S 0 0.2676,
( 0)


（3-28）
0 ，于是，
p ( 0 ) 0.6745 0.2676 ( 0 ) 0.1805 1
( 0)
ˆ ( 0 ) ( n 2)， ( n 1)，x 的预测值（Matlab 程序见后，文件名 fungry1.m）
下面结合一个具体例子加以说明 GM（1,1）模型的建立和检验过程。 【例 3-13】某市第三产业从业人数如表 3-8 所示，试建立 GM（1,1）模型，预测其 2006 年的从业人数。 表 3-8 某市从事第三产业的人数（单位：万人） 年份 人数 2000 2.97 2001 3.23 2002 3.29 2003 3.46 2004 3,59 2005 3.71
(1) ( 0)
作累加
，即作：
x (1) (i ) x ( 0 ) ( m ),
m 1
i
i 1,2, , n
（3-15）
具体地说，就是：
x (1) (1) x ( 0 ) (1) (1) ( 0) (1) x (i ) x (i ) x (i 1), i 2, , n
ˆ ( 0 ) (i ) x
2.97 3.21 3.33 3.45 3.58 3.71
ˆ ( 0 ) (i ) ( 0 ) (i ) x ( 0 ) (i ) x
0 0.02 -0.04 0.01 0.01 0
相对误差（%）
0 0.62 -1.22 0.29 0.28 0
由误差计算结果可以看到，相对误差不超过 2%，模型精度是非常高的。 2．关联度检验。关联系数定义为：
( 0)
ˆ ( 0 ) (i ), (i ) x ( 0 ) (i ) x
(i 1,2, , n ) (i 1,2, , n ) （3-22） (i 1,2, , n )
( 0 ) (i ) 相对误差： (i ) ( 0 ) 100%, x (i )
x ( 0)
然后计算残差数列
( 0)
1 n ( 0) x (i ) n i 1
的均方差 S1 。其定义为：
S1
S12 n 1
n
S12 ( 0 ) (i ) ( 0 )
i 1


2
（3-26）
( 0)
1 n ( 0) (i ) n i 1
本例中， S 0 0.2676, S1 0.0210 。由此计算方差比 c ，
从而， a 0.03653, u 3.04123 这就是说，该市第三产业从业人数的 GM（1,1）预测模型为：
ˆ (1) (i 1) 86.2281e 0.03653i 83.2581 x
下面对得到的模型进行检验。 灰色预测模型检验有残差检验、 关联度检验和后验差检验。 1． 残差检验。残差检验有两种，绝对误差和相对误差检验。 绝对误差：
一般来说，在 0.5, 0.695 时是令人满意的（具体讨论参见灰色系统理论的有关 参考书） 。 3．后验差检验。首先计算原始数列 x
( 0)
的均方差 S 0 。其定义为：
S0
n
2 S0 n 1
2 S0 x ( 0 ) (i ) x ( 0 ) i 1


2
（3-25）
1 1 1 1 1


T
(3.23,3.29,3.46,3.59,3.71) T ，于是有
4.585 7.845 4.585 7.845 11 .22 14.745 18.395 BT B 11 .22 1 1 1 1 1 14.745 18.395 764.2457 56.79 5 56.79 764.2457 56.79 ( B B) 5 56.79
x ( 0 ) (i )
x ( 0 ) (1)
x ( 0 ) ( 2)
x ( 0 ) (3)
x ( 0 ) ( 4)
( 0)
x ( 0 ) (5)
x ( 0 ) ( 6)
如果利用手工计算，则有以下步骤。首先计算累加生成数列 x 表 3-9 序号 i 1 2.97 2.97 由x 2 3.23 6.2
（i）
min((i )) max( (i )) , (i ) max( (i ))
i 1,2, , n
（3-23）
式中，
(i ) 为第 i 个数据的关联系数， 为取定的最大差百分比，一般取 0.5 ；
ˆ 而： (i ) x
( 0)
(i ) x ( 0 ) (i ) , i 1,2, n
1 GM（1,1）模型
设有预测对象的历史数据资料形成的时间序列 x 结果想知道： x
( 0) ( 0)
: x ( 0 ) (1), x ( 0 ) ( 2), x (), x ( 0 ) ( n 2),
根据灰色系统理论，可以针对这个问题建立如下的 GM（1,1）模型。模型所依据的灰色 系统理论在这里就不作介绍，下面介绍具体的模型建立方法和检验方法。首先对 x 生成，得到新的数列 x
本例中，取 0.5 ，经过计算有如表 3-11 所示的结果。