gm.txt

关于灰色系统预测GM(1,1)模型在matlab6.5中的实现

clear
%原始数据输入；
x0=[7.3017;7.5189;7.4976;7.5230;7.5973;7.6479;7.9058;8.0884;8.3062];
n=9;
%由数列x0生成数列x1;
s=0;
for i=1:n
    s=s+x0(i);
    x1(i)=s;
end
%由数列x1生成矩阵C和A；
for j=1:n-1
    C(j,1)=x1(j+1);
    A(j,1)=x1(j);
    A(j,2)=1;
end
%计算出估计值a1,a和u;
beta=inv(A'*A)*A'*C;
beta1=beta(1);
beta2=beta(2);
%求a,u值；
a=-log(beta1);
u=beta2/(1-beta1)*a;
%求出数列x1的预测值数列x2;
m=input('m=');
for k=0:(m-1)
    x2(k+1)=(x0(1)-u/a)*exp(-a*k)+u/a;
end
%求出原始数列x0的预测值数列x3;
x3(1)=x0(1);
for k=1:(m-1)
    x3(k+1)=(1-exp(a))*(x0(1)-u/a)*exp(-a*k);
end
%残差数据输入；
d0=[-0.3437   -0.2226   -0.3409   -0.2997   -0.0000];
k0=3;
%由数列d0生成数列d1;
sum=0;
for i=1:(n-k0+1)
    sum=sum+d0(i);
    d1(i)=sum;
end
%由数列d1生成矩阵C1和A1；
for j=1:(n-k0)
    C1(j,1)=d1(j+1);
    A1(j,1)=d1(j);
    A1(j,2)=1;
end
%计算出估计值b,a1和u1;
b=inv(A1'*A1)*A1'*C1;
b1=b(1);
b2=b(2);
a1=-log(b1);
u1=b2/(1-b1)*a1;
%求出数列d1的预测值数列d2;
for k=0:(n-k0)
    d2(k+1)=(d0(1)-u1/a1)*exp(-a1*k)+u1/a1;
end
%求出原始数列d0的预测值数列d3;
d3(1)=d0(1);
for k=1:(m-1)
    d3(k+1)=(1-exp(a1))*(d0(1)-u1/a1)*exp(-a1*k);
end
%求出原始数列x0的修正预测值数列x4;
x4(1)=x0(1);
for k=1:(m-1)
    if   k>=k0
         x4(k+1)=x2(k+1)-a1.*d3(k-k0+1)+u1;
    else x4(k+1)=x2(k+1);
    end
end
%求出原始数列x0的预测值数列x3;
x3(1)=x0(1);
for k=1:m-1
    x3(k+1)=(1-exp(a))*(x0(1)-u/a)*exp(-a*k);
end
x3
%残差检验；
for k=1:n
d(k)=x1(k)-x4(k);
end
d
%相对误差；
for k=1:n
    e(k)=abs(d(k)/x0(k))*100;
end
e
s5=0;
for i=1:n
    s5=s5+d(i);   
end
p=1/n*s5
%后验检验；
s1=0;s3=0;
for k=1:n
s1=s1+x0(k);
end
for k=1:n
    s3=s3+d(k);
end
mean1=1/n*s1;
mean2=1/n*s3;
s2=0;s4=0;
for k=1:n
s2=s2+(x0(k)-mean1)^2;
end
for k=1:n
s4=s4+(d(k)-mean2)^2;
end
sigma1=sqrt(1/n*s2)
sigma2=sqrt(1/n*s4)
c=sigma2/sigma1
p1=0.6745*sigma1;
sum=0;
for k=1:n
    if abs(d(k)-mean2)<p1
        sum=sum+1;
    else sum; 
    end
end
p2=sum/n[/code]