function hsyc(x0)
%灰色模型GM(1,1)matlab代码在建立灰色模型前应对原始数据进行求%%%级比%%%%重要
%输入原始时序列x0,输出各种所需参数%%%%实际使用该程序时应改变图像坐标轴%%%%
%x0代表原始数据;x1代表x0经累加后;B代表所够造的矩阵B;Y代表向量Yn;au代表系数a,u构成的矩阵
n=size(x0,2); %计算x0的列数
x1=cumsum(x0,2); %将x0进行累加得到矩阵x1
for k=1:n-1;
b(k)=-0.5*(x1(k+1)+x1(k));
y(k)=x0(k+1); %计算得到向量Y的转置矩阵y
end
a=ones(n-1,1);
B=[b',a]; %计算得到矩阵B
Y=y'; %计算得到向量Y
au=inv(B'*B)*B'*Y; %计算得到系数a,u构成的矩阵au
A=au';
k=1:n+10;
X1(k)=(x0(1)-A(2)/A(1))*exp(-A(1)*(k-1))+A(2)/A(1); %计算预测累加数列的值
X0(1)=X1(1);
k=1:n+9;
X0(k+1)=X1(k+1)-X1(k); %计算预测累加数列的还原值,即预测值
for k=2:n;
e(k)=X0(k)-x0(k); %计算残差
E(k)=(X0(k)-x0(k))/x0(k); %计算预测值与实测值的差值跟实测值的比值,即残差与实测值的比值
j(k)=x0(k-1)/x0(k); %计算实测数据的级比
end
k=1:n;
p(k)=1-(1-0.5*A(1))*j(k)/(1+0.5*A(1)); %计算级比偏差e(k)、p(k)<0.1则说明模型达到较高要求,0.2>e(k)、p(k)>0.1则说明模型达到一般要求
m=min(abs(e)); %计算残差中的最小值
M=max(abs(e)); %计算残差中的最大值
k=1:n;
g(k)=(m+0.5*M)./(abs(e(k))+0.5*M); %计算关联系数
R=sum(g')/(n-1); %计算关联度
v=[1,5,120,185];
axis(v);
grid on
plot([1996:2005],x0,'o-',[1996:2015],X0,'*:'); %画出原始数据跟时间序列的图像及预测数据跟时间序列的图像
legend('原始数据','预测数据',4)
au %输出参数a,u的值