转录组测序及其在牧草基因资源发掘中的应用前景

井赵斌，魏 琳，俞 靓，程积民

(1.西北农林科技大学动物科技学院，陕西 杨凌712100；2.西北农林科技大学资源与环境学院，陕西 杨凌 712100；3.中国科学院水利部水土保持研究所，陕西 杨凌 712100)

目前，在新一代高通量测序技术的推动下，部分动植物全基因组测序工作已经完成或正在进行中，这些基于高通量测序的相关研究逐渐深入到生物学研究的每一个角落，并成为引领生物学研究的基本手段之一。国外许多研究机构和公司正凭借其基因组测序技术、测序设备和测序成本的优势，以期利用基因资源的唯一性，抢先获得基因专利权，从而实现其“一个基因一个产业”的全球垄断局面。但新一代高通量测序技术在我国的研究和应用还处于起步阶段。

我国的牧草种质资源极为丰富，目前的研究主要集中在以全球变化与农业多样性、农业系统固碳减排、草地农业决策系统及林草生态等宏观方面[1]；同时以分子生物学为基础的微观方面也有大量研究成果，主要包括DNA水平上的遗传多样性和分子标记辅助选择育种等方面[2]。而利用新一代高通量测序技术进行牧草种质资源创新与开发的研究鲜见报道。

本研究对新一代高通量测序应用领域最为广泛的转录组测序的基本原理、实验流程、数据分析及其应用进行综述，并结合中国科学院水利部水土保持研究所林草生态研究室进行的抗旱牧草转录组测序研究思路进行浅析，旨在对转录组测序在牧草基因资源发掘及牧草作物种质创新与品种设计育种中产生的可能影响进行展望，为我国牧草资源的开发提供新的研究思路。

1 转录组测序

自1964年首个基因完整核苷酸序列发表以来，测序技术对分子生物学的发展起到了巨大的推动作用，2005年454 Life Sciences公司推出了大规模平行测序平台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)利用其新开发的可避免对靶cDNA扩增的tSMS系统，进入了转录组测序市场。

基于新一代测序技术发展了许多研究应用领域(http://www.genomics.cn/index.php)[6-7]，如：一是基于DNA水平的，主要包括全基因组重测序、目标区域捕获测序、外显子捕获测序、基因组denovo测序、SNP芯片分型；二是基于RNA水平的，主要包括转录组测序和数字基因表达谱的RNA-Seq，小RNA测序、降解组测序；三是基于蛋白水平和表观遗传学水平的，主要有用于DNA-蛋白互作的Chip-Seq，活性调控区域鉴定的DNase-Seq，拷贝数目变异的CNV-Seq，表观标记全基因组的methyl-Seq。其中用于转录组测序和数字基因表达谱研究的RNA-Seq技术是迄今应用最为广泛的领域之一[8]。

表1 目前主导市场的新一代测序技术服务商[4-5]

转录组是指特定组织或细胞在某一功能状态下转录出来的所有RNA的总和，包括mRNA和非编码RNA。转录组技术是通过新一代高通量测序技术对cDNA测序，利用统计相关reads数计算出不同mRNA的表达量，发现转录水平的SNP，新的mRNA等。转录组测序可能是迄今新一代测序技术中最复杂的应用领域之一，该技术可以从表达水平、等位基因特异性表达、RNA编辑、重要信息的融合基因转录子、差异剪接等方面展开相关研究[9]。

转录组测序技术流程主要包括样品制备和文库构建，DNA成簇扩增和高通量测序，具体实验流程如图1所示。

1.1样品分离和文库构建 RNA样品分离和制备是转录组测序工作的第一步，根据测序平台和使用的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可选择不同的商业试剂盒去除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)高保真聚合酶扩增构建成功的测序文库。为提高文库构建的质量和效率，可选择市场上现有的工作试剂盒进行文库构建，如cDNA反转录试剂盒：Clontech公司的SMARTTMPCR cDNA Synthesis kit；cDNA扩增试剂盒：Clontech公司的PCR Advantage II polymerase；去除小片段的试剂盒：Invitrogen公司的PureLinkTMPCR Purification kit。为保证测序成功，文库构建工作一般直接由测序公司完成。

1.2上机测序 新一代测序技术因其测序平台的不同，测序原理也不相同[6]，Illumina和HeliScope公司以循环可逆终止法(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级的碱基数。

实际研究中，应根据研究材料和研究目的及各测序平台的特点、成本，选择合适的测序平台和解决方案，一般对于有参比基因(即已有参考基因组序列)的物种进行转录组测序，采用Illumina 更加经济，如果没有参比基因，选择Roche GS FLX可减少后续数据处理中拼接的难度，可获得更加精确的数据。华大基因是我国高通量测序的权威研究机构与服务商，其转录组测序平台为Illumina新推出的HiSeqTM2000，该平台可同时用于有(无)参考基因组的转录组测序，通过华大自主研发的短序列组装软件SOAPdenovo获得的信息量完全可满足研究需求(个人通讯)。

1.3测序数据的生物信息学分析 转录组测序中，获得的原始数据是在测序循环中通过单个平行测序反应时，由含有碱基位信息的荧光点以图像格式输出的，对这些可达Gb单位级的原始图像数据文件处理后，可得到所测基因序列中的碱基位，然后利用相关软件进行测序数据处理可获得大量信息供后续研究利用。

图2是常见的转录组测序信息分析流程图。具体分析步骤和内容如下(以有参考基因组转录组测序为例)。(1)测序数据质量评估：用测序平台软件(如Roche 454分析软件：GS-FLX Software)去除低质量序列和衔接子区域并检测可能的样品污染，同时对得到的序列进行深度和覆盖度计算并评价测序质量。(2)序列拼接及分析：将测序序列与GeneBank等数据库同物种已发表的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进行大规模比对，选出潜在的反义链转录本，并对其进行功能注释。(6)转录本定量：对RNA-Seq结果中的基因采用 reads per kilobase of exon per million mapped sequence reads (RPKM)进行基因定量分析以推测其表达量。(7)差异基因Annotation、Synonyms及ID转换：对选出的差异表达基因进行功能注释和部分数据库的ID转换，并对各个差异基因附注对应的部分别名。(8)聚类分析：聚类分析包括层次聚类和表达模式聚类两种，层次聚类即按照样品进行聚类分析，主要分别对上调基因和下调基因以及所有差异基因进行聚类分析。差异基因聚类(又称K-Means聚类)是对各个时期的表达模式进行全面分析。(9)可变剪切预测及SSR、SNP开发。(10)基因的COG/KOG功能分类。通过对基因功能进行COG(cluster of orthologous groups of proteins) 或KOG(eukaryotic orthologous groups of proteins)分类，获得变化基因所调节的功能及研究因子对生物功能影响的情况。(11)Gene Ontology(GO)分析：分别对上调基因、下调基因及所有差异基因进行GO分析，统计每个GO term中所包括的差异基因个数，用统计检验的方法计算每个GO term中差异基因富集的显著性，并根据P-value大小判断差异基因中具有显著性意义的GO term。将差异基因相关的所有GO term用关系网络图(包括biological process、cellular component和molecular function)的形式表示，以了解各基因之间的相互层次关系和生物学功能。(12)Pathway分析：结合KEGG pathway等数据库，将上调基因、下调基因及所有差异基因进行Pathway显著性分析，统计每个Pathway中所包含的差异基因个数，用统计检验的方法计算出反映Pathway中差异基因分布富集显著性的P-value，根据P-value大小判断差异基因在生物通路中富集程度。(13)差异基因上游序列的motif分析：对差异基因进行转录因子结合位点的motif分析，从而推测差异基因中有哪些基因可能直接受该转录因子调节，并结合统计学方法进行检验。(14)差异基因的共表达网络的建立：对差异基因进行共表达网络构建，对了解基因之间可能存在的互作关系及为后续实验提供指导。

图2 转录组测序信息分析流程

实际信息分析内容取决于研究样本的复杂度和后续研究的目的，测序服务商一般提供最基础的信息分析结果，如华大基因对于有参考基因组和无参考基因组的转录组测序信息分析内容不尽相同(http://www.genomics.cn/index.php)，有参考基因组主要分析内容包括：基因表达差异、基因结构优化、新转录本预测、基因可变剪接的鉴定、基因融合检测和SNP分析等；无参考基因组主要包括：Unigene功能注释、Unigene Gene Ontology分类、Unigene 代谢通路分析、Unigene Pathway富集性分析、Unigene 差异表达分析、蛋白编码区预测(CDS)等。

转录组测序最终目的是获得大量信息的后续研究利用，信息分析内容决定了后续研究的主要方向。根据现已发表的转录组测序研究论文和测序信息的生物学分析结果来看，后续可进行的研究大致分为[10]：(1)功能注释、分类和代谢途径分析；(2)微阵列分析；(3)分子标记开发，如SNPs、SSR和EPICs等；(4)候选基因研究(如候选基因克隆、基因功能及转基因研究)；(5)发掘内含子或基因间隔区内新的转录区域；(6)反义转录区域分析；(7)利用已知外显子序列或新的转录本外显子鉴定可变剪接；(8)发掘新的融合基因转录子和其他新的RNA类型；(9)全基因组扫描。

1.4应用现状和技术特点 目前，已有许多转录组测序的研究结果发表[10-25]，如Gorodkin等[25]对猪35个不同组织和3个不同发育阶段的EST进行转录组分析，结果表明，在大脑和睾丸组织中特异性的表达基因数量最多；在不同组织之间，基因表达量也不同，通过对EST的组装，最终得到48 000个contigs和73 000个singletons。我国学者也发表了一些最新研究成果，其中Zhang等[26]用8种不同水稻(Oryzasativa)样品的不同组织不同时期混合建库，通过转录组技术分析了栽培稻的第1张转录组图谱，结果在水稻8种组织样品中检测到大约27 000个基因的表达和38 000个转录单元，证实了约9 000个基因发生可变剪接，同时鉴定出了234个由反式剪接产生的转录融合基因，表明融合基因比预期的更为普遍；Wu等[27]利用葡萄(Vitisvinifera)接种霜霉病后采集的4～8 d叶片混合样，通过Solexa技术测序获得了15 249个候选差异表达基因。

已有研究表明[9,13-16,21,27]：与以杂交方法为基础的微阵列平台测序技术、用于基因表达Gap分析和Polony多元分析的SAGE、CAGE和PMAGE技术以及基于序列标签的测序技术相比，以RNA-Seq技术为基础的转录组测序具有更多优点：(1)对已知转录本的检测不受限制，并可进一步对变异体形式鉴定、描述和量化研究；(2)可进行正确的基因功能注释，定义单个核苷酸的基因转录边界和单核苷酸多态性的表达水平；(3)“背景信号”检出率低于微阵列技术，且许多在表达水平之上转录本的动态变化可被检测到；(4)数据具有很高的重复性和再现性。

2 研究意义

黄土高原地区水土流失、荒漠化及天然草地植被退化现象严重，恢复和重建植被是保持水土、改善该地区生态环境和发展畜牧业的根本措施，而干旱是制约该地区发展栽培草地的主要限制性因素。因此，利用现有抗旱性强的优质牧草种质资源，通过植物抗旱生理与基因组研究相结合的方法，对于阐明特定条件下植物起关键作用的耐旱机制，发掘耐旱主效应基因或为多基因有效聚合提供依据，以选育高效抗旱转基因牧草或作物新品种，为解决当前黄土高原水资源缺乏、生态环境恶化和广大旱区作物高产具有重要的理论和实践意义。

我国抗旱性牧草资源丰富，目前已有大量有关于牧草抗旱性的研究结果发表。本研究在科技查新和查阅国内外关于抗旱研究选择的材料及本课题前期工作的基础上，选择国内外已有报道表明抗旱性较强和未见报道的优良抗旱性牧草，对其在干旱和复水条件下的抗旱生理和光合等指标进行测定，选出强抗旱品种，利用转录组测序和数字基因表达谱在转录水平上研究基因组序列未知的牧草在干旱胁迫处理下的基因表达和差异表达基因筛选，旨在为转基因植物育种建立一个丰富的抗旱基因资源库。

3 展望

新一代测序技术的诞生对分子生物学的深入研究发挥了巨大的促进作用，以新一代测序技术为基础的转录组测序和全基因组测序相比，成本很低，数据量大，且不易受遗传背景限制，可构建丰富的表达基因数据库，为进一步研究提供重要基础和依据。在以往的抗逆性研究中，选择的研究对象多集中在农作物中，将抗旱转基因植物研究的重点放在林草植物上更为可行，因为这方面的抗逆基因资源更为丰富，而且与一年生农作物相比，这类植物存活需求是第1位的，产量高低是第2位的，生态效益在先，只要生存下来，就有机会实现其生态经济目标[28]。我国抗逆性牧草资源丰富，利用转录组测序进行各种抗逆性牧草资源开发建立大规模的共享基因组数据库，除抗逆性牧草以外，可同时进行能源草等基因组开发研究，为转基因植物育种提供丰富的基因资源将起到巨大的推动作用。

[1]张自和.西部大开发与草业发展[J].草业科学,2010,27(8):1-4.

[2]黄文达,赵学勇,赵昕，等.分子标记在种群遗传学研究中的应用[J].草业科学，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]杨晓玲,施苏华,唐恬.新一代测序技术的发展及应用前景[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]山仑.转基因作物育种当前处于怎样的发展阶段?[N].科学时报,2010-03-26(第1版).