陳攀， 葛永剛, 孫慶敏, 梁馨月

(1.中國科學(xué)院水利部成都山地災(zāi)害與環(huán)境研究所， 成都 610041； 2.西藏大學(xué)工學(xué)院， 拉薩 850000； 3.中國科學(xué)院山地災(zāi)害與地表過程重點實驗室， 成都 610041； 4.中國科學(xué)院大學(xué)， 北京 100049)

泥石流通常認為是含有大量松散固體物質(zhì)的固液兩相流[1]。因其具有巨大的體積和極快的流速，極易損壞附近橋梁、道路和房屋建筑，威脅當?shù)鼐用竦纳踩?。近些年?shù)次重大泥石流災(zāi)害事件[2-4]給國家和人民帶來巨大的財產(chǎn)損失。進行泥石流易發(fā)性評估能識別出研究區(qū)域內(nèi)的泥石流高易發(fā)區(qū)，為工程建設(shè)和防范治理提供決策依據(jù)。

評價單元、指標因子和評價模型作為影響泥石流易發(fā)性評估結(jié)果的重要因素，不同學(xué)者已對其進行了討論研究。鄒強等[5]通過對比水文響應(yīng)單元與柵格單元，發(fā)現(xiàn)利用水文響應(yīng)單元進行泥石流易發(fā)性評價具有諸多優(yōu)點。此后，基本選擇小流域單元作為泥石流易發(fā)性評價單元[6]。易發(fā)性評價指標基本為孕災(zāi)環(huán)境因子，且主要圍繞地形地貌、構(gòu)造活動、徑流水系進行選擇[7-8]，由于區(qū)域間存在差異性，所以并沒有形成統(tǒng)一的指標因子選擇標準。泥石流易發(fā)性評價模型按評估方式不同可分為四類，包括主觀專家分析[9-10]、統(tǒng)計分析[11-15]、機器學(xué)習[16-19]和物理模型[20-21]。主觀專家分析具有較強主觀性，所得評估結(jié)果精度較低。統(tǒng)計分析較主觀專家分析更加客觀，但未反映評價指標與泥石流之間復(fù)雜的非線性關(guān)系[22]。機器學(xué)習對訓(xùn)練樣本要求較高，且易存在過擬合。物理模型需要足夠的物理力學(xué)參數(shù)以保證評估的準確性，通常適合小范圍的易發(fā)性評價。目前，主要通過耦合主觀專家分析、統(tǒng)計模型和機器學(xué)習的方式提高泥石流易發(fā)性評估精度，缺乏對指標因子選擇的討論。對此，現(xiàn)通過對比有無坡體穩(wěn)定性因子，分別建立考慮坡體失穩(wěn)破壞的試驗組指標體系和僅有孕災(zāi)環(huán)境因子的對照組指標體系，采用確定性系數(shù)模型對安寧河流域進行易發(fā)性等級劃分，不僅討論指標因子對泥石流易發(fā)性評價結(jié)果的影響，同時也耦合統(tǒng)計模型與物理模型，增加泥石流評估結(jié)果的物理意義。

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

1.1 研究區(qū)概況

安寧河流域地處四川省西南部，位于第一階梯與第二階梯過渡區(qū)內(nèi)，整個流域約11 050 km2，自北向南依次主要分布于冕寧縣、喜德縣、西昌市、德昌縣、會理縣和米易縣，海拔高度介于951 ～ 5 137 m。地貌類型包含斷陷盆地、河谷平原、構(gòu)造剝蝕中山和冰蝕高山[23]。安寧河自北向南流域內(nèi)河谷地帶分布大量第四系松散物，河谷東側(cè)主要分布沉積巖，河谷西側(cè)主要分布巖漿巖。受構(gòu)造影響安寧河流域內(nèi)發(fā)育有控制性的安寧河斷裂和則木河斷裂，滑坡、崩塌沿斷層帶分布，為泥石流形成提供充足物源。研究區(qū)屬亞熱帶季風季候，年平均降雨量達1 000 mm以上[24]，為泥石流發(fā)育提供了水源條件。此外，人類火燒和不合理工程活動降低了坡體失穩(wěn)和泥石流形成條件，使安寧河流域內(nèi)泥石流活動頻發(fā)。

1.2 數(shù)據(jù)來源

研究所需數(shù)據(jù)來源如下：①12.5 m分辨率的ASTER DEM數(shù)據(jù)提取地形地貌和水系發(fā)育信息；②1∶250 000地質(zhì)圖提取安寧河流域地層巖性和構(gòu)造斷層信息；③1∶4 000 000 安寧河流域地震峰值加速度區(qū)劃圖；④Landsat 8衛(wèi)星數(shù)據(jù)提取植被覆蓋度信息；⑤基于坡體穩(wěn)定性分析獲得潛在失穩(wěn)面積與小流域面積比數(shù)據(jù)。

2 泥石流易發(fā)性評價

2.1 小流域單元提取

小流域包含泥石流形成和運移的全過程，基于小流域獲得的評估結(jié)果更能滿足實際需求。因此，利用ArcGIS軟件中的水文分析模塊，經(jīng)過洼地填充、流向提取、流量計算、河網(wǎng)定義和集水區(qū)生成5個步驟，劃分安寧河流域范圍的小流域單元。如圖1所示，共得到491個小流域單元，小流域面積介于2.38～179.48 km2，小流域平均面積為22.52 km2。

2.2 評價指標體系構(gòu)建

2.2.1 評價指標

泥石流易發(fā)性評價因子通常圍繞物源條件和地形條件進行選擇。結(jié)合安寧河流域?qū)嶋H情況，選取地表起伏度、地表粗糙度、溝壑密度、距斷層距離、地震峰值加速度和潛在物源分布面積與流域面積比作為試驗組指標因子，并選擇與之對應(yīng)的孕災(zāi)環(huán)境因子構(gòu)建對照組評價指標體系，對照組指標因子包含植被覆蓋度、地層巖性、地表起伏度、地表粗糙度、溝壑密度、距斷層距離和地震峰值加速度(圖2)。所選評價指標可分為孕災(zāi)環(huán)境因子和坡體穩(wěn)定性因子。

圖2 泥石流易發(fā)性評價因子

孕災(zāi)環(huán)境因子植被覆蓋度從水文機制和力學(xué)強度兩個方面影響坡體穩(wěn)定性，同時植被覆蓋度對水土流失也具有重要影響。不同的地層巖性具有不同的坡體結(jié)構(gòu)強度，是常見的基礎(chǔ)地質(zhì)因子。地表起伏度和地表粗糙度均可表征地表侵蝕程度和地形的起伏變化，前者同時體現(xiàn)小流域物質(zhì)運移具有的位置勢能大小，后者則從微觀反映地表單元地勢起伏的復(fù)雜程度。溝壑密度描述了小流域范圍內(nèi)水道發(fā)育程度，反映了小流域內(nèi)受到一定規(guī)模水流沖刷的區(qū)域所占比例。距斷層距離是內(nèi)動力地質(zhì)構(gòu)造指標因子，通常斷層活動會加劇巖體破壞，影響崩塌和滑坡分布。地震峰值加速度反映了地震活動情況，表征了坡體失穩(wěn)激發(fā)因素的影響。

坡體穩(wěn)定性因子潛在物源分布面積與小流域面積比表示單位小流域面積內(nèi)的潛在失穩(wěn)面積。在保證計算精度的前提下合理簡化采樣及計算過程，通過根-土分區(qū)的方式對坡體穩(wěn)定性進行分析，具體過程見流程圖(圖3)?？紤]植被覆蓋度和土體飽和度的影響，利用安全系數(shù)計算公式[25]，計算降雨飽和狀態(tài)下流域內(nèi)坡體安全系數(shù)，以1.5作為穩(wěn)定坡體與潛在失穩(wěn)坡體閾值[26]，劃分出安寧河流域坡體潛在失穩(wěn)區(qū)。通過對比流域內(nèi)地質(zhì)災(zāi)害分布情況和坡體潛在失穩(wěn)區(qū)計算結(jié)果，發(fā)現(xiàn)從區(qū)域上看兩者具有相同的分布規(guī)律。

圖3 潛在失穩(wěn)區(qū)計算流程圖

2.2.2 多重共線性診斷

泥石流易發(fā)性評價指標體系要求評價因子間互不干擾，通常通過共線性檢驗結(jié)果進行判斷，當方差膨脹因子(variance inflation factor，VIF)大于10時，表示指標因子之間存在多重共線性，反之，則所選指標因子間相互獨立。

通過SPSS軟件對所選指標因子進行共線性診斷，容忍度(tolerance, TOL)和方差膨脹因子(VIF)結(jié)果如表1所示。由表1可知，試驗組和對照組所選各指標因子的VIF值均小于4，說明各指標因子之間基本相互獨立，可參與模型進行泥石流易發(fā)性預(yù)測。

表1 多重共線性診斷結(jié)果

2.3 確定性系數(shù)模型

確定性系數(shù)模型為概率函數(shù)，基于泥石流災(zāi)害形成條件不變的假設(shè)，可計算不同指標因子對泥石流發(fā)生的敏感性，現(xiàn)已作為泥石流易發(fā)性評價的常見方法之一[27]。其計算公式為

(1)

式(1)中：CF為確定性系數(shù)值；PPa為評價因子a中的條件概率，可用評價因子a中的泥石流災(zāi)害個數(shù)或面積與評價因子a面積比值表示；PPs為泥石流災(zāi)害在研究區(qū)的先驗概率，可用研究區(qū)的泥石流災(zāi)害個數(shù)或面積與研究區(qū)面積比值表示。

2.4 易發(fā)性評價模型計算結(jié)果

從198條泥石流溝中隨機選擇80%的樣本數(shù)據(jù)作為訓(xùn)練樣本，參考常見斷點方式和自然斷點法利用確定性系數(shù)模型計算不同指標因子的分級CF值。計算結(jié)果如表2所示。

表2 各指標因子不同分級CF值計算結(jié)果

2.5 泥石流易發(fā)性評估

基于確定性系數(shù)模型計算結(jié)果，利用ArcGIS將研究區(qū)內(nèi)各評價因子指標的CF值疊加，利用自然段點法將安寧河流域泥石流易發(fā)性劃分為高易發(fā)、較高易發(fā)、中易發(fā)、較低易發(fā)和低易發(fā)5個易發(fā)性等級區(qū)，結(jié)果如圖4所示。

圖4 安寧河流域泥石流易發(fā)性分區(qū)

2.6 泥石流易發(fā)性評估結(jié)果驗證

利用頻率比值對泥石流易發(fā)性結(jié)果進行檢驗，結(jié)果如表3所示。由表3可知,試驗組指標體系與對照組指標體系的泥石流易發(fā)性高易發(fā)區(qū)與較高易發(fā)區(qū)的頻率比值之和分別占各自總頻率比值的88.17%和85.88%，均大于85%，表明兩種指標體系均有效評價了安寧河流域泥石流易發(fā)性。

表3 泥石流易發(fā)性分區(qū)頻率比值

以泥石流災(zāi)害易發(fā)性評價中代表研究區(qū)內(nèi)各易發(fā)性等級面積累加百分比為橫軸(假陽性率)，泥石流災(zāi)害易發(fā)性評價中代表研究區(qū)各易發(fā)性等級內(nèi)真實發(fā)生泥石流災(zāi)害的累加百分比為縱軸(真陽性率)，繪制的接受者操作特征(receiver operating characteristics，ROC)曲線如圖5所示，對照組指標體系的泥石流易發(fā)性ROC曲線下面積AUC值為0.752，試驗組指標體系的泥石流易發(fā)性AUC值達0.767，與對照組指標體系相比預(yù)測精度值提升了1.99%。表明考慮坡體穩(wěn)定性分析的泥石流易發(fā)性評估結(jié)果較僅考慮孕災(zāi)環(huán)境因子的結(jié)果好。

圖5 不同指標體系ROC曲線

圖6中，中壩村泥石流溝和老洼溝泥石流溝均在試驗組指標體系泥石流易發(fā)性評價結(jié)果的較高易發(fā)區(qū)內(nèi)，表明泥石流易發(fā)性評估結(jié)果與實際泥石流發(fā)育吻合度較高。

圖6 實例驗證

3 討論

坡體穩(wěn)定性指標耦合了地形和力學(xué)強度參數(shù)，考慮了孕災(zāi)環(huán)境因子與泥石流之間非線性關(guān)系。對照組指標體系中的地表粗糙度與坡度間存在共線性，為保證評估結(jié)果的合理性同時與試驗組指標體系形成對照，綜合考慮選擇地表粗糙度作為指標因子。試驗組指標體系則將地表粗糙度(地形地貌因素)和坡度(物源穩(wěn)定性因素)在相互獨立的條件下同時引入，使構(gòu)建的指標體系更加全面。試驗組指標體系相較于傳統(tǒng)孕災(zāi)環(huán)境指標體系加入坡體穩(wěn)定性因子，增加了泥石流形成機理分析，提高了泥石流易發(fā)性評估精度。

4 結(jié)論

利用確定性系數(shù)模型，對考慮坡體穩(wěn)定性分析與僅有孕災(zāi)環(huán)境因子的不同指標體系進行泥石流易發(fā)性評價，得出如下結(jié)論。

(1)兩種指標體系所得泥石流易發(fā)性高易發(fā)區(qū)和較高易發(fā)區(qū)頻率比值之和均占總頻率比值的85%以上，表明兩種指標體系都準確有效的評價了安寧河流域泥石流易發(fā)性。

(2)考慮坡體穩(wěn)定性分析的泥石流易發(fā)性評價結(jié)果的AUC值達0.767，相較于僅有孕災(zāi)環(huán)境因子的AUC值提升了1.99%。

(3)物理強度參數(shù)通常適用于局部區(qū)域的泥石流物源區(qū)分析，綜合坡體穩(wěn)定性分析與概率模型可對區(qū)域泥石流易發(fā)性進行有效評估，為泥石流易發(fā)性評估從孕災(zāi)環(huán)境分析提升至形成機理分析提供了新思路。