郭亞雷,鄧念東,李宇新,周 陽,石 輝

(1. 西安科技大學 地質(zhì)與環(huán)境學院, 陜西 西安 710054; 2. 陜西省地質(zhì)調(diào)查院,陜西 西安 710043; 3. 陜西省水工環(huán)地質(zhì)調(diào)查中心,陜西 西安 710068)

0 引言

滑坡作為嚴重的自然災(zāi)害,對人們的生命和財產(chǎn)安全造成了巨大的威脅[1-2]?；乱装l(fā)性評價指一定區(qū)域內(nèi)滑坡發(fā)生可能性的預(yù)測,被認為是對區(qū)域滑坡進行防治措施的首要步驟以及土地空間規(guī)劃的有效借鑒。因此,科學嚴謹、高精度的滑坡易發(fā)性評價意義重大。目前,國內(nèi)外學者已進行大量滑坡易發(fā)性評價相關(guān)的研究,評價方法主要包括定性和定量2種方法,隨著評價方法研究的深入,定性方法展現(xiàn)出主觀性大的缺點,逐漸被定量方法所替代[3]。隨著計算機人工智能的飛速發(fā)展,一系列經(jīng)典的機器學習與數(shù)據(jù)發(fā)掘算法已經(jīng)應(yīng)用至滑坡易發(fā)性評價,并取得了較優(yōu)的預(yù)測結(jié)果,諸如人工神經(jīng)網(wǎng)絡(luò)[4]、支持向量機[5]、樸素貝葉斯[6]、決策樹[7]等。田述軍等[8]基于不同評價單元,論述了斜坡單元和網(wǎng)格單元對滑坡易發(fā)性評價結(jié)果的影響。李文彥等[9]對滑坡易發(fā)性不同評價模型進行對比,并驗證了各模型的精度。由于滑坡發(fā)生機理復雜,影響因子數(shù)據(jù)與標簽的關(guān)系通常呈非線性,以上機器學習算法在分析不同區(qū)域時仍存在一定適用性的差異,并且過擬合是機器學習算法進行分類時難以規(guī)避的問題。集成學習是指將多個分類器組合來共同解決分類或回歸任務(wù)的模型,分類器組成結(jié)構(gòu)主要分為同源集成與異源集成2類[10]。研究表明,集成學習能有效克服單一機器學習模型自身的局限性,具有更高的泛化能力[11]。其中Adaboost[12]、Bagging[13]、隨機子空間(random subspace, RS)[14]、隨機森林[15]等集成學習方法已應(yīng)用于國內(nèi)外滑坡易發(fā)性評價研究中。為彌補不同集成學習各自預(yù)測的缺陷,堆疊Stacking模型可將多種模型進行組合,從而提高分類精度。同時,非滑坡的選取影響著樣本數(shù)據(jù)純度,文中對評價因子進行滑坡密度分析與篩選,在滑坡分布稀疏的區(qū)域與緩沖區(qū)疊加選取負樣本。

1 研究方法

1.1 RA模型

RA模型屬于Boosting算法族,由于Boosting算法通常使用弱分類器,使得其個體學習器之間存在強依賴關(guān)系。 RA是Freund和Schapire提出的一種解決二分類問題的集成學習算法。主要目標為將弱學習器“提升”為強學習器,其核心思想是根據(jù)訓練集的權(quán)值分布來選擇各分類器所使用的子集,對權(quán)值進行調(diào)整。對集成系統(tǒng)中各個分類器的分類結(jié)果使用加權(quán)多數(shù)投票進行組合,使具有較高訓練精度的基分類器在投票時具有更大的權(quán)重,具體流程為先對每個樣本賦予相同的初始權(quán)重,每一輪學習器訓練過后都會根據(jù)其表現(xiàn)對每個樣本的權(quán)重進行調(diào)整,增加對分錯樣本的權(quán)重,從而在后續(xù)訓練中加強對分錯樣本的學習,按這樣的過程重復訓練出多個學習器,進行加權(quán)組合。最后將對應(yīng)的權(quán)值之和最大的那個類別作為分類結(jié)果(圖1)。

圖1 RA模型分類流程圖Fig. 1 Classification flow chart of RA model

1.2 DECORATE模型

DECORATE(diverse ensemble creation by oppositional relabeling of artificial training examples,DECORATE)是Melville和Mooney在2005年提出的一種集成學習算法,屬于同源集成的一類。與Bagging和Boosting算法族的一個顯著的區(qū)別是,該算法的基本思想是通過利用人工訓練樣例集來生成一些多樣化的分類器。研究表明,DECORATE比Boosting對冗余樣本數(shù)據(jù)具有更強的魯棒性,比Bagging對缺失屬性值樣本具有更好的容忍性[16]。該算法首先統(tǒng)計訓練集的分布特征,對于連續(xù)性屬性計算其均值與標準差,基于這2項指標得到一定人工樣本數(shù)據(jù);對于離散值統(tǒng)計其出現(xiàn)的頻率,根據(jù)頻率同樣隨機選取一定量數(shù)據(jù)。對上述人工樣本數(shù)據(jù)通過概率元組的方式進行類別標記,由于類標簽與集成系統(tǒng)分類結(jié)果的差異,促進了其個體分類器之間多樣性。在每次迭代中,在擴展訓練集生成一個分類器后,將擴展集成系統(tǒng)與原集成系統(tǒng)分類的訓練精度進行對比以及不斷進行篩選,來保證該模型分類精度不會降低。

1.3 RS模型

RS模型隸屬于Bagging算法族的一類,是基于對訓練集隨機抽樣的一類集成學習算法[17]。首先,從訓練集的屬性集合(A1,A2,…,An)隨機選取k個屬性(1

1.4 隨機森林模型

隨機森林主要是通過多棵決策樹中每棵樹的投票結(jié)果來獲取最優(yōu)的分類結(jié)果。每棵決策樹因有回放的方式以及隨機獲取數(shù)據(jù)特征所得到的數(shù)據(jù)集而具備更全面的輸入變量信息。通過多棵決策樹的集成來實現(xiàn)模型的魯棒性并避免過擬合。隨機森林主要特征是能夠給出相應(yīng)輸入變量的重要性排序。隨機森林中用不純度來度量最佳分割,通過環(huán)境因子k在節(jié)點分割時的基尼指數(shù)的減少值DGk實現(xiàn)基礎(chǔ)環(huán)境因子重要性計算,其中涉及平均基尼減小值占所有基礎(chǔ)環(huán)境因子平均基尼減小值綜合的百分比計算,具體如式(1):

(1)

式中:m,n,t分別為基礎(chǔ)環(huán)境因子總數(shù)、分類樹棵數(shù)和單棵樹節(jié)點數(shù);DGkhj為k個因子在第h棵樹的第i個節(jié)點上基尼指數(shù)減小值;Pk為第k個基礎(chǔ)環(huán)境因子的重要性。

1.5 Stacking模型

Stacking模型是一種異構(gòu)分類器集成的模型,這是與同源集成模型主要的區(qū)別。模型由2層框架組成:第1層由RA、DECORATE、RS多個基分類器組成;第2層為處理第1層輸出結(jié)果的元分類器。首先將數(shù)據(jù)集分成訓練集和測試集,利用訓練集訓練得到多個初級學習器,然后用初級學習器對測試集進行預(yù)測,并將輸出值作為下一階段訓練的輸入值,最終的標簽作為輸出值,用于訓練次級學習器,再讓次學習器給基分類器模型的結(jié)果分配權(quán)重,進行重復訓練后,將10倍交叉驗證后的基分類器訓練結(jié)果的概率分布,作為元分類器線性回歸模型的輸入,線性回歸模型對每個類學習了一個分類結(jié)果隸屬度的線性回歸函數(shù),歸一化后作為分類概率,最后得到集成模型的分類結(jié)果(圖2)。由于每次所使用的訓練數(shù)據(jù)不同,因此可以在一定程度上防止過擬合。

圖2 Stacking模型流程圖Fig. 2 Flow chart of Stacking model

2 研究區(qū)概況及數(shù)據(jù)源

丹鳳縣位于秦嶺東段南麓, 地處陜、 豫、 鄂三省交界之丹江通道上段。縣域東西長62.1 km, 南北寬

65.5 km,總面積為2 438 km2。屬于北亞熱帶向暖溫帶過度的季風性半濕潤山地氣候區(qū),氣候溫和,四季分明,平均氣溫13.8℃,年平均降雨量687.4 mm,年日照時間2 056 h,無霜期217 d。地勢西北高、東南低,高程介于324～2 011 m,相對高差1 687 m,縣域內(nèi)分布著自北而南的3條山脈;發(fā)育有丹江、銀花河、武關(guān)河和老君河4條河流河谷相間,呈“掌”狀地貌。通過研究區(qū)地質(zhì)災(zāi)害詳細調(diào)查,共圈定257處滑坡,通過GIS轉(zhuǎn)化為點映射至區(qū)內(nèi)(圖3)。

①優(yōu)：骨折完全愈合且功能正常，無痛、無畸形，生活自理能力同骨折前。②良：骨折愈合，無痛、無畸形，但存在10°～20°活動度之差，生活能自理。③可：骨折略有畸形愈合，無痛，功能有所恢復，生活可部分自理。④差：骨折愈合延遲，有疼痛感，功能受限，生活不能自理。

圖3 滑坡編錄圖Fig. 3 Landslide cataloging diagram

文中通過“地理空間數(shù)據(jù)云”獲取研究區(qū)DEM數(shù)字高程數(shù)據(jù)和Landsat 8遙感衛(wèi)星數(shù)據(jù);利用Bigemap地圖軟件下載1∶5萬地質(zhì)圖以及道路、水系矢量數(shù)據(jù)。根據(jù)DEM數(shù)據(jù)在ArcGIS軟件中生成坡度、坡向、曲率、地形濕度等指數(shù)因子,為模型數(shù)據(jù)庫的構(gòu)建奠定基礎(chǔ)。

3 評價模型數(shù)據(jù)預(yù)處理

3.1 評價單元與影響因子的選取

選取合適的評價單元是滑坡易發(fā)性評價的基礎(chǔ),主要劃分為柵格單元、斜坡單元、地形單元[18]。由于柵格單元具有數(shù)據(jù)結(jié)構(gòu)簡單、計算機處理高效的特點,更適用于集成學習模型的訓練與驗證。根據(jù)湯國安經(jīng)驗公式[19],文中采用30 m×30 m大小的柵格作為評價單元,將研究區(qū)共劃分為2 670 541個柵格。結(jié)合研究區(qū)地質(zhì)環(huán)境背景,初步選取高程、坡度、坡向、剖面曲率、平面曲率、TWI、NDVI、距河流距離、距道路距離、距斷層距離、地層巖性和年均降雨量共12個影響因子。對連續(xù)型因子分別采用Jenks自然間斷法和等間距法進行分級;離散型因子根據(jù)二級因子類型進行劃分(表1)。

表1 滑坡密度分析結(jié)果Table 1 Landslide density analysis results

續(xù)表

3.2 模型數(shù)據(jù)集構(gòu)建

集成學習訓練與驗證的數(shù)據(jù)庫由正負樣本組成,通常將滑坡影響因子數(shù)據(jù)作為正樣本,選取等量的非滑坡區(qū)域提取負樣本數(shù)據(jù)。因此,非滑坡區(qū)域的選取直接影響到模型擬合的效果[20]。目前,主要有以下選取方法:1)在滑坡周界外一定距離作緩沖區(qū)處理,在其余范圍進行隨機選取。該方法不足之處在于不同研究區(qū)緩沖區(qū)閾值難以確定。2)在特定的非滑坡區(qū)域進行選取,比如水系區(qū)域、坡度小于2°的區(qū)域等[21]。該方法存在負樣本選取集中、數(shù)據(jù)覆蓋不全面的問題,易造成分類模型過擬合,從而導致結(jié)果準確率降低。

為有效避免因數(shù)據(jù)冗余造成模型分類精度下降,以及非滑坡因子屬性集中造成的過擬合問題,文中通過統(tǒng)計各因子二級分類的滑坡數(shù)量與滑坡密度(表1),剔除滑坡數(shù)量與密度同時最高的二級分類范圍,選取其余區(qū)域與剔除滑坡緩沖區(qū)范圍進行疊加,選取范圍如圖4所示。在該區(qū)域隨機生成257個非滑坡點并提取因子屬性信息。將正負樣本按照7∶3隨機劃分為訓練集與驗證集,分別包括179和78個樣本。

圖4 非滑坡選取范圍Fig. 4 Selection range of non-landslide

3.3 影響因子的篩選

在R語言中利用隨機森林模型計算出各因子的權(quán)重值,對12個因子進行重要性排序,結(jié)果見圖5。篩選出的高程、坡度、坡向、剖面曲率、平面曲率、距道路距離、距斷層距離、距河流距離、TWI、NDVI、年均降雨量、地層巖性12個重要因子中,地層巖性權(quán)重值為0.029 8,遠小于其他類因子權(quán)重值,其重要性最低,因此剔除影響較弱的地層巖性因子。

圖5 各因子權(quán)重分布圖Fig. 5 Weight distribution of each factor

對模型的訓練集進行相關(guān)性和共線性分析,可以降低因數(shù)據(jù)間高度相關(guān)或共線對模型分類精度的影響。文中采用皮爾遜相關(guān)性(PCC)、方差膨脹因子(variance inflation factor, VIF)和容忍度(tolerance, TOL)進行分析,其中TOL為VIF的倒數(shù)。通常認為PCC大于0.5或VIF大于2時,數(shù)據(jù)之間存在較強烈的相關(guān)性或較嚴重的共線性,需要進行剔除[22]。同時采用相關(guān)屬性評估(correlation attribute evaluation,CAE)進行10倍交叉驗證,分析研究區(qū)滑坡發(fā)生的關(guān)聯(lián)度,其值越大代表該因子與區(qū)內(nèi)滑坡發(fā)生更密切。

由表2和表3結(jié)果可以看出,距水系距離與距道路距離的皮爾遜相關(guān)性為0.591,同時距水系距離VIF值為2.667,因此剔除距河流距離因子。根據(jù)圖6分析結(jié)果,其余10個因子對研究區(qū)滑坡均有一定作用。最終選取高程、坡度、坡向、剖面曲率、平面曲率、距道路距離、距斷層距離、TWI、NDVI、年均降雨量共10個因子作為一級指標進行評價。

表2 皮爾遜相關(guān)性指標Table 2 Pearson correlation indicators

表3 影響因子共線性分析Table 3 Collinearity analysis of impact factors

圖6 影響因子CAE分析Fig. 6 CAE analysis of impact factors

4 滑坡易發(fā)性評價

4.1 基于基分類器模型的滑坡易發(fā)性評價

文中采取的基分類器分別為RA、DECORATE和RS模型,均在WEKA3.8軟件中進行生成。對RA模型進行參數(shù)設(shè)置,其中選擇決策樹為基分類器,迭代次數(shù)為100次,收縮參數(shù)為0.1;DECORATE模型選擇J48樹算法為基分類器,其中用于剪枝的置信因子設(shè)為0.25,成員分類器數(shù)量選擇15個,迭代次數(shù)為50;RS模型采用REP樹模型作為基分類器,每一個子空間大小為0.5,迭代次數(shù)為100。將訓練集代入RA、DECORATE和RS這3種集成模型中,通過10倍交叉驗證分別得到訓練正確率為74.7%、69.9%和74.9%,代入驗證集數(shù)據(jù)得到預(yù)測率分別為76.5%、73.8%和72.6%。最終,將研究區(qū)11個因子屬性的2 670 541個柵格代入3種模型生成滑坡易發(fā)性指數(shù)(landslide susceptibility index, LSI)。根據(jù)自然間斷法,將LSI分為極低易發(fā)區(qū)、低易發(fā)區(qū)、中等易發(fā)區(qū)、高易發(fā)區(qū)和極高易發(fā)區(qū)5類[23](圖7)。

4.2 基于Stacking模型的滑坡易發(fā)性評價

同樣采用WEKA軟件構(gòu)建Stacking模型,其中基分類器使用上述RA、DECORATE以及RS模型,分類器參數(shù)選擇與單獨訓練時一致;元分類器選用線性回歸模型。通過10倍交叉驗證進行訓練,得到訓練正確率為75.5%。代入驗證集數(shù)據(jù)得到預(yù)測率為77.1%,最終生成研究區(qū)各柵格LSI值,重分類生成滑坡易發(fā)性分區(qū)圖(圖8)。

通過比較4種模型滑坡易發(fā)性分區(qū)結(jié)果,可以看出區(qū)域等級劃分趨勢基本一致。極高易發(fā)區(qū)主要集中在研究區(qū)中南部,極低易發(fā)區(qū)主要分布于研究區(qū)北至東北部。圖9對4種模型易發(fā)性分區(qū)進行統(tǒng)計,從各易發(fā)性等級對應(yīng)的滑坡密度可以看出,4種模型滑坡頻率比均隨著易發(fā)性等級提高而上升,說明分區(qū)結(jié)果與事實相符。其中,Stacking模型高至極高易發(fā)區(qū)的滑坡頻率比為2.932,高于其余模型(RA模型為2.192,DECORATE模型為2.825以及RS模型為2.821),說明Stacking模型對研究區(qū)滑坡預(yù)測更為敏感。

圖9 滑坡易發(fā)性等級分區(qū)對比Fig. 9 Comparison of landslide susceptibility grade zones

5 模型驗證與討論

文中采取受試者工作特性曲線(receiver operating characteristic,ROC)及其線下面積(area under curve,AUC)對4種模型訓練與驗證進行對比[24]。ROC曲線以敏感度(即實際為滑坡,預(yù)測為滑坡)為縱坐標、1-特異性(即實際為非滑坡,預(yù)測為滑坡)為橫坐標,通過動態(tài)分類閾值避免界限值對結(jié)果的影響[25],如圖10、圖11所示。

圖10 訓練集ROC曲線 圖11 驗證集ROC曲線Fig. 10 ROC curve of train set Fig. 11 ROC curve of validation set

從圖10、圖11中可以看出,4種模型訓練與驗證集AUC值均大于0.7,表征其預(yù)測能力均較好[26],其中集成了3種基分類器的Stacking模型AUC值高于其余單獨分類器模型,說明不同集成模型作為基分類器組合成的Stacking模型泛化能力更好。

文中分別選取滑坡點和非滑坡點對Stacking模型及其基分類器的效果進行了對比,圖12為滑坡點及其周圍模型易發(fā)性結(jié)果對比,Stacking模型預(yù)測的易發(fā)性為極高易發(fā)、高易發(fā);RA、DOCORATE及RS模型預(yù)測出的易發(fā)性為高易發(fā)、中等易發(fā);圖13為非滑坡點及其周圍模型易發(fā)性對比結(jié)果,Stacking模型預(yù)測出的易發(fā)性等級為極低和低;RA、DOCORATE及RS模型預(yù)測出的易發(fā)性等級為低和中等。結(jié)果發(fā)現(xiàn),Stacking模型在易發(fā)性預(yù)測中,與歷史滑坡分布更加吻合,說明其更適用于研究區(qū)滑坡易發(fā)性評價,其易發(fā)性分區(qū)結(jié)果可以作為相關(guān)部門進行土地規(guī)劃與滑坡防控的參考依據(jù)。

圖12 單個滑坡點易發(fā)性等級對比 圖13 單個非滑坡點易發(fā)性等級對比Fig. 12 Comparison of susceptibility grades of single landslide points Fig. 13 Comparison of susceptibility grades of single non landslide points

文中以Stacking模型進行丹鳳縣滑坡易發(fā)性評價,提供了一種新的評價方法與思路,仍存在一些問題將在后續(xù)深入研究:1)文中以集成學習的分支選取RA、DECORATE和RS模型,尚未進行基分類器數(shù)量與種類差異對堆疊效果影響的研究,來提升Stacking模型的泛化能力;2)文中基分類器參數(shù)主要依靠試驗確定,存在一定主觀性,后續(xù)可進行模型參數(shù)的優(yōu)化選取,提高模型預(yù)測精度。

6 結(jié)論

文中采用Stacking集成模型進行了丹鳳縣滑坡易發(fā)性評價,有以下結(jié)論:

1)結(jié)合相關(guān)文獻及地質(zhì)資料,文中選取高程、坡度、坡向、剖面曲率、平面曲率、TWI、NDVI、距水系距離、距道路距離、距斷層距離、地層巖性和年均降雨量共12個影響因子,CAE結(jié)果顯示其均對研究區(qū)滑坡具有關(guān)聯(lián)性,其中高程、TWI和距道路距離與滑坡發(fā)生最為密切。區(qū)內(nèi)滑坡主要集中發(fā)生在高程介于374～720 m、TWI介于8.88～12.55以及距道路距離介于3 108.58～5 364.30 m的區(qū)域中。通過VIF和皮爾遜相關(guān)性進行因子相關(guān)性分析和隨機森林計算各因子權(quán)重,剔除距水系距離、地層巖性因子后選擇剩余10個因子構(gòu)建模型數(shù)據(jù)集。

2)通過WEKA軟件分別構(gòu)建了單一的RA、DECORATE、RS模型以及基于前三者的Stacking集成模型。訓練與驗證結(jié)果表明4種模型AUC值大于0.7,4種模型均具有良好的預(yù)測能力,其中Stacking模型較其他單一模型相比,訓練成功率與驗證預(yù)測率均最高,為國內(nèi)滑坡易發(fā)性評價的模型選擇提供了新的思路。

3)4種模型的滑坡易發(fā)性分區(qū)圖劃分趨勢基本一致,研究區(qū)滑坡高至極高易發(fā)區(qū)主要分布于研究區(qū)中南部,低易發(fā)區(qū)分布于北至東北側(cè)。通過等級分區(qū)統(tǒng)計,Stacking模型的滑坡高至極高易發(fā)區(qū)滑坡頻率比達到2.932,高于3種單獨模型,表明其分區(qū)結(jié)果與歷史滑坡分布最為吻合。研究結(jié)果可作為相關(guān)部門進行滑坡防治與土地利用的參考。