槟榔（Areca catechu L.）是棕榈科槟榔属常绿乔木,已成为海南第二大热带经济作物,是海南省东部、中部和南部山区等地区农民脱贫致富的主要经济来源之一。近年来,受黄化等问题严重影响,槟榔产业发展严重受阻。黄化原因复杂,植原体、干旱、养分缺乏和除草剂药害等多种因素均会引起叶片黄化[1,2,3,4,5,6,7]。通过田间调查发现,叶片黄化症状主要有从叶缘1/4处开始发生的特异性黄化和全叶型黄化两种类型。转录组测序技术能全面揭示植物在特定组织、特定生长时期的全部基因的表达情况,且基因组信息与表型特征有很强的关联性,因此被广泛应用在模式和非模式物种转录组的研究中,以揭示植物在特定时刻、特定组织中转录水平上的调控机制[8,9,10,11]。本研究应用转录组测序手段,以表现不同类型黄化症状的叶片为实验材料,筛选出槟榔叶片发生特异性黄化症状的相关差异表达基因,并对其功能和代谢通路进行注释,分析了功能基因的表达模式,探讨了槟榔叶片表现出特异性黄化症状的分子机理,为槟榔黄化病的综合防控提供理论依据。
1 材料与方法
1.1 材料
从叶缘1/4处开始发生特异性黄化的叶片样本1、样本2和样本3分别采自屯昌、琼中、琼海,全叶型黄化的叶片样本4采自定安。健康槟榔叶片为对照组,采自中国热带农业科学院环境与植物保护研究所儋州基地。每组样品取3个生物学重复,分别取自3株槟榔树。
1.2 方法
1.2.1 总RNA提取 总RNA提取方法参照QIAGEN RNeasy Plant Mini试剂盒使用说明,分别对15个样本的总RNA进行提取后,用Nanodrop 2000测定总RNA的浓度和纯度;结合1%的琼脂糖凝胶电泳检测总RNA的质量,将所提取的RNA置于-80 ℃冰箱保存备用。
1.2.2 文库构建和测序 利用带有Oligo（dT）的磁珠与ployA进行A-T碱基配对,从总RNA中分离出mRNA,再利用Illumina HiSeq平台对短序列片段进行测序,其目的是将mRNA进行富集,随后随机将其片段化处理。以片段化mRNA为模板反转录合成cDNA,通过粘性末端补成平末端,3′末端加“A”,连接成“Y”字形的接头和PCR扩增、纯化等步骤,完成测序样本的文库构建及上机测序。委托上海美吉生物医药科技有限公司进行转录组测序。
1.3 转录组测序数据分析
1.3.1 测序数据质控 对raw data进行过滤,除去reads中的接头序列、序列末端（3′端）质量低于Q20的碱基、含未知碱基N比率超过10%的reads及长度小于30 bp的reads后,得到高质量的质控数据（clean data）,以保证后续结果分析的准确性。
1.3.3 基因表达量分析及差异表达基因（DEGs）筛选 利用软件RSEM分别对Unigenes的表达水平进行定量分析,TPM作为表达定量指标[14],TPM具体的运算公式如下：
$\text{TP}{{\text{M}}_{i}}=\frac{{{X}_{i}}}{{{l}_{i}}}\times \left( \frac{1}{\sum\nolimits_{j}{\frac{{{X}_{j}}}{{{l}_{j}}}}} \right)\times {{10}^{6}}$
式中,Xi表示比对到转录本的reads,li表示转录本长度,$\sum\nolimits_{j}{\frac{{{X}_{j}}}{{{l}_{j}}}}$表示所有能够比对到转录本上的reads按长度标准化后数值的总和。
本研究设有3个生物学重复,因此使用DESeq2软件对raw counts进行统计分析,该软件的运算基于负二项分布,设置DEGs筛选参数为：p-adjust<0.05且|log2FC|≥2[15],筛选出与槟榔叶片黄化相关的差异表达基因。将筛选出来的差异表达基因与六大数据库进行比对,阈值设为：E-value≤1×e-5,再进行功能注释和功能富集分析。
1.4 差异表达基因的qRT-PCR验证
以步骤1.2.1所提取的槟榔叶片RNA为模板,用Thermo Scientific RevertAid First Strand cDNA Synthesis Kit试剂盒反转录合成cDNA。根据转录组测序的结果,筛选出表达量在所有样本中表达均较稳定的基因,利用本地Blastx软件对这部分基因进行比对,筛选出了槟榔Actin基因作为内参基因,选取4个有代表性的差异表达基因做qRT-PCR验证,运用Primer 5.0设计荧光定量所需要的引物,引物序列见表1。用TB Green™ Premix Ex Taq™ II（Tli RNaseH Plus）配制20 µL的qRT-PCR反应体系,反应程序如下：95 ℃ 30 s;95 ℃ 3 s,60 ℃ 30 s,共40个循环;95 ℃ 15 s,60 ℃ 60 s,95 ℃ 15 s。
表1 荧光定量引物设计Tab. 1 qRT-PCR primer design |
| 引物名称 Primer name | 碱基序列（5′-3′） Primers sequences (5′-3′) |
|---|---|
| 68503F | GTGATGGGCGTGTCAAGA |
| 68503R | TGCCGAAGAGATACTCGTTC |
| 117057F | TCAAGACTGATTGGAGTAAGGC |
| 117057R | AGTGGTTCCTCACCCATCTC |
| 115124F | GCAGTTGGAATGGGAGAATC |
| 115124R | TTGTGCCTGGAAGGTTGATA |
| 116247F | GAGTATGCGAGGAGGAAAGC |
| 116247R | CCCACCCAAAGTCAGTTCTAT |
| Actin F | CGTATGAGCAAAGAGATTACTGC |
| Actin R | GAAGCCAAGATGGACCCTC |
2 结果与分析
2.1 转录组测序数据分析
15个槟榔叶片样本总共得到741 481 770个raw reads,过滤后得到730 688 448个clean reads,clean reads占raw reads的比值均在98.00%以上。从整体上看,碱基错误率均在0.03%以内,Q20为97.01%~97.54%,Q30为91.57%~92.70%,GC平均含量约为48.50%。表明整体测序质量好,数据可以用于后续生物信息学分析。组装过滤后,得到的clean reads共包含321 402个Unigene,所得Unigene的GC含量为42.09%,N50为1647 bp,E90N50为2350 bp,平均长度为857.35 bp。能满足后续做基因表达定量及差异基因表达分析基本要求。
2.2 功能注释
将组装获得的Unigenes序列与六大数据库（NR,Swiss-Prot,Pfam,COG,GO和KEGG数据库）进行注释（E-value≤1×e-5）,结果显示：至少能比对注释到一个数据库上的Unigenes序列有126 444条;在六大数据库中均能比对注释上的Unigenes序列有5693条,其比例仅有1.77%;60.66%的Unigenes在任何数据库中均无被比对注释上。
2.3 差异基因的筛选及分析
将3组实验组分别与对照组两两间进行基因的差异分析,并将DEGs筛选参数设为：p-adjust<0.05且|log2FC|≥2。结果显示,样本1特异性黄化叶片共10 619个基因显著差异表达,其中6422个上调表达,4197个下调表达;样本2特异性黄化叶片共18 433个基因显著差异表达,其中8310个上调表达,10 123个下调表达;样本3特异性黄化叶片共3788个基因显著差异表达,其中2035个上调表达,1753个下调表达;样本4全黄叶型黄化叶片共5583个基因显著差异表达,其中3085个上调表达,2498个下调表达（图1）。
根据上述筛选的不同差异分组的基因集合绘制韦恩图,发现表现为特异性黄化的样本与全叶型黄化的样本,共有654个共同显著差异表达的基因,其中314个基因在所有样本中均上调表达,309个均下调表达,还有31个基因上下调情况不一致;其中有440个基因仅在特异性黄化样本中为共同显著差异表达基因,221个基因发生共同上调表达,205个共同下调表达,14个基因在3组黄化组中上下调情况不一致（图2）。
2.4 特异性叶片黄化样本差异表达基因GO功能注释和富集分析
从整体上分析基因的表达情况,揭示槟榔叶片发生特异性黄化症状后差异表达基因参与的功能类别的DEGs富集在分子功能、生物过程、细胞组分三类GO（Level 2）条目中的情况,在生物过程得到最多的分类,DEGs主要注释在metabolic process（代谢过程）、cellular process（细胞过程）、binding（连接）、catalytic activity（催化反应）中（图3）。
根据差异基因上调或下调情况使用超几何检验分别进行GO富集分析,其中特异性黄化共同上调表达的DEGs富集最显著的条目主要是：DNA metabolic process（DNA代谢过程）、transferase activity, transferring glycosyl groups（转移酶活性,转移糖基）、transferase activity,transferring hexosyl groups（转移酶活性,转移己糖基）、amino acid transport（氨基酸转移）、GDP-mannose 3,5-epimerase activityDNA（gmp-甘露糖3,5-表位酶活性）等,具体富集情况见图4。下调表达的DEGs富集最显著的条目主要是：ADP binding（ADP结合）、post-embryonic root morphogenesis（胚后期的根系形态）、lateral root morphogenesis（侧根形态发生）、root morphogenesis（根系形态）、gibberellic acid homeostasis（赤霉酸体内平衡）等,具体富集情况见图5。
图4 特异性黄化共同上调表达的差异基因GO富集Fig. 4 GO enrichment of differentially expressed genes up-regulated in leaves with specific yellowing symptoms |
2.5 特异性叶片黄化样本差异表达基因的KEGG显著性富集分析
为了解与槟榔叶片发生特异性黄化症状相关的关键代谢通路,进行了差异基因的KEGG通路富集分析,结果显示,表现为特异性黄化样本的DEGs富集靠前的前5个通路为：betalain biosynthesis（甜菜碱生物合成）、sulfur metabolism（硫代谢）、endocytosis（内吞作用）、ascorbate and aldarate metabolism（抗坏血酸和醛酸代谢）、glutathione metabolism（谷胱甘肽代谢）（图6）。
2.6 特异性叶片黄化样本差异表达基因的代谢通路分析
2.6.1 能量代谢通路 在氧化磷酸化代谢和光合生物固碳途径中,叶片表现特异性黄化的DEGs编码琥珀酸脱氢酶上的琥珀酸脱氢酶（泛醌）细胞色素b560亚基（SDHC）的基因下调表达。在光合生物固碳通路中,编码磷酸核酮糖激酶（PRK）的基因上调表达,编码二磷酸核酮糖羧化酶（rbcL）的基因下调表达。
2.6.2 激素代谢通路 油菜素内酯信号转导和生物合成途径仅在表现特异性叶片黄化的DEGs被特异注释,在油菜素内酯信号转导途径中编码木葡聚糖内糖基转移酶的基因上调表达;在油菜素内酯生物合成途径中编码合成油菜素内酯的关键酶油菜素类固醇-6-氧化酶1的基因下调表达。
2.6.3 次级代谢通路 槟榔叶片表现出特异性黄化症状的样本特异注释到的10个次级代谢通路,分别是柠檬酸循环、抗坏血酸和醛酸代谢、醚脂类代谢、乙醛酸和二羧酸代谢、核黄素代谢、吲哚生物碱的生物合成、油菜素内酯生物合成、类黄酮生物合成、甜菜碱生物合成以及二苯乙烯、二芳基庚醇和姜辣素生物合成。其中,在抗坏血酸和醛酸代谢中,编码GDP-甘露糖3,5-差向异构酶（催化GDP-甘露糖转化为左旋GDP-L-半乳糖和GDP-L-古洛糖,以及催化GDP-L-古洛糖转化为GDP-L-半乳糖）的基因上调表达,编码L-抗坏血酸过氧化物酶的基因下调表达;在醚脂类代
谢途径中,编码磷脂酶D1/2的基因上调表达;乙醛酸和二羧酸代谢途径中,编码二磷酸核酮糖羧化酶小链（RuBisCO小亚基）的基因下调表达;在核黄素代谢途径中,编码3,4-二羟基-2-丁酮4-磷酸合成酶/GTP环化水解酶II的基因上调表达;在吲哚生物碱的生物合成途径中,编码聚尿苷-醛酯酶的基因上调表达;在甜菜碱生物合成途径中,编码4,5-多巴雌二醇双加氧酶的基因上调表达;在类黄酮生物合成、苯丙烷类生物合成和二苯乙烯、二芳基庚醇和姜辣素生物合成途径3条代谢途径中,编码莽草酸-羟基肉桂酰转移酶的基因上调表达。
2.7 qRT-PCR验证
3 讨论
油菜素内酯是一种甾醇类植物激素,对植物生长发育和抗逆等有重要的调控作用[20,21]。特异性叶片黄化的DEGs在油菜素内酯信号转导途径中编码木葡聚糖内糖基转移酶且上调表达,该酶是一种细胞壁修饰酶,该酶在植物形态建成中能够转移细胞壁中的木葡聚糖,从而影响细胞壁的形成和降解,在植物响应逆境胁迫方面发挥着重要作用[22]。在油菜素内酯合成途径中编码合成油菜素内酯的关键酶油菜素类固醇-6-氧化酶1的基因下调表达,油菜素内酯含量的降低可间接调控植物光合作用、色素合成等通路,导致叶片表现黄化症状。编码油菜素内酯信号转导和合成途径中的基因仅在特异性黄化的样本中为显著差异基因,说明上述基因可能对槟榔叶片发生特异性黄化有重要的调控作用。
本研究通过对转录组测序结果中差异基因的分析,特异性黄化样本中的差异基因主要富集在能量代谢通路、激素代谢通路和和次级代谢途径中,其中在油菜素内酯信号转导途径、类黄酮生物合成和苯丙烷类生物合成途径中显著上调表达,而在油菜素内酯合成途径中显著下调表达,可能对槟榔叶片特异性黄化症状起重要的调控作用,可能调控了植物体内的内源激素的表达,影响了细胞内的渗透等共同调控槟榔体内的代谢通路,导致槟榔叶片表现出黄化症状。
叶片表现出特异性黄化症状的样本特异注释到的10个次级代谢通路,其中编码甜菜碱生物合成、抗坏血酸和醛酸代谢、吲哚生物碱的生物合成、二苯乙烯、二芳基庚醇和姜辣素的生物合成、醚脂类代谢、核黄素代谢以及类黄酮生物合成通路的相关基因显著上调表达;编码油菜素内酯生物合成、柠檬酸循环、乙醛酸和二羧酸代谢以及抗坏血酸和醛酸代谢通路的相关基因显著下调表达,可能调控了植物体内的内源激素的表达,影响了细胞内的渗透等共同调控槟榔体内的代谢通路,导致槟榔叶片表现出相似的黄化症状。而编码类黄酮生物合成、苯丙烷类生物合成以及二苯乙烯、二芳基庚醇和姜辣素生物合成这3条次生代谢途径中的莽草酸-羟基肉桂酰转移酶的基因,且显著上调表达。黄酮类植物着色色素的生物合成途径是通过苯丙烷类代谢来完成的,对病原菌侵染的抗性和着色色素形成有密切的关系[23]。黄酮类化合物是高等植物中的一类次生代谢产物,可以通过提高过氧化物酶的活力、减少丙二醛含量以及增加脱落酸的含量等对提高植株的抗旱性起着极其关键的作用[24,25,26,27]。在发生特异性黄化的槟榔叶片样本中,莽草酸-羟基肉桂酰转移酶显著上调表达,说明调控该酶的基因上调表达可能影响了植物着色色素的合成导致叶片发生黄化症状;类黄酮合成途径中莽草酸-羟基肉桂酰转移酶的表达反而受到促进,其作用机理需进行后续验证性实验才可以确认。