基于“三生”用地變化的生態(tài)風(fēng)險時空格局分析和預(yù)測
——以粵港澳大灣區(qū)為例

聶 琪，湯江龍

（東華理工大學(xué)測繪工程學(xué)院，江西·南昌 330013）

生態(tài)風(fēng)險評價是十分重要的生態(tài)管理手段，可以分析和度量生態(tài)系統(tǒng)及其組成部分在受到外界脅迫后所產(chǎn)生不良生態(tài)后果的可能性[1-2]。城市土地利用變化是城市化的直接表現(xiàn)，隨著城市無序擴張進程加快，城市土地利用結(jié)構(gòu)趨于不合理[3-4]，生態(tài)平衡受到?jīng)_擊，人與自然面臨越來越多的考驗。而在中國“十四五”規(guī)劃階段，基于土地利用開展生態(tài)風(fēng)險研究，對形成高質(zhì)量發(fā)展的國土空間保護新格局具有重要意義。

目前國內(nèi)外生態(tài)風(fēng)險研究主要集中于風(fēng)險源與風(fēng)險受體分析、基于土地監(jiān)測數(shù)據(jù)進行生態(tài)風(fēng)險評價和生態(tài)風(fēng)險對策調(diào)控等方面[5-6]。而生態(tài)風(fēng)險評價方法基本可分為基于風(fēng)險源匯和基于土地利用景觀格局兩大類。在生態(tài)檢測資料缺乏歷史積累時，可基于土地利用數(shù)據(jù)闡述和評價土地利用方式對生態(tài)系統(tǒng)的影響[7]。土地利用對區(qū)域生態(tài)安全起決定性作用，不同的土地利用結(jié)構(gòu)對生態(tài)環(huán)境產(chǎn)生的影響具有區(qū)域性和累積性，并直接將結(jié)果反映于生態(tài)系統(tǒng)的組成和結(jié)構(gòu)上。2019 年，中央明確提出劃定“三生”空間的總體要求，基于“三生”空間的土地利用格局及其演變將土地利用變化與區(qū)域轉(zhuǎn)型發(fā)展相銜接，是研究土地利用變化問題的新視角[8]。土地利用的生態(tài)、生產(chǎn)、生活三大主導(dǎo)功能間的轉(zhuǎn)化是土地利用轉(zhuǎn)型的表現(xiàn)之一，也是地方社會經(jīng)濟轉(zhuǎn)型的必經(jīng)階段。

伴隨著聯(lián)合國千年生態(tài)系統(tǒng)評估報告發(fā)布[9]，研究視角正經(jīng)歷著由全球到區(qū)域再到地方的轉(zhuǎn)變，已有研究體現(xiàn)了基于土地利用進行生態(tài)風(fēng)險研究的優(yōu)勢，但大多側(cè)重于“從過去到現(xiàn)在”土地利用演變對生態(tài)風(fēng)險的影響的總結(jié)和驅(qū)動機制研究，而基于未來土地利用變化對生態(tài)風(fēng)險影響的研究有待深入；并且對土地生態(tài)風(fēng)險空間集聚效應(yīng)的時空演變研究相對不足，不利于國土資源和社會發(fā)展差別化政策的指定；同時，在土地預(yù)測模擬方面，多采用CLUE-S、CA-Markov、FLUS等模型。然而，F(xiàn)LUS模型[10]雖然能夠根據(jù)驅(qū)動因子與土地利用之間的定量關(guān)系有效模擬土地利用空間分布，但缺乏對土地利用需求數(shù)量的模擬。Markov 模型[11]側(cè)重時間維度分析預(yù)測土地利用變化，與FLUS結(jié)合，可實現(xiàn)土地利用在數(shù)量和空間的雙重模擬，為其他城市的土地利用生態(tài)風(fēng)險預(yù)測提供參考。

粵港澳大灣區(qū)作為中國首個正式審批建設(shè)的大灣區(qū)，其經(jīng)濟活力和城鎮(zhèn)化速度均位于中國前列，但高強度的經(jīng)濟開發(fā)和無序的城市擴張對原有生態(tài)環(huán)境的結(jié)構(gòu)和功能造成了破壞，加強大灣區(qū)生態(tài)文明建設(shè)刻不容緩。本文以粵港澳大灣區(qū)為例，引入生態(tài)風(fēng)險評價模型對區(qū)域過去、現(xiàn)在和未來的生態(tài)風(fēng)險的時空格局演變特征進行識別，借助空間統(tǒng)計分析方法、重心模型以及FLUS-Markov復(fù)合模型深入探究區(qū)域風(fēng)險的演變規(guī)律，研究結(jié)果可為促進區(qū)域“三生”空間協(xié)調(diào)發(fā)展、構(gòu)建景觀生態(tài)安全格局提供借鑒，為打造“美麗灣區(qū)”提供方向指引。

1 研究區(qū)概況與數(shù)據(jù)來源

1.1 研究區(qū)概況

粵港澳大灣區(qū)作為全球四大灣區(qū)之一，是由“九市二區(qū)”組成（圖1），即香港、澳門、肇慶市、佛山市、江門市、佛山市、廣州市、東莞市、中山市、深圳市、珠海市。該區(qū)域地處珠江下游，總面積約為5.6萬km2，呈“蝌蚪狀”，擁有約1.1萬km2“喇叭形”的平原面積，形成了“三面環(huán)山、一面臨海，三江匯合、八口分流”的特殊地貌。大灣區(qū)以全國5%的人口創(chuàng)造了約占全國1/7的GDP生產(chǎn)總值，已日漸成為全球人口和經(jīng)濟體量最大的灣區(qū)之一。但自改革開放以來，在經(jīng)濟高速發(fā)展，建設(shè)“一流灣區(qū)和世界級城市群”的同時，粵港澳大灣區(qū)土地利用程度加強，“三生”空間和生態(tài)風(fēng)險均發(fā)生了明顯的變化。

圖1 粵港澳大灣區(qū)位置及地形示意圖Fig.1 Schematic diagram of the location and topography of the Guangdong-Hong Kong-Macao Greater Bay Area

1.2 數(shù)據(jù)來源及預(yù)處理

本文主要采用了2000年、2010年以及2018年三期土地利用遙感數(shù)據(jù)、高程等地理數(shù)據(jù)、水系分布數(shù)據(jù)、全國路網(wǎng)數(shù)據(jù)和相關(guān)GDP等數(shù)據(jù)。其中，土地利用遙感數(shù)據(jù)、全國地市行政邊界、中國路網(wǎng)數(shù)據(jù)、全國居民點空間分布數(shù)據(jù)均來源于中國科學(xué)院資源環(huán)境科學(xué)與數(shù)據(jù)中心 (https://resdc.cn)，DEM數(shù)據(jù)來自于地理空間數(shù)據(jù)云(http://www.gscloud.cn/)，GDP相關(guān)數(shù)據(jù)來源于各城市的統(tǒng)計年鑒等。

利用ENVI軟件對獲取到的三期Landsat TM進行解譯，Kappa系數(shù)均高于80%，精度為30×30m，并按照《土地利用現(xiàn)狀分類標準》對用地進行分類整合。同時，從“三生”視角，參考已有研究成果[8,12]，按照國土管理行業(yè)標準將用地分為三大類和八小類形成“三生”格局，將土地利用類型基于滿足人們不同方面的需求劃分為生活、生產(chǎn)和生態(tài)功能，構(gòu)建“三生”用地分類體系以求充分體現(xiàn)粵港澳大灣區(qū)社會發(fā)展的多個維度，其具體分類如表1所示。

預(yù)測所需要的坡度、坡向數(shù)據(jù)通過ArcMap的3D Analyst tools對獲取到的DEM數(shù)據(jù)進行提取。

表1 “三生”用地功能分類Table 1 Functional classification of the production-living-ecological land use

2 研究方法

2.1 信息熵指數(shù)

信息熵能夠體現(xiàn)土地利用的有序性和均衡性，可以充分反映出土地利用系統(tǒng)的空間結(jié)構(gòu)的特征，熵值越大，土地利用系統(tǒng)就越紊亂，從而生態(tài)風(fēng)險的變化度就越大[13]。計算公式如下：

式中：P為土地利用結(jié)構(gòu)的信息熵，Ri是研究區(qū)各地類占土地總面積的比例，C為研究區(qū)土地的總面積，Ci是某一地類的面積。

2.2 區(qū)域生態(tài)風(fēng)險評價模型

生態(tài)風(fēng)險是指人為或自然因素影響下景觀格局與生態(tài)相互作用所產(chǎn)生的不利后果[14]。由于不同的生態(tài)系統(tǒng)所處的狀況發(fā)生條件對生態(tài)系統(tǒng)所造成的影響各異，可將生態(tài)風(fēng)險指數(shù)定義為各類景觀對生態(tài)所產(chǎn)生的風(fēng)險值，表示各景觀所導(dǎo)致生態(tài)風(fēng)險的貢獻度。因此參考已有研究[14-16]，將土地利用生態(tài)風(fēng)險權(quán)重指數(shù)和景觀破碎度引入生態(tài)風(fēng)險構(gòu)建區(qū)域生態(tài)風(fēng)險模型，建立土地利用類型與景觀生態(tài)風(fēng)險之間的聯(lián)系，以全面反映外界干擾下各景觀類型所代表的生態(tài)體系所受到的風(fēng)險威脅程度。景觀破碎度主要表現(xiàn)為景觀斑塊密度和景觀蔓延度。斑塊密度指數(shù)是指區(qū)域單位面積內(nèi)的斑塊數(shù)量，值越大說明景觀受干擾更強烈，景觀破碎化程度高；景觀蔓延度用來表征斑塊連通水平，值越大，景觀斑塊的連通性和完整度越好，景觀破碎化程度低。區(qū)域生態(tài)風(fēng)險計算公式如下：

式中：ERI表示土地利用生態(tài)風(fēng)險指數(shù)；Bi是第i類土地利用的面積；B代表研究區(qū)域的總面積；Qi是第i類土地利用所代表的景觀生態(tài)風(fēng)險強度指數(shù)，也就是i地類生態(tài)風(fēng)險權(quán)重；U是斑塊密度指數(shù)；V代表景觀蔓延度指數(shù)。

本文利用Matlab軟件，基于AHP分析法構(gòu)造各景觀的相對強度判斷計算各地類的生態(tài)風(fēng)險權(quán)重。但以往絕大多數(shù)的研究是通過層次分析法（AHP）采用單一方法求得權(quán)重，本研究以Matlab為平臺通過層次分析法代碼利用算術(shù)平均法、特征值法、幾何平均法求得權(quán)重，對三者進行求和平均（表2），可使權(quán)重結(jié)果愈加合理。而景觀斑塊密度和蔓延度利用Fragstats4.2軟件計算獲得。

通過計算分析，生態(tài)風(fēng)險權(quán)重一致性檢驗結(jié)果為0.029（小于0.1），表示強度指標矩陣具有滿意的一致性。各土地利用生態(tài)風(fēng)險權(quán)重如表2所示。

表2 土地利用生態(tài)風(fēng)險權(quán)重指數(shù)計算確定Table 2 Calculation of land use ecological risk weight index

2.3 空間地統(tǒng)計分析

根據(jù)景觀生態(tài)學(xué)研究[17]，單個風(fēng)險小區(qū)的面積應(yīng)是土地利用平均斑塊面積的2～6倍，利用Fragstats軟件計算可知粵港澳大灣區(qū)的平均斑塊面積約為3.02 km2，兼顧研究區(qū)域的面積和采樣工作量大小，經(jīng)過反復(fù)嘗試，決定采用5×5 km的格網(wǎng)對大灣區(qū)進行等距離采樣，共生成2595個生態(tài)風(fēng)險小區(qū)。

將計算出的2595個生態(tài)單元的生態(tài)風(fēng)險指數(shù)賦值給風(fēng)險小區(qū)的幾何中心，采用半變異函數(shù)分析方法，利用Gs+軟件對2000、2010和2018年的生態(tài)風(fēng)險指數(shù)與半變異函數(shù)進行擬合，根據(jù)擬合結(jié)果利用GIS對生態(tài)風(fēng)險小區(qū)的風(fēng)險指數(shù)進行插值，從而得到大灣區(qū)生態(tài)風(fēng)險的空間分布情況。半方差計算公式[18]為：

式中：y(h)是半方差；h為樣本距；E(xi)、E(xi+h)是在xi、xi+h的生態(tài)風(fēng)險指數(shù)；n是間距為h的樣本總數(shù)。運用Gs+軟件對三期數(shù)據(jù)進行半變異函數(shù)進行擬合后，利用GIS進行克里金空間插值?；诖鬄硡^(qū)的土地實際情況，為凸顯不同時期生態(tài)風(fēng)險變化特征將大灣區(qū)生態(tài)風(fēng)險插值結(jié)果劃分為5個等級：低風(fēng)險區(qū)（ERI＜0.03），較低風(fēng)險區(qū)（0.03＜ERI≤0.055），中風(fēng)險區(qū)（0.055＜ERI≤0.105），較高風(fēng)險區(qū)（0.105＜ERI≤0.135），高風(fēng)險區(qū)（ERI＞0.135）。

2.4 重心模型

重心模型主要基于重心位置偏離幾何重心的移動方向和距離展示出某要素的空間變化特征[19]，在探索空間要素演化軌跡上具有獨特優(yōu)勢。本研究采用幾何重心法衡量某要素的空間分布狀況，計算公式如下：

式中：X、Y為生態(tài)風(fēng)險重心的地理橫坐標和縱坐標；Pi是第i個生態(tài)風(fēng)險次級區(qū)域（每一個柵格單元）的屬性值；Xi、Yi指第i個生態(tài)風(fēng)險次級區(qū)域中心的橫、縱坐標。

2.5 FLUS-Markov復(fù)合模型

FLUS模型是基于元胞自動機（CA）模型和系統(tǒng)動力學(xué)（SD）模型，整合神經(jīng)網(wǎng)絡(luò)和輪盤賭選擇機制建立的土地利用情景預(yù)測模型，根據(jù)該模型的運行機制可以獲取基于土地利用與多個驅(qū)動因子的適宜性概率圖集和基于自適應(yīng)慣性機制的元胞自動機的土地利用情景模擬[20]。FLUS模型計算公式參見文獻[21]。

由于FLUS模型缺乏預(yù)測未來土地利用類型的數(shù)量變化，現(xiàn)耦合Markov模型，從而在數(shù)量和空間都充分挖掘土地利用動態(tài)演變信息。Markov模型的計算公式參見文獻[22]。

（1）土地利用驅(qū)動因子的選擇

綜合現(xiàn)有研究成果[14,23]、研究區(qū)域現(xiàn)狀，選取高程、坡度、坡向、到鐵路的距離、到河流的距離、到公路的距離、到居民點的距離和GDP作為土地利用驅(qū)動因子。

為驗證驅(qū)動因子的合理性，本研究利用SPSS軟件進行Logistic逐步回歸分析，檢驗土地利用與各驅(qū)動因子之間的相關(guān)性。受試者工作特征曲線 (Receiver operating characteristic curve, ROC)可用來分析各因子對地類的判斷解釋效能，當ROC值越趨近于1，各因子對地類的解釋程度越高。各土地利用類型的ROC檢驗結(jié)果表明：ROC值分別 為0.79、0.87、0.80、0.83、0.90、0.77、0.78、0.79，均大于0.73，上述所選擇的驅(qū)動因子可以用來描述大灣區(qū)土地利用的變化。

（2）精度檢驗

本文根據(jù)2010年的粵港澳大灣區(qū)土地利用模擬了2018年的土地利用情況，并采用Kappa系數(shù)以及FoM值對2018年的模擬圖層進行精度檢驗。2018年的檢驗結(jié)果顯示：Kappa系數(shù)以及FoM值分別為0.86和0.03，精度較高，該復(fù)合模型具有較好的預(yù)測效果，基于此，進一步開展2026年的土地利用格局預(yù)測。

3 結(jié)果分析

3.1 土地利用變化分析

根據(jù)“三生”格局分類標準，利用ArcMap的 “Reclass” 工具對三期土地利用圖層進行重分類得到粵港澳大灣區(qū)2000～2018年“三生”格局下土地分布圖（圖2）。整體上，大灣區(qū)地類以林地生態(tài)用地和農(nóng)業(yè)生產(chǎn)用地為主，分別占地類總面積的53.55%、22.36%，主要分布在肇慶、江門、惠州。研究區(qū)2000年、2010年和2018年土地利用信息熵分別為：1.22、1.34和1.37。大灣區(qū)總體信息熵都較高，城市化速度和對外開放的程度不斷加快，于2000～2010年信息熵變化較大，增幅達到7.14%，研究區(qū)土地地類變化明顯，土地利用系統(tǒng)的有序性降低。但進入后城市化階段，大灣區(qū)的土地利用信息熵增長變緩，2010～2018年增幅為2.11%。

圖2 研究區(qū)2000～2018年“三生”用地二級分布Fig.2 Secondary distribution of the production-living-ecological land use in 2000-2018

粵港澳大灣區(qū)地類變化趨勢主要可分為遞增、遞減、先減后增三種情況。呈現(xiàn)“階梯式”遞增的地類是城鎮(zhèn)生活用地和工礦生產(chǎn)用地，在2000～2010年期間，受產(chǎn)業(yè)轉(zhuǎn)移浪潮和“前店后廠”模式影響，大灣區(qū)城鎮(zhèn)生活用地和工礦生產(chǎn)用地面積顯著增加，分別增加了2498.71 km2和1023.16 km2，動態(tài)度也分別高達12.95%和22.06%，其增長主要集中在廣州、深圳、東莞、香港；經(jīng)歷急劇增長后，兩地類于2010～2018年出現(xiàn)增長幅度較小，動態(tài)度分別為0.97%、4.16%。呈現(xiàn)“V”型減少態(tài)勢的是牧草生態(tài)用地和農(nóng)村生活用地，兩地類在2000～2010年面積分別減少了128.39 km2、650.50 km2，在2010～2018年又呈現(xiàn)出面積回升態(tài)勢，但增長幅度較弱，僅增長了140.47 km2、2.97 km2。而伴隨著城鎮(zhèn)生活用地和工礦生產(chǎn)用地的急劇增加，大灣區(qū)社會經(jīng)濟迅速發(fā)展，農(nóng)業(yè)生產(chǎn)用地、林地生態(tài)用地、水域生態(tài)用地、其他生態(tài)用地都發(fā)生不同程度的減少。其中，農(nóng)業(yè)生產(chǎn)用地面積減少量高達1825.43 km2，其動態(tài)度為-1.26%，生態(tài)功能用地逐漸被生活功能和生產(chǎn)功能用地擠壓（表3）。

3.2 區(qū)域生態(tài)風(fēng)險變化特征

3.2.1 生態(tài)風(fēng)險時空分異

（1）生態(tài)風(fēng)險指數(shù)空間結(jié)構(gòu)分析

基于半變異函數(shù)[24]，利用GS+7.0對2595個生態(tài)單元的風(fēng)險值進行擬合，擬合結(jié)果如表4和圖3所示：由于2000至2018年的生態(tài)風(fēng)險值用球狀模型擬合效果最優(yōu)，可決系數(shù)分別為0.66、0.79、0.83，所以可基于球狀模型對粵港澳大灣區(qū)的生態(tài)風(fēng)險空間結(jié)構(gòu)進行分析。其中，塊金值呈遞增趨勢，但在研究后期增長速度放緩，表明在該研究尺度內(nèi)，生態(tài)風(fēng)險具有內(nèi)部變異性，且該變異性情況由隨機因素所引起的概率總體增大；基臺值由2000年的0.00105增加至2018年的0.00337，說明大灣區(qū)的城鎮(zhèn)化不斷發(fā)展，生態(tài)風(fēng)險空間分布不均性增強，差異性逐漸擴大。塊金效應(yīng)逐漸減小，變化范圍在20%～40%之間，變程值由2000年的93800 m增加至2018年的98500 m，這說明大灣區(qū)生態(tài)風(fēng)險的空間相關(guān)性不斷增強，生態(tài)風(fēng)險值相關(guān)性處于中等水平，相關(guān)性范圍不斷增加。

表3 粵港澳大灣區(qū)土地利用動態(tài)度變化Table 3 Changes in land use attitudes in the Guangdong-Hong Kong-Macao Greater Bay Area

表4 半變異函數(shù)擬合值Table 4 Fit values of the semi-variant functions

圖3 粵港澳大灣區(qū)生態(tài)風(fēng)險半變異函數(shù)擬合曲線Fig.3 Plot of the semi-variation function of ecological risk in the Greater Bay Area

（2）生態(tài)風(fēng)險時空格局分析

利用ArcMap地統(tǒng)計分析模塊進行克里金插值后得到三期粵港澳大灣區(qū)土地利用生態(tài)風(fēng)險等級圖（圖4）和面積變化表（表5）。結(jié)果表明：2000～2018年粵港澳大灣區(qū)的土地利用生態(tài)風(fēng)險等級面積變化較大，生態(tài)環(huán)境整體趨于惡化。較高風(fēng)險區(qū)域所占比例由1.43%增加到5.7%，高風(fēng)險區(qū)所占比例由0.39%增加到9.76%，兩區(qū)域所增加的面積基本與建設(shè)用地擴張區(qū)域重合。

自1978年改革開放以來，大灣區(qū)發(fā)揮著“試驗田”的作用，空間、人口、經(jīng)濟和社會的快速城鎮(zhèn)化，導(dǎo)致生態(tài)壓力不斷加劇，建設(shè)用地的持續(xù)擴張，景觀破碎化嚴重，生態(tài)風(fēng)險等級持續(xù)增高。2000～2010年期間，低風(fēng)險和較低風(fēng)險共減少了8968.65 km2，較高風(fēng)險和高風(fēng)險共增加6608.03 km2。2010年之后，實行了嚴格的耕地保護和節(jié)約用地制度，全力探索耕地保護嚴、建設(shè)占地少、用地效率高的科學(xué)發(fā)展之路等一系列措施，有效遏制了建設(shè)用地?zé)o效擴張，促進了土地利用結(jié)構(gòu)更加合理化。由此，在2010～2018年期間，較高風(fēng)險和高風(fēng)險面積呈現(xiàn)增加放緩趨勢，共增加了962.71 km2。

對“三生空間”土地利用生態(tài)風(fēng)險進行空間定位（圖4），粵港澳大灣區(qū)不同時期生態(tài)風(fēng)險呈現(xiàn)出明顯的空間差異性，總體呈現(xiàn)出東北、西北兩翼低，中間圈層、沿海方向半環(huán)式延伸的分布特點。

2000年，高風(fēng)險區(qū)主要呈球狀分布在人類活動密集的廣州市和深圳市，在研究期間，該高風(fēng)險區(qū)逐漸向外延伸擴張，基本與城市建設(shè)用地擴張空間重合，城市逐漸進入建設(shè)用地“存量挖掘”階段，建設(shè)用地從無序零散狀態(tài)逐漸向有序規(guī)模演變，2018年高風(fēng)險區(qū)主要分布于深圳、廣州、東莞、香港和澳門這些經(jīng)濟發(fā)達、現(xiàn)代化建設(shè)較好的中心城區(qū)，這主要因為隨著城市擴張和經(jīng)濟發(fā)展，這些區(qū)域人口急劇增多，人地矛盾逐漸惡化，需加強管制。2000年，低生態(tài)風(fēng)險區(qū)主要集中在肇慶市和惠州市，但2000～2018年期間，兩翼低生態(tài)風(fēng)險區(qū)面積減少了29.26%，主要是因為在城市發(fā)展過程中，人類行為破壞了土地的整體性，其分離度和破碎度逐漸變大。較高風(fēng)險區(qū)總體沿中部高風(fēng)險區(qū)向外延伸，分布較為集中，面積由795.70 km2增加至2018年的3163.02 km2，增加的區(qū)域主要分布在廣州、深圳、中山地區(qū)，主要由于建設(shè)用地零散分布，侵占大量農(nóng)業(yè)生產(chǎn)用地，景觀破碎化加劇，從而佐證了快速發(fā)展的城市化對生態(tài)風(fēng)險影響顯著。

圖4 2000～2018年粵港澳大灣區(qū)生態(tài)風(fēng)險等級空間分布Fig.4 Spatial distribution of ecological risk levels in the Guangdong-Hong Kong-Macao Greater Bay Area from 2000 to 2018

表5 大灣區(qū)2000～2018年土地利用生態(tài)風(fēng)險面積變化Table 5 Changes in land use ecological risk area in the Greater Bay Area from 2000 to 2018

3.2.2 生態(tài)風(fēng)險空間集聚特征分析

為探討土地利用生態(tài)風(fēng)險指數(shù)的空間相關(guān)性，利用GeoDa計算2000年、2010年和2018年的風(fēng)險值的空間相關(guān)性。三期全局Moran’s I指數(shù)分別為0.58、0.72、0.72，均為正相關(guān)并呈上升趨勢，表明粵港澳大灣區(qū)生態(tài)風(fēng)險存在集聚效應(yīng)，且有較強的空間相關(guān)性。局部空間自相關(guān)圖（圖5）表明，三期粵港澳大灣區(qū)生態(tài)風(fēng)險指數(shù)空間分布均主要以高—高和低—低聚集為主。ERI的高值聚集區(qū)主要分布在廣州、深圳，在2000～2018年期間高值聚集區(qū)范圍逐步擴大。低值聚集區(qū)在2000年主要分布在東北、西北兩翼，但自2000年后，該區(qū)域低值聚集區(qū)域面積明顯減少，于2018年，低值區(qū)主要集中在肇慶和惠州市。

圖5 2000～2018年局部空間自相關(guān)Fig.5 Local spatial autocorrelation in 2000-2018

3.2.3 城市較高和高風(fēng)險區(qū)分析

基于GIS以及EXCEL技術(shù)支持，得到粵港澳大灣區(qū)11個城市高風(fēng)險和較高風(fēng)險區(qū)的比例變化（表6），結(jié)果顯示：各城市在研究階段高風(fēng)險和較高生態(tài)風(fēng)險所占比例各異。在2000年，廣州高風(fēng)險區(qū)所占比例高達70.03%，深圳的較高風(fēng)險面積所占比例最高，達43.11%。伴隨著經(jīng)濟的高速發(fā)展，在2018年，除肇慶和香港的高風(fēng)險所占比例均為0之外，其他城市高風(fēng)險面積迅速擴大，其中廣州的面積較研究初期增加了695.72%，深圳增加了2376.32%；香港較高風(fēng)險區(qū)面積較基期變化最大，增加了266.21%，其次是廣州增加了183.48%，深圳卻減少了42.44%，東莞的較高風(fēng)險區(qū)面積比例最大，為20.5%，其次是廣州和佛山。

表6 研究區(qū)內(nèi)各城市高風(fēng)險區(qū)占比Table 6 The proportion of high risk areas in all cities in the study area

粵港澳大灣區(qū)生態(tài)高風(fēng)險和較高風(fēng)險重心軌跡見圖6，基于GIS的“Spatial Analyst”工具對克里金插值后的生態(tài)風(fēng)險指數(shù)計算高風(fēng)險和較高風(fēng)險的重心遷移軌跡，深入分析生態(tài)高風(fēng)險和較高風(fēng)險的空間演變特征。結(jié)果顯示：高風(fēng)險區(qū)重心2000～2018年均位于廣州，總體朝南移動，總位移為29206.69 m，表明大灣區(qū)南部的生態(tài)風(fēng)險變化較大，有一定程度的加劇；根據(jù)大灣區(qū)高風(fēng)險區(qū)的橢圓形狀，即扁平程度呈現(xiàn)逐漸降低的趨勢，表明長軸方向的城市生態(tài)風(fēng)險指數(shù)快速提高，推動標準差橢圓朝南部區(qū)域擴展，分布范圍逐步擴大，橢圓面積在研究期間增加了7386.01 km2。高風(fēng)險區(qū)位于珠江三角洲平原，地勢平坦，人類活動劇烈，土地破碎化嚴重，是土地利用生態(tài)風(fēng)險需要調(diào)控的重點區(qū)域。較高風(fēng)險區(qū)重心在2000年位于東莞，經(jīng)西南位移12456.44 m至廣州，再經(jīng)西北平移了2202.07 m基本與高風(fēng)險重心遷移趨勢重合；根據(jù)較高風(fēng)險區(qū)標準差橢圓形狀可知，其扁平程度呈現(xiàn)“先減后增”的趨勢，扁率自2000年的0.879降低至2010年的0.143，又在2018年增加至0.21，說明匯集在長軸方向的城市風(fēng)險指數(shù)不斷提高，推動了較高風(fēng)險標準差橢圓向西南方向發(fā)展，分布面積增加了5149.99 km2?？傮w而言，高風(fēng)險區(qū)和較高風(fēng)險區(qū)都主要位于城市中心周邊地區(qū)，由于城市的“虹吸”和“輻射”效應(yīng)，土地利用的蔓延度下降，景觀被各類生活用地切割，所以在經(jīng)濟發(fā)展轉(zhuǎn)型的同時，需對國土空間規(guī)劃加強生態(tài)指導(dǎo)，形成合理的土地利用空間布局。

圖6 較高風(fēng)險和高風(fēng)險區(qū)重心轉(zhuǎn)移軌跡Fig.6 Center of gravity shift trajectory in higher and high risk areas

3.3 大灣區(qū)土地利用生態(tài)風(fēng)險預(yù)測分析

3.3.1 土地利用格局預(yù)測

經(jīng)GeoSOS-FLUS軟件評估，2018年的Kappa系數(shù)以及FoM值分別為0.86和0.03，證明該復(fù)合模型具備較好的土地格局預(yù)測效果，基于此，本研究進一步對粵港澳大灣區(qū)2026年的土地利用格局進行預(yù)測，得到2026年的土地利用格局分布（圖7）和“三生空間”各地類的面積預(yù)測結(jié)果（表7）。結(jié)果表明：粵港澳大灣區(qū)2010～2018年的土地格局演變與2018～2026年的演變趨勢較為一致。農(nóng)業(yè)生產(chǎn)用地、林地生態(tài)用地、水域生態(tài)用地、其他生態(tài)用地面積持續(xù)下降，變化率基本與2010～2018年一致，城鎮(zhèn)生活用地、農(nóng)村生活用地、工礦生產(chǎn)用地面積均有不同程度的增加，而牧草生態(tài)用地面積下降了7.35%。其中，工礦生產(chǎn)用地變化率最大，2018～2026年的變化率為20.85%，農(nóng)村生活用地面積略微增加，2018～2026年階段的變化率為1.12%。

3.3.2 生態(tài)風(fēng)險預(yù)測

基于FLUS-Markov復(fù)合模型獲得粵港澳大灣區(qū)2026

圖7 粵港澳大灣區(qū)2026年土地利用格局預(yù)測Fig.7 Forecast of the land use pattern in the Guangdong-Hong Kong-Macao Greater Bay Area in 2026

表7 粵港澳大灣區(qū)2010～2026年景觀類型面積變化Table 7 Changes in the Guangdong-Hong Kong-Macao Greater Bay Area from 2010-2026

年生態(tài)風(fēng)險空間分布（圖8）和生態(tài)風(fēng)險等級面積占比（表8）。結(jié)果顯示，粵港澳大灣區(qū)2018～2026年生態(tài)風(fēng)險的演變規(guī)律與2010～2018年較為一致。到2026年，低風(fēng)險和較低風(fēng)險占比持續(xù)減少，主要向高風(fēng)險和中風(fēng)險轉(zhuǎn)移，高風(fēng)險區(qū)基本與中心城市重合，為廣州、深圳、佛山、東莞、香港和澳門，地區(qū)生態(tài)風(fēng)險進一步加劇，高風(fēng)險占比大幅上升至20.14%。上述現(xiàn)象反映了粵港澳大灣區(qū)2018～2026年生態(tài)風(fēng)險有明顯加劇，這主要是因為在2018～2026年生態(tài)風(fēng)險低的生態(tài)用地面積縮減，土地破碎化由于人類活動逐漸嚴重，生態(tài)風(fēng)險高的生活用地和工礦生產(chǎn)用地進一步擴張?；浉郯拇鬄硡^(qū)在建設(shè)“世界一流灣區(qū)”的同時，如果不改進現(xiàn)有的土地利用模式，城市化和工業(yè)化步伐加快，建設(shè)用地?zé)o序增加，生態(tài)風(fēng)險將進一步惡化，應(yīng)采取必要措施平衡經(jīng)濟發(fā)展和生態(tài)保護之間的關(guān)系，維護粵港澳大灣區(qū)的生態(tài)安全。

圖8 粵港澳大灣區(qū)2026年土地利用生態(tài)風(fēng)險分布預(yù)測Fig.8 Forecast of land use ecological risks in the Greater Bay Area in 2026

表8 2018～2026年不同等級生態(tài)風(fēng)險面積及比例變化Table 8 Changes in ecological risk area and proportion of different grades from 2018 to 2026

4 結(jié)論與討論

4.1 結(jié)論

本文基于2000～2018年粵港澳大灣區(qū)土地利用數(shù)據(jù)，運用區(qū)域生態(tài)風(fēng)險評價模型、重心模型以及地統(tǒng)計分析方法，對2000～2018年大灣區(qū)生態(tài)風(fēng)險時空格局分布變化進行了分析，并利用FLUS-Markov復(fù)合模型預(yù)測了2026年土地利用的動態(tài)變化及其對生態(tài)風(fēng)險變化的影響，主要結(jié)論如下：

（1）2000～2018年，大灣區(qū)土地景觀主要以農(nóng)業(yè)生產(chǎn)用地和林地生態(tài)用地為主，面積達76%以上，城鎮(zhèn)生活用地面積變化最大，增加了2841.45 km2，其次為農(nóng)村生產(chǎn)用地，減少了2058.35 km2；土地利用信息熵由2000年的1.22增加至2018年的1.37，土地利用結(jié)構(gòu)變化較大。

（2）大灣區(qū)土地利用生態(tài)風(fēng)險以中低風(fēng)險為主，但高風(fēng)險和較高風(fēng)險區(qū)面積呈上升趨勢，風(fēng)險等級空間呈“兩翼低、中間高”分布特點。生態(tài)風(fēng)險指數(shù)基臺值和塊金值均逐年增加，研究區(qū)內(nèi)生態(tài)風(fēng)險空間分布不均性增強，差異性逐漸擴大。在研究期間，低風(fēng)險和較低風(fēng)險面積分別下降了29.26%和21.74%，中風(fēng)險面積增加了45.35%；高風(fēng)險和較高風(fēng)險面積擴張迅猛，分別增加了5203.43 km2和2367.31 km2。

（3）粵港澳大灣區(qū)各城市高風(fēng)險和較高風(fēng)險區(qū)的變化各不相同，基本分為3種情況：不變、上升和下降。較高風(fēng)險和高風(fēng)險的重心均逐漸向南移動，標準差橢圓的扁率均逐漸降低，大灣區(qū)中南部的生態(tài)風(fēng)險變化較大，應(yīng)加強城市之間的生態(tài)聯(lián)系，促進粵港澳大灣區(qū)的一體化建設(shè)。

（4）粵港澳大灣區(qū)2018～2026年的土地利用和生態(tài)風(fēng)險趨勢與2010～2018年基本一致，研究區(qū)生態(tài)風(fēng)險進一步加劇，低風(fēng)險和較低風(fēng)險面積呈現(xiàn)持續(xù)減少的趨勢，分別減少3.52%和5.01%；高風(fēng)險區(qū)在空間上呈團狀集聚分布，并呈現(xiàn)出逐漸向外層擴大的趨勢，在2018～2026年期間面積增加了20.14%，需加強生態(tài)管理，平衡經(jīng)濟發(fā)展與生態(tài)保護之間的關(guān)系。

4.2 討論

本文從三生視角和建設(shè)“世界一流灣區(qū)”宏觀背景出發(fā)，基于景觀蔓延度、斑塊密度以及生態(tài)風(fēng)險權(quán)重對大灣區(qū)土地利用生態(tài)風(fēng)險進行定量分析，但土地利用生態(tài)風(fēng)險是一個極為復(fù)雜的系統(tǒng)，土壤污染和其他人類活動也會影響生態(tài)風(fēng)險，考慮到數(shù)據(jù)可獲得性及量化方法不夠成熟，本文并未將其列入。預(yù)測模型的精度與驅(qū)動因子密切相關(guān)，今后需將宏觀和微觀指標相結(jié)合，完善生態(tài)風(fēng)險評價模型，對研究區(qū)驅(qū)動因子進行合理更新，以更全方位地認識和預(yù)測土地利用生態(tài)風(fēng)險的時空特征，實現(xiàn)大灣區(qū)經(jīng)濟和生態(tài)的雙向共贏。當前，研究區(qū)存在較大的生態(tài)問題，比如土地利用結(jié)構(gòu)不合理、景觀破碎化嚴重、生態(tài)風(fēng)險持續(xù)加劇，但人們對于生態(tài)質(zhì)量的要求日益升高。為實現(xiàn)粵港澳大灣區(qū)綠色可持續(xù)發(fā)展，需對新增建設(shè)占用農(nóng)用地、生態(tài)用地的預(yù)案加強審批，嚴格遏制建設(shè)用地的無序擴張。同時，落實國土空間“三區(qū)三線”規(guī)劃工作，將由耕地資源單要素保護向山水林田湖草多要素保護轉(zhuǎn)變。本文僅對大灣區(qū)生態(tài)風(fēng)險進行定量分析，下一步研究應(yīng)討論生態(tài)風(fēng)險變化的影響因素和多尺度下的生態(tài)風(fēng)險評價，為粵港澳大灣區(qū)制定合理的土地利用決策和生態(tài)環(huán)境保護政策提供依據(jù)。