余婉風 呂科 劉洋 朱偉杰



摘要:在近幾年全國空氣質量總體有所好轉的大環境下,本文通過分析蕪湖市空氣質量指數數據,探究蕪湖市空氣質量現狀,并構建AQI短期預測模型,為蕪湖市政府控制環境污染和有效地治理提供科學的依據。2013年12月1日—2020年10月31日近8年蕪湖市空氣質量指數(AQI)數據作為研究對象,R語言為實現工具。首先,分析AQI數據曲線圖,采用非參數檢驗Kruskal-Wallis檢驗比較這8年AQI數據是否具有顯著性差異;其次,根據對AQI時間序列平穩性分析結果,選擇合理的時間序列模型—ARIMA模型,估計模型參數,建立擬合模型,并評價模型有效性;最后,利用模型預測未來幾個月AQI。
關鍵詞:AQI;非參數檢驗;時間序列分析;ARIMA預測模型; R語言
中圖分類號:TP391? ? ? 文獻標識碼:A
文章編號: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
隨著城市工業的成長壯大,空氣污染嚴重,空氣質量惡化,不僅影響到人們的正常生活,而且威脅著人們的身心健康[1]。環保作為國家戰略性重點產業,全國各級政府對本省市的環保工作高度重視,紛紛積極推動節能減排和環境治理工作,截至目前,全國幾乎所有的省市,均已出臺生態保護相關政策、資金支持或項目管理方案,為我國全面推進環保事業提供有力的支持[2-3]。蕪湖市政府以科學發展觀為指導,緊緊圍繞國家環保總局提出的七個方面整改要求,開展集中整治行動,推進環保執法,推進節能減排,推進全市經濟發展模式的轉變[4]。
同時,隨著R語言中各類工具包(package)的不斷涌現,數理統計分析工具和可視化方案的R程序包得以開發應用開來。本研究數據來自2013年-2020年日報AQI數據,非參數檢驗Kruskal-Wallis檢驗可以對多組獨立樣本的多重比較,安裝、引用R語言pgirmess程序包,調用kruskalmc函數實現Kruskal-Wallis檢驗[5];利用R語言軟件包tseries進行時間序列分析;調用軟件包forecast的auto.arima()函數進行系統定階,擬合最優ARIMA模型。本文中,運用非參數檢驗對比近幾年蕪湖市空氣質量指數有無明顯改善;應用時間序列分析建立預測模型,對模型進行參數估計,診斷和評價,確定最優模型,并檢驗模型有效性[6],最后運用模型進行預測。
1 非參數檢驗—Kruskal-Wallis檢驗
1.1抽樣數據
研究對象是蕪湖市2013年12月到2020年7月的空氣質量指數(AQI),取每個月AQI平均值。
為了更直觀地分析2013年到2020年AQI的變化趨勢,描繪數據趨勢圖。如圖1所示,2013年12月到2020年10月蕪湖市空氣質量指數(AQI)有逐年降低趨勢,但不明顯,且可能呈周期性變化。僅從圖1無法確定AQI是否有明顯的逐年降低趨勢,且考慮到AQI分布未知,樣本量少,采用非參數檢驗進行8組獨立樣本的多重比較,判斷不同年份的AQI數據是否具有統計學差異。
1.2 R語言實現Kruskal-Wallis秩和檢驗
利用抽樣技術,從2013年12月到2020年10月中抽取每個月1號的數據,即83個樣本數據,不同年份作為組別,共有8個獨立組別。安裝并加載程序包pgirmess,調用Kruskal-Wallis秩和檢驗函數計算結果p-value=0.02139,顯著性水平α=0.05,p-value小于顯著性水平,說明拒絕原假設(原假設是8組樣本數據無顯著性差異),即這8年的AQI數據有顯著性差異。結合圖1分析,可以認為近幾年蕪湖市的空氣質量指數有明顯的降低趨勢,空氣質量明顯好轉,環境治理效果顯著。為了給蕪湖市政府控制空氣污染和有效地治理提供科學的依據,分析數據特征,建立有效的AQI預測模型。
2 AQI時間序列預測模型
從2013年1月到2020年10月蕪湖市空氣質量指數(AQI)總趨勢雖無明顯降低,但AQI不是雜亂無章的,如圖2所示,分解了數據的變化趨勢、季節性和不確定性因素,分析圖中“seasonal”曲線, AQI有明顯的季節性、周期性趨勢。基于AQI的數據特征,采用時間序列分析建立AQI短期預測模型,在預測模型有效的前提下利用模型預測2020年未來幾個月的空氣質量指數,并繪制擬合圖形。
2.1 平穩性分析
(1)平穩性分析
常見的時間序列模型包括ARAM和ARIMA模型等,根據時間序列的平穩性來選擇合適的預測模型,常見的時間序列平穩性檢驗方法有PP檢驗法和ADF檢驗法。
R軟件安裝加載tseries包,進行PP檢驗和ADF檢驗。PP檢驗結果中出現警告信息,ADF檢驗p值=0.04803< 0.05,默認顯著性水平[α=0.05],拒絕原假設(原假設認為時間序列是非平穩的),可以認為AQI數據是平穩的時間序列。
(2)白噪聲檢驗
對于平穩的時間序列需要進行白噪聲檢驗,因為白噪聲是純隨機序列,對純隨機序列建模毫無意義。對AQI時間序列白噪聲檢驗結果p值=4.6637e-07小于顯著性水平0.05,拒絕原假設,可以認為AQI數據為非白噪聲時間序列,對該平穩時間序列建模有意義。
考慮到AQI呈季節性、周期性變化,選擇目前最常用的擬合平穩序列的模型ARAM(p, q)模型(自回歸移動平均模型)作為AQI時間序列預測模型。ARMA模型(自回歸移動平均模型)又可以細分為AR(p)模型、MA(q)模型和ARMA(p, q)模型三大類,確定AQI預測模型屬于哪一類這都取決于階數p和q的值,其一般準則如表1所示。
2.2 確定ARAM(p, q)模型階數p和q
對于ARAM(p, q)模型,一般通過分析ACF自相關圖和PACF偏相關圖來估計兩個未知參數p和q。在R軟件中,擬合線性回歸模型,再調用acf()和pacf()兩個函數繪制ACF自相關圖和PACF偏相關圖。
圖3中,自相關系數拖尾,且從1開始控制在置信區間之內,p= 1;圖4中,偏相關系數拖尾,且從2開始控制在置信區間之內,q= 2。根據表1確定ARAM模型準則,初步推測AQI預測最優模型為ARMA(1, 2)。
2.3 系統自動定階
為了驗證通過時間序列穩定性、ACF自相關圖和PACF偏相關圖分析確定的ARAM(1, 2)模型是否比較合理,利用R語言軟件包forecast的ARIMA模型函數,默認預測模型為ARIMA(p, d, q)模型(差分整合移動平均自回歸模型)進行系統自動定階,與通過平穩時間序列ACF自相關圖和PACF偏相關圖分析確定的ARAM(1, 2)模型進行比較,找到最優模型。安裝、加載forecast包,時間序列ARIMA建模,結果如表2,ARIMA(1,0,0)(2,0,0)12 是擬合AQI時間序列的最佳模型,其中參數p= 1,d= 0,q= 2。參數d表示差分階數,當d為0時,ARIMA模型就等同于ARMA模型,即ARIMA(1, 0, 2)模型與ARMA(1, 2)模型是等價的。由此,通過分析自相關系數圖和偏相關系數圖確定的最優模型與R語言進行系統自動定階結果一致。
在確定了最優模型ARIMA(1,0,2)模型之后,并采用極大似然估計思想進行參數估計,如表2,ARIMA(1,0,2)擬合模型:
[yt=0.6101yt-1+0.3725εt-1+0.2043εt-2+84.4825]
接下來,對ARIMA(1,0,2)模型的3個系數和1個截距進行假設檢驗,用極大似然估計思想估計的系數的絕對值除以其標準差(s.e.)得到的商與t檢驗5%的臨界值1.96比較,商的絕對值大于1.96,拒絕原假設(原假設為參數影響不顯著),否則認為參數影響顯著。t(ar1)= 6.532> 1.96, t(sar1)= 3.242> 1.96,t(sar2)= 1.625< 1.96, t(mean)= 9.621> 1.96,所以參數ar1、sar1和截距對該時間序列模型都具有顯著影響。
2.4 預測模型有效性檢驗
常用的檢驗ARIMA模型有效性的方法是純隨機性檢驗方法Box-Ljung檢驗,R語言提供了Box-Ljung檢驗的tsdiag()函數。調用tsdiag()函數,檢驗結果如圖5所示。
第二行的ACF檢驗說明殘差沒有明顯的自相關性;第三行的Box-Ljung檢驗顯示所有的p值都大于顯著性水平0.01,殘差序列不能拒絕純隨機的原假設,說明殘差為白噪聲,所以此預測模型是有效模型。利用ARIMA模型預測過去8年的AQI數據,將預測值與真實值對比,如表3所示,該預測模型在誤差允許的范圍內能反映AQI數據的變化趨勢和規律,有參考價值。
3 ARIMA模型預測未來空氣質量指數
在檢驗ARIMA(1, 0, 2) 預測模型是有效的前提下,利用預測2020年未來幾個月的空氣質量。調用forecast()函數繪制預測模型曲線,如圖6所示,兩條曲線分別表示AQI時間序列和ARIMA(1, 0, 2)模型,右邊粗線表示預測模型預測的未來幾個月AQI序列。
調用forecast()函數預測11月和12月每月平均的空氣質量指數,有:
[predict(AQI—11)=68.67];
[predict(AQI—12)=78.73]。
4 總結
時間序列模型的缺點是精確度不夠高,為了提高精確度,可以在模型中增加與空氣質量相關的變量如PM2.5、CO和N02含量等相關變量。考慮到目前蕪湖市政府正在逐漸加大環保力度,從長期來看,由于外部因素干預加強,導致模型的預測能力下降,所以此模型較適合短期內預測。
參考文獻:
[1] 史美義.淺析當前我國城市環境污染的現狀及原因[J].科技信息,2012(18):79.
[2] 國務院.國務院關于落實科學發展觀加強環境保護的決定[J].中國環境監測,2006,22(1):1-6.
[3] 安徽省人民政府.安徽省人民政府貫徹國務院關于落實科學發展觀加強環境保護決定的實施意見[J].安徽省人民政府公報,2006(21):21-26.
[4] 付偉,司紅君,盧堯,等.蕪湖市空氣質量特征及其受氣象要素的影響分析[C]//第35屆中國氣象學會年會論文集.合肥,2018:128-129.
[5] 金英良,趙華碩,孫桂香,等.基于R軟件的多組獨立樣本秩和檢驗的多重比較[J].預防醫學論壇,2016,22(11):805-806,809.
[6] 牟敬鋒,趙星,樊靜潔,等.基于ARIMA模型的深圳市空氣質量指數時間序列預測研究[J].環境衛生學雜志,2017,7(2):102-107,117.
【通聯編輯:唐一東】