余婉風(fēng) 呂科 劉洋 朱偉杰

摘要：在近幾年全國(guó)空氣質(zhì)量總體有所好轉(zhuǎn)的大環(huán)境下，本文通過(guò)分析蕪湖市空氣質(zhì)量指數(shù)數(shù)據(jù)，探究蕪湖市空氣質(zhì)量現(xiàn)狀，并構(gòu)建AQI短期預(yù)測(cè)模型，為蕪湖市政府控制環(huán)境污染和有效地治理提供科學(xué)的依據(jù)。2013年12月1日—2020年10月31日近8年蕪湖市空氣質(zhì)量指數(shù)（AQI）數(shù)據(jù)作為研究對(duì)象，R語(yǔ)言為實(shí)現(xiàn)工具。首先，分析AQI數(shù)據(jù)曲線圖，采用非參數(shù)檢驗(yàn)Kruskal-Wallis檢驗(yàn)比較這8年AQI數(shù)據(jù)是否具有顯著性差異;其次，根據(jù)對(duì)AQI時(shí)間序列平穩(wěn)性分析結(jié)果，選擇合理的時(shí)間序列模型—ARIMA模型，估計(jì)模型參數(shù)，建立擬合模型，并評(píng)價(jià)模型有效性;最后，利用模型預(yù)測(cè)未來(lái)幾個(gè)月AQI。

關(guān)鍵詞：AQI;非參數(shù)檢驗(yàn);時(shí)間序列分析;ARIMA預(yù)測(cè)模型; R語(yǔ)言

中圖分類號(hào)：TP391? ? ? 文獻(xiàn)標(biāo)識(shí)碼：A

文章編號(hào)：1009-3044（2021）11-0239-03

Study on Air Quality Index（AQI） of? Wuhu City Based on Nonparametric Test and ARIMA Model

YU Wan-feng， LV Ke， LIU Yang， ZHU Wei-jie

（College of Big Data and Artificial Intelligence， Anhui Institute of Information Technology， Wuhu 241003，China）

Abstract： In recent years， air conditions nationwide has been improving， through analyzing the air quality index data of Wuhu city， to explore the current situation of air conditions of Wuhu City， and model for AQI prediction ， which provides scientific basis for Wuhu city government to control environmental pollution effectively. From December 1， 2013 to October 31， 2020， air quality Index （AQI） data of Wuhu city in recent 8 years were taken as the research object， and R language was used as the implementation tool. Firstly， the AQI data graph was analyzed， and the kruskal-Wallis test was used to compare the significant differences of the AQI data over the past 8 years. Secondly， according to the stationary analysis of AQI time series， an reasonable time series model -- ARIMA model is selected， to estimate fitting model parameters and evaluate the effectiveness of the model. Finally， the ARIMA model is used to predict AQI in the coming months.

Key words： AQI; nonparametric tests; Time-Series analysis; ARIMA prediction model; R softwre

隨著城市工業(yè)的成長(zhǎng)壯大，空氣污染嚴(yán)重，空氣質(zhì)量惡化，不僅影響到人們的正常生活，而且威脅著人們的身心健康[1]。環(huán)保作為國(guó)家戰(zhàn)略性重點(diǎn)產(chǎn)業(yè)，全國(guó)各級(jí)政府對(duì)本省市的環(huán)保工作高度重視，紛紛積極推動(dòng)節(jié)能減排和環(huán)境治理工作，截至目前，全國(guó)幾乎所有的省市，均已出臺(tái)生態(tài)保護(hù)相關(guān)政策、資金支持或項(xiàng)目管理方案，為我國(guó)全面推進(jìn)環(huán)保事業(yè)提供有力的支持[2-3]。蕪湖市政府以科學(xué)發(fā)展觀為指導(dǎo)，緊緊圍繞國(guó)家環(huán)?？偩痔岢龅钠邆€(gè)方面整改要求，開(kāi)展集中整治行動(dòng)，推進(jìn)環(huán)保執(zhí)法，推進(jìn)節(jié)能減排，推進(jìn)全市經(jīng)濟(jì)發(fā)展模式的轉(zhuǎn)變[4]。

同時(shí)，隨著R語(yǔ)言中各類工具包（package）的不斷涌現(xiàn)，數(shù)理統(tǒng)計(jì)分析工具和可視化方案的R程序包得以開(kāi)發(fā)應(yīng)用開(kāi)來(lái)。本研究數(shù)據(jù)來(lái)自2013年-2020年日?qǐng)?bào)AQI數(shù)據(jù)，非參數(shù)檢驗(yàn)Kruskal-Wallis檢驗(yàn)可以對(duì)多組獨(dú)立樣本的多重比較，安裝、引用R語(yǔ)言pgirmess程序包，調(diào)用kruskalmc函數(shù)實(shí)現(xiàn)Kruskal-Wallis檢驗(yàn)[5];利用R語(yǔ)言軟件包tseries進(jìn)行時(shí)間序列分析;調(diào)用軟件包forecast的auto.arima（）函數(shù)進(jìn)行系統(tǒng)定階，擬合最優(yōu)ARIMA模型。本文中，運(yùn)用非參數(shù)檢驗(yàn)對(duì)比近幾年蕪湖市空氣質(zhì)量指數(shù)有無(wú)明顯改善;應(yīng)用時(shí)間序列分析建立預(yù)測(cè)模型，對(duì)模型進(jìn)行參數(shù)估計(jì)，診斷和評(píng)價(jià)，確定最優(yōu)模型，并檢驗(yàn)?zāi)Ｐ陀行訹6]，最后運(yùn)用模型進(jìn)行預(yù)測(cè)。

1 非參數(shù)檢驗(yàn)—Kruskal-Wallis檢驗(yàn)

1.1抽樣數(shù)據(jù)

研究對(duì)象是蕪湖市2013年12月到2020年7月的空氣質(zhì)量指數(shù)（AQI），取每個(gè)月AQI平均值。

為了更直觀地分析2013年到2020年AQI的變化趨勢(shì)，描繪數(shù)據(jù)趨勢(shì)圖。如圖1所示，2013年12月到2020年10月蕪湖市空氣質(zhì)量指數(shù)（AQI）有逐年降低趨勢(shì)，但不明顯，且可能呈周期性變化。僅從圖1無(wú)法確定AQI是否有明顯的逐年降低趨勢(shì)，且考慮到AQI分布未知，樣本量少，采用非參數(shù)檢驗(yàn)進(jìn)行8組獨(dú)立樣本的多重比較，判斷不同年份的AQI數(shù)據(jù)是否具有統(tǒng)計(jì)學(xué)差異。

1.2 R語(yǔ)言實(shí)現(xiàn)Kruskal-Wallis秩和檢驗(yàn)

利用抽樣技術(shù)，從2013年12月到2020年10月中抽取每個(gè)月1號(hào)的數(shù)據(jù)，即83個(gè)樣本數(shù)據(jù)，不同年份作為組別，共有8個(gè)獨(dú)立組別。安裝并加載程序包pgirmess，調(diào)用Kruskal-Wallis秩和檢驗(yàn)函數(shù)計(jì)算結(jié)果p-value=0.02139，顯著性水平α=0.05，p-value小于顯著性水平，說(shuō)明拒絕原假設(shè)（原假設(shè)是8組樣本數(shù)據(jù)無(wú)顯著性差異），即這8年的AQI數(shù)據(jù)有顯著性差異。結(jié)合圖1分析，可以認(rèn)為近幾年蕪湖市的空氣質(zhì)量指數(shù)有明顯的降低趨勢(shì)，空氣質(zhì)量明顯好轉(zhuǎn)，環(huán)境治理效果顯著。為了給蕪湖市政府控制空氣污染和有效地治理提供科學(xué)的依據(jù)，分析數(shù)據(jù)特征，建立有效的AQI預(yù)測(cè)模型。

2 AQI時(shí)間序列預(yù)測(cè)模型

從2013年1月到2020年10月蕪湖市空氣質(zhì)量指數(shù)（AQI）總趨勢(shì)雖無(wú)明顯降低，但AQI不是雜亂無(wú)章的，如圖2所示，分解了數(shù)據(jù)的變化趨勢(shì)、季節(jié)性和不確定性因素，分析圖中“seasonal”曲線， AQI有明顯的季節(jié)性、周期性趨勢(shì)。基于AQI的數(shù)據(jù)特征，采用時(shí)間序列分析建立AQI短期預(yù)測(cè)模型，在預(yù)測(cè)模型有效的前提下利用模型預(yù)測(cè)2020年未來(lái)幾個(gè)月的空氣質(zhì)量指數(shù)，并繪制擬合圖形。

2.1 平穩(wěn)性分析

（1）平穩(wěn)性分析

常見(jiàn)的時(shí)間序列模型包括ARAM和ARIMA模型等，根據(jù)時(shí)間序列的平穩(wěn)性來(lái)選擇合適的預(yù)測(cè)模型，常見(jiàn)的時(shí)間序列平穩(wěn)性檢驗(yàn)方法有PP檢驗(yàn)法和ADF檢驗(yàn)法。

R軟件安裝加載tseries包，進(jìn)行PP檢驗(yàn)和ADF檢驗(yàn)。PP檢驗(yàn)結(jié)果中出現(xiàn)警告信息，ADF檢驗(yàn)p值=0.04803< 0.05，默認(rèn)顯著性水平[α=0.05]，拒絕原假設(shè)（原假設(shè)認(rèn)為時(shí)間序列是非平穩(wěn)的），可以認(rèn)為AQI數(shù)據(jù)是平穩(wěn)的時(shí)間序列。

（2）白噪聲檢驗(yàn)

對(duì)于平穩(wěn)的時(shí)間序列需要進(jìn)行白噪聲檢驗(yàn)，因?yàn)榘自肼暿羌冸S機(jī)序列，對(duì)純隨機(jī)序列建模毫無(wú)意義。對(duì)AQI時(shí)間序列白噪聲檢驗(yàn)結(jié)果p值=4.6637e-07小于顯著性水平0.05，拒絕原假設(shè)，可以認(rèn)為AQI數(shù)據(jù)為非白噪聲時(shí)間序列，對(duì)該平穩(wěn)時(shí)間序列建模有意義。

考慮到AQI呈季節(jié)性、周期性變化，選擇目前最常用的擬合平穩(wěn)序列的模型ARAM（p， q）模型（自回歸移動(dòng)平均模型）作為AQI時(shí)間序列預(yù)測(cè)模型。ARMA模型（自回歸移動(dòng)平均模型）又可以細(xì)分為AR（p）模型、MA（q）模型和ARMA（p， q）模型三大類，確定AQI預(yù)測(cè)模型屬于哪一類這都取決于階數(shù)p和q的值，其一般準(zhǔn)則如表1所示。

2.2 確定ARAM（p， q）模型階數(shù)p和q

對(duì)于ARAM（p， q）模型，一般通過(guò)分析ACF自相關(guān)圖和PACF偏相關(guān)圖來(lái)估計(jì)兩個(gè)未知參數(shù)p和q。在R軟件中，擬合線性回歸模型，再調(diào)用acf（）和pacf（）兩個(gè)函數(shù)繪制ACF自相關(guān)圖和PACF偏相關(guān)圖。

圖3中，自相關(guān)系數(shù)拖尾，且從1開(kāi)始控制在置信區(qū)間之內(nèi)，p= 1;圖4中，偏相關(guān)系數(shù)拖尾，且從2開(kāi)始控制在置信區(qū)間之內(nèi)，q= 2。根據(jù)表1確定ARAM模型準(zhǔn)則，初步推測(cè)AQI預(yù)測(cè)最優(yōu)模型為ARMA（1， 2）。

2.3 系統(tǒng)自動(dòng)定階

為了驗(yàn)證通過(guò)時(shí)間序列穩(wěn)定性、ACF自相關(guān)圖和PACF偏相關(guān)圖分析確定的ARAM（1， 2）模型是否比較合理，利用R語(yǔ)言軟件包forecast的ARIMA模型函數(shù)，默認(rèn)預(yù)測(cè)模型為ARIMA（p， d， q）模型（差分整合移動(dòng)平均自回歸模型）進(jìn)行系統(tǒng)自動(dòng)定階，與通過(guò)平穩(wěn)時(shí)間序列ACF自相關(guān)圖和PACF偏相關(guān)圖分析確定的ARAM（1， 2）模型進(jìn)行比較，找到最優(yōu)模型。安裝、加載forecast包，時(shí)間序列ARIMA建模，結(jié)果如表2，ARIMA（1，0，0）（2，0，0）12 是擬合AQI時(shí)間序列的最佳模型，其中參數(shù)p= 1，d= 0，q= 2。參數(shù)d表示差分階數(shù)，當(dāng)d為0時(shí)，ARIMA模型就等同于ARMA模型，即ARIMA（1， 0， 2）模型與ARMA（1， 2）模型是等價(jià)的。由此，通過(guò)分析自相關(guān)系數(shù)圖和偏相關(guān)系數(shù)圖確定的最優(yōu)模型與R語(yǔ)言進(jìn)行系統(tǒng)自動(dòng)定階結(jié)果一致。

在確定了最優(yōu)模型ARIMA（1，0，2）模型之后，并采用極大似然估計(jì)思想進(jìn)行參數(shù)估計(jì)，如表2，ARIMA（1，0，2）擬合模型：

[yt=0.6101yt-1+0.3725εt-1+0.2043εt-2+84.4825]

接下來(lái)，對(duì)ARIMA（1，0，2）模型的3個(gè)系數(shù)和1個(gè)截距進(jìn)行假設(shè)檢驗(yàn)，用極大似然估計(jì)思想估計(jì)的系數(shù)的絕對(duì)值除以其標(biāo)準(zhǔn)差（s.e.）得到的商與t檢驗(yàn)5%的臨界值1.96比較，商的絕對(duì)值大于1.96，拒絕原假設(shè)（原假設(shè)為參數(shù)影響不顯著），否則認(rèn)為參數(shù)影響顯著。t（ar1）= 6.532> 1.96， t（sar1）= 3.242> 1.96，t（sar2）= 1.625< 1.96， t（mean）= 9.621> 1.96，所以參數(shù)ar1、sar1和截距對(duì)該時(shí)間序列模型都具有顯著影響。

2.4 預(yù)測(cè)模型有效性檢驗(yàn)

常用的檢驗(yàn)ARIMA模型有效性的方法是純隨機(jī)性檢驗(yàn)方法Box-Ljung檢驗(yàn)，R語(yǔ)言提供了Box-Ljung檢驗(yàn)的tsdiag（）函數(shù)。調(diào)用tsdiag（）函數(shù)，檢驗(yàn)結(jié)果如圖5所示。

第二行的ACF檢驗(yàn)說(shuō)明殘差沒(méi)有明顯的自相關(guān)性;第三行的Box-Ljung檢驗(yàn)顯示所有的p值都大于顯著性水平0.01，殘差序列不能拒絕純隨機(jī)的原假設(shè)，說(shuō)明殘差為白噪聲，所以此預(yù)測(cè)模型是有效模型。利用ARIMA模型預(yù)測(cè)過(guò)去8年的AQI數(shù)據(jù)，將預(yù)測(cè)值與真實(shí)值對(duì)比，如表3所示，該預(yù)測(cè)模型在誤差允許的范圍內(nèi)能反映AQI數(shù)據(jù)的變化趨勢(shì)和規(guī)律，有參考價(jià)值。

3 ARIMA模型預(yù)測(cè)未來(lái)空氣質(zhì)量指數(shù)

在檢驗(yàn)ARIMA（1， 0， 2） 預(yù)測(cè)模型是有效的前提下，利用預(yù)測(cè)2020年未來(lái)幾個(gè)月的空氣質(zhì)量。調(diào)用forecast（）函數(shù)繪制預(yù)測(cè)模型曲線，如圖6所示，兩條曲線分別表示AQI時(shí)間序列和ARIMA（1， 0， 2）模型，右邊粗線表示預(yù)測(cè)模型預(yù)測(cè)的未來(lái)幾個(gè)月AQI序列。

調(diào)用forecast（）函數(shù)預(yù)測(cè)11月和12月每月平均的空氣質(zhì)量指數(shù)，有：

[predict（AQI—11）=68.67];

[predict（AQI—12）=78.73]。

4 總結(jié)

時(shí)間序列模型的缺點(diǎn)是精確度不夠高，為了提高精確度，可以在模型中增加與空氣質(zhì)量相關(guān)的變量如PM2.5、CO和N02含量等相關(guān)變量。考慮到目前蕪湖市政府正在逐漸加大環(huán)保力度，從長(zhǎng)期來(lái)看，由于外部因素干預(yù)加強(qiáng)，導(dǎo)致模型的預(yù)測(cè)能力下降，所以此模型較適合短期內(nèi)預(yù)測(cè)。

參考文獻(xiàn)：

[1] 史美義.淺析當(dāng)前我國(guó)城市環(huán)境污染的現(xiàn)狀及原因[J].科技信息，2012（18）：79.

[2] 國(guó)務(wù)院.國(guó)務(wù)院關(guān)于落實(shí)科學(xué)發(fā)展觀加強(qiáng)環(huán)境保護(hù)的決定[J].中國(guó)環(huán)境監(jiān)測(cè)，2006，22（1）：1-6.

[3] 安徽省人民政府.安徽省人民政府貫徹國(guó)務(wù)院關(guān)于落實(shí)科學(xué)發(fā)展觀加強(qiáng)環(huán)境保護(hù)決定的實(shí)施意見(jiàn)[J].安徽省人民政府公報(bào)，2006（21）：21-26.

[4] 付偉，司紅君，盧堯，等.蕪湖市空氣質(zhì)量特征及其受氣象要素的影響分析[C]//第35屆中國(guó)氣象學(xué)會(huì)年會(huì)論文集.合肥，2018：128-129.

[5] 金英良，趙華碩，孫桂香，等.基于R軟件的多組獨(dú)立樣本秩和檢驗(yàn)的多重比較[J].預(yù)防醫(yī)學(xué)論壇，2016，22（11）：805-806，809.

[6] 牟敬鋒，趙星，樊靜潔，等.基于ARIMA模型的深圳市空氣質(zhì)量指數(shù)時(shí)間序列預(yù)測(cè)研究[J].環(huán)境衛(wèi)生學(xué)雜志，2017，7（2）：102-107，117.

【通聯(lián)編輯：唐一東】