苦茶全长转录组测序及基因结构分析

　　摘要为进一步探究苦茶特异性状形成的遗传基础，基于PacBio对苦茶进行Iso-Seq，最终得到Polishedconsensus序列132334个，其中有26184个为已知基因的新转录本，7115个为新基因，同时鉴定得到21106个SSR位点;基因结构分析结果中，4892个基因发生了可变剪切事件，其中内含子滞留事件占比最大;预测得到1918个新的转录因子，778个LncRNA。比较KEGG(10976+5626)、KOG(6296+1896)、GO(6221+575)三大数据库功能注释情况，发现注释到KEGG的转录本数最多，其中多个转录本注释到与茶叶品质和苦茶特异性状形成相关的代谢途径，如苦茶碱等嘌呤生物碱代谢(74+17)、氨基酸代谢(409+69)、苯丙烷生物合成(70+18)和萜类物质生物合成(152+31)等(括号内数值为未比对到基因组的转录本和新转录本数量)。这些发现将有助于苦茶特异性状相关基因标记的开发、并为其形成的机理及遗传机制的研究奠定基础。

　　关键词苦茶;全长转录组测序;LncRNA;可变剪切

　　苦茶(Camelliasinensis)是中国的一种特异茶树种质资源[1]，由其制成的茶叶感官品质表现极苦，且与其他茶类的苦味不同，因此，被命名为苦茶。苦茶多分布于云南、四川，以及广东、湖南、江西毗邻区，尤以云南以及南岭山脉两侧最多[2-3]。在苦茶的原生地，苦茶成为当地人民生活的必须品，将它作为一种药物饮用，认为长期饮用苦茶具有“退火发汗、解毒、治病”的功效[2]。

　　农业论文投稿刊物：《西南农业学报》主要刊登农牧业各个学科应用基础研究、应用技术、区域综合开发和软科学方面的研究报告、专题论文、研究简报、学科进展述评等文章，立于大西南、面向全国、系学科地、及时地、有见地反映西南地区农业科学研究新成果、前沿动态，以推动农业科研、教育和生产的发展。

　　研究已证实，苦茶嘌呤生物碱的组分与常规茶树差异较大，其以苦茶碱(1，3，7，9-四甲基尿酸)为主，其次是咖啡碱、可可碱[3-4]，苦茶碱具有镇静催眠[5]、抗抑郁[6]等药理活性，并且毒理试验鉴定苦茶碱为无毒[6]，其还可以削弱咖啡碱的兴奋作用[7]，这使得苦茶的价值不断提升，消费需求连年增加。苦茶除含有高含量的苦茶碱外，还有非儿茶素类茶多酚[8]。

　　同时，苦茶中还具有高含量的表没食子儿茶素(EGC)，EGC含量与红茶中关键品质成分茶黄素的含量呈极显著正相关[2，9]，所以苦茶可以作为选育和培育特异高级红茶茶树品种的育种材料;以往的研究表明，苦茶茶叶中存在一种具有丁香香气的特殊物质丁子香酚甙[2]。这说明苦茶还具有多种不同于常规栽培茶树的特异性状，这些性状形成的机理、性状相关基因标记的开发以及遗传机制的研究，都需要借助于苦茶资源。当前，大多数转录组测序是基于Illumina平台的第二代测序(SGS)技术生成的[10]。

　　但是，SGS技术不能产生较长的转录本，并且SGS无法获得可变剪接等信息，从而限制了该技术在转录组测序中的利用[11]。而全长转录组测序是基于第三代测序(TGS)平台进行的，该方法可以产生长读段，因此可以直接测序全长转录本[12]，在这方面全长转录测序优于短读测序;另外，它还可用于选择性剪接事件及初级-前体-成熟RNA结构的分析，以帮助更好地理解RNA的加工过程;此外，其可以在转录水平上全面分析由交替剪接和基因融合产生的结构变异，比SGS更准确地检测结构变异(SV)，使其适合于多倍体物种或具有高重复序列和高杂合性的物种的转录组分析[13]。

　　由于茶树基因组具有高的重复序列(至少包含64%的重复序列)和高的杂合度[14]，以及苦茶特异性状形成机制需进一步解析。因此，本研究以苦茶为材料，利用PacBio平台进行Iso-Seq，对所获得序列与参考基因组比对、同时对其进行功能注释、基因结构分析等，为苦茶特征形状遗传分析提供基础数据，为更好地了解和利用这一重要的特异茶树资源提供帮助。

　　1材料与方法

　　1.1试验材料

　　用于测序分析的苦茶资源‘老曼娥苦茶’来源于云南省农业科学院茶叶研究所试验基地，其6个不同组织(芽、叶、花芽、花蕾、茎和幼果)的样品采摘后用液氮冷冻，并置于-80℃冰箱，之后送至诺禾致源科技股份有限公司(北京)进行测序。

　　1.2RNA提取与构建文库

　　参照张亚真等[15]的方法，分别提取苦茶芽、叶、花芽、花蕾、茎、幼果不同部位的总RNA，检测合格后进行等量混匀。利用带有Oligo(dT)的磁珠从混匀后的RNA样品中分离出mRNA，并利用SMARTerPCRcDNASynthesisKit将其反转录为cDNA，接着用BluePippin对cDNA进行片段筛选，对全长cDNA进行损伤修复、末端修复，并连接SMRT哑铃型接头，最终使用核酸外切酶消化获得文库。

　　1.3Iso-Seq及组装分析

　　采用软件SMRTLinkv5.0对输出进行过滤和处理，参数：--minLength=200，--minRead-Score=0.75，最终得到的数据即为有效数据，校正后获得CCS序列(CircularConsensusSe-quence，参数：--minPasses=1，minPredictedAc-curacy=0.8)。

　　通过检测CCS是否包含5′-prim-er、3′-primer、poly-A对CCS进行分类找出FLNC(full-lengthnonchimera)序列，之后对FLNC进行聚类和校正，最终获得polishedcon-sensus序列用于后续分析。使用MISA(http：//pgrc.ipk-gatersleben.de/misa/misa.html)软件对单、二、三、四、五和六核苷酸的最低重复次数分别设置为10、6、5、5、5和5，进行SSR搜索。

　　2结果与分析

　　2.1Iso-Seq分析

　　2.1.1测序结果与组装全长

　　转录组测序共获得8730808个Subreads(10G)，平均Subreads的长度为1145bp，N50为1741。通过校正，CCS序列(环形一致性序列)有418424个，含有5′端引物的reads数有338306个，含有3′端引物的reads数目为368458个，含有PolyA尾的reads数目为351952个;全长序列数目为295803个，FLNC序列的数目有217701个，且FLNC的平均长度为1836bp;最终获得Pol-ishedconsensus序列132334个，用于后续分析。

　　2.1.2SSR分析

　　利用MISA对苦茶全长转录组测序获得的25735条(大于500bp)转录本序列进行搜索，共在12926条转录本序列中发现符合标准的21106个SSR位点，其中5227个转录本含有大于1个SSR位点，SSR的发生频率为50.22%，SSR出现率为82%，平均2.4kb出现1个SSR。

　　在检测到的苦茶转录组SSR中，单、二和三核苷酸重复类型的比例较大，分别占总SSR的32.49%、47.91%和16.18%;其他3种重复类型所占比例相对较少，只占总SSR的3.42%。苦茶转录组SSR重复单元的重复次数分布在5～46，其中5～10次重复的SSR位点有13625个，占位点总数的64.56%;11～20次重复的SSR位点有7299个，占总数的34.58%;20次重复以上的SSR位点有0.86%。

　　3讨论与结论

　　目前，PacBio-Iso-Seq测序手段已在多种植物中得到了广泛的应用，该方法相比于二代测序具有读长长、无需组装转录组就可以直接获得全长转录本、更准确地检测结构变异和适合于具有高重复序列和高杂合性的物种的转录组分析等优势[11-13]。本试验利用Iso-Seq技术，结合生物信息学分析，获得了苦茶全长转录组数据。

　　本研究共得到8730808个Subreads，平均Subreads的长度为1145bp，N50为1741，可以看出全长转录组测序读长长且连续性较高。全长非嵌合(FLNC)reads有217701个，通过对FLNC-reads的聚类及校正分析，最终获得pol-ishedconsensus序列132334个。庞丹丹等[17]对不同叶色的6份茶树材料进行二代测序，共获得112233个Unigene，平均长度为759bp，N50为1081bp。

　　Li等[18]采用RNA-seq对龙井43的根、茎、4个不同嫩度的叶片及花蕾和果实等13个组织进行测序分析，共获得347827条Unigene(>200bp)，平均长度为791.2bp，N50为1342bp。上述结果表明在获得的序列质量和基因数等方面，三代测序结果均优于二代测序。苦茶的全长转录组测序在很大程度上补充了现有茶树基因资源，并为发现新的或以前未被识别的蛋白质的编码基因和转录亚型提供了优势。

　　与茶树参考基因组比对，分析发现大多数PacBio转录本是已知基因的新亚型(59.32%)，有16.12%的转录本为7115个新基因提供了证据。分析表明，新的转录组数据对探究苦茶的注释具有巨大的潜力。Unmapped转录本在NR数据库中的比对结果显示，共比对上311个物种，比对到最多的物种仍是茶树，其次是葡萄，这表明比对所用到的参考基因组可能还需进一步完善。

　　比较KEGG(10976+5626)、KOG(6296+1896)、GO(6221+575)3大数据库功能统计情况，发现注释到KEGG的转录本数最多(括号内数值为未比对到基因组的转录本和新转录本数量，下同)，这些转录本中包含多个转录本与茶叶品质和苦茶特征性状形成相关的代谢途径，如苦茶碱等嘌呤生物碱代谢(74+17)、氨基酸代谢(409+69)、苯丙烷生物合成(70+18)和萜类物质的生物合成(152+31)等。Wang等[19]对3个苦茶品种进行转录组测序，通过与嘌呤生物碱代谢的KEGG通路富集分析，筛选出一批与苦茶碱和咖啡碱生物合成相关的关键基因。同样，本研究新注释到的一些转录本亦能够为后期对苦茶性状形成的研究提供基因资源。

　　可变剪切(AlternativeSplicing，AS)，即某些基因的一个mRNA前体以不同的剪接方式，形成不同的mRNA异构体[20];其能够调控基因的表达水平和促进蛋白质组的多样性，在植物生长、发育和胁迫响应中发挥着关键作用[21-23]。选择性剪接的亚型具有基于组织或时间的优先表达特征，并且它们也受环境条件的影响[24]。在茶树中，研究发现可变剪切亚型在高温、干旱、低温等不同环境胁迫条件下的表达模式存在明显差异[25-26]，同时可变剪切事件还参与了黄酮、类黄酮和花青素等重要物质的次生代谢过程中[24，27-28]。

　　因此，对基因结构进行分析，发现4892个基因发生了选择性剪接，其中内含子滞留事件占比最大，这为进一步研究可变剪切在苦茶特异成分代谢中的作用奠定了基础。以往的研究表明，转录因子(TF)参与植物多种生长代谢过程中。在茶树中，发现它们不仅参与叶片花等器官的发育[29-30]，以及胁迫反应[31]，还参与儿茶素[17]、茶氨酸[32]和花青素[33-34]等多种次生代谢过程。

　　本研究还预测得到1918个TF，在获得的TF家族中，bHLH家族有96个、C3H家族有80个、bZIP家族有70个、MYB-re-lated家族有66个、C2H2家族有66个、MYB家族有57个。这些转录因子的预测与分析，为之后研究苦茶特异性状(苦茶碱[3-4]、丁子香酚甙[2]和高EGC[9]等)相关基因表达以及详细的基因家族分析提供数据。此外，本研究还预测得到778个LncRNA，它们可能在茶树次生代谢的调控中起着特殊的作用。本研究为苦茶品质和特异性状形成机制的研究提供了基础数据，上述结果将有助于进一步开展苦茶特异性状相关基因标记的开发、并为其形成的机理及遗传机制的研究奠定基础。

　　参考文献Reference：

　　[1]张宏达.茶叶植物资源的订正[J].中山大学学报(自然科学版)，1984(1)：3-14.ZHANGHD，Aevisiononthetearesourceplants[J].ActaScientiarumNaturaliumUniversitatisSunyatseni，1984(1)：3-14.

　　[2]王新超，姚明哲，马春雷，等.我国苦茶资源主要生化成分的鉴定评价[J].中国农学通报，2008，24(6)：65-69.WANGXCH，YAOMZH，MACHL，etal.Analysisand evaluationofbiochemicalcomponentsinbitterteaplantgermplasms[J].ChineseAgriculturalScienceBulletin，2008，24(6)：65-69.

　　[3]汪云刚，刘本英，宋维希，等.云南茶组植物的分布[J].西南农业学报，2010，23(5)：1750-1753.WANGYG，LIUBY，SONGWX，etal.DistributionofSect.Thea(L.)DyerinYunnanprovince[J].SouthwestChi-naJournalofAgriculturalSciences，2010，23(5)：1750-1753.

　　[4]JOHNSONTB.Purinesintheplantkingdom：Thediscov-eryofanewpurineintea[J].JournaloftheAmericanChemicalSociety，1937，59(7)：1261-1264.

　　[5]叶创兴，林永成，苏建业，等.苦茶CamelliaassamicavarkuchaChangetWang的嘌呤生物碱[J].中山大学学报(自然科学版)，1999，38(5)：82-86.

　　作者：庞丹丹1，2，刘玉飞1，2，孙云南1，2，田易萍1，2，宋维希1，2，陈林波1，2