超声诱导对发芽花生的转录组分析及苯丙烷类合成相关基因的挖掘

解梦汐，于淼，鲁明，付欣，张良晨，石太渊

（辽宁省农业科学院食品与加工研究所，辽宁 沈阳 110034）

花生（Arachis hypogaeaL.）又名“长生果”，是我国一种重要的经济作物和油料作物，含有大量蛋白质、膳食纤维和不饱和脂肪酸。花生发芽过程中总酚含量和反式白藜芦醇等生物活性成分增加，具有高抗氧化活性，因此发芽花生被专家视为理想的新型功能芽菜品种，也称花生芽[1]。花生在发芽过程中通过调控苯丙烷代谢途径提高苯丙烷类化合物如酚酸、白藜芦醇和黄酮类物质的含量。大量研究证实了酚类物质的产生有助于植物应对多种外界生物和非生物的胁迫，而这些物质的产生涉及植物次生代谢过程信号传导途径的调控作用[2，3]。

在实际生产中，为解决花生芽生产周期长、活性成分富集消减速度慢等问题，国内外学者通过利用物理场诱导如超声处理[4-6]、紫外照射[7，8]、外源添加[9]、高压静电场[10]等手段对花生或发芽花生的生物活性成分进行富集。Sales等[8]研究表明，经过超声诱导的花生总酚含量从0.84 mg GAE/g提升至1.22～1.89 mg GAE/g；Yu等[4，5]研究表明超声诱导能够显著提高发芽花生的白藜芦醇含量，且随着超声频率和发芽时长而变化，此现象也与花生品种的差异密切相关。但是未见超声诱导对发芽花生苯丙烷类活性物质合成相关基因的表达调控的研究报道。

为了加强对发芽花生更深层次的开发利用，本研究选用超声诱导处理后的发芽花生作为研究对象进行转录组测序，并利用生物信息学方法对其测序数据进行功能注释分类、代谢通路分析，探讨超声诱导处理对发芽花生苯丙烷类合成相关基因的表达的作用。

1 材料与方法

1.1 试验材料及超声处理

材料于2020年6月采自辽宁省农业科学院沙地治理研究所花生研究基地，所选品种为阜花23号。选取大小均一且无破损的200粒花生籽粒放在纱布上，用体积分数为 1%的次氯酸钠（次氯酸钠:水=1:100）浸泡消毒20 min，再用纯净水冲洗三次，避光浸泡6 h后，连同纱布放在烧杯里，烧杯加水没过花生，放入超声仪器中进行超声处理。超声参数设置为 35 ℃、30 min、240 W（60%）、85 kHz。处理后，种子放置在托盘上，底部铺上滤纸，铺上纱布然后另一半纱布盖上进行避光发芽。经过72 h后，取空白对照组和超声诱导下的花生胚芽置于液氮中快速冷冻，于-80 ℃冰箱中保存备用，送至上海美吉生物有限公司。

1.2 总RNA提取、cDNA文库构建和测序

从组织样品中提取total RNA，利用Nanodrop 2000对所提RNA的浓度和纯度进行检测，琼脂糖凝胶电泳检测RNA完整性，Agilent 2100测定RIN值。样品检测合格后利用带有 Oligo（dT）的磁珠与 ployA 进行A-T碱基配对，加入fragmentation buffer，可以将mRNA随机断裂，通过磁珠筛选分离出300 bp左右的小片段。在逆转录酶的作用下，加入六碱基随机引物（random hexamers），以mRNA为模板反转合成一链cDNA，随后进行二链合成，形成稳定的双链结构。然后加入End Repair Mix将其补成平末端，随后在3’末端加上一“A”碱基，用于连接Y字形的接头。最后通过PCR扩增，纯化PCR产物，得到最终文库。文库质量检测合格后，不同文库按照目标下机数据量进行 pooling，利用Illumina平台进行测序（PE文库，读长2×150 bp）。

1.3 测序数据质控

使用软件SeqPrep去除reads中的接头序列，去除由于接头自连等原因导致没有插入片段的reads；将序列末端（3’端）低质量（质量值小于30）的碱基修剪掉，如剩余序列中仍然有质量值小于10的碱基，则将整条序列剔除，否则保留；去除含N比率超过10%的reads；舍弃去adapter及质量修剪后长度小于50 bp的序列得到clean data，使用TopHat2（http://tophat.cbcb.umd.edu/）软件进行序列比对分析。

1.4 差异表达基因的筛选

以样品KB组作为对照样品，与样品CS进行基因表达量的比较。将错误发现率（false discovery rate，FDR）≤0.001且差异倍数（fold change）≥2作为DEG的筛选标准，满足此筛选标准的基因即为DEG。

1.5 DEG的功能注释

通过BLAST软件将筛选到的所有DEG分别与NCBI非冗余蛋白数据库（Non-redundant Protein Sequence Database，NR）、京都基因与基因组百科全书（Kyoto Encyclopedia of Genesand Genomes，KEGG）数据库、基因本体论（GeneOntology，GO）、蛋白质直系同源簇数据库（COG）、蛋白质序列数据库（Swiss-Prot）、蛋白数据库（TrEMBL）和非冗余蛋白序列库（NR）等一系列数据库，并对基因进行功能注释、分析以及统计。

2 结果与讨论

2.1 发芽花生转录组测序数据统计分析

表1 发芽花生转录组测序统计Table 1 Transcriptomic sequencing statistics of peanut sprout

将构建好的发芽花生转录组文库通过HiSeq 2000高通量测序，得到超声处理 CS组和空白对照组 KB的三次生物学重复间的相关系数均大于0.9，ReadSum的平均值分别为46208369条和49367041条，数据过滤后超声处理CS组和空白对照组KB的平均总碱基数分别是6.98 Gb和7.45 Gb，GC含量分别为45.00%和45.13%，这些序列的Q20和Q30分别在98%和94%以上。以上数据说明全部样品转录组测序质量较高，可为后续的差异表达基因分析提供了科学的依据。

2.2 发芽花生转录组功能注释

将组装获得的所有基因和转录本与 NR、Swiss-Prot、Pfam、EggNOG、GO和KEGG六个常用数据库进行序列比对（图1），最终获得有注释信息的表达基因数目为 48660，各个数据库注释到的表达基因有较大的区别，其中NR数据库注释到47998个表达基因，占比98.64%，GO和COG分别有27880个和46411个表达基因，占比57.30%和95.38%，得到注释最少的数据库为KEGG，仅有22200个表达基因，占比为 45.62%。长度位于 300～1000 bp和长度≥1000 bp的数目分别有18310个和19233个（表2）。

表2 基因功能注释统计Table 2 Functional annotation of ref genes

2.3 花生芽响应超声诱导DEGs数目的分析

为研究发芽花生响应超声诱导的基因表达情况，对超声处理后的差异表达基因进行分析（图2）。两组样品中共筛选出1104个DEGs，其中上调DEGs共有583个，用红色点表示；下调DEGs共有521个，用绿色点表示；灰色点代表非差异表达基因。结果表明，上调DEGs的数目低于下调DEGs的数目。

2.4 DEGs的蛋白聚类分析

将DEGs和COG数据库进行比对，并对DEGs的功能进行预测和统计归类，结果显示，共有 14497条DEGs与数据库中的基因相似度高。根据功能分类，可初步将花生种子转录组中的差异基因分为21类（图3），数据分析发现，DEGs的COG功能涉及大多数的生命活动，其中种类全面，数量最多的是转录类的基因，其次是碳水化合物的运输和代谢、重复和修复及翻译后修饰、蛋白代谢、信号转导机制、氨基酸转运和代谢等类别；辅酶的运输和代谢类数量最少，只有4个；与油脂转运及代谢相关的基因有654个；此外仍有643个DEGs未知功能类别（图3）。

2.5 基因功能分类分析

对超声诱导下发芽花生的DEGs进行了GO富集分析，共分为生物进程（biological process）、分子功能（molecular functions）和细胞组分（cellular component）三个大类。在DEGs富集数目的共33个亚类中（图4），生物进程注释到13个分支，DEGs主要集中在细胞组成或生物发生和细胞进程，分别有292个和241个基因功能得到注释，说明细胞组分类的相关基因在超声诱导花生发芽的生物进程中起着非常重要的作用。富集在细胞组分的DEGs主要分为10个分支，其中细胞组分包含的基因数目最多，共有217个，其次是膜的组成部分。在分子功能的9个分支中，参与催化活性和涉及结合功能的基因数量远超过参与其他生物过程的基因数，分别有338个和245个基因。

2.6 DEGs的代谢通路富集分析

为了进一步探索DEGs的生物学功能，解析花生芽在超声诱导下的代谢调控网络，对 DEGs进行了KEGG富集分析，共有327个DEGs富集到了97条通路中。在富集显著的前20条代谢通路中（图5），包含苯丙烷生物合成（24）、植物激素信号转导（15）、胞吞作用（12）、植物与病原体的相互作用（11）、黄酮类生物合成（10）、淀粉和蔗糖代谢（8）、核糖体（8）、剪接体（7）、亚油酸代谢（6）、磷酸肌醇代谢（6）、谷胱甘肽代谢（6）、甘油磷脂代谢（6）、错配修复（6）、同源重组（6）、基因复制（6）、内质网中的蛋白质加工（5）、类胡萝卜素的生物合成（5）、氨酰基tRNA的生物合成（5）、抗坏血酸和藻酸盐代谢（5）、α-亚麻酸代谢（5）、精氨酸和脯氨酸代谢（5）、氨基糖和核苷酸糖代谢（5）、糖酵解/糖异生（5）、RNA转运（5）、核苷酸切除修复（5）等。富集分析结果说明花生芽通过代谢、合成次生代谢物等途径参与超声诱导的响应。

2.7 花生芽在超声诱导下的转录因子分析

转录因子在调控植物非生物胁迫抗性方面起着重要作用[11-13]，本研究分析超声诱导下花生芽的DEGs，共属于48个转录因子家族。其中B3家族包含293个基因，其次是MYB家族288个、bHLH家族281个，MYB_related家族261个以及ERF家族210个，说明这些家族转录因子的表达在发芽花生受超声诱导下调控效果显著。

本研究发现B3转录因子是发芽花生在超声处理下表达量变化最多的家族。研究表明已知MYB、B3、WRKY、bHLH、bZIP等转录因子表达量变化能提高植物适应逆境的能力[14]。其中B3家族是一类含有B3功能域的转录因子家族，在植物的生长发育过程中发挥重要的调控作用，其功能域是一种可以和DNA特异结合且高度保守的结构域[15]。B3家族转录因子除了有参与DNA结合B3结构域外，还有其他保守的结构域，包括AP2，生长素应答因子，生长素/吲哚-3-乙酸等。bHLH转录因子家族是植物中最大的基因家族之一，广泛的参与了植物代谢、发育及对各种胁迫的应答机制[16]。MYB转录因子是植物中最大的转录因子家族之一，广泛参与植物生理生化过程，包括植物表皮组织细胞分化、外界环境因素响应、激素应答等[17]。研究表明，MYB转录因子的功能多种多样，对植物的生长发育至关重要。例如，MYB转录因子可通过控制各个细胞分裂时期来实现对细胞周期的调控，也可通过调节植物次生代谢相关基因的表达来调节花青素、类黄酮等次生代谢产物的合成，在非生物胁迫和生物胁迫响应中也扮演了重要的角色[15]。此外，MYB转录因子还能调控植株对植物激素的响应[16]。

2.8 花生芽在超声诱导下DEGs的苯丙烷类生物合成途径

苯丙烷类生物合成途径为DEGs富集数目最多的代谢通路，分析该调控途径对研究玉米响应盐胁迫的分子机制具有重要意义（图7）。发现共有21 DEGs得到富集，其中过氧化物酶（E1.11.1.7）注释到10个，下调3个，上调7个；β-葡糖苷酶上调（EC3.2.1.21）1个；四氢大麻酸合酶上调1个；甘露醇脱氢酶上调2个；网状氧化酶样蛋白上调1个；4-香豆酸酯-CoA连接酶上调2个；亚精胺羟肉桂酰基转移酶上调1个；氰基β-葡萄糖苷酶上调1个，下调1个；咖啡酰-CoAO-甲基转移酶上调2个；肉桂酰基-CoA还原酶上调2个。

前人研究报道白藜芦醇在植物体内是以苯丙氨酸为底物，通过苯丙烷途径生成[18]，因此本研究所得的所有转录本与KEGG通路数据库中的苯丙烷合成途径（map00940）进行比对，结果发现花生芽中arahy.Tifrunner.gnm1.ann1.DXZI51、arahy.Tifrunner.gnm1.ann1.VGN2GE、arahy.Tifrunner.gnm1.ann1.Y23DM6这三个基因显著上调与白藜芦醇合成途径相关，其中arahy.Tifrunner.gnm1.ann1.DXZI51为肉桂酰基辅酶A还原酶（Cinnamoyl-CoA reductase 2）参与了肉桂酸的合成。苯丙氨酸（PHE，Phenylalanine）首先经过苯丙氨酸解氨酶（PAL，phenylalanine ammonia lyase）脱氨基作用失去氨基之后转化为肉桂酸（Cinnamic acid），肉桂酸再通过肉桂酸羟化酶（C4H，cinnamate-4-hydroxylase）发生羟基化磷脂反应得到对香豆酸（Cinnamoyl-COA），在本通路结果中arahy.Tifrunner.gnm1.ann1.VGN2GE、arahy.Tifrunner.gnm1.ann1.Y23DM6基因均为4-香豆酸酯-CoA连接酶（4-coumarate-CoA ligase），参与了由4-香豆酸辅酶A连接酶（4CL，4-coenzyme A ligase）得到对香豆辅酶A的生物合成。前人报道香豆辅酶A（4-coumaroyl-COA）最后在外界因素刺激下，开启白藜芦醇生成途径，由1分子的对香豆酰辅酶A和3分子的丙二酰辅酶A（Malonyl-COA）在芪合酶（STS，stilbene synthase）的催化作用下缩合生成1分子的反式白藜芦醇[19]。由此可以说明在花生芽在超声诱导下产生机械损伤，通过调控arahy.Tifrunner.gnm1.ann1.DXZI51、arahy.Tifrunner.gnm1.ann1.VGN2GE、arahy.Tifrunner.gnm1.ann1.Y23DM6促进白藜芦醇的富集。

本研究发现苯丙烷类生物合成途径是富集数目最显著的一条代谢通路，共有21个相关基因得到注释。在植物次生代谢途径中，所有黄酮类物质、木质素和白藜芦醇等大部分酚类物质均经苯丙氨酸代谢途径合成，研究已经表明，植物体内这些关键酶在应对外界生物及非生物逆境时起到至关重要的调控次生代谢作用[20-23]。例如，PAL活性的提高增强了大豆根部木质素的合成进；甘草黄酮类物质的提高与PAL、C4H基因及其酶活性的增加有着不可分割的联系[24，25]。上述结果表明，超声诱导提高了发芽花生苯丙烷类合成途径中基因的表达，这意味着发芽花生苯丙烷代谢途径中的下游产物也将可能被超声诱导所改变，这也进一步验证了本研究在前期试验中发现超声诱导能够显著提高发芽花生的白藜芦醇含量的结果[4-6]。殷向静[26]利用超声诱导葡萄富集白藜芦醇，发现白藜芦醇的富集受到白藜芦醇合酶的转录控制，而超声处理能够引起白藜芦醇合酶活性升高。此外，李淑莹[27]的研究表明，超声联合苯丙氨酸诱导可对花生芽中白藜芦醇的富集产生协同效应，原因是由于超声可以提高肉桂酸-4-羟基化酶、4-香豆酸-辅酶A连接酶、类黄酮3-O-葡萄糖基转移酶4种酶的活性。Ling等[28]研究超声处理（400 W，6 min）与0.4%过乙酸在20 ℃对枇杷贮藏过程中生理变化的影响，结果表明，联合超声和乙酸处理后，枇杷果实中总黄酮含量在第3 d达到最高值，且高于对照组。这可能是因为超声联合过乙酸在处理可以提高枇杷果实内多酚氧化酶和过氧化物酶的活性，而这两种酶则是参与苯丙烷途径和氧化过程的重要酶，从而产生多种具有结构和防御功能的酚类化合物[2]。此外，卞紫秀等[23]研究表明，利用超声处理可以有效地促进苦荞麦种子的萌发，从而富集苦荞芽苗中黄酮类化合物。这可能要归因于超声处理能有效激活植物种子萌发期的各种酶类的活性[24]，显著提高种子萌发率，同时诱导种子中一些生物活性成分的合成，这与本研究发现的结果相一致，说明超声诱导能够调控发芽花生苯丙烷类生物合成代谢途径的基因表达情况。

3 结论

为了揭示发芽花生如何应答超声外源场诱导，筛选在超声诱导后起重要调节作用的基因，本研究在前期研究基础上，选出适宜富集白藜芦醇的花生品种阜花23号为试验原料，是较有代表性的一个研究材料，为此次研究奠定了良好的材料基础。研究了经超声诱导后的发芽花生（CS）和未经诱导处理的发芽花生（KB）的转录组差异，证实了超声能够影响发芽花生苯丙烷类物质的生物合成，从两组样品中共筛选出1104个差异基因，其中上调583个，下调521个。通过与七大数据库的比对分析，解析了发芽花生在超声诱导下所涉及到的功能分类、代谢通路及生物进程的调控，功能注释表明差异基因主要富集在遗传信息处理、代谢途径和蛋白质转换。并挖掘到 3个（arahy.Tifrunner.gnm1.ann1.DXZI51、arahy.Tifrunner.gnm1.ann1.VGN2GE、arahy.Tifrunner.gnm1.ann1.Y23DM6）参与发芽花生苯丙烷类物质的显著基因，这些研究结果为揭示发芽花生苯丙烷类物质的生物合成途径提供一定的参考依据。但众多基因复杂的调控网络及发芽过程中酚酸、木质素及总黄酮等物质在超声诱导处理下的代谢机理需要进一步深入探讨，目的在于将来通过外援手段改变此类基因的调节功能，让发芽花生的营养和功能性成分的变化向对人类有利的发现发展。