資源描述：
《灰色模型Matlab程序》由會員上傳分享，免費在線閱讀，更多相關(guān)內(nèi)容在教育資源-天天文庫。
1、%byallen@紅嘴海鷗%灰色模型預測是在數(shù)據(jù)不呈現(xiàn)一定規(guī)律下可以采取的一種建模和預測方法，其預測數(shù)據(jù)與原始數(shù)據(jù)存在一定的規(guī)律相似性%下面程序是灰色模型GM(1,1)程序二次擬合和等維新陳代謝改進預測程序,matlab6.5,使用本程序請注明，程序存儲為gm1.m%x=[5999,5903,5848,5700,7884];gm1(x);測試數(shù)據(jù)%二次擬合預測GM(1,1)模型functiongmcal=gm1(x)sizexd2=size(x,2);%求數(shù)組長度k=0;fory1=xk=k+1;ifk>1x1(k)=x1(k-1)+x(k);%累加生
2、成z1(k-1)=-0.5*(x1(k)+x1(k-1));%z1維數(shù)減1，用于計算Byn1(k-1)=x(k);elsex1(k)=x(k);endend%x1,z1,k,yn1sizez1=size(z1,2);%size(yn1);z2=z1';z3=ones(1,sizez1)';YN=yn1';%轉(zhuǎn)置%YNB=[z2z3];au0=inv(B'*B)*B'*YN;au=au0';%B,au0,auafor=au(1);ufor=au(2);ua=au(2)./au(1);%afor,ufor,ua%輸出預測的au和u/a的值constant1
3、=x(1)-ua;afor1=-afor;x1t1='x1(t+1)';estr='exp';tstr='t';leftbra='(';rightbra=')';%constant1,afor1,x1t1,estr,tstr,leftbra,rightbrastrcat(x1t1,'=',num2str(constant1),estr,leftbra,num2str(afor1),tstr,rightbra,'+',leftbra,num2str(ua),rightbra)%輸出時間響應方程%******************************
4、************************%二次擬合k2=0;fory2=x1k2=k2+1;ifk2>kelseze1(k2)=exp(-(k2-1)*afor);endend%ze1sizeze1=size(ze1,2);z4=ones(1,sizeze1)';G=[ze1'z4];X1=x1';au20=inv(G'*G)*G'*X1;au2=au20';%z4,X1,G,au20Aval=au2(1);Bval=au2(2);%Aval,Bval%輸出預測的A,B的值strcat(x1t1,'=',num2str(Aval),estr,l
5、eftbra,num2str(afor1),tstr,rightbra,'+',leftbra,num2str(Bval),rightbra)%輸出時間響應方程nfinal=sizexd2-1+1;%決定預測的步驟數(shù)5這個步驟可以通過函數(shù)傳入%nfinal=sizexd2-1+1;%預測的步驟數(shù)1fork3=1:nfinalx3fcast(k3)=constant1*exp(afor1*k3)+ua;end%x3fcast%一次擬合累加值fork31=nfinal:-1:0ifk31>1x31fcast(k31+1)=x3fcast(k31)-x3fc
6、ast(k31-1);elseifk31>0x31fcast(k31+1)=x3fcast(k31)-x(1);elsex31fcast(k31+1)=x(1);endendendx31fcast%一次擬合預測值fork4=1:nfinalx4fcast(k4)=Aval*exp(afor1*k4)+Bval;end%x4fcastfork41=nfinal:-1:0ifk41>1x41fcast(k41+1)=x4fcast(k41)-x4fcast(k41-1);elseifk41>0x41fcast(k41+1)=x4fcast(k41)-x(1
7、);elsex41fcast(k41+1)=x(1);endendendx41fcast,x%二次擬合預測值%***精度檢驗pC************//////////////////////////////////k5=0;fory5=xk5=k5+1;ifk5>sizexd2elseerr1(k5)=x(k5)-x41fcast(k5);endend%err1%絕對誤差xavg=mean(x);%xavg%x平均值err1avg=mean(err1);%err1avg%err1平均值k5=0;s1total=0;fory5=xk5=k5+1;i
8、fk5>sizexd2elses1total=s1total+(x(k5)-xavg)^2;e
