【Matlab数学建模】灰色预测模型

一、灰色预测的概念

  1982年我国学者邓聚龙教授发表第一篇中文论文《灰色控制系统》标志着灰色系统这一学科诞生。

  白色系统是指一个系统的内部特征已知的，即系统的信息是完全充分的。黑色系统是指一个系统的内部信息对外界来说是一无所知的，只能通过它与外界的联系来加以观测研究。

  灰色系统则介于二者之间。

二、灰色关联度与优势分析

  选取参考数列：

  假设有m个比较数列：

  则称：

  为比较数列X_i对参考数列X₀在k时刻的关联系数。其中为分辨系数，一般情况下

  关联系数在每一时刻都有一个，过于分散，不利于比较，为此我们给出：

  为比较数列X_i对参考数列X₀的关联度。

  然而关联系数公式中的并不能去别处因素关联是正关联还是负关联，可采用以下方法解决该问题，记：

  则当，则X_i和X_j为正关联，否则为负关联。

  代码实现：

function [zeta, r] = func(x, x_0, rho)
	
	[m, ~] = size(x);

	if nargin < 3
		rho = 0.5;
	end
	
	if nargin < 2
		x_0 = max(x);
	end
	
	x_abs = abs(x_0 - x);
	mn = min(min(x_abs));
	mx = max(max(x_abs));
	
	zeta = (mn + rho * mx) ./ (abs(repmat(x_0, m, 1) - x) + rho * mx);
	r = mean(zeta, 2);
	
end

三、灰色生成数列

  累加生成：

  累减生成：

function res = func(a, r)
	% a为原始数列
	% 若r为正数，则res为a的+r次累加数列
	% 若r为负数，则res为a的-r次累减数列

	if r == 0
		res = a;
		return
	end
	
	n = length(a);
	
	if r > 0
		a = func(a, r - 1);
		res = cumsum(a);
		return
	end
	
	if r < 0
		a = func(a, r + 1);
		res = a - [0, a(1 : length(a) - 1)];
		return
	end
	
end

  加权邻值生成

  称为原始数列的任意一对邻值。为后邻值，为前邻值。对于常数，设有：

  将数列称为数列在权下的邻值生成数，被称为生成系数。

  特别地，若，则称：

  为均值生成数，也称等权邻值生成数。

function z = func(x, alpha)
	
	if nargin < 2
		alpha = 0.5;
	end
	
	a = alpha * x;
	b = [0, (1 - alpha) * x(1 : length(x) - 1)];
	z = a + b;
	
end

四、灰色模型GM(1, 1)

① 灰色预测的相关概念

  G代表grey(灰色)，M代表model(模型)

  设原始数列：

  的一次累加数列:

  定义的灰导数为：

  令为数列的邻值生成数列，即：

  于是定义GM(1, 1)的灰微分方程模型为：

  其中称为灰导数，a称为发展系数，称为白化背景值，a称为灰作用量。

  取由灰微分方程可得：

  引入矩阵：



  于是GM(1, 1)模型可以表示为：


  对于GM(1, 1)的灰微分方程，如果将的时刻k换成连续变量t，则对应的白微分方程即是：

② 灰色预测的步骤

第一步：级比检验

  设原始数列为：

  计算数列的级比：

  如果所有的级比都落在可覆盖区间内，则原始数列可以建立GM(1, 1)模型并进行灰色预测，否则需要对数据进行变化处理，如平移变化：

第二步：建立GM(1, 1)模型

  原始数列为：

  为它建立GM(1, 1)模型：

  对应的白化模型是：

  该微分方程的解是：

  于是便能得到预测值：

function [u, f] = func(x, alpha)
	% x是原始数列
	% alpha是生成系数，默认为0.5
	% f是预测数列的一次生成数列的函数

	if nargin < 2
		alpha = 0.5;
	end
	
	x_1 = cumsum(x);
	z = alpha * x_1(2 : end) + (1 - alpha) * x_1(1 : end - 1);
	u = [-z; ones(1, length(x) - 1)]' \ x(2 : end)';
	f = @(k) (x(1) - u(2) /u(1)) * exp(-u(1) * k + u(1)) + u(2) / u(1);
	
end

第三步：检验预测值

  ① 残差检验：

  称：

  为相对残差

  如果，则认为达到较高的要求。否则，若，则认为达到一般要求。

  ② 级比偏差值检验

  计算：

  如果，则认为达到较高的要求。否则，若，则认为达到一般要求。