基于GIS 和地理加權(quán)回歸的砂田土壤陽離子交換量空間預(yù)測①

王幼奇，張 興，趙云鵬，包維斌，白一茹*

(1 寧夏大學(xué)資源環(huán)境學(xué)院，銀川 750021；2 旱區(qū)特色資源與環(huán)境治理教育部國際合作聯(lián)合實驗室，銀川 750021)

土壤陽離子交換量 (CEC) 直接反映土壤保蓄、供應(yīng)和緩沖陽離子養(yǎng)分的能力，同時也影響著土壤有機質(zhì)、酸堿度和土壤結(jié)構(gòu)等性質(zhì)[1-2]。因此CEC 通常被作為土壤質(zhì)量評價指標之一，同時也是土壤施肥和改良等的重要依據(jù)。砂田是干旱半干旱區(qū)農(nóng)民在長期生產(chǎn)實踐中探索形成的一種獨有的保護性耕作方式，具有保溫增滲、抑制土壤蒸發(fā)和防止侵蝕等作用[3]。近年來隨著寧夏硒砂瓜產(chǎn)業(yè)的快速發(fā)展，砂田在發(fā)揮其巨大經(jīng)濟和生態(tài)效益的同時，也面臨著由于長期耕作導(dǎo)致土壤肥力下降、質(zhì)量退化等一系列嚴重問題[4]。因此精準預(yù)測砂田土壤CEC 空間分布規(guī)律對于防止砂田退化和提高土地生產(chǎn)力具有積極作用。

目前，基于地統(tǒng)計學(xué)的克里格插值方法在土壤要素空間預(yù)測中最常用[5]。普通克里格插值法(OK)對于樣點數(shù)量和樣點本身的數(shù)據(jù)質(zhì)量具有一定的依賴性，土壤自身具有的強變異性使得在復(fù)雜環(huán)境下使用OK法對土壤屬性進行空間插值已不能滿足當(dāng)前要素空間預(yù)測對于精度的要求[6]。因此利用輔助信息協(xié)助變量進行空間插值的方法如回歸克里格(RK)等得到了極大的應(yīng)用。RK 可以綜合多個變量對土壤屬性進行空間插值，有效提高了預(yù)測精度，但土壤屬性的強空間變異性使得最小二乘法(OLS)這一全局模型無法很好地體現(xiàn)土壤屬性的局部特征[7]。地理加權(quán)回歸模型(GWR)作為一種有效處理回歸分析中空間非平穩(wěn)性和空間依賴性的局部模型，近年來被廣泛應(yīng)用到科學(xué)研究中[8-9]。如江振藍等[10]探討了GWR 模型在土壤重金屬高光譜預(yù)測中的適用性及局限性；王合玲等[11]應(yīng)用GWR模型對艾比湖土壤有機質(zhì)和土壤因子響應(yīng)關(guān)系的空間非平穩(wěn)性進行了研究；袁玉蕓等[12]利用傳統(tǒng)回歸和GWR模型分析于田綠洲表層土壤鹽分及其影響因素的空間非平穩(wěn)性。以上研究均表明GWR模型對存在空間非平穩(wěn)性的數(shù)據(jù)具有更強的解釋能力和估計精度。地理加權(quán)回歸克里格法(GWRK)是將GWR 模型與RK 法相結(jié)合，綜合多個變量對土壤屬性進行局部空間插值，提高了插值精度，能更好地反映出土壤屬性的局部變異情況[13]。目前，關(guān)于寧夏砂田土壤屬性空間插值的研究多使用普通克里格插值法[14-15]，對于將GWR 模型應(yīng)用在土壤要素空間插值上的研究較少。本文以砂田土壤為研究對象，通過分析土壤CEC 及其影響因素之間的相關(guān)性，驗證和比較OK、RK 和GWRK 的空間制圖效果和插值精度，以期為土壤CEC 空間變異研究和土壤肥力管理提供科學(xué)依據(jù)。

1 材料與方法

1.1 研究區(qū)概況

研究區(qū)位于寧夏中衛(wèi)市香山鄉(xiāng)(圖1)，屬于寧夏中部干旱帶，地處 36°57′ ～ 37°07′N，104°56′ ～105°15′E，區(qū)內(nèi)氣候為溫帶大陸性氣候，多年平均氣溫6.8 ℃，年均降水量247.4 mm，降雨多集中于7—9月。全區(qū)地勢南北高中間低，平均海拔約1 740 m，土壤類型主要為淡灰鈣土。天然植被主要以旱生灌木、半灌木、耐寒的蒿屬和禾本科草類為主。土壤基本理化性質(zhì)見表1。

圖1 研究區(qū)樣點分布Fig. 1 Distribution of soil sampling sites in study area

表1 研究區(qū)土壤基本理化性質(zhì)Table 1 Soil physicochemical properties in study area

1.2 樣點布設(shè)與數(shù)據(jù)采集

以香山地區(qū)行政區(qū)劃圖為底圖，在研究區(qū)域采用1.5 km × 1.5 km 網(wǎng)格布點方式進行采樣。實際采樣過程中有些樣點落在村莊、道路等地，通常在附近200 m 內(nèi)進行調(diào)整，并用GPS 記錄調(diào)整后的坐標，共布設(shè)108 個樣點。使用土鉆采集根層土壤樣品(0 ～20 cm)。樣品帶回實驗室后經(jīng)自然風(fēng)干，剔除植物殘根及石礫等，碾磨分別過10 目和60 目篩待測。CEC的測定采用乙酸鈉-火焰光度法[1]。土壤有機質(zhì)采用重鉻酸鉀氧化-外加熱法測定[2]。土壤機械組成采用激光粒度儀(MS3000，Malvern instruments)測定，粒徑分級標準采用美國制[16]。

1.3 研究方法

1.3.1 回歸克里格法 回歸克里格法(RK)考慮到土壤屬性空間變異驅(qū)動因子的復(fù)雜性，將線性回歸與克里格插值相結(jié)合。當(dāng)目標變量與輔助變量存在相關(guān)關(guān)系時，先通過對目標變量和輔助變量的相關(guān)性分析和線性逐步回歸擬合，建立目標變量與輔助變量的多元(或一元)回歸關(guān)系，得到趨勢項。然后對所得殘差項進行半方差分析，并使用普通克里格法對殘差項進行空間插值，最后使用柵格計算器將趨勢項與殘差項兩項相加，即為回歸克里格的插值結(jié)果。公式表達為：

式中： z(s) 為目標變量在s處的預(yù)測值，m(s)為使用線性逐步回歸得到的趨勢項， ε(s)為使用普通克里格插值得到的殘差項。

1.3.2 地理加權(quán)回歸克里格法 地理加權(quán)回歸克里格法(GWRK)是通過對目標變量和輔助變量進行地理加權(quán)回歸擬合，得到局部回歸的殘差項，然后使用普通克里格插值法(OK)對所得殘差項進行插值。GWRK 將RK 中的全局回歸(式2)轉(zhuǎn)換成地理加權(quán)回歸模型(GWR)的局部回歸(式3)，能更好地體現(xiàn)土壤屬性空間變異的局部變化。

式中：yi為樣點i的因變量；xik為第i個樣點上第k個變量的觀測值；(ui,vi)為樣點i的地理空間坐標；β0為回歸的常數(shù)項；βk(ui,vi) 為第i個采樣點上的第k個回歸參數(shù)；εi為殘差項。如果βk(u,v)在空間保持不變，則模型(式3)就變?yōu)槿帜Ｐ?式2)。

1.3.3 模型精度檢驗指標 為了評價模型的預(yù)測精度，選取以下3 個指標對模型進行精度評價。平均誤差表示預(yù)測值與實測值偏差的算術(shù)平均值；均方根誤差表示預(yù)測值與實測值偏差的平方和觀測次數(shù)比值的平方根。均方根誤差對預(yù)測中特大或特小誤差反應(yīng)比較敏感；相對精度改進值(RI)是表示衡量模型模擬效果是否優(yōu)于僅對實測值取平均值的指標。

式中：ROK、RRK和RGWRK分別代表OK、RK 和GWRK擬合值與實測值的相關(guān)系數(shù)，RI 值為正值即表示RK、GWRK 較OK 的預(yù)測精度高，值越大說明預(yù)測精度提高地越多；RI 值為負值則表示RK、GWRK預(yù)測精度低于OK。

1.4 數(shù)據(jù)處理與分析 采用SPSS 18.0 和Excel 對數(shù)據(jù)進行描述統(tǒng)計分析、相關(guān)性分析和逐步回歸分析，使用GS+9.0 和ArcGIS 10.2 軟件對數(shù)據(jù)進行地統(tǒng)計分析和空間插值。為排除異常值對半方差函數(shù)穩(wěn)健性的影響，依據(jù)3σ 準則對原始取樣數(shù)據(jù)進行了異常值識別，異常值用正常的最大或最小值代替，后續(xù)分析采用處理過的原始數(shù)據(jù)進行計算。為驗證模型的預(yù)測精度，隨機選擇86 個樣點作為建模樣點進行空間插值，22 個樣點作為驗證樣點用于分析插值精度。

2 結(jié)果與分析

2.1 土壤CEC 描述統(tǒng)計分析

描述統(tǒng)計表明(表2)，建模點和驗證點CEC 均值分別為9.82、10.47 cmol/kg，含量變化范圍分別為5.75 ～ 14.87 cmol/kg 和6.15 ～ 15.83 cmol/kg。研究區(qū)土壤CEC 均值為10.145 cmol/kg。根據(jù)土壤保肥能力分級方法[17]，供肥保肥能力弱的樣本數(shù)占總樣本數(shù)的60.19%，供肥保肥能力中等的占39.81%，說明研究區(qū)土壤供肥保肥能力多處于低水平。偏度和峰度數(shù)值均較小(＜1)。從變異程度看，建模點和驗證點CEC變異系數(shù)分別為17.57% 和21.02%，屬中等變異。經(jīng)K-S 非參數(shù)檢驗，建模點和驗證點的P值均大于0.05，數(shù)據(jù)服從正態(tài)分布。

2.2 土壤CEC 與土壤理化性質(zhì)的相關(guān)性分析

影響土壤CEC 含量的因素有很多，有機質(zhì)及礦質(zhì)膠體的數(shù)量與性質(zhì)是決定土壤CEC 的重要因子[2]。相關(guān)學(xué)者亦研究發(fā)現(xiàn)CEC 與有機質(zhì)、土壤質(zhì)地存在顯著相關(guān)性，明顯影響CEC 在區(qū)域上的分布特征[18-19]。因此為提高CEC 預(yù)測精度、更好展現(xiàn)其空間分異特征，本文從CEC 與有機質(zhì)、礦質(zhì)膠體的耦合關(guān)系為出發(fā)點，選取了研究區(qū)土壤有機質(zhì)、質(zhì)地與CEC 進行相關(guān)性分析[20]。表3 可以看出，CEC 與有機質(zhì)、砂粒、粉粒和黏粒呈顯著相關(guān)關(guān)系。其中CEC 與有機質(zhì)、黏粒和粉粒呈顯著正相關(guān)，與砂粒呈顯著負相關(guān)。

表2 土壤CEC 的描述統(tǒng)計特征Table 2 Descriptive statistics of soil CEC

2.3 土壤CEC 的回歸模型與半方差模型

2.3.1 逐步回歸分析 在相關(guān)分析基礎(chǔ)上，選擇有機質(zhì)、黏粒、粉粒和砂粒作為輔助變量，采用逐步回歸篩選出對各個回歸模型具有顯著貢獻的回歸參數(shù)。逐步回歸結(jié)果表明(表4)，有機質(zhì)和砂粒是對CEC 插值的最佳自變量。方差膨脹因子(VIF)為1.12，小于7.5，說明自變量間不存在共線性。調(diào)整后的決定系數(shù)為0.567，說明OLS 模型對CEC 的方差解釋度為56.7%。擬合模型如下：y= 7.478+0.964x1-0.052x2，式中，y為CEC；x1為有機質(zhì)含量；x2為砂粒含量。為了便于對比，應(yīng)用GWR 模型時，同樣選用有機質(zhì)和砂粒進行建模。

表3 土壤CEC 與土壤理化性質(zhì)的相關(guān)性Table 3 Correlation between soil CEC and other soil physicochemical properties

表4 土壤CEC 與土壤理化性質(zhì)的逐步回歸過程Table 4 Stepwise linear regression process of soil CEC with other physicochemical properties

2.3.2 地統(tǒng)計分析 對CECOLS 殘差和GWR 殘差進行描述性統(tǒng)計分析。OLS 殘差變化范圍為 -1.49 ～4.02，均值為0.004 7；GWR 殘差變化范圍為 -1.57 ～4.44，均值為 -0.0345。經(jīng)K-S 非參數(shù)檢驗，兩種回歸方法所得殘差均服從正態(tài)分布。在GS+9.0 中對實測值、OLS 殘差和GWR 殘差進行半方差函數(shù)模型擬合(表5)。可以看出，三者均可用球狀模型擬合。塊金值(C0)代表隨機因素帶來的變異，通常由區(qū)域變量的變異性和測定誤差所造成[21]。三者的塊金值(C0)均小于0.3，表明CEC 由測定誤差或土壤性質(zhì)帶來的隨機變異較小，空間變異主要受結(jié)構(gòu)性因子影響。塊金系數(shù)(C0/(C0+C))表示由隨機因素引起的空間變異占系統(tǒng)總變異的比例。一般來說，比值＜25% 時空間變異性表現(xiàn)為強空間相關(guān)性，介于25% ～ 75% 時表現(xiàn)為中等強度空間相關(guān)性，比值＞75% 時表現(xiàn)為弱空間相關(guān)性[17]。實測值、OLS 殘差和GWR 殘差塊金系數(shù)分別為8.50%、6.36% 和7.02 %，比值均在25% 以下，說明CEC 及其殘差具有強烈空間自相關(guān)。

表5 土壤CEC、OLS 殘差和GWR 殘差的半方差模型參數(shù)Table 5 Semivariance parameters of soil CEC and regression residuals of OLS and GWR

2.4 土壤CEC 的空間分布特征

利用OK、RK 和GWRK 3 種方法進行空間插值得到的CEC 空間分布圖。根據(jù)圖2 可以看出，3 種方法得到的結(jié)果中CEC 的分布格局大致相似，即研究區(qū)西部的CEC 要大于東部，高值區(qū)分布在西北、東北和南部，CEC 存在從這些地區(qū)向中東部逐漸減小的趨勢。CEC 預(yù)測值范圍分別為5.75 ～ 14.88、5.77 ～ 13.00、5.52 ～ 12.90 cmol/kg，各方法基本一致。從制圖效果看，OK 法插值圖較為粗略，RK 法的插值結(jié)果體現(xiàn)了CEC 隨著有機質(zhì)和砂粒的變化情況，但是斑塊的邊界較為零碎，這與現(xiàn)實中CEC 的空間漸變特征不相符。GWRK 法不僅保留了回歸分析中CEC 隨土壤理化性質(zhì)變化的豐富細節(jié)，還反映出CEC 的空間漸變特征，斑塊邊界更為光滑，感官上更符合CEC 的實際分布狀況。為了更好地反映3 種插值模型對區(qū)域土壤CEC 的預(yù)測精度，利用平均誤差、均方根誤差和RI 值對插值結(jié)果進行精度評價。

2.5 不同插值方法的空間預(yù)測精度評價

表6 對OK、RK 和GWRK 方法進行了精度評價。結(jié)果表明3 種方法的平均誤差基本接近于0，表明了所構(gòu)建模型的預(yù)測具有無偏性；均方根誤差OK＞RK＞GWRK，范圍為1.084 ～ 1.853，趨近于1，說明所構(gòu)建模型的誤差有效。從相對精度改進值來看，RK 和GWRK 法相較于OK 法擬合精度分別提高了40.49%和41.50%，表明3 種插值方法中GWRK法插值精度最高，為最優(yōu)插值模型。

圖2 三種插值方法下土壤CEC 的空間分布Fig. 2 Spatial distribution of soil CEC by OK, RK and GWRK interpolation

表6 OK、RK 和GWRK 方法的精度評價指標對比Table 6 Precision evaluation of three interpolation methods

3 討論

研究區(qū)土壤CEC 均值為10.145 cmol/kg，供肥保肥弱的樣本數(shù)占總樣本數(shù)的60.19%，研究區(qū)土壤保肥供肥能力較弱。前人的相關(guān)研究也發(fā)現(xiàn)了類似情況[4,15]。其產(chǎn)生原因可能與當(dāng)?shù)赝寥滥纲|(zhì)和耕作制度有關(guān)。研究區(qū)處于我國西北內(nèi)陸腹地，土壤類型為低肥力的淡灰鈣土，自然環(huán)境條件限制了當(dāng)?shù)赝寥赖姆柿λ?。再加上研究區(qū)有覆砂的傳統(tǒng)耕作制度，農(nóng)民施肥難度大，外源性的有機肥施入量不夠，導(dǎo)致當(dāng)?shù)赝寥婪柿ζ?。砂田土壤CEC 在空間上展現(xiàn)出高度變異性，在研究區(qū)地形、氣候和母質(zhì)相對均一的情況下，土壤理化性質(zhì)可能是導(dǎo)致土壤CEC 在空間上呈現(xiàn)出較強變異性的重要原因。其中有機質(zhì)中膠體成分在土壤固相中所占比例較大，土壤膠體是CEC 的載體，隨有機質(zhì)含量增加，CEC 也隨之增加。同時在土壤顆粒組成中，黏粒膠體集中了80% 以上的負電荷，土壤黏粒含量愈高，CEC 就越高[18-20]。相關(guān)學(xué)者亦研究發(fā)現(xiàn)CEC 與有機質(zhì)、土壤質(zhì)地存在顯著相關(guān)性，明顯影響CEC 在區(qū)域上的分布特征。砂田CEC 與土壤有機質(zhì)、黏粒和粉粒呈正相關(guān)，與砂粒呈負相關(guān)，即在有機質(zhì)和黏粒含量高的地方，CEC相對充足。針對這一情況，建議在后續(xù)農(nóng)業(yè)生產(chǎn)過程中增施有機肥，推廣輪作制度以培肥土壤，從而提高作物產(chǎn)量。

空間插值技術(shù)是利用采樣點數(shù)據(jù)對研究區(qū)內(nèi)未知點進行預(yù)測的方法。常用的空間插值方法有普通克里格法、回歸克里格法等[22-23]。張慧智等[24]對比了普通克里格、泛克里格和回歸克里格在土壤溫度空間預(yù)測中的效果，發(fā)現(xiàn)回歸克里格法預(yù)測精度最高，能夠更好地表達復(fù)雜地形地區(qū)的局部變異。Hengl 等[25]采用回歸克里格對土壤有機質(zhì)含量進行空間預(yù)測，發(fā)現(xiàn)其預(yù)測精度較普通克里格表現(xiàn)更好，均方根誤差數(shù)值更低。但是，土壤理化性質(zhì)與土壤CEC 之間的相關(guān)關(guān)系因地而異，局部區(qū)域的CEC 與土壤理化性質(zhì)關(guān)系可能與做全局分析時得出的結(jié)論相反，這在回歸克里格中沒有考慮。GWR 模型能將數(shù)據(jù)的空間信息納入分析過程，通過計算回歸模型的局部參數(shù)來解決地理數(shù)據(jù)中存在的空間自相關(guān)性及空間非平穩(wěn)性問題，從而提高模型的擬合優(yōu)度及模擬效果[26]。本文對3 種插值模型的精度分析結(jié)果也證明了GWRK 法在空間制圖方面表現(xiàn)出更豐富的局部細節(jié)和更光滑的空間漸變特征。因此，在后續(xù)土壤屬性的空間制圖研究中，應(yīng)注重GWRK 法在土壤屬性空間插值中的應(yīng)用。

4 結(jié)論

1) 相關(guān)分析表明，土壤CEC 與有機質(zhì)、黏粒和粉粒含量呈顯著正相關(guān)，與砂粒含量呈顯著負相關(guān)。

2) 地統(tǒng)計分析表明，土壤CEC 及其插值殘差均可用球狀模型擬合；實測值、OLS 殘差和GWR 殘差塊金值均小于0.3，說明其空間變異主要受結(jié)構(gòu)性因子影響；塊金系數(shù)分別為8.50%、6.36% 和7.02%，表現(xiàn)出強烈空間自相關(guān)。

3)從插值效果和精度分析來看，GWRK 是3 種模型中插值效果最為理想的；RK 和GWRK 法相較于OK 法擬合精度分別提高了40.49% 和41.50%，3種插值方法中GWRK 法插值精度最高，為最優(yōu)插值模型。