Genome wide identification of chilling responsive microRNAs in Prunus persica
© Barakat et al.; licensee BioMed Central Ltd. 2012
Received: 22 June 2012
Accepted: 23 August 2012
Published: 15 September 2012
Skip to main content
© Barakat et al.; licensee BioMed Central Ltd. 2012
Received: 22 June 2012
Accepted: 23 August 2012
Published: 15 September 2012
MicroRNAs (miRNAs) are small RNAs (sRNAs) approximately 21 nucleotides in length that negatively control gene expression by cleaving or inhibiting the translation of target gene transcripts. Within this context, miRNAs and siRNAs are coming to the forefront as molecular mediators of gene regulation in plant responses to annual temperature cycling and cold stress. For this reason, we chose to identify and characterize the conserved and non-conserved miRNA component of peach (Prunus persica (L.) Batsch) focusing our efforts on both the recently released whole genome sequence of peach and sRNA transcriptome sequences from two tissues representing non-dormant leaves and dormant leaf buds. Conserved and non-conserved miRNAs, and their targets were identified. These sRNA resources were used to identify cold-responsive miRNAs whose gene targets co-localize with previously described QTLs for chilling requirement (CR).
Analysis of 21 million peach sRNA reads allowed us to identify 157 and 230 conserved and non-conserved miRNA sequences. Among the non-conserved miRNAs, we identified 205 that seem to be specific to peach. Comparative genome analysis between peach and Arabidopsis showed that conserved miRNA families, with the exception of miR5021, are similar in size. Sixteen of these conserved miRNA families are deeply rooted in land plant phylogeny as they are present in mosses and/or lycophytes. Within the other conserved miRNA families, five families (miR1446, miR473, miR479, miR3629, and miR3627) were reported only in tree species (Populus trichocarpa, Citrus trifolia, and Prunus persica). Expression analysis identified several up-regulated or down-regulated miRNAs in winter buds versus young leaves. A search of the peach proteome allowed the prediction of target genes for most of the conserved miRNAs and a large fraction of non-conserved miRNAs. A fraction of predicted targets in peach have not been previously reported in other species. Several conserved and non-conserved miRNAs and miRNA-regulated genes co-localize with Quantitative Trait Loci (QTLs) for chilling requirement (CR-QTL) and bloom date (BD-QTL).
In this work, we identified a large set of conserved and non-conserved miRNAs and describe their evolutionary footprint in angiosperm lineages. Several of these miRNAs were induced in winter buds and co-localized with QTLs for chilling requirement and bloom date thus making their gene targets potential candidates for mediating plant responses to cold stress. Several peach homologs of genes participating in the regulation of vernalization in Arabidopsis were identified as differentially expressed miRNAs targets, potentially linking these gene activities to cold responses in peach dormant buds. The non-conserved miRNAs may regulate cellular, physiological or developmental processes specific to peach and/or other tree species.
microRNAs (miRNAs) are small non-coding RNAs approximately 19 to 22 nucleotides (nt) in length that function in gene regulation at the post-transcriptional level . They are transcribed from one DNA strand as long precursors, which under the processing of various enzymes lead to mature double strand miRNAs . One strand of the miRNA is degraded; the other is incorporated in a ribonucleoprotein complex, the RNA-induced silencing complex (RISC), which is involved in RNA interference (RNAi). The mature miRNA guides the action of RISC to degrade mRNAs of its targets.
Research in recent years demonstrated that miRNAs play a significant role in the epigenetic regulation of genes controlling various cellular and developmental processes . Indeed, most miRNA targets identified to date are involved in the development of plant organs such as root, leaf, and flower [1, 3]. However, a significant fraction of miRNA target genes seems to be involved in cellular defense against abiotic stresses such as nutrient deprivation, drought, heat, and UV exposure  or biotic stresses including attack by fungi such as Phakopsora pachyrhizi and Cronartium quercuum[5, 6] , or bacteria such as Pseudomonas syringae.
Cold temperature exposure is a major abiotic stress that affects plant development and growth and significantly impacts crop production. Plants evolved various processes to cope with cold stress at the cellular and developmental levels. Flower and vegetative winter bud dormancy responses enable plants to avoid flower and leaf damage caused by fluctuating temperatures during winter and to synchronize the re-initiation of growth with the annual return of optimal temperatures in the spring. Genetics analyses of the cold response linked with studies on vernalization suggest that sRNAs play a significant role in these phenological traits. For instance, expression analysis in Arabidopsis thaliana Populus trichocarpa Brachypodium distachyon and other plant species [4, 8–10] identified several miRNAs with altered expression under cold stress. Another study showed that Arabidopsis miR402 positively regulates seed germination under dehydration or cold stress conditions  by repressing DEMETER-LIKE protein3 (DML3), which is involved in DNA demethylation . Other studies [13, 14] showed that the repression of the Flowering Locus C (FLC) by vernalization was mediated by a 24-nt siRNA.
Peach belongs to Prunoideae, one of the four subfamilies of the Rosaceae. Due to the compact genome size (227 Mb) and the extensive genomics and genetics resources available for this species (e.g. a high density integrated genetic/physical map, extensive cDNA sequence resources, availability of a doubled haploid), the peach has been highlighted as a genome model for fruiting trees in the Rosaceae. In addition, as a temperate perennial tree species with a reasonably short juvenility period (two years) and a reasonably diverse breeding germplasm, it represents an excellent model for the genetic dissection of traits important to both basic and applied tree research. As a stone fruit, peach holds a key phylogenetic position in the Rosaceae, which is comprised of over 100 genera and 3,000 species . The Rosaceae family ranks third in economic importance of the plant families in temperate regions ; it includes many species such as almond, apple, apricot, blackberry, cherry, peach, raspberry, rose and strawberry, valued for nuts, fruits and flowers. Recently the genome of the peach was completely sequenced and the assembly and gene annotation made publicly available on the Genome Database for Rosaceae[17, 18]. However, no study reporting miRNAs from peach is yet published. Identifying miRNAs from peach is very important to infer ancestral relationships from divergent pools of miRNAs in other closely related Rosaceae species [19, 20] and other more distantly related angiosperms. Moreover, the identification of miRNAs from peach whose expression is correlated with changing environmental cues provides a tool for uncovering the network of regulated gene activities associated with phenological responses in trees. Identifying miRNAs from this species enables analyses of the function of these miRNAs in Rosaceae species and the evaluation of the level of conservation of the miRNA-mediated gene regulatory networks in this botanical family.
In this work, sRNAs from young emerging leaves and chilled vegetative buds were sequenced using the SOLID Platform (Applied Biosystems, Foster City, CA) and the sRNA transcriptome was analyzed. sRNA analysis in combination with miRNA mining in the whole genome sequences of peach identified 157 and 230 unique conserved and non-conserved miRNA sequences, respectively. In silico expression analyses enabled the identification of several miRNAs displaying differential expression in winter dormant buds versus pre-dormant leaves. The prediction of miRNA targets showed that these miRNAs targets are involved in various biological processes. Several miRNAs or miRNA-regulated genes were found co-localizing with previously reported QTLs for chilling requirement (CR-QTLs) and bloom date (BD-QTLs) . Several of these genes appear to be involved in response to cold stress. Comparative and phylogenetic analysis of miRNA distribution in peach and other species provide insight into the evolution of miRNA families.
Sequencing summary of peach sRNAs
18-24 nt (reads)
More than 10 hits (reads)
Conserved miRNAs (reads)
Non-conserved sRNA sequences (reads)
Non-conserved sRNA unique sequences (reads)
Non-conserved sRNAs passed MirCheck and have miRNA* prior to visual inspection and length variant removal (reads)
All miRNAs that are conserved in Arabidopsis, rice and Populus (11) or conserved between Arabidopsis and Populus (23) were found in the peach sequences (Additional file 2: Table S2). Sixteen of these miRNA families were deeply rooted in the phylogeny of land plants; they were present in either moss or lycophytes (Additional file 2: Table S2). Among these miRNA families, five (miR1446, miR473, miR479, miR3629, and miR3627) were reported only in tree species (Populus trichocarpa, Citrus trifolia, and Prunus persica (L.) Batsch). A search for these miRNAs in the genomes of Arabidopsis and Oryza using patscan assembly with up to two substitutions identified a homolog only for miR3629 in Oryza. Twenty three conserved miRNA families were found in peach and other plant species but have not been reported in other tree species. A search of these miRNA in the Populus genome enabled the identification of homologs from 11 families (Additional file 1: Table S1). Other miRNA families (miR1030 and miR1088) were reported only in moss (Physcomitrella patens) or fern (Selaginella moellendorfii) and peach. A set of miRNA (miR2275, miR1850, miR2919, and miR5083) were found in peach and monocots but were not reported in other land plant species.
Mapping the identified conserved miRNAs on previously reported CR-QTLs and BD-QTLs  showed that 37 sequences belonging to 18 conserved miRNA families were co-localizing with CR-QTLs and BD-QTLs. Three of these families (miR5021, miR2919, and miR414) were predicted to target genes involved in response to cold stress.
MirCheck  analysis of 478,393 and 210,296 non-conserved sRNA sequences from leaves and buds, respectively, enabled the identification of 19,993 and 7326 unique miRNA sequences that satisfied the biogenesis criterion. After removing miRNA candidates that have no sequenced miRNA* and those that were represented by less than two reads in either leaves or buds datasets, as well as length variants, 230 unique non-conserved miRNAs were identified (Additional file 3: Table S3). A search of these sequences that have not been previously reported in other model plants using patscan with up to two substitutions showed that 15, 14, and 22 were found in Arabidopsis Oryza, and Populus, respectively. These non-conserved miRNAs were represented by a number of reads ranging from 1 to 5636 and from 1 to 3446 for buds and leaves respectively. Among the identified 230 non-conserved miRNAs, 31 were co-localized with CR-QTLs and BD-QTLs and five of them regulate genes involved in development processes. Six of these miRNAs (pper-3-1, pper-49-1, pper-71-1, pper-99-1, pper-161-1, and pper-164-2) target genes involved in development while two others (pper-99-1 and pper-49-1) regulate genes involved in cold response. One miRNA (pper-82-1) of these 31 was also differentially highly expressed in winter buds (see below).
A target search using the peach proteome  enabled the prediction of targets for 136 out of 157 (86.6 %) conserved miRNAs (Additional file 7: Table S7 and Additional file 8: Table S8). At least one target was identified for each miRNA family; however, several targets were identified for most miRNAs. Targets with top scores are presented in Additional file 7: Table S7. The other targets identified with scores less than four are listed in Additional file 8: Table S8. GO annotation using the Arabidopsis proteome (Figure 5) showed that 12.64 % and 23.94 %, and 21.39 % of the targets belong to the function categories “response to stimuli”, “cellular process” and “development processes” respectively. Genes from development processes are all involved in various processes of flower bud induction and flower morphogenesis and development. A few genes such as chitinase A, and AICARFT/IMPCHase are involved in regulation of plant response to cold stress. A large fraction of the predicted miRNA targets (Additional file 8: Table S8) were previously identified in other studies. However, new targets were predicted for several conserved miRNA families in peach. Confirmation of these predicted targets by analyzing the degradome  or using 5’ RACE is necessary.
At least one target gene was predicted for 166 out of 230 (72.1 %) non-conserved miRNAs (Additional file 9: Table S9 and Additional file 10: Table S10). GO annotation analysis identified several genes involved in plant development including several embryo defective proteins, Late embryogenesis abundant, PPR2, and DNA glycosylase. A few targets encode proteins involved in cold stress response such as Annexins, WD-40, phosphatidylinositol-specific phospholipase C, and cellulose synthase. Most importantly, several identified targets encode proteins involved in sRNA biogenesis and methylation such as a Dicer homolog (RNA helicase), ribonuclease III, DNA glycosylase DEMETER, histone 3 lysine 9 specific methyltransferase, and SUVH4/KYP homolog.
Several genes from QTL regions were found being targeted by conserved and non-conserved miRNAs. One of these genes (Reduced Vernalization 1 (VRN1) is a key player in regulating vernalization in other plant species. It was predicted to be targeted by a non-conserved miRNA (pper_100_1).
sRNA sequencing of over 21 million sRNA reads from predormant young leaves and dormant leaf buds enabled the identification of over 387 conserved and non-conserved miRNA sequences from peach. The conserved miRNA dataset include all miRNAs conserved in Arabidopsis Oryza and Populus or conserved between Arabidopsis and Populus suggesting that the dataset identified here is complete. Among the identified conserved miRNA set, five conserved miRNA families (miR1446, miR473, miR479, miR3629, and miR3627) and several non-conserved miRNAs were reported only in tree species (Populus trichocarpa Citrus trifolia, and Prunus persica (L.) Batsch). A search of these miRNAs in the genomes of the herbaceous species Arabidopsis and Oryza confirmed that four of them were not found in the Arabidopsis and Oryza genomes. These results confirm that a small set of miRNAs may exist only in tree species and function in biological processes specific to them. For instance, it was suggested that miR1446 and miR473 may play important roles in the growth of trees such as the formation of specialized woody tissue or in the response of trees to mechanical stress . Other miRNA families (miR1030, miR1088) were reported only in moss (Physcomitrella patens) or fern (Selaginella moellendorfii) and peach. A set of miRNAs (miR1850, miR2275, miR2919, miR2931 and miR5083) was found in peach and monocots but not in other land plant species. This pattern of gain/loss of certain miRNA families in various plant lineages was previously reported [27, 28]. The inability to detect some miRNAs in some angiosperm lineages could be due to the fact that these miRNAs were not expressed in tissues used for constructing sRNA libraries in the non-model species studied. However, the absence of these miRNAs in whole genome sequences of model species could only be explained by the fact that some miRNAs, such as the ones identified only in trees, are specific to some lineages or species. It is noteworthy that the loss of some miRNAs is not unidirectional as some miRNAs from peach were lost either in monocots (Oryza), or eudicots (Arabidopsis and Populus).
Analysis of the size of miRNA gene families showed that most of the conserved miRNA families from peach have a similar number of members in Arabidopsis but less genes than in Populus. However, one miRNA family (miR5021) was over-represented in peach compared to both Arabidopsis and Populus. The increase in gene number of this family could be associated with the biological role of these miRNAs in peach.
Since we were interested in CR, we searched for miRNAs that co-localize with CR-QTLs and BD-QTLs. We found 31 miRNAs co-localizing with seven CR-QTLs and BD-QTLs. Members of four families (miR5021, miR164, miR414, miR396, miR2919) target genes encoding proteins involved in dormancy such as embryo-defective or maternal effect embryo arrest 14 (MEE14). Others encode proteins involved in cold stress response such as annexin 5, AICARFT/IMPCHase, ATMEKK1, cellulose, and MAPK/ERK kinase. The list of targets also includes genes encoding a histone-lysine N-methyltransferase that is involved in vegetative to reproductive phase transition of meristem. Like conserved miRNAs, six non-conserved miRNAs were located in the CR-QTL and BD-QTL regions. Five of these miRNAs target genes that are known to be involved in development processes and two target genes known to respond to cold stress.
Relative expression of the conserved miRNAs varied widely, based on the number of sequences observed for each miRNA in our dataset. The ten most highly expressed miRNAs (miR156, miR157, miR159, miR164, miR167, miR172, miR393, miR396, miR414, miR2275, and miR5021) in buds and leaves are miRNAs regulating genes involved in flower and leaf development processes such as integument development, leaf morphogenesis, meristem initiation, maintenance, and growth, bilateral symmetry determination, organ morphogenesis, plant phase transition, shoot apical meristem identity, flower and fruit development, and plant architecture. All of these miRNAs are cold responsive . Unlike conserved miRNAs, only one of the highest expressed non-conserved miRNAs in pre-dormant leaves is known to play a role in a developmental process. Expression analysis using DEGseq identified several miRNAs highly expressed in winter buds. Genes targeted by these miRNAs include genes involved in apical meristem formation, pollen development, and an RNA slider that may play a role in posttranscriptional silencing. Several of these miRNAs were previously reported as induced following cold stress . The occurrence of stress-related cis elements such W-box in the regions of several of these miRNAs may play a role in their activation . Three of these miRNA genes (miR156, miR172, and miR398) were also reported as responding to cold stress in several studies [4, 8, 10, 29, 30]. miR172 was suggested to fine-tune plant development under continuously fluctuating temperature conditions . Similarly, ten non-conserved miRNAs were differentially expressed in dormant winter buds (at least a two fold higher expression) of which three were located in CR-QTLs and BD-QTLs. However, none of these miRNAs seem to regulate known genes involved in cold stress, flower or leaf development.
Target analysis enabled the prediction of target genes for all conserved and over half of the non-conserved miRNAs. miRNAs for which no target could be identified may correspond to young non-functional miRNAs. A large set of genes targeted by peach miRNAs encode genes involved in flower and leaf development. A small set of genes are involved in cold stress. Among these genes is VRN1 a known regulator in vernalization . This is in agreement with previous studies on Arabidopsis Brachypodium, and Populus showing that miRNAs are involved in response to cold stress [4, 8–10, 32]. We also identified genes involved in methylation or RNA biogenesis that may function in epigenetic regulation of plant response to cold stress. Among these genes, a Dicer homolog (RNA helicase) involved in miRNA processing , a ribonuclease III that is required for endogenous RDR2-dependent siRNA formation , a DNA glycosylase DEMETER responsible for endosperm maternal-allele-specific hypomethylation at the MEDEA (MEA) gene , a histone 3 lysine 9 specific methyltransferase involved in the maintenance of DNA methylation , and a SUVH4/KYP homolog that is involved in epigenetic control of FLC expression . This is in accordance with recent studies showing that epigenetic modifications including methylation, siRNA signatures, and histone modifications play a major role in vernalization [13, 38].
Furthermore, several genes from CR-QTLs and BD-QTLs regions that were reported as playing key roles in regulating vernalization were predicted targets by miRNAs. Several of these genes belong to the polycomb repressor 2 silencing complex that regulates the expression of the key flowering gene FLC. For instance, VRN1 is targeted by the non-conserved miRNA (pper_100_1). These two genes function in maintaining the repression of the major target of the vernalization pathway, the floral repressor FLC. However, it is currently unknown what role these two genes may play in CR regulation in perennial trees, which is a significantly different phenological character particularly since the FLC gene has not yet been reported as present in the peach genome or the genomes of its close relatives (Zhebentyayeva, personnal communication).
In conclusion, this study identified hundreds of conserved and non-conserved miRNAs and their target genes in peach, a Rosaceae model species expanding our understanding of the pattern of conservation of miRNA in land plant lineages. This study also identified several cold-responsive miRNAs and their predicted gene targets present in CR-QTLs and BD-QTLs thus, potentially connecting miRNA regulatory activities to CR and dormancy in peach. This miRNA dataset will be very useful to the scientific community working on Rosaceae and other plant families for functional analyses of genes of interest and for deciphering the gene regulatory networks of CR and bud dormancy in peach.
Young emerging leaves and chilled vegetative winter buds (600 chill-hours accumulated) of the double haploid reference genome genotype “Lovell” (Clemson University) were used for sRNA sequencing. Samples were collected from trees grown at Musser Farm, Clemson SC, frozen in liquid nitrogen and stored at −80˚. Total RNA was prepared by the method of Carra and collaborators . The quality of total RNA was assessed with an Agilent Technologies 2100 Bioanalyzer (Agilent Technologies) for possible degradation. Total RNAs from both samples used in this study have RNA Integrity Numbers (RNI) values of 7.9 and 8 indicating their high quality .
sRNA enrichment was performed using FlashPage according to the manufacturer’s instructions (Ambion/Life Technologies). sRNA libraries were prepared from enriched sRNA using the “Prepare Small RNA Libraries “protocol from the SOLID Total RNA-Seq Kit according to the manufacturer’s instructions (Applied Biosystems/ LifeTechnologies). cDNA samples were subjected to ePCR to prepare templated beads. Templated beads were deposited on a quarter of a flow cell and 50 nt sequencing was performed on a SOLID 3 Plus System (Applied Biosystems/Life Technologies) at Penn State University.
The sequences produced by SOLID sequencing, which correspond to sRNA-derived cDNAs, were processed using two approaches. The first one consisted of converting sequences from color space to nucleotide format using custom perl scripts. The second one consisted in mapping sequences in color space format to the peach genome using Bowtie software  and by tolerating up to two mismatches. Sequences obtained using both approaches were analyzed separately. sRNA sequences smaller than 18 and larger than 23 nucleotides in length were discarded. sRNA sequences that passed the adapter check and size filter were then screened using BLASTN against chloroplast and mitochondrial genomes , tRNA , rRNA , snoRNA , and all contaminating rRNA, tRNA and snoRNA sequences were removed. The cleaned sequences were then sorted by sequence identity, the size variants identified, and the relative count of each miRNA was determined. Unique sRNA sequences were queried against known miRNAs using the program Patscan  with default parameters and two mismatches in order to identify homologs of known miRNAs. Flanking sequences (150 nt at each side) of the mature miRNA were retrieved and folded using the RNA fold program . The folded sequences were then input into the MirCheck program  to check their secondary structures. Sequences that passed MirCheck were then inspected manually for miRNA features. Annotation of miRNAs was performed as described by Meyers and collaborators . Another set of miRNAs not expressed in leaf and bud datasets was identified by querying miRBase plant sequences against the peach genome using Patscan as described above. Sequences that passed MirCheck were checked for the conservation of their mature miRNA sequences and their stem loop secondary structure in other species as described previously . Sequences that satisfied this criterion were also annotated as miRNAs.
Sequences with no similarity to known miRNA sequences were used to search for non-conserved miRNAs in the peach genome using Patscan with a setting of 0 mismatches, 0 insertions, and 0 deletions. Sequences that had over than 10 matches to the peach genome may correspond to repeated sequences and were removed. Genomic regions (150 at each side) harboring sRNA sequences were then retrieved and the sequences folded using RNAfold  and the secondary structure was checked for miRNA features using MirCheck. Sequences that passed MirCheck were sorted for redundancy. When several length variants of the same miRNA were sequenced, only variants with the highest representation were considered. Only sRNAs that passed the MirCheck filter and have a miRNA* sequenced at least in one of the libraries was annotated as miRNAs.
A search for miRNA target genes was performed by querying the peach predicted coding DNA sequences (CDS)  using the scoring approach previously described . Only hits with a score less than 4 were considered good target sequences. Target sequences were then annotated using the Arabidopsis proteome (BLASTX, e-value < 1e-5.
Expression of miRNAs in leaves and buds samples were quantified using DEGseq . DEGseq uses the MA plot based method with random sampling. This method works under the assumption that RNA sequencing is random sampling and the number of reads for each gene follows a binomial distribution. Fisher exact test and likelihood ratio tests were then performed on this model to identify differentially expressed genes.
DNA complementary to RNA
Polymerase Chain Reaction
Rolling Circle Amplification
Small interfering RNA
Small nucleolar RNA
The authors thank the International Peach Genome Initiative (IPGI) for accessing the peach genome sequence. We also thank our colleague Craig Praul for performing the peach sRNA sequencing and Mike Axtell for providing us with scripts to convert the format of sRNA sequences from color-space to nucleotide space. We thank Dana Nelson and SRS’s support to Abdelali Barakat through Cooperative Agreement 11-CA-11330126-106. This work was supported by a USDA NIFA SCRI Award to Dorrie Main and Albert Abbott and the Coker Endowed Research Chair at Clemson University.
This article is published under license to BioMed Central Ltd. This is an Open Access article distributed under the terms of the Creative Commons Attribution License (http://creativecommons.org/licenses/by/2.0), which permits unrestricted use, distribution, and reproduction in any medium, provided the original work is properly cited.