基于微卫星标记探讨同域分布竹蝗间的基因渐渗

高 尚，蒋国芳*，范冰冰

(1. 泉州师范学院海洋与食品学院，福建泉州 362000；2. 南京师范大学生命科学学院，南京 210023)

渐渗是指一个种的遗传物质通过杂交与反复回交穿越种间障碍转入到另一个物种内的现象，又称渐渗杂交(Introgressive hybridization) (Anderson，1953)。自然界中10%的动物和25%的植物都存在着杂交现象(Mallet, 2005)，而杂交常常会导致基因渐渗(Bowley & Taylor, 1987)。国外对于蝗虫基因渐渗的研究已有一些报道。Arnold等(1987)使用限制性酶切技术最早发现澳洲1种蝗虫Calediacaptiva的2个亚种在1个狭窄的杂交区发生了rDNA渐渗。之后，欧洲小翅雏蝗Chorthippusparallelus2个邻接亚种在比利牛斯杂交区的核DNA渐渗(Vazquezetal., 1994)及其比利牛斯山东部杂交区的基因渐渗(Shukeretal., 2005)，澳洲Morabine蝗虫(Kawakamietal., 2011)和亚洲稻蝗Oxyahylaintricate(Lietal., 2011)的基因渐渗，瑞士1种小车蝗Oedaleusdecorus的线粒体基因的历史性渐渗(Kindleretal., 2012)以及新西兰蝗虫1个濒危属Sigaus2个种间的基因渐渗(Dowleetal., 2014)均有较为深入的研究，至2015年Pedraza-Lara等人对墨西哥具有重要经济价值的Sphenarium属蝗虫线粒体基因渐渗的情况进行了检测。在这些报道中，Kindler等(2012)和Dowle等(2014)分别使用了11个和3个微卫星位点标记数据进行分析，结果显示对渐渗现象有较好的分辨率。

竹蝗是我国竹林的重要害虫，其中著名的有黄脊竹蝗Ceracriskiangsu和青脊竹蝗C.nigricornis(图1)(陈永林，2007)。黄脊竹蝗和青脊竹蝗在我国主要分布于云南、广西、海南、广东、福建、浙江、安徽、河南、重庆、贵州、四川、陕西、甘肃和台湾等14个省，同时前者在江苏省尚有分布(郑哲民和夏凯龄，1998)。多方面证据已表明，黄脊竹蝗与青脊竹蝗并非姊妹种(易传辉，2003；许升全等，2005；和秋菊等，2011；Gaoetal., 2018)。不过，有意思的是，本课题组多年来野外采集发现，这两种竹蝗在我国南方的大片竹林中常常是同时存在，呈现出同域分布状态。

本文采用微卫星遗传标记法检测安徽金寨县青脊竹蝗和黄脊竹蝗之间是否发生基因渐渗事件，明确这两种竹蝗基因的流动方向，为今后竹蝗属物种形成和遗传分化的研究提供有效科学依据。

图1 青脊竹蝗(a)和黄脊竹蝗(b)Fig.1 Ceracris nigricornis (a) and Ceracris kiangsu (b)

1 材料与方法

1.1 实验材料

本实验所用的标本是2014年采集于安徽省金寨县张冲乡(E 116.0195, N 31.4335)(图2)的青脊竹蝗和黄脊竹蝗各16头个体(雌雄各8头)，共32头个体。所有标本采集后，送回实验室鉴定物种完毕后，使用液氮快速冷冻处理，并保存于-80℃冰箱中直到提取总DNA。

图2 竹蝗采样地点Fig.2 Sampling location for bamboo locusts

1.2 总DNA提取、引物选择、PCR产物扩增和分型

利用本课题组改良的酚-氯仿法提取研究样品的总DNA (Xuanetal., 2009)，并将DNA稀释为30～50 ng/μL的终浓度用于后续的PCR反应。从上述32头标本的后足股节肌肉中提取总DNA。本实验从Xuan等(2009)报道的12对引物中筛选出5对引物，引物详细信息参见表1。微卫星引物由上海捷瑞生物工程有限公司合成，荧光基团(TAMRA，FAM，HEX)在PCR加样操作过程中添加。

表1 本研究中使用的微卫星引物

微卫星序列的PCR反应在BIO-RAD S1000 PCR仪上进行，总反应体系为15 μL，包括包括2.5 μL 10 × PCR buffer，2 μL dNTP mixture (2.5 mM)，0.2 μL Taq polymerase (5 U/μL)，1 μL正反引物(5 μM/L)，1 μL DNA模板，加水至15 μL。扩增反应程序为：95℃预变性3 min；95℃变性30 s，退火30 s(不同的引物退火温度不同，根据引物适当调整)，72℃延伸30 s，共32个循环；最后72℃延伸10 min。

PCR扩增产物首先用琼脂糖凝胶电泳检测是否有目标产物(100～400 bp)，符合目标大小的产物可在ABI 3500 XL遗传自动分析仪上进行基因分型，片段大小参照GS 600 Tamara荧光分子量标准(ABI Inc.，USA)。等位基因分型结果使用GeneMarker V2.2.0(Soft Genetics LLC, USA)软件读取。每个结果均读取2次以上以防止人工误差，并以双栏模式(Two column format)将扩增产物片段长度作为等位基因输入Excel电子表格。使用Excel Microsatellite Tool kit程序(Park, 2001)和Convert version 1.31(Glaubitz, 2004)对结果进行整理和统计分析，转化为其他各种遗传学分析软件所需要的格式。

1.3 数据分析

使用CERVUS 2.0(Marshall, 2001)和POPGENE version 1.32(Yehetal., 1999)软件计算以下遗传多样性参数：平均等位基因数目(Na)，有效等位基因数目(Ne)，观测杂合度(Ho)和期望杂合度(He)，多态信息含量(PIC)(Botsteinetal., 1980)。

用GENEPOP version 3.4(Raymond and Rousset, 2004)通过马尔科夫链方法(Markov chain method)无偏P值估计，设置1000 dememorization number，500 batches，and 1000 iterations per batch，计算哈迪-温伯格平衡(Hardy-Weinberg Equilibrium，HWE)偏离测试以及位点间的连锁不平衡检验(Linkage Disequilibrium，LD)。用Bonferroni correction进行微卫星位点偏离HWE显著性水平的检验(Rice, 1989)。

对于遗传结构分析，本研究使用STRUCTURE 2.3.4软件(Pritchardetal., 2000)，基于贝叶斯方法(Bayesian)，使用混合祖先模型(Admixture model)和相关等位基因频率模型(Allele frequencies correlated model)对青脊竹蝗和黄脊竹蝗同域分布的地理种群32个体进行遗传聚类数K的估算，将所有个体分别纳入相对应的分组之中，并计算各推断聚类的成员系数比例(proportion of membership coefficient)。设定聚类数K=1～10，100万次MCMC重复运算，丢弃前10万次不稳定的结果，分别运行10次求平均值。根据Evanno等(2005)所提出的方法，通过计算连续K值之间Ln P (D)变化率的大小(ΔK)来做折线图选择最佳K值，即得群体遗传结构的分组数(Evannoetal., 2005)。在所有个体的分配中，依据成员最高比例(highest percentage of membership)将所有个体分配到相应的分组。

2 结果与分析

2.1 遗传多样性分析

本文采用5对微卫星引物对安徽金寨县的黄脊竹蝗和青脊竹蝗种群共32个样品进行了分析，共检测到60个等位基因，每个位点拥有等位基因10～13个，平均为12个，其中WJ625位点检测到等位基因数最少；各位点的期望杂合度(He)范围为0.711～0.885，平均预期杂合度为0.833；观测杂合度(Ho)均为1。每个位点都具有较高的多态性(WJ625最低为0.661)，平均多态信息含量(PIC)为0.799，说明这些位点都存在着丰富的遗传信息。在哈迪-温伯格平衡检测中，2个种群在所有位点上均显著偏离哈迪-温伯格平衡(以下简称为哈温平衡)(P<0.01)(表2)。在整体水平上，所有微卫星位点对之间为显著的连锁不平衡(LD)。

2.2 哈迪温伯格平衡检测

在本研究使用的5个微卫星位点，用Bonferroni校正后，根据第一种假设(H1)：杂合子缺乏，黄脊竹蝗3个位点(WJ623,WJ625和WJ630)显示出显著偏离哈温平衡(表3)。

表2 本研究中微卫星位点的遗传多样性

注：**，P<0.01。

表3 哈迪温伯格平衡检测(H1=Heterozygote deficiency)

Table 3 Results of Hardy-Weinberg Equilibrium test

A.黄脊竹蝗种群Ceracriskiangsupopulation

LocusP-valSEW&CR&HWJ6220.18360.02610.00510.0525WJ6230.01150.00130.46110.3778WJ6250.00000.00000.29290.2090WJ6260.24120.03970.07140.0166WJ6300.00000.00000.65750.5213

B.青脊竹蝗种群Ceracrisnigricornispopulation

LocusP-valSEW&CR&HWJ6221.0000-0.0690-0.0690-0.0356WJ6231.00000.0000-0.1060-0.0462WJ6250.99330.0051-0.2245-0.0325WJ6260.99800.0016-0.2162-0.0413WJ6301.00000.0000-0.3793-0.0815

注：P-val，检验的P值；SE，标准误；W&C，Weir & Cockerham’s (1984)估计；R&H，Robertson & Hill’s (1984)估计。下表同。Note:P-val，Pvalue of the test；SE，standard error；W&C，Weir & Cockerham’s (1984) estimate；R&H，Robertson & Hill’s (1984) estimate.Same below.

根据第二种假设(H2)：杂合子过量，青脊竹蝗部分位点(WJ630)显示出显著偏离哈温平衡(表4)。

通过比较上述P值，本研究发现黄脊竹蝗种群有3个微卫星位点显著偏僻哈温平衡，其中WJ625和WJ630极显著(表3)，而青脊竹蝗种群也有3个位点显著偏离哈温平衡，但只有1个位点(WJ630)极显著，说明2种竹蝗的群体结构均不稳定，且黄脊竹蝗更甚。

表4 哈迪温伯格平衡检测(H2=Heterozygote excess)

Table 4 Results of Hardy-Weinberg Equilibrium test

A.黄脊竹蝗种群Ceracriskiangsupopulation

LocusP-valSEW&CR&HWJ6220.88300.0151 0.00510.0525WJ623 0.98910.00070.46110.3778WJ6251.00000.00000.29290.2090WJ6260.73810.02980.07140.0166WJ6301.00000.00000.65750.5213

B.青脊竹蝗种群Ceracrisnigricornispopulation

LocusP-valSEW&CR&HWJ6220.34200.0538-0.0690-0.0356WJ6230.22360.0386-0.1060-0.0462WJ6250.02650.0187-0.2245-0.0325WJ6260.02340.0138-0.2162-0.0413WJ6300.00000.0000-0.3793-0.0815

2. 3 遗传结构与关系

当K=2时，L (K)和ΔK值均为最大值，表明K=2为最佳分组(图3)。

STRUCTURE的聚类分析表明(图4)，当K=2时，黄脊竹蝗种群个体聚为一组(红色)，青脊竹蝗种群聚为另一组。青脊竹蝗种群的部分样品包含黄脊竹蝗的遗传物质，而黄脊竹蝗种群的遗传背景比较一致。

图3 基于△K 最大时的最合适K值Fig.3 Optimum K value when△K is the largest

图4 基于贝叶斯分析的黄脊竹蝗种群和青脊竹蝗种群聚类结果(K=2)Fig.4 Bayesian inference of the number of clusters (K=2) for the Ceracris kiangsu and Ceracris nigricornis populations analyzed by STRUCTURE注：编码1-16为1号种群黄脊竹蝗个体，编码17-32为2号种群青脊竹蝗个体。Note: Codes 1-16 are individuals of Ceracris kiangsu and codes 17-32 are Ceracris nigricornis.

3 讨论与结论

在群体遗传学研究中，物种变异程度和群体适应性随着观测杂合度增大而增大。从本研究结果来看，黄脊竹蝗和青脊竹蝗的物种变异程度和群体适应性均很大。不过，各微卫星位点的观测杂合度均为1，处于最大值，可能与观测杂合度受样本量大小的变化而出现明显变异有很大的关系(李鸥等，2009；杨君等，2018)，因此该结果的参考价值并不大。而期望杂合度(He)则是判断遗传多样性的重要指标，它用于辨别群体中杂合度比率偏离的情况。本研究中各位点的期望杂合度(He)均大于0.70，表明这些种群属于高度多态，遗传多样性水平较高。

在假设杂合子缺乏情形下，黄脊竹蝗种群部分位点显著地偏离哈迪温伯格平衡。这可能是由于黄脊竹蝗的种群结构更加不稳定，更容易受到突变、选择、杂交和遗传漂变等因素的影响，及受近交和无效等位基因等因素的影响。造成偏离哈温平衡的原因有很多，一方面可能是存在无效等位基因(Pembertonetal., 1995)，而被误认为纯合子，因而导致偏离哈温平衡；另一方面也可能是由于本研究的样本数量不足够大(李鸥等，2009；牛东红等，2010)。同时，本研究在局域水平展示了维持或者强化遗传结构的多态性证据。本研究检测到了同域分布物种间的基因流。贝叶斯聚类分析结果显示，基因流是从黄脊竹蝗种群到青脊竹蝗种群的。也可能是因为黄脊竹蝗种群多分布于竹林郁闭度较高的区域，青脊竹蝗多分布于郁闭度较低的竹林外缘区域，而黄脊竹蝗雄性要比雌性更加具有扩散能力(Baileyetal., 2007；Ortegoetal., 2011)，所以黄脊竹蝗雄性个体由竹林深处向外缘扩散，从而可能与青脊竹蝗雌性个体发生交配。

最后，本研究聚焦于地理种群水平上的遗传差异基因渐渗事件，检测了种群内的差异。这显示同域分布的黄脊竹蝗和青脊竹蝗种群间发生了基因渐渗，但本研究目前不能确定其是不是来自遗传漂变和瓶颈期的信号。竹蝗属于植食性昆虫，影响其种群分化的原因有很多，如迁移能力和地理隔离，寄主植物的选择以及生境片段化等。一般情况下，基因流受到地理距离阻碍，导致遗传分化。一些迁飞能力较强的物种分布区较广如东亚飞蝗Locustamigratoriamanilensis，增加了基因交流，从而削弱了种群间分化，一般有比较弱的地理格局。此外，频繁、强烈的基因流能够影响遗传变异的均质化、选择、突变以及遗传漂变等(Slatkin, 1987)。然而，青脊竹蝗和黄脊竹蝗种群间遗传和种群动力的相互作用仍需要继续深入探索。对于竹蝗这样的非模式生物物种，目前尚无其全基因组数据可利用。不过，新近开发的简化基因组技术(reduced representation genome sequencing，RRGS)(魏书军等，2016)可以通过对某物种大量个体的基因组进行测序而满足此目的，并具有较高的可靠度，可用于进一步研究同域分布物种间的基因渐渗现象。

致谢：南京师范大学生命科学学院研究生李冉同学跟随本文通讯作者前往安徽金寨县采集竹蝗标本，泉州师范学院李明杰同学协助绘制插图，特此致谢！