基于外界影響及時(shí)序因素的PM2.5預(yù)測(cè)研究

楊艷梅,程宗毛

(杭州電子科技大學(xué) 理學(xué)院,浙江 杭州310037)

霧霾問(wèn)題不僅對(duì)人體健康造成極大危害，也破壞了植被的正常生理機(jī)能，對(duì)天氣及氣候的影響十分顯著。PM2.5作為霧霾的主要成分之一，是指大氣中直徑小于或等于2.5微米的顆粒物，其已成為影響城市空氣質(zhì)量的主要污染物之一[1-2]。

國(guó)內(nèi)外已開(kāi)展了關(guān)于PM2.5預(yù)測(cè)的研究，提出了一系列數(shù)據(jù)分析模型。PM2.5濃度預(yù)測(cè)方式大致分為兩類：(1)時(shí)間序列分析預(yù)測(cè)法。將歷史數(shù)據(jù)作為輸入變量，利用傳統(tǒng)統(tǒng)計(jì)模型或機(jī)器學(xué)習(xí)進(jìn)行時(shí)間序列預(yù)測(cè)。例如，文獻(xiàn)[3]利用傳統(tǒng)的自回歸滑動(dòng)模型(Autoregressive Moving Average Model，ARMA)預(yù)測(cè)長(zhǎng)沙市PM2.5濃度；文獻(xiàn)[4]則利用北京市PM2.5濃度的歷史數(shù)據(jù)訓(xùn)練機(jī)器學(xué)習(xí)中的長(zhǎng)短期記憶模型(Long Short Term Memory，LSTM)來(lái)達(dá)到預(yù)測(cè)目的；(2)外界因素分析預(yù)測(cè)法。將PM2.5濃度成因作為輸入變量，從影響因素著手提高預(yù)測(cè)模型的解釋能力。例如，文獻(xiàn)[5]利用平均溫度、地表溫度和平均風(fēng)速等氣象因素來(lái)進(jìn)行基于BP(Back Propagation，BP)神經(jīng)網(wǎng)絡(luò)的PM2.5預(yù)測(cè)；文獻(xiàn)[6]利用20種氣象因素和6種環(huán)境因素建立了基于量子粒子群算法和徑向基神經(jīng)網(wǎng)絡(luò)算法相結(jié)合的預(yù)測(cè)模型。上述兩種方式各有優(yōu)勢(shì)，也都有一定的局限性。時(shí)間序列分析預(yù)測(cè)法突出了時(shí)間因素在預(yù)測(cè)中的作用，建立了時(shí)間序列內(nèi)部關(guān)系，但其暫不考慮外界具體因素的影響，精確度有待提高。外界因素分析預(yù)測(cè)法忽略了PM2.5濃度自身的內(nèi)部聯(lián)系，且PM2.5的影響因素具有多樣性，使得該方法無(wú)法保證對(duì)外界因素?cái)?shù)據(jù)的收集及準(zhǔn)確性。

本文兼納了兩種思路，采用BP-EEMD-LSTM組合模型進(jìn)行預(yù)測(cè)。本文先依據(jù)主要外界影響因素建立基于BP神經(jīng)網(wǎng)絡(luò)的PM2.5濃度初步預(yù)測(cè)模型，此模型預(yù)測(cè)產(chǎn)生的殘差序列中包含了時(shí)序因素。針對(duì)具有典型的非線性、非平穩(wěn)性的復(fù)雜性較高的時(shí)間序列，如果直接建立預(yù)測(cè)模型，則很難達(dá)到預(yù)測(cè)精度要求[4]。為了降低時(shí)間序列復(fù)雜程度，本文先將殘差序列進(jìn)行集合經(jīng)驗(yàn)?zāi)B(tài)分解(Ensemble Empirical Mode Decomposition，EEMD)；然后完成復(fù)雜程度較低的分解分量的 LSTM建模，將各個(gè)分量的預(yù)測(cè)值累加得到最終的殘差修正值；最后，將BP模型得到的初步預(yù)測(cè)結(jié)果與EEMD-LSTM得到的殘差修正結(jié)果線性相加為最終的預(yù)測(cè)結(jié)果。本文通過(guò)實(shí)例證實(shí)了該組合模型在PM2.5日濃度預(yù)測(cè)方面有一定的優(yōu)越性。

基于BP-EEMD-LSTM組合模型的PM2.5濃度預(yù)測(cè)流程如圖1所示。

圖1 預(yù)測(cè)流程圖Figure 1. Flow chart of prediction

1 理論模型

1.1 BP神經(jīng)網(wǎng)絡(luò)

(1)

為了使得誤差函數(shù)達(dá)到最小，利用梯度下降法反復(fù)修正權(quán)值ω和閾值b。假設(shè)η、α分別為權(quán)值和閾值的學(xué)習(xí)速率，則更新計(jì)算式為式(2)。

(2)

(3)

預(yù)測(cè)時(shí)，將規(guī)范形式的數(shù)據(jù)輸入模型中，即可利用重復(fù)更新完成后的權(quán)值和閾值獲得最終的預(yù)測(cè)值。

1.2 集合經(jīng)驗(yàn)?zāi)B(tài)分解

研究人員于2009年提出了適用于分析非平穩(wěn)、非線性復(fù)雜信號(hào)的EEMD(Ensemble Empirical Mode Decomposition)[8]，其可將復(fù)雜信號(hào)分解成若干個(gè)按頻率高低排列的本征模函數(shù)(Intrinsic Mode Function，IMF)和一個(gè)趨勢(shì)項(xiàng)。該分解算法在原信號(hào)中多次添加白噪聲，利用白噪聲序列均勻分布的頻率改變信號(hào)極值點(diǎn)的特性，改善了經(jīng)驗(yàn)?zāi)B(tài)分解(Empirical Mode Decom-position，EMD)的模態(tài)混疊現(xiàn)象[9]。EEMD的算法步驟如下[10-11]：

步驟1在原始信號(hào)x(t)中添加平均值為零，方差恒定的等長(zhǎng)隨機(jī)高斯白噪聲ni(t)，得到一個(gè)新的信號(hào)xi(t)

xi(t)=x(t)+ni(t)

(4)

式中，ni(t)為第i次添加的白噪聲；xi(t)為第i的新信號(hào)；

步驟2將噪聲疊加信號(hào)xi(t)進(jìn)行EMD分解，得到n個(gè)IMF分量ci,j(t)和一個(gè)趨勢(shì)項(xiàng)ri(t)。

(5)

重復(fù)執(zhí)行步驟1和步驟2共N次，計(jì)算每個(gè)分量以及余量經(jīng)過(guò)N次實(shí)驗(yàn)后的平均值，可消除多次加入白噪聲對(duì)真實(shí)IMF分量的影響，最終得到n個(gè)IMF分量c1(t)～cn(t)和趨勢(shì)項(xiàng)r(t)。分解結(jié)果為

(6)

(7)

1.3 長(zhǎng)短期記憶神經(jīng)網(wǎng)絡(luò)

LSTM作為一種時(shí)間遞歸神經(jīng)網(wǎng)絡(luò)，適用于處理和預(yù)測(cè)時(shí)間尺度較長(zhǎng)的事件，并可避免梯度消失問(wèn)題[12]。LSTM單元包括內(nèi)部單元、輸入門、輸出門和遺忘門。長(zhǎng)短時(shí)記憶網(wǎng)絡(luò)主要通過(guò)“門”來(lái)丟棄或者增加信息[13]。LSTM單元結(jié)構(gòu)圖如圖2所示。

圖2 LSTM單元結(jié)構(gòu)圖Figure 2. Architecture of LSTM cell

LSTM網(wǎng)絡(luò)更新[14-15]為：

(1)遺忘門過(guò)程函數(shù)

ft=σ(Wf·[ht-1,xt]+bf)

(8)

式中，ft為輸出遺忘率；σ表示sigmoid函數(shù)；Wf為遺忘門權(quán)重矩陣；ht-1表示上一時(shí)刻的輸出結(jié)果；xt為該時(shí)刻的輸入信息；bf為遺忘門的偏置項(xiàng)；

(2)輸入門第1部分過(guò)程函數(shù)為

it=σ(Wi·[ht-1,xt]+bi)

(9)

式中，it為輸出結(jié)果；Wi為輸入門權(quán)重矩陣；bi為輸入門的偏置項(xiàng)。輸入門的第2部分由tanh函數(shù)產(chǎn)生1個(gè)候選向量為

(10)

(3)每個(gè)單元更新后的狀態(tài)表示為

(11)

式中，Cd為更新后的當(dāng)前時(shí)刻的單元狀態(tài)；Ct-1為上一次的單元狀態(tài)；

(4)輸出門的計(jì)算輸出門結(jié)果Ot為

Ot=σ(Wo·[ht-1,xt]+b0)

(12)

式中，Ot為輸出門；Wo為輸出門權(quán)重矩陣；bo為輸出門偏置項(xiàng)。利用Ot得出當(dāng)時(shí)刻的記憶輸出ht，如式(13)所示。

ht=Ot·tanh(Ct)

(13)

2 模型建立

本文基于BP-EEMD-LSTM組合模型進(jìn)行PM2.5日濃度預(yù)測(cè)，組合模型分為以下步驟：

步驟1BP初步預(yù)測(cè)。將訓(xùn)練集主因素作為BP的輸入變量，并將訓(xùn)練集PM2.5日濃度作為輸出變量來(lái)訓(xùn)練此模型。將主因素輸入訓(xùn)練后的BP模型，得到預(yù)測(cè)集的PM2.5日濃度BP初步預(yù)測(cè)值Yb和訓(xùn)練集的PM2.5日濃度BP預(yù)測(cè)值Ybt，并計(jì)算訓(xùn)練集的殘差時(shí)間序列Ye=Yb-Ybt；

步驟2采用EEMD分解殘差時(shí)間序列。步驟1產(chǎn)生的殘差時(shí)間序列Ye通過(guò)EEMD分解，得到n個(gè)不同頻率的IMF分量和一個(gè)趨勢(shì)項(xiàng)；

3 數(shù)據(jù)選取與預(yù)處理

本文相關(guān)數(shù)據(jù)來(lái)源于中國(guó)空氣質(zhì)量在線監(jiān)測(cè)分析平臺(tái)。選擇來(lái)源于該平臺(tái)的杭州市PM2.5、PM10、SO2、CO、NO2、O3的日濃度以及平均溫度、最高溫度、最低溫度、平均相對(duì)濕度、總降雨量/總?cè)谘┝?、平均風(fēng)速和最大持續(xù)風(fēng)速[16]共13個(gè)指標(biāo)數(shù)據(jù)作為研究對(duì)象。將2014年1月1日～2019年5月31日的PM2.5日濃度和以上12個(gè)相關(guān)數(shù)據(jù)作為實(shí)驗(yàn)數(shù)據(jù)，設(shè)置2019年5月數(shù)據(jù)作為預(yù)測(cè)集，其它月份數(shù)據(jù)作為訓(xùn)練集，各因素與PM2.5日濃度的相關(guān)性如表1所示，相關(guān)系數(shù)越大表示關(guān)系越密切。

表1 相關(guān)性分析

由表1可以看出，在參與分析的因子中，PM10、SO2、CO、NO2是影響 PM2.5變化最主要的因素[16]，相關(guān)系數(shù)均達(dá)到了0.5以上。因此，本文選取PM10、SO2、CO、NO2日濃度作為影響PM2.5日濃度的主因素，PM10、SO2、CO、NO2日濃度也將作為BP預(yù)測(cè)模型的輸入量。

4 實(shí)證分析

4.1 BP初步預(yù)測(cè)

由相關(guān)性分析可知，BP神經(jīng)網(wǎng)絡(luò)的輸入數(shù)據(jù)包含4個(gè)特征，即PM10、SO2、CO、NO2日濃度，因此設(shè)立四輸入單輸出的BP模型，具體操作為：(1)模型隱藏層設(shè)置。構(gòu)建3層隱藏層，將第1層神經(jīng)元個(gè)數(shù)設(shè)為12，將第2層神經(jīng)元個(gè)數(shù)設(shè)置為2，第3層神經(jīng)元個(gè)數(shù)設(shè)為1，輸入數(shù)據(jù)維度為4，采用tanh作為其激活函數(shù)；(2)模型參數(shù)設(shè)置。本文設(shè)置損失函數(shù)為絕對(duì)值均差(Mean Absolute Error，MAE)，選擇隨機(jī)梯度下降(Stochastic Gradient Descent，SGD)為優(yōu)化器，為防止模型訓(xùn)練出現(xiàn)過(guò)擬合現(xiàn)象，設(shè)置Dropout(0.4)；(3)模型訓(xùn)練。設(shè)置學(xué)習(xí)步長(zhǎng)為0.01；(4)模型預(yù)測(cè)。將外界成因數(shù)據(jù)輸入訓(xùn)練后的模型得到預(yù)測(cè)集的PM2.5日濃度預(yù)測(cè)值Yb以及訓(xùn)練集的PM2.5日濃度殘差值Ye。圖3為此次BP模型初步預(yù)測(cè)值的結(jié)果，圖4表述產(chǎn)生的訓(xùn)練集殘差時(shí)間序列圖。

圖3 BP預(yù)測(cè)值與真實(shí)值對(duì)比Figure 3. Comparison between BP predicted value and real value

圖4 訓(xùn)練集殘差序列Figure 4. The residual series of training set

4.2 EEMD分解殘差時(shí)間序列

將得到的訓(xùn)練集的殘差時(shí)間序列Ye進(jìn)行EEMD分解，得到10個(gè)不同頻率且相對(duì)穩(wěn)定的IMF分量c1～c10和趨勢(shì)項(xiàng)r，如圖5所示。原始?xì)埐钚蛄型ㄟ^(guò)EEMD分解后得到了變化規(guī)律更加明顯的數(shù)據(jù)，進(jìn)而降低了時(shí)間序列的復(fù)雜程度，減少了預(yù)測(cè)殘差的難度。

圖5 殘差序列的EEMD分解結(jié)果Figure 5. The decomposition results of residual series by EEMD

4.3 LSTM模型預(yù)測(cè)殘差修正值

為了適用于訓(xùn)練LSTM，需要使用滑動(dòng)窗口將各分量時(shí)間序列數(shù)據(jù)轉(zhuǎn)換為監(jiān)督學(xué)習(xí)型數(shù)據(jù)。假設(shè)LSTM模型每組輸入數(shù)據(jù)長(zhǎng)度為n，即為預(yù)測(cè)目標(biāo)之前的n天分量數(shù)據(jù)，用“過(guò)去”預(yù)測(cè)“未來(lái)”。設(shè)置預(yù)測(cè)步長(zhǎng)為1步，每次只能預(yù)測(cè)1個(gè)數(shù)據(jù)，之后將預(yù)測(cè)后的數(shù)據(jù)作為新的輸入數(shù)據(jù)，從而實(shí)現(xiàn)模型的動(dòng)態(tài)更新。輸入數(shù)據(jù)及輸出數(shù)據(jù)如表2所示。

表2 LSTM模型的輸入輸出數(shù)據(jù)

圖6 基于EEMD-LSTM模型的殘差修正Figure 6. The residual correction based on EEMD-LSTM model

4.4 結(jié)果分析

圖7 各模型的預(yù)測(cè)結(jié)果對(duì)比Figure 7. Comparison of prediction results of each model

由圖7可知，BP-EEMD-LSTM 組合模型具有良好的預(yù)測(cè)效果，能準(zhǔn)確地預(yù)測(cè)出PM2.5日濃度整體變化趨勢(shì)，并響應(yīng)PM2.5日濃度的波動(dòng)變化。BP 模型、BP-LSTM模型和BP-EMD-LSTM模型都能較好地預(yù)測(cè)PM2.5日濃度趨勢(shì)，但BP 模型整體的預(yù)測(cè)準(zhǔn)確度較差；BP-LSTM模型在5月1日～5月8日、5月19日～5月22日以及5月26日～5月30日預(yù)測(cè)存在較大誤差；BP-EMD-LSTM模型則無(wú)法準(zhǔn)確預(yù)測(cè)了5月1日～5月5日以及5月16日～5月19日的PM2.5日濃度。同時(shí)，本文采用均方根誤差(Root Mean Squared Error，RMSE)、相關(guān)系數(shù)R和平均百分比誤差(Mean Absolute Percentage Error，MAPE)來(lái)評(píng)測(cè)4種預(yù)測(cè)模型的性能，以便更直觀地定量分析預(yù)測(cè)結(jié)果。比較分析結(jié)果如表3所示。

表3 模型預(yù)測(cè)的誤差評(píng)估指標(biāo)

由表3可知，本文提出的BP-EEMD-LSTM組合模型的相關(guān)系數(shù)更接近1，說(shuō)明該組合模型能夠較好地捕捉PM2.5日濃度變化趨勢(shì)。該組合模型的均方根誤差為2.74，較BP模型提高了73.4%，較BP-LSTM組合模型提高了54.4%，較BP-EMD-LSTM組合模型提高了26.5%。此外，BP-EEMD-LSTM組合模型的平均百分比誤差為8.01，較BP模型提高了73.5%，較BP-LSTM組合模型提高了54.0%，較BP-EMD-LSTM組合模型提高了24.2%。由此可知，本文采取的預(yù)測(cè)模型比其他4種模型預(yù)測(cè)結(jié)果更加精確、有效，整體上BP-EEMD-LSTM組合預(yù)測(cè)模型的結(jié)果更接近真實(shí)值。

5 結(jié)束語(yǔ)

霧霾問(wèn)題現(xiàn)在已成為中國(guó)大部分城市的主要空氣污染問(wèn)題之一，針對(duì)霧霾主要標(biāo)志物PM2.5的濃度預(yù)測(cè)尤為重要[18-19]。影響PM2.5預(yù)測(cè)的因素較多，本文先利用相關(guān)性分析方法選擇了相關(guān)性最大的4個(gè)指標(biāo)作為預(yù)測(cè)PM2.5的主要外界因素；然后，將主要外界因素輸入BP模型得到初步預(yù)測(cè)值與歷史殘差序列，再由EEMD-LSTM模型基于歷史殘差序列做時(shí)間序列預(yù)測(cè)，預(yù)測(cè)結(jié)果為殘差修正值；最后，將BP初步預(yù)測(cè)值與殘差修正值結(jié)合從而完成對(duì)PM2.5的預(yù)測(cè)。實(shí)驗(yàn)結(jié)果證明本文提出的基于外界影響及時(shí)序因素的PM2.5預(yù)測(cè)研究即BP-EEMD-LSTM組合模型具有更好的預(yù)測(cè)效果與精度，為PM2.5的預(yù)測(cè)研究提供了新的思路。

本文主要探討了氣象因素和大氣污染物因素對(duì)PM2.5濃度預(yù)測(cè)的影響，其他因素(例如地理位置、環(huán)保政策等因素)仍有待進(jìn)一步的探索研究。