茶树(Camellia sinensis) 起源于中国西南地区，是我国重要的多年生经济作物之一。其中，乌龙茶是我国特有的一种高香型半发酵茶类，主产于福建省、广东省、中国台湾等地区，因其具有独特优雅的花果香和鲜醇浓厚口感而享誉中外^[1]。乌龙茶独特风味品质的形成与其加工工艺流程关系密切。萎凋是乌龙茶采后加工过程中的首道工序，对乌龙茶独特风味的形成具有重要作用。在乌龙茶加工过程中，需要根据茶叶的理化性状，灵活适度进行日光萎凋处理。在日光萎凋处理中，鲜叶会遭受干旱胁迫、高温胁迫和紫外线胁迫等多重环境胁迫的共同影响，但由于萎凋叶仍保持相对完整的形态结构，光合作用、蒸腾作用和植物激素信号转导等生理活动仍在叶片中正常进行^[2-3]。此外，萎凋过程中还会发生多种生物化学反应和次生代谢产物丰度的变化，来维持叶片细胞膜的通透性和代谢平衡以应对环境胁迫带来的不利影响。外界胁迫可诱导植物产生许多次生代谢产物，这些次生代谢产物具有帮助植物抵御外界胁迫的作用^[4]。此外，茶叶中的次生代谢产物还赋予了茶叶独特的风味品质^[5]。类黄酮是茶叶中最重要的次生代谢产物，与乌龙茶滋味的形成密切相关，是茶叶苦涩味形成的主要原因^[6]。挥发性化合物也是评价乌龙茶风味品质的重要因子。有报道显示，乌龙茶中花果香品质的形成与挥发性萜类化合物密切相关^[7]。同时，逆境胁迫可以触发信号分子，特别是植物激素的生物合成^[8]。这些植物激素作为信号转导的主要信使，通过调控代谢物生物合成相关基因的表达水平，影响相关代谢物含量的变化^[9]。过去的大部分研究都集中在探究茶树采前代谢物含量变化与代谢途径中关键基因表达水平之间的相关性^[10-12]，茶树采后生产相关研究仅仅局限于生产工序处理前后的次生代谢物含量测定^[13-16]。但尚未见到从转录调控角度揭示乌龙茶采后风味品质形成的报道。

高通量测序技术以较低的成本实现了较高的测序覆盖率，可以快速全面地挖掘样本中所有mRNA转录本的丰度信息。目前，转录组测序已广泛应用于茶叶风味品质形成相关基因信息的挖掘^{[11-12, 17]}。Wang等^[18]利用转录组测序，鉴定出参与儿茶素合成的12个结构基因。而萜类生物合成相关结构基因的鉴定和代谢通路绘制工作也已全部完成^[19]。随着染色体级别茶树基因组数据的公布^[20-23]，为更深入研究编码基因的上游调控机制奠定了基础。可变剪接是指一个mRNA前体通过不同的剪接方式产生多种mRNA转录本的过程，是植物调控其编码基因表达水平的重要分子机制^[24-25]。研究表明，可变剪接是植物应对逆境胁迫的重要方式，通过产生大量剪接转录本，增强转录本和编码蛋白的多样性来应对和缓解逆境胁迫带来的不利影响^[26-27]。尽管植物次生代谢产物的合成途径中存在大量可变剪接事件^[28-29]，然而可变剪接在茶叶加工过程中对风味品质形成相关基因的具体调控机制及不同剪接转录本在茶叶风味品质形成的生物学功能尚不清楚。

为揭示mRNA在乌龙茶萎凋过程中影响风味品质形成的分子机制，本研究通过对鲜叶、室内萎凋叶和日光萎凋叶进行转录组测序，挖掘差异表达基因主要富集的代谢途径，鉴定可变剪接事件和差异剪接基因，对其中参与风味品质形成及其上游调控途径中的关键差异表达基因和差异剪接基因的表达模式进行深入分析。同时，还分析了差异剪接转录本的表达丰度和相关代谢物及植物激素的含量之间的相关性，以期进一步了解转录调控及可变剪接机制在乌龙茶风味品质形成中的具体作用以及日光萎凋对乌龙茶风味品质形成的重要性奠定理论基础。

1 材料与方法 1.1 试验材料

供试材料为种植于福建农林大学教学实践茶园(E 119° 14′，N 26° 05′) 的铁观音茶树。采集长势一致的铁观音茶树叶片，采摘标准为一芽三叶。将采集后的茶树叶片等分为3部分，未进行任何处理的鲜叶作为对照组；第二部分的茶叶进行日光萎凋处理(光照强度：(40 000±1 000) lx；温度：(25±2) ℃；相对湿度：(60%±5%)，时间：45 min)；第三部分的茶叶进行室内萎凋处理(光照强度：(50±5) lx，其他参数和日光萎凋处理一致)。每个处理进行3次生物学重复，采集后的鲜叶(FL)、室内萎凋叶(IW) 和日光萎凋叶(SW) 经液氮急速冷冻后存放在-80 ℃冰箱储存备用。

1.2 RNA提取与转录组测序

使用TransZol UP试剂盒(全式金生物，北京，中国) 分别提取鲜叶、室内萎凋叶和日光萎凋叶的总RNA。提取后RNA的质量和完整性分别利用Agilent Bioanalyzer 2100 (安捷伦科技，美国) 和NanoDrop 2 000 (赛默飞世尔科技，美国) 进行检测。利用Oligo(dT)介导的mRNA富集法处理检测合格的RNA，随后加入破碎缓冲液将富集的mRNA进行随机打断。以打断后的RNA片段为模板，使用六碱基随机引物进行反转录，再合成cDNA第二链。将双链cDNA分别进行末端修复、5′末端磷酸化、连接上鼓泡接头后进行PCR扩增。PCR产物升温变性解链后，加入桥式引物将单链cDNA环化后获得cDNA文库。最后使用BGISEQ-500平台(华大基因，中国) 对cDNA文库进行高通量测序。对测序后的原始数据进行质量控制，使用SOAPnuke^[30]过滤去除接头序列、低质量序列以及含有未知碱基N比例超过5%的序列。采用HISAT^[31]将过滤后的数据比对茶树染色体水平参考基因组^[20]。

1.3 转录组测序数据分析

转录本的相对表达水平使用Bowtie2^[32]和RSEM^[33]进行计算，并采用FPKM (fragments per kilobase per million mapped reads) 值对相对表达水平进行归一化处理。

基于DEGseq算法，以|fold change|≥2且Q-value≤0.001作为筛选标准，进行差异表达基因(differentially expressed gene，DEG)的筛选^[34]。随后将差异表达基因分别进行基因本体论(gene ontology，GO) 和京都基因与基因组百科全书(kyoto encyclopedia of genes and genomes, KEGG)富集分析，以Q-value≤ 0.01为阈值，挖掘出差异表达基因发挥主要的生物学功能或参与影响的生物学通路。相关差异表达基因的表达谱使用Tbtools进行可视化^[35]。

使用rMATS检测鲜叶、室内萎凋叶和日光萎凋叶中的可变剪接事件和差异剪接基因^[36]。采用三氯化铝比色法测定茶叶样品中的类黄酮含量，使用紫外分光光度计测定茶样在415 nm处的吸光度，并参照Zhu等^[37]的研究方法进行计算。参照Tai等^[11]的研究方法，使用高效液相色谱(high performance liquid chromatography, HPLC) 系统进行儿茶素含量的测定，检测波长为278 nm，柱温为25 ℃。通过酶联免疫吸附法(enzyme-linked immunosorbent assay, ELISA) 检测样品中的赤霉素(GA) 含量，使用Infinite M200 PRO酶标仪(帝肯，瑞士) 进行吸光度测定，计算方法参照Zhu等^[9]进行，每个实验进行3次生物学重复。用于分析的相关代谢物及植物激素数据来自本课题组已报道的研究^{[9, 37]}。采用Pearson相关系数对差异剪接转录本的转录丰度和相关代谢物及植物激素的含量进行相关性分析。

1.4 差异表达基因和差异剪接基因的转录水平验证

采用试剂盒TransScript First-Strand cDNA Synthesis SuperMix (全式金生物，中国) 将提取的RNA反转录合成cDNA。使用LightCycler 480 (罗氏集团，瑞士) 对关键差异表达基因和差异剪接基因进行实时定量聚合酶链式反应(quantitative realtime polymerase chain reaction, qRT-PCR)扩增。以GAPDH和β-actin为内参基因，采用相对定量法2^–ΔΔCt算法^[38]，以检测关键差异表达基因在不同萎凋处理中的表达水平。每个实验进行3次生物学重复。使用SPSS 25软件进行差异显著性分析，相关引物信息见表 1。

2 结果与分析 2.1 转录组测序质量评估

经过高通量测序后，从所有样品中共获得662.40 Mb的原始数据(表 2)，对原始数据进行过滤后获得平均大小为66.15 Mb的有效数据。各样品中碱基数量分布范围为6.56–6.70 Gb，且测序质量参数Q30均在88.53%以上。将所有样品中的有效数据比对到参考基因组上，发现转录组数据覆盖基因组数据的比值均在80.98%以上，表明高通量测序数据覆盖度高、数据可靠，可用于后续试验分析。比对参考基因组后从9个样品中共鉴定出40 189个基因。随后，对基因的序列长度分布情况进行了分析(图 1)，其中序列长度小于1 000 bp的基因数量为16 101个(40.06%)，序列长度范围在1 000–2 000 bp内的基因数量为15 659个(38.96%)，而序列长度超过2 000 bp的基因数量为8 429个(20.98%)。进一步分析发现，序列长度在3 000 bp以上的基因数量最多，共有2 754个(6.85%)，而长度范围在2 900–3 000 bp内的基因数量最少，共有303个(0.75%)。

2.2 差异表达基因的鉴定

从3个转录组中共鉴定出10 793个差异表达基因，其中鲜叶和室内萎凋叶(FL vs. IW) 中有6 688个差异表达基因(上调基因3 703个，下调基因2 985个)；鲜叶和日光萎凋叶(FL vs. SW) 中有7 882个差异表达基因(上调基因4 415个，下调基因3 467个)；室内萎凋叶和日光萎凋叶(IW vs. SW)中有1 095个差异表达基因(上调基因558个，下调基因537个) (图 2A)。有3 991个差异表达基因为FL vs. IW和FL vs. SW 2组所共有的(图 2B)，表明萎凋处理可能引起这些基因的表达水平发生了显著变化。此外，有214个表达基因只在IW vs. SW中呈现差异表达趋势，推测这些差异表达基因可能是揭示室内萎凋和日光萎凋对茶叶风味品质形成的关键因素。

2.3 差异表达基因的GO富集分析

为初步明确差异表达基因在萎凋处理中发挥的具体生物学功能，对上述已鉴定出的10 793个差异表达基因进行了GO富集分析。结果表明共有8 789个差异表达基因注释到GO数据库。在FL vs. IW中共有5 479个差异表达基因注释到生物过程、细胞组分及分子功能这3个GO一级功能分类(图 3)，其中生物过程涉及19个GO二级功能分类，参与细胞进程、代谢过程及生物调节的差异表达基因最多；细胞组分中注释到14个GO二级功能分类，其中差异表达基因主要富集在细胞膜、细胞膜组分和细胞；而在分子功能中发现10个GO二级功能分类，主要涉及催化活性、结合及转运活性。在FL vs. SW中共有6 445个差异表达基因注释到GO数据库，生物过程中注释到21个GO二级功能分类，主要包括细胞进程、代谢过程和生物调节；细胞组分中鉴定出14个GO二级功能分类，其中差异表达基因主要富集于细胞膜、细胞膜组分和细胞；而从分子功能中发现12个GO二级功能分类，参与催化活性、结合及转运活性的差异表达基因数量最多。此后从IW vs. SW中发现共有445个差异表达基因涉及生物过程、细胞组分及分子功能。在生物过程中共鉴定出17个GO二级功能分类，主要包括细胞进程、代谢过程及生物调节；在细胞组分中注释到13个GO二级功能分类，其中差异表达基因主要富集在细胞膜、细胞及细胞膜组分；而在分子功能中发现10个GO二级功能分类，主要涉及催化活性、结合和转运活性。从上述结果中发现，萎凋处理影响了萎凋叶中代谢、合成和转运的生物学功能，这可能对茶叶风味品质的形成产生影响。

2.4 差异表达基因的KEGG富集分析

为进一步研究不同萎凋处理后差异表达基因参与的代谢通路，对所有差异基因进行KEGG富集分析。结果显示，共有4 605个差异表达基因注释到5个KEGG一级分支和19个KEGG二级分支(图 4)。在FL vs. IW中，除去全局纵览分支后，2 883个差异表达基因主要分布在次生代谢、信号转导、折叠，分类和降解、翻译及环境适应分支中。而对FL vs. SW中的3 372个差异表达基因和IW vs. SW中的436个差异表达基因进行分析后，发现在除去全局纵览分支后，两组比较中的差异表达基因都同样注释到和FL vs. IW中一致的5个KEGG二级分支中，表明萎凋处理不仅可能通过影响次生代谢中关键基因的表达，还可能通过植物激素信号转导、蛋白质翻译、折叠和降解间接影响相关代谢中结构基因生物学功能的实现，进而参与调控茶叶风味品质形成的相关次生代谢通路。

对差异基因进行KEGG富集分析，结果显示(图 5)，FL vs. IW中差异表达基因主要富集在剪接体、植物激素信号转导、类黄酮合成、生物钟和萜类化合物合成通路。在FL vs. SW中，富集排名前五的通路与FL vs. IW中的一致。而IW vs. SW中的差异表达基因主要涉及剪接体、植物激素信号转导、萜类化合物合成、类黄酮合成及卟啉与叶绿素代谢通路。综合分析上述结果发现，有4个富集通路在3组比较中均有出现，其中类黄酮合成和萜类化合物合成通路与茶叶滋味和香气的形成关系密切，而剪接体和植物激素信号转导通路在萎凋过程中可能通过间接影响次生代谢参与茶叶风味品质的形成。因此，选取这4个主要富集的KEGG代谢通路进行深入研究。

2.5 萎凋过程中茶叶风味品质形成相关基因挖掘与表达模式分析

KEGG分析发现萎凋过程中的差异表达基因主要富集于类黄酮合成、萜类化合物合成、植物激素信号转导及剪接体通路。其中，类黄酮合成通路中鉴定出的差异表达基因包括儿茶素合成基因C4H (cinnamate 4-hydroxylase)、CHI (chalcone isomerase)、F3H (flavanone 3-hydroxylase)、ANS (anthocyanidin synthase)和LAR (leucoanthocyanidin reductase)，黄烷醇合成基因FLS (flavonol synthase)以及木质素合成基因HCT (hydroxycinnamoyltransferase)。而萜类化合物合成通路中共鉴定出11个上游差异表达基因和6个下游差异表达基因。其中上游基因可以细分为MEP (methylerythritol phosphate)分支基因和MVA (mevalonate)分支基因。MEP分支基因包括DXS (1-deoxy-D-xylulose-5-phosphate synthase)、CMK (4-diphosphocytidyl-2-C-methyl-D-erythritolkinase)及HDS (4-hydroxy-3-methylbut-2-enyl-diphosphate synthase)，而MVA分支基因包括AACT (acetyl-CoA C-acetyltransferase)、HMGS (hydroxymethylglutaryl-CoA synthase)和HMGR (hydroxymethylglutaryl-CoA reductase)。萜类化合物合成中的下游差异表达基因则包括GGPPS (geranylgeranyl diphosphate synthase) 和TPS (terpene synthase)。植物激素信号转导通路中的差异表达基因包括茉莉酸信号转导基因JAR (jasmonic acid-amino synthetase)、COI1 (coronatine-insensitive 1)、JAZ (jasmonate ZIM domain)和MYC2 (myelocytomatosis2)，生长素信号转导基因GH3 (gretchen hagen3) 和SAUR (small auxin-up RNA)，脱落酸信号转导基因PYL (pyrabactin resistance-like)，油菜素内酯信号转导基因BRI1 (brassinosteroid insensitive 1)，与水杨酸信号转导基因TGA (transcription factor TGA)。而剪接体通路中的差异表达基因涉及Prp5 (pre-mRNA-splicing factor 5)、Prp22 (pre-mRNA-splicing factor 22)、Sm (spliceosomal Sm protein) 和SR (serine/arginine-rich protein)。随后，对上述基因的表达模式进行了分析。结果表明(图 6)，儿茶素合成相关基因的表达水平受到萎凋处理的显著抑制，且日光萎凋处理后相关基因的表达水平最低。黄烷醇和木质素合成相关基因在鲜叶中呈现最低的转录水平，而在日光萎凋处理后出现了大幅度升高。在植物激素信号转导通路中，茉莉酸信号转导基因JAR、COI1和MYC2的表达水平受到萎凋处理的显著诱导，而茉莉酸信号转导基因JAZ的表达则受到萎凋处理的显著抑制。在其他植物激素信号转导通路中，除油菜素内酯信号转导基因BRI1在萎凋叶中呈现低表达水平，脱落酸信号转导基因PYL、水杨酸信号转导基因TGA与生长素信号转导基因GH3和SAUR的表达水平都在萎凋处理后呈现上升趋势。在萜类化合物合成通路中相关基因都在萎凋叶中出现高表达，其中在日光萎凋叶中的转录丰度最高。相似地，剪接体通路中鉴定出的差异表达基因均在日光萎凋叶中呈现最高表达水平，在鲜叶中的表达水平最低。上述结果表明，萎凋处理不仅直接影响相关风味代谢物合成通路中关键基因的表达水平，还可能参与调节可变剪接事件，使相关风味物质代谢中结构基因形成不同的转录本，间接影响相关基因的表达水平，进而影响茶叶风味的形成。

2.6 萎凋过程中可变剪接事件分析

上述KEGG分析发现剪接体通路为主要的富集通路之一，且剪接体通路中核心基因Prp5、Prp22、Sm和SR的表达水平在不同萎凋处理后存在显著差异。有研究表明，可变剪接是转录后调控基因表达的重要手段，通过产生特异的可变剪接体，可在植物应对逆境胁迫中发挥重要作用^[39]。因此，对茶叶萎凋过程中的可变剪接事件进行了深入分析。从所有样本的11 696个基因中共鉴定出31 872个可变剪接事件(表 3)。在室内萎凋叶和日光萎凋叶中鉴定出的可变剪接事件数分别为10 334和11 446，均多于鲜叶中发现的可变剪接事件数(10 092)，表明萎凋处理促进了可变剪接事件数量的增加。进一步分析可变剪接类型，发现在3个样本中最主要的可变剪接事件类型均为内含子滞留(retained intron, RI)，这与在玉米^[40]和棉花^[41]中鉴定的主要可变剪接事件类型一致。相较于鲜叶，RI类型可变剪接事件在室内萎凋叶和日光萎凋叶中的数量出现了显著增加，表明萎凋处理特别是日光萎凋处理显著提高了RI类型可变剪接事件的发生频次。

基于可变剪接基因的相对转录丰度(FPKM)，对FL vs. IW、FL vs. SW和IW vs. SW三组比较中的差异剪接基因进行了分析鉴定(表 4)。在FL vs. IW和FL vs. SW中鉴定出的差异剪接基因数分别为303和590，均多于IW vs. SW中发现的差异剪接基因数(248)，表明萎凋处理显著诱导了大量差异剪接基因的产生。对所有差异剪接基因的可变剪接类型分析发现，差异剪接基因中发生RI型可变剪接的基因数量最多，这与鲜叶、室内萎凋叶和日光萎凋叶中的最主要可变剪接类型鉴定结果一致，表明日光萎凋不仅会影响可变剪接数量的变化，而且影响可变剪接转录本的表达水平，产生大量的差异剪接基因。为进一步明确差异剪接基因在茶叶萎凋过程中行使的生物学功能，对FL vs. IW、FL vs. SW和IW vs. SW三组比较中差异剪接基因进行KEGG富集分析，结果显示(图 7)，FL vs. IW中富集排名前五的通路分别为类黄酮合成、植物激素信号转导、果糖和甘露糖代谢、氨基酸合成及苯丙氨酸，酪氨酸和色氨酸代谢通路。而类黄酮合成和植物激素信号转导通路在FL vs. SW和IW vs. SW这两组比较中也均属于富集排名前五的通路，而萎凋过程中差异表达基因也主要活跃于上述两条通路。因此，在茶叶萎凋过程中，类黄酮合成和植物激素信号转导通路中相关结构基因的表达变化可能受到可变剪接调控机制的影响。

为更好地理解可变剪接在茶叶萎凋过程中的调控机制，对类黄酮合成和植物激素信号转导通路中鉴定出的差异剪接基因进行了进一步分析(图 8A)。在类黄酮合成通路中，2个差异剪接基因FLS各产生了5个可变剪接体。其中，2个全长转录本CsFLS1和CsFLS2分别编码313个和422个氨基酸，同时含有1个完整的黄烷醇合成结构域。而植物激素信号转导基因CsPIF (phytochrome interacting factor)和CsGID1 (gibberellin-insensitive dwarf 1) 也存在多个差异剪接体。4个CsFLS可变剪接转录本(CsFLS1-a、CsFLS1-b、CsFLS1-c和CsFLS2-d)和1个CsGID1可变剪接转录本(CsGID1-a) 结构中都发现存在一个提前终止密码子(PTC)，它在基因结构层面的显著变化可能会影响编码蛋白正常生物学功能的实现。此外，其他非PTC型剪接转录本(CsFLS1-d、CsFLS2-a、CsFLS2-b、CsFLS2-c、CsPIF-a和CsPIF-b)存在内含子滞留或5′/3′段可变剪接，可能会导致获得或失去部分保守结构域，也可能会引起可变剪接转录本产生功能分化。

随后，基于转录组检测的FPKM值，分析了上述可变剪接转录本在茶叶萎凋过程中的表达模式(图 8B)。类黄酮合成通路中的全长转录本CsFLS1在鲜叶中呈现最低的转录水平，而萎凋处理后CsFLS1的转录水平则出现了大幅度升高。CsFLS1对应的4个可变剪接转录本(CsFLS1-a、CsFLS1-b、CsFLS1-c和CsFLS1-d)也遵循类似的表达模式。而CsFLS2的4个可变剪接转录本与全长转录本却呈现截然相反的表达水平，表明不同转录本之间在茶叶萎凋过程中可能存在潜在的功能差异。植物激素信号转导基因CsGID1和CsPIF及其对应的可变剪接转录本(CsGID1-a和CsPIF-a) 的表达则受到萎凋处理的显著抑制，在室内萎凋叶中的表达水平最低。此外，PTC型剪接转录本在茶叶萎凋过程中的表达模式与其对应的全长转录本基本一致，但在各个样品中均呈现较低的表达水平。此外，CsPIF-b转录本，在鲜叶、室内萎凋叶和日光萎凋叶中均不表达，可能与内含子滞留引起的CsPIF-b基因结构发生显著改变有关。

为进一步探究差异剪接转录本在类黄酮合成和植物激素信号转导通路中的具体作用，对差异剪接转录本的转录丰度和相关代谢物(总类黄酮和儿茶素) 及植物激素(GA) 含量之间的相关性进行了关联分析(表 5)。全长转录本CsFLS1和其对应的4个可变剪接转录本(CsFLS1-a、CsFLS1-b、CsFLS1-c和CsFLS1-d) 与类黄酮和儿茶素的含量呈负相关关系。而另一个全长转录本CsFLS2与类黄酮和儿茶素的含量存在正相关，而其对应的4个可变剪接转录本(CsFLS2-a、CsFLS2-b、CsFLS2-c和CsFLS2-d) 与类黄酮和儿茶素的含量呈负相关关系。在植物激素信号转导通路中，CsGID1和CsPIF及其对应的可变剪接转录本CsGID1-a和CsPIF-a的转录水平与GA含量的变化趋势一致。此外，我们还发现PTC型剪接转录本与代谢物或植物激素含量的相关性要弱于其对应的全长转录本，推测可能与PTC型剪接转录本的基因结构不完整有关。

2.7 萎凋过程中差异表达基因和差异剪接基因的qRT-PCR分析

选择KEGG主要富集通路中的12个差异表达基因及类黄酮合成和植物激素信号转导通路中鉴定出的相关差异剪接基因进行qRT-PCR分析。结果显示(图 9)，相关差异表达基因和差异剪接基因在茶叶萎凋过程中的表达模式与转录组数据集中的结果一致，这证实了高通量测序数据及后续分析的可靠性与合理性。

3 讨论 3.1 日光萎凋叶中类黄酮合成基因的转录抑制与类黄酮含量的减少密切相关

萎凋是促进乌龙茶风味品质形成的第一道工序，大量风味品质相关代谢物的含量在萎凋前后产生了显著变化^{[3, 14]}。而在转录水平层面，相关代谢途径中结构基因在乌龙茶萎凋过程中影响风味品质形成的分子机制尚不清楚。本研究通过对鲜叶、室内萎凋叶和日光萎凋叶进行转录组测序，发现类黄酮合成通路是差异表达基因主要富集的通路之一。其中，C4H负责催化香豆酸及其衍生物，为后续合成类黄酮代谢物提供前体物质^[42-43]。有研究表明，C4H的转录水平与类黄酮和儿茶素的含量积累呈正相关^[44]。因此，日光萎凋叶中低表达的C4H，可能会减少用于后续类黄酮合成的底物含量，进而影响下游类黄酮和儿茶素的积累水平。有报道显示，葡萄中CHI的上调促进了类黄酮的积累^[45]。而萎凋处理显著抑制了CHI的正常表达，可能使萎凋叶中的类黄酮和儿茶素的合成受阻，这与Wang等^[3]在萎凋处理后发现CHI的转录水平和儿茶素含量呈现下降的现象一致。此外，类黄酮合成下游基因F3H、ANS和LAR的转录丰度也被证实与类黄酮和儿茶素的含量呈正相关^{[12, 46]}。相较于鲜叶和室内萎凋叶，日光萎凋叶中更低表达的F3H、ANS和LAR可能进一步限制了类黄酮合成代谢流进入儿茶素分支的合成方向。

在日光萎凋叶中表达大幅升高的FLS，主要承担将二氢黄酮醇转化为黄酮醇的功能，与ANS和LAR等儿茶素合成基因共同竞争同一前体物质二氢黄酮醇。有研究表明FLS的表达与儿茶素的积累呈负相关^[47]。因此，室内萎凋叶中低表达的FLS可能对儿茶素的生物合成影响较小，而在日光萎凋叶中高表达的FLS则会显著抑制儿茶素的积累。而另一个在萎凋处理后出现上调表达的HCT，隶属于木质素合成途径。与类黄酮合成途径相似，均位于苯丙烷代谢途径的下游。在拟南芥中，木质素合成和类黄酮合成共同竞争苯丙烷代谢产生的同一前体物质^[48]。Vanholme等^[49]发现HCT的表达水平直接影响木质素在植物体内的积累。而本研究中发现，日光萎凋处理显著诱导了HCT的表达。高表达的HCT可能会抑制代谢流进入类黄酮合成分支，从而降低日光萎凋叶中类黄酮合成前体物质的丰度，减少类黄酮的积累量。而在之前的研究中也发现，日光萎凋处理后类黄酮和儿茶素的积累量均出现了显著下降^{[1, 37]}。因此，我们推测日光萎凋处理后，类黄酮和儿茶素合成基因C4H、CHI、F3H、ANS和LAR的转录抑制及FLS和木质素合成基因HCT的转录激活，共同抑制了代谢流进入类黄酮和儿茶素合成分支，从而减少了日光萎凋叶中类黄酮和儿茶素的积累水平，上述结论仍需后续进一步的实验验证。此外，类黄酮和儿茶素还是评价乌龙茶滋味品质的关键代谢物。低含量的类黄酮和儿茶素有助于降低茶叶的苦涩味，适度提高茶汤的适口性^[50]。因此，日光萎凋比室内萎凋更有利于制作出滋味醇厚的乌龙茶。

3.2 日光萎凋叶中萜类化合物合成基因的转录增强促进了萜类化合物的积累

萜类化合物是茶叶中的主要香味物质，其含量的高低直接影响茶叶香气品质的形成^[51]。本研究中，萜类化合物合成也是主要的KEGG富集通路之一。在MVA分支中，以乙酰辅酶A为起始底物。而MEP的代谢底物为丙酮酸和甘油醛-3-磷酸。但这两条分支途径的最终产物均为异戊烯二磷酸(isopentene diphosphate, IPP)。同时有报道显示，IPP的生物合成量将直接调节下游萜类化合物的积累量^[52]。因此，高表达的MEP分支基因和MVA分支基因，有利于促进日光萎凋叶中积累更多的前体IPP，并增强下游萜类化合物的合成代谢通量。而位于下游的GGPPS和TPS也是参与合成萜类化合物的关键基因。有报道显示，当植物遭受逆境胁迫后，萜类化合物合成下游相关基因的表达水平均出现了上升，同时萜类化合物的含量也大幅度增加^[53-54]。在日光萎凋下，萎凋叶受到干旱、高温和紫外线等多种胁迫的协同影响，而室内萎凋叶则只受到干旱胁迫的影响。因此，日光萎凋叶相较于室内萎凋叶，在相同萎凋时间内受到逆境胁迫影响的程度更深。Zhao等^[55]发现萜类化合物还是植物遭受环境胁迫影响后产生的通讯信号，提示附近尚未受到影响的植物立即启动防御系统，减轻逆境胁迫带来的负面影响。因此，胁迫程度更深的日光萎凋处理快速激活了萜类化合物合成下游基因GGPPS和TPS的表达水平，将上游大量积累的IPP，迅速转化为具有逆境警示信号作用的萜类化合物，以缓解环境胁迫造成的不利影响，也从侧面进一步促进了日光萎凋叶中萜类化合物含量的显著上升，而在最新的报道中发现，日光萎凋处理后萜类化合物含量的大幅度增加也验证了我们上述的观点^{[1, 37]}。因此，萜类化合物合成基因的转录激活可能是日光萎凋叶相较于室内萎凋叶具有更高含量萜类化合物的关键所在，而日光萎凋叶中高含量的萜类化合物又有利于促进高品质花果香的形成，这就解释了在乌龙茶萎凋过程中，日光萎凋比室内萎凋更易生产出高品质香气的乌龙茶。

3.3 茉莉酸信号转导和可变剪接机制可能参与了乌龙茶萎凋过程中风味品质的形成

在植物界中，植物激素是参与胁迫响应的重要调节因子^[56-57]。植物激素生物学功能的实现，不仅与植物激素的合成有关，而且与下游的植物激素信号转导也密切相关。植物激素信号转导通路在整合和协调上游分子信号，诱导和激活下游转录因子方面也扮演着重要的角色^[58-59]。本研究中，我们也发现差异表达基因也主要富集于植物激素信号转导通路。其中，茉莉酸信号转导中共有4个基因JAR、COI1、JAZ和MYC2在3组样品中的表达水平存在显著差异。有报道显示，茉莉酸可以作为激活剂，通过促进挥发性萜烯类化合物的合成来改善植物的香气^[60-63]。此外，在逆境胁迫影响下，植物体内茉莉酸的积累水平会迅速升高^[64]。而萎凋处理后，萎凋叶相对于鲜叶，茉莉酸的积累量出现显著上升^[8]。而高含量的茉莉酸还需要通过植物激素信号转导途径，才能激活下游挥发性萜烯类化合物的合成通路，促进萜烯类化合物的积累。在日光萎凋叶中上调表达的JAR和COI1，有利于将激素信号继续向下游传递。而位于COI1下游的JAZ，主要负责阻碍MYC2的正常转录。在低含量茉莉酸时，JAZ阻碍茉莉酸信号的转导，是茉莉酸信号转导通路的负调控因子^[65]。而在胁迫程度更深的日光萎凋处理后，高表达的JAR和COI1可通过传递上游的茉莉酸信号，启动SCF^COI1泛素连接酶复合体进行JAZ的识别捕获，再将捕获后的JAZ通过泛素-26S蛋白酶体途径进行降解，从而抑制了日光萎凋叶中JAZ的正常表达。随后，JAZ对下游MYC2的抑制作用被解除，促进了MYC2转录水平的大幅度上升。有报道发现，过表达MYC2会增强TPS的转录水平，促进萜类化合物含量升高^[66]。这些结果表明，日光萎凋处理可能通过激活茉莉酸信号转导途径的JAR-COI1-JAZ-MYC2级联调控模式，进而促进了日光萎凋叶中萜类化合物的大量积累，为后续高香型品质的形成奠定了物质基础。而上述级联调控模式与萜类化合物含量变化的具体关系，我们将在后续实验中进一步验证。

可变剪接被认为是一种选择性的表达调控机制，在植物应对生物和非生物胁迫的反应中起着重要作用^[67]。在茶树中，已有研究发现逆境胁迫会诱导大量可变剪接事件的发生^{[26, 68]}。与之前的研究一致，在类似逆境胁迫的萎凋处理后，日光萎凋叶和室内萎凋叶中鉴定出的可变剪接事件数量相较于鲜叶均出现了大幅度增加。同时，剪接体通路组成基因Prp5、Prp22、Sm和SR的转录水平在萎凋处理后均出现了显著上升，在日光萎凋叶中的表达量到达最高值。有研究表明，剪接体通路中核心基因转录丰度的提高有利于加快可变剪接的加工效率，从而产生更多的可变剪接基因^[69-71]。因此，我们推测日光萎凋处理后上调表达的Prp5、Prp22、Sm和SR可能是日光萎凋叶中可变剪接基因数量多于鲜叶和室内萎凋叶的原因。随后，对差异剪接基因进行KEGG富集分析发现，类黄酮合成和植物激素信号转导是萎凋过程中最为富集的两条通路。进一步分析两条通路中4个差异剪接基因的不同转录本发现，PTC型剪接转录本在萎凋过程中的表达模式与其对应的全长转录本基本一致，但在各个样品中均呈现较低的表达水平，可能与PTC引起的基因关键保守结构域不完整或无义介导的mRNA降解机制^[72]有关，使发生可变剪接之后的转录本失去原有生物学功能。此外，可变剪接转录本CsFLS2-a、CsFLS2-b及CsFLS2-c与对应全长转录本CsFLS2在3个样品中呈现截然相反的表达水平。同时，CsFLS2在3个样品中表达水平无显著变化，而3个可变剪接转录本在萎凋处理后表达水平显著升高，暗示这3个可变剪接转录本可能在乌龙茶萎凋过程中替代全长转录本发挥主要调控作用，而这一现象在拟南芥生物钟相关基因的研究中也有发现^[73]。因此，日光萎凋处理后高表达的可变剪接转录本CsFLS2-a、CsFLS2-b和CsFLS2-c可能减少了进入儿茶素合成分支的代谢通量，进而减少了日光萎凋叶中儿茶素的积累量。后续将展开对特定转录本的功能验证，进一步探究萎凋处理后可变剪接机制与儿茶素合成之间的调控关系。

本研究通过对鲜叶、室内萎凋叶和日光萎凋叶进行转录组测序分析，发现日光萎凋处理后类黄酮合成基因的转录抑制和萜类化合物合成基因的转录增强，可能是促进了日光萎凋叶中苦涩类物质类黄酮含量减少与花果香物质萜类化合物含量上升的主要原因。同时，茉莉酸信号转导和可变剪接机制可能也共同参与了日光萎凋叶中高花果香和低苦涩味的风味品质形成。综上所述，本研究揭示了日光萎凋处理在乌龙茶优异风味品质形成的重要性，以及转录调控和可变剪接在风味代谢相关途径中的具体作用机制。