基于MaxEnt模拟欧亚大陆气候变化下叉子圆柏的潜在分布*

王爱君 路东晔 张国盛 黄海广,2 王 颖 呼斯楞 敖 民

(1.内蒙古农业大学林学院 呼和浩特 010019； 2.内蒙古自治区林业科学研究院 呼和浩特 010010；3.内蒙古鄂尔多斯市乌审旗林业局 鄂尔多斯 017300)

气候变化对物种分布有至关重要的影响，小到生物个体生理代谢、繁育，大到生态系统结构组成和功能发挥。近年来城市化进程加快，自然生态环境面临加速破碎化，联合国政府间气候变化专门委员会在第5次评估报告(AR5)中指出全球气候将继续变暖，至21世纪末地球的平均温度将增加0.3～4.5 ℃(相对于1986—2005的基准气候条件)，相应的降水、辐射、潜在蒸散以及其他气候变量的季节性变化也将随气候变暖而发生改变(Stockeretal.， 2014)。生态环境变化会直接影响到物种的地理分布格局及生态系统结构、功能与稳定性(Chenetal.， 2011； Gibsonetal.， 2011； Dielemanetal.， 2015)。应用特定的算法分析物种潜在地理分布，已成为区域生态学和生物地理学研究的热点之一(Bellardetal.， 2012)。

为了保护和加强植物资源管理，抑制生物多样性的减少，人们依据客观存在的种群地理位置及其环境变量信息构建物种分布模型模拟不同气候下植物潜在地理分布变化范围。所谓生态位模型(Ecological Niche Modelling，ENM)就是利用物种现有分布点与环境因子的关联性，推测物种生存所需因子的有效范围，并将运算结果投射到不同时间或地理空间中，以概率的形式来预测物种潜在分布区域(朱耿平等， 2013)。随着科学技术的发展，开发了基于不同算法的多种生态位模型，如BioClim、DoMain、GARP、GLM、MaxEnt(Elithetal.， 2006； Qiaoetal.， 2012)等，每种方法各有优劣。最大熵模型(Maximum Entropy，MaxEnt)是目前使用最为广泛的生态位模型(Ahmedetal.， 2015； Barbosaetal.， 2015； Vazetal.， 2015)，在预测史前地质时期与未来气候情景方面都有良好的应用(Bertrandetal.， 2012)，也是物种保育(Kumaretal.， 2009)、谱系地理(Wangetal.， 2015)、物种潜在分布区模拟(李璇等， 2018)等研究领域的首选模型。MaxEnt模型的优势在于其特征类的选择可以运行各种数量的样本，处理变量间复杂的相互作用(Wiszetal.， 2008； Van Proosdijetal.， 2016)。在现有的物种分布建模算法中，MaxEnt被证明表现最好，同时预测结果受个别错误数据的影响最低(Grahametal.， 2008)。管毕财等(2016)利用中国91个分布位点信息及气候数据，基于MaxEnt和ArcGIS软件模拟预测四照花(Cornuskousasubsp.chinensis)分布格局，推测其冰期避难所。郭飞龙等(2020)基于MaxEnt模型预测了胡杨(Populuseuphratica)世界范围内潜在适宜区，认为胡杨地理分布受多种环境变量影响，仅以气候变量模拟的胡杨适宜分布区与实际分布范围存在较大偏差，而多类环境变量综合模拟更能反映其现实种群的分布特征。生态位模型不仅在植物方面有广泛的应用，同时也适用于动物。朱井丽等(2018)基于MaxEnt模型对松嫩平原丹顶鹤(Grusjaponensis)秋迁期生境适宜性进行了评价。最大熵模型与地理信息科学的有机结合，可以预测不同气候情景下的物种潜在分布区，深入理解物种进化历程； 估算分布面积变化情况，探究物种迁移路线，预估未来种群可能会受到全球变化的威胁程度，也可推测拥有丰富遗传多样性的避难所。

叉子圆柏(Juniperussabina)为柏科(Cupressaceae)刺柏属(Juniperus)匍匐常绿灌木，因其发达的根系网与强大的萌蘖力，使它既可以生长在石质山坡、河谷上，也可以在固定、半固定的沙丘及滩地边扎根，是优良的水土保持与固沙造林资源树种，对维持生态环境稳定安全具有重要作用。叉子圆柏还可入药，作饲料和薪炭林，也是城市绿化、盆景园艺的珍贵树种。叉子圆柏自然地理分布区集中在欧洲南部、亚洲的中部及东部地区，是旱中生常绿灌木，属泛北极植物区系，具有极强的适应性，是重要建群种。由于气候趋干旱化、人为毁林开荒、不合理放牧等现象频繁发生，导致叉子圆柏种群退化严重，大多数为长势差、质量低的块状残林，因此叉子圆柏灌丛林的保护与修复刻不容缓。叉子圆柏在生理生化、栽培技术、生长与环境关系(孟和等， 2010)、生态学特性、遗传多样性、基因组大小(路东晔等， 2018)、谱系地理(路东晔等， 2020)等方面已开展了广泛研究。叉子圆柏在欧亚大陆广大区域内呈断续分布且各分布区生境因子差异大，存在局部地区资源分布信息不详，缺乏对叉子圆柏欧亚大陆分布区生态环境变量特征和生态过程的全面解析(郭飞龙等， 2020)。García-Cervign等(2018)选取西班牙东南部安达卢西亚自治区的4个种群，研究表明冬季降水延迟会影响叉子圆柏在夏季干旱条件下的生长。本文应用MaxEnt、BioClim、DoMain 3种预测模型和ArcGIS软件，收集已知的欧亚大陆叉子圆柏种群地理位置信息，基于气候环境变量和地形(海拔)数据，模拟其潜在地理分布区，使用AUC值(受试者工作特征(ROC)曲线下方面积)和Kappa值对3种模型进行比较分析与筛选。采用MaxEnt模型预测不同气候下适宜区的动态变化，运用ArcGIS 10.2软件“SDM Tools”估算各个时期分布区面积变化，质心总体迁移趋势。旨在探讨制约叉子圆柏地理分布的主导环境变量，明确叉子圆柏对生态环境变量的需求，为科学合理有效开展叉子圆柏种群保护修复、种质资源保育及人工栽培提供参考，对促进欧亚大陆乃至全球叉子圆柏种群和种质保育建设具有现实意义。

1 材料与方法

1.1 数据来源

叉子圆柏地理数据： 国家植物标本资源库(http:∥www.cvh.org.cn)、全球生物多样性信息平台(http:∥www.gbif.org)、中国国家标本资源平台(http:∥www.nsii.org.cn)，查阅收集植物的相关文献资料(Adamsetal.， 2017； 2018； 路东晔等， 2020)和野外采点数据。经检索，共获得标本采集记录种群信息和野外采点数据2 950条记录。

环境变量信息： 末次盛冰期(Last Glacial Maximum，LGM，距今约22 000年)、全新世中期(Middle Holocene，MID，约6 000年前)、当前(Current， 1970—2000年)、未来(Future， 2070：2061—2080年)的气候变量和地形(海拔)数据均来源于世界气候数据库(WorldClim，https:∥www.worldclim.org)。未来气候情景采用美国大气研究中心(NCAR)开发的CCSM4模型(Otto-Bliesneretal.， 2006)的推算数据，该模型对气候变量模拟与预测的准确性都比较高(李文庆等， 2019)。基于政府间气候变化专门委员会发布的第5次评估报告，选取典型浓度路径(RCPs)，RCP2.6和RCP8.5分别代表未来温室气体低和高2种排放情形(Van Vuurenetal.， 2011)，建立预测模型。即在RCP2.6的情景下，预测全球平均气温上升0.3～1.7 ℃，全球平均海平面上升0.26～0.54 m； 在RCP8.5的情景下，预测全球平均气温上升2.6～4.8 ℃,全球平均海平面上升0.45～0.81 m(Alexanderetal.， 2013)。下载的19个气候变量和地形(海拔)数据，其空间分辨率均为2.5′。

1.2 数据处理

剔除标本无详细地理位置、信息重复及人工林的数据，利用Google Earth软件(http:∥earth.google.com)核查删除部分分布位置不正确的点，确定种群位置经纬度。为降低样点的空间自相关且与环境数据范围保持一致，将世界地图分成5 km×5 km(空间分辨率为2.5′)为1个样本单元格网，删除格网内冗余点，共筛选确定267个有效样点用于建模。

应用ArcGIS10.2软件(https:∥www.Esri.com)空间分析功能，对19个气候变量进行相关性检验。以变量之间相关性值> 0.90为标准，保留有生物学意义且在初始模型中贡献率大的变量，减少由于各变量间高度相关而导致模型过度拟合，最终确定Bio1、Bio2、Bio4、Bio7、Bio8、Bio9、Bio12、Bio14、Bio15、Bio17、Bio18、Bio19和Elev(海拔)，共计13个环境变量(表1)。

表1 预测叉子圆柏地理分布的环境变量Tab.1 The environmental variables for prediction of geographical distribution of Juniperus sabina

1.3 模型筛选

将267个叉子圆柏分布数据和环境变量共同导入MaxEnt3.4.1软件中。特征类选择(feature class selection)是MaxEnt重要的参数设置(邱浩杰等， 2020)，本研究选择自动特征。设置25%为验证集(testing data)，75%为训练集(training data)，最大迭代次数为10 000，重复训练(replicates)为3，同时交叉验证(crossvalidation)(Phillips， 2010)防止人为统计化。勾选刀切法(Jackknife)检验各环境变量对叉子圆柏分布的贡献率，Logistic输出.asc。

同样将267个叉子圆柏分布数据和环境变量导入DIVA-GIS7.5软件(http:∥www.diva-gis.org)(杨芙蓉等， 2019)，设置25%为验证集(testing data)，75%为训练集(training data)，重复训练(replicates)为3，将13个环境数据转换成.grd格式并生成stack数据集，在Modeling-BioClim/DoMain模块中依次进行BioClim和DoMain模型预测。

使用受试者工作特征(receiver operating characteristic，ROC)曲线下方面积(area under curve，AUC值)和Kappa值对模型预测结果进行检验。二者取值范围均为0～1，值越接近1说明模型模拟精确性越高(Swets， 1988)。AUC值大于0.90，预测效果极好； Kappa值大于0.80，预测效果极佳(杨芙蓉等， 2019)。AUC值及其标准差共同反映模型的总体水平，AUC值越大，标准差越小，模型的稳定性越高(Elithetal.， 2011； 陈新美等， 2012)，且随机取样对模型影响小。

采用ArcGIS10.2软件，对分布区结果进行适宜性划分和可视化。基于叉子圆柏灵敏度-特异度和最大化(sensitivity-specificity sum maximization approach)阈值0.20，设定物种出现概率小于0.20的地区为叉子圆柏非适生区，出现概率0.20～0.40为叉子圆柏低适生区，概率0.40～0.60为叉子圆柏中适生区，概率大于0.60为高适生区。

将叉子圆柏欧亚大陆分布区内267个有效样点地理数据与气候、地形环境变量进行MaxEnt、BioClim、DoMain模型模拟，3种预测模型的AUC与Kappa值表明，MaxEnt和DoMain模型均可用于叉子圆柏潜在分布区模拟，且预测效果好，稳定性高(表2)。

表2 3个模型的AUC和Kappa值Tab.2 AUC and Kappa values of three models

叉子圆柏潜在适宜分布区横跨欧亚大陆(图1)，基本反映出叉子圆柏在欧亚大陆的分布特征，但存在3种模型的原理与算法不同，预测的适宜分布区域略有差异。MaxEnt模型预测亚洲的主要适宜分布区为亚洲西部(土耳其、格鲁吉亚、亚美尼亚)、亚洲中部(哈萨克斯坦、吉尔吉斯斯坦)、亚洲东部(中国)的干旱和半干旱地区； 欧洲为欧洲南部(西班牙、法国、德国、意大利、捷克、瑞士、奥地利、波兰、俄罗斯、保加利亚、克罗地亚、塞尔维亚)地区。DoMain模型预测的分布结果整体向北偏移，亚洲的主要适宜分布区为亚洲中部(哈萨克斯坦)、亚洲东部(中国、蒙古)地区； 欧洲为欧洲北部(西班牙、法国、德国、捷克、奥地利、波兰、俄罗斯、乌克兰、立陶宛、白俄罗斯、罗马尼亚)地区。BioClim模型预测的分布区进一步向北偏移，亚洲的主要适宜分布区为亚洲西部(土耳其、格鲁吉亚、亚美尼亚)、亚洲中部(哈萨克斯坦、吉尔吉斯斯坦)、亚洲东部(中国)地区； 欧洲为欧洲北部(德国、捷克、奥地利、波兰、俄罗斯、乌克兰、白俄罗斯、罗马尼亚、保加利亚、塞尔维亚、斯洛伐克)地区。MaxEnt、BioClim、DoMain模型预测高适生区面积分别为118.28×103、157.57×103、195.64×103km2，呈递增规律，这与3种模型的运算原理有关： MaxEnt模型是基于最大熵理论(马松梅等， 2010)，在已知的约束条件下，对未知事物不做任何假设，受取样点的影响最小，熵最大时事物最接近真实状态，其预测结果比较保守； BioClim模型是基于环境包络原理，假定此环境气候因素落在环境包络内，则物种能在该环境中正常生长、繁殖后代，故预测结果也相对保守； DoMain模型是根据Gower距离算法得到点-点的相似矩阵，受取样点影响较大，其预测结果与分布点关系紧密，分布非连续且均以样本点为中心向周围分散(杨芙蓉等， 2019)。基于了解环境变量对叉子圆柏地理分布贡献率为研究目标，鉴于3种模型的AUC、Kappa值、预测分布图比较相近的情况下，MaxEnt模型更具优势。

图1 MaxEnt、DoMain、BioClim模型模拟叉子圆柏潜在分布Fig. 1 Simulation of potential distribution of J. sabina by MaxEnt, DoMain, BioClim models

1.4 环境变量重要性评估

在MaxEnt3.4.1软件中基于刀切法，计算不同环境变量测试增益(training gain)，结合自动生成的环境变量贡献率与置换重要值，检验不同环境变量限制叉子圆柏地理分布的重要性大小。

1.5 适生区动态变化及迁移趋势

基于MaxEnt模型绘制末次盛冰期、全新世中期、当前、未来情景下叉子圆柏潜在地理分布区。运用Reclassfy分类工具，统计各类适宜区面积。利用ArcGIS 10.2软件中“SDM Tools”中“Distribution Changes between Binary SDMs”工具，获得末次盛冰期(LGM)—全新世中期(MID)、全新世中期(MID)—当前(Current)、当前(Current)—未来(Future)2070年叉子圆柏适生地扩张、收缩和稳定区域面积变化值。采用“SDM Tools”中“Centroid Changes(Lines)”工具，将各个时期的分布区域集合为一个中心点，创建矢量文件标注叉子圆柏分布方向及位置变化，探讨欧亚大陆叉子圆柏分布区迁移趋势。

2 结果与分析

2.1 制约叉子圆柏分布的主导环境变量

MaxEnt模型运行过程中自动生成的环境变量贡献率及置换重要值表明，基于气候和地形2类环境变量模拟，贡献率大于5%的分别为年均气温Bio1(28.80%±0.86%)、海拔Elev(29.00%±0.36%)、温度季节性变化Bio4(13.00%±0.98%)、最干月降水量Bio14(7.10%±2.14%)和最热季降水量Bio18(5.70%±0.49%)，累积贡献率为83.6%； 置换重要值排名前5的分别是海拔Elev(45.70%±5.67%)、年均气温Bio1(18.70%±2.12%)、温度季节性变化Bio4(8.00%±7.08%)、降水量季节性变化Bio15(7.70%±1.50%)和最干季平均温度Bio9(6.20%±1.08%)，累计置换重要值为86.3%(表3)。刀切法检验不同环境变量的测试增益结果表明，年均气温Bio1、温度季节性变化Bio4、海拔Elev、最干季平均温度Bio9、年温度变化范围Bio7排在前5位(图2)。统一分析环境变量贡献率、置换重要值及测试增益，采用2类环境变量模拟，地形(海拔)、大气温度(年均气温Bio1、温度季节性变化Bio4)是限制叉子圆柏分布的主导环境变量，同时也受大气降水(降水量季节性变化Bio15、最干月降水量Bio14)的调控。

表3 各环境变量贡献率Tab.3 Percent contribution of each environmental variable

图2 刀切法检验环境变量的测试增益结果Fig. 2 Jackknife test to examine training gain results of environmental variables

2.2 叉子圆柏欧亚大陆潜在适生区的变化特征

利用ArcGIS软件对MaxEnt模型模拟结果重新分类，分别统计叉子圆柏各适生区面积，即LGM时期为72.894×103km2、MID时期为72.245×103km2、当前时期为663.115×103km2、未来RCP2.6和RCP8.5时期分别为67.116×103km2和 55.355×103km2，呈现先增加后减少的变化趋势(表4)。MaxEnt 模型模拟结果表明，亚洲是叉子圆柏的主分布区，亚洲的适生区面积在LGM时期占86.9%、MID时期占87.0%、当前时期占57.8%、未来RCP2.6和RCP8.5时期分别占84.1%和79.2%。

表4 各时期叉子圆柏适宜分布区面积①Tab.4 Suitable distribution area of J. sabina in each period km2

基于MaxEnt模型模拟预测可知，叉子圆柏适生区范围及面积随着气候环境变化而变化，并呈现先扩张后退缩的特征。

从LGM—MID时期叉子圆柏适生区面积虽然下降了0.64×103km2，但在欧洲的意大利、瑞士、奥地利境内的阿尔卑斯2 000～3 000 m的山脉、俄罗斯西南部与格鲁吉亚交界的大高加索山脉，亚洲西部的格鲁吉亚与亚美尼亚交界的小高加索山脉、土耳其东北部的大阿勒山，亚洲中部的吉尔吉斯斯坦境内的帕米尔-阿赖山脉、哈萨克斯坦境内的塔尔巴哈台山，亚洲东部的中国境内的阿尔泰山、天山、祁连山脉、青海日月山和阴山的亚高山地带、青海湖、浑善达克沙地及蒙古的杭爱山和草原植被带，发生扩张，分布区面积扩增27.32×103km2。在欧洲的法国东南部、意大利中部、保加利亚、塞尔维亚、克罗地亚部分地带、俄罗斯西南部，亚洲西部的格鲁吉亚、亚美尼亚与土耳其东北部的低海拔山脉(1 000～1 500 m)，亚洲中部的吉尔吉斯斯坦境内的帕米尔-阿赖山脉、天山山脉，亚洲东部的中国境内的天山、祁连山脉、吕梁山、贺兰山、六盘山、柴达木盆地东南区、青海南山，以及欧亚大陆其他低纬度、低海拔、盆地等生境，为主要退缩区，退缩面积约为31.08×103km2。相对稳定区域的面积达39.99×103km2，多集中分布于中亚的天山和阿尔泰山山脉(表5，图3A)。

MID—当前时期，是欧亚大陆叉子圆柏种群扩张时期，适生区面积从72.245×103km2快速增加到663.115×103km2，分布区增加面积为559.71×103km2。此期间叉子圆柏分布遍及欧洲南部、亚洲中部及东部(图3B)。仅在蒙古北部高山草地区、吉尔吉斯斯坦境内的天山、中国青海湖及附近山脉有少部分的收缩，收缩面积约为17.39×103km2。

在RCP2.6的气候情景下，到2070年，亚洲中部吉尔吉斯斯坦中部天山区域及亚洲东部蒙古北部高山区域以及中国祁连山脉和青海日月山是叉子圆柏潜在适生区面积增加的主要区域(图3C)，新增面积约为22.80×103km2(表5)。除欧洲的意大利、瑞士、奥地利境内的阿尔卑斯山脉、俄罗斯西南部与格鲁吉亚交界的大高加索山脉、格鲁吉亚与亚美尼亚交界的小高加索山脉、土耳其东北部的大阿勒山，亚洲中部吉尔吉斯斯坦境内的天山山脉，亚洲东部中国境内的阿尔泰山和天山山脉、蒙古草原地带外，其余均处于收缩状态。

在RCP8.5气候情景下，叉子圆柏欧亚大陆潜在分布范围保持与在RCP2.6气候情景下相同的缩减趋势，缩减582.81×103km2(表5)。新增的范围主要分布于天山、蒙古高山区和中国境内的阿尔泰山、柴达木盆地，新增22.75×103km2。与RCP2.6气候情景下增加的范围相比较，RCP8.5气候情景下新增的区域更近一步向高纬度扩张(图3D)。

表5 不同气候时期叉子圆柏潜在分布区的动态变化Tab.5 Dynamic changes of potential distribution areas of J. sabina in different climatic periods

图3 不同气候情景下叉子圆柏分布格局动态变化Fig. 3 Distribution pattern dynamic changes of J. sabina under different climate scenariosA： 末次盛冰期(LGM)— 全新世中期(MID)的叉子圆柏适生地扩张、收缩和稳定区域； B： 全新世中期(MID)—当前(Current)的叉子圆柏适生地扩张、收缩和稳定区域； C： 当前(Current)—未来(Future)2070年RCP2.6的叉子圆柏适生地扩张、收缩和稳定区域； D： 当前(Current)—未来(Future)2070年RCP8.5的叉子圆柏适生地扩张、收缩和稳定区域。A： The expansion, contraction and unchange regions of J. sabina from LGM to MID; B： The expansion, contraction and unchange regions of J. sabina from MID to Current; C： The expansion, contraction and unchange regions of J. sabina from Current to Future 2070(RCP2.6); D： The expansion, contraction and unchange regions of J. sabina from Current to Future 2070(RCP8.5).

2.3 叉子圆柏分布质心的变化

叉子圆柏在末次盛冰期的分布中心位于中国新疆西部(25° N，100° E)，全新世中期位于中国新疆西北部(30° N，105° E)，当前的分布中心位于哈萨克斯坦西南部(28° N，65° E)，在RCP2.6(2070)气候情景下的分布中心位于中国新疆北部(33° N，110° E)，RCP8.5(2070)气候情境下的分布中心位于蒙古西部(35° N，113° E)。叉子圆柏的分布中心从末次盛冰期至全新世中期，整体向东北迁移； 从全新世中期到当前，出现一个回折现象； 但从当前至未来，又继续呈现向北延伸的趋势。

2.4 叉子圆柏适宜生境因子的阈值范围

环境变量贡献率、置换重要值表明，年均气温和海拔的贡献率相近(29%)(表3)，但是在置换重要值方面海拔明显大于年均气温，对叉子圆柏欧亚大陆潜在地理分布格局起着主导作用。将年均气温(Bio1)、海拔(Elev)、温度季节性变化(Bio4)、年降水量(Bio12)依次导入MaxEnt模型中进行单因子建模，以存在概率大于0.40为阈值计算4个主要气候环境因子的适宜范围。欧亚大陆叉子圆柏适生区年均气温适宜范围约为0～7.5 ℃，海拔适宜范围为900～2 600 m，温度季节性变化适宜范围为5.5～6.5 ℃，年均降水量的适宜范围为200～800 mm。这些数据总体可以定量地描述叉子圆柏适宜生境的气候环境特征，与孟和等(2010)的报道基本一致。

3 讨论

3.1 模型模拟结果的可靠性

样本采集范围及样本量是决定物种分布模型模拟结果可靠性的重要因素。一般而言，收集到的样点覆盖度均匀、区域广、样本量大，获取到的环境信息则越丰富，模型建立约束条件越多(庄鸿飞等， 2018)，预测物种潜在分布图就越精准，模型估计精度越高(赵晓冏等， 2018)。AUC值被公认为是模型预测准确度的最佳衡量指标(Hughesetal.， 2003)。Kappa统计量是一种阈值依赖评价指标，也被大量用作模型评价参数(Maneletal.， 2010； Seguradoetal.， 2004)。模拟过程中筛选的环境变量是否具有代表性、完整性以及种群的分布位置、环境变量的相关性，同样制约着评价结果是否可靠。本研究基于环境变量相关性分析，剔除相关系数大于0.9的变量，MaxEnt模型模拟结果较能体现叉子圆柏符合柏科分布带的特征，接近前人的描述(毛康珊， 2010)。叉子圆柏欧亚分布区的样本能代表其分布区的生境，MaxEnt模型模拟叉子圆柏在中国适宜分布区的面积(129.643×103km2)大约是王林和等(2014)报道的基于卫星遥感影像面积(0.032 7×103km2)的3 000倍以上，与郭飞龙等(2020)、刘超等(2018)和王茹琳等(2017)的研究结果类似，都存在模拟面积大于实际分布面积的情况。这可能是仅使用单一非生物环境变量，没有考虑限制植物分布的其他生物因素，如： 繁殖能力、种间竞争、人类干扰等(Dülgeroluetal.， 2019)； 也有可能是由于物种自身较弱的迁移能力和存在迁移屏障造成物种没有到达其适宜生境(Sobernetal.， 2005； 2009)。在某些地区叉子圆柏更新庇护植物乌柳(Salixcheilophila)、沙棘(Hippophaerhamnoides)等被砍伐后无法更新扩展形成叉子圆柏灌丛林(山中典和等， 2000)，此外城市扩张、土地用途变更、樵采、地表浅水过度利用等人类社会活动也是造成原始种群消失的重要因素，可见生物因子严重影响着叉子圆柏的实际分布状况。

3.2 制约叉子圆柏分布的环境变量

植物种群的分布特征是生物和非生物因素共同作用的综合结果，也是植物对环境响应的最终体现。叉子圆柏为北温带分布类型中的古北极植物种，至少起源于第三纪，曾是北半球森林植被的重要组成成分。MaxEnt模拟预测显示叉子圆柏适宜分布区与实际分布拟合程度大体一致。环境变量贡献率、置换重要值及测试增益显示，地形(海拔)、气候(温度和降水)对叉子圆柏地理分布格局贡献最大，这与石质山地是叉子圆柏的主要生境相关。温度和水分的水平分布特征决定着植物生理活动和生化过程(Wangetal.， 1998)，而随着海拔的变化对气候特征进行了重塑，在隐域植被形成中起着关键作用，水热因子在塑造植物生态适应方面发挥着重要作用(李文庆等， 2019)，此外土壤性质、水文地质等立地相关的环境变量同样制约着叉子圆柏的分布。

第四纪以来，地球气候的剧烈变化对当今物种的分布格局影响重大(Hewitt， 2004)。末次盛冰期(LGM)是距今最近的极寒冷时期，全球陆地24%被冰雪覆盖，两极冰盖充分扩张，至少覆盖至30° N中低纬度地区(王绍武， 2011)。本研究结果发现，在冰期叉子圆柏的分布区多位于30° N左右，分布中心均位于低纬度、低海拔山地和低盆地，较当前分布偏南，很可能是由于高海拔和高纬度地区冰川形成，迫使物种向低纬度、低海拔地区(1 000～1 500 m)迁移，同时受到山地地形的影响，形成特殊的微生境，为叉子圆柏及其他植物种提供了冰期避难所得以存活。

全新世中期(MID)全球呈现温暖湿润且较稳定的特征，平均年均气温比LGM时期高4～5 ℃，海平面上升到与现代相差无几的高度，降雨量增加(表6)。虽然MID时期气温和降雨量对植物生存有利，但与LGM时期相比叉子圆柏并未迅速扩张其分布面积，而分布区位置发生了迁移。通过模型预测MID时期温度升高，冰川开始融化，山地雪线普遍上升了约300～500 m，低海拔区域由陆地变成湖泊，物种开始向高海拔山地迁移(海拔1 500～2 500 m)，叉子圆柏的分布区也向北推进。MID时期气候的另一个特点是早期气候湿润，中期以后较为干旱，全球气候还处在一个波动激烈的时期，有B/A暖期，也有YD事件。在YD时期南欧比北欧更为干旱，孢粉记录表明希腊甚至比LGM时期还要干旱，年降水量可能低于150 mm(王绍武， 2011)，导致迁移过程中叉子圆柏分布区面积微缩(表4)。

表6 各时期叉子圆柏适生区的年均气温(Bio1)、海拔(Elev)、温度季节性变化(Bio4)、年降水量(Bio12)变化范围及均值Tab.6 The variation range and mean value of annual average temperature(Bio1), elevation(Elev), seasonal variation of temperature(Bio4) and annual precipitation(Bio12) in suitable growing areas of J. sabina in each period

MaxEnt模型模拟结果表明，当前叉子圆柏呈条带状分布在欧亚大陆30°～50° N，分布面积明显扩增。从全新世后期到当前，环境因素也发生了复杂的变化，基于沉错钻孔孢粉记录重建的全新世温度发现，温度最高值出现在初期(中国第四纪孢粉数据库小组， 2000)，后期(亚大西洋时期)气候趋凉爽潮湿，喜冷湿植物群落扩张，出现水青冈(Faguslongipetiolata)、云杉(Piceaasperata)、冷杉(Abiesfabri)林。叉子圆柏又作为云杉林下伴生树种，也应大量繁殖扩展其分布区。质心分析也表明当前相较于MID时期，叉子圆柏分布区向西偏南方向偏移扩大。

未来由大气CO2增加引起的全球增温与MID时期的增温效应类似(郑益群等， 2004)。随着全球气候变暖，许多动植物都有向高纬度和高海拔地区迁移的趋势(Sekerciogluetal.， 2008； 方精云等， 2018)。基于RCP2.6和RCP8.5的CO2排放量，测算到2070年叉子圆柏适宜生存区的海拔在2 500～3 000 m，低纬度、低海拔区域的叉子圆柏分布面积将大量减少，面临着消失的风险。MaxEnt模型模拟预测也表明，未来2070年(RCP2.6和RCP8.5)全球叉子圆柏主要分布在蒙古北部高山草原区，意大利、瑞士、奥地利境内的阿尔卑斯山，俄罗斯西南部与格鲁吉亚交界的大、小高加索山，土耳其东北部的大阿勒山，吉尔吉斯斯坦境内的天山，中国境内的阿尔泰山、天山、祁连山这些高山地带。分布区总面积与当前相比锐减，不及MID时期，潜在适生生境丧失或更加破碎化。胡菀等(2020)在研究观光木(Tsoongiodendronodorum)的潜在地理分布变迁中也得到了相同的结果，即未来将导致观光木的生境减少远甚于历史时期。

3.3 叉子圆柏资源保护

叉子圆柏是欧亚大陆亚高山地、石质山地及沙地的重要生态树种。20世纪以来，虽然人们的生态保护意识日益增强，如： 中国西北部地区对叉子圆柏的保护、恢复、重建工作取得了重大进展，部分地区已建立起省市级叉子圆柏自然保护区，有效地改善了当地生态环境，防止了过度放牧的破坏，但自然叉子圆柏灌丛林生境资源与发展状况仍不容乐观，天然林分结构退化没有得到根本遏制。天然叉子圆柏更新虽然拥有匍匐茎萌蘖和种子繁殖2种模式，但其种子库小、质量差、更新环境丧失等原因限制了许多种群的有性更新(王哲等， 2005； 张国盛等， 2006)。自然更新能力弱，匍匐茎萌蘖成为维持叉子圆柏种群生存的主要繁殖方式，进化潜力低，一些种群的遗传多样性下降或偏低(路东晔等， 2020)。叉子圆柏为雌雄异株异花授粉树种，自然分布区广，在历史进化过程中不同地理种群、不同发育阶段，都经受着环境变量选择压力的不同(刘建全， 2016)，在形态特征、生理特性、生态习性等方面产生了遗传差异，丰富了种内遗传变异，也形成许多变种(Adamsetal.， 2017； 2018； 路东晔等， 2018； 2020)。因此，在叉子圆柏分布区范围内，开展种质遗传资源调查、收集、管理与评价等一系列工作，并在此基础上建立原地和异地遗传资源保存林，来防止叉子圆柏多样性的丢失是十分必要的。

基于本研究的结果，意大利、奥地利、瑞士境内的阿尔卑斯山，俄罗斯与格鲁吉亚接壤的大、小高加索山脉，中亚的天山、阿尔泰山在4个时期内均存在稳定分布的叉子圆柏种群，侧面说明这些地方极有可能是叉子圆柏的冰期避难所。此外在4个时期内，叉子圆柏欧亚大陆潜在适生区的变化特征显示叉子圆柏适生分布区为山区，且在未来分布预测中，呈现向高纬度高海拔地区迁移的趋势，因此要高度重视山区保护，确定优先保护单元； 相应的沙区生境破碎化严重，面临种群消失的风险，要加强沙区人工林水分条件管理，满足其天然更新的需求。叉子圆柏在欧亚大陆分布整体而言，西部变化要大于东部，青藏高原的隆起阻挡了西风带，从而导致西部整体受季风气候影响明显。叉子圆柏多分布于山地，多为社会经济条件相对落后的牧区草场，为促进当地经济发展，过度放牧、不合理樵采、滥垦、水资源浪费等，造成天然叉子圆柏大面积锐减。因此，改善当地民生条件，创新传统农耕模式，加大科普宣传实践力度，提高保护意识，是叉子圆柏保护的基本措施。

4 结论

叉子圆柏地理分布受多种环境变量共同影响，采用气候和海拔变量进行MaxEnt模型模拟，能较完整地拟合出叉子圆柏的实际分布状况，反映出叉子圆柏真实的分布特征。基于末次盛冰期、全新世中期、当前和未来(2070年)4个时期叉子圆柏生境适宜性评价，年均气温、海拔、温度季节性变化是限制叉子圆柏分布的主导环境变量，山地是其主要生境及避难所，亚洲是其主分布区。本研究结果可为叉子圆柏资源管理与发展、种群恢复与重建等提供理论依据。