1、建模

##建立灰色模型GM（1，1）對應(yīng)的函數(shù)

##x表示原始數(shù)據(jù)數(shù)列，k表示數(shù)據(jù)個數(shù)

gm11<-function(x,k)

{

n<-length(x)

x1<-numeric(n);

for(i in 1:n)   ##一次累加

{

x1[i]<-sum(x[1:i]);

}

b<-numeric(n)

m<-n-1

for(j in 1:m)

{

b[j+1]<-(0.5*x1[j+1]+0.5*x1[j])   ##緊鄰均值生成

}

Yn=t(t(x[2:n]))                   ##構(gòu)造Yn矩陣

B<-matrix(1,nrow=n-1,ncol=2)      

B[,1]<-t(t(-b[2:n]))              ##構(gòu)造B矩陣

A<-solve(t(B)%*%B)%*%t(B)%*%Yn;   ##使用最小二乘法求得灰參數(shù)a,u

a<-A[1];

u<-A[2];

x2<-numeric(k);

x2[1]<-x[1];

for(i in 1:k-1)

{

x2[1+i]=(x[1]-u/a)*exp(-a*i)+u/a;

}

x2=c(0,x2);

y=diff(x2);                     ##累減生成，獲得預(yù)測數(shù)據(jù)數(shù)列

y

}

##x1原始數(shù)據(jù)數(shù)列，x2是預(yù)測數(shù)據(jù)數(shù)列

x1<-x

x2<-gm11(x,length(x))

##檢驗?zāi)Ｐ途?/P>

acc<-function(x1,x2)

{

     n<-length(x1);

     sum1=0;

     for(k in 2:n-1)

     {

            sum1<-sum1+(x1[k]-x1[1]);

     }

     s1<-sum1+0.5*(x1[n]-x1[1]);

     sum2=0;

     for(k in 2:n-1)

     {

            sum2<-sum2+(x2[k]-x2[1]);

     }

     s2<-sum2+0.5*(x2[n]-x2[1]);

     abs1<-abs(s1)

     abs2<-abs(s2)

     abs12<-abs(s1-s2)

     ee<-(1+abs1+abs2)/(1+abs1+abs2+abs12)

     ee

}

2、應(yīng)用：預(yù)測2013年和2014年全國規(guī)模以上快遞企業(yè)收入

##x數(shù)列是2008年-2012年全國規(guī)模以上快遞企業(yè)收入數(shù)據(jù)（來源：國家郵政局，單位：億元）

> x<-c(408.40,479.00,574.60,758.00,1055.30)   

> gm11(x,7)

[1]  408.4000  443.1355  585.3243  773.1370 1021.2131 1348.8894 1781.7069

> x1<-x

> x2<-gm11(x,length(x))

> acc(x1,x2)

[1] 0.9851449

預(yù)測結(jié)果：預(yù)計2013、2014年全國規(guī)模以上快遞企業(yè)收入分別為1349億元、1782億元

灰色絕對關(guān)聯(lián)度為0.9851449，即關(guān)聯(lián)度為一級，預(yù)測精度為優(yōu)