玉米粗缩病抗性相关miRNA的筛选

吴斌，马立平，张眉，姜珊珊，郭霞，张思聪，王升吉

（1.山东省农业科学院植物保护研究所／山东省植物病毒学重点实验室，山东济南 250100；2.平度市职业中等专业学校，山东平度 266752；3.山东省农药科学研究院，山东济南 250100）

玉米粗缩病是玉米生产中的一类重要病害，1949年首次在意大利发现后［1］，相继在亚洲、南美洲和欧洲的一些国家发生［2］。1954年，首次在我国新疆和甘肃地区发现，自此在我国各玉米产区均有不同程度的发生［3］。玉米粗缩病已报到的病原有4种，分别是水稻黑条矮缩病毒（Rice black-streaked dwarf virus，RBSDV）［4，5］、南方水稻黑条矮缩病毒（Southern rice black-streaked dwarf virus，SRBSDV）［6］、玉米粗缩病毒（Maize rough dwarf virus，MRDV）［7］和里奥夸尔托病毒（Mal de Rio Cuarto virus，MRCV）［8］，山东玉米产区粗缩病主要病原为RBSDV。

近年来，随着气候条件变化、耕作制度和种植品种的改变，玉米粗缩病在黄淮海玉米产区常年发生，严重影响玉米产量和品质，制约当地玉米产业发展。除农药防治、加强田间管理和错期播种，选育和推广种植抗病品种是防治该病害的重要措施。

microRNA（miRNA）是一类内源的、21～25 nt的小分子非编码RNA，是真核生物基因表达中的一类调控因子，通过降解mRNA、翻译抑制或通过改变靶基因的染色体结构和基因组DNA的甲基化程度，在转录和翻译水平上调控基因的表达［9，10］。近年来研究发现，miRNA调控靶基因的表达在抵御病毒入侵方面起着非常重要的作用。Wu等［11，12］发现水稻条纹病毒侵染水稻后，miR168与miR528分别调控水稻体内AGO1和抗坏血酸氧化酶的表达，进而影响对水稻条纹病毒（RSV）的抗性水平。Wang等［13］发现水稻 Os-MADS23、OsMADS27a和OsMADS57对OsRDR1转录的抑制作用受miR444负调控，从而激活抗病毒RNA沉默途径。nta-miR6019和nta-miR6020通过影响NB-LRR类抗病基因受体的表达来调控对烟草花叶病毒（TMV）的抗性［14］。目前有关玉米粗缩病抗性相关miRNA研究未见报道。

本研究通过对RBSDV侵染不同抗性玉米品种的叶片和茎秆组织样品进行miRNA高通量测序，鉴定不同品种中差异表达的miRNAs，并进行靶基因预测及功能富集分析，旨在筛选出玉米抗RBSDV相关miRNA，进而为探明玉米抗粗缩病分子机理奠定基础。

1 材料与方法

1.1 试验材料

RBSDV毒源，采自山东省济宁市鱼台县玉米植株，保存于小麦叶片上；传毒介体灰飞虱，饲养于健康的小麦苗上；玉米材料为耐病品种农大108和感病品种郑单958。

1.2 试验方法

1.2.1 人工接种 RBSDV 参照张爱红等［15］的方法对不同抗性玉米品种进行RBSDV人工接种。

1.2.2 miRNA测序 待郑单958发病后，采集两个品种接种病毒和未接种病毒处理的同时期同部位的植株叶片和茎秆组织（表1），编号后液氮速冻，干冰保存条件下送至杭州联川生物技术股份有限公司进行miRNA的文库构建和高通量测序。

表1 测序样品信息

1.2.3 miRNA统计 原始测序数据首先去除3′接头序列和碱基长度小于18 nt的序列，如果序列中有80%A或C或G或T，3N（不一定是连续的），只有 A、C没有 G、T，或只有 G、T没有 A、C，或者是连续的核苷酸二聚体、三聚体，这些序列均被过滤掉。同时，将测得的序列与mRNA、RFam（包含 rRNA、tRNA、snRNA、snoRNA等）以及 Repbase数据库进行比对过滤，得到有效序列，去除重复序列后得到单一原始序列和单一有效序列，统计miRNA的长度分布和碱基偏好性。

1.2.4 miRNA注释及预测 将测序得到的miRNA序列与miRBase中的玉米及其近缘物种成熟miRNA和miRNA前体序列进行比对，整合有效序列的基因组比对结果。对于未注释到任何信息的序列，运用mireap预测新的miRNA。

1.2.5 miRNA差异表达分析 对各样品中miRNA进行表达量的统计，并对两品种相同组织在接种病毒和未接种病毒处理间进行比对，按照表达量倍数差异｜log2（fold change）｜＞2、P-value＜0.05筛选差异表达miRNAs。

1.2.6 差异表达miRNA靶基因预测及富集分析借助Target Finder对差异表达miRNA的靶基因进行预测，同时通过与Genome、KEGG pathway、Gene Ontology、KOG等生物信息数据库进行比对，进行靶基因的GO功能和KEGG通路富集分析。

2 结果与分析

2.1 miRNA测序数据统计

对8个玉米miRNA文库中测序所得的原始序列、过滤后长度大于18 nt的序列、去除重复后的miRNA序列数量进行统计分析，见表2。

表2 miRNA测序数据统计 （个）

2.2 miRNA的长度分布及碱基偏好性

8个样品中去重前后的序列长度分布情况见图1，其丰度主要集中于18～25 nt之间，在24 nt的序列丰度最大。

图1 miRNA序列的长度分布

miRNA中不同位置具有碱基偏好性，图2为不同长度的miRNA 5′端首位碱基的出现频率。结果表明，长度18～22 nt和25 nt的miRNA 5′端碱基以U居多，长度23 nt和24 nt的miRNA 5′端碱基以A居多。

图2 miRNA首位碱基偏好性

2.3 miRNA注释与新miRNA预测

将测序得到的miRNA序列与miRBase中的玉米及近缘物种基因组进行比对，共获得150个已知miRNA前体序列、174个单一miRNA序列和474个新miRNA序列（表3）。已知miRNA分属28个 miRNA家族，其中 miR166、miR169、miR171数目最多，均含有14个成员，其次是miR156，含有12个成员。

表3 miRNA注释与新预测miRNA统计

2.4 miRNA差异表达分析

对各样品中miRNA进行表达量的统计，并用TPM算法对表达量进行归一化处理，分别对两品种相同组织在接种病毒和未接种病毒处理间进行比对，分析miRNA的表达差异性（表4）。结果表明，两品种间茎秆组织中表达趋势相反的miRNA有10个，叶片组织中表达趋势相反的miRNA有17个（表5、表6）。

表4 差异表达miRNA数目统计

表5 品种间茎秆组织中差异表达miRNA

2.5 差异表达miRNA靶基因预测及分析

对差异表达的miRNA借助Target Finder进行靶基因预测，结果表明，样品中所有差异表达的miRNA共预测到1 870个靶基因，存在1 870个靶基因结合位点。

靶基因GO富集分析（图3）表明，差异表达miRNA靶基因的生物过程主要集中在生物转录、转录调控、DNA复制等；细胞组成主要涉及细胞核内过程、质膜运输、细胞质等；分子功能主要涉及DNA结合、ATP结合、蛋白结合等。

表6 品种间叶片组织中差异表达miRNA

图3 差异表达miRNA的靶序列GO富集分析

靶基因KEGG通路分析（图4）表明，差异表达miRNA的靶基因主要富集在半胱氨酸和蛋氨酸代谢、氧化磷酸化、泛素介导的蛋白水解及苯丙氨酸、酪氨酸和色氨酸的生物合成、维生素B6代谢、吞噬体等通路。

3 讨论与结论

miRNA具有高度保守性和特定的时空表达特性，在转录调控、器官发育、细胞周期调控、信号转导、抗逆生理、病原物互作等重要过程中具有重要的功能。植物病毒的侵染可诱导寄主产生大量的miRNA，通过上调或者下调等模式，抑制防卫反应中的负调控因子或促进防卫反应中的正调控因子，进而行使功能，达到提高抗性的目的［16］。

种植抗病品种是防治玉米粗缩病的重要措施。因此，发掘玉米粗缩病抗性相关基因、miRNA或ta-siRNA等，分析其功能进而探究抗病机制，可为抗病品种的遗传改良提供理论依据和技术支持，对确保玉米安全生产具有重要意义。本研究对感染RBSDV不同抗性玉米品种进行了miRNA高通量测序，筛选出与抗性相关的miRNA，叶片组织中共获得表达趋势相反的miRNA 17个，茎秆组织中10个；靶基因GO富集分析和KEGG通路分析表明，差异表达的miRNA多与转录调控、蛋白代谢等相关，分析可能在玉米抗粗缩病中起到调控作用。

图4 差异表达miRNA的靶序列KEGG富集分析

miR156存在于多种植物中，参与植物的生长发育、代谢调节及生物胁迫等多种生物过程，是植物生长发育中重要的调控miRNA［17］。研究发现，在玉米中诱导miR156表达量升高后，玉米的分蘖增加、植株矮化［18］；在水稻中，miR156通过调控12个靶序列影响对水稻齿叶矮缩病毒（Rice ragged stunt virus，RRSV）的抗性［19］。miR171、miR166及其靶序列参与生长因子、转录因子、胁迫相关蛋白等的调控，在玉米生长发育和胁迫响应等生物过程中发挥重要的作用［20-22］。miR169家族是各种非生物胁迫和发育途径的普遍调节因子，参与作物的多种抗逆［23］。下一步，将以这些miRNA为研究对象，明确与玉米粗缩病抗性分化密切相关的miRNA及其靶基因，分析两者间的调控模式，以期阐明该调控模式的作用机制。