基于稀疏表達重構誤差模型的小區域滑坡易發性評價

孫晨昊， 鄭逸榛， 李俊斌， 霍姝涵

(中國地質大學(武漢) 地球物理與空間信息學院，湖北 武漢 430074)

滑坡是地貌演變過程中的一種重要塊體運動形式，是一種典型的重大自然地質災害，給人類生產、生活構成嚴重的威脅[1]。同時，滑坡影響因子本身的隨機性和不確定性以及滑坡發生物理過程的非線性，決定了準確預報滑坡發生時間的極度困難性[2]。

研究滑坡的方法之一是建立地區滑坡的易發性評價體系，通過構建滑坡空間預測模型，對潛在滑坡災害進行定量分析，可較準確地預測特定區域內潛在滑坡發生的空間概率。目前，滑坡易發性預測過程中廣泛使用的機器學習被認為具有比數理統計模型更好的非線性預測能力[3-4]。然而在小區域的滑坡易發性評價中，常用的全監督機器學習方法往往受到滑坡災害點的數量不足與非滑坡災害點的選取不準確的問題而影響到區域機器學習分類器的性能：滑坡災害點數量不足容易導致模型訓練出現“過擬合”現象，即模型在區域的泛化能力不足；非滑坡災害點選取不當極易使模型分類器學到非典型特征，影響分類準確度，得到了諸多錯誤的分類結果。

近年來稀疏表示技術在信號處理、圖像處理、目標識別、盲源分離等領域都有著突出的貢獻，在計算機領域的人臉識別、目標檢測、圖像檢索和恢復等方面取得巨大成功[5]。由于稀疏表示方法能有效解決高維特征空間的維數災難和小樣本問題，遙感領域的研究人員也將稀疏表示理論應用于遙感影像的恢復和重建、分類和目標識別等方面，并取得較好的效果[6]。受此啟發，本文將“稀疏表達”思想引入滑坡易發性評價模型的建立，提出了一種基于稀疏表達重構誤差的半監督滑坡易發性評價模型，用已有的滑坡災害點來構建“滑坡字典”，無需選取非滑坡像元作為樣本，充分利用已有的滑坡樣本，較好地解決了小區域、少樣本情況下的滑坡易發性評價問題，從而達到更為高效利用滑坡歷史數據信息及提高滑坡易發性評價精度的目的。

1 評價模型

1.1 稀疏表達

稀疏表達的提出源于神經生物學，可追溯至1959年，人們發現主視皮層V1區神經元的感受也能對視覺感知信息產生一種“稀疏表示”[7]，至此稀疏表示理論開始受到研究者的關注。1996年，Olshausen和Field[8]在《自然》雜志上發表開創性的文章，他們發現通過稀疏編碼學習出來的字典具有人類視覺真實感受的特征。由此稀疏編碼和稀疏表達理論開始被研究人員廣泛關注。近年來稀疏學習理論逐漸趨于成熟，在信號處理、圖像處理、人工智能等領域都有了廣泛應用。

稀疏表達的核心是：用較少的基本樣本的線性組合來表達大部分或者全部的原始樣本。這些基本樣本被稱作原子，而過完備字典則是由個數超過特征維數的原子聚集而來。稀疏表達的根本意義在于：稀疏性，即低復雜度，這是物理世界的普遍特征。而稀疏性即低復雜度的表現，稀疏表達會更加體現系統的核心特質，也是系統結構的更優化的表示。稀疏性假設則是引入低復雜度對可能的模型進行加權篩選。

1.2 稀疏表達重構誤差模型

本文基于上述稀疏表達原理，提出了稀疏表達重構誤差模型，將區域內已有的部分滑坡像元組成“滑坡字典”，用“滑坡字典”重構出的像元與原有像元間的誤差，即“重構誤差”來刻畫衡量區域內像元與滑坡特征模型的相似程度，以進行區域的滑坡易發性評價。具體步驟如下。

首先將研究區的所有像元作為數據集{x1,x2,x3…xn}。若用二維數組表示，則每一列為一個像元，每一行分別為高程、坡度、坡向等10項影響因素。緊接著用滑坡災害點像元及其特征屬性構建“滑坡字典”。找到合適的字典B，最后用盡可能少的字典原子來線性表示每一個像元，則像元的稀疏表達形式為：

xi=Bαi+e

(1)

圖1和式(1)中：d=10，為像元特征數；k為字典的詞匯量，也就是選取當作字典的滑坡數量；αi為某一像元的稀疏表達稀疏；e為對應的稀疏表達殘差。不難證明，任一像元的稀疏表達系數αi不是唯一的，要找到最稀疏的那個組合以及稀疏表達系數αi，這就遇到了多項式復雜程度的非確定性問題(Non-deterministic Polynomial)，難以求解。用數學符號表示，筆者希望找到一個αi，使其含有最多的0元素，即αi的0范數最?。?/p>

圖1 稀疏表達原理圖

(2)

s.t.xi=B·αi

(3)

這個命題可以等價轉化為如下形式，找到一個最多只有η個非零元素的稀疏系數αi，使誤差項e最小：

(4)

s.t.‖αi‖0≤η

(5)

這里可以采用正交匹配跟蹤算法(Orthoganal Matching Pursuit)[9]求解某一像元對應滑坡像元字典的最稀疏系數。正交匹配跟蹤算法的原理如圖2所示。

圖2 正交匹配跟蹤算法原理圖

最后，可以知道，選取了區域內的滑坡像元作為字典，那么對于區域內的某一像元而言，若其因子特征與滑坡像元的特征越相似，用“滑坡字典”來稀疏表達后的誤差越小；若某一像元的因子特征與滑坡發育特征大相徑庭，則用“滑坡字典”來稀疏表達后的誤差是越大的。當然，這需要“滑坡字典”是一個過完備字典，即包含有較多的滑坡災害點，記錄區域內的普遍滑坡災害點發育特征。于是就可以用稀疏表達后存在的誤差，即“重構誤差”來刻畫衡量區域內像元與滑坡特征模型的相似程度，進行區域的滑坡易發性評價。

1.3 BP神經網絡模型

神經網絡屬于動態非線性系統，它適用于處理背景知識不清楚，推理規則不明確，模糊、隨機、復雜的信息識別問題[10]，人工神經網絡能為客觀反映滑坡災害地質體穩定性程度提供一條有效途徑。

BP神經網絡模型是1986年由Rumelhart和McClelland為首的科學家小組提出[11]，是目前應用最廣泛、最成功的神經網絡模型之一。BP神經網絡模型能學習和存貯大量的輸入—輸出模式映射關系，而無需事前揭示描述這種映射關系的數學方程。它的學習規則是使用梯度下降法，通過反向傳播來不斷調整網絡的權重和閾值，使網絡的誤差平方和最小。

BP神經網絡模型拓撲結構包括輸入層(input layer)、隱藏層(hide layer)和輸出層(output layer)。輸入層中各神經元負責接收來自外界的輸入信息，并傳遞給中間層各神經元；中間層是內部信息處理層，負責信息變換，根據變換能力的需求，中間層可以設計為單隱層或者多隱層結構，其中隱藏層的節點個數經驗公式為：

S=1+[m×(n+2)]0.5

(6)

式中：S為隱藏層的節點個數；m為輸入節點數；n為輸出節點數。最后一個隱藏層傳遞到輸出層各神經元的信息，經進一步處理后，完成一次學習的正向傳播處理過程，由輸出層向外界輸出信息處理結果。

1.4 信息量模型

信息量模型以信息論為基礎，采用地質災害發生過程中熵的減少來表征地質災害事件產生的可能性，以影響地質災害發生的因素為影響因子，通過計算各因子的信息量，進而將信息量單因素或加權疊加，建立易發性評價模型，對易發性作出評價[12]。對研究區多因子共同作用下的綜合信息量的計算，可將各因子信息量進行單因素信息量疊加。如果用Ii表示影響因子xi的信息量，則有式(7)：

(7)

實際計算時可用統計頻率估計條件概率來估算信息量，則有式(8)：

(8)

式中：S為研究區評價單元總面積；N為研究區含有地質災害分布的地質災害點總數；S0為研究區內含有影響因子x1x2…xn組合的單元總面積；N0為具有相同因子組合x1x2…xn的特定類別內的地質災害點總數。上式可得到式(9)：

(9)

式中：I為某評價單元信息量預測值；Ni為影響因子xi占有的評價單元中的地質災害點總數；Si為含有影響因子xi的評價單元的面積。

2 萬州區滑坡易發性評價

2.1 研究區概況

研究區為重慶市東北部萬州區部分區域，位處三峽庫區中段長江上游，面積約為404 km2。受亞熱帶季風影響，該區氣候溫暖濕潤，降雨豐富，多年平均年降雨量為1 194.4 mm。萬州區位于四川盆地東北部，地形以山地丘陵為主，東南部地勢較高，長江兩岸地勢較低。境內的鐵峰山、方斗山等屬世界三大褶皺山系之一的川東平行嶺谷，呈現出背斜成山、向斜成谷、山谷相間、彼此平行的褶皺地貌。區域內巖層主要為侏羅系和三疊系地層，巖性以砂巖、泥巖、頁巖、灰巖為主，且人工堆積、沖積、殘積、坡積物等土層分布也較為廣泛。受構造作用影響，萬州區褶皺發育密集，但并無大的斷裂構造。區內地表水系發育，長江、數十條支流及小型溪溝構成了復雜的地表徑流網絡。正是由于上述復雜的工程地質條件，使得全區滑坡災害十分嚴重，對人民生命和財產安全威脅巨大。

圖3 萬州區地質概況圖

2.2 數據來源

本次采用的主要數據來源見表1。

表1 數據來源

2.3 指標體系與相關性分析

滑坡是內因和外因共同作用的結果，其發育受多種因素影響[13]。滑坡易發性評價影響因子的選擇通常根據研究區的滑坡發育特征而定，同時需要分析各影響因子之間的相關性，以保持評價因子的相對獨立性。本文提取了高程、坡度、坡向、地形濕度指數(TWI)、曲率、地層巖性、距構造線距離、距河流距離、距道路距離與土地利用類型共10個因子。

為確保各評價因子的相對獨立性，保證評價模型的可靠性和精度，需對影響因子進行相關性分析?；赟PSS軟件，采用Pearson(皮爾遜)相關系數分析影響因子的相關性，當系數值>0.5時，認為因子之間存在相關性，對評價結果準確性將產生影響[14]。計算結果如表2所示。

通過觀察表2以及查閱相關系數界值表可得各指標皮爾遜相關系數值最大為0.469，絕對值均<0.5，說明各指標相對獨立性較強，可用于滑坡易發性評價。

表2 因子相關系數表

采用高程、坡度、坡向、TWI、曲率、地層巖性、距構造線距離、距河流距離、距道路距離與土地利用類型共10個因子，繪制災害點在各因子數據中的分布頻率直方圖(圖4-圖13)，根據災害點頻率突變情況對各因子進行分級，各因子提取的圖層與分級如表3所示。

表3 各因子分級統計表

圖4 高程因子分級圖

圖5 坡度因子分級圖

圖6 坡向因子分級圖

圖7 地形濕度指數(TWI)分級圖

圖8 曲率因子分級圖

圖9 地層巖性因子分類圖

圖10 距構造線距離分級圖

圖11 距河流距離因子分級圖

圖12 距道路距離因子分級圖

圖13 土地利用類型因子分類圖

2.4 滑坡易發性評價

2.4.1稀疏表達重構誤差模型

基于ArcGIS10.7與MATLAB平臺，選取100 m×100 m柵格作為評價單元，研究區共劃分40 401個柵格單元。收集和現場調查滑坡樣本共計128個。具體評價流程如圖14所示。

圖14 評價流程圖

(1) 字典構建。首先將滑坡災害點數據及各評價因子數據導入。隨機選取約70%的滑坡像元來構建“滑坡字典”。這里的70%數量是經過慎重考慮的：數量過多可能導致“滑坡字典”所記錄的滑坡特征冗余，數量過少則難以構成過完備字典(即字典中的滑坡像元個數要遠大于特征行數)。之后也嘗試了將所有的滑坡像元納入字典中。

(2) 字典與數據的規范化。將字典與數據進行規范化，這是進行稀疏系數求解的必要步驟。

(3) 稀疏系數求解與重構誤差。對于圖像中的每一個規范化后像元，運用正交匹配跟蹤算法求解在“滑坡字典”D表示下的稀疏系數，同時計算其重構誤差；計算出來的重構誤差為規范化后像元每個特征的重構誤差，之后需要求其F范數作為像元總的重構誤差，將其作為像元的值。

(4) 還原圖像與精度計算。以重構誤差為評價指標，將所得的重構誤差結果還原為199×199圖像格式，并導入ArcGIS進行萬州區域滑坡易發性制圖，此處采用的分級方法為快速聚類分段法。之后繪制ROC曲線。將70%滑坡樣本作為字典，利用快速聚類法[15]得到的聚類中心為0.104，0.168，0.225，0.29，0.568，由此得到的五個區劃等級的分段點為0.136 3，0.196 8，0.262，0.433；得到的區劃等級為：極高易發區[0～0.136 3),高易發區[0.136 3,0.196 8),中易發區[0.196 8，0.262),較低易發區[0.262，0.433),低易發區[0.433,0.788]。

2.4.2信息量模型

經過因子選擇與分析，得到了影響滑坡災害發育的10個因子及其圖層，隨后利用信息量模型計算各因子的信息量值，得到每個二級指標的狀態指標的信息量。將每個因子的信息量圖層進行疊加，得到最終的結果圖層。最后采用快速聚類法[15]確定分級閾值，將結果圖層的信息量值進行五級分類處理，通過聚類得到的分級閾值為：-2.279,-1.083 8,-0.033和1.087 8。由此得到的五個分級為：低易發區[-6.13,-2.279)，較低易發區[-2.279，-1.083 8)，中易發區[-1.083 8，-0.033)，高易發區[-0.033,-1.087 8)，極高易發區[-1.087 8，2.99]。

2.4.3BP神經網絡模型

BP神經網絡(back propagation neural network，BP)具有很強的非線性映射能力[16]，可用于滑坡易發性評價。將各類滑坡影響因子作為網絡輸入層的變量，輸出層則對應危險程度的具體數值表達[17]。故筆者基于MATLAB編程，使用了BP神經網絡模型來進行萬州區的滑坡易發性評價。

首先進行數據準備，將各因子圖層從ArcGIS中導出，將導出的各因子圖層和災害點圖層導入MATLAB中；然后對數據進行預處理，將二維因子圖層轉為一維的行向量，并且對因子歸一化，這能夠避免由于各因子取值范圍的差異對BP神經網絡性能造成的影響。隨后進行標簽準備，將滑坡樣本與非滑坡樣本進行樣本標記；再選取訓練樣本，選取和已知滑坡柵格單元同等數量的非滑坡柵格單元，這里采用隨機選取非滑坡像元，將相同數量的滑坡與非滑坡像元構成訓練樣本。然后進行BP神經網絡訓練，70%為訓練集，15%為測試集，15%為驗證集。根據經驗公式，這里設置的隱藏層神經元個數為7個；最后依照分類結果進行圖像還原，采用快速聚類法[15]將滑坡易發性結果圖分為高易發區、較高易發區、中易發區、較低易發區與低易發區5類。

2.5 易發性評價結果及分析

滑坡易發性評價結果見圖15-圖18，易發性共分為5個等級，分別為低、較低、中、高、極高。從圖中可以明顯看出，歷史滑坡點基本上都位于滑坡中等及以上易發區的范圍，大多沿長江兩岸分布，遠離河流地區為滑坡的較低易發區和低易發區。

圖15 萬州區滑坡災害易發性分級圖(70%字典)

圖16 萬州區滑坡災害易發性分級圖(100%字典)

圖17 信息量易發性分級圖

圖18 BP神經網絡易發性分級圖

2.5.1精度評價

根據易發性評價結果，統計滑坡各易發性等級所占面積、滑坡各易發性等級面積占研究區總面積的比例、滑坡數量、滑坡數量占總滑坡比例、滑坡面積和滑坡面積占研究區總面積的比例，詳細數據見表4-表7。

表4 稀疏表達重構誤差模型(70%字典)不同等級區劃與滑坡分布統計表

表5 稀疏表達重構誤差模型(100%字典)不同等級區劃與滑坡分布統計表

表6 信息量模型不同等級區劃與滑坡分布統計表

表7 BP神經網絡模型不同等級區劃與滑坡分布統計表

觀察統計表格，從滑坡數量的分布來看，信息量模型和BP神經網絡模型的結果在較低易發區與低易發區占比較少，滑坡大多分布在中易發區等級以上，說明了這兩種模型有合理的預測效果，但是滑坡存在于極高易發區的比例均不是特別理想，相比之下，基于70%滑坡像元構建字典的稀疏表達重構誤差模型在區域上取得了令人滿意的結果。為了更好地比較三者的評價結果，通過計算兩種模型易發性分區圖的每個類別中存在滑坡的個數與每個類別總數的比值，即特定類別精度[18]，其公式如下：

(10)

式中：Pj為特定類別精度，j=1，2…M(M為滑坡易發性分區類別的總數)；Aj和Bj分別為第j類滑坡易發性分區中存在滑坡的斜坡單元個數和斜坡單元總數。三種模型的特定類別精度如圖19所示。由圖19可見，稀疏表達重構誤差模型在高易發區的特定類別精度遠遠高于其他兩種模型，由此說明，稀疏表達重構誤差模型預測的滑坡易發性分區中高易發區等級以上包含了更多先前調查的滑坡。

圖19 三種模型特定類別精度圖

其中基于100%滑坡像元構建字典的稀疏表達重構誤差模型所劃分的極高易發區所占面積比為20.04%，略低于70%滑坡像元構建字典的稀疏表達重構誤差模型所劃分的極高易發區所占面積比(21.73%)，但卻包括了所有的滑坡災害點。這個結果值得深思：將所有的滑坡災害點用來構建字典，則所有滑坡災害點所在像元的重構誤差往往也取到極小值，故都落入極高易發區；換而言之，此時的“滑坡字典”記錄了128個滑坡像元的所有模式，并且將區域上其他像元的特征模式與這128個滑坡像元特征比較，在一定程度上提高了字典的表達能力；然而這樣做也存在一個比較大的缺陷，現實中存在如下情況：區域內存在偶然發生的滑坡災害點，也就是說，128個滑坡災害點并非全部在極高易發區，部分災害點是偶然發生的，所記錄的特征也不具備識別滑坡的可靠性，而是一種“噪聲”，故將全部的滑坡像元充當字典是不妥的。

選用受試者工作特征(receiver operating characteristic，ROC)曲線線下面積(area under curve，AUC)、總體精度(overall accuracy，OA)和Kappa系數對三種模型進行評價。對曲線來說，曲線由圖幅面積百分比和滑坡面積百分比確定，線下面積越大，越接近1，說明滑坡預測精度越高[19]；而對OA與Kappa系數來說，指標值越大，模型精度越高。

如圖20所示，稀疏表達重構誤差模型、信息量模型及BP神經網絡模型三種模型的線下面積(AUC)分別為0.834 3、0.655 2和0.680 3，三種模型中僅有稀疏表達重構誤差模型AUC值>0.8，表明預測精度較高，較其他兩種模型性能有明顯的提升。除此以外，選用所有歷史滑坡像元與隨機選取的等量非滑坡像元進行精度評價，模型評價結果中極高易發區像元被認為是滑坡像元，其余像元則被認為是非滑坡像元，求得三種模型的混淆矩陣、OA及Kappa系數(表8)。由表8可以得知，稀疏表達重構誤差模型的OA和Kappa系數均優于信息量模型和BP神經網絡模型，表明稀疏表達重構誤差模型在小范圍的滑坡易發性評價方面有著良好的表現。上述結果說明稀疏表達重構誤差模型是有效的，對小區域滑坡預測精度的提高具有一定幫助。

圖20 模型ROC曲線圖

表8 模型精度評價表

2.5.2模型穩定性分析

本文為了檢驗稀疏表達重構誤差模型的穩定性，選取相同樣本，在保持所有參數設置相同的條件下連續進行5次實驗，得到OA、Kappa系數和AUC三個指標的穩定性分析結果(表9)。由表9可知，本文構建的稀疏表達重構誤差模型在滑坡易發性評價中平均OA為78.86%，平均Kappa系數為0.577，平均AUC為0.828，穩定性較好。

表9 模型穩定性分析表

3 結論

滑坡災害易發性研究在滑坡災害風險管理與城市規劃等方面具有非常重要的現實意義。本文結合萬州地區已有的滑坡地質災害點數據，基于GIS對其滑坡地質災害的影響因子進行了分析評價。選取高程、坡度、坡向、曲率、地形濕度指數(TWI)、距構造線距離、距道路距離、地層巖性、距河流距離以及土地利用類型10項因子參與滑坡易發性模型建立。

本文將“稀疏表達”的思想引入滑坡易發性評價，避免了非滑坡像元的選擇；小范圍隨機選取70%已有滑坡災害點構建“滑坡字典”，通過重構誤差來衡量某一像元與滑坡特征模式的相似程度，最后采用了快速聚類法進行易發性程度分段。除此之外，為了對比突出該模型的合理性，又使用了常用的信息量模型以及經典的機器學習模型——BP神經網絡，從統計分析、制圖結果、ROC曲線和模型穩定性等多角度評價結果表明，稀疏表達重構誤差模型(70%字典)是一種更行之有效的方法。AUC值高達0.834 3，遠優于信息量模型的0.655 2以及BP神經網絡模型的0.680 3，除此外OA和Kappa系數也均優于信息量模型和BP神經網絡模型，分別為78.1%和0.562，并且模型穩定性較好，制圖結果與實際更符合。

可以說，將“稀疏表達”思想引入滑坡易發性評價模型的建立，是一次大膽且創新的嘗試，取得良好結果的同時也值得進一步挖掘其背后更深層次的原理本質。未來的工作將圍繞排除字典內偶然滑坡像元，選取更具代表性的“滑坡字典”原子以及找到一個最優的“滑坡字典”詞匯量，而不是范式地選取70%滑坡像元。