陳星宇 趙筱青 楊永貴 普軍偉 王 茜 李 琴 趙祖軍

(1.云南大學(xué)地球科學(xué)學(xué)院，云南 昆明 650500；2.云南大學(xué)國際河流與生態(tài)安全研究院，云南 昆明 650500；3.云南省生態(tài)環(huán)境監(jiān)測中心，云南 昆明650034)

頻繁的人為活動導(dǎo)致自然保護(hù)區(qū)有效保護(hù)范圍縮小、土壤退化、景觀破裂程度增加、景觀多樣性改變等問題[1-2]，同時引起生態(tài)系統(tǒng)結(jié)構(gòu)和功能的變化，從而導(dǎo)致生態(tài)系統(tǒng)服務(wù)功能變化。為促進(jìn)自然保護(hù)區(qū)持續(xù)發(fā)展，給保護(hù)區(qū)內(nèi)動植物提供良好的生態(tài)環(huán)境，有必要對保護(hù)區(qū)內(nèi)景觀格局與生態(tài)系統(tǒng)服務(wù)關(guān)系進(jìn)行研究。許多學(xué)者對景觀格局和生態(tài)系統(tǒng)服務(wù)關(guān)系進(jìn)行了探討。BENNETT等[3]定量研究了生態(tài)系統(tǒng)服務(wù)受景觀異質(zhì)性的影響程度。張明陽等[4]、童晨等[5]的研究均發(fā)現(xiàn)，生態(tài)系統(tǒng)服務(wù)與一些景觀格局指數(shù)密切相關(guān)。但是對于景觀格局與生態(tài)系統(tǒng)服務(wù)的空間分異規(guī)律及兩者的空間耦合關(guān)系還未進(jìn)行深入探討。本研究擬利用GeoDA軟件來探討這一問題。

白馬雪山國家級自然保護(hù)區(qū)(以下簡稱白馬雪山保護(hù)區(qū))是全球24個重點(diǎn)生物地理區(qū)之一，也是世界上生物多樣性保護(hù)熱點(diǎn)地區(qū)之一[6]。近年來，白馬雪山保護(hù)區(qū)也出現(xiàn)了生物資源和自然環(huán)境管理不完善、農(nóng)民過度放牧和薪材消耗量大、保護(hù)區(qū)周邊基礎(chǔ)設(shè)施建設(shè)大增、旅游業(yè)快速發(fā)展等問題，導(dǎo)致整體景觀格局發(fā)生變化，生態(tài)系統(tǒng)服務(wù)受到不同程度的影響。因此，本研究定量計算了白馬雪山保護(hù)區(qū)的景觀格局指數(shù)和生態(tài)系統(tǒng)服務(wù)功能，并基于空間自相關(guān)分析方法探討景觀格局與生態(tài)系統(tǒng)服務(wù)的空間關(guān)系，以期為白馬雪山保護(hù)區(qū)生態(tài)環(huán)境保護(hù)以及區(qū)域生態(tài)紅線劃定提供參考依據(jù)。

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

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

Sentinel-Ⅱ(2019)遙感影像數(shù)據(jù)來源于美國地質(zhì)勘探局網(wǎng)站(http://glovis.usgs.com/)；歸一化植被指數(shù)(NDVI)遙感影像數(shù)據(jù)來源于2019年美國國家航空航天局(https://earthdata.nasa.gov/)的MODIS13Q1產(chǎn)品，空間分辨率為250 m，時間分辨率為16 d；數(shù)字高程模型(DEM)數(shù)據(jù)，來源于地理空間數(shù)據(jù)云(http://www.gscloud.cn/)；氣象數(shù)據(jù)(包括氣溫、降水量、日照時數(shù)、太陽輻射等)來源于中國氣象數(shù)據(jù)網(wǎng)(http://data.cma.cn/)；土壤數(shù)據(jù)(包括土壤深度，根系深度和土壤砂粒、粉粒、黏粒及有機(jī)質(zhì)質(zhì)量分?jǐn)?shù)等)來源于世界土壤數(shù)據(jù)庫(HWSD)的中國土壤數(shù)據(jù)集(v1.1)，空間分辨率為1 km。所有空間數(shù)據(jù)均統(tǒng)一轉(zhuǎn)換為WGS_1984_UTM_47N坐標(biāo)系。

1.2 遙感數(shù)據(jù)處理

1.2.1 土地利用類型數(shù)據(jù)提取

用大氣校正、波段融合及影像裁剪等技術(shù)對Sentinel-Ⅱ(2019)遙感影像數(shù)據(jù)進(jìn)行預(yù)處理，在野外考察并掌握遙感影像總體特征的基礎(chǔ)上，建立解譯標(biāo)志，通過人機(jī)交互解譯的方式提取白馬雪山保護(hù)區(qū)2019年的土地利用類型數(shù)據(jù)。參考《土地利用現(xiàn)狀分類標(biāo)準(zhǔn)》(GB/T 21010—2017)，并結(jié)合白馬雪山保護(hù)區(qū)實(shí)際情況，將土地利用類型劃分為林地、草地、耕地、水域及水利設(shè)施用地、住宅用地、交通運(yùn)輸用地、其他用地7類一級用地，并進(jìn)一步劃分為喬木林地、灌木林地、疏林地、中覆蓋草地、高覆蓋草地、河流水面、湖泊水面、冰川及永久性積雪、內(nèi)陸灘涂、農(nóng)村宅基地、公路用地、裸土地、裸巖石礫地、水田、旱地15類二級用地。然后，根據(jù)野外調(diào)查時間應(yīng)與遙感數(shù)據(jù)同期或相近原則，選取各土地利用類型中具有代表性的地塊再次進(jìn)行野外實(shí)地驗(yàn)證，進(jìn)一步完善遙感解譯結(jié)果。最終解譯結(jié)果的Kappa系數(shù)大于0.85，解譯精度為88.18%。

1.2.2 NDVI的月和年數(shù)據(jù)合成

NDVI遙感影像數(shù)據(jù)通過最大合成法消除云、大氣和太陽高度角等干擾，獲取NDVI時間序列數(shù)據(jù)集。利用MRT軟件對NDVI遙感影像數(shù)據(jù)進(jìn)行投影批處理和平滑處理，每月兩個16 d的NDVI遙感影像數(shù)據(jù)用最大合成法取最大值作為該月的數(shù)據(jù)，然后用最大合成法取每月的最大值作為該年的數(shù)據(jù)。

1.3 景觀格局指數(shù)

斑塊數(shù)、邊緣密度、蔓延度和香農(nóng)多樣性4個景觀格局指數(shù)是研究區(qū)域景觀破裂程度的常用指標(biāo)[7-8]，每個柵格(30 m×30 m)的各景觀格局指數(shù)的計算參見文獻(xiàn)[9]。

通過Fragstats4.2軟件移動窗口法來分析景觀格局指數(shù)的空間分布特征。使用移動窗口法時，窗口半徑的選擇十分重要[10]。比較了300、500、1 000 m為窗口半徑時的景觀格局指數(shù)穩(wěn)定性，發(fā)現(xiàn)當(dāng)窗口半徑為500 m時，景觀格局指數(shù)的空間變異特征最穩(wěn)定，可以有效反映景觀格局的空間變異特征，因此選擇500 m為窗口半徑。

1.4 生態(tài)系統(tǒng)服務(wù)功能

(1) 水源供給服務(wù)功能計算

(1)

式中：ETn為n月參考蒸散發(fā)量，mm；Tn為n月氣溫，℃。

AWC=SR×(54.509-0.132×SA-0.003×SA2-0.055×SI-0.060×SI2-0.738×CL+0.007×CL2-2.688×OC+0.501×OC2)

(2)

式中：SR為土壤深度或根系深度的較小值，mm；SA、SI、CL、OC分別為土壤砂粒、粉粒、黏粒和有機(jī)質(zhì)質(zhì)量分?jǐn)?shù)，%。

(2) 土壤保持服務(wù)功能計算

每個柵格的土壤保持服務(wù)功能通過InVEST模型的SDR模塊計算定量。導(dǎo)入該模塊的相關(guān)數(shù)據(jù)有降水侵蝕因子(R，通過式(3)計算得到)、土壤可蝕性因子(K，通過式(4)計算得到)，植被覆蓋、作物管理和土壤保持措施因子參考《InVEST模型用戶手冊》設(shè)置，地形因子從DEM數(shù)據(jù)中提取。

R=0.444 88·p0.969 82

(3)

式中：p為降水量，mm。

(4)

(3) 植被凈初級生產(chǎn)力(NPP)服務(wù)功能計算

每個柵格的NPP服務(wù)功能是基于改進(jìn)的CASA模型形成的植被凈初級生產(chǎn)力(NPP)插件在ENVI軟件中計算定量的。需要導(dǎo)入的相關(guān)數(shù)據(jù)有降水量、太陽輻量、氣溫、土地利用類型、NDVI等數(shù)據(jù)。通過上述數(shù)據(jù)導(dǎo)入模型求出月NPP數(shù)據(jù)，并將一年12個月數(shù)據(jù)進(jìn)行累加獲取2019年NPP服務(wù)功能。

(4) 生境支持服務(wù)功能計算

每個柵格的生境支持服務(wù)功能通過InVEST模型的Habitat Quality模塊計算定量。導(dǎo)入該模塊的相關(guān)數(shù)據(jù)有土地利用類型數(shù)據(jù)、威脅因子數(shù)據(jù)(指土地利用類型數(shù)據(jù)中的農(nóng)村宅基地、公路用地、水田、旱地數(shù)據(jù))，對威脅因子的影響距離和土地利用類型對各威脅因子的敏感度參考《InVEST模型用戶手冊》設(shè)置。

1.5 雙變量空間自相關(guān)分析

空間自相關(guān)分析能看出兩個變量是否具有聚集性(包括全局空間自相關(guān)性和局部空間自相關(guān)性)[11]。本研究運(yùn)用GeoDA軟件來分析景觀格局指數(shù)和生態(tài)系統(tǒng)服務(wù)功能的空間自相關(guān)性[12]。

2 結(jié)果與分析

2.1 景觀格局指數(shù)空間分布特征

白馬雪山保護(hù)區(qū)的景觀格局指數(shù)空間分布見圖1。

圖1 景觀格局指數(shù)空間分布

(1) 斑塊數(shù)

白馬雪山保護(hù)區(qū)景觀的破裂程度較高，斑塊數(shù)總數(shù)達(dá)到8 684個，其空間分布總體呈北高南低的趨勢，說明白馬雪山保護(hù)區(qū)的北部地區(qū)景觀破裂程度更高。

(2) 邊緣密度

白馬雪山保護(hù)區(qū)景觀的邊緣密度平均值為31.37 m/hm2，其空間分布也呈北高南低的趨勢。

(3) 蔓延度

白馬雪山保護(hù)區(qū)景觀的蔓延度平均值為72.42%，其空間分布除呈南高北低的趨勢外還表現(xiàn)出四周較高中間較低的趨勢，同樣說明北部地區(qū)景觀破裂程度更高，同時邊緣地區(qū)連通性較好，而中間地區(qū)連通性較差，說明中間存在更多的小斑塊，破碎化程度更高。

(4) 香農(nóng)多樣性

白馬雪山保護(hù)區(qū)景觀的香農(nóng)多樣性平均值為1.26，其空間分布呈北高南低的趨勢，高值區(qū)在北部及邊緣地區(qū)，景觀元素較多，景觀多樣性較大；低值區(qū)在南部地區(qū)，景觀元素較少，景觀多樣性較小。這種特殊情況可能是斑塊數(shù)造成的。

2.2 生態(tài)系統(tǒng)服務(wù)功能空間分布特征

白馬雪山保護(hù)區(qū)的生態(tài)系統(tǒng)服務(wù)功能空間分布見圖2。

圖2 生態(tài)系統(tǒng)服務(wù)功能空間分布

(1) 水源供給服務(wù)功能

白馬雪山保護(hù)區(qū)水源供給服務(wù)功能為124.94～1 094.56 mm/(hm2·a)，其空間分布呈西南高、東北低的趨勢，這是因?yàn)槲髂系貐^(qū)植被覆蓋度較高，而東北地區(qū)的土地利用類型主要為水田、旱地、交通運(yùn)輸用地、水域及水利設(shè)施用地以及其他用地等，植被覆蓋度較低，受人類活動干擾較大。

(2) 土壤保持服務(wù)功能

白馬雪山保護(hù)區(qū)土壤保持服務(wù)功能為0～3 689.15t/(hm2·a)，總體較低，其空間分布呈南高北低的趨勢，南部地區(qū)土地利用類型主要是林地和草地，植被覆蓋度較高且土壤侵蝕量較小，水土保持情況較好，而北部地區(qū)人類活動影響強(qiáng)度較大，因此水土流失比較嚴(yán)重，土壤保持能力相對較弱。

(3) NPP服務(wù)功能

白馬雪山保護(hù)區(qū)NPP服務(wù)功能為169.27～431.68 g/m2，其空間分布呈四周高中間低的趨勢，特別是西北部、東北部、東南部區(qū)域都較高，因?yàn)樗闹艿穆佣容^高，植被資源豐富。

(4) 生境支持服務(wù)功能

白馬雪山保護(hù)區(qū)生境支持服務(wù)功能為0～1，其空間分布呈北低南高的趨勢。

2.3 景觀格局與生態(tài)系統(tǒng)服務(wù)的相關(guān)性分析

2.3.1 景觀格局指數(shù)與生態(tài)系統(tǒng)服務(wù)功能的全局空間自相關(guān)性

由表1可見，白馬雪山保護(hù)區(qū)NPP服務(wù)功能、水源供給服務(wù)功能、生境支持服務(wù)功能3種生態(tài)系統(tǒng)服務(wù)功能都與蔓延度呈空間正相關(guān)關(guān)系，而與斑塊數(shù)、邊緣密度、香農(nóng)多樣性呈空間負(fù)相關(guān)關(guān)系；土壤保持服務(wù)功能與蔓延度、斑塊數(shù)、邊緣密度呈空間正相關(guān)關(guān)系，與香農(nóng)多樣性呈空間負(fù)相關(guān)關(guān)系。

表1 景觀格局指數(shù)與生態(tài)系統(tǒng)服務(wù)功能的莫蘭指數(shù)

2.3.2 景觀格局指數(shù)與生態(tài)系統(tǒng)服務(wù)功能的局部空間自相關(guān)性

白馬雪山保護(hù)區(qū)水源供給服務(wù)功能與斑塊數(shù)、邊緣密度、香農(nóng)多樣性的高-高聚集區(qū)在中部偏西北地區(qū)，說明該區(qū)域景觀斑塊類型復(fù)雜多樣且水源供給服務(wù)功能好；高-低聚集區(qū)在南部偏東地區(qū)，該區(qū)域水源供給服務(wù)較好但景觀斑塊類型較為單一；低-高聚集區(qū)集中在北部地區(qū)，該區(qū)域水源供給服務(wù)不好但景觀斑塊類型較為復(fù)雜多樣；水源供給服務(wù)功能與蔓延度在空間上聚集區(qū)較少且分散。總的來說，北部地區(qū)斑塊數(shù)、邊緣密度、香農(nóng)多樣性對水源供給產(chǎn)生正向影響；南部地區(qū)斑塊數(shù)、邊緣密度、香農(nóng)多樣性對水源供給產(chǎn)生負(fù)向影響。

NPP服務(wù)功能與斑塊數(shù)、邊緣密度、香農(nóng)多樣性的高-高聚集區(qū)分布在北部地區(qū)，景觀類型豐富，景觀邊緣形狀相對不規(guī)整，景觀破裂程度較高，說明斑塊數(shù)、邊緣密度、香農(nóng)多樣性對NPP服務(wù)功能產(chǎn)生了正向影響；高-低聚集區(qū)分布在南部地區(qū)，景觀破裂程度較小，斑塊類型較為單一而NPP服務(wù)功能則較高；NPP服務(wù)功能與蔓延度的高-高聚集區(qū)分布面積很少，主要分布于東部和東北部的部分地區(qū)?？偟膩碚f，北部地區(qū)斑塊數(shù)、邊緣密度、香農(nóng)多樣性對NPP服務(wù)產(chǎn)生正向影響；南部地區(qū)斑塊數(shù)、邊緣密度、香農(nóng)多樣性對NPP服務(wù)產(chǎn)生負(fù)向影響。

土壤保持服務(wù)功能與蔓延度、斑塊數(shù)、香農(nóng)多樣性幾乎沒有高-高聚集區(qū)；低-高聚集區(qū)分布在北部地區(qū)，這3個景觀格局指數(shù)對土壤保持服務(wù)產(chǎn)生負(fù)向影響；低-低聚集區(qū)分布在南部地區(qū)，這3個景觀格局指數(shù)對土壤保持服務(wù)產(chǎn)生正向影響；土壤保持服務(wù)功能與蔓延度的空間聚集關(guān)系不明顯?？偟膩碚f，北部地區(qū)蔓延度、斑塊數(shù)、香農(nóng)多樣性對土壤保持服務(wù)產(chǎn)生負(fù)向影響；南部地區(qū)蔓延度、斑塊數(shù)、香農(nóng)多樣性對土壤保持服務(wù)產(chǎn)生正向影響。

生境支持服務(wù)功能與景觀格局指數(shù)的高-高聚集區(qū)較少且比較分散；生境支持服務(wù)功能與蔓延度、斑塊數(shù)、香農(nóng)多樣性的低-高聚集區(qū)分布在北部地區(qū)，這3個景觀格局指數(shù)對生境支持服務(wù)產(chǎn)生負(fù)向影響；高-低聚集區(qū)分布在南部地區(qū)，生境支持服務(wù)功能較好，3個景觀格局指數(shù)對生境支持服務(wù)產(chǎn)生負(fù)向影響?？偟膩碚f，景觀格局指數(shù)對生境支持服務(wù)產(chǎn)生負(fù)向影響。

3 討 論

3.1 景觀格局與生態(tài)系統(tǒng)服務(wù)及其空間關(guān)系研究方法探討

以往研究一般采用標(biāo)準(zhǔn)法計算景觀格局指數(shù)，雖然能夠定量反映景觀結(jié)構(gòu)組成和空間配置特征的簡單量化指標(biāo)，但不能從空間上反映景觀格局的特征[13]。本研究使用Fragstats4.2軟件移動窗口法計算景觀格局指數(shù)，分析景觀格局空間分布特征，可以實(shí)現(xiàn)區(qū)域尺度下景觀格局指數(shù)空間分布的可視化，形成景觀格局指數(shù)柵格分布圖，使景觀格局的信息更加明朗[14]，并為進(jìn)一步研究景觀格局指數(shù)與生態(tài)系統(tǒng)服務(wù)功能空間關(guān)系提供數(shù)據(jù)基礎(chǔ)。對于生態(tài)系統(tǒng)服務(wù)功能計算，以往一般使用物質(zhì)量、價值量和能值量方法等，結(jié)果的空間可視化表達(dá)和分析比較困難[15]，而本研究通過空間數(shù)據(jù)量化生態(tài)系統(tǒng)服務(wù)功能，并在地圖上以可視化的形式定量表達(dá)生態(tài)系統(tǒng)服務(wù)功能。景觀格局和生態(tài)系統(tǒng)服務(wù)的空間關(guān)系分析上，目前國際上較為常用的空間關(guān)系統(tǒng)計分析軟件有R軟件和GeoDA軟件，R軟件只能以散點(diǎn)圖的形式顯示，不能以地圖的形式顯示局部聚集性[16]，本研究將GeoDA軟件運(yùn)用于生態(tài)系統(tǒng)服務(wù)及景觀格局關(guān)系研究，可以生成局部空間自相關(guān)的顯著圖和聚集圖，從而可以更加直觀的看出生態(tài)系統(tǒng)服務(wù)和景觀格局的聚集關(guān)系，在景觀生態(tài)空間分析中有很好的適用性。

3.2 景觀格局與生態(tài)系統(tǒng)服務(wù)及其空間關(guān)系結(jié)果探討

對于景觀格局指數(shù)，白馬雪山保護(hù)區(qū)北部和邊緣地區(qū)景觀格局趨于破碎化，這是因?yàn)檫@些地區(qū)受人類活動干擾較多，用地類型主要為農(nóng)村宅基地、水田、旱地，人口密度較大，人類活動頻繁；南部和中部地區(qū)景觀格局破裂程度較小，這是因?yàn)檫@些區(qū)域主要用地類型為林地和草地，人口密度較小，受人類活動干擾較少，而且保護(hù)區(qū)成立后更是加強(qiáng)了對林區(qū)的保護(hù)。

對于生態(tài)系統(tǒng)服務(wù)功能，本研究利用InVEST模型計算的土壤保持服務(wù)功能為0～3 689.15 t/(hm2·a)，水源供給服務(wù)功能為124.94～1 094.56 mm/(hm2·a)，生境支持服務(wù)功能為0～1，與孫興齊[17]用InVEST模型計算的香格里拉生態(tài)系統(tǒng)服務(wù)功能較為接近。這是因?yàn)檫@兩個地區(qū)生態(tài)環(huán)境比較類似，且都使用的是InVEST模型，同時也說明了用InVEST模型計算生態(tài)系統(tǒng)服務(wù)功能的可行性。

在景觀格局指數(shù)和生態(tài)系統(tǒng)服務(wù)功能空間關(guān)系分析上，與張金茜等[18]、王航等[19]、董敏等[20]的研究結(jié)果相似，本研究也表明生態(tài)系統(tǒng)服務(wù)與景觀格局存在相關(guān)性，研究結(jié)果合理、科學(xué)。不過，不同景觀格局指數(shù)和生態(tài)系統(tǒng)服務(wù)功能的影響比較復(fù)雜，既有正向的影響也有負(fù)向的影響[21]。因此，景觀格局與生態(tài)系統(tǒng)服務(wù)的相互作用還需要進(jìn)一步研究。

4 結(jié) 論

(1) 景觀格局指數(shù)斑塊數(shù)、邊緣密度、香農(nóng)多樣性在保護(hù)區(qū)北高南低，而蔓延度除呈南高北低的趨勢外還表現(xiàn)出四周高中間低的趨勢，景觀破裂程度高導(dǎo)致景觀多樣性增加。

(2) 生態(tài)系統(tǒng)服務(wù)功能中水源供給服務(wù)功能為124.94～1 094.56 mm/(hm2·a)，土壤保持服務(wù)功能為0～3 689.15 t/(hm2·a)，NPP服務(wù)功能為169.27～431.68 g/m2，生境支持服務(wù)功能為0～1。

(3) 白馬雪山保護(hù)區(qū)NPP服務(wù)功能、水源供給服務(wù)功能、生境支持服務(wù)功能3種生態(tài)系統(tǒng)服務(wù)功能都與蔓延度呈空間正相關(guān)關(guān)系，而與斑塊數(shù)、邊緣密度、香農(nóng)多樣性呈空間負(fù)相關(guān)關(guān)系；土壤保持服務(wù)功能與蔓延度、斑塊數(shù)、邊緣密度呈空間正相關(guān)關(guān)系，與香農(nóng)多樣性呈空間負(fù)相關(guān)關(guān)系。不同景觀格局指數(shù)和生態(tài)系統(tǒng)服務(wù)功能的影響比較復(fù)雜，既有正向的影響也有負(fù)向的影響。