高通量测序的转录组学研究范文

《生物技术通报杂志》2014年第七期

1材料与方法

1.1材料本研究中辽东栎的芽、花、叶与果（5个阶段：开花后20、40、60、80和100d均采集）均采集于北京市门头沟区小龙门国家森林公园（E11526’，N3959’），每份样品至少取自10个单株。样品采集后立即放入液氮中，后存放于-80℃冰箱备用。

1.2方法1.2.1总RNA的提取和Solexa测序总RNA的提取采用RNeasyPlantMiniKits（Qiagen，Inc.，Valen-cia，CA，USA）试剂盒按照实验手册操作步骤提取。来自芽、花、叶和果的等量RNA混合后，取10μg用于cDNA文库构建，cDNA文库的构建参考文献［17］的方法。利用IlluminaSolexaHiseq2000的paired-end测序方法进行转录组测序。

1.2.2序列的处理与拼接鉴于Solexa数据错误率对结果的影响，对原始数据进行质量预处理。先利用滑动窗口法去除低质量片段：质量阈值20（错误率=1%），窗口大小5bp，长度阈值35bp；再切除reads中含N部分序列：长度阈值35bp；最后利用软件Trinity对辽东栎有效reads进行从头拼接［18］。

1.2.3功能注释使用BLASTx程序将拼接所得的unigene与核酸、蛋白质序列数据库比对（E值<1e-5），并选取最佳注释。蛋白质数据库包括Swiss-Prot、GenBank非冗余蛋白数据库（Nr）、蛋白质直系同源簇数据库（ClustersofOrthologousGroupsofprot-eins，COGs）、京都基因和基因组百科全书（KyotoEncyclopediaofGenesandGenomes，KEGG）及GeneOntology（GO）。其中，unigene通过COG、GO和KEGG数据库的分类的参考文献［19］的方法。

1.2.4SSR位点的筛选利用MISA软件在所有unigene中搜索SSR位点，参数设置如下：二核苷酸至少重复次数为6，三核苷酸、四核苷酸、五核苷酸和六核苷酸至少重复次数均为5。

2结果

2.1IlluminaSolexa测序和序列拼接采用IlluminaSolexaHiseq2000测序技术对辽东栎的芽、花、叶和果实的混合样进行转录组测序，共获得46400862条原始序列，总长4.64Gb。经过预处理，最终得到了有效序列40621588条，数据量为3.8Gb，平均长度为94.95bp。使用软件Trinity对有效reads进行从头拼接最终得到了151339个长度大于200bp的contig，总长度约为130.92Mb，最大长度、平均长度以及N50分别为11284、865和1442bp。取每个contig下最长的转录本作为unigene，得到了95800个unigene，总长度约为73.57Mb，平均长度与N50分别为768bp和1373bp。其中，大于2000bp的序列共有7975条，占unigene总数的8.32%（图1-A）。

2.2序列的比对、功能注释及unigene的特征分析将所获得的unigene与公共数据Nr和Swiss-Prot进行Blastx比对，通过gene的相似性进行功能注释，共有57637条unigene获得了基因注释（hit），占总unigene总数的60.16%，而在其余未得到任何一个上述数据库的注释（no-hit）38163条unigene（39.84%）中，有37752条unigene（98.92%）小于或等于1000bp（图1-B）。

2.3功能分类研究为进一步研究辽东栎中unigene的功能分类，将所得到的95800条unigene在COG和GO数据库中进行比对及功能注释与分类。在COG分类中，共有36407条unigene（占unigene总数的38%）被注释到24个COG类别中。其中，“一般功能基因”是最大类别，包含8330条unigene，占被注释到unigene总数的22.88%；其次是“蛋白质翻译后修饰与转运，分子伴侣”，包含4128条unigene；而“核酸结构”是最小的类别，仅包含11条unigene。此外，有2491条unigene参与了碳水化合物运输与代谢（图2）。利用GO对获得辽东栎unigene进行功能分类，共有43766条unigene被注释到生物学过程、细胞组分和分子功能3个大类别中。其中，36372条unigene归入生物学过程，18876条unigene归入细胞组分以及40358条unigene归入分子功能。3个大的类别又被划分为45个小的类别（图3）。过程是生物学过程中的最大类别，包含31427条unigene；细胞是细胞组分中的最大类别，包含12764条unigene；蛋白结合是分子功能中的最大类别，包含28892条unigene。

2.4代谢途径分析和淀粉合成基因的筛选将辽东栎的unigene序列映射到KEGG数据库的参考代谢通路中，共有11468条unigene参与到185个代谢通路中。其中包含unigene最多的是代谢通路是剪接体（ko03010），共有1161条unigene。其次是内质网上的蛋白加工（ko04141），包含630条unigene。而参与淀粉与蔗糖的代谢通路（ko00500）的unigene共有434条（表1）。的unigene，编码9个关键酶，其中5个unigene编码α-糖苷酶；17个unigene编码β-呋喃果糖苷酶；9个unigene编码己糖激酶；10个unigene编码果糖激酶；4个unigene编码葡萄糖-6-磷酸异构酶；4个unigene编码葡萄糖磷酸变位酶；8个unigene编码葡萄糖-1-磷酸腺苷酰基转移酶；6个unigene编码淀粉合成酶及4个unigene编码1，4-α-葡聚糖分支酶（表2）。结合KEGG数据库中参考pathway的关于淀粉与蔗糖代谢通路中发掘到的unigene及在其他公共数据库中的注释，共统计筛选出67条参与淀粉合成

2.5SSR分析利用MISA软件在辽东栎的13380条unigene中共搜索到15901个SSR位点，占unigene总序列数的13.97%，平均每1.28kb出现1个SSR，其中包含有两个及两个以上SSR的unigene共有2521条。二核苷酸和三核苷酸重复类型分别占SSR总数的60.07%和38.09%；四核苷酸重复类型占SSR总数的1.57%；而五核苷酸和六核苷酸重复类型在辽东栎中转录组序列中含量较少，仅占SSR总数的0.19%和0.09%。除此之外，不同核苷酸的重复次数也有很大的变化（表3）。

3讨论

基于高通量测序技术的转录组学研究是一种非常高效、可靠的发掘功能基因手段，且在淀粉的合成及代谢中也已得到广泛应用［20，21］。目前，高通量测序技术主要是Roche公司的454测序技术、Illumina公司和ABI公司相继推出的Solexa和SOLid测序技术。其中，Solexa测序技术相对于其他两种技术在测序成本和数据量输出方面更具优势［22］。前人的研究表明，不同组织的混合取样，可在节约试验成本的基础上发掘到更多的转录本［23］。因此，本研究采用辽东栎的芽、花、叶和果实的混合样品进行转录组测序。栎属植物的基因组大小约为539-921Mb［24］，本研究共获得3.8Gb的有效数据量，约覆盖辽东栎基因组的4.13-7.05倍。同时利用在对IlluminaSolexaHiseq2000测序数据进行拼接过程中表现非常优异的Trinity软件［25］对本研究的数据进行拼接处理，共获得95800条unigene，其中有38163条unigene未在Blastx同源性搜索中得到注释，但大多数片段小于1000bp（37752条，占98.92%），因此这些片段可能是由于较短而未与公共数据库中的序列比对上，也可能是短的非编码序列或者是新的基因。COG和GO的功能分类对初步了解基因的功能起着重要作用，而KEGG数据库中的参考pathway不仅可以推测基因的功能，而且可以研究基因在不同代谢通路中所在位置及作用，三者相辅相成，成为新物种中发掘功能基因的重要手段［19］。本研究通过KEGG数据库中的pathway分析筛选与淀粉合成相关的基因，同时结合所筛选到的unigene在本研究中其他数据库中的功能注释，从而进一步确保了所获得的基因的可靠性。鉴于分布广泛和多态性较高的特性，SSR标记已被广泛应用于分子标记辅助选择育种（Molecularmarker-assistedselection，MAS）和利用分子标记通过关联分析（Associationmapping）发掘与好的农艺性状连锁紧密的相关基因［26］，尤其是对于多年生的植物，可大大缩短育种年限［27］。本研究在辽东栎中发掘到15901个SSR位点，其中二核甘酸和三核苷酸的重复占到总数的98.16%，为了保证SSR位点的潜在多态性，我们在筛选过程中对于四、五和六核苷酸的最小重复次数同样设置为5，一定程度上影响到了这3类核苷酸重复在总SSR位点中所占比例。本研究通过高通量的测序，获得了大量的辽东栎转录组序列，不仅丰富了辽东栎的基因库，而且为辽东栎及其他栎类淀粉合成基因的克隆与功能研究奠定了基础。发掘到的SSR位点可通过进一步的开发为辽东栎的进化与多样性分析及辽东栎的遗传图谱的构建和QTL（QuantitativeTraitLoci）的定位提供了数据支持。

4结论

本研究通过SolexaHiseq2000高通量测序，获得3.8Gb的辽东栎转录组序列，拼接获得95800条unigene，发掘出67条参与淀粉合成的unigene以及15901个SSR位点。

作者：刘玉林李伟张志翔单位：北京林业大学生物科学与技术学院北京林业大学自然保护区学院