基于Maxent模型的西南地区香果树地理分布预测研究

魏俊，朱坤，2，陈文德，2*，卢雨田，蔡艳琨，吴群，黄仲宣

（1.成都理工大学旅游与城乡规划学院，四川成都610059；2.四川省社会科学重点研究基地（扩展）国家公园研究中心，四川成都610059）

香果树（Emmenopterys henryiOliv.）隶属茜草科（Rubiaceae）、香果树属（EmmenopterysOliv.）落叶大乔木，是古老孑遗植物。在《中华本草》和《新华本草纲要》中记载香果树具有湿中和胃，降逆止呕之功效。香果树种子萌发力较低，主要靠根系或树桩萌条更新，繁殖比较困难，尚未发现纯林，只见其星散分布，野外资源量少，1991 年被认定为近危保护植物纳入《中国植物红皮书》，1999 年被列为第一批中国特有二级珍稀濒危保护植物，有研究认为香果树是香果树属中己知的唯一现存种［1］，具有较高的经济和科研价值［2］。

目前，关于香果树主要研究具体保护区等小区域的种群特征、群落结构、生长动态、基因组特征、繁殖能力和病害等方面且以中东部省区为主［3-8］，对西南地区香果树的地理空间范围分布情况的研究较少，其中张永华对全国香果树的亲缘地理学和景观遗传学及其生态位分化进行了研究［1］；周文娟等对贵州省境内香果树的分布现状研究发现其在贵州属广布种，但多为零星分布［9］。

最大熵模型（Maximum entropy model，Maxent）是近年来用作预测物种在地理尺度上空间分布的重要方法之一。王娟娟等利用Maxent 模型对川贝母的潜在分布和适应性进行了研究［10］；崔晋亮等利用该模型对蓝莓的潜在分布区进行了预测研究［11］；张东方等采用最大熵模型研究了当归在全球的生态适宜区和生态特征［12］；王国峥等采用4 种基于不同算法的模型（GARP、Bioclim、Domain 和Maxent）研究金钱松的潜在分布区，结果表明Maxent 模型对植物潜在地理分布预测结果比其他模型更为精确［13］。因此，本研究采用最大熵模型，以对西南地区凉山山系样地调研为基础，结合香果树标本数据，对西南地区（川、云、贵、藏、渝）香果树在历史时期、当前时期和未来时期适生区的空间分布情况进行预测分析研究。

1 数据来源与研究方法

1.1 样本数据来源

采用样地调查法对凉山山系香果树的分布情况进行实地调研，获得样本数据点44 个，在全球生物多样性信息数据库网络（GBIF，http：//www.gbif.org/）、中国数字植物标本馆（CVH，http：//www.cvh.ac.cn/）和中国自然保护区资源平台（http：//www.papc.cn/）中共查询到西南地区香果树样本点117 个，其中贵州 19 个、云南 5 个、重庆 7 个、四川 86 个，数字植物标本馆中显示西藏地区有1 份香果树标本，但是没有提供坐标信息。加上实地调研的数据一共得到样本点161 个，统计发现有5 个重复样本点，去除后共有156个样本点（图1）。

图1 西南地区香果树分布情况Fig.1 Distribution of E.henryiin Southwest China

1.2 生态因子数据来源

生态因子共36 个，其中环境因子包含19 个生物气候因子、7 类气候因子和1 个地形因子，环境因子 数 据 来 源 于 Worldclim（http：//www. worldclim.org/），此数据集为当前条件（2000～2020 年）且分辨率为30S（约1km2）的生物气候变量，这些变量主要反映温度和降水的特点及其季节性变化特征，地形因子来源于地理空间数据云（http：//www.gscloud.cn/），下载分辨率为90 m 的数字高程地图（DEM），从中提取海拔，然后在ArcGIS10.2 中重采样生成分辨率为30″的栅格数据图层［14］。土壤数据来自联合国粮食及农业组织网站中的世界土壤数据库（Har⁃monized world soil database v1.2|FAO SOILS POR⁃TAL），坐标系为WGS-1984，栅格大小约为1km2，共选取了9个土壤因子。

1.3 研究方法

1.3.1 数据预处理

表1 生态因子缩写及含义Table 1 Abbreviations and meanings of ecological factors

从国家基础地理信息系统（http：nfgis.nsdi.gov.cn/）下载获得的1∶400 000 0 的西南地区（四川、重庆、云南、贵州、西藏）标准地图作为分析底图，再利用ArcGIS 10.2 软件将西南地区的生物变量数据掩膜出来，将图层的坐标系定义为GCS－WGS－1984，并导出为ASCII格式。将香果树样本点数据文件保存为“.csv”格式。

1.3.2 生态因子筛选

由于环境变量之间具有多重共线性，会导致预测分布过度拟合［13］，用 ArcGIS 10.2 中 ArcTool⁃box 下多元分析的波段集统计功能计算生态因子的协方差和相关矩阵进行相关性检验，筛选出相关系数绝对值小于0.75 的生态因子，并依据香果树的生长环境再次筛选，最终得到4 个土壤因子：表层土砾石含量、表层土粘土含量、表层土壤有机碳含量、表层土基本饱和度，13 个环境因子：海拔、1 月和 7 月平均降水量、1 月太阳辐射、7 月最高温、1 月最低温、年均温、等温性、温度季节性变化的标准差、年均降水量、最干月降水量、降水量变异系数、最冷季度降水量，共筛选出17 个生态因子。

图2 生物气候因子相关系数图Fig.2 Correlation coefficient of bioclimatic factors

1.3.3 最大熵模型

最大熵原理（The principle of maximum entropy）起源于信息论和统计力学，是基于有限的已知信息对未知分布进行无偏推断的一种数学方法［15］。Phil⁃lips 团队基于生态位理论，考虑气候、海拔、植被等环境因子，用最大熵原理作为统计推断工具，构建了物种地理尺度上空间分布的最大熵模型［16-17］，编写了一个可以免费获取的计算软件（http：//www.cs.princeton.edu/~schapire/maxent/）［15］。该模型是借助机器学习方法建立的对最大熵原理的应用［16］，Elith等［18］在2011年对其做了统计上学的解释，使这个模型在生态学领域得到了更加广泛的应用。

最大熵模型软件的具体操作为：将已知的物种分布数据的“.csv”格式文件和筛选出来的生态因子导入Maxent 软件，定义结果输出的位置以及生态因子图层位置，勾选创建响应曲线和刀切法复选框对变量重要性进行测试，设置采取25%的样本点作为测试集，剩下75%的样本点作为训练集，其余选项采用模型的默认设定，重复运行10次进行建模。受试者工作特征曲线（Receiver operating characteristic⁃curve，ROC）分析法在物种潜在分布预测模型评价中应用比较普遍，这一诊断试验评价指标的认可度较高［12，19］。计算 ROC 曲线下方的面积（Area under curve，AUC），其中AUC 的数值范围为 0.6～0.7 时预测精度较低，0.7～0.8 时预测精度中等，0.8～0.9 时预测精度良好，大于0.9 时预测精度优秀；并利用Jackknife刀切图来评价各环境因子的权重［14］。

2 结果分析

2.1 预测精度分析

将筛选出来的17 个生态因子以及香果树样本点数据分别导入Maxent3.3.3软件中，重复运行10次上述方法建立最大熵模型，得到生态因子Maxent 模型平均训练集的AUC 值为0.956（图3），表明模型预测结果为优秀，满足进行香果树地理分布和生态因子关系模拟研究的精度要求。

图3 生态因子预测ROC曲线Fig.3 Receiver operating characteristic curve of ecological factor prediction

2.2 潜在分布区及其适生区等级划分

为了对各区域的适生程度加以区分，按照政府间气候变化专门委员会（Intergovernmental panel on climate change，IPCC）关于可能性的划分标准并结合文献研究［10］，将香果树潜在分布预测区划分为4个等级：P ＜ 0.1 为不适宜区；0.1 ≤ P ＜ 0.4 低适宜区；0.4 ≤ P ＜ 0.6 为较适宜区；P ≥ 0.6 为高适宜区。将最大熵模型模拟导出的“Oliv._avg.asc”文件导入ArcGIS 10.2，利用工具栏将“.asc”文件转为栅格文件并进行重分类，得到香果树适生区分布图（图4）。

图4 生态因子模拟香果树潜在适生区Fig.4 Ecological factors simulating the potential adaptive area of E.henryi Oliv

利用ArcGIS10.2 统计分析工具，将模拟结果的栅格数据转换为矢量数据并进行面积计算，得到香果树当前潜在分布区及其适生区面积。生态因子模拟结果中香果树分布的高适宜区面积约为3.27×104km2，中适宜区面积约为 16.78 × 104km2，低适宜区面积约为24.13 × 104km2，不适宜区面积约为198.82 × 104km2，其中适宜区总面积仅占西南地区总面积的18.19%。

2.3 主导生态因子特征分析

根据参与模型建立的生态因子对最大熵模型的贡献率（表2），可以判断影响物种分布的主要生态因子。生态因子贡献率结果显示共有8个环境因子对模拟结果的贡献率大于4%，由高到低依次为表层粘土含量（25.25%）、最干月降水量（14.15%）、海拔（13.51%）、表层土砾石含量（13.1%）、1 月太阳辐射强度（11.69%）、表层土基本饱和度（5%）、1 月最低温（4.38%）、表层土壤有机碳含量（4.1%），累积贡献率高达91.12%。

表2 生态因子贡献率Table 2 The percent contribution of ecological factors

2.3.1 香果树分布区环境特征

基于各气候因子响应曲线（图5）可知，当潜在分布概率大于0.5 时，其对应的生态因子值比较适合香果树的生存［12］，曲线中蓝色部分表示10次重复运算的波动范围，红线表示10次运算的平均值。最干月降水量小于8 mm 时香果树适生概率小于0.5；海拔在950~2 200 m 时其适生概率大于0.5；1 月太阳辐射在6 900~8 100 kJ/m2·day 时其适生概率大于0.5；1 月最低温在-4~1℃时其适生概率大于0.5，香果树在极端最低温-15℃时都可能生存下来；7 月最高温在21~32℃时其适生概率大于0.5，在35℃的高温下仍具有9%的生存概率；降水量变异系数的范围为60~95 时，适生概率大于0.5；等温性的适宜范围在25~34。

图5 香果树的潜在分布概率与环境主导影响因子的关系Fig.5 Relationship between the potential distribution probability of E.henryi Oliva environmental dominant impact factors

2.3.2 香果树分布区土壤特征

基于各土壤因子响应曲线（图6）可知适合香果树生存的土壤因子范围。表层土粘土含量在21%~40%时香果树适生概率大于0.5；表层土砾石含量在4%的时候其适生概率达到了0.5，随着表层土砾石含量升高，适生概率短暂下降，砾石含量达到5.5%之后其适生概率开始增加，当砾石含量增加到8%时适生概率最大，砾石含量适宜范围为4%~10%；表层土基本饱和度在15~22之间时香果树适生概率逐渐下降，之后香果树存在概率逐渐增加，当表层土基本饱和度达到33时香果树存在概率最大，表层土基本饱和度适宜范围为10~40；香果树适生概率随着表层土壤有机碳含量的增加逐渐增大，表层土壤有机碳含量达到1.5%的时候香果树存在概率最大，表层土壤有机碳含量适宜范围为0.5%~1.5%。

图6 香果树的潜在分布概率与土壤主导影响因子的关系Fig.6 Relationship between the potential distribution probability of E.henryi Oliva soil dominant impact factors

2.4 未来时空分布预测

由于香果树属于落叶大乔木，其群落演替时间较长，因此基于历史气候数据和未来气候数据对其近100 年的时空分布变化情况进行预测。气候数据来源于 Worldclim（http：//www.worldclim.org/），其中历史气候数据是该网站发布的1 970~2 000 年世界气候2.1 版气候数据；未来气候数据（2030s、2050s、2090s）来自于该网站 BCC-CSM2-MR 模型ssp585 共享社会经济路径，目前提供的未来气候数据的空间分辨率最高为2.5 min，提供的气候数据类型为月最低温、月最高温、月平均降水量和19 个生物气候因子。筛选出1 月和7 月平均降水量、1月和7 月最低温、1 月和7 月最高温、年均温、等温性、温度季节性变化的标准差、年均降水量、最干月降水量、降水量变异系数、最冷季度降水量，假设海拔因子没有变化，沿用前面的DEM 数据，共得到14项环境因子。

将4 个时期的数据按照前面的方式进行处理，并分别用Maxent3.3.3 软件进行模拟预测，平均训练集的 AUC 值分别为 0.942、0.939、0.933、0.930，预测结果达到了优秀水平。利用ArcGIS 10.2 进行加工处理，得到 20 世纪末、2030s、2050s和2090s 各时期香果树适生区分布情况（图7），同样利用统计分析工具计算4 个时期香果树各等级适生区面积（图8）、Maxent3.3.3 软件输出结果中包含环境变量对最大熵模型的相对贡献的估计值（图9）。

IPCC 第五次评估报告利用全球气候耦合模式，模拟发布了在温室气体排放浓度最低情景下，到本世纪中叶全球平均地表温度将上升0.4~1.6℃［20］。根据图 8 可以看出，在未来的 3 个时期，随着气候的变化，香果树在西南地区的不适宜分布区面积将会减少，低适宜区和中适宜区面积将会有所增加，高适宜区面积在2030s 和2050s 将会有所增加，但是在2090s 时期内相比于前一时期将会有所减少。由图7 可以看出，未来的三个时期内低适宜区分布范围在西藏的东南边缘和云南的东北部将会有所增加，中适宜区分布范围将主要在贵州省境内有所增加，2030s 和2050s 高适宜区分布范围将在四川的东北边缘、重庆的东南部和贵州境内有所增加，2090s 高适宜区分布范围在贵州中西部地区将会有所减少。

图7 不同时期香果树适生区的预测分布图Fig.7 Distribution map for predicting the suitable areas of E.henryi Oliv in different periods

图8 不同时期预测适生区分布面积Fig.8 Distribution acreage of adaptive area predicted in different periods

根据图9可以看出最干月降水量因子在各时期贡献率最大，尤其是2030s、2050s 这两个时期内其贡献率超过了35%；其次是海拔和等温性因子，在四个时期内贡献率都超过了10%；1 月平均降水量、温度季节性变化的标准差、7 月最高温、7 月平均降水量这几个因子的贡献率最小，不超过1%。由此可知，影响香果树适生区在西南地区分布范围的主要环境因子是最干月降水量、海拔和等温性，表明香果树具有喜湿润的特点。

图9 环境因子在各时期贡献率百分比Fig.9 The percent contribution of environmental factors in each period

3 结论

将 ArcGIS 10.2 与 Maxent3.3.3 结合，通过 17 个生态因子对西南地区香果树当前的潜在分布区进行预测，预测精度达到优秀水平（AUC=0.956），利用ArcGIS 10.2 统计分析工具计算得到西南地区香果树低适宜区面积约为21.54×104km2，占西南地区总面积的9.93%；中适宜区面积约为16.78 × 104km2，占西南地区总面积的6.91%；高适宜区面积约为3.28 × 104km2，占西南地区总面积的1.35%。香果树分布中等适宜以上区域面积占比为8.26%。分析结果表明西南地区香果树的适生区范围很小，中高适宜区主要分布在四川南部、云南东北部、贵州省和重庆市全域范围内，分布比较零散，西藏的东南边缘的唐古拉山区分布有少量中低适宜区。

在影响香果树分布的主导环境因子中，最干月降水量、海拔、1 月太阳辐射这三项的贡献率都超过了10%。适宜区最干月降水量的阈值为5~43 mm，海拔的阈值为133~3 000 m，1 月太阳辐射的阈值为5 030~9 700 kJ/m2·day，1 月最低温的阈值为-14℃~7℃。在影响香果树分布的主导土壤因子中，表层粘土含量、表层土砾石含量的贡献率超过了10%，其中表层土粘土含量的贡献率达到了25.25%；适宜区表层粘土含量的阈值为15%~56%，表层土砾石含量的阈值为0~15%，表层土基本饱和度的阈值为10%~100%，高适宜区表层土壤有机碳含量的阈值为0.2%~2.5%。

对香果树在西南地区近100 年的适生区的预测结果表明，在未来的几十年内香果树的适生区面积整体会有小幅度增加，但主要是中低适宜区取代小部分不适宜区，高适宜区面积将仅在2050s 有所增加，到2090s 则会减少至现有水平。其中最干月降水量、海拔和等温性等环境因子是主要影响因子。

通过对西南地区香果树的潜在分布及适生区分布的预测研究，为西南地区的中国二级珍稀濒危保护物种香果树的树种分布调查和保护提供了较为详细的依据；对西南地区香果树分布主导生态影响因子的分析，得出了香果树生长的生境条件；近100 年的预测结果展示了只在气候环境变化的影响下，西南地区香果树适生区范围的变化趋势，为香果树的植育保护提供了参考。