论文部分内容阅读
摘 要：為了揭示lncRNA在红苞凤梨嵌合叶片形成和生长发育过程中的调控作用机制，该文以金边红苞凤梨为材料，采用Hiseq2500测序和SMRT三代全长转录组测序联合测序分析技术，分析挖掘红苞凤梨lncRNA信息。结果表明：（1）鉴定得到6 018条lncRNA，包含3 298个基因间lncRNA，870个反义lncRNA，717个内含子lncRNA和1 109个正义lncRNA，数据量较二代测序有了极大的提高。（2）结构分析表明，红苞凤梨lncRNA的总体表达丰度低于mRNA;序列长度在 400～1 200 nt区间比例高于mRNA，而在>1 600 nt区间，lncRNA分布的比例显著小于mRNA;lncRNA中的外显子数量总体少于mRNA，开放阅读框长度总体上也短于mRNA。（3）差异表达分析表明，在全绿、全白叶片发育过程中鉴定到1 710个差异表达lncRNA。（4）靶基因预测结果表明，5 441个lncRNA通过cis作用方式预测到靶基因，1 544个lncRNA通过trans方式预测到靶基因。（5）靶基因的功能注释和富集分析显示，差异表达lncRNA的靶基因主要作为酶蛋白参与调节叶片代谢活动和信号转导等方面，与叶片的颜色形成、光合作用和生长发育密切相关。该文鉴定出的lncRNA信息以及对其结构和功能的分析，为红苞凤梨以及凤梨科其他植物的lncRNA表观遗传调控机理研究提供了数据基础，筛选出的差异表达lncRNA在金边红苞凤梨叶片嵌合性状的形成和生长发育中具有重要的调控作用。
关键词：红苞凤梨，Hiseq2500测序，SMRT全长转录组测序，lncRNA鉴定
中图分类号：Q943
文献标识码：A
文章编号：1000-3142（2021）08-1237-14
Abstract： In order to reveal the regulation function of lncRNA on the chimeric character formation and development of the leaf of Ananas comosus var. bracteatus， Hiseq2500 sequencing and SMRT the third-generation full-length transcriptome sequencing were applyed to identify lncRNA of Ananas comosus var. bracteatus. The results were as follows： （1） A total of 6 018 lncRNA were identified， containing 3 298 intergenic lncRNA， 870 antisense lncRNA， 717 intron lncRNA and 1 109 sense lncRNA， which were greatly improved compared with the second-generation information. （2） Structural analysis showed that the overall expression level of lncRNA was lower than that of mRNA. The transcript length distribution of lncRNA in the range of 400-1 200 nt was higher than that of mRNA， while in the range > 1 600 nt， the proportion of lncRNA distribution was significantly lower than that of mRNA. The number of exons in lncRNA was generally less than that of mRNA， and the open reading frame was also shorter in length than that of mRNA. （3） For analysis of differential expression， 1 710 differentially expressed lncRNA were identified during the development of complete green and complete white leaves. （4） Target gene prediction results showed that 5 441 lncRNA were predicted target genes by Cis action， and 1 544 lncRNA were predicted target genes by Trans action. （5） Functional annotation and enrichment analysis of target genes revealed that the target genes of differentially expressed lncRNA mainly act on metabolic activities and signal transduction of leaves as enzyme proteins， and were closely related to leaf color formation， photosynthesis and leaf growth. The lncRNA information identified in this paper and， as well as the analysis of its structure and functions， provide the data basis for the study of the epigenetic regulation mechanism of lncRNA in Ananas comosus var. bracteatus and other plants in Bromeliaceae. The identified differentially expression of lncRNA plays an important role in the chimeric character formation and development of leaf of Ananas comosus var. bracteatus. Key words： Ananas comosus var. bracteatus， Hiseq2500， SMRT full length transcriptome sequence， lncRNA identification
红苞凤梨 （Ananas comosus var. bracteatus）因其叶片绿白镶嵌、花果颜色艳丽且观赏期长，己成为一种重要的新型观赏植物。红苞凤梨自交不亲和，在生产中以吸芽进行繁殖，繁殖系数低、苗木整齐度差，限制了红苞凤梨的规模化应用。组织培养能快速繁殖红苞凤梨，但繁殖过程中叶片嵌合性状不稳定，再生植株叶片常失去嵌合性状而变为全绿植株（曹莉，2011）。细胞白化突变是金边嵌合叶片形成的基础，研究红苞凤梨叶片细胞白化突变的分子机理，对揭示红苞凤梨嵌合性状的形成机理，提高嵌合性状的稳定性，培育新的嵌合性状品种具有重要的理论和实践意义。
本课题组前期研究结果表明，红苞凤梨白化细胞中叶绿素含量极显著下降，但叶绿素合成代谢的结构基因表达上调（Li et al.，2017;Xue et al.，2019），说明转录后调控在红苞凤梨细胞白化失绿、金边嵌合性状形成中发挥了重要调控作用。lncRNA具有类似mRNA的结构特征，可在多个层面调控目标基因的表达。可作为信号分子、诱饵分子、引导分子以及支架分子，在表观遗传、转录调控、转录后调控等多个水平发挥功能（Zhang et al.， 2018）。目前，lncRNA在人类和动物中研究较为广泛，与人类疾病的发生及生物体的生长发育密切相关（Johnson，2012;余铖亮等，2015;王艳芳等，2018），而植物lncRNA的研究还处于起步阶段。研究结果表明lncRNA在植物的开花诱导（Csorba et al.，2014）、花粉发育（Ding et al.，2012）、逆境脅迫（Qin et al.，2017）中具有重要功能，然而其具体作用机制及调节功能等尚不清楚。由于红苞凤梨没有基因组数据，且Hiseq2500二代测序技术读长较短，使得红苞凤梨非编码RNA的挖掘具有一定的局限性。随着测序技术的发展，SMRT（single-molecule real-time）三代测序技术的出现，因其无需进行PCR扩增，大幅降低了因PCR反应引入的碱基错误，操作更为简单等优点，已得到广泛应用（Flusberg et al.，2010）。目前，SMRT三代测序技术在基因组、甲基化识别、SNP 的鉴定、基因重测序和转录组学等方面的优势越来越明显（Smith et al.，2012;Guo et al.，2018）。此外，SMRT测序技术得到几千kb的数据，读长显著增长，大大减少了测序后的 Contig数量，使得基因组和转录组的组装得到极大改善（English et al.，2012）。但是，第三代测序技术也存在一定的缺点，测序错误率普遍偏高，测序产生的错误率可高达15%（Koren et al.，2017）。因此，采取二代测序和三代测序联合分析已成为当下基因组研究的主要方法。一般有两种联合手段，可以选择以三代测序为主，进行组装，再通过二代测序得到的高质量短片段对三代的数据进行碱基纠错和矫正;也可以选择以二代测序为主，用三代测序得到的长片段Reads 进行辅助组装（马建超，2018）。最为通常的手段是利用短的但准确度高的二代数据去辅助校正长的但准确度较低的三代数据，并在此基础上进行混合组装。经验证，通过这种混合组装校正方法得到的数据准确性可达99%（马东娜等，2018）。因此，这种“2 + 3”的联合模式已被广泛认可并应用于动物、植物和微生物的基因组研究（Koren et al.，2012;Hackl et al.，2014;徐伟南等，2018）。
本研究以红苞凤梨为材料，利用二代和三代测序数据混合组装校正方法，鉴定红苞凤梨叶片中存在的lncRNA，分析在全绿和全白突变叶片生长发育过程中差异表达的lncRNA，通过靶基因的功能注释和富集分析，揭示lncRNA在红苞凤梨叶片细胞失绿白化及生长发育过程中的作用。本研究结果为红苞凤梨lncRNA表观遗传调控机理研究提供数据基础。
1 材料与方法
1.1 试验材料
本研究以金边红苞凤梨茎段为外植体，通过组织培养获得的全白和全绿植株为试验材料。选取长势均一的未展叶期、4～5叶期和10～12叶期三个发育时期的全绿及全白植株各10株的叶片作为样本提取RNA（图1）。样本取下后立即转入液氮速冻，后储存于-80 ℃用于提取RNA进行Hiseq2500二代转录组测序。
1.2 红苞凤梨RNA提取及二代Illumina Hiseq2500测序文库构建
1.2.1 红苞凤梨RNA提取及检测 取-80 ℃冻存的红苞凤梨样品，采用LABGENETM plant RNA Isolation kit多糖多酚植物RNA提取试剂盒分离总RNA，操作方法参考使用说明。高质量的RNA是实验成功的基础，为保证测序的准确性，对样品RNA进行了检测，分别采用荧光定量仪（Qubit 2.0）、微量分光光度计（Nanodrop）、生物分析仪（Aglient 2100）、电泳方法检测RNA样品的纯度、浓度、完整性和是否有基因组DNA污染等，达到质控要求的RNA样品用于后续试验。
1.2.2 Illumina Hiseq2500 cDNA文库构建及上机测序 利用epicentre Ribo-ZeroTM试剂盒去除样品中的rRNA。以rRNA-depleted RNA为模板，用六碱基随机引物（random hexamers）合成cDNA第一条链和第二条链。cDNA纯化后进行末端修复、加A尾并连接测序接头，然后用AMPure XP beads进行片段大小选择。降解含U链，通过PCR富集得到cDNA文库。文库构建后，经Qubit 2.0定量，Agilent 2100检测文库的insert size质量;釆用QPCR对文库的浓度（文库有效浓度>2nmol·L-1）进行准确定量，完成库检。库捡完成后，采用Illumina Hiseq2500平台对文库进行测序。 1.3 数据质量监控
测序所得Reads的过滤与修剪是保证分析数据可靠性的关键，结合前期SMRT测序结果（Ma et al.，2018）进行混合组装校正后，删除Raw Reads中包含adapter、ploy-N的Reads和低质量的Reads以获得Clean Reads。以菠萝基因组（Acomosus_321_v3https：//phytozome.jgi.doe.gov/pz/portal.html#！info？alias=Org_Acomosus_er）为参考基因组，利用TopHat v2.0.9（Kim et al.，2013）软件对Clean Reads进行比对。采用Scripture（beta2）（Langmead et al.，2009）和Cuffiinks（v2.1.1）软件对每个样品比对上的Reads进行装配获得转录本。
1.4 轉录本表达水平和编码潜能分析
使用Cufflinks软件的Cuffdiff组件，对转录本表达水平进行分析。根据lncRNA编码的特点进行基本筛选：选择长度≥200 bp，Exon个数≥2及FPKM≥0.1的转录本。
因lncRNA不编码蛋白，因此，通过对基本筛选得到的转录本进行编码潜能筛选，判断其是否具有编码潜能，从而可以判定该转录本是否为lncRNA。主要利用CPC分析、CNCI分析、CPAT分析、Pfam蛋白结构域四种分析方法分析lncRNA的编码能力。去掉基本筛选中具有潜在编码能力的转录本，余下的即为预测的lncRNA。
通过与已知mRNA进行比较，利用Cuff-compare分析结果中的class codes对筛选的lncRNA进行分类。
1.5 lncRNA靶基因预测及功能富集分析
基于lncRNA与其靶基因的作用方式（cis和trans），采用两种预测方法：第一种是根据lncRNA与mRNA的位置关系预测lncRNA的靶基因，即lncRNA 100 kb范围内的邻近基因为其靶基因;第二种是使用LncTar靶基因预测工具进行预测，其原理是利用lncRNA与mRNA之间碱基互补配对产生的作用来预测。
对差异表达IncRNA的靶基因利用KEGG（Kyoto Encyclopedia of Genesand Genomes）、GO（gene ontology）、NR （Non-Redundant Protein Sequences Database）、COG （clusters of orthologous groups of proteins）和 Swiss-Prot 数据库进行功能注释和富集分析，分析结果显著性用P值表示。
1.6 差异表达基因分析
以差异倍数（Fold Change）≥2且错误发现率FDR（False Discovery Rate）<0.05作为差异表达的筛选标准，使用EBseq分析6个样本间差异表达的lncRNA及mRNA。
2 结果与分析
2.1 测序数据与参考基因组比对效率分析
长链非编码数据利用率的直接体现就是比对效率，即Mapped Reads占Clean Reads的百分比。经SMRT全长转录组数据（NCBI提交号PRJNA564223）（Ma et al.，2018）修正后，6个样品的Reads与所选参考基因组的比对效率为67.74%～78.58%，比单独使用Hiseq2500二代测序数据的比对效率提高了约5%（蔺珍，2019）（表1）。说明三代数据的修正，有效提高了lncRNA测序数据的利用效率，有利于进一步深入挖掘红苞凤梨lncRNA信息。
2.2 红苞凤梨lncRNA的鉴定
对测序样品Cuffiinks的拼接结果第一步经过Cuffcompare软件分析合并，筛选出转录本长度≥200 bp，外显子个数≥2，FPKM ≥ 0.1的转录本，根据已知mRNA数据库，进行比较过滤mRNA，然后通过CNCI、CPC、CPAT 和Pfam四个软件进行蛋白编码潜能筛选，最终鉴定得到6 018条lncRNA，其中新lncRNA 5 689条（图2）。比利用Hiseq二代测序数据鉴定到的lncRNA数量增加了约70%（蔺珍，2019），极大地提高了红苞凤梨lncRNA的数据信息量，为进一步研究红苞凤梨的非编码RNA调控机理提供了数据基础。
鉴定得到的6 018个lncRNA中，包括了3 298个基因间lncRNA（intergenic lncRNA），870个反义lncRNA（antisense lncRNA），717个内含子lncRNA（intronic lncRNA）和1 109个正义lncRNA（sense lncRNA）（图2：B）。与二代测序分析鉴定结果相比，intergenic lncRNA所在比例显著提高，由17%增加到55%。而sense lncRNA所占比例显著下降，由76%下降到18.5%（蔺珍，2019）。
2.3 红苞凤梨lncRNA结构分析
为了进一步分析红苞凤梨lncRNA的结构特点，将lncRNA与蛋白质编码RNA在整体表达水平、序列长度分布、外显子数目分布以及开放阅读框长度分布情况进行了比较分析（图3）。分析结果表明，在总体表达水平上，mRNA的表达丰度高于lncRNA的表达丰度（图3：A）。在转录本的长度分布上，lncRNA在400～1 200 nt区间分布比例高于mRNA（图3：B），二代测序结果中则是表现在400～600 nt和1 400～1 600 nt区间（蔺珍，2019）;而在转录本长度>1 600 nt区间，lncRNA分布的比例显著小于mRNA，尤其在>=3 000 nt区间（图3：B）。lncRNA中的外显子数量总体少于mRNA，约82%的lncRNA只含有2个外显子（图3：C）。而二代测序分析中有41.80%的lncRNA只有2个外显子，而 31.62%的mRNA的外显子数超过 5 个（蔺珍，2019）。lncRNA的开放阅读框长度总体上也短于mRNA。约99%的lncRNA开放阅读框长度<=100 nt（图3：D），而二代测序分析中66%的lncRNA的开放阅读框在 0～100 nt之间（蔺珍，2019）。 2.4 红苞凤梨lncRNA的差异表达分析
以Fold Change≥2.0且FDR<0.05作为筛选标准，共鉴定得到了1 710个差异表达的lncRNA。对筛选出的差异表达lncRNA进行了层次聚类分析（图4：A）。在未展叶期，全绿和全白叶片中，大量差异表达的lncRNA表达丰度较高，在全白苗和全绿苗的发育到4～5叶期时，多数差异表达的lncRNA表达水平下降。而全白苗在第三发育时期有部分差异表达lncRNA表达水平显著上调。在三个发育时期全绿和全白叶片间显著差异表达的lncRNA可能是红苞凤梨嵌合性状形成的关键调控因子。
同一发育时期的全绿和全白叶片差异表达的lncRNA数量和mRNA数量见图4：B。在未展叶期，差异表达lncRNA共476个，其中在全白叶片中上调表达的lncRNA 192个，占比约40%;差异表达的mRNA共3 911个，其中上调表达的mRNA 2 152个，占比55%。在4～5叶期，差异表达lncRNA共397个，其中在全白叶片中上调表达的lncRNA 216個，占比约54%;差异表达的mRNA共2 300个，其中上调表达的1 036个，占比45%。在10～12叶期，差异表达lncRNA共594个，其中在全白叶片中上调表达的lncRNA 452个，占比约76%;差异表达的mRNA共2 100个，其中上调表达的856个，占比约41%。可以看出，随着植株的生长，lncRNA在全绿植株与全白植株间的差异表达越发显著，在全白植株中上调表达的lncRNA显著增加，lncRNA的差异表达可能在绿、白叶色的差异形成过程中起着重要的调控作用。而这个过程中，差异表达的mRNA数量在减少，上调表达的基因占比也在减少。
2.5 红苞凤梨lncRNA靶基因预测
lncRNA调控其靶基因的方式有两种，分为cis作用和trans作用。根据cis作用，我们将lncRNA 100 kb范围内的邻近蛋白质编码基因为其靶基因，经过分析预测，5 441个lncRNA通过cis作用方式预测到靶基因;trans作用是指lncRNA与mRNA由于碱基互补配对而产生作用，LncTar（Li et al.，2015）正是利用lncRNA和mRNA之间存在的互补配对关系进行预测，通过计算配对位点自由能和标准化自由能，标准化自由能阈值以下的则认为是lncRNA的靶基因。1 544个lncRNA通过trans方式预测到靶基因。靶基因的预测可以帮助理解lncRNA的功能，揭示lncRNA在红苞凤梨生长发育过程中的调控作用。
2.6 红苞凤梨差异表达lncRNA靶基因功能注释和富集分析
2.6.1 差异表达lncRNA顺式靶基因功能注释和富集分析 对差异表达lncRNA的顺式靶基因进行了COG、GO、KEGG、KOG、NR和Swiss-Prot功能富集分析，富集分析结果如表2所示。
GO数据库是一个结构化的标准生物学注释系统，在GO分析中，基因注释在三个层次上，即生物过程（biological process）、分子功能（molecular function）和细胞组分（cellular component）（图5）。在未展叶期，差异表达lncRNA的靶基因主要富集在生物过程中的生物相、节律过程和移动，细胞组分方面的细胞外基质和类核，分子功能方面的营养库活性、蛋白质结合转录因子活性、鸟苷酸交换因子活性。在4～5叶期，差异表达lncRNA的靶基因主要富集在细胞组分的类核和分子功能的蛋白质结合转录因子活性、鸟苷酸交换因子活性。在10～12叶期，差异表达lncRNA的靶基因主要富集在生物过程中的生物粘附、节律过程和移动，分子功能的营养库活性、蛋白质结合转录因子活性、鸟苷酸交换因子活性。
在生物体内，不同的基因产物相互协调以执行生物学功能，在GO分析基础上，对差异表达lncRNA反式靶基因的信号通路注释分析能够进一步解读基因的功能（图7）。在未展叶期，全绿叶片与全白叶片间的差异表达lncRNA的靶基因主要富集在核糖体、碳代谢、氧化磷酸化、淀粉和糖代谢、氨基酸代谢、脂类代谢等这些基础代谢和植物激素信号转录通路上。同时，在卟啉和叶绿素代谢途径中也富集了7个差异表达基因。在4～5叶期，差异表达lncRNA的靶基因主要富集在碳代谢、氧化磷酸化、氨基酸代谢、核糖体这些基础代谢上。同时，也富集在植物激素信号转导、嘌呤代谢、内质网蛋白过程、RNA降解、光合作用、淀粉和糖代谢方面。在10～12叶期，差异表达lncRNA的靶基因主要富集在核糖体、 碳代谢、 氨基酸代谢、氧化磷酸化、内质网蛋白过程、嘌呤代谢、RNA降解、RNA转运等方面，在卟啉和叶绿素合成代谢途径、光合作用途径中也富集了差异表达基因。与二代测序结果（蔺珍，2019）对比发现，联合测序分析增加的靶基因主要富集在碳代谢、氨基酸代谢、氧化磷酸化、RNA降解等通路上。差异表达基因的富集结果充分说明lncRNA参与了红苞凤梨叶片色素的合成、光合作用、物质代谢、生长发育调控等生理过程的调控。
2.6.2 差异表达lncRNA反式靶基因功能注释和富集分析 对差异表达lncRNA反式靶基因进行COG、GO、KEGG、KOG、NR和Swiss-Prot功能富集分析，富集分析结果如表3所示。
功能注释和富集分析结果表明，差异表达lncRNA反式靶基因功能注释到的基因数量不多，在叶片发育过程中，全绿和全白叶片间差异表达lncRNA顺式靶基因主要富集在核糖体、碳代谢、氨基酸的生物合成、氧化磷酸化、内质网中的蛋白质加工以及淀粉与蔗糖代谢等方面;而差异表达lncRNA反式靶基因主要富集在TCA循环、淀粉和糖代谢、氨基糖和核糖代谢、RNA降解、氨基酸代谢等方面。两类lncRNA在不同的生理代谢过程中发挥着调控作用。
3 讨论与结论
植物叶色嵌合体嵌合性状明显，嵌合方式多样且易于观察，是研究植物生长发育以及遗传育种的优良材料。对植物叶色嵌合体的研究，己经成为嵌合体研究的重要方向，深入研究植物嵌合体形成机理，对于植物生长发育过程中细胞间的相互作用、植物嵌合性状稳定繁殖以及植物遗传育种具有十分重要的意义。红苞凤梨叶花果具有艳丽的颜色，是重要的新型观赏植物，是研究叶色镶嵌形成机理的理想材料。叶肉细胞的白化突变是金边嵌合叶色形成的基础和关键环节，是多基因协同作用的结果。嵌合性状的形成与稳定，与基因表达的有序调控是密切相关的。前期研究表明，转录后调控在红苞凤梨金边嵌合性状形成过程中具有重要作用，而lncRNA可在表观遗传调控、转录调控、转录后调控等水平发挥功能，构建红苞凤梨lncRNA表达谱，揭示lncRNA的调控机制对于明确红苞凤梨金边嵌合性状形成机理具有重要意义。由于红苞凤梨基因组信息的缺乏，lncRNA的鉴定分析以二代转录组测序数据为基础，以菠萝基因组作为参考基因组进行分析。红苞凤梨与菠萝为同属同种不同变种植物，亲缘关系很近（Bartholomew et al.，2003），采用菠萝基因组数据作为参考基因组能有效鉴定红苞凤梨的lncRNA，为红苞凤梨的表观遗传调控研究提供基础数据。而二代测序技术由于其读长短，不能提供完整的转录本（Koren et al.，2012），很难正确预测基因的结构（Coghlan et al.，2008）。SMRT全长转录组测序克服了二代测序读长短的缺点，是研究基因结构、基因功能和比较基因组学的基本方法（Sharon et al.，2013;Luo et al.，2017）。本研究采用SMRT全长转录组测序数据和二代测序数据混合组装修正后进行lncRNA分析，以提高lncRNA分析的准确性。 经混合组装校正后的Clean Reads比对到菠萝参考基因组，比对效率达到67.74%～78.58%，比二代数据分析提高了约5%，有效地提高了长链非编码数据的利用率（蔺珍，2019）。经CPC分析、CNCI分析、CPAT分析、pfam蛋白结构域分析，共鉴定得到6 018个lncRNA，比之前鉴定得到的lncRNA数量提高了约70%，极大地丰富了红苞凤梨lncRNA数据库，为红苞凤梨表观遗传调控的研究提供了基础数据。在鉴定得到的6 018个lncRNA中，最多的是intergenic lncRNA，占比约55%，比之前结果提高了约2倍。而sense lncRNA的数量极显著下降，由之前的75%下降到18.5%（蔺珍，2019）。本次分析结果中4种lncRNA的数量分布与玉米（Zea mays）（Wang et al.，2016）基本相似，修正了之前分析结果中，sense lncRNA比例太高的情况。SMRT全长转录组测序数据与lncRNA二代转录组测序数据的联合分析，有效地提高了lncRNA测序数据的利用效率，很大程度上弥补了以近缘种菠萝基因组作为参考基因组鉴定lncRNA所带来的缺陷，鉴定得到的lncRNA数量极显著提高，丰富了红苞凤梨lncRNA数据库信息。联合分析结果中，修正了4种lncRNA数量分布的异常，提高了lncRNA鉴定结果的准确性，有效地提高了后期相关研究的可行性和可靠性。lncRNA和mRNA序列结构对比分析表明，与编码基因相比，lncRNA具有表达丰度较低、序列长度较短、外显子数目较少、ORF长度较短等结构特点。这与斑马鱼（Brachydanio rerio）（高霄霄，2017）、臭橘（Poncirus trifoliata）（Wang et al.，2017）、杨树（Populus）（田净净，2016）的分析结果一致，说明lncRNA的这些结构特征在生物中具有普遍性。
lncRNA靶基因的预测以及靶基因的功能注释和富集分析，是研究lncRNA功能的重要途径（蔺珍，2019）。本研究中GO和KEGG联合分析结果表明，在叶片发育的三个时期，差异表达lncRNA的靶基因注释到了叶片生长发育的多个方面，包括了碳代谢、氨基酸代谢、脂代谢、核糖体、淀粉和糖代谢等这些基础代谢方面，也包括了植物激素信号转导等调控机制方面。在未展叶期，差异表达lncRNA的靶基因注釋到了卟啉和叶绿素代谢途径，说明在叶片发育的早期，全绿和全白叶片的叶绿素合成代谢就已出现差异，从而导致了叶片颜色的差异。在4～5叶期和10～12叶期，光合作用途径富集了差异表达基因，说明由于叶片的失绿白化，光合作用受到抑制，相关基因差异表达。同时，由于光合作用能力的差异，全绿和全白叶片在各个生理代谢方面都表现出差异，在很多基础代谢、调控途径、核酸代谢等方面都富集了差异表达基因。这些基因表达模式说明，lncRNA对靶基因表达的调控作用可能在红苞凤梨叶片发育早期叶绿素合成代谢差异以及后期光合作用差异和其他生理代谢差异方面都起着重要的调控作用。本研究鉴定出的相关差异表达lncRNA信息，为进一步研究lncRNA对靶基因的调控作用及靶基因在红苞凤梨叶片失绿白化突变机制，以及全面揭示红苞凤梨金边嵌合体形成的机理提供了重要的数据基础。
参考文献：
BARTHOLOMEW DP， PAULL RE， ROHRBACH KG， 2003， The pineapple： Botany， production and uses [M]. Wallingford： CABI Publishing.
CAO L， 2011. A study on in vitro culture of chimera cultivars of Ananas bracteatus schultes and their stability of chimeric traits [D]. Guangzhou： South China Agricultural University： 1-57. [曹莉， 2011. 红苞凤梨嵌合体品种离体培养及其稳定性的研究 [D]. 广州： 华南农业大学： 1-57.]
CSORBA T， QUESTA JI， SUN Q， et al.， 2014. Antisense COOLAIR mediates the coordinated switching of chromatin states at FLC during vernalization [J]. Proc Natl Acad Sci USA， 111（45）： 16160-16165.
COGHLAN A， FIEDLER TJ， MCKAY SJ， et al.， 2008. nGASP—the nematode genome annotation assessment project [J]. BMC Bioinformatics， 9： 549DOI 10.1186/1471-2105-9-549.
DING JH， LU Q， OUYANG YD， et al.， 2012. A long noncoding RNA regulates photoperiod-sensitive male sterility， an essential component of hybrid rice [J]. Proc Natl Acad Sci USA， 109（7）： 2654-2659.
ENGLISH AC， RICHARDS S， HAN Y， et al.， 2012. Mind the gap： Upgrading genomes with Pacific biosciences RS long-read sequencing technology [J]. PLoS ONE， 7（11）： e47768. FLUSBERG BA， WEBSTER DR， LEE JH， et al.， 2010. Direct detection of DNA methylation during single-molecule， real-time sequencing [J]. Nat Methods， 7（6）： 461-465.
FINN RD， BATEMAN A， CLEMENTS J， et al.， 2014. Pfam： The protein families database [J]. Nucl Acid Res， 42（Database issue）： D222-D230.
GUO F， WANG D， WANG LS， 2018. Progressive approach for SNP calling and haplotype assembly using single molecular sequencing data [J]. Bioinformatics， 34（12）： 2012-2018.
GAO XX， 2017. Screening and identification of long noncoding RNAs in thepubertal female goats [D]. Hefei： Anhui Agricultural University： 1-39. [高霄霄， 2017. 初情期雌性山羊 lncRNA 的篩选与鉴定 [D]. 合肥： 安徽农业大学： 1-39.]
HACKL T， HEDRICH R， SCHULTZ J， et al.， 2014. Proovread： Large-scale high-accuracy PacBio correction through iterative short read consensus [J]. Bioinformatics， 30（21）：3004-3011.
JOHNSON R， 2012. Long non-coding RNAs in Huntington’s disease neurodegeneration [J]. Neurobiol Disease， 46（2）：245-254.
KOREN S， WALENZ BP， BERLIN K， et al.， 2017. Canu： Scalable and accurate long-read assembly via adaptive k-mer weighting and repeat separation [J]. Genome Res， 27（5）：722-736.
KOREN S， SCHATZ MC， WALENZ BP， et al.， 2012. Hybrid error correction and de novo assembly of single-molecule sequencing reads [J]. Nat Biotechnol， 30（7）： 693-700.
KIM D， PERTEA G， TRAPNELL C， et al.， 2013. TopHat2： Accurate alignment of transcriptomes in the presence of insertions， deletions and gene fusions [J]. Genome Biol， 14（4）： R36.
KONG L， ZHANG Y， YE ZQ， et al.， 2007. CPC： Assess the protein-coding potential of transcripts using sequence features and support vector machine [J]. Nucl Acid Res， 35（Web Server issue）： W345-W349.
KOREN S， SCHATZ MC， WALENZ BP， et al.， 2012. Hybrid error correction and de novo assembly of single-molecule sequencing reads [J]. Nat Biotechnol， 30（7）： 693-700 DOI 10.1038/nbt.2280.
LIN Z， 2019. Identification of A. comosus var. bracteatus lncRNAs and functional verification of lncABCG11 [D]. Ya’an： Sichuan Agricultural University： 1-113 [蔺珍， 2019. 红苞凤梨lncRNAs的鉴定及lncABCG11的功能验证 [D]. 雅安： 四川农业大学： 1-113.]
LI X， KANAKALA S ， HE YH， et al.， 2017. Physiological characterization and comparative transcriptome analysis of white and green leaves of Ananas comosus var. bracteatus [J]. PLoS ONE， 12（1）： e0169838.
LANGMEAD B， TRAPNELL C， POP M， et al.， 2009. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome [J]. Genome Biol， 10（3）： R25. LI JW， MA W， ZENG P， et al.， 2015. LncTar： A tool for predicting the RNA targets of long noncoding RNAs [J]. Brief Bioinform， 16（5）： 806-812.
LUO YH， DING N， SHI X， et al.， 2017. Generation and comparative analysis of full-length transcriptomes in sweetpotato and its putative wild ancestor I. trifida [J]. BioRxiv， https：//doi.org/10.1101/112425
MA DN， ZHANG XT， WEI LF， et al.， 2018. Benchmarking hybrid correction and assembly using short Illumina reads and long pac bio reads [J]. Genom Appl Biol， 37（4）：1547-1555. [马东娜， 张兴坦， 魏柳锋， 等， 2018. 基因组二代测序数据与三代测序数据的混合校正和组装 [J]. 基因组学与应用生物学， 37（4）： 1547-1555.]
MA J， XIANG YX， XIONG YY， et al.， 2018. SMRT sequencing analysis reveals the full-length transcripts and alternativesplicing patterns in Ananas comosus var.bracteatus [J]. Peer J， 7： e7062
MA JC， 2018. Genome sequence of a widely cultivated poplar and its lnc RNAs response to salt stress [D].Lanzhou： Lanzhou University： 1-82. [马建超， 2018. 新疆杨基因组及其lncRNA响应盐胁迫的研究 [D]. 兰州： 兰州大学： 1-82.]
QIN T， ZHAO HY， CUI P， et al.， 2017. A nucleus-localized long non-coding RNA enhances drought and salt stress tolerance [J]. Plant Physiol， 175（3）： 1321-1336.
SMITH CC， WANG Q， CHIN CS， et al.， 2012. Validation of ITD mutations in FLT3 as a therapeutic target in human acute myeloid leukaemia [J]. Nature， 485（7397）： 260-263.
SUN L， LUO HT， BU DC， et al.， 2013. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts [J]. Nucl Acid Res， 41（17）： e166.
SHARON D， TILGNER H， GRUBERT F， et al.， 2013. A single-molecule long-readsurvey of the human transcriptome [J]. Nat Biotechnol， 31：1009-1014. DOI： 10.1038/nbt.2705.
TRAPNELL C， WILLIAMS BA， PERTEA G， et al.， 2010. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation [J]. Nat Biotechnol， 28（5）： 511-515.
TIAN JJ， 2016. The application of CRISPR/Cas9 system in the zebrafish gene editing [D]. Yangzhou： Yangzhou University： 1-68. [田凈净， 2016CRISPR/Cas9 系统在斑马鱼基因编辑研究中的应用 [D]. 扬州： 扬州大学： 1-68.]
WANG YF，SU WY， ZHANG L， et al.， 2018. Advances of long non-coding RNA in plants [J]. Acta Bot Boreal-Occident Sin， （3）： 582-588. [王艳芳， 苏婉玉， 张琳， 等， 2018. 植物中lncRNAs的研究进展 [J]. 西北植物学报， （3）： 582-588.]
WANG LG， PARK HJ， DASARI S， et al.， 2013. CPAT： Coding-potential assessment tool using an alignment-free logistic regression model [J]. Nucl Acid Res， 41（6）： e74. WANG B ， TSENG E ， REGULSKI M ， et al.， 2016. Unveiling the complexity of the maize transcriptome by single-molecule long-read sequencing [J]. Nat Comm， 7： 11708.
WANG CY， LIU SR， ZHANG XY， et al.， 2017. Genome-wide screening and characterization of long non-coding RNAs involved in flowering development of trifoliate orange （Poncirus trifoliata L. Raf.） [J]. Sci Rep， 7： 43226.
XIONG YY， MA J， HE YH ， et al.， 2018. High-throughput sequencing analysis revealed the regulation patterns of small RNAs on the development of Ananas comosus var. bracteatus leaves [J]. Sci Rep， 8（1）：1947.
XIONG YY， 2019. MicroRNAs identification and screening and functional verification of key microRNAs involved in the albino of Ananas comosus var. bracteatus [D]. Ya’an： Sichuan Agricultural University： 1-95. [熊颖媛， 2019. 红苞凤梨microRNA的鉴定及白化关键microRNA的筛选与功能验证 [D]. 雅安： 四川农业大学： 1-95.]
XUE YB， MA J， HE YH， et al.， 2019. Comparative transcriptomic and proteomic analyses of the green and white parts of chimeric leaves in Ananas comosus var. bracteatus [J]. PeerJ， 7： e7261.
XU WN， HUANG RM， LIU YY， et al.， 2018. Genome sequencing and assembly strategy analyses of Flammulina filiformis [J]. Mycosystema， 37（12）： 1578-1585. [徐偉南， 黄蓉梅， 刘媛媛， 等， 2018. 金针菇基因组测序与组装策略分析 [J]. 菌物学报， 37（12）： 1578-1585.]
YU CL， LUO L， LIAO Q， 2015. Annotation and functional prediction of lncRNAs [J]. Chin J Biochem Mol Biol， （3）： 239-243. [余铖亮， 骆亮， 廖奇， 2015. lncRNAs功能注释和预测 [J]. 中国生物化学与分子生物学报， （3）： 239-243.]
ZHANG YW， TAO Y， LIAO Q， 2018. Long noncoding RNA： A crosslink in biological regulatory network [J]. Brief Bioinform， 19（5）： 930-945.
（责任编辑 李 莉）
关键词：红苞凤梨，Hiseq2500测序，SMRT全长转录组测序，lncRNA鉴定
中图分类号：Q943
文献标识码：A
文章编号：1000-3142（2021）08-1237-14
Abstract： In order to reveal the regulation function of lncRNA on the chimeric character formation and development of the leaf of Ananas comosus var. bracteatus， Hiseq2500 sequencing and SMRT the third-generation full-length transcriptome sequencing were applyed to identify lncRNA of Ananas comosus var. bracteatus. The results were as follows： （1） A total of 6 018 lncRNA were identified， containing 3 298 intergenic lncRNA， 870 antisense lncRNA， 717 intron lncRNA and 1 109 sense lncRNA， which were greatly improved compared with the second-generation information. （2） Structural analysis showed that the overall expression level of lncRNA was lower than that of mRNA. The transcript length distribution of lncRNA in the range of 400-1 200 nt was higher than that of mRNA， while in the range > 1 600 nt， the proportion of lncRNA distribution was significantly lower than that of mRNA. The number of exons in lncRNA was generally less than that of mRNA， and the open reading frame was also shorter in length than that of mRNA. （3） For analysis of differential expression， 1 710 differentially expressed lncRNA were identified during the development of complete green and complete white leaves. （4） Target gene prediction results showed that 5 441 lncRNA were predicted target genes by Cis action， and 1 544 lncRNA were predicted target genes by Trans action. （5） Functional annotation and enrichment analysis of target genes revealed that the target genes of differentially expressed lncRNA mainly act on metabolic activities and signal transduction of leaves as enzyme proteins， and were closely related to leaf color formation， photosynthesis and leaf growth. The lncRNA information identified in this paper and， as well as the analysis of its structure and functions， provide the data basis for the study of the epigenetic regulation mechanism of lncRNA in Ananas comosus var. bracteatus and other plants in Bromeliaceae. The identified differentially expression of lncRNA plays an important role in the chimeric character formation and development of leaf of Ananas comosus var. bracteatus. Key words： Ananas comosus var. bracteatus， Hiseq2500， SMRT full length transcriptome sequence， lncRNA identification
红苞凤梨 （Ananas comosus var. bracteatus）因其叶片绿白镶嵌、花果颜色艳丽且观赏期长，己成为一种重要的新型观赏植物。红苞凤梨自交不亲和，在生产中以吸芽进行繁殖，繁殖系数低、苗木整齐度差，限制了红苞凤梨的规模化应用。组织培养能快速繁殖红苞凤梨，但繁殖过程中叶片嵌合性状不稳定，再生植株叶片常失去嵌合性状而变为全绿植株（曹莉，2011）。细胞白化突变是金边嵌合叶片形成的基础，研究红苞凤梨叶片细胞白化突变的分子机理，对揭示红苞凤梨嵌合性状的形成机理，提高嵌合性状的稳定性，培育新的嵌合性状品种具有重要的理论和实践意义。
本课题组前期研究结果表明，红苞凤梨白化细胞中叶绿素含量极显著下降，但叶绿素合成代谢的结构基因表达上调（Li et al.，2017;Xue et al.，2019），说明转录后调控在红苞凤梨细胞白化失绿、金边嵌合性状形成中发挥了重要调控作用。lncRNA具有类似mRNA的结构特征，可在多个层面调控目标基因的表达。可作为信号分子、诱饵分子、引导分子以及支架分子，在表观遗传、转录调控、转录后调控等多个水平发挥功能（Zhang et al.， 2018）。目前，lncRNA在人类和动物中研究较为广泛，与人类疾病的发生及生物体的生长发育密切相关（Johnson，2012;余铖亮等，2015;王艳芳等，2018），而植物lncRNA的研究还处于起步阶段。研究结果表明lncRNA在植物的开花诱导（Csorba et al.，2014）、花粉发育（Ding et al.，2012）、逆境脅迫（Qin et al.，2017）中具有重要功能，然而其具体作用机制及调节功能等尚不清楚。由于红苞凤梨没有基因组数据，且Hiseq2500二代测序技术读长较短，使得红苞凤梨非编码RNA的挖掘具有一定的局限性。随着测序技术的发展，SMRT（single-molecule real-time）三代测序技术的出现，因其无需进行PCR扩增，大幅降低了因PCR反应引入的碱基错误，操作更为简单等优点，已得到广泛应用（Flusberg et al.，2010）。目前，SMRT三代测序技术在基因组、甲基化识别、SNP 的鉴定、基因重测序和转录组学等方面的优势越来越明显（Smith et al.，2012;Guo et al.，2018）。此外，SMRT测序技术得到几千kb的数据，读长显著增长，大大减少了测序后的 Contig数量，使得基因组和转录组的组装得到极大改善（English et al.，2012）。但是，第三代测序技术也存在一定的缺点，测序错误率普遍偏高，测序产生的错误率可高达15%（Koren et al.，2017）。因此，采取二代测序和三代测序联合分析已成为当下基因组研究的主要方法。一般有两种联合手段，可以选择以三代测序为主，进行组装，再通过二代测序得到的高质量短片段对三代的数据进行碱基纠错和矫正;也可以选择以二代测序为主，用三代测序得到的长片段Reads 进行辅助组装（马建超，2018）。最为通常的手段是利用短的但准确度高的二代数据去辅助校正长的但准确度较低的三代数据，并在此基础上进行混合组装。经验证，通过这种混合组装校正方法得到的数据准确性可达99%（马东娜等，2018）。因此，这种“2 + 3”的联合模式已被广泛认可并应用于动物、植物和微生物的基因组研究（Koren et al.，2012;Hackl et al.，2014;徐伟南等，2018）。
本研究以红苞凤梨为材料，利用二代和三代测序数据混合组装校正方法，鉴定红苞凤梨叶片中存在的lncRNA，分析在全绿和全白突变叶片生长发育过程中差异表达的lncRNA，通过靶基因的功能注释和富集分析，揭示lncRNA在红苞凤梨叶片细胞失绿白化及生长发育过程中的作用。本研究结果为红苞凤梨lncRNA表观遗传调控机理研究提供数据基础。
1 材料与方法
1.1 试验材料
本研究以金边红苞凤梨茎段为外植体，通过组织培养获得的全白和全绿植株为试验材料。选取长势均一的未展叶期、4～5叶期和10～12叶期三个发育时期的全绿及全白植株各10株的叶片作为样本提取RNA（图1）。样本取下后立即转入液氮速冻，后储存于-80 ℃用于提取RNA进行Hiseq2500二代转录组测序。
1.2 红苞凤梨RNA提取及二代Illumina Hiseq2500测序文库构建
1.2.1 红苞凤梨RNA提取及检测 取-80 ℃冻存的红苞凤梨样品，采用LABGENETM plant RNA Isolation kit多糖多酚植物RNA提取试剂盒分离总RNA，操作方法参考使用说明。高质量的RNA是实验成功的基础，为保证测序的准确性，对样品RNA进行了检测，分别采用荧光定量仪（Qubit 2.0）、微量分光光度计（Nanodrop）、生物分析仪（Aglient 2100）、电泳方法检测RNA样品的纯度、浓度、完整性和是否有基因组DNA污染等，达到质控要求的RNA样品用于后续试验。
1.2.2 Illumina Hiseq2500 cDNA文库构建及上机测序 利用epicentre Ribo-ZeroTM试剂盒去除样品中的rRNA。以rRNA-depleted RNA为模板，用六碱基随机引物（random hexamers）合成cDNA第一条链和第二条链。cDNA纯化后进行末端修复、加A尾并连接测序接头，然后用AMPure XP beads进行片段大小选择。降解含U链，通过PCR富集得到cDNA文库。文库构建后，经Qubit 2.0定量，Agilent 2100检测文库的insert size质量;釆用QPCR对文库的浓度（文库有效浓度>2nmol·L-1）进行准确定量，完成库检。库捡完成后，采用Illumina Hiseq2500平台对文库进行测序。 1.3 数据质量监控
测序所得Reads的过滤与修剪是保证分析数据可靠性的关键，结合前期SMRT测序结果（Ma et al.，2018）进行混合组装校正后，删除Raw Reads中包含adapter、ploy-N的Reads和低质量的Reads以获得Clean Reads。以菠萝基因组（Acomosus_321_v3https：//phytozome.jgi.doe.gov/pz/portal.html#！info？alias=Org_Acomosus_er）为参考基因组，利用TopHat v2.0.9（Kim et al.，2013）软件对Clean Reads进行比对。采用Scripture（beta2）（Langmead et al.，2009）和Cuffiinks（v2.1.1）软件对每个样品比对上的Reads进行装配获得转录本。
1.4 轉录本表达水平和编码潜能分析
使用Cufflinks软件的Cuffdiff组件，对转录本表达水平进行分析。根据lncRNA编码的特点进行基本筛选：选择长度≥200 bp，Exon个数≥2及FPKM≥0.1的转录本。
因lncRNA不编码蛋白，因此，通过对基本筛选得到的转录本进行编码潜能筛选，判断其是否具有编码潜能，从而可以判定该转录本是否为lncRNA。主要利用CPC分析、CNCI分析、CPAT分析、Pfam蛋白结构域四种分析方法分析lncRNA的编码能力。去掉基本筛选中具有潜在编码能力的转录本，余下的即为预测的lncRNA。
通过与已知mRNA进行比较，利用Cuff-compare分析结果中的class codes对筛选的lncRNA进行分类。
1.5 lncRNA靶基因预测及功能富集分析
基于lncRNA与其靶基因的作用方式（cis和trans），采用两种预测方法：第一种是根据lncRNA与mRNA的位置关系预测lncRNA的靶基因，即lncRNA 100 kb范围内的邻近基因为其靶基因;第二种是使用LncTar靶基因预测工具进行预测，其原理是利用lncRNA与mRNA之间碱基互补配对产生的作用来预测。
对差异表达IncRNA的靶基因利用KEGG（Kyoto Encyclopedia of Genesand Genomes）、GO（gene ontology）、NR （Non-Redundant Protein Sequences Database）、COG （clusters of orthologous groups of proteins）和 Swiss-Prot 数据库进行功能注释和富集分析，分析结果显著性用P值表示。
1.6 差异表达基因分析
以差异倍数（Fold Change）≥2且错误发现率FDR（False Discovery Rate）<0.05作为差异表达的筛选标准，使用EBseq分析6个样本间差异表达的lncRNA及mRNA。
2 结果与分析
2.1 测序数据与参考基因组比对效率分析
长链非编码数据利用率的直接体现就是比对效率，即Mapped Reads占Clean Reads的百分比。经SMRT全长转录组数据（NCBI提交号PRJNA564223）（Ma et al.，2018）修正后，6个样品的Reads与所选参考基因组的比对效率为67.74%～78.58%，比单独使用Hiseq2500二代测序数据的比对效率提高了约5%（蔺珍，2019）（表1）。说明三代数据的修正，有效提高了lncRNA测序数据的利用效率，有利于进一步深入挖掘红苞凤梨lncRNA信息。
2.2 红苞凤梨lncRNA的鉴定
对测序样品Cuffiinks的拼接结果第一步经过Cuffcompare软件分析合并，筛选出转录本长度≥200 bp，外显子个数≥2，FPKM ≥ 0.1的转录本，根据已知mRNA数据库，进行比较过滤mRNA，然后通过CNCI、CPC、CPAT 和Pfam四个软件进行蛋白编码潜能筛选，最终鉴定得到6 018条lncRNA，其中新lncRNA 5 689条（图2）。比利用Hiseq二代测序数据鉴定到的lncRNA数量增加了约70%（蔺珍，2019），极大地提高了红苞凤梨lncRNA的数据信息量，为进一步研究红苞凤梨的非编码RNA调控机理提供了数据基础。
鉴定得到的6 018个lncRNA中，包括了3 298个基因间lncRNA（intergenic lncRNA），870个反义lncRNA（antisense lncRNA），717个内含子lncRNA（intronic lncRNA）和1 109个正义lncRNA（sense lncRNA）（图2：B）。与二代测序分析鉴定结果相比，intergenic lncRNA所在比例显著提高，由17%增加到55%。而sense lncRNA所占比例显著下降，由76%下降到18.5%（蔺珍，2019）。
2.3 红苞凤梨lncRNA结构分析
为了进一步分析红苞凤梨lncRNA的结构特点，将lncRNA与蛋白质编码RNA在整体表达水平、序列长度分布、外显子数目分布以及开放阅读框长度分布情况进行了比较分析（图3）。分析结果表明，在总体表达水平上，mRNA的表达丰度高于lncRNA的表达丰度（图3：A）。在转录本的长度分布上，lncRNA在400～1 200 nt区间分布比例高于mRNA（图3：B），二代测序结果中则是表现在400～600 nt和1 400～1 600 nt区间（蔺珍，2019）;而在转录本长度>1 600 nt区间，lncRNA分布的比例显著小于mRNA，尤其在>=3 000 nt区间（图3：B）。lncRNA中的外显子数量总体少于mRNA，约82%的lncRNA只含有2个外显子（图3：C）。而二代测序分析中有41.80%的lncRNA只有2个外显子，而 31.62%的mRNA的外显子数超过 5 个（蔺珍，2019）。lncRNA的开放阅读框长度总体上也短于mRNA。约99%的lncRNA开放阅读框长度<=100 nt（图3：D），而二代测序分析中66%的lncRNA的开放阅读框在 0～100 nt之间（蔺珍，2019）。 2.4 红苞凤梨lncRNA的差异表达分析
以Fold Change≥2.0且FDR<0.05作为筛选标准，共鉴定得到了1 710个差异表达的lncRNA。对筛选出的差异表达lncRNA进行了层次聚类分析（图4：A）。在未展叶期，全绿和全白叶片中，大量差异表达的lncRNA表达丰度较高，在全白苗和全绿苗的发育到4～5叶期时，多数差异表达的lncRNA表达水平下降。而全白苗在第三发育时期有部分差异表达lncRNA表达水平显著上调。在三个发育时期全绿和全白叶片间显著差异表达的lncRNA可能是红苞凤梨嵌合性状形成的关键调控因子。
同一发育时期的全绿和全白叶片差异表达的lncRNA数量和mRNA数量见图4：B。在未展叶期，差异表达lncRNA共476个，其中在全白叶片中上调表达的lncRNA 192个，占比约40%;差异表达的mRNA共3 911个，其中上调表达的mRNA 2 152个，占比55%。在4～5叶期，差异表达lncRNA共397个，其中在全白叶片中上调表达的lncRNA 216個，占比约54%;差异表达的mRNA共2 300个，其中上调表达的1 036个，占比45%。在10～12叶期，差异表达lncRNA共594个，其中在全白叶片中上调表达的lncRNA 452个，占比约76%;差异表达的mRNA共2 100个，其中上调表达的856个，占比约41%。可以看出，随着植株的生长，lncRNA在全绿植株与全白植株间的差异表达越发显著，在全白植株中上调表达的lncRNA显著增加，lncRNA的差异表达可能在绿、白叶色的差异形成过程中起着重要的调控作用。而这个过程中，差异表达的mRNA数量在减少，上调表达的基因占比也在减少。
2.5 红苞凤梨lncRNA靶基因预测
lncRNA调控其靶基因的方式有两种，分为cis作用和trans作用。根据cis作用，我们将lncRNA 100 kb范围内的邻近蛋白质编码基因为其靶基因，经过分析预测，5 441个lncRNA通过cis作用方式预测到靶基因;trans作用是指lncRNA与mRNA由于碱基互补配对而产生作用，LncTar（Li et al.，2015）正是利用lncRNA和mRNA之间存在的互补配对关系进行预测，通过计算配对位点自由能和标准化自由能，标准化自由能阈值以下的则认为是lncRNA的靶基因。1 544个lncRNA通过trans方式预测到靶基因。靶基因的预测可以帮助理解lncRNA的功能，揭示lncRNA在红苞凤梨生长发育过程中的调控作用。
2.6 红苞凤梨差异表达lncRNA靶基因功能注释和富集分析
2.6.1 差异表达lncRNA顺式靶基因功能注释和富集分析 对差异表达lncRNA的顺式靶基因进行了COG、GO、KEGG、KOG、NR和Swiss-Prot功能富集分析，富集分析结果如表2所示。
GO数据库是一个结构化的标准生物学注释系统，在GO分析中，基因注释在三个层次上，即生物过程（biological process）、分子功能（molecular function）和细胞组分（cellular component）（图5）。在未展叶期，差异表达lncRNA的靶基因主要富集在生物过程中的生物相、节律过程和移动，细胞组分方面的细胞外基质和类核，分子功能方面的营养库活性、蛋白质结合转录因子活性、鸟苷酸交换因子活性。在4～5叶期，差异表达lncRNA的靶基因主要富集在细胞组分的类核和分子功能的蛋白质结合转录因子活性、鸟苷酸交换因子活性。在10～12叶期，差异表达lncRNA的靶基因主要富集在生物过程中的生物粘附、节律过程和移动，分子功能的营养库活性、蛋白质结合转录因子活性、鸟苷酸交换因子活性。
在生物体内，不同的基因产物相互协调以执行生物学功能，在GO分析基础上，对差异表达lncRNA反式靶基因的信号通路注释分析能够进一步解读基因的功能（图7）。在未展叶期，全绿叶片与全白叶片间的差异表达lncRNA的靶基因主要富集在核糖体、碳代谢、氧化磷酸化、淀粉和糖代谢、氨基酸代谢、脂类代谢等这些基础代谢和植物激素信号转录通路上。同时，在卟啉和叶绿素代谢途径中也富集了7个差异表达基因。在4～5叶期，差异表达lncRNA的靶基因主要富集在碳代谢、氧化磷酸化、氨基酸代谢、核糖体这些基础代谢上。同时，也富集在植物激素信号转导、嘌呤代谢、内质网蛋白过程、RNA降解、光合作用、淀粉和糖代谢方面。在10～12叶期，差异表达lncRNA的靶基因主要富集在核糖体、 碳代谢、 氨基酸代谢、氧化磷酸化、内质网蛋白过程、嘌呤代谢、RNA降解、RNA转运等方面，在卟啉和叶绿素合成代谢途径、光合作用途径中也富集了差异表达基因。与二代测序结果（蔺珍，2019）对比发现，联合测序分析增加的靶基因主要富集在碳代谢、氨基酸代谢、氧化磷酸化、RNA降解等通路上。差异表达基因的富集结果充分说明lncRNA参与了红苞凤梨叶片色素的合成、光合作用、物质代谢、生长发育调控等生理过程的调控。
2.6.2 差异表达lncRNA反式靶基因功能注释和富集分析 对差异表达lncRNA反式靶基因进行COG、GO、KEGG、KOG、NR和Swiss-Prot功能富集分析，富集分析结果如表3所示。
功能注释和富集分析结果表明，差异表达lncRNA反式靶基因功能注释到的基因数量不多，在叶片发育过程中，全绿和全白叶片间差异表达lncRNA顺式靶基因主要富集在核糖体、碳代谢、氨基酸的生物合成、氧化磷酸化、内质网中的蛋白质加工以及淀粉与蔗糖代谢等方面;而差异表达lncRNA反式靶基因主要富集在TCA循环、淀粉和糖代谢、氨基糖和核糖代谢、RNA降解、氨基酸代谢等方面。两类lncRNA在不同的生理代谢过程中发挥着调控作用。
3 讨论与结论
植物叶色嵌合体嵌合性状明显，嵌合方式多样且易于观察，是研究植物生长发育以及遗传育种的优良材料。对植物叶色嵌合体的研究，己经成为嵌合体研究的重要方向，深入研究植物嵌合体形成机理，对于植物生长发育过程中细胞间的相互作用、植物嵌合性状稳定繁殖以及植物遗传育种具有十分重要的意义。红苞凤梨叶花果具有艳丽的颜色，是重要的新型观赏植物，是研究叶色镶嵌形成机理的理想材料。叶肉细胞的白化突变是金边嵌合叶色形成的基础和关键环节，是多基因协同作用的结果。嵌合性状的形成与稳定，与基因表达的有序调控是密切相关的。前期研究表明，转录后调控在红苞凤梨金边嵌合性状形成过程中具有重要作用，而lncRNA可在表观遗传调控、转录调控、转录后调控等水平发挥功能，构建红苞凤梨lncRNA表达谱，揭示lncRNA的调控机制对于明确红苞凤梨金边嵌合性状形成机理具有重要意义。由于红苞凤梨基因组信息的缺乏，lncRNA的鉴定分析以二代转录组测序数据为基础，以菠萝基因组作为参考基因组进行分析。红苞凤梨与菠萝为同属同种不同变种植物，亲缘关系很近（Bartholomew et al.，2003），采用菠萝基因组数据作为参考基因组能有效鉴定红苞凤梨的lncRNA，为红苞凤梨的表观遗传调控研究提供基础数据。而二代测序技术由于其读长短，不能提供完整的转录本（Koren et al.，2012），很难正确预测基因的结构（Coghlan et al.，2008）。SMRT全长转录组测序克服了二代测序读长短的缺点，是研究基因结构、基因功能和比较基因组学的基本方法（Sharon et al.，2013;Luo et al.，2017）。本研究采用SMRT全长转录组测序数据和二代测序数据混合组装修正后进行lncRNA分析，以提高lncRNA分析的准确性。 经混合组装校正后的Clean Reads比对到菠萝参考基因组，比对效率达到67.74%～78.58%，比二代数据分析提高了约5%，有效地提高了长链非编码数据的利用率（蔺珍，2019）。经CPC分析、CNCI分析、CPAT分析、pfam蛋白结构域分析，共鉴定得到6 018个lncRNA，比之前鉴定得到的lncRNA数量提高了约70%，极大地丰富了红苞凤梨lncRNA数据库，为红苞凤梨表观遗传调控的研究提供了基础数据。在鉴定得到的6 018个lncRNA中，最多的是intergenic lncRNA，占比约55%，比之前结果提高了约2倍。而sense lncRNA的数量极显著下降，由之前的75%下降到18.5%（蔺珍，2019）。本次分析结果中4种lncRNA的数量分布与玉米（Zea mays）（Wang et al.，2016）基本相似，修正了之前分析结果中，sense lncRNA比例太高的情况。SMRT全长转录组测序数据与lncRNA二代转录组测序数据的联合分析，有效地提高了lncRNA测序数据的利用效率，很大程度上弥补了以近缘种菠萝基因组作为参考基因组鉴定lncRNA所带来的缺陷，鉴定得到的lncRNA数量极显著提高，丰富了红苞凤梨lncRNA数据库信息。联合分析结果中，修正了4种lncRNA数量分布的异常，提高了lncRNA鉴定结果的准确性，有效地提高了后期相关研究的可行性和可靠性。lncRNA和mRNA序列结构对比分析表明，与编码基因相比，lncRNA具有表达丰度较低、序列长度较短、外显子数目较少、ORF长度较短等结构特点。这与斑马鱼（Brachydanio rerio）（高霄霄，2017）、臭橘（Poncirus trifoliata）（Wang et al.，2017）、杨树（Populus）（田净净，2016）的分析结果一致，说明lncRNA的这些结构特征在生物中具有普遍性。
lncRNA靶基因的预测以及靶基因的功能注释和富集分析，是研究lncRNA功能的重要途径（蔺珍，2019）。本研究中GO和KEGG联合分析结果表明，在叶片发育的三个时期，差异表达lncRNA的靶基因注释到了叶片生长发育的多个方面，包括了碳代谢、氨基酸代谢、脂代谢、核糖体、淀粉和糖代谢等这些基础代谢方面，也包括了植物激素信号转导等调控机制方面。在未展叶期，差异表达lncRNA的靶基因注釋到了卟啉和叶绿素代谢途径，说明在叶片发育的早期，全绿和全白叶片的叶绿素合成代谢就已出现差异，从而导致了叶片颜色的差异。在4～5叶期和10～12叶期，光合作用途径富集了差异表达基因，说明由于叶片的失绿白化，光合作用受到抑制，相关基因差异表达。同时，由于光合作用能力的差异，全绿和全白叶片在各个生理代谢方面都表现出差异，在很多基础代谢、调控途径、核酸代谢等方面都富集了差异表达基因。这些基因表达模式说明，lncRNA对靶基因表达的调控作用可能在红苞凤梨叶片发育早期叶绿素合成代谢差异以及后期光合作用差异和其他生理代谢差异方面都起着重要的调控作用。本研究鉴定出的相关差异表达lncRNA信息，为进一步研究lncRNA对靶基因的调控作用及靶基因在红苞凤梨叶片失绿白化突变机制，以及全面揭示红苞凤梨金边嵌合体形成的机理提供了重要的数据基础。
参考文献：
BARTHOLOMEW DP， PAULL RE， ROHRBACH KG， 2003， The pineapple： Botany， production and uses [M]. Wallingford： CABI Publishing.
CAO L， 2011. A study on in vitro culture of chimera cultivars of Ananas bracteatus schultes and their stability of chimeric traits [D]. Guangzhou： South China Agricultural University： 1-57. [曹莉， 2011. 红苞凤梨嵌合体品种离体培养及其稳定性的研究 [D]. 广州： 华南农业大学： 1-57.]
CSORBA T， QUESTA JI， SUN Q， et al.， 2014. Antisense COOLAIR mediates the coordinated switching of chromatin states at FLC during vernalization [J]. Proc Natl Acad Sci USA， 111（45）： 16160-16165.
COGHLAN A， FIEDLER TJ， MCKAY SJ， et al.， 2008. nGASP—the nematode genome annotation assessment project [J]. BMC Bioinformatics， 9： 549DOI 10.1186/1471-2105-9-549.
DING JH， LU Q， OUYANG YD， et al.， 2012. A long noncoding RNA regulates photoperiod-sensitive male sterility， an essential component of hybrid rice [J]. Proc Natl Acad Sci USA， 109（7）： 2654-2659.
ENGLISH AC， RICHARDS S， HAN Y， et al.， 2012. Mind the gap： Upgrading genomes with Pacific biosciences RS long-read sequencing technology [J]. PLoS ONE， 7（11）： e47768. FLUSBERG BA， WEBSTER DR， LEE JH， et al.， 2010. Direct detection of DNA methylation during single-molecule， real-time sequencing [J]. Nat Methods， 7（6）： 461-465.
FINN RD， BATEMAN A， CLEMENTS J， et al.， 2014. Pfam： The protein families database [J]. Nucl Acid Res， 42（Database issue）： D222-D230.
GUO F， WANG D， WANG LS， 2018. Progressive approach for SNP calling and haplotype assembly using single molecular sequencing data [J]. Bioinformatics， 34（12）： 2012-2018.
GAO XX， 2017. Screening and identification of long noncoding RNAs in thepubertal female goats [D]. Hefei： Anhui Agricultural University： 1-39. [高霄霄， 2017. 初情期雌性山羊 lncRNA 的篩选与鉴定 [D]. 合肥： 安徽农业大学： 1-39.]
HACKL T， HEDRICH R， SCHULTZ J， et al.， 2014. Proovread： Large-scale high-accuracy PacBio correction through iterative short read consensus [J]. Bioinformatics， 30（21）：3004-3011.
JOHNSON R， 2012. Long non-coding RNAs in Huntington’s disease neurodegeneration [J]. Neurobiol Disease， 46（2）：245-254.
KOREN S， WALENZ BP， BERLIN K， et al.， 2017. Canu： Scalable and accurate long-read assembly via adaptive k-mer weighting and repeat separation [J]. Genome Res， 27（5）：722-736.
KOREN S， SCHATZ MC， WALENZ BP， et al.， 2012. Hybrid error correction and de novo assembly of single-molecule sequencing reads [J]. Nat Biotechnol， 30（7）： 693-700.
KIM D， PERTEA G， TRAPNELL C， et al.， 2013. TopHat2： Accurate alignment of transcriptomes in the presence of insertions， deletions and gene fusions [J]. Genome Biol， 14（4）： R36.
KONG L， ZHANG Y， YE ZQ， et al.， 2007. CPC： Assess the protein-coding potential of transcripts using sequence features and support vector machine [J]. Nucl Acid Res， 35（Web Server issue）： W345-W349.
KOREN S， SCHATZ MC， WALENZ BP， et al.， 2012. Hybrid error correction and de novo assembly of single-molecule sequencing reads [J]. Nat Biotechnol， 30（7）： 693-700 DOI 10.1038/nbt.2280.
LIN Z， 2019. Identification of A. comosus var. bracteatus lncRNAs and functional verification of lncABCG11 [D]. Ya’an： Sichuan Agricultural University： 1-113 [蔺珍， 2019. 红苞凤梨lncRNAs的鉴定及lncABCG11的功能验证 [D]. 雅安： 四川农业大学： 1-113.]
LI X， KANAKALA S ， HE YH， et al.， 2017. Physiological characterization and comparative transcriptome analysis of white and green leaves of Ananas comosus var. bracteatus [J]. PLoS ONE， 12（1）： e0169838.
LANGMEAD B， TRAPNELL C， POP M， et al.， 2009. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome [J]. Genome Biol， 10（3）： R25. LI JW， MA W， ZENG P， et al.， 2015. LncTar： A tool for predicting the RNA targets of long noncoding RNAs [J]. Brief Bioinform， 16（5）： 806-812.
LUO YH， DING N， SHI X， et al.， 2017. Generation and comparative analysis of full-length transcriptomes in sweetpotato and its putative wild ancestor I. trifida [J]. BioRxiv， https：//doi.org/10.1101/112425
MA DN， ZHANG XT， WEI LF， et al.， 2018. Benchmarking hybrid correction and assembly using short Illumina reads and long pac bio reads [J]. Genom Appl Biol， 37（4）：1547-1555. [马东娜， 张兴坦， 魏柳锋， 等， 2018. 基因组二代测序数据与三代测序数据的混合校正和组装 [J]. 基因组学与应用生物学， 37（4）： 1547-1555.]
MA J， XIANG YX， XIONG YY， et al.， 2018. SMRT sequencing analysis reveals the full-length transcripts and alternativesplicing patterns in Ananas comosus var.bracteatus [J]. Peer J， 7： e7062
MA JC， 2018. Genome sequence of a widely cultivated poplar and its lnc RNAs response to salt stress [D].Lanzhou： Lanzhou University： 1-82. [马建超， 2018. 新疆杨基因组及其lncRNA响应盐胁迫的研究 [D]. 兰州： 兰州大学： 1-82.]
QIN T， ZHAO HY， CUI P， et al.， 2017. A nucleus-localized long non-coding RNA enhances drought and salt stress tolerance [J]. Plant Physiol， 175（3）： 1321-1336.
SMITH CC， WANG Q， CHIN CS， et al.， 2012. Validation of ITD mutations in FLT3 as a therapeutic target in human acute myeloid leukaemia [J]. Nature， 485（7397）： 260-263.
SUN L， LUO HT， BU DC， et al.， 2013. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts [J]. Nucl Acid Res， 41（17）： e166.
SHARON D， TILGNER H， GRUBERT F， et al.， 2013. A single-molecule long-readsurvey of the human transcriptome [J]. Nat Biotechnol， 31：1009-1014. DOI： 10.1038/nbt.2705.
TRAPNELL C， WILLIAMS BA， PERTEA G， et al.， 2010. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation [J]. Nat Biotechnol， 28（5）： 511-515.
TIAN JJ， 2016. The application of CRISPR/Cas9 system in the zebrafish gene editing [D]. Yangzhou： Yangzhou University： 1-68. [田凈净， 2016CRISPR/Cas9 系统在斑马鱼基因编辑研究中的应用 [D]. 扬州： 扬州大学： 1-68.]
WANG YF，SU WY， ZHANG L， et al.， 2018. Advances of long non-coding RNA in plants [J]. Acta Bot Boreal-Occident Sin， （3）： 582-588. [王艳芳， 苏婉玉， 张琳， 等， 2018. 植物中lncRNAs的研究进展 [J]. 西北植物学报， （3）： 582-588.]
WANG LG， PARK HJ， DASARI S， et al.， 2013. CPAT： Coding-potential assessment tool using an alignment-free logistic regression model [J]. Nucl Acid Res， 41（6）： e74. WANG B ， TSENG E ， REGULSKI M ， et al.， 2016. Unveiling the complexity of the maize transcriptome by single-molecule long-read sequencing [J]. Nat Comm， 7： 11708.
WANG CY， LIU SR， ZHANG XY， et al.， 2017. Genome-wide screening and characterization of long non-coding RNAs involved in flowering development of trifoliate orange （Poncirus trifoliata L. Raf.） [J]. Sci Rep， 7： 43226.
XIONG YY， MA J， HE YH ， et al.， 2018. High-throughput sequencing analysis revealed the regulation patterns of small RNAs on the development of Ananas comosus var. bracteatus leaves [J]. Sci Rep， 8（1）：1947.
XIONG YY， 2019. MicroRNAs identification and screening and functional verification of key microRNAs involved in the albino of Ananas comosus var. bracteatus [D]. Ya’an： Sichuan Agricultural University： 1-95. [熊颖媛， 2019. 红苞凤梨microRNA的鉴定及白化关键microRNA的筛选与功能验证 [D]. 雅安： 四川农业大学： 1-95.]
XUE YB， MA J， HE YH， et al.， 2019. Comparative transcriptomic and proteomic analyses of the green and white parts of chimeric leaves in Ananas comosus var. bracteatus [J]. PeerJ， 7： e7261.
XU WN， HUANG RM， LIU YY， et al.， 2018. Genome sequencing and assembly strategy analyses of Flammulina filiformis [J]. Mycosystema， 37（12）： 1578-1585. [徐偉南， 黄蓉梅， 刘媛媛， 等， 2018. 金针菇基因组测序与组装策略分析 [J]. 菌物学报， 37（12）： 1578-1585.]
YU CL， LUO L， LIAO Q， 2015. Annotation and functional prediction of lncRNAs [J]. Chin J Biochem Mol Biol， （3）： 239-243. [余铖亮， 骆亮， 廖奇， 2015. lncRNAs功能注释和预测 [J]. 中国生物化学与分子生物学报， （3）： 239-243.]
ZHANG YW， TAO Y， LIAO Q， 2018. Long noncoding RNA： A crosslink in biological regulatory network [J]. Brief Bioinform， 19（5）： 930-945.
（责任编辑 李 莉）