基于地質大數據的泥石流災害易發性評價

張永宏，葛濤濤，田 偉，夏廣浩，何 靜

(1.南京信息工程大學 信息與控制學院，南京210044； 2.南京信息工程大學 計算機與軟件學院，南京 210044)(*通信作者電子郵箱tw@nuist.edu.cn)

0 引言

隨著大數據時代的到來，大數據所創造的價值滲透在當今生活的方方面面，已成為高校和企業的重點研究對象。傳統方法難以處理量大、異構、更新快速的數據[1]，而神經網絡特征提取能力強，易于整合多源信息，對動態變化的數據適應性高的特點，足夠使其成為大數據價值轉化的紐帶；同時，大數據也為神經網絡訓練供應了足夠多的訓練樣本。伴隨著計算機軟硬件技術的發展，神經網絡應用到了各個領域，如語音大數據分析[2]、醫學大數據分析[3]、農業大數據分析[4]以及電力大數據分析[5]等，大數據結合神經網絡已成為當下推動社會與科學發展的重要動力。

長期地質調查工作積累了大量的地質原始數據，通過遙感、鉆探、物探等手段獲得多源、多模態、多時態的地質數據具有大數據的典型特性，即地質大數據。泥石流是發生在山區的一種突發性地質災害，其易發性評價以地質環境為基礎，根據地質災害現狀，靜態考察泥石流災害在一定區域內發生的可能性大小，它的強弱同時受多種外在因素和內在因素的共同影響[6]。泥石流易發性評價的實施可以為當地公共設施建立、居民人生安全以及防災減災提供巨大幫助。

目前比較常用的泥石流易發性的評價方法有層次分析法[7]、模糊數學[8-9]、邏輯回歸[6]等，也有學者將已有方法結合使用，取長補短[10]。這些算法在取得不錯成績的同時也存在一些不足之處，如:層次分析法指標體系的提出需要有專家系統的支持，且泥石流指標維度較多，標度工作量太大；模糊綜合評價計算復雜，確定泥石流指標權重向量時過于主觀，造成泥石流易發性評價精度不高；邏輯回歸屬于廣義線性回歸模型，泥石流易發性評價指標中的多重共線性問題對其影響較大。

從上述傳統泥石流易發性評價方法來看，它們存在一些共同的局限性：評價主觀，難以處理多維度的大數據信息。神經網絡具有良好的非線性映射能力，能夠在海量地質數據中揭示數據內部聯系，對未來作出精準的預測，在泥石流易發性評價應用中成果顯著[11]; 然而神經網絡存在學習速度慢，難以處理非平衡數據等問題，使其作用得不到充分發揮，并且面對地質大數據，無法找出輸入指標對易發性的影響程度。針對傳統方法的不足，本文基于神經網絡，結合Borderline-SMOTE(Synthetic Minority Oversampling TEchnique)算法，遺傳算法(Genetic Algorithm, GA)和平均影響值(Mean Impact Value, MIV)算法對影響泥石流發生的地質大數據進行分析與挖掘，并建立泥石流易發性評價模型。實驗證明，在研究地處“一帶一路”重要節點之一的雅魯藏布江中上游流域泥石流孕災環境的基礎上，使用TensorFlow框架搭建的神經網絡結合優化算法可以對研究區域內山地泥石流災害易發性進行準確的評價。

1 泥石流易發性評價模型框架

為了應對區域泥石流易發性評價所面臨的挑戰，本文構建了一種能夠客觀準確地評價區域泥石流易發性的模型框架，如圖1所示。本框架一共包含5層，自下而上分別為數據層、學習框架層、易發性評價方法層、易發性評價技術層以及應用層。該泥石流易發性評價框架中，各個層次核心功能描述如下：

數據層 數據層主要提供遙感圖像，地形圖，地質圖，災害歷史事件等基礎資料。對地質大數據數據庫的管理包括定義、檢索、修改、增減等操作，并建立數據字典。由于泥石流溝的基礎數據以及各種孕災環境因子處于動態變化之中，所以該層應具有開放性，有利于隨時更新和維護數據庫中的存儲數據。

學習框架層 TensorFlow學習框架能夠在多種運行環境下使用，并且在前端支持多種主流語言調用其API(Application Programming Interface)。該框架封裝了當下流行的眾多神經網絡模型，可以根據需要構建和實現所需要的網絡模型。本文在Windows系統下調用TensorFlow框架，通過import命令導入待使用的數據庫，并使用Python編程語言搭建能夠準確評價泥石流易發程度的學習架構。

易發性評價方法層 泥石流成災條件復雜，影響其易發程度的因素隨機且模糊。BP(Back Propagation)神經網絡無須建立確定的模型，只需根據一定量的樣本來學習研究對象的內部結構，便可完成輸入集合與輸出結果之間的非線性映射。網絡采取隨機初始化的方式獲取權值和閾值，導致分類效率低下，易陷入局部最優, 因此，采用遺傳算法搜索全局最優權值和閾值，可縮減訓練時間，提升網絡的泛化性能。

易發性評價技術層 在易發性評價過程中，分類器能夠對類別分布均衡的數據作較為精準的分類; 然而，實際生活中，不同易發等級的數據極不平衡, 且雖然網絡能夠逼近任意非線性關系，但是難以表達輸入指標與模型的相關性, 因此，需要設計易發性評價方法優化技術，在解決類別間數據不平衡問題對分類器影響的同時篩選輸入指標。分別采用重采樣技術和平均影響值算法，達成上述優化目標。

應用層 泥石流易發性評價實際應用中，將已知的歷史災害數據集按訓練集、驗證集、測試集進行合理的分配。設網絡的訓練樣本數據集為{(x1,d1), (x2,d2),…,(xm,dm)}，xi為泥石流i的特征向量，di為對應的泥石流易發性等級標簽，di∈{1,2,…,r}，r為標簽的類別數，m為訓練的樣本個數。基于泥石流易發性評價方法層和技術層，實現類別不平衡數據下，根據輸入孕災環境，快速準確地評價區域泥石流易發性，并獲取影響當地泥石流發生的主要影響因子。

圖1 泥石流易發性評價模型框架

2 泥石流易發性評價方法

本文對于泥石流易發性的評價包含兩個基本過程：基于神經網絡的泥石流易發性模型搭建及其技術優化，所以將泥石流易發性評價問題劃分為以下3個步驟：

1)從采集到的泥石流災害系統的眾多影響因子中提取能夠準確反映區域特征的泥石流易發性評估指標，對其作歸一化和過采樣Borderline-SMOTE的處理；

2)構建評價泥石流易發性神經網絡模型，將過采樣后的訓練數據輸入到網絡進行訓練，調整并確定模型參數，并對已經確定結構的神經網絡采用遺傳算法進行優化，找到全局最優權值和閾值，提升模型的效率和泛化性能；

3)使用MIV算法探討指標因子與上述GA_BP模型的相關性，刪減冗余信息，最后通過網絡Softmax分類器計算模型輸出結果。

技術路線圖如圖2所示。

圖2 泥石流易發性分類技術路線

最后采用基于受試者工作特征曲線(Receiver Operating Characteristic Curve, ROC)分析的AUC(Area Under the Curve)方法評估模型的準確性和可靠性。

2.1 重采樣技術

過采樣技術在解決數據不平衡分類問題中具有廣泛的應用價值[12-13]。在山地災害研究中，往往可以發現，中易發性的地區占絕大多數，而其他易發性的地區是相對較少的。分類過程中，少數類數據稀少是影響神經網絡分類精度的主要原因之一。本文使用一種基于邊界的SMOTE算法(Borderline-SMOTE)，其具體的方法為：

首先使用K近鄰算法計算少數類樣本集T中每一個樣本在訓練集x中的K近鄰樣本集P，并分析其中屬于多數類樣本集T′的比例：若多數類樣本占少數則直接將該樣本放回少數類樣本集；若多數類樣本占多數則將該樣本作為邊界樣本，放至邊界(Borderline)集合F中。計算邊界樣本集中每個樣本Ti在T中的K近鄰，并從中隨機選取一個樣本Tj，使用Ti與Tj結合SMOTE算法，由式(1)生成新的少數類樣本Tnew，合成的新樣本Tnew加入到原始樣本集形成新的樣本集Xinput。

Tnew=Ti+random(0,1)*|Ti-Tj|

(1)

算法描述偽代碼如下：

算法1 過采樣。

InputK,m,T,T′,x；

OutputXinput。

1)

fori=1 to |T| do

2)

計算少數類樣本中Ti在訓練集x中的K近鄰集合P;

3)

if (K/2<|P∩T′|

4)

將Ti加入到邊界集合F中;

5)

計算集合F中元素在少數類樣本的K近鄰;

6)

end if

7)

end for

8)

z=|T′|-|F|;

9)

whilez≠0 do

10)

隨機選擇集合F中的Ti;

11)

diff=Ti與其近鄰Tj之間的距離;

12)

gap=0到1之間的隨機數;

13)

Tnew=Ti+gap*diff;

14)

Xinput=Tnew∪x

15)

z=z-1;

16)

end while

17)

returnXinput

2.2 易發性分類方法

人工神經網絡(Artificial Neural Network, ANN)是一種模仿生物神經系統處理信息的數學模型，BP神經網絡作為ANN中最精華的部分，通常具有3層及以上的網絡結構，包括輸入層、輸出層以及隱藏層， 其功能包括信號的前向傳播和誤差的反向傳播。

由于網絡各層的激活函數和學習規則等存在差異，所以在神經網絡的輸入端應對重采樣后的數據使用歸一化算法，定義如下：

xnew=(xij-xi_min)/(xi max-xi_min)

(2)

樣本x完成歸一化后，樣本集Xinput離開輸入層進入隱藏層，經傳遞函數(如Tanh、Relu等函數，經實驗對比，本文使用Tanh函數)激活后繼續往下一層傳遞，直至輸出層，各層輸出公式為：

(3)

其中:f(·)表示傳遞函數；θ={w,b}表示網絡參數，其中w為權重，b為閾值。網絡的輸出層采用Softmax分類器進行分類識別，其公式定義為：

(4)

式中：p(di=k|xi)表示樣本xi識別為第k個類別的概率；hλ(xi)是r行1列的向量，即r個類別分類器輸出的也是一個r維向量，所有行元素之和為1；λ1,λ2,…,λr∈Rn+1是分類器的參數，且構成參數矩陣λ。對于xi，選擇輸出向量中最大概率對應的k值作為分類結果。

模型的訓練目標是通過尋找一組最優參數θ*={w*,b*}，使預測值與實際值誤差最小，即在訓練過程中，使代價函數J(λ)達到最小，代價函數J(λ)定義為：

(5)

式中：1{·}是一個指示性函數，大括號里的值為真返回1； 否則返回0。經過實驗，利用誤差逆傳播和Adam優化算法可實現式(5)達到最小，結束訓練。

由上可知，BP算法的權重和閾值對學習結果存在很大影響，而網絡采用隨機初始化的方式使其存在波動、難以收斂或訓練速度慢的問題，最終導致學習結果不理想。遺傳算法是一種基于自然選擇和遺傳變異等自然界進化過程的計算模型，具有很強的全局搜索能力和全局優化性能[14]。遺傳算法將前述神經網絡中的初始權值和閾值作為種群中的個體，采用二進制編碼為其編碼，編碼長度由BP神經網絡的結構而定，包括網絡層數及各層節點數。其中隱層神經元個數通常由式(6)確定：

(6)

式中:l為隱層節點數，n為輸入層節點數，r為輸出層節點數，α為整數，取值范圍為1～10。一組個體的集合作為種群(pop)，種群規模應根據具體問題選擇合適的數量，初始種群規模設為S；最大遺傳代數設為Gen。遺傳算法優化神經網絡權值和閾值主要分為兩步: 第一步, 構造適應度函數用于考察個體的優劣情況，本文設定為神經網絡代價函數的倒數。輸出值與期望值之差越小，適應度值越高，個體越有優勢，反之則會在進化過程中被淘汰。 第二步，進行遺傳操作：1)根據計算的適應度，進行個體的選擇。首先通過最優策略選出適應度最高的最佳個體，直接遺傳給下一代種群，其余個體通過概率值進行選擇，公式如下：

(7)

其中fi為個體i的適應度值。2)選擇種群中的一對作為父體，以一定的交叉概率Pc進行交叉配對，產生另外兩個新個體。3)隨機選擇種群中的個體，以特定概率Pm進行變異，即通過改變個體中的某一基因位而產生新個體。4)最后進行循環操作，當遺傳代數小于最大遺傳次數時，轉至第一步，繼續運行。直至指定的進化次數時停止運行，從種群中選擇適應度最高的個體，輸出對應的權值w和閾值b。最后將這組權值和閾值應用于BP神經網絡中，當達到目標準確率(Aim_accuracy)時，結束訓練，輸出分類結果。算法偽代碼描述如下：

算法2 遺傳神經網絡(GA_BP)。

輸入Pm,Pc,Gen,Aim_accuracy,S,Xinput；

輸出 GA_BP網絡結構和參數。

1)

初始化遺傳算法的參數和種群pop;

2)

t=0;

3)

while(t

4)

基于Xinput數據集計算S組w和b的網絡誤差;

5)

計算個體適應度;

6)

根據適應度，從pop(t) 中選擇優質下一代pop(t+1);

6)

if (random(0,1)

7)

執行交叉操作;

8)

end if

9)

if (random(0,1)

10)

執行變異操作;

11)

end if

12)

end while

13)

選擇具有最佳適應度的一組w和b;

14)

accuracy=0

15)

while(accuracy

16)

使用最優w和b迭代計算輸出結果;

17)

計算準確率;

18)

最小化損失函數;

19)

更新w和b;

20)

end while

21)

return GA_BP網絡結構和參數

2.3 指標因子評價

影響泥石流易發性的因子眾多，在構建泥石流易發性模型時，應該考慮特征因子影響力的大小，但是神經網絡無法從已有的特征中確定特征所占的權重。MIV值可作為確定輸入神經元對輸出神經元影響大小的指標[15]，所以本文使用該算法確定泥石流易發性分類中各個因素與網絡的相關程度，刪減冗余信息，提升分類精度和速度，同時可為研究人員野外勘察提供參考。具體過程如下：

將訓練樣本Xinput中每一個樣本在其原有基礎上加減10%，形成兩個新的訓練樣本A1、A2。將這兩個新的訓練樣本代入已經訓練好的網絡中，得到另外兩個仿真結果B1、B2。將B1-B2的值作為改變該特征對輸出結果產生的影響變化值IV。最后將該影響值按災害觀測數求取平均值。按照上述步驟分別算出n個特征的平均影響值MIV，最后根據平均影響值絕對值的大小對特征進行從大到小排序，并構建包括排序值前1至前13個因子的13種因子組合。利用前述GA_BP算法對比計算13種組合分類效果。算法描述偽代碼如下：

算法3 指標評價與選擇。

輸入 GA_BP網絡,Xinput；

輸出 最優指標組合x*。

1)

fori=0 tondo

2)

A1=x,A2=x;

3)

A1[:,i] =A[:,i]*1.1;

4)

A2[:,i] =A[:,i]*0.9;

5)

分別計算數據集A1和A2在GA_BP的預測結果B1和B2;

6)

IV=B1-B2;

7)

MIV=總體IV的均值;

8)

end for;

9)

按照MIV值從大到小的順序為指標因子排序;

10)

根據排序構建13組輸入組合;

11)

計算各組合在GA_BP 模型的輸出準確率;

12)

根據準確率大小選擇最優組合x*作為輸入;

13)

returnx*

3 實驗與結果檢驗

本次實驗使用Windows版本的TensorFlow學習框架，采用Anaconda軟件包結合Python語言進行程序編寫，實驗平臺使用的CPU為Inter Core i7- 4790 3.6 GHz，Windows7操作系統，運行內存為8 GB。

3.1 評價因子的選取與數據源

根據《泥石流災害防治工程勘察規范》(DZ/T 0220— 2006)[16]，通過野外考察和實地勘測，收集雅魯藏布江中上游的地形地貌、地質構造、地層巖性以及泥石流發生的歷史記錄，初步選取13個影響因子作為影響泥石流易發性評價的主要評價指標，具體指標及表征意義如表1所示。

表1 指標及表征意義

研究區災害點的基礎數據，主要來源于西藏高原大氣環境科學研究所研究人員長年累月的調查與收集，總共統計出有記載的1 778個泥石流地質災害點。本文采用中國地質調查局《滑坡崩塌泥石流災害詳細調查規范》來劃分暴發頻率，如表2。

由暴發頻率圖(圖3)可以看出，雅魯藏布江中上游流域泥石流災害主要以中易發為主。

圖3 研究區域泥石流災害易發性等級統計

3.2 數據預處理

泥石流易發性受孕災環境的多重因子影響，各類別數據之間存在不平衡現象, 因此，需要彌補少數類樣本在信息分布方面的不足。本文使用Borderline-SMOTE算法為少數類樣本：不易發、低易發、高易發插值生成新的樣本，提高分類器對少數類樣本的識別精度。對評價因子進行歸一化處理既是模型建立的需要也是為防止不同指標因單位不同，數據存在偏大或偏小問題。最后將訓練樣本的目標矢量設為(1,0,0,0)、(0,1,0,0)、(0,0,1,0)以及(0,0,0,1)，作為易發性的等級，結果如表3所示, 樣本數均為1 194。

表3 泥石流易發性等級

3.3 確定BP神經網絡結構

依據神經網絡的基本原理，三層神經網絡已經能夠有較好的分類性能。本文使用兩層隱藏層，因為在實驗過程中發現多一層隱藏層能更快地實現收斂，且準確率更高。輸入層根據泥石流易發性評價因子，使用13個神經元。輸出層根據待分類類別數，使用4個神經元。實驗設置初始的權值與閾值為隨機系統值。隱藏層節點數常常使用試湊法確定，使用經驗公式(6)獲取節點數的取值范圍為5～14，使用排列組合的方式，對同一個樣本集進行訓練，最終選取高準確率對應的節點數。實驗基于雅魯藏布江中上游流域泥石流災害歷史數據集，在迭代次數相同的情況下，不同隱藏層節點數訓練后，得到對應的驗證集準確率如圖4所示。

橫坐標為第一層隱藏層節點數，圖例不同形狀表示第二層隱藏層節點數(如：5、6分別表示第二層隱藏層節點個數為5和6)，縱坐標表示分類準確率。分析圖4可知，當隱藏層結構為11- 14時，網絡準確率能夠到達最高值，因此神經網絡的結構設為13- 11- 14- 4。

圖4 不同隱節點數準確率統計

下面對比研究泥石流災害歷史數據集有無Borderline-SMOTE算法處理對模型分類準確率的影響。

如圖5(a)所示，地質數據未經預處理時，網絡驗證準確率在迭代2 000次之后趨于穩定，迭代時間為12.3 s，此時，訓練準確率為94%，但是驗證準確率僅85%左右，類別間數據不平衡導致嚴重過擬合現象。相反，從圖5(b)可以看出，經過算法處理后的數據集，進入網絡迭代至8 000次左右，網絡驗證準確率趨于穩定，訓練時間較長，為23.5 s，此時訓練準確率達96.5%，驗證準確率為94.3%，過擬合現象得到明顯緩和，可見該重采樣算法降低了不平衡數據集對神經網絡的影響，使模型更具有普適性。

3.4 遺傳算法優化神經網絡權值和閾值

GA_BP網絡在訓練時，使用與之前一樣的網絡結構，即：13- 11- 14- 4。對于GA參數，本文依據前人工作[17-18]以及多次反復實驗，選擇種初始群規模為Size=20。最大遺傳代數Gen=20；交叉概率為Pc= 0.6；變異概率為Pm=0.1。訓練過程中，設置最大訓練次數為900，驗證準確率目標為94%。圖6所示為遺傳算法優化后的準確率曲線，模型迭代到300次以后趨于穩定，訓練準確率高達96%，此時驗證準確率為94.03%，運行時間僅為3.25 s。對于神經網絡，計算復雜度最直接的性能評價參數可通過收斂速度直觀衡量[19]。實驗證明，遺傳算法的全局搜索能力使其在指定進化次數內優化了網絡的初始權值和閾值，與3.3節隨機初始化的方式相比較，訓練所需要的迭代次數和時間大幅降低。因此遺傳算法有利于降低神經網絡的計算復雜度，提升收斂速度。

圖6 GA_BP準確率曲線

3.5 計算指標MIV值

根據已經訓練好的網絡，設置訓練集每個特征數據的調節率為10%，即對其分別做乘以1.1和0.9操作，然后計算每個泥石流影響因子的MIV算法值，并按照MIV絕對值從大到小排序，計算結果如表4所示。

表4 各指標MIV

表4中x2、x4、x5、x6、x8、x9、x10、x11、x12呈現正相關，其余特征為負相關。

據泥石流易發性評價因子從大到小排序，構建包括排序值前1至前13個因子的13種因子組合。將每種組合作為模型的輸入，使用上述GA_BP網絡對其訓練，并對包含955個樣本的驗證集進行驗證，記錄模型分類準確率如圖7所示。

圖7 不同因子組合分類結果

從圖7可以看出，評價因子的不同組合對分類準確率影響較大。隨著泥石流指標因素的增加，分類準確率值總體上呈上升趨勢，在到達一定的個數時，趨于穩定。組合1、2、3的準確率相較之下明顯偏低，組合8～13的分類準確率十分接近，其中組合11的分類準確率最高(95.23%)，所以選取因子組合為11的序號，結果顯示，該方法能夠有效優化輸入維度，刪減冗余信息，提升模型精度。

3.6 模型檢驗

由于準確率評估在多類不平衡的分類評價中存在一些不足，所以其不適合作為評判算法優劣的標準。AUC方法相較于準確率具有更強的可區分性，能夠給出不同類的不同分布情況差別[20]。AUC方法定義為受試者工作特征曲線(ROC)下的面積，其值大小作為評判分類器性能的優劣，其中，ROC曲線繪制為真陽性率(True Positive Rate, TPR)與假陽性率(False Positive Rate, FPR)之間的變化關系。本文采用AUC方法，針對13種影響因子，基于未曾使用的400個泥石流災害數據集，采用Logistic、BP、GA_BP、XGBoost[21]以及本文模型進行區域泥石流易發性評價，對比測試結果如圖8和表5所示。

圖8 五種評價模型ROC測試曲線

方法輸入維度各類別比例AUC值/%測試時間/sLogistic1386∶50∶241∶2392.420.0512BP1386∶50∶241∶2395.230.0016GA_BP1386∶50∶241∶2396.900.0041XGBoost1386∶50∶241∶2397.850.0018本文方法1186∶50∶241∶2397.950.0020

選取400處泥石流災害點作為測試樣本進行模型檢驗與比較(表5)，分析數據顯示，本文模型取得97.95%的AUC值，較傳統的Logistic算法提升較大，且處理速度有較大幅度領先；BP和GA_BP模型相比，GA_BP算法有著較高AUC值，但在測試時間上稍慢于BP模型；本文在GA_BP基礎上結合MIV算法，優化了輸入維度，AUC值提高了1個百分點，災害易發性評價速度也提升了接近50%；最后與當下較為熱門的XGBoost分類模型相比，本文算法在AUC值和處理時間上與其幾乎持平。總體來說，本文基于神經網絡的區域泥石流易發性評價模型具備良好的分類性能和運行速度。

4 結語

充分利用地質大數據，基于遺傳神經網絡，結合Borderline-SMOTE算法和MIV算法建立泥石流主要指標與易發程度關系，提高泥石流易發性評價精度，從中挖掘孕災環境中的泥石流災害發生規律，以上研究表明：

1)影響研究區泥石流易發程度的原因主要是不良地質現象、主溝縱坡以及新構造影響。輸入指標存在冗余信息：松散物平均厚度與地層巖性，影響模型的分類精度和速度。

2)Borderline-SMOTE算法有選擇地對少數類數據進行過采樣，能夠有效克服已有泥石流歷史災害數據分布不平衡的困難，緩解了模型的過擬合現象；遺傳算法優化BP神經網絡的權值和閾值，避免了原有網絡隨機初始化的不足，有效提升了泥石流易發性分類的效率和泛化能力；MIV算法能夠從13個特征因子中選擇11個特征因子作為模型的輸入，降低了模型的輸入維度；最后以測試樣本進行仿真，較高的AUC值表明系統可以對泥石流易發性進行分類。本文建立的模型對雅魯藏布江中上游流域泥石流易發性進行分類，能夠根據研究區域內孕災環境的特征，對泥石流易發性作出合理的分類，避免了個人主觀因素對泥石流易發性評價的影響，對當地基礎設施(如：公路、鐵路)建設，具有重要的指導價值。與此同時，本文模型亦存在一些問題值得進一步研究：加入人類活動的相關因子對模型的影響；采用強正則模型解決過采樣可能帶來的局部過擬合，以提升化性能。