腐皮镰刀菌侵染甘薯的转录组分析

罗勤川, 唐 伟, 马居奎, 陈晶伟, 杨冬静, 高方园, 孙厚俊, 谢逸萍, 张成玲

(中国农业科学院 甘薯研究所,江苏徐淮地区徐州农业科学研究所,农业农村部甘薯生物学与遗传育种重点实验室,江苏 徐州 221131)

甘薯(IpomoaebatatasL. Lam)是旋花科(Convolvulaceae)番薯属(Ipomoea)作物[1],具有高产稳产、适应性强等特点,世界上有一百多个国家和地区分布、种植。甘薯是一种多用途作物,茎叶和块根均可食用,还是饲料、淀粉生产、食品加工、酒精发酵等的原料作物,更是保证粮食安全的底线作物[2]。

腐皮镰刀菌(Fusariumsolani)是一种土传病原菌,寄主范围广泛,可以侵染马铃薯[3]、苜蓿[4]、大豆[5]等多种作物,引起根腐病、枯萎病等病害,严重影响了作物的产量与品质。腐皮镰刀菌(F.solani)也是引起甘薯病变的重要病原菌之一,侵染甘薯后引起甘薯根腐病(sweetpotato root rot)、甘薯腐烂溃疡病(Fusariumroot rot and stem canker)等病害[6]。甘薯腐烂溃疡病在甘薯苗期、大田期和贮藏期均可发生,苗期在薯苗茎部形成不规则黑色或褐色病斑,严重时导致整株枯死;大田期在薯块上形成轮纹状病斑、深入薯块组织内部,形成蜂窝状空腔,有苦味,贮藏期发病严重时可造成烂薯烂窖,严重影响甘薯的产量与品质,造成重要的经济损失[7]。

转录组(transcriptome)是指生物体的细胞或组织在一个特定状态下转录出来的所有RNA的总和,而通常所说的转录组学研究主要是指mRNA。转录组测序可用作新基因挖掘、基因家族鉴定、转录图谱绘制、代谢途径确定以及进化分析等多个方面[8],在医学、动物、植物等多个领域得到应用和发展。主要的方法有基因芯片技术、大规模平行测序技术以及RNA测序技术(RNA-seq)等。目前,RNA-seq在主要的大田作物病原菌研究中得到极大的应用,稻瘟病菌(Magnaportheoryzae)缺氧条件下,获得了上千个表达上调基因,这些基因多与金属离子的转运、菌丝的生长发育以及甾醇的生物合成密切相关[9];对毒性差异较大的两种大丽轮枝菌(Verticilliumdahliae),通过RNA-seq及数据分析发现多个与致病相关的黑色素、黄曲霉毒素的膜蛋白基因和耐药性有关的表达差异基因,为V.dahliae致病机理研究提供了依据[10]。通过RNA-seq对玉蜀黍黑粉菌(Ustilagomaydis)、大麦坚黑粉菌(U.hordei)和丝黑穗菌(Sporisoriumreilianum)的转录组进行分析,结合已有的基因组注释,明确100多个预测蛋白可能出现了转录后修饰,为更好地研究基因功能提供了理论基础[11]。

目前,国内针对F.solani引起的甘薯腐烂溃疡病的研究报道较少,致病机理机制研究尚不明确,因此,本研究采用RNA-seq对腐皮镰刀菌(F.solani)侵染甘薯不同时期病原的转录组进行测序,对差异表达基因进行GO(Gene Ontology)功能和KEGG(Kyoto Encyclo-pedia of Genes and Genomes)通路富集分析,结合表达量和表达倍数的差异以及关键通路,挖掘关键目的基因,分析腐皮镰刀菌侵染甘薯的机制,为今后甘薯腐烂溃疡病的致病机理研究奠定理论基础,以期为甘薯腐烂溃疡病的综合防治工作提供理论依据。

1 材料与方法

1.1 材料

腐皮镰刀菌为本研究室从发病的甘薯薯块上分离鉴定并保存,采用马铃薯葡萄糖琼脂培养基(PDA)传代培养,备用。

接种用甘薯品种徐薯32薯块由本单位提供。

1.2 方法

1.2.1 病菌侵染及取样

选择大小均匀的甘薯薯块洗净晾干,75%乙醇消毒后,切成厚度为0.5 cm的圆片,放到灭菌培养皿中备用。病原菌在28 ℃、PDA上培养5 d后,制成孢子悬浮液,孢子浓度为106CFU·mL-1,吸取200 μL接种到切好的甘薯圆片上,以滴加同体积不含孢子的灭菌水为空白对照(CK),加入灭菌水湿润的脱脂棉于一侧保湿,每处理3个重复,每个重复3个甘薯圆片,置于28 ℃培养,以接种前的菌株为测序对照并于接种后6 h、24 h、3 d、5 d后分别取样,每个重复1.0 g,液氮速冻后置于-80 ℃冰箱保存,送基迪奥生物公司测序。

1.2.2 RNA提取及反转录

使用TIANGEN多糖多酚植物总RNA提取试剂盒(天根生化科技有限公司)提取总RNA,用琼脂糖凝胶电泳检测核酸样本的完整性,同时分别用Nanodrop微量分光光度计和Agilent2100检测RNA的纯度(D260/D280、D260/D230)以及RNA片段长度,用带有Oligo(dT)的磁珠富集具有polyA尾巴的真核mRNA后,用超声波把mRNA打断。以片段化的mRNA为模板,随机寡核苷酸为引物反转录,合成cDNA。

1.2.3 cDNA文库构建及质检

纯化后的双链cDNA经过末端修复、加A尾并连接测序接头,用AMPure XP beads筛选200 bp左右的cDNA,进行PCR扩增并再次使用AMPure XP beads纯化PCR产物,最终获得文库。使用安捷伦科技有限公司DNA 1000 文库质检试剂盒进行质检。

1.2.4 转录组测序与分析

利用Illumina Hiseq 4000平台对文库数据进行转录组测序。为保证数据质量,对下机数据raw reads利用fastp质控过滤掉低质量数据,去除含接头的序列、含N(不确定碱基)比例大于10%的序列、A碱基的序列以及低质量碱基含量(Q≤20)的碱基数大于50%的序列[12]。使用短序列比对工具bowtie2[13]进行核糖体比对,用HISAT2软件[14]比对测序序列和参考基因组(GCA_002215905.1),用Stringtie[15]软件进行转录本重构,计算每个样本所有基因的表达量。利用R语言计算样本的相关性。使用DESeq[16]软件分析差异表达基因,将FDR<0.05且|log2(FC)|>1的基因标注为显著性差异基因。将差异表达蛋白向GO数据库各term映射,比对检验,对差异表达基因进行GO富集分析。将测序结果与KEGG数据库比对分析其通路,并对差异表达基因进行功能注释和分类。

1.3 qRT-PCR验证

通过对不同侵染阶段的转录组测序数据分析与KEGG通路分析,随机分析氨基酸合成代谢、碳水化合物代谢等通路中10个差异基因的表达情况,利用qRT-PCR验证结果。提取各样品总RNA,用Takara反转录试剂盒反转录合成cDNA,以反转录产物为模板,tub(ncbi_9666438)基因为内参基因,采用ABclone公司SYBR Green RT-PCR试剂盒,Real-Time PCR检测系统(bio-rad, CFX96, USA)在96孔板中进行实时定量PCR反应。反应体系20 μL:正反引物各1 μL,cDNA1 μL,2×TaqRT-PCR Mix 10 μL,ddH2O 7 μL。扩增程序:95 ℃预变性5 min;94 ℃变性15 s,60 ℃退火30 s,72 ℃延伸 10 s,40个循环,每个样品3个重复,按照2-ΔΔCT定量计算相对表达量[17]。

2 结果与分析

2.1 转录组数据分析

测序共获得1 070 504 930个下机数据,经过质量控制获得1 068 861 644个有效数据。从15个样品中均获得了8.8 G以上的原始数据,质控后碱基含量均在47%以上,测序碱基质量Q20和Q30分别在96%和91%以上(表1),说明转录组测序获得的碱基质量较好,可以用于后续分析。

表1 样品测序数据及质量检验Table 1 Sample sequencing data and quality inspection

为明确不同重复和处理样本之间的关系,构建样本相关性热图,样品之间表达模式的相似度越高,相关系数(r)就越接近1,从图1中可以看出,组间相关系数均高于0.72,说明结果可靠且样本选择合适。

图1 样本相关性热图Fig.1 Sample correlation heat map

将转录组测序数据与参考基因组(GCA_002215905.1)进行比对分析,结果如表2所示,实验组和对照组转录组测序共获得的序列750 924 327 bp,占比0.45%～72.17%;比对到参考序列唯一位置上的序列长340 313 597 bp,占比0.44%～72.09%;比对到参考序列多个位置上的序列长408 135 bp,占比0～0.08%;比对到参考基因外显子区域序列有168 683 025 bp,占比90.31%～92.15%;内含子1 394 822 bp,占比在0.53%～0.90%。

表2 测序数据与参考基因比对Table 2 Comparison on sequencing data and reference genome

2.2 差异表达基因筛选

通过不同侵染时间转录组数据和对照进行比较分析,侵染后6 h时筛选出1 056个差异表达基因(differentially expressed genes,DEGs),表达上调和下调基因分别有587和469个,这些差异基因中有201个属于代谢过程(metabolism processing),占比19.03%,103个基因属于遗传信息处理过程(genetic information processing),占比9.75%,27个基因属于细胞过程(cellular processing),占比2.56%;24 h时筛选出995个差异表达基因,表达上调和下调基因分别有679和316个,其中有210个基因属于代谢过程,占比21.10%,43个基因属于细胞过程,占比4.32%,26个基因属于遗传信息处理过程,占比2.61%;3 d时筛选出737个基因,上下调基因分别有430和307个,其中有196个基因属于代谢过程,占比26.59%,90个基因属于细胞过程,占比12.21%,42个基因属于遗传信息处理过程,占比5.70%,还有7个基因属于环境信息处理(environmental information processing);5 d时筛选出935个差异表达基因,上下调基因分别有551和384个,其中有219个基因属于代谢过程,占比23.42%,122个基因属于细胞过程,占比13.05%,67个基因属于遗传信息处理,占比7.17%,还有12个基因属于环境信息处理(图2)。进一步比对筛选,有208个基因在不同的侵染阶段都被显著调控表达,其中55个基因属于代谢过程,占比26.44%,9个基因属于细胞过程,占比4.33%,还有9个基因属于遗传信息处理过程,占比4.33%。

A、B、C、D分别为对照病菌与病菌侵染6 h、24 h、3 d、5 d的差异表达基因火山图。A, B, C and D were the volcano maps of DEGs in comparisons: CK vs 6 h, CK vs 24 h, CK2 vs 3 d, CK vs 5 d, respectively.图2 差异表达基因火山图Fig.2 Volcano diagram of DEGs

2.3 差异表达基因功能注释

2.3.1 GO富集分析

GO富集分析结果显示,在病菌侵染甘薯6 h、24 h、3 d、5 d,生物过程(biological process,BP)分别富集到650、666、526、671个差异表达基因;分子功能(molecular function,MF)分别富集到577、544、427、538个差异表达基因;细胞学组分(cellular component,CC)分别富集到506、349、325、435个差异表达基因。富集最多前40个GO功能分类中,CC类别中富集较多的是细胞(cell)、细胞部分(cell part)、细胞膜(membrane)等差异表达基因较多(图3-A);BP类别中,富集较多的是代谢进程(metabolic process)、单一生物进程(single-organism process)、细胞进程(cellular process)、生物调节(biological regulation)等差异表达基因较多(图3-B);MF类别中,富集到催化活性(catalytic activity)、结合(binding)、运输活性(transporter activity)等差异表达基因较多(图3-C)。

A:a,定位;b,刺激反应;c,多有机体过程;d,信号;e,细胞成分组织或生物发生;f,生物过程正向调节;g,生物过程负向调节;h,单一生物进程;i,细胞进程;j,节律过程;k,生长;l,繁殖;m,多细胞组织过程;n,繁殖过程;o,发育过程;p,解毒;q,生物调节;r,生物进程调节;s,代谢进程;t,细胞聚集;u,运动。 B:a,结构分子活性;b,转运活性;c,分子功能调节;d,核酸结合转录因子活性;e,抗氧化活性;f,分子传感器活性;g,结合;h,信号传感器活性;i,转录因子活性、蛋白质结合;j,催化活性。C:a,细胞连接;b,大分子复合物;c,类核;d,膜部分;e,细胞;f,细胞部分;g,膜;h,超分子纤维;i,细胞器部分;j,膜封闭腔;k,细胞器。A: a, localization; b, response to stimulus; c, multi-organism process; d, signaling; e, cellular component organization or biogenesis; f, positive regulation of biological process; g, negative regulation of biological process; h, single-organism process; i,cellular process; j, rhythmic process; k, growth; l, reproduction; m, multicellular organismal process; n, reproduction process; o, develpomental process; p, detoxification; q, biological regulation; r, regulation of biological process; s, metabolicprocess; t, cell aggregation; u,locomotion. B: a, structural molecule activity; b, transporter activity; c, molecular function regulator; d, nucleic acid binding transcription factor activity; e, antioxidant activity; f, molecular transducer activity; g, binding; h, signal transducer activity; i, transcription factor activity, protein binding; j, catalytic activity.C: a, cell junction; b, macromolecular complex; c, nucleoid; d, membrane part; e, cell; f, cell part; g, membrane; h, supramolecular fiber; i, organelle part; j, membrane-enclosed lumer; k, organelle.图3 不同处理时间差异基因GO富集柱状图Fig.3 GO functional classification of differentially expressed genes at different time after infection

2.3.2 KEGG通路富集分析

进一步对不同侵染时间的差异表达基因进行KEGG通路富集分析,分别筛选出富集最显著的20条代谢通路(图4),6 h富集的主要通路有核糖体(ribosome)、氧化磷酸化(oxidative phosphorylation)、乙醛酸和二羧酸代谢(glyoxylate and dicarboxylate metabolism)相关途径;24 h富集到的主要通路有代谢途径(metabolism pathway)、次生代谢物生物合成(biosynthesis of secondary metabolites)、抗生素生物合成(biosynthesis of antibiotics)、碳代谢(carbon metabolism)、缬氨酸、亮氨酸和异亮氨酸的降解(valine, leucine and isoleucine degradation)等;3 d富集的通路主要有代谢途径(metabolic pathway)、抗生素生物合成(biosynthesis of antibiotics)、碳代谢(carbon metabolism)等;5 d富集到的通路有代谢途径(metabolic pathway)、抗生素生物合成(biosynthesis of antibiotics)、氧化磷酸化(oxidative phosphorylation)、丙酮酸代谢(pyruvate metabolism)等。这些通路多属于碳水化合物代谢、能量代谢、氨基酸代谢、核苷酸代谢,与腐皮镰刀菌的生长发育、增殖、破除寄主防御系统等密切有关。由此看出,在侵染初期,病原菌通过消耗自身的能量,合成侵染定殖寄主所需要的各种蛋白质;在侵入寄主后,寄主植物启动了防御反应,因此病原菌通过合成抗生素来消除寄主防御机制带来的不利影响,同时,通过碳代谢途径将寄主相关物质转化能量来维持自身的增殖和侵染。

A、B、C、D分别为对照病菌与病菌侵染6 h、24 h、3 d、5 d的差异表达基因的KEGG通路。A, B, C , D were the KEGG pathway of DEGs in comparisons: CK vs 6 h after infection, CK vs 24 h after infection, CK vs 3 d after infection, CK vs 5d after infection, respectively.图4 差异表达基因的 KEGG 通路归类Fig.4 KEGG pathway classification of differentially expressed genes

Y轴(左轴)为qRT-PCR验证结果(折线表示相对表达量),Y轴(右轴)为转录组测序结果(柱状图表示FPKM值)。Y-axis (left) indicates qRT-PCR data (line graphs); Y-axis (right) indicates RNA-Seq data (bar chart). Data from qRT-PCR are means of three replicates and bars represent data from RNA-seq are means of the replicates.图5 Fusarium侵染甘薯不同时间不同基因的表达差异趋势Fig.5 The expression levels of DEGs in Fusarium infecting sweetpotato

2.4 关键基因的鉴定及其功能描述

病菌侵染寄主不同阶段均差异表达显著的基因中,选择表达量高、差异表达倍数在2倍以上的差异显著的基因(read count>20,至少在一个样品中FPKM≥50),结合代谢通路,确定了20个差异表达基因(表3),主要涉及代谢、遗传信息处理过程以及细胞过程。其中15个基因与代谢相关,参与碳水化合物代谢、能量代谢、核苷酸代谢以及氨基酸代谢;3个基因与遗传信息处理过程相关,参与翻译和蛋白质加工和修饰;2个基因与细胞过程相关,参与了细胞自噬过程。

表3 关键基因及功能描述Table 3 The key genes and their functions description

2.5 表达差异基因的qRT-PCR验证

为进一步验证转录组数据的可靠性,随机选择9个表达差异和表达趋势不同的基因进行qRT-PCR验证。结果如图所示,qRT-PCR检测的9个基因的表达趋势结果与转录组统计结果基本一致,说明转录组测序数据较为可靠。

3 讨论

腐皮镰刀菌是植物上一种土传性病原菌,可在土壤和病株残体中越冬越夏,危害包括甘薯在内的多种作物的根、茎及果实,严重影响农作物的产量与品质。本研究对侵染甘薯不同时期的腐皮镰刀菌转录组测序分析发现,接种后6 h、24 h、3 d、5 d与对照相比,分别有1 056、995、737、935个差异表达基因,这些差异表达基因主要参与了氨基酸生物合成、脂质代谢、糖代谢、嘌呤代谢等通路,与细胞蛋白质代谢过程、细胞交流、细胞进程等生物过程密切相关。在GO富集中,富集到生物过程的差异基因数量明显多于细胞组分和分子功能的基因数量,这表明腐皮镰刀菌在侵染甘薯过程中与生物进程相关的基因活性较高,参与了镰刀菌侵染过程,其中包括代谢过程、刺激应答、繁殖等GO生物过程,这与纪晓贝[18]等研究白粉菌侵染橡胶树过程转录组分析GO富集分析结果一致,这些基因的差异表达可能与病原菌侵入寄主、争夺寄主植物的营养物质,从而进一步在寄主体内定殖侵染有关。

病原菌侵染寄主过程中需要各种蛋白质的参与,在破除寄主的防御系统时需要各种水解酶类来降解寄主的细胞壁,在侵入到寄主体内,在寄主体内定殖和增殖需要各种蛋白质协助夺取寄主的营养物质,这些蛋白质是病原菌维持正常生命活动的重要组成部分,也是病原菌夺取寄主营养物质的重要工具。氨基酸是蛋白质的主要组成成分,又是许多代谢循环途径的中间体,机体对蛋白质的需求本质上是对氨基酸的需求[19]。在本研究中,腐皮镰刀菌在不同侵染阶段涉及到氨基酸代谢分别有109、138、104、110个显著差异表达的基因,其中上调基因分别有34、54、40、42个,占比分别是31.2%、39.1%、38.5%和38.2%,多涉及氨基酸的生物合成与降解,例如甘氨酸、丝氨酸和苏氨酸代谢以及精氨酸生物合成。苏氨酸脱氢酶基因MoILV1缺失,可导致稻瘟病菌分生孢子形态异常,分生孢子数量减少,附着胞介导的渗透受限,毒力减弱[20],Zhang等研究证实稻瘟病菌精氨酸合成酶基因MoARG1,6和 7缺失后,病原菌的致病性下降[21];稻瘟病菌的精氨酸甲基转移酶基因MoHMT1影响了菌丝生长发育,在病菌形态发育中起着重要的调控作用,该基因缺失突变体的菌落较小,气生菌丝较薄,生长速度明显减慢,分子孢子数量仅为野生型的20%,导致其对敏感水稻致病性降低[22]。本研究转录组数据表明,甘氨酸、丝氨酸、苏氨酸代谢或合成相关基因表达差异,可能参与了菌株的致病过程。在后续腐皮镰刀菌致病基因功能研究工作中,可以以氨基酸代谢或合成相关基因为突破口,进一步开展这些基因在菌丝生长、孢子形态,致病性以及毒素分泌等方面的功能。

不同侵染时期均差异表达显著的基因中有19个与碳代谢密切相关。表达上调基因有12个,主要与糖酵解、淀粉和蔗糖代谢等相关;下调基因有7个,主要与乙醛酸和二羧酸的代谢、丙酸酯代谢、丁酸酯代谢等相关。说明碳代谢在腐皮镰刀菌侵染甘薯过程中具有重要作用。病原菌可以利用寄主植物的营养物质进行增殖,淀粉和蔗糖代谢可以将淀粉蔗糖转化为葡萄糖,而糖酵解是所有生物体进行葡萄糖分解代谢所必须经过的共同阶段。因此,可以推测腐皮镰刀菌可以利用甘薯中的营养成分实现增殖,从而实现对甘薯寄主的进一步侵染。

4 结论

腐皮镰刀菌侵染甘薯不同时期病原菌转录组数据分析表明,腐皮镰刀菌中与表达调控破除寄主防御、代谢相关的基因在侵染甘薯过程中起着重要的作用。虽然随着侵染进程阶段的不同,不同基因在侵染甘薯过程中的表达调控方式不同,相关基因的表达趋势随之改变,但关键基因涉及的途径和功能相对集中。今后可进一步筛选验证腐皮镰刀菌破除寄主细胞壁表达调控的关键基因,进行全长克隆和功能验证等方面进行深入研究。