基于高通量测序的不同生态环境东方泽泻块茎转录组分析△

韦陈彬，谷巍，2*，田荣，邱蓉丽，韩赟，李桃，季媛媛，朱瑜

1.南京中医药大学 药学院，江苏 南京 210023；

2.南京中医药大学 江苏省中药资源产业化过程协同创新中心，江苏 南京 210023；

3.南京中医药大学 附属苏州市中医医院，江苏 苏州 215009

建泽泻为植物东方泽泻Alisma orientale（Sam.）Juzep.的干燥块茎，始载于《神农本草经》，具有利水渗湿、泻热、化浊调脂的功效。其主要功效物质为原萜烷型三萜类化合物，该类成分具有利尿、降血糖、调血脂、降血压、抗癌等作用[1]。福建的建泽泻是公认的道地药材，但由于建泽泻生长周期长、生产成本高等原因，目前市场上超过90%的泽泻为来自四川的川泽泻。

转录组测序技术是功能基因组研究的一种重要手段，其能快速、准确地获取大量特定细胞或组织中表达的基因序列，可在无参考基因组的情况下，直接对药用植物进行转录组分析。转录组测序分析可为研究药用植物的基因功能和次生代谢产物代谢途径积累基础，在三七、川续断等药用植物研究中均有应用。徐若等[2]使用转录组测序分析了三七应对干旱胁迫时差异基因的变化。朱春艳等[3]通过转录组技术对川续断三萜皂苷生物合成路径上的关键酶进行挖掘。

现阶段建泽泻产量下降趋势短时间内难以逆转，为延续建泽泻道地药材的优良种质，本研究将东方泽泻同时种植于福建和四川，基于转录组学技术对两地东方泽泻的块茎进行测序，从中挖掘差异表达基因，探讨不同生态环境对东方泽泻三萜生物合成的影响。

1 材料

1.1 试药

选用本课题组前期基于内转录间隔区2（ⅠTS2）条形码鉴别的种植于福建建瓯的东方泽泻作为实验材料[4]，经育苗过程后于2020年9月同时移栽至福建省吉阳镇玉溪村（N27°8′21″，E118°8′27″）和四川省谢家镇吴堰村（N30°13′20″，E103°47′13″）继续种植，并于2020年12月泽泻药材采收期收集两地东方泽泻块茎样品，每组样品设置3 个生物学重复，样品暂存于干冰中，收集完成后转移至-80 ℃冰箱保存，以备后续分析。样品经南京中医药大学谷巍教授鉴定均为东方泽泻Alisma orientale（Sam.）Juzep.的干燥块茎。

OmniPlant RNA 提取试剂盒（北京康为世纪生物科技有限公司）；RNA 6000 Nano kit（美国Agilent公司）。

1.2 仪器

Novaseq 6000型测序平台（美国Ⅰllumina公司）；NanoDrop 2000 型超微量紫外-可见分光光度计、2100 Nano 型生物分析仪（美国Agilent 公司）；GeneAmp 9700 型聚合酶链式反应（PCR）仪（美国ABⅠ公司）。

2 方法

2.1 总RNA提取及cDNA文库构建

使用OmniPlant RNA 提取试剂盒提取东方泽泻块茎中总RNA；NanoDrop 2000 型超微量紫外分光光度计检测总RNA 浓度和纯度；琼脂糖凝胶电泳检测总RNA的完整性；Agilent 2100 Nano和RNA 6000 Nano kit 检测总RNA 的RNA 完整值（RⅠN）。使用Oligo（dT）磁珠分离mRNA，再将其断裂为300 bp左右的短片段，以其为模板合成cDNA第1链，置换合成法获得双链cDNA，对纯化后双链cDNA 进行PCR扩增，构建测序文库。

2.2 测序数据质控和组装

采用Ⅰllumina测序平台对构建的cDNA 文库进行高通量测序，测序由上海美吉生物医药科技有限公司完成。使用Fastp 0.19.5软件对原始数据进行质量控制以获得高质量的clean reads。使用Trinity v2.8.5软件对所有clean reads 进行从头组装，组装完成后使用TransRate v1.0.3软件对拼接结果进行评估，对从头组装的序列进行过滤和优化，最后使用CD-HⅠT v4.5.7软件去除冗余基因，获得非冗余基因集。

2.3 基因功能注释及差异基因筛选

使用Blast 2.9.0软件将unigenes与6个数据库进行基因功能注释：Pfam（http://pfam.xfam.org/）、SwissProt（http://web.expasy.org/docs/swiss-prot_guideline.html）、京都基因与基因组百科全书（KEGG，http://www.genome.jp/kegg/）、COG（http://www.ncbi.nlm.nih.gov/COG/）和基因本体（GO，http://www.geneontology.org）。转录因子信息通过与Plant TFDB（http://planttfdb.cbi.pku.edu.cn/）进行比对，得到同源的转录因子信息。以TPM（transcripts per million reads）作为基因表达量定量指标，采用Benjamini-Hochberg 校正方法对假设检验得到的显著性P值进行校正，并以|log2FC|＞0.58（FC 为差异倍数）、P＜0.05 为阈值筛选不同生态环境下东方泽泻的差异表达基因（different expressed genes，DEGs），并对DEGs 进行GO 和KEGG 富集分析。

3 结果

3.1 转录组测序结果

采用Ⅰllumina RNA-seq 技术进行转录组测序。在去除低质量的序列后，从6 个样品中共获得超过3.47 亿条clean reads。其错误率＜0.024 0%、质量值＞30 的碱基占比（Q30）高达95%（表1）。结果表明高通量测序数据质量高，可以进行下一步分析。

表1 东方泽泻转录组测序结果

3.2 转录组功能注释及差异基因筛选

使用Trinity v2.8.5 软件对所有有效序列进行从头组装后获得了50 652 条unigenes，其中22 689 条unigenes 被注释到6 个功能数据库中，占总数的44.79%，其中在NR 数据库中注释到22 347 条、SwissProt数据库中注释到15 399 条、Pfam 数据库中注释到16 365 条、COG 数据库注释到18 924 条、GO 数据库注释到19 017 条、KEGG 数据库注释到8770 条，得到东方泽泻转录组中的基因功能、代谢通路等数据。以|log2FC|＞0.58、P＜0.05为筛选条件，筛选不同生态环境下东方泽泻中的DEGs。结果显示福建东方泽泻较四川东方泽泻存在3878 条DEGs，其中2082 条表达上调、1796 条表达下调，得到的DEGs数量较多，说明不同生长环境对东方泽泻基因表达影响较大。

3.3 差异基因的GO富集分析

对不同生态环境东方泽泻的DEGs 进行GO 富集分析，并显示富集程度最高的前20 条GO 条目（图1），包括多糖代谢过程、细胞外区、细胞壁组织或生物发生、外部封装结构组织、细胞壁组织、化学稳态。此外还富集到大量与植物内稳态及多糖代谢相关的条目。结果表明，东方泽泻移栽至四川后能够通过调控体内稳态基因的表达提高其适应新环境的能力；此外生态环境变化还改变了东方泽泻体内代谢过程，这同样是东方泽泻适应新环境的体现。

图1 不同生态环境下东方泽泻DEGs的GO富集分析

3.4 差异基因的KEGG富集分析

对DEGs 进行KEGG 富集分析，结果显示，DEGs 共注释到126 条KEGG 通路，其中显著富集的主要为植物激素信号转导、淀粉和蔗糖代谢、苯丙素生物合成、氰基氨基酸代谢、植物-病原体相互作用、牛磺酸和亚牛磺酸代谢等通路。另外富集到萜类生物合成相关的KEGG 通路，包括萜类骨架生物合成、倍半萜和三萜生物合成、单萜生物合成和二萜生物合成。KEGG 富集分析结果表明，东方泽泻移栽至四川后其次生代谢过程发生明显变化，包括三萜类成分的生物合成过程，这一变化可能会造成泽泻药材的质量变化（图2）。

图2 不同生态环境下东方泽泻DEGs的KEGG富集分析

3.5 参与植物生理活动的DEGs分析

对DEGs 在植物生理活动中的表达情况进行分析，以探索东方泽泻在不同生态环境下生理活动的变化，DEGs 功能主要集中在防御反应、光合作用、碳水化合物代谢和糖酵解等方面（表2）。

表2 福建和四川东方泽泻中差异关键基因的名称、功能及表达情况

植物防御反应是指植物抵御外部侵染、机械损伤等细胞生理反应，能够减少外界环境、病原物等对植物组织的危害，保证植物体的生长发育，有利于次生代谢产物积累[5]。本研究对不同生态环境下东方泽泻参与防御反应的DEGs 进行分析，结果显示有28 条与防御反应相关的基因表达差异有统计学意义，其中在福建东方泽泻中上调的有21 条，下调的有7 条，主要是参与植物生物和非生物胁迫的MLO基因和参与植物-病原体相互作用过程的RPS基因。参与植物防御反应过程的基因在植物生长发育过程中发挥了重要作用，环境变化之后防御反应相关基因下调明显，这可能会影响东方泽泻的生长发育。

药用植物通过光合作用能够获得维持植物生长必需的初生代谢产物和主要活性成分来源的次生代谢产物，因此光合作用与药用植物质量形成密切相关[6]。本研究对不同生态环境下与东方泽泻光合作用相关的DEGs 进行分析，结果显示有12 条与光合作用相关基因在福建和四川东方泽泻中的表达存在差异，其中3条在福建东方泽泻中上调表达，9条下调表达。说明生态环境变化影响了东方泽泻中与光合作用相关的关键酶基因表达，影响了东方泽泻的光合效率。

碳水化合物不仅是植物生长发育过程中重要的能源提供者，而且还为次生代谢产物的生物合成提供碳骨架并具有重要的调节作用[7]。有研究表明，碳水化合物含量对块茎发育具有重要的影响，淀粉是植物体中碳水化合物的主要储存形式[8]。本研究中，在蔗糖和淀粉代谢途径中共发现38 条DEGs，福建东方泽泻中上调的有25条，下调的有13条。其中淀粉通路上的关键酶淀粉合成酶和颗粒结合淀粉合酶（granule-bound starch synthase，GBSS）基因在福建东方泽泻中高表达，而蔗糖代谢过程中的关键酶蔗糖合成酶（starch synthase，SuS）和蔗糖磷酸合成酶（sucrose phosphate synthase，SPS）基因在四川东方泽泻中高表达。此外，参与碳水化合物代谢的α-半乳糖苷酶（α-galactosidase，α-GAL）和β-半乳糖苷酶（β-galactosidase，β-GAL）基因均在福建东方泽泻中高表达。淀粉是东方泽泻块茎中的主要组成，能够为其次生代谢提供能量和碳骨架，以上结果说明东方泽泻移栽四川可能不利于其淀粉合成，影响泽泻药材的质量。

在植物体内，葡萄糖的糖酵解反应能够产生乙酰辅酶A（acetyl coenzyme A，Acetyl-CoA），Acetyl-CoA是三萜类成分生物合成的重要前体，对次生代谢产物的合成效率存在重要影响[9]。2 组东方泽泻在糖酵解途径上共存在20 个DEGs，其中12 个在福建东方泽泻中上调表达，8 个下调表达。参与糖酵解过程的关键酶磷酸果糖激酶（phosphofructokinase，PFK）基因在四川东方泽泻中高表达，丙酮酸激酶（pyruvate kinase，PK）基因在福建东方泽泻中高表达。结果表明，移栽到四川对东方泽泻的糖酵解过程存在一定的影响，可能会影响东方泽泻次生代谢过程。

3.6 差异表达转录因子分析

2 组东方泽泻共注释到568 条转录因子，其中163 个转录因子表达差异有统计学意义，在福建东方泽泻中上调表达的有67 个，下调表达的有96 个。差异转录因子被归类为24 个转录因子家族，2 种生态环境下差异表达的转录因子主要集中在AP2/ERF、MYB_superfamily、NAC、WRKY、C2C2 和bHLH家族（图3）。植物体内，调控三萜类成分生物合成的转录因子家族主要包括bHLH、MYB_superfamily、AP2/ERF 和WRKY 4 类[10]，占差异转录因子总数的50.9%。福建东方泽泻在AP2/ERF 家族上调表达20个，下调表达13 个；bHLH 家族上调表达4 个，下调表达6个；MYB家族上调表达15个，下调表达13个；WRKY家族上调表达5个，下调表达7个（图4），说明改变东方泽泻生长环境会造成与三萜生物合成相关转录因子的差异表达，影响原萜烷型三萜生物合成，造成两地东方泽泻原萜烷型三萜含量差异。

图3 不同生态环境下东方泽泻差异表达转录因子数量

图4 三萜类生物合成相关转录因子家族在福建东方泽泻中差异表达情况

3.7 参与三萜生物合成DEGs分析

东方泽泻的药效成分为原萜烷型三萜类成分，主要通过细胞质中的甲羟戊酸（mevalonate，MVA）途径合成[11]。福建和四川东方泽泻转录组中存在19个与萜类生物合成相关DEGs，其中包括5 个与原萜烷型三萜生物合成相关的关键酶基因（图5）。其中，羟甲基戊二酰CoA 合酶（hydroxymethylglutaryl-CoA synthase，HMGS）、甲羟戊酸激酶（mevalonate kinase，MK）、甲羟戊酸焦磷酸脱羧酶（mevalonate pyrophosphate decarboxylase，MPD）等位于MVA 途径上的基因表达量在福建东方泽泻中显著高于四川东方泽泻，而位于2-C-甲基-4-磷酸-D-赤藓糖醇（2-C-methyl-Derythrito-4-phosphate，MEP）途径的1-羟基-2-甲基-2-(E)-丁烯基-4-焦磷酸合酶 [1-hydroxy-2-methyl-2-(E)-butenyl 4-diphosphate synthase，HDS] 基因表达量在四川东方泽泻中显著高于福建东方泽泻，鲨烯合酶（squalene synthetase，SS）基因表达量同样在四川东方泽泻中高于福建东方泽泻。表明不同生态环境下东方泽泻原萜烷型三萜生物合成通路上关键酶基因的表达存在差异造成两地东方泽泻中原萜烷型三萜含量的变化。目前，已经在多种植物中成功鉴定了细胞色素P450（CYP450）酶基因，其中CYP51、CYP71、CYP72 和CYP85 家族基因参与三萜生物合成[12]，本研究中鉴定得到6个注释到CYP71和CYP85 家族的DEGs，均在四川东方泽泻中上调表达，该特征可能会导致2 组东方泽泻三萜组成差异形成。

图5 不同生态环境下东方泽泻原萜烷型三萜生物合成关键酶基因差异表达情况

4 讨论

植物所处生态环境改变时，植物内稳态能够提高植物对外界环境的适应能力，这一过程常涉及到植物形态和生理的变化[13]，能够直接或间接影响药用植物活性成分生物合成和药材质量。本研究开展了不同生态环境下东方泽泻的转录组测序，从分子层面探讨生态环境对东方泽泻生理活动及药效成分积累的影响。

不同生态环境下东方泽泻在防御反应、淀粉和蔗糖代谢、糖酵解/糖异生和光合作用等生理活动相关通路上存在差异表达。东方泽泻移栽至四川后参与防御反应的基因表达明显下调，说明东方泽泻在四川种植的生长过程中，对于外界环境的防御能力下降，影响了其自身的生长。淀粉能够为次生代谢产物合成和积累提供更多的底物，促进植物的次生代谢产物积累[7]。东方泽泻块茎中含有丰富的淀粉，质量分数约为25%，在福建生长的东方泽泻中与淀粉合成相关的关键酶基因表达上调，说明东方泽泻在原产地生长更有利于淀粉的生物合成与积累，有益于次生代谢产物的积累。此外，在福建东方泽泻中上调表达的α-半乳糖苷酶也在碳水化合物代谢的生理调节过程中发挥重要作用[14]。糖酵解途径相关基因上调能够为次生代谢产物的合成提供中间产物，还能为植物次生代谢提供能量，从而影响到植物的次生代谢[9]。这2 条途径上关键酶基因的差异表达可能会影响三萜的合成与积累，造成两地东方泽泻原萜烷型三萜的含量差异。

转录因子是介导外界环境因素信号转导与药用植物次生代谢的一类调控因子[15]。东方泽泻移栽至四川后转录因子表达发生显著改变，尤其是被证明能够调控三萜生物合成的转录因子家族基因，其数量在差异转录因子中占半数以上。鉴定到的差异表达转录因子主要属于AP2/ERF 和MYB 家族，占差异转录因子总数的37.4%。转录因子的活性将影响关键酶的表达，从而调控药用植物中药效成分的积累。有研究表明，AP2/ERF 转录因子家族中的AaERF1、AaERF2 和AaTAR1 能够参与调控青蒿中青蒿素的合成[16]；丹参SmMYB9b 转录因子通过促进SmDXS2、SmDXR、SmGGPPS 和SmKSL1 的表达从而提高丹参酮的含量[17]。在本研究中，不同生态型东方泽泻中bHLH、MYB_superfamily、AP2/ERF 和WRKY 等转录因子的表达存在明显差异，表明其在不同生态环境作用下参与调控东方泽泻的生理活动和药效成分积累。

MVA 途径和MEP 途径合成的焦磷酸异戊烯酯（isopentenyl diphosphate，ⅠPP）和焦磷酸二甲基烯丙酯（dimethylallyl diphosphate，DMAPP）能够通过质体膜互相交换以此调节萜类化合物的合成，因此MVA 途径、MEP 途径和下游三萜生物合成路径上的基因表达差异是造成原萜烷型三萜含量变化的主要原因。本研究中，MVA 途径上的DEGs 主要在福建东方泽泻中高表达，MEP 途径上DEGs 在四川东方泽泻高表达，说明移栽至四川改变了东方泽泻中参与原萜烷型三萜生物合成的基因表达水平，可能不利于原萜烷型三萜的生物积累。CYP450为超蛋白家族，能够催化多种代谢反应，包括氧化、去饱和、C-C 键裂解等催化反应。研究发现，CYP450 通过对萜类骨架进行复杂的修饰，可以生成不同类型的萜类化合物[12]。本研究中，注释到CYP71 和CYP85 家族的DEGs 在四川东方泽泻中显著高于福建东方泽泻，这些特征可能会促进四川东方泽泻中的23-乙酰泽泻醇B（alisol B 23-acetate）向其他三萜化合物转化[18]，导致四川东方泽泻中23-乙酰泽泻醇B 含量较福建降低，而其他原萜烷型三萜化合物含量升高。本研究丰富了东方泽泻转录组数据库并初步筛选出部分可能参与东方泽泻适应环境变化的关键基因，主要涉及到植物生理活动和次生代谢，为阐明生态环境对东方泽泻原萜烷型三萜生物合成的调控机制研究提供参考，为东方泽泻异地引种和种质保护提供理论依据。