汾河上游流域景观破碎化与生态系统服务关系研究

董敏,苏常红,王亚璐

(山西大学 黄土高原研究所,山西 太原 030006)

0 引言

生态系统服务指人类从生态系统中获得的由生态系统结构及过程所形成的多种惠益,一般分为供给服务、支持服务、调节服务和文化服务等,与人类福祉息息相关[1-2]。景观破碎化是指在自然或人文因素的干扰下,景观由简单趋向于复杂的过程[3]。景观破碎化与景观格局、功能与过程密切联系,与自然资源保护相互依存。近年来,随着人类活动的逐渐增强,流域景观格局发生改变,而景观格局的变化导致生态系统组分、结构与生态过程及生物多样性等发生改变,对生态系统服务产生影响,进而影响区域可持续发展[4-5],因此,研究景观格局变化与生态系统服务的关系是区域可持续发展的重要课题。目前,针对这一问题已有相关研究,Matthew基于生态系统服务流和供需关系,重塑景观破碎化与生态系统服务关系研究框架[6];邹月等研究发现西安市景观格局破碎化导致生态系统服务总价值减少[7];王航等从数理统计角度分析了淮河上游景观格局对生态系统服务价值的影响,结果表明生态系统服务与景观格局丰度和聚集度变化存在明显相关[8],这些研究均为景观格局与生态系统服务关系进行了积极的探索,但是对二者空间耦合关系未进行进一步探讨。

汾河上游流域作为黄土高原地区独立的流域单元,是太原市唯一水源地。近年来,人类活动及城市化的加快使得研究区景观格局发生演变,土壤保持等服务类型均受到不同程度的影响,生态状况呈现复杂的空间格局。本文以汾河上游流域为研究区,运用网格分析方法和空间自相关等方法,对流域生态系统服务时空变化及流域景观破碎化空间分异规律进行深入分析,并基于GeoDa 软件探讨流域景观破碎化与生态系统服务的相关关系,为流域景观格局优化、资源合理开发利用以及区域可持续发展提供决策依据。

1 研究区概括与数据来源

1.1 研究区域

汾河上游流域位于管涔山至太原市尖草坪区上兰村区段(图1),地理位置为东经111°21′～112°27′,北纬37°51′～38°59′,流域全长126 km,总面积5 253.56 km2。涵盖宁武、静乐、岚县、娄烦、古交、阳曲以及太原市尖草坪、杏花岭等县(市、区)、89个乡镇,总人口约70万[9]。流域属温带大陆性季风气候,四季分明,春季干旱,夏季多雨而炎热,秋季雨量相对较少,冬季干旱而寒冷;流域内地貌类型复杂多样,包含山地、丘陵和河川阶地,其土壤类型主要有褐土、钙质粗骨土等[10]。汾河上游植被类型以草地和林地为主,山地河谷河岸分布有沙棘(Hippophaerhamnoides)、小香蒲(Typhaminima)等植物群落,以宁武县为中心的管涔山林区的天然次生林林相为华北之最,西部主要以针叶林、针阔叶混交林为主,东部以阔叶林为主,是山西省著名水源涵养用材林基地[11]。

图1 汾河上游流域位置示意图Fig.1 Outline of the Upper Reach of the Fenhe River Watershed

1.2 数据来源

气象数据包括降水、温度、辐射等下载于中国气象科学数据共享服务网(http:∥data.cma.cn/)。数字高程模型(Digital Elevation Model,DEM)数据源于中国科学院数据云平台(http:∥www.gscloud.cn.),空间分辨率为30 m。研究区景观类型图数据来自Landsat ETM+遥感数据,时间为2000、2008和2015年,空间分辨率为30 m,时相为7～8月。根据土地利用/覆盖变化(Land Use/Cover Change,LUCC)分类标准和汾河上游流域的土地利用情况,将土地利用类型分为耕地、林地、草地、水域、建筑用地和废弃地共6类(图2)。

图2 土地利用类型Fig.2 Land use of the study area

2 研究方法

2.1 生态系统服务评估

选取土壤保持、产水量和NPP生产3种典型的生态系统服务进行研究。运用InVEST模型(Integrated Valuation of Ecosystem Services and Trade-offs,InVEST model)对研究区2000、2008和2015年的土壤保持和产水量服务进行评估。土壤保持服务通过模型SDR模块,采用通用土壤流失方程(USLE)对区域土壤保持量进行估算,其具体方法是基于土地利用方式分别计算栅格单元潜在土壤流失量(UKLS)和实际土壤流失量(USLE),二者的差值即为土壤保持功能(SD)[12-13]。产水服务通过InVEST模型的Hydropower/water yield模块进行评估,具体原理为基于栅格单元的降雨量减去实际蒸发量即是实际产水服务。NPP生产通过CASA模型(Carnegie-Ames-Stanford approach,CASA model)进行计算,该模型利用遥感数据、温度、降水、太阳辐射以及植被类型等数据,并基于ArcGIS和ENVI软件计算得到NPP值[14]。本研究采用网格分析方法,根据前人研究和汾河上游流域范围大小,选取8 km×8 km的网格作为研究尺度,将流域划分为154个网格[15-17]。

2.2 景观格局指数

景观格局指数能够表征景观的空间结构特征,是景观异质性的具体体现,同时又是各种生态过程在不同尺度上作用的结果[18]。本研究从景观尺度水平上选取4个景观指数:斑块数(Number of Patches,NP)、边缘密度(Edge Density,ED)、蔓延度(Contagion Index,CONTAG)和香农多样性(Shannon’s Diversity Index,SHDI)来进行分析。以8 km×8 km的网格作为研究尺度,基于FRAGSTATS软件计算不同时期各网格的景观格局指数,然后运用GS+和ArcGIS地统计分析对各景观格局指数进行普通克里格插值[19],得到汾河上游流域3个时期各景观格局指数的空间分布图。

2.3 空间自相关分析方法

运用GeoDa软件分析生态系统服务与景观格局指数的空间相关程度。

(1)构建空间权重矩阵。建立空间权重矩阵是空间自相关分析的前提,它表征不同空间对象之间的空间分布,即在一定空间尺度准则下分析生态系统服务与景观格局指数的关系[20-21]。本文选取Rook邻接准则对154个网格单元构建基于空间邻接关系的权重矩阵。

(2)空间自相关。空间自相关包括全局空间自相关和局部空间自相关(Local indicators of spatial association,LISA)。全局空间自相关用来分析流域各网格单元某项生态系统服务评估量与景观格局指数之间的空间关联性。局部空间自相关用来分析流域网格单元某项生态系统服务评估量与邻近网格单元景观格局指数的空间相关程度。结果用Moran’s I来表示,Moran’s I的取值在[-1,1]之间,当I<0时表示负相关,当I>0时表示正相关,当I=0时,表示不相关,随机分布。为了检验 Moran’s I 是否显著,在 GeoDa 中采用蒙特卡罗模拟的方法来检验[22-24]。

3 结果与分析

3.1 生态系统服务时空变化特征

土壤保持服务以网格为单位输出(图3-1)。2000年土壤保持高值区主要位于流域北部以及西南边缘。2008年土壤保持服务空间分布大体呈现北高南低的分布态势,土壤保持高值区主要集中于流域北部,低值区主要分布在流域西南部。与2008年相比,2015年流域南部土壤保持高值区范围明显增加,流域北部土壤保持服务有所降低。2000、2008、2015年3个时期的土壤保持平均值分别为:238.01 t/hm2、259.27 t/hm2、178.09 t/hm2,呈现先增加后减小的趋势,流域东南部土壤保持能力有所增强。

产水量服务也以网格为单位输出(图3-2)。2000年产水量呈现西南高东北低的分布特征,产水量高值区集中在流域西南部。2008年产水量服务空间分布与2000年相反,大体呈现出东北高西南低的分布态势,流域产水量服务高值区主要集中于流域北部。2015年流域产水量服务也呈现西南高东北低的分布态势。2000、2008、2015年3个时期的产水量服务平均值分别为:190.90 t/hm2、259.58 t/hm2、202.90 t/hm2,呈现出与土壤保持类似的先增加后减小的趋势。

NPP生产服务空间分布情况如图3-3所示。2000年NPP生产高值区主要位于流域北部、东南部以及西南边缘。2008 年NPP生产高值区范围有所增大,同时NPP生产服务较低的区域主要集中在流域西部。2015年流域NPP生产高值区范围明显增加。2000、2008、2015 年3个时期的NPP生产平均值分别为:150.89 t/hm2、190.08 t/hm2、221.31 t/hm2,呈现出持续上升的态势。

图3 汾河上游流域生态系统服务时空变化Fig.3 Spatiotemporal variation of ecosystem service of the Upper Reach of the Fenhe River Watershed

2000、2008和2015年土壤保持与产水量服务均产生不同空间分布格局的主要原因:土壤保持受降雨侵蚀力因子和土地利用类型的影响,降雨侵蚀力因子的空间分布很大程度上决定了土壤保持的空间分布,而降水侵蚀力因子由降水量计算所得,因此,降水量是导致土壤保持产生不同空间分布格局的主要因素;产水量受降水量、蒸散发和土地利用类型的影响,其中降水量和蒸散发是影响产水量空间分布格局的主要因素。此外,NPP生产呈现持续上升的态势,表明研究区植被得到一定程度的恢复。

3.2 景观格局指数的时空特征

3.2.1 景观格局指数的时间变化特征

汾河上游流域景观格局指数的时空变化如图4。2000～2015年,流域景观格局指数蔓延度呈现减小趋势,表明流域景观中不同斑块类型团聚程度下降,斑块组成趋于破碎。此外,景观格局指数斑块数(NP)、边缘密度(ED)以及香农多样性(SHDI)均呈现持续增加的趋势,其中,斑块数(NP)由5 181个增加到了10 568个,表征流域景观中斑块总数不断增多,其值越大景观破碎程度就越高;边缘密度(ED)和香农多样性(SHDI)不同程度的增大,表征整个研究区景观在边缘形状上由简单变为复杂,景观的优势度逐渐降低,多样性提高。综合表明,在2000～2015年,汾河上游流域的景观破碎化程度增大。

图4 汾河上游流域景观格局指数的时间变化Fig.4 Temporal variation of landscape pattern indices of the Upper Reach of the Fenhe River Watershed

3.2.2 景观格局指数的空间分布规律

汾河上游流域景观格局指数斑块数(NP)空间分布如图5-1。2000年流域NP高值区与低值区交错分布,高值区主要分布在流域西南部,低值区主要分布在流域西南边缘和东南部;2008年NP高值区呈面状分布,流域西南部NP高值分布面积明显增大;2015年NP低值区集中于流域西北部。

景观格局指数边缘密度(ED)的空间分布特征与斑块数(NP)的空间分布相似,高值和低值区的空间分布规律基本一致(图5-2)。

景观格局指数蔓延度的空间分布与NP和ED相反(图5-3)。2000年流域CONTAG高值区集中于流域南部,低值区分布广泛;2008年流域南部CONTAG显著提高,同时流域西南部CONTAG减小;2015年流域西北部CONTAG大幅度减小。总体来说,2000～2015年,流域CONTAG减小,尤其以流域西南部最为明显;而流域东部和北部地区CONTAG明显增大。

景观格局指数香农多样性(SHDI)空间分布如图5-4。2000年流域SHDI高值区分布呈条带状,低值区主要分布在流域西北、西南边缘;2008年SHDI高值区分布范围在流域西南部有所扩大;2015年SHDI低值区分布面积明显增大,集中于流域北部和东南部。

综合景观指数的空间分布规律,表明流域西南部景观格局趋于破碎化,其破碎化程度增强,因为该区是农耕区,人口活动频繁,景观格局容易受到人类的干扰,从而加剧景观的破碎化。此外,流域北部景观格局变化也尤为强烈,其破碎化程度减弱,主要原因是该区主要为林区,人口密度小,人类对其的影响很微弱,再加上近几年对林区保护政策和措施,使得该地的景观破碎化趋于减弱的态势。

3.3 景观破碎化与各项生态系统服务的全局空间自相关

2000、2008和2015年流域各项生态系统服务评估量与景观格局指数的空间关联性如表1。2000～2015年土壤保持服务与景观指数NP、ED和SHDI呈现显著空间正相关关系,与景观指数CONTAG相关性不显著,表征土壤保持高值区景观趋于破碎。2008年,产水量服务与景观指数CONTAG呈显著正相关。2000～2008年,NPP生产服务与景观格局指数ED和SHDI的Moran’s I值基本接近于0 (|Z|≤1.96),NPP生产服务与景观破碎化的空间相关性很弱;2015年NPP生产服务与景观格局指数NP、ED、CONTAG和SHDI的Moran’s I分别为0.267、0.108、0.138、0.104、0.101和-0.139 (|Z|≥1.96),NPP生产与景观破碎化的空间相关性增强。

图5 汾河上游流域景观格局指数空间分布Fig.5 Spatiotemporal variation of landscape pattern indices of the Upper Reach of the Fenhe River Watershed

表1 汾河上游流域2000～2015年景观格局指数与生态系统服务的Moran’s I指数

综合分析,2000～2015年流域土壤保持、产水量与景观破碎化的空间相关关系均呈现先增强后减弱的趋势,且产水量与景观破碎化的空间相关性更显著。与土壤保持和产水量服务相比,NPP生产服务与景观破碎化的空间关联性相对较弱其原因是:影响NPP生产服务的因素有多种,景观格局演变只是其中一个因素,还有温度、降雨、辐射以及地形地貌等因素交互作用,共同对NPP生产服务产生影响。

3.4 景观破碎化与各项生态系统服务的局部空间自相关

2008年汾河上游流域土壤保持与景观破碎化的空间关联性最为显著(表1),选择2008年进行研究(图6a)。在99%置信度下,土壤保持与景观指数NP、ED和SHDI的高高聚集区集中分布在流域西南部,表征该地区土壤保持服务较高的同时景观斑块类型多样且边缘形状复杂,即地区土壤保持服务与景观破碎化呈正相关关系(P<0.05);低高聚集区主要分布在流域北部及西南部分区域,该区土壤保持服务与景观破碎化呈负相关关系(P<0.05);高低聚集特征不明显且分布分散。土壤保持与景观指数CONTAG的空间正相关关系微弱(P>0.05),集聚特征不明显。综合分析,流域北部景观破碎化对土壤保持产生负向影响,西南部景观破碎化对土壤保持服务产生正向影响。

2008年汾河上游流域产水量与景观破碎化的空间关联性最为显著(表1),因此选择2008年进行研究(图6b)。在99%置信度下,产水量高值与景观指数NP、ED和SHDI高值在流域西南部形成集聚现象,表征该区产水量与景观破碎化呈正相关关系(P<0.05)。产水量与景观指数NP和ED在流域北部和东南部地区形成了高低聚集区,表明这些区域产水量与景观破碎化呈空间负相关关系(P<0.05),即景观破碎化对产水量服务产生负面影响。综合分析,流域产水量与景观破碎化空间关联性在流域北部和东南部呈负相关,这表明该区景观破碎化对产水量服务产生负向影响。除此之外,流域西南部景观破碎化对产水量服务产生正向影响。

2015年汾河上游流域NPP生产与景观破碎化的空间关联性较为显著(表1),因此选择2015年进行研究(图6c)。在99%置信度下,NPP生产量与景观指数NP、ED和SHDI的空间集聚特征基本相似,高低聚集区主要分布在流域北部,说明该区NPP生产较高的同时景观斑块类型单一、边缘形状相对规整,景观破碎程度对NPP生产服务产生负向影响;低高聚集区分布面积较少主要集中在流域西南部,表明该区NPP生产与景观破碎化呈空间负相关关系(P<0.05)。NPP生产与景观指数CONTAG的高高聚集区分布面积很小,主要分布在流域东南部,说明该区NPP生产较高的同时有连通性较高的优势斑块存在,景观破碎程度较低。综合分析,流域北部景观破碎化对NPP生产服务产生负向影响,流域西南部景观破碎化对NPP生产服务产生正向影响。

4 讨论

景观破碎化通过影响生态系统结构、功能及过程，进而影响生态系统服务,这些影响可能是正向的也可能是负向的,最终导致生态系统服务的增加或减少[25]。本研究中流域北部景观破碎化对土壤保持和产水量服务产生负向影响,而流域西南部分地区景观破碎化对土壤保持和产水量服务产生正向影响,这可能与土地利用类型相关。流域北部宁武为汾河流域源头,主要为林区,林区景观破碎度在一定程度上能够反映植被的受损情况,其景观破碎度增大必然会导致其生态环境破坏,影响其正常的生态结构、过程及功能,致使其生态系统服务减少[26];流域西南部分地区为主要农耕区和建筑用地,人口活动较大,该地区景观破碎化程度的增大可能会使其生态系统服务增加,如退耕还林还草工程的实施,部分农田用地转化为森林或草地,其景观趋于破碎化的同时生态系统功能趋于稳定,生态系统服务得到相应增大。此外,景观格局与生态系统服务的相关程度还可能受到其他因素的影响,本研究中NPP生产服务与景观破碎化的空间关联性相对较弱,主要原因是NPP生产服务受众多因素影响,除景观格局外,温度、降雨、辐射以及地形地貌等因素均对NPP生产服务产生一定程度的影响。因此,景观格局对生态系统服务的影响是一个非常复杂的过程,两者之间并非简单的线性相关,其内在的相互作用机理还需进一步的探讨。

粒度是景观生态学研究的主要问题,景观格局指数随着粒度不同而发生改变[25]。过小的粒度容易忽略总体规律,过大的粒度容易忽略大量细节。本研究根据研究区范围大小,选择8 km×8 km为粒度对景观格局与生态系统服务二者的相关关系进行研究,之后还需深入探究不同粒度对二者关系的影响。

不同景观格局指数与生态系统服务关系呈现出复杂的关系。研究发现,景观指数斑块数(NP)、边缘密度(ED)以及香农多样性(SHDI)与生态系统服务相关性显著,而蔓延度指数(CONTAG)与生态系统服务的相关性较弱,表征景观指数CONTAG对生态系统服务的敏感程度较低。不同类型景观指数对生态系统服务敏感程度不同,对生态系统服务的影响程度也不同,还需进一步研究每种景观指数对生态系统服务的敏感程度,挖掘更具生态学意义的景观格局指数,从而深入探究景观格局对生态过程及功能的影响。

图6 汾河上游流域各景观格局指数与生态系统服务的局部空间自相关聚集图Fig.6 Local spatial autocorrelation cluster map of landscape pattern indices and three ecosystem services of the Upper Reach of the Fenhe River Watershed