根据2018年哈尔滨市生态环境状况公报^[1]显示, 2018年空气质量共超标52天, 超标天数中首要污染物为PM_2.5的天数共有39天. PM_2.5又称细颗粒物, 由于其粒径小, 活性强且附带有害物质的特性, 会对呼吸道系统和心血管系统造成伤害^[2]. 空气污染问题不仅与人们的日常生活、身体健康密切相关, 还会对城市的招商引资、经济发展造成较大的影响. 因此对大气中PM_2.5浓度的预测有着重要的意义.

空气质量预测模型主要分为基于大气运动学的数值预报模型(NWP)与基于机器学习算法的统计模型^[3]. 数值预报模型是从大气内部物理规律如大气动力学、热力学等出发, 建立对应的数学物理模型, 用数理方法, 借助大型计算机的计算能力, 建立数值预测空气污染物浓度的动态分布运输和扩散模型. 例如: JM Xu等^[4]基于尺度空气质量模式系统(CMAQ)提出了CMAQ-MOS模型, 有效纠正了CMAQ中平均污染物排放清单所导致的系统性预测误差; Ming-Tung Chuang等^[5]应用结合化学的天气研究与预测模型(WRF-Chem-MADRID)对美国东南部地区进行预测, 结果显示在O₃和PM_2.5预测中表现良好. 机器学习模型则是利用统计学、概率论以及复杂的算法搭建模型, 从已知数据中挖掘数据关系, 实现精细化预测. 例如王敏、孙宝磊等^{[6, 7]}利用BP神经网络模型对PM_2.5浓度进行预测, 实验效果较好; 彭斯俊等^[8]利用ARIMA模型对于时间序列短期预测有着较高的精度的特点, 挖掘序列内部规律, 得到了不错的结果; 喻其炳等^[9]通过K-means方法筛选相似性高的数据来训练基于支持向量机的PM_2.5浓度预测模型, 提升效果明显.

虽然以上所提PM_2.5浓度预测模型各有优点, 但是其往往只根据目标站点的污染因子和气象数据, 而忽略了PM_2.5的分布也具有区域性^[10]. 针对这种问题, 本文提出了一种基于时空特征的KNN-LSTM网络模型, 通过KNN算法选择目标站点的空间相关信息, 利用LSTM具有时间记忆的特点, 从时间和空间两个维度考虑, 实现对PM_2.5浓度的预测. 同时还利用哈尔滨市多个空气质量监测站的污染物数据, 与BP神经网络模型和普通的LSTM模型进行实验对比, 结果表明本文所提KNN-LSTM模型能较好预测未来时刻的PM_2.5浓度值.

1 理论与模型 1.1 长短时记忆网络

长短时记忆网络(LSTM)主要是为了解决循环神经网络(RNN)在长序列训练过程中梯度消失的问题^{[11, 12]}. 与RNN相比, LSTM增加了3个门控: 输入门、遗忘门、输出门. 门的作用是为了控制之前的隐藏状态、当前的输入等信息, 确定哪些信息该丢弃、哪些信息该保留. LSTM的模型结构如图1所示.

图1中 $\varTheta $ 代表操作矩阵中对应的元素相乘, $ + $ 代表进行矩阵的加法.

对于LSTM当前的输入x^t和上一状态传递的h^t–1, 通过下面公式得到4个状态:

${{ }}{{{{\textit{z}}}}^i}{{ = }}\sigma ({W_{xi}}{x^t} + {W_{hi}}{h^{{{t - 1}}}} + {b_i}) $

(1)

${{{{\textit{z}}}}^f}{{ = }}\sigma ({W_{xf}}{x^t} + {W_{hf}}{h^{{{t - 1}}}} + {b_f})$

(2)

${{{{\textit{z}}}}^o}{{ = }}\sigma ({W_{xo}}{x^t} + {W_{ho}}{h^{{{t - 1}}}} + {b_o}) $

(3)

${{{\textit{z}} = \tanh}}({W_{x{\textit{z}}}}{x^t} + {W_{h{\textit{z}}}}{h^{{{t - 1}}}} + {b_{\textit{z}}})$

(4)

其中, z^f、zⁱ、z^o、z分别表示遗忘门控、输入门控、输出门控和当前输入内容, W_xi、W_xf、W_xo、W_xz分别代表输入层到输入门、遗忘门、输出门和细胞状态的权重矩阵; b_i、b_f、b_o、b_z分别为输入门、遗忘门、输出门和细胞状态的偏移量. σ为Sigmoid函数, tanh为双曲正切函数. LSTM存在两种传输状态: c^t (cell state)和h^t (hidden state), 其中LSTM中的c^t对应RNN中的h^t. t时刻在LSTM内部主要进行以下3步, 如式(5)、式(6)、式(7)所示.

(1) 经过z^f控制上一个状态c^t–1哪些需要保留或遗忘.

(2) 经过zⁱ对输入x^t选择记忆.

(3) 经过z^o控制哪些会作为当前状态的输出.

${c^t} = {{\textit{z}}^f}\varTheta {c^{t - 1}} + {{\textit{z}}^i}\varTheta {\textit{z}} $

(5)

${h^t} = {{\textit{z}}^o}\varTheta \tanh ({c^t}) $

(6)

${y^t} = \sigma (W'{h^t}) $

(7)

其中, 隐藏层到输出层的权重矩阵记为W'.

1.2 KNN算法

KNN算法^[13]是一种有监督学习的分类算法. 其实现较为简单, 训练速度较快. KNN用空间中两个点的距离来度量其相似度, 距离越小, 相似度越高. 通过最邻近的k个点归属的主要类别, 来对测试点进行分类. 常见的距离度量方式有欧氏距离、马氏距离、曼哈顿距离等. 特征向量X_i、X_j之间的欧氏距离计算公式如下:

${d_{ij}}{ = ^{}}||{X_i} - {X_j}|{|_2}^{} = \sqrt {\sum\limits_{{{m = }}1}^{{k}} {{{({x_i}({{m}}) - {x_j}({{m}}))}^2}} } $

(8)

其中, k表示特征向量的维度, x_i(m)、x_j(m)分别为X_i、X_j第m维的值, m的取值范围为1,2,3,···,k.

1.3 模型介绍 1.3.1 基于单站点的LSTM预测模型

本文首先实现了基于LSTM的单站点PM_2.5浓度预测模型, 该模型是通过目标站点空气质量六因子(NO₂、PM_2.5、PM₁₀、SO₂、CO、O₃)的历史数据, 来预测当前时刻的PM_2.5浓度. 若当前时刻为t, 其输入为t–s_t, t–s_t+1,···,t–1时刻六因子的小时数据. 输入的数据先后经由LSTM层和两层全连接层(Dense层)后得到PM_2.5的浓度预测结果. 其中损失函数采用均方误差(MSE), 计算公式如式(9)所示. 模型的优化器选用Adam优化器^[14].

${{MSE = }}\dfrac{1}{n}\sum\limits_{i = 1}^n {{{({{{X}}_{{i}}} - {X_r})}^2}} $

(9)

其中, X_i为PM_2.5浓度预测值, X_r为PM_2.5浓度真实值, n为训练集样本数.

实验选取均方根误差(RMSE)和平均绝对误差(MAE)作为评价指标. 计算公式如下:

${{RMSE = }}\sqrt {\dfrac{1}{n}\sum\limits_{i = 1}^n {{{({{{X}}_{{i}}} - {X_r})}^2}} } $

(10)

${{MAE = }}\dfrac{1}{n}\sum\limits_{i = 1}^n {|{{{X}}_{{i}}} - {X_r}|} $

(11)

1.3.2 KNN-LSTM预测模型

对于单站点LSTM模型只考虑时间特征, 而没有充分利用空间相关特征的问题, 本文利用KNN算法对邻近的空间因素进行筛选, 作为额外输入, 对其进行修正, 构建KNN-LSTM模型以实现对PM_2.5浓度更精准的预测.

KNN-LSTM模型主要可分为以下7步, 预测模型流程图如图2所示.

(1)数据预处理. 在原始数据使用前, 需要对缺失样本和异常值进行处理. 空气监测站对数据缺失的部分用–1表示, 统计后可知, 各因子的缺失部分占总数据量均不超过5%. 故采用简单方便的线性插值法进行填充. 由于其数据的收集都要经过一系列严格的审核流程, 认为其数据真实且有效, 不对数据进行去噪处理.

鉴于检测站点各个污染因子有着不同的量纲和量级, 为了取消由于量纲不同引起的误差和保证模型的高效性, 对数据进行max-min归一化处理^[15], 将数据转化到[0,1]的范围内, 并在预测结束后进行数据恢复. max-min归一化公式如下:

$x'_t = ({x_t} - {x_{\min }})/({x_{\max }} - {x_{\min }}) $

(12)

其中, ${x_t}$ 为原数据, $x'_t$ 为归一化后的数据, ${x_{\max }}$ 、 ${x_{\min }}$ 分别为数据中的最大值与最小值.

(2)采用KNN算法提取目标站点的空间相关特征,令K=1. 本文采用欧氏距离来衡量目标监测站和附近站点之间的相关程度, 距离越小, 空间相关性越强. KNN算法的选择过程如下所示.

1)根据单站点LSTM模型中的时间步长s_t构建PM_2.5状态矩阵. 假设其时间步长为s_t, 区域中共有m个站点, 为方便起见, 我们将t时刻站点i的时间步状态量 $X_t^i $ 和t时刻PM_2.5状态矩阵S_t定义如(13)、式(14)所示:

$X_t^i {{ = }} \left( {\begin{array}{*{20}{c}} {x_{t - st}^{{i}}}&{x_{t - st + 1}^i}&{...}&{x_{t - 2}^i}&{x_{t - 1}^i} \end{array}} \right)_{}^{\rm T}$

(13)

$\begin{split} & S_t^{} = \left( {\begin{array}{*{20}{c}} {X_t^1}&{X_t^2}&{...}&{X_t^{m - 1}}&{X_t^m} \end{array}} \right) = \\ & \left( {\begin{array}{*{20}{c}} {x_{t - st}^1}&{x_{t - st}^2}&{x_{t - st}^3}&{...}&{x_{t - st}^{m - 1}}&{x_{t - st}^m} \\ {x_{t - st + 1}^1}&{x_{t - st + 1}^2}&{x_{t - st + 1}^3}&{...}&{x_{t - st + 1}^{m - 1}}&{x_{t - st + 1}^m} \\ {x_{t - st + 2}^1}&{x_{t - st + 2}^2}&{x_{t - st + 2}^3}&{...}&{x_{t - st + 2}^{m - 1}}&{x_{t - st + 2}^m} \\ {...}&{...}&{...}&{...}&{...}&{...} \\ {x_{t - 1}^1}&{x_{t - 1}^2}&{x_{t - 1}^3}&{...}&{x_{t - 1}^{m - 1}}&{x_{t - 1}^m} \\ {x_t^1}&{x_t^2}&{x_t^3}&{...}&{x_t^{m - 1}}&{x_t^m} \end{array}} \right) \\ \end{split} $

(14)

其中， $x_n^u$ 为编号为u的站点在n时刻的PM_2.5浓度值(u=1, 2, 3,…, m; n=1, 2, 3,…, t).

2)通过计算目标站点与其他邻近站点(共m–1个) t时刻时间步状态量之间的欧氏距离, 并进行从小到大排序, 选取前k个对应监测站点的PM_2.5数据作为t时刻目标站点的空间相关因素, 记为X_sp:

$X_{sp}^{} {{ = }} \left( {\begin{array}{*{20}{c}} {x_{sp}^1}&{x_{sp}^2}&{..}&{x_{sp}^{k - 1}}&{x_{sp}^k} \end{array}} \right)$

(15)

其中， $x_{sp}^i$ 为t时刻与目标站点第i相关的站点PM_2.5浓度值.

(3)根据目标站点的历史污染数据, 利用LSTM网络提取时间特征. 输入仍为目标站点s_t时间步的空气质量六因子数据, t时刻输入数据 $lstm\_inpu{t_t}$ 如式(16)所示.

$\begin{array}{l} lstm\_inpu{t_t} = \left( {I_t^{\rm co}{\kern 1pt} {\kern 1pt} I_t^{{\rm n}{{\rm o}_2}}{\kern 1pt} I_t^{{\rm s}{{\rm o}_2}}{\kern 1pt} {\kern 1pt} I_t^{{{\rm o}_3}}{\kern 1pt} {\kern 1pt} I_t^{{\rm p}{{\rm m}_{10}}}{\kern 1pt} {\kern 1pt} I_t^{{\rm p}{{\rm m}_{2.5}}}} \right) = \\ \left( {\begin{array}{*{20}{c}} {i_{t - st}^{\rm co}}&{i_{t - st}^{{\rm n}{{\rm o}_2}}}&{i_{t - st}^{{\rm s}{{\rm o}_2}}}&{i_{t - st}^{{{\rm o}_3}}}&{i_{t - st}^{{\rm p}{{\rm m}_{10}}}}&{i_{t - st}^{{\rm p}{{\rm m}_{2.5}}}}\\ {i_{t - st + 1}^{\rm co}}&{i_{t - st + 1}^{{\rm n}{{\rm o}_2}}}&{i_{t - st + 1}^{{\rm s}{{\rm o}_2}}}&{i_{t - st + 1}^{{{\rm o}_3}}}&{i_{t - st}^{{\rm p}{{\rm m}_{10}}}}&{i_{t - st + 1}^{{\rm p}{{\rm m}_{2.5}}}}\\ {i_{t - st + 2}^{\rm co}}&{i_{t - st + 2}^{{\rm n}{{\rm o}_2}}}&{i_{t - st + 2}^{{\rm s}{{\rm o}_2}}}&{i_{t - st + 2}^{{{\rm o}_3}}}&{i_{t - st}^{{\rm p}{{\rm m}_{10}}}}&{i_{t - st + 2}^{{\rm p}{{\rm m}_{2.5}}}}\\ \cdots & \cdots & \cdots & \cdots & \cdots & \cdots \\ {i_{t - 2}^{\rm co}}&{i_{t - 2}^{{\rm n}{{\rm o}_2}}}&{i_{t - 2}^{{\rm s}{{\rm o}_2}}}&{i_{t - 2}^{{{\rm o}_3}}}&{i_{t - st}^{{\rm p}{{\rm m}_{10}}}}&{i_{t - 2}^{{\rm p}{{\rm m}_{2.5}}}}\\ {i_{t - 1}^{\rm co}}&{i_{t - 1}^{{\rm n}{{\rm o}_2}}}&{i_{t - 1}^{{\rm s}{{\rm o}_2}}}&{i_{t - 1}^{{{\rm o}_3}}}&{i_{t - st}^{{\rm p}{{\rm m}_{10}}}}&{i_{t - 1}^{{\rm p}{{\rm m}_{2.5}}}} \end{array}} \right) \end{array}$

(16)

其中, $i_j^{\rm {co}}$ 、 $ i_j^{\rm {n{o_2}}}$ 、 $i_j^{\rm {s{o_2}}}$ 、 $i_j^{\rm{o_3}}$ 、 $i_j^{\rm {p{m_{10}}}}$ 、 $i_j^{\rm {p{m_{2.5}}}}$ 为j时刻CO、NO₂、SO₂、O₃、PM₁₀、PM_2.5的浓度值.

将 $lstm\_inpu{t_t}$ 经过LSTM层得到的输出向量记为 $lstm\_{{out}}pu{t_t}$ .

(4)将KNN提取的空间相关特征 $x_{sp}^i$ 作为第二输入, 与LSTM层提取的时间特征 $lstm\_{{out}}pu{t_t}$ 相拼接, 而后通过全连接网络获取输出. 网络结构具体如图3所示.

(5)将数据划分成训练集与测试集, 用训练集数据进行模型训练, 用测试集数据进行验证, 并记录结果.

(6)令K=K+1, 重复步骤(3)~(5), 直到K>M为止. 其中M为邻近站点的数量.

(7)根据所记录测试集上评价指标(这里选用RMSE), 确定最优的K.

2 实验结果及分析 2.1 监测站点及数据展示

为了验证所提方法的有效性, 本文采用哈尔滨市区11个空气质量监测站4月1日至5月31日两个月共计1464条小时监测数据作为实验数据. 站点名称及坐标位置见表1.

阿城会宁空气质量监测站远离城区且与其他监测站点距离较远, 根据各站点PM_2.5时间序列的相关性可知: 阿城会宁站点与其他站点间的相关性系数均不超过0.6, 而其他站点任意两站点间的相关性系数均超过0.7, 因此在KNN进行空间站点选取时, 忽略阿城会宁站, 同时这证明了PM_2.5的浓度分布具有区域性. 图4为部分站点PM_2.5数据展示, 可以明显看出站点的数据的分布具有相似性.

2.2 实验环境及模型参数设置

实验环境及计算机配置如下: 程序设计语言为Python3.6.5; 开发环境采用Anaconda所自带的Spyder编辑器; Keras版本为2.2.4; Tensorflow版本为1.13.1; Scikit-Learn版本为0.19.1; 计算机处理器为AMD Ryzen 5 2500U, 内存为8 GB; 操作系统为Window10 x64.

对数据进行划分, 其中前1300条数据用于训练, 80%为训练集、20%为验证集, 后152条数据作为测试集. LSTM层的时间步长timestep为12, K值选取为2 (具体原因见2.3节). 具体网络参数设置如表2所示.

2.3 KNN算法中K值的选择

本文通过KNN算法来获取目标站点的空间特征, 由于K值的选择会影响到选择相关站点的数目, 进而影响实验结果. 图5给出了K值选择与评价指标RMSE变化之间的关系. 根据图5所示, 当K=2时, 预测效果最好. 故K值选择为2.

2.4 实验对比

为了验证本文所提模型的有效性, 本文选择传统BP神经网络以及传统LSTM神经网络作为对比实验, 3种模型均在相同的实验平台下进行, 各个模型预测值与实际值的对比见图6.

本文所有实验均采取均方根误差(RMSE)和平均绝对误差(MAE)作为评价指标. 为消除一次实验的偶然性, 每种实验均进行30次, 迭代50次.表3中评价指标数据均为30次实验结果的平均值(结果保留3位小数). 由表可知本文所提KNN-LSTM模型相对于BP神经网络, 平均绝对误差(MAE)、均方根误差(RMSE)分别降低了19.25%、13.23%; 相较于LSTM模型MAE、RMSE分别降低了4.29%、6.99%. 可见本文所提KNN-LSTM模型要优于BP神经网络和LSTM模型.

3 结语

本文使用哈尔滨市国控空气质量监测站空气质量六因子的小时数据进行PM_2.5浓度预测. 首先对缺失值采用线性插值的方法进行填充, 为减少由于量纲造成的误差, 对数据进行max-min归一化处理. 再利用KNN算法为目标站点选取与其空间相关的邻近站点, 得到其空间特征. 然后对目标站点搭建LSTM模型, 接着将得到的空间特征输入到LSTM模型中. 通过训练集数据进行训练, 再将训练好的模型用测试集进行评估. 通过选取最优的K值, 确定最终模型. 实验结果表明: 本文所提模型的预测结果曲线更为平滑且与真实值更加接近, 均方根误差和绝对平均误差均为最小, 预测效果相较BP神经网络和传统LSTM神经网络模型更好, 可以为PM_2.5的预警预报提供一定的参考.

由于数据收集的不充分, 仅仅利用空气质量六因子作为输入数据是相对片面的, 后续研究考虑引进更多影响因素以提升模型效果.