基于RNA-Seq技术的茶树花转录组分析

夏丽飞，朱兴正，梁名志，马 伟，宋维希，田易萍，周 萌，陈林波

(云南省农业科学院茶叶研究所/云南省茶树种质资源创新与配套栽培技术工程研究中心，云南 勐海 666201)

基于RNA-Seq技术的茶树花转录组分析

夏丽飞，朱兴正*，梁名志，马 伟，宋维希，田易萍，周 萌，陈林波**

(云南省农业科学院茶叶研究所/云南省茶树种质资源创新与配套栽培技术工程研究中心，云南 勐海 666201)

本研究利用RNA-seq技术对父本、母本、子代不育茶树花三个样本花的转录组进行测定，经组装分析获得403 469条高质量的Unigenes。将获得的Unigenes与SWISS-PROT、TREMBL、CDD、PFAM、NR和KOG库进行blast，共注释到307 291条Unigenes。KOG功能分类显示，有23 739个Unigenes被分为25类。KEGG通路分析表明，共识别出26 967个Unigenes涉及的pathway有328个。SSR查找发现，从403 469个Unigenes中找到46 440个含有SSR序列。这些信息为茶树不育基因筛选、不育机理研究以及分子标记开发奠定了基础。

茶树；转录组测序；不育基因；分子标记

茶树是重要的叶用经济作物，茶树栽培的目的是为了多采收嫩芽叶。茶树的生长过程包括营养生长和生殖生长，生殖生长期长，从当年的花芽分化和花器官发育，到下一年的茶果生长成熟，要花费1年半的时间，消耗大量的营养物质，影响了茶叶的营养生长[1-2]。茶树开花结果会争夺茶叶养分，致使茶叶产量下降、质量降低[3-4]。因此，对于以采收嫩芽叶的茶树来说，减少或阻止茶树开花结实，控制生殖生长，促进营养生长是关系到茶树产量和品质提高的关键，也是关系到茶树育种的关键，同时也是防止茶树老化的措施之一。

云南省农业科学院茶叶研究所科技人员以福鼎大白茶茶树作为父本、佛香2号茶树作为母本，通过人工杂交授粉获得一株特异茶树材料。通过多年的观测发现，该花没有花粉，经人工去雄授粉发现该花不结茶果，自然状态也不结果，对该花的可溶性糖分析发现，其可溶性糖低于正常花。因此，该花在生理上或遗传上显示出明显的不育特征[5]，是茶树花不育的天然突变体。本研究利用RNA-seq技术对母本、父本、子代不育花3样本的转录组进行测序以及生物信息学相关分析，为筛选茶树不育基因以及进一步开展功能基因组研究提供帮助。

1 材料与方法

1.1 研究材料

选用云南省农业科学院茶叶研究所实验基地的福鼎大白茶花(父本)、佛香2号花(母本)以及他们的杂交后代(不育花)为材料，11月份采摘父本花(花蕾、初开花、盛开花)、母本花(花蕾、初开花、盛开花)以及后代不育花(花蕾、初开花、盛开花)，立即放入液氮中固定并放入低温冰箱保存。

1.2 RNA的提取与检测

RNA的提取选用CTAB-LiC1方法[6-7]，分别提取佛香2号茶树(母本)花蕾、初开花、盛开花RNA，采用Nanodrop核酸检测仪及2 %的琼脂糖凝胶电泳检测总RNA的含量、稳定性和完整性，再进行等量混匀，作为母本花样本RNA，同理获得父本、不育花样本RNA。将获得的父本花RNA、母本花RNA以及不育花RNA送苏州贝斯派生物科技有限公司分别进行转录组测序。

1.3 转录组测序

1.3.1 Unigene序列组装 采用Hiseq2500的paired-end测序得到的reads，利用cutadapt软件切除reads中含N部分序列以及去除接头，再使用软件Trinity进行拼接获得Unigene。

1.3.2 序列注释、功能分类和生物学通路分析 将获得的Unigene与SWISS-PROT、TREMBL、CDD、PFAM、NR和KOG库进行比较，取相似度>30 %，且e<1e-5的注释，对注释的基因进行KOG功能分类预测以及利用KEGG数据库作参考进行Unigene的pathway注释。

1.4 SSR位点查找与分析

对拼接得到 UniGene进行SSR简单重复序列的查找，二核苷酸重复的次数在6次或6次以上，三至六核苷酸重复的次数在5次或5次以上，同时分析SSR的基元类型的特征。

图1 Unigenes序列长度分布Fig.1 The length distribution of Unigenes sequence

2 结果与分析

2.1 高通量测序与序列组装分析

通过Hiseq2500的paired-end高通量测序技术对父本、母本、子代不育花的转录组进行测序，使用软件Trinity的paired-end拼接方法进行de novo拼接，对拼接序列删除重复，最终得到403 469个长度大于200 bp的转录本，大小259 Mb。All-Unigene序列长度主要分布在250～2000 bp，大于2000 bp的有10 618条，范围在250～950 bp的较多，占全部序列数量的80 %以上见图1。

2.2 Unigenes的功能注释

将样本Unigenes与Nr、SWISS-PROT、TrEMBL、Cdd、pfam和KOG库进行比较。307 291条Unigenes被注释，有96 178的Unigenes未得到注释，说明了茶树花转录组中含有大量的未知基因。表1为注释到Nr、SWISS-PROT、TrEMBL、Cdd、pfam和KOG等6个数据库条目。其中以TrEMBL、NR和PFAM数据得到的条目最多，分别占全部条目27.42 %、27.03 %和24.22 %。

注释基因的同源序列的物种分布情况见图2，注释到葡萄(Vitisvinifera)的序列占40.78 %，其次是中粒种咖啡(Coffeacanephora)，有10.21 %，再次是可可(Theobromacacao)，有10.19 %。这是因为NR数据库中保存葡萄、咖啡和可可大量的基因信息，为本研究转录组的注释提供丰富的参考序列。

表1 Unigenes的注释统计

图2 Unigenes在NR库中物种分布Fig.2 Species classification of the Unigenes in NR bank

对Unigenes进行GO功能分类预测，Level2水平下的GO注释，主要分为生物学过程、细胞组分和分子功能三大类67个小类。其中代谢过程、细胞组成过程、细胞、细胞部分、捆绑、催化活性6类富集程度较高，这可能与茶树花发育过程分生组织不断进行细胞增殖以及花器官组织内代谢活动旺盛有关，见图3。

对预测基因进行KOG功能分类预测，共有23 739个Unigenes被注释上25种KOG分类见图4。如RNA加工与修饰、染色质的结构和动力学、能源产生和转换、细胞周期调控、氨基酸的转运和代谢等。在25种KOG分类中，注释最多的是信号转导机制(T)，其次是一般功能预测(R)。值得注意的是，有3742条被注释到信号转导机制上，为后续研究茶树花不育奠定了良好的基础。

图3 Unigenes的GO功能注释及分类统计Fig.3 GO function annotation and classification of Unigenes

图4 Unigenes的COG功能注释分布Fig.4 COG functional annotation distribution of Unigenes

编号Code代谢途径MetabolicpathwaysUnigene的条目Unigeneitems百分比(%)Percentage代谢通路IDMetabolicpathwaysID1碳代谢10613．93ko012002氨基酸的生物合成8463．14ko012303淀粉和蔗糖代谢6992．59ko005004糖酵解/糖质新生5482．03ko000105丙酮酸代谢4341．61ko00620

表3 SSR不同重复基元分布及优势碱基组成

2.3 Unigenes的代谢通路分析

生物学通路主要包括代谢通路与信号转导通路，根据KEGG pathway生物学通路的注释和预测，共识别出26 967个 Unigenes，具有对应的1081种酶(Enzy-me, EC)功能，并关联到328条生物学通路。其中与 “代谢通路”相关的 Unigenes 最多是碳代谢通路，占3.93 %(1 061条)，其次是氨基酸的生物合成、淀粉和蔗糖代谢，分别占3.14 %、2.59 %。表2列出Unigenes注释最多的5个pathway。此外，从KEGG分析中，鉴定出29个编码含MADS-box转录因子的Unigenes。这些注释为后续研究茶树花的发育以及不育机理提供有价值的信息。

2.4 SSR信息分析

对拼接得到Unigenes进行SSR简单重复序列的查找。筛选标准为二核苷酸重复的次数在6次或6次以上，三至六核苷酸重复的次数在5次或5次以上，SSR结果的统计结果见表3。分析403 469条Unigenes中挖掘出46 440 SSR位点，占11.51 %。从2个到6个碱基重复基元均存在，说明SSR的类型丰富。其中2碱基型的最多，占74.24 %，其次是3碱基型的占23.53 %，最少的是5碱基型的占0.14 %。在SSR检测中出现频率最高的优势基元有AG/TC、TGG/ACC、CACT/GTGA、TGACA/ACTGT、CGCCAT /GCGGTA见下表4。上述SSR的分析，有助于开展茶树遗传图谱的构建研究以及不育基因标记的开发。

3 讨 论

近年来，新一代测序技术凭借高通量、高灵敏度和运行成本低已成为生命科学研究的新手段，而建立在高通量测序基础上的转录组研究已逐步取代基因芯片技术成为目前从全基因组水平研究基因表达的主流方法[8-10]。RNA-seq测序技术进行的转录组分析是近年来建立的分子生物学研究分析方法，具有信息涵盖量大、分析准确、数据冗余性低、可检测到低表达的基因存在等特点[11]。本研究利用RNA-seq技术对茶树母本、父本、子代不育花三样本的转录组进行测序以及生物信息学相关分析，共获得403 469条质量较高的Unigenes，有307 291得到了注释，共识别出26 967个Unigenes对应的1081种功能酶，并关联到328条生物学通路上。这些数据的获得为下一步分析正常花与不育花之间的基因表达差异、筛选不育基因以及研究不育机制奠定基础。

MADS-box基因是一类编码转录因子的基因家族，广泛存在于动物、植物和真菌中，其编码的转录因子在真核生物的生长发育和信号转导过程中发挥着重要作用[12]，特别是在决定花器官特征方面发挥着重要的调控作用[13-15]。目前研究最为清楚的是MADS-box基因与花器官发育之间的关系，即著名的ABCDE模型[16]。本研究鉴定出29个编码含MADS-box转录因子基因片段，初步获得了一些参与决定花器官特征的基因序列信息，为下一步MADS-box转录因子基因克隆、表达模式以及在茶树花发育过程中的作用等方面的研究提供信息资料。

转录组数据是SSR和SNP开发标记的理想资源，而SSR是利用转录组数据开发最多最常用的分子标记[15-17]。本研究从3个茶树花转录组403 469条transcript中挖掘出46 440个SSR位点。这些SSR信息的获得将为下一步开发不育基因的SSR标记、遗传多样性的研究以及遗传连锁图谱的构建奠定基础。

[1]陈兴淡.茶树育种学[M].北京:中国农业出版社,1980:82-85.

[2]夏春华,束际林.茶树化学除花与茶叶增产[J].植物杂志,1979(2):6-7.

[3]杨昌云,朱永兴.茶树生殖生长的影响因素及控制方法[J].中国茶叶,1999(5):12-13.

[4]吴淑平.茶树营养生长与生殖生长的关系及调控方法[J].中国园艺文摘,2011(5):182-183.

[5]朱厂廉.太谷核不育小麦花药内游离脯氨酸的来源利用及不育花败育的关系[J].植物生理学报,1985,11(2):122-129.

[6]史成颖,宛晓春,江昌俊,等.提取高质量茶树总RNA的方法研究[J].安徽农业大学学报,2007,34(3):360-363.

[7]江昌俊,王朝霞,李叶云.茶树中提取总RNA的研究[J].茶叶科学,2000, 20(1):27-29.

[8]Sun C, Li Y, Wu Q, et al. De novo sequencing and analysis of the American ginseng root transcriptome using a GS FLX Titanium platform to discover putative genes involved in ginsenoside biosynthesis[J]. BMC Genomics, 2011,11:262.

[9]Wang W, Wang YJ, Zhang Q, et al. Global characterization of Artemisia annua glandular trichome transcriptome using 454 pyrosequencing[J]. BMC Genomics, 2009,10:465.

[10]杨 莹，高 珊，潘宝平，等.青蛤(Cyclinasinensis)IRAK-4基因的克隆及其组织间的表达分析[J].安徽农业科学,2015,43(27):33-36.

[11]Alagna F D, Agostino N, Torchia L, et al. Comparative 454 pyrosequencing of transcripts from two olive genotypes during fruit development[J]. BMC Genomics, 2009,10(1):399.

[12]Becker A, Thesissen G. The major clades of of MADS-box genes and their role in the development and evolution of flowering plants[J]. Molecular Phylogenetics and Evolution, 2003,29(3):464-489.

[13]Theissen G. Development of floral organ identity: Stories from the MADS house[J]. Curr Opin Plant Biol, 2001,4(1):75-85.

[14]Heijmans K, Morel P, Vandenbussche M. MADS-box genes and floral development: The dark side[J].J Exp Bot,2012,63(15):5397-5404.

[15]袁秀云，蒋素华，王默霏，等.蝴蝶兰MADS-Box基因克隆及植物表达载体的构建[J].南方农业学报,2014,45(3):345-351.

[16]王力娜,范术丽,宋美珍,等.植物MADS-box基因的研究进展[J].生物技术通报,2010(8):12-19.

[17]Luro F L, Costantino G, Terol J, et al. Transferability of the EST-SSRs developed on Nules clementine (CitrusclementinaHort ex Tan) to otherCitrusspecies and their effectiveness for genetic mapping[J]. BMC Genomics, 2008, 9:287.

[18]Novaes E, Drost D R, Farmerie W G, et al. High-throughput gene and SNP discovery in Eucalyptus grandis, an uncharacterized genome[J]. BMC Genomics, 2008(9): 312.

[19]袁阳阳,王青锋,陈进明.基于转录组测序信息的水生植物莕菜SSR标记开发[J].植物科学学报,2013,31(5):485-492.

(责任编辑 王家银)

Transcriptome Analysis of Flower from Tea Plant [Camelliasinensis(L.)O.Kuntz] via RNA-Seq Technology

XIA Li-fei, ZHU Xing-zheng*, LIANG Ming-zhi, MA Wei, SONG Wei-xi, TIAN Yi-ping, ZHOU Meng, CHEN Lin-bo**

(Tea Research Institute, Yunnan Academy of Agricultural Sciences/Yunnan Engineering Research Center of Tea Germplasm Innovation and Matching Cultivation, Yunnan Menghai 666201, China)

This research established a reference transcriptome sequencing and bioinformatics analysis of male parent, female parent and offspring sterility flowers fromCamelliasinensisby the RNA-Seq technology. A total of 403 469 unigenes were generated from theCamelliasinensisflowers transcriptome by using RNA-seq. A total of 307 291 unigenes were aligned to the sequences of public databases, such as Nr, TrEMBL, Cdd, pfam and the KOG database, and 23 739 unigenes were assigned at 25 KOG classifications and 26 967 unigenes at 328 KEGG pathways. The characteristic of SSR distribution showed that 46 440 SSRs loci were detected from 403 469 unigenes. These results laid the foundation for screening the sterility gene, studying sterile mechanism and developing of molecular markers.

Camelliasinensis; Transcriptome sequencing; Sterility gene; Molecular marker

1001-4829(2016)09-2058-05

10.16213/j.cnki.scjas.2016.09.009

2015-01-23

国家自然科学基金项目(31460216)；云南省重点专项(2013BB006)；云南省人才培养计划(2015HB105)；茶树生物学与资源利用国家重点实验室开放基金(SKLTOF20150105)

夏丽飞(1977-)，副研究员，研究方向为茶树生物技术，*为共同第一作者，**为通讯作者。

S571.1

A