程海云, 段家充, 張 超, 潘 昭
(河北大學(xué)生命科學(xué)學(xué)院, 生命科學(xué)與綠色發(fā)展研究院, 河北省動物系統(tǒng)學(xué)與應(yīng)用重點(diǎn)實(shí)驗(yàn)室, 河北保定 071002)
芫菁科(Meloidae)隸屬于鞘翅目(Coleoptera)擬步甲總科(Tenebrionoidea),世界已知3亞科16族133屬近3 000種,廣泛分布于除新西蘭、南極和波利尼西亞群島之外的世界各地;中國記錄2亞科8族27屬近200種(亞種),各省市自治區(qū)均有記載(Bolognaetal., 2010; Pan and Ren, 2020)。該科昆蟲是一類重要的資源昆蟲,因其復(fù)變態(tài)等獨(dú)特的生物學(xué)特性而被學(xué)者廣泛關(guān)注,尤其是斑蝥素的藥用價值更是近年芫菁應(yīng)用研究的熱點(diǎn)之一(楊玉霞和任國棟, 2005; Bolognaetal., 2010)。然而,芫菁科部分物種的外部形態(tài)特征十分近似,在形態(tài)分類研究中經(jīng)常出現(xiàn)難以區(qū)分、錯誤鑒定的情況,如北方常見物種圓點(diǎn)斑芫菁Mylabrisaulica和蘋斑芫菁Mylabriscalida(Pan and Ren, 2020)、橫紋溝芫菁Hycleussolonicus(Panetal., 2017)等。分子系統(tǒng)學(xué)的發(fā)展,特別是DNA條形碼的提出和興起,為解決以上問題提供了可行的方案。Liu等(2016)基于線粒體COI基因片段序列,對西北豆芫菁Epicautasibirica、疑豆芫菁Epicautadubia和中華豆芫菁Epicautachinensis分類研究的有效性進(jìn)行了探討,提出后兩者是西北豆芫菁的次異名,然而該研究并未使用分子物種界定方法,僅是基于分支單系性得出的結(jié)論。
為更加全面地評估分子物種界定方法在芫菁科分類研究中的可行性,本研究以中國北方常見芫菁科昆蟲為研究材料,選取線粒體COI基因片段和核CAD基因片段序列,分別采用自動條形碼間隔探索(automatic barcode gap discovery, ABGD)、廣義混合Yule溯祖模型(generalized mixed Yule coalescent,GMYC)、貝葉斯泊松樹進(jìn)程(Bayesian Poisson tree processes, bPTP)和貝葉斯系統(tǒng)發(fā)育和系統(tǒng)地理分析(Bayesian phylogenetics and phylogeography, BPP)方法進(jìn)行分子物種界定。在對比形態(tài)鑒定結(jié)果的基礎(chǔ)上,探討以上序列和分析方法界定的準(zhǔn)確性和適用性,以期為芫菁科的整合分類研究提供數(shù)據(jù)支持。
本研究以中國北方常見芫菁科昆蟲為研究對象,選取58個樣本,經(jīng)形態(tài)鑒定,共涉及2亞科4族6屬18種;選擇擬步甲科(Tenebrionidae)黃粉蟲Tenebriomolitor為外群。實(shí)驗(yàn)樣本大多來自河北大學(xué)博物館館藏,少部分樣本序列下載自GenBank,詳見表1。
選取樣本均使用Nikon SMZ1500立式顯微鏡進(jìn)行形態(tài)學(xué)鑒定。形態(tài)鑒定參考中國芫菁分類相關(guān)成果,如:Pan等(2010, 2017),潘昭等(2010a, 2010b, 2011, 2013), Li等(2020)以及Pan和Bologna(2021)。
從酒精泡制標(biāo)本中取后足肌肉組織,使用昆蟲基因組DNA提取試劑盒(倍沃醫(yī)學(xué)科技)提取總DNA。通過PCR擴(kuò)增COI和CAD基因片段,COI基因擴(kuò)增引物為TW-N1284(5′-AAACTAACARCCTTC AAAGYTGT-3′)和COI-R2329(5′-ACHGTAAAYATA TGATGGGCTCA-3′),反應(yīng)程序:94℃預(yù)變性3 min;94℃變性1 min, 48℃延伸1.5 min, 40個循環(huán); 72℃復(fù)性1.5 min;72℃最終延伸3 min。CAD基因擴(kuò)增引物為CD439F(5′-TTCAGTGTACARTTYCAYC CHGARCAYAC-3′)和CD688R(5′-TGTATACCTAG AGGATCDACRTTYTCCATRTTRCA-3′),反應(yīng)程序: 94℃預(yù)變性3 min; 94℃變性30 s, 55℃延伸30 s, 50個循環(huán); 72℃復(fù)性1 min; 72℃最終延伸3 min(Wild and Maddison, 2008)。PCR擴(kuò)增采用25 μL體系:上下游引物(10 μmol/L)各1 μL,去離子水8.5 μL, PCR Mix 12.5 μL,DNA模板2 μL。 PCR產(chǎn)物經(jīng)1%瓊脂糖凝膠電泳檢測,送通用生物系統(tǒng)(安徽)有限公司進(jìn)行雙向測序。
所得序列使用DNASTAR SeqMan v7.1.0進(jìn)行檢查和拼接,在NCBI中進(jìn)行BLAST分析,確定所得序列是否為目的序列, 使用MEGA X (Kumaretal., 2018)中的Clustal W進(jìn)行序列的比對。所有新得到的序列上傳至GenBank(登錄號見表1)。
使用MEGA X分別對COI和CAD序列數(shù)據(jù)集進(jìn)行序列堿基含量分析及種內(nèi)和種間遺傳距離分析(Kumaretal., 2018)。使用PhyloSuite v1.2.2 (Zhangetal., 2020)將COI和CAD序列數(shù)據(jù)集串聯(lián),構(gòu)建3個序列數(shù)據(jù)集,即COI,CAD和COI+CAD串聯(lián)序列數(shù)據(jù)集。使用PhyloSuite內(nèi)嵌的MrBayes 3.2(Ronquistetal., 2012),使用貝葉斯推斷法(Bayesian inference, BI),分別基于3個數(shù)據(jù)集重建系統(tǒng)發(fā)育樹,運(yùn)行2 000 000世代,采樣頻率為每100代一次,確保分裂頻率分支頻率的標(biāo)準(zhǔn)偏差小于0.01(Ronquist and Huelsenbeck, 2003),以獲得高支持率的合意樹,生成合意樹之前舍棄25%的老化樣本,最佳堿基替代模型使用ModelFinder (Kalyaanamoorthyetal., 2017)進(jìn)行評估。 最后,使用Evolview v3 (Subramanianetal., 2019)對系統(tǒng)樹進(jìn)行可視化編輯。
基于3個數(shù)據(jù)集,分別使用自動條形碼間隔探索(ABGD),廣義混合Yule溯祖模型(GMYC),貝葉斯泊松樹進(jìn)程(bPTP)和貝葉斯系統(tǒng)發(fā)育和系統(tǒng)地理分析(BPP)方法進(jìn)行物種界定分析,具體參數(shù)設(shè)置如下:
ABGD:使用ABGD在線分析(https:∥bioinfo.mnhn.fr/abi/public/abgd/abgdweb.html)(金倩等, 2017),設(shè)置參數(shù)X(relative gap width)為0.5,其他參數(shù)為默認(rèn)值,將得到的序列劃分情況與形態(tài)鑒定的結(jié)果進(jìn)行比較分析。
GMYC(金倩等, 2017):基于BEAST v1.10.4 (Drummond and Rambaut, 2007)構(gòu)建超度量二歧分支樹(有根),選擇嚴(yán)格分子鐘。MCMC 鏈為10 000萬代,每1 000代參數(shù)取樣一次,Burn-In默認(rèn)為鏈長10%,使用Tracer v1.7.2(Rambautetal., 2018)評估系統(tǒng)樹的收斂情況。將含有樹文件的*.trees.txt 文件導(dǎo)入到Tree Annotator v1.10.4生成所需系統(tǒng)樹。Burn In (as trees)根據(jù) Tracer v1.7.2分析結(jié)果(Burn-In 鏈長/采樣頻率)設(shè)置為1 000,最后生成所需最大分支可信度樹(maximum clade credibility tree),導(dǎo)入GMYC進(jìn)行分析。
bPTP(Zhangetal., 2013):本實(shí)驗(yàn)使用PTP的更新版本bPTP進(jìn)行界定分析,將GMYC分析中得到的Newick系統(tǒng)樹上傳到PTP在線網(wǎng)站(https:∥species.h-its.org/ptp/)進(jìn)行分析,設(shè)置為有根樹,輸入外群名稱,馬爾科夫鏈代數(shù)為500 000,其他參數(shù)為默認(rèn)值。
BPP(Yang, 2015):基于串聯(lián)基因數(shù)據(jù)集構(gòu)建的貝葉斯系統(tǒng)樹,在BPP v4.2.9中利用A11模板輸入指導(dǎo)樹(基于COI序列構(gòu)建的BI樹),運(yùn)行分析。
COI基因片段序列經(jīng)拼接、比對、剪切處理后的長度為561 bp。使用MEGA X軟件中的Kimura雙參數(shù)模型(Kimura 2-parameter model)進(jìn)行序列堿基含量分析及種內(nèi)、種間遺傳距離分析(表2和3),結(jié)果顯示:序列的T, C, A和G平均含量分別為33.5%, 22.3%, 26.8%和17.5%,平均A+T含量顯著高于G+C含量,變異多發(fā)生在堿基第3位;Hycleussolonicus,Mylabrisaulica,Epicautamegalocephala,Mylabrissibirica,Lyttacaraganae,Meloeproscarabaeus,Meloeviolaceus和Epicautasibirica的種內(nèi)平均遺傳距離分別為0.92%, 0.71%, 0.83%, 0.30%, 1.66%, 1.34%, 0.11%和0.24%,其余10個物種樣本量不足,無法進(jìn)行種內(nèi)遺傳距離分析;18種芫菁的種間遺傳距離在9.49%~22.10%之間,平均值為18.36%。
表2 基于COI基因片段序列計算得出芫菁科種內(nèi)平均遺傳距離和種內(nèi)遺傳距離范圍Table 2 Average intraspecific genetic distance and range of intraspecific genetic distance of Meloidaecalculated based on COI fragment sequences
CAD基因片段序列處理后長度為730 bp。使用MEGA X 軟件中的Kimura 雙參數(shù)模型進(jìn)行序列堿基含量分析及種內(nèi)種間遺傳距離分析(表4和5),結(jié)果顯示:T, C, A和G平均含量分別為33.5%, 22.3%, 26.8%和17.5%,平均A+T含量顯著高于G+C含量,沒有堿基變異偏向性,密碼子的1, 2和3位均有發(fā)生變異;Hycleussolonicus,Mylabrisaulica,Epicautamegalocephala,Epicautasibirica,Lyttacaraganae,Meloeproscarabaeus和Meloeviolaceus的種內(nèi)平均遺傳距離為1.71%, 0.65%, 0.91%, 0.73%, 0.46%, 0.27%和0.06%,Mylabrissibirica沒有變異,其余10個物種樣品量不足,無法進(jìn)行種內(nèi)遺傳距離分析; 18種芫菁的種間遺傳距離是處于0.43%~22.33%之間,平均值為13.85%。
表4 基于CAD基因片段序列計算得出芫菁科種內(nèi)平均遺傳距離和種內(nèi)遺傳距離范圍Table 4 Average intraspecific genetic distance and range of intraspecific genetic distance ofMeloidae calculated based on CAD fragment sequences
利用PhyloSuite內(nèi)嵌的ModelFinder(Kalyaanamoorthyetal., 2017)進(jìn)行最佳堿基替代模型分析, 實(shí)現(xiàn)Akaike’sInformation Criterion (AIC)選擇最適模型,3個數(shù)據(jù)集的最適模型均為GTR+F+I+G4。
對58頭芫菁樣本進(jìn)行BI分析,基于COI,CAD和COI+CAD串聯(lián)序列數(shù)據(jù)集分別建樹,顯示出高度相似的拓?fù)浣Y(jié)構(gòu)(圖1),聯(lián)合數(shù)據(jù)集分支支持度高于單基因建樹,高度支持這18個種的單系性。除Hycleussolonicus,Meloevariegatus,Meloebrevicollis和Epicautaobscurocephala這幾個分支的后驗(yàn)概率(posterior probability, PP)<1,其余種劃分PP均等于1。
本實(shí)驗(yàn)對58個芫菁樣本(不含外群)的分子物種界定結(jié)果如圖1所示,具體如下:
ABGD劃分結(jié)果包括初始劃分(initial partition)和遞歸劃分(recursive partition)兩種情況。 基于COI基因片段和COI+CAD串聯(lián)基因片段序列的劃分結(jié)果顯示,初始劃分較為穩(wěn)定,遺傳距離在0.001~0.1之間,被劃分為18個分子操作分類單元(molecular operational taxonomic unit, MOTU),支持18個形態(tài)種,且與BI樹拓?fù)浣Y(jié)構(gòu)一致。遞歸劃分先驗(yàn)值浮動較大,不夠穩(wěn)定?;贑AD基因片段序列的ABGD劃分結(jié)果顯示,屬級分類單元的劃分和形態(tài)劃分吻合,而物種階元上則將Meloeproscarabaeus和Meloeviolaceus劃分為一個MOTU,物種劃分紊亂,如圖2(A, B, C)。
基于3個數(shù)據(jù)集的GMYC分析所劃分的MOTUs結(jié)果皆與形態(tài)種相符。單閾值分析,基于COI基因片段序列劃分的似然實(shí)體數(shù)量為19,置信區(qū)間為19~22,基于CAD基因片段序列劃分的似然實(shí)體數(shù)量為19,置信區(qū)間為1~49,基于COI+CAD串聯(lián)序列劃分的似然實(shí)體數(shù)量為19,置信區(qū)間為15~21,如圖3-5所示。
bPTP分析結(jié)果表明,基于COI基因序列劃分得到21個MOTUs,將Mylabrisaulica,Meloeproscarabaeus和Lyttacaraganae各自分為了2個MOTUs;基于CAD基因序列劃分為17個MOTUs,將Meloeproscarabaeus和Meloeviolaceus歸為一個MOTU;基于COI+CAD基因串聯(lián)序列得到18個MOTUs,和形態(tài)學(xué)劃分結(jié)果一致,如表6-8所示。
表6 基于COI序列通過bPTP方法得到的芫菁科物種劃分Table 6 Species division of Meloidae based on COI sequence by bPTP method
表7 基于CAD序列通過bPTP方法得到的芫菁科物種劃分Table 7 Species division of Meloidae based on CAD sequence by bPTP method
表8 基于COI+CAD串聯(lián)序列通過bPTP方法得到的芫菁科物種劃分Table 8 Species division of Meloidae based on concatenated COI+CAD sequence by bPTP method
BPP是基于給定指導(dǎo)樹進(jìn)行分析,得到18個MOTUs的后驗(yàn)概率(PP)是0.98925,遠(yuǎn)大于其他劃分方式的PP,即支持將58個樣本劃分為18個獨(dú)立物種,支持指導(dǎo)樹給定的18個物種,與形態(tài)學(xué)劃分一致。
Liu等 (2016)首次應(yīng)用分子系統(tǒng)學(xué)方法對芫菁科昆蟲的物種劃分問題進(jìn)行了探討并提出了2個新異名,然而此工作依據(jù)的是外部形態(tài)、遺傳距離和系統(tǒng)樹拓?fù)浣Y(jié)構(gòu)的綜合結(jié)果,并未使用其他分子物種界定方法進(jìn)行分析。
本研究對中國北方常見的芫菁科6屬18個形態(tài)種進(jìn)行分子物種界定,總體而言所得結(jié)果基本與形態(tài)劃分一致,僅基于CAD序列的ABGD分析結(jié)果和基于單基因的bPTP分析結(jié)果有所差異(圖1)。該結(jié)果一方面驗(yàn)證了形態(tài)種劃分結(jié)果,例如目前使用觸角末節(jié)形狀對具有相似鞘翅斑紋的圓點(diǎn)斑芫菁和蘋斑芫菁進(jìn)行區(qū)分(Pan and Ren, 2020),該形態(tài)界定結(jié)果得到了分子物種界定結(jié)果的支持;另一方面也顯現(xiàn)出對芫菁科昆蟲相對可靠的可用于分子物種界定的DNA序列和分析方法,即多基因片段串聯(lián)序列優(yōu)于COI序列優(yōu)于CAD序列,ABGD方法最簡便、準(zhǔn)確性相對較高且穩(wěn)定,在一定程度上可減少人為界定的誤差。在進(jìn)行ABGD和bPTP分析時,發(fā)現(xiàn)基于CAD的劃分對物種劃分不清晰,因此推測是由于CAD基因進(jìn)化速率較慢,種間差異小,所以導(dǎo)致種的劃分紊亂。
此外,參考分子物種界定結(jié)果,基于COI基因種內(nèi)和種間的遺傳距離分析結(jié)果表明,在芫菁科中,當(dāng)遺傳距離在1.66%以下可認(rèn)為是同一物種的種內(nèi)變異,而大于9.49%則可認(rèn)為是不同的物種,而在1.66%~9.49%之間則需要綜合更多證據(jù)謹(jǐn)慎處理(表2和3)?;贑AD基因種內(nèi)和種間的遺傳距離和物種界定綜合分析結(jié)果表明,在芫菁科中,當(dāng)遺傳距離在1.71%以下可認(rèn)為是同一物種的種內(nèi)變異,而大于1.71%則有成為不同物種的可能,在1.71%~2.77%之間則需要綜合更多證據(jù)謹(jǐn)慎處理(表4和5)。
綜上,基于COI和CAD的獨(dú)立和串聯(lián)數(shù)據(jù)集和ABGD, GMYC, PTP和BPP 4種分析方法均適用于芫菁科昆蟲的分子物種界定,但有一定的誤差,需多基因、多方法綜合考慮。此外,形態(tài)和分子物種界定方法的結(jié)合能夠有效提高芫菁科昆蟲物種劃分的準(zhǔn)確性。