不同干擾程度下土壤有機質(zhì)空間最優(yōu)插值法研究

馬利芳,熊黑鋼,孫 迪,王 寧,葉紅云,張 芳

1 新疆大學資源與環(huán)境科學學院,綠洲生態(tài)教育部重點實驗室, 烏魯木齊 830046 2 北京聯(lián)合大學應(yīng)用文理學院, 北京 100083

有機質(zhì)是土壤必不可少的組成成分,是土壤質(zhì)量及肥力評價的基礎(chǔ),在土壤發(fā)揮功能與碳循環(huán)過程中起重要作用[1-2]。明確有機質(zhì)的空間分布特征是土壤資源與環(huán)境科學管理的依據(jù)[3-4]。近些年,人類活動所造成的干擾對地表自然狀況及生態(tài)環(huán)境的影響不斷加劇[5],不同程度的人為干擾對土壤有機質(zhì)空間分布的影響有所差異,這在有機質(zhì)十分匱乏的干旱區(qū)更加明顯。因此,對干旱區(qū)不同人為干擾程度下土壤有機質(zhì)含量進行空間插值研究,對掌握土壤有機質(zhì)的空間分布規(guī)律以及實現(xiàn)干旱地區(qū)農(nóng)業(yè)可持續(xù)發(fā)展意義重大[6]。

當前,國內(nèi)外已有大量關(guān)于土壤屬性空間插值的研究,地統(tǒng)計方法尤其是克里格插值法的應(yīng)用相對廣泛[7-8]。例如：利用普通克里格法分析了克什米爾農(nóng)業(yè)區(qū)[9]、尼羅河三角洲地區(qū)[10]、意大利耕地和牧場[11]、北京某生態(tài)功能區(qū)[12]、典型喀斯特峰叢洼地小流域[13]內(nèi)土壤屬性的空間分布特征等,但亦有研究顯示普通克里格法并不能夠很好地對土壤有機質(zhì)進行插值分析[14]。與傳統(tǒng)的克里格法相比,采用改進土地利用回歸法[15]、隨機森林法[16]、地理加權(quán)法[17]等方法或者借助輔助變量[18]對土壤屬性進行空間分布特征分析,會實現(xiàn)提高空間插值精度的目的。已有關(guān)于土壤屬性空間插值方法精度比較的研究,結(jié)論并不一致[19-20],能夠適用于所有土壤屬性,且精度在任何區(qū)域都達到最優(yōu)的插值方法并不能被確定?，F(xiàn)有研究結(jié)果最終選擇的最優(yōu)空間插值方法相差迥異。例如,有學者認為徑向基函數(shù)法能更好地表達土壤屬性的空間分布特征[21]；亦有研究采用不同插值方法對土壤屬性插值精度進行比較,發(fā)現(xiàn)反距離加權(quán)法的估算效果更佳[22]。還有相關(guān)研究表明不同土層深度的有機質(zhì)含量最優(yōu)插值方法也有所差異[23]。

現(xiàn)有成果可能是在不同土壤環(huán)境背景下的討論,故有較大的差異。因而,利用反距離加權(quán)法、徑向基函數(shù)法和局部多項式法對不同人類干擾程度下的土壤有機質(zhì)空間特征進行研究,探討不同空間插值方法對其估算精度的影響,以期尋求不同人類干擾下土壤有機質(zhì)的最優(yōu)空間插值方法,為提高空間估測精度提供價值參考。

1 材料與方法

1.1 研究區(qū)概況及樣點布設(shè)

圖1 研究區(qū)域位置及采樣點分布圖 Fig.1 Location map of the study area and distribution of the sampling sites

研究區(qū)位于天山北坡東段,準噶爾盆地南緣,地理位置為87°44′—88°46′E,43°29′—45°45′N,為典型的干旱荒漠區(qū)。該區(qū)氣候?qū)儆谥袦貛Т箨懶詺夂?冬季長干冷、嚴寒多雪,夏季短干熱,春秋季節(jié)不明顯,蒸發(fā)強烈,年均氣溫6.6℃,光照充足,熱量豐富,溫差大,降水稀少且空間分布不均,年降水量僅186mm,年蒸發(fā)潛力2064mm左右。土壤類型為灰漠土。

根據(jù)研究區(qū)土壤受人類活動的干擾程度,將其劃分為3個區(qū)域(圖1)。無人為干擾區(qū)(A區(qū))幾乎沒有受到人為活動的干擾,基本保持景觀的原有風貌,地表植被相對豐富,部分地點蓋度更高、植株較大,有梭梭、琵琶柴、鹽爪爪、白刺、紅柳及成片的鹽生雜草等,植被覆蓋度約為30%,且整個區(qū)域有大量較厚且發(fā)育良好的黑色生物結(jié)皮；人為干擾區(qū)(B區(qū))內(nèi)主要為棄耕地,地表被犁翻耕且有很明顯的駱駝踩踏痕跡。有豬毛菜及少量琵琶柴、梭梭、紅柳、鹽爪爪等植被,但植被覆蓋度相對較低約15%—20%,區(qū)內(nèi)土壤表面生物結(jié)皮較少,發(fā)育較差,部分地表無生物結(jié)皮,人為干擾較強烈；重度人為干擾區(qū)(C區(qū))包含兩片農(nóng)場,受到強烈的人為干擾。土地全部被翻耕,人工種植的榆樹林行間距均為3.5m,榆樹林株距分別為3m和1.2m,平均高度為3m左右,最高的植株約4.3m,冠幅分別為0.5m×0.5m和1.1m×1.0m。3個區(qū)域相鄰僅以欄網(wǎng)或溝渠相隔,整個研究區(qū)位于綠洲下緣的平原區(qū),其地形、土壤屬性、光照時長、熱量分布、降水量、溫度濕度等自然條件基本相似。因此,在干擾等級進行劃分時,主要考慮的是土壤受人類活動的干擾程度。

雖然各區(qū)面積大小不同,根據(jù)野外實地考察情況,使得各區(qū)采樣線間距設(shè)置有所差異,不能將其等分,但各區(qū)均采用網(wǎng)格法布設(shè)樣點(圖1)。這樣就在各區(qū)形成可以全面控制該區(qū)土壤有機質(zhì)變化情況的網(wǎng)格,以保證數(shù)據(jù)的代表性和合理性。A區(qū)位于研究區(qū)東部,由南向北布設(shè)5條采樣線,每條采樣線上分布6個樣點；C區(qū)位于A區(qū)西側(cè),由南向北布設(shè)6條采樣線,每條采樣線上設(shè)置5個樣點；B區(qū)位于C區(qū)北部,由南向北布設(shè)5條采樣線,依據(jù)當?shù)貙嶋H情況,每條采樣線上設(shè)定5—7個樣點,3個區(qū)域均設(shè)置30個采樣點,每個樣點間距保持在300—400m左右,盡可能均勻分布。

1.2 樣品采集與處理

以Landsat系列遙感影像數(shù)據(jù)為參考,于2017年10月進行野外土壤樣品采集。根據(jù)樣點設(shè)置,按照一致性、同質(zhì)性和代表性原則先選擇10m×10m范圍的植物樣方,主要記錄每個樣方內(nèi)植被種類、個體數(shù),總蓋度、高度以及土壤表層生物結(jié)皮生長情況,同時,詳細記錄每個樣點周邊環(huán)境特征。為使土壤樣品代表性更強,先測完植被樣方,然后在樣點周圍2m范圍內(nèi)以梅花樁方式采集5處土樣混合均勻后放入密封袋中,封口標記編號,并用GPS定位。

將采集的土壤樣品經(jīng)過預(yù)處理之后送至中國科學院新疆生態(tài)與地理研究所理化測試中心,進行土壤有機質(zhì)數(shù)據(jù)的測定。

1.3 數(shù)據(jù)處理方法與檢驗

采用4種方法對不同干擾程度下土壤有機質(zhì)進行空間插值分析,并對不同插值法進行交叉驗證,獲得精度評價結(jié)果,得到基于最優(yōu)插值方法的不同干擾程度下土壤有機質(zhì)含量空間分布圖。

1.3.1普通克里格(OK)空間插值法

普通克里格插值法(Ordinary Kriging)的適用條件是區(qū)域化變量具有空間相關(guān)性[24]。利用半變異函數(shù),選取球狀模型和高斯模型,來反映各區(qū)域變量的結(jié)構(gòu)性和隨機性。

1.3.2反距離加權(quán)(IDW)插值法

反距離加權(quán)插值(Inverse Distance Weighted)基于相近相似的原理：即兩個物體的性質(zhì)是否相似取決于兩者距離的遠近,離得越近性質(zhì)越相似,反之,離得遠則相似性小[25]。

1.3.3徑向基函數(shù)(RBF)插值法

徑向基函數(shù)(Radial Basis Function)插值法是一種精確的非線性計算工具[26]。每個點都能用幾個高斯函數(shù)的疊加進行逼近[27],對于距離較遠的點,其影響小于距離較近的點,這一特性排除了遠距離點的干擾,使訓練速度更快。

1.3.4局部多項式(LPI)插值法

局部多項式內(nèi)插(Local Polynomial Interpolation)基于局部加權(quán)最小二乘法對處在特定重疊鄰近區(qū)域內(nèi)的多個多項式進行擬合的一種近似插值法[28]。

1.3.5插值結(jié)果的精度檢驗

插值的精度分析采用交叉驗證的K折交叉驗證(K-fold cross-validation),它是一種能快速評價插值結(jié)果質(zhì)量的方法[29]。為更加準確地評估預(yù)測模型精度,采用最常用的十折交叉驗證法,即將樣本隨機分成10份,輪流將其中的9份用于訓練,剩余一份用于評估,循環(huán)10次后所有數(shù)據(jù)都會被驗證一次,取10次結(jié)果的均值為最終預(yù)測誤差。該方法的特點是可直接進行誤差的估算無需任何要求,操作便捷、適用性強；尤其是數(shù)據(jù)量較小時計算效率會更高,數(shù)據(jù)的循環(huán)使用能更接近原始樣本分布且不容易受到隨機因素的影響。

評價指標采用決定系數(shù)(R2)平均預(yù)測誤差(ME)和均方根預(yù)測誤差(RMSE),依據(jù)R2越接近1、ME越接近0、RMSE越小其預(yù)測精度越高的原則,尋求最優(yōu)插值方法。各指標計算公式如下：

(1)

(2)

(3)

式中,n為樣本數(shù),Xi為實測值,Xj為預(yù)測值。

2 結(jié)果與分析

2.1 不同干擾程度土壤有機質(zhì)統(tǒng)計特征

從無人為干擾、人為干擾到重度人為干擾的梯度上,土壤有機質(zhì)含量最大值、最小值、均值均逐漸減小,而標準差逐漸增大,即有機質(zhì)含量逐漸減小(表1)。根據(jù)全國第二次土壤普查有機質(zhì)分級標準[30],無人為干擾區(qū)處于四級中下等水平,人為干擾及重度人為干擾區(qū)均為五級缺乏水平。變異系數(shù)(CV)反映數(shù)據(jù)的離散程度,CV<10%為弱變異；10%≤CV≤100%為中等變異；CV>100%為強變異[31]。即無人為干擾區(qū)呈現(xiàn)弱變異,人為干擾及重度人為干擾區(qū)為中等變異。

表1 不同干擾程度土壤有機質(zhì)的描述性統(tǒng)計和K-S檢驗

N表示符合正態(tài)分布檢驗，標準差S.D.(Standard Deviation)，變異系數(shù)CV(Coefficient of Variation)，K-S檢驗(Kolmogorov-Smirnov)

2.2 不同干擾程度土壤有機質(zhì)空間變異特征

利用地統(tǒng)計學方法對研究區(qū)不同干擾程度下土壤有機質(zhì)進行模型擬合,比較各模型的決定系數(shù)(R2)和殘差平方和(RSS),選取R2更接近1 且RSS較小的最優(yōu)理論模型,得到模型相關(guān)參數(shù)值(表2)。隨著干擾程度的加劇,R2逐漸減小,由0.932下降到0.764,減小了0.168；RSS逐漸增大,由0.027增大到1.24。無人為干擾區(qū)(A區(qū))土壤有機質(zhì)最優(yōu)理論模型為球狀模型,人為干擾區(qū)(B區(qū))和重度人為干擾區(qū)(C區(qū))為高斯模型。當塊基比C0/(C0+C)<25%則以結(jié)構(gòu)性變異為主,具有強烈的空間相關(guān)性；比值>75%表示其空間變異以隨機性因素為主,空間相關(guān)性很弱；在25%—75%之間,則說明受結(jié)構(gòu)性和隨機性因素共同影響[32],空間相關(guān)性中等。無人為干擾區(qū)有機質(zhì)塊基比為12.26%,空間結(jié)構(gòu)性極強,人為干擾區(qū)塊基比為35.9%,空間相關(guān)性中等,而重度干擾區(qū)塊基比為76.21%,說明在該區(qū)內(nèi)土壤有機質(zhì)含量的空間變異受隨機性因素影響很大。

變程和分形維數(shù)(D)亦可反映區(qū)域化變量的空間自相關(guān)性大小及空間變異范圍尺度,D值越小,由空間自相關(guān)部分引起的空間變異性越弱,結(jié)構(gòu)性越好,受隨機因素影響越小[33]。從無人為干擾區(qū)(A區(qū))、人為干擾區(qū)(B區(qū))到重度干擾區(qū)(C區(qū)),土壤有機質(zhì)的變程和分形維數(shù)均在逐漸增大,說明隨著人類活動力度的加大,隨機因素對土壤有機質(zhì)含量作用越來越強,空間自相關(guān)部分引起的空間變異性越高?？梢?半方差函數(shù)充分反映了對不同程度人為干擾下土壤有機質(zhì)的空間變異特點。

表2 不同干擾程度土壤有機質(zhì)的半方差理論模型及其參數(shù)

S表示球狀模型;G表示高斯模型;決定系數(shù)R2(Coefficient of determination)，殘差平方和RSS(Residual sum of squares)，分形維數(shù)D(Dimension)

2.3 不同干擾程度土壤有機質(zhì)空間插值結(jié)果的交叉驗證

不同插值的交叉驗證結(jié)果知(表3),各方法均有誤差存在。A區(qū)土壤有機質(zhì)的4種插值方法中,R2最大、RMSE最小,且ME最接近0的是OK法,其次是IDW法,而后為RBF法,LPI法插值能力相對較弱；而對B區(qū)和C區(qū)土壤有機質(zhì)而言,各插值方法的精度分別為RBF法>IDW法>OK法>LPI法,RBF法>IDW法>LPI法>OK法,即RBF法能更精確地對B區(qū)及C區(qū)內(nèi)土壤有機質(zhì)含量進行空間插值分析。不同干擾程度下4種插值方法精度均表現(xiàn)為A區(qū)>B區(qū)>C區(qū),其中OK法的R2變化范圍較大,為0.312—0.625,且其對A區(qū)土壤有機質(zhì)的插值精度最高(RMSE為2.049)。在B區(qū)和C區(qū),OK法的插值精度相對變低,與A區(qū)相比較RMSE分別升高了0.482和1.033,R2分別下降了0.129和0.313；而RBF法和IDW法的RMSE較A區(qū)雖有所增大,R2雖有所減小,但其精度均高于OK法?？梢?人為干擾程度的加劇對各插值方法的精度都會產(chǎn)生影響,這是因為不同強度的人類活動對土壤屬性造成的干擾不同,使其理化特征發(fā)生改變,導致數(shù)據(jù)的離散程度各異,且各插值方法對不同區(qū)域內(nèi)土壤理化性質(zhì)變化的適用性有所差異。尤其是OK法的插值效果更易受隨機因素左右,在人類干擾強度大,土壤屬性變化強烈的區(qū)域進行插值是有一定局限的。

表3 不同干擾程度土壤有機質(zhì)插值方法精度檢驗

普通克里金OK，反距離加權(quán)IDW，徑向基函數(shù)RBF,局部多項式LPI，平均誤差ME (Mean Error)

2.4 基于最優(yōu)插值法土壤有機質(zhì)的空間預(yù)測

利用各區(qū)域的最適插值方法：A區(qū)采用OK法、B區(qū)和C區(qū)采用RBF法,對不同干擾程度下土壤有機質(zhì)進行空間預(yù)測,分析其空間分布特征(圖2)。A區(qū)土壤有機質(zhì)東部及東南部高,西部及東北部較低；B區(qū)中部及東部高,南部、西北部較低；C區(qū)在北部出現(xiàn)高值,中部及南部較低,且含量極值差距較大?？梢?種植、翻耕、灌溉等人為因素對土壤有機質(zhì)含量的分布格局有很大的影響,與近乎自然狀態(tài)下的A區(qū)相比,B區(qū)和C區(qū)干擾程度不同,對土壤有機質(zhì)分布特征產(chǎn)生的影響有所差異。這些分布特點與各區(qū)域土地利用類型、植被條件以及人類活動干擾程度等關(guān)系密切,A區(qū)以荒地為主,地表植被較豐富,幾乎沒有受到人為活動的干擾,從西至東有機質(zhì)變化平緩；B區(qū)以半荒地、棄耕地為主,覆蓋度相對較低,人為干擾具有隨機性,有機質(zhì)變化復(fù)雜；C區(qū)為棄耕地、人工林地,土壤全部被翻耕,人工林地種有梭梭林和榆樹林,受人為干擾最強烈。

圖2 不同干擾程度土壤有機質(zhì)含量空間分布預(yù)測Fig.2 Prediction of spatial distribution of soil organic matter in different disturbance districts

圖3 基于最適插值方法土壤有機質(zhì)預(yù)測值和實測值散點圖Fig.3 Scatter plot of predicted values and measured values of soil organic matter based on optimal interpolation method

基于各區(qū)域最優(yōu)空間預(yù)測模型,繪制出不同干擾程度下土壤有機質(zhì)預(yù)測值和實測值的散點圖(圖3),分析數(shù)據(jù)較為均勻地布局在1∶1線的兩側(cè),說明預(yù)測值和實測值總體上呈現(xiàn)出相對較好的線性關(guān)系。無人為干擾區(qū)有機質(zhì)含量值相對較高,且數(shù)據(jù)分布較集中,而隨人類干擾程度加劇,有機質(zhì)值越低,且含量變化越大,分布越分散。

3 討論

對于同一區(qū)域,采用不同的插值方法所得結(jié)果有所不同,同類插值方法在不同區(qū)域產(chǎn)生的效果亦會有差異?？臻g插值沒有一個通用的、普適的模型,而是需要根據(jù)不同特征的研究對象來選擇相對適宜的插值方法和相關(guān)參數(shù)[34]。

OK法多適用于具有強烈空間相關(guān)性的插值分析中,因為其空間插值精度一定程度上取決于待插值土壤屬性的空間變異特征,這與前人研究結(jié)果一致[35-36]。OK法雖然取得了較好的成果[37-38],但對不同干擾程度下土壤有機質(zhì)進行空間插值分析,并比較不同插值方法的精度,發(fā)現(xiàn)在樣點數(shù)量不變的條件下(各區(qū)采樣點均為30個),OK法僅對空間自相關(guān)性較強的、變異性較弱的無人為干擾區(qū)土壤有機質(zhì)插值效果最好,而在人為干擾區(qū)插值精度偏低,在重度人為干擾區(qū)最低。主要是由于其首先考慮的是空間屬性在空間位置上的變異分布,受人為干擾等隨機因素導致空間變異性較高,僅依據(jù)樣點間的地理位置信息對其進行插值必然會使其精度有所下降。因此,在不能充分滿足克里格插值前提條件且土壤屬性變化強烈的地區(qū)進行插值,會導致其精度偏低。同時,其算法較復(fù)雜,在選定變異函數(shù)時具有主觀性。

RBF法適合于需要將樣本值和擬合值保持一致的應(yīng)用中,原因是它一定程度上能克服平滑效應(yīng),是一種精確的非線性插值方法,其結(jié)果會盡可能地保留元素含量的極值信息[39]。尤其是在數(shù)據(jù)較少,土壤有機質(zhì)變異程度相對較大的情況下(人為干擾和重度人為干擾區(qū))進行空間插值效果較好。同時,其最大的特點是能在高維空間中利用高斯函數(shù)[40],可以不受任何約束地逼近任意函數(shù),而人為干擾區(qū)和重度人為干擾區(qū)土壤有機質(zhì)的空間變異特征均符合高斯模型,因此,RBF法更勝一籌。

IDW法一般應(yīng)用于對極值不要求度量的插值模型中,其算法相對簡單,較易實現(xiàn),適合分布較均勻且密集的樣點。但選擇函數(shù)冪次時具有敏感性,易受采樣集群影響造成“牛眼”現(xiàn)象[41]。與OK法和RBF法相比,它基于相近相似原理,又不能對樣本中的極大極小值進行預(yù)測,所以此法雖然能實現(xiàn)對各區(qū)域土壤有機質(zhì)插值運算,但精度不如前兩者高。而LPI法多適宜于解釋小范圍的局部變異,且平滑性較好的數(shù)據(jù)情況,因為它利用最小二乘法擬合元素含量的空間分布趨勢,趨向于得到一個平滑的表面[42],且其平滑作用比其他3種方法更加明顯。雖然它是一種非參數(shù)估計法,能消除異方差的影響[43],但其僅考慮樣本局部范圍趨勢,針對特定空間領(lǐng)域內(nèi)信息進行分析,所以不論在哪種情況下,LPI法的插值精度均不是最佳,其對各區(qū)域土壤有機質(zhì)含量的插值分析局限性較大。

影響空間插值精度和效果的因素有很多,比如：環(huán)境背景、相關(guān)參數(shù)的選擇以及采樣點分布密度狀況等。今后對土壤有機質(zhì)的空間特征進行討論時,應(yīng)首先區(qū)別其是否受到人類活動的干擾,其次應(yīng)詳細區(qū)分其受人類干擾的強度,然后再對其分析研究。這樣更有利于尋求不同干擾情況下土壤屬性的最優(yōu)插值方法,而后才能更加精準地進行空間分析。

4 結(jié)論

不同干擾程度下土壤有機質(zhì)均呈正態(tài)分布,且無人為干擾區(qū)有機質(zhì)具有弱變異性和強烈的空間自相關(guān)性,人為干擾區(qū)空間變異及相關(guān)性均呈中等強度,重度人為干擾區(qū)為中等變異,而空間相關(guān)性較弱。土壤有機質(zhì)含量受人為干擾活動等隨機因素作用越強,空間自相關(guān)部分引起的空間變異性越高。

隨著干擾程度加劇,無論采用哪種方法對土壤有機質(zhì)空間插值分析,其精度均在降低。從無人為干擾區(qū)到人為干擾區(qū)再至重度干擾區(qū),各種方法的插值精度R2由0.487—0.625降低為0.425—0.562再降至0.312—0.434。其中OK法一定程度上依賴于有機質(zhì)的空間變異特征,在空間結(jié)構(gòu)性強的無人為干擾區(qū),其插值精度最高；而在人為干擾和重度人為干擾區(qū)插值效果最好的是RBF法。