黃艷琴, 李為樂, 許 洲, 李鵬飛, 鐵永波
(1.成都理工大學(xué)地質(zhì)災(zāi)害防治與地質(zhì)環(huán)境保護(hù)國家重點(diǎn)實(shí)驗(yàn)室,四川 成都 610059;2.四川峨眉山四零三建設(shè)工程有限責(zé)任公司,四川 樂山 614200; 3.中國電建集團(tuán)貴陽勘測(cè)設(shè)計(jì)研究院有限公司,貴州 貴陽 550081;4.中國地質(zhì)調(diào)查局成都地質(zhì)調(diào)查中心,四川 成都 610081)
目前,我國公路交通建設(shè)正處于快速發(fā)展階段,山區(qū)公路建設(shè)的投入也在不斷加大[1]。山區(qū)地質(zhì)環(huán)境復(fù)雜、山高坡陡、溝谷縱橫,在地震或強(qiáng)降雨作用下滑坡災(zāi)害頻發(fā),嚴(yán)重威脅著山區(qū)公路的建設(shè)和運(yùn)營安全。近年來,隨著山區(qū)大量高速公路的建設(shè)和使用,每年都會(huì)出現(xiàn)公路滑坡災(zāi)害發(fā)生的情況,造成公路損毀、交通癱瘓甚至人員傷亡。例如: 2014年7月17日,國道213線K774+600 m處突發(fā)山體滑坡,塌方量達(dá)3 000 m3,造成國道213線約100 m的道路被掩埋阻斷[2]; 2018年7月26日,國道213線K852 +600 m處發(fā)生山體滑坡,造成交通臨時(shí)中斷[3]。
山區(qū)公路滑坡隱患一般具有高位、高隱蔽性、突發(fā)性強(qiáng)等特點(diǎn)[4],傳統(tǒng)地面調(diào)查手段無法進(jìn)行大范圍提前識(shí)別。利用光學(xué)遙感技術(shù)可大范圍識(shí)別出形態(tài)特征和變形跡象顯著的滑坡[5]。而采用合成孔徑雷達(dá)干涉測(cè)量(interferometric synthetic aperture radar,InSAR)技術(shù)可獲取滑坡地表形變,對(duì)目前正在變形的活動(dòng)滑坡進(jìn)行識(shí)別。其中,差分干涉圖疊加(Stacking-InSAR)技術(shù)相比于其他InSAR技術(shù),能夠有效抑制大氣效應(yīng)并減少DEM誤差的影響,快速獲取地表形變速率,被廣泛應(yīng)用于滑坡隱患的早期識(shí)別[5]。光學(xué)遙感和InSAR技術(shù)的結(jié)合可提高滑坡隱患識(shí)別的準(zhǔn)確性。
滑坡易發(fā)性是對(duì)滑坡空間發(fā)生概率的定量評(píng)價(jià)[6]。隨著遙感技術(shù)的發(fā)展,國內(nèi)外學(xué)者依托遙感技術(shù)對(duì)滑坡易發(fā)性評(píng)價(jià)開展了大量研究工作,評(píng)價(jià)方法也由定性評(píng)價(jià)轉(zhuǎn)為定量評(píng)價(jià)。常見滑坡易發(fā)性定量評(píng)價(jià)模型有Logistic回歸模型、隨機(jī)森林模型、支持向量機(jī)模型、信息量模型、神經(jīng)網(wǎng)絡(luò)模型等[7]。
國道213汶川—松潘段位于強(qiáng)震山區(qū),構(gòu)造活動(dòng)強(qiáng)烈,滑坡災(zāi)害頻發(fā),造成公路被堵斷或損毀,影響公路交通的正常運(yùn)營。本文以國道213汶川—松潘段為研究區(qū),結(jié)合光學(xué)遙感技術(shù)和Stacking-InSAR技術(shù)識(shí)別研究區(qū)滑坡隱患,并采用Logistic回歸模型進(jìn)行易發(fā)性評(píng)價(jià),為提前預(yù)防研究區(qū)山體滑坡,保障國道213線安全運(yùn)行,并為我國公路規(guī)劃建設(shè)、防災(zāi)減災(zāi)提供理論參考和技術(shù)支撐。
國道213汶川—松潘段位于四川強(qiáng)震山區(qū),自四川省汶川縣映秀鎮(zhèn)開始,途徑茂縣,向北延伸至松潘縣川主寺鎮(zhèn),總長約230 km(圖1)。研究區(qū)整體地勢(shì)北西高、南東低,海拔在800~4 700 m之間,相對(duì)高差500~1 200 m,屬于中-高山地貌。岷江是研究區(qū)內(nèi)的主要水系,由南向北貫穿全境,受內(nèi)部水動(dòng)力和外部地質(zhì)作用的影響,岷江河床不斷向下深切,形成了如今的高山峽谷地貌。受高海拔西北風(fēng)氣流和印度洋西南季風(fēng)的影響,研究區(qū)呈現(xiàn)青藏高原季風(fēng)氣候特點(diǎn),冬季干冷,夏季濕潤涼爽,日照量大,氣溫日差大、年差小。研究區(qū)內(nèi)地質(zhì)構(gòu)造活動(dòng)強(qiáng)烈,主要發(fā)育的斷層有龍門山斷裂帶、岷江斷裂帶、虎牙斷裂帶、松坪溝斷層、汶川—茂汶斷層、較弧場(chǎng)構(gòu)造等[8],2008年汶川Ms 8.0級(jí)地震和2017年九寨溝Ms 7.0級(jí)地震均造成研究區(qū)公路不同程度的損毀。區(qū)內(nèi)出露地層有新元古界震旦系,古生界寒武系、奧陶系、志留系、泥盆系、石炭系、二疊系,中生界三疊系以及新生界第四系,巖體質(zhì)量較差、結(jié)構(gòu)破碎。
圖1 研究區(qū)位置和影像覆蓋范圍
利用2022年2月獲取的高分2號(hào)衛(wèi)星影像和Google Earth在線三維衛(wèi)星影像數(shù)據(jù),采用人工目視解譯方法。光學(xué)遙感識(shí)別的滑坡隱患類型主要分為2類: 已發(fā)生整體失穩(wěn)的老(古)滑坡隱患和正在變形的潛在滑坡隱患[9]。不同滑坡隱患類型的識(shí)別標(biāo)志不同: 對(duì)于已經(jīng)發(fā)生過整體失穩(wěn)的古(老)滑坡,主要基于形態(tài)特征進(jìn)行識(shí)別,整體平面形態(tài)呈圈椅狀或舌形,后緣可見滑坡壁,前緣可見滑坡舌擠壓河流,滑坡體上有群居居民或耕地[10](圖2); 對(duì)于正在變形的潛在滑坡隱患,主要識(shí)別變形跡象特征明顯的區(qū)域,識(shí)別標(biāo)志是坡體裂縫的形成過程以及前緣是否出現(xiàn)局部垮塌現(xiàn)象,針對(duì)這一類型滑坡隱患的識(shí)別,相對(duì)快速有效的方法是充分收集多時(shí)相光學(xué)遙感影像,通過對(duì)比不同時(shí)期的遙感影像特征解譯滑坡隱患,從而提高遙感識(shí)別滑坡隱患的準(zhǔn)確性(圖3)。
圖2 典型老(古)滑坡隱患衛(wèi)星影像(左)與Google Earth衛(wèi)星影像(右)
(a) 2011年 (b) 2016年 (c) 2019年 (d) 2021年
InSAR技術(shù)是一種監(jiān)測(cè)地表緩慢形變的遙感技術(shù)[11],具有全天候、全天時(shí)、精度高、周期性重復(fù)觀測(cè)等優(yōu)點(diǎn)。本文選取2017年1月至2021年12月覆蓋研究區(qū)域的288景Sentinel-1A降軌數(shù)據(jù),采用Stacking-InSAR技術(shù)對(duì)研究區(qū)活動(dòng)滑坡進(jìn)行探測(cè)。
2.2.1 Stacking-InSAR技術(shù)原理
Stacking-InSAR技術(shù)由Price等[12]提出,其基本理論為: 假設(shè)地形相位誤差、大氣延遲相位誤差等具有隨機(jī)性,聯(lián)合多期解纏后的差分干涉相位圖建立形變速率和干涉相位間的線性函數(shù)模型,從而監(jiān)測(cè)時(shí)間段的線性形變速率。利用該技術(shù)可以削弱InSAR技術(shù)中存在的隨機(jī)軌道誤差、大氣延遲相位誤差和地形誤差等的影響,提高形變計(jì)算精度。Stacking-InSAR技術(shù)的計(jì)算方法可表示為[13]
。
(1)
式中:V為雷達(dá)衛(wèi)星視線方向平均形變速率,mm/a;λ為波長,m;φi為第i個(gè)解纏后的差分干涉相位;M為干涉圖數(shù)量;ti為第i個(gè)差分干涉相位的時(shí)間間隔,d。
相應(yīng)的誤差傳播公式為
。
(2)
式中:ΔVdisp為雷達(dá)衛(wèi)星視線方向形變速率誤差值,mm/a;E為假定的單幅干涉圖相位誤差;M為干涉圖數(shù)量;λ為波長,m;ti為第i個(gè)差分干涉相位的時(shí)間間隔,d。
2.2.2 Stacking-InSAR技術(shù)處理流程
通過Stacking-InSAR技術(shù)獲取研究區(qū)地表形變速率圖,可以大面積探測(cè)活動(dòng)滑坡隱患。運(yùn)用該技術(shù)首先對(duì)SAR影像數(shù)據(jù)進(jìn)行處理,之后選取主從影像并配準(zhǔn)生成干涉圖。結(jié)合外部DEM與衛(wèi)星成像參數(shù),將DEM和SAR影像精確配準(zhǔn),并將DEM數(shù)據(jù)編碼到SAR影像坐標(biāo)系,生成差分干涉圖。對(duì)差分干涉圖進(jìn)行空間濾波和相位解纏處理,解纏相位加權(quán)疊加后生成地表形變速率圖。通過地表形變速率圖篩選形變區(qū)域,從而確定滑坡隱患位置并識(shí)別滑坡隱患,數(shù)據(jù)處理流程見圖4。
結(jié)合現(xiàn)場(chǎng)調(diào)查與無人機(jī)航攝技術(shù)進(jìn)行滑坡野外調(diào)查驗(yàn)證,主要工作包括判斷室內(nèi)解譯的滑坡隱患是否為滑坡、滑坡邊界的圈定是否正確并判定滑坡類型; 查明坡體裂縫的發(fā)育情況,如位置、方向、深度、寬度等(圖5)。研究區(qū)共驗(yàn)證滑坡點(diǎn)185處,占滑坡總數(shù)的59.7%。其中有163處屬于滑坡,排除滑坡22處,野外驗(yàn)證正確率為88.1%。
(a) 滑坡隱患邊界范圍 (b) 坡體變形及裂縫發(fā)育調(diào)查
利用光學(xué)衛(wèi)星影像共識(shí)別滑坡隱患308處,基于InSAR技術(shù)共探測(cè)出活動(dòng)滑坡隱患27處?;顒?dòng)滑坡主要集中分布于汶川—茂縣段,其中25處為光學(xué)衛(wèi)星遙感與InSAR技術(shù)共同識(shí)別,另外2處(B4、B16)為InSAR技術(shù)單獨(dú)識(shí)別(圖6、圖7)。野外驗(yàn)證排除滑坡隱患22處,最終研究區(qū)共識(shí)別滑坡隱患288處(圖8)。其中,汶川—茂縣段滑坡隱患發(fā)育較多,分布較密集,該區(qū)域內(nèi)地質(zhì)構(gòu)造活動(dòng)強(qiáng)烈,地層巖性差異大,巖體結(jié)構(gòu)較破碎,受公路或河流切割明顯,臨空條件較好。在地震或惡劣天氣條件誘發(fā)作用下,發(fā)生滑坡的可能性增大。
圖6 基于Stacking-InSAR技術(shù)的活動(dòng)滑坡隱患形變探測(cè)
圖7 Stacking-InSAR技術(shù)探測(cè)滑坡隱患放大效果圖
圖8 研究區(qū)滑坡隱患分布
通常易發(fā)性評(píng)價(jià)因子體系主要由內(nèi)在因素和外在因素兩部分構(gòu)成[14]。內(nèi)在因素一般側(cè)重于斜坡本身的巖土體性質(zhì)或地形地貌條件,其中地形地貌包括斜坡坡度、坡向和坡高等,巖土體性質(zhì)包括地層巖性、斜坡及其周圍環(huán)境的地質(zhì)構(gòu)造等。外在因素則被定義為除斜坡本身成災(zāi)環(huán)境條件之外誘發(fā)滑坡發(fā)生的因素,即外部環(huán)境條件對(duì)滑坡的影響,一種是自然環(huán)境的誘發(fā)因素,如地震或降雨; 另一種是人類工程活動(dòng)的誘發(fā)因素,如交通建設(shè)、城鎮(zhèn)建設(shè)、水利水電工程建設(shè)等。
根據(jù)上述對(duì)易發(fā)性評(píng)價(jià)因子體系的分析以及對(duì)研究區(qū)滑坡特點(diǎn)、分布規(guī)律、影響因素等的研究,選取高程、坡度、坡向、地表曲率、工程地質(zhì)巖組、歸一化植被指數(shù)(normalized difference vegetation index,NDVI)、距道路距離、距河流距離以及距斷層距離共計(jì)9個(gè)評(píng)價(jià)因子,對(duì)研究區(qū)進(jìn)行了滑坡易發(fā)性評(píng)價(jià)。
評(píng)價(jià)因子數(shù)據(jù)來源主要分為4類: DEM數(shù)據(jù)來源于30 m GDEMV2數(shù)字高程數(shù)據(jù)(地理空間數(shù)據(jù)云),用于提取高程、坡度、坡向、地表曲率因子; 地質(zhì)數(shù)據(jù)來源于1∶20萬地質(zhì)圖,用于獲取巖性與斷層分布; 水文數(shù)據(jù)來源于1∶20萬地質(zhì)圖,主要提取河流與道路矢量; NDVI數(shù)據(jù)來源于2020年12月2日的全球植被覆蓋Landsat 8衛(wèi)星影像數(shù)據(jù)。
研究區(qū)采用的評(píng)價(jià)單元為柵格單元,柵格單元是滑坡易發(fā)性評(píng)價(jià)中應(yīng)用最廣泛的評(píng)價(jià)單元,把研究區(qū)按照一定尺寸劃分為規(guī)則的網(wǎng)格,具有計(jì)算機(jī)處理方便、可操作性強(qiáng)等優(yōu)點(diǎn)。柵格單元?jiǎng)澐值拇笮∮啥喾N因素決定,常用專家經(jīng)驗(yàn)公式求取[15]。同時(shí),需綜合考慮研究區(qū)地貌特征以及歸一化植被指數(shù)柵格數(shù)據(jù)的空間分辨率,因此研究區(qū)采用的柵格單元大小為30 m×30 m,共劃分出2 759 358個(gè)柵格單元。公式為
Gs=7.49+6×10-4×S-2.0×
10-6+2.9×10-15×S2。
(3)
式中:Gs為柵格單元建議大小;S為高程數(shù)據(jù)的比例尺倒數(shù),高程數(shù)據(jù)比例尺為1∶50 000。
分析各評(píng)價(jià)因子對(duì)滑坡發(fā)育的影響,統(tǒng)計(jì)各評(píng)價(jià)因子類別中滑坡發(fā)育的數(shù)量(圖9),從而對(duì)各評(píng)價(jià)因子進(jìn)行分級(jí)。
(1)高程。高程反映區(qū)域的宏觀地貌特征,表現(xiàn)為巖土體特征和地下水性質(zhì)的差異,從而間接影響滑坡的發(fā)育條件。研究區(qū)處于中-高山區(qū),海拔在828~4 694 m之間,將高程按300 m等間距分為10個(gè)區(qū)間,統(tǒng)計(jì)分析各區(qū)間內(nèi)滑坡的分布。結(jié)果表明,滑坡在[1 300,1 600) m區(qū)間分布比較集中,該區(qū)間內(nèi)滑坡較發(fā)育,易發(fā)性貢獻(xiàn)率大,而在高程<1 300 m和≥3 400 m區(qū)間滑坡發(fā)育相對(duì)最弱,易發(fā)性貢獻(xiàn)率小。根據(jù)滑坡在研究區(qū)高程的發(fā)育程度,將研究區(qū)高程分為: <1 300 m、[1 300,1 600) m、[1 600,2 500) m、[2 500,2 800) m、[2 800,3 400) m以及≥3 400 m共6個(gè)區(qū)間。
(2)坡度。斜坡的應(yīng)力分布受坡度的影響,坡度大小不同,則斜坡的應(yīng)力與穩(wěn)定性不同。研究區(qū)坡度在0°~87°之間,為了更好地分析滑坡在不同坡度范圍內(nèi)的發(fā)育數(shù)量和分布情況,將研究區(qū)的坡度因子以等間距5°劃分為11個(gè)坡度區(qū)間。分析結(jié)果顯示,滑坡在[35°,45°)區(qū)間內(nèi)分布比較集中,該坡度區(qū)間滑坡較發(fā)育,易發(fā)性貢獻(xiàn)率大,而坡度在<20°和≥50°區(qū)間內(nèi)滑坡分布較稀疏,滑坡不發(fā)育,易發(fā)性貢獻(xiàn)率小。根據(jù)滑坡在研究區(qū)坡度的分布情況,將研究區(qū)坡度分為<20°、[20°,30°)、[30°,40°)、[40°,50°)以及≥50°共5個(gè)區(qū)間。
(a) 高程-滑坡數(shù)量 (b) 坡度-滑坡數(shù)量 (c) 坡向-滑坡數(shù)量
(d) 地表曲率-滑坡數(shù)量 (e) 工程地質(zhì)巖組-滑坡數(shù)量 (f) NDVI-滑坡數(shù)量
(g) 距道路距離-滑坡數(shù)量 (h) 距河流距離-滑坡數(shù)量 (i) 距斷層距離-滑坡數(shù)量
(3)坡向。坡向是日照時(shí)數(shù)和太陽照射強(qiáng)弱的客觀反映[16]。研究區(qū)坡向在0°~360°之間,將其按方位劃分為8個(gè)方位向,并統(tǒng)計(jì)滑坡在各方位向的數(shù)量和分布情況。分析統(tǒng)計(jì)結(jié)果可知,研究區(qū)東朝向[67.5°,112.5°)、西南向[202.5°,247.5°)以及南朝向[247.5°,292.5°)滑坡分布較為密集,易發(fā)性貢獻(xiàn)率較大,而坡向[0,22.5°)與[337°,360°)區(qū)間的滑坡分布少,易發(fā)性貢獻(xiàn)率相對(duì)最低。根據(jù)滑坡在各朝向的分布情況,將研究區(qū)坡向劃分為[0°,67.5°)、[67.5°,112.5°)、[112.5°,202.5°)、[202.5°,247.5°)、[247.5°,337.5°)以及≥337.5°共6個(gè)區(qū)間。
(4)地表曲率。地表曲率的大小指示了斜坡的形態(tài)特征。在不同斜坡形態(tài)條件下,滑坡的發(fā)生具有一定的差異性。根據(jù)滑坡發(fā)育形態(tài)特征進(jìn)行分類: 地表曲率<-0.5為凹型坡,地表曲率[-0.5,0.5)的為直線型坡,地表曲率≥0.5為凸型坡。
(5)工程地質(zhì)巖組。斜坡穩(wěn)定性受自身物質(zhì)組成的影響,物質(zhì)組成不同,巖土體的物理參數(shù)不同,則巖土體的物理性質(zhì)不同。研究區(qū)巖組的分類標(biāo)準(zhǔn)參考《GB 50218—2014工程巖體分級(jí)標(biāo)準(zhǔn)》[17],共分為4個(gè)等級(jí),分別為松散巖體、較軟弱巖體、較堅(jiān)硬巖體和堅(jiān)硬巖體。統(tǒng)計(jì)分析各類巖組中的滑坡發(fā)育數(shù)量,可知松散巖體中滑坡數(shù)量最少,相對(duì)最不發(fā)育; 較堅(jiān)硬巖體和堅(jiān)硬巖體中滑坡發(fā)育最多,可能與巖體遭受風(fēng)化剝蝕有關(guān),風(fēng)化越嚴(yán)重,巖體破碎程度越明顯,斜坡越不穩(wěn)定,發(fā)生滑坡災(zāi)害的可能性越大。
(6)歸一化植被指數(shù)。植被對(duì)滑坡發(fā)生的影響主要表現(xiàn)為植被的根莖深入地底,可以起到根固的作用,還可以緩沖坡面水流流速,從而極大程度地減少山體水流對(duì)坡面的破壞,減少滑坡災(zāi)害的發(fā)生。NDVI能反映某一區(qū)域內(nèi)植被的覆蓋程度,利用Landsat8遙感影像求取NDVI,計(jì)算公式為[18]
NDVI=(IR-R)/(IR+R)。
(4)
式中:IR為近紅外波段;R為紅外波段。NDVI的取值區(qū)間為[-1,1],其值越大表示植被覆蓋程度越高,生長狀況良好。結(jié)合研究植被分布規(guī)律,本研究采用自然間斷法,將NDVI劃分為<0.17、[0.17,0.36)、[0.36,0.58)、[0.58,0.82)以及[0.82,1]共5個(gè)區(qū)間。
(7)距道路距離。在道路建設(shè)過程中,坡腳開挖、坡頂荷載等人類工程活動(dòng)會(huì)不同程度地破壞巖土體完整性和天然狀態(tài)下的自然結(jié)構(gòu)。按照100 m間距劃分出11個(gè)區(qū)間,得出距道路距離與滑坡數(shù)量的關(guān)系。在距道路100 m范圍內(nèi)滑坡分布最多,道路對(duì)滑坡發(fā)生的影響最大。根據(jù)滑坡在不同距道路距離的分布情況,將研究區(qū)劃分為<100 m、[100,200) m、[200,300) m、[300,500) m以及≥500 m共5個(gè)區(qū)間。
(8)距河流距離。河流掏蝕作用會(huì)使巖體產(chǎn)生裂縫,越靠近河流的巖土體含水率越大,抗剪強(qiáng)度降低,斜坡穩(wěn)定性隨之降低。以100 m為間隔將距河流距離劃分為11個(gè)區(qū)間,統(tǒng)計(jì)各區(qū)間的滑坡發(fā)育情況。分析表明,滑坡主要集中分布在距河流距離300 m的范圍內(nèi),而距離河流大于800 m的區(qū)域滑坡總體發(fā)育較少。通過分析滑坡發(fā)育數(shù)量及分布與河流距離的關(guān)系,將研究區(qū)距河流距離劃分為<100 m、[100,300) m、[300,500) m、[500,800) m以及≥800 m共5個(gè)區(qū)間。
(9)距斷層距離。區(qū)域構(gòu)造活動(dòng)的強(qiáng)烈程度會(huì)影響節(jié)理裂隙發(fā)育和巖土體結(jié)構(gòu),可直接或間接影響滑坡的形成和破壞。以1 km為緩沖間隔將研究區(qū)的距斷層距離劃分為11個(gè)區(qū)間,統(tǒng)計(jì)分析各緩沖區(qū)域內(nèi)滑坡的發(fā)育數(shù)量。結(jié)果表明,在[0,1) km區(qū)間內(nèi)滑坡相對(duì)最發(fā)育,而在[2,8) km區(qū)間滑坡發(fā)育數(shù)量呈現(xiàn)逐步降低趨勢(shì),當(dāng)距斷層距離≥10 km時(shí),滑坡發(fā)育數(shù)量呈增長模式,由于研究區(qū)區(qū)域的形狀較長,斷層分布跨度大,在此區(qū)間內(nèi)斷層相對(duì)集中,滑坡發(fā)育數(shù)量增多。統(tǒng)計(jì)不同斷層緩沖區(qū)內(nèi)滑坡的數(shù)量,將研究區(qū)斷層緩沖區(qū)劃分為<1 km、[1,2) km、[2,4) km、[4,10) km以及≥10 km共5個(gè)區(qū)間。
3.5.1 Logistic回歸模型
Logistic回歸模型是一種預(yù)測(cè)性建模技術(shù),其基本原理為利用回歸分析方法研究因變量(目標(biāo))與自變量(預(yù)測(cè)器)之間的關(guān)系,由一個(gè)或多個(gè)滑坡影響因素(自變量)來判定一個(gè)滑坡(自變量)是否發(fā)生(0或1)的概率[19-20]。Logistic回歸函數(shù)表達(dá)式為
Y=Logit(P)=ln[P/(1-P)]=
B0+B1X1+B2X2+…BnXn,
(5)
P=1/(1-exp(-Y))。
(6)
式中:Y為二元值(0或1),表示滑坡發(fā)生與否,發(fā)生Y=1,不發(fā)生Y=0;P為滑坡發(fā)生概率;X1,X2,…,Xn為各評(píng)價(jià)因子;B0為常量,取值-10.756;B1,B2,…,Bn為各評(píng)價(jià)因子回歸系數(shù)。
本文選取Ayalew[21]提出的評(píng)價(jià)因子量綱歸一化方法獲取歸一化值。計(jì)算公式為
,
(7)
。
(8)
表1 評(píng)價(jià)因子分級(jí)和因子歸一化值
(續(xù)表)
3.5.2 Logistic回歸模型的構(gòu)建
建立邏輯回歸模型的步驟如下:
(1)將研究區(qū)矢量滑坡分布圖轉(zhuǎn)為柵格圖層,發(fā)生滑坡的區(qū)域賦值為1,未發(fā)生滑坡的區(qū)域賦值為0。
(2)根據(jù)所求的歸一化值,分別將9個(gè)評(píng)價(jià)因子圖柵格化。
(3)在研究區(qū)域范圍內(nèi),隨機(jī)選取230個(gè)滑坡點(diǎn)(占總滑坡點(diǎn)的80%),將剩余的58個(gè)滑坡點(diǎn)(占總滑坡點(diǎn)的20%)作為本次評(píng)價(jià)的檢驗(yàn)點(diǎn)。同時(shí),借助ArcGIS軟件中的隨機(jī)點(diǎn)生成工具,在研究區(qū)范圍內(nèi),排除滑坡點(diǎn)以相同的點(diǎn)間距500 m生成隨機(jī)點(diǎn),創(chuàng)建與滑坡點(diǎn)數(shù)量相同的230個(gè)非滑坡點(diǎn),共460個(gè)點(diǎn)作為本次評(píng)價(jià)的總樣本點(diǎn)。
(4)將步驟(1)和步驟(2)中的值賦在460個(gè)樣本點(diǎn)上,得到樣本點(diǎn)屬性表,將其導(dǎo)入SPSS軟件中,利用回歸分析工具進(jìn)行總樣本點(diǎn)的邏輯回歸分析,得出回歸分析結(jié)果(表2)。
表2 Logistic回歸分析結(jié)果
將各因子回歸系數(shù)代入式(5)和式(6),建立Logistic回歸模型為
P=1/[exp(-10.756+5.423a+8.153b+
5.624c+3.334d+6.603e+1.621f+
5.025g+2.329h+4.069i)] 。
(7)
式中:P為滑坡發(fā)生概率,取值范圍為0~1;a為高程因子歸一化值;b為坡度因子歸一化值;c為坡向因子歸一化值;d為地表曲率因子歸一化值;e為工程地質(zhì)巖組因子歸一化值;f為NDVI因子歸一化值;g為距道路距離因子歸一化值;h為距河流距離因子歸一化值;i為距斷層距離因子歸一化值。
建立邏輯回歸模型后,利用ArcGIS的柵格計(jì)算器工具,將各評(píng)價(jià)因子回歸系數(shù)的柵格圖層進(jìn)行疊加運(yùn)算,得到研究區(qū)滑坡災(zāi)害易發(fā)性概率,概率區(qū)間為[0.008,0.963],采用自然斷點(diǎn)法[22],按照易發(fā)性程度將易發(fā)性概率進(jìn)行分級(jí),可分為極低易發(fā)區(qū)[0.008,0.176)、低易發(fā)區(qū)[0.176,0.323)、中等易發(fā)區(qū)[0.323,0.481)、高易發(fā)區(qū)[0.481,0.688)、極高易發(fā)區(qū)[0.688,0.963]共5級(jí),采用重分類生成研究區(qū)滑坡易發(fā)性等級(jí)分區(qū),結(jié)果如圖10所示。
圖10 研究區(qū)滑坡易發(fā)性等級(jí)分區(qū)
由圖11可知,國道213汶川—松潘段滑坡分布特征如下: ①極高易發(fā)區(qū)和高易發(fā)區(qū)主要集中分布在汶川縣映秀鎮(zhèn)至茂縣石大關(guān)鄉(xiāng)一帶,該區(qū)域地勢(shì)陡峻,斜坡坡度較大,巖體破碎,易發(fā)生滑坡災(zāi)害。其中,極高易發(fā)區(qū)面積約106.75 km2,占研究區(qū)總面積的4.31%,發(fā)育滑坡92處,滑坡點(diǎn)密度為0.862處/km2; 高易發(fā)區(qū)面積約505.25 km2,占研究區(qū)總面積的12.15%,發(fā)育滑坡54處,滑坡點(diǎn)密度為0.179處/km2; ②滑坡的中易發(fā)區(qū)面積約492.54 km2,占研究區(qū)總面積的19.83%,發(fā)育滑坡37處,滑坡點(diǎn)密度約為0.075處/km2; ③低易發(fā)區(qū)和極低易發(fā)區(qū)主要集中在茂縣石大關(guān)鄉(xiāng)以北至松潘縣城一帶。其中,低易發(fā)區(qū)面積約505.25 km2,占研究區(qū)總面積的20.34%,發(fā)育滑坡25處,滑坡點(diǎn)密度為0.049處/km2; 極低易發(fā)區(qū)面積約1 077.23 km2,占研究區(qū)總面積的43.38%,發(fā)育滑坡22處,滑坡點(diǎn)密度為0.021處/km2。低易發(fā)區(qū)和極低易發(fā)區(qū)范圍內(nèi)地勢(shì)相對(duì)平緩,斜坡坡度較小,坡度大部分約20°,發(fā)生滑坡的可能性較小。
由分析結(jié)果可知: 隨著滑坡易發(fā)性等級(jí)的升高,易發(fā)性面積占研究區(qū)總面積的比例減少,發(fā)育滑坡數(shù)量增加,因而滑坡點(diǎn)密度不斷增大,說明易發(fā)性分區(qū)是合理的。
3.7.1 顯著性檢驗(yàn)
Logistic回歸模型顯著性檢驗(yàn)指的是因變量與自變量之間相關(guān)關(guān)系的檢驗(yàn)。采用Wald檢驗(yàn)法[23],若各評(píng)價(jià)因子的顯著性均小于0.05,表示自變量通過95%的顯著水平檢驗(yàn),即邏輯回歸模型的檢驗(yàn)效果好[20-22]。從表2可知,各評(píng)價(jià)因子顯著性均小于0.05,自變量通過了95%的顯著水平檢驗(yàn),模型顯著性較高。
3.7.2 精度檢驗(yàn)
易發(fā)性評(píng)價(jià)結(jié)果精度可以采用受試者工作特征曲線(receiver operating characteristic curve,ROC)進(jìn)行檢驗(yàn)[24]。ROC曲線是依據(jù)多種不同種類的二分類方式,以敏感度為縱坐標(biāo),以1-特異度為橫坐標(biāo)繪制的曲線,表示敏感度和特異的相互關(guān)系。ROC曲線利用曲線下面積(the area under the ROC curve,AUC)來評(píng)定評(píng)價(jià)結(jié)果的精度,AUC的取值范圍為[0.5,1],其值越大,檢驗(yàn)越具有價(jià)值,評(píng)價(jià)結(jié)果的精度越高。檢驗(yàn)結(jié)果得出AUC值為0.837,表明評(píng)價(jià)結(jié)果的精度較高,說明采用邏輯回歸模型對(duì)國道213汶川—松潘段滑坡易發(fā)性進(jìn)行評(píng)價(jià)是較為客觀準(zhǔn)確的(圖11)。
圖11 Logistic回歸模型的ROC曲線
本文以國道213汶川—松潘段為研究區(qū),綜合運(yùn)用光學(xué)衛(wèi)星遙感與InSAR變形探測(cè)技術(shù),結(jié)合統(tǒng)計(jì)分析與邏輯回歸評(píng)價(jià)方法,開展研究區(qū)滑坡隱患綜合遙感識(shí)別和易發(fā)性評(píng)價(jià)研究,得出以下結(jié)論:
(1)采用光學(xué)遙感和InSAR變形觀測(cè)技術(shù)對(duì)研究區(qū)進(jìn)行滑坡隱患識(shí)別,共識(shí)別出滑坡隱患288處,其中InSAR技術(shù)探測(cè)出活動(dòng)滑坡27處,說明綜合遙感識(shí)別技術(shù)在該研究區(qū)具有較好的識(shí)別效果。
(2)極高易發(fā)區(qū)和高易發(fā)區(qū)主要集中分布在汶川—茂縣一帶,總面積約612 km2; 極低易發(fā)區(qū)和低易發(fā)區(qū)主要集中在茂縣以北至松潘一帶,總面積約1 582 km2。同時(shí),采用ROC曲線進(jìn)行易發(fā)性結(jié)果精度分析,AUC值為0.837,表明評(píng)價(jià)結(jié)果的精度較高,說明采用Logistic回歸模型對(duì)國道213汶川—松潘段滑坡隱患進(jìn)行易發(fā)性評(píng)價(jià)是較為客觀準(zhǔn)確的。
(3)分析9個(gè)評(píng)價(jià)因子的邏輯回歸系數(shù),結(jié)果表明坡度、高程、坡向、距道路距離4個(gè)評(píng)價(jià)因子對(duì)國道213汶川—松潘段滑坡易發(fā)程度的貢獻(xiàn)相對(duì)較大,其中坡度[35°,45°)、高程[1 300,1 600) m、坡向[67.5°,112.5°)、距道路距離<100 m是滑坡最容易發(fā)生的地區(qū)。