武翠翠， 王世杰*

(1.蘭州交通大學(xué)測(cè)繪與地理信息學(xué)院， 蘭州 730070； 2.地理國(guó)情監(jiān)測(cè)技術(shù)應(yīng)用國(guó)家地方聯(lián)合工程研究中心，蘭州 730070； 3.甘肅省地理國(guó)情監(jiān)測(cè)工程實(shí)驗(yàn)室， 蘭州 730070)

土壤或者土體是由于在外營(yíng)力(風(fēng)力、重力、水力、凍融)的作用下產(chǎn)生剝蝕、沖刷、堆積和遷移等現(xiàn)象即為土壤侵蝕。近年來(lái)，利用地理信息系統(tǒng)(geographic information system，GIS)技術(shù)研究流域生態(tài)敏感性已成為測(cè)繪和生態(tài)學(xué)的交叉研究熱點(diǎn)。用GIS技術(shù)單因子的疊加進(jìn)行土壤侵蝕敏感性評(píng)價(jià)是學(xué)者們普遍采用的方法。如李益敏等[1]以星云湖流域?yàn)槔?，使用GIS技術(shù)，選取體現(xiàn)流域特點(diǎn)的估算模型，統(tǒng)計(jì)分析流域的土壤侵蝕空間分布特征。徐萍[2]基于改進(jìn)的土壤侵蝕模型(the revised universal soil loss equation，RULSE)對(duì)龍江縣水土流失進(jìn)行敏感性評(píng)價(jià)。此外，相關(guān)學(xué)者以慶安市、沈陽(yáng)市、山東省等作為研究區(qū)域進(jìn)行了生態(tài)敏感性評(píng)價(jià)研究[3-5]。但小尺度流域的土壤侵蝕敏感性的評(píng)價(jià)尚且不多，如徐劍春等[6]以通渭縣牛谷河流域?yàn)槔M(jìn)行小流域尺度水土流失敏感性評(píng)價(jià)。魏偉等[7]基于遙感指數(shù)的干旱內(nèi)陸河流域土地生態(tài)敏感性時(shí)空演變特征進(jìn)行研究。有文獻(xiàn)指出，對(duì)現(xiàn)有流域生態(tài)治理前后，生態(tài)敏感性現(xiàn)狀評(píng)價(jià)多集中在某一時(shí)間節(jié)點(diǎn)上，存在一定的主觀性，缺乏系統(tǒng)性評(píng)價(jià)指數(shù)[7]。國(guó)外Wischmeier和Smith兩位學(xué)者提出傳統(tǒng)的土壤侵蝕模型(the universal soil loss equation，USLE)是進(jìn)行區(qū)域土壤侵蝕敏感性評(píng)價(jià)普遍采用的模型[8-13]，評(píng)價(jià)方法多樣。RUSLE模型是學(xué)者們?cè)赨SLE模型的基礎(chǔ)上改進(jìn)的模型，RUSLE模型比USLE所用的數(shù)據(jù)量增加，對(duì)USLE中的許多錯(cuò)誤進(jìn)行糾正，填補(bǔ)空白的原始數(shù)據(jù)，模擬不同的系統(tǒng)和增加了模型的靈活性[14]。相關(guān)研究表明[15]，用RUSLE模型評(píng)價(jià)土壤侵蝕敏感性是目前最有效的方法。但是，RUSLE模型也存在一些局限性，高原地區(qū)不適合該緩坡地模型，中國(guó)黃土高原土壤侵蝕主要以降雨侵蝕為主，而該模型是年降雨的侵蝕產(chǎn)沙模型，不適合中國(guó)許多地區(qū)使用，而且對(duì)于土壤可蝕性因子K的計(jì)算一直存在爭(zhēng)議。黃河流域土壤質(zhì)地疏松、坡地陡峭，許多風(fēng)蝕、水蝕等強(qiáng)力搬運(yùn)作用，加之人口密集，植被遭到大量破壞，夏季由于洪水決堤，土壤受到破壞，土壤侵蝕嚴(yán)重。研究區(qū)域降水時(shí)間集中，汛期導(dǎo)致土壤侵蝕更為嚴(yán)重，泥沙入河，使得災(zāi)害頻發(fā)，加之人為因素的活動(dòng)，生態(tài)系統(tǒng)極度敏感。針對(duì)上述問(wèn)題，現(xiàn)結(jié)合研究區(qū)環(huán)境現(xiàn)狀，加入評(píng)價(jià)因子權(quán)重，對(duì)RUSLE模型進(jìn)行改進(jìn)。采用研究區(qū)遙感數(shù)據(jù)，結(jié)合GIS空間分析方法和層次分析法，利用改進(jìn)后的加權(quán)RUSLE模型，對(duì)黃河流域蘭州段進(jìn)行土壤侵蝕敏感性評(píng)價(jià)。

圖1 研究區(qū)地理位置Fig.1 The geographic location of the study area

1 研究區(qū)概況

研究區(qū)為黃河流經(jīng)的蘭州市轄區(qū)，如圖1所示。區(qū)域徑流長(zhǎng)度約150.7 km，蘭州市面積為13 100 km2。東鄰定西市，北邊與武威市和白銀相鄰，南邊與臨夏回族自治州接壤。平均海拔1 530～1 580 m，年蒸發(fā)量為1 676 mm，年均降水量324.8 mm，年均氣溫10.3 ℃，是溫帶大陸性氣候，夏季降雨量多且集中暴雨。蘭州市地處黃河由青藏高原流入黃土高原的轉(zhuǎn)折點(diǎn)上，是串珠狀黃河河谷盆地，地形為南北兩山夾峙，地貌為典型的黃土地貌，是中國(guó)黃土沉積厚度最厚的地區(qū)?；意}土為土壤中的主要成分，成土母質(zhì)是黃土，有機(jī)質(zhì)含量低，植被稀疏，土壤無(wú)法保水保肥，由于水量不穩(wěn)定且雨季集中，導(dǎo)致流域土壤侵蝕嚴(yán)重。

2 研究方法

2.1 數(shù)據(jù)來(lái)源

數(shù)據(jù)源蘭州市衛(wèi)星影像下載自美國(guó)地質(zhì)調(diào)查局(united states geological survey，USGS)，軌道號(hào)為130/35和131/35，成像時(shí)間為2017年12月，云量10%以下。使用該數(shù)據(jù)提取歸一化植被指數(shù)(normalized difference vegetation index，NDVI)，在ENVI軟件中完成圖像的預(yù)處理，用來(lái)計(jì)算植被覆蓋因子。30GDEMV2數(shù)字高程模型(digital elevation model，DEM)下載自地理空間數(shù)據(jù)云，5景DEM數(shù)據(jù)拼接后覆蓋研究區(qū)域，用于計(jì)算坡度、坡長(zhǎng)因子。收集2017年黃河流域蘭州段及其周圍站點(diǎn)的月平均降水量數(shù)據(jù)，用于計(jì)算降雨侵蝕力因子。水土保持措施因子是使用中科院2017年土地利用數(shù)據(jù)計(jì)算獲取。根據(jù)世界土壤數(shù)據(jù)庫(kù)(harmonized world soil database，HWSD)數(shù)據(jù)和中科院收集的土壤類型空間分布數(shù)據(jù)，計(jì)算砂土、黏土、壤土、有機(jī)質(zhì)含量，進(jìn)而計(jì)算土壤可蝕性因子。利用GIS空間分析技術(shù)完成土壤侵蝕敏感性的評(píng)價(jià)。

2.2 改進(jìn)的土壤侵蝕敏感性評(píng)價(jià)模型

2.2.1 RUSLE模型

學(xué)術(shù)界普遍應(yīng)用的土壤侵蝕模型(RUSLE)表達(dá)式為

A=RKLSPC

(1)

式(1)中：A為單位面積年平均土壤侵蝕量，t/hm2；R為降雨侵蝕力因子，MJ·mm/(hm2·h·a)，在數(shù)值上，A與R相等；K為土壤可蝕性因子，t·h/(MJ·mm)；L為坡長(zhǎng)因子；S為坡度因子；P為水土保持措施因子；C為植被覆蓋因子。

2.2.2 改進(jìn)的加權(quán)綜合評(píng)價(jià)模型

以RUSLE 模型為基礎(chǔ)，結(jié)合研究區(qū)的實(shí)際土壤狀況，選取6個(gè)因子作為評(píng)價(jià)因子，加以權(quán)重改進(jìn)RUSLE模型。改進(jìn)后的加權(quán)綜合評(píng)價(jià)模型為

A1=PiRPiKPiLPiSPiPPiC

(2)

(3)

式中：A1為單位面積年平均土壤侵蝕量，t/hm2；Sj為土壤侵蝕敏感性指數(shù)；Pi為每個(gè)影響因子對(duì)應(yīng)的權(quán)重，綜合考慮各因子對(duì)土壤侵蝕的影響程度，應(yīng)用層次分析法確定權(quán)重。

2.2.3 因子的選取與權(quán)重確定方法

根據(jù)黃河流域蘭州段實(shí)際生態(tài)環(huán)境特點(diǎn)，構(gòu)建土壤侵蝕敏感性的合理評(píng)價(jià)指標(biāo)體系是關(guān)鍵。在各個(gè)因子的選取方法上遵循實(shí)用性、差異性、主導(dǎo)性、綜合性和所需指標(biāo)因子的可獲取性等原則。結(jié)合前人的指標(biāo)體系經(jīng)驗(yàn)，根據(jù)研究區(qū)土壤環(huán)境情況及結(jié)構(gòu)功能特性，綜合考慮自然因素和人為因素的共同影響。確定的自然因素指標(biāo)體系有坡度、坡長(zhǎng)、土壤可蝕性、水土保持措施、植被覆蓋度和降雨侵蝕力等因子，建立了相對(duì)完善的土壤侵蝕敏感性評(píng)價(jià)指標(biāo)體系，如表1、表2所示。

權(quán)重的確定采用層次分析法。提出此方法的是美國(guó)運(yùn)籌學(xué)家A L Saaty。層次分析法適用于確定權(quán)重較復(fù)雜的模糊綜合評(píng)價(jià)系統(tǒng)，既有非定量事件的定性分析，也有較精確的定量分析。步驟是建立黃河流域蘭州段土壤侵蝕敏感性評(píng)價(jià)模型層次結(jié)構(gòu)遞階關(guān)系，構(gòu)建各層次評(píng)價(jià)指標(biāo)兩兩相比的判斷矩陣，計(jì)算各層元素的權(quán)重值，對(duì)判斷矩陣進(jìn)行單層次排序與一致性檢驗(yàn)，計(jì)算各指標(biāo)因子的綜合權(quán)重值，采用算術(shù)平均值確定最終的權(quán)重。表3是判斷矩陣標(biāo)度值及其含義，ri、rj為不同指標(biāo)的影響因素，rij為因素i與因素j的重要性之比。結(jié)合5位生態(tài)行業(yè)專家打分與研究區(qū)的實(shí)際生態(tài)狀況，應(yīng)用YAAHP軟件確定的最終權(quán)重如表4所示。其中判斷矩陣一致性檢驗(yàn)結(jié)果是0.073 8<0.1，說(shuō)明權(quán)重正確。

2.3 單因子的識(shí)別和計(jì)算

2.3.1 降雨侵蝕力因子R

降雨量能夠直觀反映土壤侵蝕的強(qiáng)弱度，降雨數(shù)據(jù)與時(shí)間、雨量、強(qiáng)度、類型這些影響因素密切相關(guān)。結(jié)合黃河流域蘭州段降雨量的實(shí)際情況，選取月均降雨量數(shù)據(jù)。采用經(jīng)典簡(jiǎn)易算法來(lái)測(cè)算研究區(qū)域降雨侵蝕力因子R。簡(jiǎn)易法公式為

(4)

式(4)中：Pi為流域內(nèi)月均降雨量，mm；P為流域內(nèi)年均降雨量，mm。

黃河流域蘭州段的降雨監(jiān)測(cè)站點(diǎn)皋蘭站和榆中站分布在黃河流域兩側(cè)。根據(jù)2017年區(qū)域月均降雨量數(shù)據(jù)，應(yīng)用式(4)計(jì)算R值，并用ArcGIS進(jìn)行插值，生成如圖2(a)的R值柵格分布圖。這里根據(jù)實(shí)際情況考慮，用相對(duì)平均誤差最小的反距離權(quán)重插值法模擬降雨量。計(jì)算得出2017年均降雨侵蝕因子，再根據(jù)自然斷點(diǎn)法進(jìn)行重分類。

2.3.2 土壤可蝕性因子K

K值能夠有效地評(píng)價(jià)土壤性質(zhì)和侵蝕的關(guān)系，用來(lái)反映侵蝕營(yíng)力的破壞程度，影響著土壤的侵蝕狀態(tài)。K值的計(jì)算方法較多，但很多方法相關(guān)參數(shù)難以獲取。參考彭雙云等[16]改進(jìn)的簡(jiǎn)易估算方法進(jìn)行計(jì)算。利用蘭州市HWSD數(shù)據(jù)和土壤類型空間分布數(shù)據(jù)，計(jì)算砂粒、粉粒、黏粒和有機(jī)碳百分比含量，得到土壤可蝕性因子K值。計(jì)算公式為

表1 土壤侵蝕敏感性評(píng)價(jià)指標(biāo)體系Table 1 Evaluation index system of soil erosion sensitivity

表2 生態(tài)敏感性等級(jí)劃分Table 2 Ecological sensitivity rating

表3 判斷矩陣標(biāo)度值及其含義Table 3 Value of judgment matrix scale and its meaning

表4 黃河流域蘭州段土壤侵蝕敏感性指標(biāo)權(quán)重Table 4 The weight of soil erosion sensitivity index in Lanzhou section of the Yellow River Basin

(5)

(6)

式中：K為土壤可蝕性因子，t·h/(MJ·mm)；Sa為砂粒含量(粒徑為0.1～2.0 mm)；Si為粉粒含量(粒徑為0.002～0.100 mm)；C1為黏粒含量(粒徑<0.002 mm)；C為有機(jī)碳含量，%。

1∶100萬(wàn)的HWSD是第二次土地調(diào)查獲取的數(shù)據(jù)，用來(lái)計(jì)算有機(jī)碳含量。1∶100萬(wàn)中國(guó)土壤質(zhì)地空間分布數(shù)據(jù)可得到砂粒、粉粒、黏粒的含量。使用ArcGIS得出研究區(qū)的土壤可蝕性因子分布圖[如圖2(b)]。

2.3.3 坡度因子S

坡度因子S是重要的地形指標(biāo)，可以反映坡度與土壤侵蝕間的關(guān)系。一般情況下，徑流量大的坡度陡，坡面的沖刷力就強(qiáng)，區(qū)域內(nèi)土壤受到侵蝕作用就越強(qiáng)，故坡度因子對(duì)土壤侵蝕的影響非常大。坡度因子S采用《第4次全國(guó)土壤侵蝕普查技術(shù)規(guī)程》中規(guī)定的方法計(jì)算，該方法適用于地形起伏較大地域。計(jì)算公式為

(7)

式(7)中：θ為坡度值；S為坡度因子。

研究區(qū)DEM數(shù)據(jù)利用ENVI(the environment for visualizing images)軟件裁剪、鑲嵌等處理，再以處理后的數(shù)字高程模型為基礎(chǔ)，提取坡度專題圖。得到研究區(qū)S因子的柵格圖[圖2(c)]。其中提取坡度因子時(shí)，研究區(qū)緯度為37°，計(jì)算過(guò)程中用到的Z因子(ArcGIS軟件進(jìn)行坡度提取時(shí)一個(gè)表面Z單位中地面x，y單位的數(shù)量)是通過(guò)緯度內(nèi)插得到，Z因子為0.000 011 71。

2.3.4 坡長(zhǎng)因子L

坡長(zhǎng)因子L依據(jù)1978年Wischmeier和Smith兩位學(xué)者提出的經(jīng)典算法計(jì)算，其中，m的取值采用劉寶元等[17]提出的方法，計(jì)算公式為

(8)

式(8)中：λ為水平坡長(zhǎng)；θ為坡度值；L為坡長(zhǎng)因子；m為坡長(zhǎng)指數(shù)。

根據(jù)坡度圖，在柵格計(jì)算器中計(jì)算m值，得到流向圖?；诜直媛蕿?0 m的DEM提取研究區(qū)的坡長(zhǎng)λ：先對(duì)DEM影像進(jìn)行填洼，得到徑流方向(角度)，按照流向計(jì)算單元柵格的λ。應(yīng)用ENVI軟件中的band math工具，利用式(8)計(jì)算坡長(zhǎng)因子L。計(jì)算結(jié)果如圖2(d)所示。

2.3.5 水土保持措施因子P

P因子是指特定水土保持措施下的土壤流失量與未實(shí)施水土保持措施地塊順坡耕作時(shí)的土壤流失量之比。估算水土保持措施因子P的方法有多種。黃杰等[18]比較了各種估計(jì)方法的優(yōu)缺點(diǎn)，可以為P因子的計(jì)算提供參考。綜合考慮研究區(qū)的氣候、環(huán)境、地形等差異，因地制宜地確定P值計(jì)算方法尤為重要。運(yùn)用蘭州段的土地利用類型賦值法來(lái)確定P值，計(jì)算結(jié)果如圖2(e)所示。P因子取值范圍在0～1。一般未采取水土保持措施的土地利用類型將P取值為1，如林地、居民地、未利用土地、園地等。不發(fā)生土壤侵蝕的土地利用類型如水體P值取為0。草地起到較好地防止土壤侵蝕的作用，賦值為0.15。沒(méi)有灌溉或者灌溉不能保證的旱地取值為0.35[19-21]。不同土地利用類型P取值如表5所示。

2.3.6 植被覆蓋度因子C

植被覆蓋因子C是指具有良好的管理措施或者有綠色植被覆蓋地區(qū)與同樣條件下無(wú)植被地表裸露地區(qū)的土壤侵蝕量的比值。取值范圍是[0，1]，C值越大越可以防止土壤侵蝕，起到抑制土壤侵蝕的作用[22]。植被覆蓋可以攔截水流、涵養(yǎng)水源、調(diào)節(jié)水文狀況、固定土壤，可以有效地控制土壤侵蝕。植被覆蓋度指標(biāo)是反映植被群落覆蓋地表狀況的綜合化指標(biāo)。用地理空間數(shù)據(jù)云中的2017年Landsat8-OLI 30 m影像數(shù)據(jù)，經(jīng)過(guò)數(shù)據(jù)預(yù)處理，如輻射定標(biāo)和大氣校正等，應(yīng)用像元二分模型[23]，構(gòu)建植被覆蓋度(vegetation fraction coverage，VFC)，進(jìn)而得出歸一化植被指數(shù)(NDVI)，建立的C值模型是反映植被覆蓋度c和植被覆蓋度C因子的線性關(guān)系，應(yīng)用線性混合像元分解法來(lái)估算C因子，結(jié)果如圖2(f)所示，具體計(jì)算公式為

表5 不同土地利用類型P值Table 5 Pvalues of different land use types

(9)

式(9)中：R為紅外波段反射值；NIR為影像中近紅外波段反射值。

(10)

式(10)中：NDVI為像元的歸一化植被指數(shù)；NDVImin為完全裸土或無(wú)植被覆蓋度區(qū)域的NDVI值(研究區(qū)歸一化植被指數(shù)的最小值)；NDVImax為完全被植被覆蓋的像元NDVI值，即純植被像元NDVI值(研究區(qū)歸一化植被指數(shù)的最大值)。

(11)

式(11)中：c為植被覆蓋度，%;C為植被覆蓋度因子。

3 結(jié)果與分析

基于改進(jìn)的加權(quán)綜合土壤侵蝕評(píng)價(jià)模型，構(gòu)建土壤侵蝕敏感性評(píng)價(jià)指標(biāo)體系。使用ArcGIS軟件，進(jìn)行空間柵格運(yùn)算，綜合利用研究區(qū)遙感影像數(shù)據(jù)、DEM數(shù)據(jù)、土地利用數(shù)據(jù)、土壤類型空間分布數(shù)據(jù)、HWSD數(shù)據(jù)、月平均降雨量數(shù)據(jù)對(duì)黃河流域蘭州段做土壤侵蝕敏感性評(píng)價(jià)。評(píng)價(jià)結(jié)果如下：利用RUSLE模型計(jì)算的單位面積年平均土壤侵蝕量值分布范圍在(1.205 34，5.377 75)，如圖3(a)所示。利用改進(jìn)后加權(quán)綜合評(píng)價(jià)模型的單位面積年平均土壤侵蝕量值范圍在(1.034 699，1)，如圖3(b)所示。圖4所示為新舊模型計(jì)算的土壤侵蝕敏感性各等級(jí)對(duì)應(yīng)的土地面積分布圖，從圖4可以看出，利用改進(jìn)后評(píng)價(jià)模型的土壤侵蝕敏感性程度分布折線圖更加平緩些，更接近研究區(qū)實(shí)際土壤侵蝕和生態(tài)環(huán)境狀況。

圖3 2017年蘭州市土壤侵蝕空間分布示意圖Fig.3 Spatial distribution diagram of soil erosion in Lanzhou in 2017

圖4 兩種模型土壤侵蝕敏感性等級(jí)分布面積Fig.4 The distribution area of soil erosion sensitivity grade of the two models

(1)黃河流域蘭州段，6個(gè)因子對(duì)土壤侵蝕敏感性強(qiáng)弱影響程度如下：降雨侵蝕力因子對(duì)研究區(qū)的土壤侵蝕敏感性較弱，其中土壤侵蝕不敏感區(qū)域占區(qū)域總面積的67.70%，主要分布在永登縣、紅古區(qū)、西固區(qū)、安寧區(qū)。土壤可蝕性因子重度敏感占區(qū)域總面積的80.08%，分布廣泛，地形復(fù)雜，土壤較為貧瘠，容易遭受侵蝕。坡度、坡長(zhǎng)因子和植被覆蓋因子以輕度敏感為主，分別占區(qū)域總面積的34.15%和58.06%，因?yàn)榈匦蔚孛矎?fù)雜，海拔較高，輕度敏感性分布廣泛。水土保持措施因子敏感性大多是中度敏感，占區(qū)域總面積的60.02%，主要集中在永登縣、皋蘭縣和榆中縣，其他地區(qū)也有零星的分散分布。由于西北強(qiáng)勁風(fēng)力，加上地形地貌影響，使得黃河流域蘭州段土壤侵蝕加劇。

(2)黃河流域蘭州段，土壤侵蝕敏感性空間分布特征如下：不敏感、輕度敏感、中度敏感、重度敏感、極敏感區(qū)域面積占比分別為3.62%、8.05%、45.60%、32.10%、10.63%，以中度敏感和重度敏感為主。集中分布在研究區(qū)西北部和中部地區(qū)，有皋蘭縣、永登縣、七里河區(qū)、城關(guān)區(qū)和榆中縣，其中榆中縣50%是極敏感地區(qū)。說(shuō)明研究區(qū)近幾年生態(tài)環(huán)境受到一定破壞，生態(tài)環(huán)境質(zhì)量在逐漸下降，導(dǎo)致較高敏感區(qū)的面積在不斷擴(kuò)大，研究區(qū)土壤侵蝕敏感性整體等級(jí)呈現(xiàn)西北部向東南部遞增的趨勢(shì)。不敏感區(qū)域和輕度敏感區(qū)域的面積在減少，重度敏感區(qū)的面積增加可能與暴雨沖刷導(dǎo)致土壤可蝕性增大和植被覆蓋率降低有關(guān)。

(3)從評(píng)價(jià)結(jié)果可以看出，土壤侵蝕敏感性地區(qū)的面積占比較大。區(qū)域土壤侵蝕敏感性的空間分布和客觀實(shí)際相符，加強(qiáng)敏感區(qū)域的土地利用保護(hù)和管理，是黃河流域蘭州段生態(tài)建設(shè)的重點(diǎn)工作。

4 結(jié)論

研究對(duì)于黃河流域蘭州段土壤侵蝕敏感性進(jìn)行綜合評(píng)價(jià)，識(shí)別出研究區(qū)敏感性程度等級(jí)。對(duì)傳統(tǒng)的土壤侵蝕敏感性評(píng)價(jià)模型做了改進(jìn)。實(shí)驗(yàn)發(fā)現(xiàn)，改進(jìn)后的加權(quán)RUSLE模型的評(píng)價(jià)結(jié)果更符合研究區(qū)實(shí)際情況，可為黃河流域黃土、高寒區(qū)域的土壤侵蝕敏感性評(píng)價(jià)提供技術(shù)參照。但文中利用層次分析法確定權(quán)重時(shí)，專家打分法有一定的局限性，有待尋找更好的權(quán)重確定方法。土壤侵蝕敏感性動(dòng)態(tài)監(jiān)測(cè)也是后續(xù)研究的重要方向。