预测就是借助于对过去的探讨去推测、了解未来。灰色预测通过原始数据的处理和灰色模型的建立,发现、掌握系统发展规律,对系统的未来状态做出科学的定量预测。对于一个具体的问题,究竟选择什么样的预测模型应以充分的定性分析结论为依据。模型的选择不是一成不变的。一个模型要经过多种检验才能判定其是否合适,是否合格。只有通过检验的模型才能用来进行预测。本章将简要介绍灰数、灰色预测的概念,灰色预测模型的构造、检验、应用,最后对灾变预测的原理作了介绍。 灰色系统理论的产生和发展动态
灰色系统与模糊数学、黑箱方法的区别
2.1 R源码#灰色预测模型G(1,1)GM11<-function(x0,t){ #x0为输入学列,t为预测个数="">-function(x0,t){><-cumsum(x0) #一次累加生成序列1-ag0序列(累加序列)="">-cumsum(x0)><-numeric(length(x0)-1)# 初始化长度为length(x0)-1的整数部分个numeric类型且值为0的数据集b="">-numeric(length(x0)-1)#><-length(x0)-1# n="" 为length(x0)-1长度="" 因为需要生成mean(紧邻均值)生成序列="" 其长度少1="" for(i="" in="" 1:n){="" #生成x1的紧邻均值生成序列="" ="">-length(x0)-1#><--(x1[i]+x1[i+1]) ="" b}="" #得序列b,即为x1的紧邻均值生成序列="">--(x1[i]+x1[i+1])><-numeric(length(x0)-1)>-numeric(length(x0)-1)><-1>-1><-cbind(b,d)#作b矩阵>-cbind(b,d)#作b矩阵><-t(b)#b矩阵转置>-t(b)#b矩阵转置><-solve(bt%*%b)#求bt*b得逆>-solve(bt%*%b)#求bt*b得逆><-numeric(length(x0)-1)>-numeric(length(x0)-1)><-x0[2:length(x0)]>-x0[2:length(x0)]><-m%*%bt%*%yn #模型的最小二乘估计参数列满足alpha尖="">-m%*%bt%*%yn><-matrix(alpha,ncol=1)# 将结果变成一列="" #="" 得到方程的两个系数="">-matrix(alpha,ncol=1)#><-alpha2[1]>-alpha2[1]><-alpha2[2] ="" #="" 下面为结果输出="" #="" 输出参数估计值及模拟值="" cat('gm(1,1)参数估计值:','\n','发展系数-a=',-a,' ','灰色作用量u=',u,' \n','\n')="" #利用最小二乘法求得参数估计值a,u="">-alpha2[2]><-numeric(length(c(1:t)))# t为给定的预测个数="">-numeric(length(c(1:t)))#><-x1[1] #="" 第一个数不变="" for(w="" in="" 1:(t-1)){="" #将a,u的估计值代入时间响应序列函数计算x1拟合序列y="" ="">-x1[1]><-(x1[1]-u )*exp(-a*w)+u/a="" }="" cat('x(1)的模拟值:','\n',y,'\n')="">-(x1[1]-u><-numeric(length(y))>-numeric(length(y))><-y[1] for(o="" in="" 2:t){="" #运用后减运算还原得模型输入序列x0预测序列="" ="">-y[1]><-y[o]-y[o-1] }="" cat('x(0)的模拟值:','\n',xy,'\n','\n')="" ="" ="" ="" ="" ="" ="" ="" ="" ="" ="" ="" ="" #="" 计算残差e="">-y[o]-y[o-1]><-numeric(length(x0)) for(l="" in="" 1:length(x0)){="" ="">-numeric(length(x0))><-x0[l]-xy[l] #得残差序列(未取绝对值)="" }="" cat('绝对残差:','\n',e,'\n')="" #计算相对误差="">-x0[l]-xy[l]><-numeric(length(x0)) for(s="" in="" 1:length(x0)){="" ="">-numeric(length(x0))><-(abs(e[s]) 0[s])="" #得相对误差="" }="" cat('相对残差:','\n',e2,'\n','\n')="" cat('残差平方和=',sum(e^2),' \n')="" cat('平均相对误差=',sum(e2)/(length(e2)-1)*100,' %','\n')="" cat('相对精度=',(1-(sum(e2)/(length(e2)-1)))*100,' %','\n','\n')="" ="" #后验差比值检验="">-(abs(e[s])><><-sum((abs(e)-avge)^2);evar=esum length(e)-1);se="sqrt(evar)" #计算残差的均方差se="">-sum((abs(e)-avge)^2);evar=esum><><-sum((x0-avgx0)^2);x0var=x0sum length(x0));sx="sqrt(x0var)" #计算原序列x0的方差sx="">-sum((x0-avgx0)^2);x0var=x0sum><-se x="" #得验差比值(方差比)="" cat('后验差比值检验:','\n','c值=',cv,' \n')#对后验差比值进行检验,与一般标准进行比较判断预测结果好坏。="" #计算小残差概率="">-se><><0.6745*sx) ength(e)="" cat('小残差概率:','\n','p值=',P,' \n')="" ="" if(cv="">0.6745*sx)>< 0.35="" &&="" p="">0.95){ cat('C<0.35, p="">0.95,GM(1,1)预测精度等级为:好','\n','\n') }else{ if(cv<0.5 &&="" p="">0.80){ cat('C值属于[0.35,0.5), P>0.80,GM(1,1)模型预测精度等级为:合格','\n','\n') }else{ if(cv<0.65 &&="" p="">0.70){ cat('C值属于[0.5,0.65), P>0.70,GM(1,1)模型预测精度等级为:勉强合格','\n','\n') }else{ cat('C值>=0.65, GM(1,1)模型预测精度等级为:不合格','\n','\n') } } } #画出输入序列x0的预测序列及x0的比较图像 plot(xy,col='blue',type='b',pch=16,xlab='时间序列',ylab='值') points(x0,col='red',type='b',pch=4) legend('topleft',c('预测','原始'),title = '预测时序与原始时序对比',pch=c(16,4),lty=l,col=c('blue','red'))}0.65>0.5>0.35,> 2.2 GM11 测试> x<-c(2.67,3.13,3.25,3.36,3.56,3.72)> GM11(x,length(x)+2)GM(1,1)参数估计值: 发展系数-a= 0.04396098 灰色作用量u= 2.925617 x(1)的模拟值: 2.67 5.78087 9.031547 12.42832 15.97774 19.68668 23.56231 27.61211 x(0)的模拟值: 2.67 3.11087 3.250677 3.396768 3.549424 3.708941 3.875626 4.049803 绝对残差: 0 0.01913011 -0.0006772995 -0.03676788 0.010576 0.01105927 相对残差: 0 0.006111858 0.0002083999 0.01094282 0.002970786 0.002972921 残差平方和= 0.001952456 平均相对误差= 0.4641357 % 相对精度= 99.53586 % 后验差比值检验: C值= 0.0407599 小残差概率: P值= 1 C<0.35, p="">0.95,GM(1,1)预测精度等级为:好0.35,>-c(2.67,3.13,3.25,3.36,3.56,3.72)> |
|