劉永坤,劉 偉,王昌陽,陳嘉明,馬 進(jìn)

(1. 江蘇師范大學(xué)地理測繪與城鄉(xiāng)規(guī)劃學(xué)院，江蘇 徐州 221116； 2. 中國礦業(yè)大學(xué)環(huán)境與測繪學(xué)院，江蘇 徐州 221116； 3. 武漢大學(xué)遙感信息工程學(xué)院，湖北 武漢 430079)

土壤全鉀含量是植物生長的主要營養(yǎng)元素，同時也是土壤中影響農(nóng)作物產(chǎn)量的三要素之一[1]。研究土壤鉀含量的空間分布特征，不僅對均衡和擴(kuò)展土壤有效養(yǎng)分庫具有重要意義，還能為土地資源可持續(xù)利用和地區(qū)農(nóng)業(yè)的健康發(fā)展提供理論依據(jù)[2]。近幾年來，國外文獻(xiàn)[3—4]對土壤空間變異規(guī)律進(jìn)行了深層次的研究；國內(nèi)文獻(xiàn)[5]對土壤重金屬、土壤有機(jī)質(zhì)的空間變異進(jìn)行研究，取得了大量的研究成果。青海湖流域作為典型的復(fù)雜地貌類型區(qū)，研究土壤全鉀含量對改善地區(qū)施肥效果、提高土壤利用效率具有一定的現(xiàn)實(shí)意義[6]。

土壤屬性的連續(xù)變化是土壤資源科學(xué)管理和使用的前提，土壤各類屬性的空間插值是用來評估土壤性質(zhì)連續(xù)變化的主要方法[7]，同時也是“數(shù)字土壤”和“土壤計(jì)量學(xué)”的一個重要研究工具。空間插值精度取決于很多因素，一部分研究成果表明Kriging的模擬效果優(yōu)于反距離加權(quán)法(inverse distance weighting method,IDW)和Spline,而另一部分研究成果持相反觀點(diǎn)[8]。如文獻(xiàn)[9]通過比較3種插值方法評價土壤全氮含量的插值誤差，得出徑向基函數(shù)法插值模型優(yōu)于普通克里金插值法(ordinary kriging,OK)和IDW。修思玉等通過3種插值方法評價黑龍江省2010年32格氣象站點(diǎn)的年降雨量，采用平均誤差、平均絕對誤差、平均相對誤差等指標(biāo)進(jìn)行優(yōu)度檢驗(yàn)，得出Spline法最優(yōu)、IDW得到的降雨量的最大值與最小值與原始降雨量數(shù)據(jù)相差最大。文獻(xiàn)[10]對比了土壤Cd污染插值模型，得出不同的插值方法對極大值都存在不同程度的平滑。綜合以上研究發(fā)現(xiàn)，目前對空間屬性曲面建模的研究主要集中在單一全局插值模型。因此研究基于空間屬性的自適應(yīng)曲面建模具有重要的意義。

此外，一系列研究論證了利用地學(xué)環(huán)境變量為補(bǔ)充信息能夠有效提高插值精度。文獻(xiàn)[11]以高程、坡度、土地利用等輔助變量預(yù)測土壤精度。文獻(xiàn)[12]采用土壤電導(dǎo)率變量對不同尺度的土壤鹽分開展了空間變異特征的隨機(jī)模擬。文獻(xiàn)[13]采用NDVI、高程和距河流距離為輔助變量研究黃河三角洲地區(qū)的鹽漬化，指出結(jié)合地學(xué)環(huán)境變量的插值模擬效果更好。但是，結(jié)合地學(xué)環(huán)境變量的插值模型也存在一些不足，如很多研究使用單一、同類的地學(xué)環(huán)境變量解決不同地區(qū)的預(yù)測問題，然而不同地區(qū)受到不同因素影響，其考慮的環(huán)境變量是不同的。各個插值方法的適應(yīng)性不同，不存在絕對最優(yōu)的空間插值模型，若能根據(jù)不同的插值模型自身的適應(yīng)性，分區(qū)進(jìn)行多模型集成，理論上可以提高模擬精度。

1 研究區(qū)概況與采樣方法

1.1 青海湖自然概況

研究區(qū)(36°38′N—37°29′N,99°52′E—100°50′E)地處青海湖流域東南部。由于地殼運(yùn)動形成了復(fù)雜多樣的地貌類型，地貌以沖擊洪積扇平原為主，流域周圍大小河流40余條，屬于內(nèi)陸封閉水系，是許多野生動物的棲息繁殖區(qū)域。青海湖流域?qū)俑咴敫珊档臏貛Т箨懶詺夂?，位于凍土地帶，干旱少雨，多風(fēng)，太陽輻射強(qiáng)，溫差變化較大。研究區(qū)主要包括共和、剛察和海晏3個地區(qū)，除去青海湖，總面積約為2100 km2,海拔高度為3043～4516 m，有大量不同的土壤、植被、地貌類型，屬于典型的復(fù)雜地貌類型區(qū)。

1.2 樣品采集

研究通過空間分層組合的方式[14]進(jìn)行土壤采樣。在土壤采樣的同時記錄了采樣點(diǎn)的經(jīng)緯度、海拔、土地利用類型、土壤類型、植被類型及地質(zhì)類型等與土壤屬性密切相關(guān)的地學(xué)環(huán)境信息，每個樣點(diǎn)在0～15 cm、15～30 cm依次采樣兩次。筆者于2013年9月，在青海省環(huán)境檢測中心人員提供協(xié)助下，收集了青海湖典型地區(qū)土壤表層(0～30 cm)樣點(diǎn)數(shù)據(jù)110個。采樣結(jié)束后在實(shí)驗(yàn)室對樣品進(jìn)行風(fēng)干、研磨、過2 mm篩等流程后，取兩次土壤中全鉀元素的含量平均值作為記錄樣本值。

2 研究方法

2.1 地學(xué)環(huán)境要素篩選

已有研究表明，利用地學(xué)環(huán)境要素作為輔助變量能夠有效提高土壤全鉀含量插值精度和制圖效果。土壤全鉀空間變異性驅(qū)動因子主要包括：土地利用類型、土壤類型、植被類型、地質(zhì)類型、坡度、施肥等[15-16]。根據(jù)已有的研究結(jié)論，結(jié)合本次主要針對自然景觀類型區(qū)研究，排除施肥和土地管理措施2個驅(qū)動因子，選取土地利用類型、土壤類型、植被類型、地質(zhì)類型4個地學(xué)環(huán)境要素作為輔助變量。其中土地利用類型有灌叢、耕地和草原等6種；土壤類型有栗鈣土、流動和風(fēng)沙土等5種；植被類型有水柏枝草場型、紫花針茅雜草草場型等30種；地質(zhì)類型有沙丘、沙山和河谷平原等13種。使用SPSS計(jì)算輔助變量的土壤全鉀含量統(tǒng)計(jì)特征。為了進(jìn)一步篩選對土壤全鉀空間分布具有顯著影響的地學(xué)要素，利用SPSS軟件對土壤全鉀含量與上述4個地學(xué)要素進(jìn)行方差分析，選取具有顯著特征的地學(xué)要素作為輔助環(huán)境變量。表1的方差結(jié)果分析表明：在本文研究區(qū)域內(nèi)，植被類型的景觀破碎度太高，很多區(qū)域內(nèi)只有1個采樣點(diǎn)，一些亞類型的草地只有1～2個采樣點(diǎn)，模擬效果不理想，因此將植被類型排除在外。

2.2 土壤屬性建模方法

2.2.1 插值方法的選取

本文中對土壤鉀含量進(jìn)行空間插值的4種方法均基于ArcGIS 10.2平臺實(shí)現(xiàn)。傳統(tǒng)的插值方法有OK、IDW、回歸克里金(regression Kriging,RK)、貝葉斯克里金(Bayesian Kriging,BK)。OK在滿足二階平穩(wěn)假設(shè)和內(nèi)蘊(yùn)條件的基礎(chǔ)上，根據(jù)區(qū)域化變量的原始數(shù)據(jù)和變異函數(shù)(或協(xié)方差)的特點(diǎn)，確定待估點(diǎn)周圍已知點(diǎn)參數(shù)對待估點(diǎn)的加權(quán)值大小，再對待估點(diǎn)做出最優(yōu)估計(jì)的方法，其優(yōu)點(diǎn)是顧及各樣本點(diǎn)的空間相關(guān)性，在得到模擬結(jié)果的同時，得出估算精度的方差；IDW使用一組采樣點(diǎn)的線性權(quán)重組合來確定輸出的柵格值，權(quán)重與插值點(diǎn)和樣本點(diǎn)之間的距離成反比，是一種加權(quán)移動平均的方法；RK綜合考慮了影響鉀空間變異的多種環(huán)境因子，但由于研究區(qū)環(huán)境要素分布比較破碎，導(dǎo)致沒有足夠多的采樣點(diǎn)來估計(jì)相對準(zhǔn)確的半方差函數(shù)。對全鉀含量進(jìn)行半變異分析，得到插值時的各參數(shù)值和擬合模型，其中塊金值較小，表明自身隨機(jī)性誤差引起的變異性不大；S/N+S比值接近于1，N/S比值小于30%，表明全鉀含量具有較強(qiáng)的空間相關(guān)性，擬合的最優(yōu)模型為指數(shù)模型[15]。

表1 不同地學(xué)要素類型之間土壤全鉀的方差分析

注：*0.05顯著水平；**0.01顯著水平。

2.2.2 融合地學(xué)環(huán)境要素的自適應(yīng)曲面建模方法

2.2.2.1 融合地學(xué)環(huán)境要素的基插值模型

根據(jù)變異和空間相關(guān)理論，任何變量的空間變化都可以用下述兩個主要成分之和來表示：與局部變化有關(guān)的殘差成分、與趨勢有關(guān)的結(jié)構(gòu)成分。相應(yīng)的土壤全鉀含量的空間分布模式可以表示為

S(xi,l,k,yj,l,k)=direction(Geox,y)+residual(xi,l,k,yj,l,k)

(1)

式中,S(xi,l,k,yj,l,k)為第l種地學(xué)要素的第k種類型的土壤全鉀采樣點(diǎn)模擬值，其中(xi,yj)為采樣點(diǎn)坐標(biāo)，i和j分別表示坐標(biāo)的行和列；direction(Geox,y)是描述(xi,yj)處第l種地學(xué)要素的第k種類型S的趨勢值，其中Geo(x,y)是描述(xi,yj)處與土壤全鉀密切相關(guān)的地學(xué)環(huán)境信息；residual(xi,l,k,yj,l,k)是描述(xi,yj)處S的第l種地學(xué)要素的第k種類型中第i行、第j列柵格點(diǎn)土壤全鉀的殘差值?；诰的Ｐ偷玫节厔莺瘮?shù)，擬合描述結(jié)構(gòu)成分的趨勢面，實(shí)現(xiàn)趨勢分離；并根據(jù)土壤厚度采樣點(diǎn)的空間特征，比較選擇最佳的空間相關(guān)插值算法對殘差進(jìn)一步處理。基于以上理論知識，具體試驗(yàn)方案如圖1所示。

2.2.2.2 利用集成學(xué)習(xí)實(shí)現(xiàn)自適應(yīng)分區(qū)

根據(jù)上述建立的基插值模型Modeli對整個研究區(qū)進(jìn)行模擬生成一系列土壤全鉀插值曲面，使用采樣點(diǎn)土壤全鉀實(shí)測值減去預(yù)測值得到模擬誤差，以掃面線算法(算法步驟如圖2所示)為集成學(xué)習(xí)的分類學(xué)習(xí)器，對插值曲面進(jìn)行掃描，選擇分類正確率較高且差異較大的掃描線進(jìn)行集成，對各個插值曲面進(jìn)行多種類的自適應(yīng)分區(qū)，得到各個插值模型適用的空間范圍。

“+”和“-”分別表示滿足插值精度要求的樣點(diǎn)和不滿足插值精度要求的樣點(diǎn)。在這個過程中，使用水平和垂直的掃描線作為分類器，對土壤全鉀插值曲面進(jìn)行分區(qū)。

(1) 根據(jù)對插值曲面分區(qū)的正確率，得到一個新的樣本分布(即樣本中每個樣本點(diǎn)的權(quán)重分布)D2、一個子分類器h1。其中畫圈的樣本表示被分錯的，比較大的“+”表示對該樣本作了加權(quán)。

(2) 根據(jù)分區(qū)的正確率，得到一個新的樣本分布D3、一個子分類器h2。弱分類器h2中有3個“-”符號分類錯誤，得到分區(qū)錯誤率ε2=0.21，h2應(yīng)分配的權(quán)重α2=0.66。

(3) 得到一個子分類器h3。弱分類器h3中有2個“+”符號和一個“-”符號分類錯誤，得到分區(qū)錯誤率ε3=0.14，h3應(yīng)分配的權(quán)重α3=0.91。

(4) 整合所有子分類器，得到最終的Hfinal。經(jīng)過以上三步即可抽取所有滿足精度閾值的分區(qū)，每個區(qū)域?qū)儆谀膫€類別，由該區(qū)域所在分類器的權(quán)值綜合決定。經(jīng)過集成和NO Data區(qū)域的填充，得到土壤屬性自適應(yīng)曲面建模方法(adaptive surface modeling for soil properties,ASM-SP)模擬的空間分布圖。

2.3 插值結(jié)果精度檢驗(yàn)

獨(dú)立驗(yàn)證模型存在檢驗(yàn)精度不夠的問題，使用交叉驗(yàn)證能夠更優(yōu)評價插值方法的質(zhì)量[17-19]。在交叉驗(yàn)證中，針對110個監(jiān)測點(diǎn)中的每一個點(diǎn)位數(shù)據(jù)，假設(shè)其數(shù)據(jù)未知，基于其余109個點(diǎn)位的數(shù)據(jù)，使用不同插值方法計(jì)算出模擬值，再結(jié)合觀測值求算誤差。

交叉驗(yàn)證的評定指標(biāo)有平均估計(jì)誤差百分比(PAEE)、相對均方差(RMSE)、均方根預(yù)測誤差(RMSPE)及殘差分析等。本文采用最大誤差(max error)、最小誤差(min error)、平均誤差(ME)、均方根誤差(RMSE)。

3 模擬結(jié)果與分析

3.1 土壤全鉀含量的統(tǒng)計(jì)特征

從研究區(qū)110個表層土壤鉀含量的統(tǒng)計(jì)結(jié)果分析(如圖3所示)可以看出，土壤鉀含量的值介于1.406 5%～2.346 4%之間，平均值為1.958 8%，K-S檢驗(yàn)結(jié)果表明鉀含量總體上服從正態(tài)分布。土壤鉀含量的變異系數(shù)為10.99%，表明研究區(qū)土壤鉀含量的空間變異性處于中等水平[15]。

3.2 ASM-SP

ASM-SP由3個步驟組成：①建立一系列基插值模型，即OK-Landuse，OK-Soil，OK-Geology；②基插值曲面的自適應(yīng)分割；③基插值優(yōu)化組合。自適應(yīng)分區(qū)方法詳見2.2節(jié)。ASM-SP構(gòu)建過程如下：

3.2.1 基插值模型建立

(1) 計(jì)算每種地學(xué)要素的土壤全鉀均值并得到趨勢曲面direction(Geox,y)，根據(jù)土壤全鉀實(shí)測值，結(jié)合輔助變量，得到土壤全鉀平均值(如圖3所示)。

(2) 用實(shí)測值減去土壤全鉀含量均值，計(jì)算得到土壤全鉀的模擬殘差。對模擬殘差進(jìn)行OK插值得到殘差曲面residual(xi,l,k,yj,l,k)(如圖4所示)。

(3) 利用式(1)計(jì)算第l種地學(xué)要素的第k種類型的土壤全鉀采樣點(diǎn)模擬值S(xi,l,k,yj,l,k)，即將上述的趨勢曲面和殘差曲面相加。這是基插值曲面的構(gòu)造過程。

3.2.2 基于線性掃描算法的自適應(yīng)曲面分割

基于第2.2.2節(jié)所述的構(gòu)造誤差曲面的方法，獲得不同插值模型的誤差面(如圖5所示)，并確定每個插值模型的適用范圍。

3.2.3 插值曲面的整合

在柵格單元優(yōu)化的基礎(chǔ)上，選取具有最小誤差的柵格單元的插值結(jié)果作為最佳柵格單元進(jìn)行集成。圖6顯示了與不同插值模型相對應(yīng)的最優(yōu)分區(qū)。

3.3 不同方法預(yù)測精度分析

為了評價ASM-SP模擬土壤全鉀含量空間分布的精度，本文比較了5種插值方法的精度，即RK、BK、IDW、OK+Landuse和ASM-SP。OK、BK結(jié)合輔助變量Landuse，與ASM-SP進(jìn)行對比；RK方法結(jié)合Landuse擬合趨勢項(xiàng)。逐點(diǎn)交叉驗(yàn)證的誤差指標(biāo)(見表2),對于最大誤差值來說數(shù)值小精度高，對于最小誤差值來說數(shù)值小精度高。表明ASM-SP預(yù)測誤差范圍最小，預(yù)測值的靈敏度和反映極值的情況優(yōu)于其他4種方法。針對ME和RMSE這幾個值而言，ASM-SP分別達(dá)到-0.002 1、0.281 1，相對其他結(jié)果較優(yōu)，且精度AC值最優(yōu)，因此ASM-SP顯示出更好的性能，其ME值比傳統(tǒng)插值方法(即OK和IDW)的ME值更接近于0。這意味著結(jié)合輔助變量的插值更加具有無偏性。并且ASM-SP的精度也充分表明其回歸曲線可以更好模擬預(yù)測值與真實(shí)值的關(guān)系。綜上所述ASM-SP的插值方法最優(yōu)。

總的來說，ASM-SP插值方法優(yōu)于其他方法的主要原因有：①該方法結(jié)合了輔助變量，因此更準(zhǔn)確地刻畫土壤全鉀隨地學(xué)環(huán)境要素變化而突出的邊界；②使用線性掃描算法將最優(yōu)區(qū)分割出來，并按此方法將各個地學(xué)要素得到最優(yōu)區(qū)域整合，在很大程度上提高了插值的精度。

表2 5種方法插值有效性交叉驗(yàn)證結(jié)果

3.4 不同預(yù)測方法效果對比

本文對比了5種土壤全鉀含量的插值效果(如圖7所示)。從圖7中可以看出，RK、BK和OK刻畫的變異程度較真實(shí)值區(qū)間小，插值顯示出不同程度的弱“牛眼”效應(yīng)。IDW插值可以體現(xiàn)出土壤全鉀分布的總體格局，但具有較強(qiáng)的“牛眼”效應(yīng)，且插值精度較低。

ASM-SP能夠有效地刻畫土壤全鉀空間變異格局，生成一個適中的插值范圍(1.30～2.32)，且能體現(xiàn)土壤全鉀空間分異較多的局部變化細(xì)節(jié)。BK、RK和OK插值得到的土壤全鉀值和IDW插值無法描述土壤屬性隨地學(xué)環(huán)境變量變化導(dǎo)致的突變邊界信息。ASM-SP方法對地形復(fù)雜地貌類型區(qū)土壤屬性的空間插值具有較強(qiáng)的適應(yīng)性，能更準(zhǔn)確地描述研究區(qū)土壤屬性的空間變異格局。

4 結(jié) 語

利用地學(xué)環(huán)境變量對插值結(jié)果進(jìn)行修正，能夠解決土壤屬性空間分異突變的問題，尤其能夠解決在丘陵溝壑等復(fù)雜地貌區(qū)和不同地學(xué)要素銜接區(qū)產(chǎn)生較大突變的問題。土壤屬性的空間分異在較短水平距離內(nèi)產(chǎn)生較大的變化，導(dǎo)致單一全局插值模型精度有限。因此本文提出了一種融合地學(xué)環(huán)境要素，將不同插值模型最優(yōu)區(qū)進(jìn)行自適應(yīng)建模(ASM-SP)方法。在復(fù)雜地貌類型區(qū)，ASM-SP能夠更準(zhǔn)確地刻畫土壤屬性的空間變異，有效降低模擬誤差。此外ASM-SP結(jié)合輔助變量，使其模擬結(jié)果更符合地學(xué)規(guī)律，便于物理解釋土壤屬性的空間變異特征。

對比未使用輔助變量的空間插值模型(IDW、OK)，以及利用地形因子輔助插值的RK法的結(jié)果顯示，ASM-SP刻畫土壤全鉀含量更加符合研究區(qū)土壤屬性空間變異規(guī)律，地學(xué)要素邊緣的細(xì)節(jié)信息表現(xiàn)得更加明顯。ME和RMSE等精度評價指標(biāo)也更小，分別達(dá)到-0.002 1和0.281 1，顯示出比其他插值模型更高的模擬精度，尤其能準(zhǔn)確刻畫復(fù)雜地貌類型區(qū)土壤全鉀空間變異隨周圍地學(xué)環(huán)境要素變化而突出的邊界。ASM-SP是實(shí)現(xiàn)復(fù)雜地貌類型區(qū)土壤全鉀含量的高精度模擬的一種新方法，為以后的土壤屬性制圖研究提供了新思路和有益借鑒。