基于地統(tǒng)計(jì)學(xué)和GIS技術(shù)的北京市大興區(qū)禮賢鎮(zhèn)土壤養(yǎng)分空間變異性研究

安永龍，杜子圖，黃 勇

(1.北京市地質(zhì)勘察技術(shù)院，北京 102218；2.中國地質(zhì)調(diào)查局 發(fā)展研究中心，北京 100037)

0 引 言

土壤是由母質(zhì)經(jīng)過長時(shí)間的成土作用而形成的三維自然體，成土過程中由于受到物理、化學(xué)、生物等外界因素的不斷影響，加之人類活動對其造成的不確定性干擾，使得土壤在時(shí)間和空間上具有高度的異質(zhì)性[1-2]。單憑傳統(tǒng)的統(tǒng)計(jì)學(xué)方法無法解決空間數(shù)據(jù)與屬性數(shù)據(jù)關(guān)聯(lián)的問題，目前許多國內(nèi)外學(xué)者將地統(tǒng)計(jì)學(xué)與GIS技術(shù)相結(jié)合研究土壤中微量元素的空間變異性取得了豐富的成果。如，武婕等[3]研究了南四湖區(qū)農(nóng)田土壤有機(jī)質(zhì)和微量元素的空間分布特征及其影響因素；朱菊蘭等[4]分析了渾河太子河流域表層土壤養(yǎng)分空間分布以及地形因子、土地利用對土壤養(yǎng)分空間變異的影響；Kay等[5]通過不同土地利用方式及地形因素，對稻田土壤養(yǎng)分的空間分布進(jìn)行了預(yù)測；Roger等[6]以部分環(huán)境因子作為輔助數(shù)據(jù)，對土壤養(yǎng)分指數(shù)進(jìn)行回歸克里格插值，分析出線性回歸和克里格相結(jié)合可以提高土壤養(yǎng)分預(yù)測精度；Cordova等[7]分析了英國洛桑實(shí)驗(yàn)區(qū)土壤中氮元素的空間變異性。

雖然目前已經(jīng)有不少關(guān)于北京地區(qū)土壤養(yǎng)分指標(biāo)空間變異性的研究，但是這些研究的采樣密度較稀，趙汝東[8]在研究北京市平原地區(qū)耕地的土壤養(yǎng)分空間變異性時(shí)，共采集299件土壤樣品；董士偉等[9]在揭示北京市大興區(qū)土壤養(yǎng)分空間分布特征時(shí)，共采集1 240件土壤樣品；由于土壤的空間變異性較強(qiáng)，采樣精度的高低和不同插值方法的選擇會直接影響研究成果的準(zhǔn)確性。因此本文采用地統(tǒng)計(jì)學(xué)和GIS方法，對大興區(qū)禮賢鎮(zhèn)進(jìn)行高精度高密度采樣，共采集土壤樣品1 789件，對其中全氮、全磷、全鉀、全硼、堿解氮等10項(xiàng)養(yǎng)分指標(biāo)的空間變異性、豐缺格局及其變化特征進(jìn)行了詳細(xì)的論證，旨在為當(dāng)?shù)赝寥乐叙B(yǎng)分管理、平衡施肥和提高土地利用率提供重要依據(jù)。如此對鄉(xiāng)鎮(zhèn)級進(jìn)行高密度采樣來研究其空間變異性在全國可謂首例，同一地形單元內(nèi)的克里格插值精度隨采樣密度增加而提高[10]，對于提高鄉(xiāng)鎮(zhèn)級養(yǎng)分指標(biāo)空間克里格插值精度具有先行性作用。

1 材料與方法

1.1 研究區(qū)概況

禮賢鎮(zhèn)地處北京市大興區(qū)最南部，北鄰龐各莊鎮(zhèn)、魏善莊鎮(zhèn)、安定鎮(zhèn)，西靠榆垡鎮(zhèn)，東接河北省廊坊市，位于永定河洪積沖積平原，地勢平坦。研究區(qū)地處中緯度，屬暖帶半濕潤大陸季風(fēng)氣候，氣候特點(diǎn)是四季分明，春旱多風(fēng)，夏季炎熱多雨，冬季寒冷干燥。年平均氣溫11.6 ℃，全年無霜期平均210天左右，年平均降水量556.4 mm。研究區(qū)土壤類型主要有砂質(zhì)潮土、壤質(zhì)潮土、硫酸鹽鹽潮土、草甸沼澤土、風(fēng)沙土等。沉積物主要由永定河沖洪積而成，巖性主要為卵石、砂、黏性土及粉土。研究區(qū)總面積為92.06 km2，土地利用類型多樣，其中基本農(nóng)田保護(hù)區(qū)約占總面積的73%，建設(shè)用地區(qū)約占14%，一般農(nóng)地區(qū)約占10%，其他土地利用地區(qū)約占3%，全區(qū)農(nóng)業(yè)用地類型占比較大。

1.2 樣品采集與分析方法

本次取樣采用GPS定位，充分考慮土壤類型及土地利用方式，按照64個(gè)樣/km2的采樣密度進(jìn)行采集，共采集1 789件表層土壤樣品，如圖1所示。土壤采樣深度為0～20 cm，要求定位點(diǎn)四周輻射20 m，采集3～5個(gè)分樣點(diǎn)，等組分混合為一個(gè)組合樣，并基本保證了每個(gè)采樣單元內(nèi)采集的土樣類型的一致性。采樣工作中選擇在農(nóng)用大田、菜地、林帶地、山坡下側(cè)土層較厚處等地采樣，盡量避開水土流失嚴(yán)重或表土被破壞、明顯點(diǎn)狀污染、新近搬運(yùn)的堆積土、垃圾土、田埂等地段，以保證樣品具有代表性。

樣品由北京一零一生態(tài)地質(zhì)檢測有限公司參照國家環(huán)境質(zhì)量相關(guān)標(biāo)準(zhǔn)進(jìn)行測定。元素全量測定方法：全氮使用氧化燃燒-氣相色譜法測定，全鉀、全磷使用電感耦合等離子體光譜法測定，全硼使用發(fā)射光譜法(ES)測定，有機(jī)質(zhì)使用管式爐燃燒紅外吸收法測定，pH值使用去二氧化碳水將樣品水溶，劇烈攪拌后靜止半小時(shí)，pH計(jì)直接測定。

圖1 研究區(qū)土壤樣點(diǎn)分布Fig.1 Distribution of soil sampling sites in the study area

元素有效態(tài)測定方法：使用ICP-OES法測定有效硼，稱風(fēng)干土樣(通過2 mm篩)20.0 g于250 mL錐形瓶中，加40 mL蒸餾水，連接冷凝管，文火煮沸5 min，加4滴0.5 mol/L CaCl2溶液，離心分離得清液；使用ICP-OES法測定有效磷，稱風(fēng)干土樣(通過2 mm篩)2.5 g于浸提瓶中，加50 mL 0.5 mol/L碳酸氫鈉浸提劑，在25 ℃±1 ℃溫度中振蕩30 min后過濾；使用ICP-OES法測定速效鉀，稱風(fēng)干土樣(通過2 mm篩)5.0 g于浸提瓶中，加50 mL乙酸銨溶液，加塞振蕩30 min后過濾；稱風(fēng)干土樣(通過2 mm篩)2.0 g平鋪于螺紋絲扣密封類型的塑料擴(kuò)散皿外室，吸取3 mL 20 g/L硼酸-指示劑溶液(pH=4.5)于擴(kuò)散皿內(nèi)室，加10.0 mL 1.8 mol/L氫氧化鈉溶液于擴(kuò)散皿外室，輕搖晃動后迅速旋緊蓋子，于40 ℃保溫箱內(nèi)恒溫24 h，0.01 mol/L HCl標(biāo)準(zhǔn)溶液滴定內(nèi)室硼酸中的氨，計(jì)算堿解氮含量。

為了保證數(shù)據(jù)的準(zhǔn)確性，選取GSS1、GSS2、GSS3和GSS4共4個(gè)國家一級標(biāo)準(zhǔn)物質(zhì)以及GBW07413a、GBW07414a和GBW07415a共3個(gè)有效態(tài)標(biāo)準(zhǔn)物質(zhì)，用選定分析方法對每一個(gè)標(biāo)準(zhǔn)物質(zhì)進(jìn)行8次分析檢驗(yàn)，并參照《多目標(biāo)區(qū)域地球化學(xué)調(diào)查規(guī)范(1∶250 000)》(DZ/T 0258—2014)和《生態(tài)地球化學(xué)評價(jià)樣品分析技術(shù)要求(試行)》(DD2005—03)中關(guān)于合格率標(biāo)準(zhǔn)，分別計(jì)算被測項(xiàng)目每個(gè)樣品平均值與標(biāo)準(zhǔn)值之間的相對誤差(RE%)和相對標(biāo)準(zhǔn)偏差(RSD%)，測試結(jié)果中相對誤差的合格率全部達(dá)到100%，相對標(biāo)準(zhǔn)偏差的合格率全部達(dá)到90%以上，標(biāo)樣回收率均在93%～107%之間，全部滿足分析測試質(zhì)量的相關(guān)要求。

1.3 數(shù)據(jù)處理

正態(tài)分布性檢驗(yàn)、相關(guān)系數(shù)計(jì)算、因子分析計(jì)算均由SPSS19.0[11]完成，因子分析采用最大四次方值法，最大收斂性迭代次數(shù)為25次。只有當(dāng)數(shù)據(jù)全部服從正態(tài)分布時(shí)，才能進(jìn)行因子分析，克里格插值方法才有效，否則可能存在比例效應(yīng)[12]。判斷數(shù)據(jù)是否服從正態(tài)分布方法較多，如探索性檢驗(yàn)中柱狀圖或Q-Q圖[13]、單樣本柯爾莫哥洛夫—斯米諾夫(K-S)正態(tài)分布檢驗(yàn)[14-15]等，本文采用偏度峰度檢驗(yàn)法，結(jié)果顯示有效磷、速效鉀、堿解氮的偏度值和峰度值遠(yuǎn)大于1，而其余各指標(biāo)值基本在1附近，符合正態(tài)分布的要求，隨后對有效磷、速效鉀、堿解氮值進(jìn)行對數(shù)變換，滿足地統(tǒng)計(jì)學(xué)分析的假設(shè)條件，進(jìn)行半變異函數(shù)的計(jì)算。

半變異函數(shù)的計(jì)算和理論模型的擬合均采用地統(tǒng)計(jì)學(xué)軟件GS+9.0[16]，變異函數(shù)作為研究空間變異的重要函數(shù)，是分析土壤養(yǎng)分指標(biāo)空間變異隨機(jī)性或結(jié)構(gòu)性的基礎(chǔ)性工具。該函數(shù)如下：

(1)

式中：γ(h)表示變異函數(shù)；N表示以h為步長的觀測數(shù)據(jù)成對數(shù)目；[Z(xi)-Z(xi+h)]表示間隔為h的2個(gè)觀測點(diǎn)的實(shí)測值，也稱步長。

在本次研究中為了更加準(zhǔn)確地反映土壤養(yǎng)分指標(biāo)的空間變異性，根據(jù)半變異函數(shù)的決定系數(shù)R2和殘差RSS進(jìn)行擬合，構(gòu)建各自的理論模型[17]。本次研究所建立的半方差函數(shù)模型有指數(shù)(Exponential)和線性(Linear)兩種[18]。

指數(shù)模型：

(2)

線性模型：

γ(h)=C0+Ch/αh>0

(3)

式中：C0，塊金方差，即間距為0時(shí)的半方差，由隨機(jī)因素引起的變異；C，結(jié)構(gòu)方差，由系統(tǒng)因素引起的變異；C0+C為基臺值，即半方差函數(shù)隨間距遞增到一定程度后出現(xiàn)的平穩(wěn)值，表示系統(tǒng)內(nèi)總的變異；h，步長；α，變程，即半方差達(dá)到基臺值的樣本間距，表示隨機(jī)變量在空間上的自相關(guān)尺度。

克里格插值由ArcGIS 10.2中Geostatistics analysis模塊的普通克里格空間插值法完成[19]。通過變量在點(diǎn)x處影響范圍內(nèi)的n個(gè)有效觀測值Z(xi)的線性組合得到在該點(diǎn)x處的估計(jì)值Zx，如下：

(4)

式中：Zx為采樣點(diǎn)之間距離的估計(jì)值；Z(xi)為采樣點(diǎn)之間距離的觀測值；λi是賦予觀測值Z(xi)的權(quán)重，表示各個(gè)觀測值Z(xi)對Zx的貢獻(xiàn)。

2 結(jié)果分析

2.1 土壤養(yǎng)分的統(tǒng)計(jì)分析

研究區(qū)表層土壤養(yǎng)分指標(biāo)統(tǒng)計(jì)分析是以統(tǒng)計(jì)學(xué)為基礎(chǔ)，通過計(jì)算各樣點(diǎn)的最小值、最大值、中位數(shù)、平均值、標(biāo)準(zhǔn)差等統(tǒng)計(jì)特征來描述土壤營養(yǎng)元素的含量和分布特征(表1)。表層土壤中全量養(yǎng)分、pH值、有機(jī)碳樣品的平均值與中位數(shù)相差不大，加之偏度值分布在1附近，表明表層土壤中這些指標(biāo)的中心趨向分布呈標(biāo)準(zhǔn)正態(tài)分布，而有效養(yǎng)分的平均值與中位數(shù)相差較大，加之偏度值遠(yuǎn)超過1，表明有效養(yǎng)分的中心趨向分布可能被異常值影響而使其呈非標(biāo)準(zhǔn)正態(tài)分布。由于非正態(tài)分布的數(shù)據(jù)會增加變異函數(shù)的波動性，導(dǎo)致誤差增大，因此對有效磷、堿解氮和速效鉀進(jìn)行對數(shù)轉(zhuǎn)換，轉(zhuǎn)換后的數(shù)據(jù)服從對數(shù)正態(tài)分布，因而在計(jì)算有效磷、堿解氮和速效鉀的變異函數(shù)值時(shí)用轉(zhuǎn)換后的數(shù)據(jù)。研究區(qū)表層土壤整體為堿性環(huán)境，pH值為7.55～9.36。

大量研究表明變異系數(shù)(Cv)可以反映土壤元素的離散程度，進(jìn)而表示土壤元素分布的均勻程度。變異系數(shù)越大，元素的分布越不均勻，離散程度越大。本次研究參照測試數(shù)據(jù)范圍將變異程度分為4種類型：Cv<30%為均勻分布，30%≤Cv<60%為弱分異(中等起伏)，60%≤Cv<100%為較強(qiáng)分異(較大起伏)，Cv≥100%為強(qiáng)分異(很大起伏)。參照這一分類標(biāo)準(zhǔn)顯示有效磷、速效鉀變異系數(shù)分別為105.99%、100.43%，表明表層土壤中有效磷和速效鉀屬于強(qiáng)分異，以有效磷的分異程度最強(qiáng)。相關(guān)統(tǒng)計(jì)發(fā)現(xiàn)，土壤養(yǎng)分的來源、施用率及其在土壤當(dāng)中的流動性是影響土壤養(yǎng)分空間變異程度的主要因素。說明局部地區(qū)在施用磷肥和鉀肥之后，短期內(nèi)植物對土壤中的有效磷和速效鉀吸收不充分，導(dǎo)致土壤中有一部分殘留，但最終這些成分彼此會達(dá)到動態(tài)平衡狀態(tài)[20]；堿解氮的變異系數(shù)為80.29%，屬于較強(qiáng)分異，由此可見土壤中養(yǎng)分有效態(tài)形式的分異性相對較強(qiáng)；有效硼、有機(jī)質(zhì)、全氮的變異系數(shù)為32.46%～38.69%，說明三者的分異程度較弱；全磷、全硼、全鉀、pH值的變異系數(shù)為2.64%～21.28%，四者為均勻分布，說明這幾種養(yǎng)分指標(biāo)在土壤中較為穩(wěn)定，其中以pH值的分異程度最低，僅為2.64%，此項(xiàng)研究與國內(nèi)外部分學(xué)者的研究結(jié)果相似[21-22]，說明表層土壤酸堿性是一個(gè)較為穩(wěn)定的指標(biāo)，受施肥、灌溉等隨機(jī)性因素影響較小。

2.2 土壤養(yǎng)分含量之間的相關(guān)性分析

表層土壤10個(gè)養(yǎng)分指標(biāo)之間的相關(guān)分析(Spearman相關(guān)系數(shù))表明，大部分的指標(biāo)之間存在顯著或極顯著的相關(guān)關(guān)系(表2)。土壤中有機(jī)質(zhì)與全磷、全氮、有效硼之間呈極顯著的正相關(guān)，與有效磷、堿解氮之間呈顯著的正相關(guān)，與速效鉀、全鉀、全硼之間相關(guān)性較弱。速效鉀與全鉀呈顯著的正相關(guān)關(guān)系，而兩者與全硼之間的相關(guān)性較弱。土壤pH值與其余指標(biāo)間存在負(fù)相關(guān)關(guān)系，這與部分學(xué)者研究結(jié)果相似[23-25]。

2.3 土壤養(yǎng)分主成分分析

主成分分析法就是要在數(shù)據(jù)信息丟失最少的前提下，對高維度變量進(jìn)行降維處理[26-28]。對研究區(qū)表層土壤中9種養(yǎng)分指標(biāo)進(jìn)行主成分分析，前3個(gè)主成分特征參見表3，前3個(gè)成分的累計(jì)方差貢獻(xiàn)率為85.45%，集中了原始數(shù)據(jù)的大部分信息。第一個(gè)主成分的方差貢獻(xiàn)率為61.39%，占3種主成分中最大比例，因此可以作為衡量該研究區(qū)養(yǎng)分狀況的一個(gè)綜合性指標(biāo)。為了進(jìn)一步對主成分所包含信息進(jìn)行研究，對9種養(yǎng)分指標(biāo)的3種主成分上的因子載荷進(jìn)行統(tǒng)計(jì)(表4)。結(jié)果顯示，全磷、有效磷、全氮、堿解氮、有效硼、有機(jī)質(zhì)在第一主成分中普遍具有高載荷，而在第二主成分和第三主成分中均為低載荷；全鉀和速效鉀在第二主成分中普遍具有高載荷，而在第一主成分和第三主成分中均為低載荷；全硼在第三主成分中普遍具有高載荷，而在第二主成分和第三主成分中均為低載荷。

表1 土壤養(yǎng)分的描述性統(tǒng)計(jì)

注：全氮、全鉀和有機(jī)質(zhì)的單位為%，有效磷、堿解氮、速效鉀、全硼和有效硼的單位為10-6。

表2 土壤10個(gè)養(yǎng)分指標(biāo)的相關(guān)系數(shù)

注：*表示在 0.05 水平(雙側(cè))上顯著相關(guān)，**表示在0.01水平(雙側(cè))上顯著相關(guān)，n=1 789。

綜合相關(guān)性分析和主成分結(jié)果進(jìn)一步分析，可知全磷、有效磷、全氮、堿解氮、有效硼、有機(jī)質(zhì)來源相同，均以結(jié)構(gòu)性因素為主，全鉀和速效鉀以隨機(jī)性因素為主，全硼則以兩者混合為主。

表3土壤養(yǎng)分的主成分分析結(jié)果

Table3Principalcomponentanalysisresultofsoilnutrients

主成分特征值累計(jì)特征值方差貢獻(xiàn)率/%累計(jì)方差貢獻(xiàn)率/%主成分18.739.7361.3961.39主成分21.1510.8812.8274.21主成分31.0111.8911.2485.45

表4主成分分析旋轉(zhuǎn)后的成分載荷矩陣

Table4Rotatedcomponentmatrixofprincipalcomponentanalysis

養(yǎng)分指標(biāo)變量在各主成分上的因子載荷第一主成分第二主成分第三主成分全磷 0.819-0.203-0.048有效磷0.732-0.003-0.396全鉀 -0.0120.7350.179速效鉀0.4210.614-0.358全氮 0.8810.0110.185堿解氮0.6240.322-0.121全硼 0.2000.0710.759有效硼0.6920.2240.226有機(jī)質(zhì)0.776-0.1430.295

注：旋轉(zhuǎn)方法采用具有Kaiser標(biāo)準(zhǔn)化的四分旋轉(zhuǎn)法，旋轉(zhuǎn)在5次迭代后收斂。

2.4 土壤養(yǎng)分元素的趨勢分析

圖2 土壤微量養(yǎng)分元素趨勢分析圖Fig.2 Trend analysis of soil microelements

自然界中由于土壤易受到多種因素的共同影響，因此其空間分布具有較為明顯的趨勢效應(yīng)特征，本文運(yùn)用ArcGIS 10.2中地統(tǒng)計(jì)模塊的趨勢分析功能，探討了各項(xiàng)指標(biāo)的空間分布的趨勢效應(yīng)特征。通過異向性軸軸向自動搜索功能，最終確定各向異性的情況如下，選取無(沒有趨勢效應(yīng))、一階(區(qū)域化變量沿一定方向呈直線變化)、二階(區(qū)域化變量沿一定方向呈多項(xiàng)式變化)[29-30]。圖2為研究區(qū)中選取3類典型指標(biāo)的空間分布趨勢效應(yīng)分析圖，圖中x軸表示正東方向，y軸表示正北方向，z軸表示各樣點(diǎn)測定值的大?。蛔蠛笸队懊嫔蠝\綠色曲線表示東西向的全局性趨勢效應(yīng)變化，右后投影面上淺藍(lán)色曲線表示南北向的全局性趨勢效應(yīng)變化。圖2顯示，研究區(qū)土壤中全氮表現(xiàn)為在東西方向上先減小后增大的趨勢，屬于二階類型；有效磷表現(xiàn)為在東西方向和南北方向上呈直線分布，屬于一階類型，該類型還包括全磷、全鉀、速效鉀、全硼和有效硼；有機(jī)質(zhì)表現(xiàn)為無階效應(yīng)類型，該類型還包括堿解氮和pH值。

結(jié)合克里格插值法所造成的插值誤差進(jìn)行比較，綜合考慮各個(gè)參數(shù)選擇的合理性，評價(jià)標(biāo)準(zhǔn)為平均誤差(ME)的絕對值最接近0，如果平均標(biāo)準(zhǔn)誤差(ASE)大于均方根誤差(RMSE)，說明高估了預(yù)測值，反之則說明低估了預(yù)測值，因此兩者盡量相等。標(biāo)準(zhǔn)化平均誤差(MSE)的絕對值最接近0，標(biāo)準(zhǔn)化均方根誤差(RMSSE)值最接近1，如果標(biāo)準(zhǔn)化均方根誤差大于1，說明低估了預(yù)測值，反之為高估預(yù)測值[31]。由表5可知，全磷、有效磷、全鉀、速效鉀、全硼、有效硼的一階趨勢預(yù)測與無趨勢預(yù)測和二階趨勢預(yù)測相比，ME和MSE的絕對值更加接近0，ASE和RMSE最接近，RMSSE的值最接近1，因此克里格插值時(shí)應(yīng)選擇一階更加精準(zhǔn)，同樣堿解氮、有機(jī)質(zhì)和pH值應(yīng)該選擇無趨勢，而全氮?jiǎng)t應(yīng)該選擇二階。

2.5 土壤養(yǎng)分空間異質(zhì)性特征

半變異函數(shù)(semi-variogram)是地統(tǒng)計(jì)學(xué)中應(yīng)用最廣泛的空間格局描述的基本工具[32-33]，不但可以利用半變異函數(shù)對參數(shù)的空間分布進(jìn)行結(jié)構(gòu)性和變異性分析，而且利用這一結(jié)果結(jié)合克里格法可以對未知點(diǎn)進(jìn)行預(yù)測和模擬。因此應(yīng)用地統(tǒng)計(jì)學(xué)方法對研究區(qū)土壤養(yǎng)分指標(biāo)的空間異質(zhì)性特征進(jìn)行半變異分析。其中有3個(gè)重要的參數(shù)，塊金值(C0)代表在最小取樣距離內(nèi)參數(shù)的變異，屬于由測量誤差等隨機(jī)因素引起的隨機(jī)變異；變程(A0)也稱空間最大相關(guān)距離，反映了變量空間自相關(guān)范圍的大??；塊金系數(shù)是塊金值與基臺值的比值，當(dāng)塊金系數(shù)小于25%時(shí)說明變量具有較強(qiáng)的空間相關(guān)性，其空間變異主要由結(jié)構(gòu)性因素如種植制度、施肥措施等控制[34]。當(dāng)塊金系數(shù)介于25%～75%，說明變量具有中等強(qiáng)度空間相關(guān)性，其空間變異由隨機(jī)性因素與結(jié)構(gòu)性因素共同控制。若塊金系數(shù)大于75%，則說明研究變量的空間相關(guān)性較弱，其空間變異主要受隨機(jī)性因子控制[35-37]。隨機(jī)性因子能夠減弱土壤養(yǎng)分的空間相關(guān)性，而結(jié)構(gòu)性因素能夠增加土壤養(yǎng)分的空間相關(guān)性[38]。半變異函數(shù)模型的選擇是土壤養(yǎng)分空間變異性研究的關(guān)鍵，R2值越接近1.0，RSS值越小，表明擬合模型的效果越好。

表5 不同趨勢階數(shù)插值誤差比較

由表6可知，10項(xiàng)養(yǎng)分指標(biāo)的R2值為0.718～0.956之間，說明上述養(yǎng)分指標(biāo)的實(shí)驗(yàn)半變異與相應(yīng)的理論半變異模型非常吻合，擬合效果好。這10項(xiàng)指標(biāo)中有5項(xiàng)指標(biāo)的變程超過了8 000 m，最小的變程也在200 m以上，因此對其進(jìn)行克里格插值能得到可靠的結(jié)果[39]。

全磷的半方差函數(shù)屬于指數(shù)模型，變程為210 m，塊金系數(shù)為12.4%，有效磷的半方差函數(shù)屬于指數(shù)模型，變程為230 m，塊金系數(shù)為13.2%，全鉀的半方差函數(shù)屬于指數(shù)模型，變程為650 m，塊金系數(shù)為11.1%，全硼的半方差函數(shù)屬于指數(shù)模型，變程為270 m，塊金系數(shù)為11.8%，pH值的半方差函數(shù)屬于指數(shù)模型，變程為210 m，塊金系數(shù)為12.4%，這5項(xiàng)養(yǎng)分指標(biāo)都具強(qiáng)烈的空間相關(guān)性，說明研究區(qū)表層土壤中這5項(xiàng)指標(biāo)主要受到氣候、地形、成土母質(zhì)、土壤類型等結(jié)構(gòu)性因素的影響；全氮的半方差函數(shù)屬于線性模型，變程為8 409.81 m，塊金系數(shù)為74.2%，有效硼的半方差函數(shù)屬于線性模型，變程為8 409.81 m，塊金系數(shù)為71.9%，有機(jī)質(zhì)的半方差函數(shù)屬于線性模型，變程為8 409.81 m，塊金系數(shù)為71.9%，這3項(xiàng)養(yǎng)分指標(biāo)都具中等強(qiáng)度的空間相關(guān)性，說明研究區(qū)表層土壤中這3項(xiàng)指標(biāo)主要受到結(jié)構(gòu)性因素和隨機(jī)性因素的共同影響；速效鉀的半方差函數(shù)屬于線性模型，變程為8 409.81 m，塊金系數(shù)為86.4%，堿解氮的半方差函數(shù)屬于線性模型，變程為8 409.81 m，塊金系數(shù)為84.5%，這2項(xiàng)養(yǎng)分指標(biāo)都具較弱的空間相關(guān)性，說明研究區(qū)表層土壤中這2項(xiàng)指標(biāo)主要受到施肥、灌溉、耕作措施等隨機(jī)性因素的影響。

2.6 普通克里格空間插值法

GIS在數(shù)據(jù)統(tǒng)計(jì)模型的建立和制圖表達(dá)方面具有很強(qiáng)大的功能，它的插值圖可視化功能可以直觀地反映土壤養(yǎng)分指標(biāo)的空間變異特征?？臻g插值的方法較多，如反距離權(quán)重法[40]，泛克里格法[41]，規(guī)則樣條函數(shù)法[42]等，本次研究選取地統(tǒng)計(jì)學(xué)中最常用的普通克里格空間插值法。普通克里格空間插值法也稱空間局部估計(jì)或空間局部插值法，它是建立在半變異函數(shù)理論及結(jié)構(gòu)的分析基礎(chǔ)上，在有限區(qū)域內(nèi)對區(qū)域化變量未知點(diǎn)進(jìn)行線性無偏最優(yōu)估計(jì)的一種方法[43]。

如圖3所示，參照《北京市土壤養(yǎng)分分等定級標(biāo)準(zhǔn)》以及北京市平原區(qū)土壤元素背景值[44]，可知全氮的含量為0.05×10-2～0.27×10-2、全鉀的含量為1.37×10-2～2.36×10-2，整體較為缺乏；堿解氮的含量為47.9×10-6～241.89×10-6、全磷的含量為652.1×10-6～1 920×10-6、有機(jī)質(zhì)含量為0.55×10-2～2.66×10-2、有效硼含量為0.38×10-6～2.35×10-6，整體處于中等水平；有效磷的含量為6.42×10-6～92.37×10-6、速效鉀的含量為79.42×10-6～2 691.17×10-6和全硼含量為45.55×10-6～74.11×10-6，總體較為豐富。研究區(qū)土壤pH值范圍在8.19～8.96之間。從空間分布來看，10項(xiàng)養(yǎng)分指標(biāo)具有相對較明顯的區(qū)域性特征，低值區(qū)分布具有以下特點(diǎn)，全氮、堿解氮、全磷、有效硼、有機(jī)質(zhì)主要位于研究區(qū)中南部的2條北西向展布的低值帶，分別為王化莊—趙家莊—小馬坊—紫各莊一帶、大馬房—辛家安—祁各莊—楊各莊一帶，其余低值帶零星分布于研究區(qū)邊緣，有效磷則位于其中一條低值帶，即王化莊—趙家莊—小馬坊—紫各莊一帶。而這6項(xiàng)養(yǎng)分指標(biāo)的高值區(qū)都分布在研究區(qū)的西部和東北部邊緣地帶。全鉀和速效鉀低值區(qū)的分布與之不同，全鉀低值區(qū)分布在研究區(qū)大辛莊—良種場—楊各莊一帶、東段家務(wù)—于場—內(nèi)官莊一帶，速效鉀零星分布在研究區(qū)北部和南部邊緣。全硼則在李各莊和紫各莊等地含量相對較低。

表6 土壤養(yǎng)分含量變異函數(shù)理論模型及其相關(guān)參數(shù)

圖3 土壤養(yǎng)分指標(biāo)的空間分布Fig.3 Spatial distribution of soil nutrients

研究區(qū)中土壤pH值全部為堿性，堿性較高地區(qū)分別為王化莊—趙家莊—小馬坊—西段家務(wù)一帶、大馬房—辛家安—祁各莊—楊各莊一帶，因此可以看出，土壤堿性越高，土壤中全氮、堿解氮、全磷、有效硼、有機(jī)質(zhì)、有效磷含量越低。而全鉀、速效鉀、全硼分布則與pH值關(guān)系較弱，同時(shí)三者含量與有機(jī)質(zhì)之間無相關(guān)關(guān)系，因此得出，全鉀、速效鉀、全硼分布與土壤結(jié)構(gòu)性因素關(guān)系較弱，主要受人為活動影響較大。這與前文中相關(guān)分析、主成分分析所得的結(jié)果一致。

因此，研究區(qū)應(yīng)按照養(yǎng)分指標(biāo)空間分布格局以及所種植農(nóng)作物對養(yǎng)分吸收特點(diǎn)進(jìn)行科學(xué)有效的施肥，按照不同養(yǎng)分指標(biāo)進(jìn)行“局部圈定，分類施肥”，尤其要加大氮肥和鉀肥的投入量，調(diào)整具體肥料的用量比，實(shí)現(xiàn)合理施肥，這樣既可以避免肥料施用過量后造成的污染和浪費(fèi)，又可以保障土壤中養(yǎng)分比例得當(dāng)，增加肥料的利用效率，有效地提高農(nóng)作物的產(chǎn)量與質(zhì)量。

3 結(jié) 論

農(nóng)田區(qū)分布面積在北京市大興區(qū)禮賢鎮(zhèn)占有絕對優(yōu)勢，因此土壤養(yǎng)分指標(biāo)情況成為當(dāng)?shù)剞r(nóng)業(yè)規(guī)劃管理的關(guān)鍵性因素。通過對禮賢鎮(zhèn)表層土壤10項(xiàng)養(yǎng)分指標(biāo)空間變異及分布格局特征進(jìn)行研究，得出以下結(jié)論：

(1)經(jīng)典統(tǒng)計(jì)學(xué)分析中有效磷、速效鉀屬于強(qiáng)分異，變異系數(shù)分別為105.99%、100.43%，堿解氮屬于較強(qiáng)分異，變異系數(shù)為80.29%；其余指標(biāo)分異性相對較弱。綜合相關(guān)性分析和主成分結(jié)果進(jìn)一步分析，可知全磷、有效磷、全氮、堿解氮、有效硼、有機(jī)質(zhì)來源相同，均以結(jié)構(gòu)性因素為主，全鉀和速效鉀以隨機(jī)性因素為主，全硼則以兩者混合為主。

(2)通過異向性軸軸向自動搜索功能，結(jié)合不同趨勢階數(shù)插值誤差綜合對比，確定了全磷、有效磷、全鉀、速效鉀、全硼、有效硼應(yīng)該選擇一階，堿解氮、有機(jī)質(zhì)和pH應(yīng)該選擇無趨勢，而全氮?jiǎng)t應(yīng)該選擇二階。

(3)地統(tǒng)計(jì)學(xué)分析結(jié)果表明，全磷、有效磷、全鉀、全硼、pH值可以用指數(shù)模型較好地?cái)M合，速效鉀、全氮、堿解氮、有效硼、有機(jī)質(zhì)可以用線性模型較好地?cái)M合。其中全磷、有效磷、全鉀、全硼、pH值的塊金效應(yīng)都小于25%，因此具強(qiáng)烈的空間相關(guān)性，受到結(jié)構(gòu)性因素的影響；全氮、有效硼、有機(jī)質(zhì)的塊金效應(yīng)為25%～75%，因此具中等強(qiáng)度的空間相關(guān)性，主要受到結(jié)構(gòu)性因素和隨機(jī)性因素的共同影響；速效鉀、堿解氮的塊金效應(yīng)大于75%，具較弱的空間相關(guān)性，受到隨機(jī)性因素的影響。

(4)普通克里格插值圖直觀地反映了禮賢鎮(zhèn)土壤養(yǎng)分的分布特征。研究區(qū)發(fā)現(xiàn)2條全氮、堿解氮、全磷、有效硼、有機(jī)質(zhì)北西向展布的低值區(qū)，主要分布于中南部，即王化莊—趙家莊—小馬坊—紫各莊一帶、大馬房—辛家安—祁各莊—楊各莊一帶；有效磷的低值區(qū)歸屬于其中1條。全鉀和速效鉀的含量分布區(qū)域相似。同時(shí)發(fā)現(xiàn)研究區(qū)土壤堿性越高，土壤中全氮、堿解氮、全磷、有效硼、有機(jī)質(zhì)、有效磷含量越低，這都與經(jīng)典統(tǒng)計(jì)學(xué)結(jié)論相似。