基于MCR-ANN-CA模型的包頭市生態用地演變模擬

劉建華 張啟斌 YANG Di 岳德鵬 于 強 楊 斕

(1.北京林業大學精準林業北京市重點實驗室， 北京 100083； 2.弗羅里達大學地理系， 蓋恩斯維爾 FL 32611)

0 引言

《全國生態環境保護綱要(國發[2000]38號)》中提到了生態用地一詞，明確指出了生態用地的重要性[1]。眾多學者基于此不斷發展生態用地的概念和內涵。一種觀點認為應該從土地類型上來劃分生態用地，即提供自然生態系統服務價值的土地都可以被視為生態用地。另一種觀點認為應該從土地的主體功能來區分，主張農業用地不應作為生態用地[2-3]。根據前人理論成果，結合包頭市實地情況，本文將林地、草地、水體作為生態用地范圍。

生態用地作為自然生態系統服務的基本載體，是解決城市建設用地擴張與生態保護矛盾的綜合途徑[4]。在我國西北干旱地區，生態用地在防止土地沙漠化、水土保持等方面的生態功能更加突出，對于維持區域生態系統健康穩定具有至關重要的作用[5-6]。對生態用地的空間變化進行模擬，可以為干旱地區生態用地科學規劃和生態安全格局構建提供重要參考價值。

生態用地的演變模擬是土地利用演變模擬的一部分，當前土地利用格局的演變模擬主要包括兩種類型，一種是宏觀回顧模型，如邏輯回歸(Logistic regression)模型[7]、馬爾可夫(Markov)模型[8]、系統動力學(System dynamic)模型等[9-10]，這些模型根據土地利用的時空數據，量化土地利用在宏觀時空尺度上的變化規律，對未來情景進行預測。另一種是微觀預測模型，包括多智能體模型(Multi-agent)[11]、元胞自動機(Cellular automata, CA)[12]等，這種模型從微觀尺度入手，通過研究土地利用與生態、人類活動等因子的相互作用，預測宏觀上土地利用的未來演變。將兩種模型結合，可得到二者的混合模型，既可以考慮全局宏觀效應，又可以考慮多個因子對土地利用的微觀作用機制，相比單一模型具有較大優勢。將宏觀回顧模型與元胞自動機結合，是這種模型的一種常見形式，如CA-Markov模型[13]。在這種模型中，有兩個關鍵點對其模擬精度有較大影響，一是元胞鄰域內不同土地利用類型間轉變規則，即鄰域規則；二是地形、生態環境等因子對土地利用演變影響，即適宜性規則[14]。當前傳統元胞自動機鄰域規則制定方法較多，如支持向量機、多目標決策、Markov模型等[15]，用戶往往難以選擇合適的規則定義方法，部分線性模型也難以根據鄰域特征精確模擬生態用地演化等非線性過程；基于邏輯回歸、多標準評價模型及層次分析法等制定的適宜性規則雖然考慮了多種因素對土地利用變化的影響，但多基于單個像元利用空間疊加手段獨立評價其適宜性，不能根據生態用地演化的過程量化研究區內各景觀單元演化為生態用地的適宜性。

針對上述問題，選取西北干旱地區典型城市包頭市為研究區，基于源匯理論，利用最小累積耗費阻力(Minimum cumulative resistance, MCR)模型量化生態用地演化過程中，從“生態源地”到其他土地利用類型的適宜性，同時利用人工神經網絡(Artificial neural network, ANN)提取CA鄰域轉換規則，構建MCR-ANN-CA模型對包頭市生態用地的演變過程進行模擬，以期為區域生態用地規劃及生態建設提供理論與方法支持。

1 材料與方法

1.1 研究區概況

包頭市位于內蒙古自治區西部(東經109°13′～111°26′，北緯40°13′～42°44′)，面積27 768 km2。包頭市深處內陸，氣候為典型的溫帶干旱、半干旱大陸性氣候，冬季寒冷干燥、夏季炎熱多雨[16]，年平均氣溫2.0～7.7℃，年均降水量175～400 mm，年均蒸發量為2 100～2 700 mm。包頭市可利用地表水總量為9×108m3，地下水補給量為8.6×109m3。黃河流經包頭境內214 km，水面寬130～458 m，最大流量6 400 m3/s，年平均徑流量為2.6×1011m3，是包頭市主要用水來源[17]。包頭市海拔976～2 317 m，地勢中間高南北低，北部丘陵、中部山地、南部平原分別占總面積的14.49%、75.51%和10%[18]。干旱的氣候條件與起伏的地貌特征使得包頭市生態環境較為脆弱，面臨較高的水土流失與土地沙漠化、荒漠化風險。近年來，包頭市城市建設規模擴張迅速，房地產開發項目、工業園區等建設項目不斷涌現，導致草地、林地、濕地等生態用地遭到占用與破壞，生態風險有所升高[19]。

1.2 數據來源與處理

本研究主要數據源包括包頭市土地利用數據(2006、2011、2016年)、歸一化植被指數(Normalized difference vegetation index, NDVI)、數字地面高程模型(Digital elevation model, DEM)、各區縣旗人口密度及工業園區分布數據。其中土地利用數據由Landsat-8遙感數據解譯得到，工業園區分布數據來自《包頭市土地利用總體規劃》和《包頭市城市總體規劃》；歸一化植被指數(NDVI)、數字地面高程模型(DEM)來自地理空間數據云(http:∥www.gscloud.cn)；各區人口密度數據來自包頭市統計年鑒。

根據研究需要，從遙感影像中提取耕地、林地、草地、建設用地、水體、未利用地并通過實地驗證確保其精度，其中林地、草地、水體為本文中生態用地范圍；利用DEM數據提取高程、坡度數據。

1.3 最小累積耗費阻力模型

MCR模型由KANPPEN提出，最早應用于物種遷徙過程研究，之后在物種保護、景觀格局分析等方面取得了廣泛應用[20]。該模型主要考慮3個因素，即“源”、阻力和累積代價，通過對3個因素的分析，對“源”克服阻力向外傳播所耗費的代價或者所做的功進行描述[21]。MCR模型的一般形式為

R=f∑DijRi

(1)

式中R——最小累積阻力

f——未知負函數，表示最小累積阻力與生態適宜性的負相關關系

Dij——從源j到景觀單元i的空間距離

Ri——景觀單元i處的阻力

本文將生態用地的演化過程看作生態用地對其他景觀的競爭性控制過程，且這種演化必須通過克服阻力實現，這樣生態用地的演化過程就可以抽象為從源(現有生態用地斑塊)到匯(其他景觀單元)克服阻力做功的水平過程[22]。由于區域下墊面差異，不同空間位置的土地演化為生態用地的阻力是不同的，通過由“源”到當前像元的累積阻力可量化當前像元演化為生態用地的概率，即阻力越大，該像元演化為生態用地的概率越小，演化為其他用地類型的概率越大。本文利用MCR模型，綜合考慮土地利用、NDVI、坡度、政府規劃工業園區、人口密度、水體距離、高程因子構建累積耗費阻力面，利用該累積耗費阻力面構建CA模型的適宜性規則，對上述過程進行模擬，提高模型預測精度。

1.4 元胞自動機模型

CA模型是一種時間、空間、狀態都離散的動力學模型，具有明顯的時間與空間特征，CA模型的一般形式為[23]

Si,t+1=f(Si,t,SN,t)

(2)

式中Si,t+1——元胞i在t+1時刻的狀態

Si,t——元胞i在t時刻的狀態

SN,t——元胞i的鄰域集合在t時刻的狀態

f——轉換規則

CA模型中，中心元胞在下一時刻的狀態是其在當前時刻狀態及其鄰域集合狀態的函數，準確定義該函數對于CA模型的模擬精度具有關鍵作用[24-27]。由于生態用地的演化過程是一種非線性的復雜動力學過程，當前元胞及其鄰域狀態影響中心元胞過程很難用簡單的規則定義，因此本文采用人工神經網絡模型提取CA模型轉換規則，對土地利用模擬過程中的現有規則進行改進[28-29]，以提高CA模型的模擬精度。

本研究中，元胞形狀為正方形的柵格像元，尺寸為30 m×30 m，采用3×3經典摩爾鄰域定義鄰域空間，模型初始柵格數據中的元胞狀態定義為研究區景觀格局類型，依次為耕地、林地、草地、建設用地、水體、未利用地，根據中心元胞及鄰域元胞狀態，利用ANN模型提取元胞轉換規則，利用該規則計算元胞演化為生態用地的概率。

1.5 人工神經網絡模型

ANN模型是20世紀80年代以來人工智能領域興起的研究熱點[30]。ANN模型在復雜的非線性系統的模擬方面具有明顯優勢，其自組織、自學習、聯想以及記憶的優勢能夠有效簡化CA模型，從原始訓練數據中提取CA轉換規則，避免主觀因素影響，提高模擬精度[31]。

由于ANN模型在解決此類非線性問題中的明顯優勢，國內外學者很早就開始了相關研究，嘗試利用ANN模型提取CA模型的轉換規則[32]。此類研究往往把鄰域規則和適宜性規則統一放入ANN模型進行規則提取，然而在模擬生態用地演化的CA模型中，相比適宜性規則，鄰域規則部分更難定義也更難被人類理解，因此本文利用ANN模型提取CA模型鄰域部分轉換規則，適宜性規則采用MCR模型構建。

采用經典BP神經網絡算法，網絡結構共3層，分別是輸入層、隱含層、輸出層。模型的輸入層共7個節點，隱含層共5個節點，采用tansing激勵函數，輸出層共1個節點，采用sigmoid激勵函數，其輸出為中心元胞轉換為生態用地的概率，值域為0～1。輸入變量及其取值范圍如表1所示。

表1 ANN模型輸入變量及取值范圍Tab.1 Input variables and their ranges of ANN model

1.6 MCR-ANN-CA模型

以CA模型為基礎框架，耦合MCR模型與ANN模型，構建MCR-ANN-CA模型用于包頭市生態用地的演化模擬，技術路線如圖1所示。模型模擬的基本步驟為：

圖1 技術路線Fig.1 Technical roadmap

(1)利用MCR模型，綜合考慮多種土地適宜性因子，構建研究區范圍內各像元演化為生態用地的累積阻力面。

(2)利用歷史數據的轉換情況及鄰域特征，構建訓練數據，訓練ANN模型，提取CA模型鄰域規則。

(3)以當前土地利用為CA模型的輸入數據，針對任一元胞，利用訓練好的ANN模型，根據其鄰域狀態判斷其演化為生態用地的概率P。

圖2 阻力因子的阻力分布Fig.2 Spatial distributions of resistance values of resistance factors

(4)將P與累積阻力面中像元最大值相乘，其結果為S(取值范圍為0到累積阻力面最大值)，若S大于當前元胞所對應的累積阻力，則該元胞演化為生態用地，否則，該元胞演化為非生態用地。

(5)將模擬結果轉換為柵格圖像并輸出結果。

2 結果與分析

2.1 基于MCR模型的適宜性規則構建

對土地利用、NDVI、坡度、政府規劃工業園區、人口密度、水體距離、高程等多種因子進行空間化與標準化處理，結果如圖2所示。其中將NDVI、坡度、人口密度、高程因子數據原始值進行歸一化處理作為其阻力。通過對不同用地類型阻力賦值，并對賦值結果歸一化處理得到土地利用數據。政府規劃工業園區阻力賦值為1，其他區域賦值為0，將數據柵格化后作為一個阻力因子。

根據包頭市2016年土地利用數據，提取包頭市生態用地范圍，如圖3a所示，由提取結果可知，包頭市生態源地主要分布在北部草原帶、大青山一帶及黃河沿岸，其他區域也有零星生態用地分布，大致形成了三屏多點的格局。

利用GIS疊加分析，將圖2中各因子進行空間疊加，得到包頭市各個空間位置演化為生態用地的阻力面，結果如圖3b所示。基于該阻力面，采用MCR模型，對源地進行累積耗費阻力計算，量化由近及遠的空間范圍內各像元演化為生態用地面臨的阻力，得到如圖3c所示的累積生態阻力面。由圖3可知，包頭市南北部生態阻力較低，中部累積生態阻力較高，尤其圖3c中A區域及B區域，累積生態阻力達到包頭市最高水平。

圖3 包頭市生態用地演化累積生態阻力面構建Fig.3 Construction of ecological land transition minimum accumulative resistance surface of Baotou City

2.2 基于ANN模型的CA鄰域規則提取

基于研究區2006、2011年景觀格局分布數據，通過Matlab隨機抽取3 000個景觀像元，統計2006年景觀格局分布數據中，以3 000個像元為中心的3×3鄰域內的各輸入變量大小，形成ANN模型訓練數據庫。選其中70%進行模型訓練，15%作為驗證數據集，15%作為測試數據集。

采用量化共軛梯度法對ANN進行訓練，模型經57次迭代后均方根誤差達到最小，為0.12，如圖4所示，圖中圓圈處為最優擬合點，此時驗證集的決定系數R2為0.90，模型擬合精度較好地滿足了本研究要求。

圖5 MCR-ANN-CA與CA-Markov模擬結果對比Fig.5 Comparison of simulation results of MCR-ANN-CA model and CA-Markov model

圖4 ANN模型均方根誤差Fig.4 Root mean square error of ANN model

2.3 MCR-ANN-CA生態用地模擬

將2011年景觀格局數據輸入訓練后的ANN模型，判斷各像元在鄰域規則影響下是否會演變為生態用地，進而根據MCR模型生成的累積耗費阻力面，判斷其最終演變方向,將模型模擬結果輸出為柵格圖像。包頭市2016年實際生態用地8 409.32 km2，MCR-ANN-CA模型的模擬結果中包頭市生態用地面積為8 670.01 km2，如圖5a、5b所示。

2.4 模型模擬精度對比分析

為驗證模型模擬精度，利用CA-Markov模型對包頭市2016年生態用地進行模擬，結果如圖5c所示。以2016年生態用地實際分布為參照，利用Idrisi Selva軟件中的Cross Tab模塊對2個模型的模擬結果進行分析，定量分析模型的模擬精度，結果如表2所示。表中KIA指數為卡帕一致性指數。

表2 模型模擬精度評價Tab.2 Accuracy assessment of simulation result

由表2可知，MCR-ANN-CA模型與CA-Markov模型模擬結果中，生態用地的面積均大于實際值，但總體上與2016年生態用地的實際分布保持了較高的一致性，其中MCR-ANN-CA模型的模擬精度略高于CA-Markov模型，二者的相對誤差分別為3.10%與5.31%，KIA指數分別為0.89和0.87。

相比CA-Markov模型，MCR-ANN-CA模型通過ANN模型提取了元胞自動機鄰域內的轉換規則，同時利用MCR模型構建累積耗費阻力面，對不同空間位置元胞演化為生態用地的阻力進行了量化，因此模擬精度得到了進一步提高。

3 結論

(1)利用ANN模型提取了元胞自動機的鄰域規則，同時利用MCR模型構建累積耗費阻力面，基于MCR-ANN-CA模型對包頭市生態用地演化情況進行模擬，模擬精度較高。

(2)將MCR-ANN-CA模型與CA-Markov模型模擬結果進行對比，KIA指數分別為0.89和0.87，相對誤差分別為3.10%和5.31%，MCR-ANN-CA模型對包頭市生態用地的演化過程具有更高的模擬精度。