李 衡 陳國銀 王世晶 張 震 張志偉 靳煥煥

(東北農(nóng)業(yè)大學(xué)水利與土木工程學(xué)院, 哈爾濱 150030)

0 引言

近年來,隨著經(jīng)濟的快速發(fā)展和城鎮(zhèn)化進程的加快,土地利用/覆被變化(Land use and cover change,LUCC)成為了全球氣候與環(huán)境變化研究的熱點問題。其中,2020年中國城鎮(zhèn)化率達到63.89%,預(yù)計在2030年將達到70%[1]。城鎮(zhèn)化率的提高和土地利用的變化將大量的耕地、草地等轉(zhuǎn)變?yōu)榻ㄔO(shè)用地,這一過程將會產(chǎn)生社會、資源和環(huán)境問題,為此,減少此類問題的發(fā)生已顯得刻不容緩。當(dāng)前,國內(nèi)外學(xué)者針對土地利用/覆被變化的研究主要集中在LUCC驅(qū)動機制和演變過程[2]、土地利用變化的生態(tài)風(fēng)險及優(yōu)化[3]、LUCC的可持續(xù)性和生態(tài)效益評估以及LUCC模型與未來變化模擬[4-5]、土地利用變化帶來的碳排放時空差異和作用機理[6]等方面。然而LUCC研究的一個重要方向是未來土地利用變化的模擬,其中比較成熟的模型包括系統(tǒng)動力學(xué)模型[7]、Markov模型[8]、CA及改進模型[9-12]、智能體模型(ABM)[13]、FLUS模型[14]、PLUS模型、CLUE模型[15]以及在其基礎(chǔ)上改進的CLUE-S(Conversion of land use and its effects at small regional extent)模型[9,16-18]等。而CLUE-S模型適用于小范圍內(nèi)土地利用空間數(shù)量和空間位置變化的模擬與預(yù)測,且能夠?qū)ν恋乩米兓c其社會、經(jīng)濟和自然等驅(qū)動因子相互關(guān)系進行定量分析,應(yīng)用較廣泛[19]。

撓力河流域作為三江平原主要的糧食產(chǎn)區(qū)和國家重要商品糧基地,地理位置優(yōu)越,在受到自然環(huán)境和經(jīng)濟社會雙重影響下土地利用類型發(fā)生變化,勢必會造成土地利用適宜與不適宜,因此,研究該地區(qū)歷史年份地類時空變化、模擬未來年份地類變化格局和地類適宜性至關(guān)重要。

本文以撓力河流域作為研究區(qū),利用撓力河流域2000、2010、2020年3期土地利用數(shù)據(jù),從土地利用動態(tài)度指數(shù)和土地利用轉(zhuǎn)移矩陣兩方面分析歷史時期土地利用時空變化,利用二元邏輯回歸模型選取高程、坡度和坡向等15個驅(qū)動因子,運用CLUE-S模型對撓力河流域2020年土地利用時空變化進行模擬驗證,在此基礎(chǔ)上設(shè)置基線情景、農(nóng)業(yè)發(fā)展情景和生態(tài)保護情景預(yù)測2030年研究區(qū)土地利用空間變化并對各情景下耕地適宜性進行評價,從而為加強未來耕地質(zhì)量建設(shè)與管理、落實耕地保護目標(biāo)和科學(xué)管理土地資源提供可持續(xù)發(fā)展的科學(xué)決策。

1 研究區(qū)概況

撓力河流域位于黑龍江省三江平原腹地,地理范圍為45°42′～47°31′N、131°10′～134°09′E,行政區(qū)劃包括富錦市、饒河縣、七臺河市、雙鴨山市、集賢縣、寶清縣、友誼縣以及建三江和紅興隆農(nóng)墾分局部分地區(qū),流域總面積為2.27×104km2,占三江平原總面積的1/4[20], 如圖1所示。該流域地形主要以平原和丘陵為主,平原主要分布在流域的北部和中部的內(nèi)、外七星河及撓力河中游地區(qū),丘陵主要分布在流域的西南部和南部,地勢上呈西南高、東北低的特點[21]。撓力河流域建國以來經(jīng)歷了4次大規(guī)模的土地利用開發(fā)[22],耕地面積持續(xù)上升,呈現(xiàn)出“三山一水六分田”的土地利用特點,耕地類型主要有旱地(玉米-小麥輪作)和水田(單季稻種植為主),水澆地面積極少,是三江平原主要的糧食產(chǎn)區(qū)和國家重要商品糧生產(chǎn)基地[20-23]。

圖1 研究區(qū)域示意圖Fig.1 Schematic of study area

2 數(shù)據(jù)來源和處理

土地利用數(shù)據(jù)來源于國家基礎(chǔ)地理信息中心,本文選取2000、2010、2020年3期數(shù)據(jù),空間分辨率為30 m(2000、2010年數(shù)據(jù)總體精度為83.50%,Kappa系數(shù)為0.78;2020年數(shù)據(jù)總體精度為85.72%,Kappa系數(shù)為0.82)。該數(shù)據(jù)包括耕地、森林、草地、灌木林、濕地、水體、苔原、人造地表、裸地、冰川和永久積雪10個土地利用類型,本文選取6個地類進行研究,數(shù)據(jù)說明見表1。數(shù)字高程模型數(shù)據(jù)(DEM)來源于美國國家航空航天局(NASA)全球的全新分辨率30 m DEM數(shù)據(jù)(NASA DEM),具有質(zhì)量好、覆蓋范圍廣等優(yōu)點。人口密度、國內(nèi)生產(chǎn)總值(GDP)、中國水系矢量數(shù)據(jù)和中國道路數(shù)據(jù)均來自中國科學(xué)院資源環(huán)境科學(xué)與數(shù)據(jù)中心(http:∥www.resdc.cn)。土壤pH值數(shù)據(jù)和土壤剖面深度數(shù)據(jù)來源于國家青藏高原數(shù)據(jù)中心的中國土壤有機質(zhì)數(shù)據(jù)集。

表1 地類類型說明Tab.1 Description of earth type

利用ArcGIS 10.6軟件將柵格數(shù)據(jù)(土地利用、DEM、人口密度、GDP、土壤剖面深度和土壤pH值)經(jīng)投影、重采樣、重分類、掩膜提取得到撓力河流域100 m×100 m的柵格圖層,土壤剖面深度和土壤pH值要進行4級量化處理;將矢量數(shù)據(jù)(道路、水系)經(jīng)投影、裁剪得到撓力河流域矢量圖層。坡度、坡向柵格圖層是由數(shù)字高程數(shù)據(jù)經(jīng)空間分析—水文分析工具提取。各種道路距離是根據(jù)歐氏距離工具提取。根據(jù)CLUE-S模型需求,制備的驅(qū)動因子可以分為自然因素因子(DEM,坡度,坡向,與河流距離,與水系距離,與鄉(xiāng)道、縣道、省道距離,與城市一、二、三、四級道路距離,與鐵路距離)和社會經(jīng)濟因素因子(GDP、人口密度)。本文所用的地理數(shù)據(jù)投影統(tǒng)一采用WGS_1984_Albers,重采樣后像元大小為 100 m×100 m。

3 研究方法

3.1 土地利用動態(tài)度

選擇單一土地利用動態(tài)度和綜合土地利用動態(tài)度分析研究區(qū)土地利用變化的速率。

3.1.1單一土地利用動態(tài)度

單一土地利用動態(tài)度表征研究區(qū)內(nèi)某一地類演變的基本特征,值越大,表明該地類在研究期內(nèi)變化越活躍,值越小,表明在研究期內(nèi)變化越穩(wěn)定,計算公式為

(1)

式中K——研究期單一土地利用動態(tài)度

Sa——研究期初期單一土地利用類型面積

Sb——研究期末期單一土地利用類型面積

T——研究時段

當(dāng)T為年時,K為研究區(qū)內(nèi)土地利用類型年變化率。

3.1.2綜合土地利用動態(tài)度

綜合土地利用動態(tài)度表征研究區(qū)內(nèi)土地利用類型動態(tài)度演變的綜合特征,計算公式為

(2)

式中L——綜合土地利用動態(tài)度

LUi——研究初期第i類地類面積

ΔLUi-j——研究開始至結(jié)束時段內(nèi)第i類地類轉(zhuǎn)為其他地類的面積絕對值

n——地類總數(shù)量

3.2 土地利用轉(zhuǎn)移矩陣

土地利用轉(zhuǎn)移矩陣是Markov模型在土地利用/覆被變化研究中的應(yīng)用。利用Markov模型可以定量描述土地利用面積增加來源和減少去向[24]。轉(zhuǎn)移矩陣數(shù)學(xué)表達式為[25]

(3)

式中Aij——地類面積

i、j——研究初期和末期的地類類型

3.3 土地利用變化模擬與預(yù)測

CLUE-S模型是通過自然因素和社會經(jīng)濟兩類驅(qū)動因子相結(jié)合模擬土地利用變化及環(huán)境效應(yīng)的模型。模型輸入數(shù)據(jù)包括:根據(jù)政策紅線設(shè)置限制區(qū)域、轉(zhuǎn)移規(guī)則、土地利用類型面積需求量和空間分析。

3.3.1驅(qū)動因子的選取依據(jù)及處理

選取驅(qū)動因子時首先結(jié)合研究區(qū)實際情況,保證空間數(shù)據(jù)的可獲取性以及數(shù)據(jù)像元大小、空間范圍和空間坐標(biāo)系的時空一致性,其次是選取的數(shù)據(jù)可以實現(xiàn)空間定量化表達,最后要全面考慮社會、經(jīng)濟和自然等因素,保證模擬精度的科學(xué)準(zhǔn)確性。

結(jié)合研究區(qū)實際情況和驅(qū)動因子選取依據(jù),本文選取15個驅(qū)動因子(表2)。

表2 驅(qū)動因子處理Tab.2 Driving factor processing

3.3.2限制區(qū)域

由于研究區(qū)內(nèi)有國家級自然保護區(qū)且都劃分為人造地表用地,因此本文的限制區(qū)域為人造地表用地,即該地類不能轉(zhuǎn)移為其他地類。

3.3.3轉(zhuǎn)移規(guī)則

土地利用轉(zhuǎn)移規(guī)則有轉(zhuǎn)移矩陣和轉(zhuǎn)移彈性系數(shù)(ELAS)。轉(zhuǎn)移矩陣用1表示可以發(fā)生轉(zhuǎn)移,0表示不發(fā)生轉(zhuǎn)移,根據(jù)前文土地利用空間時空特征及研究區(qū)實際情況設(shè)定的土地利用類型轉(zhuǎn)移矩陣如表3所示。轉(zhuǎn)移彈性系數(shù)(ELAS)在0～1之間取值,值越大,表示該地類越穩(wěn)定,不易轉(zhuǎn)為其他地類。由于CLUE-S模型對參數(shù)非常敏感,設(shè)置時需要謹(jǐn)慎調(diào)節(jié),本文根據(jù)研究區(qū)的實際情況、前人的相關(guān)研究以及在模型輸入過程中反復(fù)調(diào)試,提高模擬精度,得出適合研究區(qū)各地類3種情景下(基線情景、農(nóng)業(yè)發(fā)展情景和生態(tài)保護情景)ELAS值(表4)。

表3 土地利用類型轉(zhuǎn)移矩陣Tab.3 Transfer matrix of land use type

表4 土地利用轉(zhuǎn)移彈性系數(shù)Tab.4 Coefficient of elasticity of land use transfer (ELAS)

3.3.4土地利用類型面積需求量和空間分析

假設(shè)土地利用變化幅度為線性變化,用2010—2020年現(xiàn)有土地利用數(shù)據(jù)借助線性內(nèi)插,并通過趨勢外推法得到2030年各地類面積需求量?？臻g分析中運用二元邏輯回歸分析方法分析各地類空間分布與驅(qū)動因子之間的關(guān)系,比較受試者工作特征曲線(Receiver operating characteristic, ROC),選取最優(yōu)尺度,當(dāng)ROC大于0.7時,說明回歸方程對各地類具有較好的解釋力[26],計算式為

(4)

式中pi——每個柵格可能出現(xiàn)某種地類的概率

X1,i～Xn,i——與某種地類相關(guān)的驅(qū)動因子

β0～βn——驅(qū)動因子的回歸系數(shù)

3.3.5CLUE-S模型空間模擬

模型的空間模擬是在空間分析模塊中完成,空間分析模塊是整個CLUE-S模型的核心內(nèi)容,其分析由總的分配概率決定,機制表述為

Ti=Pi,u+Iu+Eu

(5)

式中Ti——每個柵格適合于地類的總概率

Pi,u——通過二元邏輯回歸方程計算的空間分布概率

Iu——地類的迭代變量

Eu——根據(jù)轉(zhuǎn)換規(guī)則設(shè)置的參數(shù)

3.3.6CLUE-S模型精度檢驗方法

將撓力河流域2020年土地利用變化模擬圖和土地利用現(xiàn)狀圖進行對比分析驗證模擬精度。本文利用地類模擬準(zhǔn)確率和Kappa系數(shù)檢驗?zāi)M結(jié)果。

(1)地類模擬準(zhǔn)確率

基于ArcGIS 10.6中的Combine功能,將2020年土地利用變化模擬圖和土地利用現(xiàn)狀圖進行比較,將相同的柵格數(shù)除以2020年土地利用現(xiàn)狀圖柵格數(shù),即為各地類模擬準(zhǔn)確率。

(2)Kappa系數(shù)

Kappa系數(shù)是檢驗?zāi)M精度的常用手段之一,當(dāng)Kappa系數(shù)大于0.75時,模擬精度高;Kappa系數(shù)為0.4～0.75時,模擬精度一般;Kappa系數(shù)為0～0.4時,模擬精度差。

3.4 土地適宜性評價

土地適宜性評價是對某塊土地對于某種用途是否適宜以及適宜程度,做出等級評價。本文基于ArcGIS 10.6軟件,對上文3種情景下2030年耕地的適宜性進行評價,評價結(jié)果能夠為未來充分利用自然資源,開發(fā)土地潛力,實現(xiàn)作物高產(chǎn),優(yōu)化作物種植布局和為國土空間規(guī)劃提供決策支持。

通過查閱相關(guān)歷史耕地適宜性評價文獻以及研究區(qū)的實際發(fā)展情況,根據(jù)影響研究區(qū)耕地適宜性的道路因素、地形因素、土壤因素和水利因素等以及選取評價因子時要遵循差異性原則、定性定量結(jié)合原則、可操作性原則和綜合性原則[27],選取8個評價因子用于2030年耕地適宜性評價。研究區(qū)的地形因素選取DEM、坡度和坡向作為評價因子,交通便捷性因素選取與主要道路距離和與耕地距離作為評價因子,環(huán)境適宜性因素選取與水體距離作為評價因子,土壤因素選取土壤pH值、土壤剖面深度作為評價因子,建立評價指標(biāo)和評價模型,根據(jù)FAO的《土地評價綱要》將研究區(qū)劃分為高度適宜區(qū)、中度適宜區(qū)、勉強適宜區(qū)、暫時適宜區(qū)4個等級進行評價[27-29]。

3.4.1評價因子權(quán)重

確定權(quán)重的方法為層次分析法(AHP),數(shù)值越大權(quán)重越高。分析結(jié)果見表5。

表5 AHP層次分析結(jié)果Tab.5 AHP analytic hierarchy process %

3.4.2熵技術(shù)修正權(quán)重

利用熵技術(shù)對AHP確定的權(quán)重進行修正,可以使指標(biāo)體系更加優(yōu)化。修正步驟詳見文獻[30]。修正后權(quán)重見表5。

3.4.3評價模型建立

在建立評價模型之前,首先需要將8個評價因子進行量化處理,然后將各評價等級賦予一定的分值,高度適宜區(qū)、中度適宜區(qū)、勉強適宜區(qū)和暫時不適宜區(qū)分值分別為1、2、3、4。

根據(jù)設(shè)定的評價因子,將各評價因子的權(quán)重以及評價等級進行建模,模型為

(6)

式中Yi——耕地適宜性評價綜合值

m——評價因子個數(shù)

ai——第i個評價因子權(quán)重

bi——第i個評價因子等級賦分值

4 結(jié)果與分析

4.1 土地利用動態(tài)度變化

根據(jù)式(1)、(2),計算研究區(qū)研究時段內(nèi)各地類的K值和L值,見表6。從K值可以看出,2000—2010年耕地、草地和人造地表呈現(xiàn)負(fù)動態(tài),表明這3種地類轉(zhuǎn)化成其他地類,其余地類呈現(xiàn)為正動態(tài);這一時段K值最大的地類是水體,最小的是耕地。

2010—2020年耕地、草地、水體動態(tài)度呈現(xiàn)負(fù)動態(tài),森林、濕地和人造地表地類呈現(xiàn)正動態(tài)。這一時期K值最大的地類是人造地表,最小的是森林。2000—2020年K值呈現(xiàn)正動態(tài)的地類有森林、濕地、水體和人造地表,耕地和草地呈現(xiàn)負(fù)動態(tài)。這一時期K值最大的地類是水體,最小的是耕地。整個研究期間K值最大的地類是水體,達到9.84%,最小的地類是森林,為0.01%。

從L值可以看出,整個研究區(qū)20年間差別較為顯著,綜合動態(tài)度為0.15%,第1個10年時期綜合動態(tài)度最高,為0.24%,表明這10年間土地利用變化較為劇烈,主要是受到大量的未利用土地開荒,旱改水等政策的影響;2010年后用地變化逐漸趨向緩慢,綜合動態(tài)度最低的時段是2010—2020年,為0.07%?？傮w而言研究區(qū)土地利用綜合動態(tài)度呈現(xiàn)逐漸下降趨勢。

4.2 土地利用轉(zhuǎn)移矩陣變化

基于ArcGIS 10.6,利用空間分析功能計算撓力河流域2000—2010年和2010—2020年地類轉(zhuǎn)移矩陣,有助于了解研究區(qū)地類空間格局變化情況,也能直觀反映各地類之間的轉(zhuǎn)換過程和轉(zhuǎn)移方向(表7、8)。

表7 2000—2010年撓力河流域土地轉(zhuǎn)移矩陣Tab.7 Naoli River basin land transfer matrix in 2000—2010 hm2

從表7可以看出,研究區(qū)2000—2010年地類變化為:耕地、草地和人造地表面積增加,森林、濕地和水體面積減少。耕地面積凈增加3 872.93 hm2,草地面積凈增加23 675.28 hm2,人造地表面積凈增加755.42 hm2,森林面積凈減少20 056.77 hm2,濕地面積凈減少344.86 hm2,水體面積凈減少7 902 hm2。耕地的主要轉(zhuǎn)出方向是森林、草地和濕地,森林的轉(zhuǎn)出方向是耕地和草地,草地的轉(zhuǎn)出方向是耕地和森林,濕地的轉(zhuǎn)出方向是耕地和水體,水體的主要轉(zhuǎn)出方向是濕地,人造地表用地主要轉(zhuǎn)出方向是耕地和草地,這是因為隨著撓力河流域地區(qū)經(jīng)濟的不斷發(fā)展,城鄉(xiāng)建設(shè)對土地的需求進一步增多,侵占各地類的現(xiàn)象劇烈,尤其森林和水體用地面積凈減少分別達到20 056.77 hm2和7 902 hm2。

從表8可以看出,研究區(qū)2010—2020年地類變化為:耕地、草地和水體用地面積增加,森林、濕地和人造地表用地面積減少。其中,耕地面積凈增加5 002.69 hm2,草地面積凈增加2 349.15 hm2,水體面積凈增加499.38 hm2,森林面積凈減少 727.46 hm2,濕地面積凈減少1 736.8 hm2,人造地表用地面積凈減少5 377.97 hm2。耕地的主要轉(zhuǎn)出方向是森林和草地,森林的轉(zhuǎn)出方向是耕地和草地,草地的轉(zhuǎn)出方向是耕地和森林,濕地的轉(zhuǎn)出方向是耕地和水體,水體的主要轉(zhuǎn)出方向是濕地和耕地,人造地表用地主要轉(zhuǎn)出方向是耕地,這是因為在2010—2020年間隨著經(jīng)濟的發(fā)展,該地區(qū)出現(xiàn)了許多的空心房或者危舊房,國家出臺增減掛鉤政策后地塊被復(fù)墾成耕地,或者改造成草地和水塘,導(dǎo)致耕地、草地和水體面積增加,人造地表用地面積減少。

表8 2010—2020年撓力河流域土地轉(zhuǎn)移矩陣Tab.8 Naoli River basin land transfer matrix in 2010—2020 hm2

4.3 CLUE-S模型精度

由二元邏輯回歸方程檢驗結(jié)果可知,撓力河流域各地類分布概率模擬結(jié)果較好,耕地、森林、草地、濕地、水體、人造地表的ROC值分別為0.889、0.951、0.785、0.958、0.857、0.869,均大于0.7,滿足模型的回歸要求[31],所選取的影響因子可以用于土地利用變化模擬。以撓力河流域2010年土地利用數(shù)據(jù)為基礎(chǔ),將限制區(qū)域、土地利用轉(zhuǎn)移規(guī)則、土地利用需求、二元邏輯回歸方程檢驗結(jié)果等相關(guān)數(shù)據(jù)輸入CLUE-S模型中,對撓力河流域2020年土地利用變化進行模擬,如圖2所示。

圖2 撓力河流域土地利用現(xiàn)狀與CLUE-S模型預(yù)測結(jié)果Fig.2 Present situation of land use in Naoli River basin and forecast results of CLUE-S model

結(jié)果所示,模擬正確柵格數(shù)目與研究區(qū)總柵格數(shù)目的比值為0.911 8,即地類模擬準(zhǔn)確率為91.18%(表9)。從整體上看,本次研究模擬精度較高,模擬效果理想;從模擬誤差上看,土地利用類型面積大且集中連片分布的地類模擬精度高于分布相對分散的地類,同時模擬結(jié)果與土地利用類型本身的實際情況有關(guān);從本次研究來看,撓力河流域由于地類的數(shù)量和分布差異較大,且地形復(fù)雜,海拔差異大,因此低海拔平坦區(qū)域模擬精度高于高海拔地形復(fù)雜區(qū)。通過計算Kappa系數(shù)為0.894,大于模型要求的精度0.75,說明模擬效果優(yōu)且各參數(shù)值的設(shè)置比較合理,可使用CLUE-S模型和相關(guān)參數(shù)對撓力河流域土地利用變化進行模擬。

表9 各地類模擬精度Tab.9 Everywhere class simulation accuracy

4.4 土地利用變化特征及情景模擬分析

4.4.1情景模擬分析

根據(jù)《黑龍江省土地利用總體規(guī)劃(2006—2020年)》《黑龍江省國民經(jīng)濟和社會發(fā)展第十四個五年規(guī)劃和2035年遠景目標(biāo)綱要》及撓力河流域自然狀況和經(jīng)濟社會的發(fā)展需求以及前人研究的經(jīng)驗,為撓力河流域未來土地利用格局的發(fā)展,本文構(gòu)建了基線情景、農(nóng)業(yè)發(fā)展情景、生態(tài)保護情景3種情景方案,對不同情景下各地類面積變化進行分析。

(1)基線情景

在基線情景下,撓力河流域繼續(xù)保持2010—2020年的土地利用相關(guān)政策,土地利用需求變化按2010—2020年各地類的變化進行設(shè)定。如表10所示,到2030年撓力河流域耕地面積變?yōu)? 490 099 hm2,較2020年減少0.32%,人造地表地類面積由2020年的39 116 hm2增加到2030年的44 186 hm2,增加12.91%,其余地類均有相應(yīng)的增加和減少,森林和濕地較2020年增加0.17%和1.75%,草地和水體較2020年減少1.92%和3.13%。

表10 不同情景下土地利用類型面積Tab.10 Land use type area under different scenarios hm2

(2)農(nóng)業(yè)發(fā)展情景

在農(nóng)業(yè)發(fā)展情景下,撓力河流域優(yōu)先保證耕地的充足供應(yīng),保證耕地面積在不減少的情況下,優(yōu)先滿足農(nóng)業(yè)發(fā)展需求。持續(xù)推動森林地類中包括樹冠蓋度為10%～30%的疏林地,草地地類中的稀樹草原、荒漠草原以及濕地地類中有淺層積水的部分土地向耕地的轉(zhuǎn)移。設(shè)定地方政府大力開展森林、草地、濕地、水體地類的結(jié)構(gòu)調(diào)整,將撓力河流域耕地的潛力全部釋放。如表10所示,2030年撓力河流域耕地面積變?yōu)? 848 900 hm2,較2020年增加23.68%,森林、草地、濕地和水體面積變?yōu)?42 778、63 458、64 744、10 572 hm2,與2020年相比分別減少50%、50%、40.02%、29.31%,人造地表地類面積由于模型的參數(shù)設(shè)置未發(fā)生變化。

(3)生態(tài)保護情景

在生態(tài)保護情景下,撓力河流域以綠色發(fā)展為導(dǎo)向,嚴(yán)格保護森林、草地、濕地和水體地類面積的減少,大力開發(fā)潛在生態(tài)空間。設(shè)定地方政府通過退耕還林還草和植樹造林等措施,將水體地類中的部分坑塘以及耕地地類中的部分雨養(yǎng)旱地和牧草種植地轉(zhuǎn)為生態(tài)用地。如表10所示,2030年撓力河流域森林、草地、濕地和水體地類面積分別為728 146、190 508、140 572、18 064 hm2,較2020年分別增加49.94%、50%、30.23%、20.78%,在2030年撓力河流域生態(tài)用地約占研究區(qū)總面積的47.47%,耕地由2020年的1 494 961 hm2減少到2030年的1 153 158 hm2,減少22.86%,同樣人造地表地類面積由于模型的參數(shù)設(shè)置未發(fā)生變化。

4.4.2土地利用空間布局及模擬

根據(jù)3種情景方案預(yù)測撓力河流域2030年各地類的面積(表10),用2020年土地利用數(shù)據(jù)作為模擬的起始年份,繼續(xù)使用15個驅(qū)動因子、各情景下ELAS和轉(zhuǎn)移矩陣,輸入CLUE-S模型中模擬2030年撓力河流域土地利用空間布局,最終生成3種情景下土地利用類型圖(圖3)。

圖3 2030年情景模式模擬圖Fig.3 Scenario pattern simulation diagrams in 2030

(1)基線情景下,撓力河流域只有人造地表擴張趨勢程度大,森林和濕地用地擴張趨勢程度較小,其余地類均呈現(xiàn)出不同程度的縮減現(xiàn)象,對研究區(qū)的糧食和生態(tài)安全產(chǎn)生威脅。該模擬情景下,人造地表用地和森林地類擴張主要位于寶清縣的中部和西南部,占用了耕地的發(fā)展空間,對該地區(qū)農(nóng)業(yè)的發(fā)展會產(chǎn)生一定的影響。

(2)農(nóng)業(yè)發(fā)展情景下,撓力河流域只有耕地用地呈現(xiàn)擴張趨勢,耕地穩(wěn)定性較好,糧食安全得到充足的保障。盡管研究區(qū)耕地用地集中分布,但是土地-水資源需求壓力影響下,饒河縣的東南部、寶清縣的東部和西南部出現(xiàn)大量的森林用地和草地銳減情況,還出現(xiàn)了富錦市西南部濕地消失的現(xiàn)象,由此可見該情景下生態(tài)-水資源安全會對撓力河流域的發(fā)展帶來較大壓力。

(3)生態(tài)保護情景下,撓力河流域森林、草地、濕地出現(xiàn)了擴張趨勢,主要集中在寶清縣的西南部和雙鴨山的南部。由于生態(tài)空間的擴張主要是耕地面積的減少,因此該情景下糧食生產(chǎn)能力將會受到一定程度的影響,經(jīng)濟發(fā)展同樣會受到制約。

4.5 耕地適宜性評價結(jié)果

根據(jù)評價指標(biāo)體系和評價模型,利用ArcGIS空間分析功能,最后得出2030年撓力河流域模擬3種情景下的耕地適宜性評價等級以及各評價等級下的面積比,如圖4、表11所示。

表11 耕地適宜性評價面積統(tǒng)計Tab.11 Acreage statistics of cultivated land adaptability evaluation

4.5.1綜合適宜性評價

從表11和圖4可以看出,基線情景下耕地在中度適宜區(qū)是重點分布區(qū)域,面積最大,為1 341 013 hm2,占比59.09%;高度適宜區(qū)面積為450 608 hm2,分布在寶清縣、集賢縣和友誼縣部分區(qū)域;勉強適宜區(qū)面積為417 505 hm2,分布在高海拔地區(qū)、寶清縣的西南和富錦市小部分區(qū)域;暫時不適宜區(qū)同樣分布在地勢高的區(qū)域,是該情景下面積最小的區(qū)域,為 60 442 hm2,由于地勢高低不平,不利于農(nóng)作物的種植,因此被劃分為暫時不適宜區(qū)。農(nóng)業(yè)發(fā)展情景下耕地在中度適宜區(qū)面積最大,為1 197 922 hm2,和基線情景下相比面積減少143 091 hm2;高度適宜區(qū)面積為426 690 hm2,分布區(qū)域和基線情景下幾乎相同,減少的部分主要位于友誼縣和富錦市北部;勉強適宜區(qū)面積為633 974 hm2,主要位于雙鴨山市轄區(qū)南部,寶清縣的東北部、西南部,富錦市南部以及饒河縣南部區(qū)域;暫時不適宜區(qū)和基線情景下相比面積減少49 460 hm2,減少的區(qū)域演化成勉強適宜區(qū),可能原因是環(huán)境的改變和人類活動的影響造成該區(qū)域的演化。生態(tài)保護情景下耕地在中度適宜區(qū)面積最大,為1 068 453 hm2,占比為47.08%;高度適宜區(qū)面積為367 174 hm2,和前兩種情景相比分別減少83 434 hm2和59 516 hm2,減少區(qū)域分布在友誼縣的北部和撓力河下游附近;勉強適宜區(qū)面積和前兩種情景相比分別增加216 469 hm2和332 963 hm2;暫時不適宜區(qū)面積為83 473 hm2,與前兩種情景相比增加的區(qū)域位于原有基礎(chǔ)上向周圍擴張部分和富錦市北部區(qū)域。綜合來看,3種情景下耕地位于高度適宜區(qū)面積占比均達到16%以上,中度適宜區(qū)面積占比達到47%以上,勉強適宜區(qū)面積占比達到18%以上,這是因為適宜區(qū)域的土地主要分布在地勢平坦地帶,水熱條件好,土層深厚,作物產(chǎn)量高,適宜農(nóng)作。耕地位于暫時不適宜區(qū)在基線情景下面積占比為2.67%,農(nóng)業(yè)發(fā)展情景下為0.49%,生態(tài)保護情景下為3.67%,這是因為暫時不適宜區(qū)域土地主要分布在寶清縣的西南部、東南部和饒河縣的西部區(qū)域,地勢坡度較大,屬于丘陵地帶,不適宜農(nóng)作物種植。

4.5.2適宜性因素分析

影響撓力河流域2030年3種情景下耕地適宜性分布差異的重要因素是自然因素,包括土壤、地形地貌等。地形地貌對耕地適宜性評價有很大影響,表現(xiàn)在適宜區(qū)主要分布在平原地區(qū),不適宜區(qū)由于海拔的增加,導(dǎo)致耕地的水熱條件發(fā)生變化,從而影響土壤類型及土壤條件發(fā)生變化,不適宜耕作;研究區(qū)由于特色的地質(zhì)條件形成大面積的黑土、黑鈣土和草甸土,黑土和黑鈣土有機質(zhì)含量高,土壤土質(zhì)肥沃,利于農(nóng)作物生長,因此研究區(qū)絕大部分區(qū)域的土地條件對農(nóng)作物的生長十分有利。

5 結(jié)論

(1)研究區(qū)土地利用類型以耕地為主,2000—2020年森林、濕地、水體和人造地表呈增加趨勢,耕地和草地呈減少趨勢。

(2)研究區(qū)2000—2020年土地利用轉(zhuǎn)移主要發(fā)生在耕地、草地和森林之間,2000—2010年草地面積轉(zhuǎn)入最多,耕地次之,森林最少,轉(zhuǎn)出面積相反;2010—2020年森林轉(zhuǎn)入和轉(zhuǎn)出面積最多,草地次之,耕地最少。

(3)經(jīng)模擬驗證,Kappa系數(shù)為0.894,模擬總體精度為91.18%,精度較高,可以模擬未來地類變化。

(4)2030年基線情景下人造地表面積增加12.91%,占用了耕地的發(fā)展空間,將對研究區(qū)糧食安全產(chǎn)生威脅;農(nóng)業(yè)發(fā)展情景下耕地得到了擴張,穩(wěn)定性較好,糧食安全得到保障;生態(tài)保護情景下生態(tài)用地約占研究區(qū)總面積的47.47%,耕地面積減少22.86%,糧食生產(chǎn)能力受到一定程度的影響,經(jīng)濟發(fā)展同樣會受到制約。

(5)根據(jù)評價指標(biāo)體系和評價模型,分析計算出3種情景下耕地均位于適宜區(qū)的面積占比達到96%以上,不適宜區(qū)面積占比不足4%,可見研究區(qū)96%以上的土地適合耕作,土地利用基本合理。