北温带地区小叶锦鸡儿的地理分布模拟

杨光，于世永，张鑫，裴智勇，郑文胜，李贺

(1.江苏师范大学，江苏 徐州 221116；2.呼和浩特市环境监测中心站，内蒙古 呼和浩特 010030；3.二连浩特市林业综合工作站，内蒙古 二连浩特 011100；4.镶黄旗亚力盖图机械林场，内蒙古 锡林浩特 013250)

近百年来，在北温带区发生了木本植物盖度和密度增加的现象，称之为灌丛化[1-3]。灌丛化的发生一般认为由全球气候变化、过度放牧、氮沉降和火烧抑制等因素引起。灌丛化的发生对生态系统的结构和功能产生了显著影响，例如灌丛化会导致景观破碎化，原本连续的景观变成了木本植物与草本植物的镶嵌景观，进一步导致群落生物结构和多样性的改变，以及能量流动和土壤养分循环的变化，从而影响区域内人类的生存和发展[3，4-6]。早在1990年，周道玮就发现了我国内蒙古自治区锦鸡儿属植物有增加的现象[1]，后来熊小刚等[7]研究分析了锡林郭勒盟地区锦鸡儿属植物引起的草原灌丛化现象，并发现灌丛化的发生对植物群落和养分循环产生了显著影响。然而，灌丛化对生态系统是有利还是有弊，还未形成定论，例如灌丛化虽然能起到一定的防风固沙作用，但也有降低生物多样性的风险，导致草原退化[8]。当前我国有很多区域进行锦鸡儿属植物的人工种植，人为增加了灌木的分布面积，虽然取得了一定的生态和经济效益，但对于生态系统其他方面的影响还有诸多不确定性。目前，全球预计约有10%～20%的干旱、半干旱区域正在经历草原灌丛化[9]，我国内蒙古自治区有5.1×106hm2的草原发生了灌丛化[10]，如果我国温带草原地区灌丛化持续发展，将对草业和牧业的可持续发展带来严峻挑战[11-12]。未来，是遏制小叶锦鸡儿的扩张还是继续加大其种植还需要进一步研究，但无论如何，准确模拟和预测锦鸡儿属植物的地理分布格局，探明其控制因素都是非常重要的一环。此外，当前我国对于草地生态系统生物量的估算以及生态服务等方面的评价，都没有考虑灌木植物，而随着灌丛化进程的不断发展，灌木分布的动态变化以及准确评估草地系统的结构与功能都需要准确定量化模拟灌木的分布，而目前虽然对锦鸡儿属植物的地理分布已有了准确的描述，例如，赵一之[13]、周道玮等[14]的研究，但是在定量化模拟方面还缺乏相关的研究。

小叶锦鸡儿(Caraganamicrophylla)是北温带地区常见的落叶灌木[1]，也是引起我国北温带地区草原灌丛化的主要植物[15-16]，主要分布在俄罗斯贝加尔湖南部和达乌里地区，蒙古国的东部和东北部，以及我国内蒙古自治区大部分地区、松辽平原西部和华北山地等区域，一般呈斑块状散生分布[17-18]。气候是决定植被地理分布的主要因素之一，而全球气候的变化将对植物的分布带来巨大的挑战[8]，已有研究表明气候变暖将导致最低温度升高、冻害减少，降低灌木的死亡率，为灌木扩张和定居创造了条件[19]，为此，阐明小叶锦鸡儿地理分布格局及其与气候因子的关系是定量化研究小叶锦鸡儿分布的关键。在模拟过程中，模型的选择至关重要，国内外学者们基于气候因子使用模型对植物分布的模拟做了很多研究。例如，刘然[20]通过最大熵模型模拟分析了我国冷杉属植物的地理分布数据；Lutz[21]在2010年对美国加利福尼亚州中部约塞米蒂国家公园内松科、柏科和山毛榉科等共17个物种的分布随气候条件的变化进行了研究。然而，对于小叶锦鸡儿的地理分布格局模拟还存在较多的不确定性，一方面以往的研究数据量较小，并且在全球气候变化背景下，小叶锦鸡儿的分布发生了变化，需要数据的更新；另一方面现大部分物种分布模型主要以乔木植物为主，对于灌木植物的研究较少，虽然王娟和倪健[22]也曾运用BioClim模型对我国北温带5种锦鸡儿植物的地理分布进行了模拟，但是数据量的大小以及单一模型的模拟存在较大不确定性，使得通过更为详尽的数据，结合多种模型来模拟小叶锦鸡儿的分布，并进行比较筛选显得尤为必要。

为了准确模拟北温带区小叶锦鸡儿的地理分布格局，本文以小叶锦鸡儿为研究对象，利用R语言及其相关数据包为分析工具，运用多种模型探讨北温带区小叶锦鸡儿植物地理分布格局及与气候因子之间的关系，筛选定量化小叶锦鸡儿的分布模型，为干旱半干旱区灌木植物的动态变化模拟以及草原生态系统的评估、保护和管理以及可持续发展提供必要的支持。

1 材料和方法

1.1 小叶锦鸡儿地理分布数据

本研究的研究范围为北温带小叶锦鸡儿的主要分布区，地理坐标在N 34.41°～52.56°，E 101.40°～124.86°。首先对获取的小叶锦鸡儿地理分布数据进行筛选，按照史料记载，确定北温带区小叶锦鸡儿的地理分布范围边界；其次通过Google Earth软件精确检查每个数据点，去除植物园、城市绿化等非自然分布数据，不包括非自然分布，最终确定有效小叶锦鸡儿地理分布数据：(1)野外实地考察调研数据共122个；(2)全球生物多样性信息数据库(http：//www.gbif.org/)共268个；(3)文献、标本馆和史料记载数据(包括地方草原站及自然资源部门的历史数据)共146个，总计有效分布数据点536个(图1)。

图1 北温带区小叶锦鸡儿分布样点图Fig.1 Distribution of C.microphyllain North temperate zone areas

1.2 气候数据

本研究采用的气候数据来自2017年公布的WorldClim全球气候数据集Version 2(http://worldclim.org/)[3]，所选分辨率为30弧度秒，对应的栅格大小为1 km×1 km，包含了年平均气温、昼夜温差均值、等温性、温度季节性变化、年降水量等在内的19个生物气候变量(Bio1～Bio19)，并计算北温带区小叶锦鸡儿分布范围内各个气候因子的均值、变异系数等指标。

1.3 数据分析

分析19个变量之间的相关性，用主成分分析法研究北温带区小叶锦鸡儿分布与气候因子的关系，将筛选出的关键变量纳入到后续的模型分析中。结合地理分布和气候数据，用人工神经网络模型(Artificial Neural Network，ANN)、决策树模型(Classification Tree Analysis,CTA)、广义增强模型(Generalized Boosting Model,GBM)、最大熵模型(Maximum Entropy Model，MaxEnt)以及随机森林模型(Random Forest，RF)等5种模型模拟北温带区小叶锦鸡儿的地理分布格局。

采用受试者工作特征(Receiver Operator Characteristic,ROC)曲线来评价模型的模拟效果。ROC曲线是以假正率或假阳性率(在已知物种发生点中预测不适合发生点的比率)为横坐标，以灵敏度或真阳性率(在已知物种发生点中预测适合发生点的比率)为纵坐标所形成的曲线，ROC曲线与横坐标构成的封闭面积就是AUC值(Area Under Curve)，AUC值0.5～1.0，以AUC值的大小来判定模型预测的效果，AUC值在0.9～1.0表明模型模拟结果非常好[24]。此外，利用实际考察的122个小叶锦鸡儿地理分布数据点与5种模型模拟的结果进行配对T检验分析，通过平均数间差异的显著性水平以进一步判定模型模拟结果的精度。

本研究所有分析均使用R软件完成(R Core team)[25]。

2 结果与分析

2.1 小叶锦鸡儿地理分布与气候因子的关系

北温带小叶锦鸡儿分布区的年平均气温为-5.20～14.58℃，平均值为4.08℃；年降水量为94～689 mm，平均值为358.89 mm(表1)。各个气候因子的变异系数为9%～85%，其中平均日较差(Bio2)、极端最高温度(Bio5)、降水量季节性变化(Bio15)、等温性(Bio3)、年温度变化范围(Bio7)、最湿季平均温度(Bio8)等气候因子的变异系数较小，说明小叶锦鸡儿的分布对这些气候因子的依赖性较大，是小叶锦鸡儿分布的限制因子；而年平均气温(Bio1)、最干季平均温度(Bio9)、最冷季平均温度(Bio11)、最干季降水量(Bio17)等气候因子的变异系数则较大，说明小叶锦鸡儿的分布对这些气候因子的依懒性不强。总而言之，小叶锦鸡儿分布于气候较为寒冷和干燥的区域，而其分布主要受温度和降水的日、季节变化以及生长季的温湿度等的影响，而对冬季的温湿度变化不太敏感。

表1 气候变量说明及描述性统计结果

平均日较差(Bio2)和降水量季节性变化(Bio15)等限制因子和其他气候因子的相关性不高，而年温度变化范围(Bio7)、等温性(Bio3)和极端最高温(Bio5)等限制因子则与年平均气温(Bio1)、最冷季平均温度(Bio11)等有较强相关性(图2A)。19个气候因子的主成分分析，前三轴的累积解释率达到了84.26%，解释率较高(图2B，表2)。其中第1轴中年平均气温(Bio1)、温度季节性变化(Bio4)、极端最低温(Bio6)、年度温度变化范围(Bio7)、最干季平均温度(Bio9)、最冷季平均温度(Bio11)等温度相关的气候因子的载荷较大，因此第一主分量主要反映热量因素；在第2轴中，最干月降水量(Bio14)和最冷季降水量(Bio19)的载荷值较大，第二主分量主要反映极端降水情况；在第三轴中，降水量季节性变化(Bio15)的载荷较大，因此第三主分量反映季节降水情况。综合变异系数和相关性分析，上述主成分分析3轴所筛选的10个气候变量均为限制因子或与限制因子具有较强的相关性，所以这10个气候因子作为影响小叶锦鸡儿分布的关键气候因子纳入小叶锦鸡儿的分布模拟。

图2 小叶锦鸡儿地理分布与各生物气候变量的主成分分析图Fig.2 Principal component analysis of geographical distribution and bioclimatic variables of C.microphylla

表2 小叶锦鸡儿分布区气候因子的主成分分析

2.2 小叶锦鸡儿分布的模型模拟与评价

结合主成分分析筛选出的10个气候因子，选择5种模型模拟北温带区小叶锦鸡儿的地理分布格局(图3)。通过ROC曲线分析(图4)，5种模型的AUC值均在0.8以上，其中最大熵模型、随机森林模型和广义增强模型的AUC值达到0.9以上，说明这3个模型对于小叶锦鸡儿分布的模拟效果较优。对比小叶锦鸡儿的采样数据点，3种模型的模拟结果基本反映了小叶锦鸡儿的地理分布状况，小叶锦鸡儿植物主要分布在甘肃省中部、宁夏回族自治区中部和南部、陕西省北部、山西省中部和北部、河北省西北部、内蒙古自治区中南部、蒙古国东部等地区，这些区域内的小叶锦鸡儿植物分布概率达到了60%～80%，属于北温带区小叶锦鸡儿的适宜分布区。但是，3种模型的模拟结果存在一些差异，最大熵模型的模拟结果较为细腻，各个概率范围之间的替代性明显，既从核心分布区到不适宜区依次展开,随机森林模型和广义增强模型的模拟结果层次也较为明显，而决策树模型和人工神经网络模型的模拟结果各个概率的范围替代不明显。经定量化分析各分布概率下小叶锦鸡儿的生境面积(表3)，结合实际采样点与模型结果的配对T检验发现(表4)，广义增强模型、随机森林模型、决策树模型和人工神经网络模型在模拟最适生境上都有不同程度的高估现象，尤其是决策树模型模拟结果明显与实际不符。因此，本研究认为最大熵模型模拟小叶锦鸡儿植物的地理分布情况更为准确(P<0.05)。

图3 5种模型下北温带区小叶锦鸡儿当前地理分布情况Fig.3 Current geographical distribution of C.microphylla in northern temperate zone under five models注：A：最大熵模型；B：随机森林模型；C：广义增强模型；D：决策树树模型；E：人工神经网络模型

表3 小叶锦鸡儿不同分布概率下生境面积的模型模拟结果

图4 模型评估的ROC曲线及AUC值Fig.4 ROC curve and AUC value for precision evaluation of model注：A：人工神经网络模型；B：决策树模型；C：广义增强模型；D：最大熵模型；E：随机森林模型

表4 5种模型模拟结果与实际采样点的配对T检验

3 讨论

影响北温带小叶锦鸡儿分布的19个气候变量中，小叶锦鸡儿植物的分布受年平均气温(Bio1)、平均日较差(Bio2)、温度季节性变化(Bio4)、极端最低温(Bio6)、年度温度变化范围(Bio7)、最干季平均温度(Bio9)、最冷季平均温度(Bio11)、最干月降水量(Bio14)、降水量季节性变化(Bio15)以及最冷季降水量(Bio19)这10个气候因子的控制，这反映出小叶锦鸡儿植物对温度和降水的要求较高[15]，并且对热量、极端降水以及季节降水变化情况较为敏感[26]，容易受到气候变化的影响[27]。前人也揭示了一些影响小叶锦鸡儿地理分布的原因。在气候变化影响方面，最直接的原因是土壤表层水分含量的降低，导致根系分布较浅的植物(尤其是禾本科草本植物)竞争能力下降[26]，而有利于小叶锦鸡儿这类灌木植物生长，使得其广泛分布。其次，大气中二氧化碳浓度和温度的上升，更有利于灌木植物如小叶锦鸡儿而不利于一些草本植物的生长[28-30]；从植物本身分析，叶是植物散失水分的主要器官，而生活在干旱区的植物通过叶厚度、面积、被毛、叶色、革质化以及鳞片等形态变异减少水分散失[15]，使得小叶锦鸡儿这类灌木植物对干旱环境的适应能力较好，从而导致在北温带区一些木本植物盖度和密度增加的现象。

从模型运行结果来看，5种模型的AUC值都超过了0.8。决策树模型具有强大的分析和统计功能，是对特征数据实施统计建模分析的一种优秀方法[31]；人工神经网络模型是模拟生物神经网络结构组成和功能原理的模型，可用于模拟输入和输出变量之间的复杂关系[32-33]；随机森林模型在没有显著提高模型运算量的前提下提升了运算精度，它对变量的共线性不敏感，能够很好地处理数据的缺失与样本量的不足[34]；最大熵模型作为研究物种分布模型的机器语言，具有运行速度快、操作简单、受样本偏差影响小、模拟结果精度高等优点[35]；广义增强模型不同于其他传统方法的是，其使用Boosting的原理与方法，自适应联合大量相关简单回归树，以优化模型预测能力，使结果形成最优化模型[36]。但结合实际采样点数据来看，人工神经网络模型与决策树模型的结果比较粗糙，与实际情况吻合度相对较低；虽然随机森林模型以及广义增强模型结果较为细腻，符合植物的实际地理分布规律，但对小叶锦鸡儿地理分布的模拟存在不同程度高估的现象，而最大熵模型对北温带区小叶锦鸡儿地理分布模拟效果更为准确。这一结果与分布点数据的样本量有很大关系，虽然本研究尽可能地收集了小叶锦鸡儿的分布点，但仍有较大提升空间。相较于其他模型，最大熵模型具有在数据有限的情况下，预测结果精度高等特点[36]。本研究结果表明，除乔木植物之外，最大熵模型在灌木植物的地理分布模拟方面也很适用。

此外，以最大熵模型模拟结果为参考，北温带区小叶锦鸡儿植物在甘肃省中部、宁夏回族自治区中部和南部、陕西省北部、山西省中部和北部、河北省西北部、内蒙古自治区南部以及蒙古国中部和东部地区的小叶锦鸡儿植物分布概率达到了60%～80%以上，成为北温带区小叶锦鸡儿植物适宜生存的核心地带，在全球气候变化和草原灌丛化发生的大背景下，这些区域的监测需得到持续关注。

本研究所涉及的分布点均为小叶锦鸡儿的自然分布点，没有纳入非自然研究点。小叶锦鸡儿多年来是干旱和半干旱草原和荒漠区造林和生态恢复等措施的主要物种[37]。为了防风固沙和提高森林覆盖率，很多地方都种植了锦鸡儿属植物[38]，如小叶锦鸡儿、柠条锦鸡儿等灌木植物，虽然由于人为种植导致小叶锦鸡儿这一物种在我国有了一定数量的非自然分布，但是，在自然条件下小叶锦鸡儿和其它锦鸡儿属植物的分布面积和种群密度也在不断的增加，即在全球变化背景下草原灌丛化在不断加剧[2]。草原灌丛化对生态系统的影响是多方面的，其生态学效应也是多样的，但对于生态系统的影响是有利还是有弊目前在国际上并未形成定论[39]，虽然种植诸如小叶锦鸡儿等灌木在防风固沙等方面有一定的生态学效应，但在种植和管理过程中，部分地区因对其营养状况了解不充分，利用和管理措施不当，导致小叶锦鸡儿灌丛木质化严重、植被退化，使其生态和经济效益并未得到充分发挥[40]。在此条件之下，未来是继续扩大锦鸡儿属植物的人工种植还是遏制其分布还需要进一步研究。而本研究对小叶锦鸡儿地理分布的模拟及其与气候因子之间的关系的相关研究结果，作为一个基础将有助于后续研究和生态政策的制定，这也将为后续草原发生灌丛化后，草地生物量的准确计算和生态服务评价提供相关的支持。

本研究在前人相关工作的基础上，通过扩展和更新数据，采用多种方法在小叶锦鸡儿地理分布格局的定量化上进行了探索，是前人研究结果的扩展、补充和利用，后续还需继续加强。例如，本研究仅选择了气候变量作为模型模拟了北温带小叶锦鸡儿地理分布的控制因素，未对地形、土壤等因素进行考虑，并且气候因子的精确度也有待提高；另一方面，模型的模拟本身也存在一定的不确定性，需要不断地优化。此外，本研究仅研究了小叶锦鸡儿的分布，希望在此基础上能对其他锦鸡儿属植物乃至其他灌木植物地理分布进行定量化研究，并对未来气候变化下其动态变化做出预测。

4 结论

本研究通过大量的实地调查和数据收集，研究了影响北温带小叶锦鸡儿地理分布的主要气候因子，并运用多种模型定量化模拟了北温带小叶锦鸡儿的地理分布并对模型进行了评价。研究结果显示：(1)小叶锦鸡儿分布于气候较为寒冷和干燥的区域，其分布主要受温度和降水的日、季节变化以及生长季的温湿度等影响，对冬季温湿度变化不太敏感，而对热量、极端降水以及季节降水变化情况较为敏感；(2)最大熵模型、随机森林模型和广义增强模型对小叶锦鸡儿分布的模拟结果较好(AUC>0.9)，可以在灌木植物中应用，而最大熵模型的模拟结果最优；(3)模型结果准确反映出小叶锦鸡儿在甘肃省中部、宁夏回族自治区和内蒙古自治区的中部和南部以及蒙古国的中部和东部等地区分布较为集中的实际分布格局，在全球气候变化和草原灌丛化加剧的背景下，这些区域需要得到持续关注和监测。本研究将为小叶锦鸡儿分布的动态变化提供基础数据和模型支持。此外，定量化模拟小叶锦鸡儿地理分布格局的相关结果将为后期模拟精准计算灌丛化草原系统的生产力和生态服务评价提供关键参数。