杜 俊, 丁文峰, 范仲杰, 李清溪
(1.長(zhǎng)江科學(xué)院 水土保持研究所, 湖北 武漢 430010; 2.水利部 山洪地質(zhì)災(zāi)害防治工程技術(shù)研究中心, 湖北 武漢 430010)
三峽庫區(qū)是中國(guó)典型的生態(tài)環(huán)境脆弱區(qū),2003年三峽水庫蓄水以來,庫區(qū)滑坡問題愈發(fā)突出,大量典型案例的調(diào)查、觀測(cè)和試驗(yàn)表明,水庫水位的周期性升降會(huì)改變滑體靜/動(dòng)水壓力、降低巖土體抗剪強(qiáng)度,進(jìn)而促進(jìn)滑體變形、滑帶貫通以及最終的累積釋放[1-4]。防汛部門的資料顯示,早在三峽水庫蓄水以前,庫區(qū)滑坡即已十分發(fā)育,是該區(qū)山洪地質(zhì)災(zāi)害(即溪河洪水以及由其引發(fā)的滑坡、泥石流等)[5]的主要表現(xiàn)形式。顯然水庫蓄水只是外因,庫區(qū)特有的自然地理和地質(zhì)條件才是滑坡發(fā)育的實(shí)質(zhì),因此只有弄清當(dāng)?shù)鼗掳l(fā)育格局與具體根植性要素的關(guān)系,才可能從根本上認(rèn)知區(qū)域滑坡宏觀發(fā)育機(jī)制和開展更有針對(duì)性的防治工作。以往研究對(duì)于這一問題已有一定積累,鄧清祿[6]和陳劍[7]等指出庫區(qū)滑坡歷史上主要與新構(gòu)造期以來的快速抬升和間冰期極大降水有關(guān),因此許多大型滑坡在長(zhǎng)江干流和主要支流深切方向上沿程分布;白世彪等[8]基于GIS技術(shù)和頻率統(tǒng)計(jì),對(duì)三峽水庫175 m回水影響區(qū)的滑坡分布進(jìn)行了統(tǒng)計(jì),總結(jié)了滑坡分布與巖性巖組、高程、坡度、坡向、曲率等的敏感性特征;喬建平等[9]使用本底因子貢獻(xiàn)率法統(tǒng)計(jì)了地層、坡度、坡形、高差和坡向?qū)θ龒{庫區(qū)云陽—巫山段滑坡分布的敏感性特征,并認(rèn)為高差、坡度和地層的危險(xiǎn)度權(quán)重更高;Wang和Niu[10]應(yīng)用中巴地球資源衛(wèi)星遙感影像,采集植被覆蓋、坡體結(jié)構(gòu)、水庫水位、高程、工程巖組等20個(gè)指標(biāo)信息,基于決策數(shù)方法對(duì)秭歸縣郭家壩鎮(zhèn)的滑坡分布進(jìn)行了高精度預(yù)報(bào);李雪等[11]分析了三峽庫首區(qū)滑坡分布與12項(xiàng)地形因子的相關(guān)性,認(rèn)為高程、坡向、高程變異系數(shù)與滑坡的相關(guān)性較高;張俊等[12]總結(jié)的萬州區(qū)滑坡發(fā)育主要因子為地層巖性、地質(zhì)構(gòu)造、水系分布、坡度、坡向、坡體結(jié)構(gòu)及土地利用。綜觀這些工作,大部分學(xué)者都十分重視統(tǒng)計(jì)滑坡發(fā)生頻率對(duì)地質(zhì)、地形類因子的敏感性特征,而對(duì)各類因子之于滑坡分布的相對(duì)重要性缺乏興趣,更不會(huì)探究因子之間可能存在的交互影響;在因子選用方面也存在一定重復(fù)。此外,防汛部門所關(guān)注的溪洪—滑坡,主要指由于山區(qū)溪河洪水沖蝕坡岸所引發(fā)的各類急性或慢性滑坡,其收集的歷史災(zāi)害信息可能與流域水文屬性具有更好的響應(yīng),這在以往研究中也較少關(guān)注??臻g序列分析可以探究地形、土壤等短時(shí)間尺度內(nèi)不隨時(shí)間變化因子對(duì)地理現(xiàn)象的影響,是因地制宜的開展目標(biāo)對(duì)象關(guān)鍵控制因子研究的基本方法論。另一方面,地理探測(cè)器模型(Geodetector)是一種探測(cè)空間分異性并揭示其背后驅(qū)動(dòng)力的一組統(tǒng)計(jì)學(xué)方法[13],它不僅可以分析各自變量對(duì)因變量的影響程度,還可以探測(cè)自變量間交互作用于因變量的影響。據(jù)此,本文嘗試基于空間分析視角,以三峽庫區(qū)香溪河流域的溪洪—滑坡災(zāi)害為研究對(duì)象,結(jié)合經(jīng)典統(tǒng)計(jì)和Geodetector方法,系統(tǒng)分析研究對(duì)象的空間分異特征及其與主要自然因子的關(guān)系,探明以該流域?yàn)榇淼膸靺^(qū)東部川鄂褶皺山地溪洪—滑坡災(zāi)害空間分布的關(guān)鍵控制因子,對(duì)比分析地理探測(cè)器模型在本研究中所體現(xiàn)的優(yōu)勢(shì)與局限。
香溪河位于長(zhǎng)江上游三峽庫區(qū)左岸,發(fā)源于神農(nóng)架南麓,由北向南流經(jīng)湖北省興山縣和秭歸縣,于香溪鎮(zhèn)入注長(zhǎng)江,干流全長(zhǎng)約94 km,流域面積3 193 km2,是長(zhǎng)江上游一級(jí)支流中距離三峽大壩壩首最近的中尺度河流。受亞熱帶大陸性季風(fēng)氣候影響,本區(qū)四季分明,降水充沛,年均雨量1 000 mm以上,北部較南部稍高;流域地貌屬鄂西褶皺山地,為大巴山和巫山的余脈,地形以中低山地為主,山勢(shì)陡峻、河谷深切;流域植被曾于20世紀(jì)80年代遭到嚴(yán)重破壞,1989年起經(jīng)過一系列水土保持項(xiàng)目治理后得到明顯改善,現(xiàn)森林覆蓋率達(dá)56%。流域涉及的主要構(gòu)造形跡為近SN向展布的黃陵背斜和秭歸向斜,其中黃陵背斜一直處于上升趨勢(shì),每年上升2~4 mm[14];區(qū)內(nèi)斷裂活動(dòng)總體不顯著;出露巖層以碳酸鹽巖和砂、泥巖為主,主要包括三疊系中統(tǒng)嘉陵江組上部厚層白云制灰?guī)r、嘉陵江組中部中—厚層灰?guī)r、嘉陵江組下部薄層灰?guī)r,侏羅系上統(tǒng)遂寧組磚紅色泥巖與砂巖互層,侏羅系上統(tǒng)沙溪廟組紅色泥巖和石英巖等。
本研究著眼于目標(biāo)流域溪洪—滑坡災(zāi)害與主要自然因子的關(guān)系,因此在資料選取上力圖排除2002年之后水庫蓄水的影響,涉及的主要數(shù)據(jù)和來源詳見表1。
表1 研究區(qū)基礎(chǔ)數(shù)據(jù)與來源
一般認(rèn)為,溪河洪水的形成與強(qiáng)降水、大起伏地形、植被生長(zhǎng)、土壤屬性、流域形狀等自然要素密切相關(guān),但對(duì)于由溪洪引發(fā)的(及一般的)滑坡、泥石流,巖性、斷層(裂)、松散物質(zhì)分布等地質(zhì)要素的作用也非常重要。
基于學(xué)界對(duì)山洪、滑坡成因的基本認(rèn)知[8-12,14-15],結(jié)合資料收集情況,選取暴雨、土壤、巖性、斷層、歷史災(zāi)害等指標(biāo)作為分析基礎(chǔ)構(gòu)建指標(biāo)體系,各指標(biāo)圖層如圖1所示。
圖1 香溪河流域主要自然因子的過程?hào)鸥駡D
1.3.1 暴雨〔多年平均年暴雨極值(A1)與變異系數(shù)Cv(A2)〕 多年平均暴雨極值反映區(qū)域暴雨的絕對(duì)能量,暴雨極值的變異系數(shù)反映其不確定性,與溪洪、滑坡災(zāi)害均存在理論上的聯(lián)系。統(tǒng)計(jì)收集的10個(gè)站點(diǎn)[16]1971—2002多年平均年暴雨極值量和Cv值,在GIS中進(jìn)行反距離權(quán)重(IDW)插值,以結(jié)果的均方根誤差(RMSE)最小為宜,得到流域面上的A1和A2。
1.3.2 一般下墊面〔土地利用產(chǎn)匯流能力(A3)、土壤入滲性(A4)與NDVI(A5)〕 地表綜合產(chǎn)匯流能力不僅與溪洪、滑坡發(fā)育有直接關(guān)聯(lián),也不同程度影響著地下水分布及孔隙水壓力情況。研究表明黏粒占比較少、有機(jī)質(zhì)含量較高的土壤類型,或土地利用中林草覆蓋較多的類型,下墊面入滲性能較好[17-19],不易產(chǎn)流形成洪水,但理論上有增加滑坡體重力勢(shì)能的可能,參考文獻(xiàn)[15]的工作,對(duì)相關(guān)指標(biāo)進(jìn)行單向打分;NDVI直接采用2000—2002年4—9月研究區(qū)NDVI值進(jìn)行算術(shù)平均得到。
1.3.3 地形〔坡度(A6)與坡向(A7)〕 依據(jù)文獻(xiàn)[8-9]對(duì)庫區(qū)敏感性坡度的總結(jié),應(yīng)用GIS中的“slope”和“raster calculator”工具提取流域10°~45°之間的坡度,其余坡度賦予低值;坡向由GIS中“aspect”工具對(duì)DEM處理后,再對(duì)方向分級(jí)得到。
1.3.4 地質(zhì)巖性〔脆弱性(A8)、巖性軟硬程度(A9)與斷層緩沖區(qū)(A10)〕 地層巖性是滑坡發(fā)育的重要考量因素,本研究結(jié)合巖性類別、表面破碎程度、產(chǎn)狀水平和完整性,使用專家打分法得到巖性脆弱性指標(biāo)[20](表2)。此外,依據(jù)巖性軟硬分級(jí),得到軟硬程度指標(biāo);汶川震區(qū)的研究表明[21],大部分大型滑坡點(diǎn)分布在距離斷層5 km范圍內(nèi),庫區(qū)斷層活躍程度相對(duì)較弱,故緩沖區(qū)依距離斷層0.5,1,2,3,10,30 km設(shè)置。
表2 香溪河巖性脆弱性賦值
1.3.5 小流域水文屬性〔小流域面積(A11)、河網(wǎng)密度(A12)、主溝坡降(A13)和形狀系數(shù)(A14)〕 發(fā)育溪洪的小流域一般面積較小,基于DEM,利用GIS的水文分析工具,得到香溪河流域43個(gè)面積在6.5~183.6 km2的小流域,并利用GIS的地統(tǒng)計(jì)分析功能分別提取這些小流域大于2 km2的河網(wǎng)密度、主溝坡降和形狀系數(shù)。
1.3.6 歷史災(zāi)害〔歷史溪洪—滑坡災(zāi)害的綜合災(zāi)度(B)〕 由于很多歷史災(zāi)害點(diǎn)分布在小流域的邊緣地帶,不能充分反映流域特性的影響,且歷史災(zāi)害點(diǎn)本身只能表現(xiàn)已經(jīng)發(fā)生過的事件,對(duì)潛在的災(zāi)害易發(fā)范圍沒有明確指示?;凇熬嚯x歷史災(zāi)害點(diǎn)越近的地區(qū)孕災(zāi)環(huán)境越相似”這一假設(shè),使用空間插值的方法對(duì)收集的109個(gè)歷史災(zāi)害點(diǎn)進(jìn)行補(bǔ)點(diǎn)插值,具體方法是結(jié)合實(shí)地環(huán)境調(diào)查,對(duì)流域內(nèi)及周邊無明顯災(zāi)害記錄的地區(qū)補(bǔ)0點(diǎn),對(duì)發(fā)災(zāi)記錄較多地區(qū)的災(zāi)害點(diǎn)合并為高值點(diǎn),使網(wǎng)點(diǎn)在空間分布上相對(duì)均勻,最后使用IDW插值得到表征香溪河流域溪洪—滑坡災(zāi)害潛在易發(fā)程度的綜合災(zāi)度圖(圖2)。
1.4.1 逐步回歸 逐步回歸方法各類文獻(xiàn)中均有詳細(xì)介紹,這里不再累述。本研究使用GIS提取的香溪河43個(gè)小流域的綜合災(zāi)度(B)和除坡向(A7)外14項(xiàng)影響因子的平均值,以前者為因變量、后者為自變量構(gòu)筑逐步回歸模型,分析主要自然因子對(duì)溪洪—滑坡災(zāi)害的影響。
1.4.2 地理探測(cè)器模型(Geodetector model) 地理探測(cè)器模型是圍繞地理現(xiàn)象的分層(類)異質(zhì)性或空間分異性研發(fā)的一種統(tǒng)計(jì)學(xué)方法[13],它基于“兩變量空間分布越相似則關(guān)聯(lián)性越大”的假設(shè),通過對(duì)因變量的再分類計(jì)算因變量類內(nèi)方差與類間方差的關(guān)系,以此分析自變量對(duì)因變量的影響程度,其基本公式為:
(1)
除了上述基本的因子探測(cè)功能外,Geodetector還具有生態(tài)探測(cè)和交互作用探測(cè)功能,前者主要使用F檢驗(yàn)判斷兩自變量對(duì)因變量的影響是否存在顯著差異;后者的原理為比較兩自變量在GIS中疊加后生成的q值與相應(yīng)獨(dú)立自變量q值的關(guān)系,主要的交互作用方式有:
①q(A∩B) ② min(q(A),q(B)) ③q(A∩B)>max〔q(A),q(B)〕 (雙因子增強(qiáng)) ④q(A∩B)=q(A)+q(B) (獨(dú)立) ⑤q(A∩B)>q(A)+q(B) (非線性增強(qiáng)) 由于本研究自變量擬分類數(shù)較多,僅提取43個(gè)小流域的數(shù)值建立關(guān)聯(lián)可能導(dǎo)致失真,因此這里除去反映小流域統(tǒng)計(jì)屬性的A11—A14,對(duì)余下的11項(xiàng)變量的柵格圖統(tǒng)一采樣精度(100 m)并轉(zhuǎn)化為點(diǎn)矢量圖,每張柵格圖得到3 123個(gè)矢量點(diǎn),將自變量矢量點(diǎn)值分類后與因變量矢量點(diǎn)值聯(lián)立,采用因子和交互作用探測(cè)功能分析各自然因子對(duì)綜合災(zāi)度的影響。 依據(jù)上述方法,經(jīng)自然斷點(diǎn)法分類后,得到香溪河流域溪洪—滑坡綜合災(zāi)度顯示圖(圖2),綜合災(zāi)度是從宏觀格局上對(duì)區(qū)域潛在災(zāi)害易發(fā)程度的綜合反映,并不是實(shí)際發(fā)災(zāi)情況,因此不會(huì)具體到坡岸尺度。依據(jù)綜合災(zāi)度分布,香溪河流域的溪洪—滑坡災(zāi)害高易發(fā)區(qū)主要分布在興山以下的流域西南部,東部和西北也有少量分布,北部上游的大起伏山地反而發(fā)災(zāi)較少,可能與該區(qū)人口密度不高、植被覆蓋較好有關(guān)。 圖2 香溪河流域溪洪-滑坡綜合災(zāi)度分布 通過逐步回歸和地理探測(cè)器模型分析,認(rèn)為下墊面條件,特別是地層巖性與植被發(fā)育,是控制香溪河流域溪洪—滑坡綜合災(zāi)度空間分異的主要方面,降水和土壤因素的影響有限,小流域水文屬性要素的作用不明顯。逐步回歸結(jié)果顯示進(jìn)入模型的5項(xiàng)因子可以解釋因變量80.6%的變化,依據(jù)標(biāo)準(zhǔn)化回歸系數(shù)(表3),各因子的貢獻(xiàn)率依次為巖性脆弱性(A8)28.6%,斷層緩沖區(qū)(A10)19.8%,NDVI(A5)18.9%,暴雨極值(A1)16.4%,以及巖性軟硬程度(A9)16.3%,巖層多傾斜、破碎、巖性較軟弱、距離斷層較近、植被覆蓋少的地區(qū),發(fā)災(zāi)潛力更大。 表3 香溪河流域溪洪-滑坡綜合災(zāi)度與主要自然因子的逐步回歸結(jié)果 地理探測(cè)器模型分析的結(jié)果與回歸模型稍有不同(表4),依據(jù)各因子對(duì)因變量的解釋度(q值)排序,僅土壤入滲性(A4)取代了暴雨極值(A1)。巖性脆弱性(A8)不僅在因子探測(cè)分析中具有最高的解釋度,與其他4項(xiàng)因子還分別存在非線性增強(qiáng)的交互影響。實(shí)際上與地質(zhì)環(huán)境關(guān)系密切的巖性脆弱性(A8)、斷層緩沖區(qū)(A10)和巖性軟硬程度(A9)在兩類模型中均有穩(wěn)定的表現(xiàn),體現(xiàn)了地質(zhì)要素對(duì)區(qū)域滑坡分布絕對(duì)的控制力。有學(xué)者[22]認(rèn)為正是由于侏羅系黃陵背斜的褶皺隆升,引起秭歸向斜發(fā)生陸相繞曲沉積,導(dǎo)致沉積巖層產(chǎn)狀嚴(yán)重傾斜,加上沉積物質(zhì)來源不同,沉積相軟硬互層,才有了當(dāng)前香溪河中下游的滑坡集中區(qū),本文的結(jié)果在一定程度上支持這一觀點(diǎn)。此外,NDVI(A5)在兩類模型的影響力僅次于地層巖性與斷層類因子,且與巖性脆弱性(A8)在交互作用下可以解釋因變量68.8%的變化,顯示了其對(duì)控制滑坡的積極影響。 表4 香溪河流域溪洪-滑坡綜合災(zāi)度與主要 注:僅列出通過顯著性檢驗(yàn)且q值大于0.1的變量。 一些被認(rèn)為是滑坡發(fā)育基本因子的指標(biāo),如降水、地形類指標(biāo)在兩類模型中的期望低于預(yù)期,這一方面是因?yàn)榭臻g分析本身的“瓶頸效應(yīng)”,即區(qū)內(nèi)某類指標(biāo)的數(shù)值趨于均化,或早已滿足因變量變化的要求,不構(gòu)成發(fā)育瓶頸,則此類指標(biāo)的作用就很難體現(xiàn),如降水在香溪河乃至整個(gè)庫區(qū)均有非常充分的供給[16,23],不是滑坡發(fā)育的瓶頸,因此本研究中相關(guān)指標(biāo)的貢獻(xiàn)較低;地形與滑坡關(guān)系的特殊性也是原因之一,地形對(duì)滑坡發(fā)育的影響需與其他因子配合,并不是充分條件,因此Geodetector因子探測(cè)的q值較低(A6, 0.012;A7, 0.008),但交互探測(cè)顯示,坡度與巖性脆弱性的交互作用強(qiáng)度達(dá)62.0%,表明坡度對(duì)滑坡的影響需聯(lián)合地質(zhì)因素。考慮到此前也有案例反映坡度、坡向等指標(biāo)對(duì)滑坡發(fā)育的影響并不大[24-25],本研究認(rèn)為單純的坡度、坡向?qū)α饔蚧赂窬值挠绊懹邢?,它們的作用更多體現(xiàn)在協(xié)同方面。 Geodetector通過計(jì)算比較類別化了的自變量對(duì)因變量的類內(nèi)與類間方差,得出因變量與自變量間的關(guān)聯(lián)大小或相似程度,以此判斷自變量對(duì)因變量的影響程度,為定量分析非線性變量間的關(guān)聯(lián)性提供了新的視角。本例中該模型探測(cè)了各單一因子對(duì)溪洪—滑坡綜合災(zāi)度分布的影響,得到的各因子解釋度大小和排序與回歸結(jié)果相似,區(qū)別體現(xiàn)在土壤和暴雨因子的地位。直觀比較,土壤入滲性較之于暴雨極值,擁有與流域綜合災(zāi)度更好的空間對(duì)應(yīng)關(guān)系(圖1—2)。此外,即便使用Geodetector的矢量點(diǎn)數(shù)據(jù)進(jìn)行回歸,土壤入滲性的貢獻(xiàn)仍然較低(表5),Geodetector在處理全局性問題時(shí)對(duì)從屬因子擁有更高的敏感性和顯示度。然而,正是由于Geodetector對(duì)自變量的類別化處理,消除了自變量的數(shù)值屬性,使其無法探知自變量對(duì)因變量的影響方向,進(jìn)而判斷這種影響是否存在假相關(guān)。此外,雖然“生態(tài)探測(cè)”功能可以對(duì)各自變量于因變量的影響格局進(jìn)行判斷,但這種分析只能通過兩兩比較進(jìn)行,且結(jié)果僅能說明比較因子在影響格局上的相似性,而影響格局相似的因子并不一定是重復(fù)的?;貧w模型可以近似表達(dá)溪洪—滑坡綜合災(zāi)度與自然因子的定量關(guān)系和作用方向,偏相關(guān)或逐步回歸等方法也能從統(tǒng)計(jì)上剔除可能存在的假相關(guān)或重復(fù)性變量,與Geodetector具有一定的互補(bǔ)性。 表5 香溪河流域溪洪-滑坡綜合災(zāi)度與主要自然因子的矢量點(diǎn)數(shù)據(jù)逐步回歸結(jié)果 結(jié)合逐步回歸和Geodetector方法,從空間分異視角對(duì)三峽庫區(qū)香溪河流域的溪洪—滑坡災(zāi)害與其自然影響因子的關(guān)系進(jìn)行了定量研究。結(jié)果表明,進(jìn)入逐步回歸模型的5項(xiàng)因子總解釋度80.6%,貢獻(xiàn)率依次為:巖性脆弱性28.6%,斷層緩沖區(qū)19.8%,NDVI 18.9%,暴雨極值16.4%,巖性軟硬程度16.3%;地理探測(cè)中對(duì)因變量解釋度較大的5項(xiàng)因子分別為:性脆弱性24.5%,斷層緩沖區(qū)20.9%,NDVI 18.3%,土壤入滲性16.3%,巖性軟硬程度14.8%;巖性脆弱性與NDVI,坡度的交互作用強(qiáng)度分別達(dá)68.8%和62.0%;基于定量分析、定性理解的原則,認(rèn)為地質(zhì)要素該本地區(qū)溪洪—滑坡綜合災(zāi)度的空間分布處于支配地位,地層巖性越傾斜、破碎、軟弱,距離斷層越近,發(fā)災(zāi)可能性越大;植被生長(zhǎng)對(duì)綜合災(zāi)度分布也有重要影響,降水、土壤屬性的控制作用較小,小流域水文屬性與綜合災(zāi)度的關(guān)系不明顯,地形的影響體現(xiàn)在協(xié)同方面。 需要強(qiáng)調(diào)的是,本文的結(jié)果僅明確了研究區(qū)內(nèi)影響溪洪—滑坡災(zāi)害分布的瓶頸性因子,未進(jìn)入模型或貢獻(xiàn)低于預(yù)期的因子,未否認(rèn)它們對(duì)目標(biāo)現(xiàn)象的機(jī)理性影響,僅表明這些因子并非這一地區(qū)溪洪—滑坡災(zāi)害分布的瓶頸。Geodetector為定量分析非線性變量間的關(guān)聯(lián)性提供了新的視角,但它消除了自變量的數(shù)值屬性,無法探知自變量對(duì)因變量的影響方向,多變量分析時(shí)也很難判斷各自變量對(duì)因變量可能存在的重復(fù)性問題,因此建議結(jié)合經(jīng)典統(tǒng)計(jì)方法,對(duì)結(jié)果進(jìn)行綜合分析與判斷。 致謝:感謝中科院水利部成都山地災(zāi)害與環(huán)境研究所的張信寶研究員在巖性脆弱性賦值方面給予的支持與幫助;感謝長(zhǎng)江科學(xué)院重慶分院的盧陽高級(jí)工程師和秦維工程師提供的歷史災(zāi)害數(shù)據(jù)支持。2 結(jié)果與討論
2.1 香溪河流域溪洪-滑坡綜合災(zāi)度分布格局
2.2 香溪河流域溪洪-滑坡綜合災(zāi)度與主要自然因子的定量關(guān)系
2.3 Geodetector的應(yīng)用比較
3 結(jié) 論
——與非適應(yīng)性回歸分析的比較