趙明松，劉斌寅，盧宏亮，李德成，張甘霖

基于地理加權(quán)回歸的地形平緩區(qū)土壤有機(jī)質(zhì)空間建模

趙明松1,2，劉斌寅1，盧宏亮1，李德成2，張甘霖2※

（1. 安徽理工大學(xué)測(cè)繪學(xué)院，淮南 232001；2. 土壤與農(nóng)業(yè)可持續(xù)發(fā)展國(guó)家重點(diǎn)實(shí)驗(yàn)室中國(guó)科學(xué)院南京土壤研究所，南京 210008）

氣候變化效應(yīng)評(píng)估、土壤固碳潛力和肥力管理等，迫切需要詳盡的土壤有機(jī)質(zhì)（soil organic matter, SOM）空間分布信息。該文以江蘇省第二次土壤普查的1 519個(gè)典型土壤剖面的表層（0～20 cm）SOM含量為例，選擇1 217個(gè)樣本為建模集，302個(gè)為驗(yàn)證集，選取年均溫度、年均降雨、物理性黏粒和土壤pH值等因子進(jìn)行SOM的地理加權(quán)回歸（geographically weighted regression, GWR）建模。從建模集中分別隨機(jī)抽取100%（1 217個(gè)）、80%（973個(gè)）、60%（730個(gè)）、40%（486個(gè)），20%（243個(gè)）的樣點(diǎn)，對(duì)比不同樣點(diǎn)數(shù)量下GWR和傳統(tǒng)全局回歸模型的精度差異，并選擇最優(yōu)模型進(jìn)行SOM空間預(yù)測(cè)制圖。結(jié)果表明：1）江蘇省SOM含量在不同空間尺度上存在極顯著的空間自相關(guān)性。不同樣點(diǎn)數(shù)量的建模集的全局自相關(guān)性和局部空間自相關(guān)聚類圖結(jié)果相似。全局Moran’s I值介于0.25～0.61（<0.001）。SOM含量空間分布以空間聚集特征為主，“高-高”聚集區(qū)主要分布在蘇中和蘇南地區(qū)，“低-低”聚集區(qū)主要分布在蘇北地區(qū)。2）GWR建模結(jié)果均優(yōu)于傳統(tǒng)的全局回歸建模，其殘差在不同的空間尺度上均不存在空間自相關(guān)性。不同建模集的GWR的2adj較全局建模均提高0.15～0.20，其AIC和RSS均比全局模型有大幅降低，為56.08～360.19和17.40～76.67。不同建模樣本數(shù)量的GWR模型對(duì)SOM的解釋能力差異較小。3）建模樣點(diǎn)數(shù)量（除建模樣本=243）對(duì)GWR預(yù)測(cè)制圖結(jié)果的精度影響不大，RMSE介于5.56～5.75 g/kg之間，MAE介于3.87～4.05 g/kg之間，2介于0.52～0.48之間，均優(yōu)于全部建模樣點(diǎn)的普通克里格插值驗(yàn)證結(jié)果。該研究可為樣點(diǎn)數(shù)較少的省級(jí)尺度地區(qū)SOM空間建模與制圖提供借鑒。

土壤; 有機(jī)質(zhì)；模型；地理加權(quán)回歸；數(shù)字土壤制圖；地形平緩區(qū)；江蘇省

0 引 言

土壤有機(jī)質(zhì)（soil organic matter，SOM）是土壤肥力和質(zhì)量的重要指標(biāo)，是陸地生態(tài)系統(tǒng)碳儲(chǔ)量的重要組成部分，在農(nóng)業(yè)生產(chǎn)和環(huán)境保護(hù)等方面有著重要的意義[1-2]。研究SOM空間分布是評(píng)估區(qū)域土壤碳儲(chǔ)量、實(shí)現(xiàn)土壤可持續(xù)利用的前提，也是土壤質(zhì)量研究的重要內(nèi)容。20世紀(jì)90年代興起的計(jì)量土壤學(xué)和數(shù)字土壤制圖研究，主要基于土壤-景觀模型，利用土壤在地形、植被等景觀環(huán)境上的差異，采用線性回歸、廣義線性回歸、地統(tǒng)計(jì)學(xué)、模糊聚類、地理加權(quán)回歸（geographically weighted regression, GWR）、機(jī)器學(xué)習(xí)等數(shù)學(xué)方法和地理信息系統(tǒng)等技術(shù)進(jìn)行土壤屬性的空間建模，為獲取柵格化的土壤空間分布信息提供解決途徑[3-8]。建立較優(yōu)的土壤屬性空間預(yù)測(cè)模型是獲取較準(zhǔn)確的土壤空間分布信息的關(guān)鍵。但是在平原或地形平緩地區(qū)，由于地形差異減小，植被等景觀環(huán)境因素差異也變小，地形和植被等條件不能反映土壤的空間變異，因此不能有效地用于土壤屬性的空間建模。

土壤是具有高度變異的時(shí)空連續(xù)體，在形成發(fā)育過(guò)程中成土因子對(duì)土壤的作用是非線性，在較大區(qū)域中土壤屬性的主要影響因素較多，這種非線性作用更為明顯。周濤等[9]研究發(fā)現(xiàn)在<10 ℃的地區(qū)，有機(jī)碳儲(chǔ)量與年均溫的負(fù)相關(guān)性最強(qiáng)；在10～20 ℃的地區(qū)，受與降水正相關(guān)的影響，有機(jī)碳儲(chǔ)量與年均溫呈正相關(guān)；而在>20 ℃的地區(qū)，有機(jī)碳儲(chǔ)量與溫度和降水的相關(guān)性都很差。趙明松等[10]研究表明江蘇省土壤質(zhì)地對(duì)SOM的作用在不同土壤類型和地理區(qū)域間存在著差異，總體上由北至南相關(guān)性減弱。因此在較大區(qū)域尺度上，利用傳統(tǒng)的全局模型（如線性回歸，普通回歸克里格等）并不能很好地或較為準(zhǔn)確地?cái)M合環(huán)境因子與土壤屬性之間的非線性關(guān)系。

GWR模型由英國(guó)地理學(xué)家Fotheringham等[11]提出，針對(duì)自變量和因變量間的關(guān)系在不同空間上的變化進(jìn)行非參數(shù)局部空間回歸建模，模型中自變量的回歸系數(shù)隨空間位置而變化，對(duì)于空間數(shù)據(jù)具有較強(qiáng)的局部分析能力[12-13]。GWR探測(cè)空間關(guān)系的非平穩(wěn)性，廣泛用于區(qū)域經(jīng)濟(jì)[14-16]、生態(tài)環(huán)境[17-20]、土壤建模制圖[21-25]等研究。劉瓊峰等[17]利用GWR模型分析了長(zhǎng)沙城郊農(nóng)田土壤Pb、Cd含量的空間結(jié)構(gòu)及影響因素。郭龍等[26]等對(duì)比發(fā)現(xiàn)協(xié)同克里格和GWR模型在土壤預(yù)測(cè)制圖中均具有較高精度。Zeng等[25]研究表明在小流域尺度和地區(qū)尺度上GWR模型均能提高SOM建模的精度。Song等[21]、Yang等[23]研究表明GWR模型能夠提高較大復(fù)雜景觀區(qū)域的土壤有機(jī)碳和電導(dǎo)率的建模預(yù)測(cè)精度。Zhang等[24]基于GWR模型利用年均降雨、土壤類型和土地利用等因子預(yù)測(cè)了愛(ài)爾蘭全國(guó)的土壤有機(jī)碳的空間分布。上述研究表明了GWR模型在較大區(qū)域中土壤屬性空間建模的優(yōu)勢(shì)。

針對(duì)較大區(qū)域尺度、地形平緩區(qū)土壤屬性空間建模中存在的問(wèn)題，以及省級(jí)尺度上采集大量土壤樣品的困難。本研究選擇江蘇省為例，以第二次土壤普查數(shù)據(jù)為基礎(chǔ)，選擇柵格氣象數(shù)據(jù)等，采用GWR模型開(kāi)展地形平緩區(qū)SOM建模預(yù)測(cè)研究，并探討樣點(diǎn)數(shù)量對(duì)GWR建模精度的影響。研究結(jié)果可為省級(jí)尺度區(qū)域的SOM空間預(yù)測(cè)制圖提供借鑒，為江蘇省土壤肥力和質(zhì)量管理、提升區(qū)域土壤固碳能力等提供數(shù)據(jù)基礎(chǔ)。

1 材料與方法

1.1 研究區(qū)概況

江蘇?。?16°18′～121°57′ E，30°45′～35°20′ N）地處中國(guó)大陸東部，長(zhǎng)江、淮河下游，面積10.26萬(wàn)km2。全省處于亞熱帶向暖溫帶過(guò)渡地帶，年均氣溫13～16 ℃，年均降水量800～1 200 mm，年均日照2 100～2 600 h。全省平均海拔13 m，海拔低于10 m或坡度小于1%的區(qū)域占全省陸地面積的90%。按地質(zhì)地貌差異，分為低山丘陵崗地區(qū)和平原區(qū)，平原、低山丘陵和水域面積分別占69%、14%和17%。低山丘陵分布在西南和東北部；平原由北至南為徐淮黃泛平原、里下河平原、濱海平原、沿江平原和太湖平原。主要的土壤類型有水稻土、潮土、砂姜黑土、濱海鹽土、黃棕壤、黃褐土、棕壤等。水田和旱地分別占全省面積的42.46%和23.57%（2006年）。江蘇省是中國(guó)重要的糧食生產(chǎn)基地，隨著農(nóng)業(yè)現(xiàn)代化和城市化的迅猛發(fā)展，土壤保護(hù)與利用之間的矛盾比較突出，對(duì)土壤空間分布信息的需求非常迫切。

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

土壤數(shù)據(jù)來(lái)源于江蘇省第二次土壤普查資料（1979－1982年）中記錄的典型土壤剖面，主要包括《江蘇土種志》[27]和各市、縣土種志共60本。土壤剖面數(shù)據(jù)主要包括采樣點(diǎn)位置和景觀環(huán)境描述和土壤理化性質(zhì)（主要有SOM、pH值、物理性黏粒和砂粒含量等），篩選出土壤數(shù)據(jù)比較完備的1 519個(gè)樣點(diǎn)（圖1）。本研究以表層（0～20 cm）SOM含量為研究對(duì)象。對(duì)于深度大于20 cm的土層數(shù)據(jù)不做處理，對(duì)于小于20 cm的土層，以土層深度為權(quán)重取0～20 cm內(nèi)的所有土層的加權(quán)平均值。SOM含量采用重鉻酸鉀（K2Cr2O7）氧化-滴定法測(cè)定。

年均氣溫和年均降水量數(shù)據(jù)，來(lái)自中國(guó)農(nóng)業(yè)科學(xué)院農(nóng)業(yè)資源與農(nóng)業(yè)區(qū)劃研究所中國(guó)生態(tài)環(huán)境背景層面建造項(xiàng)目完成的柵格數(shù)據(jù)（1 km分辨率），為1980－1999年的逐月平均值計(jì)算合成。在ArcGIS軟件中，根據(jù)地形圖、行政區(qū)劃變更信息等資料，將剖面點(diǎn)位置空間化，從氣候柵格數(shù)據(jù)中提取各剖面點(diǎn)的年均氣溫和年均降水?dāng)?shù)據(jù)。

江蘇省地形以平原為主，研究發(fā)現(xiàn)SOM含量與海拔和坡度等地形因子的相關(guān)性不足0.01，利用逐步回歸建模分析時(shí)，無(wú)法進(jìn)入模型，可用于SOM空間建模的有效變量較少。對(duì)于省級(jí)區(qū)域尺度，氣候等條件對(duì)土壤屬性的影響較大。鑒于以上背景和原因本文選擇年均氣溫、年均降雨、物理性黏粒含量、土壤pH值等環(huán)境變量進(jìn)行SOM含量的空間建模。

圖1 土壤樣點(diǎn)分布圖

1.3 研究方法

1.3.1 地理加權(quán)回歸

GWR模型的一般形式如下

式中(u,v)是第個(gè)樣點(diǎn)的空間坐標(biāo)；為驗(yàn)證樣點(diǎn)個(gè)數(shù)；Y和X是因變量和自變量集X在空間位置(u,v)處的實(shí)測(cè)值；0(u,v)為在空間位置(u,v)處的常數(shù)項(xiàng)；系數(shù)(u,v)是連續(xù)函數(shù)(,)在點(diǎn)的值；為符合(0,2)分布的誤差項(xiàng)。通常采用加權(quán)最小二乘法進(jìn)行系數(shù)的局部估計(jì)。一般采用高斯（Gaussian）型或雙重平方（bi-square）型空間核函數(shù)，利用空間位置(u,v)與周圍一定距離（又稱帶寬）內(nèi)觀測(cè)樣點(diǎn)的距離來(lái)估算其權(quán)重。本研究采用雙重平方函數(shù)估算權(quán)重，采用黃金分割搜索（golden section search）方法、赤池信息量準(zhǔn)則（akaike information criterion，AIC）[28]，確定最優(yōu)帶寬，評(píng)價(jià)擬合模型的優(yōu)劣。根據(jù)Fotheringham等[11]提出的方法準(zhǔn)則，當(dāng)GWR模型AIC值最小時(shí)，帶寬b為最優(yōu)帶寬，此時(shí)模型擬合最優(yōu)。詳細(xì)的GWR模型擬合方法參考文獻(xiàn)[11]。如果公式（1）中系數(shù)(,)在空間中任意一點(diǎn)保持不變，則為全局模型，即傳統(tǒng)回歸模型，通常用普通最小二乘法（ordinary least square, OLS）估計(jì)參數(shù)。

1.3.2 數(shù)據(jù)處理與分析

本研究采用相關(guān)性分析、回歸分析等，研究年均氣溫、年均降雨、物理性黏粒、土壤pH值對(duì)江蘇省SOM的影響，采用容差和方差膨脹因子（VIF）對(duì)環(huán)境因子與SOM關(guān)系進(jìn)行共線性檢驗(yàn)，在SPSS 18.0中完成。在利用GWR進(jìn)行SOM空間建模時(shí)，選擇302個(gè)樣本作為驗(yàn)證集，1 217個(gè)樣本為建模集。為比較不同樣本數(shù)量下GWR建模的精度差異，從建模集中分別隨機(jī)抽取100%（1 217個(gè)）、80%（973個(gè)）、60%（730個(gè)）、40%（486個(gè)），20%（243個(gè)）的樣點(diǎn)進(jìn)行建模，并選擇最優(yōu)模型進(jìn)行SOM預(yù)測(cè)制圖。利用ArcGIS 9.3的地統(tǒng)計(jì)模塊中創(chuàng)建子集（create subsets）功能完成不同數(shù)量建模樣本的抽樣。利用GWR4.0軟件進(jìn)行GWR建模與制圖。同時(shí)，利用ArcGIS 9.3的地統(tǒng)計(jì)模塊以全部建模集（1 217個(gè)）進(jìn)行普通克里格插值（ordinary kriging，OK），繪制SOM空間分布圖與GWR預(yù)測(cè)制圖結(jié)果做對(duì)比。

在SOM空間預(yù)測(cè)制圖時(shí)，需要環(huán)境因子（建模的自變量）的空間分布。本研究中物理性黏粒含量和pH數(shù)據(jù)為點(diǎn)狀數(shù)據(jù)，利用1 519個(gè)樣點(diǎn)進(jìn)行OK插值生成其空間分布，用于SOM空間制圖。綜合考慮氣象數(shù)據(jù)的分辨率、研究區(qū)較大及數(shù)據(jù)運(yùn)算等問(wèn)題，制圖結(jié)果均采用1 km空間分辨率。

1.3.3 模型評(píng)價(jià)

選取平均誤差（mean error，ME）、平均絕對(duì)誤差（mean absolute error，MAE）、均方根誤差（root mean of squared error，RMSE）和決定系數(shù)（determination coefficient，2）等指標(biāo)評(píng)價(jià)建模結(jié)果，其計(jì)算公式如下

2 結(jié)果與分析

2.1 SOM主要影響因素分析

江蘇省SOM質(zhì)量分?jǐn)?shù)范圍較大（= 1 519），為1.3～52.4 g/kg，極差為51.1 g/kg；平均為16.55 g/kg，標(biāo)準(zhǔn)差為8.49 g/kg（表1）。變異系數(shù)為51.36%，屬于中等變異強(qiáng)度。偏度和峰度系數(shù)分別為1.12和1.75，經(jīng)自然對(duì)數(shù)變換后符合正態(tài)分布。不同樣點(diǎn)數(shù)量的建模集的SOM頻率分布均符合對(duì)數(shù)正態(tài)分布。

相關(guān)分析表明，江蘇省SOM與年均氣溫極顯著正相關(guān)，相關(guān)系數(shù)為0.47（<0.001，=1 519），以年均降雨量為控制變量時(shí)二者的偏相關(guān)系數(shù)為0.21（<0.001，= 1 519），即在溫度較高的長(zhǎng)江以南地區(qū)SOM含量較高，溫度較低的蘇北地區(qū)含量較低。這與周濤等[9]、許信旺等[29]研究結(jié)果一致。SOM與年均降雨量極顯著正相關(guān)，相關(guān)系數(shù)為0.44（<0.001，=1 519），以年均氣溫為控制變量時(shí)二者的相關(guān)性降低，偏相關(guān)系數(shù)為0.11（<0.001，=1 519）。這表明SOM與降雨量的相關(guān)性是由年均降雨與年均氣溫的相關(guān)性（=0.82，<0.001）引起的。SOM與土壤物理性黏粒和砂粒含量相關(guān)系數(shù)為0.46和?0.45（<0.001，=1 519）。SOM與pH值的相關(guān)系數(shù)為?0.30（<0.001，=1 519）。線性回歸分析表明年均氣溫、年均降雨、土壤質(zhì)地和pH值對(duì)SOM含量（=1 519）變異的綜合解釋能力為45.3%（0.001）?；貧w模型中的各因子的容差介于0.31～0.98，均大于0.1，VIF介于0.02～3.25，均不大于7.5，表明自變量之間不存在高度共線性[30]。

表1 江蘇省SOM含量統(tǒng)計(jì)值

注：括號(hào)內(nèi)為自然對(duì)數(shù)變換值。

Note: values in the brackets were logarithm transformed.

2.2 SOM空間自相關(guān)分析

利用ArcGIS軟件空間統(tǒng)計(jì)功能分析不同建模樣本集的SOM的全局和局部空間自相關(guān)性。在ArcGIS中計(jì)算全局Moran’s I值時(shí)，在0～100 km范圍內(nèi)按照10 km間隔設(shè)定不同的空間距離，計(jì)算全局Moran’s I值。圖2為不同建模集的全局Moran’s I值，結(jié)果表明不同建模樣本集的SOM在不同的空間尺度上均存在不同程度的極顯著空間自相關(guān)性，全局Moran’s I值介于0.25～0.61（<0.001）。隨著空間距離的增大，不同建模集的SOM的全局Moran’s I值均逐漸減小，表明其空間自相關(guān)性隨著空間距離的增大而遞減。在相同的空間距離上，除建模集為=486時(shí)的全局Moran’s I值稍小外，其他建模集的SOM的全局Moran’s I值總體上相差不大。

圖3為建模集為=1 217時(shí)的局部空間自相關(guān)聚類圖。由圖知，江蘇省SOM含量存在空間聚集特征，“高-高”（high-high）聚集區(qū)（SOM含量較高且空間差異較小的樣點(diǎn)聚集分布）主要分布在蘇中和蘇南地區(qū)，“低-低”（low-low）聚集區(qū)（SOM含量較低且空間差異較小的樣點(diǎn)聚集分布）主要分布在蘇北地區(qū)；“高-低”（high-low）聚集區(qū)（SOM含量較高的樣點(diǎn)被含量較低的樣點(diǎn)包圍，且含量的空間差異顯著）主要分布在蘇北地區(qū)，“低-高”（low-high）聚集區(qū)（SOM含量較低的樣點(diǎn)被含量較高的樣點(diǎn)包圍，且含量的空間差異顯著）主要分布在蘇南地區(qū)。不同建模集的局部空間自相關(guān)聚類圖結(jié)果相似，沒(méi)有一一作圖展示。以上分析表明不同建模集的SOM在不同的空間尺度均存在顯著的全局和局部空間自相關(guān)性，可以采用GWR模型進(jìn)行研究區(qū)SOM建模[17]。

注：顯著性為P<0.001。

注：圖例代表不同聚類聚集區(qū)，下同。

2.3 基于GWR的SOM建模分析

將不同建模集的SOM經(jīng)過(guò)自然對(duì)數(shù)變換后，分別進(jìn)行GWR建模和全局建模。為了比較不同自變量在建模中的貢獻(xiàn)大小，首先在建模前將各自變量進(jìn)行標(biāo)準(zhǔn)化變換到同一量綱，然后通過(guò)其系數(shù)來(lái)對(duì)比各自變量在SOM建模中的作用。

表2為不同建模集的SOM的GWR建模和全局建模結(jié)果。一般認(rèn)為GWR模型和全局模型的AIC相差大于3，即使考慮到GWR模型的復(fù)雜性，其也明顯優(yōu)于全局模型。

表2 不同樣點(diǎn)數(shù)量的GWR建模和全局建模結(jié)果比較

表2顯示對(duì)于不同建模集的GWR模型的AIC均比相應(yīng)的全局模型要小，相差56.08～360.19之間，表明GWR模型明顯優(yōu)于全局模型。不同建模集的GWR模型的調(diào)整決定系數(shù)（2adj）均大于相應(yīng)的全局模型，其2adj較相應(yīng)的全局模型增加0.15～0.20之間，表明對(duì)于相同的建模集，GWR模型對(duì)SOM變異的解釋能力大于全局模型。對(duì)于不同建模集SOM的GWR模型，2adj介于0.54～0.62，總體上隨著建模樣本數(shù)量的減少GWR模型對(duì)SOM的解釋能力稍有差異，但差異不大。這說(shuō)明對(duì)于江蘇省SOM的建模估計(jì)，采用較少樣點(diǎn)數(shù)進(jìn)行GWR建模時(shí)，也可滿足一定的精度需求（如建模樣本= 486和= 1 217的模型擬合度較接近）。對(duì)于不同建模集SOM的全局模型，其2adj介于0.39～0.44之間。

不同建模集的GWR模型的殘差平方和（residual sum of squares,RSS）均比相應(yīng)的全局模型減少17.40～76.67之間。GWR模型和全局模型殘差的方差分析結(jié)果（表3）顯示，不同建模集的GWR模型的殘差較全局模型的殘差顯著降低，其檢驗(yàn)值均達(dá)到0.01極顯著性水平，表明GWR模型的擬合效果較全局模型有顯著改善。上述這些結(jié)果表明，不同建模集的GWR建模結(jié)果均優(yōu)于全局建模，且當(dāng)建模樣點(diǎn)數(shù)較少時(shí)GWR模型也能有較好的建模效果。

表4為不同樣點(diǎn)數(shù)量的GWR模型回歸系數(shù)和全局模型回歸系數(shù)。結(jié)果顯示，在不同樣點(diǎn)數(shù)量的全局建模中（除樣點(diǎn)數(shù)量= 243），同一環(huán)境因子對(duì)SOM預(yù)測(cè)的作用相似，如各全局模型中黏粒的回歸系數(shù)介于0.24～0.25，相差較?。徊煌兞繉?duì)SOM的作用大小依次為：黏粒> 年均溫> 年均降雨> 土壤pH值。在不同樣點(diǎn)數(shù)量的GWR建模中（除樣點(diǎn)數(shù)量= 243），從回歸系數(shù)的中值和平均值來(lái)看，各環(huán)境因子對(duì)SOM的作用大小規(guī)律與全局建模一致。當(dāng)建模樣點(diǎn)數(shù)量較少（= 243時(shí)），環(huán)境因子與SOM之間的關(guān)系不能完全地表達(dá)和刻畫。GWR模型中，各環(huán)境因子的回歸系數(shù)存在不同幅度的變化范圍，同一環(huán)境因子的回歸系數(shù)變化范圍相差不大，其中黏粒的回歸系數(shù)變化范圍最穩(wěn)定且系數(shù)均為正（表4），體現(xiàn)了黏粒含量對(duì)SOM含量的正相關(guān)。

表3 GWR模型和全局建模殘差的方差分析結(jié)果

注：**為顯著性< 0.01,下同。

Note: **is significant difference at 0.01 level, the same below.

圖4為GWR模型回歸系數(shù)（建模樣點(diǎn)= 1 217），各環(huán)境因子對(duì)SOM的作用存在顯著的空間差異，即SOM與這些環(huán)境因子之間的關(guān)系存在空間非平穩(wěn)性[19]。各因子對(duì)SOM的影響程度可通過(guò)對(duì)應(yīng)的回歸系數(shù)來(lái)解釋，回歸系數(shù)的絕對(duì)值大小反映該因子對(duì)SOM影響的強(qiáng)烈程度，正反映因子對(duì)SOM的正影響或正相關(guān)，反之亦然[13, 30]。如，黏粒的回歸系數(shù)變化范圍比其他因素回歸系數(shù)的變化范圍較小，表明在研究區(qū)內(nèi)黏粒對(duì)SOM的影響相比其他因素較穩(wěn)定。黏粒對(duì)SOM的影響在江蘇省東北部和西南部最小，在中南部最大，且均表現(xiàn)為正的影響（圖4a），這與趙明松等[10]在該地區(qū)的研究結(jié)果相似。年均溫對(duì)SOM的作用在大部分地區(qū)表現(xiàn)為正相關(guān)（圖4b）；年均降雨對(duì)SOM的作用在江蘇省東部地區(qū)表現(xiàn)為負(fù)相關(guān)，在西部表現(xiàn)為正相關(guān)。雖然研究區(qū)內(nèi)年均溫和年均降雨顯著正相關(guān)，但是對(duì)于SOM的作用強(qiáng)度在空間上顯著差異，相應(yīng)系數(shù)的空間格局差異較大。土壤pH值對(duì)SOM的作用，在西北-東南方向上向南北兩側(cè)逐漸變化，總體上由負(fù)相關(guān)轉(zhuǎn)變?yōu)檎嚓P(guān)。

表4 不同樣點(diǎn)數(shù)量的GWR模型回歸系數(shù)

注：*顯著性< 0.05,下同。

Note: *, ** Siginificant difference at 0.05 level, the same below.

圖4 GWR建?；貧w系數(shù)空間分布圖（n=1 217）

2.4 基于GWR的SOM預(yù)測(cè)制圖及模型評(píng)價(jià)

圖5為江蘇省SOM含量空間分布預(yù)測(cè)圖，其中a～e為不同樣點(diǎn)數(shù)量的GWR建模預(yù)測(cè)圖，f為采用全部建模樣點(diǎn)（=1 217）的OK插值圖。

圖5 SOM含量空間分布預(yù)測(cè)圖

SOM的OK插值的變異函數(shù)采用球狀模型擬合，塊金值為0.13，偏基臺(tái)值為0.10，變程為75 km，步長(zhǎng)為10 km。結(jié)果顯示，不同樣點(diǎn)數(shù)量的GWR預(yù)測(cè)結(jié)果的范圍介于3.38～45.08 g/kg，平均值介于14.94～15.34 g/kg，標(biāo)準(zhǔn)差介于5.99～7.15 g/kg，均比建模樣點(diǎn)的相應(yīng)統(tǒng)計(jì)值稍小。不同樣點(diǎn)數(shù)量的GWR預(yù)測(cè)結(jié)果的空間分布格局基本一致：SOM含量的最高值主要分布在太湖平原和里下河平原，從蘇南到蘇北地區(qū)逐漸降低；在研究區(qū)西南部的部分地區(qū)SOM的空間分布細(xì)節(jié)存在細(xì)微差異（圖5）?？傮w上，不同樣點(diǎn)數(shù)量的GWR預(yù)測(cè)結(jié)果與采用全部建模樣點(diǎn)的OK預(yù)測(cè)結(jié)果具有較相似的空間分布格局。

表5為SOM建模的獨(dú)立驗(yàn)證結(jié)果。不同樣點(diǎn)數(shù)量（除= 243時(shí)）的GWR建模的驗(yàn)證結(jié)果差異不大：RMSE介于5.56～5.75 g/kg，MAE介于3.87～4.05 g/kg，均小于全部建模樣點(diǎn)的OK驗(yàn)證結(jié)果；2介于0.48～0.52，均高于OK插值的驗(yàn)證結(jié)果。不同樣點(diǎn)數(shù)量的GWR驗(yàn)證結(jié)果的ME介于?0.58～?0.22 g/kg，表明GWR建模結(jié)果總體上存在不同程度的高估現(xiàn)象。樣點(diǎn)數(shù)=243時(shí)GWR建模的驗(yàn)證結(jié)果（RMSE=6.31，ME=?0.22，MAE=4.47）較其他稍差，但其預(yù)測(cè)的空間分布格局與其他結(jié)果較相似。上述結(jié)果表明，在江蘇省利用GWR進(jìn)行SOM空間預(yù)測(cè)制圖時(shí)，樣點(diǎn)數(shù)量對(duì)建模精度影響不大，較少的樣點(diǎn)數(shù)量也可以獲得較好的制圖結(jié)果。與OK插值方法比較發(fā)現(xiàn)，在樣點(diǎn)數(shù)最大時(shí)，GWR建模精度提高幅度有限，2僅提高了0.06。但從另一個(gè)方面看，如果樣點(diǎn)較少時(shí)（如486個(gè)）采用GWR建模，能夠獲得和較多樣點(diǎn)的OK插值結(jié)果的相似精度，對(duì)減少制圖成本也有一定的貢獻(xiàn)。

表5 不同樣本數(shù)量的GWR建模的驗(yàn)證結(jié)果

圖6為不同樣本數(shù)量的GWR模型殘差的全局Moran’I指數(shù)。GWR模型殘差的全局Moran’s I值介于0.01～0.12（>0.1），隨著空間距離的增大其值均逐漸減小，表明在不同的空間尺度上不同樣點(diǎn)數(shù)量的GWR模型殘差均不存在顯著空間自相關(guān)性。在相同的空間距離上，GWR殘差的全局Moran’s I值總體上相差不大。圖7為GWR模型殘差（= 1 217）的局部空間自相關(guān)聚類圖。由圖知，絕大部分樣點(diǎn)的建模均不存在顯著局部空間自相關(guān)性，僅有少量樣點(diǎn)（= 45）的殘差存在著“高-高”（殘差較高且空間差異較小的樣點(diǎn)聚集分布區(qū)）、“低-低”和“高-低”聚類。不同建模集的GWR殘差的局部空間自相關(guān)聚類圖結(jié)果相似，沒(méi)有一一作圖展示。不同樣點(diǎn)數(shù)量GWR建模的殘差空間自相關(guān)性較弱，也間接表明GWR擬合的模型較優(yōu)。

注：不顯著（P=0.1）。

圖7 GWR模型殘差局部空間聚類圖（n=1 217）

3 結(jié) 論

1）江蘇省SOM含量在不同空間尺度上存在極顯著的全局和局部空間自相關(guān)性。SOM含量空間分布以空間聚集特征為主，“高-高”聚集區(qū)主要分布在蘇中和蘇南地區(qū)，“低-低”聚集區(qū)主要分布在蘇北地區(qū)。不同建模集的全局自相關(guān)性和局部空間自相關(guān)聚類圖結(jié)果相似。

2）GWR建模結(jié)果均優(yōu)于全局建模，對(duì)SOM空間變異的解釋能力更高，不同建模集的GWR的2adj較全局建模均有明顯提高，其AIC和RSS均比全局模型有大幅降低，擬合效果有較大提升。隨著建模樣本數(shù)量的減少GWR模型對(duì)SOM的解釋能力差異不大。GWR建模的殘差在不同的空間尺度上均不存在空間自相關(guān)性。

3）利用GWR模型進(jìn)行SOM空間預(yù)測(cè)制圖時(shí)，建模樣點(diǎn)數(shù)量對(duì)建模精度影響不大，不同樣點(diǎn)數(shù)量的GWR建模結(jié)果與采用全部樣點(diǎn)的OK插值結(jié)果具有相似的空間格局。樣點(diǎn)數(shù)較少時(shí)，在大區(qū)域利用GWR模型進(jìn)行SOM預(yù)測(cè)制圖，也可以達(dá)到一定的精度要求。

數(shù)字土壤制圖過(guò)程中，一般要求環(huán)境因子為空間柵格形式。本文中物理性黏粒含量和pH值為點(diǎn)狀數(shù)據(jù)，在預(yù)測(cè)制圖中首先通過(guò)插值其空間分布，再用于SOM預(yù)測(cè)制圖，可能增加二次誤差，對(duì)最終的制圖結(jié)果可能會(huì)有一定程度的影響，這些有待進(jìn)一步的深入研究。

[1] Yadav V, Malanson G. Progress in soil organic matter research: Litter decomposition, modelling, monitoring and sequestration [J]. Progress in Physical Geography, 2007, 31(2): 131－154.

[2] Davidson E A, Trumbore S E, Amundson R. Biogeochemistry: Soil warming and organic carbon content[J]. Nature, 2000, 408(6814): 789－790.

[3] McBratney A B, Santos M L M, Minasny B. On digital soil mapping [J]. Geoderma, 2003, 117: 3－52.

[4] Thompson J A, Pena-Yewtukhiw E M, Grove J H. Soil-landscape modeling across a physiographic region: Topographic patterns and model transportability [J]. Geoderma, 2006, 133: 57－70.

[5] Qi F, Zhu A X, Harrower M, et al. Fuzzy soil mapping based on prototype category theory[J]. Geoderma, 2006, 136: 774－787.

[6] Zhu A X, Qi F, Moore A, et al. Prediction of soil properties using fuzzy membership values [J]. Geoderma, 2010, 158: 199－206.

[7] de Bruin S, Stein A. Soil-landscape modelling using fuzzy c-means clustering of attribute data derived from a digital elevation model (DEM) [J]. Geoderma, 1998, 83: 17－33.

[8] 朱阿興，楊琳，樊乃卿，等. 數(shù)字土壤制圖研究綜述與展望[J]. 地理科學(xué)進(jìn)展，2018，37(1)：66－78.

Zhu Axing, Yang Lin, Fan Naiqing, et al. The review and outlook of digital soil mapping[J]. Progress in Geography, 2018, 37(1): 66－78. (in Chinese with English abstract)

[9] 周濤，史培軍，王紹強(qiáng). 氣候變化及人類活動(dòng)對(duì)中國(guó)土壤有機(jī)碳儲(chǔ)量的影響[J]. 地理學(xué)報(bào)，2003，58(5)：727－734.

Zhou Tao, Shi Peijun, Wang Shaoqiang. Impacts of climate and human activities on soil carbon storage in China[J]. Acta Geographica Sinica, 2003, 58(5): 727－734. (in Chinese with English abstract)

[10] 趙明松，張甘霖，李德成，等. 江蘇省土壤有機(jī)質(zhì)變異及其主要影響因素[J]. 生態(tài)學(xué)報(bào)，2013，33(16)：5058－5066.

Zhao Mingsong, Zhang Ganlin, Li Decheng, et al. Variability of soil organic matter and its main factors in Jiangsu Province[J]. Acta Ecologica Sinica, 2013, 33(16): 5058－5066. (in Chinese with English abstract)

[11] Fotheringham A S, Brunsdon C, Charlton M. Geographically Weighted Regression—the Analysis of Spatially Varying Relationships[M]. Chichester.UK: John Wiley & Sons Inc., 2002: 42－46.

[12] 王遠(yuǎn)飛，何洪林. 空間數(shù)據(jù)分析方法[M]. 北京：科學(xué)出版社，2007：132－141.

[13] 瞿明凱，李衛(wèi)東，張傳榮，等. 地理加權(quán)回歸及其在土壤和環(huán)境科學(xué)上的應(yīng)用前景[J]. 土壤，2014，46(1)：15－22.

Qu Mingkai, Li Weidong, Zhang Chuanrong, et al. Geographically weighted regression and its application prospect in soil and environmental sciences[J]. Soils, 2014, 46(1): 15－22. (in Chinese with English abstract)

[14] 呂萍，甄輝. 基于GWR模型的北京市住宅用地價(jià)格影響因素及其空間規(guī)律研究[J]. 經(jīng)濟(jì)地理，2010，30(3)：472－478.

Lü Ping, Zhen Hui. Affecting factors research of Beijing residential land price based on GWR model[J]. Economic Geography, 2010, 30(3): 472－478. (in Chinese with English abstract)

[15] 湯慶園，徐偉，艾福利. 基于地理加權(quán)回歸的上海市房?jī)r(jià)空間分異及其影響因子研究[J]. 經(jīng)濟(jì)地理，2012，32(2)：52－58.

Tang Qingyuan, Xu Wei, Ai Fuli. A GWR based study on spatial pattern and structural determinants of Shanghai’s housing price[J]. Economic Geography, 2012, 32(2): 52－58. (in Chinese with English abstract)

[16] 劉衛(wèi)東，劉紅光，范曉梅，等. 地區(qū)間貿(mào)易流量的產(chǎn)業(yè)：空間模型構(gòu)建與應(yīng)用[J]. 地理學(xué)報(bào)，2012，67(2)：147－156.

Liu Weidong, Liu Hongguang, Fan Xiaomei, et al. Sector-specific spatial statistic model for estimating inter-regional Trade tflows: A case study of agricultural, chemical and electronic sectors in China[J]. Acta Geographica Sinica, 2012, 67(2): 147－156. (in Chinese with English abstract)

[17] 劉瓊峰，李明德，段建南，等. 農(nóng)田土壤鉛、鎘含量影響因素地理加權(quán)回歸模型分析[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2013，29(3)：225－234.

Liu Qiongfeng, Li Mingde, Duan Jiannan, et al. Analysis on influence factors of soil Pb and Cd in agricultural soil of Changsha suburb based on geographically weighted regression model [J]. Transactions of the Chinese Society of Agricultural Engineering(Transactions of the CSAE), 2013, 29(3): 225－234. (in Chinese with English abstract)

[18] 霍霄妮，李紅，孫丹峰，等. 北京耕地土壤重金屬空間自回歸模型及影響因素[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2010，26(5)：78－82.

Huo Xiaoni, Li Hong, Sun Danfeng, et al. Spatial autogression model for heavy metals in cultivated soils of Beijing[J]. Transactions of the Chinese Society of Agricultural Engineering(Transactions of the CSAE), 2010, 26(5): 78－82. (in Chinese with English abstract)

[19] 李錦芬，瞿明凱，黃標(biāo)，等. 區(qū)域土壤CEC與相關(guān)控制因子的空間非平穩(wěn)關(guān)系評(píng)估[J]. 土壤學(xué)報(bào)，2017，54(3)：639－647.

Li Jinfen, Qu Mingkai, Huang Biao, et al. Spatially non-stationary relationships between cation exchange capacity and related control factors[J]. Acta Pedologica Sinica, 2017, 54(3): 639－647. (in Chinese with English abstract)

[20] 王景雷，康紹忠，孫景生，等. 基于PCA和GWR的作物需水量空間分布估算[J]. 科學(xué)通報(bào)，2013，58(12)：1131－1139.

Wang Jinglei, Kang Shaozhong, Sun Jingsheng, et al. Estimation of crop water requirement based on principal component analysis and geographically weighted regression[J]. China Science Bulletin, 2013, 58, 1131－1139. (in Chinese with English abstract)

[21] Song X D, Brus D J, Liu F, et al. Mapping soil organic carbon content by geographically weighted regression: A case study in the Heihe River Basin, China[J]. Geoderma, 2016, 261: 11－22.

[22] Mishra U, Lal R, Liu D S, et al. Predicting the spatial variation of the soil organic carbon pool at a regional scale[J]. Soil Science Society of America Journal, 2010, 74(3): 906－914.

[23] Yang S H, Liu F, Song X D, et al. Mapping topsoil electrical conductivity by a mixed geographically weighted regression kriging: A case study in the Heihe River Basin, northwest China[J]. Ecological Indicators, 2019, 102: 252－264.

[24] Zhang C S, Tang Y, Xu X L, et al. Towards spatial geochemical modelling: Use of geographically weighted regression for mapping soil organic carbon contents in Ireland[J]. Applied Geochemistry, 2011, 26(7): 1239－1248.

[25] Zeng C, Yang L, Zhu A X, et al. Mapping soil organic matter concentration at different scales using a mixed geographically weighted regression method[J]. Geoderma, 2016, 281: 69－82.

[26] 郭龍，張海濤，陳家贏，等. 基于協(xié)同克里格插值和地理加權(quán)回歸模型的土壤屬性空間預(yù)測(cè)比較[J]. 土壤學(xué)報(bào)，2012，49(5)：1037－1042.

Guo Long, Zhang Haitao, Chen Jiaying, et al. Comparison between cokring model and geographically weighted regression model in spatial prediction of soil attributes[J]. Acta Pedologica Sinica, 2012, 49(5): 1037－1042. (in Chinese with English abstract)

[27] 江蘇省土壤普查辦公室. 江蘇土種志[M]. 南京：江蘇科學(xué)技術(shù)出版社，1996：5－105.

[28] Akaike H. A new look at the statistical model identification[J]. IEEE Transactions on Automatic Control, 1974, 19(6): 716－723.

[29] 許信旺，潘根興，汪艷林，等. 中國(guó)農(nóng)田耕層土壤有機(jī)碳變化特征及控制因素[J]. 地理研究，2009，28(3)：601－612.

Xu Xinwang, Pan Genxing, Wang Yanlin, et al. Research of changing characteristics and control factors of farmland topsoil organic carbon in China[J]. Geographical Research, 2009, 28(3): 601－612. (in Chinese with English abstract)

[30] 許爾琪. 基于地理加權(quán)回歸的石漠化影響因子分布研究[J]. 資源科學(xué)，2017，39(10)：1975－1988.

Xu Erqi. Spatial variation in drivers of karst rocky desertification based on geographically weighted regression model[J]. Resources Science, 2017, 39(10): 1975－1988. (in Chinese with English abstract)

Spatial modeling of soil organic matter over low relief areas based on geographically weighted regression

Zhao Mingsong1,2, Liu Binyin1, Lu Hongliang1, Li Decheng2, Zhang Ganlin2※

(1.232001; 2.210008)

Accurate estimates of the spatial variability of soil organic matter (SOM) are necessary to properly evaluate climatic chagne, soil carbon sequestration potential and soil fertility. In plains and gently undulating terrains, soil spatial variability is not closely related to relief, and thus digital soil mapping (DSM) methods based on soil-landscape relationships often fail in these areas. Therefore, different predictors or methods are needed for DSM in plains. In provincial regional scale, climatic factors influence spatial distribution of soil properties. For this research, Jiangsu Province was selected as example and mean annual temperature (MAT), mean annual precipitation (MAP), physical clay content, and soil pH were selected for SOM spatial modeling using geographically weighted regression (GWR). The SOM content in the surface layer (0－20cm) of 1 519 typical soil profiles of the Second National Soil Survey in Jiangsu Province were collected. 1 217 samples were selected as the modeling set and 302 were the validation set. Fristly, 100% (1 217), 80% (973), 60% (730), 40% (486), and 20% (243) samples were randomly selected from the modeling set, and global and local spatial autocorrelation of SOM content were analyzed at different spatial scales using spatial statistics tools in ArcGIS. Secondly, comparison of the accuracy between GWR model and the global regression model under the different sampling size was conducted. Akaike information criterion (AIC), residual sum of squares (RSS) and adjustment determination coefficient (2adj) were used modeling comparison. Thirdly, the optimal model was selected for mapping SOM spatial prediction. Independent validation was used for model evaluation, using four indices: mean error (ME), mean absolute error (MAE) and root mean of squared error (RMSE), and determination coefficient (2). Results show that: 1) There was a significant spatial autocorrelation of SOM content in Jiangsu Province at different spatial scales. The clustering pattern of global and local spatial autocorrelation of modeling set with different sampling size were similar. The global Moran’s I ranged from 0.25 to 0.61 (<0.001). The spatial distribution of SOM content was mainly characterized by spatial clustering pattern. The “high-high” clustering areas were mainly distributed in the central and south of Jiangsu, and the “l(fā)ow-low” clustering areas were mainly distributed in the north of Jiangsu. 2) The modeling results of GWR were better than the global regression modeling, and the residuals had no spatial autocorrelation at different spatial scales. The2adjof GWR in different modeling sets was increased by 0.15 to 0.20 compared with the global model. The AIC and RSS were significantly lower than the global model, which were decreased by 56.08 to 360.19 and 17.40 to 76.67 respectively. There were slight difference between GWR models with different sampling size. 3) The number of modeling samples (except for the number of modeling samples was 243) had little effect on the accuracy of prediction and mapping results of GWR, the RMSE was between 5.56 and 5.75 g/kg, MAE was between 3.87 and 4.05 g/kg and2was between 0.48 and 0.52. The results were all better than the validation result of Ordinary Kriging using all modeling sampling points. This study can provide reference for SOM modeling and mapping in large and low relief areas with sparse samples.

soils; organic matter; models; geographically weighted regression; digital soil mapping; low relief areas; Jiangsu Province

趙明松，劉斌寅，盧宏亮，李德成，張甘霖. 基于地理加權(quán)回歸的地形平緩區(qū)土壤有機(jī)質(zhì)空間建模[J]. 農(nóng)業(yè)工程學(xué)報(bào)，2019，35(20)：102－110.doi：10.11975/j.issn.1002-6819.2019.20.013 http://www.tcsae.org

Zhao Mingsong, Liu Binyin, Lu Hongliang, Li Decheng, Zhang Ganlin. Spatial modeling of soil organic matter over low relief areas based on geographically weighted regression[J]. Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2019, 35(20): 102－110. (in Chinese with English abstract) doi：10.11975/j.issn.1002-6819.2019.20.013 http://www.tcsae.org

2019-03-26

2019-06-23

國(guó)家自然科學(xué)基金（41501226）；土壤與農(nóng)業(yè)可持續(xù)發(fā)展國(guó)家重點(diǎn)實(shí)驗(yàn)室開(kāi)發(fā)基金（Y412201431）；安徽省高校自然科學(xué)研究項(xiàng)目（KJ2015A034）

趙明松，博士，副教授，主要從事數(shù)字土壤制圖和空間變異研究。Email：zhaomingsonggis@163.com

張甘霖，博士，研究員，博士生導(dǎo)師，主要從事土壤發(fā)生分類、土壤地理研究。Email：glzhang@issas.ac.cn

10.11975/j.issn.1002-6819.2019.20.013

S159.9

A

1002-6819(2019)-20-0102-09