劉華偉，李朝緒，李 芬，呂朝軍，吳少英*，覃偉權(quán)*

（1.海南大學(xué)熱帶作物學(xué)院，?？?570228；2.中國(guó)熱帶農(nóng)業(yè)科學(xué)院椰子研究所，?？?571399；3.海南大學(xué)植物保護(hù)學(xué)院，?？?570228；4.海南省院士團(tuán)隊(duì)創(chuàng)新中心，?？?571339；5.海南省熱帶油料作物生物學(xué)重點(diǎn)實(shí)驗(yàn)室，?？?571399）

椰心葉甲Brontispalongissima(Gestro)屬鞘翅目Coleoptera鐵甲科Hispidae害蟲(chóng)，其成蟲(chóng)和幼蟲(chóng)均以棕櫚科植物未展開(kāi)的心葉為食，并且?guī)缀跛械淖貦翱浦参锒伎墒芷錇楹1,2]。自2002年6月首次在海南省發(fā)現(xiàn)椰心葉甲為害，至2006年僅四年的時(shí)間，已蔓延至海南全省[3]。椰心葉甲嚙小蜂TetrastichusbrontispaeFerrière是椰心葉甲的優(yōu)勢(shì)天敵寄生蜂。目前國(guó)內(nèi)外已經(jīng)對(duì)該蜂的生物生態(tài)學(xué)、室內(nèi)繁殖和野外釋放等方面進(jìn)行了相關(guān)研究[4-8]，但由于缺乏遺傳信息，椰心葉甲嚙小蜂的遺傳多樣性、基因功能等的研究還較為滯后。二代測(cè)序技術(shù)的轉(zhuǎn)錄組測(cè)序技術(shù)（RNA-seq），具有低成本、快速和高準(zhǔn)確性等優(yōu)勢(shì)，能夠在缺乏基因信息的條件下獲得物種的代謝和生長(zhǎng)規(guī)律，并揭示其基因與生物學(xué)特性內(nèi)在關(guān)聯(lián)，同時(shí)還可獲得物種大多數(shù)的轉(zhuǎn)錄產(chǎn)物[9]。本研究通過(guò)RNA-seq技術(shù)對(duì)兩個(gè)種群椰心葉甲嚙小蜂進(jìn)行轉(zhuǎn)錄組測(cè)序，并進(jìn)行了功能基因注釋和分析，旨在為椰心葉甲嚙小蜂的功能基因挖掘等分子生物學(xué)研究中提供理論基礎(chǔ)。

1 材料與方法

1.1 材料

為后續(xù)基于轉(zhuǎn)錄組數(shù)據(jù)分析連續(xù)自交是否會(huì)引起椰心葉甲嚙小蜂的種群退化，退化程度以及退化對(duì)其影響。本試驗(yàn)中所采用的兩個(gè)種群椰心葉甲嚙小蜂分別選用中國(guó)熱帶農(nóng)業(yè)科學(xué)院椰子研究所天敵工廠提供的經(jīng)過(guò)復(fù)壯（每年一次）的種群和中國(guó)熱帶農(nóng)業(yè)科學(xué)院環(huán)境與植物保護(hù)研究所天敵工廠提供的自引進(jìn)后隔離繁殖未經(jīng)復(fù)壯蜂種。兩個(gè)種群各設(shè)置3次試驗(yàn)重復(fù)，每個(gè)重復(fù)由20頭椰心葉甲嚙小蜂（雌雄蜂混合的1～2日齡成蜂）混合而成。

1.2 總RNA提取，文庫(kù)構(gòu)建和轉(zhuǎn)錄組測(cè)序

本試驗(yàn)將收集的兩個(gè)種群椰心葉甲嚙小蜂經(jīng)液氮快速冷凍后，用干冰將樣品送至上海派森諾生物科技股份有限公司進(jìn)行總RNA提取、文庫(kù)構(gòu)建和轉(zhuǎn)錄組測(cè)序。

1.3 序列組裝和轉(zhuǎn)錄本拼接

RNA-seq測(cè)序完成后，統(tǒng)計(jì)總reads的數(shù)量和長(zhǎng)度、Q30值、模糊堿基（N）所占比例以及Q20和Q30所占比例等。經(jīng)數(shù)據(jù)過(guò)濾后獲得 clean reads，繼而統(tǒng)計(jì) clean reads的數(shù)量、總長(zhǎng)度和占比。通過(guò) Trinity Software（http:///trinityrnaseq.Github.io/）軟件對(duì)clean reads進(jìn)行拼裝獲得轉(zhuǎn)錄本（Transcript），選取每個(gè)基因最長(zhǎng)的轉(zhuǎn)錄本作為代表序列（Unigene）[9,10]，并進(jìn)一步統(tǒng)計(jì)Transcript和Unigene的總長(zhǎng)度、序列總數(shù)、最大長(zhǎng)度、平均長(zhǎng)度、N50、N50%、N90、N90%以及GC含量等指標(biāo)。

1.4 Unigene功能注釋

根據(jù)基因的相似性，通過(guò) BLAST比對(duì)工具（參數(shù)設(shè)置：E≤1e-5）將椰心葉甲嚙小蜂轉(zhuǎn)錄組獲得的Unigene與 NCBI non-redundant protein sequences（NR）數(shù)據(jù)庫(kù)、Gene Ontology（GO）數(shù)據(jù)庫(kù)、Kyoto Encyclopedia of Genes and Genome（KEGG）數(shù)據(jù)庫(kù)和 evolutionary genealogy of genes: Non-supervised Orthologous Groups（eggNOG）數(shù)據(jù)庫(kù)比對(duì)，進(jìn)行功能注釋。

2 結(jié)果與分析

2.1 椰心葉甲嚙小蜂轉(zhuǎn)錄組測(cè)序數(shù)據(jù)及轉(zhuǎn)錄本拼接質(zhì)量情況

轉(zhuǎn)錄組數(shù)據(jù)集已保存在NCBI SRA數(shù)據(jù)庫(kù)中，編號(hào)：PRJNA678031。對(duì)每個(gè)樣品的原始數(shù)據(jù)（Raw data）進(jìn)行統(tǒng)計(jì)（表1），兩個(gè)種群的6個(gè)樣品分別得到48514092、43095558、44266534、49405442、47396830、47326810條reads，且每個(gè)樣本的堿基總量均在6.46 Gb以上；堿基百分比Q20均大于97.47%、堿基百分比Q30均大于93.44%，含量相近；模糊堿基占比在0.001485%～0.001512%的較低水平，通過(guò)數(shù)據(jù)過(guò)濾后統(tǒng)計(jì)發(fā)現(xiàn)（表2），Clean Reads占reads總數(shù)的93.94%～94.37%。對(duì)Transcript和Unigene序列進(jìn)行統(tǒng)計(jì)（表3），共得到78930條Transcript和29535條Unigene，Unigene的長(zhǎng)度分布統(tǒng)計(jì)圖見(jiàn)圖1，總長(zhǎng)度為51330466 bp，平均長(zhǎng)度為1737.95 bp，N50值為3547 bp。

圖1 椰心葉甲嚙小蜂轉(zhuǎn)錄組Unigene組裝長(zhǎng)度分布統(tǒng)計(jì)Fig.1 Statistics of Unigene length distribution of T.brontispae transcriptome

表1 椰心葉甲嚙小蜂轉(zhuǎn)錄組測(cè)序數(shù)據(jù)及質(zhì)量情況統(tǒng)計(jì)Table 1 Transcriptome sequencing data and quality statistics of T.brontispae

表2 椰心葉甲嚙小蜂轉(zhuǎn)錄組高質(zhì)量Reads與堿基統(tǒng)計(jì)Table 2 High quality Reads and base statistics of T.brontispae transcriptome

表3 椰心葉甲嚙小蜂轉(zhuǎn)錄組Transcript和Unigene序列統(tǒng)計(jì)表Table 3 Statistics of transcript and Unigene sequences of T.brontispae transcriptome

以上數(shù)據(jù)表明樣本的測(cè)序質(zhì)量，文庫(kù)構(gòu)建質(zhì)量和拼接的完整性都較好，可以用于后續(xù)的分析。

2.2 Unigene功能注釋情況

通過(guò)Unigene和四個(gè)公共數(shù)據(jù)庫(kù)的比對(duì)并進(jìn)行功能注釋（表4）。共有13796條Unigene被注釋，其在4個(gè)數(shù)據(jù)庫(kù)中的注釋如下：13401條基因（45.37%）在NR數(shù)據(jù)庫(kù)中獲得了注釋，3834條基因（12.98%）在GO數(shù)據(jù)庫(kù)中獲得了注釋，12707條基因（43.02%）在eggNOG數(shù)據(jù)庫(kù)中獲得了注釋，5999條基因（20.31%）在KEGG數(shù)據(jù)庫(kù)中獲得了注釋。

表4 椰心葉甲嚙小蜂轉(zhuǎn)錄組注釋結(jié)果匯總Table 4 Summary of annotated results of T.brontispae transcriptome

2.3 椰心葉甲嚙小蜂轉(zhuǎn)錄組Unigene的NR功能分類

通過(guò) Blastx將椰心葉甲嚙小蜂轉(zhuǎn)錄組中的Unigene與NR數(shù)據(jù)庫(kù)進(jìn)行比對(duì)，注釋到13401條Unigene。統(tǒng)計(jì)注釋結(jié)果并繪制物種分布圖（圖2），結(jié)果顯示椰心葉甲嚙小蜂與蠅繭蜂Diachasmaalloeum、佛羅里達(dá)弓背蟻Camponotusfloridanus、黑褐毛蟻Lasiusniger、赤眼蜂Trichogrammapretiosum、多胚跳小蜂Copidosomafloridanum、榕小蜂Ceratosolensolmsimarchali、麗蠅蛹集金小蜂Nasoniavitripennis都有同源序列分布，其中與麗蠅蛹集金小蜂相似序列最多，占50.34%，與榕小蜂、佛羅里達(dá)跳小蜂和赤眼蜂分別有13.3%、7.42%和5.95%的相似序列，與黑蟻、佛羅里達(dá)弓背蟻和蠅繭蜂中的同源序列較少（0.98%～1.57%）。

圖2 椰心葉甲嚙小蜂轉(zhuǎn)錄組NR注釋結(jié)果統(tǒng)計(jì)圖Fig.2 NR annotated species distribution of T.brontispae transcriptome

2.4 椰心葉甲嚙小蜂轉(zhuǎn)錄組Unigene的GO分類注釋

對(duì)椰心葉甲嚙小蜂轉(zhuǎn)錄組Unigene進(jìn)行GO功能分類，共有3834條（12.98%）Unigene獲得注釋。將注釋得到的Unigene劃分為三大類（生物過(guò)程、細(xì)胞成分和分子功能）67個(gè)分支（圖3），統(tǒng)計(jì)注釋到每一類的基因數(shù)量，發(fā)現(xiàn)在生物過(guò)程類中細(xì)胞過(guò)程（1527條Unigene）的占比最大（39.83%）；在細(xì)胞成分大類中膜（1324條Unigene）的占比最大（34.53%）；在分子功能大類中結(jié)合（1733條Unigene）的占比最大（45.20%）。

圖3 椰心葉甲嚙小蜂轉(zhuǎn)錄組Unigene的GO功能分類Fig.3 GO analysis of T.brontispae transcriptome Unigene

2.5 椰心葉甲嚙小蜂轉(zhuǎn)錄組Unigene的eggNOG功能注釋

將椰心葉甲嚙小蜂轉(zhuǎn)錄組的Unigene與eggNOG數(shù)據(jù)庫(kù)進(jìn)行比對(duì)，注釋到12707條Unigene，根據(jù)功能可將其劃分為26個(gè)功能區(qū)域（圖4），統(tǒng)計(jì)注釋到各類功能的基因數(shù)量，發(fā)現(xiàn)無(wú)特征基因占比最多（63.06%），其次是新陳代謝類（33.49%），細(xì)胞過(guò)程和信號(hào)傳遞，信息存儲(chǔ)與處理類分別占 20.82%和19.67%。

圖4 椰心葉甲嚙小蜂轉(zhuǎn)錄組Unigene的eggNOG功能分類Fig.4 eggNOG analysis of T.brontispae transcriptome Unigene

2.6 椰心葉甲嚙小蜂轉(zhuǎn)錄組Unigene的KEGG代謝途徑分類分析

將椰心葉甲嚙小蜂轉(zhuǎn)錄組的Unigene與KEGG數(shù)據(jù)庫(kù)進(jìn)行對(duì)比，共有5999條Unigene獲得注釋，涉及的代謝通路可歸為五個(gè)大類（新陳代謝、遺傳信息處理、環(huán)境信息處理、細(xì)胞進(jìn)程和有機(jī)系統(tǒng)）35個(gè)子類（圖 5）。統(tǒng)計(jì)注釋到各類通路的基因數(shù)量，發(fā)現(xiàn)在新陳代謝一類中獲得注釋最多的是碳水化合物代謝（299條Unigene）；遺傳信息處理一類中翻譯獲得注釋最多（361條Unigene）；環(huán)境信息處理一類中獲得注釋最多的是信號(hào)轉(zhuǎn)導(dǎo)（735條Unigene），細(xì)胞進(jìn)程和有機(jī)系統(tǒng)兩類中獲得注釋最多的分別為運(yùn)輸和分解代謝（382條Unigene）和內(nèi)分泌系統(tǒng)（388條Unigene）。

圖5 椰心葉甲嚙小蜂轉(zhuǎn)錄組KEGG注釋統(tǒng)計(jì)圖Fig.5 KEGG annotation statistical chart of T.brontispae transcriptome

3 討論

第二代高通量測(cè)序技術(shù)（Next-Generation Sequencing，NGS）因測(cè)序時(shí)間短、成本低、高準(zhǔn)確性和所獲得數(shù)據(jù)量大等優(yōu)點(diǎn)，被廣泛應(yīng)用于非模式生物分子生物學(xué)研究中[9,10]。椰心葉甲嚙小蜂是椰心葉甲的優(yōu)勢(shì)寄生蜂，但是由于其遺傳信息數(shù)據(jù)的缺乏，對(duì)其分子生物學(xué)的研究仍然較少。本研究旨在通過(guò)該技術(shù)分別對(duì)經(jīng)過(guò)復(fù)壯和未經(jīng)復(fù)壯兩個(gè)種群的椰心葉甲嚙小蜂進(jìn)行轉(zhuǎn)錄組測(cè)序和拼裝，揭示椰心葉甲嚙小蜂整體基因表達(dá)特征，并為后續(xù)基于轉(zhuǎn)錄組數(shù)據(jù)分析復(fù)壯對(duì)椰心葉甲嚙小蜂的影響，種群遺傳結(jié)構(gòu)等研究提供數(shù)據(jù)支撐。

一般來(lái)說(shuō)，Q30在80%以上，N50值不小于800 bp（N50值越大表示長(zhǎng)片段越多）就可以認(rèn)為測(cè)序質(zhì)量可靠，序列組裝的完整性較好[11-13]。本研究通過(guò)無(wú)參考基因組分析椰心葉甲嚙小蜂轉(zhuǎn)錄組的特異性，序列拼接后共獲得29535條Unigene，Q30和N50值分別為93.44%和3547 bp，可以認(rèn)為本研究中所得到的測(cè)序結(jié)果質(zhì)量和序列拼裝的完整性都較好，可以滿足后續(xù)分析的基本要求。

對(duì)組裝后所獲得的29535條Unigene在NR、GO、KEGG和eggNOG數(shù)據(jù)庫(kù)進(jìn)行基因功能注釋，共有13796條Unigenes被注釋，仍有15739條Unigenes未被注釋。這一結(jié)果在許多生物的轉(zhuǎn)錄組測(cè)序中都有存在，如植物的彩色馬蹄蓮[10,14]，動(dòng)物中的云錦杜鵑[15]和昆蟲(chóng)中的麗蠅蛹集金小蜂[16]都有出現(xiàn)，可能由于椰心葉甲嚙小蜂缺乏基因組方面研究的基礎(chǔ)資料，使得部分Unigene在數(shù)據(jù)庫(kù)中無(wú)法得以注釋[11]，也可能是由于部分Unigene片段太短或椰心葉甲嚙小蜂中存在新的功能基因而導(dǎo)致的[17]。

與Nr公共數(shù)據(jù)庫(kù)比對(duì)發(fā)現(xiàn)，有13401條Unigene獲得注釋（45.37%），其中與麗蠅蛹集金小蜂相似的序列最多（50.34%）。在椰心葉甲嚙小蜂的轉(zhuǎn)錄組數(shù)據(jù)中有19.36%的Unigene屬于其他序列，這在很多物種的轉(zhuǎn)錄組測(cè)序中也有出現(xiàn)[18-20]，可能屬于椰心葉甲嚙小蜂自身特有的與大多數(shù)物種不同的序列，也可能這些基因?qū)儆诜蔷幋aRNA或不是功能基因。此外，還有54.63%的Unigene未獲得相關(guān)注釋信息，此現(xiàn)象在其他昆蟲(chóng)中也存在，如二化螟盤(pán)絨繭蜂，菜粉蝶等[21,22]，其原因可能是Unigene片段過(guò)短，基因信息的暫時(shí)缺乏，或者椰心葉甲嚙小蜂中存在新的功能基因[23]。在G0數(shù)據(jù)庫(kù)中有3834條Unigene獲得注釋，僅占總Unigene的12.98%，相對(duì)注釋率較低，這進(jìn)一步說(shuō)明了椰心葉甲嚙小蜂基因信息的缺乏[24]，也可能與G0數(shù)據(jù)庫(kù)信息不夠完善有關(guān)[25]。eggNOG數(shù)據(jù)庫(kù)中注釋到的Unigene也較多（43.02%），除無(wú)特征序列外，注釋到與新陳代謝有關(guān)的Unigene較多（24.44%），說(shuō)明椰心葉甲嚙小蜂自身具有較強(qiáng)的代謝能力。在 KEGG數(shù)據(jù)庫(kù)中注釋到 5999條 Unigene，同樣與新陳代謝有關(guān)的代謝通路最多，有 1842條 Unigene（30.7%），與eggNOG分析結(jié)果相一致。這些注釋結(jié)果為在分子層面研究椰心葉甲嚙小蜂提供了充足的依據(jù)。

綜上，本研究利用轉(zhuǎn)錄組測(cè)序技術(shù)，對(duì)椰心葉甲嚙小蜂復(fù)壯和未復(fù)壯兩個(gè)種群進(jìn)行了轉(zhuǎn)錄組測(cè)序，并利用四大數(shù)據(jù)庫(kù)對(duì)獲取的Unigene進(jìn)行了功能注釋、分類和代謝途徑預(yù)測(cè)等，為下一步深入研究椰心葉甲嚙小蜂兩個(gè)種群的差異表達(dá)基因和其他分子生物學(xué)研究提供了參考。