中國(guó)水蝕區(qū)土壤可蝕性因子更新方法與應(yīng)用

田芷源, 梁 音?, 趙 院, 曹龍熹, 趙 艷,4, 武逸杭,4

(1.土壤與農(nóng)業(yè)可持續(xù)發(fā)展國(guó)家重點(diǎn)實(shí)驗(yàn)室,中國(guó)科學(xué)院南京土壤研究所,210008,南京;2.水利部水土保持監(jiān)測(cè)中心,100053,北京;3.成都理工大學(xué)生態(tài)環(huán)境學(xué)院,610059,成都;4.中國(guó)科學(xué)院大學(xué)南京學(xué)院,211135,南京)

土壤可蝕性因子K是土壤流失方程的關(guān)鍵參數(shù)之一。K因子是侵蝕預(yù)報(bào)模型計(jì)算和水土流失動(dòng)態(tài)監(jiān)測(cè)的基礎(chǔ)數(shù)據(jù)。K值的定義是在標(biāo)準(zhǔn)小區(qū)上單位降雨侵蝕力所引起土壤流失量的多少[1],反映土壤對(duì)侵蝕的敏感程度,其大小與土壤本身性質(zhì)有關(guān)。

我國(guó)第1次水利普查時(shí)獲得了全國(guó)范圍首張K值分布圖[2]。由于實(shí)測(cè)K值使用的徑流小區(qū)建造昂貴,花費(fèi)時(shí)間長(zhǎng)且缺少統(tǒng)一的測(cè)量規(guī)范,能收集到的K值數(shù)據(jù)十分有限,因此利用經(jīng)驗(yàn)公式對(duì)K值進(jìn)行估算[3]。然而,使用的土種志數(shù)據(jù)距今已有近40 a,難以反映目前的土壤情況。為了獲得更準(zhǔn)確的K值結(jié)果,亟需在大范圍內(nèi)對(duì)計(jì)算K值的基礎(chǔ)數(shù)據(jù)進(jìn)行更新。

依靠有限樣點(diǎn)很難滿足越來(lái)越高的制圖精度要求,因此還需探索一種能夠?qū)崿F(xiàn)高空間分辨率的K值點(diǎn)面擴(kuò)展方法。以往研究在利用樣點(diǎn)制作K值分布圖時(shí),常使用土壤圖斑鏈接法及地統(tǒng)計(jì)插值[4]。然而,圖斑鏈接法導(dǎo)致K值在土壤類型邊界發(fā)生突變,且不能反映同一土壤圖斑內(nèi)存在的K值變異性[5]。使用地統(tǒng)計(jì)學(xué)方法對(duì)K值進(jìn)行制圖時(shí),由于采樣密度對(duì)插值精度影響較大,需空間結(jié)構(gòu)因素大于隨機(jī)部分才能得到較準(zhǔn)確的結(jié)果,因此研究范圍大多局限在流域尺度[4]。隨著數(shù)字土壤制圖的發(fā)展,遙感數(shù)據(jù)已被應(yīng)用于K值空間預(yù)測(cè),例如高程作為協(xié)變量被引入克里格插值進(jìn)行K值制圖[6],高程及植被作為空間變量可以提高K值的預(yù)測(cè)精度[7]。機(jī)器學(xué)習(xí)算法中的隨機(jī)森林模型,是一種通過(guò)自助重采樣技術(shù)對(duì)多個(gè)訓(xùn)練集隨機(jī)選取特征變量,分別進(jìn)行回歸建模并平均多個(gè)預(yù)測(cè)值得到最終結(jié)果的方法[8]。該方法已被成功運(yùn)用于預(yù)測(cè)土壤有機(jī)質(zhì)含量、質(zhì)地與土層厚度等關(guān)鍵土壤屬性[9]。

為了實(shí)現(xiàn)全國(guó)土壤可蝕性因子更新目標(biāo),筆者提出一套基于近期匯編的中國(guó)土系志數(shù)據(jù),利用隨機(jī)森林模型進(jìn)行點(diǎn)面擴(kuò)展獲得全國(guó)水蝕區(qū)K值分布圖的方法,解決更新的2個(gè)關(guān)鍵問(wèn)題:第一是如何獲取全國(guó)地區(qū)最新的樣點(diǎn)K值,第二是如何將樣點(diǎn)的更新結(jié)果擴(kuò)展至區(qū)域尺度。本次更新結(jié)果將有助于改善當(dāng)前全國(guó)K值圖數(shù)據(jù)陳舊及制圖方法落后的問(wèn)題。該成果可為土壤侵蝕預(yù)報(bào)模型提供基礎(chǔ)數(shù)據(jù),為年度水土流失動(dòng)態(tài)監(jiān)測(cè)工作提供服務(wù)。

1 材料與方法

1.1 土壤理化數(shù)據(jù)來(lái)源與處理

中國(guó)土系志數(shù)據(jù)是我國(guó)自第2次土壤普查以來(lái)最新且系統(tǒng)性的土壤調(diào)查數(shù)據(jù)[10]。全國(guó)土系調(diào)查項(xiàng)目于2008—2018年基于定量標(biāo)準(zhǔn)和統(tǒng)一分類原則開(kāi)展系統(tǒng)性調(diào)查研究[11]。筆者基于各省份土系調(diào)查資料,收集并整理4 327個(gè)采樣點(diǎn)的土壤剖面數(shù)據(jù),查閱其調(diào)查地點(diǎn)(經(jīng)緯度)和理化性質(zhì)。其中有機(jī)碳質(zhì)量分?jǐn)?shù)通過(guò)除以系數(shù)0.58(采用Van Bemmelen因子,假定土壤有機(jī)質(zhì)的含碳質(zhì)量分?jǐn)?shù)為58%)統(tǒng)一轉(zhuǎn)換為有機(jī)質(zhì)質(zhì)量分?jǐn)?shù)[12]。將數(shù)據(jù)進(jìn)行預(yù)處理,主要采用加權(quán)平均的方法提取各樣點(diǎn)0～30 cm耕層土壤的有機(jī)質(zhì)和機(jī)械組成信息,當(dāng)土體厚度<30 cm時(shí)則以實(shí)際厚度進(jìn)行加權(quán)平均。

1.2 土壤可蝕性K值計(jì)算

采用全國(guó)第1次水利普查時(shí)K值的計(jì)算方法,使用通用土壤流失方程中的諾模圖(Nomo)計(jì)算K值(Nomo-K)[2]。

K=(2.1×10-4M1.14(12-O)+
3.25(S-2)+2.5(P-3))/(100×7.593);

(1)

M=N1(100-N2)

(2)

式中:K為土壤可蝕性,t·hm2·h/(MJ·mm·hm2);O為土壤有機(jī)質(zhì)質(zhì)量分?jǐn)?shù),%;S為土壤結(jié)構(gòu)系數(shù)(計(jì)算土壤粒徑幾何平均直徑后查表所得,量綱為1);P為土壤滲透性等級(jí)(根據(jù)質(zhì)地分類查表所得,量綱為1);M為質(zhì)地指數(shù),由N1和N2計(jì)算所得,%;N1為粒徑在0.002～0.100 mm之間的土壤顆粒質(zhì)量分?jǐn)?shù)比例,%;N2為粒徑< 0.002 mm的土壤顆粒質(zhì)量分?jǐn)?shù),%。

由于全國(guó)土系調(diào)查中的質(zhì)地分布僅包含砂粒(≥0.050～2.000 mm)、粉粒(≥0.002～0.050 mm)和黏粒(<0.002 mm)分級(jí),Nomo公式計(jì)算時(shí)需要把已有土壤顆粒分析結(jié)果轉(zhuǎn)化為極細(xì)砂質(zhì)量分?jǐn)?shù)(≥0.050～0.100 mm)。本研究比較了3種插值方法,分別是自然對(duì)數(shù)線性插值法、三次樣條函數(shù)法、三次樣條函數(shù)結(jié)合自然對(duì)數(shù)。計(jì)算過(guò)程如下:首先將各級(jí)別土壤顆粒質(zhì)量分?jǐn)?shù)累加計(jì)算得到小于某一粒徑的質(zhì)量分?jǐn)?shù)(<0.002 mm、<0.050 mm、<2.000 mm),并計(jì)算該粒徑對(duì)應(yīng)的自然對(duì)數(shù)(ln(0.002)、ln(0.050)和ln(2.000))。自然對(duì)數(shù)線性插值法可建立顆粒累加質(zhì)量分?jǐn)?shù)與粒級(jí)自然對(duì)數(shù)之間的線性關(guān)系;三次樣條函數(shù)法可建立顆粒累加質(zhì)量分?jǐn)?shù)與粒徑之間的三次樣條函數(shù);三次樣條函數(shù)結(jié)合自然對(duì)數(shù)法則可建立顆粒累加質(zhì)量分?jǐn)?shù)與粒徑自然對(duì)數(shù)之間的三次樣條函數(shù)?；谝陨戏椒ǖ玫搅?0.100 mm粒徑所對(duì)應(yīng)的質(zhì)量分?jǐn)?shù),并將<0.050 mm和<0.100 mm粒徑質(zhì)量分?jǐn)?shù)相減得到極細(xì)砂質(zhì)量分?jǐn)?shù),通過(guò)分析極細(xì)砂質(zhì)量分?jǐn)?shù)的合理性和比較插值方程R2確定最優(yōu)的插值方法。

當(dāng)有機(jī)質(zhì)質(zhì)量分?jǐn)?shù)>12%時(shí),使用Nomo公式可能會(huì)導(dǎo)致負(fù)值出現(xiàn),因此對(duì)有機(jī)質(zhì)>12%的樣點(diǎn)采用EPIC模型計(jì)算K值[3]

(3)

Sn1=1-Sa/100。

(4)

式中:Sa為粒徑≥0.050～2.000 mm的砂粒質(zhì)量分?jǐn)?shù),%;Si為≥0.002～0.050 mm 的土壤粉粒質(zhì)量分?jǐn)?shù),%;Cl為粒徑< 0.002 mm的土壤黏粒質(zhì)量分?jǐn)?shù),%;C為土壤有機(jī)碳質(zhì)量分?jǐn)?shù),%,通過(guò)有機(jī)質(zhì)乘以轉(zhuǎn)換系數(shù)58%所得;Sn1為土壤碳酸鈣質(zhì)量分?jǐn)?shù),%。

由于Nomo和EPIC公式的計(jì)算過(guò)程不同,采用EPIC公式得到的K值(EPIC-K)不能直接替代 Nomo-K值。因此,筆者使用線性、指數(shù)、冪函數(shù)和3階多項(xiàng)式函數(shù)分別建立2種公式計(jì)算K值的轉(zhuǎn)換關(guān)系。比較4種擬合方程的R2,除線性方程擬合精度略低以外(0.795 5),其他3種非線性方程的擬合精度十分接近(0.806 4～0.807 2)。從趨勢(shì)來(lái)看,Nomo-K和EPIC-K值之間并非簡(jiǎn)單線性關(guān)系,隨著EPIC-K值的增大,Nomo-K值出現(xiàn)先緩慢增加后快速增加的趨勢(shì)。這與指數(shù)方程所顯示的規(guī)律是一致的?？傊?當(dāng)土壤有機(jī)質(zhì)≤12%時(shí),使用Nomo公式計(jì)算K值;當(dāng)有機(jī)質(zhì)>12%時(shí),采用EPIC公式計(jì)算K值并利用其與Nomo-K值的指數(shù)關(guān)系式進(jìn)行修正。

1.3 環(huán)境因素分布圖收集與處理

沿用在太湖流域片預(yù)測(cè)K值時(shí)使用的遙感方法獲取的環(huán)境因子指標(biāo)作為預(yù)測(cè)變量,包括氣候、地表溫度、植被、光波、地形和母質(zhì)[13]。由于隨機(jī)森林模型具有避免過(guò)擬合的優(yōu)勢(shì),且模型訓(xùn)練集樣本量足夠大,沒(méi)有對(duì)環(huán)境因子進(jìn)行剔除。

各項(xiàng)環(huán)境指標(biāo)的空間分辨率從30 m到1 km不等,通過(guò)最鄰近法將所有指標(biāo)統(tǒng)一重采樣為30 m分辨率后進(jìn)行疊加運(yùn)算。由于氣候和地表溫度在較小尺度下不存在大的差異,重采樣后在1 km范圍內(nèi)氣候及地表溫度值不變符合實(shí)際情況。各樣點(diǎn)對(duì)應(yīng)的環(huán)境變量值通過(guò)坐標(biāo)位置提取。

1.4 土壤可蝕性K值制圖

采用隨機(jī)森林回歸方法建立K值預(yù)測(cè)模型。首先對(duì)樣點(diǎn)的環(huán)境要素及K值進(jìn)行訓(xùn)練,然后利用多源遙感影像將建立起K值預(yù)測(cè)模型推廣至全國(guó)地區(qū),完成K值由點(diǎn)及面的反演計(jì)算。其中模型的2個(gè)關(guān)鍵參數(shù)通過(guò)調(diào)參選擇:模型誤差隨著決策樹(shù)數(shù)量(ntree)的增加而降低,當(dāng)決策樹(shù)數(shù)量超過(guò)400以后,模型誤差趨于穩(wěn)定,因此選擇使用默認(rèn)參數(shù)ntree=500;隨機(jī)森林節(jié)點(diǎn)變量數(shù)(mtry)一般選擇環(huán)境變量數(shù)的1/3(本研究對(duì)應(yīng)mtry=20),通過(guò)對(duì)節(jié)點(diǎn)變量數(shù)在mtry=20左右逐個(gè)比較模型精度,確定最佳的mtry值為17。

K值預(yù)測(cè)模型的建模精度通過(guò)外部驗(yàn)證的方法進(jìn)行評(píng)估,采取常用的5折交叉驗(yàn)證法,重復(fù)100次取平均值作為評(píng)價(jià)結(jié)果。評(píng)價(jià)指標(biāo)包括決定系數(shù) (R2),均方根誤差 (E1)以及平均絕對(duì)誤差 (E2),計(jì)算公式如下:

(5)

(6)

(7)

式中:n為訓(xùn)練樣本個(gè)數(shù);Yi為觀測(cè)K值;i為模型預(yù)測(cè)K值;為平均觀測(cè)K值。

基于K值預(yù)測(cè)模型,以61項(xiàng)環(huán)境因子空間分布圖作為自變量,使用隨機(jī)森林方法進(jìn)行預(yù)測(cè),輸出全國(guó)尺度K因子?xùn)鸥駡D,K值更新圖的空間分辨率為30 m×30 m。

1.5 數(shù)據(jù)處理與統(tǒng)計(jì)分析

研究利用Excel 2019計(jì)算樣點(diǎn)Nomo-K值與EPIC-K值,并建立兩者的回歸關(guān)系;利用ArcGIS 10.8處理空間數(shù)據(jù),包括投影轉(zhuǎn)換、空間重采樣,提取樣點(diǎn)對(duì)應(yīng)的環(huán)境變量,統(tǒng)計(jì)柵格極值、平均值及標(biāo)準(zhǔn)差等;利用R語(yǔ)言編程實(shí)現(xiàn)土壤粒徑的批量插值轉(zhuǎn)換,同時(shí)構(gòu)建隨機(jī)森林回歸模型,進(jìn)行全國(guó)范圍K值的空間預(yù)測(cè)。在分析K值分布規(guī)律時(shí),土壤類型圖使用1995年編制的《1∶100萬(wàn)中華人民共和國(guó)土壤圖》。

2 結(jié)果與分析

2.1 土壤質(zhì)地插值轉(zhuǎn)換

插值獲得諾模圖計(jì)算所需極細(xì)砂的質(zhì)量分?jǐn)?shù),以江西省關(guān)山系樣點(diǎn)為例,使用3種插值方法得到粒徑質(zhì)量分?jǐn)?shù)結(jié)果如表1所示。其中,自然對(duì)數(shù)線性插值法擬合精度較差,獲得的<0.100 mm粒徑質(zhì)量分?jǐn)?shù)低于<0.050 mm粒徑,導(dǎo)致極細(xì)砂質(zhì)量分?jǐn)?shù)為負(fù)值;三次樣條函數(shù)法受樣本數(shù)量過(guò)少(n=3)的影響,獲得的累加極細(xì)砂粒徑(<0.100 mm)質(zhì)量分?jǐn)?shù)高于100%,明顯不合理;三次樣條函數(shù)結(jié)合自然對(duì)數(shù)法克服了函數(shù)分布不合理的缺點(diǎn),以高精度方程擬合了不同土壤粒徑的質(zhì)量分?jǐn)?shù)(R2≈1),同時(shí)避免了質(zhì)量分?jǐn)?shù)高于100%的情況(圖1)。

圓點(diǎn)為已知粒徑質(zhì)量分?jǐn)?shù),三角形為預(yù)測(cè)粒徑質(zhì)量分?jǐn)?shù)。下同。The dot represents the known particle size mass fraction, and the triangle represents the predicted particle size mass fraction. The same below. 圖1 土壤質(zhì)地插值轉(zhuǎn)換方法對(duì)比Fig.1 Comparison of soil texture interpolation conversion methods

表1 土壤質(zhì)地插值轉(zhuǎn)換結(jié)果對(duì)比(江西關(guān)山系)Tab.1 Comparison of soil texture interpolation conversion results(Guanshan series in Jiangxi)

在少數(shù)樣點(diǎn)中,土壤<0.050 mm和<2.000 mm的粒徑質(zhì)量分?jǐn)?shù)十分接近,使用三次樣條函數(shù)結(jié)合自然對(duì)數(shù)法計(jì)算所得<0.100 mm的粒徑結(jié)果可能略高于100%。針對(duì)這種情況,以山東省小島河系樣點(diǎn)為例(圖2),利用線性函數(shù)在<0.050 mm和<2.000 mm分段之間插值得到<0.100 mm的粒徑質(zhì)量分?jǐn)?shù),線性插值使<0.100 mm粒徑質(zhì)量分?jǐn)?shù)既小于<2.000 mm粒徑質(zhì)量分?jǐn)?shù),也大于<0.050 mm粒徑質(zhì)量分?jǐn)?shù),能確保結(jié)果的合理性。

實(shí)線為三次樣條函數(shù)結(jié)合自然對(duì)數(shù)法擬合,虛線為自然對(duì)數(shù)線性插值法。The solid line is cubic spline function combined with natural logarithm, and the dotted line is linear interpolation with natural logarithm. 圖2 土壤粒徑插值結(jié)果過(guò)高處理Fig.2 Processing of excessive soil particle size interpolation results

2.2 土壤可蝕性K值的計(jì)算與預(yù)測(cè)模型

采用Nomo和EPIC公式分別計(jì)算4 327個(gè)土系樣點(diǎn)K值,統(tǒng)計(jì)結(jié)果如表2所示。EPIC計(jì)算K因子的最小值比Nomo-K值大,最大值比Nomo-K值小,且標(biāo)準(zhǔn)差更小,EPIC計(jì)算K值的平均值及中位數(shù)比Nomo-K值大16%。因此,不能直接使用EPIC-K值替換高土壤有機(jī)質(zhì)條件下的Nomo-K值。

表2 全國(guó)土壤可蝕性K值計(jì)算結(jié)果統(tǒng)計(jì)Tab.2 Statistics of calculated results of soil erodibility K value nationwide

建立2種公式的轉(zhuǎn)換關(guān)系時(shí),R2均>0.79(圖3),但散點(diǎn)圖顯示隨著EPIC-K值的增大,Nomo-K值出現(xiàn)先緩慢增加后快速增加的趨勢(shì),這與指數(shù)方程的規(guī)律相符。最終選擇該方程對(duì)全國(guó)有機(jī)質(zhì)>12%的55個(gè)樣點(diǎn)進(jìn)行K值修正。這些樣點(diǎn)主要分布在東北和青藏高原的濕潤(rùn)或者滯水條件下,采樣部位在上、中、下坡和臺(tái)地均有分布,表明并非沉積導(dǎo)致K出現(xiàn)負(fù)值。根據(jù)定義,K因子為在標(biāo)準(zhǔn)小區(qū)上由單位降雨侵蝕力所引起土壤流失量的多少[1]。土壤流失量不可能為負(fù)數(shù),所以K值不應(yīng)存在負(fù)值的情況。校正前55個(gè)樣點(diǎn)的EPIC-K均值為0.034 4,校正后均值為0.026 6,彌補(bǔ)了EPIC-K值大于Nomo-K值的差異。

EPIC-K is soil erodibility calculated by EPIC model. Nomo-K is soil erodibility calculated by Nomo formula.圖3 土壤可蝕性K值計(jì)算結(jié)果擬合方程對(duì)比Fig.3 Comparison of fitting equations for the calculated results of soil erodibility K

以土系調(diào)查項(xiàng)目所建立的更新K值作為因變量,以包含氣候、地形、植被、母質(zhì)類型、地表溫度、可見(jiàn)光及近紅外波段在內(nèi)的61項(xiàng)環(huán)境因子指標(biāo)作為自變量,在mtry=17,ntree=500的模型參數(shù)下,建立的全國(guó)K因子隨機(jī)森林預(yù)測(cè)模型,運(yùn)用5折交叉驗(yàn)證方法重復(fù)100次的評(píng)估結(jié)果為:R2=0.381,E1=0.012 47,E2=0.009 63。這與太湖流域片建立K值模型的精度類似[13],但低于使用隨機(jī)森林預(yù)測(cè)土壤有機(jī)質(zhì)和全氮的空間分布[14]。這與本研究采樣密度稀疏、空間建模尺度大以及缺乏與K值強(qiáng)相關(guān)性的預(yù)測(cè)變量等原因有關(guān)。

2.3 全國(guó)土壤可蝕性K值分布規(guī)律

全國(guó)K因子分布范圍在0.005 1～0.074 5之間,這與前人在樣點(diǎn)尺度建立的我國(guó)主要土壤K因子分布范圍0.000 8～0.070 5相似[15]。分布趨勢(shì)呈現(xiàn)出華北平原和黃土高原地區(qū)高,天山地區(qū)及東北平原次之,南方山地丘陵、青藏高原及內(nèi)蒙古高原較低的規(guī)律。這與其他學(xué)者總結(jié)的我國(guó)土壤可蝕性分布規(guī)律相似[16],但該研究發(fā)現(xiàn)青藏高原的K值最高。也有研究發(fā)現(xiàn)青藏高原K值比黃土高原、川渝地區(qū)和滇東北地區(qū)的K值略低[17],這可能是由于青藏高原K值存在較大差異所造成的。

進(jìn)一步利用水土保持一級(jí)區(qū)劃分析全國(guó)更新K值的分布規(guī)律,發(fā)現(xiàn)更新K值在西北黃土高原區(qū)和北方土石山區(qū)較大,在青藏高原區(qū)最小(表3),這與自然及人為活動(dòng)等因素的綜合作用有關(guān)。其中,西北黃土高原區(qū)土壤有機(jī)質(zhì)含量普遍很低因此K值偏高[18];北方土石山區(qū)長(zhǎng)期以來(lái)受人類耕作活動(dòng)的影響導(dǎo)致K值較高;青藏高原區(qū)由于地形落差大,K值具有垂直地帶性,在少量低海拔河谷地帶K值較高,而在大面積的高海拔地區(qū)由于成土作用弱K值明顯降低[17]。此外,東北黑土區(qū)在風(fēng)砂土分布地區(qū)K值最大,在黑土集中的漫川漫崗地帶由于農(nóng)耕活動(dòng)導(dǎo)致K值較大,暗棕壤分布的大興安嶺及長(zhǎng)白山等地K值較小。

表3 土壤可蝕性更新K值在水土保持一級(jí)區(qū)劃的統(tǒng)計(jì)Tab.3 Statistics of soil erodibility update K value in the first level zoning of soil and water conservation

K值的空間分異與我國(guó)主要土壤類型的分布情況有關(guān),前人曾報(bào)道黃土高原的黃綿土和東北的黑土K值較大,南方紅壤K值較小[15,19]。本研究從主要土壤類型統(tǒng)計(jì)來(lái)看,黃綿土和潮土的更新K值較大,平均值均超過(guò)0.04,而高山土的平均K值最小,其次為栗鈣土、暗棕壤和紅壤,均不超過(guò)0.03(表4)。分析原因認(rèn)為黃綿土和潮土的粉粒含量較高,土壤更容易遭受侵蝕。此外,潮土K值的標(biāo)準(zhǔn)差最大。前人研究發(fā)現(xiàn),同一種土壤下不同母質(zhì)和植被覆蓋可以通過(guò)影響土壤的水穩(wěn)性從而改變土壤的抗蝕能力[20],因此相同土壤類型下K值的差別也可能很大。

表4 主要土壤類型下土壤可蝕性更新K值統(tǒng)計(jì)Tab.4 Statistics of updated K values of soil erodibility under main soil types

值得注意的是,由于Nomo和EPIC模型是基于美國(guó)土壤實(shí)驗(yàn)所建立的,在中國(guó)地區(qū)存在適用性問(wèn)題。張科利等[19]利用徑流小區(qū)觀測(cè)資料在我國(guó)東部地區(qū)建立了K值的校正公式,分別為K=-0.033 36+0.744 88KNomo(R=0.721),K=-0.013 83+0.515 75KEPIC(R=0.613);更進(jìn)一步說(shuō),Nomo模型在黃土高原和黑土區(qū)的應(yīng)用效果(R2=0.55～0.60)要優(yōu)于紅壤和紫色土(R2=0.11～0.37),EPIC模型除黑土區(qū)以外(R2=0.05)在其他地區(qū)的應(yīng)用效果較好(R2>0.58)[21]。直接應(yīng)用本研究結(jié)果可能會(huì)對(duì)K值產(chǎn)生不同程度的高估,因此需要在全國(guó)范圍內(nèi)進(jìn)一步收集實(shí)測(cè)資料將K值校正以后用于土壤侵蝕模數(shù)計(jì)算。

3 結(jié)論

1)本次更新進(jìn)一步完善了全國(guó)土壤可蝕性因子計(jì)算方法。本研究在美國(guó)制砂粒、粉粒、黏粒機(jī)械組成實(shí)測(cè)值的基礎(chǔ)上插值得到極細(xì)砂質(zhì)量分?jǐn)?shù),基礎(chǔ)數(shù)據(jù)更為準(zhǔn)確;通過(guò)計(jì)算土壤粒徑幾何平均直徑來(lái)獲取結(jié)構(gòu)系數(shù),更具有客觀性;建立于全國(guó)統(tǒng)一的土系基礎(chǔ)數(shù)據(jù)庫(kù),消除了不同省份之間的系統(tǒng)誤差;以30 m分辨率柵格為制圖單位,更為詳細(xì)地表達(dá)K值的空間變異情況。

2)全國(guó)土壤可蝕性更新K值的空間分布呈現(xiàn)出以下規(guī)律:黃土高原區(qū)和北方土石山區(qū)居高,西南紫色土區(qū)、南方紅壤區(qū)、西南巖溶區(qū)、東北黑土區(qū)和北方風(fēng)沙區(qū)次之,青藏高原區(qū)最小。這主要是由于不同土壤類型在各地區(qū)的分布所造成的。統(tǒng)計(jì)可知黃綿土和潮土的K值較大,紅壤、暗棕壤、栗鈣土的K值較小,高山土的K值最小。

3)本研究采用的Nomo和EPIC公式都是基于美國(guó)土壤建立的K值計(jì)算經(jīng)驗(yàn)公式,可能與實(shí)際監(jiān)測(cè)的K真實(shí)值有差距。因此需要積累我國(guó)不同地區(qū)標(biāo)準(zhǔn)徑流小區(qū)的長(zhǎng)期觀測(cè)資料,使用年平均侵蝕量和降雨侵蝕力因子進(jìn)一步率定K值,以期建立起實(shí)測(cè)值與經(jīng)驗(yàn)公式計(jì)算值之間的關(guān)系式或直接建立適用于我國(guó)的土壤可蝕性計(jì)算公式。