許譯方，周永杰 ，何煜然鄺紅艷楊大勇

1.生態(tài)環(huán)境部華南環(huán)境科學(xué)研究所，廣州 510535

2.普渡大學(xué)健康與人類科學(xué)學(xué)院，印第安納州 47906

0 前言

生態(tài)系統(tǒng)服務(wù)(Ecological Services)是指人類通過生態(tài)系統(tǒng)的結(jié)構(gòu)、功能和過程直接或間接地得到的產(chǎn)品和服務(wù)[1-2]，支撐著人類社會(huì)的發(fā)展和進(jìn)步。但是，隨著快速的工業(yè)化和城鎮(zhèn)化，土地利用結(jié)構(gòu)發(fā)生了劇烈變化，致使生態(tài)系統(tǒng)服務(wù)功能受到損害和削弱[3]。Westman 首先提出“Nature’s services”的概念及其價(jià)值評(píng)估問題，Costanza等人在《Nature》發(fā)表的文章使生態(tài)系統(tǒng)服務(wù)價(jià)值評(píng)估研究成為熱點(diǎn)，而 2001年啟動(dòng)的千禧年生態(tài)系統(tǒng)評(píng)估項(xiàng)目(Millennium Ecosystem Assessment)進(jìn)一步探討了生態(tài)系統(tǒng)及其服務(wù)的驅(qū)動(dòng)因素和對(duì)人類福祉的影響[4]，生態(tài)系統(tǒng)服務(wù)價(jià)值已成為當(dāng)前生態(tài)學(xué)、資源經(jīng)濟(jì)學(xué)、環(huán)境經(jīng)濟(jì)學(xué)以及生態(tài)經(jīng)濟(jì)學(xué)的交叉前沿領(lǐng)域[5]。隨著GIS、RS等技術(shù)的發(fā)展及其在生態(tài)系統(tǒng)服務(wù)價(jià)值分析中的廣泛應(yīng)用，生態(tài)系統(tǒng)服務(wù)價(jià)值動(dòng)態(tài)變化的研究已成為當(dāng)前研究趨勢之一[6]。因此，認(rèn)識(shí)土地利用對(duì)生態(tài)系統(tǒng)服務(wù)價(jià)值的影響，可以從土地利用規(guī)劃的角度，為決策者制定生態(tài)保護(hù)規(guī)劃與管理制度提供可靠參考和依據(jù)。

本文以粵港澳大灣區(qū)城市群為研究區(qū)域，利用遙感影像計(jì)算其 2000—2015年的土地利用類型構(gòu)成和轉(zhuǎn)移矩陣，并通過馬爾科夫鏈預(yù)測2025年的土地利用結(jié)構(gòu)，同時(shí)在此基礎(chǔ)上分析其生態(tài)服務(wù)價(jià)值的變化，預(yù)測研究區(qū) 2025年的生態(tài)系統(tǒng)服務(wù)價(jià)值變化，揭示了土地利用與生態(tài)系統(tǒng)服務(wù)價(jià)值變化的聯(lián)動(dòng)關(guān)系，為粵港澳大灣區(qū)城市群制定生態(tài)安全管控政策提供分析決策依據(jù)，以期為促進(jìn)區(qū)域土地資源合理利用和生態(tài)環(huán)境保護(hù)提供理論支撐。

1 數(shù)據(jù)與方法

1.1 研究區(qū)概況

粵港澳大灣區(qū)城市群由廣州、深圳、珠海、佛山、東莞、中山、惠州、江門、肇慶9個(gè)地級(jí)市和香港、澳門兩個(gè)特別行政區(qū)組成，地處珠江下游，位于北緯 21 °—25 °，東經(jīng) 111°—116 °之間，區(qū)域國土面積約5.6萬km2，是繼美國紐約灣區(qū)、舊金山灣區(qū)、日本東京灣區(qū)之后，世界第四大灣區(qū)。區(qū)域氣候?qū)儆谀蟻啛釒У暮Ｑ笮约撅L(fēng)氣候，年均溫 22 ℃左右，夏季多雨且臺(tái)風(fēng)頻繁，冬季冷暖變化無常。2017年末常住人口約6956.93萬人，地區(qū)生產(chǎn)總值為10.18萬億元，折合美元1.6萬億元，位列四大灣區(qū)第二位，僅次于東京灣區(qū)的1.8萬億美元。

1.2 數(shù)據(jù)來源

本研究采用的數(shù)據(jù)來源于美國國家航空航天局(https://earthdata.nasa.gov/)提供的Landsat TM/ETM+遙感數(shù)據(jù)，選擇2000、2005、2010和2015共4個(gè)時(shí)間點(diǎn)(表1)，經(jīng)過大氣校正、正射校正和幾何校正等預(yù)處理后，根據(jù) IPCC土地覆被類型的分類原則，采用面向?qū)ο蠓诸惖姆椒?，結(jié)合人工和自動(dòng)識(shí)別，實(shí)現(xiàn)土地利用類型分類。

圖1 粵港澳大灣區(qū)地理位置圖Figure 1 The location map of Guangdong-Hongkong-Macau Greater Bay Area

由于受到自然環(huán)境和波譜的相互影響，從傳感器的光譜信號(hào)提取的地表信息存在不確定性，因此遙感分類的精度評(píng)價(jià)是極其重要的環(huán)節(jié)。精度評(píng)價(jià)一般采用定量評(píng)價(jià)，評(píng)價(jià)指標(biāo)主要有總體分類精度、用戶精度和Kappa系數(shù)等。

在研究區(qū)內(nèi)隨機(jī)抽取了 300個(gè)檢查點(diǎn)，利用Google Earth不同年份的高分辨率影像進(jìn)行精度評(píng)價(jià)，總體分類精度為 93.83%，Kappa 系數(shù)為 0.9472，滿足本研究分類精度要求。其中水域的用戶精度最高，達(dá)到 100%，其次是建設(shè)用地(95.49%)，未利用地的用戶精度最低，主要被誤分為草地。

1.3 研究方法

1.3.1 馬爾科夫鏈在土地利用預(yù)測中的方法原理

馬爾科夫鏈?zhǔn)且环N特殊的隨機(jī)運(yùn)動(dòng)過程，一個(gè)運(yùn)動(dòng)系統(tǒng)在t+1時(shí)刻的狀態(tài)和t時(shí)刻的狀態(tài)有關(guān)，而與以前的狀態(tài)無關(guān)，這點(diǎn)用于土地利用格局變化是合適的[7]。

用Pij表示已知在時(shí)刻tn系統(tǒng)所處的狀態(tài)Ii的條件下，在時(shí)刻tn+1系統(tǒng)所處的狀態(tài)Ij的條件概率，稱Pij為轉(zhuǎn)移概率[8]。

現(xiàn)有土地利用類型Ii(i=1，2，…，n)的取值(面積)S(1)，S(2)，…，S(n-1)，S(n)，其中令S(n)=ai的條件下，轉(zhuǎn)變到土地利用類型Ij(j=1，2，…，n)的取值S(n+1)=aj的條件概率為:Pij=P[S(n+1)=aj/S(1)，S(2)，…，S(n-1)，S(n)=ai]=P[S(n+1)=aj/S(n)=ai]。式中，Pij則為各地類的轉(zhuǎn)變概率[9]。

假設(shè)預(yù)測的土地利用類型可能處于Ii(i=1，2，…，n)的取值Si(i=1，2，…，n)，而且每次只能處于一種狀態(tài)，那么每一種狀態(tài)都具有n個(gè)轉(zhuǎn)向(包括轉(zhuǎn)向自身)即I1→I1，I1→I2，…，I1→In其轉(zhuǎn)變概率P(I1→I1)=P(S1/S1)，P(I1→I2)=P(S2/S1)，…，P(I1→Ij)=P(Sj/S1)，…，P(I1→In)=P(Sn/S1)。令P(Ii→Ij)=Pij，則可得到以下轉(zhuǎn)變概率矩陣:

式(1)中:0≤Pij≤1，且有即式中每個(gè)元素非負(fù)且各行元素之和等于1。

用Pij表示由土地類型Ii經(jīng)過m階轉(zhuǎn)變后到達(dá)地類Ij的轉(zhuǎn)變概率，可知馬爾科夫鏈的高階土地類型面積轉(zhuǎn)變概率矩陣:

1.3.2 生態(tài)系統(tǒng)服務(wù)價(jià)值測算方法

聯(lián)合國發(fā)起的千年生態(tài)系統(tǒng)評(píng)估(MA)將生態(tài)系統(tǒng)服務(wù)歸納為供給服務(wù)、調(diào)節(jié)服務(wù)、文化服務(wù)和支持服務(wù)4個(gè)功能類別[10]。本研究采用Costanza等[11]提出的模型計(jì)算研究區(qū)的生態(tài)系統(tǒng)服務(wù)價(jià)值。

式(2)中:ESV為生態(tài)系統(tǒng)服務(wù)價(jià)值；Ak為土地利用類型k面積(hm2)；VCk為土地利用類型k單位面積的生態(tài)價(jià)值系數(shù)。

2 土地利用變化和預(yù)測

2.1 土地利用變化

通過統(tǒng)計(jì)解譯后的遙感影像數(shù)據(jù)，粵港澳大灣區(qū)城市群土地利用類型以林地為主。2015年，林地面積約為31250.17 km2，占區(qū)域總面積的56.82%，16年間減少1505.24 km2，主要分布在羅浮山、鼎湖山和古兜山構(gòu)成的綠色屏障區(qū)；耕地在迅速減少，2000年面積約為 10436.72 km2，占區(qū)域總面積的18.98%，2015年面積約為 9298.03 km2，16年間減少1138.69 km2，主要分布在珠三角沖積平原地區(qū)；水域也占有較大面積，并且呈逐漸減少的趨勢，2015年水域面積為4888.30 km2，占區(qū)域總面積的6.89%；建設(shè)用地主要分布在各城市的主城區(qū)，其增長較快，由 2000年 11.68%增至 2015年 16.61%，面積增長2711.25 km2；草地與未利用地分布較少。

2.2 土地利用類型轉(zhuǎn)換特征

從 2000—2015年粵港澳大灣區(qū)土地利用轉(zhuǎn)移矩陣來看，林地在過去16年中發(fā)生了較大變化，其中有5.46%的面積轉(zhuǎn)入其它土地利用類型，面積最大的為建設(shè)用地，達(dá)到 957.39 km2，其次是轉(zhuǎn)入耕地，面積為394.81km2；從水域轉(zhuǎn)出情況來看，30.92%的原濕地面積流失，主要轉(zhuǎn)移方向是耕地和建設(shè)用地，其中轉(zhuǎn)入耕地面積為587.05 km2，占原水域面積的13.86%，轉(zhuǎn)入建設(shè)用地面積為634.83 km2，占原水域面積14.99%；16年間研究區(qū)耕地一直處于下降趨勢，有14.44%的耕地被轉(zhuǎn)變?yōu)榻ㄔO(shè)用地，面積為2087.68 km2，說明耕地的流失主要是人類的建設(shè)活動(dòng)所造成的。

圖2 2000—2015年粵港澳大灣區(qū)土地利用類型分布圖Figure 2 The maps of land use types in Guangdong-Hongkong-Macau Greater Bay Area from 2000 to 2015

2.3 土地利用預(yù)測

運(yùn)用馬爾柯夫鏈進(jìn)行預(yù)測的關(guān)鍵在于轉(zhuǎn)移概率的確定[12]。為了檢驗(yàn)?zāi)Ｐ皖A(yù)測的精度，本文以2015年的土地利用數(shù)據(jù)作為樣本，并利用MPE(平均相對(duì)誤差)對(duì)其進(jìn)行驗(yàn)證。結(jié)果表明，利用馬爾科夫鏈預(yù)測的各土地利用類型面積與實(shí)際各土地利用類型面積的平均相對(duì)誤差為4.08%，結(jié)果較為可靠。

本文以 2010—2015年的各土地利用類型轉(zhuǎn)換概率作為初始概率矩陣，利用2000—2015年各地類的年均轉(zhuǎn)換面積得到年均轉(zhuǎn)移概率矩陣(表3)，并用初始概率矩陣與年均轉(zhuǎn)移概率矩陣相乘，得到最終狀態(tài)矩陣(表4)，最后計(jì)算出粵港澳大灣區(qū)2025年各地類面積，其中耕地面積為10860.61 km2，林地面積為 28309.98 km2，草地面積為 1395.22 km2，水域面積為 3024.18 km2，建設(shè)用地面積為 11387.10 km2，未利用地面積為3.63 km2。

表2 2000—2015年粵港澳大灣區(qū)土地利用轉(zhuǎn)移矩陣Table 2 The transition matrix of ecological system in Guangdong-Hong Kong-Macau Greater Bay Area from 2000 to 2015 km2

表3 年均轉(zhuǎn)移概率矩陣Table 3 The transition probability matrix of land use types per year %

表4 最終轉(zhuǎn)移概率矩陣Table 4 The ultimate transition probability matrix of land use types %

3 生態(tài)系統(tǒng)服務(wù)價(jià)值變化和預(yù)測

3.1 以耕地為基準(zhǔn)的生態(tài)服務(wù)價(jià)值當(dāng)量修訂

生態(tài)系統(tǒng)服務(wù)價(jià)值核算方法主要有兩種，即基于單位服務(wù)功能價(jià)格的方法和基于單位面積價(jià)值當(dāng)量因子的方法[13]。本研究根據(jù)粵港澳大灣區(qū)的實(shí)際情況，采用謝高地[14]提出的生態(tài)服務(wù)價(jià)值當(dāng)量表進(jìn)行核算，并根據(jù)徐麗芬等[15]提出的方法以耕地為基準(zhǔn)對(duì)當(dāng)量表進(jìn)行修訂。

2015年粵港澳大灣區(qū)和全國的耕地單位面積糧食產(chǎn)量分別為5.25噸每公頃和5.48噸每公頃，得到生態(tài)服務(wù)當(dāng)量的地區(qū)修訂系數(shù)為0.96。根據(jù)廣東華南糧食交易中心華南糧網(wǎng)統(tǒng)計(jì)的2015年末糧食采購價(jià)格為 2939.7元每噸，在沒有人力投入的自然生態(tài)系統(tǒng)提供的經(jīng)濟(jì)價(jià)值是現(xiàn)有單位面積耕地提供的食物生產(chǎn)服務(wù)經(jīng)濟(jì)價(jià)值的 1/7[16]，計(jì)算出粵港澳大灣區(qū) 1個(gè)生態(tài)服務(wù)價(jià)值量因子的經(jīng)濟(jì)價(jià)值量為 2204.00元·hm-2，最終得到研究區(qū)不同生態(tài)系統(tǒng)單位面積的生態(tài)服務(wù)價(jià)值(表5)。

3.2 生態(tài)系統(tǒng)服務(wù)價(jià)值變化分析

從 2000—2025年粵港澳大灣區(qū)生態(tài)系統(tǒng)服務(wù)價(jià)值(圖3)來看，整體呈現(xiàn)逐漸下降趨勢，36年間共減少755.98億元，其中，時(shí)間序列2015—2020年的下降速度最快，共減少487.37億元。從各個(gè)生態(tài)系統(tǒng)類型的服務(wù)價(jià)值來看，林地和濕地的價(jià)值占總價(jià)值的 95%以上，其次是耕地；林地和濕地的服務(wù)價(jià)值損失明顯，26年間分別減少215.86和560.36億元，主要原因是林地退化和濕地萎縮較為顯著。

表5 粵港澳大灣區(qū)不同生態(tài)系統(tǒng)單位面積的生態(tài)服務(wù)價(jià)值Table 5 The ecological service values of the different ecological system unit area in Guangdong-Hongkong-Macau Greater Bay Area in 2015 元/(hm2·a)

圖3 2000—2025年各生態(tài)系統(tǒng)類型服務(wù)價(jià)值Figure 3 The service values of different ecosystem types from 2000 to 2025

從空間分布格局上看，粵港澳大灣區(qū)生態(tài)系統(tǒng)服務(wù)價(jià)值高值區(qū)主要分布在惠州、肇慶和江門，其變化相對(duì)較小；佛山和中山部分地區(qū)因水域面積較大而成為生態(tài)系統(tǒng)服務(wù)價(jià)值次高值區(qū)；而珠江口其他區(qū)域因建設(shè)用地比例較高而成為生態(tài)系統(tǒng)服務(wù)價(jià)值低值區(qū)，其縮減量也較為明顯。

4 結(jié)論與討論

本研究以土地利用變化為基礎(chǔ)，對(duì)粵港澳大灣區(qū)的生態(tài)系統(tǒng)服務(wù)價(jià)值進(jìn)行研究，主要結(jié)論如下:

(1)粵港澳大灣區(qū)在城市化的進(jìn)程中，主要表現(xiàn)為城鎮(zhèn)建設(shè)用地的迅速增長，大量的耕地、林地轉(zhuǎn)化為建設(shè)用地，從現(xiàn)狀來看，時(shí)間序列 2010—2015年更為顯著。根據(jù)馬爾科夫鏈預(yù)測的結(jié)果，2025年城鎮(zhèn)建設(shè)用地達(dá)到11387.10 km2。

(2)粵港澳大灣區(qū)的生態(tài)系統(tǒng)服務(wù)價(jià)值主要來源于林地、水域和耕地，但是由于城市發(fā)展的需要，建設(shè)用地的擴(kuò)張擠壓了其它地類空間，造成大量自然、半自然生態(tài)系統(tǒng)面積的喪失，導(dǎo)致 2000—2015年間生態(tài)系統(tǒng)服務(wù)價(jià)值整體減少 145.54億元，到2025年整體減少755.98億元。

綜上所述，盡管城市化的過程會(huì)導(dǎo)致自然生態(tài)系統(tǒng)的面積驟減，相應(yīng)的生態(tài)系統(tǒng)服務(wù)下降，但是通過土地利用優(yōu)化與景觀格局優(yōu)化，可以作為提升生態(tài)系統(tǒng)服務(wù)價(jià)值的有效優(yōu)化平臺(tái)和實(shí)際操作工具[4]。從地類開發(fā)的角度，對(duì)林地和水域進(jìn)行科學(xué)合理的旅游景區(qū)開發(fā)，可以提升其文化服務(wù)價(jià)值；從地類結(jié)構(gòu)的角度，將不適宜耕作的區(qū)域恢復(fù)成林地或?yàn)┩?，不僅可以有效防止水土流失，同時(shí)也有利于提升生態(tài)系統(tǒng)服務(wù)價(jià)值。