气候是决定物种地理分布的主要环境因子,物种分布格局的变化是对气候变化的直接反映[1]。人口膨胀、温室气体排放使气候变暖、植物栖息地破碎化,进而导致全球生态系统退化和物种丧失[2]。气候变暖对物种的生活史和分布会产生严重威胁[3],尤其是对一些濒危或分布范围狭窄的国家重点保护植物。全球变暖等气候变化过程将导致物种原生境无法继续提供稳定适宜的生存条件,物种为适应新环境条件不得不改变其原有的地理分布范围,若物种无法完成适应过程,可能会面临灭绝风险[4]。因此,运用生态位模型预测当前和未来不同气候情景下物种的潜在分布区变化及可能的灭绝风险,可以为濒危物种保护措施的制定及其可持续发展提供理论依据。
生态位模型利用物种已知的分布数据和环境变量,分析物种分布点与环境因子的相关性,探究物种的生态需求,并将运算结果投射到不同时间和空间中以预测物种的潜在适生区[5]。近年来,利用生态位模型预测未来气候情景下的物种分布范围变化已成为植物保护的研究热点[5-6]。常用的生态位模型较多,但各模型都具有一定的偏好性。Elith等[7]和Phillips等[8]运用广义线性模型、广义加性模型、生物气候分析方法、最大熵模型(MaxEnt)等16种模型预测物种适生区,结果表明这些模型均能较好地模拟物种的地理分布,但MaxEnt模型的预测准确性最高且稳定程度最优。MaxEnt模型对样本量没有严格要求,样本个数大于5就可以获得精确度较高的预测结果且不受分布信息缺失的约束,是预测物种潜在分布区的首选模型[9],被广泛用于预测植物适生区。
华南五针松(Pinus kwangtungensis)为松科(Pinaceae)松属(Pinus)常绿乔木,是中国特有的国家二级保护野生植物。相对于双维管束植物,华南五针松作为单维管束植物,是一种较原始、较古老的物种,对裸子植物系统分类学研究具有重要价值[10]。华南五针松喜温暖湿润气候,主要分布于广东、广西、贵州、湖南和海南等地,虽然其生态适应性较强,现代地理分布范围跨度大,但呈现出分布零散、数量少等特点[11]。目前,对华南五针松的研究多集中在种间关系[12]、种群动态[13]、群落特征[14]和生理特征[15]等方面,其在气候变化情景下在我国的地理分布变迁情况有待研究。本文运用地理信息系统(GIS)和MaxEnt模型,研究华南五针松地理分布对气候变量的响应关系,探究影响华南五针松地理分布的主导气候因子,揭示当前与未来气候情景下华南五针松适生区的变化趋势,以期确定华南五针松在中国的核心适生区,为华南五针松的保护管理及可持续利用提供科学依据。
1 数据与方法
1.1 物种分布数据
本研究所需华南五针松自然分布数据主要来源于实地调查、全球物种多样性信息库(http://www.gbif.org)和中国数字植物标本馆(http://www.cvh.org.cn)。共收集华南五针松种群在中国的自然分布记录263条,去除重复记录和无法精确到乡镇行政区采集地点的分布数据,消除预测中群集效应可能造成的误差,每2.5'×2.5'栅格内仅选择1个点[16],最终得到64个分布点用于后续研究(图1)。对每个分布点进行经纬度确认,并转化为CSV格式进行保存,用于MaxEnt模型的建立。
1.2 环境变量来源与筛选
研究所用当前时期(1970—2000年)和未来2个时期(2041—2060年、2061—2080年)的气候数据均来源于WorldClim数据库(http://www.worldclim.org),包括19个生物气候变量(Bio1~Bio19,表1),各因子的空间分辨率均为2.5'。未来模拟气候数据基于美国大气研究中心开发的CCSM 4模型,该模型对气候变量的模拟和预测准确性较高[17]。选择未来气候系统情景模式下的3种不同温室气体排放情景RCP 2.6、RCP 4.5和RCP 8.5。RCP是21世纪人类活动影响下气候变化预测模型的输入参数,用于描述未来人口、社会经济、科学技术、能源消耗和土地利用等方面发生变化时,温室气体、反应性气体、气溶胶的排放量以及大气成分的浓度。RCP 2.6表示在政府干预下温室气体排放浓度较低,总辐射强迫在2100年后稳定在2.6 W/m2;RCP 4.5表示总辐射强迫在2100年后稳定在4.5 W/m2;RCP 8.5为无气候变化政策干预下的高浓度温室气体排放情景[18]。从自然资源部标准地图服务系统(http://www.bzdt.ch.mnr.gov.cn)下载1∶2 000万标准中国地图,用于ArcGIS软件的分析底图,审图号为GS(2019)1822。
表1 本研究的环境变量Tab.1 Environmental variables in this study |
| 变量 | 变量描述 | 单位 |
|---|---|---|
| Bio1 | 年平均气温 | ℃ |
| Bio2 | 平均气温日较差 | ℃ |
| Bio3 | 等温性 | |
| Bio4 | 温度季节变化方差 | |
| Bio5 | 最暖月最高气温 | ℃ |
| Bio6 | 最冷月最低气温 | ℃ |
| Bio7 | 气温年较差 | ℃ |
| Bio8 | 最湿季度平均气温 | ℃ |
| Bio9 | 最干季度平均气温 | ℃ |
| Bio10 | 最暖季度平均气温 | ℃ |
| Bio11 | 最冷季度平均气温 | ℃ |
| Bio12 | 年均降水量 | mm |
| Bio13 | 最湿月降水量 | mm |
| Bio14 | 最干月降水量 | mm |
| Bio15 | 降水量季节变化方差 | |
| Bio16 | 最湿季度降水量 | mm |
| Bio17 | 最干季度降水量 | mm |
| Bio18 | 最暖季度降水量 | mm |
| Bio19 | 最冷季度降水量 | mm |
因最大熵模型的构建过程中存在自相关和多重线性重复等问题,会影响拟合结果的准确性,本研究利用R 3.6.3中的“corrplot”程序包对所有环境因子进行相关性分析,保留相关系数小于0.7的气候因子用于建模,当2个变量的相关性大于0.7时,仅保留对预测概率贡献较大的变量用于模型运算[19]。经相关性分析,最终筛选出年平均气温(Bio1)、平均气温日较差(Bio2)、等温性(Bio3)、温度季节变化方差(Bio4)、气温年较差(Bio7)、年均降水量(Bio12)、最干月降水量(Bio14)、降水量季节变化方差(Bio15)、最干季度降水量(Bio17)、最暖季度降水量(Bio18)等10个气候因子,并转化为ASC格式,用于华南五针松在中国的适生区预测。
1.3 物种分布模型的建立与预测结果准确性评价
将华南五针松的分布数据(CSV格式)和筛选后的10个气候因子(ASC格式)分别导入MaxEnt 3.4.1中,重复10次建模预测,将25%的华南五针松分布点设为测试集,剩下75%为训练集,采用刀切法检验权重,设置受试者操作特征曲线(receiver operating characteristic curve,ROC),模型结果以Logistic格式和ASC类型文件输出,其余参数采用模型默认设置。采用受试者操作特征曲线下面积值(area under curve,AUC)作为评价模型模拟结果精度的依据[20]。当AUC≤0.6时,表示模型预测失败;当0.6<AUC≤0.7时,表示预测结果较差,可信度低;当0.7<AUC≤0.8时,表示预测精度一般;当0.8<AUC≤0.9时,表示预测结果较准确;当0.9<AUC≤1.0时,表示预测结果非常准确,可信度高[21]。
1.4 适生区等级划分与环境因子重要性评估
通过MaxEnt模型输出华南五针松的潜在适生区,运用ArcGIS软件中的Re-classify重分类工具与自然间断点分级法将适宜生境划分为4个等级:非适生区(0~0.1)、低适生区(0.1~0.3)、中适生区(0.3~0.5)和高适生区(0.5~0.8)。计算不同时期各适生区的面积,比较华南五针松在不同气候变化情景下的适生区分布范围变化。依据MaxEnt模型的刀切法结果,筛选出主导气候因子,运用主导气候因子的响应曲线分析其对华南五针松适生区的影响。
2 结果
2.1 模型预测的精确性
将筛选后的10个气候因子导入MaxEnt模型,得到平均训练集与测试集的AUC值。在当前和未来不同气候情景下,平均AUC值均在0.995以上(图2),表明模型预测水平优秀,即本研究所构建的华南五针松适生区预测模型结果精确、可靠。
2.2 影响华南五针松分布的主导气候因子
MaxEnt模型运行结果表明,对华南五针松适生区分布影响较大的气候因子为最暖季度降水量和温度季节变化方差,二者的贡献率分别为60.2%和17.6%;置换重要值较高的依次为等温性、年均温、气温年较差,其重要值分别为47.22%、22.08%、17.44%。分析刀切法检验结果,发现仅使用单一环境变量时,正则化训练增益值、测试增益值和AUC值排在前3位的变量均为最暖季度降水量、年均降水量和年平均气温(表2)。影响华南五针松分布的气候因子主要为2种降水因子(最暖季度降水量、年均降水量)和3种气温因子(等温性、年平均气温、温度季节变化方差),其中最暖季度降水量是影响其分布的关键因子。
表2 主要气候因子参数Tab.2 The parameters of main climatic factors |
| 气候因子 | 贡献率/% | 置换重要值/% | 正则化训练增益值 | 测试增益值 | AUC值 |
|---|---|---|---|---|---|
| 最暖季度降水量 | 60.2 | 2.54 | 2.60 | 2.73 | 0.97 |
| 温度季节变化方差 | 17.6 | 0.80 | 1.25 | 1.36 | 0.91 |
| 降水量季节变化方差 | 6.7 | 4.11 | 1.02 | 1.14 | 0.88 |
| 最干月降水量 | 6.7 | 0.17 | 1.47 | 1.64 | 0.93 |
| 平均气温日较差 | 5.3 | 2.04 | 1.05 | 1.27 | 0.90 |
| 年均降水量 | 1.4 | 0.99 | 2.12 | 2.27 | 0.96 |
| 年平均气温 | 0.9 | 22.08 | 1.90 | 2.03 | 0.95 |
| 气温年较差 | 0.6 | 17.44 | 1.18 | 1.36 | 0.91 |
| 等温性 | 0.5 | 47.22 | 1.09 | 1.20 | 0.89 |
| 最干季度降水量 | 0.1 | 2.61 | 1.40 | 1.55 | 0.92 |
为找到华南五针松的适生区间(即响应曲线存在概率大于0.5的区间),将5个主导气候因子运用MaxEnt进行建模,得到主导气候因子与华南五针松的响应曲线图,如图3所示。
对于降水因子,最暖季度降水量(Bio18)为300~610 mm时,华南五针松的存在概率呈上升趋势;最暖季度降水量为610 mm时,华南五针松的存在概率最大,约为66%,此后存在概率逐渐降低;华南五针松的最暖季度降水量最适分布区间为560~820 mm。年均降水量(Bio12)在700~1 700 mm时,华南五针松的存在概率呈上升趋势;年均降水量为1 700 mm时,华南五针松的存在概率最大,约为70%;年均降水量大于1 700 mm时,华南五针松的存在概率逐渐降低,其年均降水量的最适分布区间为1 400~1 950 mm。
对于气温因子,温度季节变化方差(Bio4)为570~750时,华南五针松的存在概率均大于50%;温度季节变化方差为700时,华南五针松的存在概率最大。等温性(Bio3)数值为27~33时,华南五针松的存在概率均大于50%;等温性数值为29时,华南五针松的存在概率最大。年平均气温(Bio1)为14~20 ℃时,华南五针松的存在概率均大于50%;年平均气温为16 ℃时,华南五针松的存在概率最大。
2.3 当前时期华南五针松的适生区模拟
由图4可知,华南五针松的适生区主要集中在中国南部,总适生区面积为162.27万km2,约占中国国土面积的17%。其中,高适生区面积为35.36万km2,占中国国土面积的3.6%,主要位于广西、广东北部、福建东北部、贵州与广西交界处,另外在重庆与湖南交界处以及贵州与四川交界处也有少量高适生区。中适生区面积为61.48万km2,占中国国土面积的6.3%,集中在贵州东北部与南部、重庆、湖南西部、广东南部、福建西部与南部、江西东部。低适生区面积为65.43万km2,占中国国土面积的6.7%,主要分布在台湾、海南、湖南东部、湖北南部、浙江,另外安徽南部、上海南部、云南、四川等地区也有少量低适生区。
2.4 未来气候情景下华南五针松的适生区模拟
不同二氧化碳浓度气候情景下,未来不同时期预测的华南五针松总适生区面积和高适生区面积有显著差异(图5)。在RCP 2.6气候情景下,2061—2080年的总适生区面积与高适生区面积均显著大于2041—2060年(P<0.01),2041—2060年的总适生区面积及2061—2080年的高适生区面积显著大于当前(P<0.01)。在RCP 4.5气候情景下,当前的总适生区面积显著大于2061—2080年预测的总适生区面积(P<0.01),但当前高适生区面积显著低于2041—2060年预测的高适生区面积(P<0.05)。在RCP 8.5气候情景下,2041—2060年的总适生区面积与高适生区面积均显著大于当前及2061—2080年(P<0.01),同时,当前的总适生区面积与高适生区面积明显大于2061—2080年。综上,随二氧化碳浓度气候情景的变化,不同时期预测的华南五针松总适生区面积与高适生区面积出现明显波动,表明未来不同二氧化碳浓度可能会影响华南五针松的适生区范围。
与当前气候相比,未来3种气候情景下,华南五针松的适生区表现出不同的变化格局。分析华南五针松在2041—2060年的适生区发现:在RCP 2.6气候情景下(图6a),高适生区增加区域主要集中在广西西南部、广东西南部与东北部,重庆东南部从高适生区零星分布变为中适生区;在RCP 4.5气候情景下(图6b),高适生区仅广东西南部有所增加,中适生区面积大幅降低;在RCP 8.5气候情景下(图6c),高适生区增加集中在四川东南部、重庆西南与东南部、贵州北部等区域,而广西南部、广东东北与西南部均从高适生区转为中适生区。分析华南五针松在2061—2080年的适生区发现:在RCP 2.6气候情景下(图6d),高适生区增加区域为重庆东南部,高适生区减少区域为福建、四川地区,而广东与广西南部由高适生区破碎化分布转为中适生区;在RCP 4.5(图6e)与RCP 8.5(图6f)气候情景下,虽然部分地区的适生区面积有所增加,但增加面积小于减少面积,致使适生区面积整体呈下降趋势。适生区面积减少区域主要集中在广东西南部和广西南部,增加区域集中在重庆东南部和贵州东南部。值得注意的是,广东与广西南部之前一直是中、高适生区,但在2061—2080年的3种不同大气二氧化碳浓度下均有低适生区出现。总体来看,与当前气候条件相比,未来2个时期华南五针松适生区的扩张与缩减趋势较为相似,其高适生区都有不同程度的“北移”。
由华南五针松适生区的空间变换格局(图7)可以看出,在未来2个时期的3种气候情景下,华南五针松适生区的集中性增强,破碎化程度降低。在2041—2060年RCP 4.5及2061—2080年RCP 4.5、RCP 8.5情景下,适生区变化情况较为相似,其丧失生境均集中在台湾、湖南与湖北中部以及安徽与江苏部分区域,增加生境仅有四川北部地区。在2041—2060年RCP 2.6情景下,华南五针松的分布丧失区主要在台湾、云南东部、重庆北部、拉萨北部地区,适宜生境增加区域为湖北中部、海南西部以及浙江、安徽、江苏一带。在2041—2060年RCP 8.5与2061—2080年RCP 2.6情景下,华南五针松的增加生境都集中在云南与贵州交界处、重庆北部、四川北部与中部以及浙江、安徽、江苏部分地区,而台湾、湖南中部、湖北中部地区的部分生境将丢失。总体来看,在未来气候情景下华南五针松的分布丧失区主要集中在湖南、湖北、台湾、浙江和江苏等省份,而增加区主要出现在云南、贵州、广西三省交界处以及四川中部和大巴山部分地区。
2.5 华南五针松适生区的质心迁移方向
当前时期华南五针松适生区的质心在湖南省邵东市。在RCP 2.6气候情景下,华南五针松质心在2041—2060年向西北迁移至湖南洞口县,在2061—2080年又向东部迁移至湖南邵东市。在RCP 4.5气候情景下,华南五针松质心在2041—2060年向西部迁移至贵州三穗县,到2061—2080年又回迁至湖南祁东县。在RCP 8.5气候情景下,华南五针松质心在2041—2060年主要向北迁移,在2061—2080年又向西南方向迁移至湖南绥宁县。由此可知,华南五针松在不同气候情景下其质心迁移方向具有一定差异,但总体呈现2041—2060年向西北地区迁移、2061—2080年又向东部地区迁移的趋势(图8)。
3 讨论
3.1 制约华南五针松适生区的主要气候因子
植物地理分布是环境因子综合作用的结果,温度和降水对植物的生长和分布起主导作用[22],但二者的作用在不同物种间有较大差异。例如,胡忠俊等[23]运用MaxEnt模型对紫花针茅(Stipa purpurea)的地理分布进行研究,结果表明温度对其分布的影响比降水量大。对红松(Pinus koraiensis)的潜在地理分布进行研究,发现水分对红松种群分布的影响高于其他因子[24]。本研究采用5种评估途径对10个气候因子进行位次排序,结果显示最暖季度降水量、等温性、年平均气温、温度季节变化方差、年均降水量是影响华南五针松分布的主要环境因子。研究结果表明,降水和温度对华南五针松的地理分布均有影响,其中最暖季度降水量是影响华南五针松分布的关键因子。韩洁[25]和戴旻峻[26]分别对巴山松(Pinus henryi)和柳杉(Cryptomeria japonica var. sinensis)的种群分布进行研究,结果也发现最暖季度降水量是影响2个种群分布的主要气候因子。
植物对生长环境的需求在一定程度上决定了物种的分布范围,有研究表明华南五针松喜好温暖湿润的生境条件[11]。黄蕴凯等[27]研究发现降水量与华南五针松径向生长显著正相关,表明华南五针松生长需要充足的水分。夏季高温会加快水分蒸发,导致土壤和树体内的含水量降低,从而抑制树木生长,此时较多的降水会缓解植物水分缺失压力[28-29],因此华南五针松分布主要受最暖季度降水量的影响,而中国的最暖季度降水量较高地区为南方地区[30],这与本研究中华南五针松适生区的研究结果一致。陶翠等[11]研究认为温度因子对华南五针松的影响大于降水因子。这与本研究结果存在差异,其可能原因是所选取的环境因子不同,本研究选取了19个气候因子,相对更为全面;此外,分析方法也有差异,前者运用主成分分析来确定关键气候因子,本文运用当下预测效果最好的MaxEnt模型对华南五针松在不同时期的潜在分布区进行模拟预测,预测结果更直观、可靠[7]。
3.2 华南五针松在未来气候变化下的分布格局与质心迁移
华南五针松适生区质心在2041—2060年不同气候情景下均向西北方向迁移,在2041—2060年RCP 8.5 情景下其质心北移现象尤其明显。政府间气候变化专门委员会发布的报告显示,由于温室气体的排放,全球气候变暖加速,到21世纪末在RCP 2.6气候情景下的气温上升不超过2 ℃,而在RCP 8.5情景下全球温度上升可达4.8 ℃[34]。Wu等[35]根据区域气候模式,对未来温度带动态进行了预估,提出未来温度带普遍北移,最显著的将是亚热带、暖温带和高原温带北界。在对长苞铁杉(Nothotsuga longibracteata)和桫椤(Alsophila spinulosa)的潜在适生区进行研究中也发现,2个物种的适生区在未来气候模式下均有不同程度的北移[36-37]。在2041—2060年的不同大气二氧化碳浓度排放情景下,随着二氧化碳浓度的增加,华南五针松高适生区面积呈先减少后增加的趋势,在RCP 8.5排放情景下,其增加幅度最大,总适生区面积约增加4.4万km2,适生区增加区域集中在贵州南部、重庆南部和四川东南部。高浓度二氧化碳排放情景下的降水量较高,能缓解降水因子对物种分布的限制,使适生区得以扩张[38],这可能是2041—2060年RCP 8.5排放情景下华南五针松适生区面积增加的主要原因。然而,华南五针松在2061—2080年RCP 8.5情景下的适生区面积缩减,这可能是由于随着年份的增长,未来全球温度的持续升高导致环境过度干旱,可供物种吸收的有效水分减少,增加的降水量不能降低或缓解降水因子对物种分布的限制[39],因此该情景下华南五针松适生区面积减少。
根据本文模拟结果,未来2个时期的华南五针松适生区有部分增加,但物种的实际分布不仅与降水和温度有关,还受紫外线、生物、人为活动和物种扩散等因素的制约[40],物种适生区实际增加面积未必能达到模型预测值[41],且龚粤宁等[42]研究发现目前华南五针松的自然演替处于衰退之中,加之全球气候变暖将加剧华南五针松资源的减少,为保护这一珍稀树种,建议将广西、广东西北部、贵州南部边缘等现有核心分布区作为华南五针松应对未来气候变化的重点保护区域。另一方面,华南五针松在广西南部边缘、广东南部与东南部地区的高度适宜区在未来2个时段显著减少,揭示了该种在这些区域的适生区变化规律,能为气候变化背景下华南五针松的保护提供科学依据,也可为准确预测珍稀植物的地理分布范围提供方法参考。但本研究仍有不足,在今后研究中可综合地形、植被、海拔和人为活动等复合因子,提高预测的准确性,以更准确地揭示华南五针松的适生区分布格局,为其可持续保护管理提供更科学的依据。
4 结论
华南五针松在当前气候条件下适生区主要集中在中国南部,涵盖广西、广东北部、福建东北部、贵州与广西交界处等地区。最暖季度降水量、年均降水量、温度季节变化方差、等温性和年均温是影响华南五针松分布的主要环境因子,华南五针松在最暖季度降水量560~820 mm、年均降水量1 400~1 950 mm、温度季节变化方差570~750、等温性27~33、年均温14~20 ℃的地区出现生存概率。其中,最暖季度降水量是影响华南五针松分布的关键因子。未来,气温升高及降水增加到一定程度时,有利于华南五针松原有适生区的保留及新适生区的扩增,尤其在2041—2060年RCP 8.5气候情景下,华南五针松的扩增区域面积达到最大,重庆北部、四川北部与中部、云南与贵州交界处将成为华南五针松的新增适生区。但RCP 8.5气候情景并不利于华南五针松的长期生长,在2061—2080年RCP 8.5情景下,华南五针松的丧失生境增加。不同温室气体情境下,华南五针松质心整体呈现2041—2060年向西北地区迁移,2061—2080年又向东部地区迁移的趋势。