数学建模--灰色预测

灰色预测

一、简介

目前常用的一些预测方法(如回归分析等)，需要较大的样本，若样本较小，常造成较大误差，使预测目标失效。灰色预测模型对时间序列短、统计数据少、信息不完全系统的分析与建模，具有独特的功效，因此得到了广泛的应用，是处理小样本预测问题的有效工具。

简单介绍

灰色预测通过鉴别系统因素之间发展趋势的相异程度，即进行关联分析，并对原始数据进行灰色生成来寻找系统变动的规律，生成有较强规律性的数据序列，然后建立相应的微分方程模型，从而预测事物未来发展趋势的状况。

四种类型

(1)数列预测

用观察到的反映预测对象特征的时间序列来构造灰色预测模型，预测未来某一时刻的特征量，或达到某一特征量的时间。

(2)拓扑预测

原始数据作曲线，在曲线上按定值寻找该定值发生的所有时点，并以该定值为框架构成时点数列，然后建立模型预测该定值所发生的时点。

(3)系统预测

通过对系统行为特征指标建立一组相互关联的灰色预测模型，预测系统中众多变量间的相互协调关系的变化。

(4)灾变与异常值预测

通过灰色模型预测异常值出现的时刻，预测异常值什么时候出现在特定时区内。

今天主要介绍的是数列预测

二、灰色生成

定义

在灰色预测中数据预处理叫做灰色生成

级比生成

先判断数据是否有关联，级比生成就是用来判断是否有关联

不重要：将原始数据列中的数据，按某种要求作数据处理称为生成。客观世界尽管复杂，表述其行为的数据可能是杂乱无章的，然而它必然是有序的，都存在着某种内在规律，不过这些规律被纷繁复杂的现象所掩盖，人们很难直接从原始数据中找到某种内在的规律。对原始数据的生成就是企图从杂乱无章的现象中去发现内在规律

灰色生成方式

（1）累加生成

通过数列间各时刻数据的依次累加以得到新的数据与数列。累加前的数列称原始数列x^{(0)==，累加后的数列称为生成数列==x}(1)。

将非负的摆动数列转化为递增数列

如果有负数，所有值加一个数让最小值为0，再累加

（2）累减生成

数列相邻两数据的差，累减生成是累加生成的逆运算，累减生成可将累加生成还原为非生成数列。

（3）均值生成

将序列中的前后相邻的两数取平均数，以获得生成新的序列。

三、GM（1,1）模型

概念

对灰色系统建立的模型称为灰色模型。灰色模型就是利用离散随机数经过生成变为随机性被显著削弱而且较有规律的生成数，建立起的微分方程形式的模型，这样便于对其变化过程进行研究和描述。灰色预测模型又称为GM模型，揭示了系统内部事物连续发展变化的过程。

现在只介绍GM(1,1)模型

基本步骤

数据处理(1)(2),建立模型(3)(4)(5)(6)

(1)对原始白色信息和典型曲线拟合信息进行关联分析

数据不一定都在这个范围，如果不在说明模型不是特别好

(2)由原始数据序列x^{(0)计算一次累加序列x}(1)

(3)建立矩阵B,y

(4)求逆矩阵( B^(T) B)^(-1)

(5)根据U =(B^(T)B)(-1) B^(T) y，求参数估计值a和b

(6)得出响应函数，计算x^(1) (t) ，进行累减还原，得到x^(0) (t)的拟合值

(7)精度检验与预测

四、例题和作业

例题

clc;
clear;
%% 
%%建立符号变量a(发展系数)和b(灰作用量)
syms a b;
c=[a b]';  %转置矩阵，相当于原文的U

%原始数列A  
A=[2.874,3.278,3.337,3.39,3.679];
n=length(A);

%对原始数列A做累加得到数列D
D=cumsum(A);

%对数列D做均值生成
for i =2:n
    C(i)=(D(i)+D(i-1))/2;  %生成累加矩阵
end
C(1)=[];

%% 
%构造数据矩阵
B=[-C',ones(n-1,1)];
Y=A; Y(1)=[]; Y=Y';


%使用最小二乘法算参数a(发展系数)和b(灰作用量)
c=(inv(B'*B))*B'*Y;
c=c';
a=c(1);b=c(2);%带入到F(i) = (A(1)-b/a)/exp(a*(i-1))+b/a

%% 
%预测后续数据

F=[];F(1) = A(1);
for i=2:(n+1) %%改变数据
    F(i) = (A(1)-b/a)/exp(a*(i-1))+b/a;
end

%对数列F累减还原，得到预测出的数据
G=[];G(1)=A(1);
%预测往后一年内的数据，所以n+1
for i=2:(n+1)
    G(i)=F(i)-F(i-1); %得到预测出来的数据
end
disp('预测数据为：');
%% 
%模型检验

H=G(1:n);
%计算残差序列
epsilon=A - H;

%法一：相对残差检验法
%计算相对误差序列
delta=abs(epsilon./A);
%计算相对误差Q
disp('相对残差检验')
Q=mean(delta)

%法二：后残差检验法
%计算方法比C
disp('方差比C检验：')
C=std(epsilon,1)/std(A,1)

%计算小误差概率P
S1=std(A,1);
tmp=find(abs(epsilon-mean(epsilon))<0.6745*S1);
disp('小误差概率P检验：')
P=length(tmp)/n

%% 
%绘制曲线图
t1=2:6;%原始数据横轴数据个数
t2=2:7;%预测值横轴数据个数

plot(t1,A,'ro');hold on;
plot(t2,G,'g-');
xlabel('年份');ylabel('公司营业额');
legend('实际销售额','预测销售额');
title('企业销售额变化曲线');
grid on;

作业