郭 楊，丁振軍

(遼寧省生態(tài)環(huán)境監(jiān)測中心，遼寧 沈陽 110161)

1 引言

水資源是人類生存不可或缺的物質(zhì)資源，也是生態(tài)環(huán)境的重要組成部分，水資源水體的質(zhì)量一直受到人們的廣泛關(guān)注。遼河流域是全國七大流域之一，是國家重點治理的“三河三湖”之一。遼河流域的水環(huán)境質(zhì)量的好壞直接影響流域人們的健康。隨著遼寧省水污染治理力度的加大，遼河流域水質(zhì)明顯改善。根據(jù)2020年中國環(huán)境公報顯示，遼河流域為輕度污染，主要污染指標(biāo)為化學(xué)需氧量、高錳酸鹽指數(shù)和五日生化需氧量[1]。在水環(huán)境質(zhì)量管理和控制的工作中，需要了解和掌握水環(huán)境質(zhì)量在未來的變化趨勢，遼河流域水質(zhì)能否繼續(xù)保持良好態(tài)勢，是否存在超標(biāo)風(fēng)險，對全省未來水環(huán)境形勢進行預(yù)判對未來水污染防治工作方向可以提供參考。

2 預(yù)測方法

目前較為常用的水質(zhì)預(yù)測方法有灰色系統(tǒng)理論法[2～4]、人工神經(jīng)網(wǎng)絡(luò)法[5，6]和時間序列分析法[7]等。其中時間序列法是指基于歷史數(shù)據(jù)利用統(tǒng)計方法建立模型預(yù)估未來水質(zhì)情況。時間序列預(yù)測方法分為自回歸模型(AR)、移動平均模型(MA)、自回歸移動平均模型(ARIMA)等，尤其以ARIMA模型最為經(jīng)典，它融合了時間序列和回歸分析優(yōu)點，并且適用于平穩(wěn)時間序列，在預(yù)測分析中應(yīng)用最為廣泛。童俊以氨氮為例，構(gòu)建SVR-ARIMA組合模型預(yù)測了金澤水庫取水口的水質(zhì)[8]。Parmar等使用自回歸整合移動平均時間序列模型預(yù)估了未來的水質(zhì)參數(shù)值，為水質(zhì)管理決策提供基礎(chǔ)[9]，羅學(xué)科等采用ARIMA-SVR組合方法預(yù)測了巢湖水域pH值和溶解氧濃度[10]。顧杰等基于ARIMA模型與BP神經(jīng)網(wǎng)絡(luò)算法對嘉興市河道水質(zhì)進行了預(yù)測[11]。李娜等利用ARIMA預(yù)測了太湖水體Chl-a濃度[12]，江南等采用ARIMA乘積季節(jié)模型擬合序列，建立氟化物時間序列預(yù)測模型，預(yù)測了氟化物月均濃度[13]。謝建輝等利用AMIMA模型預(yù)測了化學(xué)需氧量和氨氮兩個指標(biāo)未來1天的濃度[14]。杜鑫運用ARMA模型對遼河流域東陵大橋斷面COD水質(zhì)變化趨勢進行了預(yù)測，結(jié)果顯示COD呈現(xiàn)增長態(tài)勢[15]?！笆濉逼陂g，遼河流域治理取得較大成就，水生態(tài)環(huán)境質(zhì)量改善明顯。但隨著經(jīng)濟結(jié)構(gòu)的轉(zhuǎn)型升級，水資源開發(fā)利用率高，遼河流域人均水資源量不足全國1/3，水污染防治形式依然十分嚴(yán)峻。為更好掌握未來水質(zhì)變化情況，有必要基于歷史監(jiān)測數(shù)據(jù)對遼河流域水質(zhì)狀況進行預(yù)測，為水污染治理規(guī)劃提供引導(dǎo)方向，減少盲目性，以實現(xiàn)精準(zhǔn)治污。

3 預(yù)測過程

建立ARIMA模型的步驟一般包括確定時序是平穩(wěn)的→找到一個(或幾個)合理的模型(即選定可能的p值和q值)擬合模型→從統(tǒng)計假設(shè)和預(yù)測準(zhǔn)確性等角度評估模型→預(yù)測。本次預(yù)測選取遼河流域部分水質(zhì)類別為Ⅳ類以上的斷面作為研究對象，基于2016～2020年監(jiān)測數(shù)據(jù)，對COD、氨氮、總磷等主要指標(biāo)建立時間序列進行預(yù)測。選取預(yù)測一個斷面的氨氮濃度為例來說明預(yù)測過程。

3.1 驗證序列的平穩(wěn)性

將數(shù)據(jù)轉(zhuǎn)化為時間序列，畫出時序圖判斷平穩(wěn)性，如原始數(shù)據(jù)平穩(wěn)無需對數(shù)據(jù)進行變換，原始數(shù)據(jù)不平穩(wěn)對數(shù)據(jù)進行差分，觀察差分后數(shù)據(jù)是否為平穩(wěn)。在R語言中加載forecast包和tidyverse包：library(forecast)library(tidyverse)，讀取斷面數(shù)據(jù)：water <-read.csv("斷面.csv",sep = ",",header = T)，將數(shù)據(jù)轉(zhuǎn)化為時間序列格式：waterts <-ts(water$氨氮,start=c(2015,1),end=c(2020,12),frequency=12)，查看氨氮數(shù)據(jù)趨勢圖：plot.ts(waterts)(圖1、圖2)。

圖1 原始數(shù)據(jù)時序圖

圖2 數(shù)據(jù)一階差分后時序圖

原始數(shù)據(jù)呈現(xiàn)出兩個先升后降的趨勢。需對原始數(shù)據(jù)進行差分，差分階數(shù)為1：waterdiff <-diff(waterts, differences=1)，查看差分后數(shù)據(jù)趨勢：plot.ts(waterdiff)。原始數(shù)據(jù)一階差分后數(shù)據(jù)點似乎呈現(xiàn)出平穩(wěn)分布狀態(tài)。

3.2 選擇模型

平穩(wěn)后查看ACF和PACF的結(jié)果：acf(waterdiff,lag.max = 30)；pacf(waterdiff,lag.max = 30)(圖3、圖4)。

圖3 ACF

圖4 PACF

一般結(jié)合ACF和PACF圖選擇參數(shù)p和q。我們通過程序包中的auto.arima()函數(shù)自動計算合適的(p，d，q)值：auto.arima(waterts,trace=T)。根據(jù)auto.arima()計算結(jié)果，選用(1,1,1)非季節(jié)性模型進行參數(shù)估計：waterarima1 <-arima(waterts,order=c(1,1,1),method="ML")。

3.3 模型檢驗

對殘差序列進行白噪聲檢驗：Box.test(waterarima1$residuals,type="Box-Pierce",lag = 1)；Box-Pierce test；data:waterarima1$residuals；X-squared=0.17719,df=1,p-value=0.6738，得出p=0.6738>0.05,殘差序列白噪聲檢驗說明，模型顯著成立，ARIMA(1,1,1)模型對該時間序列擬合成功。

3.4 預(yù)測

采用forecast()函數(shù)預(yù)測1年后數(shù)據(jù)：waterforecast <-forecast(waterarima1,h=12)，并將預(yù)測結(jié)果作圖：plot(waterforecast)，預(yù)測顯示為藍線，淺灰色和深灰色區(qū)域分別代表95%置信區(qū)間河80%置信區(qū)間(圖5)。

圖5 氨氮預(yù)測結(jié)果

4 預(yù)測結(jié)果

ARIMA(1,1,1)能夠用于氨氮濃度預(yù)測?；谝陨蠗l件對選取的全部斷面的其他主要指標(biāo)濃度進行了預(yù)測。結(jié)果顯示，在經(jīng)濟環(huán)境、自然條件狀況及環(huán)保治理力度不變的情況下，未來一年預(yù)測的斷面中，COD、BOD5、高錳酸鹽指數(shù)濃度總體保持穩(wěn)定，個別斷面氨氮、總磷濃度存在升高的趨勢，在“十四五”期間應(yīng)重點關(guān)注主要污染物濃度存在上升趨勢的斷面，提前做好相關(guān)管控計劃，避免水質(zhì)反彈惡化，實現(xiàn)精準(zhǔn)施策。

5 結(jié)語

本次研究采用ARIMA模型對遼河流域部分?jǐn)嗝娴闹饕笜?biāo)未來濃度進行了預(yù)測，取得較好結(jié)果。水環(huán)境防治的重要依據(jù)就是對水質(zhì)的污染情況進行預(yù)測并制定相應(yīng)的治理措施，達到提前對水環(huán)境進行變化趨勢分析的效果。ARIMA是對單因素短期預(yù)測效果較好，但在實際中，一方面，河流水質(zhì)受水體內(nèi)部各物質(zhì)之間存在相互影響；另一方面，遼河流域水質(zhì)也受降水等自然情況及水資源開發(fā)利用、環(huán)保治理力度等人為因素的影響，因此在河流水質(zhì)情況預(yù)測時要充分考慮各因素的影響來篩選出合適的預(yù)測方法。