基于GIS和Logistic模型的白龍江流域泥石流危險(xiǎn)性分析

田 菲，馬金輝，彌沛峰，屈 創(chuàng)

(蘭州大學(xué)資源環(huán)境學(xué)院，甘肅蘭州730000)

白龍江流域位于甘肅南部，屬秦嶺西段，地處青藏高原、黃土高原和四川盆地交匯處，是我國(guó)泥石流危害最為嚴(yán)重的地區(qū)之一［1］。白龍江中下游的泥石流分布面積和爆發(fā)頻率在省內(nèi)均居首位，泥石流災(zāi)害頻發(fā)給當(dāng)?shù)厝嗣裆?cái)產(chǎn)安全造成了極大威脅，嚴(yán)重阻礙了當(dāng)?shù)厣鐣?huì)經(jīng)濟(jì)的發(fā)展。開展該區(qū)域的泥石流災(zāi)害危險(xiǎn)性防治工作成為該地區(qū)的一項(xiàng)重要基礎(chǔ)工作?？茖W(xué)合理地防治泥石流災(zāi)害首先要對(duì)災(zāi)害發(fā)生的可能性、危險(xiǎn)性、危害范圍和程度有一個(gè)基本的認(rèn)識(shí)和評(píng)價(jià)［2］。泥石流危險(xiǎn)性評(píng)價(jià)是根據(jù)對(duì)山地斜坡賦存的環(huán)境條件和泥石流活動(dòng)歷史記錄的分析，劃分出各區(qū)域泥石流危險(xiǎn)等級(jí)并評(píng)估對(duì)人類和財(cái)產(chǎn)可能造成的災(zāi)害損失嚴(yán)重程度［3］。20世紀(jì)90年代以來，隨著GIS技術(shù)的發(fā)展，國(guó)外基于GIS技術(shù)的區(qū)域地質(zhì)災(zāi)害危險(xiǎn)性評(píng)價(jià)研究發(fā)展較快［4-5］，同時(shí)國(guó)內(nèi)也有不少學(xué)者針對(duì)地質(zhì)災(zāi)害危險(xiǎn)性評(píng)價(jià)進(jìn)行了深入研究［6-9］。當(dāng)前，在地質(zhì)災(zāi)害危險(xiǎn)性評(píng)價(jià)中采用的主要研究方法為專家打分法、層次分析法、人工神經(jīng)網(wǎng)絡(luò)法、信息量法、概率模型、判別分析等。其中，Logistic回歸模型是由標(biāo)準(zhǔn)回歸方程變換得來的一種非線性模型，原理簡(jiǎn)單、實(shí)現(xiàn)方便，結(jié)合GIS的空間分析功能在泥石流危險(xiǎn)性區(qū)劃方面應(yīng)用較多。國(guó)外Gorsevski等提出將Logistic回歸模型和GIS技術(shù)應(yīng)用于地質(zhì)災(zāi)害的空間預(yù)測(cè)［10］。Ohlmacher等在美國(guó)東北地區(qū)對(duì)地質(zhì)災(zāi)害進(jìn)行預(yù)測(cè)時(shí)，也用到了Logistic回歸模型［11］。Ayalew等將Logistic回歸模型應(yīng)用在滑坡易發(fā)性制圖方面［12］。國(guó)內(nèi)潘赟等［13］、叢威青等［14］結(jié)合 GIS 技術(shù)將 Logistic 模型應(yīng)用于地質(zhì)災(zāi)害危險(xiǎn)性區(qū)劃中。苗雨應(yīng)用GIS與Logistic回歸對(duì)延安市寶塔區(qū)做出滑坡災(zāi)害危險(xiǎn)性區(qū)劃［15］。

由此可見，GIS技術(shù)和數(shù)學(xué)模型結(jié)合應(yīng)用越來越廣泛，Logistic模型憑借其出色的非線性擬合優(yōu)勢(shì)也逐漸受到諸多學(xué)者的重視，在泥石流危險(xiǎn)性評(píng)價(jià)過程中應(yīng)用逐漸增多。該研究通過對(duì)白龍江流域武都段的野外實(shí)地調(diào)查和歷史資料收集，在結(jié)合研究區(qū)泥石流活動(dòng)特點(diǎn)、發(fā)育特征的基礎(chǔ)上，借鑒以往研究的經(jīng)驗(yàn)，并綜合專家學(xué)者意見，最終確定泥石流危險(xiǎn)性的主要影響因子。并在GIS的支持下，將影響因子根據(jù)特性分類，并與災(zāi)害發(fā)生與否(0，1)在統(tǒng)計(jì)軟件中進(jìn)行二元Logistic回歸分析，將擬合的模型運(yùn)用于整個(gè)研究區(qū)，最終得到了白龍江流域武都段泥石流災(zāi)害危險(xiǎn)性分區(qū)結(jié)果。

1 研究對(duì)象與研究方法

1.1 研究區(qū)概況 研究范圍為白龍江流域武都段，屬于白龍江中下游區(qū)，包括武都區(qū)、舟曲縣東部、宕昌縣南部和文縣北部山區(qū)，為泥石流連續(xù)分布區(qū)［16］。研究區(qū)范圍是33°02＇～33°50＇N 和104°25＇～ 105°20＇E，流域面積約 3 397.3 km2。該區(qū)地形陡峻、滑坡發(fā)育，降雨集中且多暴雨，泥石流的分布面積和爆發(fā)頻率在省內(nèi)均居首位并以粘性泥石流為主。經(jīng)統(tǒng)計(jì)，全研究區(qū)具有災(zāi)害性泥石流溝309條，占溝谷總面積的70%，區(qū)域上主要密集分布于武都區(qū)白龍江干流兩岸。圖1為研究區(qū)泥石流溝分布圖。

1.2 研究方法 通過收集研究區(qū)歷史資料和野外實(shí)地調(diào)查，借鑒以往研究經(jīng)驗(yàn)，并征求本地相關(guān)專家意見，選取了高程、坡度、地層巖性、距斷層距離、滑坡點(diǎn)密度、土地利用和NDVI共7個(gè)因子。在GIS支持下根據(jù)其性質(zhì)分別劃分為若干類。將整個(gè)研究區(qū)按300 m×300 m方格劃分為37 425個(gè)網(wǎng)格單元，其中已發(fā)生泥石流單元數(shù)為18 126個(gè)。將災(zāi)害發(fā)生與否(0，1)數(shù)據(jù)和分組處理后的因子數(shù)據(jù)在統(tǒng)計(jì)軟件中進(jìn)行二元Logistic回歸分析，計(jì)算泥石流發(fā)生概率。最終生成白龍江流域武都段泥石流災(zāi)害危險(xiǎn)性分區(qū)圖。

圖1 研究區(qū)泥石流溝分布示意

1.2.1 評(píng)價(jià)指標(biāo)的選擇。泥石流的形成除了受豐富的松散固體物質(zhì)、充足的水源和陡峭的地形3個(gè)必備條件之外還受控于多種因素。劉希林等選取地形指標(biāo)、巖性指標(biāo)、構(gòu)造指標(biāo)、地震指標(biāo)等提出了區(qū)域泥石流危險(xiǎn)度的多因子綜合評(píng)價(jià)模型［17］。邢釗從地層巖性、地質(zhì)構(gòu)造、地形地貌、地震活動(dòng)、水文氣象因子及人類活動(dòng)因子6個(gè)方面確定因子，進(jìn)行了白龍江流域的泥石流危險(xiǎn)性評(píng)價(jià)［18］。該研究通過分析當(dāng)?shù)氐哪嗍鳉v史災(zāi)害特征和地貌地質(zhì)環(huán)境條件，綜合考慮因子的完整性的基礎(chǔ)上，結(jié)合研究區(qū)實(shí)際工作資料，選取高程、坡度、地層巖性、距斷層距離、滑坡點(diǎn)密度、NDVI和土地利用7個(gè)因子作為參評(píng)因子。由于降水資料缺乏，未能選擇降雨因子。

1.2.2 影響因子分析。

1.2.2.1 高程。通過研究區(qū)泥石流分布和高程的關(guān)系分析可以得出，高程對(duì)泥石流的分布有重要影響。據(jù)統(tǒng)計(jì)，白龍江干流及其二級(jí)支流沿線多為河谷地帶，高程范圍在1 600～2 500 m和2 500～3 100 m內(nèi)屬于泥石流高發(fā)區(qū)。相反，高程值較大的區(qū)域由于地表多為基巖或有高山苔原植被覆蓋，不容易產(chǎn)生松散堆積物，泥石流不易發(fā)生。

1.2.2.2 坡度。地形坡度是控制著斜坡上松散物質(zhì)的堆積厚度以及斜坡體穩(wěn)定程度的因素。隨著坡度的增大，坡體穩(wěn)定性降低，坡面物源增多，相應(yīng)泥石流發(fā)生的概率也增大。如圖2所示，在坡度0～35°范圍內(nèi)，泥石流分布較少，而在＞25°范圍內(nèi)較活躍。

圖2 研究區(qū)坡度示意

1.2.2.3 地層巖性。研究區(qū)的地層巖性比較復(fù)雜，從志留系到第四系均有出露，而以志留系、泥盆系、石炭系、二疊系、三疊系分布較為廣泛。根據(jù)區(qū)域地質(zhì)資料將研究區(qū)巖性分為7類巖組(表1)，通過分析及實(shí)地考察驗(yàn)證，泥石流溝在D和H的巖組中分布較多，說明這兩類巖組巖性比較脆弱，極易風(fēng)化成泥石流的物源。

1.2.2.4 距斷層距離。通過計(jì)算任一網(wǎng)格到斷層的最近距離，發(fā)現(xiàn)距離斷層中心線越近，泥石流發(fā)生的危險(xiǎn)性就越大。原因是地層巖石的破碎程度越高，所產(chǎn)生的松散固體物質(zhì)就越多，為泥石流提供了充足的物源條件。

表1 影響因子分類

1.2.2.5 滑坡點(diǎn)密度?；率悄嗍鳟a(chǎn)生的又一重要物源，滑坡主要沿白龍江干流及其二級(jí)支流的河谷兩側(cè)發(fā)育，尤以武都段密度最高。結(jié)合研究區(qū)歷史災(zāi)害情況，采用密度分析，發(fā)現(xiàn)滑坡點(diǎn)密度的分布與泥石流溝分布規(guī)律基本相一致(圖3)。

圖3 研究區(qū)滑坡密度示意

1.2.2.6 NDVI。植被分布是影響泥石流是否形成的重要因素之一。NDVI常常用來表示植被的覆蓋程度。統(tǒng)計(jì)發(fā)現(xiàn)，研究區(qū)內(nèi)NDVI＜0.3，是泥石流災(zāi)害密集區(qū)和高發(fā)區(qū)，說明在此范圍內(nèi)，植被覆蓋度較低，植被固土保水能力差，巖土體失穩(wěn)，易于泥石流的發(fā)生。研究區(qū)內(nèi)NDVI＞0.60的區(qū)域面積較小，此范圍內(nèi)分布的泥石流面積也較小。

1.2.2.7 土地利用類型。研究區(qū)的土地利用類型主要為耕地、林地、草地、水域、工礦居民用地和未利用土地。通過實(shí)地走訪調(diào)查發(fā)現(xiàn)該地區(qū)城鄉(xiāng)工礦及居民點(diǎn)主要建設(shè)在溝谷里、溝口泥石流沖積扇上，或古滑坡體上。加之長(zhǎng)期以來不合理的人類活動(dòng)改變了原有的地貌形態(tài)、坡度與下墊面條件，使得地表土層結(jié)構(gòu)遭到嚴(yán)重破壞，形成了大量泥石流物源，在強(qiáng)降雨條件下容易爆發(fā)泥石流。

1.2.3 Logistic回歸模型。Logistic回歸，指因變量為二值(0-1)分類變量的回歸分析(Binary Logistic Regression Analysis)。以災(zāi)害發(fā)生與否為因變量，影響因子為自變量，設(shè)P為災(zāi)害發(fā)生的概率，取值范圍［0，1］(發(fā)生為1，否則為0)并服從二項(xiàng)分布，通過Logit變換，即有Logistic回歸模型:

式中，X1，X2，…，Xm為影響 P 的變量，β0為常數(shù)，β1，β2…βm表示邏輯回歸系數(shù)，如果 β0，β1，β2…βm確定，能夠根據(jù)不同的指標(biāo)值Xi計(jì)算某一區(qū)域發(fā)生災(zāi)害的概率。根據(jù)P值的大小，可以劃分災(zāi)害發(fā)生可能性等級(jí)。

2 泥石流災(zāi)害危險(xiǎn)性分析

2.1 回歸模型的建立 Logistic回歸分析要求因變量是分類變量，自變量可以是數(shù)值型連續(xù)變量，也可以是順序型分類變量，如果是名義變量，可以轉(zhuǎn)換成啞變量來處理。關(guān)于處理變量，國(guó)內(nèi)許多研究人員采用了不同方法，李雪平等采用所有因子分組的辦法評(píng)價(jià)了斜坡穩(wěn)定性［19］。王衛(wèi)東等單個(gè)災(zāi)害因子分級(jí)，并在GIS環(huán)境下結(jié)合已發(fā)生的地質(zhì)災(zāi)害進(jìn)行二元空間統(tǒng)計(jì)分析，得到單因子各級(jí)別下災(zāi)害面積和實(shí)現(xiàn)各因子量綱統(tǒng)一［20］。該研究收集的7個(gè)因子數(shù)據(jù)，其中有坡度、高程等連續(xù)變量，也有巖性、土地利用等名義變量。根據(jù)數(shù)據(jù)情況，對(duì)名義變量進(jìn)行了分類，結(jié)果見表1。

在GIS中將數(shù)據(jù)按以上的方法處理完畢后，再將數(shù)據(jù)表導(dǎo)入統(tǒng)計(jì)分析軟件中進(jìn)行Logistic回歸，可得到結(jié)果見表2。

表2 回歸系數(shù)分析

通過表2可以看出，回歸系數(shù)的檢測(cè)結(jié)果是顯著的，有統(tǒng)計(jì)學(xué)意義，可以得出研究區(qū)的泥石流易發(fā)性模型:

Risk=-1.119 6+2.591 Slide+0.000 3Elevation+0.000 1 Near-Dis+0.015 5Lithology-0.062 2LUC+0.017 8Slope-0.204 9NDVI

與線性回歸一樣，擬合Logistic回歸模型時(shí)也對(duì)因子的多重共線性很敏感，王濟(jì)川等提出采用線性回歸模型的形式檢驗(yàn)多重共線性［21］。表3為對(duì)7個(gè)因子的多重共線性診斷結(jié)果。

表3 共線性診斷結(jié)果

VIF為方差膨脹因子，容差為VIF的倒數(shù)。當(dāng)VIF＞5時(shí)表明該因子與其他因子存在多重共線性，表3中VIF的最大值為1.468，可以認(rèn)為因子不存在明顯多重共線性，即上述模型多重共線性對(duì)回歸系數(shù)的影響有限。

2.2 泥石流危險(xiǎn)性分區(qū)及檢測(cè) 為了檢驗(yàn)?zāi)Ｐ偷男Ч?，采用受試者工作特征曲線(Receiver Operating Characteristic Curve，簡(jiǎn)稱ROC曲線)來反映(圖4)。ROC曲線是以虛報(bào)概率為橫軸，擊中概率為縱軸，由被測(cè)試的樣本在特定閾值下根據(jù)不同判斷標(biāo)準(zhǔn)得出的不同結(jié)果繪出曲線。AUC(Area under the Curve)為ROC曲線之下面積，可以用作衡量模型精度的標(biāo)準(zhǔn)，其值越接近于1越好。

圖4 ROC曲線

上述模型的ROC曲線AUC面積約為0.728，表明模型預(yù)測(cè)準(zhǔn)確率基本準(zhǔn)確。

對(duì)整個(gè)研究區(qū)使用上述模型，采用自然間斷點(diǎn)分級(jí)法將滑坡危險(xiǎn)度分為高風(fēng)險(xiǎn)、中風(fēng)險(xiǎn)、較低風(fēng)險(xiǎn)、低風(fēng)險(xiǎn)4個(gè)級(jí)別，得到結(jié)果見圖5。

由圖5可以分析出高危險(xiǎn)區(qū)和中危險(xiǎn)區(qū)主要分布在白龍江干流和二級(jí)支流的沿岸。與實(shí)地基本狀況相符，其中以舟曲縣南部所占比例較高。較低危險(xiǎn)區(qū)分布比較散亂，主要分布于高度危險(xiǎn)區(qū)的周圍;輕度危險(xiǎn)區(qū)主要分布在武都南部及文縣西北部區(qū)域;其中高危險(xiǎn)區(qū)達(dá)到研究區(qū)總面積的13.7%。中危險(xiǎn)區(qū)達(dá)到24.4%，真實(shí)發(fā)生泥石流落入高危險(xiǎn)區(qū)和中危險(xiǎn)區(qū)的面積占到已發(fā)生泥石流的溝谷面積的73.3%。說明應(yīng)用Logistic回歸模型進(jìn)行泥石流危險(xiǎn)度分區(qū)可以達(dá)到較好的效果。

圖5 白龍江流域武都段危險(xiǎn)度分區(qū)示意

3 結(jié)論與討論

隴南地區(qū)由于其特殊的地形、降雨等因素，是我國(guó)泥石流、滑坡、崩塌等地質(zhì)災(zāi)害四大高發(fā)區(qū)之一。泥石流災(zāi)害的風(fēng)險(xiǎn)區(qū)劃能夠?yàn)楫?dāng)?shù)卣姆罏?zāi)減災(zāi)工作提供科學(xué)決策的依據(jù)。從模型的評(píng)價(jià)結(jié)果可以看出，該研究所選取的致災(zāi)因子能夠較好地反映研究區(qū)的地形地貌特征、地質(zhì)環(huán)境條件。通過Logistic回歸方法擬合了泥石流災(zāi)害與各種致災(zāi)因子的關(guān)系，獲取的評(píng)估結(jié)果與實(shí)際情況基本相符，能夠?yàn)楫?dāng)?shù)胤乐文嗍鳛?zāi)害提供科學(xué)的依據(jù)。但由于降雨數(shù)據(jù)的缺乏，在研究中未能引入降雨因子，使得模型的解釋效果并不完美，因此在今后的研究中，為了建立更完整、解釋能力更好的模型，還需要引入其他致災(zāi)因子，并且從泥石流的成因機(jī)理方面進(jìn)行深入研究。

［1］孟興民，陳冠，郭鵬，等.白龍江流域滑坡泥石流災(zāi)害研究進(jìn)展與展望［J］.海洋地質(zhì)與第四紀(jì)地質(zhì)，2013(4):1-15.

［2］侯蘭功，崔鵬.單溝泥石流災(zāi)害危險(xiǎn)性評(píng)價(jià)研究［J］.水土保持研究，2004(2):125-128.

［3］唐川，朱大奎.基于GIS技術(shù)的泥石流風(fēng)險(xiǎn)評(píng)價(jià)研究［J］.地理科學(xué)，2002(3):300-304.

［4］VAN WESTEN C J，TERLIEN M J T.An approach towards deterministic landslide hazard analysis in GIS.A case study from Manizales(Colombia)［J］.Earth Surface Processes and Landforms，1996，21(9):853-868.

［5］LEE S，CHOI J，MIN K.Probabilistic landslide hazard mapping using GIS and remote sensing data at Boun，Korea［J］.International Journal of Remote Sensing，2004，25(11):2037-2052.

［6］唐川，朱靜.瀾滄江中下游滑坡泥石流分布規(guī)律與危險(xiǎn)區(qū)劃［J］.地理學(xué)報(bào)，1999(S1):84-92.

［7］沈芳，黃潤(rùn)秋，苗放，等.區(qū)域地質(zhì)環(huán)境評(píng)價(jià)與災(zāi)害預(yù)測(cè)的GIS技術(shù)［J］.山地學(xué)報(bào)，1999(4):338-342.

［8］李雪平.基于GIS的區(qū)域斜坡穩(wěn)定性評(píng)價(jià)Logistic回歸模型研究［D］.武漢:中國(guó)地質(zhì)大學(xué)，2005.

［9］李軍，周成虎.基于柵格GIS滑坡風(fēng)險(xiǎn)評(píng)價(jià)方法中格網(wǎng)大小選取分析［J］.遙感學(xué)報(bào)，2003(2):86-92，161.

［10］GORSEVSKI P V，GESSLER P L，F(xiàn)OLTZ R B.Spatial prediction of landslide hazard using logistic regression and GIS［C］//4th International Conference on Integrating GIS and Environmental Modeling:Proceedings Propects and Research Needs.Banff，Alberta，Canada，2000.

［11］OHLMACHER G C，DAVIS J C.Using multiple logistic regression and GIS technology to predict landslide hazard in northeast Kansas，USA［J］.Engineering Geology，2003，69(3):331-343.

［12］AYALEW L，YAMAGISHI H.The application of GIS-based logistic regression for landslide susceptibility mapping in the Kakuda-Yahiko Mountains，Central Japan［J］.Geomorphology，2005，65(1):15-31.

［13］潘赟，叢威青，潘懋.基于GIS的遼寧省岫巖縣泥石流災(zāi)害危險(xiǎn)性區(qū)劃研究［J］.北京大學(xué)學(xué)報(bào):自然科學(xué)版，2010(4):601-606.

［14］叢威青，潘懋，李鐵鋒，等.基于GIS的滑坡、泥石流災(zāi)害危險(xiǎn)性區(qū)劃關(guān)鍵問題研究［J］.地學(xué)前緣，2006(1):185-190.

［15］苗雨.基于GIS和Logistic模型的地質(zhì)災(zāi)害危險(xiǎn)性區(qū)劃研究［D］.西安:長(zhǎng)安大學(xué)，2010.

［16］李鴻璉，曾思偉.甘肅泥石流［M］.北京:人民交通出版杜，1982.

［17］劉希林，張松林，唐川.溝谷泥石流危險(xiǎn)度評(píng)價(jià)研究［J］.水土保持學(xué)報(bào)，1993(2):20-25.

［18］邢釗.基于信息熵與AHP模型的白龍江流域泥石流危險(xiǎn)性評(píng)價(jià)［D］.蘭州:蘭州大學(xué)，2012.

［19］李雪平，唐輝明.基于GIS的分組數(shù)據(jù)Logistic模型在斜坡穩(wěn)定性評(píng)價(jià)中的應(yīng)用［J］.吉林大學(xué)學(xué)報(bào):地球科學(xué)版，2005(3):361-365.

［20］王衛(wèi)東，鐘晟.基于GIS的Logistic回歸模型在地質(zhì)災(zāi)害危險(xiǎn)性區(qū)劃中的應(yīng)用［J］.工程勘察，2009(11):5-10.

［21］王濟(jì)川，郭志剛.Logistic回歸模型:方法與應(yīng)用［M］.北京:高等教育出版社，2001.