不同气候环境的我国白桦潜在分布区预测1)

2023-02-22高明龙萨如拉铁牛张丽丽

东北林业大学学报 2023年2期

高明龙萨如拉铁牛张丽丽

(内蒙古农业大学，呼和浩特，010019)(内蒙古自治区林业科学研究院)(内蒙古阿尔山森林工业有限公司

全球气候变暖是当今世界面临的主要气候问题之一，从上世纪初至今全球增温1.2 ℃[1]。根据IPCC第六次气候变化评估报告，在本世纪末期，全球平均温将继续增加1.0～5.7 ℃[2]。全球气候变暖将对自然生态系统造成严重的影响[3-4]。研究表明，植物物种分布主要受到气候变化的影响，气候的持续变暖将使温带森林植被向高纬度和高纬度地区迁移[5-6]。近年来，物种分布模型逐渐成为分析气候变化与物种分布关系的一个研究热点[7]。物种分布模型依据生态位理论，通过物种实际地理分布数据和环境因子对该物种的潜在分布区进行模拟[8]。目前，主要的物种分布模型有最大熵模型(MaxEnt)、分类与回归树模型(CTA)、广义线性模型(GLM)、生境模型(HABITAT)、遗传算法模型(GARP)和生物气候模型(BIOCLIM)等模型。其中，最大熵模型(MaxEnt)相较其他同类模型具有预测结果易读性高、准确度高和灵活性强的优点[9]。因此，最大熵模型在酸枣(Ziziphusjujubavar.spinosa)、沙冬青(Ammopiptanthusmongolicus)和祁连圆柏(Juniperusprzewalskii)等众多物种的适生区预测中已有较为广泛应用[10-12]。白桦(BetulaplatyphyllaSuk.)为我国常见桦木属(Betula)落叶乔木，是我国温带地区主要次生林先锋树种之一，白桦形成的白桦次生林在保育土壤、涵养水源及固碳释氧方面具有较高的生态价值，其木材广泛用于建筑、器具和地板等日常生产中，同时白桦汁可用于多种药品和保健饮品的加工生产，具有较高的社会、生态及经济价值[13-15]；但以白桦为主体的白桦次生林在林分结构、林木生长、林地生境和生态功能等多个方面均与当地原始林分有较大差异，多表现为生物多样性、稳定性和抗逆性等不同程度的减弱[16]。目前，对白桦的研究多为林分生产力特征、生理生化特点和遗传多样性分析等，而关于白桦潜在分布区的研究尚未见报道[17-19]。因此，本研究通过优化后的最大熵模型模拟我国白桦不同气候环境下的潜在分布区变化，分析影响白桦分布的主要限制因子，对我国未来白桦资源的进一步开发与保护及低效白桦次生林改造具有重要理论指导意义。

1 材料与方法

白桦地理分布数据的收集：在2017—2021年，通过野外实地调查，在内蒙古东部和黑龙江北部的白桦分布区域共获取99个白桦分布点数据；再结合国家标本馆(NSII)、中国数字标本馆(CVH)、中科院昆明植物研究所标本馆(KUN)及其他相关文献公开资料，获取白桦101个分布点数据。为避免分布点密集所导致的模拟结果出现偏差，本研究在每10 km网格中仅保留1个白桦分布点[12]，最终得到166个分布点(见图1)。

图1 中国白桦地理分布点

1.1 环境因子及筛选

本研究选用的34个环境因子是空间分辨率均为30″的19个生物气候因子、14个土壤因子和1个地形因子。气候数据来源于世界气候数据库(www.worldclim.org)，包括末次间冰期、末次冰盛期、全新世中期、当前(2 000 s)、2050 s和2090 s共6个时期的最暖季平均气温、月平均昼夜气温差、最热月最高气温、最湿月降水量、年均气温、最湿季平均气温、最干季平均气温、气温变异系数、气温年较差、年降水量、最暖季降水量、等温性、最冷季平均气温、最干季降水量、最冷季降水量、最冷月最低气温、降水量变异系数、最湿季降水量和最干月降水量。土壤和地形数据来源于联合国粮农组织世界土壤数据库(HWSD，https://www.fao.org/soils-portal/soil-survey/soil-maps-and-databases/harmonized-world-soil-database-v12/en/)，包括表层土砾石体积分数、阳离子交换量、土壤质地类型、砂砾质量分数、黏粒质量分数、粉砂质量分数、黏性成分阳离子交换量、土壤密度、pH值、交换性盐基总和、盐基饱和度、可交换钠盐百分比、电导率、有机碳质量分数和海拔。

采用社会经济共享浓度(Shared Socio-economic Pathways,SSP)中提供的过去、当前和未来气候数据。其中以第二代国家(北京)气候中心中等分辨率气候系统模式(BCC-CSM2-MR)作为未来气候系统模式，在未来土壤和地形因子无显著变化的前提下，选择低水平温室气体排放量(SSP126)、中水平温室气体排放量(SSP245)和高水平温室气体排放量(SSP585)三个环境对白桦未来适宜分布区进行预测[20]。

为避免各因子间相似性过高导致模型过度拟合，通过方差膨胀因子(VIF)和皮尔逊(Person)相关性分析，结合MaxEnt模型预实验中各个因子贡献率，选择相关系数小于0.8且VIF值小于10的因子中对白桦分布具有较高贡献率和生态学意义的因子[21]。最终，保留5个气候因子、3个土壤因子和1个地形因子进行建模(表1)。

表1 参与MaxEnt模型运算的环境因子

1.2 模型建立

基于爪哇(Java)语言的MaxEnt模型是一种结合机器学习和统计模型模拟预测物种分布概率的自学习模型，将白桦166个分布点数据和9个环境因子数据导入MaxEnt3.4.1软件，随机选取75%的样点作为机器学习的训练数据集，其余25%的样点作为数据验证集，重复模式为自举法(Bootstrap)，进行10次以独立重复建模，建模结果以逻辑斯蒂(Logistic)形式输出，按照各环境因子贡献率和置换重要值确定主导环境因子。

根据Jackknife法分析评判各环境因子的贡献率，MaxEnt模型准确度可通过受试者工作特征曲线(ROC)进行检验，利用曲线下面积(AUC)评估所建模型精度，AUC值在0.5～1.0，其值越大预测越精确。其中，AUC值在0.5～0.7时表示预测结果较差，0.7～0.8时表示预测结果一般，0.8～0.9时表示预测结果较好，0.9～1.0时表示预测结果非常好[22]。

1.3 模型优化

本研究模型优化方法采用棋盘格2法(Checkerboard2)，调整调控倍频(RM)和特征组合(FC)2个参数改变模型正则化水平。通过改变特征数量和组合方式以及RM的数值，将MaxEnt模型中的5种特征：线性特征(L)、二次型特征(Q)、片段化特征(H)、乘积型特征(P)和阈值性特征(T)，与16个RM值(0.5～8.0，以0.5为间隔梯度)排列组合出148种组合方式。利用R语言中的ENMeval程序包对其进行优化测试，结果中各组合更正的赤池信息准则值(delta.AICc)和10%测试遗漏率越低，表示模型预测精准度越高[23-24]。

1.4 数据处理

运用地理信息系统10.4.1软件(ArcGIS 10.4.1)将模型运行后生成的数据进行可视化，根据MaxEnt模型生成的白桦适宜性阈值，对自然断点法划分的适宜性等级进行修正，将白桦生境适宜性划分为4个等级，依次为不适生区(0～0.1]、较不适生区(0.1～0.3]、一般适生区(0.3～0.5]和高度适生区(0.5～1.0][25]。在ArcGIS软件中，根据上述划分的结果，对比不同时期白桦的适生区和非适生区的地理空间变化，计算并绘制各时期气候变化背景下白桦未来空间分布格局变化图；在ArcGIS软件中，使用物种分布模型工具箱2.0(SDMtoolbox2.0)工具包计算白桦在当前和未来2个时期3种气候环境下的适生区质心位置及其随时间迁移方向，并计算其质心迁移距离。

2 结果与分析

2.1 模型优化结果与评价

根据166个分布点和9个环境因子使用MaxEnt模型对白桦的潜在分布区域进行模拟预测。对比148类参数组合方式，当MaxEnt模型参数设定调控倍频(RM)为2，特征组合(FC)为片段化特征(H)、乘积型特征(P)和阈值性特征(T)时，赤池信息准则值(delta.AICc)为0，10%训练遗漏率相较默认值降低18.5%(表2)。因此，选取该参数作为模型最终参数，使用该参数的10次模拟训练AUC均值为0.933(图2)，说明预测结果准确度非常好。

表2 不同参数MaxEnt数模型优化结果

图2 MaxEnt模型的ROC响应曲线

2.2 当前时期我国白桦分布

通过MaxEnt模型推算出中国白桦当前时期的适生分布区(图3)，当前总适生区面积(一般适生区面积与高度适生区面积之和)为1 418 524.3 km2，占我国陆地面积约14.7%，主要分布于大兴安岭地区、小兴安岭地区、长白山地区、祁连山地区以及横断山脉地区，在河北、陕西、山西、甘肃和宁夏等省份也有小面积分布。高度适生区占总适生区面积约43.5%，其主要位于内蒙古东北部、黑龙江西北部和东南部、吉林省东北部。在河北北部和中部、山西中部、甘肃南部、青海东部、四川西部、西藏东南部、云南西北部等地区呈破碎状分布。

图3 当前我国白桦潜在分布区

2.3 过去时期白桦地理分布

在过去3个时期中，白桦在东北大兴安岭地区的分布状况最为稳定，与当前白桦适生区分布基本一致。西南方向白桦适生区在3个时期内均向东发生不同程度的收缩。在末次间冰期，白桦在我国具有较为广泛的分布，其适生区主要位于我国青藏高原和东北大兴安岭地区，其中高度适生区占总适生区面积的67.3%(图4，表3)。在末次冰盛期，白桦总适生区面积相较末次间冰期仅有0.1%的小幅度增加，但总适生区在我国分布位置出现较大变化，在东北大兴安岭、小兴安岭、长白山地区具有大面积的新增适生区出现，同时西藏西部地区的适生区大面积消失(图5)，仅在青海、西藏、四川3省交界处有较大面积留存。在全新世中期，白桦总适生区面积缩减较为明显，仅为末次冰盛期面积的66.8%，分布范围与当前时期最为相似(图4)，青藏高原大部分地区已不适合白桦生长，西南地区适生区开始呈现出破碎化趋势，东北地区白桦适生区与当前气候条件下基本一致。

2.4 未来气候环境下白桦地理分布

根据表3，在未来时期3类气候条件下，白桦总适生区面积均为缩小趋势，2050 s和2090 s总适生区面积相较前一时期均有不同程度的缩减。在SSP126未来气候环境下，总适生区面积损失最小，2050 s总适生区面积相较当前总适生区缩小23.0%，约合面积为325 486 km2；2090 s总适生区面积相较2050 s总适生区小幅度缩小0.8%，损失面积仅为9 201 km2。在SSP245环境下，2050 s总适生区面积相较当前总适生区缩小23.7%，约合面积为336 076 km2；2090 s总适生区面积相较2050 s缩小15.8%，损失面积约为170 711.9 km2。在SSP585环境下，总适生区面积损失最大，2050 s总适生区面积相较当前总适生区缩小32.9%，约合面积为467 239.3 km2；2090 s总适生区面积相较2050 s缩小37.4%，约合面积为355 746.5 km2，此时总适生区缩减比例最高。

表3 不同时期白桦适生区空间变化

通过未来各时期适生区变化可以看出，SSP245和SSP585环境下白桦适生区对气候变化的响应最敏感(图4)。在SSP245环境下，2050 s高度适生区面积相较当前减少33.9%，面积约为208 854.2 km2；在2090 s白桦高度适生区面积的变化幅度最大，高度适生区面积相较当前减少64.7%，面积约为399 323 km2。在SSP126环境下，2050 s高度适生区面积的变化幅度最小，新增加适生区增加率最低，高度适生区面积相较当前减少28.3%，面积约为174 635.5 km2；2090 s高度适生区面积将持续减小，相较当前时期减少34.9%，减少面积约为215 312.5 km2，但此时期新增适生区面积的增加率最高。在SSP585环境下，2050 s高度适生区面积相较当前时期减少55.0%，面积约为339 444.5 km2；2090 s高度适生区面积为当前时期适生区面积的36.3%，减少面积约为393 020.9 km2。

LIG表示为末次间冰期，LGM表示为末次冰盛期，MH表示为全新世中期，Current表示为当前时期。

2.5 质心迁移趋势

在空间格局方面，过去时期白桦适生区的质心总体向东北迁移，而未来不同气候环境下质心总体向西南迁移，迁移距离随气候变化加剧有进一步扩大的趋势(图6)。当前时期白桦适生区的质心在河北省张家口市沽源县(115°33′36″E，41°35′24″N)。当气候环境为SSP126-2090s时，适生区质心向西迁移至内蒙古自治区乌兰察布市商都县(113°58′12″E，41°40′12″N)，迁移距离为133 121 m；当气候环境为SSP585-2090s时，质心向西南迁移至内蒙古自治区巴彦淖尔市乌拉特中旗(109°38′24″E，40°38′24″N)，距离为506 837 m。

图5 不同时期白桦空间分布格局变化

图6 不同气候变化下白桦适生区质心变化

2.6 环境因子分析

根据Jackknife法分析结果(表4)，对白桦潜在地理分布贡献率超过5%的因子分别为最暖季度降水量、最热月最高气温、年均温、年降水量、最冷月最低气温和海拔,这6个环境因子的总贡献率高达89.8%，同时这6个环境因子的总置换重要值为94.4%。

表4 参与建模的环境因子贡献率及置换重要值

2.7 未来分布区生态特征变化

影响白桦分布的主要环境因子与白桦生境适宜度关系见表5。在SSP126气候环境下，2050 s和2090 s白桦生境适宜度分别比当前降低了0.07和0.09；在SSP245气候环境下，2050 s和2090 s白桦生境适宜度分别比当前低0.09和0.15；在SSP585气候环境下，2050 s白桦生境适宜度比当前低0.13，2090 s白桦生境适宜度相较当前大幅降低47.2%，仅为0.28。

表5 不同模式情境下白桦适生区环境因子与白桦生境适宜度的变化

166个白桦分布点的年平均气温与白桦生境适宜度的变化呈相反趋势。在SSP126、SSP245和SSP585气候环境下，2050 s的年平均气温相较当前时期分别增加158.5%、205.7%和264.2%；2090 s相较2050 s年平均气温分别增加6.9%、29.3%和62.7%。年降水量与年平均气温相似，均与白桦生境适宜度变化相反。在SSP126、SSP245和SSP585环境下，2050 s的年降水量相较当前分别增加14.67、23.38和41.80 mm；2090 s与2050 s相比，年降水量分别增加19.97、23.91和47.9 mm。

3 讨论

3.1 优化后模型对白桦适生区预测

优化后的MaxEnt模型预测结果与白桦实际地理分布高度一致，表明经过优化后的MaxEnt模型用于白桦的分布区预测时，预测结果较为准确、可靠。当前白桦主要分布于内蒙古东部、黑龙江、吉林、河北、河南、陕西、宁夏、甘肃、青海、四川、云南和西藏东南部等地区[26]，通过对比分析各参数组合预测结果和实际分布区域的拟合程度以及对地理预测图进行视觉检查来判断模型预测精度[27]，优化后的MaxEnt模型对白桦适生区的预测结果也恰好位于上述地区。

3.2 环境因子对白桦地理分布制约

水热因子是影响白桦分布的决定性因子。在本研究中，白桦高度适生区的平均最热月最高气温为22.7 ℃和最冷月最低气温为-23.3 ℃，与我国白桦连续分布面积最大的大兴安岭北部地区最热月最高气温16.0～17.9 ℃和最冷月最低气温-25.4～-30.0 ℃相比较略高[28]，这是由于本研究在分析白桦最适环境因子时将各地白桦高度适生区对应环境因子数值按分布区面积进行加权平均，所以导致整体平均值有所提高。同时，降水作为影响白桦地理分布的主要环境因子之一，白桦高度适生区最暖季度降水量占年降水量的69.1%，这说明白桦喜生于夏季较为湿润的立地环境，但我国东南沿海地区白桦分布较少的现象也说明过量降水也会对白桦的分布存在一定的限制。

本研究白桦分布高度适生区平均海拔为1 537.5 m，与先前研究得出白桦主要分布在400～4 100 m的结论存在一定差异[26]，这是因为本研究中白桦在东北海拔1 000 m以下地区和华北海拔1 300～1 900 m地区的高度适生区面积占比最大，而西南高海拔地区白桦高度适生区呈小面积破碎状分布，这种分布格局导致白桦分布高度适生区平均海拔相较先前研究有所降低[29-30]。参与建模的土壤因子对白桦分布的影响相对较小，说明白桦对土壤不敏感，这也是白桦在我国有较为广泛分布的一个重要原因。

3.3 各时期气候变化对中国白桦分布影响

通过对白桦历史时期分布区模拟表明，在末次间冰期，白桦作为一种耐寒树种，广泛分布于此时较为寒冷的青藏高原和大兴安岭地区；在末次盛冰期，青藏高原地区的白桦适生区向东较暖地区进行了大幅度迁移；在全新世中期，由于气候与当前时期类似，白桦适生区也呈现出与当前适生区相似的分布情况[31-32]。古植被研究表明，在末次间冰期，温带树种在北半球广泛分布，具有良好的生态适应性；到末次冰盛期，其在高纬度地区分布范围发生收缩；在全新世中期，温带树种适宜生境有一定程度的扩大[33-35]。白桦过去时期适生区变化基本印证了这一结论。在3个历史时期中，大兴安岭地区均保留有大面积白桦适生区，因此本研究推断白桦在冰川期的气候避难所为大兴安岭地区。

在未来SSP126、SSP245和SSP585气候变化下，白桦适生区面积均有不同程度的缩小，并且对SSP245和SSP585气候变化环境下的影响最敏感。白桦在我国东北适生区分布范围向北收缩尤为明显，与其混交林常见树种兴安落叶松向北迁移趋势基本一致[36]。在西南部的适生区面积也随高温高湿的加剧而收缩，但西南适生区收缩速率明显低于东北适生区。在温室气体排放浓度相同时，除SSP126环境下的2090 s白桦适生区质心点向西北小幅度迁移外，在SSP245和SSP585两种较高温室气体排放浓度下，由于东北地区的白桦适生区面积相较西南地区白桦适生区收缩幅度更大，2090 s和2050 s白桦适生区质心均向西南方向迁移，并且迁移距离随温室效应加剧而进一步加大。这与我国温带森林树种地理分布在未来全球变暖下向高纬度区域和高海拔地区迁移的趋势基本相符[37]。