單 捷, 邱 琳, 田 苗, 王志明, 王晶晶, 盧必慧, 黃曉軍
(江蘇省農業科學院農業信息研究所,江蘇 南京 210014)
耕地是人類繁衍生息和生存發展的基礎性資源,也是農業生產最基本的物質條件,更是糧食安全的關鍵保障[1]。隨著中國城鎮化建設進程的加快,部分耕地轉變為建設用地、工礦用地,全國耕地總面積呈現出逐年減少的態勢。《第三次全國國土調查主要數據公報》顯示,截至2019年底,中國的耕地面積為1.279×108hm2,比2009年減少了7.523×107hm2。為有效加強耕地資源保護,促進耕地質量持續提升,國家相關部門先后出臺了《中華人民共和國國民經濟和社會發展第十三個五年規劃綱要》、《中共中央、國務院關于做好2022年全面推進鄉村振興重點工作的意見》等一系列政策文件,明確提出“堅持最嚴格的耕地保護制度”、“落實‘長牙齒’的耕地保護硬措施”、“嚴守18億畝耕地紅線”。政府對耕地資源保護的高度重視以及耕地資源保護面臨的現實困境,都需要加強對耕地資源利用與保護的科學研究,這對于合理規劃、科學管理和有效利用耕地資源具有重要意義。
圍繞耕地資源的利用與保護,近20年來國內外學者也開展了一系列研究,重點聚焦在耕地的時空變化與驅動力[2-5]、耕地景觀與生態安全[6-8]、耕地集約利用與整理[9-11]等方面。近年來,景觀生態學方法和地理空間分析方法被廣泛應用于耕地資源時空格局變化研究中,任平等[12]采用核密度估算和空間自相關等研究方法開展耕地空間分布格局及其變化特征研究。李黎等[13]運用核密度估算和景觀指數揭示了都江堰市耕地的時空演變特征。張揚等[14]結合核密度估算和空間自相關等方法研究喀斯特山區耕地分布與時空演變規律及其驅動因素。目前,已有的關于江蘇省耕地時空分布特征的研究中,所用的耕地數據多為各類統計數據[15-18],但統計數據并不能體現耕地地塊的形狀特征和空間分布特征,而利用遙感數據進行江蘇省耕地空間分布特征分析的研究并不多見。基于此,本研究擬以江蘇省耕地為研究對象,利用遙感數據提取高精度耕地地塊,綜合運用耕地指數、景觀指數、核密度估算、空間自相關等方法,揭示江蘇省耕地的空間分布特征,以期為江蘇省耕地的空間格局優化、耕地資源的合理利用以及耕地管理政策的制定提供理論基礎。
江蘇省地處中國大陸東部沿海地區,北緯30°45′~35°08′,東經116°21′~121°56′,是長江三角洲地區的重要組成部分。江蘇省土地面積為1.072×107hm2,其中耕地面積為4.099×106hm2,土地資源以平原為主,土層深厚,肥力中上,農業生產條件得天獨厚,適宜種植水稻、小麥等糧食作物,被稱為“魚米之鄉”。
本研究以2012年覆蓋江蘇全省的RapidEye衛星影像為數據源,各期影像質量完好。首先對各期影像進行預處理,然后根據野外實地調查資料確定耕地和其他地類的解譯標志,最后采用目視解譯的方法對耕地地塊進行人工勾繪,得到2012年江蘇省耕地數據。在全省范圍內隨機建立220個500 m×500 m地面樣方,采用亞米級差分全球定位系統(GPS)對地面樣方進行實地測量,以實測結果對解譯的耕地數據進行精度檢驗,解譯精度高于98%。
本研究所用的江蘇省行政區數據來源于全國地理信息資源目錄服務系統中的1∶250 000全國基礎地理數據庫。考慮到研究區耕地數據的統一性,本研究分別對江蘇省各設區市內的市轄區進行合并,最終得到56個研究單元(圖1)。
為了準確揭示研究區耕地的空間分布特征,本研究采用耕地指數(Ri)[18]分析耕地的空間分布特征,計算公式如下:
Ri=ci/Si
(1)
式中,ci表示第i個研究單元的耕地面積,Si表示第i個研究單元的總面積。

圖1 研究區研究單元示意Fig.1 Schematic diagram of research units in study area
本研究參考相關方法[19-21],并結合研究區耕地的實際特點,在斑塊類型尺度水平上,選取平均斑塊面積、邊界密度、斑塊密度、面積加權平均形狀指數和面積加權平均分維數等5個景觀指數,對江蘇省耕地地塊的形狀特征進行分析,各景觀指數及生態學意義見表1。
空間自相關常被用來分析某一變量在空間上的分布特征,該方法通過判斷變量的變化是否取決于其相鄰位置的變化,從而確定該變化是否具有空間自相關性[22]。空間自相關方法按功能分為2類:全局空間自相關和局部空間自相關。本研究借助該方法對研究區耕地的相關指數的空間分布特征進行分析。
2.3.1 全局空間自相關 莫蘭指數(Moran’sI)常被用于全局空間自相關分析,其表達式為:
(2)

表1 景觀指數及其生態學意義

通常還要將I標準化為Z,即采用Z檢驗對其結果進行統計檢驗,進一步判斷變量空間相關的正負性。Z的計算公式如下:
(3)
式中,E[I]為期望值,V[I]為方差。當Z為正值且大于1.96時,表明存在正的空間自相關,呈集聚分布;當Z值為負值且小于-1.96時,表明存在負的空間自相關,呈離散分布;當Z取值在[-1.96,1.96]時,空間自相關不明顯,呈隨機分布[27]。
2.3.2 局部空間自相關 由于全局Moran’sI僅能描述變量的整體分布狀況,判斷變量在空間是否有集聚特征,但其并不能確切指出集聚在哪些地區。因此,本研究選取局部空間自相關指數Local Moran’sI分析空間集聚區域,公式為:
(4)
(5)

核密度估算(Kernel density estimation,KDE)利用核函數計算各樣點xi在以h為半徑的圓內的各柵格單元中心點的密度貢獻值[28],估算模型為:
(6)
式中,h為搜索半徑或帶寬,n為帶寬內樣點的數量,k()為核函數,(x-xi)為估計點x到樣本點xi的距離。
2.4.1 粒度的選擇 核密度估算前,需要將矢量數據轉換為柵格數據,柵格數據粒度的大小會影響估算結果[28]。本研究參考相關研究[28-29],選擇矢量數據與柵格數據面積差最小時的粒度進行核密度估算。
2.4.2 帶寬的確定 核密度估算時,h增大,估計點的密度會變得平滑,但會掩蓋密度的結構,h減小時,估計點密度變化突兀不平[30]。相關研究[28]選用基于“Silverman經驗規則”帶寬估算公式確定搜索半徑(SR),計算公式為:
(7)
式中,SD為標準距離,Dm為中值距離,n為帶寬內樣點的數量。本研究根據研究區耕地的特點,并參考相關方法[21],經過對不同帶寬的多次試驗并與公式(7)進行對比,確定最優帶寬然后進行核密度估算。
3.1.1 耕地指數空間分布分析 計算江蘇省56個研究單元的耕地指數,并采用自然斷點法對其進行分級(圖2):低比重區(0.096 7~0.179 1)、中低比重區(0.179 2~0.364 6)、一般比重區(0.364 7~0.477 5)、中高比重區(0.477 6~0.560 8)、高比重區(0.560 9~0.660 5)。江蘇省耕地指數空間分布呈北高南低的特征,中高比重區和高比重區主要分布于耕地資源豐富的蘇北5市(徐州市、宿遷市、連云港市、淮安市和鹽城市),一般比重區分布于蘇中3市(揚州市、泰州市和南通市),低比重區和中低比重區分布于經濟發展速度較快、城市化進程較快的蘇南5市(南京市、鎮江市、常州市、無錫市和蘇州市)。全省56個研究單元中,處于低比重區的研究單元有3個,處于中低比重區的有13個,處于一般比重區的有14個,處于中高比重區的有12個,處于高比重區的有14個。可見,全省56個研究單元在除低比重區以外的其他4個等級區分布的數量都較為平均。
3.1.2 耕地指數的全局空間自相關分析 利用GeoDa 1.20計算耕地指數的全局Moran’sI,進行全局空間自相關分析。江蘇省耕地指數全局Moran’sI為0.657,P值為0.001,Z值為7.814 5,表明江蘇省耕地空間分布具有很強的空間正相關,呈現顯著的集聚狀態。
3.1.3 耕地指數的局部空間自相關分析 在GeoDa 1.20中計算耕地指數的局部Moran’sI,分析其局部空間集聚特征。高-高型表示與耕地指數高值研究單元相鄰的研究單元都為高值;低-低型表示與耕地指數低值研究單元相鄰的研究單元都為低值;低-高型表示與耕地指數低值研究單元相鄰的研究單元都為高值;高-低型表示與耕地指數高值研究單元相鄰的研究單元都為低值。
圖3顯示,耕地指數局部正相關類型中,高-高型和低-低型的研究單元各有11個;高-高型集中分布于蘇北的鹽城市和連云港市;低-低型集中分布于蘇南的蘇錫常3市;局部負相關類型中,高-低型和低-高型各有1個,分別是蘇南的丹陽市和蘇北的連云港市區;不顯著型有32個,這些地區的耕地指數在空間上呈隨機分布。可見江蘇省耕地高比重集聚區主要分布在蘇北地區,耕地低比重集聚區分布在蘇南地區。
3.2.1 耕地景觀指數的空間分布分析 首先,在Arcgis 10.2中將研究區的耕地矢量數據轉換為空間分辨率為10 m×10 m的柵格數據,再利用Fragstats 4.2計算景觀指數并分級(圖4)。由圖4可以看出,江蘇省耕地平均斑塊面積指數的分布特征是蘇北>蘇中>蘇南,指數較高的地區主要分布于蘇北5市,

圖2 江蘇省耕地指數空間分布Fig.2 Spatial distribution of cultivated land index in Jiangsu province

圖3 江蘇省耕地指數集聚區分布Fig.3 Spatial aggregation distribution of cultivated land index in Jiangsu province

圖4 景觀指數空間分布Fig.4 Spatial distribution of landscape indices in Jiangsu province
指數較低的地區則集中分布在蘇南5市。斑塊密度指數和邊界密度指數的分布特征均為蘇中>蘇北>蘇南,除南通市整體較高以外,其他各市都有高有低。面積加權平均形狀指數和面積加權平均分維數的分布特征相似,均為蘇南>蘇中>蘇北,指數較高的地區分布在蘇南的南京市和常州市以及蘇中的南通市,蘇北和蘇中的指數均較低。可見,蘇南耕地地塊形狀復雜度最大,破碎度最高,其次是蘇中,蘇北最小。
3.2.2 耕地景觀指數的全局空間自相關分析 利用GeoDa 1.20計算耕地景觀指數的全局Moran’sI指數,對其進行全局空間自相關分析(表2)。5個景觀指數的全局Moran’sI指數均大于0,邊界密度指數和面積加權平均分維數的Z值均大于1.960 0,說明這2個指數在空間分布上具有顯著的正相關,集聚特征較為顯著;平均斑塊面積指數、斑塊密度指數和面積加權平均形狀指數的Z值均小于1.960 0,說明這3個指數的空間自相關性不明顯,沒有明顯的集聚特征。

表2 全局空間自相關分析
3.2.3 耕地景觀指數的局部空間自相關分析 利用GeoDa 1.20計算各景觀指數的局部Moran’sI(圖5)。
平均斑塊面積指數的局部正相關類型有3個,均為低-低型(南通市區、常熟市和高郵市);局部負相關類型有3個,其中高-低型1個(溧陽市),低-高型2個(泗陽縣和盱眙縣);其他50個均為不顯著型。
斑塊密度指數的局部正相關類型有4個,其中高-高型1個(如東縣),低-低型3個(宜興市、無錫市區和蘇州市區);局部負相關類型有1個,為低-高型(高郵市);其他51個研究單元均為不顯著型。
邊界密度指數的局部正相關類型有9個,其中高-高型4個(集中分布在南通市),低-低型5個(集中分布在無錫市和蘇州市);局部負相關類型有1個,為低-高型(高郵市);其他46個均為不顯著型。
面積加權平均形狀指數的局部正相關類型有6個,其中高-高型3個(南京市、鎮江市區和揚中市),低-低型3個(東海縣、新沂市和沭陽縣);局部負相關類型有3個,其中高-低型1個(徐州市區),低-高型2個(常州市區和金壇市);其他47個均為不顯著型。
面積加權平均分維數的局部正相關類型有10個,其中,高-高型4個(集中分布在南京市和鎮江市),低-低型6個(集中分布在徐州市、宿遷市和連云港市);局部負相關類型有3個,其中高-低型1個(徐州市區),低-高型2個(常州市區和金壇市);其他43個均為不顯著型。
可見,無論是局部正相關類型還是局部負相關類型,邊界密度指數和面積加權平均分維數所包含的研究單元個數均高于其他3個景觀指數,這與全局Moran’sI指數的結果一致。
3.3.1 耕地核密度估算結果 本研究采用5 km×5 km的網格對研究區進行劃分,計算每個網格里耕地比重并分級,據此建立耕地面積點狀空間分布圖(圖6a),從而進行核密度估算。最終得到耕地密度變化范圍為1 km20~84.037 2點,并將耕地密度分為低密度區(1 km20~9.263 9點)、中低密度區(1 km29.264 0~27.791 8點)、中密度區(1 km227.791 9~45.658 0點)、中高密度區(1 km245.658 1~61.208 2點)和高密度區(1 km261.208 3~84.037 2點)(圖6b)。
由圖6可知,比重較高的耕地集中分布在中高密度區到高密度區,如蘇北5市和蘇中3市;比重較低的耕地分布在低密度區到中密度區,如蘇南5市。全省0.52%的耕地分布在低密度區,5.36%的耕地分布在中低密度區,19.60%的耕地分布在中密度區,37.93%的耕地分布在中高密度區,36.60%的耕地分布在高密度區。可見,江蘇省70%以上的耕地集中分布在中高密度到高密度區。

圖5 江蘇省耕地景觀指數集聚區分布Fig.5 Spatial aggregation distribution of cultivated landscape indices in Jiangsu province

a:耕地比重;b:核密度估算結果。圖6 耕地比重與核密度估算結果分布Fig.6 Spatial distribution of cultivated land proportion and kernel density estimation results
分別對各耕地密度等級區分布的研究單元個數進行統計,由表3可知,隨著耕地密度等級的提高,各密度區的分布范圍呈先增大后減小的趨勢。中密度區分布范圍最廣,分布在全部56個研究單元中,低密度區分布的研究單元個數最少,只有37個。
對各個研究單元包含的耕地密度等級個數進行統計,由圖7可知,全省56個研究單元中,25個研究單元包含5個等級(低密度區至高密度區),約占所有研究單元的44.6%,集中分布在蘇北5市和蘇中的南通市;22個研究單元包含4個等級,占比約39.3%,大部分集中分布在蘇南5市(低密度區至中高密度區),少數分布于蘇北和蘇中(中低密度區至高密度區);9個研究單元包含3個耕地密度等級,占比約16.1%,集中分布在蘇南的鎮江市(中低密度區至中高密度區)、蘇州市南部(低密度區至中密度區)。可見,蘇北的耕地等級結構復雜且耕地密度等級較高,其次是蘇中,蘇南的耕地等級結構簡單且耕地密度等級較低。

表3 不同耕地密度等級區分布的研究單元個數

a:各研究單元包含的耕地密度等級個數;b:各研究單元包含的耕地密度區間。圖7 耕地密度等級空間分布Fig.7 Spatial distribution of cultivated land density grade
3.3.2 不同耕地密度區的景觀指數變化特征分析 分別對各研究單元內不同耕地密度等級區內的景觀指數求平均值,由表4可知,隨著耕地密度等級的提高,各密度區的平均斑塊面積指數的均值呈上升趨勢,在高密度區達到峰值,其他4個景觀指數的均值都是先上升后下降,在中高密度區達到峰值。

表4 不同耕地密度區的景觀指數的平均值統計
對不同耕地密度區內各研究單元的景觀指數進行統計,分析各景觀指數在不同密度區的變化特征(圖8)。
(1)平均斑塊面積指數:34個研究單元的平均斑塊面積指數隨著耕地密度的增加而增加,集中分布在徐州市、宿遷市、淮安市、揚州市、泰州市、南京市、鎮江市和常州市等;7個研究單元的平均斑塊面積指數呈先上升再下降的趨勢,集中分布于泰州市北部及其相鄰的張家港市;7個研究單元的平均斑塊面積指數呈先上升再下降最后再上升的趨勢,主要分布于南通市中部和鹽城市北部;8個研究單元的平均斑塊面積指數呈先下降再上升的趨勢,分散分布于蘇北、蘇中和蘇南。可見,全省60%地區耕地的平均斑塊面積指數隨著耕地密度由低到高逐漸增大。

圖8 景觀指數隨耕地密度等級變化趨勢空間分布Fig.8 Spatial distribution of landscape indices changing with cultivated land density grade
(2)斑塊密度指數:14個研究單元的斑塊密度指數隨著耕地密度的增加而增加,集中分布在沿海地區的連云港市和鹽城市;41個研究單元的斑塊密度指數呈先上升再下降的趨勢,集中分布在徐州市、宿遷市、淮安市、蘇中3市和蘇南5市;僅有揚中市的斑塊密度指數下降。可見,全省73%地區耕地的斑塊密度指數隨著耕地密度由低到高逐漸增大再減小。
(3)邊界密度指數:13個研究單元的邊界密度指數隨著耕地密度的增加而增加,集中分布在鹽城市;42個研究單元的邊界密度指數呈先上升再下降的趨勢,集中分布在徐州市、宿遷市、淮安市、蘇中3市和蘇南5市;僅有揚中市的邊界密度指數下降。可見,全省75%地區耕地的邊界密度指數隨著耕地密度由低到高逐漸增大再減小。
(4)面積加權平均形狀指數:30個研究單元的面積加權平均形狀指數隨著耕地密度的增加而增加,集中分布于鹽城市、淮安市、徐州市與宿遷市的交界地區、泰州市、鎮江市、常州市和蘇州市;10個研究單元的面積加權平均形狀指數呈先上升再下降的趨勢,集中分布在南通市和連云港市北部;2個研究單元(金湖縣和高郵市)的面積加權平均形狀指數呈下降趨勢;5個研究單元的面積加權平均形狀指數先上升再下降最后再上升,分散分布于徐州市、宿遷市和鹽城市;7個研究單元的面積加權平均形狀指數先下降再上升,集中分布于宿遷市、連云港市和淮安市的交界地區;2個研究單元(徐州市區和南京市區)的面積加權平均形狀指數先下降再上升最后再下降。可見,全省54%地區耕地的面積加權平均形狀指數隨著耕地密度由低到高逐漸增大。
(5)面積加權平均分維數:30個研究單元的面積加權平均分維數隨著耕地密度的增加而增加,集中分布于鹽城市、淮安市、徐州市與宿遷市的交界地區、連云港市南部、泰州市、鎮江市、常州市和蘇州市;10個研究單元的面積加權平均分維數先上升再下降,集中分布在南通市;3個研究單元(金湖縣、高郵市和東海縣)的面積加權平均分維數呈下降趨勢;4個研究單元(邳州市、泗洪縣、響水縣和濱海縣)的面積加權平均分維數先上升再下降最后再上升;7個研究單元的面積加權平均分維數先下降再上升,集中分布于徐州市、宿遷市和淮安市交界處;2個研究單元(徐州市區和南京市區)的面積加權平均分維數先下降再上升最后再下降。可見,全省54%地區耕地的面積加權平均分維數隨著耕地密度由低到高逐漸增大。
本研究以江蘇省為研究區,利用高空間分辨率遙感影像提取耕地數據,在此基礎上采用耕地指數、景觀指數和多種空間分析方法,探討江蘇省耕地地塊的空間分布特征,結論如下:
(1)耕地指數的空間分布特征:江蘇省耕地指數空間分布呈現北高南低的特征,耕地高比重區集中分布于耕地資源豐富的蘇中和蘇北地區,耕地低比重區集中分布于經濟發展速度較快、城市化進程較快的蘇南地區。
(2)耕地指數的空間自相關特征:江蘇省耕地指數的空間分布整體上呈現很強的空間正相關,呈顯著的集聚狀態;耕地高比重集聚區分布在蘇北地區,耕地低比重集聚區分布在蘇南地區。
(3)耕地景觀指數的空間分布特征:平均斑塊面積指數總體分布特征是蘇北>蘇中>蘇南;斑塊密度指數和邊界密度指數總體分布特征為蘇中>蘇北>蘇南;面積加權平均形狀指數和面積加權平均分維數總體分布特征為蘇南>蘇中>蘇北;蘇南耕地地塊形狀復雜度最大,其次是蘇中,蘇北最小。
(4)耕地景觀指數的空間自相關特征:邊界密度指數和面積加權平均分維數存在明顯的空間自相關,在空間分布上具有顯著的正相關,集聚特征較為顯著;平均斑塊面積指數、斑塊密度指數和面積加權平均形狀指數的空間自相關不明顯,沒有明顯的集聚特征。
(5)不同耕地密度區的耕地景觀指數空間分布特征:江蘇省70%以上的耕地集中分布在中高密度到高密度區;隨著耕地密度等級的提高,各密度區的分布范圍呈先增加后減少的趨勢;全省50%以上地區耕地的平均斑塊面積指數、面積加權平均形狀指數和面積加權平均分維數隨著耕地密度的增大而升高;全省70%以上地區耕地的斑塊密度指數和邊界密度指數隨著耕地密度的增加先上升再下降。
同時,本研究還存在一定的不足,如未對形成耕地空間分布特征差異的原因進行分析,所以在后續的研究中將運用相關分析方法并結合蘇北、蘇中、蘇南的自然、社會、經濟等因素對其進行深入分析和探討。另外,在今后的研究中,將對江蘇省耕地地塊數據進行實時更新,并運用如耕地重心模型、土地利用動態度指數模型、土地利用轉移矩陣、地理加權回歸模型等地理空間分析方法,對江蘇省耕地資源的時空分布格局、演變特征及驅動力進行探索和分析,為耕地資源的可持續利用、合理規劃和管理以及耕地保護政策的制定等提供科學依據。