灰色预测专设工⑼他QA—叫吋)为原始数列.其1次累❖加生成数列为恥=妙①曲⑵，…卅何),其中X° 仇)二工* ° (0.址=1=2= -:n5-1卷定义卫的决导数为d(k) = *町(上)=x 叫咼-x cl)(Jt-l).令为数列工①的邻值生成数列.即却(去)=^(*) + (1- a)x0)(t-lX于是定义GM (L 1)的灰微分方程模型为d(k)-血⑴住)=K即或严＞(£) + “尹⑻=人⑴在式(1)中』。

＞(灼称为灰导数，我称为发展系数, 弧称为白化背景值，b称为灰作用量乜将时刻表殳二2「3「/代入(1)式有V!1「—ay=代⑶ B =Ib*- :X闵0)-Z，：］(K)1于是G\I <1»1）複至可表示为Y = Bu.現在问题归结为求sb 在值。

用一元线性回归・即最小二秦法求它们的活计值 为注二实陌上回归分析中求估计值是用软件计尊的・有标准程序求解，iOmaClab 等。

GM <1» 1>的白化晏対于G\I <1> 1）的灰微分方程（1） >如果将灰导数打（Q 的时刻 视为连绫变里"则x°）视为时问（函数卅⑺,于是*〉（Q 対血于导数里级 心2 >白化背臬值申的对应于导数卅⑴。

于是G\I （1，1）的坝徽 分方樂対应于的白微分方程为内・则数堀列X©可以塗互G\I <19 1） 且可以进行页色预测。

否朋，対数摄做适当的克换处理■如平移叢换:取C 使得鞍据列严伙）=工⑴伙）+ G 上=1,2,…,的级比都華住可吝禎盖内。

心⑴⑴ + o?i> (r)二◎ dr<2)GM mi ）质色预测的步骤1 •教摇的枪绘与处連为了ftilGAl （1,1）建複方法的可行性，亲要为已知期S 做必要的检蛉处理。

设原始教据列为了 逛=（乂°（1）*6（2）严炉00； >计算数列的级比如果所有的级比都落在可容覆盖区间 • fc =A-2,3"・如果対所有的|p 伙)|<0・1 -则认为达到较高的要求,否则 若旳所有的|。

(鋼<0・2 >则认为送到一般要求。

2.逢立G\I C1F 1）樓型不妨设少=(由(1)/®(2〉,…¥切)彌足上面的要 求，以它为埶据列建立G 、I CU 1)複塑J°)仏)十血⑴(上)=b.用回归分析求得场b 的佔计值，于是相应的白化僅型为絹为dthk *）© = （2 ⑴-？上 YE +三 a a于是得到预则值护住十1)=(少(1)-2)穴+色,上=12…1,aa从而相应地得到预测值：x <0) (^ +1) = x (1)(^ +l)-x o )(t)5 k = 1,2,…,M-1,s.检验预测值(1)死差检验二计尊相对残差如果对所有的 ^<0.1，认为达到技高的要求,否则，若 対所有的|f(/t)|<02 >则认为达到一骰要求。

(2)级比偏差值检验；计算皿)=1上0% l+0.5a 律）,•:•例北方某城市1986〜1992年道路交通噪声平均 芦级数据见表6 第一步:级比检验 建立交通噪声平均声级数粥时间序列如F ：x (0) = (2),…工⑼⑺)-(71.1, 72A 724, 72L71A 72.0,71.6)⑴求级比/住)心少(上-1) • 2 = (A(2)/(3),…,班7» ‘丿 =(0.982.1.1.0042.1.0098. 0.9917.1.0056) (2)级比判断由于所有的入(k)曰0.982,1.0098], k=2,3、 ，7,故可以用x(0)作满意厉GM (b 1)建模。

第二步:GM (1, 1)建模(1) 对原始数据X®作一次累加，即 X 。

)=(71.1443.5215.9.28&359.4.431.4.503) (2) 构造数据矩阵B 及数据向量Y•:♦表6市近年来交通噪声数据[dB(A)]序号年份旳z1 1996 71.12 1997 72.43 1988 72.4 4 1989 72.15 1990 71.46 1991 72.07 1992 71.6-駅“弋 6 ) + x (l < 7 )(3) ifMif= ^r ^)-1B r y =(0-0023)<4*®7V72.6573dr+ 0.0023 丿〉=72.6573 dt❖求解得 沙⑴= (x"°)(l)-0w + £.aar^(Jt +l) = (x^Xl)--)^ +?=-30929&7陀％ +31O QO(5)求生成教列值: i w (上+1)令k=HS4J6由上面的南冏五鈕数可算得，其中取cX 由 丘⑴(1) =沁玖得？ 11讣(上)=尹(上)-炉(上-1),G1-1,72.4: 72.1 72.1, 719： 71.7, 71 6)曲= 3(1)，旳(2)，…，旳(7)>=•:・第三步:模型检验•:•模型的各种检验指标値的计算结果见裘7. •:・表7 G 、I(1,1濃型检验表・：・序号年份原始値模型值残差相对误差级比偏差•:• 1 1986 71.1 71.1 0 0❖ 21987 72.4 72.4 -0.0057 0.01%0.0023=b r 1 12> 0 3>•:・3198872.472.20.16380.23%0.0203• 4198972.172.10.03290.05%-0.0018❖5199071.471.9-0.49840.7%-0.0074•:- 6199172.071.70.26990.37%0.0107❖7199271.671.60.0378 0.05%-0.0032•:•经验证，该模型的精度较高，可进行预测和预报。

计算的M 耐程序如下■蚁和预报lamda = x0(l ：n^l)./id )(2：n) % 计算级比 range = minnax( lamda*) %计算锁比的范国 xl =cumsum(xO)； % 累加运算B = f -0^5* (xl(l ;n - 1) +xl(2：n))t ones(n-lj)]； V = xO(2：n) * u = B\Y % 拟合 ^»u(l)^atU (2)=bx=dsolve(^Dx + a*x = b\f x(Q)=池')< 绻戏徴裁方程的挣号餐X = subfi(xJ ,a ,Jb ,/xO ,| J u( 1)+ u(2) ,xO( I) I )； % 代入牯计赛魏值和初始值 yucel = Bubs (x 「t/* [0： n - 1 ]): %求巳他第据的预測值y = vpa(x t 6) %其中曲6表示显示&位教字 yuce = [xO(l)f diff(yucel)] * 差分运勲还原擬据 delta - abs( epsilon % 讣算祀对课產rbo = l-(l-0.5*u(l))41 心 F ⑴厂“京 Y^rhuJst 模型作用：求累加数列、求ab 的值、求预测方程、求残差clc %清屏，以使结果独立显示x=[71.1 72.4 72.4 72.1 71.4 72.0 71.6]; format long ; if len gth(x(:,1))==1 x=x'; end% 设置计算精度 %对输入矩阵进行判断，如不是一维列矩阵，进行转置变换n=len gth(x); z=0; for i=1: n z=z+x(i,:); be(i,:)=z; %取输入数据的样本量%计算累加值，并将值赋予矩阵 beendfor i=2: n %对原始数列平行移位y(i-1,:)=x(i,:);endfor i=1: n-1%计算数据矩阵 B 的第一列数据 c(i,:)=-0.5*(be(i,:)+be(i+1,:));xO=[71 J 72,4 72 A 72 Jn =length(xO)；71 A 72 .0 71 #6]\%注意疑里为列向epsilon - xO" — yuce % 计算残差 屈丫 t 计鼻奴比備垦扰・讥1endfor j=1:n-1 e(j,:)=1;%计算数据矩阵B 的第二列数据endfor i=1:n-1B(i,1)=c(i,:);B(i,2)=e(i,:);%构造数据矩阵Bendalpha=inv(B'*B)*B'*y; for i=1:n+1 值%计算参数矩阵即a b 的值%计算数据估计值的累加数列，如改为n+1 为n+m 可预测后m-1 个ago(i,:)=(x(1,:)-alpha(2,:)/alpha(1,:))*exp(-alpha(1,:)*(i-1))+alpha( 2,:)/alpha(1,:); %显示输出预测值的累加数列endvar(1,:)=ago(1,:) for i=1:n %显示输出预测值%如改n 为n+m-1 ，可预测后m-1 个值var(i+1,:)=ago(i+1,:)-ago(i,:); %估计值的累加数列的还原，并计算出下一预测值endfor i=1:nerror(i,:)=x(i,:)-var(i,:); %计算残差endc=std(error)/std(x);c%调用统计工具箱的标准差函数计算后验差的比值Lxago alpha var %显示输出预测值的累加数列%显示输出参数数列%显示输出预测值error %显示输出误差c %显示后验差的比值作用：数据处理判断是否可以用灰色预测、求级比、求累加数列、求a b 的值、求预测方程clc,clearx0=[71.1 72.4 72.4 72.1 71.4 72.0 71.6]'; %注意这里为列向量n=length(x0);lamda=x0(1:n-1)./x0(2:n) range=minmax(lamda')%计算级比%计算级比的范围x1=cumsum(x0) %累加运算B=[-0.5*(x1(1:n-1)+x1(2:n)),ones(n-1,1)];Y=x0(2:n);u=B\Y %拟合参数u(1)=a,u(2)=bx=dsolve( 'Dx+a*x=b' , 'x(0)=x0' ); %求微分方程的符号解x=subs(x,{ 'a' , 'b' , 'x0' },{u(1),u(2),x0(1)}) %代入估计参数值和初始值yuce1=subs(x, 't' ,[0:n-1]); %求已知数据的预测值y=vpa(x,6) %其中的6 表示显示6 位数字yuce=[x0(1),diff(yuce1)] %差分运算，还原数据。