韭菜遲眼蕈蚊雌雄成蟲差異基因表達(dá)的轉(zhuǎn)錄組分析

胡靜榮，史彩華

長(zhǎng)江大學(xué)農(nóng)學(xué)院，湖北 荊州 434025

韭菜遲眼蕈蚊(BradysiacellarumFrey)俗稱韭蛆，隸屬雙翅目眼蕈蚊科[1]。該蟲是我國蔬菜作物的重要地下害蟲，其寄主范圍廣，可為害大蒜、大蔥、甘藍(lán)、蘿卜、甜瓜、芹菜等7科30多種蔬菜、瓜果和食用菌，尤其喜歡取食百合科的韭菜(AlliumtuberosumRottlerexSpreng)[2，3]。在無任何防治措施的情況下，可造成韭菜減產(chǎn)50%左右，嚴(yán)重時(shí)毀種絕收[4]。

防治韭菜遲眼蕈蚊的方法有臭氧水、昆蟲病原線蟲、粘蟲板、防蟲網(wǎng)和生物菌劑等[5-8]。然而，這些方法均存在一定的局限性，如防治效果差、成本高、操作復(fù)雜等。因此，韭菜遲眼蕈蚊仍然是我國韭菜種植業(yè)的重大生物災(zāi)害。目前，防治韭菜遲眼蕈蚊的首選方法仍然是澆灌化學(xué)農(nóng)藥，然而農(nóng)藥的過度使用存在抗藥性、食品污染和環(huán)境污染等風(fēng)險(xiǎn)[9]。為有效防治韭菜遲眼蕈蚊，研究新的、可持續(xù)性的防控措施勢(shì)在必行。

研究表明，將雌蚊轉(zhuǎn)化為雄蚊是當(dāng)前蚊蟲防控工作的重要突破口。然而，通過生物技術(shù)調(diào)控昆蟲性別必須建立在昆蟲性別決定與性別分化的基礎(chǔ)之上。目前，為揭示蚊蟲的性別差異，KOUTSOS等[10]研究了岡比亞按蚊(Anophelesgambiae)在胚胎、5個(gè)齡期的幼蟲、蛹、成蚊等8個(gè)連續(xù)生命周期中的基因表達(dá)譜，結(jié)果表明其雌蚊表達(dá)量較高的基因主要與免疫有關(guān)，雄蚊表達(dá)量較高的基因主要與代謝過程有關(guān)。WHITTLE等[11]研究表明，埃及伊蚊(Aedesaegypti)卵巢高表達(dá)的基因有2927個(gè)，精巢高表達(dá)的基因有2013個(gè)。然而，關(guān)于韭菜遲眼蕈蚊雌雄差異表達(dá)基因并不清楚。為明確韭菜遲眼蕈蚊雌雄個(gè)體在生殖中的作用與貢獻(xiàn)，筆者采用PacBio Sequel平臺(tái)進(jìn)行第三代轉(zhuǎn)錄組測(cè)序，從基因水平比較雌雄成蟲之間的差異，以期為揭示韭菜遲眼蕈蚊性別決定與性別分化的分子機(jī)理奠定基礎(chǔ)。

1 材料與方法

1.1 材料與試劑

韭菜遲眼蕈蚊由中國農(nóng)業(yè)科學(xué)院蔬菜花卉研究所提供，并在長(zhǎng)江大學(xué)農(nóng)學(xué)院實(shí)驗(yàn)室恒溫飼養(yǎng)3代，設(shè)置溫度(25 ± 1)℃、濕度(70 ± 5)%、光照14h∶10h(晝∶夜)。用吸蟲管分別將未交配的雌雄成蟲收集到不同的1.5mL無菌離心管中，作好標(biāo)記，放液氮中徹底冷凍、備用。

RNA提取試劑Trizol購置于Invitrogen公司(Invitrogen, Carlsbad, CA, USA)；cDNA文庫構(gòu)建試劑盒購置于Illumina公司(TruseqTM RNA Sample prep Kit)；PCR產(chǎn)物純化試劑盒購置于Qiagen公司(Qiagen, CA, USA)；cDNA合成試劑盒購置于Thermo公司(Revert AidTM First Strand cDNA Synthesis Kit)；dNTP購置于TIANGEN。

1.2 方法

1.2.1 總RNA提取、測(cè)序、數(shù)據(jù)校正

將樣品送北京諾禾致源科技股份有限公司測(cè)序，采用Trizol法提取韭菜遲眼蕈蚊雌雄成蟲的總RNA。用帶有Oligo(dT)的磁珠，通過A-T互補(bǔ)配對(duì)的方法與mRNA的ployA結(jié)合富集mRNA，再加入fragmentation buffer打斷mRNA，然后用六堿基隨機(jī)引物(random hexamers)合成cDNA，再加dNTPs、DNA polymerase I合成雙鏈cDNA，并純化。

將純化后的雙鏈cDNA進(jìn)行末端修復(fù)，加A尾和測(cè)序接頭，用AMPure XP beads進(jìn)行片段大小選擇，最后采用PCR富集得到最終的cDNA文庫。對(duì)檢測(cè)合格的文庫進(jìn)行上機(jī)測(cè)序。

用PacBio軟件處理測(cè)序后的原始數(shù)據(jù)獲得Subreads序列，再通過CCS (Circular Consensus Sequence) 算法進(jìn)行校正，獲得環(huán)形一致的CCS序列。將CCS序列分為全長(zhǎng)序列(FLNC)與非全長(zhǎng)序列(nFL)，再將同一轉(zhuǎn)錄本的FLNC序列進(jìn)行hierarchical n*log(n)算法聚類，獲得cluster consensus序列。然后對(duì)全長(zhǎng)序列進(jìn)行polish，獲得高質(zhì)量的Polished consensus序列，用二代轉(zhuǎn)錄組數(shù)據(jù)校正Polished consensus序列，并用CD-HIT軟件去冗余，得到非冗余轉(zhuǎn)錄本，包括CDS序列和unigene序列的數(shù)量。

1.2.2 基因功能注釋

分別用NCBI blast、Blast2GO和KAAS軟件，將1.2.1中拼接得到的unigenes序列與Nr非冗余蛋白質(zhì)序列數(shù)據(jù)庫、Nt核酸序列數(shù)據(jù)庫、KOG同源蛋白簇?cái)?shù)據(jù)庫、GO 基因本體論數(shù)據(jù)庫，以及KEGG日本京都基因和基因組百科全書數(shù)據(jù)庫(e-value<0.00001)進(jìn)行比對(duì)，獲得All Unigene 在每個(gè)數(shù)據(jù)庫中的功能注釋信息。

1.2.3 基因表達(dá)差異分析

按照基因的表達(dá)量(FPKM值)計(jì)算基因在韭菜遲眼蕈蚊雌雄樣本間的差異表達(dá)倍數(shù)，并結(jié)合差異檢驗(yàn)的FDR(false discovery rate)值篩選雌雄樣本間的差異表達(dá)基因。

1.2.4 基因差異富集分析

用KEGG注釋信息得到Unigene基因的pathway注釋，選擇q-value≤0.05的pathway為差異表達(dá)基因中的顯著富集條目。

2 結(jié)果與分析

2.1 韭菜遲眼蕈蚊雌雄成蟲基因功能注釋

2.1.1 Nr數(shù)據(jù)庫注釋

將韭菜遲眼蕈蚊雌雄成蟲基因序列與Nr數(shù)據(jù)庫進(jìn)行比對(duì)，結(jié)果表明，雌成蟲注釋到11950個(gè)基因，對(duì)應(yīng)312個(gè)物種；雄成蟲注釋到16571個(gè)基因，對(duì)應(yīng)353個(gè)物種。同時(shí)，分別對(duì)雌雄成蟲注釋到的物種按基因相似度進(jìn)行排序，取相似度排前20位的物種作圖，其中，雌雄成蟲共有的高相似物種占95%，詳見圖1。

圖1 韭菜遲眼蕈蚊雌雄成蟲基因序列與Nr數(shù)據(jù)庫物種比對(duì)的注釋信息Fig.1 Annotation information on gene sequences of male and female adults of Bradysia cellarum compared with species in Nr database

圖2 韭菜遲眼蕈蚊雌雄成蟲的GO分類 Fig.2 GO taxonomy of male and female adults of Bradysia cellarum

2.1.2 GO分類

將韭菜遲眼蕈蚊雌雄成蟲基因序列進(jìn)行GO注釋，結(jié)果表明，雌蟲注釋上的基因數(shù)為9705個(gè)，雄蟲注釋上的基因數(shù)為12168個(gè)。將基因功能按照細(xì)胞組分、分子功能和生物學(xué)過程進(jìn)行分類，雌蟲分別占11114、11910、19600個(gè)，雄蟲分別占13473、14394、25332個(gè)，詳見圖2。

2.1.3 KOG分類

將韭菜遲眼蕈蚊雌雄成蟲基因序列與KOG數(shù)據(jù)庫進(jìn)行比對(duì)分類，結(jié)果表明，雌蟲比對(duì)上10781個(gè)基因，雄蟲比對(duì)上14520個(gè)基因，詳見圖3。

圖3 韭菜遲眼蕈蚊雌雄成蟲的KOG分類圖Fig.3 KOG taxonomic map of male and female adults of Bradysia cellarum

2.1.4 KEGG分類

圖4 韭菜遲眼蕈蚊雌雄成蟲基因功能注釋KEGG代謝通路分類Fig.4 Classification of gene function annotation to KEGG metabolic pathways in male and female adults of Bradysia cellarum

將韭菜遲眼蕈蚊雌雄成蟲基因序列與KEGG數(shù)據(jù)庫進(jìn)行比對(duì)，結(jié)果表明，從雌雄成蟲轉(zhuǎn)錄組數(shù)據(jù)中分別獲得13064和23488個(gè)轉(zhuǎn)錄本。對(duì)雌雄成蟲的基因分別按6大通路進(jìn)行注釋，結(jié)果表明，雌蟲中參與細(xì)胞過程(cellular processes)、環(huán)境信息處理(environment information processing)、遺傳信息處理(genetic information processing)、人類疾病(human disease)、新陳代謝(metabolism)和機(jī)體系統(tǒng)(organismal systems)的基因通路分別為1560、1237、1363、2643、1701、2160條；雄蟲中參與細(xì)胞過程、環(huán)境信息處理、遺傳信息處理、人類疾病、新陳代謝和機(jī)體系統(tǒng)的基因通路分別為228、2055、1276、4132、3301、3865條，詳見圖4。

2.1.5 數(shù)據(jù)庫之間的維恩圖分析

韭菜遲眼蕈蚊雌成蟲注釋到Nr、Nt、KOG、KEGG、GO等5個(gè)數(shù)據(jù)庫中的基因分別為11950、4091、9562、10012、9704條，其中共同注釋到的基因有3209條，詳見圖5(a)。雄成蟲注釋到Nr、Nt、KOG、KEGG、GO等5個(gè)數(shù)據(jù)庫中的基因分別為16571、6046、12873、16203、12168條，其中共同注釋到的基因有4402條，詳見圖5(b)。

2.2 差異基因表達(dá)分析

2.2.1 韭菜遲眼蕈蚊雌雄成蟲差異基因篩選

根據(jù)韭菜遲眼蕈蚊雌雄成蟲基因表達(dá)倍數(shù)比表明，上調(diào)的基因有9794個(gè)，下調(diào)的基因有8881個(gè)，無顯著差異的基因有5472個(gè)，詳見圖6。

2.2.2 韭菜遲眼蕈蚊雌雄成蟲差異基因維恩圖

雄性成蟲特有基因6833個(gè)，雌性成蟲特有基因3898個(gè)，雌雄成蟲共同擁有的基因17763個(gè)，詳見圖7。

2.3 韭菜遲眼蕈蚊雌雄成蟲差異基因KEGG富集分析

對(duì)雌雄成蟲的差異基因進(jìn)行KEGG富集，結(jié)果表明，上調(diào)基因的富集通路282個(gè)，平均q-value為0.00088；下調(diào)基因的富集通路285個(gè)，平均q-value為0.0048。挑選富集最顯著的20條pathway作圖，其中包含上調(diào)基因1821個(gè)，下調(diào)基因1351個(gè)。差異表達(dá)基因明顯富集在有機(jī)酸代謝、肽代謝過程等功能類別上，參與核糖體生物合成、RNA轉(zhuǎn)運(yùn)等，詳見圖8。

注：橫坐標(biāo)代表基因在雌雄成蟲樣品中表達(dá)倍數(shù)，縱坐標(biāo)代表基因表達(dá)量在統(tǒng)計(jì)學(xué)上的顯著程度。其中，校正后的p-value越小，-log10越大，即差異越顯著。圖中散點(diǎn)代表各個(gè)基因，藍(lán)色圓點(diǎn)表示無顯著差異的基因，紅色圓點(diǎn)表示顯著上調(diào)的基因，綠色圓點(diǎn)表示顯著下調(diào)的基因。圖6 韭菜遲眼蕈蚊雌雄基因差異表達(dá)分析Fig.6 Differential expression analysis of male and female genes of Bradysia cellarum

注：每個(gè)大圓中數(shù)字的和代表某種樣品差異基因總數(shù)，圓圈交疊處的數(shù)字表示雌雄成蟲共有的差異基因。圖7 韭菜遲眼蕈蚊雌雄成蟲差異基因維恩圖Fig.7 Venn map of differential genes of male and female adults of Bradysia cellarum

3 討論與結(jié)論

轉(zhuǎn)錄組通常指物種或組織中所有 RNA 轉(zhuǎn)錄產(chǎn)物的集合[12]。不同物種或同一物種的不同組織或同一組織的不同時(shí)期，它們的轉(zhuǎn)錄組表達(dá)情況均存在很大差異[13]。分析生物雌雄樣本的轉(zhuǎn)錄組數(shù)據(jù)，有利于挖掘系統(tǒng)發(fā)育基因、性別關(guān)聯(lián)基因、特異遺傳性狀基因等[14]。

注：縱軸表示pathway名稱，橫軸表示pathway對(duì)應(yīng)的富集因子。qvalue值越小則顏色越接近紅色，每個(gè)pathway下包含的差異基因的多少用點(diǎn)的大小來表示。圖8 韭菜遲眼蕈蚊雌雄成蟲差異基因KEGG富集散點(diǎn)圖Fig.8 Scatter plot of differential gene KEGG enrichment of male and female adults of Bradysia cellarum

該研究分別構(gòu)建了韭菜遲眼蕈蚊雌雄成蟲的2個(gè)cDNA文庫。其中，在雌性cDNA文庫中獲得13064條unigenes，同時(shí)被Nr、Nt、KOG、KEGG和GO等數(shù)據(jù)庫注釋的基因有3209條；在雄性cDNA文庫中獲得23488條unigenes，同時(shí)被Nr、Nt、KOG、KEGG和GO等數(shù)據(jù)庫注釋的基因有4402條。雄性成蟲高表達(dá)的基因數(shù)顯著多于雌性成蟲，可能與雄性成蟲活動(dòng)能力略強(qiáng)于雌性成蟲有關(guān)[15]。

韭菜遲眼蕈蚊雌雄成蟲表現(xiàn)出高度的性別二態(tài)性，比如雌性個(gè)體大，雄性個(gè)體小；雌蟲腹部粗大，末端細(xì)而尖，雄蟲有1對(duì)抱握器等[16]。轉(zhuǎn)錄組數(shù)據(jù)分析有助于闡明性別二態(tài)性的分子機(jī)理，對(duì)揭示性別分化和性別決定具有重要意義。研究結(jié)果表明，韭菜遲眼蕈蚊偏向雄性的基因多于雌性，其中雄性特有的基因6833條，雌性特有的基因3898條，與EADS等[17]研究水蚤的結(jié)果一致。然而，赤擬谷盜(Triboliumcastaneum)[18]、岡比亞蚊[19]和埃及伊蚊[11]的特有基因數(shù)量偏向雌性，從而暗示物種不同，性別基因的偏向并不相同。KOUTSOS等[10]研究表明，岡比亞按蚊雌蚊的高表達(dá)基因主要與免疫有關(guān)，雄蚊高表達(dá)基因主要與代謝有關(guān)。韭菜遲眼蕈蚊的高表達(dá)基因之所以偏向雄性，可能與雄性超強(qiáng)的交配能力和強(qiáng)于雌性的飛行能力需要超強(qiáng)的代謝功能有關(guān)[20]，但具體原因有待進(jìn)一步研究。

轉(zhuǎn)錄組序列的拼接算法包括有參考基因組序列的拼接方法和無參考基因組序列的拼接方法：參考基因組序列的拼接算法主要有Cufflinks[21]和Scripture[22]等開發(fā)軟件；無參考基因組序列的拼接算法主要有ABySS[23]、SOAP denovo[24]、Trinity[25]等開發(fā)軟件。笗者利用新一代高通量測(cè)序技術(shù)對(duì)韭菜遲眼蕈蚊雌雄成蟲進(jìn)行轉(zhuǎn)錄組測(cè)序，再采用無參考基因組序列的拼接算法對(duì)韭菜遲眼蕈蚊序列進(jìn)行拼接，與序列表達(dá)標(biāo)簽(EST)技術(shù)進(jìn)行比較，具有測(cè)序數(shù)據(jù)覆蓋度高的優(yōu)點(diǎn)。該研究為揭示韭菜遲眼蕈蚊性別決定、生殖生長(zhǎng)、營養(yǎng)代謝等分子機(jī)制奠定了基礎(chǔ)，同時(shí)，豐富了昆蟲遺傳信息庫，有助于韭菜遲眼蕈蚊綠色防控新技術(shù)研發(fā)。