烟草脉斑驳病毒侵染本氏烟的转录组研究

张荠丹，孙振琪，岳建英，王海娟，赵倩倩，杨向东，赵明敏，1b*

（1.内蒙古农业大学 a.园艺与植物保护学院，b.生物农药创制与资源利用内蒙古自治区高等学校重点实验室，呼和浩特 010019；2.吉林省农业科学院 a.农产品加工研究所，b.吉林省农业生物技术重点实验室，长春 130033）

烟草种植常常受到不同病毒病的影响［1-3］。马铃薯Y 病毒属（Potyvirus）病毒是目前已知的最大的植物病毒病原体类群，是所有烟草传染性病毒中破坏性最大、传播最广的病毒［4］。马铃薯Y 病毒属病毒约占植物病毒总数的30%，该属病毒寄主范围广，可通过蚜虫以非持久方式传播。因此，普通的种子处理、清洁田园、轮作和化学防治等疾病防治措施对其无效。马铃薯Y 病毒属病毒的基因组RNA 为正义单链RNA，长为10 kb。基因组两端分别为5'-非翻译区（untranslated region，UTR）和3'-UTR。它的基因组包含一个大的开放阅读框，编码一个由3 000～3 300 个氨基酸组成的多聚蛋白，然后被自身的蛋白酶加工成10 个成熟的蛋白［5］。1972 年，Gooding 等［1］报道了美国北卡罗来纳州一种发生在白肋烟上的新病毒——烟草脉斑驳病毒（tobacco vein mottling virus，TVMV），且发现该病毒在自然界中通过特定的媒介——蚜虫传播。其侵染症状表现为花叶及系统斑驳，叶脉两侧组织呈深绿色带状斑，称为“脉斑驳”。TVMV 的侵染可显著降低烟草的产量及品质，改变烟草内部化学成分［6］。病毒侵染造成的产量损失取决于多种因素，包括病毒株系的致病性强弱以及侵染时间等［7-8］。

与马铃薯Y病毒属家族的其他RNA病毒一样，TVMV 基因组最初被翻译成一个大的多聚蛋白，然后被3种病毒蛋白酶加工成10 个单个蛋白［9］。2000年，Pirone 等［10］研究发现，TVMV 的基因组由大约9 400 个碱基的单链正义RNA 组成，包括205 nt 的5'-UTR、340 kD 的多聚蛋白和250 nt 的3'-UTR［不包括poly（A）尾］。其中，多聚蛋白在蛋白水解过程中生成病毒外壳蛋白（coat protein，CP）和至少6 个非结构蛋白［11］。

高通量测序又称第二代测序技术或者深度测序（deep-sequencing），其最突出的特点是单次运行可获得的序列数据量极大［12-13］，能一次对多达几百万条的DNA 分子进行序列测定。目前，第二代测序技术己被广泛应用于生物基因组及表达研究的各个领域，例如，物种全基因组的从头测序（denovo sequencing）、RNA 测序（RNA-seq）技术等，为挖掘基因表达水平和调控网络奠定了基础。高通量测序以及新转录元的注释、差异基因的表达分析等技术的推广应用，为病毒与寄主植物之间的相互作用关系和致病机制研究提供了依据。目前，利用第二代测序技术对病毒侵染寄主后转录组进行分析的报道已有很多，其不仅在模式植物烟草或者拟南芥中进行了研究，还用于一些农作物和经济作物［14-21］。2019 年，许跃语［22］以甘蔗抗病品种ROC24 和感病品种CP72-1210 为供试材料，研究了两个甘蔗品种在黄叶病毒侵染后转录组的差异表达情况，并对抗病基因进行了进一步的分析研究。2021 年，雷阳等［23］报道了辣椒在黄瓜花叶病毒侵染下的分子响应机制，并发掘了抗病相关基因，对接种黄瓜花叶病毒的叶片进行了高通量转录组测序，利用生物信息学方法对基因表达和功能进行了研究。2023 年，苏伟华等［24］为了探究大豆抵抗大豆花叶病毒（soybean mosaic virus，SMV）侵染的分子调控机制，挖掘抗病相关基因，对侵染病毒后的大豆品种进行了转录组测序。

本研究将侵染了TVMV 的本氏烟（Nicotianabenthamiana）叶片进行了转录组测序，挖掘了其差异表达基因，对差异表达基因进行了基因本体论（gene ontology，GO）和京都基因基因组百科全书（Kyoto encyclopedia of genes and genomes，KEGG）基于途径的通路富集分析，揭示了与抗TVMV 基因相关的代谢及信号通路，为进一步开展与抗TVMV 相关基因的克隆、表达及功能研究奠定了重要的基础。

1 材料与方法

1.1 研究方法

本氏烟在温度为23℃、光照16 h 和黑暗8 h 的温室中进行培养，挑选5～6 叶期、大小一致的植物用于试验。试验所用病毒为TVMV。

采用机械摩擦接种法接种病毒。将冷冻保存的TVMV 侵染的本氏烟叶片与磷酸缓冲液（pH ＝7.2）混匀，得到病毒粗提液。每棵植株接种2片叶，接种前先撒少许二氧化硅，再将20 μL 病毒粗提液滴在接种叶片上，左手托住叶片背面，用右手食指在接种叶片上轻轻摩擦。摩擦时注意力度，使叶片表皮细胞造成微伤口。摩擦完成后在23℃条件下培养，10 d 后观察其发病情况，拍照记录发病症状。14 d 时取上部系统发病叶片，以未接种TVMV 的健康本氏烟为阴性对照，各设置3 个生物学重复，使用2 mL离心管保存于实验室超低温冰箱（-80℃）。

利用Western blot 检测病毒积累量。将烟草植株发病叶片采集后，研磨成粉末，对其叶片进行蛋白质提取，然后用Western blot 法检测病毒含量。具体步骤如下：使用包涵体蛋白提取液（含尿素）［inclusion Body Protein Buffer（with Urea），UREA］缓冲液提取蛋白质，并在12%十二烷基硫酸钠聚丙烯酰胺凝胶电泳（sodium dodecyl sulphate-polyacrylamide gel electrophoresis，SDS/PAGE）中进行电泳，然后转移至免疫印迹（polyvinylidene fluoride，PVDF）膜中。通过丽春红染色10 min，检查蛋白质的转移情况，然后通过蒸馏水洗涤将膜完全脱色。将膜放在封闭溶液中，在室温下摇床上封闭2 h。用通用型Western blot 抗体稀释液（Tris Buffer Saline，TBS）短暂冲洗膜3 次，每次10 min，然后添加一抗，在4℃过夜孵育。用TBS缓冲液洗涤3次，每次10 min后，在室温下将膜在二抗中孵育1 h。用TBS 缓冲液洗涤3 次，每次10 min，最后加入增强型化学发光试剂（enhanced chemiluminescence，ECL）显色液，在近红外双色激光和化学发光双功能成像系统上观察蛋白条带。

1.2 总RNA提取、cDNA文库构建及高通量测序

以侵染TVMV 的本氏烟叶片和健康叶片为材料，通过Trizol 法提取总RNA。先进行匀浆处理。植物组织经过液氮速冻处理后，称取大约0.1 g，并加入1 mL Trizol，将其旋涡充分混匀，室温放置5 min；然后4℃，12 000 r/min 离心10 min，取上清，弃沉淀（需在冰上操作），加入200 μL 氯仿，剧烈震荡15 s，室温放置3 min，以除去蛋白；接着4℃，12 000 r/min 离心10～15 min。取水相（约500 μL）转移至新的1.5 mL 无酶离心管中，加入等体积的异丙醇混匀，在冰上放置25 min；4℃，12 000 r/min 离心10 min，弃上清（得到胶状沉淀），然后加入1 mL 75%乙醇（需用无酶水配制）；4℃，5 000 r/min 离心3 min，倒出收集管中液体，注意不要倒出沉淀，剩余的少量液体短暂离心，然后用枪头吸出，注意不要吸到沉淀，室温放置2 min 以晾干乙醇，加50 μL无酶水用枪头反复吹打混匀，充分溶解RNA，最后保存于-80℃冰箱中。

用Oligo-dT 磁珠富集带有polyA 的mRNA。随后将完整的mRNA 进行片段化，即测序数据中的Reads。对RNA 片段进行转录组测序文库构建，使用Illumina Novaseq 6000测序平台进行高通量测序，最终获取RNA-seq数据。

1.3 主要试剂和仪器

主要试剂：异丙醇、氯仿、无水乙醇等购自园益鑫科研用品经销部；Trizol RNA 提取试剂购自天根生化科技有限公司。

主要仪器：超低温冰箱（海尔、DW-86LJ）、台式高速离心机（德国Eppendorf公司）、高速冷冻离心机（德国Eppendorf 公司5424R）、移液器（德国Eppendorf 公司）、超微量核酸蛋白浓度测定仪（中国赛默飞世尔科技有限公司Nano Drop One/One C）。

1.4 生物信息学分析

首先，使用Cutadapt 软件［25］以及本地perl 脚本去除低质量序列、污染序列以及测序仪接头序列，得到CleanData。然后，使用FastQC软件（http://www.bioinformatics.babraham.ac.uk/projects/ fastqc/）对CleanData 进行质控，使用HISAT2 序列比对软件［26］的默认参数将Reads 比对到参考基因组上（拉丁文：Nicotiana benthamiana，基因组版本：V1.0.1）；使用exomePeak 序列分析软件［27］和ChIPseeker 生信数据分析R 包［28］进行Peak calling 分析和Peak 注释。最后，使用MEME［29］和HOMER 软件［30］对富集区域进行motif 分析。基因定量软件为StringTie［31］，归一化方式为每千个碱基转录每百万映射读取的片段（fragments per kilobase of exon model per million mapped fragments，FPKM）。基因差异分析采用R 包edgeR［32］。采用的数据库信息如表1 所示。采用的数据分析软件如表2所示。

表1 数据库信息Tab.1 Database information

表2 数据分析软件Tab.2 Data analysis software

2 结果与分析

2.1 TVMV侵染本氏烟症状及其病毒检测

对生长至5～6 叶龄的本氏烟植株摩擦接种TVMV。接种一周后，对本氏烟进行症状观察，发现与健康对照相比，接种TVMV 的植株叶片呈花叶以及系统斑驳，叶脉两侧组织呈深绿色带状斑。采用TVMV 的α-CP 特异性抗体对侵染TVMV 样品进行Western blot 检测，结果显示，接种TVMV 的本氏烟叶片样品中均检测到TVMV的CP蛋白（图1），表明TVMV成功侵染本氏烟。

图1 TVMV侵染本氏烟症状及其病毒检测Fig.1 Symptom of TVMV infection in N.benthamiana and virus detection

2.2 转录组测序数据结果统计

对测序产生的原始数据（raw data）进行预处理，使用cutadapt 过滤掉不合格的序列后得到有效数据（clean data），再进行下一步分析。具体处理步骤如下：1）去除带接头（adaptor）的Reads；2）去除含有N（N 表示无法确定碱基信息）的比例大于5%的Reads；3）去除低质量Reads（质量值Q≤10的碱基数占整个Reads 的20%以上）；4）统计原始测序量、有效测序量、Q20、Q30、GC含量，并进行综合评价。

由表3 可以看出，测序序列有效Reads 的占比在健康样本1 和2 中分别为98.11%和97.76%，在侵染TVMV 样本1 和样本2 中分别为98.15%和97.96%，都在97.00%以上，表明参考数据真实可靠，可用于进一步分析。

表3 测序序列统计与质控Tab.3 Sequence statistics and quality control

2.3 参考基因组比对

使用比对软件HISAT2 或tophat 对预处理后的有效数据（valid data）进行参考基因组比对，同时根据基因组注释文件gtf 指定的基因位置信息对其分别进行统计，包括测序数据与参考基因组比对Reads 统计、Reads 在外显子和内含子比例统计以及Reads分布等信息。

2.4 差异表达基因分析

为了分析TVMV 侵染本氏烟后内源的抗病相关基因的表达差异，本研究通过转录组测序检测了TVMV 侵染的本氏烟与健康植物的基因表达情况，并绘制了火山图。根据log2（fold change）以及-lg（Pvalue）值判断出TVMV 侵染可显著影响寄主基因表达水平。结果表明，与健康本氏烟相比，TVMV侵染的样品中有4 593个差异表达基因，其中3 564个基因上调表达，1 029个基因下调表达（图2）。

图2 TVMV侵染的本氏烟样本差异表达基因火山图Fig.2 Volcano picture of differentially expressed genes in TVMV infected N.benthamiana

进一步分析与抗病相关的基因以及富集于代谢途径的主要差异表达基因（表4），共筛选出6 个上调差异表达基因，在侵染TVMV 植株中高表达。基因Niben101Scf06275g03011 的差异倍数（fold change，FC）为12.09，具有防御素蛋白相关功能；基因Niben101Scf00436g05003 的FC 为10.05，具有超氧化物歧化酶［Cu-Zn］相关功能；基因Niben101Scf09876g00011 的FC 为6.83，具有组蛋白H2A.1 相关功能；基因Niben101Scf09741g00010的FC 为3.77，具有蛋白酶抑制剂相关功能；基因Niben101Ctg11463g00002 的FC 为2.85，具有30S核糖体（ribosome）蛋白S4 相关功能；基因Niben-101Scf09268g00007 的FC 为2.73，具有富含甘氨酸的RNA 结合蛋白3 相关功能。筛选出4 个下调差异表达基因，在健康植株样本中高表达。基因Niben101Scf18347g00012 的FC 为0.46，属未知蛋白；基因Niben101Ctg06879g00004 的FC 为0.41，具有叶绿素a-b 结合蛋白相关功能；基因Niben-101Scf05692g00006 的FC 为0.28，具有PG1 蛋白相关功能；基因Niben101Scf00390g01010 的FC 为0.48，具有S-腺苷甲硫代脱羧酶蛋白酶相关功能。

表4 TVMV侵染的本氏烟样本主要差异表达基因Tab.4 Main differentially expressed genes after TVMV infection of N.benthamiana

综合上述分析可知：TVMV 的侵染极大地促进了本氏烟防御素蛋白的表达，FC 达12.09；同时，其抑制了S-腺苷甲硫代脱羧酶蛋白酶的表达，FC 达0.48（表4）。

2.5 聚类热图分析

比对差异表达基因绘制热图，结果显示：共有60 个差异表达基因在侵染TVMV 的植株样本中高表达，而在健康植株样本中低表达；有40 个差异表达基因在侵染TVMV 的植株样本中低表达，而在健康植株样本中高表达。基因Niben-101Scf05813g03001、Niben101Scf03849g00014、Niben101Scf06686g00001 在样本TVMV1 中表达水平最高，基因Niben101Scf04918g02009、Niben-101Scf02909g02010、Niben101Scf03834g00005在样本TVMV2 中表达水平最高，基因Niben-101Scf02132g00027 在样本Healthy1 中表达水平最高，基因Niben101Scf03839g05003、Niben101Ctg13152g00001、Niben101Scf07188g02011在样本Healthy2中表达水平最高（图3）。

图3 TVMV侵染的本氏烟样本差异表达基因热图Fig.3 Differential expression gene heatmap of TVMV infection of N.benthamiana

2.6 基因的GO功能注释

对差异基因进行GO 分析，结果发现，在生物学过程大类中，关于翻译所涉及的差异基因数最多最显著（257 个差异基因）。在细胞组分大类中，关于细胞内部的差异基因数最多且最显著（266 个差异基因）。在分子功能大类中，蛋白质结合（501 个差异基因）（P＞0.05）以及ATP 结合（453 个差异基因）、DNA 结合（293 个差异基因）所涉及的差异基因数最多且最显著（图4）。

图4 TVMV侵染的本氏烟样本差异表达基因GO 富集柱状图Fig.4 GO enrichment histogram of differentially expressed genes in TVMV infected N.benthamiana

上调差异基因Niben101Scf06275g03011 主要富集在防御反应过程，隶属于生物学过程大类；基因Niben101Scf00436g05003 主要富集在氧化减少过程、超氧化物歧化酶活性、超氧化物代谢过程以及金属离子结合过程，其中，氧化减少过程、超氧化物代谢过程隶属于生物学过程大类，超氧化物歧化酶活性、金属离子结合过程隶属于分子功能大类；基因Niben101Scf09876g00011 主要富集在核小体、DNA 结合、核以及蛋白质异构化活性过程，其中，核小体、核隶属于细胞组分大类，DNA结合、蛋白质异构化活性隶属于分子功能大类；基因Niben101Scf09741g00010 主要富集在丝氨酸型内肽酶抑制剂活性以及对伤害的反应过程，其中，丝氨酸型内肽酶抑制剂活性隶属于分子功能大类，对伤害的反应隶属于生物学过程大类；基因Niben101Ctg11463g00002 主要富集在RNA 结合、核糖体的结构成分、细胞内部、翻译、小核糖体亚基以及rRNA 结合过程中，其中，RNA 结合、核糖体的结构成分以及rRNA 结合隶属于分子功能大类，细胞内部、小核糖体亚基隶属于细胞组分大类，翻译隶属于生物学过程大类；基因Niben-101Scf09268g00007 主要富集在核苷酸结合以及核酸结合，其中，核苷酸结合和核酸结合都隶属于分子功能大类（表5）。

表5 TVMV侵染的本氏烟样本中上调差异基因富集分析数据统计Tab.5 Statistics of up-regulation differential gene enrichment data in TVMV infected N.benthamiana

下调差异基因Niben101Ctg06879g00004 主要富集在光合作用（photosynthesis）、收获光以及膜过程，其中，光合作用、收获光隶属于生物学过程大类，膜隶属于细胞组分大类；基因Niben-101Scf00390g01010 主要富集在腺苷甲硫代脱羧酶活性、精胺生物合成过程以及亚精胺生物合成过程，其中，腺苷甲硫代脱羧酶活性隶属于分子功能大类，精胺生物合成过程、亚精胺生物合成过程隶属于生物学过程大类（表6）。

表6 TVMV侵染的本氏烟样本中下调差异基因富集分析数据统计Tab.6 Statistics of down-regulation differential gene enrichment data in TVMV infected N.benthamiana

2.7 基因的KEGG功能注释

对差异基因进行KEGG 功能注释，可以确定差异表达基因参与的主要生化代谢途径或信号转导途径，将有助于我们更加系统地分析和破译基因的功能。以显著性富集或富集差异基因数量最多为条件，列出前20 个Pathway 进一步分析。结果如图5所示，TVMV 侵染的本氏烟植株中，差异基因富集在芪类化合物、二芳基庚烷与姜酚生物合成（stilbenoid，diarylheptanoid and gingerol biosynthesis）、淀粉和蔗糖代谢（starch and sucrose metabolism）、真核生物核糖体的生物发生（ribosome biogenesis in eukaryotes）、核糖体、嘌呤代谢（purine metabolism）、植物激素信号转导（plant hormone signal transduction）、光合作用、磷脂酰肌醇信号系统（phosphatidylinositol signaling system）、戊糖和葡萄糖醛酸相互转化（pentose and glucuronate interconversions）、核苷酸切除修复（nucleotide excision repair）、错配修复（mismatch repair）、亚油酸代谢（linoleic acid metabolism）、肌醇磷酸代谢（inositol phosphate metabolism）、同源重组（homologous recombination）、葡萄糖苷生物合成（glucosinolate biosynthe-sis）、类黄酮生物合成（flavonoid biosynthesis）、胞吞作用（endocytosis）、DNA复制、甜菜素生物合成（betalain biosynthesis）以及角质、木栓素和蜡的生物合成（cutin,suberine and wax biosynthesis）。

图5 TVMV侵染的本氏烟样本差异表达基因KEGG富集性散点图Fig.5 KEGG enrichment of differentially expressed genes in TVMV infected N.benthamiana

此外，本研究还注释到：核糖体途径的差异基因富集程度最高，差异表达基因数为307 个（Pvalue=0）；真核生物核糖体生物反应的差异表达基因数为92（Pvalue=0）；注释到“DNA 复制”途径的差异表达基因数为55（Pvalue=0）。核糖体代谢途径的富集因子（rich factor）最高，接近0.3，表明富集此通路中的基因多数差异表达。

3 讨论

近年来，关于TVMV 的研究主要集中在核内蛋白酶等方面，尚未有关于其侵染本氏烟转录组学的研究报道。本文采用Illumina Novaseq 6000 高通量测序平台对侵染了TVMV 的本氏烟进行转录组测序和分析，获得了高质量的转录组数据，为挖掘本氏烟抗TVMV 相关基因以及信号和代谢通路筛选出具有重要生物学意义的差异基因及信号和代谢通路提供了依据。本研究共筛选出10 个具有抗病功能的显著差异表达基因，并且分析了其GO 和KEGG富集通路，结果发现，蛋白质结合、ATP结合、核糖体、真核生物核糖体的生物反应、DNA 复制等通路可能在调控本氏烟抗TVMV 病毒中起着重要作用。

通过对发现的抗病相关差异表达基因进行分析，2008 年，赵朴等［35］研究发现，防御素是一类由粒性白细胞和上皮细胞产生、富含半胱氨酸的内源阳离子抗微生物肽。近年来，越来越多的证据表明，哺乳动物防御素在抗病毒免疫中不仅能直接抗病毒，而且能调节宿主免疫反应。2009 年，王猛等［36］发现，防御素具有较高的抗菌活性，并且已被证实对病毒具有杀伤作用，这些防御素对单疱疹病毒、流感病毒等衣壳病毒有较强的杀伤力。另外，防御素不仅作为抵御微生物入侵的首道屏障发挥重要作用，而且也是机体免疫反应的重要组成部分。孙勇如等［37］将兔子的防御素基因NP-1基因构建成植物表达载体质粒，将该质粒导入植物后，可以在植物体内表达，在植物体内产生防御素，从而赋予植物抗病毒的能力。1989 年，袁勤生［38］研究发现，超氧化物歧化酶（superoxide dismutase，SOD）是含有Cu、Zn 的金属酶，它能够催化超氧阴离子的歧化反应，并且其具有抗病毒的作用。2014 年，曾秀存等［39］研究发现，SOD 是一种逆境条件下清除细胞内活性氧的关键酶，而铜锌超氧化物歧化酶（Cu/Zn-SOD）是该酶系中最主要的活性氧清除剂，与植物的抗逆性密切相关。2020 年，张明达等［40］研究发现，组蛋白可以通过特殊的沉积和清除机制参与基因的表达调控，其在特定的条件下可通过调控细胞衰老和代谢抑制肿瘤增殖。2003 年，曲晓华等［41］研究发现，蛋白酶抑制剂可以抑制昆虫的生长和发育。有研究表明，蛋白酶抑制剂与植物的抗病性也有一定的关系。例如：2003 年，杨静等［42］发现，蛋白酶抑制剂泛指具有蛋白酶抑制作用的一类物质，它既可以是大分子的蛋白质、寡核苷酸，又可以是小分子的有机物，目前其应用于抗炎、抗感染、抗病毒感染、抗肿瘤等方面；2020 年，马宁等［43］克隆了一个番茄蛋白酶抑制剂基因，研究发现，其在抗土传病害和叶部病害方面具有重要作用；2010 年，卢秀萍等［44］克隆烟草富含甘氨酸RNA 结合蛋白基因，并进行原核表达，为制备抗体和研究烟草抗逆性分子机理奠定了基础；2020 年，吴廷全等［45］首次公开了黄瓜捕光叶绿素a/b 结合蛋白在抗瓜类疫病中的应用，特别是对侵染瓜类的瓜类疫霉具有抗病功能。

本研究通过对侵染TVMV 的本氏烟转录组分析揭示了与其相关联的基因和信号通路，为今后研究主要关联基因及信号通路的功能，阐明抗TVMV机理奠定了重要的基础。