千针万线草根转录组测序及黄酮类化合物合成相关基因挖掘

孙诗瑶，王晓丽，曹子林 ，张博林，郭紫微，赵思懿

（1.西南林业大学生态与环境学院，云南 昆明 650224；2.西南林业大学林学院，云南 昆明 650224）

0 引言

【研究意义】千针万线草（Stellaria yunnanensis）为石竹科（Caryophyllaceae）繁缕属（Stellaria）的一种多年生草本植物，其根具有补气健脾、养肝活血的药用价值[1]。据研究，繁缕属植物普遍含有大量的黄酮类化合物，从中分离出的化合物主要有黄酮类、异黄酮类和黄酮醇类物质3类，具广谱抗病毒药物的功效[2]。因此，研究千针万线草根的转录组和黄酮类代谢，对其开发利用具重要意义。【前人研究进展】谭瑞璞等[3-4]优化了千针万线草中总黄酮的提取工艺，并测定不同产地千针万线草中牡荆素的含量。在模式植物中，与类黄酮合成相关的功能基因已被克隆并进行广泛的研究[5]。目前，黄酮类化合物的生物合成路径已较为明确[6]。而就缺乏基因组信息的物种而言，对于植物次生代谢中类黄酮合成途径的研究，多停留在红花、中麻黄、半枝莲、金线莲、铁皮石斛、藤茶等药用植物的基因挖掘和鉴定上[7-12]。【本研究切入点】千针万线草基因组和转录组数据缺乏，制约了该物种种质资源鉴定、遗传多样性分析和分子机制等研究的开展。因此，对千针万线草根进行高通量转录组测序及基因功能注释，可为后续研究奠定基础。【拟解决的关键问题】利用Illumina Novaseq 6000高通量测序平台对千针万线草根进行转录组测序，应用公共数据库对其基因组进行功能注释及分类，借助KEGG挖掘出该物种黄酮类化合物合成相关基因，为开展黄酮类化合物合成的调控机制研究提供依据。

1 材料与方法

1.1 样品采集

千针万线草（Stellaria yunnanensis）根采自会泽县迤车镇箐口村山坡上，埋于土中带回学校，种于花盆放在苗圃，每天正常浇水萌出根蘖苗后，经西南林业大学曾觉民教授鉴定为千针万线草。于2021年7月取其根部以无菌水洗净，滤纸吸去水分后装入冻存管，液氮速冻、备用。

1.2 RNA提取、文库的构建及测序

提取根的total RNA后，检测RNA的浓度、纯度和完整性。参照文献[13]的方法，将带有Oligo（dT）的磁珠和ployA配对A-T碱基。mRNA从total RNA中分离后，富集得到完整的RNA序列，添加fragmentation buffer随机断裂mRNA，利用磁珠筛选分离获得约300 bp小片段；加入六碱基随机引物，在逆转录酶存在下，以mRNA作为模板，反转录合成一链cDNA，再进行二链合成，进而形成稳定的双链结构。连接adaptor，对目的条带大小合适的片段用2%的琼脂糖胶回收后，再用PCR技术进行扩增；最后，对短序列片段用Illumina Novaseq 6000进行测序。

1.3 测序原始数据过滤及组装

测序得到原始数据后，用软件Seq Prep将数据中带有的大量接头序列、低质量序列末端、adapter以及N的比例高于10%的碱基去除，获得高质量的clean reads。使用Trinity软件[14]组装denove，用TransRate[15]评估从头拼接结果中普遍存在的错误，获取每一条contig的质量评分，并将这些contigs评分整合后，得到整个组装结果的分数，筛选优化转录组获得的序列。利用CD-HIT软件[16]对序列比对聚类，剔除冗余、相似序列，最终获得非冗余的序列。用BUSCO软件[17]评价基因组或者转录组的组装完整性。

1.4 基因功能注释

利用NCBI将Unigene序列比对到Nr、Pfam、SwissProt和 egg NOG数据库中（evalue＜0.000 01），用BLAST2 GO软件对Unigene序列进行GO注释，利用KOBAST软件对Unigene序列进行KEGG注释[18]，并分类统计注释结果。用TF数据库对预测的蛋白序列进行hmmscan比对并得到转录因子家族。

2 结果与分析

2.1 千针万线草根转录组组装质量分析

对千针万线草根进行转录组测序后，经过严格的质量控制及数据过滤得到55 400 486个高质量且干净的序列。用Trinity软件将这些高质量数据进行de novo组装，共获得37 331 463条有效序列片段，N50为1 799 bp，平均长度为1 093.58 bp，GC核苷酸的含量（GC%）为41.02%，千针万线草根Unigene34 137条（不含N的组装片段）长度分布于201～14 683 bp（表1）。用长度分布柱状图分析组装后的基因如图1所示，长度为200～500 bp、501～1 000 bp、1 001～1 500 bp、1 501～2 000 bp、2 001～2 500 bp、2 501～3 000 bp、3 001～3 500 bp、3 501～4 000 bp、4 001～4 500 bp、大于4 500 bp的unigene分别占40.2%、20.2%、13.3%、10.3%、6.7%、3.8%、2.2%、1.3%、0.8%、1.2%。

图1 千针万线草根转录组Unigenes的长度分布Fig.1 Distribution of length of unigenes from S.yunnanensis roots

表1 Unigene基本信息表Table 1 Basic information on unigenes

2.2 千针万线草根基因的总体注释情况

将组装得到的基因在Nr、Swiss-Prot、Pfam、egg NOG、GO、KEGG数据库中逐一进行比对，所得到注释的基因条数分别为21 510、17 208、17 702、19 414、17 191、10 285条。至少有22 369条unigenes在一个数据库中进行了匹配和功能注释，有11 768条unigenes未比对到数据库，这表明功能未知的基因还较多。相比之下，注释到Nr数据库的基因数最多有21 510条，占比为63.01%；注释到KEGG数据库的基因最少，共有10 285条，占比为30.13%（表2）。

表2 Unigene注释统计Table 2 Annotation of unigenes

2.3 千针万线草根基因的NR功能注释

共有21 510条unigenes比对到NR数据库，与其他物种比对后的结果如图2所示。注释序列物种中注释到基因数量最多的是甜菜（Beta vulgaris），共有6 473条，占比为30.09%；其余依次为藜麦（Chenopodium quinoa）、菠菜（Spinacia oleracea），分别有 6 124、3 954条unigenes，占比分别为28.47%、18.38%，这3种植物占比高达76.94%；其他297种物共有4 959条unigenes，占比23.06%。

图2 千针万线草根基因比对物种分类Fig.2 Species classification of unigenes

2.4 egg NOG功能注释

将千针万线草根基因比对到egg NOG数据库中，得到注释到的Unigenes有19 980条，占总Unigenes的58.53%，被分为23类（图3）。有11 394个未知的功能注释信息，其准确的生物学功能未确定，占所有功能注释信息的57.03%；翻译后修饰、蛋白开关和分子伴侣注释到1 264条unigenes，所占比例为6.33%；转录注释到1 144条unigenes，所占比例为5.73%；信号传导机制注释到895条unigenes，占4.48%；胞内运输、分泌和囊泡运输注释到702条unigenes，所占比例为3.51%；碳水化合物运输和代谢注释到632条unigenes，所占比例为3.16%；翻译、核糖体结构与生物发生注释到615条unigenes，所占比例为3.08%；最少的是细胞运动和细胞核结构的unigenes，分别仅有6条（0.03%）和4条（0.02%）。这些结果表明，千针万线草在翻译后修饰、蛋白开关和分子伴侣、信号转导机制、胞内运输、分泌和囊泡运输、碳水化合物运输和代谢、翻译、核糖体结构与生物发生等基因表达丰度较高。

图3 千针万线草unigene的egg NOG分类Fig.3 Egg NOG function classification of S.yunnanensis

2.5 GO功能注释

注释到GO数据库中的unigenes一共有17 191条，获得69 356个功能注释，共分为生物过程、细胞组分和分子功能三类，依次得到21 481、27 214和20 661个功能注释，分别占比为30.97%、39.24%、29.79%。3个本体细分为53个功能亚类，分子功能类包括16个功能亚类，结合（8 993个注释）所占比例最多，其次是催化活性（8 393个注释）；生物学过程类占其中23个功能亚类，细胞进程（6 958个注释）和代谢过程（5 972个注释）占比较大，其次是生物调控（2 581个注释）；细胞组分类包括14个亚类，“细胞”所占比例最多（8 164个注释），其次是“膜”（5 897个注释）和“细胞器”（4 792个注释）。只有少数转录本被注释参与细胞杀伤、毒素活性、碳利用、介导分子活动、细胞增殖及发信号。这些结果表明，以结合、催化活性、细胞、细胞过程、代谢过程相关的基因较多（图4）。

2.6 KEGG功能注释

在KEGG数据库共注释到6 505个unigenes，可归为6类，与代谢有关的通路共100条、注释到的unigenes有3 190条，所占比例最高（49.04%）；与遗传信息处理有关的通路有21条，注释到的unigenes有2 115条（占32.51%）；环境信息处理相关的通路有4条，注释到的unigenes有409条（占6.29%）；细胞过程相关的通路有6条，注释到的unigenes有501条（占7.70%）；生物体系统相关的通路有2条，注释到的unigenes有259条（占3.98%）；而人类疾病相关的通路有2条，注释到的unigenes有31条（占0.48%）。代谢相关的通路细致分成10个亚类，其中排名前三的是碳水化合物代谢、氨基酸代谢、脂类物质代谢，分别占比为11.79%、7.29%、6.41%，遗传信息加工分为4个亚类，环境信息处理和人类疾病各分为2个亚类，细胞过程和生物体系统的通路各分为1个亚类。其中“翻译”“折叠、分类和降解”“运输和代谢”“能量代谢”以及“信号转导”分别占14.05%、9.99%、7.70%、5.70%、5.21%，其余12个亚类所占比例均小于5%（图5）。

2.7 千针万线草根中黄酮类物质生物合成相关基因挖掘

基于KEGG信号通路富集分析结果如表3所示，在千针万线草根的转录组数据中共筛选获得80条与黄酮类生物代谢相关的基因。其中 18条unigene 映射到了苯丙烷代谢通路（ko00940），45条unigene映射到类黄酮生物合成代谢通路上（ko00941），10条unigene映射到黄酮醇生物合成代谢通路（ko00944），7条unigene映射到异黄酮生物合成代谢通路上（ko00943）。

表3 千针万线草根转录组中黄酮类化合物合成相关基因Table 3 Flavonoid biosynthesis-related gene in transcriptome of S.yunnanensis roots

在苯丙烷代谢通路中，18条unigene分别编码苯丙烷途径中的3个关键酶：PAL、C4H、4CL。其中5条unigene编码该途径中第1个催化酶PAL，5条unigene编码第2个催化反应酶C4H，该酶属于细胞色素P450（CY450）家族。8条unigene编码该途径的第3个酶4CL，其控制苯丙烷类代谢走向不同的代谢途径，是一个重要分支点。在类黄酮生物合成代谢通路中，45条unigene共编码10种该通路下的关键酶，包括查尔酮合成酶（CHS）、查耳酮异构酶（CHI）、黄烷酮-3-羟化酶（F3H）、肉桂酸-4-单加氧酶（CYP73A）、类黄酮-3′单加氧酶（CYP75B1）、根皮苷合酶（PGT1）、花青素还原酶（ANR）、莽草酸羟基肉桂转移酶（HCT）、香豆酰脂3′羟化酶（C3′H）、咖啡酰辅酶A-O甲基转移酶（CCOAOMT）。其中，CHS是催化黄酮类生物合成的第一个限速步骤，是植物聚酮合成酶超家族的成员之一，是植物类黄酮合成途径中的第一个关键结构。CHI制约合成黄酮醇的含量。F3H是黄烷酮分支点的一个核心酶。在黄酮醇生物合成代谢通路中，分别编码该通路下的2种关键酶：类黄酮-3′单加氧酶（CYP75B1）和黄酮醇-3-0-葡萄糖苷葡萄糖基转移酶（FG3）。在异黄酮生物合成代谢通路中，7 条 unigene均编码异黄酮2′-羟化酶（CYP81E）。

2.8 转录组Unigene的转录因子预测

千针万线草转录组共预测到724个转录因子，可分为33个转录因子大家族。其中，MYB_superfamily（116个，16.02%）属于最大家族，其次是ERF（59个，8.15%），接着为bHLH（57个，7.87%）、C2C2（53个，7.32%）、WRKY（52个，7.18%）、bZIP（47个，6.49%）、NAC（39个，5.39%）、C2H2（35个，4.83%）、B3-superfamily（31个，4.28%）、GRAS（31个，4.28%）、LBD（30个，4.14%）、C3H（29个，4.01%），其他21个家族共有基因145个，占20.03%（图6）。

图6 千针万线草根Unigene的转录因子预测Fig.6 Predicted transcription factors of unigenes of S.yunnanensis roots

3 讨论

开展药用植物转录组研究是综合诠释植物基因组信息、发掘植物次生代谢生物合成途径中关键酶以及植物次生代谢成分生物合成途径与调控的新途径[19]。本研究首次对药用植物千针万线草根进行转录组测序，测序产生clean reads的Q20值为98.42%，Q30值为95%以上，说明测序质量合格。经denove拼接组装，去除冗余后，共得到34 137条Unigene，N50为1 799 bp，平均长度为1 093.58 bp。组装结果与杭白芷（Angelica dahurica）根（N50为1 703 bp，平均长度1 164 bp）较为接近[20]。N50是评价组装序列完整性的重要指标，本研究中N50＞800 bp，表明组装片段完整，所得Unigenes 数量基本涵盖全部转录信息，序列信息量庞大，可以满足转录组分析的要求。

甜菜（Beta vulgaris）、藜麦（Chenopodium quinoa）、菠菜（Spinacia oleracea）同属黎科，在Nr数据库注释中，注释到基因占比高达76.94%，可能黎科与石竹科在亲缘关系上比较接近。而石竹科蝇子草属的滇白前注释到黎科这三种植物基因占比达64.06%[21]，千针万线草相比于滇白前，在与黎科亲缘关系上可能较近。在GO功能注释方面，可以分为3个大类和53个亚类，其中结合、催化活性、细胞过程及代谢过程的Unigene数目较多。通过egg NOG功能分类可知，参与翻译后修饰、蛋白开关和分子伴侣、转录、信号传导机制的Unigene较多，但有11 394条未知功能基因，占57.03%。这是由于egg NOG更新比较快，收录来源很广泛，数据库收录的未知序列占比就高，所以注释的结果未知序列占比较高。转录因子是能从DNA上读取并解释遗传“蓝图”的蛋白质组之一，它和DNA相结合，有助于启动调节基因转录的程序，对植物应对非生物逆境胁迫具有重要影响[22]。研究表明，MYB、AP2/ERF、bHLH、bZIP、NAC、WRKY和C2H2转录因子对植物应对干旱、高盐和寒冷胁迫均有显著影响[21]。对这七类转录因子进行统计，千针万线草根中有405个（55.94%），可能与其抗逆适应相关。在野外观察发现，千针万线草在向阳干旱贫瘠的陡坡也能正常生长，这是其抗逆性强的一个方面，是否能适应其他逆境，还有待进一步研究。

通过KEGG数据库和通路分析发现，涉及代谢途径的Unigene（3 190个）最多，通路数量最多的也是涉及代谢途径的，共100条。从本研究获得的基因注释结果来看，有5个生物合成代谢途径80个unigene涉及到药用成分的合成，其中包括黄酮、黄酮醇、异黄酮等物质。通过进一步筛选，这些基因分别编码16种关键酶，包括苯丙烷代谢途径所需的3种酶（PAL、C4H、4CL），类黄酮代谢所需的10种关键酶（CHS、CHI、F3H、CYP73A、CYP75B1、PGT1、ANR、HCT、C3′H、CCOAOMT）和黄酮醇代谢途径所需的2种关键酶（CYP75B1、FG3）和异黄酮代谢途径关键酶CYP81E。千针万线草根的药用成分主要是黄酮类物质，通过对千针万线草根转录组测序、注释、分析以及黄酮类化合物合成相关基因的挖掘，丰富了千针万线草根的遗传信息，可以为进一步研究该物种药用成分合成的关键基因及调控机制、克隆关键基因及功能分析奠定基础，也为千针万线草品质的形成提供理论依据。