Response of the suitable distribution areas of Prunus pedunculata to Climate Change

本研究中459个长柄扁桃的地理数据来自实地调查（2013—2019）和网络数据平台。实地调查时利用GPS采集分布点地理数据（实地调查获得433个地理数据，全部分布在中国）。网络平台的数据（共26个，19个位于蒙古，7个位于俄罗斯）来源于全球生物多样性信息平台 （GBIF）（http://www.Gbif.org/）。

本研究采用19个生物气候变量来预测物种的适宜分布，生物气候变量分别为：年均温（bio1）、平均气温日较差（bio2）、等温性（bio3）、温度季节性（bio4）、最温暖月份的最高温度（bio5）、最冷月份的最低温度（bio6）、温度年变化范围（bio7）、最潮湿季度的平均温度（bio8）、最干燥季度的平均温度（bio9）、最温暖季度的平均温度（bio10）、最冷季度的平均温度（bio11）、年降水量（bio12）、最潮湿月份的降水量（bio13）、最干燥月份的降水量（bio14）、降水季节性变化（bio15）、最潮湿季度的降水量（bio16）、最干燥季度的降水量（bio17）、最温暖季度的降水量（bio18）、最冷季度的降水量（bio19）。所用的气候情景包括基准气候情景和未来气候情景，所有情景均包含有19个生物气候变量。基准气候情景的时间跨度为1970—2000年，每个生物气候变量取值为31年间的平均值。未来气候情景有4个，依据CO₂排放体量分为低排放情景（SSP1-2.6），中等排放情景（SSP2-4.5），中高排放情景（SSP3-7.0）和高排放情景（SSP5-8.5）^[15-16]，时间跨度为未来80 a（2021—2100），有4个分段2021—2040、2041—2060、2061—2080和2081—2100，每个情景均包括4个时间分段，每个时间分段内生物气候变量取值为20 a间的平均值。以上所有情景数据精度为2.5弧分（2.5 km），均来自世界气候数据库（http://www.worldclim.org/）^[17]。为了便于记录，未来气候情景中每20 a的时间尺度仅用起始年份标记，例：2021—2040 仅标记为 2021。

用MaxEnt（版本3.4.4）预测长柄扁桃的适宜分布区，将分布数据以25%的占比作为检验集，75%的占比作为训练集输入MaxEnt模型中，将模型运行10次以计算平均值，其他参数是软件默认值。模型输出数据包括预测结果、受试者工作特征（ROC）曲线、刀切法（Jackknife）检验结果和变量贡献率统计表等。

为避免变量间的高相关性导致模型过度拟合，在运行全变量预测后对气候变量进行筛选^[18-19]。使用筛选后剩余变量，重新编译MaxEnt模型，模型的计算效率和预测精度都将得到提高。筛选过程如下：

（1）在SPSS（版本26.0）软件中用皮尔逊相关系数（Pearson correlation coefficient）计算19个气候变量两两之间的相关性。

（2）删除预测过程中贡献率百分比过小的变量（贡献率百分比<1%）。在剩余具有高相关性的变量中（相关系数绝对值大于0.8），留下贡献率最高的变量。将筛选后的变量用于重新编译模型。

使用 ArcGIS （版本10.2）中的To Raster工具将 MaxEnt 导出的 ASCII 文件转换为栅格图层（长柄扁桃适宜分布区被包含在该图层中）。模型预测的适宜度值（也称物种存在概率“P”）为连续栅格数据，数值范围为 0～1。使用ArcGIS中的Reclassify工具和人工分级方法，将适宜度值分为4个等级：不适宜（0 ≤ P < 0.15）、低适宜（0.15 ≤ P < 0.33）、中等适宜（0.33 ≤ P < 0.66）和高适宜（P > 0.66），以地图形式导出。

由于适宜分布区域的边缘不规则，用适宜分布区几何中心的运动来表示其变化更为合适，也更具代表性。对长柄扁桃适宜分布区几何中心的提取涉及ArcGIS中的3个工具：Raster Calculator工具、Raster Domain工具和Mean Center工具。首先，使用Raster Calculator工具来切取栅格图层，只保留适宜分布区这部分栅格数据。然后，使用Raster Domain工具将适宜分布区域栅格数据转换为平面图形。最后，使用Mean Center工具提取出平面图形的几何中心，即长柄扁桃适宜分布区的几何中心。

预测结果表明：基准气候情景（1970—2000）下，长柄扁桃的适宜分布区集中在中国、蒙古和俄罗斯三国境内（图1，即长柄扁桃集中适宜分布区，分布范围为97.23°～122.62° E和35.49°～54.91° N）。在我国，长柄扁桃总适宜面积为94.54 × 10⁴ km²，主要分布在甘肃、青海、宁夏、内蒙古、陕西、山西和河北等地，中、高适宜分布区主要分布在内蒙古，陕西和山西北部有少量高适宜分布区。长柄扁桃在蒙古国的总适宜面积为88.35 × 10⁴ km²，适宜分布区主要分布在库苏古尔、后杭爱和巴彦洪戈尔往东的地区。在俄罗斯，长柄扁桃的总适宜面积为13.33 × 10⁴ km²，集中在伊尔库兹克州、布里亚特共和国和赤塔州南部。

Figure 1.  Global suitable distribution (A) and concentrated suitable distribution region (B, C) of Prunus pedunculata under the reference climate scenario

图2为未来气候情景中（2021—2100）长柄扁桃的适宜分布区面积的变化。2041年前，长柄扁桃集中适宜分布区面积在未来4个气候情景中都呈现增加趋势；2041年后，情景SSP1-2.6中的适宜面积开始减少，情景SSP2-4.5中适宜面积基本保持稳定，情景SSP3-7.0和SSP5-8.5中适宜面积持续增加，并且情景SSP3-7.0中的适宜面积在2061年前超过情景SSP5-8.5，随后持续增加。

Figure 2.  Changes of the concentrated suitable distribution area of Prunus pedunculata under different climate scenarios

在我国境内，4个情景中长柄扁桃适宜面积皆呈波浪状，这表明长柄扁桃适宜面积在我国并不是单一的增加或者减少，应对气候的变化表现出多样性。

在蒙古境内，2021年前，情景SSP1-2.6中长柄扁桃的适宜面积呈增加趋势，2021年后基本保持稳定，变化较小；其它情景中，适宜面积表现为增长趋势，在2061年前，排放量越大的情景适宜面积增加量也越高，而在2061年后，情景SSP3-7.0中的适宜面积反超情景SSP5-8.5。

在俄罗斯，情景SSP2-4.5、SSP3-7.0和SSP5-8.5中，长柄扁桃适宜面积在2061年前，都呈增加趋势；2061年后，情景SSP3-7.0中的适宜面积超过情景SSP5-8.5，而情景SSP2-4.5中适宜面积开始减少。与其他情景不同，SSP1-2.6中，适宜面积于2021年前呈增加趋势，随后减少，2041年前后，先减少后增加，而2061年前后类似于2021年，是另外一个峰值点。

19个气候变量经过两两之间的相关性筛选，最终确定5个气候变量（bio1、bio7、bio13、bio15和bio17）用于编译MaxEnt模型的运行。

刀切法（Jackknife）测试用于分析气候变量对预测结果的影响，以确定每个气候因素的重要性。当单独使用气候变量时，增益最大变量是降水季节性变化（bio15）（图3），表明该变量是预测过程中限制长柄扁桃分布的最关键变量。

Figure 3.  Jackknife test for evaluating the relative importance of bioclimatic variables

图4显示了基准气候情景和未来4种气候情景下长柄扁桃集中适宜分布区的几何中心。在情景SSP1-2.6中，2021年、2061年和2081年的几何中心位于基准气候情景几何中心的西北方，2041年几何中心位于基准气候情景几何中心的东北方。在情景SSP3-7.0中，2021年几何中心位于基准气候情景几何中心的东北方，其余时间段的几何中心都位于基准气候情景几何中心的西北方。在情景SSP2-4.5和SSP5-8.5中，所有时间段的几何中心全部位于基准气候情景几何中心的西北方。结果表明，未来长柄扁桃集中适宜分区的几何中心整体有往西北方移动的趋势。

Figure 4.  Geometric centers of the concentrated suitable distribution regions of Prunus pedunculata under different climate scenarios

将几何中心的迁移进行量化，以基准气候情景下的几何中心为原点，距离原点最远的几何中心的迁移率表示为1，其余几何中心和原点间距离与最远几何中心和原点间距离之比就是这些几何中心的迁移率（表1）。几何中心迁移率最高的气候情景为2061和2081年的SSP5-8.5情景，迁移率为1。几何中心迁移率最低的气候情景为2021年的SSP3-7.0情景，迁移率为0.31。4种气候情景的平均迁移率分别为0.39（SSP1-2.6）、0.46（SSP2-4.5）、0.54（SSP3-7.0）和0.79（SSP5-8.5）。随着排放情景的升级，几何中心的平均迁移率也在增加，表明高排放情景下长柄扁桃集中适宜分布区的变化比低排放情景更活跃。

MaxEnT模型运用受试者工作曲线（ROC）对预测结果进行精确性验证。在该方法中，使用曲线下的面积，即AUC（area under curve）值（范围为0～1）来判断模型的预测精度。当AUC ≤ 0.7时，表明预测精度较差；当0.8 ≤ AUC < 0.9时，表明预测精度较好；当AUC > 0.9时，表明预测精度很好^{[11, 14]}。本研究中，模型的AUC值为0.982，表明预测结果精度高，预测效果好。

内蒙古和蒙古国的中东部是长柄扁桃集中适宜分布区。适宜分布区纬度跨度很大，从35.49°N至54.91° N，总体上长柄扁桃属于高纬度经济林树种。随着纬度的增加往往伴随着无霜期短、年均温度低、降水量少等低温、干旱的环境，导致大部分经济林树种无法发展^[20-21]。而长柄扁桃在高纬度的分布，为高纬度地区提供了一个适宜的树种。

降水季节性变化（bio15）是限制长柄扁桃分布的最关键变量。温度和降水是限制植物生长的2个关键因素。水对植物尤其重要，它不仅对植物的光合作用有重要影响，而且还决定着植物的生长状况，是限制植物分布的一个重要因素^[22]。使用MaxEnt模型预测唇萼薄荷（Mentha pulegium L.）时，结果表明最冷季度的降水量（bio19）是限制其分布的最关键变量^[23]；对五叶木通（Akebia quinata （Thunb.） Decne.）的分布预测表明，限制其分布的最关键变量是年降水量（bio12）^[15]。降水季节性变化（bio15）、最冷季度的降水量（bio19）和年降水量（bio12）都是与降水相关的变量，表明降水条件对长柄扁桃、薄荷和五叶木通等物种的分布限制优先级高于温度条件。而温度变量最冷季平均温度（bio11）是限制三叶木通（Akebia trifoliata （Thunb.） Koidz.）的最关键变量^[16]，限制竹叶花椒（Zanthoxylum armatum DC.）分布的关键环境因子的是1月份最高气温^[24]，对于三叶木通和竹叶花椒等物种来说，温度条件对分布的限制要高于降水条件。由此可知，尽管限制不同物种分布区的变量不同，但都或与温度，或与降水息息相关。

不同物种针对气候变化做出的响应不同，尤其是随着CO₂排放体量增加的情况下，不同的物种适宜面积既有增加也有减少。本研究中，长柄扁桃适宜面积总体呈增加趋势，具体情况是：2041年之前，长柄扁桃集中适宜分布区面积在未来4个情景中都为增加趋势；2041年之后，情景SSP1-2.6中的适宜面积开始减少，情景SSP2-4.5中适宜面积基本保持稳定，情景SSP3-7.0和SSP5-8.5中适宜面积持续增加，并且情景SSP3-7.0中的适宜面积在2061年之前就将超过情景SSP5-8.5，随后继续增加。同样，紫斑牡丹（Paeonia rockii （S. G. Haw & Lauener） T. Hong & J. J. Li）的适宜面积在排放情景RCP2.6和排放情景RCP8.5中都有所增加^[22]；而对松茸（Tricholoma matsutake （lto et lmai） Singer）的研究表明，在RCP8.5情景下，其适宜面积将大大减少，适宜分布区甚至会支离破碎^[25]。由于气候变化，物种的适宜分布区也将发生变化^[26-28]。物种分布区的扩张或收缩、物种迁移、甚至是分布区的破碎化是物种分布区几何中心变化的原因所在^{[15-16, 18]}。未来气候情景中，长柄扁桃集中适宜分布区的几何中心有向西北移动的趋势，而五叶木通在东亚的适宜分布区的几何中心随着排放情景的升级而持续向东北移动^[15]，三叶木通适宜分布区的几何中心将向北移动^[16]。随着CO₂排放体量的增加，全球气候的变化不管是提升还是降低适宜分布区的面积，都会使物种适宜分布区的几何中心迁移，也就是原本适宜的地区不再适宜原生物种，这在很大程度上会威胁物种的生存。因此，控制碳排放，降低人类活动对气候的影响是非常必要的。

长柄扁桃适宜分布区主要集中在中国、蒙古和俄罗斯，范围为97.23°～122.62° E和35.49°～54.91° N。降水季节性变化（bio15）是限制长柄扁桃分布的最关键气候变量。在未来的4种气候情景下，长柄扁桃集中适宜分布区有向西北方迁移的趋势。随着排放情景的升级，几何中心的平均迁移率也在增加，高排放情景下集中适宜分布区的变化比低排放情景更活跃。本研究预测长柄扁桃分布区响应气候变化的规律，对其保护和开发利用具有重要意义。