基于BP神經(jīng)網(wǎng)絡(luò)和SVR的Fund?o尾礦壩排水數(shù)據(jù)預測對比研究*

戴健非，楊 鵬,2，王昕宇

(1.北京聯(lián)合大學 北京市信息服務(wù)工程重點實驗室，北京 100101； 2.北京科技大學 土木與資源工程學院，北京 100083)

0 引言

尾礦庫是維持礦山企業(yè)正常運營的必要設(shè)施，由1個或多個尾礦壩堆筑攔截谷口或圍地所構(gòu)成的礦山生產(chǎn)設(shè)施，是具有高勢能的人造危險源。其穩(wěn)定性與安全性對礦山企業(yè)的正常生產(chǎn)和區(qū)域安全至關(guān)重要[1-3]。因此，研究尾礦壩潰壩失效的早期預警是很有必要的。在Fund?o尾礦壩倒塌1 年后，超過4 300萬m3的鐵礦石尾礦仍然造成環(huán)境破壞，污染了從Doce河流域到大西洋的668 km水道。污染物的數(shù)量和受影響的生態(tài)系統(tǒng)的范圍占據(jù)了前所未有的比例，涉及巴西大西洋森林河口，沿海和海洋環(huán)境[4]。對潰壩原因調(diào)查后發(fā)現(xiàn)[5]，尾礦壩由于一系列“施工缺陷”導致了砂土液化。其中，基礎(chǔ)排水溝的施工缺陷最為嚴重，排水系統(tǒng)有缺陷，導致庫內(nèi)尾礦飽和度過高，小型地震誘發(fā)潰壩。由此可見，排水能力是否正常在尾礦壩安全領(lǐng)域起著重要的指示作用。本文根據(jù)專家小組的數(shù)據(jù)支持預測排水數(shù)據(jù)flow即每秒排水量(L/s)。國內(nèi)外學者在尾礦庫安全預警方面做了大量研究。在尾礦壩數(shù)據(jù)預測方面， Dong等[6]通過考慮潛水線預測方程，極限平衡狀態(tài)參數(shù)，水庫水位和降雨量建立數(shù)值模擬模型，通過遠程實時預警獲得尾礦壩穩(wěn)定或危險的預警信號；王英博等[7]采用修正型果蠅算法優(yōu)化的BP神經(jīng)網(wǎng)絡(luò)對灘頂高度、庫水位、浸潤線、干灘長度和安全超高5個指標擬合,顯示出較高的預測精度；李娟等[8]使用SVM預測尾礦壩浸潤線,實現(xiàn)了小樣本情況下的高精度預測； Tongle等[9]基于遺傳算法和BP神經(jīng)網(wǎng)絡(luò)預測尾礦浸潤線的未來狀態(tài)。在尾礦壩安全理論方面，何學秋等[10]經(jīng)過多次實驗得出尾礦壩變形需經(jīng)歷衰減、穩(wěn)定、加速3個階段的結(jié)論，基于流變-突變理論,分別制定預警準則；王昆等[11]綜合分析了國內(nèi)外尾礦庫潰壩災害防控現(xiàn)狀及發(fā)展情況；張興凱等[12]采用物理模型試驗方法建立尾礦庫漫頂潰壩演化模型，證明了尾礦庫洪水漫頂潰壩位移與壩體飽和程度有關(guān)，壩體浸潤線越高,尾礦庫潰壩時滑動位移越大，潰口破壞程度取決于溢流對壩體的沖刷侵蝕作用。在傳感器穩(wěn)定采集數(shù)據(jù)方面，陳凱等[13]研究了礦山監(jiān)測系統(tǒng)的防護、供電能力、網(wǎng)絡(luò)通信等技術(shù),來確保極端天氣環(huán)境下的數(shù)據(jù)獲取。

為實現(xiàn)尾礦庫排水流量數(shù)據(jù)的實時預測，本文建立基于BP神經(jīng)網(wǎng)絡(luò)(Back Propagation Neural Network)和支持向量回歸(Support Vector Regression， SVR)的巴西Fund?o尾礦壩流量數(shù)據(jù)的預測回歸模型。分析庫水位、降雨量和干灘長度對排水數(shù)據(jù)的影響。本文所述模型可為Fund?o壩的潰壩預警提供決策支持，對防止?jié)问鹿屎捅Ｗo尾礦安全具有參考價值。

1 數(shù)據(jù)來源與傳感器

本文所用數(shù)據(jù)均來自于文獻[5]。

巴西Fund?o尾礦壩的流量監(jiān)測是采用位于海拔826 m處的流量傳感器采集得到。在826 m處安裝有27個排洪隧洞進行排水，材料為混凝土材質(zhì)。2011年的ITRB(Independent Tailings Review Board)報告指出，排水系統(tǒng)朝向右側(cè)和左側(cè)基臺，但是，只能從右側(cè)基臺附近的水槽進行測量[5]。

排水數(shù)據(jù)監(jiān)測有每小時排水立方米和每秒排水升數(shù)2種計量模式。本文僅采用每秒排水升數(shù)作為模型的輸出特征，并記為flow。庫水位和降雨量來自于小組審查的Samarco現(xiàn)場監(jiān)測數(shù)據(jù)。安裝的監(jiān)測儀器類型是振弦式壓力計和水位指示器(立管)。主要數(shù)據(jù)來源是Samarco的月度儀表報告和Samarco電子數(shù)據(jù)表格。干灘長度數(shù)據(jù)每月更新1次。數(shù)據(jù)來源共有4種：月度尾礦排放報告、使用無人機地形測量估算干灘長度及上升率電子表格、每月巖土監(jiān)測報告和地形測量數(shù)據(jù)。選取的潰壩災害數(shù)據(jù)包含3個特征屬性：庫水位(m)、降雨量(mm)、干灘長度(m)。本文提取了2011年7月至2015年11月的巴西Fund?o壩的上述3種潰壩災害特征數(shù)據(jù)和排水數(shù)據(jù)共62組數(shù)據(jù)作為數(shù)據(jù)源。其中，2011年7月至 2015年3月的數(shù)據(jù)作為訓練集，部分歷史數(shù)據(jù)如表1所示。

表1 Fund?o壩部分歷史數(shù)據(jù)Table 1 Partial historical datas of the Fund?o

2 基于BP神經(jīng)網(wǎng)絡(luò)預測模型

人工神經(jīng)網(wǎng)絡(luò)技術(shù)是對大腦功能的模擬、抽象和簡化，具有學習、聯(lián)想和記憶功能，以及高度的非線性預測能力[14]。本文采用基于Levenberg-Marquardt算法的BP雙隱層神經(jīng)網(wǎng)絡(luò)對巴西Fund?o壩的排水數(shù)據(jù)進行預測建模。該雙隱層網(wǎng)絡(luò)包含1個輸入層，1個輸出層，2個隱藏層。其中輸入層有3個神經(jīng)元節(jié)點，輸出層包含1個神經(jīng)元節(jié)點。為了全面分析巴西Fund?o壩的排水流量規(guī)律，本文提取了2011年7月至2015年3月的庫水位、降雨量和干灘長度數(shù)據(jù)做訓練集的輸入，而對應的排水流量則為輸出特征。在BP神經(jīng)網(wǎng)絡(luò)模型中，回歸系數(shù)R用于評估模型的有效性，R越接近1，模型的預測能力越精準。

2.1 Levenberg-Marquardt算法

Levenberg-Marquardt算法結(jié)合了高斯牛頓法和梯度下降法的優(yōu)點，是目前使用廣泛的優(yōu)化算法。設(shè)x(i)為當前狀態(tài)的權(quán)重向量，x(i+1)為BP神經(jīng)網(wǎng)絡(luò)中更新后的權(quán)重向量，Δx為權(quán)重更新值，有：

x(i+1)=x(i)+Δx

(1)

Δx=-[2E(x)]-1E(x)

(2)

(3)

式中：E(x)為損失函數(shù)，選用均方誤差；ei為排水數(shù)據(jù)真實值與預測值之差；2E(x)為損失函數(shù)的Hessian矩陣，即JT(x)J(x)；E(x)為其梯度；N為樣本個數(shù)。這時對于BP神經(jīng)網(wǎng)絡(luò)：

x(i+1)=x(i)-[JT(xi)J(xi)+μiI]-1JT(xi)ei

(4)

即

Δx=-[JT(xi)J(xi)+μiI]-1JT(xi)ei

(5)

式中：μ為大于0的比例系數(shù)；I為單位矩陣；J(x)為損失函數(shù)E(x)的雅可比矩陣，即誤差向量與更新權(quán)重向量的一階微分矩陣。

此時的Levenberg-Marquardt算法由μ來控制，μ越大，越接近梯度下降法，當μ=0時，Levenberg-Marquardt變?yōu)楦咚古ｎD法。

2.2 神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)和模型訓練

多層前向神經(jīng)網(wǎng)絡(luò)具有理論上可逼近任意非線性連續(xù)映射的能力，因而非常適合于非線性系統(tǒng)的建模及控制[15]。根據(jù)專家小組提供的數(shù)據(jù)和Fund?o大壩排洪隧洞遭到破壞的實際情況，影響尾礦排水能力的應變因素主要包括庫水位(Res)、降雨量(Rai)和干灘長度(Bre)3個因素[13]。將這3個因素作為BP網(wǎng)絡(luò)的輸入層神經(jīng)元，排水數(shù)據(jù)flow作為輸出層的神經(jīng)元。而隱含層神經(jīng)元個數(shù)則根據(jù)式(6)初步確定[16]：

(6)

式中：m為隱含層節(jié)點數(shù)；n為輸入層節(jié)點數(shù)；l為輸出層節(jié)點數(shù)；α為1～10之間的常數(shù)。網(wǎng)絡(luò)的初步結(jié)構(gòu)設(shè)定為3-3-3-1。

網(wǎng)絡(luò)隱含層的傳遞函數(shù)選用雙曲正切S型傳遞函數(shù)tansig，輸出層則使用線性傳遞函數(shù) purelin。模型的訓練選用trainlm，即Levenberg-Marquardt優(yōu)化算法。憑此設(shè)置，使用MATLAB工具建立BP神經(jīng)網(wǎng)絡(luò)模型，經(jīng)

過多次調(diào)整確定網(wǎng)絡(luò)結(jié)構(gòu)為3-4-4-1，網(wǎng)絡(luò)結(jié)構(gòu)如圖1。此時的回歸系數(shù)是0.985 6，說明網(wǎng)絡(luò)的訓練效果良好，訓練效果如圖2。其中，最大迭代次數(shù)為5 000次，網(wǎng)絡(luò)的學習率為0.01，目標誤差為0.65×10-3。

圖1 BP雙隱層神經(jīng)網(wǎng)絡(luò)結(jié)構(gòu)Fig.1 Structure of BP double hidden layer neural network

圖2 BP神經(jīng)網(wǎng)絡(luò)訓練效果Fig.2 BP neural network training effect

2.3 模型預測效果驗證

為了驗證模型的效果好壞，取10組訓練樣本之外的數(shù)據(jù)即2015年3月-11月的數(shù)據(jù)輸入模型中進行預測，其中，只保留庫水位、降雨量和干灘長度。預測的結(jié)果與實際真值對比如圖3、表2所示。

圖3 BP神經(jīng)網(wǎng)絡(luò)預測結(jié)果Fig.3 BP neural network prediction results

序號日期庫水位/m降雨量/mm干灘長度/m樣本真值/(L·s-1)樣本預測值/(L·s-1)相對誤差/%12015-3-1884.988.66209.1370.5169.131.9622015-4-1885.4033.82194.7170.5169.151.9332015-5-1886.9011.67197.1164.6066.162.4142015-6-1887.542.82290.8656.2857.672.4752015-6-5887.9713.17304.0861.7559.653.4062015-6-10888.180.00295.8759.7857.184.3572015-7-1889.462.86168.2665.4764.681.2182015-9-1890.952.90290.8657.8156.622.0792015-10-1891.590.00216.3570.9569.192.48102015-11-1892.8716.26240.3869.8569.081.10

由圖3可知，神經(jīng)網(wǎng)絡(luò)的預測結(jié)果與實際結(jié)果之間相差不大，變換趨勢基本相同，表示該模型訓練成功。由表2可知，排水數(shù)據(jù)flow(L/s)對應的10組驗證數(shù)據(jù)的相對誤差范圍在1.10%至4.35%之間，該BP神經(jīng)網(wǎng)絡(luò)的運行時間為8 s。可見，BP神經(jīng)網(wǎng)絡(luò)預測模型的誤差范圍相對較小，準確度較高。

3 基于SVR算法的回歸預測模型

3.1 SVR回歸預測模型的建立

SVR回歸即找到1個高維的回歸平面，讓數(shù)據(jù)集內(nèi)的所有樣本到該平面距離最近。SVR也是小樣本數(shù)據(jù)集下性能可觀的預測模型。

已知數(shù)據(jù)樣本是三維向量，訓練集可表示為{(xi,yi),…,(xk,yk)}，其中xi∈R3是輸入變量即庫水位、降雨量和干灘長度，yi∈R是排水數(shù)據(jù)實際值，i=1,2,…,k是訓練樣本的個數(shù)。此時構(gòu)建回歸超平面：

f(x)=ω·x+b

(7)

式中：f(x)為排水數(shù)據(jù)的預測值；ω為權(quán)值向量，ω·x做內(nèi)積運算；b為偏移值。

(8)

式中：ε為損失界限，控制了模型對數(shù)據(jù)的誤差寬度。引入拉格朗日(Lagrange)乘子a,a*將凸二次優(yōu)化問題轉(zhuǎn)變?yōu)閷ε紗栴}：

(9)

求解參數(shù)，符合決策變量的數(shù)據(jù)成為支持向量回歸，此時的目標函數(shù)為：

(10)

因為輸入特征與輸出特征屬于非線性關(guān)系，本模型采用徑向基核函數(shù)，函數(shù)如式(11)：

(11)

式中：y為核函數(shù)中心；σ為核函數(shù)的寬度參數(shù)， 控制了函數(shù)的徑向作用范圍；x為樣本點。高斯核函數(shù)將數(shù)據(jù)從其特征空間映射到無窮維空間，解決了非線性不可分的問題，擬合出1條高維的曲線達到預測效果[18]。最后模型的目標函數(shù)為：

(12)

模擬過程如圖4所示。

圖4 SVR預測模型流程Fig.4 SVR prediction model flow chart

1)回歸預測模型是基于特征屬性搭建的，將2011年7月至2015年3月Fund?o尾礦大壩的庫水位(Res)、降雨量(Rai)和干灘長度(Bre)作為特征輸入，排水數(shù)據(jù)作為輸出。

2)采用和BP神經(jīng)網(wǎng)絡(luò)預測排水數(shù)據(jù)一樣的測試集。調(diào)用sklearn中的svm庫，核函數(shù)設(shè)置為rbf徑向基函數(shù)，擬合SVR回歸曲線，采用試湊法經(jīng)過多次試驗找到最優(yōu)參數(shù)確定為gamma=26,C=1。

3)完成模型回歸，得到預測結(jié)果，與真值對比計算誤差。

3.2 結(jié)果驗證

借助Anaconda+Pycharm集成開發(fā)環(huán)境，調(diào)用sklearn庫中的svm模塊，對數(shù)據(jù)集進行支持向量回歸預測。目前，常用的最優(yōu)參數(shù)選取方法有K折交叉驗證法、留一法和試湊法。留一法和交叉驗證法不能固定訓練樣本，無法與BP神經(jīng)網(wǎng)絡(luò)方法做結(jié)果對比，因此選用試湊法。試湊法的基本原理為：首先為函數(shù)固有的參數(shù)賦初始值，然后開始實驗測試，根據(jù)測試精度重復調(diào)整參數(shù)值，直至得到滿意的測試精度為止[19]。選取默認參數(shù)gamma=1，C=1，固定gamma選用不同C值確定最低均方根誤差(Root Mean Squard Error, RMSE)，此時gamma=1，C=5，RMSE=3.61。再固定C值，選用不同gamma確定最低RMSE值，此時gamma=25，C=5，RMSE=2.76。重復上述步驟直至RMSE值不再降低。經(jīng)過多次調(diào)試，最優(yōu)參數(shù)設(shè)置為gamma=26,C=1，此時RMSE=2.74。調(diào)試過程見圖5。表3為SVR模型預測結(jié)果，可知，預測排水數(shù)據(jù)相對誤差范圍在1.64%～9.21%之間，最大誤差小于10%，運行時間為0.08 s，所以認為該模型的預測效果良好。

圖5 SVR模型參數(shù)調(diào)試Fig.5 SVR model parameter modulation

序號日期樣本真值/(L·s-1)樣本預測值/(L·s-1)相對誤差/%12015-3-170.5171.671.6522015-4-170.5171.681.6532015-5-164.6067.405.0142015-6-156.2859.535.7852015-6-561.7562.921.8862015-6-1059.7865.299.2172015-7-165.4761.496.0882015-9-157.8158.972.0192015-10-170.9572.121.64102015-11-169.8571.011.67

綜合對比BP神經(jīng)網(wǎng)絡(luò)與SVR回歸預測模型，二者的預測機制不同，預測效果也不一樣。前者通過反向傳播誤差更新權(quán)值來預測數(shù)據(jù)，后者通過在高維空間擬合超平面使樣本點距離平面最近預測數(shù)據(jù)。2個模型對比結(jié)果如表4所示，SVR回歸預測模型的平均相對誤差為3.66%，而BP神經(jīng)網(wǎng)絡(luò)模型的平均相對誤差則只有2.34%。從時間上看，SVR預測模型要更快。

表4 模型對比結(jié)果Table 4 Results of model comparison

4 結(jié)論

1)以庫水位、降雨量和干灘長度為輸入，運用BP神經(jīng)網(wǎng)絡(luò)預測模型和SVR回歸預測模型2種方法預測排水數(shù)據(jù)，2種模型的最大誤差分別為4.35%和9.21%，均小于10%，說明模型的建立和預測結(jié)果是可行的。

2)BP神經(jīng)網(wǎng)絡(luò)和SVR預測模型克服了排水數(shù)據(jù)復雜且非線性的特點。前者的預測精度高，平均誤差為2.34%，后者的運行速度快，運行時間為0.08 s。其中，核函數(shù)為徑向基函數(shù)，參數(shù)gamma=26，C=1時，SVR模型預測結(jié)果最優(yōu)。

3)借助BP神經(jīng)網(wǎng)絡(luò)預測模型和SVR回歸預測模型，礦山安全管理部門可根據(jù)模型預測結(jié)果預測排水數(shù)據(jù)，對于尾礦庫安全管理具有一定參考價值。