张 敏,王 梅,2,3,王红霞,2,3,王俊玲,2,3,穆国俊,2,3
(1.河北农业大学 华北作物种质资源教育部重点实验室/河北种质资源实验室,河北 保定 071001;2.河北省藜麦产业技术研究院,河北 张家口 075000;3.中国农村专业技术协会河北张北藜麦科技小院,河北 张家口 075000)
花青素属于类黄酮物质,广泛存在于植物的叶片、果实等组织器官中,是一种可溶于水、甲醇等溶剂的强抗氧化剂[1]。Liu Y 等[2]对紫色大白菜的结构基因进行研究表明,PAL(苯丙氨酸解氨酶,Phenylalanine ammonia-lyase)、C4H(肉桂酸羟化酶,Cinnamate 4-hydroxylase)、4CL(4-香豆酸辅酶A 连接酶,4-coumarate COA ligase)、CHS(查尔酮合成酶,Chalcone synthase)、CHI(查尔酮异构酶,Chalcone isomerase)、F3H(黄烷酮3-羟化酶,Flavanone 3-hydroxylase)、F3'H(类黄酮3'-羟化酶,Flavonoid 3 '-hydroxylase)、FLS(黄酮醇合酶,Flavonol synthase)、DFR(氢黄酮醇还原酶,Dihydroflavonol reductas)和ANS(花青素合成酶,Anthocyanidin synthase)在子叶期、苗期和大叶期均上调表达。向奕等[3]对茶树紫芽中花青素合成相关基因差异表达的分析结果表明,PAL、C4H、CHS、CHI、F3H、F3'H、F3'-5'H(类黄酮3'-5' 羟化 酶,Flavonoid 3'-5' hydroxylase)、DFR 和ANS关键酶基因均为上调表达。上述调控花青素生物合成的基因均为结构基因,此外还存在一些调控花青素合成途径的结构基因,以转录因子(Transcription factor,TF)的形式激活结构基因的表达,间接调控花青素的合成。Nesi N[4]以及Baudry A 等[5]的研究结果均表明,MYB(v-Myb avian myeloblastosis viral oncogene homolog)转录因子TESTA 2(TT2)通过与因子TT8(bHLH42,Basic helix-loop-helix)和TTG1(WD40)结合,调控晚期生物合成基因的表达,从而使种皮着色。在矮牵牛花瓣细胞中,花青素的生物合成也受MYB-bHLH-WD40 复合体的调控,但一些MYB 如AtMYB12 和VvMYBF1 可以在没有bHLH 和WD40 辅助因子的情况下调节类黄酮的生物合成[6-7]。此外bZIP 转录因子家族也可促进花青素含量的积累。An JP 等[8]研究发现bZIP 转录因子MdbZIP44 可调控脱落酸,促进苹果花青素积累。Liu CC 等[9]发现bZIP 转录因子HY5 介导CRY1a诱导番茄花青素的生物合成。
藜麦(Chenopodium quinoaWilld.)是一年生苋科藜属双子叶植物,藜麦可作为粮食作物,保健食品[10],因其叶片可呈现粉、绿、黄、紫、红、橙等多种色彩所以具有较高的观赏价值。目前为止,藜麦在种子[11]方面研究较多,关于藜麦叶色分子调控机理的相关研究在国内外尚未见报道。本研究选取在河北省张家口种植的翠绿色和粉红色藜麦叶片的种质材料,通过转录组-代谢组联合分析探索藜麦不同色泽叶片的花青素生物合成的差异表达基因(Differentially expressed genes,DEGs)及差异代谢物(Differential accumulated metabolites,DAMs),以期为深入探讨藜麦叶片花青素生物合成的分子机理提供重要参考。
供试藜麦材料M202(翠绿色)和M146(粉红色)由河北省科技创新服务中心提供。试验于2021 年5月25 号和2022 年5 月26 号将M202 和M146 覆膜种植于河北农业大学张北实验站(海拔1 440.6 m,东经114°50'13'',北纬41°10'51''),行距40 cm、株距18 cm。分别在生长到70 d 和110 d 取M202和M146 的叶片,每个样品量约1 ~3 g,用塑封袋密封,并加入硅胶5 ~10 粒,-80 ℃保存,3 次生物学重复。
采用色差仪(CR-10plus,日本)测定M202和M146 的叶色,3 次重复。L 表示黑白值,L 值为+表示偏白,L 值为-表示偏暗;a 表示红绿值,a 值为+表示偏红,a 值为-表示偏绿;b 表示黄蓝值,b 值为+时表示偏黄,b 值为-时表示偏蓝。并以M202 和M146 的70 d 为标准,利用SPSS26进 行 差 异 计 算(ΔE=[(ΔL)2+(Δa)2+(Δb)2]1/2)。当ΔE≤2 时,人眼几乎无法辨别颜色差异;当2<ΔE≤6 时,人眼可察觉颜色差异,但不是很明显;当ΔE>6 时,人眼很容易辨别颜色差异。
花青素含量的测定参考Serrano[12]的测定方法。将M202 70 d(N1)、M202 110 d(N2)、M146 70 d(F1)和M146 110 d(F2)4 份叶片样品,每份样品取0.5 g 在0.1 mol/L 的盐酸甲醇溶液中研磨至粉末状,用锡箔纸包好,4 ℃避光过夜,待样品组织变白,10 000 r /min 离 心 机 离 心10 min,在530 和657 nm 处测定上清液的吸光度。利用公式A530-0.25×A657计算花青素在530 nm 处的吸光度差值,再用吸光度差值/鲜重表示总的花青素含量。
参 照Rio DC[13]的Trizol 沉 淀 法 提 取 样 品RNA。选用RIN ≈ 10;28S/18S ≥1.5;1.7 < OD260/OD280< 2.0 的样品以保障使用合格的样品进行转录组测序,转录组测序由广州基迪奥生物有限公司提供技术支持,3 次生物学重复。将|log2(FC)|≥1 和错误发现率(False discovery rate,FDR)< 0.05 作为DEGs 筛选标准,其中FC(差异倍数,Fold change)表示基因相对表达水平的差异。当log2(FC)>1 时DEGs 为 上 调,当log2(FC)<1时DEGs 为下调。利用GO 数据库(http://www.geneontology.org/)和KEGG 数 据 库(Kyoto Encyclopedia of Genes and Genomes)对差异基因做GO 和KEGG 富集分析。DEGs 的GO 和KEGG显著富集通路分析以q_value (Benjamini- Hochbergcorrected q_value)<0.05 为筛选标准。
利用液相色谱串联质谱(LC-MS/MS)对代谢物进行定性和定量分析。根据质谱公共数据库HMDB(http://www.hmdb.ca/)、Massbank(http://www.massbank.jp/)、LipidMaps(http://www.lipidmaps.org)、METLIN(http://metlin.scripps.edu/index.php)等公共质谱数据库确定代谢物的注释。变量对模型的重要性(Variable importance in the projection,VIP)描 述 了 每 一 个 变 量 对模型的总体贡献,通常设定阈值为VIP >1。以|log2(FC)|≥1,P< 0.05 为标准鉴别DAMs。
将DEGs 和DAMs 注释到KEGG 数据库中,以获得它们共同的KEGG 路径。本研究对从转录组和代谢组中筛选出的DEGs 和DAMs 计算Pearson 相关系数(PCC)。以PCC 阈值≥0.6,P<0.05 为标准,对代谢组和转录组数据进行联合分析。PCC ≥0.6表示DEGs 和DAMs 之间存在显著相关关系。利用Cytoscape 软件对DEGs 的表达量和DAMs 丰度的pearson(|r|> 0.8)相关性系数,绘制相关网络图。
对筛选出的差异表达转录因子进行分析并进一步分析出与藜麦花青素生物合成相关的TF 家族及相关基因。根据转录组数据中各转录因子基因的FPKM 值,分析差异表达转录因子的表达模式。
提取样品RNA,经反转录后合成相应的cDNA进行qRT-PCR 验证。反应程序为:95 ℃变性10 s,60 ℃退火20 s,72 ℃延伸15 s,40 个循环。以编码ACT7 基因作为内参基因,依据目的基因所需序列库对应的CDS 序列,用Premier5.0 software 设计引物,采用2-△△Ct法对qRT-PCR 结果进行相对表达量计算,采用t-test 检验进行相对表达量数据显著性分析,3 次生物学重复。
M202 和M146 在70 d 和110 d 时叶片颜色发生了明显的差异变化(图1A、图1B)。花青素总含量测定结果表明,随着藜麦的生长,藜麦叶片的花青素总含量逐渐升高。N1 时期叶片的花青素总含量是1.57 OG/g,N2 时期花青素总含量是2.11 OG/g;F1 时期的叶片花青素总含量是1.52 OG/g,F2 时期的花青素总含量是2.50 OG/g(图1C)。
图1 两个藜麦品种叶片颜色、吸光值和总花青素含量Fig.1 Leaf color, anthocyanin content and light absorption value of two quinoa varieties
2.2.1 差异表达基因筛选 DEGs 以|log2(FC)| ≥1及FDR<0.05 为筛选依据,结果表明,N1_vs_N2、F1_vs_F2、N1_vs_F1 和N2_vs_F2 4 个比较组差异基因总数分别为12 563、9 479、8 475 和 7 519,其中上调基因分别为7 620、5 054、3 591 和3 480,下调基因分别为4 943、4 425、4 884 和4 039(图2A)。基因重叠分析结果表明,N1_vs_N2、F1_vs_F2、N1_vs_F1 和N2_vs_F2 4 个比较组中独有的基因分别有3 339、2 026、1 896 和1 373 个,共有的差异基因为1 089 个(图2B)。
图2 4 个比较组中的差异表达基因的MA 图和Venn 图Fig.2 MA and Venn diagrams of differentially expressed genes in the 4 comparison groups
2.2.2 GO 和KEGG 富集分析 GO 富集分析结果表 明,DEGs 在N1_vs_N2、F1_vs_F2、N1_vs_F1和N2_vs_F2 中富集到的GO Terms 分别为2 107、1 857、1 815 和1 705,其中生物过程(Biological process,BP)涉及1 397、1 251、1 235 和1 136,分 子 功 能(Molecular function,MF)涉 及504、431、410 和403,细胞组分(Cellular component,CC)涉及206、175、170 和166 条。注释分析结果表明,7 条GO Terms 与花青素生物合成密切相关。被注释到BP 中,与花青素生物合成密切相关的GO Terms 为GO:0009812(类黄酮代谢过程,Flavonoid metabolic process)、GO:0010468(基因表达调控,Regulation of gene expression)、GO:0051553(黄酮生物合成过程,Flavone biosynthetic process)、GO:0009813(类黄酮生物合成过程,Flavonoid biosynthetic process)和GO:0009411(应对紫外线,response to UV);被注释到MF 中,与花青素生物合成密切相关的GO Terms 为GO:0043169(阳离子绑定,cation binding)和GO:0016703(氧化还原酶活动,Oxidoreductase activity)(图3B)。
图3 各比较组富集到的GO terms 数量统计图和藜麦叶片GO 差异表达基因富集分析弦图Fig.3 Statistical diagram of GO terms enriched by each comparison group and chordal gragh of differentially expressed genes enriched on GO analysis of leaves in quinoa
N1_vs_N2、F1_vs_F2、N1_vs_F1 和N2_vs_F2 的KEGG 富集分析结果表明,存在富集显著差异的代谢通路条数分别为133、132、128 和131,与花青素生物合成显著相关的代谢途径6 条,分别为苯丙氨酸代谢(Phenylalanine metabolism,ko00360)、苯丙素的生物合成(Phenylpropanoid biosynthesis,ko00940)、类黄酮生物合成(Flavonoid biosynthesis,ko00941)、植物激素信号转导(Plant hormone signal transduction,ko04075)、黄酮和黄酮醇生物合成(Flavone and flavonol biosynthesis,ko00944)以及昼夜节律-植物(Circadian rhythm-plant,ko04712),其中苯丙素的生物合成(Phenylpropanoid biosynthesis,ko00940)富集得到的基因最多,分别为120、114、74 和111 个;类黄酮生物合成(Flavonoid biosynthesis,ko00941)的富集因子最高,分别为0.222、0.222、0.167 和0.444(图4)。
图4 KEGG 差异表达基因富集分析气泡图Fig.4 Bubble diagram of differentially expressed genes based on KEGG enrichment analysis
代谢组学分析结果表明,各比较组中与花青素生物合成相关的DAMs 个数分别为6、6、6 和5 个(图5)。对DAMs 的log2(FC)值进行分析发现:在N1_vs_N2 比 较 组 的6 种DAMs 中,N2 时 期 的天竺葵色素3-O-(6-咖啡酰基-β-D-葡萄糖苷)(Pelargonidin 3-O-(6-caffeoyl-beta-D-glucoside))、飞燕草素3-(6-对香豆酰基)葡萄糖苷[Delphinidin 3-(6-p-coumaroyl) glucoside] 和矢车菊素3- 芸香糖苷5-葡萄糖苷(Cyanidin 3-rutinoside 5-glucoside)含量均高于N1 时期;在F1_vs_F2 比较组的6 种DAMs 中,F2 时期的天竺葵色素3-O-(6- 咖啡酰基-β-D-葡萄糖苷)、飞燕草素3-(6-对香豆酰基)葡萄糖苷、矢车菊素3-芸香糖苷5-葡萄糖苷和天竺葵素3-葡萄糖苷5-咖啡基葡萄糖苷(Pelargonidin 3-glucoside 5-caffeoyl glucoside)含量均高于F1 时期;在N1_vs_F1 比 较 组 的6 种DAMs 中,F1 时期的无色矢车菊素(Leucocyanidin)、矢车菊素3-O- 芸香糖苷(Cyanidin 3-O-rutinoside)和矢车菊素3-O-3",6"-O- 二丙二基葡萄糖苷(Cyanidin 3-O-3",6"-O-dimalonyl glucoside)含量均高于N1时期;在N2_vs_F2 比较组的5 种DAMs 中,F2 时期的天竺葵色素3-O-(6-咖啡酰基-β-D-葡萄糖苷)、飞燕草素3-(6-对香豆酰基)葡萄糖苷、矢车菊素3-O-芸香糖苷、矢车菊素3-O-3",6"-O-二丙二基葡萄糖苷和柚皮素(Naringenin)含量均高于N2 时期。
图5 不同比较组的DAMs 表达调控Fig.5 Expression regulation of DAMs in different comparison groups
转录组和代谢组学联合分析结果表明,4 个比较组的DEGs 与DAMs 在苯丙素的生物合成(Phenylpropanoid biosynthesis,ko00940)、类黄酮生物合成途径(Flavonoid biosynthesis,ko00941)和花青素生物合成途径(Anthocyanin biosynthesis,ko00942)富集通路中存在相关关系(图6)。DEGs 的表达量和DAMs 丰度的相关网络图(图9)显示,在N1_vs_N2 比较组存在15 个正相关2 个负相关;F1_vs_F2 比较组存在12 个正相关4 个负相关;N1_vs_F1 比较组存在8 个正相关2 个负相关;N2_vs_F2 比较组存在11 个正相关1个负相关。
图6 与花青素生物合成相关的差异表达基因和差异代谢物的相关性网络图Fig.6 Correlation network diagram of differentially expressed genes and differentially expressed metabolites related to anthocyanin biosynthesis
在转录因子分析中,2 645 个DEGs 被鉴定,这些DEGs 属于55 个转录因子家族,其中定位到与花青素调控相关的转录因子家族有bHLH、bZIP 和WRKY。bHLH 和bZIP 家族是花青素转录调控的重要转录因子家族,分别富集到201 和91 个DEGs;WRKY 家族富集到100 个DEGs(图7)。
图7 前20 个转录因子家族数目统计图Fig.7 The number of the top 20 transcription factor families
富集到的相关转录因子主要为bHLH 家族的MYC2、bHLH14;bZIP 家族的HY5、TGA;WRKY家族的WRKY22、WRKY24 和WRKY46,富集于ko04712、ko04075、MAPK 信号通路-植物(MAPK Signaling pathway-plant,ko04016)和植物病原菌互作(Plant-pathogen interaction,ko04626)4 条 代 谢途径中。这些转录因子在藜麦叶片生长的3 个时期都有不同程度的表达(表1)。
表1 转录因子在各样品中的FPKM 值Table 1 FPKM values of transcription factors in each sample
利用皮尔逊系数法测定筛选出的20 个差异表达转录因子与12个DEGs间的相关性并绘制热图(图8)。
图8 花青素生物合成途径以及TF 家族基因与DEGs 之间相关性热图Fig.8 heat map of correlation between TF family genes and DEGs
结果发现:UGT79B1 与WRKY22 呈显著正相关;ANR、CHI、F3H 和bHLH14 呈显著正相关;CHS与TGA 呈显著正相关;PAL 与bHLH14 和WRKY24呈显著负相关;CYP73A 和WRKY24 呈显著负相关;HCT、C3'H、FG3 与WRKY24 呈显著负相关;C3'H与MYC2 和WRKY24 呈显著负相关。转录因子调控机制较为复杂,其具体的调控机制在接下来的工作中我们会进一步进行研究,并主要对bHLH 和WRKY家族的作用模式进行深入发掘。
在富集到的6 条通路中选择10 个DEGs 进行qRT-PCR 验 证,在N1_vs_N2 中 有PAL、CHI、CYP75B1、F3H、FG3和CYP73A6 个基因被验证;在F1_vs_F2 中 有PAL、CHS、4CL、FG3和C3’H5 个基因被验证;在N1_vs_F1 中有PAL和4CL2个基因被验证以及在N2_vs_F2 中有PAL、4CL、CHI、F3H、FG3、HCT、C3'H和CYP73A8 个基因被验证,验证结果与转录组结果趋势一致(图9)。
图9 10 个关键基因的qRT-PCR 检测Fig.9 10 key genes were detected by qRT-PCR
通过表型观察和花青素含量的变化进行分析,各时期粉红色藜麦叶片比翠绿色藜麦叶片的花青素含量都高,在红肉猕猴桃和绿肉猕猴桃中[14]紫色和绿色紫苏叶[15]中都得到了相同结果,表明表型颜色与花青素含量有着密切关系。对色差值数据进行分析发现N1_vs_N2 和F1_vs_F2 的△E 大于6,肉眼可明显观察到颜色的变化。根据花青素含量的测定结果发现随着藜麦叶片变色之后的不断生长,到35 d 花青素的含量呈现上升趋势,颜色稍微加深。花青素含量变化与色差值结果一致。
4CL 具有高度趋异的底物偏好性及底物特异性。植物中CHI是黄酮类物质合成途径的第二个关键酶,上调CHI 的表达可以增加植物体内黄酮类次生代谢物质的合成。过表达的CHI基因的甘草毛状根中黄酮类化合物的含量明显提高[16]。红花CHI基因可以参与拟南芥中黄酮的合成代谢,提高黄酮类含量[17]。在本研究中4CL的表达量在各时期均呈现上升趋势,与总花青素含量测定结果变化趋势一致,表明4CL与花青素含量呈现正相关。CHI的表达量在各时期呈现下降趋势,TGA与CHI呈显著正相关关系,且表达量趋势一致,推测在藜麦叶片变色生长过程中CHI受其TGA的调控,共同参与花青素的合成。
CHS 是合成黄酮类化合物的关键酶,也是其限速酶。Hinderer 等[18]在研究胡萝卜CHS 活性时表明胡萝卜CHS受到柚皮素和柚皮素查尔酮的抑制。Fukusaki 等[19]利用RNAi 技术,使蓝色夏堇(Torenia lybrida)中的CHS基因沉默后,获得了开白色和灰白色花的转基因植株。此外成功的例子在非洲菊(Gerbera hybrida)[20]中也被报道。其现象表明CHS的表达与花青素含量呈正相关,CHS表达量高其花青素含量越高。在本研究中CHS的表达量在各时期均表现为上调,在N2 和F2 时期CHS的表达量显著上升,其变化趋势与总花青素含量变化一致。推测CHS参与花青素含量的累计,与花青素含量呈正相关,与前人研究结果一致。
本研究中在M202 品种中PAL的表达量差异逐渐下降,M146 品种中PAL的表达量差异逐渐上升,且M146 的表达量显著高于M202,推测PAL的高表达可能是粉红色藜麦叶片的颜色形成的关键。李云萍等[21]研究表明紫心甘薯的花青素含量与PAL的活性呈正相关;燕雪芬[22]通过多克隆抗体的制备也表明PAL的活性与黄酮的积累呈正相关性。本研究中M146 的PAL表达量高于M202 的PAL表达量,且粉红色藜麦叶片的花青素含量高于翠绿色藜麦叶片,与花青素含量变化趋势一致,研究结果与前人研究结果一致,表明PAL的花表达量与花青素含量呈正相关。
F3H是使类黄酮羟基化的重要催化物是花青素生物合成过程中的重要途径。本研究中,通过转录组和qRT-PCR 共同验证表明在同一品种不同时期的藜麦叶片中F3H的表达量逐渐下降,对苜蓿和康乃馨进行研究表明[23-24],F3H的低表达会导致颜色更深,此结果与藜麦叶片色差值和花青素含量测定结果一致,F3H 的低表达可促进藜麦叶片花青素含量的积累。
本研究发现羟基肉桂酰基转移酶(Hydroxy cinnamoyl transferase,HCT)参与调控藜麦叶片花青素的生物合成,催化4-香豆酰辅酶a 形成咖啡酰辅酶a。Liu Y[25]通过同时改变花青素、叶绿素和类胡萝卜素含量使粉红叶观赏羽衣甘蓝产生绿色斑纹揭示了花青素生物合成的分支途径。该途径是通过对4-香豆酰辅酶a 经过HCT作用生成咖啡酰辅酶a,然后通过CHS、F3H、DFR 和ANS等酶的催化作用合成花青素,其研究的花青素合成过程与本研究发现的合成过程一致。植物中HCT的沉默抑制了木质素的生物合成,导致CHS的激活和几种黄酮醇苷和酰基化花青素的积累[26]。在李佳伟等人[27]的研究中推测HCT具有底物偏好性。本研究的转录组数据表明HCT在M146 中的表达量显著高于M202,表明HCT可能促进粉红色藜麦叶片的合成,可能具有底物偏好性。
柚皮素是以柚皮素查耳酮(Naringenin chalcone)为前体经CHI作用后合成的,是花青素合成途径中重要的代谢产物。在本研究中柚皮素在N2_vs_F2 比较组中首次出现了上调差异表达,推测柚皮素是藜麦叶片变为粉红色的关键代谢物。有研究表明柚皮素是植物中一种重要的多酚类类黄酮物质,具有一定的药理作用[28],但其在植物颜色变化领域却鲜有报道。本研究中在N1_vs_F1 和N2_vs_F2 比较组中矢车菊素3-O-3",6"-O-二丙二基葡萄糖 苷(Cyanidin 3-O-3",6"-O-dimalonyl glucoside)出现了差异上调表达,表明矢车菊素3-O-3",6"-O-二丙二基葡萄糖苷是粉红色藜麦叶片形成的关键代谢物。Yan Yifan 等[29]研究表明蓝莓表面颜色的变化与花青素浓度的变化有关。在本研究中藜麦叶片中定位到的各代谢物含量差距较大,推测藜麦叶片颜色深浅的变化与各代谢物含量比例或花青素的总含量的有一定的关系。
综上所述,通过转录组和代谢组的联合分析,确定了藜麦叶片的花青素合成途径及相关基因及代谢物。1)代谢组学分析确定了在藜麦叶片中矢车菊素3-O-3",6"-O-二丙二基葡萄糖苷(Cyanidin 3-O-3",6"-O-dimalonyl glucoside)、柚 皮 素(Naringenin)是粉红色藜麦叶片形成的DAMs。2)PAL、4CL、CHS、CHI、HCT、CYP75B1、F3H和UGT79B1的差异表达对花青素积累起重要作 用。4CL基 因(ncbi_110730923)、F3H基 因(ncbi_110724781)、HCT基因(ncbi_110713661)可能是控制藜麦叶片呈现粉红色的关键结构基因。3)对转录组中富集到的转录因子分析表明bHLH14、WRKY24、TGA和MYC2调控某些关键结构基因共同参与藜麦叶片花青素的生物合成。本文综合分析了藜麦叶片花青素代谢产物及其相关基因的表达,确定了藜麦叶片中花青素的合成路径(图10),为揭示藜麦叶片花青素积累和颜色变化的调控机制提供参考。
图10 藜麦叶片的花青素生物合成途径Fig.10 Anthocyanin biosynthesis pathway of quinoa leaves