基于GIS與分區(qū)Kriging的采煤沉陷區(qū)土壤有機碳含量空間預(yù)測

徐占軍，張 媛，張紹良，李樂樂，余明成

（1. 山西農(nóng)業(yè)大學(xué)資源環(huán)境學(xué)院；晉中 030801；2. 中國礦業(yè)大學(xué)環(huán)境與測繪學(xué)院；徐州 221008）

0 引 言1

2014年，中國煤炭產(chǎn)量為2.6×109t，占世界煤炭總產(chǎn)量的 46.9%[1]。煤炭開采過程中的一系列環(huán)節(jié)占用甚至損毀了大量的土地資源，導(dǎo)致地上植被被移除，土壤動物和其他動物因生境喪失而消亡，礦區(qū)的農(nóng)作物減產(chǎn)，生態(tài)環(huán)境遭到嚴重破壞[2]。而且，煤礦開采塌陷的土地大部分是平原農(nóng)業(yè)區(qū)的可耕地，即華北平原、東北平原和長江中下游平原的農(nóng)業(yè)區(qū)[3]。相關(guān)研究表明，中國年塌陷的耕地面積約為200 km2[4]。人類煤炭開采活動對農(nóng)業(yè)生態(tài)系統(tǒng)危害嚴重[5]，對農(nóng)田土壤有機碳碳庫擾動十分劇烈[6]。

由于農(nóng)田的土壤有機碳庫是減少陸地生態(tài)系統(tǒng)碳排放的最大潛在因素[7-8]，中國以及世界上的其他煤炭開采大國必須在區(qū)域尺度上就煤炭開采對區(qū)域農(nóng)田碳儲量的影響進行定量的研究，才可以更好地對煤炭開采區(qū)的土壤有機碳庫進行科學(xué)管理，同時也實現(xiàn)煤炭低碳開采。而適宜于煤炭開采沉陷區(qū)的區(qū)域土壤有機碳含量空間預(yù)測模型是必須要解決的科學(xué)問題。但是由于礦區(qū)受到人類影響較大，區(qū)域內(nèi)部的土壤都不同程度地出現(xiàn)了地表破壞、土壤侵蝕、地面沉陷、植被毀損等問題[9-13]，無論哪一方面的變化都會對土壤碳庫產(chǎn)生影響，使礦區(qū)土壤有機碳庫通常存在比較強的空間變異性[14-19]，因此，有必要結(jié)合煤炭開采活動引起的擾動影響因素對礦區(qū)土壤有機碳庫進行估算，以更好地了解煤炭開采活動帶來的區(qū)域尺度上土壤有機碳庫含量的變化。

當(dāng)前關(guān)于區(qū)域土壤含量預(yù)測主要有 2種思路，一種是基于遙感數(shù)據(jù)對區(qū)域土壤有機碳庫預(yù)測，另一種是運用地統(tǒng)計學(xué)方法進行預(yù)測。但目前，基于遙感的方法存在著一個主要的缺陷，即尚未建立航天遙感與土壤有機碳庫之間的定量關(guān)系，因此地統(tǒng)計學(xué)思路就成為預(yù)測區(qū)域土壤有機碳庫的主要方法。

Kriging空間插值法是應(yīng)用最廣泛的土壤屬性空間預(yù)測方法[20-23]。但是普通Kriging比較適用于地質(zhì)及人文環(huán)境等各方面條件相對比較均質(zhì)的區(qū)域，在空間區(qū)域較小但地形復(fù)雜，土壤屬性受人類活動擾動、變化特別強烈的地區(qū)，其應(yīng)用效果并不理想。

隨著計算機技術(shù)的不斷進步，結(jié)合輔助信息的Kriging方法已經(jīng)得到了廣泛的應(yīng)用，該方法可以提高預(yù)測精度。已有研究表明，許多變量可以作為Kriging插值法的輔助變量，提高土壤有機碳庫空間預(yù)測的精度[24]。例如Chai等[25]以北京平谷區(qū)為研究區(qū)，分別利用基于約束性極大似然模型的最佳線性無偏估計法和回歸Kriging法，對區(qū)域土壤有機質(zhì)含量進行了空間預(yù)測，并比較了二者的預(yù)測精度；Mishra 等[20]以印第安納州為研究區(qū)，分別使用剖面深度分布函數(shù)和普通Kriging模型，預(yù)測了表層土壤有機碳的含量，并比較了 2種預(yù)測模型的預(yù)測精度。Aladamat R.等[26]建立了全球土壤有機碳建模系統(tǒng)，并結(jié)合研究區(qū)土地利用的變化情況，預(yù)測了約旦全國尺度上2000年到2030年的土壤有機碳的儲量。

本文提出在煤炭開采劇烈擾動的開采沉陷區(qū)，利用結(jié)合輔助變量的 Kriging法來預(yù)測土壤有機碳含量的思路。在高潛水位煤礦區(qū)，煤炭開采沉陷形成了不同的積水區(qū)，根據(jù)不同的積水情況人們采取了不同的土地利用措施，導(dǎo)致不同區(qū)域之間土壤有機碳含量呈現(xiàn)出很大的差異。本文以徐州九里煤炭開采沉陷區(qū)為研究區(qū)，采用分區(qū)Kriging法，以沉陷積水情況為分區(qū)輔助變量，進行了土壤有機碳含量預(yù)測；同時也利用普通Kriging法對沉陷區(qū)土壤有機碳庫含量進行了預(yù)測。通過對 2種方法的預(yù)測精度的比較，選出精度更好的礦區(qū)土壤有機碳庫儲量估算模型，作為礦區(qū)煤炭開采沉陷區(qū)內(nèi)土壤有機碳庫含量空間預(yù)測的方法。

1 材料與方法

1.1 研究區(qū)概況

九里采煤塌陷區(qū)（34°13'39''–34°26'16''N, 117°06'21''–117°12'16''E）位于江蘇省徐州市九里區(qū)，面積約為42.15 km2，是徐州市張小樓煤礦、龐莊煤礦、夾河煤礦的煤炭開采沉陷區(qū)。3礦區(qū)的井田邊界以及沉陷區(qū)域如圖1所示（紅線范圍內(nèi)）。研究區(qū)煤礦的開采方式為井工開采，選取長壁采煤的技術(shù)。這種采煤技術(shù)造成了地面塌陷，并伴隨著一系列的生態(tài)環(huán)境問題，主要包括：（1）煤炭開采導(dǎo)致農(nóng)田大面積塌陷；（2）煤炭開采活動導(dǎo)致一些位于高潛水位的農(nóng)田變成了沉陷積水水域。

圖1 研究區(qū)位置以及采樣點分布圖Fig.1 Location of study area and distribution of soil sampling sites

研究區(qū)屬于中國華東平原區(qū)，平均海拔在 35～44 m之間，地下水位較淺，屬于溫暖季風(fēng)氣候區(qū)，年平均氣溫14 ℃，年平均無霜期為200～220 d，年平均降水量為800～-930 mm。該區(qū)是古黃河沖積平原，水稻土、砂姜土為該區(qū)主要土壤類型。農(nóng)業(yè)生產(chǎn)按照冬小麥/夏玉米的方式進行輪作，屬于中國東部的糧食生產(chǎn)區(qū)。

1.2 土壤樣品采集與處理

采樣點包括預(yù)測樣點和驗證樣點 2部分，預(yù)測樣點采樣采用的方法是正方形格網(wǎng)法。根據(jù)研究區(qū)的實際面積，格網(wǎng)大小為1000 m 1000 m×，共采集了54個預(yù)測樣點。本研究只采用了每個土樣的0～20 cm土壤有機碳含量數(shù)據(jù)。在采樣過程中，用手持GPS機記錄下每一個采樣點的地理坐標(biāo)，其中沉陷濕地區(qū)內(nèi)采集8個土壤樣點，季節(jié)性積水區(qū)內(nèi)采集13個土壤樣點，沉陷未積水區(qū)內(nèi)采集33個土壤樣點。另外，在整個沉陷區(qū)內(nèi)再隨機采集18個樣點作為驗證樣點，與 2種估算模型的預(yù)測值對比，對模型所得出的結(jié)果進行精度評價。外業(yè)采樣完成后對土壤有機碳含量采用重鉻酸鉀氧化法進行測定[27]。本次采樣工作是在2010年10月份研究區(qū)域內(nèi)所有的農(nóng)作物收割完成以后進行的。

1.3 煤炭開采沉陷區(qū)土壤有機碳含量空間預(yù)測方法

對區(qū)域土壤理化性質(zhì)進行空間預(yù)測常用的方法是Kriging空間插值法。它是基于變異函數(shù)模型在有限區(qū)域內(nèi)對區(qū)域化變量的取值進行無偏最優(yōu)估計的空間插值方法。其基本原理和方法在許多文獻中均有詳細描述[23,28-29]，是基于區(qū)域化變量理論，利用原始數(shù)據(jù)和半方差函數(shù)的結(jié)構(gòu)性，計算可獲取變量的線性加權(quán)組合，從而對待估點預(yù)測值進行無偏的最佳估計，而且這種方法能夠計算估計值的可靠性程度（估計值方差）。

研究區(qū)包括多個煤礦，且處于高潛水位區(qū)。由于每個礦區(qū)對煤炭資源的開采程度不同，導(dǎo)致各個礦區(qū)內(nèi)部的地表沉陷程度不同，進而形成了不同程度的沉陷積水區(qū)（包括沉陷未積水區(qū)、季節(jié)性積水區(qū)和常年積水濕地區(qū)），人們針對礦區(qū)不同的積水情況采取了不同的土地利用方式，以及農(nóng)業(yè)耕作措施，因而土壤有機碳含量受到影響也會產(chǎn)生一定的差異[30]。在高潛水位礦區(qū)由開采沉陷導(dǎo)致的沉陷積水情況，是影響沉陷區(qū)土壤有機碳庫的主導(dǎo)因素，因此根據(jù)不同的積水情況來分區(qū)進行Kriging插值，將樣本中土壤有機碳庫的含量分為 2部分，即各個積水區(qū)的均值和殘差。其計算公式如下：

式中 Z （ Xkj）表示各個采樣點的有機碳含量， Xkj表示的是土樣所在的位置， u（ tk）表示各個積水區(qū)有機碳含量的均值， tk表示的是土樣所在沉陷區(qū)的積水情況：主要分為常年積水濕地區(qū)、季節(jié)性積水區(qū)以及未積水區(qū) 3個部分。其中，k表示不同的區(qū)域積水情況，j表示不同的土壤樣點。r(xkj)表示殘差，即各個采樣點的有機碳減去有機碳均值的差。其中，均值間的差異反映的是沉陷積水情況不同區(qū)域中有機碳含量的變異性，而殘差之間的差異反映的是沉陷積水情況相同的地區(qū)有機碳含量的變異性。在這種方法中，將殘差看作是新的區(qū)域變量來進行Kriging插值，其變異函數(shù) γr（ h）及待估點Xkj的預(yù)測公式分別為式（2）和式（3）。每一個待估點的土壤有機碳庫含量預(yù)測值 Z*（xkj）為所屬沉陷積水情況區(qū)的均值u(tk)與

所估計的殘差r*（xkj）之和（式4）。

式中h為2個樣本間的點空間分隔距離，N(h)是分隔距離為h時的樣本點對總數(shù)。

式中kjλ為權(quán)重系數(shù)。

在本次研究中，首先用SPSS來對土壤有機碳庫數(shù)據(jù)及其殘差進行統(tǒng)計分析，然后利用ArcGIS 10.0軟件中的Geostatistical Analyst 模塊，檢驗數(shù)據(jù)的質(zhì)量和分布，本次的采樣數(shù)據(jù)滿足正態(tài)分布和隨機性檢驗，用該軟件進行空間插值并生成空間插值圖。

1.4 預(yù)測結(jié)果準(zhǔn)確性檢驗

本次對于研究結(jié)果準(zhǔn)確性的檢驗是通過比較18個驗證樣點位置處的預(yù)測值與采樣實測值之間的相關(guān)系數(shù)（r）[31]，以及Isaaks和Shrivastava定義的3個驗證指標(biāo)，即平均誤差（mean error，ME）、平均絕對誤差（mean absolute error，MAE）和均方根誤差（root mean square error，RMSE）來評估不同插值方法的性能[32]。其中，相關(guān)系數(shù)（r）這一指標(biāo)是用以反映變量之間相關(guān)關(guān)系密切程度的。均方根誤差RMSE實際上是觀測值與真值偏差的平方和再與觀測次數(shù)N來求比值，然后再求平方根來計算得到的。

式（5）～（7）中，N表示驗證點的總數(shù)量，本次研究中N為18，xoi表示驗證點的實測值，xpi表示根據(jù)估算模型得到的驗證點預(yù)測值。

在實際的測量過程中，由于對樣點進行的觀測次數(shù)總是有限的，因此只能用最可信賴的值來代替真值，而對一組測量中的特大或特小誤差，標(biāo)準(zhǔn)誤差的反映非常敏感，所以，可以用標(biāo)準(zhǔn)誤差來很好地反映出預(yù)測的精確度。r值越大、RMSE值越小，說明所預(yù)測精度越高, 反之則說明預(yù)測精度越低。

2 結(jié)果與分析

2.1 土壤有機碳的統(tǒng)計分析

首先對試驗區(qū)所采集的土壤樣本進行了描述性統(tǒng)計分析（表 1）。分析顯示，區(qū)域內(nèi)采樣點的土壤有機碳含量符合正態(tài)分布規(guī)律，不同積水區(qū)內(nèi)土壤有機碳庫含量的差異較為顯著，因此能夠進行Kriging插值；同時也為下一步對土壤有機碳庫含量進行半方差函數(shù)結(jié)構(gòu)分析以及相應(yīng)的分區(qū)Kriging插值方案的確定提供了依據(jù)，并且減少了統(tǒng)計分析過程中所產(chǎn)生的不確定性。根據(jù)表 2可知，3個區(qū)域內(nèi)土壤有機碳含量從小到大依次為未積水區(qū)、季節(jié)性積水區(qū)、濕地，這說明土壤有機碳庫的含量在不同的積水區(qū)域內(nèi)有很大的不同。而且積水區(qū)內(nèi)土壤有機碳庫的變異系數(shù)為 0.07，比未積水區(qū)內(nèi)土壤有機碳庫的變異系數(shù)0.13要小，尤其在季節(jié)性積水區(qū)內(nèi)其標(biāo)準(zhǔn)差為 0.86，也較小，說明積水區(qū)內(nèi)的土壤有機碳庫含量離散程度較小，比較集中，這也從另一方面驗證了將沉陷區(qū)內(nèi)積水量的變化作為輔助變量的分區(qū)Kriging插值法的可行性。

表1 沉陷區(qū)域外業(yè)采樣樣本土壤有機碳庫含量描述性統(tǒng)計結(jié)果Table 1 Descriptive statistical results of soil organic carbon content of samples in subsidence area

2.2 地統(tǒng)計分析及Kriging插值

利用ArcGIS 10.0，對采樣點樣本土壤有機碳含量的原始數(shù)據(jù)、去除區(qū)域積水分區(qū)內(nèi)有機碳含量均值以后的殘差半方差函數(shù)及其擬合參數(shù)進行分析[33-36]，其結(jié)果見表3。

表2 沉陷區(qū)不同積水情況下土壤樣本有機碳描述性統(tǒng)計結(jié)果Table 2 Descriptive statistics results of soil organic carbon samples with different water accumulation in subsidence area

表3 土壤有機碳含量及其殘差半方差函數(shù)模型和參數(shù)Table 3 Soil organic carbon content and its residual semi variance function model and parameters

表 3顯示了樣本有機碳含量原始數(shù)據(jù)和土壤有機碳殘差各自的半方差函數(shù)擬合指數(shù)和球狀模型的參數(shù)。根據(jù)塊金值的相關(guān)定義，理論上當(dāng)采樣點的距離為 0時，半變異函數(shù)值應(yīng)為 0，但由于存在測量誤差以及空間變異，使得 2個采樣點非常接近時，它們的半變異函數(shù)值不為0，即存在塊金值，因此塊金值反映的是最小抽樣尺度以下變量的變異性及測量誤差。結(jié)合對表 3的分析，可以看出采樣點土壤有機碳含量與殘差的塊金值分別為0.91和 0.01，土壤有機碳的空間異質(zhì)性較大，殘差的空間異質(zhì)性比較小，說明該地區(qū)不同積水沉陷區(qū)域內(nèi)部由于試驗誤差和一些隨機因素引起的空間變異都很小。

C是結(jié)構(gòu)方差，表示系統(tǒng)屬性或區(qū)域變量最大空間變異，是由土壤母質(zhì)、地形、氣候等結(jié)構(gòu)性因素引起的變異。由于該沉陷區(qū)域開采前后氣候并未發(fā)生變化，土壤有機碳庫的空間變異性基本是由于煤炭開采造成的地表形變和沉陷積水等結(jié)構(gòu)因素引起的變異，該沉陷區(qū)域土壤有機碳結(jié)構(gòu)方差C是15.17，這種結(jié)構(gòu)性因素是該區(qū)域土壤有機碳庫變異的主導(dǎo)因素。

C/(C0+C)即塊金值與基臺值的比值，被稱為空間相關(guān)度，表示可度量空間自相關(guān)的變異所占的比例，該值的大小表示系統(tǒng)變量的空間相關(guān)性的程度。由表3可知2種方法得出的 C/(C0+C)分別為 0.998與 0.943，均大于75%，這說明其相關(guān)性很強，尤其是樣本的土壤有機碳含量的相關(guān)性極強，達到了0.998。

根據(jù)表 3可知，與原始數(shù)據(jù)相比，土壤有機碳殘差擬合函數(shù)的塊金值很小，基臺值降低，變程也降低，這主要是因為，殘差是土壤有機碳含量減去各個積水區(qū)內(nèi)的有機碳含量均值得到的，而不同積水區(qū)域內(nèi)部的土壤有機碳含量的變化比較小，相對比較均勻，這也證明了對沉陷區(qū)內(nèi)土壤有機碳含量進行預(yù)測時考慮不同積水沉陷區(qū)的合理性，間接地說明了以分區(qū)對沉陷區(qū)土壤碳庫儲量進行Kriging插值預(yù)測的有效性和可行性。

2.3 兩種方法預(yù)測精度的對比分析

采用研究區(qū)18個驗證點處的土壤有機碳含量實測值與模型的預(yù)測值進行對比，對分區(qū)Kriging模型的預(yù)測精度進行評價。圖2是預(yù)測值與實測值的散點分布圖，由2種方法所得到的線性回歸方程顯示，以區(qū)域積水情況為輔助變量進行的Kriging插值得到的預(yù)測值與實測值的決定系數(shù)（0.756 4）明顯高于普通的Kriging插值法（0.508 6）。

表4為2種預(yù)測方法的驗證點處土壤有機碳含量預(yù)測的準(zhǔn)確性評估指標(biāo)。ME表示預(yù)測的平均偏差，從表4中可以看出，普通Kriging（0.020 2 g/kg）和結(jié)合積水情況的Kriging（0.021 8 g/kg）的ME均接近0，說明2種方法預(yù)測的平均偏差都比較小，總體預(yù)測精度都比較高。平均絕對誤差（MEA）反映了預(yù)測值誤差的實際情況，由表4可知，普通Kriging的MAE為1.851 1 g/kg，結(jié)合積水情況的Kriging的MAE為1.287 8 g/kg，因此，結(jié)合積水情況的Kriging預(yù)測實際誤差更小。作為平均值和方差偏差的聯(lián)合度量指標(biāo)，如表4所示，直接進行Kriging插值對土壤有機碳含量的預(yù)測RMSE為0.55，明顯高于分區(qū)Kriging插值（0.35）。綜上所述，考慮區(qū)域內(nèi)土壤積水情況的Kriging插值法對土壤有機碳含量的預(yù)測值精度更高。

圖2 兩種方法得到的土壤有機碳預(yù)測值與實測值散點圖Fig.2 Scatter diagram of predicted and measured values of soil organic carbon by 2 methods

表4 兩種方法預(yù)測土壤有機碳含量的準(zhǔn)確性評估指標(biāo)Table 4 Accuracy assessment indices of estimated soil organic carbon content with 2 methods

2.4 土壤有機碳的空間分布特征

通過 2種方法預(yù)測得到的土壤碳庫儲量空間分布如圖3所示，可以發(fā)現(xiàn)沉陷區(qū)2種土壤有機碳含量空間預(yù)測方法（普通Kriging和分區(qū)Kriging）得到的土壤有機碳殘差的空間分布格局，而分區(qū)Kriging預(yù)測的土壤有機碳含量圖的空間分布特征及空間遞變規(guī)律更加明顯，普通Kriging得到的土壤有機碳含量范圍為：9.34～16.252 g/kg，而分區(qū)Kriging插值的結(jié)果為9.338～18.058 g/kg，2種方法的估計范圍大致相同，總體上看研究區(qū)土壤有機碳分布的空間格局基本一致，有機碳含量由中部向四周逐漸減少，中部的土壤有機碳含量最高，據(jù)實地調(diào)查該區(qū)域?qū)儆跐竦貐^(qū)，與其他土地利用相比，濕地土壤長期飽和度較低，其凋落物和有機質(zhì)分解速度較慢，這增加了濕地土壤中的土壤有機碳含量[37]。

圖3 兩種方法得到的土壤有機碳含量空間分布Fig.3 Spatial distribution of soil organic carbon content by 2 methods

西南部的土壤有機碳含量較少，該區(qū)域?qū)儆谖捶e水區(qū)，由于礦區(qū)煤炭開采對農(nóng)田土壤和植被擾動影響嚴重，造成植被NPP (net primary productivity) 降低，從而減少了土壤有機碳的一個重要補給來源，而且開采沉陷形成很多沉陷坡面，土壤有機碳由于土壤侵蝕的原因流失嚴重[38-39]。從整體上來說 2種空間插值方法都能夠反映出整個研究區(qū)的土壤有機碳分布情況。

但是在局部分布上，2種預(yù)測方法之間存在著一定的差異。對比分析發(fā)現(xiàn)，直接進行Kriging插值得到的土壤有機碳空間分布的圖斑比較連續(xù)規(guī)整，這明顯跟研究區(qū)的實際情況不太相符，說明平滑效應(yīng)使不同積水情況區(qū)域之間的差異性降低，其差值結(jié)果只能反映研究區(qū)土壤有機碳的大致分布情況，卻不能更為精確地反映研究區(qū)的實際情況，比如圖3a局部S1區(qū)以及S2區(qū)域由于個別預(yù)測樣點的土壤有機碳含量偏高，造成區(qū)域內(nèi)土壤有機碳含量空間插值結(jié)果異常偏高。普通Kriging插值是基于采樣點的空間分布特征以及樣點與樣點之間的空間位置關(guān)系來進行空間插值，事先未對插值樣點進行分類，所以無法消除鄰近的不同類別空間預(yù)測樣點對空間插值結(jié)果的影響。而分區(qū)Kriging由于按照影響沉陷區(qū)土壤有機碳含量的最大因素——沉陷積水情況進行分區(qū)，就可以有效消除這方面的誤差影響。而且，與分區(qū)Kriging的模擬結(jié)果相比，直接Kriging法模擬的礦區(qū)土壤有機碳分布特征規(guī)律性不明顯，忽略了土壤有機碳空間平穩(wěn)過渡。因此，分區(qū)Kriging模擬更能反應(yīng)土壤有機碳的空間遞變特點，有利于分析不同因素對土壤有機碳空間分布的影響，能更好地反映研究區(qū)域內(nèi)部土壤有機碳含量的分布情況。

分區(qū)Kriging空間插值的思想在相關(guān)領(lǐng)域也得到了學(xué)者的認可，Zhang等[31]通過對比在土壤現(xiàn)場采樣過程中分配樣點的 4種不同模式，即不分類的網(wǎng)格模式，基于土壤類型的分配模式，基于土地利用類型的分配模式以及基于土地利用模式-土壤類型的分配模式下得到的土壤有機碳的變異系數(shù)，得出評估土壤有機碳空間分布的最有效方法是基于土地利用模式-土壤類型分配模式。同樣地，Sandeep等[40]在研究美國賓夕法尼亞州的土壤有機碳時，采用結(jié)合環(huán)境變量的地理加權(quán)回歸克里金法，與原來的回歸克里金法進行對比發(fā)現(xiàn)，因為前者考慮了空間非平穩(wěn)性以及殘差的空間自相關(guān)性，其精度會提高。所以在受人類煤炭開采活動擾動劇烈的沉陷區(qū)土壤有機碳含量的預(yù)測應(yīng)用中，其預(yù)測精度也高于回歸克里金法。

3 結(jié) 論

本文提出了利用分區(qū) Kriging法來對煤炭開采沉陷區(qū)土壤有機碳含量進行空間插值預(yù)測的方法，并將該方法與傳統(tǒng)的直接Kriging方法得到的預(yù)測結(jié)果進行對比，發(fā)現(xiàn)：

1）結(jié)合區(qū)域內(nèi)部積水情況來進行的分區(qū) Kriging法得到的研究區(qū)土壤有機碳含量的范圍為：9.338～18.058 g/kg，而直接Kriging方法得到的范圍為：9.34～16.252 g/kg。

2）從區(qū)域預(yù)測精度上分析，與直接進行 Kriging插值相比，分區(qū)Kriging預(yù)測方法精度更高。研究區(qū)內(nèi)由于煤炭開采造成了不同的沉陷積水區(qū)，土地利用方式也各不相同，導(dǎo)致區(qū)域內(nèi)土壤有機碳含量出現(xiàn)差異。結(jié)合區(qū)域積水情況的分區(qū)Kriging，可以消除鄰近不同積水情況區(qū)的采樣點土壤有機碳含量對其空間插值預(yù)測精度的影響，提高土壤有機碳含量的預(yù)測精度。

綜上所述，可以選擇結(jié)合區(qū)域積水情況的分區(qū)Kriging空間插值作為煤炭開采沉陷區(qū)土壤有機碳含量空間預(yù)測的預(yù)測方法。在煤礦區(qū)，優(yōu)化的土壤有機碳空間預(yù)測方法，也有利于更精確地進行土壤有機碳動態(tài)模擬和理解本地土壤有機碳庫的時空演化，為礦區(qū)土地低碳復(fù)墾乃至區(qū)域內(nèi)的土地資源的低碳利用提供更加科學(xué)的依據(jù)。

[1] BP. The BP Statistical Review of World Energy 2014[EB/OL]. http://www.bp.com/statisticalreview, 2018-03-12.

[2] 原野, 趙中秋, 白中科, 等. 露天煤礦復(fù)墾生物多樣性恢復(fù)技術(shù)體系與方法：以平朔礦排土場為例[J]. 中國礦業(yè),2017, 26(8): 93－98.Yuan Ye, Zhao Zhongqiu, Bai Zhongke, et al. Technology system and method of biodiversity restoration for the reclamation of opencast coal mine:a case study from the dumps in Pingshuo mine[J]. China Mining Magazine, 2017,26(8): 93－98. (in Chinese with English abstract)

[3] 羅愛武. 淮北市采煤塌陷區(qū)土地復(fù)墾研究[J]. 安徽師范大學(xué)學(xué)報（自然科學(xué)版）, 2002, 25(3): 286－289.Luo Aiwu. The land restoration in the sunk areas of coal extraction in HuaiBei City[J]. Journal of AnHui Normal University (Nature Science), 2002, 25(3): 286－289. (in Chinese with English abstract)

[4] 彭蘇萍, 趙建慶. 中國西部煤礦區(qū)生態(tài)環(huán)境控制及改善研究[C]//中國科協(xié)2000年學(xué)術(shù)年會. 2000.

[5] 侯湖平, 徐占軍, 張紹良, 等. 煤炭開采對區(qū)域農(nóng)田植被碳庫儲量的影響評價[J]. 農(nóng)業(yè)工程學(xué)報, 2014, 30(5): 1－9.Hou Huping, Xu Zhanjun, Zhang Shaoliang, et al. Effect evaluation on vegetation carbon pool of region agro-ecosystem by coal mining in mining area[J].Transactions of the Chinese Society of Agricultural Engineering (Transactions of the CSAE), 2014, 30(5): 1－9.(in Chinese with English abstract)

[6] 王坤. 高潛水位采煤塌陷區(qū)充填復(fù)墾土壤碳動態(tài)研究[D].徐州：中國礦業(yè)大學(xué), 2016.Wang kun. Study on Carbon Dynamics of Reclaimed Soil in Coal Mining Subsidence Area with High Groundwater Level[D]. Xuzhou: China University of Mining and Technology, 2016.

[7] Miller B A, Koszinski S, Wehrhan M, et al. Comparison of spatial association approaches for landscape mapping of soil organic carbon stocks[J]. Soil, 2015(1): 217－233.

[8] Zhang G S, Ni Z W. Winter tillage impacts on soil organic carbon, aggregation and CO2emission in a rainfed vegetable cropping system of the mid–Yunnan plateau, China[J]. Soil& Tillage Research, 2017, 165: 294－301.

[9] Lan C Y, Shu W S, Wong M H. Revegetation of lead/zinc mine tailings at shaoguan, guandong province, China:Phytotoxicity of the tailings[J]. Studies in Environmental Science, 1997, 66(97): 119－130.

[10] Yan Chu. The influence of coal mining on environmental quality of mining area and the research progress of repairing technology[J]. Advances in Environmental Protection, 2016,6(1): 1－6.

[11] Anthony B. Restoration of mined land-using natural processes[J]. Ecological Engineering, 1997, 8(4): 255－269.

[12] Zhang Ling, Wang Jinman, Bai Zhongke, et al. Effects of vegetation on runoff and soil erosion on reclaimed land in an opencast coal-mine dump in a loess area[J]. Catena, 2015,128(5): 44－53.

[13] Huang Yi, Tian Feng, Wang Yunjia, et al. Effect of coal mining on vegetation disturbance and associated carbon loss[J]. Environmental Earth Sciences, 2013, 73(5): 1－14.

[14] álvaro-Fuentes J, Easter M, Paustian K. Climate change effects on organic carbon storage in agricultural soils of northeastern Spain[J]. Agriculture Ecosystems &Environment, 2012, 155(155): 87－94.

[15] Yang X M, Drury C F, Wander M M, et al. Evaluating the effect of tillage on carbon sequestration using the minimum detectable difference concept[J]. Pedosphere, 2008, 18(4):421－430.

[16] Mishra U, Ussiri David A N, Lal R. Tillage effects on soil organic carbon storage and dynamics in Corn Belt of Ohio USA[J]. Soil & Tillage Research, 2010, 107(2): 88－96.

[17] Rosemary F, Vitharana U W A, Indraratne S P, et al.Exploring the spatial variability of soil properties in an Alfisol soil catena[J]. Catena, 2017, 150(3): 53－61.

[18] Teng Mingjun, Zeng Lixiong, Xiao Wenfa, et al. Spatial variability of soil organic carbon in Three Gorges Reservoir area, China[J]. Science of the Total Environment, 2017,599(6): 1308－1316.

[19] Wu Lizhi, Li Long, Yao Yunfeng, et al. Spatial distribution of soil organic carbon and its influencing factors at different soil depths in a semiarid region of China[J]. Environmental Earth Sciences, 2017, 76 (19): 654.

[20] Mishra U, Lal R, Slater B, et al. Predicting soil organic carbon stock using profile depth distribution functions and ordinary kriging[J]. Soil Science Society of America Journal,2009, 73(2): 614－621.

[21] Liu Yaolin, Guo Long, Jiang Qinghu, et al. Comparing geospatial techniques to predict SOC stocks[J]. Soil &Tillage Research, 2015, 148(3): 46－58.

[22] 杜挺, 楊聯(lián)安, 張泉, 等. 縣域土壤養(yǎng)分協(xié)同克里格和普通克里格空間插值預(yù)測比較——以陜西省藍田縣為例[J].陜西師范大學(xué)學(xué)報(自然科學(xué)版), 2013, 41(4): 85－89.Du Ting, Yang Lian'an, Zhang Quan, et al. Spatial prediction comparison of soil nutrient between ordinary Kriging and Cokring at county scale-A case study in Lantian county of Shaanxi province[J]. Journal of Shaanxi Normal University(Natural Science Edition), 2013, 41(4): 85－89. (in Chinese with English abstract)

[23] 施周, 閆杭召, 畢晨, 等. 基于地統(tǒng)計學(xué)——克里格插值法的村鎮(zhèn)地表水體水質(zhì)監(jiān)測[J]. 環(huán)境工程學(xué)報, 2017(4):2607－2613.Shi Zhou, Yan Hangzhao, Bi Chen, et al. Surface water quality monitoring in town based on Geotatistics-Kriging interpolation[J]. Chinese Journal of Environmental Engineering, 2017(4): 2607－2613. (in Chinese with English abstract)

[24] Dai Fuqiang, Zhou Qigang, Lv Zhiqiang, et al. Spatial prediction of soil organic matter content integrating artificial neural network and ordinary kriging in Tibetan Plateau[J].Ecological Indicators, 2014, 45(5): 184－194.

[25] Chai Xurong, Shen Chongyang, Yuan Xiaoyong, et al.Spatial prediction of soil organic matter in the presence of different external trends with REML-EBLUP[J]. Geoderma,2008, 148(2): 159－166.

[26] Aladamat R, Rawajfih Z, Easter M, et al. Predicted soil organic carbon stocks and changes in Jordan between 2000 and 2030 made using the GEF SOC Modelling System[J]. Agriculture Ecosystems & Environment, 2007, 122(1): 35－45.

[27] Fang Yu, Yan Zhilei, Chen Jichen, et al. Effect of chemical fertilization and green manure on the abundance and community structure of ammonia oxidizers in a paddy soil[J]. Chilean Journal of Agricultural Research, 2015, 75(4): 487－496.

[28] Wang L X, Okin G S, Caylor K K, et al. Spatial heterogeneity and sources of soil carbon in southern African savannas[J]. Geoderma, 2009, 149(3-4): 402－408.

[29] Sun Zhili, Wang Jian, Li Rui, et al. LIF: A new Kriging based learning function and its application to structural reliability analysis[J]. Reliability Engineering & System Safety, 2017, 157(4): 152－165.

[30] Shrestha R K, Lal R. Changes in physical and chemical properties of soil after surface mining and reclamation[J].Geoderma, 2011, 161(3–4): 168－176.

[31] Zhang Z Q, Yu D S, Shi X Z, et al. Effect of sampling classification patterns on SOC variability in the red soil region, China[J]. Soil & Tillage Research, 2010, 110(1): 2－7.

[32] Isaaks E H, Srivastava R M. An Introduction to Applied Geostatistics[M]. Oxford University Press, 1989, 33(33): 483－485.

[33] Lévesque J, King D J. Airborne digital camera image semivariance for evaluation of forest structural damage at an acid mine site[J]. Remote Sensing of Environment, 1999,68(2): 112－124.

[34] Li Xinrong. Influence of variation of soil spatial heterogeneity on vegetation restoration[J]. China Science:Earth Science, 2005, 48(11): 2020－2031.

[35] 郭凌俐, 王金滿, 白中科, 等. 黃土區(qū)露天煤礦排土場復(fù)墾初期土壤顆粒組成空間變異分析[J]. 中國礦業(yè), 2015,24(2): 52－59.Guo Lingli, Wang Jinman, Bai Zhongke, et al. Analysis of spatial variability of soil granules in early stage of reclamation at opencast coal minedump in loess area[J].China Mining Magazine, 2015, 24(2): 52－59. (in Chinese with English abstract)

[36] 代富強, 周啟剛, 劉剛才. 基于回歸克里格和遙感的紫色土區(qū)土壤有機質(zhì)含量空間預(yù)測[J]. 土壤通報, 2014, 45(3):562－567.Dai Fuqiang, Zhou Qigang, Liu Gangcai. Spatial prediction of soil organic matter contents in a purplish soil region with regression Kriging and remote sensing[J]. Chinese Journal of Soil Science, 2014, 45(3): 562－567. (in Chinese with English abstract)

[37] Qualls R G, Richardson C J. Decomposition of Litter and Peat in the Everglades: The Influence of P Concentrations[M].Springer New York, 2008.

[38] Hu Zhenqi, Xiao Wu. Optimization of concurrent mining and reclamation plans for single coal seam: A case study in northern Anhui, China[J]. Environmental Earth Sciences,2013, 68(5): 1247－1254.

[39] Indorante S J, Jansen I J, Boast C W. Surface mining and reclamation: Initial changes in soil character[J]. Journal of Soil & Water Conservation, 1981, 36(6): 347－351.

[40] Sandeep K, Rattan L, Liu Desheng. A geographically weighted regression kriging approach for mapping soil organic carbon stock[J]. Geoderma, 2012, 189–190 (6): 627－634.