眼優(yōu)角蚱兩性成體轉(zhuǎn)錄組差異分析

覃英燦，容萬韜，李月妹，陸桂紅，李曉東

(河池學(xué)院化學(xué)與生物工程學(xué)院，廣西宜州 546300)

蚱類昆蟲在分類上屬于直翅目Orthoptera蝗亞目Caelifera蚱總科Tetrigoidea，由于其背板呈菱形，舊稱菱蝗。蚱類昆蟲為不完全變態(tài)昆蟲，在整個昆蟲綱中屬于較為原始的類群，約在2.2億年前，由一個共同祖先進化而來(Songetal.,2015)。大多數(shù)種類分布于熱帶和亞熱帶，是生物多樣性和動物區(qū)系的重要組成部分，世界已知 9科 270屬1 880多種，我國已知7科57屬788種，占世界三分之一多(鄧維安,2016)。由于蚱類昆蟲既不是益蟲也不是害蟲，所以國內(nèi)外對蚱類昆蟲的研究不多，主要集中在分類學(xué)方面，其它方面較少。但近年來發(fā)現(xiàn)，由于蚱類昆蟲多以苔蘚、地衣或腐殖質(zhì)為食，屬地棲性昆蟲，運動能力不強，與環(huán)境特征聯(lián)系緊密，可以作為重要的生態(tài)指示物種(鄧維安等,2007;Tsuruietal.,2010;Wennerstenetal.,2012)。因此有必要采用新的研究技術(shù)，對這一古老的類群進行多方面的研究，以挖掘該類群的科學(xué)和應(yīng)用方面的研究價值。

眼優(yōu)角蚱Eucriotettixoculatus屬于刺翼蚱科Scelimenidae優(yōu)角蚱屬Eucriotettix，廣泛分布于中國南部各省及南亞地區(qū)，屬于蚱類中的廣布和常見種類(鄭哲民,2005)，具有代表性。本研究基于轉(zhuǎn)錄組高通量測序技術(shù)，對眼優(yōu)角蚱E.oculatus雌雄兩性蟲體之間的轉(zhuǎn)錄差異進行研究，以從轉(zhuǎn)錄組層面了解蚱總科昆蟲雌雄分化的分子基礎(chǔ)，為昆蟲分子生物學(xué)、發(fā)育生物學(xué)和進化生物學(xué)提供一定的基礎(chǔ)資料。

1 材料與方法

1.1 材料準(zhǔn)備

為了盡可能避免不同環(huán)境因素和發(fā)育時期引起的雌、雄蟲之間轉(zhuǎn)錄組的差異，本研究所用的眼優(yōu)角蚱是2019年6-7月分別在廣西河池市宜州區(qū)和南丹縣不同地點采集的雌性和雄性成蟲帶回實驗室交配產(chǎn)卵后，取其卵塊在昆蟲飼養(yǎng)室中孵化并在穩(wěn)定的環(huán)境條件下從1齡飼養(yǎng)到成蟲后收集的進入成蟲期時間相近的雌、雄成蟲樣品。雌蟲樣品組(編號AF)由9個生物學(xué)重復(fù)樣品(編號為AF1-AF9)組成，且每個生物學(xué)重復(fù)樣品均由3頭雌蟲組成。雄蟲樣品組(編號AM)也是由9個生物學(xué)重復(fù)樣品(編號為AM1-AM9)組成，每個樣品也由3頭雄蟲組成。所有樣品均為活蟲經(jīng)48 h饑餓后用液氮凍存。

1.2 轉(zhuǎn)錄組測序

使用QIAGEN動物組織RNA提取試劑盒提取各樣品總RNA并構(gòu)建cDNA文庫，然后基于Illumina技術(shù)測序平臺，采用雙末端測序(Pair-End)法，完成18個樣品二代轉(zhuǎn)錄組測序分析。

1.3 基因差異表達分析

本研究測定的所有二代轉(zhuǎn)錄組數(shù)據(jù)與本實驗室前期已測定的眼優(yōu)角蚱三代轉(zhuǎn)錄組數(shù)據(jù)(由卵和各齡期雌、雄蟲所組成的混合樣品經(jīng)Pacbio Sequel II技術(shù)測序平臺利用單分子實時測序法獲得的三代數(shù)據(jù)，編號為EO)聯(lián)合分析。首先使用二代轉(zhuǎn)錄組數(shù)據(jù)對三代轉(zhuǎn)錄組數(shù)據(jù)進行轉(zhuǎn)錄本的校正，然后利用CD-HIT軟件(Fuetal.,2012)對校正后的consensus序列進行去冗余，接著對去冗余的序列基于Nr(Lietal.,2002)，Nt，Pfam(Finnetal.,2008)，KOG/COG(Tatusovetal.,2003)，Swiss-prot(Amosetal.,2000)，KEGG(Kanehisaetal.,2004)，GO(Ashburneretal.,2000)共7種數(shù)據(jù)庫進行基因功能注釋。同時以去冗余后得到的轉(zhuǎn)錄本作為基因的參考序列(ref)，將Illumina測序得到的二代轉(zhuǎn)錄組測序數(shù)據(jù)中每個樣品的clean reads比對到ref上。利用RSEM軟件(Lietal.,2011)對比對結(jié)果進行統(tǒng)計，進一步得到每個樣品比對到每個基因上的readcount值，并對其進行FPKM轉(zhuǎn)換。基于基因表達的RPKM值，使用DEseq差異表達基因鑒定程序計算AF和AM樣品組之間表達量不同的基因。差異基因篩選標(biāo)準(zhǔn)為|log2 (FoldChange)|>0和qvalue<0.05。

1.4 差異基因表達的GO和KEGG分析

采用GOseq(Youngetal.,2010)和KOBAS軟件(Chenetal.,2011)對差異基因集進行GO功能富集分析，KEGG通路富集分析。富集分析基于超幾何分布原理，其中差異基因集為差異顯著分析所得差異基因并注釋到GO或KEGG數(shù)據(jù)庫的基因集，背景基因集為所有進行差異顯著分析的基因并注釋到GO或KEGG數(shù)據(jù)庫的基因集。富集分析結(jié)果是對差異比較組合的所有差異基因集、上調(diào)基因集、下調(diào)基因集進行富集。

1.5 實時熒光定量PCR驗證

從眼優(yōu)角蚱雌性和雄性轉(zhuǎn)錄組中篩選出的差異表達基因中隨機挑取9個，以β-actin作為內(nèi)參基因，采用實時熒光定量PCR (qRT-PCR)的方法對轉(zhuǎn)錄組的準(zhǔn)確性進行驗證(驗證樣品與轉(zhuǎn)錄組測序樣品一致)，驗證基因和對應(yīng)的引物序列見表1。使用反轉(zhuǎn)錄試劑盒將總RNA反轉(zhuǎn)錄合成cDNA，反轉(zhuǎn)錄體系為40 μL(冰上配制)：5×PrimeScriptRTMaster Mix (Perfect Real-time) 8 μL和RNA模板2 000 ng，最后加入RNase-free dH2O補全至40 μL。反應(yīng)條件為：37°C 15 min，85°C 5 s，置于4℃冰箱保存。采用ChamQTM SYBR?qPCR Master Mix染料試劑盒進行實時熒光定量PCR擴增，PCR反應(yīng)體系20 μL：SYBR Premix Ex Taq (2×) 10 μL，10 μmol/L上、下游引物各1 μL，cDNA模板2 μL，ddH2O補全至20 μL。PCR反應(yīng)程序為：95℃預(yù)變性30 s；95℃變性5 s，60℃退火34 s，40個循環(huán)。

表1 實時熒光定量引物Table 1 Primer sequences for Real-time quantitative PCR

2 結(jié)果與分析

2.1 轉(zhuǎn)錄本校正和去冗余

通過本研究測定的二代轉(zhuǎn)錄組數(shù)據(jù)對三代轉(zhuǎn)錄組數(shù)據(jù)校正后，三代consensus的堿基數(shù)目有所提升，其它指標(biāo)基本無變化(表2)，說明眼優(yōu)角蚱轉(zhuǎn)錄組二代和三代數(shù)據(jù)整體質(zhì)量較高。進一步按照序列之間95%的相似性對校正后的轉(zhuǎn)錄本序列進行聚類去冗余，最終得到13 432條transcripts，總長度26 977 000 bp，平均長度2 009 bp，最大長度8 297 bp，N50長度遠大于1 000 bp(表3)，說明測序效果較好，可以進行后續(xù)分析。

表2 轉(zhuǎn)錄本校正前后長度分布統(tǒng)計表Table 2 Statistical table of length distribution of transcripts before and after correction

表3 去冗余后序列長度分布統(tǒng)計表Table 3 Statistical table of sequence length distribution after de-redundancy

2.2 基因功能注釋

對去冗余之后的序列進行基因功能注釋，通過七大數(shù)據(jù)庫共有12 074條transcripts被注釋，占整個transcripts的90%，未成功注釋的transcripts僅占比10%，說明注釋成功率比較高(圖1)。對Nr數(shù)據(jù)庫所獲得的基因注釋與已公布的基因組信息的物種進行匹配，其中被匹配最多的物種有20個(圖2)，分別是內(nèi)華達古白蟻Zootermopsisnevadensis、亞利桑那沫蟬Clastopteraarizonana、飛蝗Locustamigratoria、麥莖蜂Cephuscinctus、赤擬谷盜Triboliumcastaneum等皆為昆蟲綱物種，因此可證明本研究轉(zhuǎn)錄組注釋結(jié)果的可靠性。

圖1 七大數(shù)據(jù)庫注釋結(jié)果統(tǒng)計圖Fig.1 Statistical chart of annotation results of seven databases

圖2 Nr數(shù)據(jù)庫物種注釋統(tǒng)計圖Fig.2 Statistical chart of species annotation in NR database

2.3 差異表達基因篩選

為了獲得雌性和雄性眼優(yōu)角蚱的表達量不同的基因，運用FPKM法對雌性組(AF)和雄性組(AM)的transcripts進行分析，共獲得兩者之間的差異表達基因3 597個，其中雌性組(AF)相對于雄性組(AM)，共獲得1 295個顯著上調(diào)基因，獲得2 302個顯著下調(diào)基因(圖3)。

圖3 眼優(yōu)角蚱雌性與雄性間差異表達基因的比較Fig.3 Comparison of differential expressed genes between male and female Eucriotettix oculatus

2.4 差異表達基因GO富集結(jié)果

對差異表達基因進行GO功能分析可以了解差異表達基因的生物學(xué)意義。差異表達基因GO功能分類顯示，共有2 344條transcripts被GO注釋，并被顯著富集到生物學(xué)過程(Biological process)、細胞組分(Cellular component)和分子功能(Molecular function)三大類的主要35個亞類中。其中細胞成分組織或合成(Cellular component organization or biogenesis)、細胞成分組織(Cellular component organization)、細胞器組織(Organelle organization)、細胞(Cell)、細胞成分(Cell part)、細胞內(nèi)(Intracellular)、細胞內(nèi)成分(Intracellular part)、有機環(huán)狀化合物結(jié)合(Organic cyclic compound binding)、雜環(huán)化合物結(jié)合(Heterocyclic compound binding)、核酸結(jié)合(Nucleic acid binding)分別是三大類中顯著富集差異基因數(shù)量最多的條目(圖4)。

圖4 差異表達基因GO功能分類Fig.4 Go functional classification of differentially expressed genes

雌性相對于雄性，表達上調(diào)的1 295個基因中，被GO注釋的有922個，但沒有顯著富集的亞類。表達下調(diào)的2 302個基因中，被GO注釋的有1 422個，顯著富集到三大類的47個亞類。其中富集差異基因數(shù)目最多的是分子功能(Molecular function)大類中的蛋白質(zhì)結(jié)合(Protein binding)亞類的468個。其次是細胞組分(Cellular component)大類中的細胞器組分(Organelle part)和細胞內(nèi)細胞器組成(Intracellular organelle part)，分別是157和155個。再次是生物過程(Biological process)中的細胞成分組織或生物發(fā)生(Cellular component organization or biogenesis)和細胞成分組織(Cellular component organization)，分別是140個和133個(圖5)。其余亞類富集到的基因數(shù)據(jù)均在100以下。

圖5 下調(diào)表達基因GO功能分類Fig.5 Go functional classification of down regulated expressed genes

2.5 差異表達基因KEGG富集結(jié)果

雌性相對于雄性，表達上調(diào)的1 295個基因經(jīng)過KEGG注釋，共有1 188個基因被注釋到248條通路，富集顯著的前20條通路(圖6)，其中神經(jīng)活性配體受體相互作用(Neuroactive ligand-receptor interaction)、DNA復(fù)制(DNA replication)、蛋白質(zhì)消化吸收(Protein digestion and absorption)、錯配修復(fù)(Mismatch repair)、嘧啶代謝(Pyrimidine metabolism)和胰腺分泌(Pancreatic secretion)通路是富集最顯著的6條通路(qvalue<0.05)。表達下調(diào)的2 303個基因中共有1 886個基因被注釋到264條通路，富集顯著的前20條通路(圖7)，其中蛋白酶體(Proteasome)、氮代謝(Nitrogen metabolism)、泛素介導(dǎo)的蛋白水解(Ubiquitin mediated proteolysis)、糖酵解/糖異生(Glycolysis/Gluconeogenesis)、丙酮酸代謝(Pyruvate metabolism)和RNA轉(zhuǎn)運(RNA transport)通路是富集最顯著的6條通路(qvalue<0.05)。

圖6 雌性相對于雄性上調(diào)基因KEGG富集散點圖Fig.6 Scatter plot of KEGG enrichment in females relative to males of up-regulated gene

圖7 雌性相對于雄性下調(diào)基因KEGG富集散點圖Fig.7 Scatter plot of KEGG enrichment in females relative to males of down-regulated gene

2.6 實時熒光定量PCR驗證結(jié)果

qRT-PCR驗證結(jié)果顯示，上調(diào)表達的5個基因(VgR、VG2、PRSS1_2_3、TUBA、PLIN2)在轉(zhuǎn)錄組測序中的表達均是呈上調(diào)趨勢，下調(diào)表達的4個基因(gdhA、PKLR、UBE2O、CA)在轉(zhuǎn)錄組測序中的表達均是呈下調(diào)趨勢，其中g(shù)dhA和PKLR在qRT-PCR驗證中均是下調(diào)1倍(圖8)。

圖8 眼優(yōu)角蚱雌性與雄性差異表達基因qRT-PCR與RNA-Seq分析結(jié)果Fig.8 Results of qRT-PCR and RNA-Seq analysis of differential expressed genes between male and female Eucriotettix oculatus

3 結(jié)論與討論

動物的雌雄分化是進化的產(chǎn)物。昆蟲作為生物多樣性最豐富的動物類群，大多數(shù)種類在形態(tài)上都表現(xiàn)出了性別二態(tài)性，即同種昆蟲的雌、雄個體除生殖器官的結(jié)構(gòu)差異外，在大小、顏色、結(jié)構(gòu)等方面也常有明顯差異(王孟卿等,2005)。但是，目前對昆蟲雌雄性別在轉(zhuǎn)錄組層面的差異了解還不充分。

眼優(yōu)角蚱雌、雄成蟲在外部形態(tài)上也具有明顯的差異，雌蟲腹部末端具有細長的產(chǎn)卵瓣，體型較雄蟲粗大，前胸背板背面較雄性粗糙，其它構(gòu)造與體色與雄蟲相似(鄧維安等,2007)。作為昆蟲中較古老的物種，對眼優(yōu)角蚱雌性和雄性成蟲的轉(zhuǎn)錄差異開展研究，可為揭示昆蟲雌雄差異的分子機制提供一定的基礎(chǔ)數(shù)據(jù)。本研究測定了2組18個樣品的二代轉(zhuǎn)錄組，并與本實驗室已測定的三代轉(zhuǎn)錄本進行校正后，最終得到13 432條轉(zhuǎn)錄本，平均長度2 009 bp，最大長度8 297 bp，N50長度遠大于1 000 bp，說明本研究所得到的轉(zhuǎn)錄本質(zhì)量較較高。本研究僅有10%的轉(zhuǎn)錄本沒有注釋成功，主要是因為數(shù)據(jù)庫中有關(guān)直翅目昆蟲的基因信息匱乏。

本研究對眼優(yōu)角蚱雌性和雄性的轉(zhuǎn)錄組進行了基因差異表達分析，共獲得1 295個顯著上調(diào)和2 302個顯著下調(diào)基因。使用qRT-PCR方法驗證了9個基因的表達情況，雖然qRT-PCR結(jié)果與mRNA測序結(jié)果在倍數(shù)之間有稍許差異，但是每個基因有同樣的上下調(diào)趨勢，可以說明轉(zhuǎn)錄組數(shù)據(jù)分析的可靠性。在上調(diào)的基因中有與雌性發(fā)育有關(guān)的卵黃蛋白原受體VgR、卵黃蛋白原2VG2、胰蛋白酶PRSS1_2_3、圍脂滴蛋白2PLIN2等基因，在下調(diào)的基因中有與雄性發(fā)育有關(guān)的碳酸酐酶CA等基因，說明了差異表達基因與性別緊密聯(lián)系。為了解差異基因的生物學(xué)意義，通過差異基因的GO富集可以發(fā)現(xiàn)，差異基因顯著富集的GO條目大多是與細胞的形成有關(guān)，特別是雄性，其相對于雌性來說，結(jié)合(Binding)、細胞組分(Cellular component)和細胞過程(Cellular processes)等條目相關(guān)基因顯著高表達，說明雌、雄個體在成蟲后都有大量的細胞增殖活動，雄性較雌性更加劇烈，這應(yīng)該是與雄性進入成蟲期后要形成更多的精子為交配繁殖做好準(zhǔn)備有關(guān)(劉杰等,2020)。同時，基于整個基因組背景，應(yīng)用超幾何檢驗對差異基因進行了的KEGG Pathway富集分析，發(fā)現(xiàn)在雌性成蟲體內(nèi)高表達的基因顯著富集在主要的6條通路中。其中神經(jīng)活性配體受體相互作用(Neuroactive ligand-receptor interaction)通路屬于生殖相關(guān)通路，該通路的基因能影響促性腺素合成和釋放，在性腺發(fā)育過程中主要是分別作用于腦和卵巢(儲明星等,2009;Xuetal.,2015;Zhangetal.,2019)。蛋白質(zhì)消化吸收(Protein digestion and absorption)和胰腺分泌(Pancreatic secretion)通路屬于消化系統(tǒng)，說明雌性眼優(yōu)角蚱相對雄性消化能力更強，可能是由于雌性個體較大，需要消化吸收更多的營養(yǎng)物質(zhì)，或者是雌性要儲備更多營養(yǎng)物質(zhì)以備后期生殖需要。而DNA復(fù)制(DNA replication)、錯配修復(fù)(Mismatch repair)和嘧啶代謝(Pyrimidine metabolism)通路均與遺傳物質(zhì)的復(fù)制和修復(fù)有關(guān)。雄性體內(nèi)高表達的基因也是顯著富集在主要的6條通路。其中蛋白酶體(Proteasome)通路可能與雄性眼優(yōu)角蚱的生精過程有關(guān)，因為已有研究發(fā)現(xiàn)敲除蛋白酶體REGγ基因可導(dǎo)致精原干細胞數(shù)目減少，使雄性小鼠生殖能力降低(Alietal.,2013;李寧等,2016)。氮源是昆蟲生長發(fā)育和繁殖所必需的營養(yǎng)物質(zhì)之一。由于食物中氮含量與蛋白質(zhì)含量直接相關(guān)，食物中的氮能直接影響昆蟲的交配、取食速率和繁殖等(張碧堯,2019)，研究發(fā)現(xiàn)與氮代謝有關(guān)的精氨酸對精子發(fā)生發(fā)揮關(guān)鍵作用(Moralesetal.,2003;O’Flahertyetal.,2004;徐學(xué)玉等,2017)，因此雄性眼優(yōu)角蚱氮代謝通路被富集。泛素介導(dǎo)的蛋白水解(Ubiquitin mediated proteolysis)屬于信號通路，參與卵發(fā)生、精子發(fā)生和性腺成熟等一系列生殖過程的調(diào)控活動(Mengetal.,2015;Pengetal.,2015;Yangetal.,2016;董忠典等,2021)。與能量代謝相關(guān)的糖酵解(Glycolysis / Gluconeogenesis)、丙酮酸代謝(Pyruvate metabolism)通路在雄性中被顯著富集，可能是雄性眼優(yōu)角蚱在交配前較雌性活動多，需要消耗更多的能量(肖舒晴等,2019)。RNA轉(zhuǎn)運(RNA transport)通路顯著富集，可能是由于蛋白質(zhì)合成與RNA轉(zhuǎn)運密切相關(guān)，雄性需要大量合成與精子形成有關(guān)的蛋白，以完成精子儲備，這與雄性按蚊轉(zhuǎn)錄表達相似(劉杰等,2020)。綜上所述，眼優(yōu)角蚱雌性和雄性成蟲在轉(zhuǎn)錄水平的差異主要還是集中在與生殖細胞形成和生殖活動有關(guān)的基因表達方面。本研究對進一步深入解析昆蟲雌雄分化的分子調(diào)控差異提供了一定的理論基礎(chǔ)。