基于RNA-seq技術的甘藍型油菜新轉錄本的鑒定

徐 文,安素妨,王 艷,賈琳琳,魯?shù)さ?張瑩瑩,劉建豐,李保全*

(1.河南省農業(yè)科學院 作物設計中心，河南 鄭州 450002; 2.河南省農業(yè)科學院 農業(yè)經濟與信息研究所，河南 鄭州 450002)

轉錄組是特定時間特定組織內由基因轉錄產生的全部RNA轉錄本，包括編碼蛋白質的RNA和非編碼RNA。轉錄組測序(RNA-seq)是借助于近年來發(fā)展迅速的高通量測序技術對特定細胞在某一功能狀態(tài)下幾乎所有轉錄出來的RNA進行測序的一種技術方法[1]。隨著二代測序技術通量高、成本低、準確度和靈敏度高等優(yōu)勢，RNA-seq近年來成為揭示基因轉錄和表達調控規(guī)律的有效手段，在植物發(fā)育、抗病抗逆、分子標記開發(fā)和系統(tǒng)發(fā)生等方面得到了廣泛的應用[2-6]。同時RNA-seq在基因結構變異、低豐度轉錄本和未知轉錄本的檢測方面也發(fā)揮了重要的作用[7-13]。

甘藍型油菜(Brasscianapus)是世界四大油料作物之一，是我國重要的食用油來源和植物蛋白飼料資源之一，在國民經濟和人們日常生活中占有重要地位。甘藍型油菜基因組(AACC) 是白菜(B.rap)基因組(AA)與甘藍(B.oleracea)基因組(CC)雜交后經過染色體加倍形成的[14]，其基因組參考序列于2014年公開發(fā)布，為油菜功能基因組學的研究奠定了良好的基礎[15]。測序發(fā)現(xiàn)，油菜基因組大小為849.7 Mb，由101 040個基因組成，但是由于其基因組極其復雜和注釋方法的局限性，油菜基因組還存在基因組組裝不完整、注釋信息不完善、新轉錄本的遺漏等問題。目前，利用RNA-seq研究甘藍型油菜發(fā)育、抗逆等的報道較多[16-18]，但是甘藍型油菜中新轉錄本的鑒定尚缺乏深入的研究。鑒于RNA-seq技術的優(yōu)點和甘藍型油菜基因組注釋的不完整性，利用RNA-seq和生物信息分析技術預測甘藍型油菜基因組中尚未注釋的轉錄本，并通過RT-PCR(Reverse transcription-PCR)和克隆測序的方法來驗證，以期進一步完善油菜基因組和轉錄組注釋信息，為進一步挖掘油菜基因組中優(yōu)良的功能基因提供數(shù)據基礎。

1 材料和方法

1.1 試驗材料

本研究以甘藍型油菜W系為試驗材料。選取試驗基地生長3周的油菜幼苗，取葉片于液氮中,置于-80 ℃ 冰箱中保存。

1.2 試驗方法

1.2.1 RNA提取及檢測 油菜葉片RNA的提取采用Trizol法，按照說明書嚴格操作。利用Nanodrop檢測RNA的純度，Agilent 2100檢測RNA的完整性。

1.2.2 cDNA文庫構建和測序 樣品RNA檢測合格后，進行cDNA文庫構建，主要流程如下：用帶有Oligo(dT)的beads富集mRNA；加入fragmentation buffer將mRNA打斷成短片段；以mRNA為模板合成第一鏈cDNA，后加入buffer、dNTPs和DNA Polymerase Ⅰ合成第二鏈cDNA；隨后利用AMPure XP beads純化雙鏈cDNA；純化的雙鏈cDNA再進行末端修復、加A尾并連接測序接頭；然后用AMPure XP beads進行片段大小選擇；最后通過PCR擴增富集得到最終的cDNA文庫。

應用Illumina HiSeq 2500 高通量測序平臺對構建好的文庫進行測序，測序讀長為雙端125 bp。

1.2.3 測序數(shù)據的質控及比對 對測序得到的raw data，去除含有adapter的reads、去除N堿基比例大于5%的reads、去除低質量的reads得到clean data。從http://plants.ensembl.org/index.html下載已公布的甘藍型油菜參考基因組序列，用Bowtie 2軟件對油菜參考基因組序列建立索引[19]，用Tophat 2軟件將得到的clean data與參考基因組序列進行mapping[20],mapped reads以bam格式文件輸出。

1.2.4 轉錄本的表達豐度計算和新轉錄本的鑒定 轉錄本的豐度通過計算FPKM(Fragments per kilobase per million mapped reads)值來度量。通過Cufflinks軟件包中的cuffdiff命令來實現(xiàn)對每個樣品表達豐度的計算[21]。

將每個樣品比對得到的bam文件，借助Cufflinks軟件與甘藍型油菜參考基因組的注釋文件進行位置信息的比較和整合，初步建立每個樣品的轉錄本文庫；通過cuffmerge命令的合并功能對得到的每個樣品的轉錄本文庫進行整合，得到1個完整的基因組注釋文件merged.gtf；然后用cuffcompare命令與甘藍型油菜已知的轉錄本信息進行比較，從而鑒定候選的新轉錄本。

1.2.5 新轉錄本的RT-PCR和測序驗證 從鑒定的候選轉錄本中選取一部分，以Primer Premier 5.0根據其序列設計引物(表1)；提取甘藍型油菜葉片RNA，反轉錄成cDNA，進行PCR擴增，將擴增產物用瓊脂糖電泳檢測，然后純化回收，連接T載體，轉化DH5α感受態(tài)細胞，挑選陽性克隆，送公司測序。

表1 甘藍型油菜新轉錄本及其PCR擴增引物

續(xù)表1 甘藍型油菜新轉錄本及其PCR擴增引物

2 結果與分析

2.1 甘藍型油菜RNA質量的檢測

經檢測，17個甘藍型油菜葉片RNA均符合cDNA文庫構建的標準：OD260/OD280為1.8～2.2、28S rRNA/18S rRNA≥1.5、OD260/OD230≥1.96、RIN值≥9.7，表明RNA純度和完整性較好。

2.2 甘藍型油菜RNA測序數(shù)據處理和分析

經過質量控制后總計得到848 866 766條clean reads，與甘藍型油菜參考基因組序列進行比對。從表2可以看出，各樣品的reads與其參考基因組序列的比對率在65.00%～85.40%，說明測序數(shù)據的比對率正常。

表2 clean data與參考基因組序列比對結果

2.3 甘藍型油菜新轉錄本的識別

通過Tophat/Cufflinks一系列流程的分析，將比對成功的序列進行組裝整合。然后與已知甘藍型油菜參考基因組注釋轉錄本信息進行比較，最終獲得了由612 085個外顯子與467 743個內含子組成的103 310個基因位點(包含158 004個mRNA)，其中137 756個轉錄本是有多個外顯子組成的。其中，甘藍型油菜已知注釋的101 040個基因位點全部包含在內。新鑒定的外顯子有33 811個，新鑒定的內含子有26 839個。

對于鑒定到的158 004個轉錄本，可劃分為6類(表3):與內含子鏈匹配的轉錄本; 潛在的新轉錄本; 與已知外顯子重疊的轉錄本;內含子與反義鏈上已知內含子重疊的轉錄本; 未知的基因間隔區(qū)轉錄本;外顯子與反義鏈上已知轉錄本重疊的轉錄本。 其中，未獲得注釋的7 720個轉錄本被劃分為未知的基因間隔區(qū)轉錄本。

表3 新鑒定的甘藍型油菜轉錄本和已知轉錄本的比較結果

2.4 甘藍型油菜新轉錄本的鑒定及驗證

根據轉RNA-seq預測的轉錄本單元及其表達豐度，本研究選取了18個平均表達豐度在100以上的新鑒定的轉錄本進行RT-PCR擴增并測序驗證，結果發(fā)現(xiàn)15個新轉錄本擴增條帶很清晰(圖1)。

為了驗證結果的準確性和可靠性，將15個轉錄本PCR擴增產物進行克隆,挑選轉化子送公司測序，結果發(fā)現(xiàn)，測序得到的序列和轉錄組數(shù)據分析組裝出來的序列是一致的(圖2)。

M：Marker；1：XLOC_028977； 2：XLOC_056607；3：XLOC_009511；4：XLOC_042556；5：XLOC_097486；6：XLOC_095366；7:XLOC_097244；8：XLOC_017634；9：XLOC_015102；10：XLOC_066937；11：XLOC_092505；12：XLOC_073270；13：XLOC_046341；14：XLOC_100972；15：XLOC_041903

圖1新轉錄本RT-PCR擴增電泳檢測結果

Sequencing表示PCR產物測序序列，XLOC_095366表示通過轉錄組分析預測的轉錄本序列圖2 XLOC_095366轉錄本比對結果

3 結論與討論

本研究利用RNA-seq和生物信息學技術對甘藍型油菜葉片轉錄組文庫中的新轉錄本進行了挖掘和鑒定，并采用RT-PCR和克隆測序對鑒定到的新轉錄本進行了驗證。研究結果表明，在甘藍型油菜基因組已知基因的間隔區(qū)發(fā)現(xiàn)了7 720個新轉錄本，通過RT-PCR擴增驗證了15個新鑒定的轉錄本，為油菜基因組提供了補充性的注釋信息，為進一步挖掘甘藍型油菜基因組中優(yōu)良的功能基因提供數(shù)據基礎。

同一個基因通過可變剪接后形成多種mRNA成熟體即不同的轉錄本,經翻譯產生結構和功能不同的蛋白質或者具有調控功能的非編碼RNA。從PCR擴增檢測的電泳圖中可以看出，新鑒定的轉錄本有部分擴增的條帶不是單一的，說明這些轉錄本可能存在不同的可變剪接體。那么需要進一步深入研究這些轉錄本的特點，區(qū)分每個轉錄本是編碼蛋白質的基因還是具有調控功能的非編碼RNA。對于具有編碼功能的轉錄本，要深入研究其功能并挖掘可能存在的不同剪接體。對于在植物生長發(fā)育、逆境脅迫等方面發(fā)揮著重要作用的非編碼RNA，要克隆這些非編碼RNA并研究它們對功能基因的調控機制。因此，深入挖掘和研究這些新轉錄本及其作用機制對完善甘藍型油菜基因組注釋、為分子育種提供豐富的基因信息具有重要的理論意義。

參考文獻:

[1] Wang Z,Gerstein M,Snyder M.RNA-Seq:A revolutionary tool for transcriptomics[J].Nature Reviews Genetics,2009,10(1):57-63.

[2] Schlueter U,Denton A K,Braeutigam A.Understanding metabolite transport and metabolism in C4 plants through RNA-seq[J].Current Opinion in Plant Biology,2016,31:83-90.

[3] Martin L B,Fei Z,Giovannoni J J,etal.Catalyzing plant science research with RNA-seq[J].Frontiers in Plant Science,2013,4:66.

[4] Kakumanu A,Ambavaram M M,Klumas C A,etal.Effects of drought on gene expression in maize reproductive and leaf meristem tissue revealed by RNA-Seq[J].Plant Physiology,2012,160(2):846-867.

[5] Wei L,Li S,Liu S,etal.Transcriptome analysis ofHouttuyniacordataThunb.by Illumina paired-end RNA sequencing and SSR marker discovery[J].PLoS One,2014,9(1):e84105.

[6] Wei Z,Sun Z,Cui B,etal.Transcriptome analysis of colored calla lily (ZantedeschiarehmanniiEngl.) by Illumina sequencing:De novo assembly,annotation and EST-SSR marker development[J].Peer J,2016,4:e2378.

[7] Chettoor A M,Givan S A,Cole R A,etal.Discovery of novel transcripts and gametophytic functions via RNA-seq analysis of maize gametophytic transcriptomes[J].Genome Biology,2014,15:414.

[8] Alasoo K,Martinez F O,Hale C,etal.Transcriptional profiling of macrophages derived from monocytes and iPS cells identifies a conserved response to LPS and novel alternative transcription[J].Scientific Reports,2015,5:12524.

[9] Du Z Q,Eisley C J,Onteru S K,etal.Identification of species-specific novel transcripts in pig reproductive tissues using RNA-seq[J].Animal Genetics,2014,45(2):198-204.

[10] Roberts A,Pimentel H,Trapnell C,etal.Identification of novel transcripts in annotated genomes using RNA-Seq[J].Bioinformatics,2011,27(17):2325-2329.

[11] Weirick T,Militello G,Mueller R,etal.The identification and characterization of novel transcripts from RNA-seq data[J].Briefings in Bioinformatics,2016,17(4):678-685.

[12] Zhao J,Song X,Wang K.lncScore:Alignment-free identification of long noncoding RNA from assembled novel transcripts[J].Scientific Reports,2016,6:34838.

[13] Weikard R,Hadlich F,Kuehn C.Identification of novel transcripts and noncoding RNAs in bovine skin by deep next generation sequencing[J].BMC Genomics,2013,14:789.

[14] Allender C J,King G J.Origins of the amphiploid speciesBrassicanapusL.investigated by chloroplast and nuclear molecular markers[J].BMC Plant Biology,2010,10:54.

[15] Chalhoub B,Denoeud F,Liu S,etal.Early allopolyploid evolution in the post-neolithicBrassicanapusoilseed genome[J].Science,2014,345:950-953.

[16] Yan X,Dong C,Yu J,etal.Transcriptome profile analysis of young floral buds of fertile and sterile plants from the self-pollinated offspring of the hybrid between novel restorer line NR1 and Nsa CMS line inBrassicanapus[J].BMC Genomics,2013,14:26.

[17] Gill R A,Ali B,Cui P,etal.Comparative transcriptome profiling of twoBrassicanapuscultivars under chromium toxicity and its alleviation by reduced glutathione[J].BMC Genomics,2016,17:885.

[18] Lu K,Peng L,Zhang C,etal.Genome-wide association and transcriptome analyses reveal candidate genes underlying yield-determining traits inBrassicanapus[J].Frontiers in Plant Science,2017,8:206.

[19] Langmead B,Salzberg S L.Fast gapped-read alignment with Bowtie 2[J].Nature Methods,2012,9(4):357-359.

[20] Kim D,Pertea G,Trapnell C,etal.TopHat2:Accurate alignment of transcriptomes in the presence of insertions,deletions and gene fusions[J].Genome Biology,2013,14(4):R36.

[21] Trapnell C,Roberts A,Goff L,etal.Differential gene and transcript expression analysis of RNA-seq experiments with TopHat and Cufflinks[J].Nature Protocols,2012,7(3):562-578.