陳羽璇， 楊勤科,2?， 劉寶元， 黃晨璐， 王春梅,2， 龐國偉,2

(1.西北大學(xué)城市與環(huán)境學(xué)院，710127，西安；2.旱區(qū)生態(tài)水文與災(zāi)害防治國家林業(yè)局重點(diǎn)實(shí)驗(yàn)室，710127，西安；3.中國科學(xué)院 水利部 水土保持研究所，712100，陜西楊凌；4.北京師范大學(xué)地理學(xué)院，100875，北京)

土壤侵蝕已嚴(yán)重威脅世界許多地區(qū)的社會(huì)經(jīng)濟(jì)發(fā)展[1]，對(duì)其治理受到廣泛關(guān)注[2]，土壤侵蝕的治理需以區(qū)域土壤侵蝕調(diào)查制圖為數(shù)據(jù)基礎(chǔ)。自20世紀(jì)90年代以來，國內(nèi)外學(xué)者對(duì)區(qū)域土壤侵蝕評(píng)價(jià)進(jìn)行一系列研究[2]。Lu等[3]和Teng等[4]先后用RUSLE對(duì)澳洲土壤侵蝕進(jìn)行評(píng)價(jià)，Borrelli等[5]在全球尺度上進(jìn)行土壤侵蝕定量分析，劉寶元等[6]基于抽樣調(diào)查方法進(jìn)行全國土壤水力侵蝕普查制圖，楊志成等[7]在抽樣調(diào)查基礎(chǔ)上進(jìn)行縣域和省域土壤侵蝕評(píng)價(jià)。目前已初步形成對(duì)坡面和小流域尺度土壤侵蝕速率的精度評(píng)價(jià)方法和對(duì)區(qū)域尺度侵蝕危險(xiǎn)性評(píng)價(jià)方法，但區(qū)域尺度上實(shí)現(xiàn)土壤侵蝕速率的調(diào)查和制圖方法有待研究。

珠江是中國第三大河流，珠江三角洲是我國乃至世界最具活力的經(jīng)濟(jì)區(qū)之一。由于地處南亞熱帶，降水量大、土層薄，加上人類活動(dòng)干擾，珠江流域土壤侵蝕愈發(fā)嚴(yán)峻[8-9]。自1990年以來，研究者對(duì)珠江流域土壤侵蝕及石漠化等進(jìn)行一些研究[9-13]，但目前的研究對(duì)流域土壤侵蝕空間分異特征和主控因子的認(rèn)識(shí)不足，與珠江三角洲繁榮的經(jīng)濟(jì)不相適應(yīng)，不能有效支撐珠江流域水土保持規(guī)劃和決策。為此，迫切需要利用現(xiàn)代科學(xué)技術(shù)手段，對(duì)該區(qū)土壤侵蝕進(jìn)行系統(tǒng)研究。

筆者基于抽樣調(diào)查單元獲取的相關(guān)信息和區(qū)域土壤侵蝕因子數(shù)據(jù)，結(jié)合CSLE模型[14]與GIS技術(shù)，采用地圖代數(shù)和空間插值2種方法對(duì)研究區(qū)土壤侵蝕速率制圖并進(jìn)行土壤侵蝕分析，以期為該地區(qū)生態(tài)環(huán)境治理及水土保持工作實(shí)施提供可靠的科學(xué)依據(jù)，同時(shí)促進(jìn)我國亞熱帶地區(qū)土壤侵蝕的研究。

1 研究區(qū)概況

研究區(qū)包括珠江流域及海南省等相鄰區(qū)(圖1)，總面積約59.14萬km2。珠江流域地貌類型復(fù)雜，山地、丘陵占總面積94.4%，平原面積小且分散；西江中上游地區(qū)以石灰?guī)r山地為主，東江流域以花崗巖丘陵為主，北江中上游地區(qū)干流西側(cè)為石灰?guī)r山地、東側(cè)為紅巖盆地、下游以花崗巖丘陵為主。流域地處南亞熱帶，年平均氣溫在14～22 ℃之間；流域內(nèi)年平均降水量1 470 mm，自然植被以常綠闊葉林為主，主要土地利用類型為林地、草地和耕地。主要土壤類型有赤紅壤、紅壤、黃壤和水稻土等。

圖1 研究區(qū)位置圖Fig.1 Location map of the study area

本研究的數(shù)據(jù)主要包括2部分，一是土壤侵蝕抽樣調(diào)查數(shù)據(jù)，研究區(qū)(包括陸地部分55 km緩沖區(qū))共布設(shè)抽樣調(diào)查單元1 695個(gè)(圖2)，該數(shù)據(jù)基于高分辨率影像解譯而來，主要包含土地利用類型和水土保持措施2部分解譯內(nèi)容，精度優(yōu)于1∶5 000[15]；二是用于計(jì)算區(qū)域土壤侵蝕因子專題數(shù)據(jù)。數(shù)據(jù)及其來源詳見表1。

圖2 抽樣單元分布圖Fig.2 Distribution of sampling units

表1 基礎(chǔ)數(shù)據(jù)一覽表Tab.1 List of basic data

2 研究方法

2.1 土壤侵蝕因子計(jì)算

筆者使用CSLE模型[14]計(jì)算土壤侵蝕速率。所用參數(shù)為全球降雨侵蝕力因子(R因子)、中國土壤可蝕性因子(K因子)、全球L坡長因子與S坡度因子(LS因子)。這3個(gè)因子的計(jì)算方法詳見第一次全國水利普查[17]和相關(guān)文獻(xiàn)，本文主要介紹水土保持措施因子的計(jì)算。

植被覆蓋與生物措施因子(B因子)：基于Borrelli等的[5]C因子計(jì)算方法，根據(jù)研究區(qū)實(shí)際情況改進(jìn)。對(duì)農(nóng)業(yè)用地C因子賦值為1；對(duì)非農(nóng)業(yè)用地C因子，結(jié)合GLC30土地利用分類，定義相應(yīng)土地利用C因子值的約束范圍。

CNA=MINCF+(MAXCF-MINCF)(1-TC)

(1)

式中：CNA為林地、灌木覆蓋C因子值；MINCF=0.000 1、MAXCF=0.003；TC為林地覆蓋的比例。

CP=MINC+(MAXC-MINC)·(1-NVS)

(2)

式中：Cp為非林地覆蓋(草地覆蓋等)C因子值；MINC=0.01、MAXC=0.05；NVS為非林地(草地等)覆蓋的比例。

裸土按照式(2)計(jì)算，MINC、MAXC取值為0.1、0.5，其計(jì)算結(jié)果分辨率為250 m，通過重采樣得到30 m分辨率數(shù)據(jù)。最后基于GLC30土地利用分類合并CNA、Cp、裸土C和耕地C，得到B因子圖層。

工程措施因子(E因子)：利用抽樣單元土地利用數(shù)據(jù)，根據(jù)《區(qū)域水土流失動(dòng)態(tài)監(jiān)測(cè)技術(shù)規(guī)定》中計(jì)算方法編寫程序計(jì)算得到抽樣單元E因子均值，再對(duì)E因子均值的點(diǎn)數(shù)據(jù)進(jìn)行反距離權(quán)重空間插值得到E因子圖層。

耕作措施因子(T因子)：根據(jù)《區(qū)域水土流失動(dòng)態(tài)監(jiān)測(cè)技術(shù)規(guī)定》，對(duì)抽樣單元內(nèi)耕地的輪作措施屬性字段賦值，非耕地賦值1，并基于Borrelli等[5]對(duì)農(nóng)業(yè)中水田的處理方法，將水田賦值為0.15，最終得到T因子圖層。

2.2 流域土壤侵蝕速率圖制作

筆者采用2種方法：空間插值法和地圖代數(shù)法。空間插值法是基于抽樣單元的侵蝕速率圖[15]，用地統(tǒng)計(jì)協(xié)同克里金插值方法，以R、K、B因子作為協(xié)變量，對(duì)耕、林、草、灌木及裸地分別插值生成不同地類圖，并將水體、濕地、不透水面設(shè)置為0，最終得到土壤侵蝕速率專題圖。地圖代數(shù)法是根據(jù)準(zhǔn)備好的各因子圖層數(shù)據(jù)，直接進(jìn)行各因子相乘得到土壤侵蝕圖，考慮到LS、K因子分辨率均為30 m，本研究的運(yùn)算將分辨率設(shè)置為30 m。

2.3 主控性因子分析方法

基于地理探測(cè)器進(jìn)行土壤侵蝕主控因子分析。地理探測(cè)器是一種通過探測(cè)要素的空間分層異質(zhì)性及揭示背后驅(qū)動(dòng)力的一種統(tǒng)計(jì)學(xué)方法[18]，其大小由q來衡量，q表達(dá)式如下：

(3)

(4)

3 結(jié)果與分析

3.1 土壤侵蝕影響因子分布特征

影響土壤侵蝕速率各因子中(圖3)，R因子總體上呈從東至西遞減的規(guī)律，其值域?yàn)? 786.81～16 324.18 MJ·mm/(hm2·h·a)，平均值為6 813.56 MJ·mm/(hm2·h·a)。K因子值域?yàn)?～0.012 5 t·hm2·h/(MJ·mm·hm2)，均值為0.004 8 t·hm2·h/(MJ·mm·hm2)，有從東到西遞增趨勢(shì)。LS因子值域?yàn)?.04～21.24，平均值為8.70，高值見于云貴高原、東部山地和海南五指山等地。R、LS因子值普遍較高，是土壤侵蝕的主要誘發(fā)因子。B因子、E因子和T因子，值域分布范圍依次為0～1、0.15～1、0.15～1，其均值依次為0.21、0.91、0.87。B因子值較低，而E和T因子值較高，研究區(qū)植被覆蓋度較好，是土壤侵蝕主要抑制因子。

圖3 CSLE各因子空間分布圖Fig.3 Spatial distribution of CSLE factors

3.2 2種計(jì)算方法下土壤侵蝕的對(duì)比分析

3.2.1 空間分布特征對(duì)比分析 在水利部批準(zhǔn)SL190—2007《土壤侵蝕分類分級(jí)標(biāo)準(zhǔn)》[19]的基礎(chǔ)上進(jìn)行細(xì)分，編制了土壤侵蝕分級(jí)圖(圖4和圖5)，研究區(qū)土壤侵蝕均主要集中分布在研究?jī)?nèi)的貴州省及云南省、廣西中部和廣東省沿海區(qū)，其中強(qiáng)烈和極強(qiáng)烈侵蝕區(qū)域范圍較小且分布在較零散的坡耕地上。2張圖上土壤侵蝕速率較高的地區(qū)與土地利用圖的耕地分布較吻合，而林、草地的侵蝕則相對(duì)較弱。從微觀角度看，2圖在細(xì)節(jié)紋理結(jié)構(gòu)上略有不同?；诘貓D代數(shù)法的結(jié)果局地變化比較明顯，空間插值結(jié)果整體較為平滑，主要原因在于插值所用的侵蝕速率為抽樣單元均值，且插值過程中有平滑作用。

圖4 地圖代數(shù)法土壤侵蝕速率圖Fig.4 Soil erosion rate map based on map algebra method

圖5 空間插值法土壤侵蝕速率圖Fig.5 Soil erosion rate map based on spatial interpolation method

3.2.2 統(tǒng)計(jì)特征對(duì)比分析 采用地圖代數(shù)和空間插值法制圖的土壤侵蝕速率均值分別為791.78 t/(km2·a)和615.37 t/(km2·a)，兩者值較接近。地圖代數(shù)法計(jì)算的土壤侵蝕速率均值較高，可能與水保措施難以全面納入計(jì)算有關(guān)。筆者計(jì)算的土壤侵蝕速率值與前人研究[20-22]對(duì)比(表2)，比較接近，可以認(rèn)為本研究的結(jié)果較為可信。對(duì)比研究區(qū)2種方法的主要土地利用類型的土壤侵蝕速率結(jié)果(表3)發(fā)現(xiàn)，土壤侵蝕速率值較為接近。空間插值法制圖結(jié)果，耕地、林地的土壤侵蝕速率均值略小，這是由于地圖代數(shù)法是從整個(gè)研究區(qū)的計(jì)算，分辨率相對(duì)較粗，難以獲取詳細(xì)的水土保持措施，空間插值法從抽樣單元小流域可計(jì)算得到準(zhǔn)確的水土保持措施值，計(jì)算的土壤侵蝕速率值較小。

表2 土壤侵蝕速率對(duì)比參照表Tab.2 Soil erosion rates of this study and extracted from references t/(km2·a)

表3 主要土地利用土壤侵蝕速率對(duì)比表Tab.3 Average soil erosion rates of land use types t/(km2·a)

3.2.3 土壤侵蝕強(qiáng)度對(duì)比分析 根據(jù)SL190—2007《土壤侵蝕分類分級(jí)標(biāo)準(zhǔn)》[19]對(duì)土壤侵蝕強(qiáng)度進(jìn)行統(tǒng)計(jì)(表4)，2種方法計(jì)算的水土流失面積(輕度及其以上侵蝕)差值為10.17%，空間插值法略大；土壤侵蝕均主要集中于微、輕度侵蝕；但其在微、輕度侵蝕頻率分布略有差異，這可能與空間插值時(shí)耕地土壤侵蝕速率值較集中有關(guān)；在中度及以上的侵蝕強(qiáng)度上，土壤侵蝕速率結(jié)果頻率相差不大。

表4 空間插值法與地圖代數(shù)法土壤侵蝕強(qiáng)度統(tǒng)計(jì)Tab.4 Soil erosion intensity statistics by spatialinterpolation and map algebra %

總體上看，2種方法計(jì)算的土壤侵蝕速率結(jié)果均可表達(dá)研究區(qū)土壤侵蝕狀況，考慮到土壤侵蝕防控的需要，優(yōu)先選用空間插值法。

3.3 土壤侵蝕主控因子分析

基于地理探測(cè)器的分析結(jié)果(表5)表明，整個(gè)研究區(qū)和水土保持區(qū)劃的各個(gè)區(qū)，各因子影響程度從大到小的排序基本為土地利用、水保措施(特別是B因子)、地形、土壤和降水。除浙閔山地丘陵區(qū)外，土地利用方式為土壤侵蝕的主控因子，影響因子q值均達(dá)到47%以上；其次，B因子對(duì)土壤侵蝕速率及其空間分布的影響較大；T因子的影響排第三，華南沿海丘陵臺(tái)地區(qū)和海南及南海諸島丘陵區(qū)的q值達(dá)到40%以上，可能與該區(qū)復(fù)雜的耕作方式有關(guān)，其原因有待進(jìn)一步探究分析。K因子和R因子的影響較小，這與研究區(qū)以林草植被為主、植被覆蓋普遍較高有關(guān)。

表5 各影響因子q值統(tǒng)計(jì)表Tab.5 Statistical table of q values of each impact factor

4 結(jié)論

1)研究區(qū)土壤侵蝕特征環(huán)境是降雨侵蝕力較大、坡度較陡及植被覆蓋度較高，前2項(xiàng)是土壤侵蝕的主要誘發(fā)因子，后者則是土壤侵蝕的主要抑制因子。

2)本文通過地圖代數(shù)和空間插值2種方法制圖得到土壤侵蝕速率均值分別為791.78 t/(km2·a)、615.37 t/(km2·a)，研究區(qū)土壤侵蝕廣泛分布，坡耕地是土壤侵蝕主要的發(fā)生部位。

3)土地利用類型為影響土壤侵蝕速率的主控因素，B因子對(duì)土壤侵蝕速率及其空間分布影響較大，總體上講，由于研究區(qū)植被覆蓋度普遍較高，因而土壤和降雨的影響較小。調(diào)整土地利用結(jié)構(gòu)、不斷優(yōu)化植被的水土保持生態(tài)功能是值得注意的問題。

4)地圖代數(shù)和空間插值兩種方法制圖結(jié)果相比，空間插值法的結(jié)果更接近土壤侵蝕的實(shí)際情況，可作為區(qū)域土壤侵蝕制圖首選方法，但兩者之間的差異和相互結(jié)合有待討論。

中國科學(xué)院水利部水土保持研究所張曉萍研究員、水利部水土保持監(jiān)測(cè)中心王愛娟高工、珠江水利委員會(huì)水土保持監(jiān)測(cè)中心金平偉主任對(duì)本研究給予多方支持，一并致謝。