基于輔助環(huán)境變量的土壤有機碳空間插值
——以黃土丘陵區(qū)小流域為例

文 雯，周寶同，汪亞峰 ，黃 勇

(1.西南大學(xué)地理科學(xué)學(xué)院，重慶 400715;2.中國科學(xué)院生態(tài)環(huán)境研究中心城市與區(qū)域生態(tài)國家重點實驗室，北京 100085;3.東北林業(yè)大學(xué)林學(xué)院，哈爾濱 150040)

土壤碳庫是陸地生態(tài)系統(tǒng)中最大的碳庫，準(zhǔn)確估算土壤碳庫儲量對正確評價土壤在陸地生態(tài)系統(tǒng)碳循環(huán)、全球碳循環(huán)以及全球環(huán)境變化中的作用有重要意義［1-3］。土壤有機碳是土壤中較為活躍的部分，其含量和動態(tài)在土壤質(zhì)量演變和全球碳循環(huán)中起著十分重要的作用，也是估算土壤碳庫儲量的關(guān)鍵因素［4-5］。對土壤有機碳含量的空間差異進(jìn)行研究，獲得區(qū)域土壤有機碳含量的連續(xù)空間分布，是精準(zhǔn)農(nóng)業(yè)技術(shù)體系中管理和分析土壤背景數(shù)據(jù)的一項基礎(chǔ)性工作，然而大規(guī)模的土壤有機碳空間分布信息的采集需要花費大量的人力、物力和財力。因此，研究如何利用有限的樣本數(shù)據(jù)來獲得更為詳盡的土壤有機碳空間分布信息的方法具有重要意義。

空間插值模型是實現(xiàn)土壤有機碳含量從離散的點狀信息向面狀連續(xù)信息轉(zhuǎn)換的有力工具，是表征土壤有機碳空間分布特征的重要手段，可以將有限的采樣點用來預(yù)測整個研究區(qū)域的土壤屬性空間信息［6-9］。空間插值方法包括確定性方法和地理統(tǒng)計方法［10］:確定性方法是使用數(shù)學(xué)函數(shù)進(jìn)行插值，以研究區(qū)域內(nèi)部的相似性(如反距離加權(quán)法)，或以平滑度為基礎(chǔ)(如徑向基函數(shù)法)，由已知樣點來創(chuàng)建預(yù)測表面的插值方法;地理統(tǒng)計方法是利用已知樣點的統(tǒng)計屬性，在半變異函數(shù)理論分析基礎(chǔ)上，對區(qū)域化變量的取值實現(xiàn)無偏最優(yōu)估計的一種方法，如普通克里格法。從學(xué)者們的研究中不難發(fā)現(xiàn)，普通克里格插值法是使用頻率最高的空間預(yù)測方法，但在地形復(fù)雜的區(qū)域，普通克里格方法的應(yīng)用有一定的局限性［11］，因而一些結(jié)合地貌類型、土壤類型、土地利用類型和海拔高度等信息的克里格插值模型相繼提出，通過輔助信息可以提高目標(biāo)變量的預(yù)測精度［10-12］。

土地利用和管理水平是影響土壤有機碳含量最普遍、最直接的因素之一，很大程度上影響著土壤有機碳變化的程度和方向［13］，國內(nèi)外學(xué)者已對不同土地利用方式下土壤碳庫特征進(jìn)行了對比分析［14-16］。黃土丘陵地區(qū)，地形復(fù)雜，土地利用類型多樣，土地利用類型是土壤有機碳含量空間分布差異的重要影響因素之一，有必要在土壤有機碳的空間插值中考慮土地利用因素。對土壤屬性空間預(yù)測的插值方法的對比研究很多，然而對于適宜黃土丘陵地形復(fù)雜地區(qū)插值方法的對比研究還較少，黃土丘陵地區(qū)基于土地利用類型修正的土壤有機碳含量空間插值研究更不多見。本研究針對黃土丘陵地區(qū)現(xiàn)有的土壤有機碳研究中較少考慮的不同空間插值方法對預(yù)測精度的影響、土地利用類型修正對插值方法精度的影響等問題，以羊圈溝流域作為研究區(qū)域，分析比較3種常用的普通克里格法(OK)、反距離加權(quán)法(IDW)、徑向基函數(shù)法(RBF)，及經(jīng)過土地利用類型修正的普通克里格法(OK_LU)進(jìn)行精度對比分析，揭示了最適宜地形復(fù)雜的黃土丘陵溝壑區(qū)土壤有機碳空間插值研究的方法，以期為類似區(qū)域的土壤屬性空間分布信息研究提供參考依據(jù)。

1 材料與方法

1.1 流域自然概況

羊圈溝流域位于延安市寶塔區(qū)東北方向14 km處，屬寶塔區(qū)李渠鎮(zhèn)轄區(qū)，為延河左岸的二級支溝(36°42'N，109°31'E)，該流域?qū)儆谀肭f溝流域的一級支流。流域內(nèi)存在一系列早期建造的淤地壩，自1999年流域內(nèi)實施了退耕還林政策。該區(qū)域?qū)冱S土丘陵溝壑區(qū)第二副區(qū)，年太陽總輻射量5800 kJ/cm2，年日照時數(shù)2528.8 h，一般年份能完全滿足小麥、玉米、谷類作物的生長要求。流域植被在區(qū)劃上屬于森林草原過渡帶。流域內(nèi)由于人為活動的干擾，自然植被破壞殆盡，多為人工種植而形成的次生植被，植物種類主要有刺槐(Robiniapseudoacacia)、柳樹(Salix spp.)和楊樹(Poplar spp.)等［17］。

2011年羊圈溝流域的八種主要土地利用類型為林地、灌木林地、草地、果園、梯田、壩地、水域和建設(shè)用地(表1)。主要的分布格局為:林地、灌木林地、草地遍布在整個流域;果園分布在流域中部;水域、居住用地和梯田主要分布在流域南部;壩地沿溝分布(圖1)。

1.2 樣品采集與分析

2011年分別在延安羊圈溝流域進(jìn)行了兩次野外全流域采樣。本研究采用地形剖面線法進(jìn)行坡面采樣，采樣時按土地利用類型和土壤斷面布設(shè)采樣點，同時考慮微地貌，每個土壤斷面取3至6個點，取樣器為荷蘭Eijkekamp公司產(chǎn)直徑為6 cm的半圓鑿型土鉆。采樣深度為0—20 cm，深度間隔為0—5 cm，5—10 cm和10—20 cm。為體現(xiàn)采樣點分布的均勻性，根據(jù)各地類所占比例的大小，分地類進(jìn)行樣點的采集，以保證樣點在不同土地利用類型中的合理分布，如圖1、表1所示。將野外采集的土壤樣品帶回實驗室后，先經(jīng)過風(fēng)干、研磨、過篩、稱重和裝盒等過程后，再進(jìn)入實驗室分析階段。有機碳含量的測定使用K2Cr2O4外加熱法。

1.3 數(shù)據(jù)處理

在Mapgis軟件中通過GPS將野外均勻布點采集樣點的經(jīng)緯度坐標(biāo)與羊圈溝流域的空間分布坐標(biāo)相連接，將測得的數(shù)據(jù)導(dǎo)入羊圈溝流域的點位表中，添加到土地利用類型圖中，得到研究區(qū)域采樣點的分布圖(圖1)。將研究區(qū)域內(nèi)采樣點0—5 cm，5—10 cm，10—20 cm測得的土壤有機碳含量值綜合，得到0—20 cm土層厚度的土壤有機碳含量值。本研究主要采用0—20 cm土層的土壤有機碳含量值作為研究對象。通過對比分析3種常用的插值方法與基于土地利用類型修正的方法的預(yù)測精度，分析本研究區(qū)域的最優(yōu)插值方法。

1.3.1 插值方法的選擇

本研究中土壤有機碳的空間插值的4種方法均基于ArcGIS10.0中的ArcMap平臺實現(xiàn)。

圖1 研究區(qū)域樣點分布圖Fig．1 Map of sample distribution

表1 黃土丘陵區(qū)羊圈溝流域土地利用面積、采樣數(shù)及其所占比例Table 1 Land-use type area，sample size and their proportion in Yangjuangou watershed of hilly region，Loess Plateau

(1)普通克里格法(Ordinary Kriging，簡稱OK)

普通克里格法是利用區(qū)域化變量的原始數(shù)據(jù)和變異函數(shù)的結(jié)構(gòu)特點，對未采樣點的區(qū)域化變量的取值進(jìn)行線性無偏最優(yōu)估計的一種方法。使用公式表示為:

式中，Z為待估計的土壤有機碳含量柵格值，n為用于插值的采樣點數(shù)目，Z(xi)為采樣點土壤有機碳含量值，λi為賦予土壤有機碳含量數(shù)據(jù)的一組權(quán)重系數(shù)，并且其和等于1，即:=1，選取的λi使Z(x0)的估計無偏，并且使方差小于任意觀測值線性組合的方差［18］。

(2)逆距離權(quán)重法(Inverse Distance Weighting，簡稱IDW)

逆距離權(quán)重法［19］是對采樣點進(jìn)行線性的加權(quán)來決定輸出的柵格值，加權(quán)與距離成反比，輸入點離輸出柵格越遠(yuǎn)，它對輸出柵格的影響越小。使用公式表示為:

式中，Z為待估計的土壤有機碳含量值柵格值，z(xi)為第i(i=1，2，…，n)個采樣點的土壤有機碳含量值，n為用于土壤有機碳含量值插值的采樣點個數(shù)，di為插值點到第i個采樣點的距離，p為距離的冪。在本文中，p取為2，即反距離平方插值。

(3)徑向基函數(shù)插值

徑向基函數(shù)插值［18］是假設(shè)有一個實值函數(shù)F=F(x)，要逼近其在某些點x處的函數(shù)值。用徑向基函數(shù)來完成逼近的方法是構(gòu)造一個具有以下形式的逼近函數(shù)S(x):

其中，φ(t)，t≥0，是一個確定的實值函數(shù)，即徑向基函數(shù)，||·||表示歐幾里得距離，ai，i=1，…，N 是待定系數(shù)。點 xi，i=1，…，N為徑向基函數(shù)逼近的節(jié)點。相比于其他需要復(fù)雜處理的方法，徑向基函數(shù)曲面的簡單結(jié)構(gòu)使得徑向基函數(shù)逼近方法更加具有吸引力。

(2)結(jié)合土地利用類型修正的克里格法(OK_LU)

以采樣點分布圖為底圖，導(dǎo)出樣點土壤有機碳含量與土地利用類型的關(guān)系數(shù)據(jù)表，并利用SPSS18.0軟件對數(shù)據(jù)進(jìn)行統(tǒng)計分析，分別計算平均值、標(biāo)準(zhǔn)差、變異系數(shù)和相關(guān)系數(shù)，并利用單因素方差分析中的LSD進(jìn)行多重比較。明確土地利用類型對土壤有機碳含量分布的影響后，在普通克里格插值方法的基礎(chǔ)上，結(jié)合土地利用類型修正的克里格法即結(jié)合土地利用類型因子進(jìn)行修正，具體分兩步進(jìn)行:

第一步，先將所有預(yù)測點按照土地利用類型分為6類，求出各土地利用類型的平均土壤有機碳含量值。將每一個樣本的土壤有機碳含量值Z(Xkj)分為相同類型均值μ(Tk)和離差R(Xkj)之和。用公式表示為:

式中，Z(Xkj)為樣品的土壤有機碳含量值;μ(Tk)為相同土地利用類型樣本的土壤有機碳含量均值;R(Xkj)為樣本土壤有機碳含量與其相同土地利用類型樣品的均值之差，稱為“離差”。研究將離差作為一個新的區(qū)域變量R(Xkj)進(jìn)行普通克里格插值，空間插值利用ArcGIS10.0軟件中的地統(tǒng)計分析模塊完成。

第二步，將插值過后的離差值圖層A按土地利用類型的不同，分別加上其地類SOC均值，得到SOC初步預(yù)測值圖層B。各待估點的土壤有機碳含量預(yù)測值為各土地利用類型均值與離差估計值之和。由此得到研究區(qū)域土壤有機碳含量空間插值結(jié)果。

1.3.2 檢驗方法

為檢驗不同插值方法在研究區(qū)域的預(yù)測精度，將188個實測樣本在Arcgis10.0地統(tǒng)計工具下的Subset中隨機構(gòu)建兩個子集，75%的樣本(N=141)做預(yù)測子集，25%(N=47)的樣本做驗證子集(圖1)。本研究對驗證樣點的預(yù)測值與相同位置的實測值進(jìn)行相關(guān)性分析，對SOC的空間插值結(jié)果進(jìn)行精度檢驗。用預(yù)測子集中的樣點進(jìn)行插值預(yù)測，得到驗證子集中樣點的SOC值，然后通過計算所有驗證樣點的預(yù)測值和觀測值的Pearson相關(guān)系數(shù)(R)，平均絕對誤差(MAE)，均方根誤差(RMSE)，準(zhǔn)確度(AC)作為檢驗插值精度的標(biāo)準(zhǔn)，具體公式如下［20］:

式中，Voi是驗證點土壤有機碳含量的的實測值，Vpi是驗證點土壤有機碳含量的值，n是驗證點數(shù)量。

其中

式中，n為樣本數(shù);PE為可能的誤差變化(Potential error variance)，Pi和Oi分別為觀測值和預(yù)測值;O為實測值的平均值。AC取值范圍為0到1，1代表預(yù)測值和實測值完全一致，而0代表預(yù)測值和實測值完全不同。

2 研究結(jié)果

2.1 描述性統(tǒng)計分析

由表2可知，羊圈溝流域各土地利用類型中土壤有機碳含量平均值的順序為草地＞灌木林地＞林地＞壩地＞果園＞梯田。不同的土地利用方式下，土壤有機碳含量值存在顯著性差異，在探究土壤有機碳含量信息時，有必要考慮土地利用類型為輔助因子，分地類進(jìn)行插值。

表2 不同土地利用方式下土壤有機碳含量狀況Table 2 Soil organic carbon content in different land-use type

2.2 SOC空間結(jié)構(gòu)分析

空間自相關(guān)理論是地理統(tǒng)計插值法進(jìn)行預(yù)測的理論基礎(chǔ)［21］，所謂空間自相關(guān)是指在一定區(qū)域單元上，某種地理現(xiàn)象或某一屬性值與臨近區(qū)域單元上同一現(xiàn)象或?qū)傩灾档南嚓P(guān)性，它用于衡量空間變量的變異對鄰近區(qū)域相同變量的變異的依賴程度［22-23］，半變異函數(shù)是地統(tǒng)計學(xué)中研究土壤變異性的關(guān)鍵函數(shù)，是用來描述土壤性質(zhì)的空間連續(xù)變異的一個連續(xù)函數(shù)，反映土壤性質(zhì)的不同距離觀測值之間的變化。存在空間依賴性的數(shù)據(jù)點對，距離越近，其值越相近;反之越遠(yuǎn)，其值越不同［24］。半變異函數(shù)中，N表示其塊金值，S表示基臺值，本文應(yīng)用N/S和S/(N+S)量化原始的SOC數(shù)據(jù)(應(yīng)用于OK)及SOC離差數(shù)據(jù)(應(yīng)用于OK_LU)的空間結(jié)構(gòu)強度。N/S比值按空間結(jié)構(gòu)強度的大小分為三類(N/S＞0.6，0.3＜N/S＜0.6，N/S＜0.6):N/S＜30%，表示強烈的空間相關(guān)性;0.3＜N/S＜0.6，表示空間相關(guān)性較溫和;N/S＞0.6，表示較弱空間相關(guān)性［25］。S/(N+S)比率越接近于1，顯示其具有越強的空間相關(guān)性，反之，比率越接近于0顯示其空間相關(guān)性越弱［26］。OK和OK_LU的半變異函數(shù)模型如圖2所示。

2.3 插值結(jié)果的驗證與分析

圖3中別為運用OK、IDW、RBF和OK_LU插值法得到的羊圈溝流域SOC空間分布圖。從整體上看，四圖均反映出在流域的西北部SOC含量較高，流域南部含量低。然而，各插值方法的結(jié)果圖在局部上存在很大差異。由圖3可以看出，OK的預(yù)測值較連續(xù)，有機碳含量值分布很均勻，但是沒有體現(xiàn)出不同土地利用類型的差異;圖3中，IDW和RBF插值的結(jié)果圖出現(xiàn)了“牛眼”現(xiàn)象，這是由于樣點空間分布不均勻?qū)е碌模?0］，反距離權(quán)重法的插值表面不平滑，“牛眼”現(xiàn)象較為嚴(yán)重;樣條函數(shù)法的插值牛眼現(xiàn)象較反距離權(quán)重法有所減輕，但表面仍不光滑;結(jié)合圖3和圖2的2011年羊圈溝流域的土地利用現(xiàn)狀圖可以發(fā)現(xiàn)，OK_LU插值表面比較連續(xù)、光滑，在林地、灌木林地和草地地區(qū)體現(xiàn)了較高的土壤有機碳含量值，而在梯田、壩地、裸地表現(xiàn)出的土壤有機碳含量值較低，更詳細(xì)具體地體現(xiàn)了不同土地利用方式對于土壤有機碳空間分布的影響，對于流域尺度的土壤有機碳的空間變化，是一種比較理想的插值方法。

圖2 黃土丘陵區(qū)羊圈溝流域土壤有機碳空間插值半變異函數(shù)圖Fig．2 Semi-variagram map of spatial interpolation of SOC in YangJuanGou watershed of hilly region，Loess Plateau

如表3所示，SOC的N/S比值(0.666)大于0.6，OK_LU離差的N/S比值(0.485)在0.3—0.6之間，SOC離差的S/(N+S)值(0.8175)大于SOC的S/(N+S)值(0.7997)，即SOC和SOC離差均具有較強的空間相關(guān)性，其中SOC離差的空間相關(guān)性更強，表明OK和OK_LU兩種插值方法在此研究區(qū)域內(nèi)較為適宜。

表3 半變異函數(shù)模型Table 3 Semi-variagram models

從表4可以看出，4種方法的實測值與預(yù)測值均具有一定的相關(guān)性，IDW的R=0.3146(P=0.0397)和RBF的R=0.3012(P=0.0385)是顯著相關(guān);OK的R=0.4026(P=0.0050)和OK_LU的R=0.5589(P=0.0000)是極顯著相關(guān)。同時，通過均方根誤差(MAE)、平均絕對誤差(RSME)和準(zhǔn)確度(AC)來衡量4種預(yù)測方法的系統(tǒng)誤差及預(yù)測結(jié)果的準(zhǔn)確性。OK_LU平均絕對誤差和均方根誤差(MAE=0.7930，RSME=0.2136)均小于 OK(MAE=0.8665，RSME=0.2927)、IDW(MAE=0.8665，RSME=0.2927)、RBF(MAE=0.8665，RSME=0.2927);而 OK_LU(R=0.5589，AC=0.9505)相關(guān)性和準(zhǔn)確度均高于 OK(R=0.4026，AC=0.9081)、IDW(R=0.3146，AC=0.0397)、RBF(R=0.3012，AC=0.0385)。可以看出，在對土壤有機碳進(jìn)行空間插值時(特別是在黃土丘陵溝壑地區(qū))結(jié)合土地利用類型，能使插值結(jié)果的精度得到很大的提高，是獲得精確的土壤有機碳空間分布圖是關(guān)鍵因素。

表4 驗證樣點實測值與估算值相關(guān)性檢驗Table 4 The correlation test of the measured value and the estimated value

3 討論

3.1 OK_LU的優(yōu)越性

地統(tǒng)計法在地形變化不十分劇烈、比較平坦的平原地區(qū)可以得到較好的探測和模擬，其不僅考慮了待預(yù)測點與鄰近樣點數(shù)據(jù)間的空間距離關(guān)系，還考慮了各參與預(yù)測的樣點間的位置關(guān)系，充分利用了各樣點數(shù)據(jù)的空間分布結(jié)構(gòu)特征，地理統(tǒng)計方法的這種優(yōu)越性對于常常出現(xiàn)樣本不充分、采點未經(jīng)優(yōu)化等問題的土壤屬性數(shù)據(jù)的空間插值來說尤為重要［22］。本研究區(qū)的土壤有機碳數(shù)據(jù)來自隨機采樣，OK使其插值結(jié)果比IDW、RBF等傳統(tǒng)方法更精確，更有效避免了系統(tǒng)誤差的出現(xiàn)。

圖3 黃土丘陵區(qū)羊圈溝流域土壤有機碳空間分布圖Fig．3 Distribution of SOC in Yangjuangou watershed of hilly region，Loess Plateau

然而在黃土丘陵溝壑區(qū)，OK的應(yīng)用卻有一定的局限性［10］。由本文分析可知，在研究區(qū)內(nèi)土地利用類型是影響SOC空間分布差異的重要因素，因此結(jié)合土地利用類型對SOC空間插值修正是非常必要的。該區(qū)地形破碎、土地利用類型多樣，基于半變異函數(shù)的克里格插值難以將復(fù)雜的地形參數(shù)用于研究SOC空間變化。本文通過對普通克里格插值方法進(jìn)行土地利用類型修正，在黃土丘陵溝壑區(qū)對SOC空間格局的模擬得到了較好效果。范勝龍［11］等以福建省龍海市為例，設(shè)計了結(jié)合地貌類型、土壤類型和土地利用類型等信息的克里格插值模型，結(jié)果表明設(shè)計的插值模型預(yù)測精度均高于普通克里格法，進(jìn)一步驗證了OK_LU的優(yōu)越性。

3.2 土壤屬性數(shù)據(jù)插值方法的選擇

由于受到采樣點密度、地表的復(fù)雜程度等影響，不同的插值方法生成的結(jié)果可能大相徑庭［22］，沒有任何一種插值方法適合所有地區(qū)，應(yīng)依據(jù)對數(shù)據(jù)的定性分析、對研究區(qū)的地形分析，及對數(shù)據(jù)的空間探索分析，經(jīng)過反復(fù)實驗，選擇最優(yōu)的空間插值方法［26-27］。

若研究目的僅在于對研究區(qū)域內(nèi)幾個粗略劃分的亞區(qū)進(jìn)行比較，如對黃土丘陵區(qū)最高值和最低值的區(qū)域進(jìn)行大致定位，則此時較小尺度的空間變異可以被忽略，只需把握住整體特征和全局趨勢，以實施一系列插值優(yōu)化措施為前提，本研究選用的4種方法均可符合要求。若SOC數(shù)據(jù)樣點密度足夠大、研究區(qū)域的地形變化不大時，可以直接采用OK進(jìn)行插值，插值方法簡便快捷，預(yù)測結(jié)果準(zhǔn)確性較高;但使用OK插值相對其他地統(tǒng)計方法，步驟繁多，參數(shù)設(shè)置要求較高，耗時耗力，在研究目的許可的情況下，簡單、靈活的IDW、RBF也是合適的選擇。如需利用插值結(jié)果進(jìn)行更細(xì)致，更高級的空間分析，如圖形疊加、構(gòu)建生態(tài)系統(tǒng)模型等，而數(shù)據(jù)具有某種變化趨勢，則必須選擇插值精度最高、在小空間尺度表現(xiàn)最好的插值方法，可以考慮將某種因素對克里格插值進(jìn)行修正，剔除趨勢再進(jìn)行預(yù)測［28］，如本研究中的OK_LU，因為誤差會發(fā)生上行傳遞和放大，造成更深遠(yuǎn)的影響。

4 結(jié)論

在不同研究區(qū)域應(yīng)根據(jù)不同研究目的選擇不同的插值方法，不同插值方法各有利弊。當(dāng)常規(guī)方法不能滿足預(yù)測精度時，可通過輔助信息提高目標(biāo)變量的預(yù)測精度。在本研究區(qū)內(nèi)，普通克里格法是3種常規(guī)方法中插值效果最好的方法，徑向基函數(shù)法一般，反距離權(quán)重法插值效果最差。黃土丘陵溝壑區(qū)生態(tài)環(huán)境脆弱，地形較復(fù)雜，土地利用類型多樣，不同土地利用類型的土壤有機碳含量平均值差異明顯，基于土地利用類型修正的普通克里格插值法對SOC含量空間分布的預(yù)測精度有一定程度的提高，是地形復(fù)雜的黃土丘陵溝壑區(qū)土壤有機碳空間插值較適宜的方法。

［1］ Post W M，Peng T H，Emanuel W R，King A W，Dale V H，De Angelis D L.The global carbon cycle.American Scientist，1990，78:310-326.

［2］ Batjes N H.Total carbon and nitrogen in the soils of the world.European Journal of Soil Science，1996，47(2):151-163.

［3］ Lal R.Soil carbon sequestration impacts on global climate change and food security.Science，2004，304(5677):1623-1627.

［4］ Post W M，Kwon K C.Soil carbon sequestration and land-use change:processes and potential.Global Change Biology，2000，6(3):317-327.

［5］ Lal R，Kimble J，Levine E，Whitman C.World soils and greenhouse effect:an overview//Lal R，Kimble J，Levine E，Stewart B A，eds.Soils and Global Changes.Boca Raton:CRC Lewis Publishers，1995:1-25.

［6］ Hu K L，Li B G，Lin Q M，Li G T，Chen D L.Spatial variability of soil nutrient in wheat field.Transactions of the Chinese Society of Agricultural Engineering，1999，15(3):33-38.

［7］ Wang J，Qiu Y.Review on spatial variability and scale effects of land quality.Progress in Geography，2005，24(4):28-35.

［8］ Zhang S J，F(xiàn)ang H，He Y.Sampling strategies of field information on precision agriculture.Transactions of the Chinese Society for Agricultural Machinery，2004，35(4):88-92.

［9］ Sauer T J，Meek D W.The spatial variation of plant-available phosphorus in pastures with contrasting management.Soil Science Society of America Journal，2003，67(3):826-836.

［10］ Yao X L，F(xiàn)u B J，Lü Y H，Wang S，Liu M.Comparison of four spatial interpolation methods for estimating soil moisture in a complex terrain catchment.PLoS One，2013，1(8):1-13.

［11］ Zhang T C，Chang Q R，Liu J.Comparison of spatial interpolation methods for soil nutrient elements A case study of Yuyang County，Shaanxi Province.Agricultural Research in the Arid Areas，2010，28(2):177-182.

［12］ Liu X N，Huang F，Wang P.Principles and Methods of GIS Spatial Analysis.Beijing:Science Press，2008:185-223.

［13］ Liu J，Cai G X，Liu H B.Research of the interpolation methods for determining soil organic matter content in the hilly regions of southwestern China.Journal of Southwest University:Natural Science Edition，2008，30(3):107-112.

［14］ Fan S L，Huang Y H，Lin J S.The optimized interpolation models and its relationship with soil sampling density on detecting spatial variability of farmland soil organic carbon A case study in Longhai city，F(xiàn)ujian province.Research of Soil and Water Conservation，2011，18(6):1-5.

［15］ Allmaras R R，Schomberg H H，Douglas C L，Dao T H Jr.Soil organic carbon sequestration potential of adopting conservation tillage in U.S.croplands.Journal of Soil and Water Conservation，2000，55(3):365-373.

［16］ Chen J Y，Sun B，Li Z P，Li H Y，Pan J J.Pool size of soil organic carbon and dynamics under different land use.Journal of Soil and Water Conservation，2008，22(1):91-95，179-179.

［17］ Wang Y F，F(xiàn)u B J，Chen L D，Lü Y H，Luo C Y.Effects of land use change on soil erosion intensity in small watershed of Loess Hilly Region:A quantitative evaluation with 137-Cesium tracer.Chinese Journal of Applied Ecology，2009，20(7):1571-1576.

［18］ Fang X，Xu G L，Hong Y，Jing W F.Study on soil organic carbon density under different land use types in central hilly area of Hu'nan province.Journal of Soil and Water Conservation，2009，23(1):54-58，63-63.

［19］ Zhu H Y，Liu S L，Jia S F.Problems of the spatial interpolation of physical geographical elements.Geographical Research，2004，23(4):425-432.

［20］ Yang L，Zhu A X，Qin C Z，Li B L，Pei T，Liu B Y.Soil property mapping using fuzzy membership A case study of a study area in Heshan farm of Heilongjiang province.Acta Pedologica Sinica，2009，46(1):9-15.

［21］ Lin L，Li C H，Dai M，Cai W G，Lin Q.Optimization of the spatial interpolation for marine phytoplankton abundance.Acta Ecologica Sinica，2007，27(7):2880-2888.

［22］ Gringarten E，Deutsch C V.Teacher's aide variogram interpretation and modeling.Mathematical Geology，2001，33(4):507-534.

［23］ Li X，Cheng G D，Lu L.Comparison study of spatial interpolation methods of air temperature over Qinghai-Xizang plateau.Plateau Meteorology，2003，22(6):565-573.

［24］ Leenaers H，Okx J P，Burrough P A.Comparison of spatial prediction methods for mapping floodplain soil pollution.Catena，1990，17(6):535-550.

［25］ Cambardella C A，Karlen D L.Spatial analysis of soil fertility parameters.Precision Agriculture，1999，1(1):5-14.

［26］ Lam N S N.Spatial interpolation methods:A review.American Cartographer，1983，10(2):129-149.

［27］ Liu J S，Chen H，Yang B Y，Wang W，Xiang Y，Zhao C.Comparison of interpolation methods on annual mean precipitation in Hebei Province.Acta Ecologica Sinica，2009，29(7):3493-3500.

［28］ Lado L R，Polya D，Winkel L，Berg M，Hegan A.Modelling arsenic hazard in Cambodia:A geostatistical approach using ancillary data.Applied Geochemistry，2008，23(11):3010-3018.

參考文獻(xiàn):

［6］ 胡克林，李保國，林啟美，李貴桐，陳德立.農(nóng)田土壤養(yǎng)分的空間變異性特征.農(nóng)業(yè)工程學(xué)報，1999，15(3):33-38.

［7］ 王軍，邱揚.土地質(zhì)量的空間變異與尺度效應(yīng)研究進(jìn)展.地理科學(xué)進(jìn)展，2005，24(4):28-35.

［8］ 張淑娟，方慧，何勇.精細(xì)農(nóng)業(yè)田間信息采樣策略.農(nóng)業(yè)機械學(xué)報，2004，35(4):88-92.

［11］ 張鐵嬋，常慶瑞，劉京.土壤養(yǎng)分元素空間分布不同插值方法研究——以榆林市榆陽區(qū)為例.干旱地區(qū)農(nóng)業(yè)研究，2010，28(2):177-182.

［12］ 劉湘南，黃方，王平.GIS空間分析原理與方法.北京:科學(xué)出版社，2008:185-223.

［13］ 劉靜，蔡國學(xué)，劉洪斌.西南丘陵地區(qū)土壤有機質(zhì)含量的空間插值法研究.西南大學(xué)學(xué)報:自然科學(xué)版，2008，30(3):107-112.

［14］ 范勝龍，黃炎和，林金石.表征土壤有機碳區(qū)域分布的優(yōu)化空間插值模型研究——以福建省龍海市為例.水土保持研究，2011，18(6):1-5.

［16］ 陳錦盈，孫波，李忠佩，李海鷹，潘劍君.不同土地利用類型土壤有機碳各庫大小及分解動態(tài).水土保持學(xué)報，2008，22(1):91-95，179-179.

［17］ 汪亞峰，傅伯杰，陳利頂，呂一河，羅春燕.黃土丘陵小流域土地利用變化的土壤侵蝕效應(yīng):基于137Cs示蹤的定量評價.應(yīng)用生態(tài)學(xué)報，2009，20(7):1571-1576.

［18］ 方晰，徐桂林，洪瑜，金文芬.湘中丘陵區(qū)不同土地利用方式下土壤有機碳密度.水土保持學(xué)報，2009，23(1):54-58，63-63.

［19］ 朱會義，劉述林，賈紹鳳.自然地理要素空間插值的幾個問題.地理研究，2004，23(4):425-432.

［20］ 楊琳，朱阿興，秦承志，李寶林，裴韜，劉寶元.運用模糊隸屬度進(jìn)行土壤屬性制圖的研究——以黑龍江鶴山農(nóng)場研究區(qū)為例.土壤學(xué)報，2009，46(1):9-15.

［21］ 林琳，李純厚，戴明，蔡文貴，林欽.海洋浮游植物豐度的空間插值優(yōu)化.生態(tài)學(xué)報，2007，27(7):2880-2888.

［23］ 李新，程國棟，盧玲.青藏高原氣溫分布的空間插值方法比較.高原氣象，2003，22(6):565-573.

［27］ 劉勁松，陳輝，楊彬云，王衛(wèi)，相云，趙超.河北省年均降水量插值方法比較.生態(tài)學(xué)報，2009，29(7):3493-3500.