水田土壤有機碳時空演變下的最優(yōu)插值方法

喬 婷, 姚彩燕, 于東升, 史學正, 邢世和, 張黎明

(1.福建農(nóng)林大學資源與環(huán)境學院，福建 福州 350002;2.土壤生態(tài)系統(tǒng)健康與調控福建省高校重點實驗室，福建 福州 350002;3.中國科學院南京土壤研究所土壤與農(nóng)業(yè)可持續(xù)發(fā)展國家重點實驗室，江蘇 南京 210008)

土壤在全球碳循環(huán)中扮演著重要的角色，全球土壤有機碳(soil organic carbon, SOC)僅變化0.1%就會給大氣二氧化碳濃度造成劇烈波動[1].農(nóng)業(yè)土壤碳庫是陸地生態(tài)系統(tǒng)的重要組成部分，具有易受人為干擾、短時間內可快速調節(jié)的特點，其固碳效應可有效延緩全球氣候變化[2].農(nóng)業(yè)土壤中的水田由于特殊的淹水條件，固碳效應明顯大于旱地[3].我國水田面積超過330 000 km2[4]，占全國糧食作物播種面積的27.4%，占世界水田總面積的23%.太湖地區(qū)水稻土分布廣泛且剖面發(fā)育好，是長江中下游水稻土區(qū)中的一個典型[2]，明確該土地利用類型的SOC動態(tài)變化對于合理制定全國農(nóng)業(yè)“碳匯”管理政策尤為重要[4-5].由于受土壤類型、母質、地形、氣候、生物和農(nóng)業(yè)管理方式的影響，SOC存在明顯的時空變異.有研究表明，通過空間插值模型可有效定量化連續(xù)土壤屬性的時空變異[6].按照不同理論可將常用空間插值法歸為非地統(tǒng)計、地統(tǒng)計、互結合3大類，根據(jù)數(shù)學原理又可分為確定性插值和地統(tǒng)計插值兩大類[7].非地統(tǒng)計主要包括徑向基函數(shù)、反距離權重、樣條函數(shù)和趨勢面分析等確定性方法，以及最近鄰函數(shù)、TIN三角網(wǎng)和回歸樹等方法；地統(tǒng)計主要包括普通克里金、簡單克里金、泛克里金、析取克里金和指示克里金等；而互結合是基于前兩者的拓展方法，典型方法有回歸克里金和趨勢面分析結合克里金方法等.

目前，很多學者基于不同插值方法對土壤有機碳等屬性進行了研究[8-12].如Bhunia et al[8]采用反距離加權、局部多項式、徑向基函數(shù)、普通克里金和經(jīng)驗貝葉斯克里金5種插值方法，分析了印度梅迪尼浦爾地區(qū)同時期3個土壤深度(0～20、20～40和40～100 cm)下SOC的空間變化，表明普通克里金對于表層土是一種均方根誤差最小和決定系數(shù)最高的較優(yōu)方法.Addis et al[9]使用普通克里金、反距離權重和徑向基函數(shù)插值方法，發(fā)現(xiàn)在埃塞俄比亞塔納流域農(nóng)業(yè)區(qū)普通克里金下SOC的插值圖效果最好，且認為沒有一種特定的方法明顯適用于pH、質地和有效磷等指標.Li et al[10]對中國西部九臺縣同年度的土壤有機質、總氮、有效磷、有效鉀進行普通克里金、反距離權重、回歸克里金插值，發(fā)現(xiàn)普通克里金對有機質和有效鉀插值效果好.Long et al[11]基于廣泛使用的12種插值方法對2008年農(nóng)業(yè)農(nóng)村部配方施肥項目的235 309個采樣點的SOC進行了福建省復雜地形區(qū)的適宜性研究，結果發(fā)現(xiàn)簡單克里金對復雜地形、丘陵山區(qū)和河谷盆地地區(qū)插值效果好，反距離權重適用于平原平臺地區(qū).馬利芳等[12]運用普通克里金、反距離權重、徑向基函數(shù)和局部多項式法預測了不同人為干擾程度下干旱區(qū)土壤有機質的空間分布，結果表明普通克里金在無人為干擾的地區(qū)效果最好，而徑向基函數(shù)在重度人為干擾區(qū)效果最好.

從以上研究可以看出，各插值法在不同區(qū)域的適用性有很大差異[8-13]；其次，這些研究大多局限在單一時間尺度下，較少考慮不同時期土壤屬性差異對插值精度的影響；最后，地表條件和人為活動的差異將導致SOC產(chǎn)生明顯的時空變異，影響最優(yōu)插值模型的篩選[14].因此，為了準確估計SOC含量，基于樣點的適宜插值方法選擇時必須要考慮空間和時間維度.基于此，本研究以我國太湖地區(qū)23 200 km2水田土壤為研究對象，利用1982年第2次土壤普查和2000年土壤質量演變與持續(xù)利用“973”項目的實測樣點數(shù)據(jù)，分析目前土壤學常用的11種插值方法對不同時期水田SOC預測精度影響，篩選出整個地區(qū)、不同水稻土亞類和行政轄區(qū)的最優(yōu)插值模型，旨在揭示土壤屬性變異對最優(yōu)插值模型影響的內在機理，為我國南方水田SOC時空演變下的最優(yōu)插值模型篩選和合理制定碳匯動態(tài)管理措施提供理論依據(jù).

1 材料與方法

1.1 研究區(qū)概況

太湖地區(qū)(118°50′—121°54′E，29°56′—32°16′N)位于長江三角洲南緣，屬于亞熱帶季風氣候，氣候溫和，水熱條件優(yōu)越；年太陽輻射量高達4 686 MJ·m-2，年均氣溫15～17 ℃，年均降水量1 177 mm，農(nóng)田管理主要為水稻和冬小麥輪作制度，是中國傳統(tǒng)糧食高產(chǎn)區(qū)[5].該地區(qū)總面積約36 900 km2，其中66%被水稻土覆蓋[15].

本研究區(qū)分屬兩省一市，包括浙江省、江蘇省和上海市3個轄區(qū)，共覆蓋37個縣(市、區(qū)).根據(jù)成土母質、土地利用方式、耕作制度和農(nóng)業(yè)利用分為4個主要土區(qū)(圖1a)，分別為以湖(河、海)相沉積物、長江沖積物形成的太湖平原土區(qū)(面積占比17%)、沖積平原土區(qū)(面積占比26%)、各種巖石殘坡積物及下蜀黃土形成的丘陵低山土區(qū)(面積占比31%)和沼澤土形成的低洼圩田土區(qū)(面積占比27%)，并逐漸發(fā)育為潴育、滲育、脫潛、漂洗、淹育和潛育水稻土等6個水田土壤亞類[16](圖1b).本研究選取太湖地區(qū)水田土壤面積最大的3個主要亞類作為研究對象，分別為潴育水稻土、脫潛水稻土和滲育水稻土，面積占比分別53%、18%和16%.

1.2 數(shù)據(jù)來源

1982年土壤屬性數(shù)據(jù)來自浙江省、江蘇省和上海市第2次土壤普查的各縣(市、區(qū))土壤志，共計1 096個水田土壤剖面樣點；2000年土壤屬性數(shù)據(jù)來自中國科學院南京土壤研究所主持的“973土壤質量研究”項目，共計1 370個水田土壤表層樣點.根據(jù)Bemmelen轉換系數(shù)計算得到SOC含量后，為保證歷史數(shù)據(jù)的真實性和最大限度被利用，通過5S法(均值加減五倍標準差)剔除異常值，當異常值比總體均值大(小)5倍標準差以上(下)且比近鄰8個樣點均值大(小)3倍標準差以上(下)時，使用該規(guī)則下的最大最小值代替[17].此外，基于ArcGIS地統(tǒng)計分析Create Subset模塊按照4∶1的比例將兩期樣點集隨機分成預測集與驗證集，得到太湖整個地區(qū)1982年預測集樣點共877個，驗證集219個，2000年預測集樣點共1 096個，驗證集274個(圖2)[18].

1.3 數(shù)據(jù)探索與空間插值模型擬合

基于ArcGIS 10.2 Geostatistacal Analyst模塊下常用的7種確定性插值方法[反距離權重(inverse distance weighting、IDW)、全局多項式(global polynomial interpolation、GPI)、局部多項式(local polynomial interpolation、LPI)、規(guī)則樣條函數(shù)(completely regularized spline、CRS)、張力樣條函數(shù)(spline with tension、SWT)、高次曲面函數(shù)(multiquadric function、MF)、反高次曲面函數(shù)(inverse multiquadric function、IMS)]，以及4種地統(tǒng)計方法[普通克里金(ordinary Kriging、OK)、簡單克里金(simple Kriging、SK)、泛克里金(universal Kriging、UK)和析取克里金(disjunctive Kriging、DK)]進行插值[18].IDW是基于相近相似原理，以插值點與樣本點之間的距離為權重進行加權平均的精確插值方法；GPI和LPI是分別利用1個低階和多個的多項式進行擬合的不精確插值方法；同屬于徑向基函數(shù)的CRS、SWT、MF、IMS，基本原理是經(jīng)過每一個已知觀測點、形成總曲率最小的橡膠薄膜狀預測表面，而不同的方法則對應膜擬合點集的不同方式[18].地統(tǒng)計的OK、SK、UK和DK是以變異函數(shù)理論和結構分析為基礎，在有限區(qū)域內對區(qū)域化變量進行無偏最優(yōu)估計，不同方法對應點集適用的條件不同[19].

為滿足地統(tǒng)計插值方法正態(tài)分布的前提與平穩(wěn)性理論假設，在SPSS、ArcGIS和GS+9.0中進行正態(tài)檢驗及空間自相關分析[18-21].由于本研究各分類預測集均為超過100的大樣本，利用效果較理想的偏度峰度值聯(lián)合K-S檢驗對各預測集進行正態(tài)分布性檢驗，若K-S方法中Z值的相伴概率值大于當前顯著性水平，認為樣本來自的總體與指定正態(tài)分布無顯著差異；否則需進行轉換.其次，由于土壤樣品采樣點離散分布會掩蓋空間分布相關性，需要進行離散度與相關性分析以避免出現(xiàn)純塊金效應，造成克里金插值估值僅為簡單的數(shù)學平均值[19].全局莫蘭I指數(shù)可驗證區(qū)域整體空間模式和樣本空間自相關性，指數(shù)I=0表示隨機分布，I>0表示存在空間正相關，且關系顯著時表示變量存在明顯空間聚集現(xiàn)象，數(shù)值越大關系越密切，性質越相似；I<0表示空間異常，數(shù)值絕對值越大說明差異越大越分散[20-21].最后，在GS+9.0模擬得到半變異方差擬合結果[18-19].半變異擬合模型參數(shù)平方殘差和RSS越小越好，模型決定系數(shù)R2越大越好，且優(yōu)先考慮RSS.

1.4 預測結果精度檢驗

通過預測集數(shù)據(jù)估測自相關模型的交叉驗證方法確定某一種插值方法的最優(yōu)參數(shù)設置；然后通過外部驗證對最優(yōu)參數(shù)設置的不同插值模型進行精度評價，基于“平均預測誤差(ME)越趨近于0、均方根誤差(RMSE)越小且優(yōu)先考慮后者”的原則，得到1982年和2000年太湖地區(qū)水田SOC最適宜插值方法[7,18,22].

(1)

(2)

(3)

式中，N為采樣點數(shù)量，Xoi為采樣點SOC含量實測值，Xpi為采樣點SOC含量預測值，R為實測值和預測值之間的Pearson相關系數(shù)，其越大表明模型對各樣點集含量的解釋能力越強.

2 結果與分析

2.1 太湖地區(qū)水田SOC樣點基本統(tǒng)計特征

由表1可得，1982年和2000年太湖水田SOC含量(0～20 cm)范圍分別在1.91～37.39和1.91～44.59 g·kg-1之間，均值為(15.95±5.32)和(18.50±5.31) g·kg-1.從變異系數(shù)來看，兩個時期太湖水田SOC含量的離散程度相差不大，均在10%～100%之間，為中等變異[23].從整個太湖地區(qū)和不同水稻土亞類、行政轄區(qū)各預測集樣點SOC的偏度、峰度和K-S檢驗值來看，1982年和2000年研究區(qū)各預測集均符合正態(tài)分布，全局莫蘭I指數(shù)的范圍為0.146～0.448，對應Z得分范圍是6.45～35.43，Z高于2.58且顯著性均小于0.01，表明兩期預測集的SOC含量有明顯的空間聚集現(xiàn)象，具有空間正相關模式，符合地統(tǒng)計學插值前提[21].

表1 太湖地區(qū)水田SOC含量描述性統(tǒng)計特征1)

第2次土壤普查以來，太湖地區(qū)整個區(qū)域水田SOC的含量平均升高了2.55 g·kg-1，這與Liu et al[22]分析的1982—2000年太湖地區(qū)水稻土SOC含量呈上升趨勢的結果一致.太湖地區(qū)水稻土SOC含量上升的主要驅動因素是肥料的大量施用和降雨量的明顯增加.據(jù)統(tǒng)計，1982年太湖地區(qū)農(nóng)家肥和氮肥的平均施用量分別為231和243 kg·hm-2，而2000年為270和345 kg·hm-2，增幅分別達到19.86%和37.66%[24].其次，氣候因子也在太湖地區(qū)水田SOC提升中起著重要作用.據(jù)統(tǒng)計，該地區(qū)的年均降雨量由1982年的1 102 mm提高到2000年的1 175 mm，降雨量的增加有利于SOC積累[24].進一步統(tǒng)計表明，1982—2000年間整個地區(qū)水田SOC的變異系數(shù)增幅為-13.89%，而3個主要亞類和3個行政轄區(qū)變異系數(shù)增幅有正有負(潴育、滲育和脫潛水稻土分別為-13.56%、-4.15%和2.69%，江蘇、上海和浙江分別為-27.90%、-23.74%和7.26%)，可見過去19年間太湖地區(qū)水稻土SOC整體空間變異性變化緩和，但局域變異變化差異較大，這與第2次普查時期統(tǒng)一經(jīng)營轉變?yōu)榧彝ヂ?lián)產(chǎn)承包責任制后的農(nóng)田管理方式有關[14].

2.2 太湖地區(qū)水田SOC半變異擬合

1982年和2000年太湖地區(qū)水田SOC半變異函數(shù)模型決定系數(shù)R2均高于0.824(表2)，殘差平方和RSS均低于0.375，說明半變異擬合效果能夠充分反映各指標的空間結構特征.同時，半變異模型符合常見的指數(shù)、高斯和球狀模型，其中指數(shù)模型更適用于SOC含量擬合[18-19].

由表2可得，太湖整個水田SOC的塊金值(C0)在1982年為0.015，2000年為0.002，均小于對應年份各土壤亞類和各行政轄區(qū)的C0，且C0越大RSS誤差越大，說明SOC的C0與空間尺度為負相關，與空間分布隨機性為正相關[25].基底效應范圍為2.44%～51.81%，表示太湖地區(qū)SOC具有中強空間自相關性，同時受到施肥、耕作、種植制度等人為隨機因素和母質、土壤類型、氣候、地形等結構性因素影響.基底效應值越大說明空間自相關性越被削弱，空間變異越受到隨機效應的影響[26-27].同時有研究證明，在人為因素影響下土壤屬性變異表現(xiàn)出強烈的隨機性，其結構性和相關性被削弱，空間分布將朝均一化方向發(fā)展[10,28].

表2 兩期太湖區(qū)域不同土壤亞類和行政轄區(qū)水田SOC半方差模型1)

從變程變化來看，2000年太湖地區(qū)各水稻土亞類比1982年減少1.14×104～8.25×104m，行政轄區(qū)減少4.24×104～7.24×104m，而整個區(qū)域僅減少0.21×104m，表明過去19年間隨著SOC的積累，不同水稻土亞類和行政轄區(qū)的局部空間變異緩和，但對整個區(qū)域自相關范圍影響不大.

2.3 整個地區(qū)水田SOC的最優(yōu)插值方法

由圖3可得，1982年太湖整個地區(qū)不同插值方法下SOC的分布圖較為相似，總體上SOC含量均表現(xiàn)為中部和西南部區(qū)域含量較高，西北部、東南部和北部沿江區(qū)域含量較低，其他地區(qū)含量為中等水平；絕大部分地區(qū)處于12～14 g·kg-1及以上含量.隨著時間變換，2000年各對應插值方法的SOC空間分布格局與1982年相似，但極值差異縮小，高值區(qū)域分布增多，絕大部分地區(qū)處于14～16 g·kg-1及以上含量.從插值效果來看，7種確定性插值方法SOC分布差異很大，而4種地統(tǒng)計插值方法空間分布相似性高.其中，LPI和GPI呈現(xiàn)明顯的階梯狀變化，形成高值區(qū)到低值區(qū)的緩沖帶，GPI的空間分布圖在研究區(qū)西北邊角發(fā)生了扭曲變形.同屬于徑向基函數(shù)系列的IMS、MF、SWT和CRS大區(qū)域空間分布相似，但在小區(qū)域不同方法表現(xiàn)出不同程度的高低波動.IDW插值的空間分布圖像表面不平滑，“牛眼”現(xiàn)象嚴重.地統(tǒng)計插值OK、SK、UK和DK的空間分布連續(xù)且表面平滑，但局部的插值結果有明顯不同.

從不同確定性插值方法的精度評價來看(表3)，1982年和2000年太湖地區(qū)水田SOC均是IDW精度最高，而MF的預測值與實測值差異最大，說明IDW能較好地模擬研究區(qū)SOC的空間變異.進一步從插值函數(shù)適用性來看，由于1982年和2000年整個太湖地區(qū)SOC含量極值分別高達35.48和42.58 g·kg-1，基于人工神經(jīng)網(wǎng)絡的MF允許模擬值超出實測極值，短距離鄰域的SOC出現(xiàn)劇烈變異時易產(chǎn)生較高均值誤差，從而降低預測精度；但基于相近相似原理的IDW則是周圍實測點數(shù)量越多、距離越近，預測點位越受影響、值越相似，因此預測值更符合太湖地區(qū)SOC含量值[29].張海平等[30]在評估IDW和徑向基函數(shù)適用性時也認為在預測極值與采樣點一致的應用中IDW效果更好，這與本研究的結果相一致.

從不同地統(tǒng)計插值方法的精度評價來看(表3)，1982年太湖地區(qū)水田SOC最優(yōu)地統(tǒng)計方法為OK，2000年為DK.進一步從插值函數(shù)適用性來看，由于太湖地區(qū)總體地勢平坦，沖積平原和太湖平原占研究區(qū)總面積的55.17%，OK能較好地預測平坦地區(qū)SOC的細部變異.Addis et al[9]同樣認為最優(yōu)無偏估值的OK對農(nóng)業(yè)流域塔納盆地的SOC擬合效果最好.然而，與1982年相比，2000年太湖地區(qū)水稻土SOC的預測集樣點增加19.98%，而基于非線性連續(xù)型數(shù)據(jù)轉化的DK可充分利用大量樣點詳盡信息，能給出基于現(xiàn)有信息更精確的SOC概率分布估計值. Li et al[7]的研究也表明在低樣點密度下DK可能增加預測面平滑度從而顯著高估或低估預測誤差，而高密度的采樣有利于提高SOC的預測精度.

表3 太湖地區(qū)水田土壤整個區(qū)域插值模型預測精度評價1)

總體來看，1982年和2000年太湖地區(qū)水田SOC最優(yōu)插值方法分別是OK和DK.兩期中，地統(tǒng)計插值均為最優(yōu)的主要原因在于，此方法不僅考慮了SOC預測點位與實測點的距離，而且兼顧了SOC空間分布與空間方位關系，適用于預測結構性與隨機性并存的連續(xù)變量[25].同樣，在空間插值方法在不同領域適用性研究中，李海濤等[29]也提出克里金插值算法可適用于土壤樣本數(shù)據(jù)存在隨機性和結構性特征的場景.從整個太湖地區(qū)水田SOC含量變異來看，2000年的自相關程度和隨機因素影響程度比1982年有所提高，表現(xiàn)在2000年基底效應升高0.52%并且塊金值降低6.64倍，說明除植被、地形、氣候和母質等自然因素影響外，長期施肥等隨機人為因素引起的變異導致兩期最優(yōu)插值方法發(fā)生變化[24,26,28].

2.4 不同時期各水田土壤亞類SOC的最優(yōu)插值方法

對于面積占整個太湖地區(qū)水田53%的潴育水稻土SOC來說，1982年最優(yōu)插值方法為確定性插值的LPI，而2000年為地統(tǒng)計插值的SK(表4和表5).1982年77.06 %的樣點密集分布在低洼圩區(qū)和沖積平原土區(qū)交界處，SOC變異系數(shù)高達30.17%，而全局莫蘭I指數(shù)僅0.207，該時期的潴育水稻土與相鄰近區(qū)域SOC含量相似程度低且局部變異明顯，LPI可捕捉數(shù)據(jù)集的短程變異，因此更為適用[8].2000年太湖地區(qū)潴育水稻土的采樣點在各土區(qū)分布均勻，全局莫蘭I指數(shù)提高了0.127，變程則在各亞類中減少最多(8.25×104m) (表2)，盡管確定性插值IDW對自相關性密切且鄰值近似的SOC插值效果較好，但1982—2000年期間由于較高年降水量(1 219 mm)、年均溫度(16.6 ℃)和大量秸稈還田等因素的長期影響，潴育水稻土SOC的變異系數(shù)、塊金值和基底效應比1982年分別降低了4.09%、95.36%和33.05%；在這種情況下，兼顧平穩(wěn)性和空間關系的SK更適用于該亞類SOC的預測[6,24].史舟等[19]也認為線性SK比OK在二階平穩(wěn)的隨機過程中預測精度高，而長期人為影響造成水田SOC趨向均一，因此有利于提高簡單克里金的預測精度[28].

滲育水稻土占整個太湖地區(qū)水稻土總面積的16%，1982年和2000年該亞類SOC含量的最優(yōu)預測插值方法均為確定性插值，分別是GPI和LPI(表4和表5).這是因為1982年該亞類分布均勻，同時全局莫蘭I指數(shù)在所有土壤亞類中最高，表明該時期滲育水稻土SOC受局部影響弱且臨近區(qū)域相似度高，而GPI可忽略局部異常值而對整個預測表面進行整體平滑，較適用于分布范圍廣而屬性變化緩慢的樣點預測[26].但2000年滲育水稻土55.87%的采樣點分布在太湖地區(qū)北部沖積平原土區(qū)，SOC含量鄰域相似度弱且多局域存在較大差異，此時局部插值的LPI更適用；此外，1982—2000年間隨著經(jīng)濟的發(fā)展和農(nóng)業(yè)投入的增加，滲育水稻土亞類的化肥和有機肥年均施用量分別達到324和240 kg·hm-2，較高的肥料施用量一方面會使該亞類SOC含量升高，另一方面也會導致不同微地貌下SOC的空間變異增大，在這種情況下，LPI較其他插值方法能更準確地預測SOC含量[24,31].Bogunovic et al[32]對人為干擾強烈的農(nóng)村森林交界區(qū)土壤有機質預測時也得到LPI是優(yōu)于IDW、OK、SK的最好插值方法，該研究中基底效應高達59.56%，但變程僅是采樣間隔的8.64倍，一定程度上也證明了LPI適用于隨機干擾大、短距離自相關性強且空間屬性變異劇烈的農(nóng)耕區(qū).

表4 兩期太湖地區(qū)水田土壤主要亞類SOC含量預測誤差1)

表5 兩期太湖地區(qū)整個區(qū)域和主要亞類最優(yōu)插值方法預測精度對比

脫潛水稻土占整個太湖地區(qū)水稻土總面積的18%，介于潴育和滲育水稻土之間.1982年和2000年該亞類SOC最優(yōu)預測插值方法分別為地統(tǒng)計的DK和確定性插值的IDW(表4和表5).這主要是因為1982年滲育水稻土在所有亞類中基臺值最高，高比例隨機影響下，變程內的相鄰樣點之間可能無法滿足非簡單相關關系；在這種情況下非線性插值DK可根據(jù)概率分布估計SOC含量，提高預測精度[7].2000年研究區(qū)脫潛水稻土55.46%樣點分布在低洼圩區(qū)、33.61%分布在太湖平原與沖積平原土區(qū)，基底效應和全局莫蘭I指數(shù)大而變程小，對于樣點自相關性密切且鄰值近似、局部屬性弱波動的SOC，IDW預測精度較高.Long et al也證明在預測地勢平緩的平原臺地區(qū)SOC時，IDW的精度要高于地統(tǒng)計插值方法[11].

總體來看，1982年和2000年太湖地區(qū)潴育、滲育、脫潛3個主要水稻土亞類的SOC含量最優(yōu)空間預測方法隨著時間的推移發(fā)生了明顯變化，并且與整個太湖地區(qū)最優(yōu)方法不一致，這可能是因為亞類在不同時期的施肥量和氣候差異導致土壤屬性變異程度不同，進而使得適宜的插值模型產(chǎn)生了改變.

2.5 不同時期各行政區(qū)水田SOC的最優(yōu)插值方法

太湖地區(qū)浙江轄區(qū)水田1982年和2000年的SOC最優(yōu)插值方法分別為SK和CRS(表6和表7).這主要是因為1982年56.50%的采樣點位于低洼圩土區(qū)(海拔低于4 m)，其余樣點位于低山丘陵土區(qū)和沖積平原土區(qū)(海拔4.5～7.0 m)，鄰域地勢多變且SOC含量空間差異大；SK在SOC空間預測中易忽略細部極值影響而產(chǎn)生平緩預測面，適用于地形復雜區(qū)[11].2000年太湖地區(qū)浙江轄區(qū)50%的采樣點集中于平原區(qū)域，低山丘陵土區(qū)的樣點大幅減少，再加上大量施肥等隨機因素影響，基底效應升高8.78%，變程降低42 400 m；而CRS作為局部收斂的精確性插值器，更適用于樣點聚集度極高、局部變異大的SOC預測[31].有研究表明，與克里金和IDW模型相比，徑向基函數(shù)模型可以更真實地反映分析空間變量的結構[11].

表6 兩期太湖地區(qū)水田土壤行政轄區(qū)SOC含量預測誤差1)

表7 兩期太湖地區(qū)整個區(qū)域和行政轄區(qū)最優(yōu)插值方法預測精度對比

上海轄區(qū)水田1982年和2000年SOC最優(yōu)插值方法分別為OK和SK(表6和表7).這主要是因為1982年該地區(qū)59.05%樣點分布在東北部平原土區(qū)，其余分布于西南低洼圩區(qū)，區(qū)域內極值適中，OK對于地勢平坦的縣域土壤屬性插值效果更好.龍軍等[18]也認為OK方法可以較好預測平原縣域耕地SOC.但在2000年，屬于上海轄區(qū)的沖積平原土區(qū)采樣點比例比1982年提高了5.31%，尤其城市擴張造成崇明島和上海市區(qū)耕地大量減少形成更大范圍無樣點分布的空白區(qū)，SK可利用經(jīng)驗模型平滑預測面進行非平穩(wěn)結構的線性預測，尤其適用于這類微地形等結構性因素影響更多的情況.

江蘇轄區(qū)水田1982年和2000年SOC最優(yōu)插值方法均為OK(表6和表7).該轄區(qū)兩期樣點平原區(qū)的占比較高，分別達50.98%和56.77%.在不同確定性方法中，1982年江蘇轄區(qū)水田SOC的全局莫蘭I指數(shù)在各行政轄區(qū)中最高達到0.433，表明空間鄰域上SOC相似度高，但同時潴育和脫潛水稻土SOC含量極值差異大造成區(qū)域內部變異大，CRS可通過最小化表面曲率對平原區(qū)為主的區(qū)域實現(xiàn)較高精度的預測，因此為1982年最優(yōu)確定性插值方法[18].2000年，一方面江蘇轄區(qū)潴育水稻土樣點比例增加到了60.23%，另一方面SOC含量處于中高水平的潴育水稻土和脫潛水稻土一般分布于該轄區(qū)中部，而SOC含量處于低水平的滲育水稻土分布于該轄區(qū)邊緣部分；LPI在區(qū)域特征值聚集明顯情況下的預測精度較高，因此該方法在2000年確定性插值中最優(yōu).OK在弱化細節(jié)變異、有效地避免邊緣效應上優(yōu)于所有確定性插值方法，因此江蘇轄區(qū)水稻土兩時期SOC含量的最優(yōu)預測方法為OK[7].

總的來看，1982年和2000年江蘇、上海、浙江轄區(qū)水田SOC含量的預測地統(tǒng)計學方法的插值精度最高，適用性最好.很多研究也發(fā)現(xiàn)，地統(tǒng)計方法尤其是OK在SOC含量數(shù)字制圖中能夠平滑極值保證預測精度，同時也可展示較完善的細部變異情況[9-10].從時間變換來看，浙江轄區(qū)的水田SOC含量受1982年和2000年間樣點分布土區(qū)、施肥情況等因素影響，最優(yōu)預測方法從確定性插值轉變?yōu)榈亟y(tǒng)計學插值，而上海和江蘇轄區(qū)兩期最優(yōu)預測插值方法皆為地統(tǒng)計，原因可能與水稻土亞類屬性等結構性因素和耕作方式等隨機因素影響相似有關.從空間變換來看，兩期均方誤差排序整體均為浙江>上海>江蘇，插值精度與3個行政轄區(qū)的塊金值趨勢一致，說明行政轄區(qū)水田SOC插值誤差受微觀效應引起的空間變異顯著[16].

3 討論

明確SOC最優(yōu)插值方法的時空變化對動態(tài)制定我國農(nóng)田SOC“碳匯”管理政策具有重要意義.本研究結果表明，1982年和2000年整個太湖地區(qū)水田SOC最優(yōu)插值方法分別為OK和DK，這前人研究得到的地統(tǒng)計方法中的OK為農(nóng)業(yè)流域SOC預測的最佳插值方法[9]和高密度采樣時DK適用性更高[7]的結論一致.從不同亞類來看，1982年和2000年各個水稻土亞類LPI、GPI和IDW 3種確定性插值的SOC預測精度高于常用的OK、SK、UK和DK 4種地統(tǒng)計方法.較少學者討論過不同亞類土壤屬性插值方法的適用性，但在基于亞類的插值預測效果上，Biswas et al[33]提出按土壤類型布設樣點進行數(shù)字土壤制圖的合理性，張忠啟[13]也認為土壤亞類樣點可以更詳盡地反映SOC地形地貌地方性特點的空間變異，從側面證明了本研究基于土壤亞類的插值預測可以充分反應SOC含量的空間變異特征.另外，本研究的結果表明各土壤亞類的SOC最優(yōu)預測方法除了受不同土壤分區(qū)微地形因素影響外，長期不同農(nóng)業(yè)管理方式也會造成SOC含量變異程度的差異，從而改變其最優(yōu)插值方法.對于此，馬利芳等[12]的研究同樣表明人類干擾強度大的地帶，確定性插值優(yōu)于地統(tǒng)計插值精度，與本研究結果一致.從不同行政轄區(qū)來看，1982年和2000年大部分地區(qū)的最優(yōu)插值方法是地統(tǒng)計，其中OK適用性比SK、CRS等插值方法更廣.同樣，學者們在縣域行政區(qū)進行土壤屬性空間分布預測時，結果也表明OK在較為平坦的行政區(qū)內預測面平滑且插值效果較好[8,10].

總之，SOC最優(yōu)插值方法的選擇與結構性、隨機性不均衡影響下的時空變異有關，建議在進行SOC空間特征、定量預測和適宜樣點數(shù)研究時，應充分考慮樣點分布、人為活動和環(huán)境因素等造成時空變異的影響因子.隨著時間變化應關注最適宜插值方法的變化，且在土壤亞類插值研究時應考慮確定性方法的使用，這樣更有利于精準選擇適宜插值方法進行空間預測和定量分析.

盡管本研究基于1982年和2000年兩期大樣點實測數(shù)據(jù)，從行政轄區(qū)和土壤亞類的樣點分布、插值原理、土區(qū)、氣象及人為因素等多方面分析了太湖地區(qū)水稻土SOC最優(yōu)插值方法時空變化的機理，但SOC最優(yōu)插值方法也受到其他土壤屬性，如成土母質和堿解氮等的影響，有待于進一步研究.