全球变化格局下重要传粉昆虫大蜜蜂的潜在适生区变化

徐春阳,刘秀嶶,贺春玲,高 洁,彭艳琼

(1.中国科学院西双版纳热带植物园,云南勐腊 666303;2.中国科学院大学,北京 100049;3.云南省农业科学院农产品加工研究所，昆明 650221;4.河南科技大学林学院,河南洛阳 471023)

已知的传粉昆虫占了全球传粉媒介的80%～85%，能为全球约90%的野生植物和75%的粮食作物传粉(Kearnsetal.,1998;Smithetal.,2015;Ollerton,2017)。其中，膜翅目蜜蜂总科(Apoidea)是最重要的传粉者(Ollertonetal.,2011)，最新的全球蜜蜂分布研究显示在中纬度地区蜜蜂的丰富度最高，干旱、温带地区的蜜蜂物种比热带地区的蜜蜂物种更丰富，这与相关研究主要集中在温带和地中海热点地区有关(Michener,2007;Orretal.,2021)。在热带地区，研究多集中在巴西、澳大利亚北部等地区(Gianninietal.,2012)，而在其他热带地区的研究很少。当前，全球气候变化、土地利用改变、人类活动加剧等正极大地威胁着野生传粉昆虫的多样性及分布(Pottsetal.,2016)；气候变化影响一些野生传粉昆虫的季节性活动、丰度和活动范围(Bartomeusetal.,2011;Martinsetal.,2015)；在全球变化背景下植物-传粉者互作网络的结构和功能被证实在过去的120年间发生了退化，到2010年已有50%的蜜蜂物种灭绝，导致传粉服务的数量和质量均下降(Burkleetal.,2013)。

蜜蜂属Apis物种全球有11种，除西方蜜蜂外均为亚洲起源(Loetal.,2010)。该属中的大蜜蜂A.dorsata是南亚和东南亚特有的一种野生蜜蜂，仅分布于24°N以南的湿润热带地区(Michener,2007)。大蜜蜂是热带地区重要的传粉昆虫，筑巢高度可达32.6 m，能为林冠和林下的诸多植物传粉；该蜂具有较强的迁徙能力和广泛的觅食种类，其传粉、活动的生境非常多样(Momoseetal.,1998;Nagiretal.,2016;杨培等,2018)。森林是大蜜蜂的主要栖息地和活动场所，在印度南部的卡纳塔克邦(Karnataka)，有74科共252种植物与大蜜蜂的觅食有关，最多的为乔木(49.3%)，其次是草本(23.5%)、灌木(21.7%)和攀缘植物(5.5%)(Raghunandan and Basavarajappa,2014)。在半自然栖息地和农田，大蜜蜂能为油菜、洋葱、芒果等重要经济作物传粉，已被证实是一些重要经济作物最丰富和最有效的传粉者(Alietal.,2011;Robinson,2012;Sajjadetal.,2017)。大蜜蜂还具有更为灵敏的视觉，可为某些仅在夜晚开花的植物传粉(Somanathanetal.,2009)。尽管大蜜蜂是泛化的传粉者，但有一些植物专性地依赖大蜜蜂传粉，例如金虎尾科风筝果Hipyagebenghalensis仅依赖大蜜蜂传粉进行繁殖(Renetal.,2013)。在南亚和东南亚，大蜜蜂还能为当地土著人提供多种服务，包括蜂产品和文化服务(Matiasetal.,2019)。然而，该地区正面临人类干扰增强、森林覆盖率减少(Achardetal.,2002;Hansenetal.,2013)，大蜜蜂的栖息环境遭受到不可逆转的影响。因此，本研究以大蜜蜂为对象，探讨其如何响应未来全球变化，有利于保护这一重要的传粉昆虫资源。

蜜蜂属被认为最早出现在更新世，东南亚的历次冰期和伴随的气候动荡可能严重影响蜜蜂属的物种形成和适应性辐射，也影响具体物种分布区的扩张退缩以及种群大小的变化(Hepburn and Radloff,2011)。上新世至更新世的气候变化引起了海平面的变化，冰期海平面下降导致巽他古陆(Sundaland)的陆架暴露，更新世时期巽他古陆反复与亚洲大陆相连，使得生物能从大陆迁移到东南亚的群岛上(Sodhietal.,2004)。冰期低海平面时期，气候更为干燥，暴露在海平面以上的陆地面积增加而适合生物生存的地理区域减少，此期东南亚地区成为了大多数生物类群的避难所。根据地理隔离和形态差异，大蜜蜂被分为3个亚种：广布种A.dorsatadorsata，分布于苏拉威西岛的A.dorsatabinghami及分布于菲律宾的A.dorsatabreviligula(Maa,1953;Loetal.,2010)。大蜜蜂3个亚种的分布显示大蜜蜂的扩散可能受到了第四纪冰期的影响，东南亚的南部岛屿可能是其避难所。因此，可以从第四纪更新世背景下考虑大蜜蜂的起源、扩散，并分析大蜜蜂如何在冰期、间冰期应对连续的气候和地质波动。

生态位模型在保护生物学中已被广泛应用，基于最大熵原理的MaxEnt模型，因为便于使用和较高的预测精度，已经成为物种分布和环境生态位模拟最流行的工具之一(Phillipsetal.,2017;祝梓杰等,2017)。模型通过输入物种的分布位置，也被称为仅出现数据(presence-only data)，以及一组环境预测因子，基于最大熵原理估计物种的分布并给出此时的环境参数(Elithetal.,2011;Merowetal.,2013)。模型模拟的复杂性会对结果产生重要影响，要求在使用过程中慎重选择环境变量，根据模拟物种优化模型参数，才能较为合理地反映物种对环境因子的响应并准确地模拟出该物种的潜在分布区(朱耿平和乔慧捷,2016)。本研究根据已报道的大蜜蜂物种分布数据，补充了馆藏标本信息及我们在野外调查的分布数据，以5 km的空间距离空间稀释后，得到了1960-2020年间的313个大蜜蜂物种分布数据；使用13个环境变量，通过MaxEnt模型模拟了大蜜蜂当前的潜在适生区；再使用9个气候变量并结合公共地球系统模型(CCSM4)模拟了大蜜蜂过去、当前和未来的潜在适生区。主要回答以下问题：(1)基于当前13个环境变量条件下大蜜蜂的潜在适生区；(2)筛选影响大蜜蜂分布的主要环境因子；(3)模拟过去、当前和未来气候变化下大蜜蜂潜在适生区的空间格局变化。为确定大蜜蜂的重点保护地区，有效保护大蜜蜂提供科学依据。

1 材料与方法

1.1 大蜜蜂物种分布数据

从全球生物多样性信息网络数据库(Global Biodiversity Information Facility)(https:∥www.gbif.org/)和公开发表的文献中获取全球大蜜蜂分布点数据，截至2020年5月13日。中国地区的大蜜蜂分布数据来源于中国科学院动物研究所国家博物馆、昆明动物博物馆，以及我们野外调查的分布信息(图1)。

馆藏标本记录有些缺乏精确的纬度和经度信息，因此筛选中位置精确到县及县以下的记录，并只保留一个同一级别的重复记录，最后用Google Earth定位该记录的经纬度。为了与气候数据(1960-1990)相匹配，在所有的分布点中，我们保留了1960-2020年间的物种分布数据。对存在空间上聚集的分布点，在Arcgis10.7中以5 km的分辨率对分布点进行空间过滤，保证在2.5 arc-minute分辨率的环境图层中每个网格单元只有一个分布位点，以防止模型过度拟合，提高模型的性能(Kramer-Schadtetal.,2013;Boriaetal.,2014)。再将分布位点描绘至地图上，删去落于亚洲大陆以外的及有明显错误的分布点。最后，保留了1960-2020年间的313个大蜜蜂物种分布数据用于模型模拟(图2:A)。

1.2 环境变量

物种分布主要受到非生物和生物因素的影响，大蜜蜂广泛分布于南亚和东南亚的湿润热带地区(Michener,2007)，在大尺度空间上气候、生物地理及人类活动因子可能影响大蜜蜂的潜在适生区(Pearson and Dawson,2003;Hortaletal.,2010;朱耿平等,2013)。其次，森林是大蜜蜂的主要栖息地和活动场所，森林覆盖和林冠高度作为大蜜蜂的栖息地因子(Raghunandan and Basavarajappa,2014;Nagiretal.,2016)。综合考虑这些因素后，从世界气候数据库WorldClim version 1.4(https:∥www.worldclim.org/data/v1.4/worldclim14.html)中，获取空间分辨率为2.5 arc-minute的19个标准生物气候变量的全球气候栅格数据(1960-1990)和全球海拔数据(Hijmansetal.,2005)；获取2000年全球森林覆盖图层(http:∥earthenginepartners.appspot.com/science-2013-global-forest/download_v1.7.html)(Hansenetal.,2013)；获取2005年全球林冠层高度图层(https:∥webmap.ornl.gov/wcsdown/dataset.jsp?ds_id=10023)(Simardetal.,2011)；从社会经济数据与应用中心(Socioeconomic Data and Application Center)(https:∥sedac.ciesin.columbia.edu/)数据库中的全球人类影响指数数据集(The Global Human Influence Index Dataset)中，获取空间分布率为30 arc-second(约1 km)的人类影响指数(human influence index)数据(1995-2004)。人类影响数据收集自野生动物保护协会(Wildlife Conservation Society,WCS)和哥伦比亚大学国际地球科学信息网络中心(Columbia University Center for International Earth Science Information Network,CIESIN)，通过9个包括人口压力、人为土地利用、人造基础设施、人造通道的环境图层叠加而成。最后，我们将模拟结果匹配至世界野生动物基金会(World Wildlife Fund,WWF)所划分的867个陆地生态区(Olsonetal.,2001;Dinersteinetal.,2017)，以探究不同生态区地理隔离对大蜜蜂潜在适生区的影响。

对19个标准生物气候变量进行Pearson相关性(Pearson correlation)分析，以避免环境变量之间相关性较高而导致共线性，提高模型的预测精度(Merowetal.,2013)，剔除高度相关的变量(即相关系数>0.8的变量)，最后保留9个气候变量。最终，使用9个气候变量、1个地形变量、1个人类影响变量以及2个栖息地变量，共13个环境变量来模拟大蜜蜂在当前的潜在适生区；模型中使用的环境变量及其代码见表1。

表1 本研究中使用的环境变量Table 1 Environmental variables used in this study

为了研究全球气候变化对大蜜蜂物种分布格局的影响，本研究使用9个气候变量和WorldClim中的CCSM4对大蜜蜂在过去、当前和未来3个不同时期气候下的潜在适生区进行模拟。在世界气候数据库(WorldClim version 1.4)中下载空间分辨率为2.5 arc-minute的CCSM4，该公共地球系统模型是IPCC第5次评估报告中模拟全球气候变化的“耦合模式比较计划”(Coupled Model Intercomparison Project,CMIP5)中基于1961-2005年的观测资料模拟的，能覆盖本研究的所有时间段。过去的气候情景选择第四纪的末次冰盛期(Last Glacial Maximum)和全新世中期(Mid Holocene)，末次冰盛期是距今最近的极寒冷时期，气候较当前偏冷偏干，而全新世中期是距今最近的一个暖期，气候较现在偏暖偏湿(Bartleinetal.,2011)。未来的气候情景选择2070年RCP2.6和RCP8.5两种典型浓度路径(representative concentration pathways,RCPs)，RCP2.6情景代表全球表面平均温度有可能比1986-2005年期间上升0.3～1.7℃，RCP8.5情景则代表全球表面平均温度可能升高2.6～4.8℃(IPCC,2013)。

1.3 模型模拟参数

MaxEnt模型的复杂度对物种的潜在分布预测有较大的影响(朱耿平和乔慧捷,2016)，其中调控倍频(regularization multiplier,RM)和特征组合(feature combination,FC)是影响模型复杂度的主要参数(Elithetal.,2011)。调控倍频默认参数为1，该参数是对模型进行平滑处理的参数，避免模型过度拟合，设置的值越大则模型预测的泛化能力越强，高热点地区更为突出。特征组合有5种，分别为：线性(linear,L)、二次型(quadratic,Q)、片段化(hinge,H)、乘积型(product,P)和阈值性(threshold,T)，在默认情况下模型通过物种分布点的样本量来选择特征组合的方式。为了选取最优的模型，我们将调控倍频设置为0.5-5，间隔为0.5，特征组合设置为L,LQ,H,LQH,LQHP和LQHPT，使用ENMeval包中的区块法(block)在R3.6.1中获得以上组合测试返回的修正AIC值(delta.AICc)，选取修正AIC值最小时对应的参数组合进行模型模拟(Muscarellaetal.,2014)。

MaxEnt模型要求输入的环境变量图层拥有一致的空间分辨率、投影坐标系和图层范围，因此在Arcgis10.7中将环境图层数据统一为2.5 arc-minute的空间分辨率，再对图层进行相同的投影。为了提高模型的预测准确性，将环境图层统一裁剪至稍大于大蜜蜂所有分布范围的地区；除调控倍频和特征组合外，设置随机测试数据集(random test percentage)为30%，重复运行(replicates)5次，并设置运行模型自带的刀切法(Jackknife)和响应曲线(response curve)来评估环境变量数据；其他参数均为模型默认设置。将313个大蜜蜂物种分布点和13个环境变量数据导入MaxEnt模型(3.4.1版本)(Phillipsetal.,2017)，模拟大蜜蜂在当前的潜在适生区；将大蜜蜂物种分布点、9个气候变量数据和CCSM4导入MaxEnt模型，模拟大蜜蜂在过去(末次冰盛期、全新世中期)、当前和未来(2070年)RCP2.6、RCP8.5情景下的潜在适生区。

1.4 模型评估和输出

AUC(area under the curve,AUC)值即受试者工作特征(receiver operating characteristic,ROC)曲线下的面积，是一种用于评估模型准确性的方法，它能独立于阈值而预测模型的精度，但也存在局限性(朱耿平等,2017)。本研究采用PartialROC方案计算AUC比率(AUC ratio)来评估MaxEnt模型的准确性，弥补了AUC值评估的不足。由于本研究中大多数据来自公共数据库，因此判定90%的测试数据是可靠的，使用ENMGadgets包计算AUC比率(Petersonetal.,2008;朱耿平等,2017)。AUC比率>1时，说明模型相对随机预测较好，反之则预测较差。

MaxEnt模型自动输出5次模拟结果的平均值，输出结果为0-1的大蜜蜂适生逻辑指数，值越高则大蜜蜂的生境适宜指数越高。以10%训练存在阈值(10 percentile training presence threshold)作为划分大蜜蜂适生区及非适生区的阈值，再按逻辑值对适生区进行等值划分。最低训练存在阈值(minimum training presence threshold)是拥有校准位置的任何像素预测的最小值，而10%训练存在阈值则是将排除具有最低预测值10%位置的值设置为阈值，已被广泛地应用于生态位模型中(Nazerietal.,2012)。为了更直观地显示不同时期、不同情境下大蜜蜂适生区的变化，使用Arcgis 10.7中的分区几何统计工具，计算并显示大蜜蜂潜在适生区的几何质心(centroid)位置变化，探究大蜜蜂适生区的空间格局变化。

使用模型输出的环境变量贡献率(percent contribution)、置换重要值(permutation importance)、刀切法检验以及响应曲线对参与模拟的环境变量进行评估。贡献率表示的是在模型迭代过程中，对模拟位点拟合结果的贡献，贡献率越高，对模拟位点拟合结果的贡献越大。置换重要值来自于训练位点和背景数据随机采集变量的排序，排序重要性值越大，位点采集到的环境变量比例越大。刀切法检验描述的是环境变量的贡献值，当单独考虑某一环境因子时刀切法的值越高则该环境变量越重要。响应曲线选取不考虑其他变量影响下，单个环境变量的变化，横坐标表示该环境变量的变化范围，纵坐标表示适生指数。

2 结果

2.1 模型准确性评估

使用13个环境变量时，ENMeval包模拟结果显示：修正AIC值最小时对应的模型参数优化设置是：特征组合为线性和二次型(L和Q)，调控倍率为0.5，运行后模型的AUC比率为1.31。使用9个气候变量时，模型模拟修正AIC值最小时对应的模型参数优化设置是：特征组合为线性、二次型、片段化、乘积型和阈值性(L,Q,H,P和T)，调控倍率为1.5，运行后模型的AUC比率为1.33。AUC比率显示MaxEnt模型对大蜜蜂的潜在适生区模拟结果具有较高的准确性。

2.2 大蜜蜂当前的实际分布点和潜在适生区

大蜜蜂分布于整个东南亚地区，并向西延至印度河(Hepburn and Radloff,2011)。本研究中涉及的大蜜蜂物种分布数据西起印度半岛，向东沿着喜马拉雅山脉到中国南部，并向南延至中南半岛和马来群岛，物种数据分布点已覆盖了大蜜蜂的自然分布地区(图2:A)。

图2 大蜜蜂的实际分布点与MaxEnt模型预测的当前潜在适生区Fig.2 Presence data of Apis dorsata and its current potential suitable distribution areas predicted by MaxEntA:用于MaxEnt模型中3个不同来源的大蜜蜂物种空间稀释后的分布点Filtered presence data of A. dorsata spatially from three different sources used for MaxEnt;B:基于9个气候变量的大蜜蜂当前潜在适生区Current potential suitable distribution areas of A. dorsata modeled by MaxEnt based on nine climate variables;C:基于13个环境变量的大蜜蜂当前潜在适生区Current potential suitable distribution areas of A. dorsata modeled by MaxEnt based on 13 environmental variables。审图号Map Approval No.:GS(2021)5543.

9个气候变量的模拟结果显示大蜜蜂当前的潜在适生区为：印度半岛，包括印度中南部、印度北部沿喜马拉雅山脉经过尼泊尔和孟加拉国，到中国南部包括云南、贵州、广西和海南，中南半岛除缅甸部分地区、泰国东部以外的全部区域，马来群岛除婆罗洲中北部、苏门答腊岛西部的全部地区(图2:B)。而中高潜在适生区则在印度半岛南部的西高止山脉(Western Ghats)、东高止山脉(Eastern Ghats)，印度北部的喜马拉雅山脉(Himalayas)，中国、缅甸、泰国和老挝的交界处的横断山脉(Hengduan Mountains)的南部地区，越南东南部的长山山脉(Annamite Mountains)地区，缅甸、泰国、柬埔寨的南部沿海地区，多山的菲律宾群岛、爪哇岛和苏拉威西岛。当加入人类影响、海拔、森林覆盖和林冠高度数据，整合13个环境变量后分析显示，大蜜蜂当前的潜在适生区显著收缩，且呈破碎化趋势(图2:C)。

2.3 影响大蜜蜂分布的关键因子

模型模拟过程中显示了13个环境变量的贡献率和置换重要性，从表2中可以看出13个环境变量中重要的变量为人类影响(human influence index,55.8%)、温度季节性变化(Bio4,24.5%)、最冷季均温(Bio11,4.6%)、等温性(Bio3,3.1%)及海拔(elevation,3.3%)，它们的累积贡献值达91.3%。这些因子也具有较高的置换重要性，累计贡献值达83.1%(表2)。

表2 基于MaxEnt模型的13个环境变量的贡献值和排序重要性Table 2 Percent contribution and permutation importance modeled by MaxEnt with 13 environmental variables

13个环境变量的刀切法检验结果显示：当使用某个单独环境变量时，AUC值显示重要的环境变量按重要性排前3的为：人类影响(human influence index)、等温性(Bio3)和温度季节性变化(Bio4)(图3)。综合模拟结果表明：人类影响是最主要的环境因子，刀切法的结果也显示排除人类影响后对结果的影响很大；温度季节性变化、等温性、最冷季均温是影响大蜜蜂潜在适生区的3个最主要气候因子，而海拔也对大蜜蜂的潜在适生区结果有一定影响。森林覆盖和林冠高度这两个栖息地因子对大蜜蜂潜在适生区的模拟结果影响较小。

图3 基于刀切法检验评价环境变量重要性Fig.3 Assessment of the importance of environmental variables based on Jackknife test

不考虑其他变量影响的情况下，随着单个环境变量的变化，影响大蜜蜂潜在适生区的4个重要环境变量的响应曲线显示大蜜蜂的潜在适生区要求：等温性在35%～100%之间，在70%左右适生指数最高，等温性是量化昼夜温度相对于年度温度的波动幅度，低于100%则表示日温度变化对于年温度变化幅度较小。温度季节性变化在0～6500%之间，在1600%左右适生指数最高，温度季节性变化是一年中温度变化的度量，百分比越大，则温度的变异性就越大。最冷季均温较高，最冷季均温是一年中最冷的3个月的平均温度。海拔小于3 000 m，在1 000 m左右适生指数最高(图4)。综上所述，大蜜蜂最适宜的气候条件为：全年高温，最冷季的均温高，日温差、年温差波动小，此外大蜜蜂的分布限于海拔3 000 m以下的地区，而在海拔1 000 m左右其出现的概率最大。

图4 MaxEnt模型预测影响大蜜蜂潜在适生区的4个重要环境变量的响应曲线Fig.4 Response curves of four important environmental variables affecting the potential suitable distribution areas of Apis dorsata predicted by MaxEntA:等温性Isothermality;B:温度季节性变化Temperature seasonality;C:最冷季均温Mean temperature of the coldest quarter;D:海拔Elevation.

2.4 过去、当前和未来气候变化下大蜜蜂潜在适生区的变化

因为人类影响、海拔、森林覆盖和林冠高度这些环境变量没有基于过去和未来的模型模拟数据，为了保持模型的一致，只用了9个气候因子。比较过去、当前和未来不同时期大蜜蜂的潜在适生区面积(图5)，并以不同时期的低、中、高潜在适生区的平均适生指数和面积计算不同时期的大蜜蜂潜在适生指数，结果显示：在末次冰盛期，大蜜蜂的潜在适生区比当前扩大了17.8%，适生指数扩大24.7%；在全新世中期，大蜜蜂的潜在适生区比当前缩小了8.9%，适生指数减少12.0%；在2070年的RCP2.6情景下，大蜜蜂的潜在适生区比当前缩小了1.6%，适生指数扩大6.1%；在2070年的RCP8.5情景下，大蜜蜂的潜在适生区比当前扩张了2.1%，适生指数扩大17.6%。不同时期大蜜蜂的潜在适生区及其提取的质心结果表明：末次冰盛期大蜜蜂的潜在适生区相比全新世中期偏向西北方向，过去的两个时期其潜在适生区相比当前更偏向马来群岛和中南半岛；在未来，2070年的RCP2.6情景下，大蜜蜂的潜在适生区与当前接近，稍偏向于马来群岛和中南半岛；而在RCP8.5情景下，潜在适生区面积与当前接近，但显著偏向印度半岛。

图5 MaxEnt模型预测大蜜蜂过去、当前和未来的潜在适生区Fig.5 Potential suitable distribution areas of Apis dorsata in the past,present,and future predicted by MaxEntA:末次冰盛期Last Glacial Maximum;B:全新世中期 Mid Holocene;C:当前 Present;D:RCP2.6;E:RCP8.5;F:不同时期潜在适宜分布区质心的变化Shift of distribution centroid under different scenarios.审图号Map Approval No.:GS(2021)5543.

2.5 大蜜蜂潜在适生区的生态区

将大蜜蜂当前潜在适生区匹配至陆地生态区(图6)，结果显示：大蜜蜂的中高潜在适生区主要位于南亚和东南亚地区湿润的热带雨林、热带季节性雨林和低地雨林，并通过与这些地区相邻的常绿森林连接、扩张其潜在适生区，这些地区的森林组成通常为常绿和半常绿落叶乔木，组成以龙脑香科(Dipterocarp)植物为主(Scrivenetal.,2015)。大蜜蜂的非潜在适生区是一些干旱林，如中南半岛干旱林、缅甸中部的伊洛瓦底干旱林、印度南部的德干高原荆棘灌木林等。在人类影响下，大蜜蜂的中高潜在适生区有向山区和大片连续的湿润常绿森林区收缩的趋势。

图6 基于陆地生态区MaxEnt模型预测的大蜜蜂潜在适生区分布Fig.6 Potential suitable distribution areas of Apis dorsata predicted by MaxEnt based on the terrestrial ecoregionA:研究地区的陆地生态区Terrestrial ecoregion of our study region;B:基于9个气候变量的大蜜蜂潜在适生区Potential suitable distribution areas of A. dorsata based on nine climate variables;C:基于13个环境变量的大蜜蜂潜在适生区Potential suitable distribution areas of A. dorsata based on 13 environmental variables.审图号Map Approval No.:GS(2021)5543.

3 讨论

3.1 影响大蜜蜂分布的关键因子

本研究筛选出影响大蜜蜂分布的关键因子，人类影响最大，其次是气候因子等温性、温度季节性变化、最冷季均温以及海拔；而森林覆盖和林冠高度在本研究中对大蜜蜂潜在适生区的影响较小，这与Soberón等(2005)提出的观点一致，在大尺度上主导物种分布格局的是以气候因子为主的非生物因子和生物地理分布区等，而物种间相互作用则被认为在较小的空间尺度上作用于物种分布(Soberón and Peterson,2005)。模拟结果显示温差和最冷季均温是限制大蜜蜂潜在分布的最主要气候因子，大蜜蜂潜在适生区要求全年温差不大，最冷季均温高。大蜜蜂潜在适生区覆盖的地区为南亚和东南亚地区的热带雨林和热带季风气候区，其特点均为全年高温，温差小。然而，热带季风气候相比热带雨林气候有明显的干旱季之分，最冷季均温较低(Kotteketal.,2006)；结合基于气候因子模拟的大蜜蜂潜在适生区在印度半岛和中南半岛的西部、北部边界均与热带季风气候带的边界相近，这可能表示大蜜蜂在内陆地区无法跨越的障碍是上述的气候条件。

在加入人类影响后，模拟结果显示大蜜蜂的中高潜在适生区显著减少，并有向山区收缩的趋势。Matias等(2019)的研究表明，在菲律宾的巴拉望岛上由于受到人类活动影响，大蜜蜂基本都迁移到了山上，很少能在低地见到(Matiasetal.,2019)。过去的研究显示大蜜蜂主要分布在海拔1 000 m以下的地区，部分分布至1 000-3 000 m之间的地区，显示大蜜蜂在海拔梯度上有一定适应性(Hepburn and Radloff,2011)。本研究模拟结果显示单一考虑海拔因素时，大蜜蜂的分布偏向于海拔1 000 m左右的区域，1 000 m以上潜在适生区的适生指数增加，但是大蜜蜂向高处迁移的区域有限；加之人类对大蜜蜂的影响加剧，其适生区可能将急剧减少，因此海拔3 000 m以下的南亚和东南亚的低地森林将是大蜜蜂的重点保护地区。

Orr等(2021)对全球蜜蜂的分布研究显示：全球蜜蜂的丰富度与森林呈负相关，热带地区不是蜜蜂分布最丰富的地区，也是研究相对薄弱的地区(Orretal.,2021)。大蜜蜂是南亚和东南亚的热带生态系统中重要的、泛化的传粉昆虫(Raghunandan and Basavarajappa,2014;Matiasetal.,2019)。森林是大蜜蜂的主要栖息地，本研究也考虑了森林覆盖和林冠高度这两个生物因素对大蜜蜂潜在适生区的影响，结果表明这两个环境因素对大蜜蜂的潜在适生区影响较小。可能反映了栖息地变量需要在小尺度区域性上研究，才能更好地解释；进一步可结合蜜粉源植物、土地利用等开展深入研究。

3.2 历史气候下大蜜蜂分布格局的变化

关于蜜蜂属起源的研究相对较少，早期的研究提出蜜蜂族起源于中国的喜马拉雅和横断山脉地区然后向南迁移(匡邦郁和李有泉,1988)。Cao等(2012)的系统发育研究显示：分布于中国大陆的大蜜蜂种群属于广泛分布的A.dorsatadorsata，菲律宾巴拉望岛和东、西马来西亚的大蜜蜂种群都属于中国大陆的种群，同时东、西马来西亚的大蜜蜂种群不仅和中国大陆的种群聚一起，还和尼泊尔、泰国的种群聚一起，具有高度的遗传多样性；推测大蜜蜂可能在第四纪冰川时期通过陆桥从中国南部向南迁移至东南亚地区，东南亚可能是第四纪冰期大蜜蜂的避难所(Caoetal.,2012)。

婆罗洲和中南半岛也被证实至少在过去的5 Ma期间一直是巽他古陆动植物类群的主要避难所(De Bruynetal.,2014)。本研究模拟的结果也显示在末次冰盛期，东南亚的气候适合大蜜蜂生存，尤其是中南半岛和马来群岛。第四纪中晚期巽他古陆海平面重建结果显示，当时最多的情景为海平面低于当前40-50 m时，常绿热带雨林从婆罗洲一直延伸到苏门答腊岛；其次的情景为类似于末次冰盛期的海平面非常低的时期，热带季节性雨林分布广泛，末次冰盛期可能还存在从中南半岛到爪哇岛的半常绿森林扩散走廊；最后，与当前的海平面情景相似时，常绿热带雨林从克拉地峡(Kra Isthmus)一直延伸到爪哇岛西部，包括婆罗洲地区(De Bruynetal.,2014)。因此，第四纪中晚期有适宜大蜜蜂生存的常绿和半常绿森林，大蜜蜂能从中南半岛一直向爪哇岛扩散。综上所述，冰期的寒冷气候使得大蜜蜂随着低海平面时陆架上的常绿雨林、热带季节性雨林而迁移到气候适宜的东南亚地区，并沿着克拉地峡向南一直到达爪哇岛，再向东通过陆架到达菲律宾群岛和苏拉威西岛。因为大蜜蜂的迁徙有限，无法到达巽他古陆右侧没有陆架的地区(Hepburn and Radloff,2011)。第四纪冰期和间冰期可能奠定了大蜜蜂在东南亚的主要分布格局，在末次盛冰期之后，大蜜蜂的分布随着气候的动荡发生了较小的收缩，至当前时期的气候条件下，又发生了显著的扩张。就蜜蜂属最早出现在更新世来说，横断山南部和喜马拉雅山区可能是大蜜蜂的起源地，但这需要更多的系统发育和分子证据来证明大蜜蜂的起源时间、地点及辐射演化。

3.3 未来气候下大蜜蜂潜在适生区的改变

有研究表明分布于温带的熊蜂Bombusbellicosus，因为偏好寒冷气候，随着气候的变暖而向高海拔地区、高纬度地区迁移(Martinsetal.,2015;Pykeetal.,2016)。对于分布于新热带巴西的一种兰花蜂Eulaemanigrita，能为兰花和另外几种植物传粉，在未来气候变暖下潜在适生区发生了扩张(Silvaetal.,2015)。横跨温带、热带、亚热带，适应干旱、分布广泛及取食泛化的澳大利亚木蜂Ceratinaaustralensis能适应气候变暖，潜在栖息地也有所扩张(Dewetal.,2019)。这些结果显示，泛化传粉的不同蜜蜂种类均面对全球气候的变化有积极响应。我们的模拟结果表明未来气候变暖、干旱的情况下，大蜜蜂能积极地应对未来气候变暖，其潜在适生区与当前相近，适生指数增加。除了有泛化的蜜粉源植物可依赖，与大蜜蜂的生物学特性也息息相关，大蜜蜂营露天生活，其工蜂花费大量的时间和精力来维持巢上大蜜蜂帘，让群体适应天气的变化(Seeley,1982)；此外大蜜蜂随季节和蜜粉源的变化具有迁徙特性，适应环境能力强(董诗浩等,2014)。因此，大蜜蜂对气候变化不仅积极响应，还能调整、适应。

然而，根据加入人类影响变量后大蜜蜂当前适生区的剧烈变化，可以推测大蜜蜂在未来的潜在适生区的处境并不乐观。目前在南亚和东南亚地区人类活动加剧、土地利用改变、生境退化，加之气候改变正对生物多样性构成极大的威胁(Devy and Davidar,2003;Sodhietal.,2004;Sodhietal.,2010；Corlett,2011)。此外，多个区域性调查证据也显示极端高温高湿环境对大蜜蜂的生存并不利，例如在巴基斯坦的研究显示大蜜蜂的月丰度和相对湿度呈显著的负相关(Sajjadetal.,2017)；在印度西北部干旱地区的观察则发现高温干燥了花蜜，会直接影响采蜜的大蜜蜂(Sihag,2014)；而在苏拉威西岛上高温高湿的条件下，大蜜蜂活动减少(Nagiretal.,2016)。因此，本研究模拟结果显示：尽管在2070年随着全球变暖，大蜜蜂的潜在适生区有所扩大，如果考虑到人类影响的作用，大蜜蜂将仍然面临栖息地减少、分布区改变的危机。

3.4 生态位模型在大蜜蜂保护中的应用

影响物种分布主要有3个因素：环境条件、物种相互作用与物种迁移能力。在没有物种作用情况下的物种能够生存的环境条件为物种的基础生态位；而包括了物种之间的相互作用下的物种能够生存且繁殖的条件为物种的现实生态位。有学者认为，MaxEnt模型模拟的是物种的基础生态位，而不是现实生态位，因为生物的相互作用没有被考虑在模型中(Stockwell,1999)；而又有学者认为用于建立模型的物种数据已经包括了在其原始环境下物种相互作用的信息，因此模拟出的应该就是物种的现实生态位(Austinetal.,2002;Guisan and Thuiller,2005)。而实际上应用于物种保护目的生态位模拟趋向于模拟物种的现实生态位，更能反映出物种在现实环境下集合了物种相互作用变量后的实际分布趋势。因此，在我们的研究中，考虑了除气候变量之外的栖息地变量以及人类影响变量，这对于探究影响大蜜蜂生态位分布的主要因素具有重要作用，更能综合地反映出大蜜蜂未来分布受威胁的主要因素。此外，生态位模型的模拟是建立在物种生态位保守的理论基础上，即同一物种在不同时间上基础生态位保持一致和稳定(Petersonetal.,1999)。因而在我们的研究中我们默认基于气候条件的大蜜蜂过去和未来潜在适生区的基础生态位与当前的基础生态位一致，并将其映射至不同的空间和时间上。尽管实际上物种是趋于进化的，由于我们预测的未来时间为2070年，在一个比较短的时期内，物种的进化变异有限，因此我们的模拟预测依然具有较高的准确性。

本研究通过比较仅基于13个环境因子模拟出的大蜜蜂的基础生态位与考虑栖息地变量和人类活动影响的现实生态位的结果显示：人类活动对大蜜蜂影响较大。目前东南亚大规模的人为栖息地改造、森林火灾和野生动物的过度猎取等正对当地的生物多样性构成极大的威胁(Sodhietal.,2004)；而南亚地区在20世纪初就已经有许多低地热带雨林消失(Devy and Davidar,2003)，这意味着尽管未来南亚、东南亚的气候条件适宜大蜜蜂生存，但是适宜大蜜蜂生存的生境正日益减少。在2070年RCP2.6情景下时，大蜜蜂的潜在适生区在印度半岛南部稍有减少，而适生指数在马来群岛显著增加。在RRCP8.5情景下，大蜜蜂的潜在适生区在印度半岛南部增加，中南半岛南部和马来群岛的适生指数增加显著。因此这些地区都是未来气候环境适宜大蜜蜂生存的地区，需要特别关注。此外，缅甸也是值得关注的地区，尽管当前的模拟结果显示缅甸的大蜜蜂潜在适生区较少，这可能是因为使用的缅甸地区大蜜蜂分布点太少，这也是MaxEnt模型使用物种出现数据的缺点，因为模拟结果会偏向包含更多物种分布点的区域(Kramer-Schadtetal.,2013)。从生态区和气候带来看缅甸的大部分地区都符合大蜜蜂的适生区要求，需要加强在缅甸地区对大蜜蜂的调查和研究。