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

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

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

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

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

1 材料與方法

1.1 研究區(qū)概況

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

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

1.2 數(shù)據(jù)來(lái)源

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

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

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

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

1.4 預(yù)測(cè)結(jié)果精度檢驗(yàn)

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

(1)

(2)

(3)

式中，N為采樣點(diǎn)數(shù)量，Xoi為采樣點(diǎn)SOC含量實(shí)測(cè)值，Xpi為采樣點(diǎn)SOC含量預(yù)測(cè)值，R為實(shí)測(cè)值和預(yù)測(cè)值之間的Pearson相關(guān)系數(shù)，其越大表明模型對(duì)各樣點(diǎn)集含量的解釋能力越強(qiáng).

2 結(jié)果與分析

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

由表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ù)來(lái)看，兩個(gè)時(shí)期太湖水田SOC含量的離散程度相差不大，均在10%～100%之間，為中等變異[23].從整個(gè)太湖地區(qū)和不同水稻土亞類、行政轄區(qū)各預(yù)測(cè)集樣點(diǎn)SOC的偏度、峰度和K-S檢驗(yàn)值來(lái)看，1982年和2000年研究區(qū)各預(yù)測(cè)集均符合正態(tài)分布，全局莫蘭I指數(shù)的范圍為0.146～0.448，對(duì)應(yīng)Z得分范圍是6.45～35.43，Z高于2.58且顯著性均小于0.01，表明兩期預(yù)測(cè)集的SOC含量有明顯的空間聚集現(xiàn)象，具有空間正相關(guān)模式，符合地統(tǒng)計(jì)學(xué)插值前提[21].

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

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

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

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

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

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

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

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

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

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

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

表3 太湖地區(qū)水田土壤整個(gè)區(qū)域插值模型預(yù)測(cè)精度評(píng)價(jià)1)

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

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

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

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

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

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

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

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

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

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

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

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

上海轄區(qū)水田1982年和2000年SOC最優(yōu)插值方法分別為OK和SK(表6和表7).這主要是因?yàn)?982年該地區(qū)59.05%樣點(diǎn)分布在東北部平原土區(qū)，其余分布于西南低洼圩區(qū)，區(qū)域內(nèi)極值適中，OK對(duì)于地勢(shì)平坦的縣域土壤屬性插值效果更好.龍軍等[18]也認(rèn)為OK方法可以較好預(yù)測(cè)平原縣域耕地SOC.但在2000年，屬于上海轄區(qū)的沖積平原土區(qū)采樣點(diǎn)比例比1982年提高了5.31%，尤其城市擴(kuò)張?jiān)斐沙缑鲘u和上海市區(qū)耕地大量減少形成更大范圍無(wú)樣點(diǎn)分布的空白區(qū)，SK可利用經(jīng)驗(yàn)?zāi)Ｐ推交A(yù)測(cè)面進(jìn)行非平穩(wěn)結(jié)構(gòu)的線性預(yù)測(cè)，尤其適用于這類微地形等結(jié)構(gòu)性因素影響更多的情況.

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

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

3 討論

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

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

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