气候变化对燕山地区北苍术适生区分布的影响

王冶 ，薛忠财 ，杨越 ，常佳宁 ，赵春颖 ，李云峰 ，曾凡然

1.河北民族师范学院资源与环境科学学院，河北 承德 067000；2.河北省山区地质环境重点实验室，河北 承德 067000；3.冀北山区土壤与特色植物营养创新研究中心，河北 承德 067000；4.承德医学院河北省中药研究与开发重点实验室，河北 承德 067000；5.承德笔架山生态农业科技发展有限公司，河北 承德 067000

物种分布区特征是物种重要的生态学和进化学特征[1]，对制定物种保育策略具有重要参考价值[2]，诸多因子均可能影响物种分布区特征[3-5]，其中气候变化已成为影响生态系统和生物类群的最主要因素之一[6]。气候变化导致我国一些野生药用植物资源面临巨大威胁[7-8]，预测未来气候变化对物种地理分布的影响，有利于指导濒危物种保育策略的制定。目前用于预测物种地理分布的常见模型已较为成熟[9]，其中最大熵（MaxEnt）模型基于最大熵原理，对物种的潜在分布进行预测，因其预测效果好、预测覆盖面积大等特点[10]，被广泛用于分析未来气候对物种分布的影响、野生动植物保护、物种适宜性等[11]。北苍术Atractylodes chinensis(DC.) Koidz.是菊科苍术属多年生草本植物[12]，其根茎是中药苍术的主要来源[13]。燕山地区所产北苍术“朱砂点密，香气浓郁”，因此作为该药材的道地产区[14]。燕山地区地形地貌、生态系统相对复杂，气候变化对北苍术适生区的潜在影响尚不明确。因此，研究气候变化影响下的北苍术适生区分布规律十分必要。

本研究基于实地调查和查阅文献获取燕山地区北苍术地理分布信息，应用MaxEnt模型和地理信息系统软件ArcGIS，分析当前燕山地区北苍术适生区分布，预测不同气候情景下未来不同时期潜在适生区变化，为北苍术引种驯化、资源保护和产业发展提供理论依据。

1 资料与方法

1.1 研究区概况与样点分布

本研究区域为广义燕山地区，位于坝上高原以南、河北平原以北、白河谷地以东、辽河平原以西的山地（93°41＇～43°19＇N、116°60＇～112°13＇E），调查区域总面积1.3×105km2，包含北京、张家口、唐山、承德、秦皇岛、赤峰、葫芦岛、朝阳、锦州等市的全部或部分区域。燕山地区处于温带大陆性季风气候区，年均温8.7～12.3 ℃，年降水量487.6 mm左右，作为京津冀地区的天然屏障和水源涵养地，海拔梯度较大，涵盖森林、灌丛、草地、农田、湿地等多样的生境类型。

野生北苍术地理分布信息收集：①本团队于2020年8月和2021年7月赴燕山地区进行实地考察，基于主要环境因子在燕山地区的空间分布特征进行布点，收集北苍术地理分布信息共36个；②查阅中国数字植物标本馆（https://www.cvh.ac.cn/）、中国国家标本资源平台（http://www.nsii.org.cn/）及相关文献[15-17]，经过数据筛选与核查，收集记录信息河北12 个、内蒙古5个、辽宁4个，利用ArcGIS10.0软件，基于燕山地区范围进行分布点信息校对，除去无效记录，剩余13个。最终确定燕山地区野生北苍术地理分布信息共49 个（见图1）。本文涉及的地图来自自然资源部标准地图服务系统（http://bzdt.ch.mnr.gov.cn/）。

图1 燕山地区北苍术样点分布

1.2 环境因子数据

选择生物气候变量、地形地貌变量及土壤属性变量等30个环境因子数据（见表1），19个温度、降水相关的生物气候变量数据[18]来自世界气候数据库（WorldClim，Version2.0，http://www.worldclim.org），包括当前时期（1970-2000 年）及未来3 个时期（2021-2040年、2041-2060年、2061-2080年），并基于联合国政府间气候变化专门委员会第六次评估报告，选择全球气候模型预测的温室气体低排放情景（SSP1-2.6）和高排放情景（SSP5-8.5）气候数据。地形地貌变量中数字高程模型（DEM）数据来源于地理空间数据云平台（http://www.gscloud.cn），坡度和坡向数据基于DEM，采用ArcGIS10.0软件表面分析工具计算获取。土壤属性变量（0～30 cm）数据来自世界土壤属性数据库（https://geodata.pku.edu.cn/）。以上数据空间分辨率均统一为2.5 km。

表1 30个环境因子

1.3 模型构建与精度评价

采用MaxEnt3.4.1 软件将北苍术样点分布信息（CSV格式）与不同时期的环境因子（ASCⅡ格式）数据输入软件中进行建模运算，刀切法计算各环境因子的贡献率，随机选择75%数据作为训练集，25%数据作为测试集，设置10次重复，同时设定最大迭代次数为10 000，其他参数为默认值[19]。利用ArcGIS10.0软件进行作图。采用受试者工作特征（ROC）曲线下面积（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表明预测结果优秀。AUC取值范围为0～1，越接近1表明模型预测越准确[21]。

1.4 适生区划分

基于模型结果，采用ArcGIS10.0软件重分类工具对适生区进行等级划分。根据概率P将适生区划分为4个等级：高适生区（P≥0.8）、中适生区（0.6≤P＜0.8）、低适生区（0.4≤P＜0.6）和非适生区（P＜0.4）[22]。

1.5 时空变化分析

采用重心迁移模型分析不同时期环境因子变化影响下，燕山地区北苍术适生区重心迁移轨迹，从而反映其时空格局变化趋势。重心迁移模型广泛应用在物种迁移、景观格局演变等空间过程与格局的研究中[7，23-24]，是反映对象空间演变过程的有效方法，其计算公式为：

式中，zi为第i个平面空间单元的AUC值，(xi,yi)为第i个平面空间单元的投影坐标值，点(xˉ,yˉ)为北苍术适生区地理分布重心的投影坐标值。重心移动的距离计算公式为：

式中，dm表示地理分布重心移动的距离；yt+m和yt分别表示在时期t+m和t的地理分布重心投影y坐标；xt+m和xt分别表示在时期t+m和t的地理分布重心投影x坐标。

2 结果与分析

2.1 模型评价

当前时期北苍术MaxEnt模型AUC平均值为0.913；未来3个时期（2021-2040年、2041-2060年、2061-2080年）不同温室气体排放情景下，MaxEnt模型AUC平均值均大于0.9，表明模型预测结果优秀，可用于北苍术适生区分布的预测研究。见图2。

图2 不同时期北苍术生态分布预测模型AUC比较

2.2 环境因子筛选

影响北苍术生长的主要环境因子为温度季节性变化（46.6%），其次为高程、降水量季节性变化、最湿季度降水量、年均降水量，累积贡献率达73.4%。主要环境因子适宜值和贡献率见表2。

表2 影响北苍术分布主要环境因子适宜值与贡献率

2.3 适生区及其动态变化

2.3.1 当前时期适生区分布

当前时期燕山地区北苍术适生区分布见图3。总适生区面积约7.54×104km2，占燕山地区总面积的59.10%，其中高适生区主要集中在承德西南部、张家口东部、秦皇岛西北部和唐山北部，总面积约1.06×104km2，占燕山地区总面积的8.29%；中适生区主要集中在承德东部、秦皇岛东部和赤峰西南部，总面积约2.12×104km2，占燕山地区总面积的16.61%；低适生区主要集中在朝阳西部、葫芦岛西部和北京北部，总面积约4.36×104km2，占燕山地区总面积的34.19%。综上表明，燕山地区北苍术适生区主要以承德、秦皇岛和张家口为核心，高、中、低适生区相连并呈发散状向东北方向分布。本研究北苍术地理分布信息虽无位于燕山地区西北部的分布点，但根据模型预测结果可以看出，该区域也适合北苍术生长。

图3 当前时期燕山地区北苍术适生区分布

2.3.2 未来气候变化对适生区分布的影响

未来不同时期2种气候情景下燕山地区北苍术适生区分布见图4，适生区面积变化见表3。结果表明，在SSP1-2.6情景下，随着时间变化，适生区面积增加的区域主要集中在朝阳中、西部和葫芦岛南部，而承德北部与赤峰南部适生区面积减少。2021-2040年适生区面积减少2 358.40 km2，占4.16%；2041-2060年适生区面积减少3 849.21 km2，占6.79%；2061-2080年适生区面积减少5 599.63 km2，占9.88%。在SSP5-8.5情景下，承德中北部与赤峰南部适生区面积随时间增加，而朝阳、葫芦岛、秦皇岛、唐山地区适生区面积均减少。2021-2040年适生区面积增加1 056.04 km2，占1.86%；2041-2060年适生区面积增加2 887.32 km2，占5.09%；2061-2080年适生区面积增加5 610.13 km2，占9.90%。综上表明，SSP1-2.6情景导致的气候变化使燕山地区北苍术适生区面积呈退化趋势，主要体现为高、中适生区向低、非适生区转化，而SSP5-8.5情景导致的气候变化使燕山地区北苍术适生区呈扩张趋势，主要体现为低、非适生区向高、中适生区转化。

表3 未来不同气候情景下燕山地区北苍术适生区面积转移矩阵

图4 未来不同气候情景下燕山地区北苍术适生区变化

2.4 气候变化对适生区空间迁移的影响

2.4.1 适生区水平分布变化

未来2种气候情景下，燕山地区北苍术重心迁移路线见图5。结果表明，当前适生区重心位于承德中部，SSP1-2.6情景下，适生区重心总体向燕山地区中、东部迁移，该区域为内蒙古高原与辽河平原的过渡带，最终迁至承德与朝阳交界处，各时期迁移距离分别为21.02、26.34、36.16 km，总体迁移距离为83.52 km。SSP5-8.5情景下，适生区重心总体向燕山地区北部迁移，该区域为山区与内蒙古高原，最终迁至赤峰南部地区，各时期迁移距离分别为26.27、27.27、46.94 km，总体迁移距离为100.44 km。

图5 不同气候情景下燕山地区北苍术适生区重心迁移路线

2.4.2 适生区垂向分布变化

提取2种气候情景下，不同时期的燕山地区北苍术适生区海拔数据，进行高度为300 m的频率分析。结果表明，各时期适生区面积垂向分布规律一致，随海拔升高，适生区面积呈先增加后减少的特征，主要分布海拔区间为300～900 m。见图6。

图6 不同气候情景下燕山地区北苍术适生区垂向分布变化

对燕山地区北苍术适生区垂向分布变化统计后发现：到2061-2080年，SSP1-2.6情景下，适生区面积在海拔区间0～300 m 累积增加30.12%，在海拔区间300～600、600～900 m，适生区面积分别累积减少19.06%、7.56%，其他海拔区间变化不显著；SSP5-8.5情景下，在海拔区间0～300 m，适生区面积累积减少24.96%，在海拔区间300～600、600～900 m，适生区面积累积增加10.14%、14.57%，其他海拔区间变化不显著。总体表明，在SSP1-2.6情景下，燕山地区北苍术适生区垂向分布呈现由高海拔向低海拔区域迁移的趋势，平均海拔下降88 m；在SSP5-8.5情景下，呈现由低海拔向高海拔区域迁移的趋势，平均海拔上升119 m。

3 讨论

根据环境因子贡献率结果，对燕山地区北苍术影响最大的气候变量为温度季节性变化，地形地貌变量主要为高程，二者累积贡献率达58.6%。温度季节性变化在1 100～1 350为最适宜条件，即全年4个季度温度标准偏差在11～13.5 ℃最有利于北苍术生长，与标本资料、实地调查获取及相关文献[25-26]基本一致；北苍术最适海拔为300～900 m，此结果与中国数字植物标本馆统计数据及《中国植物志》[12]记载较为一致。燕山地区属温带大陆性季风气候区，具有鲜明的季度气温差异，且由于坡向、海拔等地理因素[27]，气温呈现显著的东西差异，西部地区温度季节性变化为1 050～1 457，东部地区温度季节性变化为890～980，表明西部地区具有更大的季节性气温变化。且燕山地区的东、西部地形地貌差异显著，西部地区以高原山地为主，海拔在200～2 100 m，东部以丘陵平原为主，海拔30～180 m，总体呈现西高东低的地形特征。北苍术适生区集中在燕山西部地区，主要环境因子特征与燕山气温变化及地形特征高度吻合，表明模型计算的主要环境因子结果合理。

当前条件下，燕山地区北苍术适生区主要集中在承德、秦皇岛和张家口东部地区，与文献记载[15，28]一致。SSP1-2.6情景下，适生区分布有所退化，地理分布重心向东迁移，总体东向迁移至海拔较低的内蒙古高原与辽河平原过渡带；SSP5-8.5情景下，适生区分布有所扩张，地理分布重心向东北迁移，总体向更高海拔的内蒙古高原方向迁移。可能是温室气体低排放情景下，低海拔地区温度下降，达到了适生条件，导致北苍术向海拔较低的平原方向迁移；温室气体高排放情景下，气温升高使北苍术原适生区物种丰富度增加，出现生态位重叠，进而产生种间竞争，且高海拔山区温度有所升高，导致北苍术向海拔更高的高原方向迁移。赵文龙等[7]研究表明，气候变暖使红花绿绒蒿低海拔适生区出现退化，高海拔适生区面积扩大；郭彦龙等[29]研究发现，气候变暖使桃儿七适生区的平均海拔逐渐升高。

本研究显示，温度季节性变化与高程是影响燕山地区北苍术生长的主要环境因子，当前总适生区面积约7.54×104km2，适生区较为广泛。未来随着气候变化，北苍术分布范围将出现向不同海拔地区迁移的现象，因此在北苍术人工驯化和引种区域的海拔选择上，应注意气候变化，防止适生区退化。本研究结果有助于了解气候变化对燕山地区北苍术适生区分布的影响，为北苍术资源可持续利用及产业发展提供参考。