CN112508244A - 一种用户级综合能源系统多元负荷预测方法 - Google Patents

Abstract

本发明公开了一种用户级综合能源系统多元负荷预测方法。运用自适应噪声的完全集合经验模态分解分别对电负荷、冷负荷、热负荷进行本征模态分解，对分解得到的强非平稳分量运用变分模态分解再次进行分解；基于核主成分分析对天气日历规则特征集提取主成分实现数据降维；将分解得到的非平稳、平稳序列分量分别与天气日历规则特征集提取的主成分组成输入集，分别用DBiLSTM神经网络、MLR进行预测，最后将预测结果进行重构得到最终预测结果。本发明能够深度挖掘用户级综合能源系统中多元负荷的耦合特性，有效提升了用户级综合能源系统多元负荷预测精度，在面对随机性、波动性大的多元负荷时也能有很好的预测效果。

Description

一种用户级综合能源系统多元负荷预测方法

技术领域

本发明涉及综合能源系统多元负荷预测领域，具体涉及一种用户级综合能源系统多元负荷预测方法。

背景技术

为解决能源短缺，实现能量精确管理，综合能源系统(Integrated EnergySystem，IES)应运而生。多种能源在IES中进行能源供给、转换，使得电、冷、热负荷之间存在一定的耦合关系。伴随电力市场的发展、能量管理系统的应用，对负荷预测提出了更高的要求。

传统的负荷预测方法为统计学方法、机器学习法，包括卡尔曼滤波、多元线性回归(Multiple Linear Regression，MLR)，支持向量回归、随机森林等。统计学方法模型简单，预测速度快，但只适用于平稳序列，难以对非线性序列准确拟合；机器学习在面对非线性序列时预测性能有所提升，但忽略了样本的时序关系。

近年来深度学习在负荷预测中应用广泛，最为典型的代表就是长短期记忆(LongShort-Term Memory，LSTM)神经网络及其变体，此外，深度置信网络等也有所应用。尽管深度学习能够挖掘序列之间存在的非线性关系，可以学到其时序特征，但在面对非连续序列时，预测精度也有待提高。因此，有学者提出了组合预测方法，通过结合不同模型的各自优势来提高预测精度。

但是，目前大部分研究方法只对电力负荷进行预测，单一负荷预测手段忽视了多元负荷之间的耦合特性，无法对IES多元负荷进行准确预测。目前IES多元负荷预测的研究还较少，并且目前研究大多针对区域级IES，对于用户级IES来说，多元负荷随机、波动性更强，负荷预测面临更大的挑战。

发明内容

有鉴于此，本发明的目的在于提供一种基于数据预处理与DBiLSTM-MLR的用户级综合能源系统多元负荷预测方法，针对IES多元负荷随机、波动性强的特点，对历史数据进行预处理，通过运用二次模态分解(Quadratic Modal Decomposition，QMD)将多元负荷分解为不同的本征模态函数(Intrinsic Mode Functions，IMF)，并运用核主成分分析(Kernel Principal Component Analysis，KPCA)对天气日历规则特征集数据进行降维，进而降低预测难度。紧接着运用深度双向长短期记忆(Deep Didirectional Long Short-Term Memory，DBiLSTM)神经网络深度挖掘多元负荷的耦合特性，结合MLR能够对分解后的平稳序列准确、快速预测的特点，基于DBiLSTM神经网络和MLR分别对分解后的非平稳、平稳分量进行预测，最后将预测结果进行重构得到最终预测结果，在数据预处理下结合DBiLSTM神经网络和MLR各自优势，有效提高了预测精度。

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

一种用户级综合能源系统多元负荷预测方法，包括如下步骤：

步骤1：运用自适应噪声的完全集合经验模态分解分别对电负荷、冷负荷、热负荷进行本征模态分解，对分解得到的强非平稳分量运用变分模态分解再次进行分解；

步骤2：构建天气日历规则特征数据集，通过KPCA算法对天气日历规则特征集进行降维得到降维后天气日历规则主成分数据集；

步骤3：将二次分解得到的原始电负荷平稳序列分量、二次分解得到的原始冷负荷平稳序列分量、二次分解得到的原始热负荷平稳序列分量、降维后天气日历规则特征数据集组成MLR输入数据集，通过优化训练MLR模型得到优化后MLR模型；

步骤4：将二次分解得到的原始电负荷非平稳序列分量、二次分解得到的原始冷负荷非平稳序列分量、二次分解得到的原始热负荷非平稳序列分量、降维后天气日历规则特征数据集组成神经网络输入数据集，通过优化训练神经网络得到优化后神经网络模型；

步骤5：将原始电负荷平稳序列分量、原始冷负荷平稳序列分量、原始热负荷平稳序列分量与降维后天气日历规则特征数据集组成输入数据集，通过优化后MLR模型预测得到待预测的电负荷平稳序列分量、待预测的冷负荷平稳序列分量、待预测的热负荷平稳序列分量；将原始电负荷非平稳序列分量、原始冷负荷非平稳序列分量、原始热负荷非平稳序列分量与降维后天气日历规则特征数据集组成输入数据集，通过优化后神经网络模型预测得到待预测的电负荷非平稳序列分量、待预测的冷负荷非平稳序列分量、待预测的热负荷非平稳序列分量，将优化后MLR模型预测得到的所有待预测的负荷平稳序列分量与优化后神经网络模型预测得到的所有待预测的负荷非平稳序列分量相加得到最终的待预测电负荷值、待预测冷负荷值、待预测热负荷值。

作为优选，所述步骤1：

步骤1.1：引入原始负荷序列；

所述原始负荷序列为：x(t),t∈[1,T]；

其中，x(t)为原始负荷序列中t时刻的原始负荷，T为原始负荷的数量；

步骤1.2：向原始负荷序列加入正白噪声序列得到加入正白噪声后的负荷序列；

所述加入正白噪声后的负荷序列为：

m⁺(t)＝x(t)+n⁺(t),t∈[1,T]

其中，n⁺(t)为正白噪声序列中t时刻引入的正白噪声，T为原始负荷的数量即正白噪声的数量，m⁺(t)为加入正白噪声后的负荷序列第t时刻的负荷；

所述正白噪声序列为：n⁺(t),t∈[1,T]；

所述加入正白噪声后的负荷序列为：m⁺(t),t∈[1,T]；

步骤1.3：向原始负荷序列加入负白噪声序列得到加入负白噪声后的负荷序列；

所述加入负白噪声后的负荷序列为：

m^-(t)＝x(t)+n^-(t),t∈[1,T]

式中，n^-(t)为负白噪声序列中t时刻引入的负白噪声，T为原始负荷的数量即负白噪声的数量，m^-(t)为加入负白噪声后的负荷序列第t时刻的负荷；

所述负白噪声序列为：n^-(t),t∈[1,T]；

所述加入负白噪声后的负荷序列为：m^-(t),t∈[1,T]；

步骤1.4：运用经验模态分解法对加入正白噪声后的负荷序列进行模态分解得到经验模态分解后对应的一组本征模态负荷序列组F⁺(u,t),(u∈[1,U],t∈[1,T])，运用经验模态分解法对加入负白噪声后的负荷序列进行模态分解得到经验模态分解后对应的一组本征模态负荷序列组F^-(v,t),(v∈[1,V],t∈[1,T])，其中F⁺(u,t)表示加入正白噪声后的负荷序列进行经验模态分解后第u个序列分量第t时刻的负荷，F^-(v,t)表示加入负白噪声后的负荷序列进行经验模态分解后第v个序列分量第t时刻的负荷，U和V均表示分解后的序列分量数且U＝V；

所述加入正白噪声后的负荷序列为：[m⁺(1),m⁺(2),...,m⁺(T)]^T；

所述加入负白噪声后的负荷序列为：[m^-(1),m^-(2),...,m^-(T)]^T；

所述加入正白噪声负荷序列经验模态分解后的一组本征模态负荷序列中第u个序列分量为：[F⁺(u,1),F⁺(u,2),...,F⁺(u,T)]^T,u∈[1,U]；

所述加入负白噪声负荷序列经验模态分解后的一组本征模态负荷序列中第v个序列分量为：[F^-(v,1),F^-(v,2),...,F^-(v,T)]^T,v∈[1,V]；

步骤1.5：重复执行步骤1.1至步骤1.4多次，将多次得到的经验模态分解后的正白噪声负荷序列对应的本征模态负荷序列组通过均值计算得到经验模态分解后正白噪声负荷序列对应的一组本征模态平均负荷序列组，将多次得到的经验模态分解后的负白噪声负荷序列对应的本征模态负荷序列组通过均值计算得到经验模态分解后负白噪声负荷序列对应的一组本征模态平均负荷序列组；

所述通过均值计算得到经验模态分解后的正白噪声对应的一组本征模态平均负荷序列中第u个序列分量第t时刻平均负荷为：

式中，

为第i次经验模态分解后的正白噪声负荷序列对应的本征模态负荷序列组中第u个序列分量第t时刻的负荷，

为经验模态分解后的正白噪声负荷序列对应的一组本征模态平均负荷序列中第u个序列分量第t时刻平均负荷，N为重复执行的次数；

所述经验模态分解后的正白噪声负荷序列对应的一组本征模态平均负荷序列组为：

所述通过均值计算得到经验模态分解后的负白噪声对应的一组本征模态平均负荷序列中第v个序列分量第t时刻平均负荷为：

式中，

为第i次经验模态分解后的负白噪声负荷序列对应的本征模态负荷序列组中第v个序列分量第t时刻的负荷，

为经验模态分解后的负白噪声负荷序列对应的一组本征模态平均负荷序列中第v个序列分量第t时刻平均负荷，N为重复执行的次数；

所述经验模态分解后的负白噪声负荷序列对应的一组本征模态平均负荷序列组为：

步骤1.6：将经验模态分解后的正白噪声负荷序列对应的一组本征模态平均负荷序列组、经验模态分解后的负白噪声负荷序列对应的一组本征模态平均负荷序列组通过均值计算得到最终一组负荷分解序列组；

所述通过均值计算得到最终负荷分解序列组为：

式中，load(k,t)为最终负荷分解序列组中第k个序列分量第t时刻负荷；

步骤1.7：将最终负荷分解序列组中每个序列分量通过近似熵法得到每个序列分量的平稳性强度，根据平稳性强度从左到右由弱至强的顺序，对最终负荷分解序列组中每个序列分量进行排序，得到平稳性强度排序后序列分量；

所述最终负荷分解序列组中每个序列分量为：

其中，load表示原始负荷序列经CEEMDAN最终分解的序列组，G_k表示最终分解序列组中第k个本征模态分量序列，G_k(t)表示最终分解序列组中第k个本征模态分量序列第t时刻负荷，U表示分解后的序列分量数，T表示一个序列分量内包含样本的数量；

所述每个序列分量的平稳性强度为：

load_stat＝[ξ₁,...,ξ_k,...,ξ_U],k∈[1,U]

式中，load_stat为每个序列分量的平稳性强度集；ξ_k,k∈[1,U]表示第k个序列分量的平稳性强度，U表示分解后的序列分量数；

所述平稳性强度排序后序列分量为：

其中，load_SE表示平稳性强度排序后序列分量，g_γ表示平稳性强度排序后序列分量中第γ个本征模态分量序列，g_γ(t)表示平稳性强度排序后序列分量中第γ个本征模态分量序列第t时刻负荷，U表示分解后的序列分量数，T表示一个序列分量内包含样本的数量。

load_SE进一步简化为：

load_SE＝[g₁,...,g_γ,...,g_U],γ∈[1,U]，g_γ＝[g_γ(1),...,g_γ(t),...,g_γ(T)]^T,t∈[1,T]

所述平稳性强度排序后序列分量组load_SE对应的每个序列分量的平稳性强度更新为：

式中，

为平稳性强度排序后序列分量组内每个序列分量的平稳性强度集；χ_γ,γ∈[1,U]表示第γ个序列分量的平稳性强度，U表示分解后的序列分量数；

以τ₁作为平稳与非平稳的分界值，τ₂作为非平稳与强非平稳的分界值，对每个序列的平稳性强度

进行分析；

对任意γ∈[1,U]，若χ_γ＜τ₁则该序列分量归为平稳序列分量，若τ₁≤χ_γ≤τ₂则该序列分量归为非平稳序列分量，若χ_γ＞τ₂则该序列分量归为强非平稳序列分量；

U1个强非平稳分量序列构建的强非平稳分量集合即：

load_sns＝[g₁,g₂,...,g_q,...,g_U1],q∈[1,U1]；

U2个非平稳分量序列构建的非平稳分量序列集合即：

load_ns＝[g_U1+1,g_U1+2,...,g_U1+z,...,g_U1+U2],z∈[1,U2]

U3个平稳分量序列构建的平稳分量序列集合即：

load_s＝[g_U1+U2+1,g_U1+U2+2,...,g_U1+U2+b,...,g_U1+U2+U3],b∈[1,U3]

其中：

U1+U2+U3＝U

g_q＝[g_q(1),g_q(2),...,g_q(T)]^T,q∈[1,U1]

g_U1+z＝[g_U1+z(1),g_U1+z(2),...,g_U1+z(T)]^T,z∈[1,U2]

g_U1+U2+b＝[g_U1+U2+b(1),g_U1+U2+b(2),...,g_U1+U2+b(T)]^T,b∈[1,U3]

步骤1.8：将每个强非平稳分量序列通过变分模态分解进行进一步分解，得到强非平稳分量序列经变分模态分解后的平稳分量序列与非平稳分量序列组；

所述强非平稳分量序列通过VMD进行进一步分解得到的强非平稳分量序列经VMD分解所得序列组为：

简写

u_r＝[u_r(1),...,u_r(t),...,u_r(T)]^T,t∈[1,T]，u_r表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列；u_r(t)表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列第t时刻负荷；R表示分解得到的序列分量数；T表示一个序列分量内包含样本的数量；

将load_sns＝[g₁,g₂,...,g_q,...,g_U1],q∈[1,U1]中g_q＝[g_q(1),g_q(2),...,g_q(T)]^T,q∈[1,U1]中每个强非平稳分量序列通过VMD进行进一步分解为：

所述VMD方法的原理为：寻找R个本征模态函数分量集合{u₁,...,u_R}，每个本征模态函数分量对应的中心频率为{ω₁,...,ω_R}，以R个本征模态函数分量之和等于待分解序列g_q为约束条件，使得分解后每个u_r,r∈[1,R]的估计带宽之和最小；

所述VMD方法的具体步骤如下：

构造约束变分最优问题：

式中，g_q(t)为load_sns中第q个序列第t时刻负荷，{u₁,...,u_R}为g_q分解后的模态集合，其中u₁(t)指VMD分解后的第1个模态的第t时刻负荷，u_R(t)指VMD分解后的第R个模态的第t时刻负荷；{ω₁,...,ω_R}为{u₁,...,u_R}中各序列分量对应的中心频率，其中ω₁为负荷序列[u₁(1),u₁(2),...,u₁(T)]^T的中心频率，ω_R为负荷序列[u_R(1),u_R(2),...,u_R(T)]^T的中心频率；

表示对t求偏导；δ(t)为Dirac分布；

*为卷积运算符；

表示求二范数的平方。

利用二次惩罚因子α和Lagrange乘子λ(t)构造Lagrange函数，将上述约束问题转化为无约束问题，转变后的扩展的Lagrange表达式为：

式中，L(·)表示Lagrange函数；λ(t)为Lagrange乘子，用来严格保持约束条件；α为二次惩罚因子，保证重构精度；＜·＞表示求内积。

利用交替方向乘子算法寻找最优解，不断更新u_r、ω_r，公式如下：

式中，m为迭代次数，m≤M，M为最大迭代次数；^表示傅里叶变换；

为当前信号的维纳滤波，为

的傅里叶变换，表示分解后第r个模态第m+1次迭代后的频域值；

为对应序列

的功率谱重心；

为对应g_q(t)经傅里叶变换后的频域值；

为对应u(t)经傅里叶变换后的频域值；

为对应λ(t)经傅里叶变换后的频域值；

所述将强非平稳分量序列通过VMD进一步分解得到的强非平稳分量序列经变分模态分解所得的分量序列组为：

简写

u_r＝[u_r(1),...,u_r(t),...,u_r(T)]^T,t∈[1,T]，u_r表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列；u_r(t)表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列第t时刻负荷；R表示分解得到的序列分量数；T表示一个序列分量内包含样本的数量；

将强非平稳分量序列通过VMD进一步分解得到的分量序列组中每个序列分量通过近似熵法得到每个序列分量的平稳性强度：

以τ₁作为平稳与非平稳的分界值，τ₂作为非平稳与强非平稳的分界值，对每个序列的平稳性强度

进行分析。对任意r∈[1,R]，若κ_r＜τ₁则该序列分量归为平稳序列分量，若τ₁≤κ_r≤τ₂则该序列分量归为非平稳序列分量；

步骤1.9：构建原始负荷序列经过二次分解最终得到的平稳序列分量组、原始负荷序列经过二次分解最终得到的非平稳序列分量组；

将原始负荷序列经CEEMDAN分解得到的平稳序列组以及原始负荷序列经CEEMDAN分解得到的强平稳序列组再次经VMD进行分解后所得的强平稳序列组经VMD再次进行分解所得的平稳序列组进行组合，最终得到I个平稳序列分量数，并将原始负荷序列经二次分解最终得到的平稳序列组定义为Load_s；

将原始负荷序列经CEEMDAN分解得到的非平稳序列组以及原始负荷序列经CEEMDAN分解得到的强平稳序列组再次经VMD进行分解后所得的强平稳序列组经VMD再次进行分解所得的非平稳序列组进行组合，最终得到J个非平稳序列分量数，并将原始负荷序列经二次分解最终得到的非平稳序列组定义为Load_ns。其中，Load_ns中每个序列分量已通过近似熵法得到每个序列分量的平稳性强度，并根据平稳性强度从左到右由弱至强的顺序进行排序。

所述原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组Load_s为：

简写L_s：Load_s＝[A₁,...,A_α,...,A_I],α∈[1,I]，A_α＝[A_α(1),...,A_α(t),...,A_α(T)]^T,t∈[1,T]，A_α表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；A_α(t)表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I表示经过二次分解后最终的平稳序列分量数；T表示每个序列分量包含的样本数。

所述原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组Load_ns为：

简写Load_ns：Load_ns＝[B₁,...,B_β,...,B_J],β∈[1,J]，B_β＝[B_β(1),...,B_β(t),...,B_β(T)]^T,t∈[1,T]，B_β表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；B_β(t)表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J表示经过二次分解后最终的非平稳序列分量数；T表示每个序列分量包含的样本数。

由步骤1.1至步骤1.9，可以分别得到原始电、冷、热负荷序列分别经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的原始电、冷、热负荷序列各自对应的平稳、非平稳序列分量组为：

所述原始电负荷序列为：x^e(t),t∈[1,T]；

其中，x^e(t)为原始电负荷序列中t时刻的原始电负荷，T为原始电负荷的数量。

所述原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后，得到原始电负荷平稳序列分量组

为：

简写

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I表示经过二次分解后最终的平稳序列分量数；T表示每个序列分量包含的样本数。

所述原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始电负荷非平稳序列分量组

为：

简写

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J表示经过二次分解后最终的非平稳序列分量数；T表示每个序列分量包含的样本数。

所述原始冷负荷序列为：x^c(t),t∈[1,T]；

其中，x^c(t)为原始冷负荷序列中t时刻的原始冷负荷，T为原始冷负荷的数量。

所述原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始冷负荷平稳序列分量组

为：

简写

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I表示经过二次分解后最终的平稳序列分量数；T表示每个序列分量包含的样本数。

所述原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始冷负荷非平稳序列分量组

为：

简写

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J表示经过二次分解后最终的非平稳序列分量数；T表示每个序列分量包含的样本数。

所述原始热负荷序列为：x^h(t),t∈[1,T]；

其中，x^h(t)为原始热负荷序列中t时刻的原始热负荷，T为原始热负荷的数量。

所述原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始热负荷平稳序列分量组

为：

简写

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I表示经过二次分解后最终的平稳序列分量数；T表示每个序列分量包含的样本数。

所述原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始热负荷非平稳序列分量组

为：

简写

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J表示经过二次分解后最终的非平稳序列分量数；T表示每个序列分量包含的样本数。

作为优选，步骤2所述天气日历规则特征数据集为：

data_h＝(month_h,day_h,week_h,hour_h,holiday_h,d_h,1,d_h,2,...,d_h,Q)

h∈[1,T]

其中，data_h为天气日历规则特征数据集中第h个样本，h表示第h个小时，T为天气日历规则特征数据集中样本的数量即小时的数量，month_h为天气日历规则特征数据集中第h个样本中月份序号标识，day_h为天气日历规则特征数据集中第h个样本中天序号标识，week_h为天气日历规则特征数据集中第h个样本中星期序号标识，hour_h为天气日历规则特征数据集中第h个样本中小时序号标识，holiday_h为天气日历规则特征数据集中第h个样本中节假日序号标识，d_h,q*为天气日历规则特征数据集中第h个样本中第q^*种类型的天气因素，q^*∈[1,Q]，Q为天气因素类型的数量；

步骤2所述通过KPCA算法对天气日历规则特征数据集进行降维得到降维后天气日历规则的主成分数据集为：

data_h中天气日历规则特征数据集中样本的数量即小时的数量T在降维中保持不变，对data_h中Q+5个元素即Q+5个维度通过KPCA算法进行降维；

步骤2所述降维后天气日历规则主成分数据集为：

data^* _h＝(W_h,1,W_h,2,...,W_h,H)

h∈[1,T]

其中，data^* _h为降维后天气日历规则特征数据集中第h个样本，h表示第h个小时，T为降维后天气日历规则特征数据集中样本的数量即小时的数量，

为降维后天气日历规则因素特征数据集中第h个样本中第p^*个元素，H为降维后天气日历规则特征数据集中每个样本的维度；

作为优选，步骤3所述二次分解得到的原始电负荷平稳序列分量为：

其中，

表示二次分解得到的原始电负荷平稳序列分量组内第α个分量；

步骤3所述二次分解得到的原始冷负荷平稳序列分量为：

其中，

表示二次分解得到的原始冷负荷平稳序列分量组内第α个分量；

步骤3所述二次分解得到的原始热负荷平稳序列分量为：

其中，

表示二次分解得到的原始热负荷平稳序列分量组内第α个分量；

步骤3所述通过优化训练MLR模型得到优化后MLR模型为：

步骤3.1，建立MLR回归模型；

模型内待优化参数即各回归系数为：

其中，

表示第k^*个电负荷序列回归系数，

表示第k^*个冷负荷序列回归系数，

表示第k^*个热负荷序列回归系数，k^*表示参与回归的变量共有k^*个，常数项、

的维度、

的维度、

的维度、降维后天气日历规则特征数据集中H个主成分的维度构成k^*；

步骤3.2，将二次分解得到的原始电负荷平稳序列分量、二次分解得到的原始冷负荷平稳序列分量、二次分解得到的原始热负荷平稳序列分量、降维后天气日历规则特征数据集组成输入集输入MLR进行回归；

步骤3.3，采用最小二乘估计求得电负荷序列回归系数、冷负荷序列回归系数、热负荷序列回归系数，用于构建优化后MLR回归模型；

作为优选，步骤4所述二次分解得到的原始电负荷非平稳序列分量为：

其中，

为二次分解得到的原始电负荷非平稳序列分量组内第β个分量；

步骤4所述为二次分解得到的原始冷负荷非平稳序列分量：

其中，

为二次分解得到的原始冷负荷非平稳序列分量组内第β个分量；

步骤4所述二次分解得到的原始热负荷非平稳序列分量为：

其中，

为二次分解得到的原始热负荷非平稳序列分量组内第β个分量；

步骤4.1，所述神经网络为包含多个隐含层的DBiLSTM神经网络；

各隐含层神经元个数分别为：

θ表示隐含层数量；

每个隐含层中待优化参数为：

第Γ个隐含层的遗忘门权重矩阵即

与第Γ个隐含层的遗忘门偏置项即

第Γ个隐含层的输入门第一个权重矩阵

与第Γ个隐含层的输入门第一个偏置项

第Γ个隐含层的输入门第二个权重矩阵

与第Γ个隐含层的输入门第二个偏置项

第Γ个隐含层的输出门权重矩阵

与第Γ个隐含层的输出门偏置项

Γ∈[1,θ]；

步骤4.2，将

以及降维后天气日历规则特征数据集组成神经网络输入数据集输入至神经网络进行训练；

步骤4.3，建立模型训练的损失函数，运用优化算法优化神经网络内部每个隐含层中待优化参数即

所述神经网络的损失函数表达式为：

式中，loss为损失函数值；S为参加计算损失值的预测点个数；

为第ζ时刻电负荷的真实值；

为第ζ时刻冷负荷的真实值；

为第ζ时刻热负荷的真实值；

为第ζ时刻电负荷的预测值；

为第ζ时刻冷负荷的预测值；

为第ζ时刻热负荷的预测值。

步骤4.4，结合神经网络的损失函数，将将二次分解得到的原始电负荷非平稳序列分量、二次分解得到的原始冷负荷非平稳序列分量、二次分解得到的原始热负荷非平稳序列分量、降维后天气日历规则特征数据集组成神经网络输入数据集，通过Adam算法优化训练得到优化后神经网络模型。

本发明采用了QMD二次分解方法，即首先运用CEEMDAN对多元负荷进行一次模态分解后再次运用VMD进行二次模态分解，将多元负荷强非平稳序列分解为多个更为平稳的模态分量，降低预测难度；

本发明采用KPCA方法提取输入特征集的主成分，保留各特征变量之间非线性关系下提取出主要成分代替原始多维特征，降低模型输入维度，简化了模型；

本发明采用的DBiLSTM神经网络相比常规LSTM神经网络能够更加全面地学习历史数据的整体信息以及多元负荷之间的耦合关系，具有更高的学习能力；

本发明结合了MLR在面对平稳序列时预测精度高、速度快的特点，将平稳序列由MLR进行预测，提升预测精度的同时大大降低了时间成本。

附图说明

图1：为本发明基于数据预处理与DBiLSTM-MLR的用户级综合能源系统多元负荷预测流程图；

图2：为本发明LSTM单元结构图；

图3：为本发明DBiLSTM神经网络结构图；

图4：为本发明采用的实验数据波形图；

图5：为本发明采用CEEMDAN对电负荷进行分解结果图；

图6：为本发明采用CEEMDAN对冷负荷进行分解结果图；

图7：为本发明采用CEEMDAN对热负荷进行分解结果图；

图8：为发明采用VMD对电负荷IMF1分量再次进行分解结果图；

图9：为发明采用VMD对冷负荷IMF1分量再次进行分解结果图；

图10：为发明采用VMD对热负荷IMF1分量再次进行分解结果图；

图11：为本发明采用KPCA对特征集提取主成分结果图；

图12：为本发明采用模态分解前后电、冷、热负荷预测结果对比分析图；

图13：为本发明采用不同方法与本发明方法进行的电、冷、热负荷预测结果对比分析图。

具体实施方式

下面将结合本发明实施例中的附图，对本发明实施例中的技术方案进行清楚、完整地描述，显然，所描述的实施例仅是本发明一部分实施例，而不是全部的实施例。基于本发明中的实施例，本领域普通技术人员在没有做出创造性劳动前提下所获得的所有其他实施例，都属于本发明保护的范围。

下面结合附图1-附图13介绍本发明的具体实施方式。

如附图1-附图3所示，本发明的一种用户级综合能源系统多元负荷预测方法，具体实现包括以下步骤：

步骤1：运用自适应噪声的完全集合经验模态分解分别对电负荷、冷负荷、热负荷进行本征模态分解，对分解得到的强非平稳分量运用变分模态分解再次进行分解；

进一步的，所述步骤1：

步骤1.1：引入原始负荷序列；

所述原始负荷序列为：x(t),t∈[1,T]；

其中，x(t)为原始负荷序列中t时刻的原始负荷，T＝17520为原始负荷的数量；

步骤1.2：向原始负荷序列加入正白噪声序列得到加入正白噪声后的负荷序列；

所述加入正白噪声后的负荷序列为：

m⁺(t)＝x(t)+n⁺(t),t∈[1,T]

其中，n⁺(t)为正白噪声序列中t时刻引入的正白噪声，T＝17520为原始负荷的数量即正白噪声的数量，m⁺(t)为加入正白噪声后的负荷序列第t时刻的负荷；

所述正白噪声序列为：n⁺(t),t∈[1,T]；

所述加入正白噪声后的负荷序列为：m⁺(t),t∈[1,T]；

步骤1.3：向原始负荷序列加入负白噪声序列得到加入负白噪声后的负荷序列；

所述加入负白噪声后的负荷序列为：

m^-(t)＝x(t)+n^-(t),t∈[1,T]

式中，n^-(t)为负白噪声序列中t时刻引入的负白噪声，T＝17520为原始负荷的数量即负白噪声的数量，m^-(t)为加入负白噪声后的负荷序列第t时刻的负荷；

所述负白噪声序列为：n^-(t),t∈[1,T]；

所述加入负白噪声后的负荷序列为：m^-(t),t∈[1,T]；

步骤1.4：运用经验模态分解法对加入正白噪声后的负荷序列进行模态分解得到经验模态分解后对应的一组本征模态负荷序列组F⁺(u,t),(u∈[1,U],t∈[1,T])，运用经验模态分解法对加入负白噪声后的负荷序列进行模态分解得到经验模态分解后对应的一组本征模态负荷序列组F^-(v,t),(v∈[1,V],t∈[1,T])，其中F⁺(u,t)表示加入正白噪声后的负荷序列进行经验模态分解后第u个序列分量第t时刻的负荷，F^-(v,t)表示加入负白噪声后的负荷序列进行经验模态分解后第v个序列分量第t时刻的负荷，U和V均表示分解后的序列分量数且U＝V＝12；

所述加入正白噪声后的负荷序列为：[m⁺(1),m⁺(2),...,m⁺(T)]^T；

所述加入负白噪声后的负荷序列为：[m^-(1),m^-(2),...,m^-(T)]^T；

所述加入正白噪声负荷序列经验模态分解后的一组本征模态负荷序列中第u个序列分量为：[F⁺(u,1),F⁺(u,2),...,F⁺(u,T)]^T,u∈[1,U]；

所述加入负白噪声负荷序列经验模态分解后的一组本征模态负荷序列中第v个序列分量为：[F^-(v,1),F^-(v,2),...,F^-(v,T)]^T,v∈[1,V]；

步骤1.5：重复执行步骤1.1至步骤1.4多次，将多次得到的经验模态分解后的正白噪声负荷序列对应的本征模态负荷序列组通过均值计算得到经验模态分解后正白噪声负荷序列对应的一组本征模态平均负荷序列组，将多次得到的经验模态分解后的负白噪声负荷序列对应的本征模态负荷序列组通过均值计算得到经验模态分解后负白噪声负荷序列对应的一组本征模态平均负荷序列组；

所述通过均值计算得到经验模态分解后的正白噪声对应的一组本征模态平均负荷序列中第u个序列分量第t时刻平均负荷为：

式中，

为第i次经验模态分解后的正白噪声负荷序列对应的本征模态负荷序列组中第u个序列分量第t时刻的负荷，

为经验模态分解后的正白噪声负荷序列对应的一组本征模态平均负荷序列中第u个序列分量第t时刻平均负荷，N＝10为重复执行的次数；

所述经验模态分解后的正白噪声负荷序列对应的一组本征模态平均负荷序列组为：

所述通过均值计算得到经验模态分解后的负白噪声对应的一组本征模态平均负荷序列中第v个序列分量第t时刻平均负荷为：

式中，

为第i次经验模态分解后的负白噪声负荷序列对应的本征模态负荷序列组中第v个序列分量第t时刻的负荷，

为经验模态分解后的负白噪声负荷序列对应的一组本征模态平均负荷序列中第v个序列分量第t时刻平均负荷，N＝10为重复执行的次数；

所述经验模态分解后的负白噪声负荷序列对应的一组本征模态平均负荷序列组为：

步骤1.6：将经验模态分解后的正白噪声负荷序列对应的一组本征模态平均负荷序列组、经验模态分解后的负白噪声负荷序列对应的一组本征模态平均负荷序列组通过均值计算得到最终一组负荷分解序列组；

所述通过均值计算得到最终负荷分解序列组为：

式中，load(k,t)为最终负荷分解序列组中第k个序列分量第t时刻负荷；

步骤1.7：将最终负荷分解序列组中每个序列分量通过近似熵法得到每个序列分量的平稳性强度，根据平稳性强度从左到右由弱至强的顺序，对最终负荷分解序列组中每个序列分量进行排序，得到平稳性强度排序后序列分量；

所述最终负荷分解序列组中每个序列分量为：

其中，load表示原始负荷序列经CEEMDAN最终分解的序列组，G_k表示最终分解序列组中第k个本征模态分量序列，G_k(t)表示最终分解序列组中第k个本征模态分量序列第t时刻负荷，U＝12表示分解后的序列分量数，T＝17520表示一个序列分量内包含样本的数量；

所述每个序列分量的平稳性强度为：

load_stat＝[ξ₁,...,ξ_k,...,ξ_U],k∈[1,U]

式中，load_stat为每个序列分量的平稳性强度集；ξ_k,k∈[1,U]表示第k个序列分量的平稳性强度，U＝12表示分解后的序列分量数；

所述平稳性强度排序后序列分量为：

其中，load_SE表示平稳性强度排序后序列分量，g_γ表示平稳性强度排序后序列分量中第γ个本征模态分量序列，g_γ(t)表示平稳性强度排序后序列分量中第γ个本征模态分量序列第t时刻负荷，U＝12表示分解后的序列分量数，T＝17520表示一个序列分量内包含样本的数量。

load_SE进一步简化为：

load_SE＝[g₁,...,g_γ,...,g_U],γ∈[1,U]，g_γ＝[g_γ(1),...,g_γ(t),...,g_γ(T)]^T,t∈[1,T]

所述平稳性强度排序后序列分量组load_SE对应的每个序列分量的平稳性强度更新为：

式中，

为平稳性强度排序后序列分量组内每个序列分量的平稳性强度集；χ_γ,γ∈[1,U]表示第γ个序列分量的平稳性强度，U＝12表示分解后的序列分量数；

以τ₁＝1作为平稳与非平稳的分界值，τ₂＝2作为非平稳与强非平稳的分界值，对每个序列的平稳性强度

进行分析；

对任意γ∈[1,U]，若χ_γ＜τ₁则该序列分量归为平稳序列分量，若τ₁≤χ_γ≤τ₂则该序列分量归为非平稳序列分量，若χ_γ＞τ₂则该序列分量归为强非平稳序列分量；

U1＝1个强非平稳分量序列构建的强非平稳分量集合即：

load_sns＝[g₁,g₂,...,g_q,...,g_U1],q∈[1,U1]；

U2＝2个非平稳分量序列构建的非平稳分量序列集合即：

load_ns＝[g_U1+1,g_U1+2,...,g_U1+z,...,g_U1+U2],z∈[1,U2]

U3＝9个平稳分量序列构建的平稳分量序列集合即：

load_s＝[g_U1+U2+1,g_U1+U2+2,...,g_U1+U2+b,...,g_U1+U2+U3],b∈[1,U3]

其中：

U1+U2+U3＝U

g_q＝[g_q(1),g_q(2),...,g_q(T)]^T,q∈[1,U1]

g_U1+z＝[g_U1+z(1),g_U1+z(2),...,g_U1+z(T)]^T,z∈[1,U2]

g_U1+U2+b＝[g_U1+U2+b(1),g_U1+U2+b(2),...,g_U1+U2+b(T)]^T,b∈[1,U3]

步骤1.8：将每个强非平稳分量序列通过变分模态分解进行进一步分解，得到强非平稳分量序列经变分模态分解后的平稳分量序列与非平稳分量序列组；

所述强非平稳分量序列通过VMD进行进一步分解得到的强非平稳分量序列经VMD分解所得序列组为：

简写

u_r＝[u_r(1),...,u_r(t),...,u_r(T)]^T,t∈[1,T]，u_r表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列；u_r(t)表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列第t时刻负荷；R＝80表示分解得到的序列分量数；T＝17520表示一个序列分量内包含样本的数量；

将load_sns＝[g₁,g₂,...,g_q,...,g_U1],q∈[1,U1]中g_q＝[g_q(1),g_q(2),...,g_q(T)]^T,q∈[1,U1]中每个强非平稳分量序列通过VMD进行进一步分解为：

所述VMD方法的原理为：寻找R个本征模态函数分量集合{u₁,...,u_R}，每个本征模态函数分量对应的中心频率为{ω₁,...,ω_R}，以R个本征模态函数分量之和等于待分解序列g_q为约束条件，使得分解后每个u_r,r∈[1,R]的估计带宽之和最小；

所述VMD方法的具体步骤如下：

构造约束变分最优问题：

式中，g_q(t)为load_sns中第q个序列第t时刻负荷，{u₁,...,u_R}为g_q分解后的模态集合，其中u₁(t)指VMD分解后的第1个模态的第t时刻负荷，u_R(t)指VMD分解后的第R个模态的第t时刻负荷；{ω₁,...,ω_R}为{u₁,...,u_R}中各序列分量对应的中心频率，其中ω₁为负荷序列[u₁(1),u₁(2),...,u₁(T)]^T的中心频率，ω_R为负荷序列[u_R(1),u_R(2),...,u_R(T)]^T的中心频率；

表示对t求偏导；δ(t)为Dirac分布；

*为卷积运算符；

表示求二范数的平方。

利用二次惩罚因子α和Lagrange乘子λ(t)构造Lagrange函数，将上述约束问题转化为无约束问题，转变后的扩展的Lagrange表达式为：

式中，L(·)表示Lagrange函数；λ(t)为Lagrange乘子，用来严格保持约束条件；α为二次惩罚因子，保证重构精度；＜·＞表示求内积。

利用交替方向乘子算法寻找最优解，不断更新u_r、ω_r，公式如下：

式中，m为迭代次数，m≤M，M为最大迭代次数；^表示傅里叶变换；

为当前信号的维纳滤波，为

的傅里叶变换，表示分解后第r个模态第m+1次迭代后的频域值；

为对应序列

的功率谱重心；

为对应g_q(t)经傅里叶变换后的频域值；

为对应u(t)经傅里叶变换后的频域值；

为对应λ(t)经傅里叶变换后的频域值；

所述将强非平稳分量序列通过VMD进一步分解得到的强非平稳分量序列经变分模态分解所得的分量序列组为：

简写

u_r＝[u_r(1),...,u_r(t),...,u_r(T)]^T,t∈[1,T]，u_r表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列；u_r(t)表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列第t时刻负荷；R＝80表示分解得到的序列分量数；T＝17520表示一个序列分量内包含样本的数量；

将强非平稳分量序列通过VMD进一步分解得到的分量序列组中每个序列分量通过近似熵法得到每个序列分量的平稳性强度：

以τ₁＝1作为平稳与非平稳的分界值，τ₂＝2作为非平稳与强非平稳的分界值，对每个序列的平稳性强度

进行分析。对任意r∈[1,R]，若κ_r＜τ₁则该序列分量归为平稳序列分量，若τ₁≤κ_r≤τ₂则该序列分量归为非平稳序列分量；

步骤1.9：构建原始负荷序列经过二次分解最终得到的平稳序列分量组、原始负荷序列经过二次分解最终得到的非平稳序列分量组；

将原始负荷序列经CEEMDAN分解得到的平稳序列组以及原始负荷序列经CEEMDAN分解得到的强平稳序列组再次经VMD进行分解后所得的强平稳序列组经VMD再次进行分解所得的平稳序列组进行组合，最终得到I＝89个平稳序列分量数，并将原始负荷序列经二次分解最终得到的平稳序列组定义为Load_s；

将原始负荷序列经CEEMDAN分解得到的非平稳序列组以及原始负荷序列经CEEMDAN分解得到的强平稳序列组再次经VMD进行分解后所得的强平稳序列组经VMD再次进行分解所得的非平稳序列组进行组合，最终得到J＝2个非平稳序列分量数，并将原始负荷序列经二次分解最终得到的非平稳序列组定义为Load_ns。其中，Load_ns中每个序列分量已通过近似熵法得到每个序列分量的平稳性强度，并根据平稳性强度从左到右由弱至强的顺序进行排序。

所述原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组Load_s为：

简写L_s：Load_s＝[A₁,...,A_α,...,A_I],α∈[1,I]，A_α＝[A_α(1),...,A_α(t),...,A_α(T)]^T,t∈[1,T]，A_α表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；A_α(t)表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I＝89表示经过二次分解后最终的平稳序列分量数；T＝17520表示每个序列分量包含的样本数。

所述原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组Load_ns为：

简写Load_ns：Load_ns＝[B₁,...,B_β,...,B_J],β∈[1,J]，B_β＝[B_β(1),...,B_β(t),...,B_β(T)]^T,t∈[1,T]，B_β表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；B_β(t)表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J＝2表示经过二次分解后最终的非平稳序列分量数；T＝17520表示每个序列分量包含的样本数。

由步骤1.1至步骤1.9，可以分别得到原始电、冷、热负荷序列分别经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的原始电、冷、热负荷序列各自对应的平稳、非平稳序列分量组为：

所述原始电负荷序列为：x^e(t),t∈[1,T]；

其中，x^e(t)为原始电负荷序列中t时刻的原始电负荷，T＝17520为原始电负荷的数量。

所述原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后，得到原始电负荷平稳序列分量组

为：

简写

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I＝89表示经过二次分解后最终的平稳序列分量数；T＝17520表示每个序列分量包含的样本数。

所述原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始电负荷非平稳序列分量组

为：

简写

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J＝2表示经过二次分解后最终的非平稳序列分量数；T＝17520表示每个序列分量包含的样本数。

所述原始冷负荷序列为：x^c(t),t∈[1,T]；

其中，x^c(t)为原始冷负荷序列中t时刻的原始冷负荷，T＝17520为原始冷负荷的数量。

所述原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始冷负荷平稳序列分量组

为：

简写

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I＝89表示经过二次分解后最终的平稳序列分量数；T＝17520表示每个序列分量包含的样本数。

所述原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始冷负荷非平稳序列分量组

为：

简写

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J＝2表示经过二次分解后最终的非平稳序列分量数；T＝17520表示每个序列分量包含的样本数。

所述原始热负荷序列为：x^h(t),t∈[1,T]；

其中，x^h(t)为原始热负荷序列中t时刻的原始热负荷，T＝17520为原始热负荷的数量。

所述原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始热负荷平稳序列分量组

为：

简写

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I＝89表示经过二次分解后最终的平稳序列分量数；T＝17520表示每个序列分量包含的样本数。

所述原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始热负荷非平稳序列分量组

为：

简写

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J＝2表示经过二次分解后最终的非平稳序列分量数；T＝17520表示每个序列分量包含的样本数。

步骤2：构建天气日历规则特征数据集，通过KPCA算法对天气日历规则特征集进行降维得到降维后天气日历规则主成分数据集；

步骤2所述天气日历规则特征数据集为：

data_h＝(month_h,day_h,week_h,hour_h,holiday_h,d_h,1,d_h,2,...,d_h,Q)

h∈[1,T]

其中，data_h为天气日历规则特征数据集中第h个样本，h表示第h个小时，T＝17520为天气日历规则特征数据集中样本的数量即小时的数量，month_h为天气日历规则特征数据集中第h个样本中月份序号标识，day_h为天气日历规则特征数据集中第h个样本中天序号标识，week_h为天气日历规则特征数据集中第h个样本中星期序号标识，hour_h为天气日历规则特征数据集中第h个样本中小时序号标识，holiday_h为天气日历规则特征数据集中第h个样本中节假日序号标识，

为天气日历规则特征数据集中第h个样本中第q^*种类型的天气因素，q^*∈[1,Q]，Q＝7为天气因素类型的数量；

步骤2所述通过KPCA算法对天气日历规则特征数据集进行降维得到降维后天气日历规则的主成分数据集为：

data_h中天气日历规则特征数据集中样本的数量即小时的数量T＝17520在降维中保持不变，对data_h中Q+5个元素即Q+5个维度通过KPCA算法进行降维；

步骤2所述降维后天气日历规则主成分数据集为：

data^* _h＝(W_h,1,W_h,2,...,W_h,H)

h∈[1,T]

其中，data^* _h为降维后天气日历规则特征数据集中第h个样本，h表示第h个小时，T＝17520为降维后天气日历规则特征数据集中样本的数量即小时的数量，W_h,p*为降维后天气日历规则因素特征数据集中第h个样本中第p^*个元素，H＝2为降维后天气日历规则特征数据集中每个样本的维度；

步骤3：将二次分解得到的原始电负荷平稳序列分量、二次分解得到的原始冷负荷平稳序列分量、二次分解得到的原始热负荷平稳序列分量、降维后天气日历规则特征数据集组成MLR输入数据集，通过优化训练MLR模型得到优化后MLR模型；

步骤3所述二次分解得到的原始电负荷平稳序列分量为：

其中，

表示二次分解得到的原始电负荷平稳序列分量组内第α个分量；

步骤3所述二次分解得到的原始冷负荷平稳序列分量为：

其中，

表示二次分解得到的原始冷负荷平稳序列分量组内第α个分量；

步骤3所述二次分解得到的原始热负荷平稳序列分量为：

其中，

表示二次分解得到的原始热负荷平稳序列分量组内第α个分量；

步骤3所述通过优化训练MLR模型得到优化后MLR模型为：

步骤3.1，建立MLR回归模型；

模型内待优化参数即各回归系数为：

其中，

表示第k^*个电负荷序列回归系数，

表示第k^*个冷负荷序列回归系数，

表示第k^*个热负荷序列回归系数，k^*＝6表示参与回归的变量共有k^*个，常数项、

的维度、

的维度、

的维度、降维后天气日历规则特征数据集中H＝2个主成分的维度构成k^*＝6；

步骤3.2，将二次分解得到的原始电负荷平稳序列分量、二次分解得到的原始冷负荷平稳序列分量、二次分解得到的原始热负荷平稳序列分量、降维后天气日历规则特征数据集组成输入集输入MLR进行回归；

步骤3.3，采用最小二乘估计求得电负荷序列回归系数、冷负荷序列回归系数、热负荷序列回归系数，用于构建优化后MLR回归模型；

步骤4：将二次分解得到的原始电负荷非平稳序列分量、二次分解得到的原始冷负荷非平稳序列分量、二次分解得到的原始热负荷非平稳序列分量、降维后天气日历规则特征数据集组成神经网络输入数据集，通过优化训练神经网络得到优化后神经网络模型；

步骤4所述二次分解得到的原始电负荷非平稳序列分量为：

其中，

为二次分解得到的原始电负荷非平稳序列分量组内第β个分量；

步骤4所述为二次分解得到的原始冷负荷非平稳序列分量：

其中，

为二次分解得到的原始冷负荷非平稳序列分量组内第β个分量；

步骤4所述二次分解得到的原始热负荷非平稳序列分量为：

其中，

为二次分解得到的原始热负荷非平稳序列分量组内第β个分量；

步骤4.1，所述神经网络为包含多个隐含层的DBiLSTM神经网络；

各隐含层神经元个数分别为：

θ＝2表示隐含层数量；

每个隐含层中待优化参数为：

第Γ个隐含层的遗忘门权重矩阵

与第Γ个隐含层的遗忘门偏置项

第Γ个隐含层的输入门第一个权重矩阵

与第Γ个隐含层的输入门第一个偏置项

第Γ个隐含层的输入门第二个权重矩阵

与第Γ个隐含层的输入门第二个偏置项

第Γ个隐含层的输出门权重矩阵

与第Γ个隐含层的输出门偏置项

Γ∈[1,θ]；

步骤4.2，将

以及降维后天气日历规则特征数据集组成神经网络输入数据集输入至神经网络进行训练；

步骤4.3，建立模型训练的损失函数，运用优化算法优化神经网络内部每个隐含层中待优化参数即

所述神经网络的损失函数表达式为：

式中，loss为损失函数值；S＝672为参加计算损失值的预测点个数；

为第ζ时刻电负荷的真实值；

为第ζ时刻冷负荷的真实值；

为第ζ时刻热负荷的真实值；

为第ζ时刻电负荷的预测值；

为第ζ时刻冷负荷的预测值；

为第ζ时刻热负荷的预测值。

步骤4.4，结合神经网络的损失函数，将将二次分解得到的原始电负荷非平稳序列分量、二次分解得到的原始冷负荷非平稳序列分量、二次分解得到的原始热负荷非平稳序列分量、降维后天气日历规则特征数据集组成神经网络输入数据集，通过Adam算法优化训练得到优化后神经网络模型；

步骤5：将原始电负荷平稳序列分量、原始冷负荷平稳序列分量、原始热负荷平稳序列分量与降维后天气日历规则特征数据集组成输入数据集，通过优化后MLR模型预测得到待预测的电负荷平稳序列分量、待预测的冷负荷平稳序列分量、待预测的热负荷平稳序列分量；将原始电负荷非平稳序列分量、原始冷负荷非平稳序列分量、原始热负荷非平稳序列分量与降维后天气日历规则特征数据集组成输入数据集，通过优化后神经网络模型预测得到待预测的电负荷非平稳序列分量、待预测的冷负荷非平稳序列分量、待预测的热负荷非平稳序列分量，将优化后MLR模型预测得到的所有待预测的负荷平稳序列分量与优化后神经网络模型预测得到的所有待预测的负荷非平稳序列分量相加得到最终的待预测电负荷值、待预测冷负荷值、待预测热负荷值。

实施例

本实施例以美国亚利桑那州立大学坦佩校区IES中电、冷、热多元负荷数据作为实验数据，数据信息如附图4所示，该IES属于用户级IES。其中，输入特征集包括当地天气与日历规则两种类型。天气数据考虑温度、风速、湿度、太阳垂直辐射、太阳水平辐射、露点和大气压，日历规则考虑月数、周数、天数、小时数和节假日，特征集数据共12维。选择2018年1月1日至2019年12月31日，按8：1：1划分训练集、验证集、测试集，以1h为步长进行预测。

首先运用QMD对多元负荷数据进行模态分解。分别对三种负荷进行CEEMDAN模态分解，得到的电、冷、热一次模态分解结果如附图5-附图7所示。电、冷、热负荷经过CEEMDAN分解后产生的IMF1分量均为强非平稳分量，因此运用VMD分别对其进行二次模态分解。考虑到分解产生的IMF分量越多，重构后与分解前的IMF1越相近，同时由于VMD二次模态分解产生的IMF分量均为平稳分量，可以由MLR进行预测，所需耗费的时间很短，因此将三种负荷的IMF1分量均由VMD再次分解为80个IMF分量，电、冷、热负荷二次模态分解结果见附图8-附图10。

本发明采用的特征集数据共12维，若与多元负荷数据共同输入模型就多达15维，增加了模型的复杂度。因此运用KPCA对特征集进行数据降维。表1为经KPCA提取的前四个主成分，可以看出，第一个主成分的方差贡献率为93.32％，前两个主成分的累积方差贡献率已高达98.30％，已经包含原来12维特征集的大部分信息，因此选择前两个主成分代替原来12维特征集，提取后的前两个主成分波形见附图11。

表1 KPCA提取主成分结果

模型的输入为预测时刻前6小时的电、冷、热负荷数据及特征集降维后的两个主成分共5维，输出为待预测时刻电、冷、热负荷。

DBiLSTM神经网络层数2层，隐含神经元分别为50和100个；学习率0.01；优化算法为Adam；迭代次数600代。

为了分析模态分解对多元负荷预测精度的影响，以不进行模态分解、进行CEEMDAN一次模态分解和进行CEEMDAN、VMD二次模态分解进行实验对比，结果如附图12和表2所示。

表2模态分解前后预测精度结果

可以看出，不进行模态分解的预测精度最差，电、冷、热负荷的预测误差E_MAPE分别为2.49％、4.68％、4.07％；一次模态分解的预测精度次优，E_MAPE分别为1.48％、2.36％、2.51％，精度提升了将近1倍；而进行二次模态分解的预测精度最好，E_MAPE分别为0.47％、1.05％、0.95％，精度提升了4倍左右，其中电负荷预测精度提升了5倍。由于冷负荷曲线相对比较平稳，三种方法都能得到较好的预测结果；而电、热负荷曲线随机、波动性强，直接预测产生的误差很大。虽然进行一次模态分解能够一定程度上提高预测精度，但分解产生的强非平稳分量无法进行准确预测，而这部分分量正是包含原始负荷曲线里的随机波动分量，因此一次模态分解只是提升了对负荷大致趋势的预测，对于负荷的小范围波动则无法准确预测；而进行二次模态分解后，VMD将CEEMDAN分解产生的强非平稳分量再次分解为多个平稳分量，大大提升了这部分强非平稳分量的预测精度，使最终的预测精度有了大幅度的提升。

为验证本发明所提出的模型能够有效提高多元负荷的预测精度，将本发明所提的预测模型(简称KQDM)与几种常见组合模型进行对比，对照模型为：反向传播神经网络组合高斯过程回归(简称BP-GPR)；KPCA组合深度置信网络(简称KPCA-DBN)；LSTM神经网络组合MLR(简称LSTM-MLR)；集合经验模态分解组合门控循环单元神经网络和MLR(简称EEMD-GRU-MLR)，实验结果如附图13和表3所示。

表3不同模型预测精度结果

结合实验结果可知，本发明所提出的预测模型在预测精度上远远高于其他四种模型，电、冷、热负荷的预测误差E_MAPE仅为0.47％、1.05％、0.95％。其中，EEMD-GRU-MLR组合模型表现最差，E_MAPE分别达到了3.36％、10.47％、4.46％，比本发明的模型分别高出7倍、10倍、4.5倍；浅层神经网络结合统计学方法的BP-GPR在面对随机波动性大的时间序列表现也不佳；另外，近年来比较流行的深度置信网络结合KPCA相较于上述两种模型虽然预测精度有所提高，但E_MAPE也比本发明模型分别高出5倍、6倍、4倍左右；而相较于上述三种模型而言，LSTM-MLR在预测精度上有了进一步提高，但还达不到本发明所能达到的预测精度。

对上述结果进一步分析：

1)EEMD-GRU-MLR组合模型虽然采用了门控循环单元神经网络，能够在时间序列预测上取得较好的结果，但由于集合经验模态分解法分解时，引入了白噪声叠加在原始序列上，导致所有IMF重构后不等于原始序列，预测精度受到较大的影响；

2)BP-GPR组合模型中，反向传播神经网络在训练时容易陷入局部最优，并且容易出现梯度消失等问题，而高斯过程回归也只适用于预测平稳序列，对此类非平稳序列难以准确预测；

3)KPCA-DBN组合模型中，深度置信网络是由底层多层受限玻尔茨曼基和顶层反向传播神经网络组成，结构内部缺乏类似于LSTM神经网络的记忆结构，性能上略差于LSTM神经网络；

4)LSTM-MLR组合模型采用了LSTM神经网络，在预测精度上有了很大的提高，但由于只对历史数据进行单向训练，没有学习到历史数据的全局信息，忽略了前后时间的关联性，无法训练得到精确的预测模型，同时由于用户级IES多元负荷随机、波动性强，运用LSTM-MLR组合模型直接预测也无法得到较高的预测精度。

而相对上述四种模型而言，本发明所提的KQDM模型，运用DBiLSTM神经网络，通过对历史数据进行双向训练，更加详细地学到历史数据包含的全局信息，并结合MLR在平稳序列预测上的优势，在前期KPCA和QMD数据预处理下取得很好的预测结果。

以上所述仅为本发明的较佳实施例，凡依本发明申请专利范围所做的均等变化与修饰，皆应属本发明的涵盖范围。

进一步的，本发明对预测结果的评价指标为：

1)平均绝对百分误差E_MAPE：

2)均方根误差E_RMSE：

式中，N＝168为参加计算误差值的预测点个数；Ψ(t)为t时刻负荷真实值；

为t时刻负荷预测值。

应当理解的是，上述针对较佳实施例的描述较为详细，并不能因此而认为是对本发明专利保护范围的限制，本领域的普通技术人员在本发明的启示下，在不脱离本发明权利要求所保护的范围情况下，还可以做出替换或变形，均落入本发明的保护范围之内，本发明的请求保护范围应以所附权利要求为准。

Claims (5)

1.一种用户级综合能源系统多元负荷预测方法，其特征在于，包括如下步骤：

步骤1：运用自适应噪声的完全集合经验模态分解分别对电负荷、冷负荷、热负荷进行本征模态分解，对分解得到的强非平稳分量运用变分模态分解再次进行分解；

步骤2：构建天气日历规则特征数据集，通过KPCA算法对天气日历规则特征集进行降维得到降维后天气日历规则主成分数据集；

步骤3：将二次分解得到的原始电负荷平稳序列分量、二次分解得到的原始冷负荷平稳序列分量、二次分解得到的原始热负荷平稳序列分量、降维后天气日历规则特征数据集组成MLR输入数据集，通过优化训练MLR模型得到优化后MLR模型；

步骤4：将二次分解得到的原始电负荷非平稳序列分量、二次分解得到的原始冷负荷非平稳序列分量、二次分解得到的原始热负荷非平稳序列分量、降维后天气日历规则特征数据集组成神经网络输入数据集，通过优化训练神经网络得到优化后神经网络模型；

步骤5：将原始电负荷平稳序列分量、原始冷负荷平稳序列分量、原始热负荷平稳序列分量与降维后天气日历规则特征数据集组成输入数据集，通过优化后MLR模型预测得到待预测的电负荷平稳序列分量、待预测的冷负荷平稳序列分量、待预测的热负荷平稳序列分量；将原始电负荷非平稳序列分量、原始冷负荷非平稳序列分量、原始热负荷非平稳序列分量与降维后天气日历规则特征数据集组成输入数据集，通过优化后神经网络模型预测得到待预测的电负荷非平稳序列分量、待预测的冷负荷非平稳序列分量、待预测的热负荷非平稳序列分量，将优化后MLR模型预测得到的所有待预测的负荷平稳序列分量与优化后神经网络模型预测得到的所有待预测的负荷非平稳序列分量相加得到最终的待预测电负荷值、待预测冷负荷值、待预测热负荷值。

2.根据权利要求1所述的用户级综合能源系统多元负荷预测方法，其特征在于：

所述步骤1：

步骤1.1：引入原始负荷序列；

所述原始负荷序列为：x(t),t∈[1,T]；

其中，x(t)为原始负荷序列中t时刻的原始负荷，T为原始负荷的数量；

步骤1.2：向原始负荷序列加入正白噪声序列得到加入正白噪声后的负荷序列；

所述加入正白噪声后的负荷序列为：

m⁺(t)＝x(t)+n⁺(t),t∈[1,T]

其中，n⁺(t)为正白噪声序列中t时刻引入的正白噪声，T为原始负荷的数量即正白噪声的数量，m⁺(t)为加入正白噪声后的负荷序列第t时刻的负荷；

所述正白噪声序列为：n⁺(t),t∈[1,T]；

所述加入正白噪声后的负荷序列为：m⁺(t),t∈[1,T]；

步骤1.3：向原始负荷序列加入负白噪声序列得到加入负白噪声后的负荷序列；

所述加入负白噪声后的负荷序列为：

m^-(t)＝x(t)+n^-(t),t∈[1,T]

式中，n^-(t)为负白噪声序列中t时刻引入的负白噪声，T为原始负荷的数量即负白噪声的数量，m^-(t)为加入负白噪声后的负荷序列第t时刻的负荷；

所述负白噪声序列为：n^-(t),t∈[1,T]；

所述加入负白噪声后的负荷序列为：m^-(t),t∈[1,T]；

步骤1.4：运用经验模态分解法对加入正白噪声后的负荷序列进行模态分解得到经验模态分解后对应的一组本征模态负荷序列组F⁺(u,t),(u∈[1,U],t∈[1,T])，运用经验模态分解法对加入负白噪声后的负荷序列进行模态分解得到经验模态分解后对应的一组本征模态负荷序列组F^-(v,t),(v∈[1,V],t∈[1,T])，其中F⁺(u,t)表示加入正白噪声后的负荷序列进行经验模态分解后第u个序列分量第t时刻的负荷，F^-(v,t)表示加入负白噪声后的负荷序列进行经验模态分解后第v个序列分量第t时刻的负荷，U和V均表示分解后的序列分量数且U＝V；

所述加入正白噪声后的负荷序列为：[m⁺(1),m⁺(2),...,m⁺(T)]^T；

所述加入负白噪声后的负荷序列为：[m^-(1),m^-(2),...,m^-(T)]^T；

所述加入正白噪声负荷序列经验模态分解后的一组本征模态负荷序列中第u个序列分量为：[F⁺(u,1),F⁺(u,2),...,F⁺(u,T)]^T,u∈[1,U]；

所述加入负白噪声负荷序列经验模态分解后的一组本征模态负荷序列中第v个序列分量为：[F^-(v,1),F^-(v,2),...,F^-(v,T)]^T,v∈[1,V]；

步骤1.5：重复执行步骤1.1至步骤1.4多次，将多次得到的经验模态分解后的正白噪声负荷序列对应的本征模态负荷序列组通过均值计算得到经验模态分解后正白噪声负荷序列对应的一组本征模态平均负荷序列组，将多次得到的经验模态分解后的负白噪声负荷序列对应的本征模态负荷序列组通过均值计算得到经验模态分解后负白噪声负荷序列对应的一组本征模态平均负荷序列组；

所述通过均值计算得到经验模态分解后的正白噪声对应的一组本征模态平均负荷序列中第u个序列分量第t时刻平均负荷为：

式中，F_i ⁺(u,t)为第i次经验模态分解后的正白噪声负荷序列对应的本征模态负荷序列组中第u个序列分量第t时刻的负荷，

为经验模态分解后的正白噪声负荷序列对应的一组本征模态平均负荷序列中第u个序列分量第t时刻平均负荷，N为重复执行的次数；

所述经验模态分解后的正白噪声负荷序列对应的一组本征模态平均负荷序列组为：

所述通过均值计算得到经验模态分解后的负白噪声对应的一组本征模态平均负荷序列中第v个序列分量第t时刻平均负荷为：

式中，F_i ^-(v,t)为第i次经验模态分解后的负白噪声负荷序列对应的本征模态负荷序列组中第v个序列分量第t时刻的负荷，

为经验模态分解后的负白噪声负荷序列对应的一组本征模态平均负荷序列中第v个序列分量第t时刻平均负荷，N为重复执行的次数；

所述经验模态分解后的负白噪声负荷序列对应的一组本征模态平均负荷序列组为：

步骤1.6：将经验模态分解后的正白噪声负荷序列对应的一组本征模态平均负荷序列组、经验模态分解后的负白噪声负荷序列对应的一组本征模态平均负荷序列组通过均值计算得到最终一组负荷分解序列组；

所述通过均值计算得到最终负荷分解序列组为：

式中，load(k,t)为最终负荷分解序列组中第k个序列分量第t时刻负荷；

步骤1.7：将最终负荷分解序列组中每个序列分量通过近似熵法得到每个序列分量的平稳性强度，根据平稳性强度从左到右由弱至强的顺序，对最终负荷分解序列组中每个序列分量进行排序，得到平稳性强度排序后序列分量；

所述最终负荷分解序列组中每个序列分量为：

其中，load表示原始负荷序列经CEEMDAN最终分解的序列组，G_k表示最终分解序列组中第k个本征模态分量序列，G_k(t)表示最终分解序列组中第k个本征模态分量序列第t时刻负荷，U表示分解后的序列分量数，T表示一个序列分量内包含样本的数量；

所述每个序列分量的平稳性强度为：

load_stat＝[ξ₁,...,ξ_k,...,ξ_U],k∈[1,U]

式中，load_stat为每个序列分量的平稳性强度集；ξ_k,k∈[1,U]表示第k个序列分量的平稳性强度，U表示分解后的序列分量数；

所述平稳性强度排序后序列分量为：

其中，load_SE表示平稳性强度排序后序列分量，g_γ表示平稳性强度排序后序列分量中第γ个本征模态分量序列，g_γ(t)表示平稳性强度排序后序列分量中第γ个本征模态分量序列第t时刻负荷，U表示分解后的序列分量数，T表示一个序列分量内包含样本的数量；

load_SE进一步简化为：

load_SE＝[g₁,...,g_γ,...,g_U],γ∈[1,U]，g_γ＝[g_γ(1),...,g_γ(t),...,g_γ(T)]^T,t∈[1,T]

所述平稳性强度排序后序列分量组load_SE对应的每个序列分量的平稳性强度更新为：

式中，

为平稳性强度排序后序列分量组内每个序列分量的平稳性强度集；χ_γ,γ∈[1,U]表示第γ个序列分量的平稳性强度，U表示分解后的序列分量数；

以τ₁作为平稳与非平稳的分界值，τ₂作为非平稳与强非平稳的分界值，对每个序列的平稳性强度

进行分析；

对任意γ∈[1,U]，若χ_γ＜τ₁则该序列分量归为平稳序列分量，若τ₁≤χ_γ≤τ₂则该序列分量归为非平稳序列分量，若χ_γ＞τ₂则该序列分量归为强非平稳序列分量；

U1个强非平稳分量序列构建的强非平稳分量集合即：

load_sns＝[g₁,g₂,...,g_q,...,g_U1],q∈[1,U1]；

U2个非平稳分量序列构建的非平稳分量序列集合即：

load_ns＝[g_U1+1,g_U1+2,...,g_U1+z,...,g_U1+U2],z∈[1,U2]

U3个平稳分量序列构建的平稳分量序列集合即：

load_s＝[g_U1+U2+1,g_U1+U2+2,...,g_U1+U2+b,...,g_U1+U2+U3],b∈[1,U3]

其中：

U1+U2+U3＝U

g_q＝[g_q(1),g_q(2),...,g_q(T)]^T,q∈[1,U1]

g_U1+z＝[g_U1+z(1),g_U1+z(2),...,g_U1+z(T)]^T,z∈[1,U2]

g_U1+U2+b＝[g_U1+U2+b(1),g_U1+U2+b(2),...,g_U1+U2+b(T)]^T,b∈[1,U3]

步骤1.8：将每个强非平稳分量序列通过变分模态分解进行进一步分解，得到强非平稳分量序列经变分模态分解后的平稳分量序列与非平稳分量序列组；

所述强非平稳分量序列通过VMD进行进一步分解得到的强非平稳分量序列经VMD分解所得序列组为：

简写

u_r＝[u_r(1),...,u_r(t),...,u_r(T)]^T,t∈[1,T]，u_r表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列；u_r(t)表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列第t时刻负荷；R表示分解得到的序列分量数；T表示一个序列分量内包含样本的数量；

将load_sns＝[g₁,g₂,...,g_q,...,g_U1],q∈[1,U1]中g_q＝[g_q(1),g_q(2),...,g_q(T)]^T,q∈[1,U1]中每个强非平稳分量序列通过VMD进行进一步分解为：

所述VMD方法的原理为：寻找R个本征模态函数分量集合{u₁,...,u_R}，每个本征模态函数分量对应的中心频率为{ω₁,...,ω_R}，以R个本征模态函数分量之和等于待分解序列g_q为约束条件，使得分解后每个u_r,r∈[1,R]的估计带宽之和最小；

所述VMD方法的具体步骤如下：

构造约束变分最优问题：

式中，g_q(t)为load_sns中第q个序列第t时刻负荷，{u₁,...,u_R}为g_q分解后的模态集合，其中u₁(t)指VMD分解后的第1个模态的第t时刻负荷，u_R(t)指VMD分解后的第R个模态的第t时刻负荷；{ω₁,...,ω_R}为{u₁,...,u_R}中各序列分量对应的中心频率，其中ω₁为负荷序列[u₁(1),u₁(2),...,u₁(T)]^T的中心频率，ω_R为负荷序列[u_R(1),u_R(2),...,u_R(T)]^T的中心频率；

表示对t求偏导；δ(t)为Dirac分布；

*为卷积运算符；

表示求二范数的平方；

利用二次惩罚因子α和Lagrange乘子λ(t)构造Lagrange函数，将上述约束问题转化为无约束问题，转变后的扩展的Lagrange表达式为：

式中，L(·)表示Lagrange函数；λ(t)为Lagrange乘子，用来严格保持约束条件；α为二次惩罚因子，保证重构精度；＜·＞表示求内积；

利用交替方向乘子算法寻找最优解，不断更新u_r、ω_r，公式如下：

式中，m为迭代次数，m≤M，M为最大迭代次数；^表示傅里叶变换；

为当前信号的维纳滤波，为

的傅里叶变换，表示分解后第r个模态第m+1次迭代后的频域值；

为对应序列

的功率谱重心；

为对应g_q(t)经傅里叶变换后的频域值；

为对应u(t)经傅里叶变换后的频域值；

为对应λ(t)经傅里叶变换后的频域值；

所述将强非平稳分量序列通过VMD进一步分解得到的强非平稳分量序列经变分模态分解所得的分量序列组为：

简写

u_r＝[u_r(1),...,u_r(t),...,u_r(T)]^T,t∈[1,T]，u_r表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列；u_r(t)表示强非平稳分量序列通过VMD进一步分解后序列分量中第r个本征模态分量序列第t时刻负荷；R表示分解得到的序列分量数；T表示一个序列分量内包含样本的数量；

将强非平稳分量序列通过VMD进一步分解得到的分量序列组中每个序列分量通过近似熵法得到每个序列分量的平稳性强度：

以τ₁作为平稳与非平稳的分界值，τ₂作为非平稳与强非平稳的分界值，对每个序列的平稳性强度

进行分析；对任意r∈[1,R]，若κ_r＜τ₁则该序列分量归为平稳序列分量，若τ₁≤κ_r≤τ₂则该序列分量归为非平稳序列分量；

步骤1.9：构建原始负荷序列经过二次分解最终得到的平稳序列分量组、原始负荷序列经过二次分解最终得到的非平稳序列分量组；

将原始负荷序列经CEEMDAN分解得到的平稳序列组以及原始负荷序列经CEEMDAN分解得到的强平稳序列组再次经VMD进行分解后所得的强平稳序列组经VMD再次进行分解所得的平稳序列组进行组合，最终得到I个平稳序列分量数，并将原始负荷序列经二次分解最终得到的平稳序列组定义为Load_s；

将原始负荷序列经CEEMDAN分解得到的非平稳序列组以及原始负荷序列经CEEMDAN分解得到的强平稳序列组再次经VMD进行分解后所得的强平稳序列组经VMD再次进行分解所得的非平稳序列组进行组合，最终得到J个非平稳序列分量数，并将原始负荷序列经二次分解最终得到的非平稳序列组定义为Load_ns；其中，Load_ns中每个序列分量已通过近似熵法得到每个序列分量的平稳性强度，并根据平稳性强度从左到右由弱至强的顺序进行排序；

所述原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组Load_s为：

简写L_s：Load_s＝[A₁,...,A_α,...,A_I],α∈[1,I]，A_α＝[A_α(1),...,A_α(t),...,A_α(T)]^T,t∈[1,T]，A_α表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；A_α(t)表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I表示经过二次分解后最终的平稳序列分量数；T表示每个序列分量包含的样本数；

所述原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组Load_ns为：

简写Load_ns：Load_ns＝[B₁,...,B_β,...,B_J],β∈[1,J]，B_β＝[B_β(1),...,B_β(t),...,B_β(T)]^T,t∈[1,T]，B_β表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；B_β(t)表示原始负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J表示经过二次分解后最终的非平稳序列分量数；T表示每个序列分量包含的样本数；

由步骤1.1至步骤1.9，可以分别得到原始电、冷、热负荷序列分别经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的原始电、冷、热负荷序列各自对应的平稳、非平稳序列分量组为：

所述原始电负荷序列为：x^e(t),t∈[1,T]；

其中，x^e(t)为原始电负荷序列中t时刻的原始电负荷，T为原始电负荷的数量；

所述原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后，得到原始电负荷平稳序列分量组

为：

简写

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I表示经过二次分解后最终的平稳序列分量数；T表示每个序列分量包含的样本数；

所述原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始电负荷非平稳序列分量组

为：

简写

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；

表示原始电负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J表示经过二次分解后最终的非平稳序列分量数；T表示每个序列分量包含的样本数；

所述原始冷负荷序列为：x^c(t),t∈[1,T]；

其中，x^c(t)为原始冷负荷序列中t时刻的原始冷负荷，T为原始冷负荷的数量；

所述原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始冷负荷平稳序列分量组

为：

简写

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I表示经过二次分解后最终的平稳序列分量数；T表示每个序列分量包含的样本数；

所述原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始冷负荷非平稳序列分量组

为：

简写

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；

表示原始冷负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J表示经过二次分解后最终的非平稳序列分量数；T表示每个序列分量包含的样本数；

所述原始热负荷序列为：x^h(t),t∈[1,T]；

其中，x^h(t)为原始热负荷序列中t时刻的原始热负荷，T为原始热负荷的数量；

所述原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始热负荷平稳序列分量组

为：

简写

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列；

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的平稳序列分量组中第α个本征模态分量序列第t时刻负荷；I表示经过二次分解后最终的平稳序列分量数；T表示每个序列分量包含的样本数；

所述原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到原始热负荷非平稳序列分量组

为：

简写

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列；

表示原始热负荷序列经过二次分解并按照平稳性强度从左到右由弱至强的顺序进行排序后最终得到的非平稳序列分量组中第β个本征模态分量序列第t时刻负荷；J表示经过二次分解后最终的非平稳序列分量数；T表示每个序列分量包含的样本数。

3.根据权利要求1所述的用户级综合能源系统多元负荷预测方法，其特征在于：

步骤2所述天气日历规则特征数据集为：

data_h＝(month_h,day_h,week_h,hour_h,holiday_h,d_h,1,d_h,2,...,d_h,Q)

h∈[1,T]

其中，data_h为天气日历规则特征数据集中第h个样本，h表示第h个小时，T为天气日历规则特征数据集中样本的数量即小时的数量，month_h为天气日历规则特征数据集中第h个样本中月份序号标识，day_h为天气日历规则特征数据集中第h个样本中天序号标识，week_h为天气日历规则特征数据集中第h个样本中星期序号标识，hour_h为天气日历规则特征数据集中第h个样本中小时序号标识，holiday_h为天气日历规则特征数据集中第h个样本中节假日序号标识，

为天气日历规则特征数据集中第h个样本中第q^*种类型的天气因素，q^*∈[1,Q]，Q为天气因素类型的数量；

步骤2所述通过KPCA算法对天气日历规则特征数据集进行降维得到降维后天气日历规则的主成分数据集为：

data_h中天气日历规则特征数据集中样本的数量即小时的数量T在降维中保持不变，对data_h中Q+5个元素即Q+5个维度通过KPCA算法进行降维；

步骤2所述降维后天气日历规则主成分数据集为：

data^* _h＝(W_h,1,W_h,2,...,W_h,H)

h∈[1,T]

其中，data^* _h为降维后天气日历规则特征数据集中第h个样本，h表示第h个小时，T为降维后天气日历规则特征数据集中样本的数量即小时的数量，

为降维后天气日历规则因素特征数据集中第h个样本中第p^*个元素，H为降维后天气日历规则特征数据集中每个样本的维度。

4.根据权利要求1所述的用户级综合能源系统多元负荷预测方法，其特征在于：

步骤3所述二次分解得到的原始电负荷平稳序列分量为：

其中，

表示二次分解得到的原始电负荷平稳序列分量组内第α个分量；

步骤3所述二次分解得到的原始冷负荷平稳序列分量为：

其中，

表示二次分解得到的原始冷负荷平稳序列分量组内第α个分量；

步骤3所述二次分解得到的原始热负荷平稳序列分量为：

其中，

表示二次分解得到的原始热负荷平稳序列分量组内第α个分量；

步骤3所述通过优化训练MLR模型得到优化后MLR模型为：

步骤3.1，建立MLR回归模型；

模型内待优化参数即各回归系数为：

其中，

表示第k^*个电负荷序列回归系数，

表示第k^*个冷负荷序列回归系数，

表示第k^*个热负荷序列回归系数，k^*表示参与回归的变量共有k^*个，常数项、

的维度、

的维度、

的维度、降维后天气日历规则特征数据集中H个主成分的维度构成k^*；

步骤3.2，将二次分解得到的原始电负荷平稳序列分量、二次分解得到的原始冷负荷平稳序列分量、二次分解得到的原始热负荷平稳序列分量、降维后天气日历规则特征数据集组成输入集输入MLR进行回归；

步骤3.3，采用最小二乘估计求得电负荷序列回归系数、冷负荷序列回归系数、热负荷序列回归系数，用于构建优化后MLR回归模型。

5.根据权利要求1所述的用户级综合能源系统多元负荷预测方法，其特征在于：

步骤4所述二次分解得到的原始电负荷非平稳序列分量为：

其中，

为二次分解得到的原始电负荷非平稳序列分量组内第β个分量；

步骤4所述为二次分解得到的原始冷负荷非平稳序列分量：

其中，

为二次分解得到的原始冷负荷非平稳序列分量组内第β个分量；

步骤4所述二次分解得到的原始热负荷非平稳序列分量为：

其中，

为二次分解得到的原始热负荷非平稳序列分量组内第β个分量；

步骤4.1，所述神经网络为包含多个隐含层的DBiLSTM神经网络；

各隐含层神经元个数分别为：

θ表示隐含层数量；

每个隐含层中待优化参数为：

第Γ个隐含层的遗忘门权重矩阵即

与第Γ个隐含层的遗忘门偏置项即

第Γ个隐含层的输入门第一个权重矩阵

与第Γ个隐含层的输入门第一个偏置项

第Γ个隐含层的输入门第二个权重矩阵

与第Γ个隐含层的输入门第二个偏置项

第Γ个隐含层的输出门权重矩阵

与第Γ个隐含层的输出门偏置项

Γ∈[1,θ]；

步骤4.2，将

以及降维后天气日历规则特征数据集组成神经网络输入数据集输入至神经网络进行训练；

步骤4.3，建立模型训练的损失函数，运用优化算法优化神经网络内部每个隐含层中待优化参数即

所述神经网络的损失函数表达式为：

式中，loss为损失函数值；S为参加计算损失值的预测点个数；

为第ζ时刻电负荷的真实值；

为第ζ时刻冷负荷的真实值；

为第ζ时刻热负荷的真实值；

为第ζ时刻电负荷的预测值；

为第ζ时刻冷负荷的预测值；

为第ζ时刻热负荷的预测值；

步骤4.4，结合神经网络的损失函数，将将二次分解得到的原始电负荷非平稳序列分量、二次分解得到的原始冷负荷非平稳序列分量、二次分解得到的原始热负荷非平稳序列分量、降维后天气日历规则特征数据集组成神经网络输入数据集，通过Adam算法优化训练得到优化后神经网络模型。

Images