寧夏沿黃生態經濟帶耕地的時空演變及驅動力

舒 瑞, 孫鴻睿, 鄒業斌, 侯 琴

(1.寧夏自然資源勘測調查院 寧夏 銀川 750002； 2.寧夏大學 土木與水利工程學院 寧夏 銀川 750021)

耕地是人類依存度最高的復合生態系統，具有極大的經濟價值、社會價值及生態價值[1]。在耕地的開發利用過程中，從基礎的食物生產到就業保障、從生態服務到社會安定維護[2]，其功能不斷外延。但受經濟建設活動的驅動，大量優質耕地非農化[2-3]，并且隨著城市化、工業化進程的加速推進，人地矛盾將日益突出[4]。因此厘清耕地的演變過程，找到耕地變化的影響因素，對端穩“中國飯碗”具有重要意義。

近年來國內外研究者主要圍繞構建適合不同區域的模型來解釋和分析典型區域耕地演變與驅動力機制間的關系。在研究尺度方面，形成了微觀(縣區級)、中觀(省市級、流域)和宏觀(全球及國家級)尺度等一系列研究成果[4-7]，但以中觀較為多見。在驅動因子選擇方面，主要集中在社會經濟、自然條件、國家政策、農業科技水平等方面，如耿藝偉等[5]采用了自然地理基底、經濟發展水平、社會生活狀況3個方面10項指標對河南省耕地演變進行了驅動力分析；周晨晴等[6]從地形影響因素、社會經濟影響因素、空間距離影響方面對廣河縣土地利用變化進行了分析。在研究方法方面，早期多以定性分析和因子分析法[7]、主成分分析法[8]、多元線性回歸分析法等[9]定量研究為主，此類方法無法分析其空間異質性[10]，而Logistic不僅能處理耕地變化過程中的空間變量，還能很好的解決多元線性回歸因變量必須為連續變量、因子分析不能解釋因變量的發生概率等局限，目前已逐步應用到土地變化預測模擬、驅動力分析等方面。

自2010年間，中共中央國務院《關于深入實施西部大開發戰略的若干意見》《國家主體功能區規劃》《國家“十二五”規劃綱要》等相關文件相繼出臺，寧夏“沿黃經濟區”上升為國家戰略以來，該區域土地利用變化頻率較快，經濟社會快速發展產生的用地需求與耕地保護的矛盾更為突出。因此，本文以寧夏回族自治區(以下簡稱“寧夏”)沿黃生態經濟帶為研究對象，分析2010—2018年耕地演變，并選取影響研究區耕地變化的自然環境和社會經濟等變量，運用二元Logistic回歸模型，揭示該區域耕地的演變規律和驅動機制，以期為政府管理部門制定合理的耕地資源保護政策與措施、優化耕地保護格局提供一定的參考借鑒依據。

1 研究區概況

自20世紀90年代初寧夏首次提出“黃河經濟”發展戰略，到2017年寧夏第十二次黨代會提出全力打造生態優先、綠色發展、產城融合、人水和諧的“沿黃生態經濟帶”以來，經過近30 a的探索與發展，寧夏沿黃生態經濟帶已成為寧夏經濟發展的核心和精華地帶，并躍升為國務院確定的18個重點開發區之一。寧夏沿黃生態經濟帶位于寧夏北部，地處東經104°17′—106°57′，北緯37°25′—39°23′之間，海拔956～3 542 m，包括銀川市、石嘴山市全域以及吳忠市利通區、青銅峽市和中衛市沙坡頭區、中寧縣共13個縣(市、區)，面積226.92 km2，占寧夏總面積的43.68%。截止2018年底，常住人口4.53×106人，占寧夏全區65.91%(城鎮人口占全區78.22%)；地區生產總值3.21×1011元，占寧夏全區86.67%；耕地面積占寧夏全區33.81%；糧食產量為2.22×109kg，占寧夏全區的56.42%；耕地質量等別為9.09等，遠高于寧夏全區平均11.07等，略高于全國平均9.96等(2015年度)。

2 研究方法與數據處理

2.1 數據來源

主要數據來源： ①本文土地利用現狀數據來源于寧夏各縣(市、區)土地利用現狀變更調查成果； ②DEM數據來源于地理空間數據云(http:∥www.gscloud.cn)，分辨率為30 m，同時利用ArcGIS空間表面分析功能提取坡度； ③農村居民點、城市中心點、道路、黃河流域數據來源于寧夏基礎地理數據庫； ④城鎮化率、農村人均可支配收入、第三產業產值數據來源于《寧夏統計年簽》(2010—2018)； ⑤人口及GDPkm格網數據來源于中國科學院資源環境數據云平臺(http:∥www.resdc.cn)，時間為2010，2015年。

耕地的空間分布及變化特征主要受區域自然環境和社會經濟因素的影響[11]。為了對寧夏沿黃生態經濟帶耕地演變驅動力進行全面的分析，同時考慮到因子的代表性、可量化性、區域差異性和資料的可獲得性等因素，選取了7個自然因素因子，5個社會經濟因素因子(詳見表1)。其中：①溝渠密度從土地利用現狀數據中提取，并以行政村為單位統計計算；②各距離因子(到最近道路、農村居民點、城市、黃河及其主要支流的距離)由ArcGIS中的近鄰分析功能計算得到；③相比自然因素，社會經濟因素相對活躍，因此采用研究時段內的變化量作為因子進行驅動力分析。

表1 寧夏沿黃生態經濟帶耕地變化及影響因子變量

2.2 研究方法

(1) 核密度分析。核密度分析是通過樣本數據來計算和估計數據聚集情況，并通過既定的距離衰減函數來度量研究要素密度的變化情況，以此來探索空間區域中的熱點分布和變化特征[12-13]。本文通過核密度分析方法實現耕地分布整體性和連續性的空間表達。核密度計算公式為：

(1)

式中：f(x)為耕地地塊分布核密度估計值(即單位面積上的耕地面積)；k(x)稱為核函數；h為帶寬，帶寬的選擇決定了生成密度圖形的光滑性和準確性；x-xi表示估計點到樣本點xi處的距離。

(2) 空間自相關性分析。空間自相關是指某空間單元與其周圍單元間就某種特征值進行空間自相關性程度計算，以分析這些空間單元在空間上的分布特征[14]。本文通過全局莫蘭指數來判斷研究區耕地在空間上是否具有集聚性及聚集程度，計算公式如下：

(2)

通過局部莫蘭指數來反映研究區內耕地聚集地區的空間分布，計算公式如下：

(3)

式中：Ii為i空間單元的局部莫蘭指數，其他變量含義與公式(2)相同。

(3) 土地利用轉移矩陣。土地利用轉移矩陣可以全面、定量地反映研究時段內各土地利用方式轉換的流向、數量和速率，其數學模型為：

(4)

式中：Sij為研究期初與期末的土地利用狀態；n為土地利用類型數量。

(4) 耕地變化動態度。耕地變化動態度大小可以反映耕地變化的劇烈程度，計算公式如下：

(5)

式中：K為耕地動態度；Ua為期初耕地面積；Ub為期末耕地面積；T為研究時段長度。

(5) 二元Logistic回歸模型。

① 模型原理。當因變量為二分類變量時，常使用二元Logistic回歸模型，該模型是基于抽樣數據得到各自變量的回歸系數，并通過這些系數來討論模型中因變量與自變量的關系。設p為事件發生概率，值域為0～1，其概率可以用Logistic函數計算，表達式如下：

(6)

Logistic函數是協變的非線性函數，為求回歸系數，首先求得事件發生與不發生的概率之比p/(1-p)，記為Odds，然后對Odds進行對數轉換，得到Logistic回歸模型的線性模式：

(7)

式中：x1,x2,…,xn為耕地變化驅動因子；β0,β1,…,βn為待求的回歸系數。

本文利用SPSS 22.0計算出回歸系數β及統計量Waldχ2，估計的顯著性水平P和發生比率OR。其中，β表示對應驅動因子單位變化導致的耕地變化情況，變化程度可以用OR來衡量，當β>0且統計顯著(95%置信區間時，p<0.05)，表示其他因子不變的情況下，OR隨對應因子的增加而增加，相反，則表示OR隨對應因子的增加而減少；Waldχ2表示各驅動因子對耕地變化的解釋程度，其值越大，代表解釋程度越高[15]。

② 抽樣過程。研究區共有100 m×100 m的格網2 659 458個，為保證樣本數量至少應是自變量20～30倍[16]，且避免數據的空間自相關，本文對耕地增加和耕地減少分別采用分層抽樣法，并使用ArcGIS的創建隨機點功能，隨機抽取耕地增加、減少、未發生變化樣本點各5 000個，保證了耕地變化與未變化大致相同的預測精度。根據得到的樣本點矢量圖層分別提取自變量的值，并利用Z-score法對樣點數據庫進行標準化處理后作為Logistic回歸模型的基礎數據。

③ 模型檢驗。模型檢驗包含自變量之間的多重共線性檢驗和模型運算結果的有效性檢驗。自變量間的共線性診斷常用的統計量有方差膨脹因子VIF或容限TOL，一般若VIF>10且TOL<0.1，自變量存在共線性[17-18]。通過檢驗的自變量可以進入Logistic回歸模型，模型運算后常用皮爾遜χ2(Pearsonχ2) 、偏差D(deviance，D)和Homsmer-Lemeshow(HL)等指標進行擬合優度檢驗。其中，皮爾遜χ2和偏差D僅在分類變量或協變量類型的數據較少時適用，當自變量數量增加，尤其是連續自變量納入模型之后，HL是廣為接受的擬合優度指標。當HL指標統計不顯著表示模型擬合好，即顯著性水平取0.05時，若p>0.05則反而說明顯著性較好。

3 結果與分析

3.1 耕地分布格局特征

3.1.1 耕地核密度分析 首先使用ArcGIS的Feature to point工具將耕地圖斑轉換成矢量點，然后利用Kernel Density工具，以耕地圖斑面積為統計分析字段，通過多次試驗，確定帶寬為5 km，分別就2010年和2018年耕地進行核密度估算〔公式(1)〕，最后按照自然斷點法將耕地核密度分為低、中低、中、中高、高密度區5類等級(見封3附圖1)，并統計各等級比例(表2)。

表2 2010-2018年寧夏沿黃生態經濟帶耕地核密度等級比例

整體來看，研究區耕地的空間分布呈現出明顯差異性，其中：低密度區面積最大，約占43%，主要分布在遠離城市的邊緣地區和城市郊區，以沙坡頭區、中寧縣、大武口區、靈武市分布最多；高密度區面積最小，僅占約7%，主要以吳忠市中心城區為中心團狀分布及賀蘭縣、平羅縣北部條狀分布；中低、中及中高密度區占比均在15%～17%，主要分布在近郊區和遠郊區之間。

從耕地密度的時空分布來看，2010，2018年研究區耕地密度分別為28.94，28.97 hm2/km2，耕地密度變化不大。各等級耕地密度中，低密度區和高密度區均略有下降，向中低、中及中高密度集中，主要是由于城市向外擴張占用了周邊密度較高地區耕地，以及土地整治、占補平衡等項目實施，使低密度地區耕地更加集中。

3.1.2 空間分布關聯性分析 將研究區劃分為1 049個5 km×5 km的格網(其中703個格網含有耕地)，作為耕地面積空間自相關分析的空間單元，并統計格網內的耕地面積，然后分別使用ArcGIS的Spatial Autocorrelation工具，Cluster and Outlier Analysis工具進行全局和局部空間自相關分析。

全局自相關分析結果顯示，2018年研究區內耕地面積的全局Moran’sI為0.643且檢驗結果高度顯著(p<0.001)，表明研究區耕地面積不是隨機分布的，具有明顯的空間集聚特征。局部自相關分析結果顯示(見封3附圖2)，除不顯著區域外，耕地面積的高—高與低—低正相關類型共有295個網格，占96.72%，負相關類型占3.28%。其中，高—高型主要分布在銀川平原，包括惠農區南部、平羅縣、賀蘭縣、永寧縣中南部、青銅峽東北部、利通區北部、靈武市西部，少量分布在寧衛平原沙坡頭區東北部，這些地區灌溉水源充足、地形平坦、土壤肥沃，是寧夏的主要糧食產地；低—低型主要分布在靈武市東部、沙坡頭區中西部，這些地區屬于生態脆弱區，耕地分布數量少、質量等別低。

3.2 耕地演變格局特征

3.2.1 土地利用轉移矩陣 對研究期初(2010年)和期末(2018年)土地利用現狀數據進行空間疊置分析，得到2010—2018年寧夏沿黃生態經濟帶土地利用類型轉移矩陣(表3)。

表3 2010-2018年寧夏沿黃生態經濟帶土地利用類型轉移矩陣

從表3中可以看出，耕地的雙向流轉是以未利用地—耕地—建設用地為主的生態轉換系統。耕地主要流出至建設用地，面積為13 686.38 hm2，占耕地流出總面積的86.78%。城市擴張、工礦企業增加是耕地流向建設用地的主要原因，其次農村道路、溝渠等農業基礎設施不斷完善，占耕地流出總面積的7.68%，流出至其他地類的比例均小于3%。相應的耕地流入來源主要是其他土地(主要為其他草地、鹽堿地)、牧草地和園地，分別占流入總面積的55.49%，14.19%，9.17%，主要是由于永久基本農田政策全面鋪開，生態移民與城鄉建設用地增減掛鉤政策相繼出臺使耕地面積減少、建設用地擴張態勢得以有效控制，同時“十二五”“十三五”期間的土地整治、中北部重大工程、高標準農田建設等項目的實施進一步增加了耕地面積。

3.2.2 耕地變化動態度分析 根據公式(5)，分別以城鎮群、城市規模和鄉鎮為分析單元，進行耕地動態度變化分析(表4)。由表4可知，研究區耕地動態度為2.52%，耕地數量穩中有增，耕地保護基本實現動態平衡。從城市規模來看，城市規模與耕地動態度成反比。大城市耕地動態變化劇烈，耕地數量大幅減少，動態度最小，為-7.74%；城市規模最小的Ⅱ型小城市，耕地數量增加幅度最大，耕地動態度為4.29%。

表4 2010-2018年寧夏沿黃生態經濟帶耕地變化動態度

從空間分布來看，耕地動態度的空間差異性較大(見封3附圖3)。其中，耕地動態度最小的鄉鎮低達-87.38%，最大的高達57.59%，有27個(占鄉鎮總數的20%)鄉鎮耕地動態度小于-11.37%，有14個鄉鎮動態度高于12.89%。同時耕地動態度變化有明顯空間聚類現象(圖3，表4)，靠近城市周邊的鄉鎮耕地動態度低，遠離城市周邊的鄉鎮動態度高；核心城鎮群耕地減少聚集明顯，動態度為-0.35%，北翼城鎮群和南翼城鎮群耕地增加聚集明顯，動態度分別為3.43%，8.27%。

3.3 耕地變化驅動力分析

自變量多重共線性診斷結果顯示，自變量的方差膨脹因子VIF在1.35～4.5之間，容限TOL在0.22～0.74之間，表明自變量間多重共線性不嚴重，都可以納入Logistic回歸模型。在耕地變化Logistic回歸模型中，坡度用4個虛擬變量分別代表坡度Ⅱ(2°～6°)，Ⅲ(6°～15°)，Ⅳ(15°～25°)，Ⅴ(>25°)，坡度Ⅰ(<2°)作為它們的參照對象。然后將其他自變量代入到Binary logistic模塊中，通過計算分別得到耕地減少和耕地增加的logistic回歸模型(表5—6)。

表5 2010-2018年寧夏沿黃生態經濟帶耕地減少的Logistic回歸模型相關系數

表6 2010-2018年寧夏沿黃生態經濟帶耕地增加的Logistic回歸模型相關系數

3.3.1 耕地減少的驅動力分析 通過計算，耕地減少驅動力模型預測準確率達75.9%，高—低指標為7.067，顯著性為0.529，統計檢驗不顯著，即模型擬合程度較好。通過Waldχ2可以得出，寧夏沿黃生態經濟帶耕地減少的驅動因子影響效應從大到小分別為：距最近城市的距離>溝渠密度>距最近道路距離>GDP變化量>城鎮化率變化量>人口密度變化量>地均第三產業產值變化量>高程>距最近農村居民點的距離>農村人均可支配收入變化量。

經濟發展是耕地減少的重要原因。從回歸結果可以看出，GDP變化量、城鎮化率變化量、地均第三產業產值變化量、農村人均可支配收入變化量每增加一個單位，耕地減少的概率分別約增加1.52，1.57，1.34倍和1.11倍。作為寧夏社會經濟最發達、產業最為集中、城鎮化率最高的區域，寧夏沿黃生態經濟帶“虹吸效應”明顯，吸引了大量產業和人口的聚集，加劇了城市承載壓力，勢必導致城市擴張。而存量建設用地數量少、挖潛難度大、成本高，難以滿足需求時，耕地資源不得不受到擠占。經統計2010—2018年，寧夏沿黃生態經濟帶新增建設用地4.01×104hm2，占寧夏整體的78%，可以從側面反映出耕地的減少情況。

優越的區位條件增加了耕地減少概率，其中距離城市和道路的距離每近一個1個單位，耕地減少的概率將分別增加3.03倍(發生比率OR=0.33)，1.69倍(OR=0.59)。城市和道路作為建設用地擴展的中心和極軸，其周邊耕地有著優越的區位條件，同時相比農業生產，變為建設用地更能迅速帶來的巨大經濟收益，這些都吸引著耕地迅速向建設用地流轉。值得注意的是，距離農村居民點近的耕地容易被占用，距離居民點每近1個單位的距離，耕地減少的概率就增加1.15倍(OR=0.87)。

3.3.2 耕地增加的驅動力分析 通過計算，耕地增加驅動力模型預測準確率達70.9%，高—低指標為10.169，顯著性為0.253，模型擬合程度較好。通過Waldχ2統計量可知，溝渠密度、高程、距最近農村居民點距離、距最近城市距離等自然資源因素是耕地增加的主要驅動因素，貢獻力達86.49%，其中溝渠密度和高程貢獻最大，分別為37.43%，24.47%。在高程低、灌溉條件好的區域，土壤、水分、養分流失量小，耕地開墾費用低，土地更加適宜耕作，因此耕地增加的概率大。距離最近居民點、城市的距離每遠一個單位，耕地增加的概率分別增加1.54，1.67倍，主要是由于距離城市和道路近的區域能夠開墾為耕地的后備資源本就寥寥無幾，只能向遠離城市和道路的邊緣地區挖掘耕地潛力。

人口密度變化、農村人均可支配收入變化等社會經濟因素對耕地增加貢獻力為13.51%。人口密度變化與耕地增加呈正比：①因為大量生態移民從寧南山區遷移到沿黃生態經濟帶的興慶區、永寧縣等地方，外來移民進入農村從事農業生產，為滿足生產與生活的需求，他們會進行土地開發，促進整個地區耕地數量的增加；②外來人口涌入城市，會促進城市發展，由于城市的運作會消耗大量的資源，其中對于農產品的需求增加會間接引起耕地數量增加；③隨著第三產業產值、農村人均可支配收入的變化值每增加一個單位，耕地增加的概率分別降低47%(OR=0.68)，56%(OR=0.64)，主要是而隨著農村人均可支配收入的提高，農民收入日益多元化，耕地收入占比不斷下降，并且農業收入低、穩定性差，因此造成耕地非農化和非糧化(據寧夏第三次國土調查初步統計顯示，研究區非糧化率高達約25%)。

4 討論與結論

本文以寧夏沿黃生態經濟帶2010—2018年耕地為研究對象，借助核密度分析、空間相關性分析等方法，摸清了研究區耕地分布特征及演變規律，并采用Logistic回歸模型，揭示了耕地變化的驅動機制。

(1) 研究區耕地的空間分布呈正相性，且具有明顯的差異性。高—高聚類區域，耕地密度以高、中高和中密度區為主，主要集中在平羅縣、惠農區南部、賀蘭縣、永寧縣、青銅峽東北部、利通區北部。低—低聚類區域，耕地密度以低、中低、中密度區為主，主要分布在靈武市東部、沙坡頭區中西部等生態相對脆弱區。

(2) 耕地動態度變化有明顯的空間聚類現象。耕地動態度與距離城市的距離成正比，與城市規模成反比，這與寧夏沿黃經濟區城市帶“一核、兩軸、兩翼”的發展規劃(圖3)密切相關：①自2010年該規劃實施以來，沿黃生態經濟區地區生產總值(以2010年可比價為準)由1.37×1011元增加至2.81×1011元，年均增速為13.25%，占全區的比重也由82.81%上升至86.67%，經濟高速發展使城市建設用地沿“兩軸”鏈式擴展，沿“一核、兩翼”多中心輻射擴展，占用了大量耕地。同時，規模大的城市“虹吸效應”愈加明顯，吸引了大量產業和人口聚集，加劇了城市承載壓力，勢必導致城市擴張占用耕地；②受占補平衡政策約束，建設占用耕地后只能向遠離中心城區的城市邊緣補充耕地，使城市邊緣耕地大量增加，耕地動態度高。

(3) 耕地增加是自然資源稟賦、區位可達性及社會經濟條件等綜合作用下的區位與效益擇優過程。一般來說，自然稟賦是耕地增加的根本制約因子，在高程低、坡度小、灌溉條件好的區域，土壤、水分、養分流失量小，耕地開墾費用低，土地更加適宜耕作，因此耕地增加的概率大。而距離城市、道路、居民點距離近，這些本該有利于耕地增加的區域，反而會因為區位可達性高、經濟收益潛力大，非耕地會優先轉變為建設用地以實現收益的快速變現，同時這些區域能夠開墾為耕地的后備資源本就寥寥無幾，只能向遠離城市和道路的邊緣地區挖掘耕地潛力。這種耕地的邊緣化、集聚發展趨勢，不僅減少了因人類劇烈活動使土地生態環境破壞、灌溉水源和土壤被污染帶來的隱形耕地數量減少，也能使大量零碎、分散的土地得到適當歸并，繼而提升耕地聚集程度，為后期加快農用土地向種田大戶、家庭農場流轉奠定基礎。

(4) 耕地減少主要受社會經濟因素影響較大，GDP，城鎮化率，地均第三產業產值，農村人均可支配收入的增加都增大了耕地減少的概率。同時城市和道路作為建設用地擴展的中心和極軸，其周邊耕地有著優越的區位條件，相比農業生產，變為建設用地更能迅速帶來的巨大經濟收益，這些都吸引著耕地迅速向建設用地流轉。

(5) 距離農村居民點近的耕地更容易被占用。①由于土地產出效益低下，《土地管理法》宣傳貫徹不到位，農民珍惜土地、節約用地意識淡薄；同時，亂占耕地違法成本低，在高額利益的驅使下，一些人會鋌而走險占用耕地搞開發建設；②隨著農村居民生活水平的提高，希望獲取更大居住面積，同時，新增人口住房剛性需求加之老宅基地缺乏有效退出機制，以及部分市縣多年未新批宅基地等一系列綜合因素作用下，使農村居民點就近無序擴張，占用了耕地。“大棚房”“違建別墅”“農村亂耕耕地建房”等專項整治都很好地印證了這一點。因此要通過“疏、堵、促”結合的方式減少耕地的非農化。“疏”指通過統一規劃、統一計劃住宅用地，保障農民的合理住宅需求；“堵”指對于違法占用耕地的，要以“零容忍”態度進行處置，同時要借助遙感影像、人工智能、互聯網+等先進技術，將耕地占用情況能夠發現在初始，解決在萌芽，嚴防增量；“促”指繼續加大耕地保護宣傳力度和農業種植補貼，提高農民保護耕地意識、調動農民種糧積極性。

本文采用的是土地利用現狀變更調查數據，能夠真實地反映實際用地情況，但也存在農用地、自然保留地內部轉化，而未及時變更的問題。同時，耕地變化是一個復雜的過程，其變化不僅受到自然條件、人口變化、經濟發展的影響，還受到土地政策制度、價值觀念等難以量化因素的影響。因此，下一步研究應注重應用最新遙感影像解釋作為補充資料，提高土地利用現狀數據時效性、準確性；同時，要運用系統思維，深入探討政策制度等難以量化的因素對耕地變化的作用機制。