CN103729501A - 基于灰色理论的短期电力负荷预测方法 - Google Patents

Abstract

本发明公开了一种基于灰色理论的短期电力负荷预测方法，包括以下步骤：辨识负荷数据中的缺陷数据，并对有缺损、突变的负荷数据补齐和修正；根据上述补齐和修正的负荷数据，对GM(1,1)模型的建模可行性进行验证，并利用对数处理方法对不合格的序列进行修正；对上述步骤二构建的GM(1，1)模型中的

参数进行修正；从不同角度选取数据序列构利用步骤三修正

参数后的GM(1，1)模型成不同预测方案进行预测，对预测日分段并计算各时段内各方案关联系数的均值，选择关联系数均值最大的方案为该时段预测方案；利用后验差校验法对GM(1,1)模型精度进行检验。达到预测使用普遍性，且准确度高的目的。

Description

基于灰色理论的短期电力负荷预测方法

技术领域

本发明涉及电力领域的数据处理，具体地，涉及一种基于灰色理论的短期电力负荷预测方法。

背景技术

在很难正确选择影响因素的前提下，要用正确的数学模型来描述电力负荷预测十分困难。建立在各种数学模型基础上的算法在国内的研究虽然非常广泛，但是不结合影响因素分析，不结合正确的数学模型，单靠改进算法是不能得到正确的负荷预测结果。短期负荷预测的突出特点是以日为周期呈现变化的相似性，同时明显受天气因素的影响和特别事件影响。且电力系统的负荷从本质上是不可控的，因此预测未来负荷变化最有效的方法就是观察负荷的历史纪录数据，针对实际情况和现有的资料找到适合实际系统的负荷预测方法。

短期负荷预测的三个关键环节是：①正确考虑影响因素、②选择正确的数学模型、③正确的预测算法。负荷预测受很多不确定因素的影响，例如：国民经济的发展情况、国家的政策、人口的增长、住宅区的发展趋势、家用电器的使用情况、气候环境等。很多预测方法对这些因素的考虑往往不够全面，甚至完全忽略。显然这样所得到的结果并不适应不确定因素的影响，与实际情况有较大差别。负荷预测的核心是负荷预测方法。国内外对于负荷预测方法的研究主要经历了以下几个阶段:

第一阶段

为负荷预测传统方法研究，如回归分析法、线性外推法、时间序列法、最小二乘法等，研究重点放在负荷序列本身的规律上。此类方法的优点是计算简单，要求的历史数据少，但计算精度差。负荷预测传统方法的研究比较早，所以无论是机理研究还是实际应用，此类方法都比较成熟。但是随着社会经济的快速发展以及各种新机制的引入，负荷预测传统方法的缺陷和不足越来越明显，己经不能满足电力系统负荷预测的要求。

第二阶段

为负荷预测现代方法，如灰色系统、神经网络、小波分析、模糊理论等，重点放在以智能化的新技术代替传统方法。负荷顶测的现代方法也称为不确定性方法，不确定性方法由于数学理论和算法的不断进步，其优越性得到充分显示，当然由于电力系统负荷预测的难度很大，各种预测的适用性也各不相同。

公开号为：CN102831488A，基于改进的灰色预测法的电力负荷预测方法，本发明的技术方案客服了传统GM灰色预测方法的缺点，有效处理气温突变对于负荷的影响，并预防了可能产生的负荷曲线偏移和畸变预测结果，预测精度较传统GM方法而言大幅提升。缺点是没有对GM(1,1)模型进行可行性检验，只是很盲目的应用GM(1,1)模型。但现有的预测模型存在不能适应不同地区的预测需求，预测结果不够完备和准确的缺陷。

发明内容

本发明的目的在于，针对上述问题，提出一种基于灰色理论的短期电力负荷预测方法，以实现具有普遍适用性，且准确度高的优点。

为实现上述目的，本发明采用的技术方案是：

一种基于灰色理论的短期电力负荷预测方法，包括以下步骤：

步骤1、辨识负荷数据中的缺陷数据，并对有缺损、突变的负荷数据补齐和修正；

步骤2、根据上述补齐和修正的负荷数据，对GM(1,1)模型的建模可行性进行验证，并利用对数处理方法对不合格的序列进行修正；

步骤3、对上述步骤二构建的GM(1，1)模型中的α参数进行修正；

步骤4、从不同角度选取数据序列构利用步骤三修正α参数后的GM(1，1)模型成不同预测方案进行预测，对预测日分段并计算各时段内各方案关联系数的均值，选择关联系数均值最大的方案为该时段预测方案；

步骤5、利用后验差校验法对GM(1,1)模型精度进行检验。

根据本发明的优选实施例，上述步骤1中的辨识负荷数据中的缺陷数据，并对有缺损、突变的负荷数据补齐和修正具体如下：

设某天的负荷为：L＝{l(1),l(2),...,l(96)}，对某天的96点负荷进行判断，若：|l(i)-l(i-1)|＜λ,(i＝2,3,...96λ为常数)，假如96点负荷相邻两点间的负荷全都落入λ阈内,将该日数据标为异常，并在预测取历史数据时直接予以剔除；

若有连续2个小时以上的数据点落入λ阈内，则按如下原则处理：

按照日期类型相似的原理寻找相似日，并剔除已被判断为异常的相似日，找到数据异常日之前的N-1个负荷正常数据源，并将待修正日的数据源标记为N，N为常数；已知第N个数据源待补足时段96-k个点以外的k个点的历史负荷数据，k为常数，这k个点数据正常，则：

$\left\{\begin{array}{c}{P}_{\mathrm{AVE}}\left(t\right)=\frac{1}{N-1}\underset{i-1}{\overset{N-1}{\mathrm{\Σ}}}{P}_R(i,t),(t=\mathrm{1,2},...,96)\\ {\stackrel{\‾}{P}}_k=\frac{1}{k}\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{P}_R(N,t),(t=\mathrm{1,2},...,k)\\ {\stackrel{\‾}{P}}_{\mathrm{AVEK}}=\frac{1}{k}\underset{t=1}{\overset{k}{\mathrm{\Σ}}}{P}_{\mathrm{AVE}}\left(t\right),(t=\mathrm{1,2},...k)\\ P_R(N,t)=P_{\mathrm{AVE}}\left(t\right)-({\stackrel{\‾}{P}}_{\mathrm{AVEK}}-{\stackrel{\‾}{P}}_k),(t=k+1,...96)\end{array}\right.$

式中P_R(i,t)是第i个数据源观测值；P_AVE(t)是N-1个数据源观测值平均；是第N个数据源k个正常点数据平均；P_AVEK是N-1个数据源在待修正日已知负荷时刻平均值。

根据本发明的优选实施例，上述步骤2中根据上述补齐和修正的负荷数据，对GM(1,1)模型的建模可行性进行验证，并利用对数处理方法对不合格的序列进行修正具体如下：

对于给定序列x⁽⁰⁾，可否建立精度较高的GM(1,1)，用x⁽⁰⁾的级比σ⁽⁰⁾(k)的大小与所属区间判断；

令x⁽⁰⁾为原始序列x⁽⁰⁾＝{x⁽⁰⁾(1),x⁽⁰⁾(2),x⁽⁰⁾(3),...,x⁽⁰⁾(n)}，对任意x⁽⁰⁾(k)∈x⁽⁰⁾k＝1,2,3,...,n，令σ⁽⁰⁾(k)为x⁽⁰⁾的级比，则

则当σ⁽⁰⁾(k)满足σ⁽⁰⁾(k)∈(0.1353,7.389)时，x⁽⁰⁾可作非畸形的GM(1,1)建模；

但建立精确度高的GM(1,1)模型，级比σ⁽⁰⁾(k)应落于靠近1的一个子区间(1-ε,1+ε)，即(1-ε,1+ε)

(0.1353,7.389)；

令x⁽⁰⁾为GM(1,1)建模序列x⁽⁰⁾＝{x⁽⁰⁾(1),x⁽⁰⁾(2),x⁽⁰⁾(3),...,x⁽⁰⁾(n)}，则x⁽⁰⁾的GM(1,1)模型发展系数a满足

建模序列x⁽⁰⁾的级比σ⁽⁰⁾(k)若满足

${\mathrm{\σ}}^{\left(0\right)}\left(k\right)\∈(e^{\frac{-2}{n+1}},e^{\frac{2}{n+1}})$

则x⁽⁰⁾是可以作GM(1,1)建模的序列数据。

根据本发明的优选实施例，对于上述级比σ⁽⁰⁾(k)检验不合格的序列，根据公式

σ_y(k)表示序列级比，Δ_y(k)为差异信息，y(k)变换后的数据序列元素，利用对数变换得到符合要求的序列y_m，

具体为：

令x为原始序列，y_m为x的m次对数序列，x表示为

x(x(1),x(2),...,x(n)),

y_m＝(y_m(1),y_m(2),...y_m(n)),

则，y_m(k)＝ln^mx(k)＝ln(ln(...(ln x(k))...))，

对m进行选择即可得到符合要求的序列y_m。

根据本发明的优选实施例，所述步骤3中的α参数修正包括以下步骤：

步骤301、设变量为x⁽⁰⁾的原始数据序列：x⁽⁰⁾＝[x⁽⁰⁾(1),x⁽⁰⁾(2),...x⁽⁰⁾(n)]用1-AGO累加生成，形成一阶累加生成序列：x⁽¹⁾＝[x⁽¹⁾(1),x⁽¹⁾(2),...x⁽¹⁾(n)]，其中： $x^{\left(1\right)}\left(k\right)=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{x}^{\left(0\right)}\left(i\right),k=\mathrm{1,2},...n;$

令z⁽¹⁾为x⁽¹⁾的均值序列：z⁽¹⁾＝(z⁽¹⁾(2),z⁽¹⁾(3),...z⁽¹⁾(n))其中：z⁽¹⁾(k+1)＝αx⁽¹⁾(k+1)+(1-α)x⁽¹⁾(k)第一次计算时，取α＝0.5，构造GM(1,1)的白化微分方程：

用最小二乘法求参数记为 $P=\left[\begin{array}{c}\hat{a}\\ \hat{u}\end{array}\right],$ 在最小二乘准则下Y_n＝BP的解为

$P=\left[\begin{array}{c}\hat{a}\\ \hat{u}\end{array}\right]={\left(B^{T}B\right)}^{-1}{B}^{T}Y$

上式为GM(1,1)参数

的矩阵辨识算式，式中(B^TB)^-1B^T是数据矩阵B的广义逆矩阵；其中： $B=\left[\begin{array}{cc}{-z}^{\left(1\right)}\left(2\right)& 1\\ {-z}^{\left(1\right)}\left(3\right)& 1\\ \·& \·\\ \·& \·\\ \·& \·\\ {-z}^{\left(1\right)}\left(n\right)& 1\end{array}\right],$ Y_n＝[x⁽⁰⁾(2),x⁽⁰⁾(3),...x⁽⁰⁾(n)]^T；

步骤302：

将上述步骤301中解得的参数

重新计算记为

将

与上一次计算所用的进行比较，给定任意小的正整数ε，如果

转入上述步骤301继续计算，将

代入z⁽¹⁾(k+1)＝αx⁽¹⁾(k+1)+(1-α)x⁽¹⁾(k)计算背景值z⁽¹⁾(k+1)，再次进行GM(1,1)建模和预测运算；

当

时，迭代结束，转入步骤303；

步骤303：

建立GM(1,1)的预测模型为： ${\hat{x}}^{\left(1\right)}(k+1)=(x^{\left(0\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})e^{-\hat{a}k}+\frac{\hat{u}}{\hat{a}},k=\mathrm{1,2},...n,$ 对x⁽¹⁾(k)进行累减还原，得x⁽⁰⁾(k)的预测值，即GM(1,1)预测模型为：

${\hat{x}}^{\left(0\right)}(k+1)={\hat{x}}^{\left(1\right)}(k+1)-{\hat{x}}^{\left(1\right)}\left(k\right)=(x^{\left(0\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})(1-e^{\hat{a}})e^{-\hat{a}k},k=\mathrm{1,2},...n$

输出运行结果。

本发明的技术方案具有以下有益效果：

本发明的技术方案，通过GM(1,1)模型以及α参数修正，计算速度快、迭代次数少，迭代一次既可满足精度要求，解决负荷增长率较大即α参数较大时的预测精度差的问题，达到预测使用普遍性，且准确度高的目的。

下面通过附图和实施例，对本发明的技术方案做进一步的详细描述。

附图说明

图1为本发明实施例所述的基于灰色理论的短期电力负荷预测方法的流程图；

图2为本发明实施例所述的α参数修正算法流程图；

图3为本发明实施例所述的基于普通GM(1,1)预测模型的方案一预测示意图；

图4为本发明实施例所述的基于普通GM(1,1)预测模型的方案二预测示意图；

图5为本发明实施例所述的基于普通GM(1,1)预测模型的方案三预测示意图；

图6为本发明实施例所述的基于灰关联分段优选组合GM(1,1)模型的方案预测示意图。

具体实施方式

以下结合附图对本发明的优选实施例进行说明，应当理解，此处所描述的优选实施例仅用于说明和解释本发明，并不用于限定本发明。

如图1所示，一种基于灰色理论的短期电力负荷预测方法，包括以下步骤：

步骤1、辨识负荷数据中的缺陷数据，并对有缺损、突变的负荷数据补齐和修正；

步骤2、根据上述补齐和修正的负荷数据，对GM(1,1)模型的建模可行性进行验证，并利用对数处理方法对不合格的序列进行修正；

步骤3、对上述步骤二构建的GM(1，1)模型中的α参数进行修正；

步骤4、从不同角度选取数据序列构利用步骤三修正α参数后的GM(1，1)模型成不同预测方案进行预测，对预测日分段并计算各时段内各方案关联系数的均值，选择关联系数均值最大的方案为该时段预测方案；

步骤5、利用后验差校验法对GM(1,1)模型精度进行检验。

其中，步骤1中的辨识负荷数据中的缺陷数据，并对有缺损、突变的负荷数据补齐和修正具体如下：

设某天的负荷为：L＝{l(1),l(2),...,l(96)}，对某天的96点负荷进行判断，若：|l(i)-l(i-1)|＜λ,(i＝2,3,...96λ为常数)，假如96点负荷相邻两点间的负荷全都落入λ阈内,将该日数据标为异常，并在预测取历史数据时直接予以剔除；

若有连续2个小时以上的数据点落入λ阈内，则按如下原则处理：

按照日期类型相似的原理寻找相似日，并剔除已被判断为异常的相似日，找到数据异常日之前的N-1个负荷正常数据源，并将待修正日的数据源标记为N，N为常数；已知第N个数据源待补足时段96-k个点以外的k个点的历史负荷数据，k为常数，这k个点数据正常，则：

$\left\{\begin{array}{c}{P}_{\mathrm{AVE}}\left(t\right)=\frac{1}{N-1}\underset{i-1}{\overset{N-1}{\mathrm{\Σ}}}{P}_R(i,t),(t=\mathrm{1,2},...,96)\\ {\stackrel{\‾}{P}}_k=\frac{1}{k}\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{P}_R(N,t),(t=\mathrm{1,2},...,k)\\ {\stackrel{\‾}{P}}_{\mathrm{AVEK}}=\frac{1}{k}\underset{t=1}{\overset{k}{\mathrm{\Σ}}}{P}_{\mathrm{AVE}}\left(t\right),(t=\mathrm{1,2},...k)\\ P_R(N,t)=P_{\mathrm{AVE}}\left(t\right)-({\stackrel{\‾}{P}}_{\mathrm{AVEK}}-{\stackrel{\‾}{P}}_k),(t=k+1,...96)\end{array}\right.$

式中P_R(i,t)是第i个数据源观测值；P_AVE(t)是N-1个数据源观测值平均；

是第N个数据源k个正常点数据平均；P_AVEK是N-1个数据源在待修正日已知负荷时刻平均值。

步骤2中根据上述补齐和修正的负荷数据，对GM(1,1)模型的建模可行性进行验证，并利用对数处理方法对不合格的序列进行修正具体如下：

对于给定序列x⁽⁰⁾，可否建立精度较高的GM(1,1)，用x⁽⁰⁾的级比σ⁽⁰⁾(k)的大小与所属区间判断；

令x⁽⁰⁾为原始序列x⁽⁰⁾＝{x⁽⁰⁾(1),x⁽⁰⁾(2),x⁽⁰⁾(3),...,x⁽⁰⁾(n)}，对任意x⁽⁰⁾(k)∈x⁽⁰⁾k＝1,2,3,...,n，令σ⁽⁰⁾(k)为x⁽⁰⁾的级比，则则当σ⁽⁰⁾(k)满足σ⁽⁰⁾(k)∈(0.1353,7.389)时，x⁽⁰⁾可作非畸形的GM(1,1)建模；

但建立精确度高的GM(1,1)模型，级比σ⁽⁰⁾(k)应落于靠近1的一个子区间(1-ε,1+ε)，即(1-ε,1+ε)

(0.1353,7.389)；

令x⁽⁰⁾为GM(1,1)建模序列x⁽⁰⁾＝{x⁽⁰⁾(1),x⁽⁰⁾(2),x⁽⁰⁾(3),...,x⁽⁰⁾(n)}，则x⁽⁰⁾的GM(1,1)模型发展系数a满足

建模序列x⁽⁰⁾的级比σ⁽⁰⁾(k)若满足

${\mathrm{\σ}}^{\left(0\right)}\left(k\right)\∈(e^{\frac{-2}{n+1}},e^{\frac{2}{n+1}})$

则x⁽⁰⁾是可以作GM(1,1)建模的序列数据。

上述级比σ⁽⁰⁾(k)检验不合格的序列，根据公式

σ_y(k)表示序列级比，Δ_y(k)为差异信息，y(k)变换后的数据序列元素，利用对数变换得到符合要求的序列y_m，

具体为：

令x为原始序列，y_m为x的m次对数序列，x表示为

x(x(1),x(2),...,x(n)),

y_m＝(y_m(1),y_m(2),...y_m(n)),

则，y_m(k)＝ln^mx(k)＝ln(ln(...(ln x(k))...))，

对m进行选择即可得到符合要求的序列y_m。

步骤3中的α参数修正包括以下步骤：

步骤301、设变量为x⁽⁰⁾的原始数据序列：x⁽⁰⁾＝[x⁽⁰⁾(1),x⁽⁰⁾(2),...x⁽⁰⁾(n)]用1-AGO累加生成，形成一阶累加生成序列：x⁽¹⁾＝[x⁽¹⁾(1),x⁽¹⁾(2),...x⁽¹⁾(n)]，其中： $x^{\left(1\right)}\left(k\right)=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{x}^{\left(0\right)}\left(i\right),k=\mathrm{1,2},...n;$

令z⁽¹⁾为x⁽¹⁾的均值序列：z⁽¹⁾＝(z⁽¹⁾(2),z⁽¹⁾(3),...z⁽¹⁾(n))其中：z⁽¹⁾(k+1)＝αx⁽¹⁾(k+1)+(1-α)x⁽¹⁾(k)第一次计算时，取α＝0.5，构造GM(1,1)的白化微分方程：

用最小二乘法求参数

记为 $P=\left[\begin{array}{c}\hat{a}\\ \hat{u}\end{array}\right],$ 在最小二乘准则下Y_n＝BP的解为

$P=\left[\begin{array}{c}\hat{a}\\ \hat{u}\end{array}\right]={\left(B^{T}B\right)}^{-1}{B}^{T}Y$

上式为GM(1,1)参数

的矩阵辨识算式，式中(B^TB)^-1B^T是数据矩阵B的广义逆矩阵；其中： $B=\left[\begin{array}{cc}{-z}^{\left(1\right)}\left(2\right)& 1\\ {-z}^{\left(1\right)}\left(3\right)& 1\\ \·& \·\\ \·& \·\\ \·& \·\\ {-z}^{\left(1\right)}\left(n\right)& 1\end{array}\right],$ Y_n＝[x⁽⁰⁾(2),x⁽⁰⁾(3),...x⁽⁰⁾(n)]^T；

步骤302：

将上述步骤301中解得的参数重新计算记为

将

与上一次计算所用的进行比较，给定任意小的正整数ε，如果

转入上述步骤301继续计算，将

代入z⁽¹⁾(k+1)＝αx⁽¹⁾(k+1)+(1-α)x⁽¹⁾(k)计算背景值z⁽¹⁾(k+1)，再次进行GM(1,1)建模和预测运算；

当时，迭代结束，转入步骤303；

步骤303：

建立GM(1,1)的预测模型为： ${\hat{x}}^{\left(1\right)}(k+1)=(x^{\left(0\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})e^{-\hat{a}k}+\frac{\hat{u}}{\hat{a}},k=\mathrm{1,2},...n,$ 对x⁽¹⁾(k)进行累减还原，得x⁽⁰⁾(k)的预测值，即GM(1,1)预测模型为：

${\hat{x}}^{\left(0\right)}(k+1)={\hat{x}}^{\left(1\right)}(k+1)-{\hat{x}}^{\left(1\right)}\left(k\right)=(x^{\left(0\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})(1-e^{\hat{a}})e^{-\hat{a}k},k=\mathrm{1,2},...n$

输出运行结果。

1、灰色建模过程

GM(1,1)模型是最常用的一种灰色模型，它是一个只包含单变量的一阶微分方程构成的模型，在电力负荷预测中应用极为广泛。GM(1,1)模型的实质是对原始数据序列做一次累加生成，使生成的数列呈一定规律，其相应趋势可以用典型曲线逼近，然后用逼近的曲线作为模型，用以对系统进行预测，建立GM(1,1)模型只需一个数据列{x⁽⁰⁾}。

设有变量为x(0)的原始数据序列

x⁽⁰⁾＝[x⁽⁰⁾(1),x⁽⁰⁾(2),...,x⁽⁰⁾(n)]   (1)

用1-AGO累加生成，形成一阶累加生成序列

x⁽¹⁾＝[x⁽¹⁾(1),x⁽¹⁾(2),...x⁽¹⁾(n)]

其中 $x^{\left(1\right)}\left(k\right)=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{x}^{\left(0\right)}\left(i\right)---\left(2\right)$

由于序列x⁽¹⁾(k)具有指数增长规律，而一阶微分方程的解恰好是指数增长形式，因此可以认为x⁽¹⁾序列满足下述一阶微分方程模型，也称为白化微分方程，建立相应的微分方程为

$\frac{{\mathrm{dx}}^{\left(1\right)}}{\mathrm{dt}}+{\mathrm{ax}}^{\left(1\right)}=u---\left(3\right)$

根据导数的定义，有

$\frac{\mathrm{dx}}{\mathrm{dt}}={\lim }_{\mathrm{\Δt}\→0}\frac{x^{\left(1\right)}(t+\mathrm{\Δt})-x^{\left(1\right)}\left(t\right)}{\mathrm{\Δt}}.$

若以离散形式表示，微分项可以写成

$\frac{\mathrm{\Δx}}{\mathrm{\Δt}}=\frac{x^{\left(1\right)}(k+1)-x^{\left(1\right)}k}{k+1-k}=x^{\left(1\right)}(k+1)-x^{\left(1\right)}\left(k\right)={\mathrm{\α}}^{\left(1\right)}\[x^{\left(1\right)}(k+1)\]$

其中x⁽¹⁾(k)的值只能取k和k+1的平均值，即

因此，式(3）可改写成 ${\mathrm{\α}}^{\left(1\right)}\[x^{\left(1\right)}(k+1)\]=\frac{1}{2}\mathrm{\α}\[x^{\left(1\right)}(k+1)+x^{\left(1\right)}\left(k\right)\]=u.$

为此可推出

$\begin{array}{c}k=1,x^{\left(0\right)}\left(2\right)+\frac{1}{2}a\[x^{\left(1\right)}\left(1\right)+x^{\left(1\right)}\left(2\right)\]=u\\ k=2,x^{\left(0\right)}\left(3\right)+\frac{1}{2}a\[x^{\left(1\right)}\left(2\right)+x^{\left(1\right)}\left(3\right)\]=u\\ \·\\ \·\\ \·\\ k=n-1,x^{\left(0\right)}\left(n\right)+\frac{1}{2}a\[x^{\left(1\right)}\left(n\right)+x^{\left(1\right)}(n+1)\]=u\end{array}---\left(4\right)$

将上述结果(4)写成矩阵形式有

$\left[\begin{array}{c}{x}^{\left(0\right)}\left(2\right)\\ x^{\left(0\right)}\left(3\right)\\ \·\\ \·\\ \·\\ x^{\left(0\right)}\left(n\right)\end{array}\right]=\left[\begin{array}{cc}-\frac{1}{2}\[x^{\left(1\right)}\left(1\right)+x^{\left(1\right)}\left(2\right)\]& 1\\ -\frac{1}{2}\[x^{\left(1\right)}\left(2\right)+x^{\left(1\right)}\left(3\right)\]& 1\\ \·& \·\\ \·& \·\\ \·& \·\\ -\frac{1}{2}\[x^{\left(1\right)}(n-1)+x^{\left(1\right)}\left(n\right)\]& 1\end{array}\right]\left[\begin{array}{c}a\\ u\end{array}\right]---\left(5\right)$

简记为Y_n＝BA, $Y_n=\left[\begin{array}{c}{x}^{\left(0\right)}\left(2\right)\\ x^{\left(0\right)}\left(3\right)\\ \·\\ \·\\ \·\\ x^{\left(0\right)}\left(n\right)\end{array}\right]$ $A=\left[\begin{array}{c}a\\ u\end{array}\right]$

$B=\left[\begin{array}{cc}-\frac{1}{2}\[x^{\left(1\right)}\left(1\right)+x^{\left(1\right)}\left(2\right)\]& 1\\ -\frac{1}{2}\[x^{\left(1\right)}\left(2\right)+x^{\left(1\right)}\left(3\right)\]& 1\\ \·& \·\\ \·& \·\\ \·& \·\\ -\frac{1}{2}\[x^{\left(1\right)}(n-1)+x^{\left(1\right)}\left(n\right)\]& 1\end{array}\right]---\left(6\right)$

上述方程组中，Y_n和B为已知量，A为待定参数。由于变量只有a和u两个，而方程个数却有(n-1)个，而(n-1)＞2，故方程组无解。但可用最小二乘法得到最小二乘近似解。因此上式可改写为:其中E为误差项

欲使 $\min {\left|\right|Y_n-B\hat{A}\left|\right|}^2=\min {(Y_n-B\hat{A})}^T(Y_n-B\hat{A}),$ 利用矩阵求导公式，可得 $A={\left(B^{T}B\right)}^{-1}{B}^T{Y}_n=\left[\begin{array}{c}\hat{a}\\ \hat{u}\end{array}\right]$ 将所求得的

代回原来的微分方程，有 $\frac{{\mathrm{dx}}^{\left(1\right)}}{\mathrm{dt}}+\hat{a}{x}^{\left(1\right)}=\hat{u}$ 解之可得

$x^{\left(1\right)}\left(t\right)=\[x^{\left(1\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}}\]e^{-\hat{a}t}+\frac{\hat{u}}{\hat{a}}.$ 写成离散形式(令x⁽¹⁾(1)＝x⁽⁰⁾(1))，得下式：

$x^{\left(1\right)}(k+1)=\[x^{\left(1\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})e^{-\hat{a}k}+\frac{\hat{u}}{\hat{a}},(k=\mathrm{0,1,2},...)---\left(7\right)$

式(7)称为GM(1,1)模型得时间响应函数模型，它是GM(1,1)模型灰色预测的具体计算公式，对此式再做累减还原，得原始数列x⁽⁰⁾的灰色预测模型为

${\hat{x}}^{\left(0\right)}(k+1)={\hat{x}}^{\left(1\right)}(k+1)-{\hat{x}}^{\left(1\right)}\left(k\right)=(e^{-a}-1)(x^{\left(0\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})e^{-\hat{a}k},(k=\mathrm{0,1,2},...)---\left(8\right)$

灰色模型的建模优劣精度通常用后验差检验方法进行分析。

数据预处理

在进行负荷预测之前，首先要对这些异常数据进行处理。本技术方案采用了如下的数据预处理方法：

数据缺损处理：

设某天的负荷为：L＝{l(1),l(2),...,l(96)}，对某天的96点负荷进行判断，若：|l(i)-l(i-1)|＜λ,(i＝2,3,...96λ为常数)，假如96点负荷相邻两点间的负荷全都落入λ阈内,将该日数据标为异常，并在预测取历史数据时直接予以剔除。若有连续2个小时以上的数据点落入λ阈内，则按如下原则处理：

按照日期类型相似的原理寻找相似日，并剔除已被判断为异常的相似日，找到数据异常日之前的N-1个负荷正常数据源，并将待修正日的数据源标记为N（N为常数）。现已知第N个数据源待补足时段96-k个点以外的k个点（这k个点数据正常）的历史负荷数据，则：

$\left\{\begin{array}{c}{P}_{\mathrm{AVE}}\left(t\right)=\frac{1}{N-1}\underset{i-1}{\overset{N-1}{\mathrm{\Σ}}}{P}_R(i,t),(t=\mathrm{1,2},...,96)\\ {\stackrel{\‾}{P}}_k=\frac{1}{k}\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{P}_R(N,t),(t=\mathrm{1,2},...,k)\\ {\stackrel{\‾}{P}}_{\mathrm{AVEK}}=\frac{1}{k}\underset{t=1}{\overset{k}{\mathrm{\Σ}}}{P}_{\mathrm{AVE}}\left(t\right),(t=\mathrm{1,2},...k)\\ P_R(N,t)=P_{\mathrm{AVE}}\left(t\right)-({\stackrel{\‾}{P}}_{\mathrm{AVEK}}-{\stackrel{\‾}{P}}_k),(t=k+1,...96)\end{array}\right.---\left(9\right)$

上式中P_R(i,t)是第i个数据源观测值；P_AVE(t)是N-1个数据源观测值平均；

是第N个数据源k个正常点数据平均；P_AVEK是N-1个数据源在待修正日已知负荷时刻平均值.

GM(1,1)建模可行性检验及修正

GM(1,1)建模可行性分析

对于给定序列x⁽⁰⁾，可否建立精度较高的GM(1,1)，一般用x⁽⁰⁾的级比σ⁽⁰⁾(k)的大小与所属区间来判断。

令x⁽⁰⁾为原始序列x⁽⁰⁾＝{x⁽⁰⁾(1),x⁽⁰⁾(2),x⁽⁰⁾(3),...,x⁽⁰⁾(n)}，对任意x⁽⁰⁾(k)∈x⁽⁰⁾k＝1,2,3,...,n，令σ⁽⁰⁾(k)为x⁽⁰⁾的级比，则

则当σ⁽⁰⁾(k)满足σ⁽⁰⁾(k)∈(0.1353,7.389)时，x⁽⁰⁾可作非畸形的GM(1,1)建模。

因此只根据原始序列x⁽⁰⁾的级比σ⁽⁰⁾(k)的大小，就可以判断GM(1,1)建模的可行性，这就是灰建模GM(1,1)建模可行性的级比判断。

级比可容区（0.1353,7.389）是GM(1,1)建模的基本条件，然而不是实用条件。也就是说要想建立满意有效的GM(1,1)模型，级比σ⁽⁰⁾(k)应落于靠近1的一个子区间(1-ε,1+ε)，即(1-ε,1+ε)

(0.1353,7.389)。

令x⁽⁰⁾为GM(1,1)建模序列x⁽⁰⁾＝{x⁽⁰⁾(1),x⁽⁰⁾(2),x⁽⁰⁾(3),...,x⁽⁰⁾(n)}，则x⁽⁰⁾的GM(1,1)模型发展系数a满足

可行性分析的准则是：建模序列x⁽⁰⁾的级比σ⁽⁰⁾(k)若满足

${\mathrm{\σ}}^{\left(0\right)}\left(k\right)\∈(e^{\frac{-2}{n+1}},e^{\frac{2}{n+1}})$

则认为x⁽⁰⁾是可作GM(1,1)建模的。

数据处理原则与机理

数据处理原则

对于级比检验不合格的序列，必须做数据变换处理，使其变换后的序列，其级比落于可容区中。灰建模序列x的级比σ(k)必须落在可行域ItG中，ItG=（0.1353，7.389），才能作GM(1,1)建模，而为了获得精度高的GM(1,1)模型，级比σ(k)被限制在ItG中靠近1的子区间ItGM中，即ItGM

ItG,ItGM＝(1-ε,1+ε)，ε是指定的足够小的实数，因此灰建模数据处理的原则是：经过处理后的序列级比σ_y(k)应尽量靠近1，也就是σ_y(k)应尽可能小。

数据处理机理

在数据处理原则中已指出，数据处理原则是尽量减少级比偏差σ_y(k)，但由于

因此数据处理的机理是：选择合适的处理序列y，使差异信息Δ_y(k)与变换数据y(k)。通常变换处理途径有：对数变换、方根变换、平移变换，本技术方案使用的是对数变换，具体来说就是，在对数处理中，通过选取合适的对数阶次来达到。

对数处理方法

令x为原始序列，y_m为x的m次对数序列

x(x(1),x(2),...,x(n)),

y_m＝(y_m(1),y_m(2),...y_m(n)),

y_m(k)＝ln^mx(k)＝ln(ln(...(ln x(k))...))   (10)

(例如，m＝2,y_m(k)＝y₂(k)＝ln(ln x(k)))

又记Δ(k)＝|x(k)-x(k-1)|,Δ_m(k)＝\ln^mx(k),ln^mx(k-1)|，则有： $\mathrm{\δ}\left(k\right)=\frac{\mathrm{\Δ}\left(k\right)}{x\left(k\right)},$ ${\mathrm{\δ}}_m\left(k\right)=\frac{{\mathrm{\Δ}}_m\left(k\right)}{{\ln }^{m}x\left(k\right)}.$

对于指定的小于1的正实数ε，可通过合适的m使得δ_m(k)＜ε，此时，y_m为对数变换的满意序列，但必须指出：并不是m越大越好，因为过分大的m，将导致ln^mx(k)过分减小，这样，基于关系

在过分小的ln^mx(k)下，有可能导致δ_m(k)增大。

采用α参数修正提高灰色预测精度

通过对大量数据列进行预测发现，当数据列的发展速率比较快（即|a|的值较大）时，用GM(1,1)模型进行预测，预测精度较差，针对此种情况，本文对该方法进行了分析，提出了α参数修正法。这种方法是用：z⁽¹⁾(k+1)＝αx⁽¹⁾(k)+(1-α)x⁽¹⁾(k+1)来取代传统的

减少因|α|较大带来的预测误差，从而提高预测精度，并扩大预测范围。

有关文献指出背景值z⁽¹⁾(k)的精确计算公式应当是：设X⁽⁰⁾＝{x⁽⁰⁾(1),x⁽⁰⁾(2),...x⁽⁰⁾(n)}是一个非负序列，X⁽¹⁾是一次累加序列，给定z⁽¹⁾(k+1)＝αx⁽¹⁾(k)+(1-α)x⁽¹⁾(k+1)，则a与α具有以下关系：

利用罗必塔法则可以证明，a→0时，α的极限值为0.5。当|a|较小时，α非常接近0.5，当|a|较大时，α偏离0.5较大。然而目前GM(1,1)建模传统计算方法背景值生成采用：，只是简单地取α＝0.5不考虑a这种盲目采用α＝0.5的简单计算方法是导致在|a|较大时预测失效的重要原因。因此，可以根据不同的a值，来选择不同的α进行背景值z⁽¹⁾(k+1)的计算，从而有可能解决在|a|较大时的预测精度问题。

利用以上结论，本文提出了GM(1,1)模型的α参数修正预测方法。算法如下：

1：GM(1,1)的建模和预测运算。设有变量为x⁽⁰⁾的原始数据序列：x⁽⁰⁾＝[x⁽⁰⁾(1),x⁽⁰⁾(2),...x⁽⁰⁾(n)]用1-AGO累加生成，形成一阶累加生成序列：x⁽¹⁾＝[x⁽¹⁾(1),x⁽¹⁾(2),...x⁽¹⁾(n)]，其中：

令z⁽¹⁾为x⁽¹⁾的均值序列：z⁽¹⁾＝(z⁽¹⁾(2),z⁽¹⁾(3),...z⁽¹⁾(n))其中：z⁽¹⁾(k+1)＝αx⁽¹⁾(k+1)+(1-α)x⁽¹⁾(k)第一次计算时，可取α＝0.5。构造GM(1,1)的白化微分方程：用最小二乘法求参数记为 $P=\left[\begin{array}{c}\hat{a}\\ \hat{u}\end{array}\right],$ 在最小二乘准则下Y_n＝BP的解为

$P=\left[\begin{array}{c}\hat{a}\\ \hat{u}\end{array}\right]={\left(B^{T}B\right)}^{-1}{B}^{T}Y---\left(11\right)$

式(11)为GM(1,1)参数

的矩阵辨识算式。式中(B^TB)^-1B^T事实上是数据矩阵B的广义逆矩阵。其中： $B=\left[\begin{array}{cc}{-z}^{\left(1\right)}\left(2\right)& 1\\ {-z}^{\left(1\right)}\left(3\right)& 1\\ \·& \·\\ \·& \·\\ \·& \·\\ {-z}^{\left(1\right)}\left(n\right)& 1\end{array}\right],$ Y_n＝[x⁽⁰⁾(2),x⁽⁰⁾(3),...x⁽⁰⁾(n)]^T

2：将解得的参数

代入

中，重新计算

记为将

与上一次计算所用的

进行比较。给定任意小的正整数ε，如果

表明还可能大幅度提高预测精度，应转入1中，将

代入z⁽¹⁾(k+1)＝αx⁽¹⁾(k+1)+(1-α)x⁽¹⁾(k)计算背景值z⁽¹⁾(k+1)，再次进行GM(1,1)建模和预测运算。

当

时，迭代结束，转入下面的3。

3：建立GM(1,1)的预测模型为： ${\hat{x}}^{\left(1\right)}(k+1)=(x^{\left(0\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})e^{-\hat{a}k}+\frac{\hat{u}}{\hat{a}},k=\mathrm{1,2},...n,$ 对x⁽¹⁾(k)进行累减还原，得x⁽⁰⁾(k)的预测值，即GM(1,1)预测模型为：

${\hat{x}}^{\left(0\right)}(k+1)={\hat{x}}^{\left(1\right)}(k+1)-{\hat{x}}^{\left(1\right)}\left(k\right)=(x^{\left(0\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})(1-e^{\hat{a}})e^{-\hat{a}k},k=\mathrm{1,2},...n---\left(12\right)$

最后输出运行结果。

该方法原理简单、计算速度快、迭代次数少，一般迭代一次既可满足精度要求，可基本解决负荷增长率较大即|a|较大时的预测精度差的问题。算法流程见图2：

模型之间的比较

为提高灰色预测的精度和适用范围，本技术方案结合工程应用提出了一种内外部优化模式相结合的方法：①提出基于关联系数分析的分段多方案优选组合预测算法，这是模型外部环境优化法；②增加建模可行性验证、以及为提高预测精度而进行的α参数修正二个关键处理算法，这是模型内部优化法；③为了保证预测用的原始数据都是真实可靠的，本技术方案还增加了负荷数据预处理算法，整个模型的基本流程如下：

辨识缺陷数据，并对有缺损、突变的负荷数据补齐和修正；

对GM(1,1)模型的建模可行性进行验证，并利用对数处理方法对不合格的序列x⁽⁰⁾进行修正；

对GM(1,1)模型中的α参数进行修正，以提高模型预测精度；

从不同角度选取数据序列构成不同预测方案，对预测日分段并计算各时段内各方案关联系数的均值，选择关联系数均值最大的方案为该时段预测方案；

最后利用后验差校验法对GM(1,1)模型精度进行检验；

分别从某电网数据中中选取三类不同的原始数据序列构建基于普通GM(1,1)模型的三种预测方案：

方案一：代表前5天同时刻数据，反映日同时刻负荷变化规律所得到的预测有功功率；

方案二：代表前5周同天同时刻数据，反映周同类型日同时刻负荷变化规律(7天)所得到的预测有功功率；

方案三：选取前一天24小时负荷，反映日平均负荷水平规律所得到的预测有功功率。

将预测日分为若干时段，通过比较以上三种预测方案关联系数的大小来确定实际采用的负荷预测方案并进行预测。三种普通GM(1,1)预测模型方案及改进GM(1,1)预测模型方案的预测图和结果分析如图3至图6所示：

表1列出了基于本技术方案提出的灰关联分段优选组合GM(1,1)模型的平均误差和基于普通GM（1,1）预测模型的三种方案预测结果的平均误差。

表1、采用不同预测方案预测结果对照表：

模型	预测精度
		基于灰关连分段优选组合GM(1，1)模型的预测方案	3.30％
基于普通GM(1，1)模型的预测方案一	4.35％
		基于普通GM(1，1)模型的预测方案二	3.89％
基于普通GM(1，1)模型的预测方案三	4.48％

通过表1的结果表明采用基于灰关联分段优选组合GM(1,1)模型的预测方案平均误差3.30%，优于基于普通GM(1,1)模型的预测方案。

最后应说明的是：以上所述仅为本发明的优选实施例而已，并不用于限制本发明，尽管参照前述实施例对本发明进行了详细的说明，对于本领域的技术人员来说，其依然可以对前述各实施例所记载的技术方案进行修改，或者对其中部分技术特征进行等同替换。凡在本发明的精神和原则之内，所作的任何修改、等同替换、改进等，均应包含在本发明的保护范围之内。

Claims (5)

1.一种基于灰色理论的短期电力负荷预测方法，其特征在于，包括以下步骤：

步骤1、辨识负荷数据中的缺陷数据，并对有缺损、突变的负荷数据补齐和修正；

步骤2、根据上述补齐和修正的负荷数据，对GM(1,1)模型的建模可行性进行验证，并利用对数处理方法对不合格的序列进行修正；

步骤3、对上述步骤二构建的GM(1，1)模型中的α参数进行修正；

步骤4、从不同角度选取数据序列构利用步骤三修正α参数后的GM(1，1)模型成不同预测方案进行预测，对预测日分段并计算各时段内各方案关联系数的均值，选择关联系数均值最大的方案为该时段预测方案；

步骤5、利用后验差校验法对GM(1,1)模型精度进行检验。

2.根据权利要求1所述的基于灰色理论的短期电力负荷预测方法，其特征在于，上述步骤1中的辨识负荷数据中的缺陷数据，并对有缺损、突变的负荷数据补齐和修正具体如下：

设某天的负荷为：L＝{l(1),l(2),...,l(96)}，对某天的96点负荷进行判断，若：|l(i)-l(i-1)|＜λ,(i＝2,3,...96λ为常数)，假如96点负荷相邻两点间的负荷全都落入λ阈内,将该日数据标为异常，并在预测取历史数据时直接予以剔除；

若有连续2个小时以上的数据点落入λ阈内，则按如下原则处理：

寻找相似日，并剔除已被判断为异常的相似日，找到数据异常日之前的N-1个负荷正常数据源，并将待修正日的数据源标记为N，N为常数；已知第N个数据源待补足时段96-k个点以外的k个点的历史负荷数据，k为常数，这k个点数据正常，则：

$\left\{\begin{array}{c}{P}_{\mathrm{AVE}}\left(t\right)=\frac{1}{N-1}\underset{i-1}{\overset{N-1}{\mathrm{\Σ}}}{P}_R(i,t),(t=\mathrm{1,2},...,96)\\ {\stackrel{\‾}{P}}_k=\frac{1}{k}\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{P}_R(N,t),(t=\mathrm{1,2},...,k)\\ {\stackrel{\‾}{P}}_{\mathrm{AVEK}}=\frac{1}{k}\underset{t=1}{\overset{k}{\mathrm{\Σ}}}{P}_{\mathrm{AVE}}\left(t\right),(t=\mathrm{1,2},...k)\\ P_R(N,t)=P_{\mathrm{AVE}}\left(t\right)-({\stackrel{\‾}{P}}_{\mathrm{AVEK}}-{\stackrel{\‾}{P}}_k),(t=k+1,...96)\end{array}\right.$

式中P_R(i,t)是第i个数据源观测值；P_AVE(t)是N-1个数据源观测值平均；

是第N个数据源k个正常点数据平均；P_AVEK是N-1个数据源在待修正日已知负荷时刻平均值。

3.根据权利要求1所述的基于灰色理论的短期电力负荷预测方法，其特征在于，上述步骤2中根据上述补齐和修正的负荷数据，对GM(1,1)模型的建模可行性进行验证，并利用对数处理方法对不合格的序列进行修正具体如下：

对于给定序列x⁽⁰⁾，可否建立精度较高的GM(1,1)，用x⁽⁰⁾的级比σ⁽⁰⁾(k)的大小与所属区间判断；

令x⁽⁰⁾为原始序列x⁽⁰⁾＝{x⁽⁰⁾(1),x⁽⁰⁾(2),x⁽⁰⁾(3),...,x⁽⁰⁾(n)}，对任意x⁽⁰⁾(k)∈x⁽⁰⁾k＝1,2,3,...,n，令σ⁽⁰⁾(k)为x⁽⁰⁾的级比，则

则当σ⁽⁰⁾(k)满足σ⁽⁰⁾(k)∈(0.1353,7.389)时，x⁽⁰⁾可作非畸形的GM(1,1)建模；

但建立精确度高的GM(1,1)模型，级比σ⁽⁰⁾(k)应落于靠近1的一个子区间(1-ε,1+ε)，即(1-ε,1+ε)

(0.1353,7.389)；

令x⁽⁰⁾为GM(1,1)建模序列x⁽⁰⁾＝{x⁽⁰⁾(1),x⁽⁰⁾(2),x⁽⁰⁾(3),...,x⁽⁰⁾(n)}，则x⁽⁰⁾的GM(1,1)模型发展系数a满足

建模序列x⁽⁰⁾的级比σ⁽⁰⁾(k)若满足

${\mathrm{\σ}}^{\left(0\right)}\left(k\right)\∈(e^{\frac{-2}{n+1}},e^{\frac{2}{n+1}})$

则x⁽⁰⁾是可以作GM(1,1)建模的序列数据。

4.根据权利要求3所述的基于灰色理论的短期电力负荷预测方法，其特征在于，对于上述级比σ⁽⁰⁾(k)检验不合格的序列，根据公式

σ_y(k)表示序列级比，Δ_y(k)为差异信息，y(k)变换后的数据序列元素，利用对数变换得到符合要求的序列y_m，

具体为：

令x为原始序列，y_m为x的m次对数序列，x表示为

x(x(1),x(2),...,x(n)),

y_m＝(y_m(1),y_m(2),...y_m(n)),

则，y_m(k)＝ln^mx(k)＝ln(ln(...(lnx(k))...))，

对m进行选择即可得到符合要求的序列y_m。

5.根据权利要求1至4任一所述的基于灰色理论的短期电力负荷预测方法，其特征在于，所述步骤3中的α参数修正包括以下步骤：

步骤301、设变量为x⁽⁰⁾的原始数据序列：x⁽⁰⁾＝[x⁽⁰⁾(1),x⁽⁰⁾(2),...x⁽⁰⁾(n)]用1-AGO累加生成，形成一阶累加生成序列：x⁽¹⁾＝[x⁽¹⁾(1),x⁽¹⁾(2),...x⁽¹⁾(n)]，其中： $x^{\left(1\right)}\left(k\right)=\underset{i=1}{\overset{k}{\mathrm{\Σ}}}{x}^{\left(0\right)}\left(i\right),k=\mathrm{1,2},...n;$

令z⁽¹⁾为x⁽¹⁾的均值序列：z⁽¹⁾＝(z⁽¹⁾(2),z⁽¹⁾(3),...z⁽¹⁾(n))其中：z⁽¹⁾(k+1)＝αx⁽¹⁾(k+1)+(1-α)x⁽¹⁾(k)第一次计算时，取α＝0.5，构造GM(1,1)的白化微分方程：

用最小二乘法求参数

记为 $P=\left[\begin{array}{c}\hat{a}\\ \hat{u}\end{array}\right],$ 在最小二乘准则下Y_n＝BP的解为

$P=\left[\begin{array}{c}\hat{a}\\ \hat{u}\end{array}\right]={\left(B^{T}B\right)}^{-1}{B}^{T}Y$

上式为GM(1,1)参数

的矩阵辨识算式，式中(B^TB)^-1B^T是数据矩阵B的广义逆矩阵；其中： $B=\left[\begin{array}{cc}{-z}^{\left(1\right)}\left(2\right)& 1\\ {-z}^{\left(1\right)}\left(3\right)& 1\\ \·& \·\\ \·& \·\\ \·& \·\\ {-z}^{\left(1\right)}\left(n\right)& 1\end{array}\right],$ Y_n＝[x⁽⁰⁾(2),x⁽⁰⁾(3),...x⁽⁰⁾(n)]^T；

步骤302：

将上述步骤301中解得的参数

重新计算

记为将与上一次计算所用的

进行比较，给定任意小的正整数ε，如果

转入上述步骤301继续计算，将

代入z⁽¹⁾(k+1)＝αx⁽¹⁾(k+1)+(1-α)x⁽¹⁾(k)计算背景值z⁽¹⁾(k+1)，再次进行GM(1,1)建模和预测运算；

当

时，迭代结束，转入步骤303；

步骤303：

建立GM(1,1)的预测模型为： ${\hat{x}}^{\left(1\right)}(k+1)=(x^{\left(0\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})e^{-\hat{a}k}+\frac{\hat{u}}{\hat{a}},k=\mathrm{1,2},...n,$ 对x⁽¹⁾(k)进行累减还原，得x⁽⁰⁾(k)的预测值，即GM(1,1)预测模型为：

${\hat{x}}^{\left(0\right)}(k+1)={\hat{x}}^{\left(1\right)}(k+1)-{\hat{x}}^{\left(1\right)}\left(k\right)=(x^{\left(0\right)}\left(1\right)-\frac{\hat{u}}{\hat{a}})(1-e^{\hat{a}})e^{-\hat{a}k},k=\mathrm{1,2},...n$

输出运行结果。

Images