基于GIS和Maxent模型的白唇鹿(Cervus albirostris)潜在适宜生境及保护GAP分析

唐中海,罗华林,王建华,刘家林,游章强,*

1 绵阳师范学院生态安全与保护四川省重点实验室,绵阳 621000

2 四川省甘孜州白玉县林业和草原局白玉县自然保护区服务中心,甘孜 627150

生物多样性保护与气候变化的应对是当今两大全球性热点环境问题[1]。气候被认为是在区域和全球尺度影响生物多样性的重要因子,气候变暖已经导致全球许多野生动物面临严重的生存危机[2]。据联合国政府间气候变化委员会(Inter-Government Panel on Climate Change,IPCC)评估分析,在20世纪,全球平均气温上升了0.85℃,并且还将持续变暖[3]。为适应全球气候变暖,许多野生动物向高海拔与或高纬度区域迁徙[4—5],而仅少部分物种通过调节自身物候停留在原分布区内[6]。通过模型模拟预测野生动物在全球气候变暖背景下未来潜在适宜生境分布格局与变化趋势,了解物种对全球气候变暖的适应对策,对制定生物多样性保护与管理策略具有非常重要的作用与意义[7]。

物种分布模型(Species Distribution Models SDMs)可以根据生物气候及其他环境因子模拟预测物种的潜在地理分布[8],因此在野生动植物生境分布研究中得到广泛应用。其中最大熵模型(maximum entropy modeling, Maxent)是基于环境变量及物种出现点位,通过机器学习和根据最大熵原理对物种地理分布进行模拟预测的物种分布模型[9]。由于该模型可以在物种分布位点和环境变量之间相关不明确的情况下模拟预测物种潜在适宜生境分布,并以概率的形式反映物种对生境的偏好程度等优点,该模型自推出至今,已经广泛应用于潜在适宜生境分布的模拟预测[10—12]。

保护空白分析(Protection gap analysis,GAP)是将物种的栖息地评估、潜在适宜生境范围与保护区网络或优先保护区域等空间信息进行交集分析,提供研究区域的物种组成、物种分布、当前的保护状态及保护空缺[13]。GAP分析强调每一物种在已有的保护区网络中都得到保护,而未出现的物种分布区就是空缺区域[14]。该方法广泛应用于物种保护和自然保护区建设规划。

白唇鹿(Cervusalbirostris),别名“白鼻鹿”、“黄鹿”、“岩鹿”、“红鹿”,是中国特有鹿种之一,仅分布于青藏高原及其周边区域的甘肃、青海、四川、云南和西藏境内3500—5500m海拔区域[15]。白唇鹿属国家I级重点保护野生动物,在中国生物多样性红色名录中将其列为濒危(EN)物种[16],被世界自然保护联盟(IUCN)列为易危(VU)物种[17]。由于白唇鹿种群数量较小,分布区域环境特殊,研究难度很大。目前对白唇鹿的研究主要集中在形态特征、繁殖、生境选择、食性及肠道菌群特征等[15—18]。基于生物气候(WorldClim data v2.1)对白唇鹿潜在适宜生境分布模拟预测并进行相应的保护评价还未见相关报道。本研究基于最大熵模型模拟预测在SSP1- 2.6和SSP5- 8.5两种气候变化情况下,白唇鹿在2040(2021—2040年)年和2100(2081—2100年)年潜在适宜生境,探讨全球气候变化对白唇鹿分布格局的影响,揭示白唇鹿对全球气候变暖的适应对策,并开展白唇鹿潜在适宜生境的保护评价,对在全球气候变暖背景下白唇鹿的管理与保护提供科学依据。本研究旨在解决以下问题：(1)分析影响白唇鹿潜在适宜生境分布的主要环境因子；(2)模拟预测在SSP1- 2.6和SSP5- 8.5两种气候变化情况下,白唇鹿在2040和2100年的潜在适宜生境分布格局及变化趋势；(3)分析全球气候变暖情况下白唇鹿潜在适宜生境的保护状态。通过本研究以期了解该物种对全球气候变暖的适应对策及潜在适宜生境的保护情况,为该物种的管理与保护提供科学依据。

1 数据来源及分析方法

1.1 研究区域概况

青藏高原被誉为“世界第三极”,是生物多样性的敏感区域[19],植被类型以高寒草甸和高寒草原为主[20]。近年来,研究发现气候变化和人类活动对青藏高原的影响越来越明显,随着全球气候变化的加剧,1960—2006年期间,青藏高原平均气温升高约0.03—0.05℃,年平均降水变化不显著[21]。由于青藏高原的生态脆弱性、对气候变化的敏感性以及不断加剧的人类活动,该地区面临严峻的退化演替形势[22]。祁连山位于青藏高原东北缘,横跨甘肃和青海两省,是青藏、内蒙古、黄土高原的分界线,在自然气候分区上起着非常重要的作用。祁连山脉位于季风区、非季风区和青藏高原区三大气候区的交汇处,降水量呈自东向西的递减趋势,山区多于平原,年均降水量30—600mm左右[23]。祁连山的中、东段属大陆性半干旱高山草原气候,西段属大陆性干旱荒漠气候。独特的地理环境和气候条件孕育了森林、草原、湿地、荒漠等生态系统。白唇鹿主要分布于青藏高原东缘及其毗邻区域和祁连山中西段区域。这两个区域都是全球气候变化敏感区域,开展潜在适宜生境模拟分析,能为在全球气候变化对生物多样性影响与日俱增的情况下,深入探究白唇鹿及其他野生动物对高原环境的适应提供重要理论依据。

1.2 白唇鹿分布数据

本研究共获得491个白唇鹿分布位点数据。其中来源于文献记录的分布点372个,文献资料包括《中国哺乳动物分布》[24]、《中国白唇鹿》[25],《中国哺乳动物多样性及地理分布》[26]；野外调查记录的分布点119个,为本研究组在2009—2020年期间在西藏工布江达、四川察青松多白唇鹿国家级自然保护区和甘肃祁连山国家级自然保护区的野外调查记录。根据Maxent模型对物种分布数据格式要求,将分布位点数据保存为CSV格式以满足Maxent模型运算需求。

1.3 环境数据获取

研究表明植被类型、人为干扰和地形因子是影响白唇鹿夏季生境选择的关键因子[27],因此本研究选择了包含以上环境因子和生物气候因子变量共26个,其中生物气候变量19个,地形变量3个,到河流距离变量1个,到湖泊距离变量1个,归一化植被指数变量1个,植被类型变量1个。

生物气候数据来自WorldClim数据库(WorldClim v2.1)(https://www.worldclim.org/)。生物气候数据2.1版中的社会共享经济路径(SSPs)比2.0版中浓度路径(representative concentration pathways, RCPs)的起始点更高,预测情景更加平缓,更接近未来气候变化真实值[28]。生物气候数据2.1版是基于5种社会共享经济路径(SSP1- 5)模拟了全球未来80年内气候变化情况,其中SSP1是模拟基于全球所有国家于2020年加入以可持续发展为目标的全球经济发展策略(SPAs: shared policy assumptions)条件下全球气候变化情况,该模式下社会经济发展对全球气候变化影响最小；SSP5则相反,是模拟全球社会经济发展模式未能在2020年实现以可持续发展为目标的情况下全球气候变化情况,该模式下社会经济发展对全球气候变化影响最明显。本研究中我们选择了社会经济发展模式对全球气候变化影响最小(SSP1)和最明显(SSP5)的生物气候数据(GCM: BCC-CSM2-MR, 分辨率为2.5min)用于模拟白唇鹿潜在适宜生境分布格局与变化。选择2040(2021—2040)年期间的模拟情况作为当前白唇鹿潜在适宜生境分布格局,选择2100年的模拟作为未来白唇鹿潜在适宜生境分布格局,旨在模拟两种不同社会共享经济模式下,比较分析在未来80年期间全球气候变化对白唇鹿潜在适宜生境的影响及保护状况。地形数据包括海拔(Elevation)、坡度(Slope)和坡向(Aspect),三个地形变量都是利用ArcGIS(v 10.8)从中国数字高程(DEM)提取。中国数字高程数据来自NASA (http://srtm.csi.cgiar.org/, SRTM30, version 2.1)。中国河流、湖泊和植被类型数据来自资源环境科学与数据中心(http://www.resdc.cn/)。人类足迹指数数据来自社会经济数据和应用中心(Socioeconomic Data and Applications Center, SEDAC)(https://sedac.ciesin.columbia.edu/)。中国国家生态屏障区数据来自(http://geodoi.ac.cn/),中国自然保护区数据来自(https://www.resdc.cn/)。

1.4 数据处理

利用ArcGIS(v 10.8)将所有环境数据的地理坐标系和投影坐标分别统一设置为GCS_WGS_1984和WGS_1984_UTM_Zone_48N,并以四川、甘肃、青海、西藏和云南行政区域范围(数据来源：http://www.resdc.cn/)为掩膜进行提取,并对提取后所有环境数据统一重采样(分辨率设置30″)。利用ArcPy脚本对每个生物气候数据图层19波段提取分析(将每个波段提取为单个数据图层),然后对所有环境图层进行多值到点数据提取,获得在SSP1- 2.6和SSP5- 8.5两种情况下2040和2100年期间4个模型(SSP1- 2.6 2040; SSP1- 2.6 2100; SSP5- 8.5 2040和SSP5.8.5 2100)中白唇鹿每个分布位点的环境变量值。利用Pearson相关分析和主成分分析分别对每个模型中白唇鹿分布位点的环境变量进行共线性分析。同一模型中,若两个环境变量的相关系数绝对值大于或等于0.8,则只保留特征向量值较大的变量,经筛选后获得4个模型运行所需要的环境变量共11个(表1)。

表1 环境环境因子筛选

1.5 预测模型参数设置及GAP分析

将白唇鹿分布位点数据(CSV格式)和环境数据(ASCII格式)导入Maxent模型(v 3.4.1)[29]构建白唇鹿分布模型。在程序主页面中勾选构建反应曲线(Create response curves),刀切法检测环境变化的重要性(Do Jacknife to measure variable importance)和逻辑斯蒂输出(Logistic)三个选项。在设置菜单中勾选随机种子选项,并将随机测试百分比设置为0.25,模型重复运行次数设置为10,重复运行类型选择Bootstrap,阀值规则选择最大训练敏感度和特异度法(Maximum training sensitivity plus specificity, MaxSS),其余参数选择默认设置。利用受试者工作特征曲线下面积(Area under the receiving operator curve, AUC)评价模型优劣,采用刀切法(Jacknife)检验环境因子的重要性。

选择测试AUC(TEST AUC)值最大的模型结果用于白唇鹿潜在适宜生境分布分析。将所选择模型结果的ASC格式进行栅格转化并定义坐标系统后(地理坐标统一设置为：GCS_WGS_1984和投影坐标统一设置为：WGS_1984_UTM_Zone_48N),利用ArcGIS重分类工具对白唇鹿潜在适宜生境进行重分类分析。将白唇鹿潜在适宜生境按0—0.25,0.25—0.50,0.50—0.75,0.75—1进行分类。在重分类过程中我们将0—0.25设为NoData,并定义不适宜生境类别,0.25—0.50定义为低适宜生境类别,0.50—0.75定义为中适宜生境类别,0.75—1.00定义为高适宜生境类别。利用ArcGIS空间分析工具计算低、中、高适宜等级潜在生境的面积与质心,并以低、中、高适宜生境图层为掩膜,提取中国行政范围的数字高程图层到适宜生境范围,计算适宜生境的平均分布海拔。

潜在适应生境保护空白(GAP)面积的计算采用如下方法。首先,利用ArcGIS(v 10.8)对中国国家生态屏障区数据和中国自然保护区数据(国家级自然保护区和地方级自然保护区)图层进行矢量化后分别计算面积。其次,将每个模型模拟的白唇鹿潜在适宜生境图层重分类后的数据进行矢量化后分别计算面积。然后,以每个模型模拟的白唇鹿低、中和高潜在适宜生境图层为“输入要素”,分别以中国国家生态屏障区数据和中国国家级自然保护区数据和地方级自然保护区数据为“输出要素”进行裁剪,裁剪结果图层的面积分别为三个等级潜在适宜生境范围在三个“输出要素类”图层中的面积。三个等级潜在适宜生境面积减去裁剪结果图层的面积即为白唇鹿潜在适宜生境的保护空白面积。

2 结果

2.1 模型检验评价与环境因子重要性分析

Maxent模型结果表明,模型的AUC平均值为0.95±0.01(Mean±SD)。根据Phillips等对模型自身评价依据[29],本研究中模型运行结果中的AUC值较为接近1,表明模型拟合程度较高。刀切法检验结果表明,在SSP1- 2.6和SSP5- 8.5气候变化情况下,植被类型和海拔都是影响白唇鹿在2040和2100年的潜在适宜生境分布的关键环境因子。在SSP1- 2.6气候变化情况下,在2040年(2021—2040年)期间,影响白唇鹿潜在适宜生境分布的关键生物气候因子是最暖季均温(Mean temperature of the warmest quarter, Bio_10)；在2100年(2081—2100年)期间,影响白唇鹿潜在适宜生境分布的关键气候因子是最热月最高气温(Max. temperature of the warmest month, Bio_5)。在SSP5- 8.5气候变化情况下,在2040年对白唇鹿潜在适宜生境影响最大的气候因子是最热月最高气温(Max. temperature of the warmest month, Bio_5)；在2100年对白唇鹿潜在适宜生境影响最大的气候因子是最暖季均温(Bio_10)。

2.2 潜在适宜生境分布格局与面积

模型模拟结果表明,在SSP1- 2.6气候变化情况下,从2040到2100年,白唇鹿潜在适宜生境分布格局变化较为明显的区域主要集中在青海中东部和西藏南部(图1),低、中和高适宜生境面积均在不同程度上有所下降,其中高适宜生境面积下降了26.65%(表2)。

在SSP5- 8.5气候变化情况下,从2040年和2100年,白唇鹿潜在适宜生境分布格局变化明显的区域主要集中在西藏南部和青海、四川与青海三省交界区域(图1),低、中和高适宜生境面积大小均在不同程度上有所下降,其中高适宜生境面积下降了16.75%(表2)。

图1 不同气候变化情况下白唇鹿总潜在适宜生境分布格局

2.3 潜在适宜生境质心与平均海拔

为了探讨白唇鹿对全球气候变暖的适应对策,分析了在未来80年内,两种气候变化情况下,白唇鹿潜在低、中、高适宜生境在水平和垂直方向的变化情况。结果表明,在两种共享社会经济模式影响下,白唇鹿低、中、高适宜潜在生境的质心向高纬度移动(图2),平均海拔变化有升有降(表2)。白唇鹿低适宜生境质心向西北移动了0.28km,平均海拔下降了236m；中适宜生境质心向西北移动了0.60km,平均海拔下降了181m；高适宜生境质心向西北移动了0.78km,平均海拔上升了59m(表2)。在SSP5- 8.5气候变化情况下,从2040到2100年,白唇鹿低适宜生境质心向东南移动了0.13km,平均海拔变化不明显；中适宜生境质心向东北移动了0.33km,平均海拔上升了54m；高适宜生境质心向东北移动了1.05km,平均海拔变化不明显(表2)。

图2 白唇鹿潜在适宜生境质心

表2 白唇鹿潜在适宜生境面积、平均海拔与GAP分析结果

2.4 潜在适宜生境GAP分析

为了评价在全球气候变暖情况下潜在适宜生境的保护状态,我们分别分析了白唇鹿潜在适宜生境位于现有保护地(中国生态屏障区、中国国家级和地方级自然保护区)的情况。总体而言,白唇鹿潜在适宜生境的保护空白面积比例较大(表2)。在SSP1- 2.6气候变化情况下,白唇鹿2040年的低、中、高潜在适宜生境位于中国国家生态屏障区之外的面积比例平均值为45.15%,位于国家级自然保护区之外的面积比例平均值为83.14%,位于地方级自然保护区外的面积比例平均值为95.60%；与2040年比较,白唇鹿2100年的低、中、高潜在适宜生境位于中国国家生态屏障区、国家级自然保护区和地方级自然保护区的面积比例都明显下降。在SSP5- 8.5气候变化情况下,白唇鹿2040年的低、中、高潜在适宜生境位于中国国家生态屏障区之外的面积比例平均值为53.92%,位于国家级自然保护区之外的面积比例平均值为85.25%,位于地方级自然保护区外面积比例平均值为95.92%；与2040年比较,白唇鹿2100年的低、中、高潜在适宜生境位于中国国家生态屏障区、国家级自然保护区和地方级自然保护区的面积比例略有上升。

3 讨论

3.1 影响潜在适宜生境的环境因素

本研究结果表明植被类型、海拔、最暖季均温和最热月最高温是影响白唇鹿潜在适宜生境的关键环境因子,这与游章强等[27]发现植被类型和海拔是影响白唇鹿夏季生境选择的关键环境因子的研究结果相似。研究表明植被类型是影响鹿科动物潜在适应生境的关键环境因素,李鑫泽等研究发现植被类型是影响马鹿(Cervuselaphus)、东北梅花鹿(Cervusnippon)和驼鹿(Alcesalces)潜在适宜生境的关键环境因素[30]。白唇鹿的喜食植物主要为禾本科和莎草科植物[15,27],这些食物资源主要分布在高寒草甸和高山灌丛生境,这两类生境有明显的海拔分布,因此,食物资源的特殊类型和分布格局导致植被类型和海拔成为影响白唇鹿潜在适宜生境的关键因子。温度是影响生物生存和分布的重要环境因子,其直接影响有机体的体温,体温的高低又决定了动物的新陈代谢过程的强度和特点以及有机体的生长和发育速度、繁殖、行为、数量和分布等[31]。在最适温度范围内,动物生存压力较小,而在最适温度之外动物的生存压力增加,故动物趋向于最适温度范围,而远离最高温度范围。最暖季均温体现了动物对最适温度的适应,而最热月最高温则体现了动物对最高温度的耐受性。因此最暖季均温和最热月最高温是影响白唇鹿潜在适宜生境的关键环境因子。

在所有生物气候因子中,最暖季均温和最热月最高温是影响白唇鹿潜在适宜生境的主要气候因子,这可能与青藏高原气候特点相关。青藏高原被誉为世界“第三极”,平均海拔4000m,青藏高原地区是近百年来全球增暖最显著的区域[32],大部分地区年均温度在0℃以下[33],最热月平均温度5.5—13.6℃[34]。因此,在20世纪,尽管全球平均气温上升了0.85℃[3],青藏高原冷季气温上升仍然不会超出其年平均温度,对该区域的野生动物影响不明显。然而,在暖季,尤其是在最热月期间,受全球气候变暖的影响,最高温继续升高将超出该区野生动物的高温耐受范围,从而威胁到该区域内的生物多样性。

3.2 潜在适宜生境分布格局与面积

向高海拔或高纬度迁徙是野生动物对全球气候变暖适应的主要策略[7],然而大多野生动物在向高海拔或高纬度区域迁徙过程中,都面临生境丧失的危机[35]。本研究中,在SSP1- 2.6和SSP5- 8.5两种气候变化情况下,白唇鹿潜在适宜生境的分布格局从2040到2100年均有明显的变化,且主要发生在三江(金沙江、澜沧江和怒江)并流区域以及青海、四川和甘肃交界区域,认为这可能是白唇鹿的迁徙受三江和人类活动干扰影响所致。三江不但河面宽且流量大,受地型陡峭的影响流速快[36],成为野生动物迁徙中不可逾越的障碍[37]。在适应全球气候变暖的迁徙中,受三江水体的阻碍,白唇鹿在三江区域及其邻近区域的迁徙将必然受阻,呈现出往返迁徙的现象。在青海、四川和甘肃交界带区域内,白唇鹿潜在适宜生境变化明显可能是受交通道路和城市分布相对集中,海拔相对较低,植被分布连续性较低等因素的影响。

在SSP1和SSP5两种气候变化情况下,从2040年到2100年期,白唇鹿潜在适宜生境面积下降的最高比例达到26.65%。白唇鹿潜在适宜生境面积下降可能与受全球气候变暖情况以及白唇鹿向高海拔和高纬度区域迁徙适应策略有关。研究表明,在全球山地区域中,青藏高原属“倒金字塔(Inverse pyramid)”山型,当海拔超过3000米以后,山体陆地表面积将急速下降[35]。如在察青松多保护区内,地形以高山峡谷为主[38],受该区域的山形特征影响,随着海拔上升,其陆地表面下降,白唇鹿可用生境面积随之下降。另一方面,在全球气候变暖背景下,动物的迁徙速度快于植被的扩散速度[39],白唇鹿向高海拔高纬度区域迁徙适应过程中,虽然能回避全球气候变暖的影响,但受植被扩散滞后的影响,适宜生境面积也将下降。

3.3 潜在适宜生境保护GAP分析

利用物种模型模拟预测物种潜在适宜生境分布现状,能为物种管理与保护提供科学依据,然而不能对物种潜在适宜生境保护状态进行评估。本研究结果表明在全球气候变暖影响下,白唇鹿潜在适宜生境处于现有保护地范围的比例仍然较低。在SSP1- 2.6和SSP5- 8.5气候变化情况下,白唇鹿潜在适宜生境位于现有保护地的范围中的比例都变化(表2),表明在全球气候变暖和人类活动的干扰影响的情况下,现有保护措施在白唇鹿适应全球气候变暖和人类活动干扰起到了积极的作用。然而,在SSP1- 2.6气候变化情况下,超过一半的潜在适宜生境位于中国国家生态屏障区内,尽管在SSP5- 8.5气候变化情况下略有下降。截至2018年,我国保护区面积已经覆盖了国土面积的18%,超过了2020爱知生物多样性的目标11所拟定的17%的覆盖目标[40]。然而,在SSP1- 2.6和SSP5- 8.5两种气候变化影响下,白唇鹿潜在适宜生境位于国家级自然保护区的面积比例都低于2020爱知生物多样性目标11所拟定17%。这表明,我国在生物多样性保护中,尽管在生态屏障区和保护地建设方面提前实现了爱知生物多样性目标2020中关于保护地建设目标,在某些物种的保护力度仍然还有待提高。