基于組合模型對(duì)云南省昭通市空氣質(zhì)量指數(shù)的預(yù)測(cè)

趙玉鳳，何淑仙，王彭德

（大理大學(xué)數(shù)學(xué)與計(jì)算機(jī)學(xué)院，云南大理 671003）

隨著城市化建設(shè)的加速，大氣污染嚴(yán)重，對(duì)空氣質(zhì)量的預(yù)測(cè)和評(píng)價(jià)愈顯重要。目前在空氣質(zhì)量評(píng)價(jià)上主要采用空氣質(zhì)量指數(shù)（air quality index，AQI），AQI是以環(huán)境空氣質(zhì)量標(biāo)準(zhǔn)為基礎(chǔ)，根據(jù)各種污染物、生活環(huán)境和生態(tài)的影響，將對(duì)空氣污染物濃度進(jìn)行綜合監(jiān)測(cè)的一系列數(shù)據(jù)化為單一的概念性數(shù)值的形式〔1〕，該指數(shù)值越小意味著空氣污染較輕，如果指數(shù)值較大，則說明空氣污染嚴(yán)重。空氣質(zhì)量是根據(jù)空氣中污染物的濃度反映空氣污染程度的指標(biāo)，主要有細(xì)顆粒物（PM2.5）、可吸入顆粒物（PM10）、二氧化氮（NO2）、二氧化硫（SO2）、一氧化碳（CO）和臭氧（O3）。當(dāng)前有越來越多的研究者關(guān)注空氣質(zhì)量預(yù)測(cè)問題，提出了不同的預(yù)測(cè)方法，進(jìn)而為空氣污染防治提供了科學(xué)的依據(jù)。如：于萍〔2〕利用空氣質(zhì)量數(shù)據(jù)建立了ARIMA（1，1）模型，并對(duì)短期內(nèi)的空氣質(zhì)量指數(shù)進(jìn)行了預(yù)測(cè)和判斷。許允之等〔3〕針對(duì)徐州市霧霾形勢(shì)，采用多元線性回歸方法建立模型，在該模型的基礎(chǔ)上，剔除不重要的變量和樣本異常值。經(jīng)過兩次改進(jìn)，將九元線性模型簡(jiǎn)化為四元線性模型，通過檢驗(yàn)和診斷，驗(yàn)證了四元線性回歸模型的準(zhǔn)確性和實(shí)用性。王娜娜〔4〕利用特定時(shí)間范圍內(nèi)的空氣質(zhì)量監(jiān)測(cè)數(shù)據(jù)，建立隨機(jī)森林模型對(duì)北京市未來一段時(shí)間的塵埃污染濃度進(jìn)行了預(yù)測(cè)分析，進(jìn)而為北京市的大氣污染治理提供了參考建議。郭利進(jìn)等〔5〕利用ARIMA模型，用卡爾曼濾波建立數(shù)學(xué)模型，將徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)融合到卡爾曼濾波中，提出了一種新的融合算法，并對(duì)AQI進(jìn)行了預(yù)測(cè)，改進(jìn)了單方法的延遲估計(jì)，提高了預(yù)測(cè)精度。呂霽洲〔6〕在數(shù)據(jù)訓(xùn)練后，驗(yàn)證了神經(jīng)網(wǎng)絡(luò)算法對(duì)空氣質(zhì)量預(yù)測(cè)的可行性。由于采用單一模型對(duì)空氣質(zhì)量進(jìn)行預(yù)測(cè)，存在一定的局限性。為提高預(yù)測(cè)的精準(zhǔn)度，以下將利用2018年1月1日至2020年8月31日的云南省昭通市AQI歷史數(shù)據(jù)，建立多元線性回歸、時(shí)間序列、隨機(jī)森林、回歸樹模型對(duì)昭通市AQI進(jìn)行預(yù)測(cè)，選擇預(yù)測(cè)精度較高的兩種模型進(jìn)行組合，進(jìn)一步改進(jìn)空氣質(zhì)量預(yù)測(cè)的方法。

1 方法及原理

1.1 多元線性回歸 多元線性回歸有廣泛的應(yīng)用，可以對(duì)不同變量之間的關(guān)系進(jìn)行分析。Y表示因變量，對(duì)其產(chǎn)生影響的自變量表示為X1，X2，…，Xp，數(shù)量為p，其關(guān)系式為：

ε為隨機(jī)誤差項(xiàng)，β0代表常數(shù)項(xiàng)，βj(j=1,2,…,p)代表Y對(duì)于Xj的回歸系數(shù)參數(shù)值，通常情況下βj代表在其余自變量穩(wěn)定的條件下，Xj增大一個(gè)單位之后Y的變化均值。

1.2 ARIMA（p，d，q）模型 ARIMA（p，d，q）模型為yt=θ0+φ1yt-1+φ2yt-2+…+φpyt-p+εt-θ1εt-1+θ2εt-2+…+θqεt-q，p、q對(duì)應(yīng)著自回歸多項(xiàng)式、移動(dòng)平均多項(xiàng)式的整數(shù)；yt代表樣本值；φ和θ為模型參數(shù)；ε是隨機(jī)誤差，其均值為0。該模型實(shí)際上是在ARMA模型的基礎(chǔ)上進(jìn)行擴(kuò)展得到，針對(duì)非平穩(wěn)時(shí)間序列實(shí)現(xiàn)d階逐期差分，從而提升yt的平穩(wěn)性，然后通過ARMA（p，q）模型進(jìn)行建模、逆變換，通過這種方式可以獲取到原序列。

1.3 隨機(jī)森林 隨機(jī)森林在統(tǒng)計(jì)分析中有廣泛應(yīng)用，其分類器為基礎(chǔ)分類器主要是h(x,θk)。在樣本預(yù)測(cè)過程中需要通過隨機(jī)方式形成一個(gè)森林，在森林中的決策樹結(jié)構(gòu)較多，但是彼此之間的聯(lián)系很少。通過回歸樹的擴(kuò)展可以得到隨機(jī)森林模型，訓(xùn)練集通過隨機(jī)向量分布獲取。已知訓(xùn)練集為代表隨機(jī)的序列，通過隨機(jī)抽樣的方式獲取到T的子集Ti，其中含有的樣本數(shù)目為k，而表示各個(gè)子集對(duì)應(yīng)的模型。如果已知滿足獨(dú)立同分布條件，則通過k輪訓(xùn)練后得到序列，因此隨機(jī)森林預(yù)測(cè)模型為

1.4 回歸樹 分類和回歸樹（CART）算法是將樣本集劃分為兩大部分，也就是對(duì)應(yīng)著兩個(gè)子集，在此過程中采用的是二分遞歸分割技術(shù)，此時(shí)各個(gè)非葉子節(jié)點(diǎn)中均有兩種決策樹方法?；贑ART算法得到的決策樹屬于基本的二叉樹結(jié)構(gòu)。

針對(duì)最佳的切分點(diǎn)、切分變量進(jìn)行確定，可以得到滿足區(qū)間平方差最小的值a1、a2。求解公式為〔7〕：

然后對(duì)全部的輸入變量進(jìn)行遍歷，獲取到最佳的切分變量，對(duì)輸入空間進(jìn)行劃分，可以得到如下兩部分：

按照相同的劃分方式繼續(xù)針對(duì)各個(gè)區(qū)域進(jìn)行劃分，在達(dá)到結(jié)束條件時(shí)輸出最終的結(jié)果，由此得到回歸樹。

1.5 模型評(píng)價(jià)指標(biāo) 在模型評(píng)價(jià)過程中需要利用合適的評(píng)價(jià)指標(biāo)，以下選用平均相對(duì)誤差（MRE）、均方根誤差（RMSE）兩個(gè)指標(biāo)，其公式如下：為觀測(cè)值，為預(yù)測(cè)值，n為樣本數(shù)量。

2 數(shù)據(jù)來源

數(shù)據(jù)源于天氣網(wǎng)中的云南省昭通市2018年1月1日至2020年8月31日AQI，共975個(gè)數(shù)據(jù)，對(duì)AQI的影響因素有多種類型，在收集資料時(shí)無法考慮所有因素，主要篩選了以下指標(biāo)：AQI（Y），PM2.5（X1）、PM10（X2）、SO2（X3）、NO2（X4）、CO（X5）、O3（X6）。

3 預(yù)測(cè)模型的建立

3.1 多元線性回歸 利用最小二乘法建立多元線性回歸模型，同時(shí)需對(duì)參與的影響因素進(jìn)行分析。結(jié)合之前的討論，AQI與PM2.5，PM10，SO2，NO2，CO，O3等多個(gè)因素有關(guān)，這些因素屬于自變量，其中AQI為因變量。選取昭通市2018年1月1日至2020年8月1日共945個(gè)空氣質(zhì)量數(shù)據(jù)作為訓(xùn)練集，2020年8月2日至31日共30個(gè)數(shù)據(jù)為測(cè)試集，使用R語(yǔ)言中的lm（）函數(shù)建立模型，summary（）輸出模型，得到模型為：Y=β0+β1X1+β2X2+…+βP X P+ε，結(jié)果見表1。

表1 多元線性回歸模型

從表1發(fā)現(xiàn)X3和X5不顯著，而其余變量比較顯著，通過逐步回歸法將變量進(jìn)行剔除，輸出模型得到表2。

表2 調(diào)整后的多元線性回歸模型

由此得多元線性回歸方程：

從回歸方程中發(fā)現(xiàn)AQI與PM2.5，PM10，NO2，O3呈正相關(guān)關(guān)系，且O3對(duì)AQI的影響最大，擬合優(yōu)度R2=0.813，即認(rèn)為多元線性回歸模型是有效的，用該模型對(duì)昭通市2020年8月的30個(gè)數(shù)據(jù)進(jìn)行預(yù)測(cè)。

將預(yù)測(cè)值與2020年8月2日至31日的空氣質(zhì)量實(shí)際數(shù)據(jù)進(jìn)行比較，二者分別用實(shí)線、虛線表示，真實(shí)值與預(yù)測(cè)值的擬合圖見圖1。發(fā)現(xiàn)預(yù)測(cè)的AQI波動(dòng)大，結(jié)果誤差大，有可能是因?yàn)橛绊懸蛩氐亩鄻有院湍承┮蜃拥牟豢蓽y(cè)性，使得回歸分析在某些情況下受到限制，真實(shí)值與預(yù)測(cè)值的擬合不太好，但總體趨勢(shì)變化相近，預(yù)測(cè)結(jié)果的平均相對(duì)誤差為9.83%、均方根誤差為5.02。

圖1 多元線性回歸擬合

3.2 時(shí)間序列 繪制出昭通市2018年1月1日至2020年8月31日逐日AQI時(shí)間序列，可知并不存在顯著的特征，因此可以將其認(rèn)為是一個(gè)平穩(wěn)序列，接著通過自相關(guān)圖進(jìn)行分析。見圖2。自相關(guān)圖盡管并未控制在2倍標(biāo)準(zhǔn)差以內(nèi)，但是保持了一定的波動(dòng)性，該特征與平穩(wěn)時(shí)間序列的隨機(jī)性特征是一致的。

圖2 自相關(guān)圖

根據(jù)純隨機(jī)性檢驗(yàn)得到LB檢驗(yàn)統(tǒng)計(jì)量的P=2.2e-16（＜0.05）比較小，據(jù)此可以推斷屬于非白噪聲序列。圖3偏相關(guān)圖顯示拖尾，圖2自相關(guān)圖未快速趨于0，則需進(jìn)一步差分處理，差分后自相關(guān)圖顯示拖尾，可以嘗試擬合ARIMA模型。研究發(fā)現(xiàn)ARIMA（4，1，7）模型的赤道信息準(zhǔn)則（AIC）、貝葉斯準(zhǔn)則（BIC）最低，具體數(shù)值為7 364.51、7 363.11。殘差檢驗(yàn)結(jié)果顯示P=0.956 8＞0.05，可將殘差序列作為白噪聲序列，由此驗(yàn)證了該擬合模型的應(yīng)用效果。

圖3 偏相關(guān)圖

ARIMA模型對(duì)AQI預(yù)測(cè)結(jié)果見圖4，平均相對(duì)誤差為2.20%、均方根誤差14.66。該模型預(yù)測(cè)的AQI難以達(dá)到較高的精度，原因可能是模型沒有抓住原始數(shù)據(jù)“季節(jié)性”的特點(diǎn)，且實(shí)際中有很多因素都會(huì)影響到AQI，包括天氣條件變化等，特別是這些不確定性因素產(chǎn)生的影響具有隨機(jī)性，最終影響了模型的擬合精度。

圖4 時(shí)間序列預(yù)測(cè)圖

3.3 隨機(jī)森林 在隨機(jī)森林建模過程中，利用R語(yǔ)言中的RandomForest程序包，模型如下：r_model=randomForest,（aqi～.,data=dat,ntree=1000,importance=TRUE），ntree代表隨機(jī)森林樹的數(shù)目，這是對(duì)模型精度產(chǎn)生影響的關(guān)鍵因素，importance函數(shù)描述變量的重要性，通過精度平均較少值進(jìn)行評(píng)價(jià)。

模型應(yīng)用效果的影響因素，主要與備選變量（mtry）、隨機(jī)森林樹數(shù)目（ntree）有關(guān)，各個(gè)參數(shù)應(yīng)該保持合適的取值，前者一般取值為森林樹數(shù)目的1∕3±1，如果變量個(gè)數(shù)小于3時(shí)，則將其設(shè)為1，該模型mtry取值3。

通過訓(xùn)練集完成訓(xùn)練，結(jié)果見圖5，橫、縱軸分別對(duì)應(yīng)著隨機(jī)森林樹數(shù)目、均方誤差，均方誤差隨著隨機(jī)森林樹數(shù)目的增大而減少，二者存在一定的負(fù)相關(guān)性，但是減小到一定程度后基本保持穩(wěn)定，此時(shí)隨機(jī)森林樹數(shù)目=400。

圖5 樹目變化情況下誤差變化

利用訓(xùn)練集訓(xùn)練模型，幾個(gè)指標(biāo)在模型中的重要程度見圖6，左邊為各因子相對(duì)重要性圖，右邊為節(jié)點(diǎn)純度圖，數(shù)值越大，表明在模型中的重要性就越明顯，也就是對(duì)AQI的預(yù)測(cè)結(jié)果影響越大。根據(jù)圖中的信息可知，O3濃度和PM10對(duì)AQI預(yù)測(cè)結(jié)果影響最大，NO2對(duì)AQI預(yù)測(cè)結(jié)果影響最小。

圖6 變量重要性

輸入測(cè)試集數(shù)據(jù)，將預(yù)測(cè)值與2020年8月的空氣質(zhì)量實(shí)際數(shù)據(jù)進(jìn)行比較，見圖7，可看到隨機(jī)森林模型的預(yù)測(cè)結(jié)果準(zhǔn)確性較高，模型預(yù)測(cè)值與實(shí)際值基本吻合，在高峰值和低峰值的點(diǎn)，實(shí)際值比預(yù)測(cè)值誤差略大一點(diǎn)。模型擬合優(yōu)度達(dá)到96.95%。預(yù)測(cè)結(jié)果的平均相對(duì)誤差1.50%、均方根誤差0.80。

圖7 隨機(jī)森林?jǐn)M合

3.4 回歸樹 將空氣質(zhì)量數(shù)據(jù)前945個(gè)作為訓(xùn)練集，其余30個(gè)作為測(cè)試集，選用PM2.5，PM10，SO2，NO2，CO，O3作為模型輸入。得到AQI預(yù)報(bào)值，由回歸樹結(jié)果圖展現(xiàn)了空氣質(zhì)量指標(biāo)重要性的排序和分類，可知該回歸樹的根節(jié)點(diǎn)是基于O3來劃分的，O3對(duì)AQI的影響最大。對(duì)根節(jié)點(diǎn)最優(yōu)的切分變量是O3和切分點(diǎn)是106.5，即以O(shè)3≥106.5和O3＜106.5劃分，最后按照此規(guī)律對(duì)下面的子節(jié)點(diǎn)進(jìn)行劃分，直到葉節(jié)點(diǎn)為止。為了達(dá)到模型最優(yōu)化的目標(biāo)，通過廣義交叉驗(yàn)證模型的預(yù)測(cè)效果，選擇了合適的最優(yōu)回歸樹，確保誤差處于最低的要求。利用回歸樹模型預(yù)測(cè)的30組測(cè)試數(shù)據(jù)中有1組達(dá)到了零誤差預(yù)測(cè)，也就是預(yù)測(cè)值與實(shí)際值完全相同；最大絕對(duì)誤差僅第17組樣本，最大誤差值為3.286；其他測(cè)試樣本的誤差值都比較小，即回歸樹模型對(duì)測(cè)試集擬合效果很好，其預(yù)測(cè)精度高且泛化能力也強(qiáng)。

擬合結(jié)果見圖8，最終平均相對(duì)誤差1.40%、均方根誤差為0.85。結(jié)果表明：回歸樹模型較好地預(yù)測(cè)AQI，在最低點(diǎn)、最高點(diǎn)均可以達(dá)到精度的要求，回歸樹的預(yù)測(cè)精度高，穩(wěn)定性好。

圖8 回歸樹擬合

3.5 隨機(jī)森林-回歸樹組合模型 利用單一預(yù)測(cè)模型對(duì)昭通市2020年8月2日至31日的空氣質(zhì)量進(jìn)行預(yù)測(cè)，各自都存在一定的不足，為了增強(qiáng)模型的預(yù)測(cè)性能，在此基礎(chǔ)上設(shè)計(jì)組合模型，提升對(duì)于空氣質(zhì)量預(yù)測(cè)的效果。通過數(shù)據(jù)建立單一預(yù)測(cè)模型發(fā)現(xiàn)多元線性回歸和ARIMA模型預(yù)測(cè)誤差較大，而隨機(jī)森林、回歸樹模型精度較高，故用后者構(gòu)建組合模型預(yù)測(cè)昭通市AQI。運(yùn)用標(biāo)準(zhǔn)差法計(jì)算兩個(gè)預(yù)測(cè)模型的權(quán)重，由隨機(jī)森林模型與回歸樹模型的標(biāo)準(zhǔn)差分別為0.80，0.85，得到權(quán)重分別是w1=0.515、w2=0.485，組合預(yù)測(cè)模型為

Z1為隨機(jī)森林的預(yù)測(cè)值，Z2為回歸樹的預(yù)測(cè)值。

由此組合模型計(jì)算出8月2日至31日AQI的預(yù)測(cè)值，平均相對(duì)誤差為1.27%，均方根誤差為0.42，預(yù)測(cè)的AQI擬合效果見圖9，可以看出兩條線幾乎重合，預(yù)測(cè)效果顯著。

圖9 組合模型預(yù)測(cè)擬合

用平均相對(duì)誤差、均方根誤差檢驗(yàn)單一模型和組合模型的預(yù)測(cè)效果，其結(jié)果見表3。

表3 預(yù)測(cè)模型誤差對(duì)比

可知，相對(duì)于單一的回歸樹以及隨機(jī)森林模型，采用組合模型的預(yù)測(cè)效果更佳，各個(gè)指標(biāo)處于更優(yōu)的水平，均方根誤差低，表明預(yù)測(cè)結(jié)果的可靠性較高；而平均相對(duì)誤差較低，說明預(yù)測(cè)結(jié)果的偏離性最低，預(yù)測(cè)結(jié)果的準(zhǔn)確性更高。因此組合預(yù)測(cè)模型在預(yù)測(cè)性能上能夠達(dá)到更高的要求，相對(duì)于其他單一模型的預(yù)測(cè)效果更佳。

圖10可看出，時(shí)間序列誤差曲線有多個(gè)高峰值，誤差最大，線性回歸誤差比時(shí)間序列更平穩(wěn)一些，但是誤差在0到5之間變化，隨機(jī)森林和回歸樹的誤差很接近，少數(shù)幾個(gè)點(diǎn)預(yù)測(cè)差距明顯，而相對(duì)于單一模型而言，組合模型的預(yù)測(cè)誤差穩(wěn)定性相對(duì)較高，波動(dòng)性不明顯，且誤差在0到2以內(nèi)，預(yù)測(cè)結(jié)果與實(shí)測(cè)值的一致性更高，在實(shí)際預(yù)測(cè)中，能夠?yàn)榭諝赓|(zhì)量預(yù)測(cè)提供較為準(zhǔn)確的結(jié)果。

圖10 模型誤差對(duì)比圖

4 結(jié)論

利用云南省昭通市2018年1月至2020年8月的空氣質(zhì)量監(jiān)測(cè)數(shù)據(jù)構(gòu)建了多元線性回歸、隨機(jī)森林等單一預(yù)測(cè)模型，通過平均相對(duì)誤差、均方根誤差以及數(shù)據(jù)擬合驗(yàn)證了模型的預(yù)測(cè)效果。針對(duì)單一預(yù)測(cè)模型的不足構(gòu)建了組合模型，即隨機(jī)森林-回歸樹組合預(yù)測(cè)模型，該模型融合了兩種單一模型的優(yōu)點(diǎn)，在預(yù)測(cè)時(shí)的精度更高。然后對(duì)2020年8月的AQI進(jìn)行預(yù)測(cè)，組合模型預(yù)測(cè)結(jié)果與實(shí)際值一致性較好，效果遠(yuǎn)遠(yuǎn)優(yōu)于單一模型，因此可以將該組合模型應(yīng)用到實(shí)際空氣質(zhì)量預(yù)測(cè)中。同時(shí)隨機(jī)森林具備了一定的診斷功能，發(fā)現(xiàn)影響昭通市空氣質(zhì)量的主要因素是PM10和O3，這些顆粒物主要來自工業(yè)生產(chǎn)領(lǐng)域中，特別是在化石能源燃燒過程中容易形成大量此類污染物。因此應(yīng)采取有效的空氣污染防治措施，提高空氣的質(zhì)量，從而為居民創(chuàng)造適宜的生活環(huán)境。

致謝：本研究是在2020年云南省教育廳首批線下一流課程《概率論與數(shù)理統(tǒng)計(jì)》課題支持下完成的，在此表示感謝。