翟文華,王小東,吳明堂,吳曉亮,李倩倩
(1. 華北水利水電大學(xué) 地球科學(xué)與工程學(xué)院,河南 鄭州 450000; 2. 中國電建集團(tuán)華東勘測設(shè)計研究院有限公司,浙江 杭州 310000)
中國是一個地質(zhì)災(zāi)害多發(fā)的國家,依據(jù)自然資源部發(fā)布的數(shù)據(jù),2021年全國發(fā)生地質(zhì)災(zāi)害4772起,造成大量的人員傷亡和經(jīng)濟(jì)損失。而在我國南部地區(qū)由于受雨水浸泡、河流沖刷和人類工程活動的影響,多數(shù)山區(qū)城鎮(zhèn)面臨滑坡、崩塌此類地質(zhì)災(zāi)害的威脅,需進(jìn)行區(qū)域性地質(zhì)災(zāi)害易發(fā)性評價。地質(zhì)災(zāi)害易發(fā)性評價[1]是區(qū)域預(yù)警評價前期重要的準(zhǔn)備工作,主要研究在一定空間下由于環(huán)境因素導(dǎo)致地質(zhì)災(zāi)害發(fā)生概率的大小,其對區(qū)域防災(zāi)減災(zāi)有著重要意義。
目前地質(zhì)災(zāi)害易發(fā)性評價模型眾多,以知識驅(qū)動模型和數(shù)據(jù)驅(qū)動模型為主[2]。知識驅(qū)動模型以專家經(jīng)驗為主,評估過程簡單,容易理解,但其主觀性過強,導(dǎo)致評價結(jié)果無法進(jìn)行對比分析[3],常用的知識驅(qū)動模型有模糊邏輯[4]、層次分析[5-6]和專家打分法[7]。數(shù)據(jù)驅(qū)動以機器學(xué)習(xí)為主,能夠高效地解決地質(zhì)災(zāi)害與環(huán)境因子間復(fù)雜的非線性關(guān)系,常用的數(shù)據(jù)驅(qū)動模型有支持向量機[8-9]、邏輯回歸[10]、隨機森林[11]和神經(jīng)網(wǎng)絡(luò)[12]等。如黃發(fā)明等[13]將聚類分析與支持向量機相結(jié)合運用到重慶市地質(zhì)災(zāi)害多發(fā)區(qū)段的萬州區(qū)進(jìn)行滑坡易發(fā)性評價,結(jié)果表明聚類分析-支持向量機模型比支持向量機模型具有更好的預(yù)測效果。何書等[11]將自組織特征映射網(wǎng)絡(luò)模型和隨機森林模型結(jié)合對江西省贛州市西南部地區(qū)進(jìn)行了滑坡易發(fā)性評價模型,結(jié)果顯示自組織特征映射網(wǎng)絡(luò)與隨機森林相結(jié)合的模型比隨機森林模型具有更高的預(yù)測精度。黃立鑫等[12]將徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)模型與信息量模型進(jìn)行耦合,進(jìn)而對甘肅省岷縣進(jìn)行滑坡災(zāi)害易發(fā)性評價,結(jié)果表明徑向基函數(shù)神經(jīng)網(wǎng)絡(luò)與信息量耦合模型比單一的神經(jīng)網(wǎng)絡(luò)及信息量模型具有更好的評價精度。綜上所述,數(shù)據(jù)驅(qū)動模型在地質(zhì)災(zāi)害易發(fā)性評價領(lǐng)域具有廣泛的應(yīng)用。為了對模型進(jìn)行有效訓(xùn)練和測試,需選取與已知地質(zhì)災(zāi)害單元量相匹配的非地質(zhì)災(zāi)害單元量,而非地質(zhì)災(zāi)害單元的選取目前只是主觀推測或隨機選取,可能導(dǎo)致選取的非災(zāi)害單元位于地質(zhì)災(zāi)害單元附近的潛在地質(zhì)災(zāi)害單元之上[14],并不能保證所選的非地質(zhì)災(zāi)害單元具有較低的易發(fā)性。為減少此類情況的發(fā)生,提高地質(zhì)災(zāi)害單元與非地質(zhì)災(zāi)害單元的辨識度,采用從頻率比模型結(jié)果中選取非地質(zhì)災(zāi)害單元。在此基礎(chǔ)上,構(gòu)建頻率比支持下的隨機森林模型進(jìn)行地質(zhì)災(zāi)害易發(fā)性評價。
近年來,浙江省長興縣李家巷地區(qū)經(jīng)濟(jì)快速發(fā)展,工程活動不斷加劇,造成該地區(qū)地質(zhì)環(huán)境不斷惡化,多地區(qū)發(fā)生滑坡、崩塌,對當(dāng)?shù)鼐用竦纳柏敭a(chǎn)造成了威脅。因此,本文以李家巷為研究區(qū),利用頻率比-隨機森林模型對李家巷進(jìn)行地質(zhì)災(zāi)害易發(fā)性評價,隨后將頻率比-隨機森林模型與隨機森林模型進(jìn)行比較,探討本文提出方法的可行性。
李家巷鎮(zhèn)位于浙江省長興縣,總面積約53.13km2,地勢總體西低東高如圖1所示,研究區(qū)內(nèi)水系眾多,多年平均氣溫為23℃,多年平均降水量達(dá)1701 mm,屬亞熱帶季風(fēng)氣候。地層巖性主要以砂巖、碎屑巖、礫巖和碳酸鹽巖為主。近年來,李家巷人類工程活動加劇,如城市建設(shè)、礦產(chǎn)開采等,都會對環(huán)境造成一定程度的危害,進(jìn)而引發(fā)地質(zhì)災(zāi)害的發(fā)生。以上復(fù)雜的環(huán)境條件,共同孕育了李家巷地區(qū)頻發(fā)的地質(zhì)災(zāi)害。研究區(qū)典型的地質(zhì)災(zāi)害有第二石礦礦區(qū)南側(cè)滑坡如圖2所示。
圖1 研究區(qū)位置Fig. 1 Location of the study area
圖2 研究區(qū)滑坡現(xiàn)場照片F(xiàn)ig. 2 Photo of landslide site in the study area
本研究采用的數(shù)據(jù)源于歷史地質(zhì)災(zāi)害數(shù)據(jù)庫、全球地表覆蓋數(shù)據(jù)庫、地理空間數(shù)據(jù)云平臺、李家巷地質(zhì)災(zāi)害風(fēng)險評價報告,主要包括:LandsatTM8影像1景(2019-05-23,軌道號119/038),用于獲取歸一化植被指數(shù)(normalized difference vegetation index, NDVI); 數(shù)字高程模型(digital elevation model, DEM)數(shù)據(jù)用于提取坡向、坡度和曲率等信息;通過1∶50000地質(zhì)圖提取斷層,巖組等信息;地質(zhì)災(zāi)害編錄數(shù)據(jù)和野外調(diào)查資料等。研究區(qū)地質(zhì)災(zāi)害分布的總面積達(dá)975209 m2,本文采用10 m×10 m的柵格,研究區(qū)共劃分為156051個地質(zhì)單元柵格。
易發(fā)性評價的基礎(chǔ)是對某一環(huán)境指標(biāo)因子進(jìn)行分類并計算各分類狀態(tài)對災(zāi)害的影響程度,常采用頻率比[15]來實現(xiàn)該過程以提高分類的準(zhǔn)確性。頻率比(frequency ratio, FR)可以歸納為某因子分類區(qū)間內(nèi)災(zāi)害柵格數(shù)占所有災(zāi)害柵格百分比與該分類區(qū)間柵格數(shù)占研究區(qū)總柵格數(shù)百分比的比值。計算公式為:
(1)
式中:Fj為某因子在分類區(qū)間內(nèi)出現(xiàn)地質(zhì)災(zāi)害的柵格數(shù);F為區(qū)間內(nèi)所有地質(zhì)災(zāi)害柵格總數(shù);Cj為某因子在分類區(qū)間內(nèi)的柵格數(shù);C為研究區(qū)柵格總數(shù)。
FR表明了評價因子各分類區(qū)間對于地質(zhì)災(zāi)害發(fā)生的影響程度:FR>1表明該評價因子分類區(qū)間對地質(zhì)災(zāi)害的發(fā)生具有較強的影響。FR≤1說明該評價因子分類區(qū)間對災(zāi)害發(fā)生影響不大。
為了彌補單棵決策樹的不足,本文采用隨機森林(random forest, RF)算法。隨機森林是一種新型集成分類器[16]。使用Bootstrap重采樣技術(shù),從樣本集T中有放回地隨機抽取n個樣本(占總樣本的2/3)作為訓(xùn)練樣本集,每條訓(xùn)練樣本數(shù)據(jù)用于訓(xùn)練一棵樹,共生成n個樹組成的森林,預(yù)測數(shù)據(jù)時需根據(jù)n棵樹的結(jié)果,取眾數(shù)或平均值作為結(jié)果,流程如圖3所示。隨機森林參數(shù)的選擇將會直接干擾模型泛化能力和預(yù)測效果,本文利用網(wǎng)格搜索算法確定模型的基本參數(shù)(森林中的樹木數(shù)量、樹的最大深度、特征數(shù)量等)。
圖3 隨機森林算法流程Fig. 3 Flow of random forest algorithm
頻率比-隨機森林(FR-RF)模型主要包含4步:①對評價因子進(jìn)行相關(guān)性分析及利用頻率比計算各評價因子對地質(zhì)災(zāi)害的影響程度;②根據(jù)頻率比模型結(jié)果進(jìn)行非地質(zhì)災(zāi)害單元的選取;③將地質(zhì)災(zāi)害單元、選取的非地質(zhì)災(zāi)害單元和評價因子輸入到的RF模型中,進(jìn)行易發(fā)性評價計算;④對RF模型和FR-RF模型進(jìn)行評估和對比分析,具體流程如圖4所示。
圖4 頻率比-隨機森林模型流程Fig. 4 Frequency ratio-random forest model process
評價因子的選取和評價單元的確定是進(jìn)行災(zāi)害易發(fā)性評價的重要前提條件。地質(zhì)災(zāi)害的發(fā)生不僅由環(huán)境因素引起還與人類活動有關(guān),在結(jié)合相關(guān)文獻(xiàn)[1-13]中相似地區(qū)的地質(zhì)條件、災(zāi)害發(fā)育特征及形成背景的基礎(chǔ)上,從本文1.2節(jié)數(shù)據(jù)源中獲取評價指標(biāo),并進(jìn)行相關(guān)性分析(表1)以確保環(huán)境因子之間具有較弱的相關(guān)性。最終選取9個環(huán)境因子:高程(DEM)、坡度、坡向、巖組、剖面曲率、距離斷層的距離、距離水系的距離、土地利用和歸一化植被指數(shù)(NDVI),如圖5所示。地質(zhì)災(zāi)害易發(fā)性評價的基本評價單元主要包含柵格單元、斜坡單元、行政區(qū)劃單元等[17-18],柵格單元具有結(jié)構(gòu)簡單、方便計算等優(yōu)勢,因此本文以柵格作為基本評估單元。
表1 評價因子間的相關(guān)系數(shù)Table 1 Correlation coefficients among evaluation factors
圖5 評價因子專題圖Fig. 5 Evaluation factors
這些數(shù)據(jù)包含連續(xù)型因子和離散型因子,結(jié)合各因子分布規(guī)律利用等間隔或自然斷點法進(jìn)行劃分,劃分結(jié)果如表2所示。其中巖組中NT為以黏性土為主的巖組,Qg為堅硬塊狀以花崗巖為主的酸性巖巖組,Sc為堅硬較堅硬以砂巖、礫巖為主的粗碎屑巖巖組,Sf為軟硬不均較堅硬中層以粉砂巖、泥巖為主的細(xì)碎屑巖巖組,Tc為堅硬中至厚層狀以碳酸鹽巖為主的巖組,Tcc為較堅硬薄至中層狀以碳酸鹽巖、碎屑巖為主的巖組。本文使用頻率比計算環(huán)境因子的各分類區(qū)間對地質(zhì)災(zāi)害發(fā)生的影響程度。由表2可知,高程在17~302 m范圍內(nèi)、坡度大于10°的坡面、坡向為北、南和西南、NDVI在0.0~0.23和大于0.51、斷層距離小于600 m、水系距離小于700 m、土地利用類型為草地林地以及巖組為Sc,Sf和Tc分別在各所屬因子里面占據(jù)較高的頻率比值,易導(dǎo)致地質(zhì)災(zāi)害的發(fā)生。
表2 各評價因子頻率比值Table 2 Frequency ratio of each evaluation factor
利用FR對環(huán)境因子的各分類區(qū)間進(jìn)行計算,將9個環(huán)境因子的各分類區(qū)間按頻率比值高低進(jìn)行重分類,將重分類后的柵格數(shù)據(jù)進(jìn)行疊加分析。結(jié)合地質(zhì)災(zāi)害易發(fā)性指數(shù)分布規(guī)律,將李家巷鎮(zhèn)劃分為5個易發(fā)性等級區(qū)間:極高[0.9,1)、高[0.7,0.9)、中等[0.3,0.7)、低[0.1,0.3)和極低[0,0.1),5個易發(fā)區(qū)由高至低分別占總面積的2.36%、15.85%、40.41%、41.21%、0.16%,FR模型計算得到的李家巷地質(zhì)災(zāi)害易發(fā)性結(jié)果如圖6所示。由表3可知,極高和高易發(fā)區(qū)約占地質(zhì)災(zāi)害總面積的93.34%,而極低易發(fā)區(qū)約占地質(zhì)災(zāi)害總面積的0%,在研究區(qū)的中部地區(qū)地形平坦災(zāi)害發(fā)生概率幾乎很小,但在該模型分析結(jié)果中,有大量的中易發(fā)區(qū)在中部地區(qū),可見FR模型預(yù)測出的易發(fā)性結(jié)果具有偏差。
表3 基于頻率比模型地質(zhì)災(zāi)害易發(fā)性等級的頻率比Table 3 Frequency ratio of geological hazard susceptibility grades based on frequency ratio model
圖6 疊加分析計算得到的易發(fā)性分布及非地質(zhì)災(zāi)害單元選取Fig. 6 Distribution of susceptibility calculated by superposition analysis and the selection of non-geological hazard units圖7 隨機森林計算得到的易發(fā)性分布及非地質(zhì)災(zāi)害單元選取Fig. 7 Susceptibility distribution by random forest calculation and the selection of non-geological hazard units
由于地質(zhì)災(zāi)害往往發(fā)生于暫未發(fā)生地質(zhì)災(zāi)害的區(qū)域,為避免隨機選取的非地質(zhì)災(zāi)害單元落到潛在地質(zhì)災(zāi)害上,在FR模型結(jié)果中的低和極低易發(fā)區(qū)隨機選取與已知地質(zhì)災(zāi)害單元樣本量匹配的非地質(zhì)災(zāi)害單元,非地質(zhì)災(zāi)害單元的選取結(jié)果如圖6所示。隨機選取的非地質(zhì)災(zāi)害單元主要分布在平原區(qū),部分點分布在坡度較緩的山區(qū),表明經(jīng)頻率比模型篩選非地質(zhì)災(zāi)害單元的空間分布合理。
在研究區(qū)隨機選取與地質(zhì)災(zāi)害單元樣本量匹配的非地質(zhì)災(zāi)害單元如圖7所示,并將地質(zhì)災(zāi)害單元和非地質(zhì)災(zāi)害單元的集合分為兩部分:70%用于模型訓(xùn)練,30%用于模型測試。將地質(zhì)災(zāi)害單元設(shè)為1,非地質(zhì)災(zāi)害單元設(shè)為0。并將所有環(huán)境因子都進(jìn)行特征標(biāo)準(zhǔn)化處理作為模型的輸入變量。采用網(wǎng)格搜索算法和交叉驗證對RF模型進(jìn)行參數(shù)優(yōu)化和檢驗。計算的易發(fā)性結(jié)果采用同樣的等級分區(qū),5個易發(fā)區(qū)由高至低分別占總面積的1.53%、1.99%、8.37%、10.66%、77.44%,RF模型計算得到的李家巷地質(zhì)災(zāi)害易發(fā)性如圖7所示。由表4可知,極高和高易發(fā)區(qū)約占地質(zhì)災(zāi)害總面積的55.92%,而極低易發(fā)區(qū)約占地質(zhì)災(zāi)害總面積的7.68%,可見RF模型預(yù)測的地質(zhì)災(zāi)害易發(fā)性結(jié)果也存在一定的偏差。
表4 基于隨機森林模型地質(zhì)災(zāi)害易發(fā)性等級的頻率比Table 4 Frequency ratio of geological hazard susceptibility grades based on random forest model
利用4.1節(jié)中選取的非地質(zhì)災(zāi)害單元,并結(jié)合地質(zhì)災(zāi)害單元作為RF模型的輸入變量,采用網(wǎng)格搜索算法和交叉驗證對RF模型進(jìn)行參數(shù)優(yōu)化和檢驗。為確保不同模型的易發(fā)性評價結(jié)果具有可對比性,FR-RF模型的易發(fā)性結(jié)果同樣將研究區(qū)劃分為5個等級,5個易發(fā)區(qū)由高至低分別占總面積的9.43%、4.40%、9.60%、4.03%、72.53%。FR-RF模型計算得到的李家巷地質(zhì)災(zāi)害易發(fā)性如圖8所示。由表5中可知,極高和高易發(fā)區(qū)約占地質(zhì)災(zāi)害總面積的81.90%,而極低易發(fā)區(qū)約占地質(zhì)災(zāi)害總面積的0%,可見FR-RF模型具有較高的準(zhǔn)確度。
圖8 頻率比-隨機森林計算得到的易發(fā)性分布Fig. 8 Susceptibility distribution by frequency radio-random forest calculation圖9 評價指標(biāo)因子的貢獻(xiàn)程度雷達(dá)圖Fig. 9 Contribution degree radar diagram of evaluation index factors
表5 基于頻率比-隨機森林模型地質(zhì)災(zāi)害易發(fā)性等級的頻率比Table 5 Frequency ratio of geological hazard susceptibility grade based on frequency ratio-random forest model
RF模型和FR-RF模型兩者分析結(jié)果都表明了在研究區(qū)的東北部和東南部地區(qū)是地質(zhì)災(zāi)害的主要發(fā)生地區(qū),中部、北部以及南部地區(qū)幾乎沒有地質(zhì)災(zāi)害發(fā)生。FR-RF模型在極高和高易發(fā)區(qū)所占地質(zhì)災(zāi)害總面積數(shù)比RF模型高25.98%,在極低和低易發(fā)區(qū)所占地質(zhì)災(zāi)害總面積數(shù)比RF模型低16.7%。由圖8可知,FR-RF模型的極高易發(fā)區(qū)在東北部以及南部邊界地帶以及中部偏東地區(qū)明顯高于RF模型,且根據(jù)實際調(diào)查結(jié)果,表明FR-RF模型與實際災(zāi)害分布更為吻合,劃分效果更好。利用FR-RF模型,得到研究區(qū)地質(zhì)災(zāi)害易發(fā)性評價因子的貢獻(xiàn)程度雷達(dá)圖如圖9所示,結(jié)果表明,坡度和工程巖組是該研究區(qū)地質(zhì)災(zāi)害的主控因素,相對重要性占比分別為37.7%和28.0%。
本文通過接受者操作特性曲線(receiver operating characteristic curve, ROC)下面積(area under curve,AUC)、準(zhǔn)確率(Accuracy)、特異性(Specificity)和敏感性(Recall)對模型性能進(jìn)行評估,其中AUC通常大于0.5小于1,AUC值越大,模型性能越好,采用測試樣本繪制ROC曲線,如圖10所示。Accuracy是正確預(yù)測地質(zhì)災(zāi)害單元和非地質(zhì)災(zāi)害單元的樣本占所有樣本比例,值介于0和1之間;Recall為能正確預(yù)測地質(zhì)災(zāi)害單元的能力;Specificity能正確預(yù)測非地質(zhì)災(zāi)害單元的能力。RF和FR-RF模型在測試集上的評估結(jié)果如表6所示,FR-RF模型具有較高的AUC、Accuracy、Specificity和Recall,FR-RF模型比RF模型的AUC、Accuracy、Specificity和Recall分別高出6%、11.66%、9.51%和13.71%。FR-RF模型的Recall和Specificity分別為98.33%和93.33%,這表明FR-RF模型對地質(zhì)災(zāi)害單元和非地質(zhì)災(zāi)害單元的預(yù)測都具有較好的性能。由此可見,FR-RF具有比隨機選取災(zāi)害單元的RF模型更高的預(yù)測率,能更加有效地反映李家巷地質(zhì)災(zāi)害的分布特征。
圖10 模型計算的災(zāi)害易發(fā)性預(yù)測率曲線Fig. 10 Disaster susceptibility prediction rate curves calculated by the model表6 隨機森林和頻率比-隨機森林模型性能對比Table 6 Performance comparison of random forest and frequency ratio-random forest models %評估指標(biāo)RFFR-RFAUC90.0096.00Accuracy84.1795.83Specificity83.8293.33Recall84.6298.33
本文以浙江省長興縣李家巷鎮(zhèn)作為研究區(qū),分別采取FR模型、RF模型及FR-RF模型,對研究區(qū)進(jìn)行地質(zhì)災(zāi)害易發(fā)性評價,得出如下結(jié)論:
1) 本文采用高程、坡度、巖組、剖面曲率、坡向、歸一化植被指數(shù)、土地利用、距水系距離和距斷層距離9個評價因子進(jìn)行地質(zhì)災(zāi)害易發(fā)性評價,經(jīng)模型重要性分析表明,坡度、工程巖組是該研究區(qū)地質(zhì)災(zāi)害的主控因素,相對重要性占比分別為37.7%和28.0%。
2) 根據(jù)浙江省長興縣李家巷鎮(zhèn)地質(zhì)災(zāi)害發(fā)育與分布特征,本文利用FR與RF耦合的方法,確定的非地質(zhì)災(zāi)害單元的選取,建立了FR-RF模型,并對研究區(qū)進(jìn)行了地質(zhì)災(zāi)害易發(fā)性評價。FR-RF模型比單一的RF模型AUC、Accuracy、Specificity和Recall分別高出6%、11.66%、9.51%和13.71%。表明FR-RF模型不僅能夠有效選取非地質(zhì)災(zāi)害單元,而且對地質(zhì)災(zāi)害單元和非地質(zhì)災(zāi)害單元預(yù)測都具有較好的性能,并且在地質(zhì)災(zāi)害易發(fā)性預(yù)測方面具有較高的精度。
3) 基于FR-RF模型對研究區(qū)進(jìn)行地質(zhì)災(zāi)害易發(fā)性評價結(jié)果與實際地質(zhì)災(zāi)害分布相對應(yīng),極高易發(fā)區(qū)和高易發(fā)區(qū)主要分布于存在碎屑巖和坡度較大的區(qū)域。隨著易發(fā)性級別的增大,各等級內(nèi)地質(zhì)災(zāi)害柵格數(shù)量隨之增大,說明采用FR-RF模型確定的易發(fā)性分區(qū)與實際地質(zhì)災(zāi)害發(fā)生情況相吻合。