【数学建模】预测模型笔记之灰色预测

学习来源：清风老师

灰色预测

灰色系统

GM(1,1)模型: Grey(Gray) Model

GM(1,1)是使用原始的离散非负数据列，通过一次累加生成削弱随机性的较有规律的新的离散数据列，然后通过建立微分方程模型，得到在离散点处的解经过累减生成的原始数据的近似估计值，从而预测原始数据的后续发展。
GM(1,1)模型，第一个‘1’表示微分方程是一阶的，后面的‘1’表示只有一个变量

GM(1,1)原理介绍

最小二乘法

OLS原理介绍

完全多重共线性问题再探究

一阶微分方程

一阶齐次线性微分方程

一阶非齐次线性微分方程

准指数规律的检验

发展系数与预测情形的探究

GM(1,1)模型的评价

GM(1,1)模型的拓展

什么时候用灰色预测？

1. 数据是以年份度量的非负数据（如果是月份或者季度数据一定要用我们上一讲学过的时间序列模型）
2. 数据能经过准指数规律的检验（除了前两期外，后面至少90%的期数的光滑比要低于0.5）
3. 数据的期数较短且和其他数据之间的关联性不强（小于等于10，也不能太短了，比如只有3期数据），要是数据期数较长，一般用传统的时间序列模型比较合适。

灰色预测的例题

预测的题目的一些小套路

1. 看到数据后先画时间序列图并简单的分析下趋势（例如：我们上一讲学过的时间序列分解）
2. 将数据分为训练组和试验组，尝试使用不同的模型对训练组进行建模，并利用试验组的数据判断哪种模型的预测效果最好（比如我们可以使用SSE这个指标来挑选模型，常见的模型有指数平滑、ARIMA、灰色预测、神经网络等）
3. 选择上一步骤中得到的预测误差最小的那个模型，并利用全部数据来重新建模，并对未来的数据进行预测
4. 画出预测后的数据和原来数据的时序图，看看预测的未来趋势是否合理

GM(1,1)模型代码

1. 画出原始数据的时间序列图，并判断原始数据中是否有负数或期数是否低于4期，如果是的话则报错，否则执行下一步
2. 对一次累加后的数据进行准指数规律检验，返回两个指标：
   指标1：光滑比小于0.5的数据占比（一般要大于60%）
   指标2：除去前两个时期外，光滑比小于0.5的数据占比（一般大于90%）
   并让用户决定数据是否满足准指数规律，满足则输入1，不满足则输入0
3. 如果上一步用户输入0，则程序停止；如果输入1，则继续下面的步骤。
4. 让用户输入需要预测的后续期数，并判断原始数据的期数：
   4.1 数据期数为4：
   分别计算出传统的GM(1,1)模型、新信息GM(1,1)模型和新陈代谢GM(1,1)模型对于未来期数的预测结果，为了保证结果的稳健性，对三个结果求平均值作为预测值。
   4.2 数据期数为5,6或7：
   取最后两期为试验组，前面的n-2期为训练组；用训练组的数据分别训练三种GM模型，并将训练出来的模型分别用于预测试验组的两期数据；利用试验组两期的真实数据和预测出来的两期数据，可分别计算出三个模型的SSE；选择SSE最小的模型作为我们建模的模型。
   4.3 数据期数大于7：
   取最后三期为试验组，其他的过程和4.2类似。
5. 输出并绘制图形显示预测结果，并进行残差检验和级比偏差检验。

main

%%  输入原始数据并做出时间序列图
clear;clc
year =[1995:1:2004]';  % 横坐标表示年份，写成列向量的形式（加'就表示转置）
x0 = [174,179,183,189,207,234,220.5,256,270,285]';  %原始数据序列，写成列向量的形式（加'就表示转置）
%测试数据
% year = [2009:2015]; % 其实本程序写成了行向量也可以，因为我怕你们真的这么写了，所以在后面会有判断。
% x0 = [730, 679, 632, 599, 589, 532, 511];
% year = [2010:2017]';   % 该数据很特殊，可以通过准指数规律检验，但是预测效果却很差
% x0 = [1.321,0.387,0.651,0.985,1.235,0.987,0.854,1.021]';
% year = [2014:2017]';
% x0 = [2.874,3.278,3.337,3.390]';

% 画出原始数据的时间序列图
figure(1); % 因为我们的图形不止一个，因此要设置编号
plot(year,x0,'o-'); grid on;  % 原式数据的时间序列图
set(gca,'xtick',year(1:1:end))  % 设置x轴横坐标的间隔为1
%改到自己论文，改标签
xlabel('年份');  ylabel('排污总量');  % 给坐标轴加上标签


%% 因为我们要使用GM(1,1)模型，其适用于数据期数较短的非负时间序列
ERROR = 0;  % 建立一个错误指标，一旦出错就指定为1
% 判断是否有负数元素
if sum(x0<0) > 0  % x0<0返回一个逻辑数组(0-1组成)，如果有数据小于0，则所在位置为1，如果原始数据均为非负数，那么这个逻辑数组中全为0，求和后也是0~
    disp('灰色预测的时间序列中不能有负数')
    ERROR = 1;
end

% 判断数据量是否太少
n = length(x0);  % 计算原始数据的长度
disp(strcat('原始数据的长度为',num2str(n)))    % strcat()是连接字符串的函数，第一讲学了，可别忘了哦
if n<=3
    disp('数据量太小')
    ERROR = 1;
end

% 数据太多时提示可考虑使用其他方法（不报错）
if n>10
    disp('这么多数据量，一定要考虑使用其他的方法，例如ARIMA，指数平滑等')
end

% 判断数据是否为列向量，如果输入的是行向量则转置为列向量
if size(x0,1) == 1
    x0 = x0';
end
if size(year,1) == 1
    year = year';
end


%% 对一次累加后的数据进行准指数规律的检验(注意，这个检验有时候即使能通过，也不一定能保证预测结果非常好，例如上面的第三组数据)
if ERROR == 0   % 如果上述错误均没有发生时，才能执行下面的操作步骤
    disp('------------------------------------------------------------')
    disp('准指数规律检验')
    x1 = cumsum(x0);   % 生成1-AGO序列，cumsum是累加函数哦~    注意：1.0e+03 *0.1740的意思是科学计数法,即10^3*0.1740 = 174
    rho = x0(2:end) ./ x1(1:end-1) ;   % 计算光滑度rho(k) = x0(k)/x1(k-1)
    
    % 画出光滑度的图形，并画上0.5的直线，表示临界值
    figure(2)
    plot(year(2:end),rho,'o-',[year(2),year(end)],[0.5,0.5],'-'); grid on;
    text(<