软腐病菌侵染下菜心叶片的转录组分析

曾小玲 赵瑞丽 钟开勤 朱朝辉 谢鑫鑫 温庆放

摘  要：為探讨软腐病菌侵染下菜心转录组功能基因信息，采用BGISEQ-500高通量测序技术对软腐病菌侵染下菜心叶片进行转录组测序，平均每个样获得47.27 M clean reads，Q20值均大于95%，共检测到表达的Unigene为36 760个，其中已知的有35 327个，预测新Unigene有1 433个;共检测出19 549个新转录本，其长度主要集中在300～2000 nt。功能注释结果显示，有32 047个Unigene在Nr数据库获得注释，其中注释到芸薹属白菜上的Unigene最多;GO功能共注释到12 588个Unigene;KEGG数据库注释到18 583个Unigene，涉及到137个代谢通路。共获得21 776个组间差异表达基因（differential expression genes， DEGs），其中对照与发病前期组间的DEGs有15 007个，6 137个上调表达，8 870个下调表达;发病前期与发病中期组间DEGs有13 118个，10 278个上调表达，2 840个下调表达;发病中期与发病后期组间DEGs有11 293个，1 790个上调表达，9 503个下调表达。DEGs的Ven图分析显示5 110个基因在每组间均存在差异，DEGs功能分析显示DEGs参与泛素蛋白降解途径、过氧化物酶体途径、光合碳固定途径、糖酵解途径等多种生命活动。本研究结果为深入开展菜心抗软腐病基因组学和分子生物学研究奠定基础。

关键词：菜心;软腐病;转录组;差异表达基因

中图分类号：S436.34      文献标识码：A

Abstract： To explore the functional gene information in flowering Chinese cabbage （Brassica campestris L. ssp. Chi-nensis （L.） var. utilis Tsen et Lee） inoculated by soft rot disease， the transcriptome of flowering Chinese cabbage was sequenced by BGISEQ-500 high-throughput sequencing technology. The results showed that average obtained 47.27 M clean reads and the Q20 base was more than 95%. A total of 36 760 Unigene transcripts were detected， of which 35 327 were known and 1 433 were predicted. A total of 19 549 new transcripts were detected， the length of which was mainly between 300 and 2000 nt. Functional annotation results showed that 32 047 Unigene were annotated in NR database， among which Unigenes were annotated most on Brassica rapa ssp. Pekinensis. GO function was annotated with 12 588 Unigene. Kegg database was annotated with 18 583 Unigene， involving 137 metabolic pathways. A total of 21 776 differential expression genes （DEGs） were obtained， including 15 007 DEGs， 6 137 up-regulated genes and 8 870 down- regulated genes between the control and premorbid groups. There were 13 118 DEGs， 10 278 up-regulated genes and 2840 down-regulated genes between pre-onset and mid-onset. There were 11 293 DEGs， 1 790 up-regulated genes and 9503 down-regulated genes. Analysis of the DEGs Ven map showed that 5 110 genes were different in each group. Functional analysis of DEGs showed that the differentially expressed genes were involved in various life activities， such as ubiquitin degradation pathway， peroxisome pathway， photosynthetic carbon fixation pathway， glycolytic process pathway and so on. The results of this study would lay a foundation for further studies on the ge-nome and molecular biology of resistance to soft rot of flowering Chinese cabbage.

Keywords： flowering Chinese cabbage; soft rot; transcriptome; differentially expressed genes

DOI： 10.3969/j.issn.1000-2561.2021.05.032

菜心（Brassica campestris L. ssp. Chinensis （L.） var. utilis Tsen et Lee）为十字花科芸薹属白菜亚种的一个变种，又名菜薹，是中国南方的特产蔬菜之一。近年来，随着菜心栽培面积的扩大和重茬现象的出现，菜心病害逐年加重，其中细菌性软腐病是危害菜心的主要病害之一。蔬菜细菌性软腐病是由软腐欧文氏菌胡萝卜亚种（Erwinia carotovora ssp. carotovora， Ecc）引起的一种流行性死体营养型细菌病害，其病原菌可存于植物表面和土壤中，通过伤口或自然孔口进入植物体内，主要依靠Types II蛋白分泌系统合成和分泌水解酶降解植物细胞壁，使植物组织呈浸渍状[1-2]，主要侵染大白菜[3]、魔芋[4]、西葫芦[5]等蔬菜。该病害发病较重时可造成蔬菜严重减产，而且在储藏运输过程中易造成腐烂，严重影响品质。目前，现有的蔬菜育种材料中具有抗软腐病的材料较少，且关于该病害的研究主要集中于物理防治[3-4， 6]、转基因[7-8]和分子标记[9]等方面，而对于该病害的发病机理尚不清楚。因此，开展与软腐病抗性相关的功能基因的研究，对提高现有育种材料的抗软腐病能力、培育具有抗软腐病蔬菜新品种有重要的意义。

随着测序技术的发展，利用转录组分析植物功能基因的研究成果层出不穷。练冬梅等[10]分析了冰菜盐胁迫下的转录组数据，筛选出8个与抗盐性相关的差异表达基因;王华等[11]在转录组分析的基础上研究不同水分条件下杜鹃花转录因子的表达发现，杜鹃花在不同水分胁迫下的转录调控主要是通过3种转录因子（ERF、bHLH和MYB）协同完成的;戎利勤等[12]完成了小花草玉梅的转录组测序，并从中筛选出12个与花发育紧密相关的MADS基因;张尉欣等[13]通过对镉胁迫下的菜心进行转录组测序，分析了菜心响应镉胁迫的主要通路及基因;陈静芳[14]通过分析热胁迫下菜心转录组数据，发现4份菜心材料的42个差异表达基因构成菜心共同热胁迫响应机制，2份强耐热性菜心材料的132个差异表达基因构成了特有的热胁迫响应机制。由此可见，转录组测序技术已经成为功能基因研究的高效手段。而关于软腐病菌侵染菜心的相关研究尚未见报道，其分子调控机制仍不清楚。本研究采用BGISEQ-500平台的测序技术，对软腐病菌侵染下的菜心叶片进行转录组测序，分析菜心软腐病发病过程中基因的表达变化，探究菜心对软腐病应答的分子机理，为菜心抗软腐病基因组学和分子生物学研究提供依据。

1  材料与方法

1.1  材料

以福州市蔬菜科学研究所自主研发的菜心新品种‘绿星为材料，种子经消毒点播于32孔的穴盘中，置于人工气候箱中培养。待幼苗长至3～4片新叶时，在每株植株上部第2片叶片基部的中脉注射OD值为0.5软腐病菌菌液15 μL，套透明薄袋保湿，观察其发病情况。选取发病前期、中期、后期的叶片作为材料，以未侵染的叶片作为对照，每个处理3次重复，每个重复的叶片混合均匀后置于冷冻管中用液氮速冻，于-80 ℃冰箱保存，样品委托华大基因科技股份有限公司进行转录组测序。

1.2  测序数据分析

将获得的数据用SOAPnuke软件进行过滤分析，采用HISAT软件（Hierarchical Indexing for Spliced Alignment of Transcripts）将获得的clean reads与白菜基因组数据进行序列比对，使用Bowtie2将clean reads比对到基因序列，采用FPKM（Fragments Per Kilobase of transcript per Million fragments mapped）作为衡量转录本或基因表达水平的指标，用perl脚本对获得的新转录本进行长度统计。

1.3  序列注释和功能分类

采用Blast比对工具，将菜心转录组获得的Unigene与公共数据库进行比对，将E值设置为≤ E-5，根據基因的相似性进行功能注释，得到与菜心Unigene具有最高序列相似性的蛋白，从而得到该Unigene的蛋白功能注释信息。公共数据库包括非冗余蛋白数据库（non-redundant protein database， Nr）、基因本体论数据库（gene ontology， GO）、东京基因与基金组百科全书（kyoto encyclopedia of genes and genomes， KEGG）。

1.4  差异基因分析

用DEGseq方法[15]对差异表达的基因进行检测，以Fold Change≥2且Q-value≤0.001作为显著差异表达基因筛选的标准，同时分析差异基因的表达和功能。

2  结果与分析

2.1  测序产出数据与质量评估

使用BGISEQ-500平台共测12个样品，每个样品平均产出7.09 Gb数据，Q20的值均大于95%，Q30的值均大于86%，说明测序数据质量满足后续要求。将获得的clean reads与白菜基因组进行比对，平均比对率为78.25%，比对到基因的平均比对率为59.66%。共检测到表达的Unigene数为36 760个，其中已知的Unigene有35 327个，预测新Unigene有1 433个;共检测出19 549个新转录本，其中13 279个属于已知蛋白编码基因的新的可变剪接亚型，1 466个属于新的蛋白编码基因的转录本，剩下的4 804个属于长链非编码RNA，其中新转录本长度主要集中在300～2000 nt，占78.58%（图1）。

2.2  基因注释

2.2.1  Nr功能注释  将获得的Unigene在Nr数据库中进行比对，共获得32 047条注释，注释到芸薹属的Unigene占93.08%，其中注释到白菜的Unigene最多，占70.77%，其次为油菜（17.71%），注释到甘蓝的Unigene为4.57%，注释到除此以外的其他科属植物的Unigene占比较小，每种植物均小于3%（表1）。

2.2.2  GO功能分类  GO功能共注释到12 588条菜心Unigene，可划分为分子功能、细胞组分和生物学过程3大类，又分为50个亚类。分子功能大类中注释到结合活性（10 031个）和催化活性（8 838个）的Unigene数量最多;细胞组分功能大类中注释到细胞（8 727个）、膜（7 373个）、细胞器（6 457个）、膜部分（6 423个）的Unigene数量较多;生物学过程大类中以细胞过程（7 554个）和代谢过程（6 790个）的Unigene较多（表2）。

2.2.3  KEGG代谢通路分析  菜心18 583个Unigene在KEGG数据库获得注释，共注释到137个代谢通路，其中注释到植物激素信号转导（1211个）、RNA转运（1 204个）、植物-病原菌相互作用（1 179个）、MAPK信号转导途径（1 111个）4条代谢通路的Unigene数量较多，其他代谢通路的Unigene相对较少。此结果说明，菜心的代谢活动较为活跃（表3）。

2.3  FPKM值分析

对FPKM值统计分析表明，菜心低表达量的单基因簇集中在发病后期，高表达量的单基因簇集中在发病中期。菜心单基因簇FPKM值低于0.1的有2547个，说明能检测到极低水平的基因表达，检测的灵敏度较高;FPKM值大于300的有29个，参与信号转导、氨基酸代谢、光合作用等多种生物途径。

2.4  差异基因分析

通过组间差异表达基因（differential expression genes， DEGs）分析（表4），共获得21 776个DEGs，其中对照与发病前期组间的DEGs有15 007个，6137个上调表达，8870个下调表达;发病前期与发病中期组间DEGs有13 118个，10 278个上调表达，2840个下调表达;发病中期与发病后期组间DEGs有11 293个，1790个上调表达，9503个下调表达。将DEGs作Ven图分析发现，有5110个基因在每组间均存在差异（图2）。

经基因的相似度和同源性分析，从菜心4个不同感病时期转录组数据中发掘大量与软腐病病程相关的基因。其中泛素蛋白降解途径相关基因有E1、E2、E3和UCH-L等，过氧化物酶体途径有CAT、XDH和SOD等，光合碳固定途径相关基因有PEPCK、ALDO和RubisCO等，糖酵解途径相关基因有PGK、ALDH3I1、PFK和PDH E1等，以及其他抗病相关基因WRKY、NBS-LRR、NSP和AP2等。

进一步分析与抗氧化相关的DEGs，发现参与到过氧物酶体途径的DEGs有112个，其中2-羟基酸氧化酶有8个，超氧化物歧化酶有6个，长链酰基辅酶A合成酶和脂肪酰辅酶A还原酶均有5个;参与到抗坏血酸代谢途径的DEGs有87个，其中醛脱氢酶有7个，肌醇磷酸酶有6个，UDP-糖基转移酶有5个（表5）。

3  讨论

利用高通量测序技术分析植物转录组数据，挖掘植物代谢途径、调控机制等信息已成为研究热点。利用该技术，紫背天葵[16]、党参[17]、青篱柴[18]、黄秋葵[19]、番茄[20]等物种已完成了转录组测序和功能注释，并从中发掘出一些重要功能基因加以利用。本研究利用高通量测序技术构建了软腐病菌侵染下菜心叶片的转录组数据库，对获得的序列在相关数据库进行比对和功能注释，发现植物激素信号转导、蛋白代谢、糖代谢和能量代谢等多种途径相关基因的表达均发生改变。

WRKY转录因子编码的蛋白能够与目标元件发生特异性结合，抑制或激活目标基因的表达，从而参与植物的胁迫应答、衰亡等一系列生理活动，在植物防卫反应机制中具有重要的调控作用[21]。赵雪等[22]、范荣伟等[23]研究大白菜BrWRKY33基因上游调控序列功能及其在软腐病抗性方面的作用，认为BrWRKY33在拟南芥抗软腐病SA/ET介导的防御过程中正调控下游基因的表达。菜心中WRKY基因序列BraA03g006670.3C、BraA07g 030410.3C在软腐病发病过程中均差异上调表达，正调控下游基因的表达，从而调控植株抗病能力;BraA09g035200.3C序列在软腐病发病前期差异上调表达，发病中期和后期下调表达;BraA02g 019440.3C和BraA06g041930.3C序列在软腐病发病前期和中期差异上調表达，发病后期下调表达;BraA01g006320.3C、BraA07g042700.3C和BraA09g029190.3C序列在软腐病发病过程中差异下调表达，负调控下游基因的表达，改善植株抗病能力。由此可见，菜心WRKY基因在软腐病发病过程中调控功能的复杂性。

油菜素类固醇（brassinosteroids， BRs）是一种新型植物生长激素，在调节细胞增殖、形态建成、顶端优势、叶片和叶绿体衰老、基因表达中起关键作用，其中油菜素内酯（brassinolide， BL）是生物活性最强的一种BRs[24]，外源施用可提高植株对根腐病、病毒病、白粉病的抗性[25-26]。BR合成基因编码的蛋白大多属于细胞色素P450超级家族，如CYP90B1/DWF4、CPD和DWF1等[27]。菜心细胞色素P450 85A2家族的不同基因调控软腐病的方式不同，序列BraA02g039080.3C在软腐病发病前期和中期差异上调表达，发病后期差异下调表达，序列BraA06g036240.3C在软腐病发病中期差异上调表达，前期和发病后期差异下调表达;P450 90A1（BraA10g028480.3C）在软腐病发病前期差异上调表达，发病中期和后期差异下调表达。BRI1是位于质膜上的富含亮氨酸重复序列（leucine-rich repeat repeat， LRR）的受体样蛋白激酶，是BR信号转导起始的受体[28]。BR信号通路打开的第一步是BRI1与BR类物质（如BL）特异性结合，此时BRI1的负调控蛋白BRI1 associated kinase1（BAK1）从细胞膜上解离而解除对BRI1的抑制，使BR信号向下传递[29]。菜心BAK1（BraA10g010250.3C）和BRI1（BraA06g 044360.3C）在软腐病发病中期均差异上调表达，前期和发病后期均差异下调表达。