基于4种生态位模型的金钱松潜在适生区预测

王国峥,耿其芳,肖孟阳,张孟源,张云燕,王中生

南京大学生命科学学院,南京 210023

金钱松(Pseudolarixamabilis(Nelson) Rehd.)是我国特有单属种松科植物,也是著名的第三纪孑遗植物,现主要分布于江苏南部、安徽南部、福建北部、浙江、江西、湖南、湖北利川至四川万县交界地区[1]。最早的化石记录出现于西伯利亚东北部的晚白垩纪地层中,始新世时出现于挪威斯匹次卑尔根岛西部,美国、欧洲、亚洲中部以及日本在第三纪的不同时期都曾发现有该属化石分布[2]。据化石资料推测,金钱松属可能起源于白垩纪中世纪晚期,第三纪时曾发展为北半球中高纬度广布属,其中在日本和欧洲一直维持至更新世第一次冰期降临[3]。在新生代末期冰期,它从欧洲发生地向我国南方迁移,最后定居于长江流域,形成现今的地理分布格局[4]。

由于金钱松特殊的进化及迁移历史,其在松科植物系统发育研究中备受关注,同样对古生态以及古气候的研究具有重要意义[5]。从其化石记录可发现金钱松的分布与地质变化过程密切相关,地质变化导致气候剧烈变化,金钱松分布范围随之发生改变。一些学者已从金钱松植物学特征[5]、种群生态学[6]、以及遗传多样性[7- 8]等方面开展了相关研究,其中高燕会[7]和刘俊[8]等利用分子标记发现金钱松天然种群具有较高的遗传多样性,认为遗传因素不是金钱松的致濒因子,其呈狭域分布的主要原因可能来自外界自然因素。

近年来生态位模型在预测、解析物种分布范围及其变化等方面应用广泛,利用物种已知的分布数据和环境变量,根据一定的算法运算来构建模型,并将运算结果投射到不同的时间和空间中以预测物种的潜在适生区[9]。目前较为常用的生态位模型较多,每种模型都可以独立预测出潜在适生区,但各模型都存在一定的偏好性[10],本文参考Miguel B. Araújo[10]的方法,采用了4种基于不同算法的模型(GARP、Bioclim、Domain和Maxent)预测金钱松的全国潜在适生区,利用集合预测系统思想,综合多种模型(3种以上)预测结果,尽可能降低由经验性选择一种模型导致的假阴性或假阳性影响,同时一种模型的缺陷可能被另一种模型所弥补,从而提高预测的科学性。GARP是基于遗传算法的规则组合进行局域环境空间建模[11]；Bioclim基于矩形框架模型,每一种环境因子被抽象成一个矩形框,绝大多数分布点(95%以上点)位于这些框内,多个矩形框构成限制范围,若某点位于这个范围内,便可以认为其适宜该物种生存[12]；Domain模型基于Gower算法,通过计算所有已知分布点之间的最大相似度,来评估预测点对于环境变量的适宜程度[13]；Maxent模型通过物种的分布数据和环境数据,找出物种分布规律的最大熵,从而对物种的分布进行预测[14- 20]。

1 数据获取和研究方法

1.1 金钱松样本分布数据

金钱松的地理分布位点主要通过中国数字植物标本馆(http://www.cvh.org.cn/)以及中国科学院北京植物所标本馆(http://pe.ibcas.ac.cn/)查阅,共获取44个分布记录数据,准确到县。用Google Earth 转换出分布点的经纬度坐标,输入Excel保存成.csv格式。

1.2 环境变量数据及预处理

环境变量共22个,包括19个气候因子和3个地形因子,19个气候因子来源于Worldclim (http://www.worldclim.org/),包括末次盛冰期(Last Glacial Maximum,LGM)、当前(Current)和2070年的数据,分辨率2.5 arc-minutes,将环境因子数据用ArcGIS 10.2统一转化为ASCII格式以便使用。

由于环境变量之间具有多重共线性,会导致预测分布过度拟合[21],因此本文参考张琴[15]对环境共线性诊断的方法,对环境因子进行Spearman秩相关分析,筛选出相关性较低的环境变量(Spearman系数<0.75),并在Spearman系数大于0.75的变量中选出具有生态学重要意义的环境变量,共得到10个变量(表1)。

表1 环境变量

1.3 相关软件

4种生态位模型所用软件：基于GARP模型预测的Desktopgarp(Version1.1.6),基于Bioclim和Domain模型预测DIVA-GIS(Version 7.5)；基于最大熵模型预测软件Maxent(Version 3.3.3)。

数据处理软件有：Excel用于分布数据记录以及数据格式转换,DIVA-GIS(Version 7.5)用于图层格式转换,ArcGIS(Version 10.2)用于图层的处理和数据的格式转换。

1.4 模型预测与检验方法

1.4.1GARP模型预测

将物种已知分布点数据转换成Desktopgarp模型软件支持的格式,通过模型软件中 Upload Data Points功能加载,选择70%的分布数据作为训练数据,其余数据作为检验数据[14]。环境数据.asc文件通过Desktopgarp中的Dataset Manager处理转换成Desktopgarp可识别的格式.raw,以数据集的形式加载到软件中。

1.4.2Maxent模型预测

将已知分布点数据文件格式在Excel中转为.csv格式文件导出,通过Browse加载到Maxent软件中。再将10个环境变量ASCⅡ文件的环境数据通过Browse加载到Maxent软件中。并进一步利用Maxent软件对金钱松在末次盛冰期(LGM)气候和未来气候(2070年)条件下进行模拟分布[22],再用ArcGIS 10.2 SDMs Tool模拟金钱松从末次盛冰期到当前的分布迁移变化[23]。

1.4.3Bioclim和Domain模型分析

在DIVA-GIS中首先添加.shp格式的训练数据集,再将10个环境变量文件ASCⅡ转换成.grd格式,生成stack数据集。在Modeling-Bioclim/Domain模块中添加stack格式的环境数据集,依次进行Bioclim和Domain模型的预测[15]。

1.5 模型预测结果检验

受试者工作特征曲线(Receiver Operating Characteristic,ROC)和Kappa统计量广泛用于模型精度的检验。本文参考张琴[15]的方法,利用DIVA-GIS软件,将10组不同的训练和测试数据用4种模型预测得到的40个预测分布图生成stack数据集,再导入软件创建的验证点.shp文件,最后得出各模型的AUC值和Kappa值。

受试者工作特征曲线下的面积即为AUC值(Area under recriver operating characteristic curve,AUC),AUC的数值范围为0.5—1.0,“0.5”说明预测结果为随机分布,预测值越大表示该模型预测结果越精确。一致性检验统计量(Kappa)是对两种方法测定结果的一致部分进行检验,取值范围是[-1,1],“-1”说明两次判断的结果完全一致,“1”说明两次判断的结果完全不一致,值越大一致性越好,通用的经验法则是Kappa大于0.75表示一致性好,小于0.4表示一致性差[24]。

1.6 影响适生分布的环境因子分析

利用Maxent软件中的“刀切法”(Do jackknife to measure variable importance)可以显示出环各个境变量对预测金钱松分布的贡献。其他3种模型(GARP、Bioclim、Domain)预测软件不具备分析环境因子对预测分布贡献的功能,因此将GARP、Bioclim、Domain等3种模型基于10个环境因子预测的适生区作为对照组,分别以“刀切法”分析影响较大的三个环境因子为单因子,或作为共同因子来预测金钱松的适生区,进行预测结果的比较。

2 结果与分析

2.1 4种模型预测结果

从GARP模型预测图中(图1)可以看出金钱松在我国的适生区集中在江苏、浙江、安徽、江西、河南、湖北、湖南、贵州、重庆、四川东部、陕西南部以及云南北部地区,集中分布在24.43°—33.35°N和106.41°—123.42°E之间,而搜集的样本数据在27.27°—34.05°N和108.41°—121.80°E之间。

Maxent预测分布图(图1)显示金钱松高适区主要集中在江苏南部、安徽南部、湖北南部、江西北部以及浙江西北部。绿色代表中度适生区,依次扩散开来,最大范围向北达江苏省潥阳市,向南至湖南衡山磨镜台,最西可达重庆市万州区,最东可达浙江东天目山。

通过DIVA-GIS的Modeling功能运算得出的Bioclim模型预测图(图1)和Domain模型预测图(图1),可观察到代表最适区的红色区域在这两种模型中占较小比例,Bioclim预测分布范围较小,显示金钱松可在浙江、湖南、安徽南部、湖北南部以及江西西北部高度适生。Domain模型预测结果分布范围较大,总体与GARP模型预测分布范围接近,高适与最适分布范围与Maxent相似,但在河南新乡、南阳以及云南昆明出现最适分布区。

2.2 不同模型预测准确性评估

本文采用10组训练数据和测试数据对4个模型进行ROC曲线分析和Kappa一致性检验。GARP模型AUC平均值为0.922,Maxent模型AUC平均值为0.980,Bioclim模型AUC平均值为0.957,Domain模型的AUC平均值为0.940(见表2),均远远超过随机模型(AUC=0.5)。并且4种模型的Kappa均值都大于0.75,表明4种模型预测的一致性均较显著,预测精度较高。4种模型AUC值的标准偏差排序为Bioclim > Domain > Maxent >GARP。

表2 4种模型的AUC值和Kappa值

2.3 影响适生区的环境因子分析

各个境变量对预测金钱松分布的贡献如图2,在10个环境因子中,年均温(Bio1)对于金钱松分布的贡献最大,年降水量(Bio12)、最湿月降水量(Bio13)和最热月极高温度(Bio5)对金钱松分布的影响也比较大,温度日较差(Bio2)和海拔对金钱松分布的影响比较小。

图1 4种生态位模型的金钱松适生区预测图Fig.1 Potential distribution of P.amabilis based on four ecological niche models

图2 刀切法检测生物气候变量对分布的重要性Fig.2 Import of environmental variables for prediction based on jackknife test

从单因子预测结果来看(图3),3种模型基于年均温(Bio1)的模拟结果与各自预测适生区最为接近,但范围较适生区均偏大。横向对比3种模型分别基于Bio1、Bio12、Bio13等3种单因子的模拟结果,基于年均温(Bio1)的预测结果相对较为接近,说明年均温(Bio1)对这3种模型的预测结果也都有显著的影响。基于年降水量(Bio12)与最湿月降水量(Bio13)两种因子的预测分布图都与适生区差异较大,且三种模型之间对比预测分布图差异也较大,说明这两种因子不是最主要的预测分布影响因子。将Bio1、Bio12、Bio13等作为共同影响因子,3种模型的预测分布图与其对应的适生区都十分接近,说明这3种因子是预测分布结果的确存在显著影响,尤其是Bioclim模型,二者基本完全重叠。此外,基于年均温(Bio1)的预测分布与3因子联合预测分布最为接近,进一步证明年均温(Bio1)对金钱松预测分布起“框架”作用,而其余环境因子对于预测分布起到“修饰”作用。

图3 三种环境因子模拟结果Fig.3 The results of the other three models图a1—a5分别为GARP模型基于年均温(Bio1)、年降水量(Bio12)、最湿月降水量(Bio13)、三因子(Bio1&12&13)以及对照组结果,图b1—b5分别为Domain模型基于年均温(Bio1)、年降水量(Bio12)、最湿月降水量(Bio13)、三因子(Bio1&12&13)以及对照组结果,图c1—c5分别为Bioclim模型基于年均温(Bio1)、年降水量(Bio12)、最湿月降水量(Bio13)、三因子(Bio1&12&13)以及对照组结果

2.4 金钱松不同时期分布格局变化模拟

金钱松从末次盛冰期到当前,分布范围主要向南迁移,分布扩张面积达到178510.08 km2,北部分布范围有部分收缩,达到122865.39 km2,净分布范围面积扩大了55644.69 km2。从当前到2070年分布范围呈现收缩趋势,主要是南部适生区大面积收缩,达到246396.07 km2,占现生分布范围的14.45%,而向北扩张68311.77 km2,仅为收缩面积的27.72%,分布范围收缩了178084.30 km2(图4)。

图4 基于Maxent和ArcGIS金钱松不同时期的分布格局变化Fig.4 Comparison of the distribution changes of P.amabilis based on ArcGIS and Maxent

3 讨论

3.1 金钱松适生范围预测

近年来物种分布模型发展迅速,多个可用的方法模型相继被开发[21]。本文通过模型检验得出四个模型都可以较为准确地预测出金钱松的潜在适生区(AUC>0.9),但GARP预测结果的AUC均值最小,精度较差,这可能是金钱松样本数据较少导致。杨会枫[19]用不同样本容量检验GARP预测结果精度,发现受到样本容量干扰,只有当样本容量超过一定阈值,AUC均值才会趋于稳定(大于0.9)。Domain和Bioclim预测结果相比其他模型适生区范围较小,这可能与模型受样点信息有关,金钱松分布数据主要来自于标本馆,样本数据存在一定偏好性,科研工作者一般根据自己研究所需进行采样或者单纯为标本馆收集标本,样本信息较为零散,缺乏系统性和代表性；其次金钱松样本数据较少可能导致生态位空间缩小的现象[26],并且由于某些地区金钱松种群遗传多样性较为丰富[5- 6],其物种在不同区域的生态位也可能发生漂移[26]。

Maxent模型预测结果比其他模型更为精确(AUC平均值为0.980),相关文献也证明Maxent在大、小样本数据条件下均能很好地预测物种适生区[14- 20]。但不同模型预测结果可以相互补充[21,27- 29],如在Maxent预测分布图上,河南郑州、信阳大别山以及云南昆明均属于金钱松的低适生区,而Domain模型中这一区域则表现为高适区,经查阅文献发现河南郑州、信阳大别山及云南昆明等地区均有人工成功引种金钱松的记录[28- 29],可确认是金钱松适生区。

3.2 金钱松分布范围变化

金钱松为高大落叶乔木,为深根系树种[5],其一年生幼苗生长旺盛期为4—5月和8—9月,这两个阶段的月均温和月降水量接近,为金钱松幼苗生长的最适条件[29]。本文发现年均温(Bio1)与年降水量(Bio12)是金钱松潜在适生区预测的主要影响因子,也决定了金钱松分布格局的形成以及迁移方向。末次冰期气候条件恶劣[23,30],而在金钱松分布区的浙江天目山、大明山、安徽黄山以及湖南、江西、湖北三省交界处的罗霄山脉是其良好的避难所[31],这些地区山峰地形多变,小生境水热条件优越,有效的降低了冰期气候对植物分布的影响。从末次间冰期(LGM)到当前气候条件下金钱松预测分布范围变化可知,因冰期之后气温回暖,金钱松的分布格局逐渐向温热的南方扩张。从当前到2070年金钱松的预测分布范围变化可知,未来其分布会一定程度向北扩张,但在南方却大面积收缩,总缩小面积达到178084.30km2。随着全球气候变化加剧[32],温室效应将大大影响植物分布格局以及其物候期[33-35],而植物分布变化是个较为缓慢的过程,全球气温上升速度远远超过金钱松自身迁移的能力,因此必须适时考虑人为辅助金钱松的北向迁移[33]。

3.3 金钱松潜在适生区与金钱松保护

金钱松作为濒危物种,自然分布范围有限,且生境受到各种威胁,建立金钱松自然保护区与种子园可以有效保护金钱松[5],自然保护区可以最大程度防止人为破坏,种子园可以有效保障种质资源和遗传水平,而二者的建立与物种适宜生境密切相关[36- 38]。本文结合多种生态位模型对金钱松的潜在适生区进行预测及其分布格局变化分析,对金钱松保护提出以下建议：除现有浙江西天目山金钱松自然保护区外,还可在江西铜鼓县、湖南张家界和衡阳这几处当前金钱松高适分布区域内建立自然保护区与种子园,保护与育苗、培种相结合,最大程度保护现有种群；在目前金钱松中等适宜区以及未来气候条件下高适分布区,如安徽北部、河南南部、湖北东南部等地区可以考虑人工引种,进行原土栽植协助金钱松北移实验,未来可将其作为园林树种,发展为园林绿植,而4种模型都预测为低适区或不适区的地方则不适于选择种植金钱松。