鸭茅bZIP 基因家族鉴定及非生物胁迫表达模式

马茜茜，张 欢，朱 杰，金献垚，陈雪松，陈 涛，邬 晓，王小珊，黄琳凯

（四川农业大学草业科技学院, 四川 成都 611130）

转录因子(transcription factor，TF)别称为反式作用因子，在生物体内，它是一类调节基因表达的蛋白，可与真核基因上游启动子区域中的特殊DNA序列(顺式作用元件)特异结合，从而激活或抑制相关靶基因的转录[1]。碱性亮氨酸拉链(basic leucine zipper，bZIP)转录因子是广泛分布于真核生物中的保守调节蛋白，在植物、动物、人、微生物中均有发现[2]。其由约20 个氨基酸的碱性结构域(basic region)与亮氨酸拉链二聚体结构域(leucine zipper)组成，前者包含N-x7-R/K 结构，可结合特异DNA 序列；后者参与寡聚化并与碱性结构域紧密结合，亮氨酸分布在每7 个氨基酸的第7 位的特定位置[3]。目前许多物种的bZIP 转录因子已被鉴定和分析，如拟南芥(Arabidopsis thaliana)中鉴定出75 个bZIP 转录因子[3]，水稻(Oryza sativa)中有89 个[4]，玉米(Zea mays)中有125 个[5]，高粱(Sorghum bicolor)中有92 个[6]，大豆(Glycine max)中有131 个[7]等。

鸭茅(Dactylis glomerata)又称果园草或鸡脚草，为禾本科鸭茅属，起源于欧洲、北非等地[8]，广泛分布于世界各地[9]。作为世界四大重要经济禾本科牧草之一，不仅叶量丰富、产量高，而且营养价值高、适口性好，还具有较强的耐阴、耐旱和越冬能力，可用于放牧、调制青贮饲料和干草等[10-11]。基于我国现有土地资源和现存生态环境问题，通常将牧草种植于边际土地，而边际土地容易受到贫瘠、干旱、涝害、寒冷等不良环境因素的影响[12]，致使牧草品质和产量降低。因此，牧草抗逆方面的研究至关重要。

植物在生长发育过程中会受到许多不利外界环境条件的影响，包括高温、干旱、盐害、寒冷等。bZIP转录因子基因家族具有强大的生物学功能，不仅参与植物生长发育过程，而且可以通过CACGTC (C 盒)、TACGTA (A 盒)、CACGTG (G 盒)等顺式作用元件与受脱落酸(abscisic acid，ABA)诱导的靶基因结合[13]，从而激活下游相关基因转录使植物更好地抵御不利外界环境条件。研究表明，过表达OsbZIP71的水稻植株可以提高其抵抗干旱和盐胁迫的能力[14]；在拟南芥中过表达小麦(Triticum aestivum)TabZIP60基因，植株更耐热、抗旱[15-16]；Gao 等[17]也发现，bZIP类转录因子可能与植物响应高温胁迫有关。但目前尚缺乏关于鸭茅bZIP基因家族的研究，因此对鸭茅bZIP基因家族进行鉴定和分析，同时进一步探究其在干旱和热胁迫下的表达模式，初步了解其在抗逆性方面的调控作用，为深入研究鸭茅bZIP基因的功能提供一定的理论和实践基础，对于其他植物的抗逆研究也具有一定的参考价值。

1 材料与方法

1.1 数据来源

拟南芥和鸭茅基因组数据分别从Ensembl Plants (https://plants.ensembl.org/index.html)下载和四川农业大学牧草课题组测序[18]获得。

1.2 鸭茅基因组中bZIP 基因家族成员的鉴定和染色体定位

使用HMMER 3.3.2 软件(http://www.hmmer.org/)，从Pfam (https://pfam.xfam.org/)数据库中查询bZIP保守域的Pfam 号(PF00170、PF07716)，并下载其隐马尔可夫模型，初步鉴定得到鸭茅bZIP基因，利用SMART 网站(http://smart.embl-heidelberg.de/)验证鉴定出的bZIP基因是否含有bZIP 保守结构域，经筛选得到候选基因。基于鸭茅的基因组注释信息，使用在线软件MapGene2Chrom (http://mg2c.iask.in/mg 2c_v2.0/)对候选DgbZIP基因的染色体定位图进行绘制。

1.3 鸭茅bZIP 基因家族的系统发育分析

使用MEGA 7.0[19]软件中的ClustalW 工具对拟南芥(二倍体，2n= 10)、水稻(二倍体，2n= 24)和鸭茅(二倍体，2n= 14) bZIP 蛋白的氨基酸序列进行多序列比对，参数设置为默认，生成.meg 文件；利用Phylogeny 工具采用邻接法(neighbor joining，NJ)构建系统进化树，评估方法为自举法(bootstrap)，数值设为1 000，其他参数使用系统默认值。并利用iTOL(https://itol.embl.de/)在线网站进行进化树的美化修饰。

1.4 鸭茅bZIP 基因家族的理化性质、保守基序和基因结构分析

使用在线工具Protparam (http://web.expasy.org/protparam/)预测候选bZIP 蛋白的物理化学性质。通过MEME (https://meme-suite.org/meme/tools/meme)在线工具查询候选bZIP 蛋白的保守氨基酸基序，重复次数设为任何，motif 数目设置为15，其余参数保持默认值。最后通过TBtools v1.098696 软件[20]将鸭茅DgbZIP基因的系统进化树、保守基序、基因结构可视化并绘制三合一图。

1.5 鸭茅bZIP 基因的顺式作用元件分析

使用TBtools 软件提取鸭茅bZIP基因CDS 序列上游2 000 bp 的启动子序列，并利用Plant CARE(http://bioinformatics.psb.ugent.be/webtools/plantcare/ht ml/)预测启动子区域潜在的顺式作用元件。

1.6 非生物胁迫下鸭茅bZIP 基因的表达模式分析

从四川农业大学草学系实验室牧草课题组前期研究中收集了DgbZIP基因在不同组织和不同非生物胁迫下的转录组数据，包括‘宝兴’鸭茅叶片和根干旱处理18 d[21]和‘宝兴’(耐热)和‘01998’(热敏感)鸭茅热胁迫10 和26 d[22]后DgbZIP基因的表达量数据，从而分析bZIP基因在不同鸭茅品种、不同组织、胁迫处理不同时间下差异基因的表达情况并推测其相应功能。若log2(FC) ＞ 1 或 ＜ -1，则认为该基因显著上调或下调并绘制DgbZIP基因的表达热图(P＜ 0.05)。

1.7 qRT-PCR 验证bZIP 基因胁迫诱导下的表达模式

‘宝兴’鸭茅种子来自四川农业大学草业科技学院草学系实验室，遗传背景一致。

将栽种有‘宝兴’鸭茅种子的白色塑料方盆置于植物生长室中培养，生长环境设为白天23 ℃，光照16 h，夜晚18 ℃，黑暗8 h。待鸭茅长到三叶期时，选取长势一致的幼苗进行胁迫处理。利用20%的PEG-6000 浸根处理模拟干旱胁迫[23]，以1/2 霍格兰营养液为对照；40 ℃高温模拟热胁迫。分别在0、2、4、8、12、24、48 h 进行取样，每个处理3 个重复，放于-80 ℃冰箱中保存备用。

使用植物总RNA 提取试剂盒(美基生物科技有限公司，中国广州)提取RNA，并用NanoDrop2000微量分光光度计和琼脂糖凝胶电泳检测RNA 纯度和完整性。然后用反转录试剂盒(莫纳生物科技有限公司，中国武汉)反转录得到cDNA。将合成的cDNA 用ddH2O 稀释10 倍，保存于-20 ℃冰箱中备用。利用Primer Premier 5 (http://www.premierbiosoft.com/primerdesign/)软件设计荧光定量PCR 引物，扩增长度介于80～300 bp，并用本地BLAST 检测引物特异性，具体引物序列如表1 所列。

表1 荧光定量PCR 所用引物Table 1 Primers used for qRT-PCR analysis

以反转录获得的cDNA 为模板，GAPDH为内参基 因，Bio-RAD CFX Connet 仪 器 进 行qRT-PCR 试验，设置3 次技术重复。总反应体系为10 μL，包含1 μL cDNA，0.5 μL 上游引物，0.5 μL 下游引物，5 μL SYBR 酶，其余用ddH2O 补齐。PCR 反应程序为95℃，30 s；95 ℃，10 s，60 ℃，30 s，40 个循环。采用2-ΔΔCt法评估基因相对表达水平的倍数差异，基因上调或下调1.5 倍以上被认为是表达差异显著(P＜ 0.05)。

1.8 数据处理

使用Excel 2016 整理数据，SPSS 23 进行单因素方差分析，运用GraphPad Prism 8.0 进行绘图。

2 结果与分析

2.1 基因家族成员鉴定和染色体定位

使用HMMER 3.3.2 软件和结构域预测网站SMART，共筛选了103 个鸭茅bZIP基因。按照染色体位置顺序编号，将鸭茅bZIP基因依次命名为DgbZIP1～DgbZIP103。

101 个bZIP基因全部定位在鸭茅的7 条染色体上，另外两个(DgbZIP1、DgbZIP2)定位在两个不同的Scaffolds 上(图1)。总的来说，101 个bZIP基因相对均匀地分布在鸭茅7 条染色体上，其中分布DgbZIP基因最多和最少的染色体分别是第4 条和第1 条染色体，分别有18 个bZIP基因(17.8%)和10 个bZIP基因(9.9%)。此外，不同DgbZIP基因在同一条染色体上分布的具体位置也有差异，如DgbZIP基因大多分布在2、5、6、7 号染色体的两端，中部分布较少。

图1 鸭茅bZIP 基因的染色体定位Figure 1 Chromosomal locations of DgbZIP genes

2.2 系统发育分析

参照拟南芥[3]和水稻[4]bZIP基因家族的分类标准，将103 个鸭茅bZIP基因分为11 个亚家族，即A～I、N 和S，每个亚家族包含的DgbZIP基因数目分别为12、3、6、14、6、3、10、3、11、8 和27 个(图2)。

图2 鸭茅、拟南芥和水稻bZIP 基因系统发育分析Figure 2 Phylogenetic analysis of DgbZIPs, AtbZIPs, and OsbZIPs

2.3 理化性质分析

鸭 茅 bZIP 蛋 白 的 氨 基 酸 数 目 从 129 个(DgbZIP65)到655 个(DgbZIP84)不等；分子量大小为14.19 (DgbZIP65)～70.32 kDa (DgbZIP84)，等电点范围介于4.47 (DgbZIP81)～11.49 (DgbZIP56)。

2.4 保守基序和基因结构分析

使用MEME 5.3.0 程序预测候选bZIP 蛋白中的保守氨基酸基序，共搜索到15 个motif。大部分基因所编码的蛋白中都有motif 1 或motif 2 (图3)，90.3%的bZIP 蛋白含有motif 1，82.5%的bZIP 蛋白含有motif 2。同一亚族的DgbZIP 蛋白含有相似的motif 组成和位置，不同亚族间存在差异。如I 亚族只含motif 1、motif 2 和motif 12；motif 1、motif 2、motif 13、motif 15分布于A 亚族；motif 1、motif 3、motif 4、motif 6、motif 8、motif 10、motif 14 分布于D 亚族；motif 2、motif 5、motif 7motif 9、motif 11 分布于N 亚族。某些保守基序具有特异性，如motif 15 特异性分布在A 亚族，只有D 亚族成员含有motif 4、motif 6、motif 8、motif 10、motif 14，只有I 亚族成员含有motif 12。

图3 鸭茅bZIP 基因进化树、motif、基因结构三合一图Figure 3 Phylogenetic relationships, motif pattern, and gene structure of DgbZIP genes in Dactylis glomerata

鸭茅bZIP基因的外显子数量及其长度变化较大，数量从1～13 个(DgbZIP16、DgbZIP76)不等，23 个bZIP基因(22.3%)无内含子且大多集中分布于S 亚族。同一亚族内的基因内含子—外显子结构相似，如C 亚族bZIP基因均含有6 个外显子、5 个内含子；除DgbZIP28、DgbZIP39外，E 和I 亚族bZIP基因均含有4 个外显子、3 个内含子；D、G 亚族成员的外显子和内含子数量较多；S 亚族成员大多含有1～3个外显子、无内含子，同时也说明这些亚族较为保守。

2.5 顺式作用元件分析

本研究利用Plant CARE 在线网站对103 个bZIP基因家族成员转录起始位点上游2 000 bp 的启动子序列进行分析，按照功能可将其分为光响应元件、胁迫响应元件以及与植物生长发育相关的元件三大类(图4)。其中，光响应元件主要有Sp1、G-box、G-Box、Box-4、GT1-motif 等；与胁迫相关的顺式作用元件有ABRE (ABA 响应元件)、ARE (厌氧诱导的反应元件)、MBS (干旱诱导MYB 结合位点)、TCrich repeats (应激反应元件)、DRE (干旱、低温、盐胁迫响应元件)、CGTCA motif 和TGACG motif (茉莉酸甲酯响应元件)等。88%的基因含有ABRE 元件，只有DgbZIP12基因含有DRE 元件。88%的bZIP基因含有CGTCA motif 和TGACG motif，且这两种茉莉酸甲酯响应元件在鸭茅bZIP基因家族成员中的分布数量相同。

图4 顺式作用元件分析Figure 4 Cis-acting elements in the promoters of DgbZIP genes

2.6 非生物胁迫下鸭茅bZIP 基因的表达模式分析

2.6.1 干旱胁迫下鸭茅bZIP基因在不同组织中的表达模式分析

从Ji 等[21]关于结合小RNA、RNA 和降解组测序揭示鸭茅叶片和根系适应干旱胁迫的microRNAmRNA 表达模式的研究中收集了‘宝兴’鸭茅叶片和根组织中干旱胁迫处理18 d 后的转录组数据，绘制了DgbZIP基因的表达热图(图5a)。干旱处理前2/3 的基因在根中的表达量高于其在叶中的表达量。经干旱胁迫处理后，在叶和根中分别有34 和49 个bZIP基因差异表达。同一bZIP基因在根、叶不同组织中的表达趋势不同，如DgbZIP27(I 亚族)在叶中表达上调，在根中表达下调，DgbZIP21(H 亚族)、36(S 亚族)的表达趋势与之相反，且与对照相比，DgbZIP27和DgbZIP21胁迫后分别在叶和根中表达增长了几十倍。S 亚族成员DgbZIP35和DgbZIP92在干旱处理前后的表达量都很高，分别在叶和根中上调和下调表达。

图5 bZIP 基因在鸭茅不同组织(a)和品种(b)中的表达模式热图Figure 5 Heatmap showing the expression pattern of bZIP genes under different (a) tissues and (b) varieties

2.6.2 热胁迫下不同鸭茅品种中bZIP基因的表达模式分析

从Huang 等[22]关于通过转录组分析鉴定高温胁迫10 和26 d 后两种基因型鸭茅(耐热品种‘宝兴’和热敏感品种‘01998’)的差异表达基因并开发分子标记的研究中收集相关转录组数据。如图5b所示，在两个品种中，热胁迫处理10 d 后差异表达的bZIP基因数量多于热胁迫处理26 d 的，上调表达的基因数量多于下调表达的。DgbZIP基因的表达具有品种特异性，如DgbZIP11只在‘宝兴’鸭茅中显著下调表达，DgbZIP101只在‘01998’鸭茅中显著下调表达；热胁迫处理10 d 后，DgbZIP14在‘宝兴’鸭茅中下调表达，在‘01998’鸭茅中上调表达。同一鸭茅品种DgbZIP基因在热胁迫不同时间表达情况有异，如热胁迫处理10 d 后，‘宝兴’鸭茅中DgbZIP13、DgbZIP77和‘01998’鸭茅中DgbZIP35、DgbZIP69显著上调表达，‘01998’鸭茅中DgbZIP92显著下调表达；热胁迫处理26 d 后，DgbZIP73在两个材料中均显著上调表达。热胁迫处理后，DgbZIP12、DgbZIP36在两个材料中均下调表达。上述差异表达的基因大多属于A、S 亚家族，推测这两个亚族基因成员可能是鸭茅耐热的关键基因。

2.6.3 qRT-PCR 验证bZIP基因在非生物胁迫诱导下的表达模式

基于干旱胁迫下鸭茅bZIP基因在根、叶中的转录组数据和热胁迫下‘宝兴’和‘01998’品种鸭茅bZIP基因的转录组数据，筛选出9 个差异表达较为显著的DgbZIP基因，利用qRT-PCR 分析其在干旱和热胁迫下的表达模式。对9 个DgbZIP基因相对于内参基因(GAPDH)进行归一化处理。

DgbZIP6在干旱处理后各个阶段表达明显下调(图6)，其相对表达量介于0.01 (干旱4 h)～0.10 (干旱48 h)。DgbZIP76、DgbZIP101在干旱处理后各个阶段上调表达，均在干旱处理4 h 时表达量最高，分别 是CK 的136 和59 倍。除DgbZIP6、DgbZIP13、DgbZIP21、DgbZIP58外，其余基因经热处理后在各个阶段均上调表达，大部分基因在热处理24 h 时表达量最高，其中DgbZIP84和DgbZIP101在24 h 时的表达量分别是其CK 的128 和185 倍。总的来说，各个基因在干旱胁迫处理下，胁迫前期(2～8 h)上调表达，胁迫后期(12～48 h)表达量下降。热胁迫处理下总体上表现为上调表达，尤其在胁迫后期(12～48 h)表现明显。

图6 鸭茅DgbZIP 基因在干旱和热胁迫下的qRT-PCR 分析Figure 6 Expression profiles of DgbZIP under drought and heat stresses as assessed by qRT-PCR

3 讨论与结论

bZIP 转录因子是一类广泛存在于真核生物中的调节蛋白，对植物生长发育、抵御逆境胁迫等至关重要。随着鸭茅高质量参考基因组的公布[18]，使用生物信息学方法对鸭茅bZIP基因家族进行鉴定和分析成为可能。

本研究共鉴定出103 个基因家族成员。按照各个基因在7 条染色体上分布的位置，依次将其命名为DgbZIP1～DgbZIP103。对拟南芥、水稻和鸭茅bZIP基因家族进行系统发育分析，可将103 个DgbZIP基因分为11 个亚家族。其中N 亚族是本研究中新命名的亚家族，该亚族未与拟南芥和水稻的bZIP基因聚类，这可能是由于其进化方向与拟南芥和水稻bZIP基因不同，但具体的基因结构与功能还有待进一步探索[24]。S 亚族包含成员数量最多，这与拟南芥的聚类情况相同[25]。内含子的插入和缺失在植物进化过程中非常普遍。鸭茅DgbZIP基因内含子数量在0～14 个，而D 亚族和G 亚族的内含子数量较多，柳枝稷(Panicum virgatum)PvbZIP基因家族中也发现了同样的结果[26]。顺式作用元件分析发现众多胁迫响应元件，其中88%的基因含有ABRE元件。ABRE 是ABA 的应答元件，保守序列为PYCGTGGC，在植物体内水分不足时，其内源ABA含量大量增加，相关干旱胁迫应答基因受到ABA诱导而表达[27-28]。所有DgbZIP基因都至少含有一种胁迫相关顺式作用元件，表明DgbZIP基因在响应非生物胁迫过程中发挥重要作用[29]。

研究表明，同一亚族内成员往往行使相似的生物学功能[30]。拟南芥bZIP基因家族中，A 亚族成员ABF 转录因子和S 亚族成员AtbZIP2可以响应干旱、ABA、热、低温、高盐等多种胁迫[31-32]，由此推测DgbZIP基因家族中A、S 亚族基因可能是鸭茅抗逆的关键基因。基于干旱和热胁迫下DgbZIP基因的RNA-seq 数据，发现DgbZIP基因的表达具有组织、品种、时间特异性，这与黄瓜(Cucumis sativus)CsbZIP的研究结果一致[33]。选取9 个差异表达较为显著的DgbZIP基因进行qRT-PCR 检测分析，干旱胁迫后DgbZIP6均下调表达，DgbZIP76、DgbZIP101均上调表达；热胁迫后，6 个基因均上调表达。两种胁迫下大多数DgbZIP基因都呈现出先上升后下降的表达趋势，大部分DgbZIP基因在干旱处理后2～8 h 表达量最高，在热处理后24 h 表达最高。综合分析表明，鸭茅bZIP基因家族S 亚族成员(DgbZIP12,DgbZIP13,DgbZIP14,DgbZIP35,DgbZIP36,DgbZIP92), A 亚 族 成 员(DgbZIP84)、 C 亚 族 成 员(DgbZIP57,DgbZIP101)在干旱和热胁迫下能被强烈诱导表达，这与先前有关A、C、S 亚族成员可能在抵御逆境胁迫中发挥作用相吻合[34]。但还需通过功能验证进一步分析其功能。

综上所述，本研究对鸭茅bZIP基因家族成员的生物信息学特征及不同组织、品种、胁迫条件下的表达模式进行了分析。A、C、S 亚族成员在不同的时间点较强烈地响应干旱和热胁迫并参与表达调控，这可能是一个复杂的过程，可为后续鸭茅bZIP基因功能的深入研究奠定理论和实践基础，也可为鸭茅抗性育种提供参考。