基于地形分區IDW 的復墾排土場土壤養分空間分布

趙艷玲，趙彬程，王 鑫，肖 武，劉慧芳

(1.中國礦業大學(北京)土地復墾與生態重建研究所，北京 100083；2.浙江大學 公共管理學院，浙江 杭州 310058；3.中國恩菲工程技術有限公司，北京 100083)

露天礦的開采過程中，表土層的剝離堆存和回填不可避免地產生土壤養分流失，導致復墾土壤貧瘠，植被恢復緩慢，使復墾區難以及時覆蓋植被，加重土壤侵蝕和退化，生態系統循環難以建立。王金滿等[1]通過對露天煤礦排土場的研究認為，土壤養分對植物生長起到促進作用；肖禮等[2]通過對露天煤礦排土場土壤與植被間關系研究發現，土壤養分對于植被恢復有著重要意義；為準確掌握土壤養分的空間分布，王黨朝等[3]利用傳統克里金插值法對勝利露天礦進行空間插值，發現排土場邊坡與平臺土壤的物理性質具有明顯的空間差異；郭凌俐等[4]根據傳統克里金插值法對排土場土壤進行研究，發現傳統克里金插值法能直觀揭示復墾區不同層次土壤顆粒的空間分布特征；王子良[5]通過對銅陵礦區土壤元素的插值，發現在最優情況下，IDW插值方法具有良好穩定的插值效果。但由于研究區不同，無法對插值方法做出明確的優劣性判斷，同時，對于地形復雜、人類活動干擾較多且變化強烈的地區(如排土場)，一些學者認為全局插值方式并不能直觀有效地對區域信息進行預測，而應將研究區域通過特征變量進行分區，并對其分別插值，此種方式能夠提高空間預測的精度。胡剛等[6]通過分區插值方法，發現不同地形使用不同的插值方法可以有效提高預測精度；U.Mishra 等[7]以印第安納州為研究區，將剖面深度作為輔助變量和傳統克里金插值結合使用，有效預測了表層土壤有機碳的含量。

排土場是由人工堆積而成，地形復雜，且邊坡和平臺的過渡較自然地形更為明顯。基于以往研究，本文提出在礦區排土場中結合地形將其分為邊坡與平臺2 類地形，采用分區反距離權重法對排土場土壤養分進行插值，同時與傳統克里金插值法和反距離權重法進行比較。通過對3 種插值方式預測精度的比較，選出最適宜的估算模型，作為排土場土壤養分含量預測的方法，為排土場土壤評價奠定基礎。

1 材料與方法

1.1 研究區概況

研究區位于內蒙古自治區錫林浩特市北郊，地理坐標為東經115°30′～116°26′，北緯 43°57′～44°14′。研究區屬半干旱草原氣候，年溫差較大，1 月平均氣溫-18.8℃，7 月平均氣溫21.2℃，年平均氣溫2.6℃。受季風影響，降雨多集中在夏季。年平均降水量294.74 mm，年平均蒸發量1 794.64 mm。由于研究區降雨少蒸發強，導致植被多為耐旱性植物且植被覆蓋率低。隨著礦區不斷開發，煤炭大量開采，礦區植被受到人為破壞，生態環境不斷惡化，因此，需進行土地復墾與植被重建。勝利煤田一號露天礦北排土場于2006 年開始治理，治理時間為3 a，治理面積107 萬m2，排土場高度為60 m，共分為4 個平臺(平臺1—平臺4)，3 個邊坡(邊坡1—邊坡3)，綠化面積為101 萬m2，其中坡面35.5 萬m2、平臺面積為65.5 萬m2。研究區位置及地形如圖1 所示。

1.2 樣品采集與測定

圖1 研究區位置及地形分區Fig.1 Study area location and terrain division

2017 年6 月，對勝利一號露天礦北排土場進行實地調查及土壤采樣，采樣點分布如圖1 所示。利用網格法布設采樣點，并考慮北排土場的地形因素，增加邊坡采樣點，共獲得117 個土壤樣品，舍棄17個不合格樣品，最終采用100 個樣品。所有樣品均來自0～10 cm 的表土層，利用環刀采集土樣后放入鋁盒，裝入自封袋，采樣同時記錄樣點編號、取樣坐標等信息。帶回實驗室后，測定土壤有機質、速效氮、有效磷、速效鉀含量等。

土壤有機質(SOM)含量通過在加熱條件下，用一定量的重鉻酸鉀-硫酸溶液氧化土壤中的有機碳，并使用油浴鍋分析獲得；速效氮含量使用半微量凱氏定氮法獲得；有效磷含量使用Olsen 法測定；速效鉀含量通過乙酸銨浸提-火焰光度法測定[8-10]。

1.3 插值方法

1.3.1 傳統克里金插值法

克里金空間插值法是基于變異函數模型在有限區域內對區域化變量的取值進行無偏最優估計的空間插值方法[11-13]。半方差函數是地統計學中用以研究土壤空間變異性的關鍵函數之一，能夠定量反映土壤性質在不同距離上觀察值之間的變化。在一維條件下，區域化變量Z(x)在x軸方向上的變異函數，被定義為區域變量Z(x)在點x和x+h處的值Z(x)與Z(x+h)差的方差的一半，表示為y(h)，計算公式如下：

根據實測數據及半方差函數，計算獲取變量的線性加權組合，以此對預測值進行無偏的最佳估計。本研究區為排土場，根據各個采樣點獲取的實測數據并利用傳統克里金插值法進行插值，在結合半方差函數y(h)的基礎上，使用若干實測點Z(xi)來推求出未知點的Z0(x0)，通過采用線性加權法表達未測點x0的預測值，具體見式(2)：

式中：y(h)為半方差函數；h為樣本間距；N(h)表示間距為矢量h的所有觀測點個數；Z(xi)、Z(xi+h)為Z(x)在xi和(xi+h)位置上的實測值；Z0(x0)為未知點處x0的線性、無偏、最優估計值；n為實測值個數；Wi為各實測點權重。

1.3.2 反距離權重插值法(IDW)

IDW 即使用若干個實測點Z(xi,yi)來推求未知點Z0(x0,y0)的一種方法，一般來說，首先要計算出預測點周圍每個實測點的權重，然后根據線性加權法計算出Z0的預測值[14]，具體見式(3)、式(4)：

式中：di為預測點周圍各實測點到預測點之間的距離。

1.3.3 分區反距離權重插值法(PIDW)

結合排土場的地形特征，將排土場分為7 個子區，利用反距離權重法在每個分區分別進行插值，獲取各個分區的空間分布圖，并整合為整個研究區的空間分布預測圖。

1.3.4 模型驗證方案

為檢驗插值方法預測土壤養分空間分布的效果，以平均誤差ME、相對誤差RE、均方根誤差RMSE、平均絕對誤差MAE 以及R2這5 個指標來比較不同插值方法的插值精度。多種誤差指標的綜合使用能更好地反映預測值誤差的實際情況[15-16]。

Shanghai-Kunming high-speed rail XU Sheng LI Jia-yu ZHANG Hui-gang et al.(35)

2 結果分析

2.1 土壤養分指標的描述性統計分析

研究區土壤養分的描述性統計結果見表1。根據變異系數(CV)等級劃分，CV＜10%屬于弱變異性，若10%≤CV≤100%則為中等變異，當CV＞100%是強變異性。表1 顯示，速效氮、有效磷、速效鉀及有機質的變異系數均大于10%卻小于100%。由此，4 種養分元素都屬于中等變異程度。

描述性統計分析僅能表示各種土壤養分含量的特征，不能反映整個研究區各土壤養分的空間分布情況及邊坡與平臺之間的空間分布差異。因此，需繪制土壤養分空間分布圖。

表1 土壤養分基本信息Table 1 Basic information of soil nutrients

2.2 不同插值方法的結果

在進行克里金插值時，要求插值數據符合正態分布。本文對不符合正態分布的4 種養分數據進行對數轉換后再進行K-S 檢驗，結果顯示，有效磷、速效氮、速效鉀、有機質4 種養分均符合正態分布。

本次選取反距離權重法與傳統克里金插值法進行插值，同時由于邊坡與平臺土壤養分含量差異較大，再選取反距離權重法進行分區插值，共3 種空間插值方法。通過ME、MAE、RMSE、RE、R2共5 個指標的驗證結果，得出最適合研究區的插值方法。表2 為3 種插值方法交叉驗證結果。

據表2 可知，4 種養分的空間分布插值結果中，插值效果最好的方法均為分區反距離權重法。對有效磷和有機質而言，反距離權重法次之，傳統克里金插值法精度最差。對速效氮而言，傳統克里金插值法次之，反距離權重法精度最差。對速效鉀而言，各種指標均為3 種方法的最小值。

根據圖2 顯示，4 種土壤養分元素采用分區反距離權重法插值的R2最高，其中速效氮與速效鉀的R2提升較為明顯，達到0.30 以上。相較來說，反距離權重法對土壤養分的R2整體上也要優于傳統克里金插值法。以土壤N 為例，分區反距離權重法得到的研究區土壤速效氮含量為8.09～63.73 mg/kg，傳統克里金方法得到的含量為11.22～55.62 mg/kg，而反距離權重法得到的含量為11.12～63.68 mg/kg；相較于反距離權重法和傳統克里金插值法，分區反距離權重法對土壤N 含量的預測范圍更加準確，對土壤N 的預測范圍與實測值相差不大，能較好地預計土壤養分的極值點(較大值和較小值)，反距離權重法次之，傳統克里金插值結果最差，其他土壤元素(P、K 和有機質)也有相同的結果。這是因為克里金插值對于預測結果有平滑效果，使其對極值點的預測效果較差，但較小值對于排土場而言屬于土壤貧瘠區，不利于植被的恢復，應當受到重視。綜合比較不同空間插值方法的交叉驗證結果顯示，針對本研究區，分區反距離權重法插值的精度相對較高。

表2 不同插值法交叉驗證結果Table 2 Cross-validation results for different interpolation methods

圖2 不同插值法驗證點散點圖Fig.2 Scatter plots of verification points for different interpolation methods

3 討論

圖3 北排土場土壤養分空間分布Fig.3 Spatial distribution of soil nutrients in the north dump

據圖3 可知，排土場的速效氮含量整體呈現出由北至南逐漸減少的規律；而有效磷含量均偏低，僅中心地區含量有所提升；整體來看排土場的速效鉀含量偏低，但東北地區含量處于水平值之上；土壤有機質整體含量偏高，中心地區含量相對偏低。平臺與邊坡的養分分布空間差異較大，但這種空間差異并沒有表現出規律性[17]。這是由于排土場的土壤是由人為堆積而成，與自然土壤的分布規律不同，表現出明顯的空間異質性。文獻[18]研究發現邊坡區域對水流的保持作用比平臺差，且水土流失、土壤受侵蝕狀況與坡度成正比，使邊坡土壤養分低于平臺。但本研究并無此類結論，這可能是由于研究區排土場復墾年限較短，且降水較少，多以風蝕為主，致使邊坡平臺之間養分分布并無規律。

分區反距離權重插值方法基于排土場的地形在常規反距離插值方法上進行了改進，但是在每個分區上的反距離插值方法仍需遵循“地理學第一定律”[19]，否則插值效果不會得到改進，因此，分區反距離權重法仍然具有一定的局限性。交叉驗證結果顯示，速效氮和速效鉀的分區反距離權重插值結果明顯好于反距離權重法和傳統克里金插值法，但有效磷的插值結果顯示3 種插值方法沒有明顯差異，此結果表明，有效磷在分區后，插值點仍然是空間不平穩的，因此，不適應本文采取的所有插值方法。

本次基于地形分區，使插值點在相同地形條件下進行插值，屬于局部插值方法，這一插值方法同時也被應用于類似的研究中。基于局部插值思想，徐占軍等[11]根據積水狀態對采煤沉陷區土壤有機碳含量進行分區克里金插值并與傳統克里金插值進行比較，得出分區克里金模擬更能反映土壤有機碳的空間遞變特點，有利于分析不同因素對土壤有機碳空間分布的影響。吳子豪等[20]利用3 種結合土地利用的克里金插值模型對農田土壤有機碳含量進行預測，發現中位數中心化克里金法結合土地利用類型可生成精度更高的土壤有機碳密度空間分布圖。

4 結論

a.排土場的速效氮含量整體呈現出由北至南逐漸減少的規律；有效磷含量均偏低，僅中心地區含量有所提升；排土場的速效鉀含量偏低，東北部含量偏高；土壤有機質整體含量偏高，中心地區含量相對偏低。

b.根據分區反距離權重法得到排土場土壤速效氮含量為 8.09～63.73 mg/kg，速效鉀含量為37.19～488.26 mg/kg，有效磷含量為0.89～6.76 mg/kg，有機質含量為2.08～21.3 g/kg，與實測數據更接近，分區反距離權重法對養分含量的預測范圍較為精確，對極值預測精度較高。

c.平臺與邊坡的土壤養分空間差異較大，表現出明顯的空間異質性，但并無明顯規律。由于排土場土壤由人為堆積而成，與自然土壤分布規律不同，且研究區排土場復墾年限較短，降雨較少，土壤侵蝕多以風蝕為主，致使邊坡平臺間養分分布并無規律。

d.從預測精度上來看，北排土場最優插值方法為分區反距離權重法。這是由于研究區地形復雜，坡度起伏較大且土壤分布受人為影響導致土壤中養分含量分布復雜，結合基于地形的分區反距離權重法插值，可以有效消除不同地形采樣點土壤養分含量對空間插值預測精度的影響，從而提高土壤養分含量的預測精度。

請聽作者語音介紹創新技術成果等信息，歡迎與作者進行交流

OSID 碼