涉及Calyx持久性的Korla芳香梨花菌脱落区mRNA-SEQ和miRNA-SEQ的综合分析

背景

本研究旨在通过转录组和小RNA测序研究库尔勒香梨花萼持久性的分子机制。用50 mg/L赤霉素对花期后期的第1、5和9天花的脱落带组织进行了测定_3.（calyx persistence treatment, C_1, C_5, C_9) or 500 mg/L PP_333.（Calyx脱落处理，收集T_1，T_5，T_9），并同时进行转录组和小RNA测序。

结果

通过转录组与小RNA测序的关联分析，进行mRNA-miRNA网络。与同期花萼脱落组比较，C_1与T_1组、C_5与T_5组、C_9与T_9组分别获得145、56和150对mRNA-miRNA;C_1与C_5和C_9比较，分别筛选出90和506对mRNA-miRNA, C_5与C_9比较得到255对mRNA-miRNA;T_1与T_5、T_9比较分别得到206、796对mRNA-miRNA, T_5与T_9比较得到383对mRNA-miRNA。这些mrna在miRNA-mRNA对中显著富集到萜类主干生物合成、光合-天线蛋白、卟啉和叶绿素代谢、类胡萝卜素生物合成、玉米素生物合成和植物激素信号转导中。此外，我们从miRNA-mRNA对中获得了一些可能与花萼脱落相关的关键基因，包括蛋白磷酸酶2C (pci - mir394a -HAB1)、受体样蛋白激酶(psi-miR396a-5p-Herk1.），纤维素合酶 - 样蛋白D3（PSI-miR827-CSLD3），β-​​半乳糖苷酶（psi-mir858b-β-半乳糖苷酶），SPL8 ' -羟化酶1 (psi-miR396a-5p-CYP707A1）和生长素响应因子（Psi-mir160a-3p-ARF6psi-miR167d -ARF18psi-miR167a-5p -ARF25）， 等等。

结论

通过对mRNA和miRNA的综合分析，进一步了解库尔勒香梨花萼脱落相关的重要基因及其调控途径。我们还建立了miRNA- mrna对网络，以了解miRNA对花萼脱落的精确调控。

库尔勒香梨(Pyrus sinkiangensis玉），是新疆维吾尔族自治区，中国一个古老的地方优质品种。在库尔勒香梨，花萼是一些花落叶但在其他持久性。一个宿存花萼是畸形果的库尔勒香梨的主要原因。这可以梨形和质量产生负面影响，并直接影响到库尔勒香梨[经济回报1］．

研究人员对库尔勒香梨花萼持久性与树木形态的关系进行了研究[2]，植物生长调节剂〔3.，4]、砧木类型[5]，授粉[6],修剪[7[光明[8］．Ma等观察了PP喷施后花萼发育阶段花萼管显微结构的差异_333.或遗传算法_3.库尔勒香梨的花期。经GA处理后，花萼管组织维管束平均面积增大，逐渐出现大量筛管细胞和异母细胞_3.为花萼提供营养和水分，阻止花萼筒脱落层的形成，最终形成宿存花萼果。用聚丙烯处理树木_333.，幼果花萼管中只有维管束。花萼筒发育后期幼果出现脱落层，最终花萼筒脱落，幼果成为无花萼果[9］．花萼持久性与GA密切相关_3.梨花萼含量。高的遗传算法_3.花萼中GA含量低则有利于花萼的存留_3.花萼含量有利于花萼脱落[10.，11.，12.］．前人研究表明，在同一花序中，随着花萼排列顺序的增加，花萼脱落能力增强。第1位至第4位花萼存留率逐渐降低，花萼脱落率逐渐增加[2］．花萼持久性也与树木活力有关。树势旺盛的梨树萼存果较多，树势较弱的梨树萼存果较少[13.］．此外，一些研究还分析了萼持久性的分子机制。据建议my例如基因(kfpMYB)涉及花萼持久性[14.，15.，16.］．Sun等。确定了三个SPL通过高通量测序与萼片发育相关的基因[17.］．Tian等人分离并克隆了果胶酸裂解酶PSPL.基因在库尔勒香梨[18.］．Qi等通过数字转录丰度测量，鉴定了库尔勒香梨花萼脱落带中参与花萼脱落过程的候选基因[19.］．Pei等人使用RNA-seq和数字基因表达(DGE)技术鉴定了50多个与花萼持久性相关的基因[20.］．周等人。使用小RNA测序确定了几种涉及Calyx持久性的miRNA [21］．大部分研究对库尔勒香梨花萼和子房的脱落带进行了研究，而对其脱落带的研究较少。此外，这些关于花萼脱落分子机制的研究是通过转录组和小RNA测序的分离应用进行的，但结合这两种组学的关联分析尚未见报道。

microrna (mirna)是一种短的非编码rna，长度约为19-24个核苷酸，与较长的靶标rna具有可变的序列互补[22］．然而，目前miRNA靶标的确定仍然是一个重大的挑战，因为miRNA的靶标序列可能不是保守的[23］．而mirna与其靶基因之间的相互作用可能是一对多或多对一，而不是严格的一对一，从而产生大量潜在的调控效应[24］．此外，估计基于生物信息学软件预测的靶mRNA的功能的miRNA的调节功能是不准确的[25］．有一种可能更可靠的方法来预测miRNA的靶基因:首先,小RNA序列用于识别和预测的目标基因差异表达microrna,然后RNA-seq技术用于筛选差异表达mrna,最后重叠目标基因和基因之间的差异表达mrna选择确定microrna的监管关系和mrna (26］．

本研究以库尔勒香梨花萼脱落带为研究对象，利用mRNA-seq和miRNA-seq分析花萼保留组和花萼脱落组之间的miRNA和mRNA表达谱，以阐明花萼保留的分子机制。此外，我们将差异表达的mrna和差异表达的mirna结合在一起，确定了库尔勒香梨花萼脱落区参与花萼保留的关键miRNA-mRNA网络和通路。这是第一次对库尔勒香梨的mRNA-seq和miRNA-seq进行综合分析，从而对花萼持久性的分子机制有了更深入的了解。

植物材料和治疗

在这项研究中所使用的植物材料在春天2018年Shayidong园艺场，库尔勒，新疆获得。Two uniform twenty-year-old Korla fragrant pear trees were selected and treated with either 50 mg/L GA_3.(花萼持久性处理)或500 mg/L PP_333.（Calyx脱落处理）在盛开的阶段。在晚绽放阶段的第一个，第五和第九天收集鲜花。在树丛中开放的第一朵花有一个持久的花萼。四棵树在树木中开放的花朵有一个落叶花萼。因此，GA后收集群集上的第一朵花_3.治疗，PP后群体上的第四朵花_333.治疗。采集后，摘除花瓣和萼片，人工从花萼管标本中分离出花萼离断带组织，其中包括花萼管和离断线近端几层离断带细胞及相邻细胞。用GA采集花期后期的第1、5、9天的花_3.分别称为C_1, C_5, C_9。用PP法采集花期后期第1、5、9天的花_333.治疗分别称为T_1，T_5，T_9。将Calyx脱落区（AZ）组织立即在液氮中冷冻并储存在-80℃直至使用。含有花萼持久性和花萼脱落的Korla芳香梨在附加档案中显示1．花萼AZ组织样品示于附加文件2．

RNA提取和质量评估

使用EASYspin Plant microRNA Kit (Aidlab, Beijing, China)按照制造商的说明从花萼AZ组织中提取总RNA。将同一组6个不同AZ组织的总RNA样本按等量RNA进行汇总。在1%琼脂糖凝胶上监测RNA降解和污染。使用纳米光度计®分光光度计(IMPLEN, CA, USA)检测RNA纯度。使用Qubit®RNA Assay Kit在Qubit®2.0荧光计(Life Technologies, CA, USA)中测量RNA浓度。使用安捷伦生物分析仪2100系统(美国安捷伦技术公司，CA，美国)的RNA Nano 6000检测试剂盒对RNA完整性进行进一步评估。

转录组测序

使用mRNA-seq的技术库尔勒香梨花的AZ组织的mRNA（成绩单）进行了分析。For six cDNA library constructions, a total amount of 1.5 μg RNA per group was used as input material for the RNA sample preparations. The library for sequencing was generated using a NEBNext® Ultra™ RNA Library Prep Kit for Illumina® (NEB, USA) and index codes were added to attribute sequences to each sample. Transcriptome sequencing was carried out on an Illumina HiSeq 2500 platform and paired-end raw reads were generated. After removing reads containing adapter, reads containing ploy-N, and low quality reads, the remaining clean reads were assembled using Trinity softwareas described for de novo transcriptome assembly without a reference genome [27］．Illumina测序在中国北京诺二烯序列进行。

要理解基因的功能，基于以下数据库的基因函数：Swissprot（手动注释和审查的蛋白质序列数据库）;NR（NCBI非冗余蛋白序列）;NT（NCBI非冗余核苷酸序列）;PFAM（蛋白质）;KO（kegg ortholog数据库）;去（基因本体论）;和Kog / Cog（蛋白质蛋白质群/真核癌的簇）。

差异表达基因(DEGs)分析

每个转录本的表达水平被测量为映射到其参考序列的干净读的数量。用RSEM将映射的干净读数归一化为RPKM(每千万个碱基对测序的转录本序列每千碱基的预期片段数)[28］．采用DEGSeq R包对两组进行差异表达分析。Q值被用来调节p值(29］．设定q≤0.005 & |log2 (fold change)| > 1为DEGs选择的阈值。DEGs进一步应用于KEGG (Kyoto Encyclopedia of Genes and Genomes)通路。我们使用KOBAS [30.]软件检测KEGG通路中DEGs的统计学富集程度。

小分子rna序列

A total amount of 3 μg total RNA per group was used as input material for the small RNA library. Sequencing libraries were generated using NEBNext® Multiplex Small RNA Library Prep Set for Illumina® (NEB, USA.) following manufacturer’s recommended protocol and index codes were added to attribute sequences to each sample. Then the libraries were sequenced by Illumina Hiseq 2500/2000 platform and 50 bp single-end reads were generated. Raw data of fastq format were firstly processed through custom perl and python scripts. Clean data were obtained by removing reads containing ploy-N, with 5′ adapter contaminants, without 3′ adapter or the insert tag, containing ploy A or T or G or C and low quality reads from raw data. Then sRNAs were mapped to thePyrus sinkiangensisYu transtripome使用领结法[31没有不匹配)。来自蛋白质编码基因、重复序列、rRNA、tRNA、snRNA和snoRNA的标签，小RNA标签被映射到RepeatMasker。从干净的数据中删除了rafm数据库。

随后，通过改进的mirdeep2软件将长度为18-30 nt的独特序列定位到miRBase 22.0中的特定物种前体[32]和srna-tools-clito（http://srna-tools.cmp.uea.ac.uk/)来识别已知的microrna。可用的软件miREvo [33]和mirdeep2结合，通过探索前面步骤中未注释的小RNA标签的二级结构、Dicer裂解位点和最小自由能来预测新的miRNAs。

差异表达mirna (dem)分析

通过TPM将mirna的表达水平归一化[34:标准化表达式=映射的readcount/总读数× 1,000,000。采用DEGSeq (2010) R软件包对两组之间的dem进行分析[35］．本试验设置显著性阈值为q≤0.01 & |log2 (fold change)| > 1。

目标预测和的DEM的功能注释

利用psRNATarget (http://plantgrn.noble.org/psRNATarget/?function = 3)(最大期望值= 5.0)。为便于描述，将“dem的靶基因”简称为“候选靶基因”。

DEGs和dem的综合分析

为了确定所有可能的miRNA-mRNA相互作用，包括miRNA与mRNA表达之间的正或负关系，我们使用Cytoscape 3.2.0构建miRNA-mRNA调控网络。miRNA-seq与mRNA-seq的整合是通过整合候选靶基因和DEGs完成的。为便于描述，将“候选靶基因”与“DEG”的交集称为“差异靶基因”。KEGG富集分析有助于阐明差异靶基因的生物学功能和关键信号通路。

实时定量PCR (qRT-PCR)验证

为验证测序数据，采用qRT-PCR检测各组间dem和DEGs的表达模式。验证了8个成熟dem和10个DEGs在miRNA-mRNA互作网络中的表达谱。总RNA样本与小RNA测序样本相同。然后将总RNA聚腺苷化，使用Mir-X miRNA First-Strand Synthesis Kit (Takara, Beijing, China)提取2 μL总RNA生成cDNA。利用primer Premier 5.0设计基因特异性引物、mirna特异性正向引物和通用反向引物。采用美国Bio-Rad cfx96 Touch RT-PCR检测系统，SYBR Green Real-time PCR Master Mix (Toyobo, Osaka, Japan)进行QRT-PCR。mRNAs和miRNA qRT-PCR方法参照前人研究[20.，21］．U6和Actin分别作为miRNAs和mRNAs的内参对照。表达量以2计算^-ΔΔct方法 [36］．计算各miRNA/mRNA表达水平和mRNA- seq /miRNA- seq数据的Pearson相关系数。额外的文件3.显示上述引物的序列。

转录组测序分析

表示萼持久性基团（C_1，C_5，C_9）和那些花萼脱落基团（T_1，T_5，T_9）在库尔勒香梨的花萼AZ组织六cDNA文库用总RNA构建并进行Illumina测序。测序的概观和装配结果花萼持久性基团和萼持久组示于表1．去除低质量的原料后读取，308085228清洁读取依然存在。Through the Trinity de novo assembly method, 158,918 non-redundant genes were obtained, and 180,542 transcripts were achieved with an N50 of 1835 bp and an N90 of 566 bp (Table2)．The length distribution of genes and transcripts larger than 300 bp are shown in Fig.1．

功能辅助和分类

158,918个组装的基因在NR、NT、KO、SwissPort、PFAM、GO和KOG数据库中均有显著匹配，基因数量分别为126,102(79.35%)、144,025(90.62%)、49,033(30.85%)、92,788(58.38%)、85,811(53.99%)、85,952(54.08%)和34650(21.8%)(图1)。2)．

在氧化石墨烯生物过程分析中，将85952个基因分为55个亚类(图1)。3.)．涉及“细胞过程”（49212），“单生物体过程”（36400）和“代谢过程”的基因（45402）基团是高度在生物过程类别表示。在这些细胞成分，“细胞”（26268）是最常见的代表，其次是“细胞部分”（26265）和“细胞器”（17538）。在分子功能类别，聚类的显著比例被分配到“结合”（51256）和“催化活性”（39876）。

为了对外贸基因产品进行分类，34,650个基因被细分为25个Kog分类。其中，“后期修饰，蛋白质周转，伴侣”（4944）的集群代表了最大的基团，其次是“一般函数预测”（4480），“翻译，核糖体结构和生物发生”（3143）和“信号”转导机制“（2644）。“细胞外结构”（38）和“细胞运动”（19）是最小的组（图。4)．

KEGG分析发现，49,033个基因分为19个途径类别。最具显着富集的途径是“碳水化合物代谢”（4226），“翻译”（4102），“折叠，分选和降解”（3682），“概述”（3163）和“氨基酸代谢”（2703）（图．5)．

花萼保留组与脱落组之间的差异

DEGs分析显示，保留花萼组与脱落花萼组的AZ组织差异有统计学意义。将花萼保留组与花萼脱落组、同一时期的2个处理、不同时期的相同处理进行比较，共进行9对比较，包括(i) C_1 vs T_1;(ii) C_5 vs T_5;(iii) C_9 vs T_9;(iv) C_1 vs C_5;(v) C_5 vs C_9;(vi) C_1 vs C_9;(vii) T_1 vs T_5;(viii) T_5 vs T_9;和(ix) T_1 vs T_9。 When calyx persistence groups compared with calyx abscission groups during the same stage, 2124, 1280 and 1845 DGEs were identified in C_1 vs T_1, C_5 vs T_5 and C_9 vs T_9, respectively; When C_1 compared with C_5 and C_9, 1835 and 6366 DEGs were screened respectively, and 3991 DEGs were obtained from the comparison between C_5 and C_9; When T_1 compared with the T_5 and T_9, respectively, 2549 and 7577 DEGs were obtained, and 4386 DEGs were obtained from the comparison between T_5 and T_9 (Fig.6一个)。

二氢酶的途径富集分析

为了进一步研究DEGs在花萼脱落过程中的作用，我们分析了明显富集的KEGG通路。以保留花萼组为对照，分别与脱落花萼组比较筛选出的DEGs主要富集于α -亚麻酸代谢、植物激素信号转导和光合作用。C_1 vs C_5、C_1 vs C_9、C_5 vs C_9中DEGs主要富集于植物激素信号转导，同时富集于半乳糖代谢和类胡萝卜素生物合成。T_1对T_5、T_1对T_9、T_5对T_9的DEGs主要富集于核糖体、光合作用和半胱氨酸、蛋氨酸代谢(图)。7)．这些结果意味着这些途径中涉及的基因可能在Calyx脱落中起到至关重要的作用。

miRNA测序分析

六个图书馆共获得95,888,963份原始阅读资料。去除低质量读和适配器后，C_1、C_5、C_9、T_1、T_5、T_9样本中分别得到14,816,133、13,952,977、14,223,817、14,122,942、19,289,704、14,150,546个clean reads。小RNA分类注释结果统计概况见表3.．选取长度为18 - 30 nt的Clean reads进行进一步分析。数字8显示了6个文库中srna的长度分布。干净读的长度分布显示，大部分读的长度在23 ~ 25 nt之间，读计数以24 nt最多。

保守和新颖mirna的鉴定

总的来说，我们发现在6个小RNA文库属于24个miRNA家族48倍保守的miRNA和84倍预测的新的miRNA（附加文件4:表S2)。关于保守miRNA家族成员数目的详细信息在附加文件中进行了总结5：表S3。共有18个保守的miRNA家庭包含多个成员。

花萼保留组与花萼脱落组之间的组织学差异

在mirna测序研究中，我们分别在C_1和T_1、C_5与T_5、C_9与T_9中获得了49个、30个和49个差异表达的mirna。在C_1与C_5、C_1与C_9、C_5与C_9中分别鉴定出37、57和43个dem。T_1对T_5、T_1对T_9、T_5对T_9的dem分别为51、69、62(图4)。6b).这一趋势与deg的数量相似。

候选靶基因的途径分析

为了更好地了解所鉴定的dem的功能，我们预测了这些mirna的潜在靶基因。预测了20744个靶基因(附加文件)6:表S4)。通过KEGG富集分析，我们发现有45条KEGG通路显著(q值≤0.05)与miRNAs靶向基因相关(Additional file)7:表S5)。碳水化合物代谢和氨基酸代谢涉及多种途径。

DEGs和dem的综合分析

为了在AZ组织中探索miRNA和mRNA调节网络，将miRNA和mRNA的表达谱组合用于进一步分析。6个处理组共2587对miRNA-mRNA，存在正相关和负相关。总共通过综合分析获得1206个差异靶基因和114 miRNA。与同期花萼脱落组比较，C_1与T_1组、C_5与T_5组、C_9与T_9组分别获得145、56和150对mRNA-miRNA;C_1与C_5和C_9比较，分别筛选出90和506对mRNA-miRNA, C_5与C_9比较得到255对mRNA-miRNA;当T_1与T_5和T_9相比，分别获得206和796 mRNA-miRNA对，并且从T_5和T_9之间的比较获得383 mRNA-miRNA对（附加文件8)．大多数mRNA-miRNA对呈现出负相关的表达模式。大多数miRNA具有多于一种可能的靶基因，而不同的miRNA可以调节相同的目标。例如，在C_1 VS T_1中，PSI-MIR827是Cluster-9706.71346，Cluster-9706.109338和Cluster-9706.121318的调节器，而Psi-MiR171A-3P和PSI-MiR171B-3P可以调节Cluster-9706.61759的表达。

差异靶基因KEGG富集分析

为了更好地了解这些差异靶基因的生物学功能，我们进行了KEGG通路分析。C_1和T_1的差异靶基因主要富集于光合-天线蛋白、半胱氨酸和蛋氨酸代谢以及苯丙氨酸、酪氨酸和色氨酸生物合成;C_5和T_5的差异靶基因主要富集于植物激素信号转导、卟啉和叶绿素代谢、半胱氨酸和蛋氨酸代谢;在C_9和T_9中，除了光合生物的α -亚麻酸代谢和碳固定外，差异靶基因主要富集于戊糖磷酸途径。C_1与C_5、C_1与C_9、C_5与C_9的差异靶基因主要富集于光合-天线蛋白、戊糖磷酸途径、卟啉和叶绿素代谢。T_1对T_5、T_1对T_9、T_5对T_9的差异靶基因主要富集于光合-天线蛋白、类胡萝卜素生物合成和戊糖磷酸途径(图)。9额外的文件9)．

为了进一步说明参与花萼脱落的可能机制，的miRNA表达分析缩小到关键基因，并通过KEGG分析鉴定显著丰富的途径。结合Degs的Kegg浓缩分析，我们得出结论，一些重要的代谢途径可能与Calyx脱落，包括萜类骨架生物合成，光合作用 - 天线蛋白，光合作用，卟啉和叶绿素代谢，类胡萝卜素生物合成，Zeepin生物合成和植物激素信号转导。．此外，我们从MiRNA-mRNA对中获得了一些关键基因，以根据它们的注释及其与脱落响应MiRNA的潜在关系，包括蛋白磷酸酶2C（Psi-MiR394a-HAB1），纤维素合酶 - 样蛋白D3（PSI-miR827-CSLD3)、β-半乳糖苷酶(psi-miR858b-β-半乳糖苷酶)、锌指蛋白(psi-miR319a-ZAT12），脱落酸8'-羟化酶1（psi-mir396a-5p-CYP707A1), Laccase-7 (psi-miR397a -LAC7），聚半乳糖酶（Psi-mir396b-3p-GSVIVT00026920001)、受体样蛋白激酶(psi-miR396a-5p-Herk1.）和生长素响应因子（Psi-mir160a-3p-ARF6psi-miR167d -ARF18psi-miR167a-5p -ARF25）， 等等。

通过qRT-PCR确认DEGs和dem

8米的DEM（PSI-miR394a，PSI-miR858b，PSI-miR397a，PSI-miR396a-5P，PSI-miR156j，PSI-miR160a-3P，PSI-miR167d，PSI-miR167a-5P）和10度的视角的表达谱{HAB1（簇-9706.44181），β-半乳糖苷酶（簇-9706.74090），Herk1.(集群- 9706 - 30698),LAC7(集群- 9706.1522),ZAT12(集群- 9706.14884),SPL13A(集群- 9706.35292),ARF6(集群- 9706.66959),ARF18(集群- 9706.105236),ARF25（Cluster-9706.68003），CYP707A1（使用QRT-PCR进一步验证mRNA-miRNA相互作用网络中的簇-9706.109952）}（图。10.)．qRT-PCR结果显示，这些mRNA / mirna大部分与mRNA- seq /miRNA-Seq数据具有相似的表达趋势，Pearson相关性也显示，大多数mirna /mRNA表达水平与mRNA- seq /miRNA-Seq数据具有很强的相关性。这可能在一定程度上验证了我们序列数据和本研究结果的可靠性。

在本研究中，我们同时分析mRNA和miRNA谱并构建miRNA-mRNA监管网络，以提高我们对Calyx持续性分子机制的理解。这是关于Calla脱脂梨的Calyx脱落的平行mRNA和miRNA表达的第一个详细信息。通过综合分析，我们获得了与Calyx脱落，它们的相互作用和监管途径相关的一组完整的miRNA / mRNA。

库尔勒香梨花发育对花萼持久性的影响

无论是在C_1 / C_5，C_5 / C_9，C_1 / C_9还是在T_1 / T_5，T_5 / T_9，T_1 / T_9中，逐渐增加了DEG的数量。他。发现最多果实的Calyx是否落后或不在盛开后的8天内掉落，即观察到Calyx管脱落症状（在脱落区的黄色环）（附加文件1：图S1）[37］．在图中。6， T1/T9和T5/T9之间的DEGs数量远高于T1/T5之间的DEGs数量，可以认为，大量与花萼脱落相关的基因在PP处理后花期后期的5-9天表达_333.．

Calyx脱落有关的差异靶基因的途径分析

反式玉米蛋白核苷(ZR)是一种天然的细胞分裂素。植物内源激素如吲哚-3-乙酸(IAA)、脱落酸(ABA)、赤霉素(GA)等_3.Zr参与梨中花萼持久性过程的调节。高含量的IAA，GA_3.而花萼中的ZR和ABA含量低是梨果实经久不衰的原因之一[12.］．Terpenoids是植物中非常重要的二级代谢物。一些三萜类化合物，如胃肠杆菌素，吲哚脲和其他植物激素，在植物生长和发育中发挥着重要作用。类胡萝卜素和叶绿素参与植物光合作用[38］．在植物中，类胡萝卜素主要存在于花和果实的叶绿体和有色体中。它参与光合作用中的光吸收和光形态建成[39，40，也参与植物对外界刺激的反应[41］．类胡萝卜素也可用作植物激素（例如脱钙酸）的前体，以参与其合成代谢过程[42］．光合作用的减少在花和果实的脱落过程中起着重要的作用[43］．我们的研究发现光合作用 - 天线蛋白和光合作用途径被显著丰富。有关光合作用例如卟啉和叶绿素代谢和类胡萝卜素生物合成途径等也显著丰富，表明光合作用起着花萼脱落的过程中起重要作用。植物激素能调节植物生长和发育，以及各种生理代谢过程。几乎所有的激素可以影响器官脱落，如乙烯，植物生长素，脱落酸，赤霉素和茉莉酸。有研究证实，乙烯和生长素是占主导地位的激素调节器官的衰老和脱落[44］．在本研究中，玉米素生物合成、萜类主干生物合成、类胡萝卜素生物合成、光合作用和植物激素信号转导途径显著富集。提示这些途径在花萼脱落中起重要作用。

miRNA-mRNA对的功能分析

我们通过组合小RNA测序和转录数据获得1206个差动靶基因和114种miRNA。这些差分靶基因代表与靶基因组差异表达花萼脱落/持久治疗后AZ组织通过的DEM调节。该基因组是基于miRNA和mRNA表达之间的调控性关系构造，反映在miRNA的表达所造成的变化花萼脱落/持久性治疗后的基因的异常改变。因此，这些的miRNA表达对非常有可能参与花萼脱落。在这里，我们专注于已知miRNA及其相应的差的靶基因。我们通过整合分析获得517个差动靶基因和41倍已知的miRNA（附加文件10.)．

漆酶（LAC)基因家族已被报道参与木质素生物合成。在拟南芥蒂利亚纳，类黄酮的氧化聚合和木质素的生物合成已被证明是由漆酶-15催化的[45］．木质素是一种重要的大分子有机材料。木质素的含量仅比该内植物纤维素的低级[46］．木质素和纤维素，半纤维素是植物骨架的主要成分。与地层和木质素的降解的基因已经被报道参与器官脱落[47］．4-香豆酸-辅酶a连接酶在木质素的合成中起重要作用[48］．综合分析发现漆酶-7 (LAC-7psi-miR397a在C_1 vs C_9、C_9 vs T_9、T_5 vs T_9、T_1 vs T_9中表达差异显著。在Populus Trichocarpa.miR397a -LAC对参与木质素的形成[49］．QRT-PCR也提示miR397a和LAC7在C_9和T_9中表现出相反的结果，并且有显著的差异表达(图。10.)．4-香豆酸-辅酶a连接酶样5 (4CLL5.)被psi-miR396b-5p和psi-miR396a-5p靶向。根据之前的研究，miR396参与花、叶和果实的发育，因此psi-miR396可能起到调节作用4CLL5.促进木质素合成的增加抑制了Calyx脱落[37］．纤维素酶是细胞壁水解酶，在器官脱落过程中起重要作用。它能降解纤维素和半纤维素，是细胞壁的主要成分[50］．纤维素合酶样蛋白D3（CSLD3）预计预计是由PSI-MiR827的目标。PSI-MIR827的表达水平在C_9 VS T_9中调节，而CSLD3被下调，这表明的下调CSLD3在调控持续的花萼重要的作用。E3泛素 - 蛋白连接酶（桶)已被证明在植物激素调节、光形态建成、花发育、衰老和病原体防御等方面发挥重要作用[51］．在这项研究中,桶psi-miR171a-3p和psi-miR397a在C_1 vs T_1、C_9 vs T_9、C_5 vs C_9、C_1 vs C_9中表达上调。

研究还发现，psi-miR396b-5p的表达与其靶标ABA 8 ' -羟化酶1 (CYP707A1)．CYP707A1是ABA氧化分解的关键酶，这在ABA的积累中发挥了负面调节作用52］．过度的CYP707A1减少ABA水平和显示出ABA-缺陷表型[53］．激素如ABA，细胞分裂素和山梨糖醇内酯苹果脱落发挥关键作用[54］．我们的研究结果表明，psi-miR396b-5p随着其靶标的下调而下调CYP707A1在T_5和T_9中上调。这表明psi-mir396b-5p-CYP707A1在下调的psi-miR396b-5p可能上调的情况下，对与花萼脱落有关CYP707A1抑制ABA表达的增加，促进花萼脱落。经qRT-PCR证实，psi-miR396b-5p和CYP707A1参与花萼滞留(图。10.)．

认为细胞分离与由几种降解酶，尤其是多肢乳突酶（PG）引起的果胶降解有关。PGS是涉及植物细胞壁和细胞间层中果胶降解的细胞壁水解酶系列之一。它涉及种子萌发，裂开，花粉成熟，花粉囊裂纹，果实成熟和器官脱落[55，56］．在本研究中，psi-miR396b-3p在T_9相对于T_1表达下调，而其靶蛋白PG表达上调。这说明减少psi-miR396b-3p的表达可增强PG的表达，从而调控花萼脱落。

许多研究证明，蛋白磷酸酶2C参与植物信号转导途径、抗病等多种生命活动。近年来，一些蛋白磷酸酶也被报道参与植物器官脱落的信号转导[57］．Corbacho等人[58发现蛋白质磷酸酶可以参与甜瓜果实脱落的调节。Meir等人。［59]检测了番茄花AZ组织的转录组变化，发现一个蛋白磷酸酶基因可以影响脱落过程。提示蛋白磷酸酶2C可能参与了器官脱落的调控。通过综合分析，预测psi-miR394a在C_9 vs T_9中靶向1个蛋白磷酸酶2C基因和3个可能的蛋白磷酸酶2C基因，在T_1 vs T_9中靶向1个蛋白磷酸酶2C基因和6个可能的蛋白磷酸酶2C基因。psi-miR394a在T_5和T_9中的表达水平下调，而蛋白磷酸酶2C (HAB1)的表达上调，提示HAB1起着萼脱落（图重要的作用。10.)．

目前的研究表明，植物受体样激酶在调节植物生长发育、免疫应答和细胞死亡等方面发挥着重要作用[60，61］．一些受体样激酶，如HAESA/HSL2，EVR.，SERK1和CSJ已被报道可控制花器官的脱落[62，63，64，65］．艾达，一种可能编码短肽配体的基因，影响着花器官中AZ细胞分离的最后阶段拟南芥蒂利亚纳［66］．富含亮氨酸的重复受体样激酶HAESA在拟南芥中艾达受体，其中也起着花器官脱落了重要的作用。在拟南芥蒂利亚纳的淘汰赛HAESA及其同源基因HAESA样2（HSL2)可显著延迟花器官的脱落[67］．在差异的靶基因中，我们发现了大量与受体样蛋白激酶和可能受体样蛋白激酶相关的基因。有一个基因编码cysteine-rich受体蛋白激酶10 psi-miR160a-3p的目标和psi-miR171a-3p c₁和T_1和c₁vs C_9,分别两个基因编码受体胞质丝氨酸/ threonine-protein激酶在C_5、C_9 RBK2 psi-miR399b的目标,四个受体蛋白激酶基因c₁vs C_9,其中1个受体蛋白激酶基因在T_1和T_5中被psi-miR396a-3p靶向，16个受体蛋白激酶基因在T_5和T_9中被靶向。QRT-PCR显示，随着PP的增加，psi-miR396a-5p的表达逐渐降低_333.治疗，而其靶受体样蛋白激酶Herk1表达逐渐增加（图。10.)．我们推测这些mRNA-miRNA对可能与花萼脱落有关。半乳糖苷酶参与细胞壁降解。吴和伯恩斯[68在柑桔脱落过程中寻找β-半乳糖苷酶的差异表达基因。Southern blot分析表明，至少有两个密切相关的b-半乳糖苷酶基因与果实脱落有关。在花期的雄蕊、花瓣和幼果中均有β -半乳糖苷酶的表达，表明β -半乳糖苷酶基因可能在花和幼果的脱落和早期生长发育过程中发挥作用。我们在不同的靶基因中发现了多个编码β -半乳糖苷酶的基因。C_1和C_5中有1个β -半乳糖苷酶基因，C_1和C_9中有2个β -半乳糖苷酶基因，T_1和T_9中有3个β -半乳糖苷酶基因。这些细胞均被pci - mir858b靶向，其可能通过调节β-半乳糖苷酶的表达来促进花萼脱落。

生长素具有抑制器官脱落的效果。在拟南芥蒂利亚纳，脱落层的IAA信号通路对调控器官脱落是必要的[69］．一种营养响应因素ARF1，ARF2，ARF7， 和ARF19基因的拟南芥蒂利亚纳参与了花器官脱落的调节[70，71］．我们的结果表明ARF6，ARF18， 和ARF25分别被psi-miR160a-3p、psi-miR156j和psi-miR167d/167a-5p靶向。qRT-PCR和测序结果均显示，这些mRNA-miRNA对在花萼保留组织和花萼脱落组织之间表达有显著差异(图。10.)．在茶树(茶树） 和芦笋， ARF也通过5'RLM-RACE被鉴定为miR160的靶基因，并表现出负相关[72］．三叶橙(枳壳trifoliata），ARF被鉴定为MIR167至5'RLM-rame的靶基因[73］．miR160的核苷酸序列中的茶树和芦笋与psi-miR160a-3p的核苷酸序列相同。三叶橙中miR167的核苷酸序列与psi-miR167a-5p的核苷酸序列相同。这些都可以证明我们结果的可靠性。因此，psi-miR160a-3p、psi-miR167d、psi-miR167a-5p和psi-miR156j可能通过调控促进花萼脱落东盟地区论坛表达。

通过综合分析，我们还获得了其他可能与花萼脱落相关的miRNA-mRNA对。我们发现psi-miR319a靶向锌指蛋白(ZAT12）在C_1 VS C_9和T_1 VS T_9，PSI-miR396b-3P在T_1靶向B框锌指蛋白32 VS T_5和T_1 VS T_9。Cho等人。［74]报道过表达的青指蛋白(AtZFP2)能延缓花器官的脱落拟南芥蒂利亚纳．SPL可以由PSI-miR156j / PSI-miR157d来调节，并且因此与器官脱落相关联75，76］．NAC100在T_1 vs T_5和T_1 vs T_9中，预测psi-miR164a靶向。脱落酸诱导的基因表达受NAC调控，在橄榄分离区高度表达[77]和苹果[78］．三叶橙(枳壳trifoliata)， SPL转录因子基因也通过5'RLM-RACE验证为ptr-miR156的靶点[73］．Jeyaraj等人[72]通过5'RLM-RACE证实了csn-miR164a的靶向作用NAC100呈现出负相关。李等人[79)表明,NAC100可以通过Zm-miR164的5 '端配对的第10和第11个碱基之间的结合区域的切割来调控，这与Jeyaraj的报道一致。这些研究支持我们的目标预测分析的有效性。在不同的靶基因中，我们获得了4个SPLC_5和T_5的基因，2SPL基因C_1和T_1，四SPL在C_1 VS C_5，四个SPLC_5和C_9的基因，2SPLC_1和C_9的基因，1SPL基因在T_1和T_5，和7SPLT_1和T_9基因。

应该注意到这项研究的几个不足之处。首先，基于转录组测序和小rna测序进行生物信息学分析，并通过qRT-PCR验证多个mRNA和miRNA的表达。预测的miRNA靶点的切割位点将在未来通过5'RLM-RACE进行验证。其次，鉴于现有研究发现的mirna下调，本研究重点探索miRNA-mRNA负调控对。然而，我们的研究也发现了大量的阳性miRNA-mRNA对。mirna对基因的正调控还有待进一步研究。第三，我们通过小RNA测序鉴定了84个新的mirna，许多新的mirna在整合分析中也存在关联。这些新型mirna的真实性有待进一步鉴定。

综上所述，通过对mRNA和miRNA的综合分析，我们建立了miRNA-mRNA对的网络，以了解miRNA对花萼脱落的精确调控。6个处理组共2587对miRNA-mRNA，存在正相关和负相关。这些差异靶基因主要参与萜类主链的生物合成、光合-天线蛋白、卟啉和叶绿素的代谢、类胡萝卜素的生物合成、玉米素的生物合成和植物激素信号转导。此外，一些关键的miRNA-mRNA对可能与花萼脱落有关，包括蛋白磷酸酶2C (psi-miR394a-HAB1)、受体样蛋白激酶(psi-miR396a-5p-Herk1.），纤维素合酶 - 样蛋白D3（PSI-miR827-CSLD3），β-​​半乳糖苷酶（psi-mir858b-β-半乳糖苷酶），SPL8 ' -羟化酶1 (psi-miR396a-5p-CYP707A1）和生长素响应因子（Psi-mir160a-3p-ARF6psi-miR167d -ARF18psi-miR167a-5p -ARF25）等，通过综合分析获得。虽然这项研究不能完全考虑到Korla香梨中的Calyx脱落，但本研究中透露的MRNA和MiRNA将有助于理解Calyx脱落所涉及的可能机制。

阿坝:

脱落酸

c₁:

用GA采集花期晚的第一天的花_3.治疗

C_5:

用GA3处理在花期后期第5天采花

C_9:

收集的花在晚绽放舞台的第九天与ga_3.治疗

DEG：

差异表达基因

DEM：

微分microrna表达

遗传算法_3.：

赤霉素

去：

基因本体论

IAA：

Indole-3-acetic酸

Kegg：

Kyoto基因和基因组的百科全书

ko：

KEGG数据库直接同源

Kog / Cog：

直际蛋白质群体

NR：

NCBI非冗余蛋白序列

NT：

NCBI非冗余核苷酸序列

PFAM：

蛋白家族

PG：

polygalacturonase.

页_333.：

摘要

QRT-PCR：

定量实时聚合酶链反应

SWISS-PROT：

一个手工注释和审查的蛋白质序列数据库

T_1:

与PP晚花期的第一天收集的花朵_333.治疗

T_5:

用PP采集花期后期第5天的花_333.治疗

T_9：

与PP晚花期的第九天收集的花朵_333.治疗

zr：

Trans-Zeatin-riboside

我们感谢维明秦（Shayidong园艺场，新疆库尔勒市）为我们提供了植物材料。

资金

国家自然科学基金项目(no . 31360474);高等学校博士学科点专项科研基金项目(no . 2013651810002);国际科技合作项目(no . 2009YD32);资助方不参与研究设计、数据分析和解释、文稿撰写，只提供资金。

数据和材料的可用性

支持数据包含在文章和其他文件中。

从属关系

1. 石河子大学农学院园艺系，新疆石河子832003

   马丽，周丽，权绍文，徐航，杨杰平，牛建新

2. 新疆生产特种水果和蔬菜栽培生理和种质资源利用建设兵团重点实验室，石河子，832003新疆，中国

   马丽，周丽，权绍文，徐航，杨杰平，牛建新

贡献

JXN领导和协调项目，JXN和LM设计研究。LM、SWQ和JPY收集植物材料，分离RNA。LM和JPY进行实时定量PCR。LM、LZ、HX进行生物信息学分析。LM写了论文。所有作者已阅读并同意最终手稿。JXN是通讯作者，负责所有的联系和通信。所有作者阅读并批准了最终的手稿。

通讯作者

对应到牛建新．

伦理批准和同意参与

不适用。

同意出版物

不适用。

利益争夺

两位作者宣称他们没有相互竞争的利益。

出版商的注意事项

施普林格《自然》杂志对已出版的地图和机构附属机构的管辖权要求保持中立。

开放获取本文遵循知识共享署名4.0国际许可协议(http://creativecommons.org/licenses/by/4.0/）如果您向原始作者和源给出适当的信用，则允许在任何介质中进行不受限制的使用，分发和再现，提供指向Creative Commons许可证的链接，并指示是否进行了更改。Creative Commons公共领域奉献豁免（http://creativecommons.org/publicdomain/zero/1.0/）除非另有说明，否则适用于本文中提供的数据。

重印和权限