春季低温(倒春寒)一直是黄淮及长江中下游流域冬小麦生产的主要气象灾害之一。近年来，倒春寒频繁发生[1]。从1981-2019年的28 a间，中国主要小麦产区共发生了14次倒春寒，频率高达50%。低温几乎影响植物生理生化的所有方面[2]。雌雄蕊原基分化期是小麦低温敏感时期[3]，该发育时期如遇倒春寒会直接影响小麦幼穗发育，造成穗数、穗粒数严重下降，对产量造成不可估量的损失[4-5]。随着分子生物学的发展，小麦对低温胁迫的耐受机制研究逐渐从表型[3，6-8]和生理生化水平[3，9-10]转到分子水平[11-16]。目前，小麦抗低温耐受机制的研究主要集中在越冬前的苗期，但研究表明小麦苗期抗冻性与春季生殖生长阶段抗冻性没有显著相关性[17]。频繁的倒春寒灾害迫切需要开展小麦抗倒春寒的相关机制研究。倒春寒对小麦的影响相关研究目前主要是低温胁迫下表型与生理生化机制的研究[3，8，18]，对于倒春寒胁迫响应基因的克隆、功能鉴定以及抗倒春寒基因表达调控分子机制等方面研究较少。
microRNA(miRNA)是一种典型的小非编码RNA，可指导转录后基因调控，在各种生物过程中发挥关键作用，包括发育、激素反应和应激适应[19-21]。自Sunkar等[22]首次证明了miRNA参与低温应激的调节以来，近年来许多研究也证实了植物miRNA在寒冷胁迫反应中的作用[23-26]。拟南芥中共鉴定出16种miRNA与低温胁迫有关[22，27-28]。小麦幼苗中miR159、miR164、miR169、miR319、miR398、miR1029和miR1126具有低温响应性[24]。虽然对于miRNA调控小麦抗低温研究已有报道，但对miRNA在小麦应对低温胁迫中的作用认识仍然有限，对于miRNA调控倒春寒的作用机制仍不明确。鉴于此，前期筛选的2个抗倒春寒能力有显著差异的小麦品种为材料[3]，在雌雄蕊原基分化期进行低温胁迫处理，利用RNA-seq技术，获得抗倒春寒能力不同的小麦品种幼穗miRNA表达谱，采用psRobot(v1.2)进行miRNA靶基因预测，解析差异表达miRNA在低温胁迫中的作用，为深入研究小麦抗倒春寒胁迫的分子机制提供理论基础。
1 材料和方法
1.1 试验材料
供试小麦材料为矮抗58(AK58，抗倒春寒)、郑麦366(ZM366，不抗倒春寒)，均由河南科技学院小麦中心提供。
1.2 试验方法
1.2.1 试验设计 试验于2017-2018年在河南省新乡市河南科技学院试验田进行。采用盆栽方法种植，10月10日播种，播种后将盆埋入大田，盆内土壤与盆外大田齐平。三叶期定苗，每盆留苗10株。2个品种分别设对照和低温处理2组，每组15盆。同时，另外播种5盆用于取样镜检跟踪幼穗发育进程观察。低温处理参照张自阳等[3]的方法进行，从小麦越冬起身之后，每隔2 d选取主茎和分蘖，调查AK58、ZM366幼穗发育进程，待2个小麦品种发育到颖片原基分化期后，将材料移至人工气候室中生长(气候室设置温度为白天22 ℃，夜晚13 ℃)。小麦生长发育至雌雄蕊原基分化期后，分别将低温处理的AK58(AKYJC)、ZM366(ZMYJC)移入人工气候箱中(人工气候箱白天温度0 ℃、湿度70%、时间12 h，光照强度30 000 lx；夜晚温度0 ℃，湿度70%，时间12 h，光照强度0 lx)进行低温处理72 h。以人工气候室中的AK58(AKYJCK)、ZM366(ZMYJCK)为对照。处理结束后，将低温处理的材料分为两部分，一部分移至人工气候室继续生长，用于调查结实情况。另一部分取幼穗用于miRNA测序，幼穗样品按0.5 g分装于冻存管中，迅速用液氮冷冻，-80 ℃冰箱贮存备用。结实率=低温处理结实粒数/正常对照结实粒数×100%。
1.2.2 总RNA提取、RNA样品的完整性和浓度测定 将冻存在-80 ℃的雌雄蕊原基分化期的幼穗组织取出，用TRIzol reagent(Invitrogen，CA，USA)试剂盒提取总RNA，琼脂糖电泳检测RNA样品质量，凯奥K5500分光光度计检测RNA纯度，用安捷伦2100检测样品RNA的完整性和浓度。
1.2.3 小RNA文库构建和Illumina测序 委托安诺优达基因科技(北京)有限公司完成小RNA文库构建和Illumina测序。总RNA样本检测合格后，对其进行PCR扩增建立测序文库。采用SE50测序策略[29]对检验合格的测序文库进行Illumina HiSeq高通量测序。
1.2.4 数据处理与生物信息学分析 首先将原始序列过滤得到的Clean reads通过Bowtie(1.1.2版)与中国春基因组(http：//www. ensembl. org/index. html)进行比对。根据miRBase数据库(Release 21)中的基因组注释信息，将Clean reads在基因组中的比对信息和miRNA的定位信息进行匹配(100%的位置重叠)[30]，鉴定已知miRNA，并计算RPM(Reads per million)值。利用软件DEGseq对AK58 与ZM366差异表达的miRNA进行分析，筛选条件为padj<0.05且|log2(Fold_change)|≥1[31]，找出差异表达的miRNA。采用psRobot(v1.2)[32]对差异表达的miRNA进行靶基因预测，并对靶基因进行GO(Gene ontology，http：//geneontology.org/)和KEGG功能注释分析。以P<0.05为阈值，筛选显著富集的GO term，以P<0.01筛选显著富集代谢通路(Pathway)。
1.2.5 差异表达miRNA的RT-PCR质量检测 对AK58和ZM366对照和低温处理幼穗中差异表达miRNA进行实时定量PCR。以U6为内参，利用miRNA cDNA synthesis kit(天根生化科技有限公司，北京)对样品进行反转录。反应体系包括Nuclease-free water(6.8 μL)、上游引物(0.4 μL)、下游引物(0.4 μL)、qPCR Master Mix(10 μL)、样品2 μL、miRNA第一链cDNA 0.4 μL，构成总体积为20 μL反应体系。扩增程序为95 ℃预变性40 s；95 ℃变性10 s，60 ℃退火35 s，72 ℃ 25 s，46个PCR 循环；采用默认设置自动生成Ct值。RT-PCR验证的miRNA及其引物列于表1中，每个样品3个重复。
表1 荧光定量引物列表
Tab.1 The list of primers of qRT-PCR
品种VarietymiRNA引物序列(5′-3′)Primersequences矮抗58tae-miR9672bAAGTACCACGACTGTCATTAAGCAAK58tae-miR9776TTGGACGAGGATGTGCAACtae-miR171aTGATTGAGCCGTGCCAATATtae-miR531ATCGCTCGCCGGAGCAGCtae-miR9774ATACAAGATATTGGGTATTTCTGTC郑麦366tae-miR1127aATCCTTCCGTTCGGAATTACZM366tae-miR1120b-3pATGTTCTTATATTGTGGGACAGAGtae-miR160TGCCTGGCTCCCTGTATGtae-miR531AATCGCTCGCCGGAGCAGtae-miR156ATTGACAGAAGAGAGTGAGCACA
2 结果与分析
2.1 低温胁迫对矮抗58、郑麦366穗子表型特征及结实率的影响
雌雄蕊原基低温胁迫72 h，AK58穗形受到轻微影响(图1-A，AKYJC)，而ZM366的穗形受影响较大，出现结实率极低的穗子(图1-A，ZMYJC-1)、光秆穗(图1-A，ZMYJC-2)。低温处理后，AK58和ZM366的结实率分别比对照下降12%,89%(图1-B)，说明雌蕊原基分化期低温胁迫对2个小麦品种均有影响，但对ZM366的影响较大。
2.2 矮抗58、郑麦366小麦幼穗高通量测序数据过滤分析
对低温处理72 h的小麦幼穗测序文库进行Illumina测序并进行数据处理(表2)，AK58共得到66 327 682条Raw reads，经过滤后得到45 301 441条Clean reads，ZM366共得到59 337 012条Raw reads，过滤后得到45 396 985条Clean reads。过滤掉17 bp以下的reads，大多数reads长度分布在23～25 nt(图2)，AK58、ZM366雌雄蕊原基分化期低温处理幼穗中24 nt reads分别占总reads的71.48%，72.96%，AK58、ZM366雌雄蕊分化期对照幼穗中24 nt reads分别占总reads的61.88%，68.64%。
不同小写字母表示同一品种不同处理之间差异显著(P<0.05)。
Different lowercase letters indicate significant differences between different treatments of the same variety(P<0.05).
图1 低温胁迫对2个小麦品种幼穗表型及结实率的影响
Fig.1 Effects of low temperature stress on the young ear phenotype and seed setting rate of two wheat varieties
表2 矮抗、郑麦366小麦幼穗测序原始数据过滤情况
Tab.2 Filtering of the original data of young spike sequencing of AK58 and ZM366
处理Treatment片段总数Totalrawreads总的有效序列TotalCleanreads总的有效序列占总片段比例/%CleanreadsrateQ20(Cleanreadsphred数值大于20的碱基占总体碱基的比例)/%Q20(ProportionofCleanreadsphredvaluegreaterthan20basesintotalbases)Q30(Cleanreadsphred数值大于30的碱基占总体碱基的比例)/%Q30(TheproportionofbaseswithCleanreadsphredvaluegreaterthan30tototalbases)AKYJC315917382262867465.3097.6395.05AKYJCK347359442267276765.2797.5995.42ZMYJC294522812314103078.5798.0395.58ZMYJCK298847312225595574.4797.5296.20
2.3 矮抗58、郑麦366 幼穗Clean reads与参考序列比对分析
将过滤后得到的17～30 nt Clean reads通过基因组比对分析软件与参考序列进行比对(表3)，AKYJC、AKYJCK、ZMYJC、ZMYJCK分别有16 454 955，17 878 373，17 544 621，18 919 501条序列可以完全匹配到基因组上，分别占Clean reads的80.69%，80.01%，72.66%，85.74%。
图2 不同处理下矮抗58、郑麦366 miRNA片段长度分布
Fig.2 Length distribution of AK58 and ZM366 sequenced miRNAs under different treatments
表3 矮抗58、郑麦366低温胁迫与对照幼穗中 Clean reads与参考序列比对
Tab.3 Comparison of Clean reads with reference sequence in young panicles of AK58 and ZM366 under low temperature stress and control
处理Treatmentreads总数Totalreads完全匹配的总reads数Perfectmatchreads不完全匹配比率/%Notmatchrate完全匹配的比率/%ReadsperfectmatchrateAKYJC203930781645495519.3180.69AKYJCK223459821787837319.9980.01ZMYJC228871831754462127.3472.66ZMYJCK220673651891950114.2685.74
2.4 已知miRNA鉴定
根据miRBase数据库(Release 21)中的基因组注释信息，共鉴定到112类miRNA，分属于38个不同的MIRNA家族。一般来说，miRNAs由保守和物种特异的miRNAs组成。保守miRNAs在植物发育和逆境反应中具有重要作用。以短柄草、水稻、玉米、拟南芥为参照物种进行同源搜索，识别小麦中MIRNA家族同源性，发现16个MIRNA家族在小麦中特异表达。12个MIRNA家族在以上4个物种中保守，3个MIRNA家族在以上3个物种中保守，1个MIRNA家族在以上2个物种中保守，6个MIRNA家族在1个物种中保守。此外，还发现了一些已知但不保守的MIRNA家族，如MIR9654、MIR9657、MIR9662、MIR9666、MIR3522，表明他们可能参与小麦幼穗的形成。
2.5 矮抗58、郑麦366低温胁迫与对照幼穗组织中差异表达miRNA分析
利用软件DEGseq将获得的已知miRNA进行差异表达分析，与对照相比，低温处理后AK58共鉴定出85个差异表达miRNA，ZM366共鉴定出88个差异表达miRNA(图3)。设定差异显著表达miRNA的筛选条件为padj<0.05、|log2(Fold_change)|≥1，AK58共鉴定出27个miRNA表达差异显著(图4，AK58)，其中23个显著上调表达，4个显著下调表达，差异表达极显著的miRNA为2个，分别是tae-miR9672b、tae-miR9672a-3p。ZM366共鉴定出48个差异显著表达的miRNA(图4，ZM366)，13个上调表达，35个下调表达。差异表达极显著的miRNA为4个，分别为tae-miR9672b、tae-miR9672a-3p、tae-miR6201、tae-miR9674b-5p。
图3 矮抗58、郑麦366低温胁迫与对照幼穗中miRNA差异表达火山图
Fig.3 Volcano map of miRNA differential expression in young ears of AK58 and ZM366 under low temperature stress
热图表达值是对miRNAs FPKM取log2对数值。
The heatmap expression value is log2 pair value for miRNAs FPKM.
图4 矮抗58、郑麦366低温胁迫与对照幼穗中miRNA差异显著表达热图
Fig.4 Heat map of miRNA differential expression in young ears of AK58，ZM366 under low temperature stress and control
低温胁迫后tae-miR9777、tae-miR9653b、tae-miR531在AK58和ZM366幼穗中均显著上调表达。tae-miR9666a-3p、tae-miR5062-5p、tae-miR171a、tae-miR9668-5p、tae-miR9670-3p、tae-miR9772、tae-miR9677a、tae-miR9657b-3p、tae-miR9657c-3p、tae-miR9663-5p、tae-miR9661-5p、tae-miR1131在低温胁迫的AK58幼穗中上调表达，但在低温胁迫的ZM366幼穗中下调表达。tae-miR9672b、tae-miR9672a-3p、tae-miR9774在低温胁迫的AK58幼穗中下调表达，而在低温胁迫的ZM366幼穗中上调表达。tae-miR156、tae-miR398在低温胁迫的ZM366幼穗中显著上调表达，而在低温胁迫的AK58幼穗中表达无差异或下降。与对照相比，低温胁迫后AK58有9条miRNA(图4，AK58)，ZM366有30条miRNA分别在各自幼穗中特异显著表达(图4，ZM366)。
2.6 低温胁迫下矮抗58、郑麦366差异表达miRNA靶基因预测及GO功能分析
miRNA通过与靶基因上的特定区域互补配对，结合特定蛋白质从而对靶基因的转录或者翻译过程进行调节。利用psRobot(v1.2)对低温胁迫下AK58、ZM366差异表达的miRNA进行靶基因预测，AK58幼穗中差异表达的85个miRNA对应848个靶基因，每个miRNA对应的靶基因3～144个不等；ZM366幼穗中88个差异表达miRNA对应6 537个靶基因。每个miRNA对应的靶基因3～913个不等，靶基因最多的为tae-miR531，共预测到了136个靶基因。
为了更好地理解这些差异表达显著miRNA的功能，对低温胁迫下AK58幼穗中848个靶基因，ZM366幼穗中5 376个靶基因进行GO功能富集分析，结果显示，AK58差异表达miRNA靶基因有20个生物学过程、6个细胞成分、6个分子功能类别极显著富集(FDR<0.01)；ZM366有10个生物学过程、14个细胞成分、19个分子功能类别显著富集(FDR<0.01)(图5)。可见，差异miRNA靶基因广泛涉及GO功能分类中的分子功能、细胞组分、生物学过程3大类别，AK58共涉及32个小类别，ZM366共涉及43个小类别。AK58中涉及较多的条目分别是细胞发育过程[GO：0048869](44个基因)、细胞分化[GO：0030154](42个基因)、细胞分裂[GO：0051301](31个基因)、细胞分裂[GO：0048468](27个基因)、蛋白质复合物[GO：0043234](78个基因)、转移酶复合物[GO：1990234](27个基因)、大分子络合物[GO：0032991](92个基因)、DNA结合[GO：0003677](86个基因)；ZM366涉及的较多的条目是固有免疫应答[GO：0002376](160个基因)、程序性细胞死亡[GO：0012501](131个基因)、防御反应[GO：0006952](387个基因)、细胞对胁迫的反应[GO：0033554](256个基因)、线粒体(GO：0005739)(222个基因)、细胞器腔[GO：0043233](213个基因)、细胞核部件[GO：0044428](237个基因)、小分子结合[GO：0036094](799个基因)、碳水化合物衍生物结合[GO：0097367](705个基因)、嘌呤核苷酸结合[GO：0017076](691个基因)、核糖核苷酸结合[GO：0032553](698个基因)、阴离子结合[GO：0043168](780个基因)。可以看出，低温对2个小麦品种的影响不同。低温胁迫下，AK58幼穗中差异表达miRNA靶基因大多数富集在生物学过程，而ZM366差异表达miRNA靶基因富集在分子功能，说明低温主要影响AK58生物学过程，影响ZM366幼穗抵御低温的分子功能。
进一步对2个小麦品种幼穗中差异表达miRNA进行KEGG代谢通路(Pathway)富集分析，结果显示，AK58差异表达miRNA靶基因显著富集(P<0.05)到萜类、哌啶和吡啶生物碱的生物合成，RNA聚合酶，核糖体，嘧啶代谢，嘌呤代谢，核苷酸切除修复，烟酸和烟酰胺代谢，错配修复，DNA复制，原核生物固碳途径，自噬11个Pathway(图6)；ZM366差异表达miRNA靶基因显著富集(P<0.05)在植物-病原互作，植物昼夜节律，缬氨酸、亮氨酸和异亮氨酸降解，细胞周期，苯丙氨酸代谢，苯丙烷生物合成，RNA聚合酶，植物激素信号转导，氰胺酸代谢，次生代谢产物的生物合成，光合生物中碳固定，丙酮酸代谢，类胡萝卜素生物合成13个代谢通路(图6)。
图5 低温胁迫下矮抗58(A)与郑麦366(B)GO 功能富集
Fig.5 GO analysis of miRNA target genes in AK58 and ZM366
miRNA调控mRNA参与小麦应答低温胁迫植物昼夜节律代谢途径中，miRNA调控光形态建成转导途径和开花过程的靶基因表达模式不同。在光形态建成的信号转导途径中，ZM366幼穗中tae-miR1127a、tae-miR1130b-3p作用于PRR5基因，上调了PRR5基因的表达，tae-miR1122c-3p、tae-miR1128作用于CRY1和CRY2基因，上调了CRY1和CRY2基因的表达。同时，tae-miR6197-5p、tae-miR1128、tae-miR1127a、tae-miR1120b-3p、tae-miR5200分别上调RFWD2、SPA1、HY5、CSNK2B、FT的基因表达。而在AK58中无变化。2个小麦品种miRNA参与低温应激反应phyA和phyB介导的远红光和红光信号转导途径，AK58和ZM366均上调了PRR5基因的表达，但CSNK2B基因在ZM366中上调表达，在AK58中变化不明显(图7)。
2.7 差异表达miRNA的qRT-PCR验证
为了验证测序结果，进行了qRT-PCR，分析miRNA的表达。基于高通量测序结果，使用10个代表不同测序数据表达水平miRNA(包括5个AK58的miRNA和5个ZM366的miRNA)。对RNA-seq数据与qRT-PCR数据进行相关性分析，相关系数分别为0.978 4和0.940 7，说明RNA-seq数据与qRT-PCR数据呈正相关关系(图8)。miRNA表达趋势与RNA-seq结果相似，说明miRNA测序数据是有价值的。
图6 矮抗58与郑麦366差异显著表达miRNA靶基因的KEGG富集分析
Fig.6 KEGG enrichment analysis of miRNA target gene with significantly different expression in AK58 and ZM366
3 结论与讨论
miRNA在植物对生物和非生物胁迫的适应反应中具有重要的调节功能[21-33]。miRNA作为转录后的调节因子，已被检测出单独或与多种miRNA协同对环境胁迫做出反应[21]。大量研究证明，miRNA在不同植物体内参与对低温胁迫的响应[34-35]。低温胁迫后大量的miRNA在葡萄[36]、水稻[37]、小麦[24，38]、大麦[39]、棉花[40]、烟草[41]等植物中差异表达。说明miRNA可能在低温胁迫耐受调节方面发挥着特殊作用。本研究用低温胁迫后2个抗倒春寒能力不同的小麦品种构建了4个miRNA文库，试图找到小麦参与倒春寒胁迫响应的miRNAs。结果表明，低温胁迫下2个小麦品种幼穗中均有大量miRNA差异表达，这与Song等[38]的研究结果相似，AK58鉴定出85个差异表达miRNA，27个miRNA表达出现显著差异；ZM366鉴定出88个差异表达miRNA，48个miRNA表达差异显著，说明差异表达的miRNA介导的分子途径响应了AK58和ZM366的春季低温胁迫。
同一植物品种的不同基因型对低温胁迫的反应能力也可能不同，从低温处理后AK58和ZM366的穗子表型和结实率可以看出，AK58抗倒春寒能力强于ZM366。2个小麦品种雌雄蕊原基分化期低温处理后，tae-miR9666a-3p、tae-miR5062-5p、tae-miR171a、tae-miR9668-5p、tae-miR9670-3p、tae-miR9772、tae-miR9677a、tae-miR9657b-3p、tae-miR9657c-3p、tae-miR9663-5p、tae-miR9661-5p、tae-miR113在AK58幼穗中上调表达，但在ZM366的幼穗中下调表达。tae-miR9672b、tae-miR9672a-3p、tae-miR9774在低温胁迫的AK58幼穗中下调表达，而在低温胁迫的ZM366幼穗中上调表达。说明相同miRNAs对低温胁迫的反应在同种植物不同基因型间不一致。AK58中上调和下调的miRNA可能与低温胁迫响应关系密切，可能通过调节代谢途径对低温胁迫做出反应。值得注意的是，miRNA171a在抗倒春寒能力不同的小麦品种幼穗中出现相反的表达模式。miRNA171具有高度保守性[42]。miRNA171a是miRNA171的同源miRNA，属于MIRNA171家族。在已有研究中，miRNA171参与调控番茄的株高[43]、大豆的育性[44]、小麦的抗盐胁迫[45]。miRNA171在抗低温的茶树品种中上调表达，在低温敏感的茶树品种中下调表达[46]。本研究中，tae-miRNA171a在抗倒春寒能力强的AK58中显著上调表达，而在抗倒春寒能力差的ZM366中显著下调表达，与miRNA在低温胁迫茶树中的表达特征相似[46]，提示MIR171家族成员在低温胁迫下可能具有不同的功能，miRNA171a高表达很可能增强了小麦幼穗的抗低温能力。
虚线方框内表示miRNA和靶基因上调表达，实线方框内的miRNA表示下调表达。
The dotted box indicates up-regulated miRNA and target gene expression，while miRNAs in the solid line box indicates down-regulated expression.
图7 矮抗58和郑麦366 miRNA调控的植物昼夜节律代谢途径
Fig.7 Circadian metabolic pathways regulated by miRNA in ZM366 and AK58
图8 对照和低温处理下矮抗58(A)和郑麦366(B)miRNA的qPCR验证
Fig.8 qPCR validation of miRNA the groups of contrast and low temperature treatment in AK58(A)and ZM366(B)
值得注意的是，低温胁迫下miRNA156在ZM366幼穗中显著上调表达。研究表明，miR156在植物穗部生长发育过程中发挥重要的调控作用。它通过负向调控TaSPL基因的表达水平来调节小麦穗部生长发育[47]。OsmiR156b和OsmiR156h在水稻中超表达导致开花推迟，穗分化枝数和小穗数减少[48]；小麦中超表达miRNA156会严重阻碍了穗部的生长发育，平均穗节长度和小穗数降低[49]。本研究表明，低温胁迫后，miRNA156在ZM366中显著上调表达，在AK58中表达变化不显著。同时，ZM366低温处理后出现光秆穗、空穗和结实率极低的穗。说明雌雄原基分化期低温胁迫导致miRNA156在ZM366幼穗中的大量表达，miRNA156的超表达可能导致ZM366小花发育异常，进而造成穗长和穗粒数严重下降。
miR398是ZM366幼穗中显著上调表达的miRNA。已有研究表明，miR398在植物逆境应答中发挥了重要作用[50]。miR398对编码Cu-Zn型超氧化物歧化酶(SOD)靶基因CSD剪切实现负调控。在小麦中通过抑制miRNA398表达增加CSDs水平清除体内过多的ROS[45]。miR398与其靶基因表达之间的负调控关系在其他植物中也得到了证实[34，36]。当植物遭遇低温胁迫时，为了保护自身机体，体内的SOD会升高，以便清除过量ROS来保护自己。本研究表明，雌雄蕊原基分化期0 ℃低温处理下，ZM366幼穗中miR398表达量显著升高，而AK58中miR398表达量下降。miR398的大量表达可能会造成CSD基因(编码SOD)表达量下降，SOD减少降低了ZM366的抗倒春寒能力。AK58幼穗中miR398表达量下降，负调控CSD基因产生大量的SOD来抵抗低温胁迫。
miRNAs对低温胁迫没有直接的调节作用，而是通过与靶基因mRNA进行反向互补配对介导靶基因的剪切，进而在转录水平上充当靶基因表达的调控因子。确定参与低温反应的靶基因对于揭示miRNAs的调控功能非常重要。miRNA的上调表达与靶基因表达的减少有关，反之亦然。通过对低温胁迫下AK58和ZM366差异表达的miRNA靶基因进行预测及功能分析发现，低温胁迫后AK58和ZM366幼穗中miRNA靶基因存在独有的显著富集代谢通路，AK58差异表达miRNA靶基因显著富集在萜类、哌啶和吡啶生物碱的生物合成，核糖体，嘧啶代谢，嘌呤代谢，核苷酸切除修复，烟酸和烟酰胺代谢，错配修复，DNA复制，原核生物固碳途径，自噬等代谢通路；ZM366差异表达miRNA靶基因显著富集在植物-病原互作，植物昼夜节律，缬氨酸、亮氨酸和异亮氨酸降解，细胞周期，苯丙氨酸代谢，苯丙烷生物合成，植物激素信号转导，氰胺酸代谢，次生代谢产物的生物合成，光合生物中碳固定，丙酮酸代谢，类胡萝卜素生物合成等代谢通路。靶基因富集的不同代谢通路表明，2个小麦品种在应对倒春寒胁迫时可能采用了不同的调控方式。
植物昼夜节律钟有助于耐旱和渗透胁迫，优化其功能是一种潜在的作物改良策略。当植物受到非生物胁迫时，它们的昼夜节律会发生变化[51]，使其能够更好地应对环境的变化。当昼夜节律功能异常时，植物对干旱、渗透胁迫、盐度和低温的耐受性发生改变[52]。昼夜节律通路关键节点的基因表达也会随着冷热的变化而变化[53]。昼夜节律振幅的变化和对低温的反应导致数千个基因表达的变化[54]。本研究表明，miRNA参与调控2个小麦品种应答低温胁迫，2个品种幼穗差异表达miRNA靶基因大量富集在植物昼夜节律代谢途径。此外，植物昼夜节律途径中，信号转导的关键节点基因表达模式也不同。植物中的光敏色素和隐花色素相互作用调控其光形态建成及光周期反应，进而调控植物的生长和发育过程[10，30-32]。CRY1、CRY2、RFWD2、SPA1、HY5、CSNK2B、FT是植物昼夜节律光形态建成的信号转导途径中光形态建成及光周期反应的关键节点基因，低温胁迫后CRY1、CRY2、RFWD2、SPA1、HY5、CSNK2B、FT在ZM366幼穗中均上调表达，而在AK58幼穗中变化不明显。说明低温胁迫下，ZM366幼穗中miRNA调控了光形态建成及光周期反应的关键节点基因的表达，造成光形态建成信号转导途径及开花调节出现异常，影响了小麦幼穗的生长发育，这可能是2个小麦品种抗倒春寒能力不同的原因所在。
应用高通量测序技术和生物信息学的方法对抗倒春寒能力不同的2个小麦品种幼穗的差异表达miRNA进行比较和分析，同时对差异表达miRNA的靶基因功能进行注释，发现小麦抗倒春寒通过miRNA171a、miRNA156、miRNA398等多个miRNA的差异表达及多个与抗倒春寒相关的代谢通路协同作用抵抗倒春寒。研究结果提供了倒春寒胁迫下抗倒春寒能力不同的2个小麦品种幼穗在雌雄蕊原基分化期的miRNA表达谱及差异表达的信息，为阐明小麦抗倒春寒的分子机制奠定了基础。该miRNA数据是小麦miRNA数据库的重要补充，有助于了解小麦抗春季低温胁迫的基因调控网络。
[1] 赵虹，王西成，胡卫国，曹廷杰，李博.黄淮南片麦区小麦倒春寒冻害成因及预防措施[J].河南农业科学，2014，43(8)：34-38.doi：10.3969/j.issn.1004-3268.2014.08.008.
Zhao H，Wang X C，Hu W G，Cao T J，Li B.Genetic analysis and countermeasures of wheat late-spring-coldness injury in South Huang-huai wheat region[J].Journal of Henan Agricultural Sciences，2014，43(8)：34-38.
[2] Sanghera G S，Wani S H，Hussain W，Singh N B.Engineering cold stress tolerance in crop plants[J].Current Genomics，2011，12(1)：30-43.doi：10.2174/138920211794520178.
[3] 张自阳，王智煜，王斌，王志伟，朱启迪，霍云凤，茹振钢，刘明久.春季穗分化阶段低温处理对不同小麦品种幼穗结实性及生理特性的影响[J].华北农学报，2019，34(4)：130-139.doi：10.7668/hbnxb.20190321.
Zhang Z Y，Wang Z Y，Wang B，Wang Z W，Zhu Q D，Huo Y F，Ru Z G，Liu M J.Effects of low temperature treatment at spring spike differentiation stage on young ear fruiting and physiological characteristics of different wheat varieties[J].Acta Agriculturae Boreali-Sinica，2019，34(4)：130-139.
[4] Fowler D B，Limin A E.Interactions among factors regulating phenological development and acclimation rate determine low-temperature tolerance in wheat[J].Annals of Botany，2004，94(5)：717-724.doi：10.1093/aob/mch196.
[5] Paulsen G M，Heyne E G.Grain production of winter wheat after spring freeze injury[J].Agronomy Journal，1983，75(4)：705-707.doi:10.2134/agronj 1983.00021962007500040031x.
[6] Maqbool A，Shafiq S，Lake L.Radiant frost tolerance in pulse crops-a review[J].Euphytica，2010，172(1)：1-12.doi：10.1007/s10681-009-0031-4.
[7] 高艳，唐建卫，殷贵鸿，韩玉林，黄峰，王丽娜，于海飞，李楠楠，张倩，杨光宇，李新平.倒春寒发生时期和次数对冬小麦产量性状的影响[J].麦类作物学报，2015，35(5)：687-692.doi：10.7606/j.issn.1009-1041.2015.05.17
Gao Y，Tang J W，Yin G H，Han Y L，Hung F，Wang L N，Yu H F，Li N N，Zhang Q，Yang G Y,Li X P.Effect of different periods and frequency of late spring coldness on winter wheat yield related traits[J].Journal of Triticeae Crops，2015，35(5)：687-692.
[8] 刘璇，王瑞丽，周伟，方保停，郑宏远，张艳林，詹克慧.春季低温对冬小麦穗部发育和粒重的影响[J].河南农业大学学报，2013,47(4)：373-380.doi：10.3969/j.issn.1000-2340.2013.04.001
Liu X，Wang R L，Zhou W，Fang B T，Zheng H Y，Zhang Y L，Zhan K H.Effect of spring low temperature on ear development and grain weight of winter wheat[J].Journal of Henan Agricultural University，2013，47(4)：373-380.
[9] Rinalducci S，Egidi M G，Karimzadeh G，Jazii F R，Zolla L.Proteomic analysis of a spring wheat cultivar in response to prolonged cold stress[J].Electrophoresis，2011，32(14)：1807-1818.doi：10.1002/elps.201000663.
[10] 吴青霞，杨林，邵慧，冉从福，杨子博，余静，李立群，李学军.药隔期低温胁迫对小麦生理及产量的影响[J].麦类作物学报，2013，33(4)：752-757.doi：10.7606/j.issn.1009-1041.2013.04.023
Wu Q X，Yang L，Shao H，Ran C F，Yang Z B，Yu J，Li L Q，Li X J.Effects of low temperature stress at anther connective formation phase on physiological characteristics and yield of wheat[J].Journal of Triticeae Crops，2013，33(4)：752-757.
[11] Vítámvás P，Saalbach G，Prášil I 
á V，Opatrná J，Ahmed J.WCS120 protein family and proteins soluble upon boiling in cold-acclimated winter wheat[J].Journal of Plant Physiology，2007，164(9)：1197-1207.doi：10.1016/j.jplph.2006.06.011.
[12] Zhang S J，Song G Q，Li Y L，Gao J，Wang J，Chen G J，Li H S,Li G Y, Zhao Z D.Comparative proteomic analysis of cold responsive proteins in two wheat cultivars with different tolerance to spring radiation frost[J].Frontiers of Agricultural Science and Engineering，2014，1(1)：37-45.doi：10.15302/J-FASE-2014008.
[13] Monroy A F，Dryanova A，Malette B，Oren D H，Farajalla M R，Liu W C，Danyluk J，Ubayasena L W C，Kane K，Scoles G J,Sarhan F, Gulick P J.Regulatory gene candidates and gene expression analysis of cold acclimation in winter and spring wheat[J].Plant Molecular Biology，2007，64(4)：409-423.doi：10.1007/s11103-007-9161-z.
[14] Winfield M O，Lu C G，Wilson I D，Coghill J A，Edwards K J.Plant responses to cold：Transcriptome analysis of wheat[J].Plant Biotechnology Journal 2010，8(7)：749-771.doi：10.1111/j.1467-7652.2010.00536.x.
[15] Herman E M，Rotter K，Premakumar R，Elwinger G，Bae R，Ehler-King L，Chen S X，Livingston D P.Additional freeze hardiness in wheat acquired by exposure to-3 ℃ is associated with extensive physiological，morphological，and molecular changes[J].Journal of Experimental Botany，2006，57(14)：3601-3618.doi：10.1093/jxb/erl111.
[16] Kang G Z,Li G Z,Yang W P ,Han Q X,Ma H Z,Wang Y H,Ren J P, Zhu Y J, Guo T C.Transcriptional profile of the spring freeze response in the leaves of bread wheat (Triticum aestivum L.)[J].Acta Physiologiae Plantarum，2013，35(2)：575-587.doi：10.1007/s11738-012-1099-3.
[17] Fuller M P，Fuller A M，Kaniouras S，Christophers J，Fredericks T.The freezing characteristics of wheat at ear emergence[J].European Journal of Agronomy，2007，26(4)：435-441.doi：10.1016/j.eja.2007.01.001.
[18] 王瑞霞，闫长生，张秀英，孙果忠，钱兆国，亓晓蕾，牟秋焕，肖世和.春季低温对小麦产量和光合特性的影响[J].作物学报，2018，44(2)：288-296.doi：10.3724/SP.J.1006.2018.00288.
Wang R X，Yan C S，Zhang X Y，Sun G Z，Qian Z G，Yuan X L，Mou Q H，Xiao S H.Effect of low temperature in spring on yield and photosynthetic characteristics of wheat[J].Acta Agronomica Sinica，2018，44(2)：288-296.
[19] Voinnet O.Origin，biogenesis，and activity of plant MicroRNAs[J].Cell，2009，136(4)：669-687.doi：10.1016/j.cell.2009.01.046.
[20] Budak H，Kantar M.Harnessing NGS and big data optimally：comparison of miRNA prediction from assembled versus non-assembled sequencing data：the case of the grass Aegilops tauschii complex genome[J].Omics：A Journal of Integrative Biology，2015，19(7)：407-415.doi：10.1089/omi.2015.0038.
[21] Budak H，Kantar M，Bulut R，Akpinar B A.Stress responsive miRNAs and isomiRs in cereals[J].Plant Science，2015，235：1-13.doi：10.1016/j.plantsci.2015.02.008.
[22] Sunkar R，Zhu J K.Novel and stress-regulated microRNAs and other small RNAs from Arabidopsis[J].Plant Cell，2004，16(8)：2001-2019.doi：10.1105/tpc.104.022830.
[23] Zhang J Y，Xu Y Y，Huan Q，Chong K.Deep sequencing of Brachypodium small RNAs at the global genome level identifies microRNAs involved in cold stress response[J].BMC Genomics，2009，10(1)：449.doi：10.1186/1471-2164-10-449.
[24] Tang Z H，Zhang L P，Xu C G，Yuan S H，Zhang F T，Zheng Y L，Zhao C P.Uncovering small RNA-mediated responses to cold stress in a wheat thermosensitive genic male-sterile line by deep sequencing[J].Plant Physiology，2012，159(2)：721-738.doi：10.1104/pp.112.196048.
[25] Shu Y J，Liu Y，Li W，Song L L，Zhang J，Guo C H.Genome-wide investigation of microRNAs and their targets in response to freezing stress in Medicago sativa L.，based on high-throughput sequencing[J].G3(Bethesda)，2016，6(3)：755-765.doi：10.1534/g3.115.025981.
[26] Yang C H，Li D Y，Mao D H，Liu X，Ji C J，Li X B，Zhao X F，Cheng Z K，Chen C Y，Zhu L H.Overexpression of microRNA319 impacts leaf morphogenesis and leads to enhanced cold tolerance in rice(Oryza sativa L.)[J].Plant Cell & Environment，2013，36(12)：2207-2218.doi：10.1111/pce.12130.
[27] Liu H H，Tian X，Li Y J，Wu C A，Zheng C C.Microarray-based analysis of stress-regulated microRNAs in Arabidopsis thaliana[J].RNA，2008，14(5)：836-843.doi：10.1261/rna.895308.
[28] Zhou X F，Wang G D，Sutoh K，Zhu J K，Zhang W X.Identification of cold-inducible microRNAs in plants by transcriptome analysis[J].Biochimica et Biophysica Acta，2008，1779(11)：780-788.doi：10.1016/j.bbagrm.2008.04.005.
[29] Hafner M，Landgraf P，Ludwig J，Rice A，Ojo T，Lin C，Holoch D，Lim C，Tuschl T.Identification of microRNAs and other small regulatory RNAs using cDNA library sequencing[J].Methods，2008，44(1)：3-12.doi：10.1016/j.ymeth.2007.09.009.
[30] Quinlan A R，Hall I M.BEDTools：a flexible suite of utilities for comparing genomic features[J].Bioinformatics，2010，26(6)：841-842.doi：10.1093/bioinformatics/btq033.
[31] Wang L K，Feng Z X，Wang X，Wang X W，Zhang X G.DEGseq：an R package for identifying differentially expressed genes from RNA-seq data[J].Bioinformatics，2010，26(1)：136-138.doi：10.1093/bioinformatics/btp612.
[32] Wu H J，Ma Y K，Chen T，Wang M，Wang X J.PsRobot：a web-based plant small RNA meta-analysis toolbox[J].Nucleic Acids Research，2012，40(W1)：22-28.doi：10.1093/nar/gks554.
[33] Budak H,Khan Z,Kantar M.History and current status of wheat miRNAs using next-generation sequencing and their roles in development and stress[J].Briefings in Functional Genomics，2015，14(3)：189-198.doi：10.1093/bfgp/elu021.
[34] Cao X，Wu Z，Jiang F L，Zhou R，Yang Z.Identification of chilling stress-responsive tomato microRNAs and their target genes by high-throughput sequencing and degradome analysis[J].BMC Genomics，2014，15(1)：1130.doi：10.1186/1471-2164-15-1130.
[35] Xu S C，Liu N，Mao W H，Hu Q Z，Wang G F，Gong Y M.Identification of chilling-responsive microRNAs and their targets in vegetable soybean(Glycine max L.)[J].Scientific Reports，2016，6：26619.doi：10.1038/srep26619.
[36] Sun X M，Fan G T，Su L Y，Wang W J，Liang Z C，Li S H，Xin H P.Identification of cold-inducible microRNAs in grapevine[J].Frontiers in Plant Science，2015，46：595.doi：10.3389/fpls.2015.00595.
[37] Lü D K，Bai X，Li Y，Ding X D，Ge Y，Cai H，Ji W，Wu N，Zhu Y M.Profiling of cold-stress-responsive miRNAs in rice by microarrays[J].Gene，2010，459(1-2)：39-47.doi：10.1016/j.gene.2010.03.011.
[38] Song G Q,Zhang R Z, Zhang S J, Li Y L, Gao J, Han X D,Chen M L,Wang J,Li W, Li G Y.Response of microRNAs to cold treatment in the young spikes of common wheat[J].BMC Genomics，2017，18(1)：212.doi：10.1186/s12864-017-3556-2.
[39] Kruszka K，Pacak A，Swida-Barteczka A，Nuc P，Alaba S，Wroblewska Z，Karlowski W，Jarmolowski A，Szweykowska-Kulinska Z.Transcriptionally and post-transcriptionally regulated microRNAs in heat stress response in barley[J].Journal of Experimental Botany，2014，65(20)：6123-6135.doi：10.1093/jxb/eru353.
[40] Wang Q S，Liu N，Yang X Y，Tu L L，Zhang X L.Small RNA-mediated responses to low-and high-temperature stresses in cotton[J].Scientific Reports，2016，6(1)：35558.doi：10.1038/srep35558.
[41] Chen L，Luan Y S，Zhai J M.Sp-miR396a-5p acts as a stress-responsive genes regulator by conferring tolerance to abiotic stresses and susceptibility to Phytophthora nicotianae infection in transgenic tobacco[J].Plant Cell Reports，2015，34(12)：2013-2025.doi：10.1007/s00299-015-1847-0.
[42] Llave C，Xie Z X，Kasschau K D，Carrington J C. Cleavage of Scarecrow-like mRNA targets directed by a class of Arabidopsis miRNA[J].Science，2002，297(5589)：2053-2056.doi：10.1126/science.1076311.
[43] Huang W，Peng S Y，Xian Z Q，Lin D B，Hu G J，Yang L，Ren M Z，Li Z G.Overexpression of a tomato miR171 target gene SlGRAS24 impacts multiple agronomical traits via regulating gibberellin and auxin homeostasis[J].Plant Biotechnology Journal，2017，15(4)：472-488.doi：10.1111/pbi.12646.
[44] Ding X L，Li J J，Zhang H，He T T，Han S H，Li Y W，Yang S P，Gai J Y.Identification of miRNAs and their targets by high-throughput sequencing and degradome analysis in cytoplasmic male-sterile line NJCMS1A and its maintainer NJCMS1B of soybean[J].BMC Genomics，2016，17(1)：24.doi：10.1186/s12864-015-2352-0.
[45] Wang B，Sun Y F，Song N，Wei J P，Wang X J，Feng H，Yin Z Y，Kang Z S.MicroRNAs involving in cold，wounding and salt stresses in Triticum aestivum L.[J].Plant Physiol Biochem，2014，80：90-96.doi：10.1016/j.plaphy.2014.03.020.
[46] Zhang Y，Zhu X J，Chen X，Song C N，Zou Z W，Wang Y H，Wang M L，Fang W P，Li X H.Identification and characterization of cold-responsive microRNAs in tea plant(Camellia sinensis)and their targets using high-throughput sequencing and degradome analysis[J].BMC Plant Biology，2014，14(1)：271.doi：10.1186/s12870-014-0271-x.
[47] 曹杰.小麦穗发育及过氧化氢胁迫相关miRNA的克隆与鉴定[D].北京：中国农业大学，2019.
Cao J.Characteriztion of miRNAs associated with spike development and H2O2response in wheat[D].Beijing：China Agricultural University，2019.
[48] Xie K B，Wu C Q，Xiong L Z.Genomic organization，differential expression，and interaction of SQUAMOSA promoter-Binding-like transcription factors and microRNA 156 in rice[J].Plant Physiology，2006，142(1)：280-293.doi：10.1104/pp.106.084475.
[49] Liu J，Cheng X L，Liu P，Sun J Q.miR156-targeted SBP-Box transcription factors interact with DWARF53 to regulate TEOSINTE BRANCHED1 and BARREN STALK1 expression in bread wheat[J].Plant Physiology，2017，174(3)：1931-1948.doi：10.1104/pp.17.00445.
[50] Sunkar R, Kapoor A，Zhu J K.Posttranscriptional induction of two Cu/Zn superoxide dismutase genes in Arabidopsis is mediated by downregulation of miR398 and important for oxidative stress tolerance[J].The Plant Cell，2006，18(8)：2051-2065.doi：10.1105/tpc.106.041673.
[51] Seo P J，Mas P.Stressing the role of the plant circadian clock[J].Trends in Plant Science，2015，20(4)：230-237.doi：10.1016/j.tplants.2015.01.001.
[52] Kant P, Gordon M, Kant S, Zolla G, Davydov O, Heimer Y M, Chalifa-Caspi V, Shaked R, Barak S.Functional-genomics-based identification of genes that regulate Arabidopsis responses to multiple abiotic stresses[J].Plant Cell & Environment，2008，31(6)：697-714.doi：10.1111/j.1365-3040.2008.01779.x.
[53] Gould P D, Locke J C W, Larue C, Southern M M, Davis S J, Hanano S, Moyle R, Milich R, Putterill J, Millar A J, Hall A.The molecular basis of temperature compensation in the Arabidopsis circadian clock[J].The Plant Cell，2006，18(5)：1177-1187.doi：10.1105/tpc.105.039990.
[54] Bieniawska Z，Espinoza C，Schlereth A，Sulpice R，Hincha D K，Hannah M A.Disruption of the Arabidopsis circadian clock is responsible for extensive variation in the cold-responsive transcriptome[J].Plant Physiology，2008，147(1)：263-279.doi：10.1104/pp.108.118059.
