菜豆象Acanthoscelides obtectus（Say），属鞘翅目Coleoptera豆象科Bruchidae，是世界范围内重要的仓储害虫，每年因豆类粮食被侵蚀而造成严重经济损失（Padı́n et al.，2002）。由于菜豆象成虫和幼虫在豆粒内部蛀食，该虫极易随寄主豆粒传播，而全球化贸易进一步加速了该虫的跨境传播和扩散（王辉等，2013）。口岸截获数据显示，2021年中国大陆口岸菜豆象截获量同比增长37%（海关总署数据）。部分研究者基于MaxEnt模型对菜豆象在中国的适生区展开预测，结果表明菜豆象在中国的适生区呈现向北扩张的趋势（易山青等，2023）。面对日益严峻的菜豆象扩散疫情，研究者亟需从种群遗传学角度阐明该虫的种群遗传结构及其演变动态，为区域性综合治理防控策略制定提供理论依据。

种群遗传结构解析是揭示物种适应性进化与扩散机制的核心手段。线粒体细胞色素氧化酶Ⅰ（COI）基因作为分子标记，因其母系遗传特性、进化速率适中及高度保守性，被广泛应用于昆虫系统发育重建和种群遗传结构分析（Hebert et al.，2003；Tembrock et al.，2023）。近年来，基于COI的DNA条形码技术已成功应用于多种害虫的种群遗传学研究，旨在揭示害虫种群结构和扩散路径（游丽斌等，2019；Tembrock et al.，2023）。已有研究借助COI序列分子表明菜豆象起源于中美洲，伴随着欧美的殖民和贸易活动，先后传入南美洲和欧洲、非洲和亚洲，并且显示菜豆象非洲与南美洲种群间存在显著遗传分化（Alvarez et al.，2005；Oliveira et al.，2013）。但是，有限的种群样本不能全面反映国内外菜豆象种群特征，尤其是我国菜豆象种群遗传特征。因此，全球尤其是亚洲菜豆象种群遗传分化、种群结构以及潜在入侵路径有待系统阐明。

本研究以南美洲（智利）、非洲（9个国家）、欧洲（德国和芬兰）和亚洲（中国）的菜豆象种群为对象，通过测序获取COI基因条形码序列并结合公开数据库中COI序列，旨在解决以下科学问题：（1）菜豆象不同地理种群的遗传多样性如何？（2）不同地理种群间遗传分化程度如何？（3）不同地理种群间遗传分化是否与地理距离相关联？（4）不同地理种群规模如何变化？研究结果将为揭示菜豆象的适应性进化潜力，为构建基于遗传信息的区域性检疫策略、不同来源的种群检测及可持续防控体系提供理论支撑。

1.   材料与方法

1.1   样本采集

本研究共采集到161个菜豆象样本（75个为国内样本，86个为国外样本），涉及21个采样点（9个为国内采样点，12个为国外采样点），详见表 1。采集到的样本，经过编号后置于对应收集管中备用，样本原始记录详见网络版增强出版材料附表 1。

1.2   DNA提取和PCR扩增

对于采集的样本，所有样品单头提取DNA，使用凯杰DNeasy Blood & Tissue Kit试剂盒（凯杰，上海）按说明书步骤操作，纯化后的DNA用80~100 µL AE缓冲液洗脱。本研究使用已经发表的文章中COI基因通用扩增引物对目标序列进行扩增，详细引物信息见下表 2（Hebert et al.，2003）。PCR扩增反应使用湖南艾科瑞生物工程公司的PCR试剂，PCR反应总体系为30 µL，包括10×PCR Buffer（含Mg²⁺）3 µL，dNTPs 2.4 µL，上下游引物（20 pmol/µL）各0.5 µL，1 µL DNA模板（10~100 ng/µL），22.4 µL ddH₂O，0.2 µL Taq酶。PCR反应程序为95℃预变性5 min，94℃变性1 min，45℃退火1.5 min，72℃延伸1.5 min，以上3步循环5次；再94℃变性1 min，50℃退火1.5 min，72℃延伸1 min，共36个循环，72℃终延伸10 min。PCR扩增产物直接送往测序公司测序。

1.3   序列拼接和对齐

对于COI基因测序结果，使用DNAStar中SeqMan模块进行正负链拼接。在删除测序结果中首末端不准确碱基，同时剔除拼接质量不佳的序列后，获得153条一致性序列，长度为645 bp。使用测序的COI基因序列比对NCBI数据库，从比对结果中筛选出片段大小合适、有明确采样地点的数据，共7条序列。合并测序序列，共获得160条一致性序列。将序列导入MEGA v10.1.6软件，使用内置的Muscle方法进行序列对齐排列，手动删除首末端未平齐部分，获得平齐的一致性序列，长度为642 bp（Kumar et al.，2018）。

1.4   进化分析

使用jModelTest2软件导入对齐后的一致性序列，选择AIC准则运行软件，得到构建进化树的最佳模型即HKY模型（Darriba et al.，2012）。使用MEGA v10.1.6软件导入对齐后的一致性序列，选择NJ（遗传距离计算模型：p-distance）和ML（核苷酸替换模型：HKY）两种模型构建系统发育树，试验方法为自展法（bootstrap），迭代1 000次（Felsenstein，1985）。使用MrBayes v3.1.2进行Bayesian建树分析（BI）（Ronquist and Huelsenbeck, 2003）。马尔科夫链（MCMC）运算总代数为1 000 000代，分析是从一个随机树开始，每500代抽样一次，核苷酸替换模型为HKY。分裂频率的平均标准偏差小于0.01，说明后验分布的抽样是充分的。使用潜在规模缩减系数（PSRF）和分裂频率的平均标准偏差验收收敛性。

1.5   单倍型分析

使用DnaSP 6.12.03软件导入对齐后的菜豆象COI一致性序列，生成单倍型文件（Rozas et al.，2017）。结合Python脚本处理单倍型文件，使用PopART软件导入标准NEXUS文件，选择Median-Joining Network，生成单倍型网络（Leigh and Bryant, 2015）。使用Python脚本处理单倍型文件。使用RStudio导入单倍型文件，载入R包（ggplot2、reshape2），编写并且运行R脚本绘制金字塔图。

1.6   种群动态分析

使用DnaSP 6.12.03软件导入按照采集地分组后的菜豆象COI一致性序列，计算单倍型多样性、核苷酸多样性和Tajima's D值。使用Alrequin 3.5.2.2软件导入经过DnaSP软件预处理后的标准输入文件，计算固定系数Fst值，同时根据公式Nm = (1 - Fst) / (4 * Fst)计算基因流Nm（Wright，1931；Wright，1951；Excoffier and Lischer, 2010）。

2.   结果与分析

2.1   菜豆象种群的单倍型地理分布和丰度

分析结果表明，菜豆象COI序列总共有15种单倍型H1~H15（表 3），其中单倍型Hap_2分布广泛，遍及欧洲、亚洲和非洲。总体上，国内外单倍型分布不均衡，显示出地理偏好性，其中国内云南和贵州仅发现为Hap_1、Hap_2、Hap_3、Hap_4型，且Hap_1和Hap_2为优势型（表 3），而国外12个样点发现有Hap_2、Hap_5、Hap_6、Hap_7、Hap_8、Hap_9、Hap_10、Hap_11、Hap_12、Hap_13、Hap_14、Hap_15，且Hap_2为优势型（表 3、图 1）。

图  1  国内外菜豆象COI单倍型丰度差异

注：纵坐标代表单倍型类型。横坐标代表样本数量。图中数字代表具体样本数，红色代表国内，绿色代表国外。

Fig.  1  Differences in COI haplotype abundance of Acanthoscelides obtectus (Say) between domestic and foreign

Note: The ordinate represented haplotype type. The abscissa represented the number of sample. The number in the figure represented the specific sample number. Red represented the domestic. Green represented the foreign.

下载: 全尺寸图片

2.2   菜豆象COI序列单倍型网络和进化树

按照地理种群对COI单倍型类型和数量进行分类统计，统计结果表明不同地理种群间单倍型多样性具有差异性（表 3），其中刚果（CG）、喀麦隆（CMR）和埃塞俄比亚（ETH）的单倍型多样性最高，各有4种单倍型，而贵州省毕节市纳雍县雍熙（GZBJNY）、贵州省毕节市织金县绮陌街道（GZZJ）、尼日利亚（NGR）、卢旺达（RWA）和埃及（EGY）的单倍型多样性最低，仅各有1种单倍型。

单倍型网络结果显示，菜豆象单倍型网络结构为多中心点辐射模式，其单倍型形成三个线粒体群组，分别以Hap_2、Hap_5和Hap_6为祖先单倍型（图 2）（Oliveira et al.，2013）。Hap_1、Hap_4、Hap_8、Hap_11、Hap_13、Hap_14由Hap_2经过突变后形成。Hap_9由Hap_6经过突变后形成。Hap_3、Hap_7、Hap_10、Hap_12、Hap_15由Hap_5突变形成。进化树结果显示，菜豆象COI单倍型序列聚集成两个Clade，形成三个亚类群。CladeⅠ是Hap_2和Hap_6线粒体群组，CladeⅡ是Hap_5线粒体群组。CladeⅠ形成两个亚类群，其中Hap_1、Hap_2、Hap_4、Hap_8、Hap_11、Hap_13、Hap_14聚成一簇即亚群Ⅰ，而Hap_6、Hap_9聚成另一簇即亚群Ⅱ。CladeⅡ形成一个亚类群，Hap_3、Hap_5、Hap_7、Hap_10、Hap_12、Hap_15聚成一簇即亚群Ⅲ。进化树聚类结果与单倍型网络结果一致（图 3、详见网络版增强出版材料附图 1和附图 2）。

图  2  COI序列的单倍型网络

注：彩色饼图的大小代表样本数量，而每个采样地点用不同颜色进行标注。Hap_后缀的文本代表每种单倍型。右边图示中空心圆形黑框代表样本数量的梯度。右边图示中每个彩色实心圆和黑色字母代表每个采样地点。

Fig.  2  Haplotypes network of COI sequences

Note: The size of the color pie chart represented sample number, and each sampling site was marked with a different color. The text suffixed with 'Hap_' represented each haplotype. The hollow black box on the right represented the gradient of the number of samples. Each colored solid circle and black letter in the right figure represented each sampling site.

下载: 全尺寸图片

图  3  菜豆象COI单倍型ML进化树

注：分类群名代表种群和单倍型。比例尺代表单位遗传距离。枝长代表遗传距离。进化枝旁的数字代表自展值。不同颜色区分进化枝聚类成的簇，其中蓝色枝为亚群Ⅰ、红色枝为亚群Ⅱ、绿色枝为亚群Ⅲ。建树序列来源信息见附表 3。

Fig.  3  Maximum likelihood evolutionary tree of COI haplotype of Acanthoscelides obtectus (Say)

Note: Taxa names represent haplotype type. Scale bars represented unit genetic distance. Branch length represented genetic distance. Numbers next to clades represented bootstrap values. Different colors were used to distinguish the clusters of clades. The blue clade was subgroup Ⅰ. The red clade was subgroup Ⅱ. The green clade was subgroup Ⅲ. The information of data source of sequence used to construct the evolutionary tree was listed in supplementary Table 3.

下载: 全尺寸图片

2.3   菜豆象的遗传多样性和种群历史动态

菜豆象种群COI单倍型多样性、核苷酸多样性及Tajima's D中性检验结果如表 4所示。结果表明，菜豆象21个地理种群的总体单倍型多样性指数Hd为0.665，种群间核苷酸多样性指数Pi为0.00291。种群内单倍型多样性为0.000~1.000，平均为0.422，其中单倍型多样性最高的种群为FIN（芬兰瓦萨市），单倍型多样性最低的种群为EGY（埃及）、GZBJNY（贵州省毕节市纳雍县雍熙）、NGR（尼日利亚）和RWA（卢旺达）。种群内核苷酸多样性为0.00000~0.00935，平均为0.00219，其中最高的种群为CHI（智利），最低的种群为EGY（埃及）、GZBJNY（贵州省毕节市纳雍县雍熙）、NGR（尼日利亚）和RWA（卢旺达）。结果表明EGY（埃及）、GZBJNY（贵州省毕节市纳雍县雍熙）、NGR（尼日利亚）和RWA（卢旺达）的遗传多样性最低。

中性检验是一种检测自然选择或遗传漂变的有效方法，Tajima's D值＞0，表明群体经历平衡选择，=0表明群体没有经历选择，符合中性进化，而＜0表明群体经历定向选择，可能近期经历种群扩张（Krutovsky and Neale, 2005）。Tajima's D中性检验（表 4）结果显示，除UGA种群外，中性检验结果均不显著，表明该20个菜豆象地理种群在过去没有出现群体扩张和持续增长模式，群体大小保持相对稳定状态。UGA种群的结果为负值（Tajima's D=-1.79107，P < 0.05），表明该种群可能在近期出现了种群扩张（Krutovsky and Neale, 2005）。

此外，通过对21个菜豆象地理种群的遗传进化和基因流分析，研究发现不同地理种群间的Fst值差异较大，取值范围从0~1（详见网络版增强出版材料附表 2），表明不同地理种群间出现不同程度的遗传分化。基因流Nm分析表明，不同地理种群间基因交流程度存在较大差异，比如YNES和GZBJNY间基因交流活跃（Nm=15.0032），而NGR和GZBJNY间无基因交流（Nm=0.0000）。同时，从各个洲种群间的Fst和Nm的尺度看，Fst值和Nm的差异性依旧存在，比如欧洲-非洲（Fst=0.0851、Nm=2.6881）的遗传分化程度最低，基因交流最为活跃，而南美洲-亚洲（Fst=0.7690、Nm=0.0751）的遗传分化程度最高，基因交流最为消极，但整体上比地理种群间的差异性更小（Fst的取值范围0.0851~0.7690）（表 5）。

进一步通过Mantel检验，计算成对地理种群间的遗传距离Fst和各采集地间地理距离的相关性，得到相关系数r=0.2102（P=0.031 < 0.05），表明遗传分化与地理距离间基本不存在相关性。

3.   结论与讨论

本研究基于线粒体COI条形码序列解析国内外菜豆象种群结构和种群遗传分化现象，结果发现菜豆象种群广泛分布于世界各地，遍及亚洲、欧洲、非洲、南美洲，总计21个地点，共有15种单倍型。经统计发现，不同地理种群间单倍型类型和数量具有差异，其中Hap_1、3、4仅在亚洲（中国）发现，Hap_7、8、9、10、12、14、15仅在非洲发现，Hap_2在亚洲、欧洲和非洲有发现，Hap_5在非洲、欧洲和南美洲有发现，Hap_6在非洲和欧洲有发现，Hap_11在非洲和南美洲有发现，Hap_13仅在欧洲有发现，揭示菜豆象世界范围的仓储害虫的特征和地理偏好性。特别是，非洲保留有很多特有的单倍型（Hap_7、8、9、10、12、14、15）。同时，Tajima's D中性检验结果显示仅有非洲的UGA（乌干达）种群的结果为负值，表明该种群经历定向选择，并且可能在近期出现了种群扩张（Krutovsky and Neale, 2005）。此外，根据以往研究报道，菜豆象起源于中美洲，伴随着欧美的殖民和贸易活动，先传入南美洲，之后再传入欧洲、非洲和亚洲（Oliveira et al.，2013）。这些特有单倍型产生的原因可能是因为菜豆象传入非洲之后，为了适应当地高温、干旱环境而产生适应性进化，最终形成非洲特有的单倍型。

单倍型网络和进化分析结果表明菜豆象种群已经出现一定程度的分化，形成三个亚群，并且不同洲的种群也会同属一个亚群。跨越地理隔离的种群聚类可能和人类历史上的贸易航线息息相关。15世纪末哥伦布开辟新航路之后，欧洲殖民者通过三角贸易（欧洲→非洲→美洲）大规模迁移人口和物种，而菜豆象可能伴随三角贸易，漂洋过海，扩散种群（Jordan，2016）。

种群遗传多样性分析发现，菜豆象种群总体单倍型多样性指数Hd为0.665，总体核苷酸多样性指数Pi为0.00291，而种群内单倍型多样性平均为0.422，种群内核苷酸多样性平均为0.00219，整体水平的Pi值小于0.005，表明整体上国内外菜豆象mtDNA的COI基因多态性水平较低，为低核苷酸多态类群（Grant and Bowen, 1998）。菜豆象蛀食行为需要能量供给，而COI基因参与线粒体能量代谢的相关信号通路（Richter-Dennerlein et al.，2016）。本研究推测COI基因的保守性可能是维持能量代谢功能的关键，从而间接地帮助菜豆象保持蛀食豆粒的特性。目前，已经有研究发现将自然产生mtDNA突变的雌性菜豆象与雄性菜豆象进行交配可造成雄性繁殖力降低。未来可能利用该项技术减少菜豆象种群规模，实现菜豆象防控防治（Lea et al，2024）。

群体分化聚焦于时空尺度下不同种群间遗传分化和基因交流。固定系数Fst（Fixation index）是种群和进化遗传学中使用最广泛的描述性统计，为进化提供了重要的见解影响种群内部和种群之间遗传变异结构的过程（Holsinger and Weir, 2009）。基因流Nm可以揭示出群体间可能的基因渗透及影响遗传分化的遗传现象，Nm值大于1说明群体间的基因交流活跃，小于1则说明群体的遗传分化很可能由遗传漂变引起（Wright，1931）。遗传分化和基因流分析发现菜豆象不同地理种群间的Fst值和Nm值存在较大差异，表明不同地理种群间具有一定程度的遗传分化，基因交流程度存在较大差异。通过分析具有代表性种群（种群数量大于3、Fst值不为负），研究发现菜豆象中国种群和欧洲种群之间的遗传分化程度差异较小，整体上遗传分化程度较小，基因交流消极（Fst取值范围为0.2443~0.8477，Nm取值范围0.0449~0.7732），而中国种群和非洲种群之间遗传分化程度差异大，基因交流程度差异大（Fst取值范围为0.0153~1.0000，Nm取值范围为0.0000~16.0792）。这种现象的产生，一方面可能是不同国家之间的贸易措施和检疫政策影响菜豆象种群的扩散。另外一方面，由于非洲种群单倍型种类多，而且存在非洲特有单倍型，某些非洲种群和中国种群之间可能已经出现生殖隔离，从而导致遗传分化和基因交流的差异程度大。此外，通过将各个种群按照洲进行合并数据后计算各个洲种群间的Fst和Nm，结果发现洲-洲种群间的Fst和Nm值均有差异，其中欧洲-非洲（Fst=0.0851、Nm=2.6881）的遗传分化程度最低，基因交流最为活跃，而南美洲-亚洲（Fst=0.7690、Nm=0.0751）的遗传分化程度最高，基因交流最为消极。这种现象可能源自于欧洲在非洲众多殖民地的人员交往密集，而南美洲与亚洲之间因传统贸易和人员交往相对偏少有关。

Mantel检验显示遗传距离与地理距离基本无相关性（r = 0.2102，P = 0.031 < 0.05）。这种现象可能是由于近几百年人类商贸活动引发跨板块基因交流，特别是工业革命后的全球化贸易，通过集装箱运输等方式形成"跳板扩散"，导致种群出现跨大陆单倍型混杂现象，从而掩盖了自然地理隔离的遗传信号（Hulme，2009；Estoup and Guillemaud, 2010；Seebens et al.，2017）。

由于菜豆象的世界性分布、寄主的分布和扩散能力等诸多因素的制约，本研究无法对世界各地的菜豆象进行系统采样，同时受限于采集到的样本量，无法系统描述菜豆象种群间的进化路线，但是本研究的结果为进一步探讨国内外菜豆象的种群遗传结构、入侵来源和入侵路径奠定了基础。

本研究基于COI基因条形码序列，通过人工采样测序并结合NCBI公共数据，对菜豆象21个地理种群共160个样本序列进行了系统的种群结构分析，得到如下结论：

1）地理分布分析发现，菜豆象种群广泛分布于世界各地，遍及亚洲、欧洲、非洲、南美洲，并且菜豆象的单倍型种群分布不均衡，显示出地理偏好性。

2）遗传多样性分析发现，菜豆象mtDNA COI基因片段在整体水平上具有较高的单倍型多样性，但是具有低核苷酸多样性。同时，种群间的单倍型多样性和核苷酸多样性具有较大差异，显示一定的种群分布。

3）单倍型网络和进化分析发现，菜豆象单倍型网络结构为中心点辐射模式。菜豆象种群间出现一定程度的分化，具有一定程度的聚类，最终形成三个亚群。

4）遗传分化和基因流分析发现，21个菜豆象地理种群间的Fst值和Nm值存在较大差异，表明不同地理种群间具有一定程度的遗传分化，基因交流程度具有较大差异。同时，Mantel测试结果显示种群间遗传距离与地理距离之间基本不存在相关性，表明地理隔离对遗传分化的基本无影响。

5）中性检验结果显示除UGA种群外，其它20个菜豆象地理种群在过去没有出现群体扩张和持续增长模式，群体大小保持相对稳定状态。

附录：附表 1 菜豆象原始采样记录

附表 2 种群间地理距离、固定系数Fst和基因流Nm

附表 3 用于构建进化树的序列的数据源

附图 1 菜豆象COI单倍型NJ进化树

附图 2 菜豆象COI单倍型贝叶斯树

详细数据见网络版增强出版材料（http://hjkcxb.alljournals.net/）

附图 1  菜豆象COⅠ单倍型NJ进化树

Appendix Fig. 1  Evolutionary tree of the COⅠ haplotype of Acanthoscelides obtectus

下载: 全尺寸图片

附图 2  菜豆象COⅠ单倍型贝叶斯树

Appendix Fig. 2  Bayesian tree of COⅠ haplotypes of Acanthoscelides obtectus

下载: 全尺寸图片