林興雨, 宋 南

(河南農(nóng)業(yè)大學(xué)植物保護(hù)學(xué)院,河南 鄭州 450002)

蝽總科(Pentatomoidea)是半翅目(Hemiptera)異翅亞目(Heteroptera)蝽次目(Pentatomomorpha)中物種數(shù)量最多的一個類群。目前,全世界有18個科1 410個屬超過8 000個已知種被描述[1]。蝽總科的大多數(shù)昆蟲以糧食作物或野生植物為食,有些是農(nóng)業(yè)和林業(yè)上重要的經(jīng)濟(jì)害蟲[2];但有些屬于天敵昆蟲,在農(nóng)業(yè)生產(chǎn)上被用于害蟲的生物防治[3]。

斯氏珀蝽 (PlautiastaliScott, 1874)屬于蝽總科蝽科(Pentatomidae),是農(nóng)業(yè)生產(chǎn)上危害較嚴(yán)重的一種昆蟲。若蟲有5個齡期,整個若蟲期約1個月,成蟲壽命較長,一般為半個月至1個月,有一定的趨光性,一年發(fā)生2代[4-5]。目前,有關(guān)斯氏珀蝽的研究較少,主要集中在生物學(xué)和昆蟲-微生物腸道共生關(guān)系方面[5-6]。青革土蝽(MacroscytusjaponensisScott, 1874)屬于蝽總科土蝽科(Cydnidae)。有關(guān)該蟲的研究很少,僅Zhu et al[7]曾對青革土蝽進(jìn)行了簡要描述。由上可知,有關(guān)斯氏珀蝽和青革土蝽的分子結(jié)構(gòu)和系統(tǒng)發(fā)育關(guān)系尚不明確。

昆蟲線粒體基因組呈環(huán)狀雙鏈結(jié)構(gòu),由37個基因組成,即13個蛋白質(zhì)編碼基因(protein-coding genes, PCGs)、2個核糖體RNA基因(ribosomal RNA genes, rRNAs)、22個轉(zhuǎn)運(yùn)RNA基因(transfer RNA genes, tRNAs),還有非編碼控制區(qū)[8-10]。通過研究昆蟲線粒體基因組可獲得其分子進(jìn)化的關(guān)鍵信息,如核苷酸組成、密碼子的使用頻率和進(jìn)化率等[11]。高通量測序技術(shù)的快速發(fā)展,使得昆蟲線粒體基因組較易通過組裝而得到,且因其具有保守的基因排列和母系遺傳等優(yōu)點(diǎn)被廣泛用于昆蟲的分子進(jìn)化、生物地理學(xué)和系統(tǒng)發(fā)育等研究[12-18]。

本研究利用高通量測序技術(shù)獲得斯氏珀蝽和青革土蝽完整的線粒體基因組數(shù)據(jù),并在此基礎(chǔ)上構(gòu)建蝽總科的系統(tǒng)發(fā)育樹,為更好地理解蝽總科下科級和種級水平的系統(tǒng)發(fā)育關(guān)系提供依據(jù)。

1 材料與方法

1.1 標(biāo)本采集與DNA提取

斯氏珀蝽和青革土蝽于2019年7月采集于河南省鄭州市人民公園(34°76′N、113°66′E)。將采集到的斯氏珀蝽和青革土蝽浸泡在無水乙醇中,并置于-20 ℃超低溫冰箱中備用。使用天根生化科技有限公司的動物基因組DNA提取試劑盒(TIANamp Genomic DNA Kit, China)提取斯氏珀蝽和青革土蝽樣本胸部肌肉的總DNA,按照說明書操作。使用NanoDrop 2000分光光度計和1.5%(質(zhì)量分?jǐn)?shù))瓊脂糖凝膠電泳分別檢測總DNA的質(zhì)量和濃度。

1.2 高通量測序與線粒體基因組組裝

將質(zhì)量合格且濃度大于20 ng·mL-1的DNA送至上海派森諾生物科技有限公司,利用全基因組鳥槍法(whole genome shotgun, WGS)構(gòu)建350 bp的小片段文庫,并使用Illumina HiSeq 2500測序平臺對構(gòu)建的文庫進(jìn)行雙末端(2×150 bp PE)測序。將下機(jī)的原始數(shù)據(jù)去除接頭,進(jìn)行質(zhì)量控制,獲得clean reads。利用GetOrganelle v1.7.5.2軟件[19]對雙末端測序獲得的clean reads進(jìn)行組裝和拼接,得到高質(zhì)量的線粒體基因組數(shù)據(jù)。

1.3 線粒體基因組注釋和分析

先在MITOS網(wǎng)站[20]上對斯氏珀蝽和青革土蝽的線粒體基因組數(shù)據(jù)進(jìn)行初步的功能注釋和tRNAs二級結(jié)構(gòu)預(yù)測,具體參數(shù)設(shè)置:Genetic Code為05-inverterbrate, Reference為RefSeq 63 Metazoa,其余均為默認(rèn)參數(shù)。然后,將斯氏珀蝽和青革土蝽與近緣種的線粒體基因組序列進(jìn)行比對驗(yàn)證,手動檢查、矯正13個PCGs、22個tRNAs、2個rRNAs的基因邊界。將斯氏珀蝽和青革土蝽的線粒體基因組數(shù)據(jù)注釋后上傳至GenBank(GenBank序列號分別為OP919327和OP930890)。利用MEGA 11.0軟件[21]計算和分析斯氏珀蝽和青革土蝽線粒體基因組中所有PCGs的核苷酸組成和密碼子使用頻率的情況。根據(jù)AT-skew=(A-T)/(A+T)和GC-skew=(G-C)/(G+C)分別計算它們的AT偏斜和GC偏斜[22]。

1.4 多序列比對

先使用MAFFT軟件[23]中的E-INS-I策略對所有PCGs進(jìn)行序列比對。然后,利用trimAl v1.4軟件[24]對質(zhì)量較差的序列片段進(jìn)行修剪、刪除和對齊。最后,使用FASconCAT-G v1.04腳本[25]串聯(lián)數(shù)據(jù)獲得數(shù)據(jù)矩陣。用于系統(tǒng)發(fā)育分析的數(shù)據(jù)矩陣包括:(1)13個PCGs的核苷酸序列構(gòu)建的數(shù)據(jù)矩陣;(2)13個PCGs的第1、2密碼子的核苷酸序列構(gòu)建的數(shù)據(jù)矩陣,因?yàn)槿コ齈CGs的第3密碼子可以有效減小異質(zhì)性飽和度對系統(tǒng)發(fā)育的不利影響[18];(3)13個PCGs的氨基酸序列構(gòu)建的數(shù)據(jù)矩陣。

1.5 系統(tǒng)發(fā)育分析

以已經(jīng)公布線粒體基因組數(shù)據(jù)的蝽總科9個科[蝽科、龜蝽科(Plataspidae)、盾蝽科(Scutelleridae)、兜蝽科(Dinidoridae)、荔蝽科(Tessaratomidae)、土蝽科、同蝽科(Acanthosomatidae)、異蝽科(Urostylididae)和Megarididae]的113個物種和本研究的斯氏珀蝽和青革土蝽作為內(nèi)群,2個紅蝽總科(Pyrrhocoroidea)物種(Dysdercusevanescens和Pyrrhocoristibialis)作為外群,分別采用最大似然法(maximum likelihood)和貝葉斯法(Bayesian inference)構(gòu)建系統(tǒng)發(fā)育樹。利用 IQ-TREE v2.0.6軟件[26]對13個PCGs的核苷酸序列構(gòu)建的數(shù)據(jù)矩陣和13個PCGs的第1、2密碼子的核苷酸序列構(gòu)建的數(shù)據(jù)矩陣進(jìn)行最大似然系統(tǒng)發(fā)育分析:由軟件自動篩選的最適模型為GTR+F+I+G4,系統(tǒng)發(fā)育樹的節(jié)點(diǎn)支持率通過自舉檢驗(yàn)置信度進(jìn)行評估(運(yùn)行次數(shù)設(shè)置為10 000次)。使用MrBayes v3.2軟件[27]對13個PCGs的氨基酸序列構(gòu)建的數(shù)據(jù)矩陣進(jìn)行貝葉斯系統(tǒng)發(fā)育分析:執(zhí)行4條蒙特卡羅馬爾柯夫鏈(Markov chain Monte Carlo, MCMC),并運(yùn)行10 000 000代,每運(yùn)行1 000代抽樣一次,每次均丟棄25%的老化樹,最后獲得多數(shù)一致樹。

2 結(jié)果與分析

2.1 線粒體基因組分析

利用OGDRAW v1.3.1軟件[28]繪制斯氏珀蝽和青革土蝽的線粒體基因組結(jié)構(gòu)圖(圖1-2)。斯氏珀蝽和青革土蝽的線粒體基因組全長分別為15 346、14 632 bp。斯氏珀蝽的A、T、G、C含量分別為42.20%、32.60%、10.60%、14.60%;青革土蝽的A、T、G、C含量分別為41.10%、31.60%、11.20%、16.00%。由上可知,斯氏珀蝽和青革土蝽的線粒體基因組均具有較高的AT含量。它們的線粒體基因組均包含37個基因(13個PCGs、22個tRNAs、2個rRNAs)和1個非編碼控制區(qū)(表1-2)。

表1 斯氏珀蝽線粒體基因組注釋1)Table 1 Annotation of P.stali mitochondrial genome

圖1 斯氏珀蝽線粒體基因組結(jié)構(gòu)Fig.1 Structure of P.stali mitochondrial genome

在斯氏珀蝽和青革土蝽的線粒體基因組中,H鏈(heavy strand)和L鏈(light strand)的編碼基因不同;H鏈共編碼23個基因,包含9個PCGs和14個tRNAs;L鏈共編碼14個基因,包含4個PCGs、8個tRNAs以及2個rRNAs。斯氏珀蝽的非編碼控制區(qū)處于rrnS和trnI之間,全長為689 bp,A、T、G、C含量分別為34.54%、34.69%、11.32%、19.45%;青革土蝽的非編碼控制區(qū)在nad4l和trnT之間,長度較短,為68 bp,A、T、G、C含量分別為54.41%、39.71%、1.47%、4.41%。斯氏珀蝽的整個線粒體基因組AT含量為74.80%,AT偏斜和GC偏斜分別為0.128和-0.157(表3);青革土蝽的整個線粒體基因組AT含量為72.70%,AT偏斜和GC偏斜分別為0.130和-0.177(表4)。

2.2 PCGs分析

斯氏珀蝽的13個PCGs中,只有nad2、cox1、nad6和nad1利用ATC或TTG作為起始密碼子,其余9個PCGs以ATA、ATT或ATG作為起始密碼子;cox1具有不完整的終止密碼子T,cox2以TTG作為終止密碼子,其余11個PCGs均具有完整的終止密碼子TAA或TAG(表1)。青革土蝽的13個PCGs中,只有nad1以GTG作為起始密碼子,其余12個PCGs以ATT、ATA或ATG作為起始密碼子;5個PCGs(cox1、cox2、atp6、cox3和nad5)具有不完整的終止密碼子T或TA,其余8個PCGs均具有完整的終止密碼子TAA或TAG(表2)。

表2 青革土蝽線粒體基因組注釋1)Table 2 Annotation of M.japonensis mitochondrial genome

表3 斯氏珀蝽線粒體基因組的核苷酸組成和偏斜1)Table 3 Nucleotide composition and skewness of P.stali mitochondrial genome

表4 青革土蝽線粒體基因組的核苷酸組成和偏斜1)Table 4 Nucleotide composition and skewness of M.japonensis mitochondrial genome

斯氏珀蝽和青革土蝽所有PCGs的氨基酸中,丙氨酸(Ala)占42.21%和41.09%,半胱氨酸(Cys)占14.56%和16.04%,甘氨酸(Gly)占10.60%和11.21%,蘇氨酸(Thr)占32.61%和31.64%,未發(fā)現(xiàn)其他氨基酸。斯氏珀蝽和青革土蝽的線粒體基因組的相對密碼子使用頻率見表5-6。

表5 斯氏珀蝽線粒體基因組的相對密碼子使用頻率(RSCU)1)Table 5 Relative synonymous codon usage (RSCU) of P.stali mitochondrial genome

表6 青革土蝽線粒體基因組的相對密碼子使用頻率(RSCU)1)Table 6 Relative synonymous codon usage (RSCU) of M.japonensis mitochondrial genome

2.3 tRNAs和rRNAs分析

斯氏珀蝽的22個tRNAs的序列長度為63～74 bp。其中:最長的基因是trnK,為74 bp;trnC和trnP最短,均為63 bp(表1)。青革土蝽的22個tRNAs的序列長度為62～72 bp。其中:trnK最長,為72 bp;trnR最短,為62 bp(表2)。斯氏珀蝽和青革土蝽的24個RNA基因中,10個RNA基因位于L鏈,14個RNA基因位于H鏈。這兩種昆蟲的trnS1基因都因缺少二氫尿嘧啶(DHU)臂而無法形成完整的三葉草結(jié)構(gòu),其余tRNAs均能形成完整的三葉草結(jié)構(gòu)(圖3-4)。此外,斯氏珀蝽和青革土蝽的trnS1基因的反密碼子不是常見的GCU,而是UCG。

圖3 斯氏珀蝽線粒體基因組中22個tRNAs的二級結(jié)構(gòu)Fig.3 Secondary structure of 22 tRNAs in P.stali mitochondrial genome

圖4 青革土蝽線粒體基因組中22個tRNAs的二級結(jié)構(gòu)Fig.4 Secondary structure of 22 tRNAs in M.japonensis mitochondrial genome

斯氏珀蝽線粒體基因組中rrnL的全長為1 235 bp(表1),AT含量為77.70%(表3),GC含量為22.30%;rrnS的全長為793 bp(表1),AT含量為76.30%(表3),GC含量為23.70%。青革土蝽線粒體基因組中rrnL的全長為1 271 bp(表2),AT含量為75.80%(表4),GC含量為24.20%;rrnS的全長為802 bp(表2),AT含量為73.60%(表4),GC含量為26.40%。

2.4 系統(tǒng)發(fā)育分析

兩種系統(tǒng)發(fā)育分析方法構(gòu)建的系統(tǒng)發(fā)育樹(圖5-7)一致支持蝽科、龜蝽科、盾蝽科、兜蝽科、荔蝽科、土蝽科、同蝽科、異蝽科和Megarididae均為單系群?；?3個PCGs的核苷酸序列構(gòu)建的系統(tǒng)發(fā)育樹(圖5)顯示:蝽科與龜蝽科為姐妹群,但節(jié)點(diǎn)支持率較低(BS=63);土蝽科與(荔蝽科+兜蝽科)為姐妹群,節(jié)點(diǎn)支持率相對較高(BS=93)?；?3個PCGs的第1、2密碼子的核苷酸序列構(gòu)建的系統(tǒng)發(fā)育樹(圖6)顯示:蝽科與[盾蝽科+(荔蝽科+兜蝽科)+(同蝽科+土蝽科)]關(guān)系較近,但節(jié)點(diǎn)支持率較低(BS=58);土蝽科與同蝽科為姐妹群,但節(jié)點(diǎn)支持率較低(BS=43)?；?3個PCGs的氨基酸序列構(gòu)建的系統(tǒng)發(fā)育樹(圖7)顯示:蝽科與[(荔蝽科+兜蝽科)+(盾蝽科+土蝽科)]親緣關(guān)系較近,且節(jié)點(diǎn)支持率相對較高(PP=0.97);土蝽科與盾蝽科的姐妹群關(guān)系被強(qiáng)烈支持(PP=1)。

*表示本研究中的兩種昆蟲。

*表示本研究中的兩種昆蟲。

*表示本研究中的兩種昆蟲。

此外,3個系統(tǒng)發(fā)育樹均強(qiáng)烈支持:斯氏珀蝽與珀蝽(P.fimbriata)為姐妹群(BS=1,PP=1),青革土蝽與M.subaeneus為姐妹群(BS=1,PP=1)。

3 討論

本研究使用高通量測序技術(shù)獲得斯氏珀蝽和青革土蝽的完整線粒體基因組數(shù)據(jù)。在兩種昆蟲的線粒體基因組中,大多數(shù)PCGs以ATA、ATG或ATT作為起始密碼子,但也有一些PCGs未以這些密碼子開頭,這種現(xiàn)象在其他半翅目昆蟲[18,29-31]中也有。同時,大多數(shù)PCGs以TAA或TAG作為終止密碼子;而斯氏珀蝽的cox1以及青革土蝽的cox1、cox2、atp6、cox3和nad5具有不完整的終止密碼子T或TA,這種情況通常被認(rèn)為是mRNA通過轉(zhuǎn)錄后,由開頭端的堿基將終止密碼子的結(jié)尾補(bǔ)充完整。trnS1基因由于缺失DHU臂而無法形成完整的三葉草結(jié)構(gòu),這種現(xiàn)象在其他半翅目昆蟲[32]中也有。此外,斯氏珀蝽和青革土蝽的rrnL和rrnS基因都展示出較高的AT含量(73.60%～77.70%)。

近年來,蝽總科的單系性得到了分子生物學(xué)和形態(tài)學(xué)研究[33-35]的支持,然而它的科級水平的系統(tǒng)發(fā)育關(guān)系仍然存在爭議。例如:Lis et al[36]基于rrnL和rrnS基因通過貝葉斯法分析發(fā)現(xiàn),兜蝽科為單系群,且荔蝽科與盾蝽科為姐妹群;而Wu et al[37]利用貝葉斯法、最大似然法和最大簡約法分析發(fā)現(xiàn),蝽科與(同蝽科+澳絲蝽科Lestoniidae)為姐妹群,土蝽科與(荔蝽科+兜蝽科)有較近的親緣關(guān)系。異蝽科與同蝽科的姐妹群關(guān)系得到了許多研究[3,18,38-40]的支持。例如,Xu et al[18]基于線粒體的37個基因以及PCGs的第1、2密碼子加24個RNA基因的研究結(jié)果均支持蝽科與(異蝽科+同蝽科)為姐妹群。

本研究基于測序獲得的斯氏珀蝽和青革土蝽的線粒體基因組數(shù)據(jù),以及蝽總科9個科的113個物種(內(nèi)群)和紅蝽總科2個物種(外群)的13個PCGs的核苷酸序列,第1、2密碼子的核苷酸序列以及氨基酸序列,分別利用最大似然法和貝葉斯法構(gòu)建了系統(tǒng)發(fā)育樹。這兩種分析結(jié)果均強(qiáng)烈支持荔蝽科與兜蝽科親緣關(guān)系較近,與Wu et al[37]和Xu et al[18]的研究結(jié)果一致。然而,兩種分析結(jié)果均不支持異蝽科與同蝽科的姐妹群關(guān)系,且蝽科與土蝽科的系統(tǒng)發(fā)育關(guān)系也不一致。因此,下一步需要通過選取更多的類群以及使用核基因、全基因組等方法更深入地研究蝽總科的系統(tǒng)發(fā)育關(guān)系。此外,兩種分析結(jié)果均表明,斯氏珀蝽與珀蝽為姐妹群,青革土蝽與M.subaeneus為姐妹群。