欧阳泽怡，李志辉，欧阳硕龙，程 勇，周志春，吴际友
（1.中南林业科技大学 林学院，湖南 长沙 410004；
2.湖南省植物园，湖南 长沙 410116；
3.湖南省林业科学院，湖南 长沙 410004；
4.中国林业科学研究院 亚热带林业研究所，浙江 富阳 311400）
赤皮青冈Cyclobalanopsis gilva，又名红椆，隶属于壳斗科Fagaceae青冈属Cyclobalanopsis常绿乔木[1]。其木材坚硬美丽，材质优良，能培育成大径材，可用于运动器材、高档家具和地板的制作，果实可食用，具有潜在的开发前景[2-3]。然而，赤皮青冈的野生资源较为稀缺，加之市场需求不断提升造成过度采伐，赤皮青冈资源日益减少，故被列为濒危树种[4]。在加强天然林保护的同时加大人工林扩繁力度，扩大赤皮青冈分布面积，符合《国家储备林建设规划》实施需要[5]。
目前，针对赤皮青冈的研究主要集中在种群生态学与群落特征、生长与生理特性、繁育技术以及分类与鉴定方面[6]，虽然对赤皮青冈的资源分布做了初步的调查与整理，但对其大尺度潜在自然分布区与生境因子的研究尚未见报道。
目前，研究物种分布状态与生境条件多采用野外调查的方式，该方法虽直观且准确，但历时一般较长，人力成本较高，也难以精确反映物种的所有分布区。近年来，基于物种分布信息和环境变量来预测物种适生区的分布已成为生物学研究中的热点之一，由此也开发出了众多物种分布模型，如BIOCLIM（Bioclimatic prediction system）、DOMAIN（Domain model）、GLMs（Generalized linear models）等[7-9]。其中，MaxEnt因其在较小样本量的情况下仍能进行无偏判断，得出较精准的预测结果，因而在众多模型中脱颖而出，成为物种适生区模拟预测的热点工具之一[10]。本研究在MaxEnt模型的基础上，利用R语言对其进行优化，同时与地理信息系统ArcGIS相结合，大尺度预测和模拟赤皮青冈当前的潜在分布区并评估其适生范围，有助于揭示赤皮青冈自然分布特征，探究影响赤皮青冈自然分布的主要环境因子，进行赤皮青冈天然林自然保护区规划，更好地维持生态多样性。对赤皮青冈在中国未来不同时期的适生区进行模拟，得到未来不同气候情境下的适生区分布及其动态变化预测，相关结果可应用于中长期造林工程规划，为该树种培育引种地选择提供参考，对未来进一步开展赤皮青冈人工林扩繁与研究具有重要意义。
1.1 赤皮青冈地理分布数据及处理
赤皮青冈的地理分布数据主要有以下3个来源：1）对赤皮青冈进行实地采样所记录的经、纬度信息；
2）国内标本数据库，包括中国国家标本资源共享平台（NSII, http://nsii.org.cn）、中国教学标本资源共享平台（http://mnh.scu.edu.cn）以及中国植物图像库（http://ppbc.iplant.cn/）；
3）全球多样性信息网络（GBIF, https://www.gbif.org/）。通过以上方式共获得赤皮青冈地理分布数据共193条。对所获得的分布数据进行筛选，剔除重复、无效（位置记录模糊且年代久远）的分布点，利用GPS定位系统对没有经、纬度记录的分布点进行经、纬度信息提取，最终获得赤皮青冈分布数据81条。
1.2 环境变量数据及处理
环境变量数据包括气候和土壤变量。气候变量为世界气候数据库（http://worldclim.org）提供的19个生物气候因子。当前时期气候变量的分辨率为30 s，约1 km；
未来气候数据为基于联合国政府间气候变化专门委员会在第6次气候变化评估报告中发表的在共享社会途径（SSPs）下的气候模型，包括2020—2040和2040—2060两个时间段，每个时间段分别使用4种升温模型：SSPs126、SSPs370、SSPs245和SSPs585，图层分辨率为2.5 min，约4 km。
为降低环境因子之间的共线性导致模型过度拟合的风险，增加模型预测的准确性，先对赤皮青冈进行预模拟试验。同时利用ArcGIS将赤皮青冈的分布点与环境因子进行数值采样，对采样得到的数值进行“Pearson”相关性分析。综合比较预模拟试验以及相关性分析的结果，去除预模拟试验中贡献率为0以及相关性系数|r|＞0.8的2个变量中贡献率较小的变量，最终选取7个变量用于赤皮青冈的适生区预测（图1和表1）。
表1 环境变量对赤皮青冈适生区预测的贡献率及置换重要值Table 1 The contribution and permutation importance of environmental variables to the predicted distribution of C.gilva
图1 环境变量之间的相关性分析Fig.1 Correlation analysis of environmental variables
1.3 模型建立与精度评价
将赤皮青冈分布数据以及筛选出的环境变量导入MaxEnt 3.4.1模型，随机选择75%赤皮青冈的分布数据作为模型训练集，用于适生区建模；
剩余25%的分布数据作为测试集，用于适生区模型精确度的检验。调用R语言程序包sdmtune对模型进行优化，选择“刀切法”分析各环境因子对赤皮青冈分布的影响，并绘制相应的环境变量响应曲线。模型迭代运行10次，取其计算的平均值作为模型预测的最终结果。采用受试者工作曲线以及AUC值对模型的精确度进行检验，AUC值越接近1，预测结果越准确。
1.4 适生区等级划分
利用ArcGIS 10.2软件中的“Conversion tools-ASCII to Raster”工具将预测结果中的“asc”格式文件转换成栅格文件，以“自然间断点分级法（Jenks）”对赤皮青冈当前及未来的适生区进行等级划分：非适生区（适生概率P＜0.3）、低适生区（0.3≤P＜0.5）、中适生区（0.5≤P＜0.7）以及高适生区（P≥0.7），最终获得赤皮青冈在不同时期、不同气候情境下的适生区分布图。以预测结果中最大检测灵敏度加特异性累积阈值（Maximum test sensitivity plus specificity cumulative threshold，0.360 4）为分类阈值，将各时期的适生区分布图分为2类，概率P≥0.360 4为适生区，以“1”表示；
P＜0.360 4为非适生区，以“0”表示，最终得到不同时期赤皮青冈分布的二进制图层，用于各时期各情境下赤皮青冈适生区面积以及质心变化的运算。
1.5 赤皮青冈适生区的动态分析
将赤皮青冈当前时期和未来时期各情境下的二进制图层分别导入ArcGIS中，利用SDMtoolbox工具包中的“Distribution changes between binary SDMs”工具对其进行适生区面积变化以及质心迁移（Centroid changes）的运算，从而获得赤皮青冈在未来时期不同气候情境下的适生区动态变化。
2.1 模型精度评价
基于R语言程序包sdmtune优化过后的MaxEnt模型对赤皮青冈的适生区进行模拟，采用ROC曲线和AUC值检验模型结果的精确度。其中，当0.5≤AUC＜0.6时，模拟失败；
0.6≤AUC＜0.7，模拟效果较差；
0.7≤AUC＜0.8，模拟效果一般；
0.8≤AUC＜0.9，模拟效果良好；
0.9≤AUC＜1，模拟效果极好。经过10次的重复运行，得到赤皮青冈当前适生区模型的平均AUC值为0.89，表明模型预测结果良好。
图2 MaxEnt预测赤皮青冈分布模型的ROC曲线Fig.2 ROC curves of the MaxEnt model for C.gilva
2.2 主导环境变量
通过MaxEnt模型中的“刀切法”对影响赤皮青冈适生区分布的环境因子进行贡献率和置换重要值分析。在参与适生区模拟的7个环境变量中，水生态因子对赤皮青冈分布的影响最大，总贡献率达到64%。其中，最干月降水量（Bio_14）和年平均降水量（Bio_12）的贡献值较高，分别达到了46.2%和14.6%，而降水量变异系数（Bio_15）的贡献率较低，只有3.2%；
除水生态因子外，温度因子对赤皮青冈生长的总贡献率为36%，其中温度年较差（Bio_07）和最冷季度平均温度（Bio_11）的共享率分别达到了11.9和11.8。置换重要值体现模型对环境变量依赖性的强弱，最干月降水量（26.1%）、温度年较差（22.6%）以及最冷季度平均温度（24.9%）的置换重要值都较高，可见赤皮青冈的适生区主要受最干月降水量、年平均降水量、温度年较差以及最冷季度平均温度影响。
2.3 当前时期的分布格局及适生区分布
我国赤皮青冈的分布区属于青冈属的中国—日本分布区[12]，其西界在滇黔桂地区。从图3A可以看出，赤皮青冈的实际分布格局东起浙江和台湾省，西至贵州西部，与罗艳等[12]对该物种的分布描述大致相符。此外，赤皮青冈在湖南省分布最广，其次是浙江贵州等省。从图3B来看，赤皮青冈现有的分布点均在模型所预测的范围内，由此表明模型的预测结果与其实际分布基本吻合。
基于MaxEnt预测的结果，对赤皮青冈的适生区进行适宜等级划分。按照分布概率P确定赤皮青冈适生区等级划分：P＜0.3为非适生区，0.3≤P＜0.5为低适生区，0.5≤P＜0.7为中适生区，P≥0.7为高适生区。经过MaxEnt模型预测分析，得到赤皮青冈在我国境内的潜在适生区分布结果（图3）。
图3 赤皮青冈当前分布格局及其适生区预测Fig 3.Distribution status and potentially suitable areas for C.gilva in the current and future periods
结果显示，赤皮青冈在我国的适生区（P≥0.3）面积高达674 566 km2，其中赤皮青冈在我国的高适生区总面积为31 206 km2，约占总适生区面积的4.6%。高适生区分布面积最多的是湖南省（14 838 km2），高适生区占比为47.5%，赤皮青冈在湖南省已被列入优先研究的珍贵乡土用材树种名单[13]；
其次为台湾省（7 169 km2），高适生区占比为23%；
再次为浙江省（4 373 km2）和贵州省（3 452 km2），高适生区占比分别为14%和11%；
而福建、江苏、广东以及重庆的高适生区面积均较少，其面积共占高适生区总面积的4%。中等适生区以及低等适生区沿高等适生区周围分布，主要集中在浙江省、福建省、江西省、广东省、广西壮族自治区、湖南省、湖北省、贵州省、重庆市、四川省以及台湾省。其中，中适生区的面积为232 510 km2，低适生区的面积为410 850 km2，分别占总面积的34.4%和60.9%。低适生区分布最为广泛，较高、中适生区分布地有明显向东、向南扩散的趋势。赤皮青冈在各省区的不同适生区面积预测结果见表2。
表2 当前时期赤皮青冈在中国各省区的适生区面积Table 2 Predicted areas for C.gilva in China under current climate conditions
2.4 未来不同气候情境下的适生区分布及其动态变化
基于未来不同时期、不同气候情境对赤皮青冈的潜在适生区进行模拟。结果显示，2021—2040年，赤皮青冈在4种不同的气候情境下的适生区面积相较于当前时期均有所增加，而随着升温途径的升高，高适生区面积呈逐步扩增趋势。在2041—2060年时间段内，相较于前20年，赤皮青冈的潜在分布区面积有所减少，其高适生区的面积也呈缩减趋势。其中，在SSPs126和SSPs245的气候情境下，赤皮青冈在我国西南至东南地区仍能维持较大的高适生区面积；
而在SSPs370和SSPs585情境下高适生区面积大大减少；
在SSPS585情景下，赤皮青冈高等适生区面积收缩明显（图4）。可见湖南省范围内高适生区面积减少一半以上，江西省出现较为明显的大块适生区空白，但与赤皮青冈当前的潜在分布区相比较，高适生区面积整体仍然呈增加趋势。
图4 赤皮青冈在不同时期以及不同情境下的适生区预测Fig.4 The predicted areas for C.gilva in different stages under different climate conditions
从赤皮青冈核心分布区的质心迁移结果来看，该树种在不同时期以及不同气候情境下的适生区范围总体呈向高纬度地区迁移的趋势（图5）。当前气候条件下，赤皮青冈的分布质心位于江西省吉安市永新县龙门镇六团村竹山（114°8′13″E，27°8′20″N）。
2021—2040年（图5A），赤皮青冈分布质心在SSP126情景下，将向北移动至湖南省浏阳市七星岭（114°7′55″E，28°0′26″N）；
在SSP245情景下，将向高纬度方向移动至江西省万载县白水乡（114°12′52″E，28°16′23″N）；
在SSP370情景下，赤皮青冈分布质心将向低纬度方向移动至（114°8′46″E，27°41′48″N）；
在SSP585情景下，将向东北移动至（114°25′57″E，28°17′37″N），该坐标和在SSP245情景下一样位于江西省万载县。
续图4Continuation of Fig.4
图5 赤皮青冈不同时期以及不同情境下的质心迁移Fig.5 The centroid changes of C.gilva in different stages under different climate conditions
2041—2060年（图5B）,赤皮青冈分布质心在SSP126情景下，经度无明显变化，将向北移动至江西省宜春市袁州区（27°48′18″N，114°32′36″E）；
在SSP245情景下，将向高纬度方向移动至湖南省平江县余段乡（114°21′54″E，28°49′41″N）；
在SSP370情景下，赤皮青冈分布质心和在2021—2040年模拟的结果一样，将向低纬度方向移动，来到湖南省浏阳市张坊镇（113°10′10″E，28°22′12″N）；
在SSP585情 景下，将向东北移动至江西省九江市修水县山口镇（114°28′27″E，28°48′52″N）。
在2021—2040年的4种气候变化背景下，赤皮青冈分布质心的移动幅度比2041—2060年的4种气候变化背景下更小，在文中所示不同时期以及不同情境下的赤皮青冈质心迁移总的来说幅度较小，经向移动幅度远小于纬向移动幅度。
植物的生长与分布受多个环境因素的综合影响。从刀切法分析的环境因子结果来看，赤皮青冈的适生区分布主要受到水热条件的限制。在参与建模的环境因子中，最干月降水量、年平均降水量、温度年较差以及最冷季度平均温度对赤皮青冈的影响较大。吴丽君等[14]的研究表明，水分胁迫在赤皮青冈幼苗生长过程中起着重要作用，在短期干旱下幼苗成熟、叶片变薄、栅栏组织和海绵组织厚度降低，气孔变小以适应群落的季节性干旱，而在长期干旱条件下则不利于幼苗的生长[15-16]。因此，在赤皮青冈的栽培引种过程中，需要加强对该树种在旱季的水分管理，以提高幼苗的存活率。
本研究利用MaxEnt模型预测赤皮青冈当前的潜在适生区，并对其适生区进行了适宜性划分。结果显示，赤皮青冈在我国的适生区（P≥0.3）面积高达674 566 km2，其分布区包括了台湾、浙江、福建、江西、湖南和贵州等多个省份，赤皮青冈高适生区主要分布于湖南、台湾、川西以及川黔交界处等赤皮青冈实际分布较多的区域。理论上赤皮青冈在我国应具有较大面积的适生区，但受种子雨密度、幼龄个体间的竞争强度影响，赤皮青冈在幼苗时期的死亡率较高，使得种群增长速度受阻[17]。此外，种间竞争抑制了赤皮青冈的种群增长，其在群落中的优势度并不明显，因此赤皮青冈在我国中部地区虽然有较高的适生度，但实际分布却相对较少，仅湖南通道、湖南靖州、湖南绥宁、福建建瓯和江西婺源等地有少量天然群落分布[17-20]。虽然受上述因素影响，赤皮青冈在群落中的优势度不高，但赤皮青冈环境适应能力较强，仍能发展至顶级群落[21]。
面对气候的变化，植物的适生区分布会发生不同程度的响应。已经被证实的是，在当前全球变暖的情况下，植物将会逐渐向高海拔或高纬度地区迁移[22-23]。面对气温的逐渐升高，北方地区一些难以耐受高温胁迫的物种将面临适生区缩小和遗传多样性减少的风险，而一些喜热物种，随着全球变暖其适生区则会进一步扩展[24-26]。模型结果显示在2021—2040年和2041—2060年这2个时间段内，赤皮青冈的适生区面积均会呈现不同程度的扩增。相较于当前时期，赤皮青冈的高适生区面积有所增加，其核心分布区也呈现向高纬度地区迁移的趋势。基于其主导环境因子可知，赤皮青冈还受温度年较差以及最冷季度平均温度影响，随着气温增加，这2类温度因子对其分布的影响也将逐渐减弱，由此推测适度的增温有助于赤皮青冈的生长。
基于本研究结果与讨论可知：1）赤皮青冈主要受水热条件影响，尤其对水分要求较高；
在栽培引种过程中应当加强对赤皮青冈幼苗的水分管理，避免长时期的过度干旱，以提高幼苗的存活率；
2）赤皮青冈当前在我国的适生区面积较高，在未来时期还将进一步增加，可见珍贵树种赤皮青冈的人工林营造仍有较大发展空间。赤皮青冈的繁育栽培技术已较为成熟，可在不影响物种多样性的前提下尽量在模型预测的中高适生区内开展人工林扩繁工作，提升造林效率的同时获得更好的生态及经济效益；
3）模型结果显示，温度对赤皮青冈的生长与分布至关重要，然而尚未见针对赤皮青冈温度胁迫的研究报道，因此还需加强温度对赤皮青冈生长影响方面的研究。