王圣偉, 李 萍, 婁天瀧, 綻玉林, 李鴻鴻
(西北師范大學(xué) 計(jì)算機(jī)科學(xué)與工程學(xué)院,甘肅 蘭州 730070)
隨著我國(guó)經(jīng)濟(jì)的迅速發(fā)展、工業(yè)化和城市化的進(jìn)程不斷加快,在消耗水資源的同時(shí)也排放大量的污水和其他污染物,其中對(duì)環(huán)境的污染中重金屬污染尤為突出。因此精準(zhǔn)預(yù)測(cè)流域重金屬污染程度有著重要的研究意義。法國(guó)的學(xué)者Ruelle與荷蘭的學(xué)者Takens[1]在他們共同發(fā)表的《論湍流的本質(zhì)》中首次對(duì)湍流使用混沌理論來(lái)進(jìn)行描述,并在研究過(guò)程中通過(guò)數(shù)學(xué)論證的方法證明了動(dòng)力系統(tǒng)中“奇怪吸引子”的存在。1973年,日本學(xué)者發(fā)現(xiàn)一種雜亂的振動(dòng)形態(tài)—Ueda吸引子[2]。1975年,李天巖和Yorke共同給出了在閉區(qū)間上連續(xù)的自映射的混沌定義[3]。而人工神經(jīng)網(wǎng)絡(luò)(artificial neural network,ANN)預(yù)測(cè)土壤重金屬含量與混沌相空間相結(jié)合研究較為鮮見,在過(guò)去幾十年,ANN被廣泛應(yīng)用于不同環(huán)境工程領(lǐng)域[4~6]。實(shí)際應(yīng)用中流域重金屬含量的高低受到溫度、降水、pH值等多因素的共同影響,單一依靠重金屬的含量去預(yù)測(cè)重金屬含量的方法準(zhǔn)確性較低,所以本文根據(jù)影響重金屬含量的多種因素構(gòu)建一個(gè)多元時(shí)間序列的混沌系統(tǒng)。通過(guò)構(gòu)建多變量混沌時(shí)間序列,改進(jìn)神經(jīng)網(wǎng)絡(luò)長(zhǎng)短期記憶(long short-term memory,LSTM)模型,準(zhǔn)確挖掘變量間內(nèi)在聯(lián)系是進(jìn)行污染物預(yù)測(cè)的研究方向。
本文通過(guò)多元混沌時(shí)間序列的重構(gòu)相空間構(gòu)建,利用相空間(phase space,PS)重構(gòu)和增加窺視孔連接的LSTM循環(huán)神經(jīng)網(wǎng)絡(luò)LSTM時(shí)間序列預(yù)測(cè)模型相結(jié)合,提出多元混沌時(shí)間序列的PS-LSTM預(yù)測(cè)模型,實(shí)現(xiàn)對(duì)流域污染物重金屬含量的精準(zhǔn)預(yù)測(cè)。
大夏河是黃河上游的重要支流,其發(fā)源于青海省的同仁縣,流經(jīng)夏河縣、臨夏縣等地,最后注入劉家峽水庫(kù),整個(gè)流域總面積達(dá)到7 154 km2。大夏河整個(gè)流域包含了30個(gè)鄉(xiāng)鎮(zhèn)268個(gè)村,流域內(nèi)生活的人口達(dá)到59.96萬(wàn)人。多年平均流量在大夏河縣城段為9.44 m3/s,多年平均徑流量達(dá)到3.12億m3。在降雨方面,流域內(nèi)年平均降水量為44.4 cm,且分布不均勻在時(shí)間和空間范圍上,降水主要集中在7~9月,大部分以連續(xù)降雨或高強(qiáng)度暴雨的形式出現(xiàn)。大夏河流域獨(dú)特的地形和分布不均的降水在一定程度上影響著土壤中重金屬的含量。
1.2.1 多元變量的選取
在各類流域的復(fù)雜生態(tài)系統(tǒng)預(yù)測(cè)研究中,其影響因素具有多樣性。Lorenz于1936年在《確定性的非周期流》[7]中指出:在三階非線性自治系統(tǒng)中可能會(huì)出現(xiàn)混亂解,并因此提出了著名的Lorenz混沌方程,由此產(chǎn)生三變量的混沌時(shí)間序列。因此,本文結(jié)合重金屬含量的相空間經(jīng)過(guò)反復(fù)實(shí)驗(yàn),最終選擇三個(gè)影響因素共同構(gòu)建多元混沌的相空間。對(duì)于影響因素的選擇本文使用較廣泛的相關(guān)系數(shù)法,通過(guò)Spss軟件計(jì)算各影響因素與重金屬含量的相關(guān)性,最終選擇相關(guān)性最高的三個(gè)影響因素溫度、日徑流和重金屬含量。
本文以大夏河2016年11月到2017年12月的數(shù)據(jù)計(jì)算溫度、pH值、降水、日徑流,四種數(shù)據(jù)每星期求一平均值與周數(shù)據(jù)的重金屬含量的相關(guān)性,結(jié)果如表1所示。
表1 各時(shí)間序列與當(dāng)?shù)仄骄亟饘俸康南嚓P(guān)性
如表1所示可以明顯看出溫度、降水以及日徑流數(shù)據(jù)與平均重金屬含量相關(guān)性較高。后續(xù)研究中表明降水?dāng)?shù)據(jù)不具備混沌特性,因此,選擇溫度、日徑流和重金屬含量共同重構(gòu)相空間,其中重金屬含量和溫度、日徑流分別作為預(yù)測(cè)模型的因子。
1.2.2 延遲時(shí)間和嵌入維數(shù)的確定
混沌時(shí)間序列分析的基礎(chǔ)是重構(gòu)相空間,Takens F等人[8]使用延遲坐標(biāo)的方法對(duì)混沌時(shí)間序列進(jìn)行相空間重構(gòu)。
X={Xi|Xi=[xi,xi+τ,…,xi+(m-1)τ]T,i=1,2,…,M}
(1)
式中m為嵌入維數(shù),τ為時(shí)間延遲,M=N-(m-1)τ為相空間中實(shí)際點(diǎn)數(shù)。Kugiumtzis D提出的嵌入窗法,明確指出兩者相關(guān)且依靠嵌入窗τw=(m-1)τ相互關(guān)聯(lián)[9]。本文采用陸振波改進(jìn)的C—C方法[10]在MATLAB中實(shí)現(xiàn)。得出各時(shí)間序列的時(shí)間延遲τ以及嵌入窗τw,根據(jù)公式τw=(m-1)τ計(jì)算三個(gè)時(shí)間序列相關(guān)參數(shù),結(jié)果如表2所示。
表2 重構(gòu)相空間中各影響因素的參數(shù)
1.2.3 最大Lyapunov指數(shù)
本文通過(guò)最大Lyapunov指數(shù)驗(yàn)證該序列混沌特性,采用小數(shù)據(jù)量法計(jì)算該時(shí)間序列的Lyapunov指數(shù)進(jìn)行混沌判別,計(jì)算原理如下:
首先尋找相空間中Yj點(diǎn)的近鄰點(diǎn)并對(duì)其進(jìn)行分離,距離計(jì)算方法如下
(2)
其次,計(jì)算相空間中每個(gè)點(diǎn)Yj對(duì)應(yīng)i個(gè)里斯按時(shí)間步后的距離dj(i)為
dj(i)=|Yj+1-Yi+i|
(3)
對(duì)上述兩步每個(gè)i求得的j的lndj(i)平均y(i),通過(guò)最小二乘法計(jì)算回歸直線的斜率,其斜率為此時(shí)間序列的最大Lyapunov指數(shù)。當(dāng)最大Lyapunov指數(shù)λ大于0時(shí),可以對(duì)該時(shí)間序列的混沌特性進(jìn)行判定,反之亦然。
經(jīng)過(guò)MATLAB計(jì)算可得,重金屬含量、溫度和日徑流的時(shí)間序列的最大Lyapunov指數(shù)分別為0.060 5,0.021 4和0.054 9,結(jié)果均大于0,說(shuō)明以上三個(gè)分量均具有混沌特性。降水時(shí)間序列雖然與重金屬含量的相關(guān)性較高,但由于大夏河流域降水概率較小甚至為0,所以該時(shí)間序列不具有混沌特性,在此予以剔除,不參與相空間的重構(gòu)。
1.2.4 多元混沌相空間
在此簡(jiǎn)單介紹多元相空間建立的基本步驟:設(shè)有Q個(gè)變量時(shí)間序列X1,X2,…,XQ且Xq={Xq,i,Xq,2,…,Xq,n},q=1,2,…,Q,N為觀測(cè)分量時(shí)間序列長(zhǎng)度,對(duì)應(yīng)變量的時(shí)間延遲為τ1,τ2,…τQ,嵌入維數(shù)為m1,m2,…,mQ,此時(shí)多元相空間產(chǎn)生的相點(diǎn)數(shù)L=N-max[(mq-1)τq]。相空間的P時(shí)刻的狀態(tài)坐標(biāo)表示為
Vp=[x1,p,x1,p-τ1,…,x1,p-(m1-1)τ1;x2,p,x2,p-τ2,…,
x2,p-(m2-1)τ2;…;xQ,p,xQ,p-τQ,…,xQ,p-(mQ-1)τQ]
(4)
多變量時(shí)間序列重構(gòu)的相空間矩陣為
(5)
P=N0,N0+1,…,N
(6)
1.2.5 重建關(guān)于重金屬含量的多元混沌相空間
根據(jù)2016年到2017年兩年共117個(gè)數(shù)據(jù)選擇前100個(gè)數(shù)據(jù)建立相空間,后17個(gè)數(shù)據(jù)作為預(yù)測(cè)模型的測(cè)試集。由于其中重金屬含量、溫度和日徑流數(shù)據(jù)間數(shù)值差距較大,因此,對(duì)其進(jìn)行歸一化從而降低最終預(yù)測(cè)誤差。根據(jù)表2給出的m和τ兩個(gè)參數(shù)及式(5)進(jìn)行相空間重構(gòu),得到一個(gè)8×88的流域污染系統(tǒng)相空間矩陣式如下
(7)
N0=13,相空間總維數(shù)M=2+4+2=8。
矩陣中x1,x2,x3分別對(duì)應(yīng)的是重金屬含量、溫度和日徑流指標(biāo)原始時(shí)間序列,PS-LSTM神經(jīng)網(wǎng)絡(luò)的輸入向量和輸出向量為相空間矩陣中的每一行,即將每一個(gè)狀態(tài)點(diǎn)作為上一神經(jīng)網(wǎng)絡(luò)訓(xùn)練集樣本輸出和下一神經(jīng)網(wǎng)絡(luò)訓(xùn)練集樣本的輸入,則訓(xùn)練集樣本輸入和輸出均為8×88。
LSTM神經(jīng)網(wǎng)絡(luò)主要通過(guò)記憶單元來(lái)實(shí)現(xiàn)節(jié)點(diǎn)的不斷更新,從而有效解決了傳統(tǒng)徑向基神經(jīng)網(wǎng)絡(luò)中對(duì)于遠(yuǎn)距離數(shù)據(jù)的依賴特性,LSTM循環(huán)神經(jīng)網(wǎng)絡(luò)作為改進(jìn)的RNN結(jié)構(gòu)[11],在RNN模型的基礎(chǔ)上引入了“記憶細(xì)胞”的概念,2000年,Gers和Schmidhuber提出窺視孔連接,彌補(bǔ)忘記門的一個(gè)缺點(diǎn):目前已有的“記憶細(xì)胞”不能夠影響到“輸入門”和“輸出門”在下一時(shí)刻的輸出,對(duì)上一時(shí)間的序列處理會(huì)使整個(gè)“記憶細(xì)胞”丟失了部分信息,其基本原理是首先將上一個(gè)神經(jīng)元輸出的Ct與此神經(jīng)元的輸入共同輸入到“輸入門”和“忘記門”,然后將其輸入到“記憶細(xì)胞”當(dāng)中;由“記憶細(xì)胞”輸出的數(shù)據(jù)在輸入到“輸出門”的同時(shí)也輸入到下一時(shí)刻的“輸入門”和“忘記門”;最后,“忘記門”輸出的數(shù)據(jù)與激活后的數(shù)據(jù)一起作為整個(gè)“記憶細(xì)胞”的輸出。增加窺視孔后“記憶細(xì)胞”改變的運(yùn)行公式如下
it=σ(Wixxt+Wihht-1+WicCt-1)
(8)
ft=σ(Wfxxt+Wfhht-1+WfcCt-1)
(9)
ot=σ(Woxxt+Wohht-1+WocCt-1)
(10)
式中x為輸入向量,h為輸出向量,Ct為t時(shí)刻在記憶細(xì)胞中的更新狀態(tài),it,ft,ot,Ct,ht分別為輸入門、遺忘門、輸出門、記憶細(xì)胞以及在t時(shí)刻的輸出;xt則為t時(shí)刻的輸入,ht-1,Ct-1分別為隱含層及記憶細(xì)胞在t-1時(shí)刻的輸出。
在PS-LSTM時(shí)間序列預(yù)測(cè)模型中使用多門協(xié)作的方式避免了梯度彌散并且能夠最大程度上模擬輸入因子的共同預(yù)測(cè)。PS-LSTM時(shí)間序列預(yù)測(cè)模型通過(guò)Python編碼實(shí)現(xiàn),以大夏河流域重金屬的多元混沌相空間的狀態(tài)坐標(biāo)為多變量輸入,本次模型將按照測(cè)試集85 %的比例選取,所以前75個(gè)數(shù)據(jù)作為訓(xùn)練樣本,重金屬含量中的最后13個(gè)數(shù)據(jù)作為測(cè)試樣本進(jìn)行實(shí)驗(yàn)。使用單層隱含層的PS-LSTM模型,總的訓(xùn)練迭代次數(shù)為300。本文在訓(xùn)練模型時(shí)將其一次訓(xùn)練所選取的樣本數(shù)設(shè)置為45,時(shí)間步長(zhǎng)為3,隱含層節(jié)點(diǎn)為100,輸出節(jié)點(diǎn)為1,利用平均絕對(duì)誤差(mean absolute error,MAE)做損失函數(shù),Adam的隨機(jī)梯度下降做優(yōu)化,根據(jù)均方根誤差(root mean square error,RMSE)進(jìn)行模型評(píng)估。當(dāng)0迭代到50次時(shí),測(cè)試損失與訓(xùn)練損失下降較快,迭代到第100次時(shí),模型損失達(dá)到最低并趨于平緩,如圖1所示可知測(cè)試損失在絕大多數(shù)的情況下高于訓(xùn)練損失,證明模型擬合較為成功,適合作為重金屬含量的預(yù)測(cè)模型。
圖1 不同訓(xùn)練時(shí)期訓(xùn)練集和測(cè)試集損失
Volterra級(jí)數(shù)在收斂區(qū)域的判定上有著嚴(yán)格的要求,增加了它的使用難度[12]。本文以重金屬含量為輸入數(shù)據(jù)使用陸振波的MATLAB混沌工具箱中的Volterra級(jí)數(shù)一步預(yù)測(cè)模型進(jìn)行預(yù)測(cè),其中只使用2016年11月至2017年12月的117個(gè)重金屬含量數(shù)據(jù),并沒有建立關(guān)于重金屬的相空間。從圖2所示的預(yù)測(cè)結(jié)果中可以看出開始預(yù)測(cè)的4個(gè)數(shù)據(jù)較為準(zhǔn)確,但是隨后的13個(gè)預(yù)測(cè)數(shù)據(jù)相比真實(shí)值有較大的差距。
圖2 Volterra級(jí)數(shù)預(yù)測(cè)結(jié)果
徑向基函數(shù)(radial basis function,RBF)是前向型神經(jīng)網(wǎng)絡(luò),主要由輸入層、隱含層、輸出層三部分組成[13]。RBF神經(jīng)網(wǎng)絡(luò)依靠其簡(jiǎn)單的網(wǎng)絡(luò)結(jié)構(gòu)、較快的收斂速度和能夠逼近任意非線性函數(shù)的特征使其在時(shí)間序列分析和圖形識(shí)別處理等領(lǐng)域獲得較為廣泛的應(yīng)用[14]。
其在關(guān)于有噪聲的重金屬含量時(shí)間序列預(yù)測(cè)上能夠?qū)崿F(xiàn)非線性的局部預(yù)測(cè)。本文使用陸振波的MATLAB混沌工具箱中的RBF神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)模型實(shí)現(xiàn)對(duì)重金屬含量的預(yù)測(cè),預(yù)測(cè)結(jié)果如圖3所示。
圖3 RBF神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)結(jié)果
支持向量回歸(support vector regression,SVR) 在函數(shù)回歸預(yù)測(cè)方面的應(yīng)用依據(jù)支持向量機(jī)(support vector machine,SVM)原理,由SVM在其分類問題上所推廣而來(lái),其是以預(yù)測(cè)的相關(guān)性原則作為基礎(chǔ),確定預(yù)測(cè)目標(biāo)的各特征與預(yù)測(cè)目標(biāo)之間的數(shù)學(xué)模型,利用該模型對(duì)其特征的變化值進(jìn)行目標(biāo)預(yù)測(cè)。
大夏河流域重金屬含量在MATLAB中使用SVR方法進(jìn)行了單因子的預(yù)測(cè),由于重金屬含量屬于小數(shù)據(jù)的預(yù)測(cè),因此,SVR的小數(shù)據(jù)量非線性預(yù)測(cè)能夠具有較強(qiáng)的對(duì)比性。由圖4的預(yù)測(cè)結(jié)果圖所示前5個(gè)重金屬含量的預(yù)測(cè)準(zhǔn)確程度相對(duì)較為理想,但后面的預(yù)測(cè)準(zhǔn)確程度較差。SVR預(yù)測(cè)準(zhǔn)確程度相比其他兩個(gè)對(duì)比實(shí)驗(yàn)具有較好的效果,但經(jīng)過(guò)改進(jìn)的LSTM預(yù)測(cè)模型在預(yù)測(cè)準(zhǔn)確程度上更準(zhǔn)確。
圖4 基于SVR回歸預(yù)測(cè)結(jié)果圖
由表3可得,在平均絕對(duì)誤差和均方根誤差指標(biāo)方面,本文的預(yù)測(cè)模型的評(píng)價(jià)指標(biāo)值RMSE均小于Volterra級(jí)數(shù),RBF神經(jīng)網(wǎng)絡(luò),SVR測(cè)試模型,而決定系數(shù)比其他三種預(yù)測(cè)模型高,說(shuō)明PS-LSTM預(yù)測(cè)模型的擬合程度比其他三種模型優(yōu)度高,因此證明建立關(guān)于重金屬含量的多元混沌相空間是有一定的實(shí)際意義,并且LSTM預(yù)測(cè)模型對(duì)于有噪聲的混沌時(shí)間序列預(yù)測(cè)精度上具有優(yōu)勢(shì)。
表3 各預(yù)測(cè)模型的結(jié)果評(píng)價(jià)
本文利用較為完善的增加窺視孔的LSTM—RNN模型,LSTM對(duì)有噪聲的重金屬含量預(yù)測(cè)相較Volterra級(jí)數(shù)預(yù)測(cè)、RBF神經(jīng)網(wǎng)絡(luò)預(yù)測(cè)以及基于支持向量機(jī)回歸預(yù)測(cè)模型在預(yù)測(cè)精度上具有一定優(yōu)勢(shì)。因此本文的預(yù)測(cè)模型具有一定的意義。采用相空間與門控循環(huán)單元(GRU)進(jìn)一步的改進(jìn)相結(jié)合的方法可以實(shí)現(xiàn)更加精準(zhǔn)預(yù)測(cè),這將是下一步的研究重點(diǎn)。