怀玉山三叶青2 个栽培种块根的转录组分析

洪森荣，黄丹丹，黄诗慧，黄夏梦，蒋 椰，李万平，蔡 红，陈荣华

（1. 上饶师范学院生命科学学院，江西 上饶 334001；2. 上饶市药食同源植物资源保护与利用重点实验室，江西上饶 334001；3. 上饶市三叶青保育与利用技术创新中心，江西 上饶 334001；4. 上饶农业技术创新研究院，江西 上饶 334001；5. 上饶市红日农业开发有限公司，江西 上饶 334700）

0 引言

【研究意义】三叶青（Tetrastigma hemsleyanumDiels et Gilg）为葡萄科（Vitaceae）崖爬藤属（Tetrastigma）多年生蔓生藤本植物，为我国特有珍稀药用植物[1]，主要分布在长江以南的大部分地区，一般以块根入药，具清热解毒、消肿散结、活血止痛之功效，多用于抗肿瘤、调节免疫力以及治疗小儿热性惊厥、黄疸、急性和慢性肝炎等[2−3]，被誉为“植物抗生素”[4−6]。中药材产区一般与其质量密切相关，不同产区的三叶青外观性状、药理活性及质量和产量均有差异[7]。研究表明[8−11]，怀玉山三叶青不但能减慢癌细胞生长、抑制癌细胞增殖、促进癌细胞凋亡，还可作为癌症预防治疗的功能性食品。由于怀玉山三叶青临床疗效显著，市场需求剧增，野生资源已无法满足其入药需求，为此，上饶市红日农业开发有限公司积极开展品种选育，筛选出两个优良栽培品种怀玉1 号和怀玉2 号。这两个品种与普通种的块根相比具有显著差异。怀玉1 号形成块根早，个头偏小，块根纺锤形，叶薄而小；怀玉2 号形成块根稍晚，个头偏大，块根长柱状，叶厚且大。因此，对怀玉山三叶青2 个栽培种的块根进行转录组分析具有重要的现实意义。【前人研究进展】转录组指植物体内所有表达基因的身份以及其转录水平，广义上指所有信使RNA 、核糖体 RNA 、转运RNA 及非编码RNA，狭义上特指所有信使RNA[12]。转录组测序有助于植物新基因的发现、基因功能的注释、基因的差异表达及其分子标记发展[13−15]。目前，对三叶青的研究主要集中在种质多样性、主效活性成分分析、栽培管理和种苗繁育等方面[16]，在转录组水平上研究三叶青野生种或栽培种与黄酮类化合物相关差异基因的表达研究未见报道。【本研究切入点】药用成分分析表明，怀玉三叶青怀玉1 号和怀玉2 号两者在黄酮类物质含量方面存在较大差异[17]。其分子机制尚无相关报道。怀玉三叶青2 个栽培种在转录组水平上的分子差异还有待深入研究。【拟解决的关键问题】本研究利用Illumina HiSeq 高通量测序技术对怀玉山三叶青怀玉1 号和怀玉2 号块根的转录组进行研究，结合生物信息学分析方法对差异表达基因进行GO 功能富集分析和KEGG 代谢通路富集分析，研究怀玉山三叶青2 个栽培种差异基因的表达，挖掘与黄酮类化合物合成相关的基因，为怀玉山三叶青相关基因的克隆和功能验证以及怀玉山三叶青后续新品种选育等研究奠定理论基础。

1 材料与方法

1.1 材料

怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2 组）块根（同一大棚同一栽培基质同一栽培时间），由上饶市红日农业开发有限公司提供（鉴定人：上饶师范学院生命科学学院王艾平教授，标本存于上饶市红日农业开发有限公司怀玉山基地三叶青种质圃）。

1.2 试验方法

1.2.1 RNA 提取和测序 从怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2 组）块根中提取total RNA，利用Nanodrop2000 对所提RNA 的含量和纯度进行检测。利用上海美吉生物医药科技有限公司Illumina Novaseq 6000 测序平台对完整 cDNA 文库进行测序、序列拼接和组装[17]。

1.2.2 数据分析 转录组分析采用3 次生物学重复。使用 Trinity 软件对clean reads 进行de novo 拼接，形成一个不能在两端扩展的 Unigene 序列。再将序列注释 到 Swiss-Protein、KEGG 和 GO 数 据 库 中 获 得 基因功能信息。根据错误发现率（false discovery rate，FDR≤0.001）和表达量倍数变化（|foldchange|＞1.5）的阈值筛选出两个样本间差异表达的基因（Degs）。对显著性差异表达基因进行 KEGG Pathway 显著性富集分析。

1.2.3 与黄酮类化合物合成相关差异表达基因的筛选通过KEGG 途径分析找出参与了黄酮类化合物生物合成以及黄酮和黄酮醇生物合成的相关差异基因。

1.2.4 查尔酮合酶和黄酮醇合酶的荧光定量 PCR 验证 在1.2.3 的基础上，找出黄酮类化合物生物合成以及黄酮和黄酮醇生物合成的差异基因查尔酮合酶基因和黄酮醇合酶基因，以GAPDH 为内参基因，设计基因引物，查尔酮合酶的引物（F: CACGAGTCCC ACCTCGACTC；R: AGGGGACGCTCCACGGACA），黄酮醇合酶的引物（F: ATATGTACCCACCATGCC CACA；R: CGGCGATCCAGTTATCGTCTT），进行荧光定量 PCR 验证。PCR 反应程序：预变性95 ℃10 min，变性95 ℃10 s，退火延伸60 ℃ 34 s，95 ℃ 15 s，40个循环。

2 结果与分析

2.1 测序数据产出统计

经过测序质量控制， HY1 组Raw reads（统计原始序列数据的条数）为42 869 358，Raw bases（Raw reads 的条数乘以长度，得到测序的总碱基数）为6 473 273 058，Clean reads（计算方法同 Raw Reads，统计的文件为过滤后的测序数据）为42 311 662，Clean bases（Clean reads 的条数乘以长度，得到质控后的总碱基数）为6 293 446 534，Error rate（测序错误率）为0.023 4%，Q20 和Q30（分别计算 Phred 数值大于20、30 的碱基数占总碱基的百分比）分别为98.66%和95.75%，GC content（G、C 碱基的总数量占总碱基数量的百分比）为47.67%；

HY2 组Raw reads 为41 989 034，Raw bases 为6 340 344 134，Clean reads 为41 411 202，Clean bases为6 179 128 919，Error rate 为0.023 2%，Q20 和Q30分别为98.74%和96.02%，GC content 为47.54%。表明，此次转录组测序数据质量较好，可以满足转录组分析的基本要求。

2.2 测序数据与组装结果比对

HY1 组Clean reads（pair reads，过滤后测序数据的条数）为42 311 662，Mapped reads（能比对到组装转 录 本 上 的Clean reads 数）为33 375 480，Mapped ratio（能定位到组装转录本上的Clean reads 所占百分比）为78.88%；HY2 组Clean reads 为41 411 202，Mapped reads 为33 321 456，Mapped ratio 为80.46%。

2.3 转录因子分析

图1 横坐标为不同转录因子家族；纵坐标为落到该转录因子家族的unigene /transcript 数目。从图1可知，HY1 组和HY2 组的转录因子家族多为MYBsuperfamily、bHLH、AP2/ERF、NAC、C2C2、WRKY 等。

2.4 基因表达量分布

获得每个样本中基因的表达量，是后续进行样本间关系分析、表达量差异分析以及实验验证的基础。为了从整体来看某一样本的表达量水平和比较不同样本表达量的整体分布情况，通常以不同类型的图形直观地展示。从图2 的盒形图（A）和密度图（B）可知，HY1 组和HY2 组表达量密度在0-2。从盒形图和密度图的结果来看，怀玉1 号转录物的整体表达水平较高，怀玉2 号转录物的整体表达水平较低。

图 1 HY1 组和HY2 组转录因子家族统计Fig. 1 Statistics on transcription factor families of HY1 and HY2

图 2 HY1 组和HY2 组表达量分布盒形图（A）和密度图（B）Fig. 2 Diagrams of box (A) and density (B) on quantitative distribution of gene expressions in HY1 and HY2

2.5 表达量样本间Venn 分析

通过样本间Venn 分析可以获得样本间或组间的共表达和特表达基因/转录本，能够帮助迅速锁定目标基因。表达量样本间Venn 分析一般以样本间Venn图表示：不同颜色的圆圈代表一组样本中表达的unigene/transcript 数目，圆圈的交叉区域代表各组共有的 unigene/transcript 数目。从图3 可知，HY1 组和HY2 组表达的共有基因数为22367，HY1 组单独表达的基因数为18196，HY2 组单独表达的基因数为8137。HY1 组单独表达的基因数显著高于HY2 组单独表达的基因数。

2.6 表达量样本间相关性分析

表达量样本间相关性分析是取全部基因的每一个基因，根据这些基因在任意两两比较样本的表达量，计算每两个样本之间的相关系数，再把这些相关系数以热图形式直观地展示任意两个样本之间的相关性的过程。该分析有助于理解样本间特别是生物学重复间的相关性。从图4 可知，HY1 组和HY2组表达量的相关系数为0.913，表明基因在样本间的表达量相似度越高，即样本间相关性越好。同时也表明这2 个样本（HY1 组和HY2 组）具有相对相似的表达模式，保证了测序和采样结果的可靠性。

2.7 差异表达分析

HY1 组和 HY2 组共产生差异表达基因（differentially expressed gene，DEG）12199 个。与HY1组比较，HY2 组上调基因数为3 551，下调基因数为8 648。从图5 和表1 可知，HY1 组和HY2 组的差异表达基因大多与光合作用碳同化有关。

图 5 HY1 组和HY2 组表达量差异统计柱状图（A）以及表达量差异火山图（B）和表达量差异散点图（C）Fig. 5 Statistic histogram on expression difference (A), volcano map (B), and scatter map (C) of HY1 and HY2

2.8 差异表达基因GO 富集分析

从图6 和表2 可知，通过富集，光合作用光系统I 光捕获、光合作用捕获、叶绿素代谢过程、蛋白质发色团连锁、前体代谢产物和能量的产生、叶绿素生物合成过程、氧化应激反应、α-氨基酸代谢过程、光合作用、质体小叶、光系统I、光系统II、质体类核仁、光系统、叶绿素结合、单加氧酶活性、铁离子结合、血红素结合、裂解酶活性等GO term富集显著。这些功能类别大多与光合作用相关。

2.9 差异表达基因KEGG 注释

从图7 和表3 可知，HY1 组和HY2 组差异表达基因富集的KEGG 通路主要有光合作用-天线蛋白、核糖体、乙醛酸和二元酸代谢、苯丙酸生物合成、二苯乙烯类、二芳基庚烷类和姜辣素的生物合成、类黄酮生物合成、光合作用、光合生物的固碳作用、甘氨酸、丝氨酸和苏氨酸代谢、植物激素信号转导、谷胱甘肽代谢、丙酮酸代谢、苯丙氨酸代谢、植物昼夜节律、黄酮和黄酮醇的生物合成、半胱氨酸与蛋氨酸代谢、氰胺酸代谢、类胡萝卜素生物合成、α-亚麻酸代谢、卟啉与叶绿素代谢等。其中核糖体（基因数306）、植物激素信号转导（基因数132）、苯丙酸生物合成（基因数100）、丙酮酸代谢（基因数92）、乙醛酸和二元酸代谢（基因数90）、半胱氨酸与蛋氨酸代谢（基因数89）、光合生物的固碳作用（基因数82）注释的基因数均超过80，这可能是HY1 组和HY2 组两个栽培种的内在差异。

2.10 与黄酮类化合物合成相关差异表达基因的筛选

经江西省分析测试中心和南昌大学食品学院检测，怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2 组）块 根 的 总 黄 酮 含 量 分 别 为189 mg·hg−1和286 mg·hg−1，其中儿茶素含量分别为19.54 mg·g−1DW 和24.54 mg·g−1DW、芦丁含量分别为38.21 mg·g−1DW 和46.47 mg·g−1DW、异槲皮素含量分别为25.64 mg·g−1DW 和38.39 mg·g−1DW、槲皮素含量分别为1.59 mg·g−1DW 和3.49 mg·g−1DW、山柰酚含量分别为1.96 mg·g−1DW和5.38 mg·g−1DW、绿原酸含量分别为0.68 mg·g−1DW和1.64 mg·g−1DW。由表4 可知，与黄酮类化合物相关差异表达基因如芪合酶（stilbene synthase） 、无色花色素双加氧酶（leucoanthocyanidin dioxygenase）、查尔酮异构酶蛋白（CHI protein）、查尔酮合酶2（chalcone synthase 2）、黄烷酮3-羟化酶（flavanone 3-hydroxylase）、无色花色素还原酶（1leucoanthocyanidin reductase 1）、类黄酮3′-羟化酶（flavonoid 3′-hydroxylase）基因在怀玉2 号（HY2组）块根中上调，而查尔酮合酶（chalcone synthase）、黄酮醇合酶（flavonol synthase）、类黄酮3′，5′-甲基转移酶（flavonoid 3′,5′-methyltransferase）基因在怀玉2 号（HY2 组）块根中下调，导致怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2 组）块根总黄酮含量的差异。

表 2 HY1 组和HY2 组部分差异表达基因GO 富集结果Table 2 GO enrichment on certain differentially expressed genes in HY1 and HY2

2.11 查尔酮合酶和黄酮醇合酶的荧光定量 PCR 验证结果

查尔酮合酶和黄酮醇合酶的荧光定量 PCR 分析表明，查尔酮合酶基因在怀玉山三叶青怀玉1 号（HY1组）和怀玉2 号（HY2 组）的RQ 值为1 和0.00051，表达下调；黄酮醇合酶在怀玉山三叶青怀玉1 号（HY1组）和怀玉2 号（HY2 组）的RQ 值为1 和0.450961，表达下调。这个结果与转录组的结果是一致的，表明本次转录组分析结果是有效的。

3 讨论与结论

表 3 HY1 组和HY2 组差异表达基因KEGG 富集分析（前20）Table 3 KEGG enrichment on top 20 differentially expressed genes in HY1 and HY2

目前利用转录组分析来研究植物不同品种差异基因的报道很少。差异基因的KEGG 富集分析有助于更进一步了解差异基因的生物学功能[18]。4 个大豆栽培种多荚与少荚品种的花序转录水平差异主要集中在物质代谢、胁迫应答、刺激应答、病害抗性等方面，与脱落直接相关的基因表达不相关[19]。2 个紫花苜蓿品种（WL319HQ 和准格尔）转录组富集的KEGG 途径主要有核糖体、植物激素信号转导等，其中核糖体途径注释的DEGs 最多[20]。3 个不同品种百香果（黄果、满天星和紫果）果皮和果肉的KEGG富集分析表明，黄果与满天星的差异在于细胞胞吞作用，黄果与紫果的差异在于玉米素的生物合成和糖酵解途径，满天星与紫果的差异在于缬氨酸、亮氨酸和异亮氨酸的生物合成[21]。在本实验中，怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2组）块根差异表达基因富集的KEGG 通路主要有核糖体、植物激素信号转导、苯丙酸生物合成、丙酮酸代谢、乙醛酸和二元酸代谢、半胱氨酸与蛋氨酸代谢、光合生物的固碳作用等，这可能是怀玉山三叶青怀玉1 号和怀玉2 号两个栽培种的内在分子差异。

黄酮类化合物是植物体内一类重要的保护性物质，可抵御病原菌侵害、防止紫外线伤害、促进色素形成等[22]，具有软化血管、抗癌、清除体内自由基等功能，已成为植物次生代谢物应用的研究热点[23]。黄酮类化合物主要包括白藜芦醇、黄酮、异黄酮、黄烷酮、查耳酮、花色素苷、原花色素等[24]。

白藜芦醇（resveratrol）是一种植保素，可抵御病菌侵染、增强抗病性[25]。研究表明，一些芪类物质还具有抗炎抗氧化[26]、防癌抗癌[27]等药理活性。许多植物中存在合成白藜芦醇的底物，但缺乏合成白藜芦醇必需的酶——芪合酶[28−29]。吴凤颖等[30]对自然状态下转基因植株中白藜芦醇合成相关基因进行实时定量 PCR 分析，发现过表达芪合酶基因时白藜芦醇积累量大幅提高。在本实验中，芪合酶在怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2组）上调和下调，导致怀玉山三叶青怀玉1 号（HY1组）白藜芦醇合成增加，而怀玉2 号（HY2 组）白藜芦醇合成下降。

表 4 与次生代谢物合成相关差异基因的表达Table 4 Expressions of differentially expressed genes related to secondary metabolites synthesis

无色花色素双加氧酶是植物花色素苷生物合成途径中的关键酶，属于双加氧酶家族，可将无色原花色素催化生成有色花色素，对花色素苷的合成至关重要[31−32]。在山葡萄果实着色过程中，无色花色素双加氧酶基因在不同时期的果皮中均有表达，在转色期的叶片、茎、果肉中也均有表达，且表达量相近[31]。在本实验中，无色花色素双加氧酶在怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2 组）下调和上调，导致怀玉山三叶青怀玉1 号（HY1 组）花色素苷合成下降，而怀玉2 号（HY2 组）花色素苷合成增加。

无色花色素还原酶基因是黄酮类化合物生物合成下游的一个重要基因。Tanner 等[33]通过构建LAR基因的重组载体进行体外表达，发现无色花色素还原酶基因原核表达蛋白具有催化合成儿茶酸、黄烷-3-醇、儿茶素等黄酮类化合物的活性。黄烷酮 3 -羟化酶与花色素合成相关，调控着黄酮与花青素苷产物的合成[34]。前人研究表明，在烟草、拟南芥等模式植物中过量表达外源黄烷酮3-羟化酶基因，能提高转基因植株的类黄酮含量和抗逆能力。在本实验中，黄烷酮3-羟化酶和无色花色素还原酶在怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2 组）下调和上调，导致怀玉山三叶青怀玉1 号（HY1 组）黄酮类化合物合成下降，而怀玉2 号（HY2 组）黄酮类化合物合成增加。

查尔酮合酶是将苯丙烷类代谢途径引向黄酮类化合物合成的第1 个关键酶，催化3 分子的丙二酰辅酶A 与1 分子的香豆酰CoA 生成查尔酮[35]。查尔酮异构酶（CHI）催化查尔酮异构化，形成生物活性二羟基黄烷酮，是植物黄酮合成途径的关键酶，推动下游黄酮类化合物的合成[36]。在本试验中，查尔酮合酶在怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2 组）上调和下调，而查尔酮合酶2 和查尔酮异构酶（CHI）在怀玉山三叶青怀玉1 号（HY1组）和怀玉2 号（HY2 组）下调和上调，这可能是由于查尔酮合酶与查尔酮合酶2 基因虽属于同一基因家族但执行不同功能，其具体的功能差异还有待进一步的分析，查尔酮异构酶（CHI）可能促进了怀玉山三叶青怀玉2 号（HY2 组）二羟基黄烷酮的合成。

类黄酮3′-羟化酶为细胞色素 P450 亚家族成员之一，影响着黄酮类化合物活性和稳定性的羟基化模式，被认为是黄酮类化合物生物合成途径中关键的限速酶，其表达水平决定了黄酮类化合物含量，是黄酮类化合物生物合成的重要调控点，常与查尔酮异构酶（CHI）、查尔酮合酶 等共同作用合成各类黄酮类化合物[37]。在本实验中，类黄酮3′-羟化酶在怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2 组）下调和上调，导致怀玉山三叶青怀玉1 号（HY1 组）各类黄酮类化合物合成下降，而怀玉2 号（HY2 组）各类黄酮类化合物合成增加。

类黄酮是次生代谢形成的最主要功能性成分之一。甲基化是类黄酮基本结构形成后由类黄酮3′，5′-甲基转移酶催化的最主要的修饰反应之一，甲基化与酰基化、糖基化等基本修饰反应导致了类黄酮种类及功能的多样性，可降低类黄酮活性基团的化学反应活性，增加其脂亲和性，扩大其在细胞内的分布范围，同时提高其抗菌性[38]。在本试验中，类黄酮3′，5′-甲基转移酶在怀玉山三叶青怀玉1 号（HY1 组）和怀玉2 号（HY2 组）上调和下调，导致怀玉山三叶青怀玉1 号（HY1 组）类黄酮种类及功能的多样性提高，而和怀玉2 号（HY2 组）类黄酮种类及功能的多样性下降。