基于补充泄漏积分回声状态网络的空气质量预测方法

1.本发明属于时间序列预测领域，涉及一种基于补充泄漏积分回声状态网络的空气质量预测方法。


背景技术：

2.随着城市化和工业化的快速推进，导致空气质量严重恶化，这对生活环境的质量产生了负面影响。特别是工业污染产生的可吸入颗粒物粒径小、扩散面积大、活动性强，可通过呼吸道进入人体，对人体健康产生不利影响。因此，大气污染物的预测对于环境污染的预警和控制具有至关重要的作用。对空气质量指数(如pm2.5浓度)进行建模和预测已成为防治空气污染的有效途径，也为制定有效措施提供了科学依据。这一理念的实施可以有效降低空气污染对健康的危害，从而实现预警和合理规划。
3.由于空气污染的成因非常复杂，分析空气质量指数的主要污染物和影响变量将为建立预测模型奠定基础。比如zhang等人在论文“zhang c,ni z,ni l.multifractal detrended cross-correlation analysis between pm2.5 and meteorological factors[j].physica a.2015,438:114
–
123.”利用多重分形去趋势互相关分析方法研究了pm2.5与4个气象因子之间的相互关系，揭示了气象变量对pm2.5浓度的影响。chen等人在论文“chen z,xie x,cai j,chen d,gao b,he b,cheng n,xu b.understanding meteorological influences on pm2.5 concentrations across china:a temporal and spatial perspective[j].atmospheric chemistry and physics.2018,18:5343
–
5358.”中利用收敛交叉映射方法分析了京津冀地区气象因子与pm2.5浓度之间的因果关系，得到了定量的因果分析结果。然而，在以往的研究中，大多数预测模型的建立并没有考虑到影响因素与pm2.5之间的相关性，空气污染的原因是非常复杂的。以pm2.5浓度为例，其浓度不仅受no2、co、o3、so2等空气污染物的影响，还受温度、气压、湿度、风速、风向等气象变量的影响，若能从各种影响变量中筛选出有效信息并加以利用有助于pm2.5浓度的准确预测，能够为未来决策提供辅助。


技术实现要素：

[0004]
本发明要解决的技术问题是，针对上述现有技术未能消除不相关和冗余变量的问题，提出一种基于补充泄漏积分回声状态网络的空气质量预测方法。
[0005]
为了达到上述目的，本发明采用的技术方案如下：
[0006]
一种基于补充泄露积分回声状态网络的空气质量预测方法。首先，收集空气质量序列，采用最小冗余最大相关性方法对空气质量序列进行特征选择，确定输入特征的排序结果，进而选择出最优的输入特征子集；其次，对空气质量序列中最优的输入特征子集进行相空间重构，旨在充分挖掘输入特征子集中有用信息，形成新的特征输入集；再次，对新的特征输入集按照比例分成训练集和测试集，在训练集上使用岭回归算法训练补充漏积分器回声状态网络模型，训练完毕得到最佳预测模型，在测试集上进行pm2.5时间序列进行预测
输出；最后，利用评价指标度量模型的预测精度。具体步骤如下：
[0007]
步骤1：从现实世界中收集空气质量序列t，包括pm2.5、pm10、no2及so2等t维特征，即t＝[x1,x2,
…
,x
t
]，也被称为多维空气质量序列；之后采用最小冗余最大相关性方法对空气质量序列中的输入特征进行特征选择，得出输入特征的排序结果，进而选择出最优的输入特征子集。具体如下：
[0008]
首先，计算最大相关度d和最小冗余度r：
[0009][0010][0011]
其中，xi和xk分别表示空气质量序列t中的第i个和第k个特征；c表示目标特征，即pm2.5；i(xi；c)表示第i个特征和目标特征c的相关函数；i(xi；xk)表示第i个特征和第k个特征的相关函数。相关函数i(
·
；
·
)通过以下公式计算：
[0012][0013]
其中，u和u
′
表示两个变量；p(
·
,
·
)表示概率。
[0014]
联立公式(1)-(3)，确定最小冗余最大相关性目标函数，定义如下：
[0015][0016]
之后，采用前向选择算法求解上述目标函数，即可得到输入特征的排序结果，以上过程称为最小冗余最大相关性方法。最后，通过交叉验证得到空气质量序列的最优t
′
维输入特征子集，其中t
′
＜t。
[0017]
步骤2：对步骤1中获得的最优t
′
维输入特征子集进行相空间重构，形成新的输入特征集作为空气质量预测模型的输入。具体如下：
[0018]
以空气质量序列中最优t
′
维输入特征子集的第j个特征xj为例，采用坐标延迟法对其进行相空间重构，计算公式如下：
[0019]
xj(n)＝[xj(n),xj(n+τj),
…
,xj(n+(m
j-1)τj)]
ꢀꢀꢀꢀ
(5)
[0020]
其中，xj(n)为第n时刻第j个特征的输入特征集；mj和τj分别为第j个特征的嵌入维数和延迟时间。同理，最优t
′
维输入特征子集的其他特征均采用如上方法。最终，形成新的高维输入特征集u，即u＝[u(1),u(2),
…
,u(n)]，其中n为总样本数或总时刻；在第n时刻，u(n)＝[x1(n),x2(n),
…
,x
t
′
(n)]。
[0021]
步骤3：将步骤2形成的高维输入特征集u按照比例划分成训练集和测试集。在训练集上训练所提补充漏积分器回声状态网络模型，该模型考虑前一状态对当前状态的影响，采用改进的储备池状态更新形式进行训练，从而得到空气质量预测最佳模型；之后，在测试集上应用最佳模型，对pm2.5时间序列进行预测输出。具体如下：
[0022]
首先，将高维输入特征集u的前75％作为训练集，后25％作为测试集。
[0023]
其次，在训练集上训练所提出的补充泄漏积分回声状态网络模型，确定空气质量预测最佳模型。所提模型的状态更新公式如下所示：
[0024]
s(n+1)＝a
×
tanh(w
in
×
[1；u(n+1)]+w
×
s(n))+b
×
s(n)+(1-a-b)
×
s(n-1)
ꢀꢀ
(6)
[0025]
其中，s(n+1)和s(n-1)分别表示回声状态网络中储备池第n+1和n-1时刻的状态，a和b是衰减参数，且a和b均大于0，w
in
表示输入层与储备池之间的权重，w表示储备池的权重。
[0026]
最后，基于空气质量预测最佳模型，在测试集上计算pm2.5时间序列的预测输出如下所示：
[0027][0028]
其中，w
out
为储备池与输出层的权重，采用岭回归对其求解。具体言之，岭回归是一种用于共线性数据分析的有偏估计回归方法，能够克服最小二乘解中的不适定问题，避免过拟合问题，其表达式如下:
[0029]wout
＝(s
t
s+ki)-1styꢀꢀꢀꢀ
(8)
[0030]
其中，k是正则化参数，i是单位矩阵，s和y分别是由s(n)和y(n)组成的矩阵。
[0031]
步骤4：利用评价指标度量所提模型的预测精度。
[0032]
进一步的，所述步骤4采用5个预测指标评价pm2.5时间序列的预测结果，即均方根误差(rmse)、归一化均方根误差(nrmse)、平均绝对误差(mae)、对称平均绝对百分比误差(smape)和皮尔逊相关系数(r)。
[0033]
与现有技术相比，本发明具有以下有益效果：
[0034]
(1)本发明从现实角度出发，为处理空气质量序列中的多维特征，采用最小冗余最大相关性方法确定空气质量序列中最优的输入特征子集，在保证最大相关性的前提下去除冗余特征，既能选择最相关的信息，又能减少冗余信息；
[0035]
(2)本发明利用相空间重构对最优的输入特征子集进行处理，有利于充分提取空气质量序列中的有用信息。
[0036]
(3)本发明提出补充泄漏积分回声状态网络模型，该模型考虑了前一状态对当前状态的影响，对储备池的状态更新公式进行了改进，使其在面对空气质量序列强烈的复杂性时仍可获得较高的预测性能。
附图说明
[0037]
图1为本发明空气质量预测模型流程图。
[0038]
图2(a)为不同因素与pm2.5的相关性。
[0039]
图2(b)为最优子集选择结果。
[0040]
图3(a)为本发明对pm2.5浓度一步(1小时)的预测曲线。
[0041]
图3(b)为本发明对pm2.5浓度一步(1小时)的误差曲线。
[0042]
图4为本发明对pm2.5浓度一步(1小时)的散点图。
[0043]
图5(a)为本发明对pm2.5浓度五步(5小时)的预测曲线。
[0044]
图5(b)为本发明对pm2.5浓度五步(5小时)的误差曲线。
[0045]
图6(a)为本发明对pm2.5浓度十步(10小时)预测曲线。
[0046]
图6(b)为本发明对pm2.5浓度十步(10小时)的误差曲线。
具体实施方式
[0047]
以下将结合附图及具体实施例对本发明做进一步说明。
[0048]
本发明所使用硬件设备包括pc机器一台。
[0049]
如图1所示，基于补充泄漏积分回声状态网络的空气质量预测方法，具体实施措施如下：
[0050]
步骤1：获取北京市海淀区2016年1月至12月的逐时空气质量序列(数据集)，共计8759组样本。数据集来自美国大使馆(哈佛大学地理分析中心数据集)，包括每小时pm2.5、pm10、no2、co、o3和so2的平均浓度，以及每小时的温度(t)、气压(p)、湿度(h)、风速(ws)和风向(wd)。
[0051]
之后，采用最小冗余最大相关性方法对空气质量序列中11维输入特征(pm2.5、pm10、no2、co、o3、so2、t、p、h、ws和wd)进行特征选择，得出排序结果，进而选择出最优的输入特征子集。具体如下：
[0052]
首先将pm2.5时间序列作为目标特征，其他10维特征作为比较变量，基于最小冗余最大相关性方法对空气质量序列进行定量分析，并绘制不同因素对pm2.5时间序列的相关性，如图2(a)所示，得到排序结果为pm2.5、co、ws、pm10、h、wd、so2、no2、p、t、o3。然后，采用交叉验证得到的空气质量序列的最优输入特征子集，如图2(b)所示。根据图2(b)的预测结果，当t
′
＝5维时，预测误差最小，为9.199。最终，确定最优5维输入特征子集为pm2.5、co、ws、pm10和h。
[0053]
步骤2：对步骤1中获得的最优5维输入特征子集进行相空间重构，利用坐标延迟法计算得出5维特征的延迟时间τ和嵌入维数m，结果如表1所示。
[0054]
表1变量的相空间重构参数
[0055][0056]
由表1的结果可以推出，对于pm2.5、pm10、co、h和ws，最优5维输入特征子集的延迟时间为[8,8,6,4,4]，嵌入维数为[2,2,2,4,4]，重构后形成14维新的输入特征集，将其用于预测模型的输入。
[0057]
步骤3：将步骤2获得的14维新的输入特征集的前75％划为训练集，后25％划为测试集。利用公式(6)训练所提的补充泄露积分回声状态网络模型，确定得出空气质量预测最佳模型；之后基于最佳模型，利用公式(7)和公式(8)在测试集上对pm2.5时间序列进行预测输出。
[0058]
步骤4：本发明采用5个评价指标评估所提模型对pm2.5时间序列的预测结果，即均方根误差(rmse)、归一化均方根误差(nrmse)、平均绝对误差(mae)、对称平均绝对百分比误差(smape)和皮尔逊相关系数(r)，计算公式如下所示：
[0059][0060]
其中，n是样本数，是预测输出，y(n)是真实值，和分别为真实值和预测输出的平均值。
[0061]
在上述评价指标中，rmse、nrmse、mae、smape的估计指标值越小，模型的预测结果越好。r＝1表示和y(n)是线性相关，r＝0表示不相关。当r值在0到1之间时，值越大，线性相关性越强。
[0062]
本发明选择回声状态网络(esn)、泄漏积分回声状态网络(li-esn)、极限学习机(elm)、分层极限学习机(h-elm)和堆叠自编码器(sae)作为比较方法。具体pm2.5时间序列(浓度)一步(1小时)预测结果见表2。由表2可以看出，本发明在1小时的单步预测中取得了较好的预测结果。pm2.5浓度一步(1小时)预测结果如图3所示。图4绘制了预测值与实际数据的散点图。从图中可以看出，预测值与实际值有很好的线性关系。
[0063]
表2pm2.5浓度一步(1小时)预测结果对比
[0064][0065]
同时，图5给出了pm2.5浓度五步(5小时)预测的模拟结果。预测曲线能较好地跟踪原始输入，中期预测效果也较好。表3给出了pm2.5浓度五步(5小时)预测结果。本发明的pm2.5浓度十步(10小时)预测结果如图6所示。如图6所示，在部分峰值处，预测曲线仍能大致拟合原始数据的波动趋势，这正是由于本发明充分利用了历史信息的作用，与其他算法的比较见表4。本发明除smape外，其余4个参数指标均达到最优值，充分证明了本发明在长期预测中的有效性。
[0066]
表3pm2.5浓度五步(5小时)预测结果对比
[0067][0068]
表4十步(10小时)预测结果对比
[0069][0070]
为了进一步证明本发明的性能，各比较方法的运行时间见表5。结果表明，基于回声状态网络和elm模型的训练时间比基于深度学习模型的训练时间要短得多。这主要是因为深度学习的训练过程耗费了大量的时间，而这种耗时的过程是其他神经网络不需要的。此外，本发明能够在可接受的时间范围内完成训练和测试。结果表明，无论是在预测精度上还是在耗时上，本发明都能获得较好的预测结果。
[0071]
表5pm2.5浓度十步(10小时)预测结果对比
[0072][0073]
结合表1～表4可以看出，本发明在各个评价指标上获得了优秀的预测性能，不仅增强了储备池对空气质量序列的特征提取和记忆能力，而且改善了历史状态对当前状态的影响。
[0074]
最后应说明的是：以上实例仅用以表达本发明的实施方式，应理解实例仅用于说明本发明而不用于限制本发明的范围，在阅读了本发明之后，对于本领域的技术人员来说，在不脱离本发明构思的前提下，还可以做出若干变形和改进，对本发明的各种等价形式的修改均落于本技术所附权利要求所限定的范围，均属于本发明的保护范围。

技术特征：
1.一种基于补充泄漏积分回声状态网络的空气质量预测方法，其特征在于，首先，收集空气质量序列，采用最小冗余最大相关性方法对空气质量序列进行特征选择，确定输入特征的排序结果，进而选择出最优的输入特征子集；其次，对空气质量序列中最优的输入特征子集进行相空间重构，旨在充分挖掘输入特征子集中有用信息，形成新的特征输入集；再次，对新的特征输入集按照比例分成训练集和测试集，在训练集上使用岭回归算法训练补充漏积分器回声状态网络模型，训练完毕得到最佳预测模型，在测试集上进行pm2.5时间序列进行预测输出。2.根据权利要求1所述的一种基于补充泄漏积分回声状态网络的空气质量预测方法，其特征在于，具体步骤如下：步骤1：收集多维空气质量序列t，t＝[x1,x2,
…
,x
t
]，并采用最小冗余最大相关性方法对空气质量序列中的输入特征进行特征选择，得出输入特征的排序结果，进而选择出最优的输入特征子集；具体如下：首先，计算最大相关度d和最小冗余度r：首先，计算最大相关度d和最小冗余度r：其中，x
i
和x
k
分别表示空气质量序列t中的第i个和第k个特征；c表示目标特征，即pm2.5；i(x
i
；c)表示第i个特征和目标特征c的相关函数；i(x
i
；x
k
)表示第i个特征和第k个特征的相关函数；相关函数i(
·
；
·
)通过以下公式计算：其中，u和u
′
表示两个变量；p(
·
,
·
)表示概率；联立公式(1)-(3)，确定最小冗余最大相关性目标函数，定义如下：之后，采用前向选择算法求解上述目标函数，即可得到输入特征的排序结果，以上过程称为最小冗余最大相关性方法；最后，通过交叉验证得到空气质量序列的最优t
′
维输入特征子集，其中t
′
＜t；步骤2：对步骤1中获得的最优t
′
维输入特征子集进行相空间重构，形成新的输入特征集作为空气质量预测模型的输入；具体如下：以空气质量序列中最优t
′
维输入特征子集的第j个特征x
j
进行说明，采用坐标延迟法对其进行相空间重构，计算公式如下：x
j
(n)＝[x
j
(n),x
j
(n+τ
j
),
…
,x
j
(n+(m
j-1)τ
j
)]
ꢀꢀꢀꢀ
(5)其中，x
j
(n)为第n时刻第j个特征的输入特征集；m
j
和τ
j
分别为第j个特征的嵌入维数和延迟时间；同理，最优t
′
维输入特征子集的其他特征均采用如上方法；最终，形成新的高维输入特征集u，即u＝[u(1),u(2),
…
,u(n)]，其中n为总样本数或总时刻；在第n时刻，u(n)＝
[x1(n),x2(n),
…
,x
t
′
(n)]；步骤3：将步骤2形成的高维输入特征集u按照比例划分成训练集和测试集；在训练集上训练所提补充漏积分器回声状态网络模型，该模型考虑前一状态对当前状态的影响，采用改进的储备池状态更新形式进行训练，从而得到空气质量预测最佳模型；之后，在测试集上应用最佳模型，对pm2.5时间序列进行预测输出；具体如下：在训练集上训练所提出的补充泄漏积分回声状态网络模型，确定空气质量预测最佳模型；所提模型的状态更新公式如下所示：s(n+1)＝a
×
tanh(w
in
×
[1；u(n+1)]+w
×
s(n))+b
×
s(n)+(1-a-b)
×
s(n-1)
ꢀꢀꢀꢀ
(6)其中，s(n+1)和s(n-1)分别表示回声状态网络中储备池第n+1和n-1时刻的状态，a和b是衰减参数，且a和b均大于0，w
in
表示输入层与储备池之间的权重，w表示储备池的权重；基于空气质量预测最佳模型，在测试集上计算pm2.5时间序列的预测输出如下所示：其中，w
out
为储备池与输出层的权重，采用岭回归对其求解；具体言之，岭回归是一种用于共线性数据分析的有偏估计回归方法，能够克服最小二乘解中的不适定问题，避免过拟合问题，其表达式如下:w
out
＝(s
t
s+ki)-1
s
t
y
ꢀꢀꢀꢀꢀꢀ
(8)其中，k是正则化参数，i是单位矩阵，s和y分别是由s(n)和y(n)组成的矩阵。3.根据权利要求2所述的一种基于补充泄漏积分回声状态网络的空气质量预测方法，其特征在于，所述的步骤3中，将高维输入特征集u的前75％作为训练集，后25％作为测试集。

技术总结
一种基于补充泄漏积分回声状态网络的空气质量预测方法，属于时间序列预测领域。首先，收集空气质量序列，采用最小冗余最大相关性方法对空气质量序列进行特征选择，确定输入特征的排序结果，选择出最优的输入特征子集；其次，对最优的输入特征子集进行相空间重构，形成新的特征输入集；再次，将特征输入集分成训练集和测试集，在训练集上使用岭回归算法训练补充漏积分器回声状态网络模型，训练完毕得到最佳预测模型，在测试集上进行PM2.5时间序列进行预测输出；最后，利用评价指标度量模型的预测精度。本发明不仅增强了储备池对多元时间序列的特征提取和记忆能力，而且改善了历史演化状态对当前状态的影响。在实际应用中，准确的序列预测有助于监测空气质量，对大气污染防治做出合理、科学的决策。科学的决策。科学的决策。


技术研发人员：许星晗 夏慧娟 胡磊 任伟杰 张成坤
受保护的技术使用者：大连理工大学
技术研发日：2023.05.10
技术公布日：2023/8/9