垂柳因姿态优美、易成活、生长快等特点，广泛应用于园林绿化当中。近年来，垂柳丛枝病不断发生且呈上升趋势。丛枝病一般多发生于阔叶树种、针叶树种和竹类上。枝条感染丛枝病后，病原物会抑制顶芽的生长，侧芽受到刺激信号提前萌发生成小枝，整株生长十分缓慢，并且在病原物的反复刺激下，顶芽得不到生长，侧芽不断抽出小枝，最终枝条变为丛生状。丛枝病主要由类菌原体和真菌侵染所致。前者是系统性丛枝病，病害由个别枝条开始，逐渐扩及全株；后者所致病害只作局部扩展，丛枝症状仅表现在直接受侵染的个别枝条上(耿站军和刘伟,2020,现代园艺,43(3):177-178)。目前，关于植原体病害的研究多集中在枣树(李继东等,2019;李玲等,2020)、泡桐(曹亚兵等,2017)、橡胶树(He et al.,2020)等具有较高经济价值的林木上，而对于垂柳丛枝病的研究却很少。对新疆阿拉尔市垂柳花变叶进行鉴定，将其划分到AY植原体组16Sr I-B亚组；加拿大发生的垂柳丛枝病与三叶草增殖组植原体相关。随后，关于垂柳丛枝病的研究，多集中于防治、研究进展等方面(Khadhair et al.,1997;耿站军和刘伟,2020,现代园艺,43(3):177-178)。

转录组受内源和外源因子的共同调控，可通过高通量测序研究植物不同时期组织生理状况下的基因差异表达水平，筛选特异性生理功能层次的高表达基因(Landers et al.,2001)。转录组学作为转录组的延伸学科，它不同于之前研究单个基因的模式，不光可以获知全部基因的表达代谢调节系统，还可以研究蛋白质的功能。基因组学研究在整体的研究方法基础上被带入了快速发展的阶段(薛建江和邱景富,2007,河北北方学院学报:医学版,24(5):63-66)。本研究利用Illumina测序平台，通过测定其RNA序列，分析健康垂柳与感病垂柳在表达定量、功能注释、SNP及SSR等方面的差异性，为垂柳丛枝病致病机理研究提供一定的参考。

1、结果与分析

1.1转录组分析

试验共采集垂柳枝叶样品6份，其中健康植株(LS＿CK) 3份、感病植株(LS＿CZ) 3份。通过原始数据过滤、测序错误率检查、GC含量分布检查等，获得样本测序数据质量汇总。测序结果显示，6份样品clean bases共计39.94 Gb，其中3个健康植株clean bases总计19.57 Gb、平均6.52 Gb,raw reads平均为23 115 859条、clean reads平均为21 750 610条；3个感病植株clean bases总计20.37 Gb、平均6.65 Gb,raw reads平均为23 670 022条、clean reads平均为21 750 610条。感病植株比健康植株总clean bases多了0.8 Gb、平均多了0.27 Gb。健康植株GC平均含量为43.84%，感病植株GC平均含量为43.73%，健康植株比感病植株高0.11%。数据整体测序错误率低于0.03%，且Q20值高于97%、Q30值高于92%(表1)，说明此次转录组测序质量较高，达到了后续分析要求，能够进行下一步信息分析。

1.2基因表达水平分析

以Trinity拼接转录组为参考序列(Ref)，用Ref做mapping。采用RSEM软件统计各样品比对到基因上的数目，根据FPKM基因表达水平开展转化分析。结果表明，从FPKM密度分布(图1A)和盒行图(图1B)可以看出，在FPKM为0时感病植株的密度明显高于健康植株，且感病植株的总体表达量明显高于健康植株，说明二者具有一定的差异性。样品间相关系数热图显示(图1C)，健康植株与感病植株最高为0.692，最低为0.482，表明二者基因的相关性比较差。相同样品间，最高为0.899，最低为0.737，说明二者之间表达模式的相似度较好。

1.3差异性表达基因的筛选

对健康植株(图2A)和感病植株(图2B)的基因表达数据进行统计学分析，筛选不同样本间显著差异基因，共检测差异性基因(DEGs) 18 536个，其中感病植株相对于健康植株表达上调的DEGS有10 316个，表达下调的DEGS有8 220个(图2D)。说明两者表达基因的种类和水平不同，感病植株通过增强或者抑制部分基因的表达达到侵染目的。差异基因韦恩图结果显示(图2C)，健康植株与感病植株有相同的差异基因51 093个，感病植株有19 017个差异性基因与健康的不同，健康植株有4 974个差异性基因与感病植不同，说明健康植株与感病植株间差异较大，其主要原因与发生垂柳丛枝病关系密切。

1.4 GO功能富集分析

对筛选的18 536个差异性基因进行GO富集性分析(图3)，结果显示，共有11 247个差异基因获得功能注释，其中上调基因6 729个，下调基因4 518个。差异基因在代谢过程和生物过程中被注释得最多；在细胞组分本体中，被注释得最多次数的条目是细胞器、蛋白复合物和细胞内细胞器；注释到分子功能本体中最多的条目为催化活性分子和氧化还原酶活性。针对差异基因的分析，在生物过程本体中下调表达的差异基因主要富集在氧化还原过程(GO:0055114)、碳水化合物代谢过程(GO:0005975) 2个方面；上调表达的差异基因主要富集在基因表达(GO:0010467)、RNA聚合酶Ⅱ转录(GO:0006366) 2个方面。在细胞组分本体中，下调表达的差异基因主要富集在类囊体(GO:0009579)、光合膜(GO:0034357) 2个方面；上调表达的差异基因主要富集在细胞器(GO:0043226)、细胞器膜结合细胞器(GO:0043227)、转录因子TFIIA复合物(GO:0005672) 3个方面。在分子功能本体中下调表达的差异基因主要富集在催化活性(GO:0003824)、ADP结合(GO:0043531) 2个方面；上调表达的差异基因主要富集在结构分子活性(GO:0005198)、转录调节活性(GO:0140110)、核糖体的结构成分(GO:0003735)、肽酶活性(GO:0008233)4个方面。

图1基因表达水平分析   

图2差异基因统计分析  

图3差异基因GO富集

1.5 KEGG通路富集分析

KEGG富集分析结果显示，5 013个差异性表达基因，上调2 951个，下调2 062个，共获得120个通路的注释，其中11个通路显著富集(P<0.05)，包括光合作用、亚油酸的新陈代谢、光合作用-天线蛋白、卟啉和叶绿素代谢、抗坏血酸和醛达酸代谢、光合作用生物的碳固定作用、乙醛酸和二羧酸的代谢、DNA复制、苯丙氨酸、酪氨酸和色氨酸的生物合成、酪氨酸代谢、异喹啉生物碱生物合成(图4)。其中，下调差异基因显著富集在植物-病原体相互作用、淀粉和蔗糖代谢和氨基糖和核苷酸糖代谢3个通路上；上调差异基因显著富集在核糖体、氧化磷酸化、真核生物核糖体生物发生、内分泌细胞增生4个通路上。注释到KEGG通路上的基因数用点的大小来表示，富集的显著性大小用从红到紫的颜色表示。

1.6垂柳丛枝病发生相关通路及候选基因筛选分析

通过对差异基因的GO和KEGG富集结果进行分析，筛选出光合作用、植物-病原相互作用、苯丙氨酸、酪氨酸和色氨酸生物合成、光合作用-天线蛋白、卟啉与叶绿素代谢、亚油酸、氨基糖和核苷酸糖代谢和α-亚麻酸代谢共计8条通路与垂柳丛枝病发生密切相关(表2)。

图4差异表达基因KEGG分类  

注:1:酪氨酸代谢;2:卟啉和叶绿素代谢;3:光合作用-天线蛋白;4:光合作用;5:苯丙氨酸,酪氨酸和色氨酸生物合成;6:吞噬体;7:戊糖和葡萄糖醛酸相互转化;8:氮代谢;9:亚油酸代谢;10:异喹啡生物碱生物合成;11:乙醛酸和二羧酸代谢;12:果糖和甘露糖代谢;13:基因复制;14:二菇生物合成;15:类胡萝卜素生物合成;16:光合生物中的碳固定;17:抗坏血酸和醛糖代谢;18:精氨酸和脯氨酸代谢;19:氨基糖和核苷酸糖代谢;20:α-亚麻酸代谢

进一步对富集的植物-病原相互作用通路中的差异基因进行分析，筛选出K02183 (CALM)、K13424(WRKY33)、K13448 (CML)、K13459 (RPS2)、K13457(RPM1,RPS3)、K13412 (CPK) 6个垂柳丛枝病发生的关键基因(图5)。其中，K13424 (WRKY33)是一个转录因子，其转录活性依赖于MPK3/MPK6介导的磷酸化修饰，MPK3/MPK6通过其SIM结构域特异识别SUMO化修饰的WRKY33，诱导WRKY33的磷酸化修饰，激活WRKY33转录活性，诱导植物抗病响应(Tao et al.,2022)，证实了筛选出与垂柳丛枝病发生相关的候选基因是可靠的。

2、讨论

垂柳丛枝病是由植原体引起的病害，感染植原体病害后，光合作用、代谢进程、生物合成等被干扰，对各类物质的活性影响较大，蛋白质、酶、激素等活性物质发生改变，导致垂柳正常生理紊乱，造成叶片发黄、花器返祖、丛枝等现象。本研究采用转录组测序技术，对健康垂柳和感丛枝病垂柳之间的差异基因表达及相关通路进行了研究。通过分析，共检测到上调表达基因10 316个，下调表达基因8 220个。KEGG富集分析光合作用、植物-病原相互作用、苯丙氨酸、酪氨酸和色氨酸生物合成、光合作用-天线蛋白、卟啉与叶绿素代谢、亚油酸、氨基糖和核苷酸糖代谢和α-亚麻酸代谢8条通路与垂柳丛枝病发生密切相关，进一步分析筛选出K02183 (CALM)、K13424 (WRKY33)、K13448 (CML)、K13459 (RPS2)、K13457 (RPS3)、K13412 (CPK) 6个垂柳丛枝病发生相关的候选基因。

枣树感染丛枝病后，内源激素发生明显变化，其中IAA含量下降，ZT、GA、ABA含量升高(Du et al.,2013)。橡胶患丛枝病后，富集淀粉、糖类和氨基酸的途径的关键基因影响内源激素的变化(He et al.,2020)。本研究中，垂柳感丛枝病后，编码苯丙氨酸的代谢、酪氨酸和色氨酸的生物合成、酪氨酸代谢的相关基因显著表达，推测其调控患病垂柳体内的内源激素的表达，改变了垂柳体内细胞内的渗透等，导致垂柳表现出丛枝症状。这与前人研究相一致。本研究也发现，垂柳感染丛枝病后，显著富集到卟啉和叶绿素代谢通路，其中46个卟啉和叶绿素代谢差异性基因中下调的有35个，抑制了卟啉和叶绿素，导致光合作用受到限制，说明垂柳感染后光合作用受到严重的抑制。

图5 6个候选基因热图和信息   

3、材料与方法

3.1样本采集与处理

垂柳供试材料品系为W1665，于2020年7月采自甘肃省武威市林业科学研究院综合试验站(38.02'28"N,102.42'16"E)，海拔1 445 m，采样所在基地水肥管理条件良好，选择长势基本一致的植株，健康植株(LS＿CK)采取健康垂柳样本的嫩枝条和树叶，感病植株(LS＿CZ)采取有花变叶、丛枝、增殖等症状的柳树样本，用RNase free水擦拭干净后，将健康垂柳的枝条和患丛枝病垂柳的枝条慢慢剥离树皮，用单面刀片切成两组厚度一致的薄片，放入液氮中保存。

3.2 RNA提取、文库构建和RNA序列

将不同时期采集的花药送往北京诺禾致源科技股份有限公司进行RNA序列测定。RNA样品的纯度和完整性通过纳米光度计分光光度计和Agilent2100生物分析仪测定。使用1%琼脂糖凝胶电泳进一步确定RNA完整性和可能的污染。每个样品总共3μg，总RNA用于RNA序列文库制备(每个阶段有3个生物复制)。使用Illumina Hi Seq 4000测序平台进行测序，具有150个碱基对(bp)配对末端(PE)读数。FPKM (每百万碱基对测序的每千碱基转录序列片段)用于表示基因表达水平(Trapnell et al.,2010)。具有log2(折叠变化)≥1和调整后的P值<0.01被认为是差异表达的。Venn图和热图是使用Novellagic服务器1生成的。GOseq和KOBAS 2.0 R软件包分别使用了2 010和1 275个来自HR和HS的差异表达基因(DEG)来分析GO术语和KEGG途径(Young et al.,2010;Xie et al.,2011)。校正P值<0.05的GO术语和KEGG通路被认为是显著的，并被选择用于进一步分析。

3.3数据处理

FASTq格式的原始数据使用Perl脚本进行处理，以确定数据的质量，包括GC含量的百分比，Q20和Q30。过滤并消除低质量读数，以获得用于后续分析的干净读数。DESeq软件用于通过负二项分布、Benjamini和Hochberg方法计算错误发现率(FDR)以及调整后的P值来确定差异基因表达。调整后的P值<0.05被用作样本间DEG识别的阈值。

作者贡献：高静涛是本研究的实验设计者和实验研究执行人，完成数据分析、论文初稿的写作；殷稳娜参与实验设计和实验结果分析；高静涛和叶春秀是项目的构思者和负责人，指导实验设计、论文写作与修改。全体作者都阅读并同意最终的文本。

参考文献:李玲,闫旭宇,张昊栗现芳,王延峰,2020,基于高通量测序的枣疯病转录组分析分子植物育种,18(11):3537 -3543.

文章来源:高静涛,殷稳娜,叶春秀.基于转录组水平解析垂柳丛枝病发病机理[J].分子植物育种,2023,21(17):5663-5670.DOI:10.13271/j.mpb.021.005663.

基金资助:武威市市级科技计划项目(WW160226)资助;