杜仲内生拮抗细菌DZSY21诱导玉米抗病基因表达变化的转录组学研究

王 其，陈小洁，顾双月，张欣悦，杭天露，丁 婷

(安徽农业大学 植物保护学院，安徽 合肥 230036)

植物内生细菌是指植物内生菌中的细菌，绝大部分内生细菌可与寄主植物建立和谐共生的关系且不会对寄主植物造成实质性危害[1]。相关报道表明，在植物病虫害防治方面应用植物内生细菌可有效增强植物的抗病性[2-4]。目前，有关内生细菌增强玉米抗病性的研究主要集中在病害胁迫下玉米的生理特性如苯丙氨酸解氨酶(PAL)、过氧化物酶(POD)和多酚氧化酶(PPO)等[5-6]方面，利用分子生物学手段来研究玉米应答内生细菌的抗病相关转录调控鲜有报道。

转录组学研究作为一种高通量的研究方法，能够从整体上分析植物在应答逆境胁迫过程中所有mRNA的变化情况。随着这项技术的不断发展完善，目前已通过转录组测序技术从拟南芥[7-8]、水稻[9]、玉米[10]、油菜[11]等植物中鉴定了大量的逆境胁迫应答基因，并对其抗逆作用机理进行了探讨。前期研究发现，杜仲内生拮抗细菌DZSY21可在玉米中稳定定殖并增强玉米植株的抗病能力[3]，在此基础上，本研究运用转录组测序技术，比较分析玉米在喷施杜仲内生拮抗细菌后不同时间段内基因的表达情况，并对其进行生物信息学分析，为进一步寻找玉米抗病基因及其抗病机理提供理论依据。

1 材料与方法

1.1 菌种来源、菌悬液制备

试验菌株为前期研究中获得的一株在离体条件下对玉米小斑病菌、玉米纹枯病菌具有较好抑菌活性的杜仲内生拮抗细菌[3，6]，编号为DZSY21 (BacillussubtilisDZSY21)，NCBI登录号为KP777560，已保藏于中国微生物菌种保藏管理委员会普通微生物中心(保藏编号：CGMCCNO.11749)。

菌株DZSY21接种于牛肉膏蛋白胨(NA)液体培养基中，摇床转速为160 r·min-1，于37 ℃培养8～10 h，调节发酵液中DZSY21的菌体浓度至1.0×106CFU·mL-1，获得DZSY21菌悬液。

1.2 供试玉米品种与栽培方式

供试玉米品种为自交系昌7-2，选取籽粒饱满、大小一致的种子，催芽后播种在安徽农业大学教学实习基地中，种植方式为垄栽，垄面宽70 cm，垄距30 cm，株距30 cm，行距60 cm，每垄栽植20粒健康玉米种子，共栽植2垄。待玉米长至大喇叭口期时选取长势相当的玉米植株，在健康叶片上喷洒DZSY21菌悬液(含菌量为1.0×106CFU·mL-1)，每个叶片表面喷洒5 mL，每个处理5株玉米，3个重复。在DZSY21菌悬液处理0、12和24 h时选取玉米叶片样品，经液氮速冻后置于-70 ℃保存，用于RNA的提取。DZSY21菌悬液处理0、12、24 h时的玉米样品分别记为T01、T02和T03。

1.3 叶片总RNA提取、文库构建及转录组测序

玉米总RNA提取方法参照Invitrogen公司的Trizol Reagent说明书进行。使用Nanodrop检测RNA的纯度(D260/D280)和浓度。RNA样品检测合格后，进行文库构建。cDNA文库构建方法参照RNA-Seq文库构建试剂盒(Illumina公司，美国)说明书进行。文库构建完成后，由北京百迈克生物科技有限公司完成测序，测序平台为Illumina HiSeq。

1.4 测序数据处理与注释

测序产出大量高质量的原始数据(raw data)，为保证测序分析的质量，对原始数据文件进行处理，去除含有接头的Reads和低质量的Reads(包括去除N的比例大于10%的Reads；去除质量值Q≤10的碱基数占整条Read的50%以上的Reads)之后，获得高质量的Clean Data。利用TopHat2[12]将Clean Reads与玉米B73参考基因组进行序列比对，获取在玉米B73参考基因组或基因上的位置信息，以及测序样品特有的序列特征信息。利用Cufflinks软件的Cuffquant和Cuffnorm组件，通过Mapped Reads(比对到指定的参考基因组上的Reads)在基因上的位置信息，对转录本和基因的表达水平进行定量，使用FPKM[13]作为基因表达水平的衡量指标。

1.5 差异表达基因筛选

使用DESeq[14]软件进行样品组间的差异表达分析，符合错误发现率(false discovery rate，FDR)小于0.01且差异倍数(fold change)在2倍及以上的基因被定义为样品间的差异表达基因(DEG)。

1.6 抗病显著差异表达基因筛选及分析

基于筛选的差异表达基因，使用BLAST[15]软件将差异表达的基因与NR[16]、Swiss-Prot[17]、GO[18]和KEGG[19]等数据库进行序列比对，获得差异表达基因的注释信息。根据注释信息筛选抗病相关基因，再根据其注释信息进行GO和KEGG等相关分析。

1.7 qRT-PCR分析

对于保留的另一份样品进行RNA提取，挑选10个(表1)DEG进行验证。所选反转录试剂为Prime ScriptTMRT reagent Kit RT with gDNA Eraser (Perfect Real Time)[宝日医生物技术(北京)有限公司]，荧光定量试剂盒为SYBR®Premix ExTaqTMII(Tli RNaseH Plus)[宝日医生物技术(北京)有限公司]，qRT-PCR体系(20 μL)：SYBR Green混合物12.5 μL，模板cDNA 2 μL，上下游引物各1 μL，ddH2O 4 μL。反应程序为95 ℃预变性30 s，95 ℃ 5 s，60 ℃ 30 s，共40个循环。ZmActin为内参基因，引物设计采用软件Primer Premier 6.0。实验包含两次生物学重复，每次生物学重复进行3次技术重复。结果使用ABIPRISM7300 系统相对定量分析软件进行分析，采用2-ΔΔCt法计算。

2 结果与分析

2.1 测序数据结果

DZSY21菌悬液处理0 h(T01)、12 h(T02)、24 h(T03)的3个样品的转录组测序结果如表2所示，上述3个处理分别获得25 143 751、35 383 050、35 383 050条Raw Reads，经过对原始数据的整理和筛选分析，总共获得21.93 Gb的Clean Data，其中各样品Clean Reads分别是24 188 231、34 070 148和29 095 437；各样品的Clean Data均达到6.08 Gb。 GC含量分别为58.00%、56.11%和56.93%；Q30碱基百分比均在87.03%及以上，满足后续分析要求。

2.2 差异表达基因筛选

根据Fold change≥2且FDR<0.01的筛选标准，以DZSY21菌悬液处理0 h(T01)样品为参照，比对发现DZSY21菌悬液处理玉米叶片12 h(T02)时，T01与T02之间有2 413个差异表达基因，其中1 278个DEGs为上调表达，1 135个DEGs下调表达；DZSY21菌悬液处理玉米叶片24 h(T03)时，T01与T03之间有737个差异表达基因，其中538个DEGs为上调表达，199个DEGs为下调表达，如图1所示，在DZSY21菌悬液处理后，玉米叶片中上调基因数量明显多于下调基因数量。

表1qRT-PCR所用引物序列

Table1Primers used in qRT-PCR

表2测序数据统计

Table2Sequencing data statistics

样品SamplesRawReadsCleanReadsGC含量GC content/%≥Q30/%T01251437512418823158.0087.03T02353830503407014856.1187.32T03353830502909543756.9387.08

Raw Reads，Raw Data中pair-end Reads总数；Clean Reads，Clean Data中pair-end Reads总数；GC含量，Clean Data中G和C两种碱基占总碱基的百分比；≥Q30，Clean Data质量值大于或等于30的碱基所占的百分比。

Raw Reads， The number of pair-end Reads in Raw Data; Clean Reads， The number of pair-end Reads in Clean Data; GC content， The percentage of G and C in the total base in Clean Data; ≥Q30， The percentage of the base whose quality value was more than or equal to 30 in Clean Data.

2.3 抗病差异表达基因筛选

在获得差异表达基因(DEGs)的基础上，依据差异表达基因(DEGs)的功能注释，结合已有的差异表达基因(DEGs)的功能报道，以抗病功能为筛选指标，从差异表达基因中筛选与抗病相关的基因，在处理12 h的样品中，筛选得到218个与抗病相关的基因；处理24 h的样品中，获得71个基因；其中，两个不同时间段获得的抗病相关基因中有重复基因23个，最终，通过对12 h和24 h两个样品的分析，共获得267个抗病相关基因，这些基因主要涉及32类抗病相关途径，结果如表3显示。

图1 DZSY21处理12 h和24 h后玉米叶片中差异表达基因的统计Fig.1 Statistic of DEGs in maize leaves at 12 h and 24 h in DZSY21 treatments

以DZSY21菌悬液处理0 h为对照，结果显示：处理12 h时与抗性相关的差异表达基因共有218个基因，其中130个上调表达，88个下调表达，相关基因主要涉及编码脂质转移蛋白(bifunctional inhibitor/lipid-transfer protein)、MATE转运蛋白(MATE efflux family protein)及LysM受体蛋白激酶(LysM domain receptor-like kinase)等26类抗病相关途径。处理24 h时，与抗性相关的差异表达基因共有71个基因，其中55个上调表达基因，16个下调表达基因，涉及编码异黄酮还原酶(isoflavone reductase)、纤维素合成酶(probable cellulose synthase)、苯丙氨酸解氨酶(phenylalanine ammonia-lyase)、L-抗坏血酸过氧化物酶(L-ascorbate peroxidase)、GLK转录因子(probable transcription factor GLK)及ACC氧化酶(ACC oxidase)等30类抗病途径。此外，分析发现，在DZSY21菌悬液处理玉米叶片后12 h和24 h，出现重叠抗病相关差异表达基因23个，具体如表4所示。

2.4 抗病差异表达基因分析

对玉米叶片中抗病相关的差异表达基因进行GO分析，结果(图2)显示：在生物学过程(biological process)分类中主要聚集在应激反应(response to stimulus)、代谢过程(metabolic process)和单一生物过程(single-organism process)；在细胞组分(cellular component)分类中主要聚集在细胞器(organelle)、细胞(cell)和细胞部分(cell part)；在分子功能(molecular function)分类中主要聚集在核酸结合转录因子活性(nucleic acid binding transcription factor activity)、催化活性(catalytic activity)和结合(binding)。

对玉米叶片中抗病相关的差异表达基因进行KEGG分析，结果(图3)显示：T01-VS-T02中有30个抗病相关差异表达基因注释到9个代谢通路分支中，T01-VS-T03中有18个抗病相关差异表达基因注释到10个代谢通路分支中。其中主要注释为谷胱甘肽代谢(glutathione metabolism)、苯丙氨酸代谢(phenylalanine metabolism)和苯丙素生物合成(phenylpropanoid biosynthesis)代谢通路中。

A，复制；B，生长；C，免疫系统进化；D，信号；E，多生物体进程；F，生殖过程；G，细胞组件生物起源；H，定位；I，多细胞生物体进程；J，发育进程；K，生物调节；L，响应刺激；M，代谢过程；N，单生物体进程；O，细胞进程；P，大分子复合物；Q，细胞连接；R，胞外区；S，膜组分；T，细胞器组分；U，膜；V，细胞器；W，细胞；X，细胞部分；Y，转运蛋白活性；Z，受体活性；AB，分子转导蛋白活性；AC，核酸结合转录因子活性；AD，催化活性；AE，连接；AF，抗氧化活性。A， Reproduction; B， Growth; C， Immune system process; D， Signaling; E， Multi-organism process; F， Reproductive process; G， Cellular component organization or biogenesis; H， Localization; I， Multicellular organismal process; J， Developmental process; K， Biological regulation; L， Response to stimulus; M， Metabolic process; N， Single-organism process; O， Cellular process; P， Macromolecular complex; Q， Cell junction; R， Extracellular region; S， Membrane part; T， Organelle part; U， Membrane; V， Organelle; W， Cell; X， Cell part; Y， Transporter activity; Z， Receptor activity; AB， Molecular transducer activity; AC， Nucleic acid binding transcription factor activity; AD， Catalytic activity; AE， Binding; AF， Antioxidant activity.图2 抗病相关差异表达基因GO注释图Fig.2 GO annotations of disease-associated DEGs

A，自噬调节；B，植物激素信号转导；C，RNA运输；D，氨基酸和核苷酸糖代谢；E，谷胱甘肽代谢；F，苯丙氨酸代谢；G，苯丙氨酸生物合成；H，植物-病原互作；I，植物昼夜节律；J，ABC转运蛋白；K，半胱氨酸和蛋氨酸代谢；L，赖氨酸降解；M，抗坏血酸和醛糖代谢。A， Regulation of autophagy; B， Plant hormone signal transduction; C， RNA transport; D， Amino sugar and nucleotide sugar metabolism; E， Glutathione metabolism; F， Phenylalanine metabolism; G， Phenylpropanoid biosynthesis; H， Plant-pathogen interaction; I， Circadian rhythm-plant; J， ABC transporters; K， Cysteine and methionine metabolism; L， Lysine degradation; M， Ascorbate and aldarate metabolism.图3 抗病相关差异表达基因KEGG注释图Fig.3 KEGG annotations of disease-associated DEGs

2.5 qRT-PCR验证

通过上述抗病基因的筛选，共获得267个抗病相关基因。任意选取其中的10个基因进行qRT-PCR验证，结果如图4所示。将差异表达基因表达趋势与转录组数据的表达趋势进行比较，发现两者趋势相同，说明了转录组数据的准确性。

3 讨论

基因调控在植物应答病害胁迫的过程中发挥重要的作用，在病害胁迫下植物中抗性基因呈现出不同的应答反应模式和表达变化水平。本研究结果表明，内生菌DZSY21引入玉米叶片的不同时间段内可以诱导或抑制不同抗性基因的表达变化。

研究发现，内生菌DZSY21处理玉米叶片12 h时，涉及编码脂质相关蛋白、LysM受体蛋白激酶、几丁质酶、过氧化物酶和VQ-motif家族蛋白等5类蛋白的相关基因仅上调表达；处理玉米叶片24 h时，涉及编码纤维素合成酶、索马甜功能域蛋白家族、丝氨酸/苏氨酸蛋白激酶、肉桂酰辅酶A还原酶、凝集素蛋白激酶、抗病蛋白、VQ-motif家族蛋白、几丁质酶、过氧化物酶、UDP-糖基转移酶、MYB转录因子、L-抗坏血酸过氧化物酶、GLK转录因子等17类蛋白的相关基因仅上调表达。对比分析发现，随着内生菌DZSY21处理时间的延长，编码脂质相关蛋白和LysM受体蛋白激酶的基因表达逐渐消失，与此同时，编码植物抗病相关蛋白及酶类的基因则呈显著上调表达趋势，推测，编码抗病相关蛋白及酶类基因的显著上调表达可能抑制了编码脂质相关蛋白和LysM受体蛋白激酶等基因的表达，从而出现基因表达消失的现象。此外，研究中发现，在内生菌DZSY21引入玉米的时间延长24 h，玉米植株中编码丝氨酸/苏氨酸蛋白激酶和凝集素蛋白激酶等基因，以及GLK和MYB等转录因子家族的转录因子均呈上调表达。丝氨酸/苏氨酸蛋白激酶和凝集素蛋白激酶均属于植物类受体蛋白激酶(RLK)且定位在细胞质膜上，主要通过胞外信号分子与其胞外区结构域特异结合，结合后激活胞内激酶域而完成跨膜信号的转导；它可被病原菌的相关分子模式PAMP或MAMP识别，而激活下游的PTI抗性反应[20]，在植物体内参与信号转导、抗逆反应和病原反应等途径[21]。本研究中，多种RLK相关基因在DZSY21处理24 h的玉米中上调表达，说明内生菌DZSY21引入玉米，可诱导RLK相关基因的表达，使玉米产生PTI免疫反应，增强其抗病性。GLK作为GARP家族的一类，在水稻、玉米、辣椒等植物中广泛存在[22]，已有研究表明，GLK可增强拟南芥对禾谷镰刀菌等病原菌的抗性[23-24]。而MYB类转录因子已被证实广泛参与类黄酮代谢途径的调控，而类黄酮可在植物受到病原侵染时被诱导合成，是具有广谱抗菌活性的植物保护素[25]。本试验中，在生物因子DZSY21的胁迫下，诱导玉米GLK和MYB转录因子上调表达，试验结果与已有的研究报道相一致，可知，GLK和MYB转录因子可通过增强玉米对病原菌的抗性以及增加植株体内类黄酮含量等方式提高玉米植株抗病性。

图4 差异表达基因的qRT-PCR验证Fig.4 Validation of DEGs data by qRT-PCR

表3抗病相关差异表达基因数量统计信息

Table3Statistic information of disease related DEGs

基因功能Gene functionT01-VS-T02上调Up-regulated下调Down-regulatedT01-VS-T03上调Up-regulated下调Down-regulated脂质转移蛋白Bifunctional inhibitor/lipid-transfer protein2000异黄酮还原酶Isoflavone reductase0001乙烯响应转录因子Ethylene-responsive transcription factor5101纤维素合成酶Probable cellulose synthase0010索马甜功能域蛋白家族Putative thaumatin domain family protein2110丝氨酸/苏氨酸蛋白激酶Serine/threonine-protein kinase31120肉桂酰辅酶A还原酶Cinnamoyl-CoA reductase2130凝集素蛋白激酶Lectin protein kinase family protein5510枯草杆菌蛋白酶样蛋白酶Putative subtilase family protein1521抗病蛋白Putative disease resistance protein10180MATE转运蛋白MATE efflux family protein2100几丁质酶Chitinase4020过氧化物酶Peroxidase9020谷胱甘肽S-转移酶Glutathione S-transferase1511富含亮氨酸重复的蛋白激酶家族蛋白2712104Putative leucine-rich repeat receptor-like protein kinase family protein富含半胱氨酸的类受体激酶Cysteine-rich RLK (RECEPTOR-like protein kinase)7331苯丙氨酸解氨酶Phenylalanine ammonia-lyase0001WRKY转录因子WRKY transcription factor7311VQ-motif家族蛋白 VQ motif family protein3010UDP-糖基转移酶UDP-glycosyltransferase4620MYB转录因子MYB transcription factor91140L-抗坏血酸过氧化物酶L-ascorbate peroxidase0020lysM受体蛋白激酶LysM domain receptor-like kinase4000

续表3

表4DZSY21处理玉米不同时间段的重叠抗病相关差异表达基因

Table4Statistic information of common disease releated DEGs after DZSY21 treatment

基因编号Gene ID差异倍数FC(12 h/0 h)差异倍数FC(24 h/0 h)基因功能Gene functionMaize_newGene_25824.214.99抗病蛋白Putative disease resistance proteinGRMZM5G8913715.626.77枯草杆菌蛋白酶样蛋白酶Putative subtilase family proteinGRMZM5G8700772.16-2.11F-box蛋白F-box domain containing proteinGRMZM2G4686103.523.44富含半胱氨酸的类受体激酶Cysteine-rich RLK (RECEPTOR-like protein kinase)GRMZM2G4559092.776.12抗病蛋白Putative disease resistance proteinGRMZM2G4553215.299.20抗病蛋白Putative disease resistance proteinGRMZM2G441031-7.55-3.55WRKY转录因子WRKY transcription factorGRMZM2G4163222.062.04富含半胱氨酸的类受体激酶Cysteine-rich RLK (RECEPTOR-like protein kinase)GRMZM2G4026315.293.54抗病蛋白Putative disease resistance proteinGRMZM2G3612563.082.70ABC转运蛋白ABC transporter G family memberGRMZM2G3122264.595.17几丁质酶ChitinaseGRMZM2G3017383.433.38凝集素蛋白激酶Lectin protein kinase family proteinGRMZM2G1591194.345.79MYB转录因子MYB transcription factorGRMZM2G1309044.493.83过氧化物酶PeroxidaseGRMZM2G1233143.212.95富含亮氨酸重复的蛋白激酶家族蛋白Putative leucine-rich repeat receptor-like protein kinase family proteinGRMZM2G0563352.422.42UDP-糖基转移酶 UDP-glycosyltransferaseGRMZM2G055684-4.04-2.25枯草杆菌蛋白酶样蛋白酶Putative subtilase family proteinGRMZM2G032856-2.72-2.02谷胱甘肽S-转移酶Glutathione S-transferaseGRMZM2G0131705.949.24抗病蛋白Putative disease resistance proteinGRMZM2G0129262.382.14抗病蛋白Putative disease resistance proteinGRMZM2G0128612.982.83富含亮氨酸重复的蛋白激酶家族蛋白Putative leucine-rich repeat receptor-like protein kinase family proteinAC234526.1_FG0052.642.70肉桂酰辅酶A还原酶Cinnamoyl-CoA reductaseAC204711.3_FG002-2.31-2.66富含亮氨酸重复的蛋白激酶家族蛋白Putative leucine-rich repeat receptor-like protein kinase family protein

本试验发现，DZSY21菌株处理24 h，玉米中AP2/EREBP转录因子、乙烯响应转录因子以及苯丙氨酸解氨酶均呈现负调控现象。AP2/EREBP转录因子与植物对逆境胁迫的应答有着密切的关系，根据AP2/EREBP结构域的数量，可以把AP2/EREBP转录因子家族分成2个亚族：AP2亚族(具有1个AP2/ERF结构域)和EREBP亚族(具有2个AP2/ERF结构域)。AP2亚族转录因子有调控花、胚珠和种子发育的功能， 而EREBP亚族转录因子的主要功能是调节植物对激素(乙烯和ABA等)、病原和胁迫(低温、干旱及高盐)等的应答反应[26]。而乙烯响应转录因子大部分属于ERFs(ethylene responsive factors)转录因子家族。ERF转录因子在转录水平主要通过激活或抑制下游抗病相关基因表达，提高作物抗病性[27]。本研究中的相关AP2/EREBP转录因子以及乙烯响应转录因子的表达呈下调，说明拮抗菌引入玉米植株，不能诱导寄主植物形成乙烯信号途径，此结果再次证实了本课题组前期的研究发现[3]，即拮抗菌株DZSY21引入玉米叶片，可诱导玉米中的水杨酸(SA)和茉莉酸(JA)信号通路。此外，苯丙氨酸解氨酶(PAL)作为苯丙烷类代谢过程的关键酶，催化苯丙氨酸向肉桂酸的转化，促使产生木质素、黄酮类、植物抗毒素等多种次生代谢产物[28]。相关研究表明，PAL的表达和某些功能蛋白的表达呈正相关[28]，如拟南芥中的KFB蛋白能和PAL在植物体内、体外相互作用，KFB蛋白表达的下调直接影响PAL，使其表达下调，从而增加植株体内苯丙素类中间次生代谢物的积累，间接增强植物抗性[29-30]。由此可知，本研究中编码苯丙氨酸解氨酶的相关基因受到某种蛋白表达的影响而呈下调趋势，造成苯丙氨酸解氨酶的数量减少，进而增加苯丙素类中间次生代谢物的积累达到抗病的目的。

植物基因数目众多，因此其调控方式表现出多样性。从本试验测序数据分析的结果可以看出，受外源内生拮抗菌DZSY21诱导的玉米在不同时间段，其编码不同种类的蛋白、酶以及转录因子相关的基因上调和下调表达的数目或表达变化程度不同；有的基因家族中既有上调表达也有下调表达，究其原因，这种表达的差异可能是受外源拮抗菌的诱导，植物重新组织和调节生理生化活动，以提高或改善某些代谢途径来适应生存环境的变化，减少伤害[31]。此外，同一玉米品种在不同时间段，所表达的基因也有部分相同，表明其在外界生存环境变化时，在生长的不同阶段植株有某种相同的反应机制，这种机制可能与基础抗性有关。