基于MCE-CA-Markov模型的佛山市土地利用景觀格局演變及模擬預測

莫惠萍，黃宇斌，胡欣雨，林美玲，陸東芳

(福建農林大學 園林學院，福州 350002)

隨著城市化進程的加快，人類對于自然資源的開發和利用不斷加大，同時也帶來了氣候、環境、土地利用變化等各方面的問題.

景觀格局變化是人為和自然因素共同作用下生態環境變化的綜合反映[1]，它能夠反映區域內的生態狀況及空間變化[2].土地是人們生存和發展的基本場所[3]，它的變化一直受到人類活動的影響，但人類對于土地的不斷開發和對土地覆被的改造都是人類發展過程中不可避免的.對土地利用景觀格局變化與預測進行研究，能夠幫助人們有效了解其成因與機制，對合理的土地利用計劃和管理具有重大意義.

目前，國內外對土地利用景觀格局的研究主要集中在時間和空間的動態格局演變、各種類型的數量轉移、時間和空間的模擬和預測等方面.所運用到的模型主要有馬爾可夫鏈(Markov)模型[4-5]、元胞自動機(CA)模型[6-8]、CA-Markov 模型[9-11]、Ann-CA-Markov 模型[12]、FLUS 模型[13-14]等.其中CA-Markov 模型主要用于模擬各種用地類型間的相互轉化的趨勢[15]，它很好地將CA模型和Markov 模型兩者的優點有效結合，在模擬土地變化的時空格局和避免受到人為因素的干擾方面具有較好的表現.本研究所采用的MCE-CA-Markov 模型由多準則評價(MCE)方法與CA-Markov 模型結合，它量化了自然和社會經濟等因素，補充了CA 模型的轉換規則，使模型更加優化[16-17]，具有科學性和實用性.

珠三角地區是我國改革開放的前沿陣地.21世紀以來，隨著工業化和城市化進程的大大加快，區域和生態系統土地利用結構頻繁調整[18]，城鎮建設與生態環境的協調成為珠江三角洲城市群需迫切解決的難題[19].佛山市地處珠江三角洲腹地，是我國工業化、城鎮化快速發展的典型代表[20].改革開放以來，佛山市始終以工業產業為主導，導致區域內生態環境遭到破壞，嚴重制約了佛山市的可持續發展.如今隨著廣佛都市圈和粵港澳大灣區戰略的大力推進，城市擴張加劇、人地矛盾劇烈等問題仍然存在，生態保護與經濟社會發展之間的矛盾愈發突出，佛山市未來可持續發展面臨嚴峻挑戰.

因此，本文選取佛山市為研究區，以佛山市2005，2010 及2015 年的Landsat 遙感影像作為數據源，對佛山市2005—2015 年土地利用景觀類型的分布格局和空間演變進行分析，并考慮自然和社會2 大影響因素，基于MCE-CA-Markov 模型對佛山市未來的土地利用景觀格局進行模擬和預測，揭示研究區域的景觀格局時空變化特征及規律；結合佛山市發展特點提出區域內可持續發展的對策，以期為同類型區域的景觀格局演變和模擬預測提供參考，為研究區內未來土地利用計劃和生態保護提供理論支持.

1 研究區概況

佛山市位于廣東省的中南部(112°2 2′～113°23′E、22°28′～23°34′N)，東接廣州，西靠中山，目前是粵港澳大灣區的極點城市以及廣佛都市圈的重要組成部分；全市土地總面積3 797.72 km2，下轄禪城、南海、順德、高明和三水5 區(見圖1)；境內地勢平緩，坡度多＜15°.改革開放以來，佛山市一直是廣東省經濟最發達的城市之一，其經濟總量在中國大陸城市中居于前列.然而，其經濟社會的快速發展也導致了區域內城市的快速蔓延、農業用地不斷減少、生態環境惡化等問題的不斷涌現.

圖1 研究區概況

1.1 數據來源及預處理

本研究所用數據主要包括佛山市2005 和2010 年Landsat5 TM 及2015 年Landsat8 OLI 遙感影像數據、ASTER 30 m 數字高程數據和2019年全國路網數據.其中，Landsat 遙感數據主要來源于USGS 平臺(https://earthexplorer.usgs.gov/)；ASTER 30 m 數字高程數據來自地理空間數據云平臺(http://www.gscloud.cn/)；路網數據源自OpenStreetMap.遙感數據的處理主要通過ENVI完成，以2005，2010 和2015 年的遙感影像為數據源，分別對其進行幾何校正、輻射定標、大氣校正和監督分類，參照《土地利用現狀分類標準》[21](GB/T21010—2017)和研究區域的實際情況，將各年份的景觀類型劃分為耕地、林地、草地、水體、建設用地和未利用地6 種用地類型，并通過精度驗證Kappa 系數使其達到0.8 以上，最終獲得2005，2010 和2015 年3 期景觀格局現狀數據，如圖2 所示.

圖2 佛山市2005，2010 和2015 年景觀格局現狀

1.2 研究方法

1.2.1 土地利用動態度

土地利用動態度可以表現研究區的土地利用動態變化情況.本文主要運用了單一土地利用動態度和綜合土地利用動態度2 個指標來反映研究區土地利用面積變化的劇烈程度.其中，單一土地利用動態度一般用以表達某一地類在一定時間范圍內的數量變化情況[22]；綜合土地利用動態度則用來描述區域內總的土地利用變化速度[23].其具體公式如下：

其中，k為單一土地利用動態度；Ua和Ub分別為研究初期和末期某一地類的面積；T為研究初期與末期的間隔時長，單位一般為a；本研究k為某一地類的年變化率，年變化量大小與k成正比，k值越大，則表示該地類的年變化量越大；R表示綜合土地利用動態度；Ui為第i類土地利用類型研究初期的面積.

1.2.2 土地利用轉移矩陣

土地利用轉移矩陣一般用來表示某時期中，期初與期末各土地利用類型之間的轉換關系[16]，它能夠通過定量描述的方式，更為直觀地將各種土地利用類型的變化方向全面、具體地表現出來[24].其表達式如下：

其中，S為土地利用類型的面積；i為期初的土地利用類型；j為期末的土地利用類型；n為研究區域各土地利用類型.

1.2.3 景觀格局指數

景觀格局分析是進一步研究景觀功能及其變化的基礎[20]，而景觀格局指數能夠定量描述景觀的差異及其變化.本研究通過對同類型研究中與本研究地塊相似的景觀格局指數的選取進行綜合分析，并根據各指標的生態學意義分別從斑塊類型水平和景觀水平2 方面選取更為符合本研究區域的指標共11 個.其中，斑塊類型水平方面主要選取了斑塊類型面積、斑塊數量、景觀聚集度指數、最大斑塊面積指數、景觀破碎化指數、散布與并列指數；景觀水平方面選擇了斑塊密度、面積加權平均形狀指數、蔓延度指數、香農多樣性指數和香農均勻度指數.各指標生態學意義見表1.

表1 景觀格局特征指標及其生態學意義

1.2.4 CA-Markov 模型

CA-Markov 模型是CA 模型和Markov 模型的組合模型[25]，該模型可以定量并有效地模擬土地利用格局在時間上的時空動態變化[26].其模擬主要通過以下步驟實現.

1)設置轉化規則.本研究區域面積較小，為了使模擬結果更符合實際，將佛山市景觀類型柵格圖像的像元大小設置為30 m×30 m，然后以2005 和2010 年分別作為初期和末期，通過IDRISI軟件的Markov 模塊，計算佛山市2005—2010 年的土地利用轉移面積矩陣和轉移概率矩陣.

2)利用多標準評價模塊(MCE)制作適宜性圖集.首先，根據研究區已出臺政策中對生態的保護以及未來各用地類型的管控方向，設置空間分布概率圖中的限制因素和影響因素；其次，通過層次分析法(AHP)確定各因子的權重，結合限制因素和影響因素對景觀類型轉化的影響，根據MCE 和COLLETION EDIT 模塊將各景觀類型的空間分布概率圖集合，生成適宜性圖集.

3)確定濾波器、起始年份和CA 循環次數.使用5×5 濾波器，以佛山市2010 年為起始年，將CA 循環次數設為5，結合2005—2010 年的土地利用轉移概率矩陣和適宜性圖集，對研究區2015 年的景觀類型進行模擬.

4)檢驗模擬精度.第3)步的模擬完成后，將2015 年的模擬結果與2015 年景觀格局現狀進行對比，檢驗該模型的模擬精度，精度≥85%則表示模型的模擬精度較高.

5)預測研究區未來的景觀類型.以2015 年為預測起始年，CA 循環次數設定為15，結合2010—2015 年的土地利用轉移概率矩陣和適宜性圖集，預測佛山市2030 年的景觀類型分布.

2 結果與分析

2.1 土地利用結構變化及面積流轉

佛山市2005—2015 年土地利用結構變化及動態度變化情況分別如表2、表3 和表4 所示.從土地利用結構變化進行分析(見表2)：2005—2015 年，耕地、林地和建設用地始終處于主導地位，其中林地面積所占比重逐漸下降，耕地和建設用地的面積則呈上升趨勢，若進一步分析2015年土地利用情況，發現該態勢依舊持續；水域的面積則呈逐年減少的趨勢，其中2005—2015 年其占比由25.96%持續減少到了13.26%.

表2 佛山市2005，2010 和2015 年土地利用結構變化

從土地利用動態度分析(見表3 和表4)可知，2005—2015 年間各土地利用類型的動態度變化大小依次為：水域＞耕地＞建設用地＞未利用地＞草地＞林地.其中，在2005—2010 年間，佛山市土地利用各類型中，水域的變化最明顯，其面積減少了470.542 5 km2，動態度為-7.95%；然后依次是未利用地、耕地和建設用地；草地、林地的動態度較小，基本處于穩定狀態.2010—2015 年間，耕地較2005—2010 年動態度明顯減小，這說明耕地流失現象依然嚴峻；水域在該研究時段由-7.95%上升到了-0.39%，這說明水域面積大量減少的情況得到了有效控制；建設用地在該研究時段內較2005—2010 年動態度也從3.35%下降到了1.03%，這說明建設用地的擴展速度在一定程度上得到有效控制；林地由于基數大，動態度變化雖然減弱，面積卻也減少了16.403 4 km2，這說明林地的保護措施尚待加強；草地和未利用地的動態度雖然較高，但由于本身初期面積較小，數量上仍然沒有大的變化，這說明草地和未利用地都相對穩定.

表3 佛山市2005—2010 年土地利用動態變化情況

表4 佛山市2010—2015 年土地利用動態變化情況

佛山市2005—2015 年的土地利用流轉趨勢如表5 所示.不考慮自身轉化即未變化的部分，可以發現，在該時段內，耕地轉化為建設用地量較大，轉移了159.719 4 km2；林地同樣是轉化為建設用地的概率較大，達65.9%；草地較為穩定，轉出量較少，但也有51.9%轉為了建設用地；就水域而言，大多數轉化為了耕地，占轉出量的72.8%；未利用地大部分轉成了建設用地；建設用地中有69.0%轉為了耕地，但建設用地總體是在不斷擴張的，并且占用了大量的耕地以及部分林地和草地，這也是城鎮化的必然趨勢.

表5 佛山市2005—2015 年土地利用類型轉移矩陣 km2

2.2 景觀格局變化分析

2.2.1 基于類型水平的景觀格局指數分析

佛山市2005—2015 年基于類型水平的景觀格局指數如表6 所示.從CA、NP和AI這3 個指數分析，研究區2005—2015 年間，耕地的斑塊面積最大，且呈現持續變大的趨勢，斑塊數量則呈先增后減的變化規律，AI指數變化規律相同，以上變化說明該景觀類型受到人類活動的干擾，開始朝破碎化方向發展，在后期趨于穩定；林地的CA、NP指數逐漸減少，AI指數變化幅度較小，說明該景觀類型破碎化程度較低，但各斑塊間的連通性不強；水域在該時段內的變化也十分明顯，CA、NP和AI這3 個指數均減少，表明水域的破碎化程度較高；建設用地則表現為斑塊數量減少，斑塊面積和聚集度指數逐漸增加的變化規律，說明該時段內建設用地逐漸增多且連片發展；草地的NP、CA指數在研究時段內呈現前5 年不斷減少，后5 年開始增加的變化趨勢；未利用地整體變化幅度不大.

表6 佛山市2005—2015 年類型水平的景觀格局指數

從最大斑塊面積指數(LPI)分析，2005—2015年間，耕地的LPI指數呈現先減少后增加的變化規律，該景觀類型的優勢度先降低后凸顯；林地和水域的LPI指數逐年降低，說明林地和水域的優勢度在逐漸減弱；而建設用地在初始階段的LPI值雖然偏小，但在2005—2015 年間呈現先減少后增加的態勢，優勢度已經逐漸高于林地和水域，說明建設用地的優勢度也愈來愈明顯；草地和未利用地的LPI指數較低且變化幅度較小，這主要是由于二者本身的面積基數小.

從散布與并列指數(IJJ)方面分析可知，在研究時段內，耕地的IJJ指數逐年減少，說明該景觀類型聚集程度降低，分離程度增加；林地和草地的IJJ指數均逐年增加，說明2 類景觀類型的斑塊聚集程度增加；水域的IJJ指數則呈現先減少后增加的趨勢，說明該景觀類型斑塊間近年來受人類活動干擾強度降低，聚集度增加；建設用地的散布與并列指數最小，說明其聚集程度最低，容易受人類活動的影響，且與之相鄰的斑塊，大多也容易受人類活動的影響，排列比較復雜.

2.2.2 基于景觀水平的景觀格局指數分析

佛山市2005—2015 年基于景觀水平的景觀格局指數如表7 所示.由表7 可知，佛山市2005—2015 年整體的景觀格局變化相對穩定.其中，研究區的斑塊密度出現先減少后增加的變化，說明該時段內景觀破碎化程度逐漸增加.SHAPE_AM指數具有相同特征，說明區域內景觀的空間復雜程度增加，景觀斑塊的形狀逐漸變得不規則.而SHDI和SHEI指數在2005—2015 年間一直處于減少的狀態，說明研究區由于外界的干涉，導致了自然生物多樣性減少，景觀分布不均勻.CONTAG指數則呈持續上升趨勢，說明研究區的各斑塊類型的連接性較好.

表7 佛山市2005—2015 年景觀水平的景觀格局指數

2.3 土地利用景觀格局變化模擬及預測

2.3.1 土地利用轉移概率矩陣

佛山市2005—2010 年、2010—2015 年土地利用類型單位轉移概率矩陣由IDRISI 軟件中的Markov 模塊計算得到，結果如表8 和表9 所示.

表8 佛山市2005—2010 年土地利用變化轉移概率矩陣

表9 佛山市2010—2015 年土地利用變化轉移概率矩陣

2.3.2 適宜性圖集制作

1)限制因素.限制因素主要選取了水域和建設用地.通過IDRISI 中的reclass 工具對提取的水域和建設用地進行重分類，將水域和建設用地的部分賦值為0，其他背景值賦值為1，生成二值化圖像.

2)影響因素.影響因素主要選擇了高程、坡度、距道路距離、距城鎮用地距離和距水域距離5 個指標.首先，利用IDRISI 中的FUZZY 模塊中的“S”型和“J”型函數對各因子進行標準化處理，得到各地類的標準化圖層；其次，按照各自的轉換規則，分別將制定的限制條件和標準化后的影響因子圖層導入MCE 模塊，選擇加權線性合并法，并賦予各因子權重值，得到各景觀類型的空間分布概率圖；最后，根據 MCE 和COLLETION EDIT 模塊將各景觀類型的空間分布概率圖集合，生成所有景觀類型的適宜性圖集.

2.3.3 佛山市2015 年土地利用景觀格局演變模擬

基于以上數據的準備，以佛山市2005 年為初始年，利用IDRISI 中的CA-Markov 模塊，輸入土地利用變化轉移概率矩陣和各景觀類型的適宜性圖集，模擬循環次數設置為10，實現對佛山市2015 年土地利用景觀格局的模擬，結果如圖3(b)所示.為了驗證該模型的模擬精度，利用IDRISI中的CROSSTAB 模塊對2015 年模擬預測結果與2015 年現狀進行了對比(見圖3)，得到的Kappa系數為0.935 2，說明該模型模擬效果較好，與實際較符合，可行性較高.

圖3 2015 年佛山市景觀類型現狀與模擬結果對比

2.3.4 佛山市2030 年土地利用景觀格局預測

以佛山市2015 年的土地利用為基期，利用Markov 模型得到佛山市2010—2015 年的土地利用轉移概率矩陣，以及基于MCE 模塊構建的適宜性圖集，將循環次數設置為 15，重復CA-Markov 操作，實現了對佛山市2030 年景觀土地利用景觀格局的模擬，結果如圖4 所示.

圖4 2015 年景觀類型現狀與2030 年預測結果對比

為了直觀分析佛山市2015—2030 年各景觀類型的數量變化，本文對佛山市2030 年預測結果與2015 年現狀中各景觀類型的面積進行了統計對比，結果如表10 所示.由表10 可知，2030 年，佛山市各景觀類型變化繼續呈現出耕地與林地的適當減少，建設用地的迅速增加，以及水域、草地與未利用地不明顯變化等特點，這種變化趨勢是經濟發展較好的區域普遍存在的.其中，耕地面積較2015 年減少了171.31 km2，占比減少了4.51 %.為防止未來耕地持續減少，建議政府應加大政策管控，嚴格控制耕地保護紅線，限制盲目開墾.林地面積減少了102.45 km2，占比減少了2.70%，其變化強度要大于2005—2015 年該類型的變化，且林地作為研究區域內面積大小僅次于耕地和建設用地的景觀類型，其生態地位十分重要，為了實現區域內生態環境和社會經濟的可持續發展，政府應采取措施嚴格控制林地的減少量.草地面積所占百分比由 0.27%增長到了0.35%，面積增加了3.177 km2，占比增加了0.08%，變化幅度較小，但草地具有涵養水源、保持水土等生態功能，其地位也不能忽視.水域面積所占百分比由13.26%上升到了13.59%，占比增加了0.33 %，說明佛山市2015—2030 年間水環境得到改善，這與佛山市人民政府2015 年提出的恢復河流水生態系統等舉措有關.建設用地面積由2015 年的31.20%上升到了2030 年的38.00%，面積增加了258.15 km2，雖增加強度較2005—2015年有所減弱，但建設用地的轉換具有不可逆性，在規劃時，應合理控制建設用地的增加量.未利用地面積整體變化量較小，尚有1.395 9 km2處于閑置，可見，未來佛山市的城市建設應將未利用地充分利用.

表10 佛山市2015 與2030 年土地利用景觀格局面積變化對比

3 結論與討論

3.1 結論

1)從土地利用動態變化來看，佛山市2005—2015 年間，各用地類型的動態度變化依次為：水域＞耕地＞建設用地＞未利用地＞草地＞林地；土地利用流轉趨勢上，耕地、林地、草地、未利用地轉化為建設用地的概率較大，水域大部分轉化為耕地，小部分也轉為建設用地，建設用地則主要轉為耕地.

2)從景觀格局變化來看，2005—2015 年佛山市整體景觀格局變化相對穩定，景觀破碎化程度逐漸升高，空間復雜程度增加，斑塊形狀趨于不規則，景觀分布越發不均勻.耕地斑塊由于受到人類活動的干擾，存在一定破碎化的情況但總體趨于穩定；林地和水域的優勢度逐漸減弱，破碎化程度較高；建設用地逐漸增多且連片發展，破碎化程度較低，但聚集程度最低.

3)通過MCE-CA-Markov 模型對佛山市2030年景觀格局的模擬預測結果看，佛山市林地和耕地的面積持續減少，水域和建設用地面積出現增加，草地和未利用地基數小，總體變化不明顯.

3.2 討論

在利用MCE-CA-Markov 模型模擬佛山市未來的景觀格局時，還存在一定的不足.如在選擇限制因素制作模擬所需的適宜性圖集時考慮的因素不夠全面，地理、社會、經濟、資源環境等客觀因素都應考慮進去，這樣模擬結果才會更精準；其次是對尺度的考慮不夠，本文主要以30 m×30 m 的柵格大小進行模擬，但柵格大小對于模擬精度會有一定的影響，所以導致模擬結果與實際情況產生一定的差距；另外，影響模擬精度的還有時間間隔長短以及濾波器大小的設置.基于此，今后的研究應進一步探究不同柵格大小、時間間隔長短和濾波器大小對模擬精度的影響，以選擇最適合的尺度進行模擬，進而提高模擬精度.