井趙斌，魏 琳，俞 靚，程積民

(1.西北農(nóng)林科技大學動物科技學院，陜西 楊凌712100；2.西北農(nóng)林科技大學資源與環(huán)境學院，陜西 楊凌 712100；3.中國科學院水利部水土保持研究所，陜西 楊凌 712100)

目前，在新一代高通量測序技術的推動下，部分動植物全基因組測序工作已經(jīng)完成或正在進行中，這些基于高通量測序的相關研究逐漸深入到生物學研究的每一個角落，并成為引領生物學研究的基本手段之一。國外許多研究機構和公司正憑借其基因組測序技術、測序設備和測序成本的優(yōu)勢，以期利用基因資源的唯一性，搶先獲得基因專利權，從而實現(xiàn)其“一個基因一個產(chǎn)業(yè)”的全球壟斷局面。但新一代高通量測序技術在我國的研究和應用還處于起步階段。

我國的牧草種質資源極為豐富，目前的研究主要集中在以全球變化與農(nóng)業(yè)多樣性、農(nóng)業(yè)系統(tǒng)固碳減排、草地農(nóng)業(yè)決策系統(tǒng)及林草生態(tài)等宏觀方面[1]；同時以分子生物學為基礎的微觀方面也有大量研究成果，主要包括DNA水平上的遺傳多樣性和分子標記輔助選擇育種等方面[2]。而利用新一代高通量測序技術進行牧草種質資源創(chuàng)新與開發(fā)的研究鮮見報道。

本研究對新一代高通量測序應用領域最為廣泛的轉錄組測序的基本原理、實驗流程、數(shù)據(jù)分析及其應用進行綜述，并結合中國科學院水利部水土保持研究所林草生態(tài)研究室進行的抗旱牧草轉錄組測序研究思路進行淺析，旨在對轉錄組測序在牧草基因資源發(fā)掘及牧草作物種質創(chuàng)新與品種設計育種中產(chǎn)生的可能影響進行展望，為我國牧草資源的開發(fā)提供新的研究思路。

1 轉錄組測序

自1964年首個基因完整核苷酸序列發(fā)表以來，測序技術對分子生物學的發(fā)展起到了巨大的推動作用，2005年454 Life Sciences公司推出了大規(guī)模平行測序平臺GS20后(2007年該公司被美國Roche收購)，使測序技術邁入了新一代測序技術(next generation sequencing，NGS)的新紀元。目前，新一代測序技術平臺的市場，主要為Roche(http://www.454.com)、Illumina (http://www.illumina.com)(2006年推出Solexa Genome Analyzer platform)和ABI (http://www3.appliedbiosystems.com/AB_Home/index.htm)(2007年推出SOLiD sequencing)3家公司所占領，這3家公司所用測序技術信息見表1，其測序原理和方法已有較多報道[3-5]。最近，Helicos Biosciences公司(http://www.helicosbio.com)利用其新開發(fā)的可避免對靶cDNA擴增的tSMS系統(tǒng)，進入了轉錄組測序市場。

基于新一代測序技術發(fā)展了許多研究應用領域(http://www.genomics.cn/index.php)[6-7]，如：一是基于DNA水平的，主要包括全基因組重測序、目標區(qū)域捕獲測序、外顯子捕獲測序、基因組denovo測序、SNP芯片分型；二是基于RNA水平的，主要包括轉錄組測序和數(shù)字基因表達譜的RNA-Seq，小RNA測序、降解組測序；三是基于蛋白水平和表觀遺傳學水平的，主要有用于DNA-蛋白互作的Chip-Seq，活性調控區(qū)域鑒定的DNase-Seq，拷貝數(shù)目變異的CNV-Seq，表觀標記全基因組的methyl-Seq。其中用于轉錄組測序和數(shù)字基因表達譜研究的RNA-Seq技術是迄今應用最為廣泛的領域之一[8]。

表1 目前主導市場的新一代測序技術服務商[4-5]

轉錄組是指特定組織或細胞在某一功能狀態(tài)下轉錄出來的所有RNA的總和，包括mRNA和非編碼RNA。轉錄組技術是通過新一代高通量測序技術對cDNA測序，利用統(tǒng)計相關reads數(shù)計算出不同mRNA的表達量，發(fā)現(xiàn)轉錄水平的SNP，新的mRNA等。轉錄組測序可能是迄今新一代測序技術中最復雜的應用領域之一，該技術可以從表達水平、等位基因特異性表達、RNA編輯、重要信息的融合基因轉錄子、差異剪接等方面展開相關研究[9]。

轉錄組測序技術流程主要包括樣品制備和文庫構建，DNA成簇擴增和高通量測序，具體實驗流程如圖1所示。

1.1樣品分離和文庫構建 RNA樣品分離和制備是轉錄組測序工作的第一步，根據(jù)測序平臺和使用的Priming方法對RNA的質量和濃度要求不同。如Roche GS FLX測序平臺對樣品濃度的要求是總RNA質量濃度≥400 ng/μL，總量≥15 μg，而Illumina Solexa平臺對樣品濃度的要求是總RNA質量濃度≥300 ng/μL，總量≥6 μg。在樣品總RNA中，大量rRNA(>90%)是影響RNA-Seq測序結果精確性的關鍵因素，為獲得高濃度高質量的mRNA可選擇不同的商業(yè)試劑盒去除rRNA(如Qiagen公司的Oligotex mRNA kit和Invitrogen公司的RiboMinus)。

圖1 轉錄組測序實驗流程

文庫構建質量直接影響RNA-Seq測序的成功與否。文庫構建過程因測序平臺略有差異，最基本的文庫構建過程包括下面5個主要步驟[6]：(1)使用oligo dT微珠純化mRNA及mRNA片段化處理；(2)反轉錄反應合成雙鏈cDNA；(3)雙鏈DNA末端修復及3′末端加‘A’；(4)使用特定的測序接頭連接DNA片段兩端；(5)高保真聚合酶擴增構建成功的測序文庫。為提高文庫構建的質量和效率，可選擇市場上現(xiàn)有的工作試劑盒進行文庫構建，如cDNA反轉錄試劑盒：Clontech公司的SMARTTMPCR cDNA Synthesis kit；cDNA擴增試劑盒：Clontech公司的PCR Advantage II polymerase；去除小片段的試劑盒：Invitrogen公司的PureLinkTMPCR Purification kit。為保證測序成功，文庫構建工作一般直接由測序公司完成。

1.2上機測序 新一代測序技術因其測序平臺的不同，測序原理也不相同[6]，Illumina和HeliScope公司以循環(huán)可逆終止法(cyclic reversible termination,CRT)為主，雖然Illumina和HeliScope測序平臺使用的方法相同，但是Illumina使用的是四色CRT法，而HeliScope是單色(Cye5 dye)CRT法；SOLiD平臺是應用連接法測序，同時利用了獨特的雙堿基編碼原理，Roche 454主要應用基于Pyrosequencing原理的非電泳生物熒光法。雖然各個公司推出的測序平臺原理不同，但也有許多共同的特點：(1)采用連接法測序(sequencing by ligation)或邊合成邊測序(synthesize by sequencing)；(2)將待測基因組序列打斷成較短的片段，對該片段進行測序，然后對測序的結果進行拼接得到待測基因組的序列；(3)一次測序可以得到Gb級的堿基數(shù)。

實際研究中，應根據(jù)研究材料和研究目的及各測序平臺的特點、成本，選擇合適的測序平臺和解決方案，一般對于有參比基因(即已有參考基因組序列)的物種進行轉錄組測序，采用Illumina 更加經(jīng)濟，如果沒有參比基因，選擇Roche GS FLX可減少后續(xù)數(shù)據(jù)處理中拼接的難度，可獲得更加精確的數(shù)據(jù)。華大基因是我國高通量測序的權威研究機構與服務商，其轉錄組測序平臺為Illumina新推出的HiSeqTM2000，該平臺可同時用于有(無)參考基因組的轉錄組測序，通過華大自主研發(fā)的短序列組裝軟件SOAPdenovo獲得的信息量完全可滿足研究需求(個人通訊)。

1.3測序數(shù)據(jù)的生物信息學分析 轉錄組測序中，獲得的原始數(shù)據(jù)是在測序循環(huán)中通過單個平行測序反應時，由含有堿基位信息的熒光點以圖像格式輸出的，對這些可達Gb單位級的原始圖像數(shù)據(jù)文件處理后，可得到所測基因序列中的堿基位，然后利用相關軟件進行測序數(shù)據(jù)處理可獲得大量信息供后續(xù)研究利用。

圖2是常見的轉錄組測序信息分析流程圖。具體分析步驟和內容如下(以有參考基因組轉錄組測序為例)。(1)測序數(shù)據(jù)質量評估：用測序平臺軟件(如Roche 454分析軟件：GS-FLX Software)去除低質量序列和銜接子區(qū)域并檢測可能的樣品污染，同時對得到的序列進行深度和覆蓋度計算并評價測序質量。(2)序列拼接及分析：將測序序列與GeneBank等數(shù)據(jù)庫同物種已發(fā)表的EST序列合并(有參考基因組序列)，采用denovo拼接方法進行序列拼接，并分析Conting和Unigene長度分布等(如Roche 454用：GS De Novo Assembler Software)。(3)reads比對至基因組：將拼接結果與參考基因組進行比對，挑選出Unigene map的所有reads用于peak分析，并進行reads的全基因組分布和基因定位分析。(4)新轉錄本的尋找：將能mapping至基因組上的reads，通過這些reads和基因組序列的比對，利用同源預測或者denovo預測新的轉錄本。(5)檢測反義鏈轉錄本并注釋：通過對符合轉錄本標準的reads進行大規(guī)模比對，選出潛在的反義鏈轉錄本，并對其進行功能注釋。(6)轉錄本定量：對RNA-Seq結果中的基因采用 reads per kilobase of exon per million mapped sequence reads (RPKM)進行基因定量分析以推測其表達量。(7)差異基因Annotation、Synonyms及ID轉換：對選出的差異表達基因進行功能注釋和部分數(shù)據(jù)庫的ID轉換，并對各個差異基因附注對應的部分別名。(8)聚類分析：聚類分析包括層次聚類和表達模式聚類兩種，層次聚類即按照樣品進行聚類分析，主要分別對上調基因和下調基因以及所有差異基因進行聚類分析。差異基因聚類(又稱K-Means聚類)是對各個時期的表達模式進行全面分析。(9)可變剪切預測及SSR、SNP開發(fā)。(10)基因的COG/KOG功能分類。通過對基因功能進行COG(cluster of orthologous groups of proteins) 或KOG(eukaryotic orthologous groups of proteins)分類，獲得變化基因所調節(jié)的功能及研究因子對生物功能影響的情況。(11)Gene Ontology(GO)分析：分別對上調基因、下調基因及所有差異基因進行GO分析，統(tǒng)計每個GO term中所包括的差異基因個數(shù)，用統(tǒng)計檢驗的方法計算每個GO term中差異基因富集的顯著性，并根據(jù)P-value大小判斷差異基因中具有顯著性意義的GO term。將差異基因相關的所有GO term用關系網(wǎng)絡圖(包括biological process、cellular component和molecular function)的形式表示，以了解各基因之間的相互層次關系和生物學功能。(12)Pathway分析：結合KEGG pathway等數(shù)據(jù)庫，將上調基因、下調基因及所有差異基因進行Pathway顯著性分析，統(tǒng)計每個Pathway中所包含的差異基因個數(shù)，用統(tǒng)計檢驗的方法計算出反映Pathway中差異基因分布富集顯著性的P-value，根據(jù)P-value大小判斷差異基因在生物通路中富集程度。(13)差異基因上游序列的motif分析：對差異基因進行轉錄因子結合位點的motif分析，從而推測差異基因中有哪些基因可能直接受該轉錄因子調節(jié)，并結合統(tǒng)計學方法進行檢驗。(14)差異基因的共表達網(wǎng)絡的建立：對差異基因進行共表達網(wǎng)絡構建，對了解基因之間可能存在的互作關系及為后續(xù)實驗提供指導。

圖2 轉錄組測序信息分析流程

實際信息分析內容取決于研究樣本的復雜度和后續(xù)研究的目的，測序服務商一般提供最基礎的信息分析結果，如華大基因對于有參考基因組和無參考基因組的轉錄組測序信息分析內容不盡相同(http://www.genomics.cn/index.php)，有參考基因組主要分析內容包括：基因表達差異、基因結構優(yōu)化、新轉錄本預測、基因可變剪接的鑒定、基因融合檢測和SNP分析等；無參考基因組主要包括：Unigene功能注釋、Unigene Gene Ontology分類、Unigene 代謝通路分析、Unigene Pathway富集性分析、Unigene 差異表達分析、蛋白編碼區(qū)預測(CDS)等。

轉錄組測序最終目的是獲得大量信息的后續(xù)研究利用，信息分析內容決定了后續(xù)研究的主要方向。根據(jù)現(xiàn)已發(fā)表的轉錄組測序研究論文和測序信息的生物學分析結果來看，后續(xù)可進行的研究大致分為[10]：(1)功能注釋、分類和代謝途徑分析；(2)微陣列分析；(3)分子標記開發(fā)，如SNPs、SSR和EPICs等；(4)候選基因研究(如候選基因克隆、基因功能及轉基因研究)；(5)發(fā)掘內含子或基因間隔區(qū)內新的轉錄區(qū)域；(6)反義轉錄區(qū)域分析；(7)利用已知外顯子序列或新的轉錄本外顯子鑒定可變剪接；(8)發(fā)掘新的融合基因轉錄子和其他新的RNA類型；(9)全基因組掃描。

1.4應用現(xiàn)狀和技術特點 目前，已有許多轉錄組測序的研究結果發(fā)表[10-25]，如Gorodkin等[25]對豬35個不同組織和3個不同發(fā)育階段的EST進行轉錄組分析，結果表明，在大腦和睪丸組織中特異性的表達基因數(shù)量最多；在不同組織之間，基因表達量也不同，通過對EST的組裝，最終得到48 000個contigs和73 000個singletons。我國學者也發(fā)表了一些最新研究成果，其中Zhang等[26]用8種不同水稻(Oryzasativa)樣品的不同組織不同時期混合建庫，通過轉錄組技術分析了栽培稻的第1張轉錄組圖譜，結果在水稻8種組織樣品中檢測到大約27 000個基因的表達和38 000個轉錄單元，證實了約9 000個基因發(fā)生可變剪接，同時鑒定出了234個由反式剪接產(chǎn)生的轉錄融合基因，表明融合基因比預期的更為普遍；Wu等[27]利用葡萄(Vitisvinifera)接種霜霉病后采集的4～8 d葉片混合樣，通過Solexa技術測序獲得了15 249個候選差異表達基因。

已有研究表明[9,13-16,21,27]：與以雜交方法為基礎的微陣列平臺測序技術、用于基因表達Gap分析和Polony多元分析的SAGE、CAGE和PMAGE技術以及基于序列標簽的測序技術相比，以RNA-Seq技術為基礎的轉錄組測序具有更多優(yōu)點：(1)對已知轉錄本的檢測不受限制，并可進一步對變異體形式鑒定、描述和量化研究；(2)可進行正確的基因功能注釋，定義單個核苷酸的基因轉錄邊界和單核苷酸多態(tài)性的表達水平；(3)“背景信號”檢出率低于微陣列技術，且許多在表達水平之上轉錄本的動態(tài)變化可被檢測到；(4)數(shù)據(jù)具有很高的重復性和再現(xiàn)性。

2 研究意義

黃土高原地區(qū)水土流失、荒漠化及天然草地植被退化現(xiàn)象嚴重，恢復和重建植被是保持水土、改善該地區(qū)生態(tài)環(huán)境和發(fā)展畜牧業(yè)的根本措施，而干旱是制約該地區(qū)發(fā)展栽培草地的主要限制性因素。因此，利用現(xiàn)有抗旱性強的優(yōu)質牧草種質資源，通過植物抗旱生理與基因組研究相結合的方法，對于闡明特定條件下植物起關鍵作用的耐旱機制，發(fā)掘耐旱主效應基因或為多基因有效聚合提供依據(jù)，以選育高效抗旱轉基因牧草或作物新品種，為解決當前黃土高原水資源缺乏、生態(tài)環(huán)境惡化和廣大旱區(qū)作物高產(chǎn)具有重要的理論和實踐意義。

我國抗旱性牧草資源豐富，目前已有大量有關于牧草抗旱性的研究結果發(fā)表。本研究在科技查新和查閱國內外關于抗旱研究選擇的材料及本課題前期工作的基礎上，選擇國內外已有報道表明抗旱性較強和未見報道的優(yōu)良抗旱性牧草，對其在干旱和復水條件下的抗旱生理和光合等指標進行測定，選出強抗旱品種，利用轉錄組測序和數(shù)字基因表達譜在轉錄水平上研究基因組序列未知的牧草在干旱脅迫處理下的基因表達和差異表達基因篩選，旨在為轉基因植物育種建立一個豐富的抗旱基因資源庫。

3 展望

新一代測序技術的誕生對分子生物學的深入研究發(fā)揮了巨大的促進作用，以新一代測序技術為基礎的轉錄組測序和全基因組測序相比，成本很低，數(shù)據(jù)量大，且不易受遺傳背景限制，可構建豐富的表達基因數(shù)據(jù)庫，為進一步研究提供重要基礎和依據(jù)。在以往的抗逆性研究中，選擇的研究對象多集中在農(nóng)作物中，將抗旱轉基因植物研究的重點放在林草植物上更為可行，因為這方面的抗逆基因資源更為豐富，而且與一年生農(nóng)作物相比，這類植物存活需求是第1位的，產(chǎn)量高低是第2位的，生態(tài)效益在先，只要生存下來，就有機會實現(xiàn)其生態(tài)經(jīng)濟目標[28]。我國抗逆性牧草資源豐富，利用轉錄組測序進行各種抗逆性牧草資源開發(fā)建立大規(guī)模的共享基因組數(shù)據(jù)庫，除抗逆性牧草以外，可同時進行能源草等基因組開發(fā)研究，為轉基因植物育種提供豐富的基因資源將起到巨大的推動作用。

[1]張自和.西部大開發(fā)與草業(yè)發(fā)展[J].草業(yè)科學,2010,27(8):1-4.

[2]黃文達,趙學勇,趙昕，等.分子標記在種群遺傳學研究中的應用[J].草業(yè)科學，2010,27(11)：115-120.

[3]Mardis E R.Next-generation DNA sequencing methods[J].Annual Review of Genomics and Human Genetics,2008,9:387-402.

[4]Shendure J,Ji H.Next-generation DNA sequencing[J].Nature Biotechnology,2008,11:1135-1145.

[5]Metzker M L.Sequencing technologies the next generation[J].Nature Reviews Genetics,2010,11:31-36.

[6]Valerio C,Claudia A,Italia D F,etal.Uncovering the complexity of transcriptomes with RNA-Seq[J].Journal of Biomedicine and Biotechnology,2010，10:1-20.

[7]楊曉玲,施蘇華,唐恬.新一代測序技術的發(fā)展及應用前景[J].生物技術通報,2010,10:76-81.

[8]Wilhelm B T,Landry J R.RNA-Seq quantitative measurement of expression through massively parallel RNA Sequencing[J].Methods,2009,48:249-257.

[9]Velculescu V E,Zhang L,Zhou W,etal.Characterization of the yeast transcriptome[J].Cell,1997,88:243-251.

[10]Samuel M,Jurg B.RNA-seq:from technology to biology[J].Cellular and Molecular Life Sciences,2010,67:569-579.

[11]He Y,Vogelstein B,Velculescu V E,etal.The antisense transcriptomes of human cells[J].Science,2008,322:1855-1857.

[12]Maher C A,Kumar S C,Cao X,etal.Transcriptome sequencing to detect gene fusions in cancer[J].Nature,2009,458:97-101.

[13]Wang E T,Sandberg R,Luo S,etal.Alternative isoform regulation in human tissue transcriptomes[J].Nature,2008,456:470-476.

[14]Pan Q,Shai O,Lee L J,etal.Deep surveying of alternative splicing complexity in the human transcriptome by high-throughput sequencing[J].Nature Genetics,2008,40:1413-1415.

[15]Ingolia N T,Ghaemmaghami S,Newman J R S,etal.Genome-wide analysis in vivo of translation with nucleotide resolution using ribosome profiling[J].Science,2009,324:218-223.

[16]Hahn D A,Ragland G J,Shoemaker D D,etal.Gene discovery using massively parallel pyrosequencing to develop ESTs for the flesh fly Sarcophaga crassipalpis[J].BMC Genomics,2009,10:234-241.

[17]Vera J C,Wheat C W,Fescemyer H W,etal.Rapid transcriptome characterization for a nonmodel organism using 454 pyrosequencing[J].Molecular Ecology,2008,17:1636-1647.

[18]Lister R,O’Malley R C,Tonti F J.Highly integrated single-base resolution maps of the epigenome inArabidopsis[J].Cell,2008,133:523-536.

[19]Emrich S J,Barbazuk W B,Li L,etal.Gene discovery and annotation using LCM-454 transcriptome sequencing[J].Genome Research,2008,17:69-73.

[20]Barbazuk W B,Emrich S J,Chen H D,etal.SNP discovery via 454 transcriptome sequencing[J].Plant,2007,51:910-918.

[21]Trick M,Long Y,Meng J,etal.Single nucleotide polymorphism (SNP) discovery in the polyploidBrassicanapususing Solexa transcriptome sequencing[J].Plant Biotechnology,2009,7:334-346.

[22]Mortazavi A,Williams B A,McCue K,etal.Mapping and quantifying transcriptomes by RNA-Seq[J].Nature Methods,2008,5:621-628.

[23]Cloonan N.Stem cell transcriptome profiling via massive-scale mRNA sequencing[J].Nature Methods,2008,5:613-619.

[24]Marioni J C,Mason C E,Mane S M,etal.RNA-seq:an assessment of technical reproducibility and comparison with gene expression arrays[J].Genome Research,2008,18:1509-1517.

[25]Gorodkin J,Cirera S,Hedegaard J,etal.Porcine transcriptome analysis based on 97 non-normalized cDNA libraries and assembly of 1,021,891 expressed sequence tags[J].Genome Biology,2007,8:R45.1-R45.16.

[26]Zhang G J,Guo G W,Hu X D,etal.Deep RNA sequencing at single base-pair resolution reveals high complexity of the rice transcriptome[J].Genome Research，2010,20(5):646-654.

[27]Wu J,Zhang Y L,Zhang H Q,etal.Whole genome wide expression profiles ofVitisamurensisgrape responding to downy mildew by using Solexa sequencing technology[J].BMC Plant Biology,2010,10:234.

[28]山侖.轉基因作物育種當前處于怎樣的發(fā)展階段?[N].科學時報,2010-03-26(第1版).