Expansion of Capsicum annum fruit is linked to dynamic tissue-specific differential expression of miRNA and siRNA profiles

PLOS ONE, Jul 2018

Small regulatory RNAs, such as microRNAs (miRNAs) and small interfering RNAs (siRNAs) have emerged as important transcriptional and post-transcriptional regulators controlling a wide variety of physiological processes including fruit development. Data are, however, limited for their potential roles in developmental processes determining economically important traits of crops. The current study aimed to discover and characterize differentially expressed miRNAs and siRNAs in sweet pepper (Capsicum annuum) during fruit expansion. High-throughput sequencing was employed to determine the small regulatory RNA expression profiles in various fruit tissues, such as placenta, seed, and flesh at 28 and 40 days after anthesis. Comparative differential expression analyses of conserved, already described and our newly predicted pepper-specific miRNAs revealed that fruit expansion is accompanied by an increasing level of miRNA-mediated regulation of gene expression. Accordingly, ARGONAUTE1 protein, the primary executor of miRNA-mediated regulation, continuously accumulated to an extremely high level in the flesh. We also identified numerous pepper-specific, heterochromatin-associated 24-nt siRNAs (hetsiRNAs) which were extremely abundant in the seeds, as well as 21-nt and 24-nt phased siRNAs (phasiRNAs) that were expressed mainly in the placenta and the seeds. This work provides comprehensive tissue-specific miRNA and siRNA expression landscape for a developing pepper fruit. We identified several novel, abundantly expressing tissue- and pepper-specific small regulatory RNA species. Our data show that fruit expansion is associated with extensive changes in sRNA abundance, raising the possibility that manipulation of sRNA pathways may be employed to improve the quality and quantity of the pepper fruit.

A PDF file should load here. If you do not see its contents the file may be temporarily unavailable at the journal website or you do not have a PDF plug-in installed and enabled in your browser.

Alternatively, you can download the file locally and open with any standalone PDF reader:


Expansion of Capsicum annum fruit is linked to dynamic tissue-specific differential expression of miRNA and siRNA profiles

July Expansion of Capsicum annum fruit is linked to dynamic tissue-specific differential expression of miRNA and siRNA profiles De nes Taller 0 1 2 Jeannette Ba lint 0 1 2 Pe ter Gyula 0 1 2 Tibor Nagy 0 1 2 Endre Barta 0 1 2 Ivett Baksa 0 1 2 GyoÈ rgy Szittya 0 1 2 Ja nos Taller 0 2 Zolta n Havelda 0 1 2 0 Data Availability Statement: The datasets generated and analysed during the current study are available in the National Center for Biotechnology Informations (NCBI) BioProject database under Accession number PRJNA388388 (https:// 1 National Agricultural Research and Innovation Centre, Agricultural Biotechnology Institute, GoÈdoÈllő, Hungary, 2 Department of Plant Science and Biotechnology, Georgikon Faculty, University of Pannonia , Keszthely , Hungary 2 Editor: Turgay Unver , Dokuz Eylul Universitesi , TURKEY Small regulatory RNAs, such as microRNAs (miRNAs) and small interfering RNAs (siRNAs) have emerged as important transcriptional and post-transcriptional regulators controlling a wide variety of physiological processes including fruit development. Data are, however, limited for their potential roles in developmental processes determining economically important traits of crops. The current study aimed to discover and characterize differentially expressed miRNAs and siRNAs in sweet pepper (Capsicum annuum) during fruit expansion. High-throughput sequencing was employed to determine the small regulatory RNA expression profiles in various fruit tissues, such as placenta, seed, and flesh at 28 and 40 days after anthesis. Comparative differential expression analyses of conserved, already described and our newly predicted pepper-specific miRNAs revealed that fruit expansion is accompanied by an increasing level of miRNA-mediated regulation of gene expression. Accordingly, ARGONAUTE1 protein, the primary executor of miRNA-mediated regulation, continuously accumulated to an extremely high level in the flesh. We also identified numerous pepper-specific, heterochromatin-associated 24-nt siRNAs (hetsiRNAs) which were extremely abundant in the seeds, as well as 21-nt and 24-nt phased siRNAs (phasiRNAs) that were expressed mainly in the placenta and the seeds. This work provides comprehensive tissue-specific miRNA and siRNA expression landscape for a developing pepper fruit. We identified several novel, abundantly expressing tissue- and pepper-specific small regulatory RNA species. Our data show that fruit expansion is associated with extensive changes in sRNA abundance, raising the possibility that manipulation of sRNA pathways may be employed to improve the quality and quantity of the pepper fruit. Introduction The fruit, this highly specialized complex organ, is the presumed crucial factor of the reproduction and dispersal success of angiosperm plants during evolution. Various types of tissues compose the complex structure of the fruit requiring sophisticated, multi-layered developmental processes during fruit formation. Understanding the molecular bases of fruit development is very NKFI-K119701 was granted to GS). The funders had no role in study design, data collection and analysis, decision to publish, or preparation of the manuscript. important from a practical point of view as well since different types of fruits provide essential nutrients and minerals for healthy and balanced diet. Eukaryotic organisms, including plants, employ several RNA interference (RNAi) pathways. Some of these pathways work through small RNAs (sRNAs), 20 to 24 nucleotide in length, which mediate the negative regulation of their target RNAs in a sequence-specific manner. In plants, RNAi pathways play indispensable roles in development and biotic and abiotic stress responses. Distinct RNAi pathways generate different classes of small regulatory RNAs such as microRNAs (miRNAs) and small interfering RNAs (siRNAs) [1±3]. MiRNAs are crucial components in the fine regulation of gene expression during development. They modulate gene expression by regulating the degradation or translational repression of their target mRNAs [ 4 ]. Fully functional mature miRNAs are generated from their precursors by some phylogenetically conserved protein families [ 2 ] in the following manner: the nuclear genome-coded miRNA genes (MIR) are transcribed by RNA polymerase II (Pol II) to produce primary miRNAs (pri-miRNAs) possessing a specific and stable hairpin-like secondary stem-loop structure. In the nucleus, an RNase III enzyme, DICER-LIKE1 (DCL1) cleaves the hairpin out of the pri-miRNA resulting in the production of 80±300-nt-long miRNA precursors (pre-miRNAs). In plants, the same DCL1 enzyme then cleaves the precursor at sites determined by structural features to generate small double-stranded miRNA intermediates, the miRNA: miRNA duplexes. These duplexes are 2'-O-methylated at their 3' termini by the HUA ENHANCER 1 (HEN1) methyltransferase and then transported to the cytoplasm where they are loaded into the RNA-induced silencing complex (RISC). This complex displaces one of the strands (miRNA ) while the mature strand (miRNA) remains bound. The miRNA-containing RISC mediates sequence-specific repression via mRNA-cleavage or translational repression. The main executors of RISCs are the ARGONAUTE (AGO) proteins [ 5 ]. The Arabidopsis thaliana genome encodes ten AGO proteins which are specialized for various RNAi pathways but often display functional redundancies. AGO1 is the central component of the RISC playing a dominant role in miRNA-mediated regulation, and its expression is regulated by feedback mechanism-mediated by miR168 [ 6 ]. The generation of siRNAs is different from that of miRNAs because it is dependent on an RNA-dependent RNA polymerase (RDRP) [ 3 ]. The perfect double-stranded (ds)RNA intermediates generated by specific RDRPs can be recognized and cleaved by Dicer-like enzymes other than DCL1 to produce various classes of siRNAs [ 3 ]. HetsiRNAs are mainly 24-nt in length and generated by DCL3 via processing of Pol IV/RDR2-dependent dsRNA species originating from intergenic or repetitive regions of the genome. These 24-nt-long siRNAs are loaded into the RNA-induced transcriptional silencing (RITS) complex, which is functionally different from the RISC. They recognize Pol V-dependent short transcripts emerging from repetitive and intergenic regions and recruit chromatin-modifying complexes. HetsiRNAs play fundamental roles in maintaining genome integrity, by initiating de novo DNA methylation on transposable and repetitive elements mainly during reproductive transitions (meiosis, gametogenesis, and embryogenesis), especially in plants with larger, repetitive genomes. Phased (pha)siRNAs, however, are derived from Pol II-dependent mRNAs which are cleaved by a special class of miRNAs/ siRNAs of 22-nt length [ 7 ]. This cleavage renders one of the cleavage products to be a substrate of RDR6/SUPPRESSOR OF GENE SILENCING 3 (SGS3). The resulting dsRNA is then cleaved by DCL4 endonuclease following a distinctive, regular (phased) pattern starting from the cleavage site. The phasiRNAs are preferentially loaded into a distinctive RISC containing AGO5. A subset of phasiRNAs, the trans-acting siRNAs (tasiRNAs) act mainly at the post-transcriptional level by negatively regulating their target transcripts including gene families which encode disease resistance proteins or transcription factors [ 8 ]. The continuously increasing number of experimental studies demonstrate that specific sRNAs regulate plant tissue and organ differentiation. Scientific experiments proved the importance of miRNAs and tasiRNAs in many aspects of plant development, including leaf, 2 / 20 shoot, flower and root development, phase-change, auxin signaling and response to environmental/biotic stresses [9±12]. In addition to the well-known evolutionarily conserved miRNA families, there are non-conserved, lineage-specific miRNAs that are only present in closely related species [ 13, 14 ]. The identification of new miRNAs and the genome-wide, tissue-specific analyses of miRNA expression can help to specify miRNA-regulated developmental processes substantially affecting the quality of economically important traits of crop plants [15]. The role of miRNAs in fruit development has been extensively investigated by Next Generation Sequencing (NGS) technology [16±19]. Several studies analyzed the expression of miRNAs at a genome-wide level during fruit development in tomato, a Solanaceae model plant [20±24]. Fruit development in another Solanaceae species, the pepper (Capsicum annuum), is less extensively investigated. Conserved and potentially new, species-specific miRNAs have been described as a result of an investigation of different pepper organs including the developing fruit as a whole [25±28]. However, since the activity of miRNAs is spatiotemporally highly regulated [29±31], it is inevitable to study the expression pattern of miRNAs in a tissue-specific manner during the fruit development. Moreover, little is known about the contribution of sRNA-mediated regulation to the fleshy fruit development at the phase of cell expansion where the economically valuable fruit volume is often to reach a two magnitude fold increase. To get high-resolution spatiotemporal information about dynamic sRNA expression pattern associated with the expansion phase of the pepper fruit development, we dissected the fruit. We collected flesh, seed, and placenta samples at 28 and 40 days after anthesis (DAA), and quantified their small RNA populations by high-throughput sequencing. Using the gained datasets, we generated the comparative tissue-specific differential expression profiles of the conserved, the already described, and the newly predicted pepper-specific miRNAs. We found that differential expression activity of miRNAs are intriguingly high at 40 DAA and especially the expansion of flesh is associated with profound changes in the miRNA expression profile. The pronounced role of miRNA-mediated regulation in flesh development was also reflected by the high level of AGO1 protein present in flesh samples. Moreover, the analyses of nonmiRNA-like sRNAs revealed an ample of pepper-specific hetsiRNAs and phasiRNAs playing a role mainly in the development of the seed and the placenta. As a summary, our work demonstrates that the cell expansion phase of the economically important sweet pepper is associated with extensive changes in sRNA abundance implying that breeding technologies capable of modifying the expression of key sRNA regulators can be potent tools to enhance the economic value of pepper. Materials and methods Plant material We used a Hungarian sweet cultivar of Capsicum annuum called FeheÂroÈzoÈn for the experiments that we obtained from ZKI-Vetőmag Kft. (H-6000 KecskemeÂt, MeÂszoÈly Gyula u. 6). We germinated the seeds on plates between wet filter papers in a growth chamber (stable 21ÊC, in darkness). We planted the seedlings into pots and grow them in a light room under long day (16 h light/8 h dark) conditions. We collected the fruits at 28 and 40 days after anthesis (DAA) and dissected to separate the seed, the placenta and the flesh tissues. RNA extraction and small RNA library construction We isolated total RNA using Trizolate Reagent (UD Genomed Ltd., Debrecen, Hungary) from the different tissues. We determined the RNA concentration and 260/280 ratio with NanoDrop spectrophotometer (Thermo Scientific, Wilmington, DE) and verified the RNA integrity by agarose gel electrophoresis. We separated 30 μg of total RNA per sample on an 8% 3 / 20 denaturing polyacrylamide/urea gel along with RNA size markers, and then isolated the sRNA range. The pellets were disolved in 20 μl RNase-free ultrapure water, of which we used 2.5 μl for library preparation. We constructed the sRNA libraries using the Truseq Small RNA Sample Prep Kit (Illumina, CA, US), according to the manufacturer's instructions. We submitted the pools of eight PCR-amplified, gel-purified and bar-coded cDNAs to UDGenomed Ltd. (Debrecen, Hungary) for sequencing on Illumina HiSeq 2000 platform. Bioinformatic analysis We processed the the raw sequences using cutadapt 1.9.1. [ 32 ]. We performed a quality check before and after processing with FastQC 0.11.5 [ 33 ]. For miRNA analysis, we aligned the clean reads to the reference genome of Capsicum annuum cv. CM334 (version 1.55). We used MirCat [ 34 ] and miRDeep-P [ 35 ] in parallel to find new miRNAs; then we combined and filtered the two datasets by custom Python scripts. We used three subsets: miRNAs from the papers [ 25, 26 ] describing pepper genome and miRNA sequencing, miRNAs from MirProf results and miRNAs from mirCat and miRDeep-P predictions. We determined the raw expression values by Patman [36], and we used custom Python scripts to create a raw, and also a normalized abundance matrix (Tables A-C in S2 File). We performed the differential expression analysis with DESeq2 [ 37 ] on the raw abundance matrix. We got the basemean values and the shrunk foldchanges (in log2 form) from DESeq2. We considered those miRNAs differentially expressed, whose basemean value was at least 10 and shrunk fold-change was higher than 6. For the visualization of the results of the differential expression analysis, we created heat maps and MA-plots. We carried out the clustering of miRNAs into families with cd-hit EST [ 38 ] by setting the sequence identity cut-off to 0.8. For multiple alignments of the sequences of miRNA families, we used Clustal Omega, Color Align Conservation [ 39 ] and WebLogo [ 40 ] applications. We identified the sRNA producing genomic clusters with ShortStack 3.4 [ 41 ]. During the alignment, we allowed one mismatch and a maximum of 1000 possible placements for multimappable reads. To allow reporting the predominant RNA size in all loci, we set the parameters called dicermin and dicermax to 18 and 35, respectively. We kept all the other parameters at the default value. To identify differentially expressed loci between tissues and developmental stages, we applied DESeq2 with the filtering parameters used for miRNAs. We carried out the annotation of the loci with Patman using the Zunla RNA and TE sequences. Western blot analysis We homogenized tissue samples in an ice-cold mortar in 400 μl 2× Laemli buffer. Protein samples were resolved on 8% SDS-polyacrylamide gel and blotted to PVDF Blotting Membrane (AmershamTM HybondTM) with Trans-Blot Turbo (BIO-RAD) and subjected to Western blot analysis. We blocked the membrane with 5% non-fat milk powder in PBS containing Tween 20 (PBST) for 1 hour. We probed the membrane with anti-Arabidopsis thaliana AGO1 (Agrisera) in PBST containing 5% non-fat milk powder at room temperature for three hours. After washing the membrane in PBST, we added the secondary peroxidase-conjugated goat antirabbit IgG (H&L). The signals were visualized by chemiluminescence (ClarityTM Western ECL substrate; BIO-RAD) according to the manufacturer's instructions. We stained the membrane with Ponceau reagent to check even protein loading. Small RNA Northern blot analysis For the sRNA analysis, we separated 5 μg of total RNA on denaturing 12% polyacrylamide gel containing 8 M urea and transferred to HybondTM-NX membrane (Amersham) with Fast Blot (Biometra). RNA was chemically crosslinked to the membrane (2008 Pall and Hamilton). We 4 / 20 hybridized the membranes with radioactively labeled LNA or DNA oligonucleotide probes to detect sRNAs [ 42 ]. We end-labeled 10 pmol of each oligonucleotide probe with [γ-32P]ATP and T4 polynucleotide kinase. The hybridization was carried out at 37ÊC (DNA) or 50ÊC (LNA) depending on the type of oligonucleotides. We washed the membrane for 5 min two times with washing solution containing 0.1% SDS and 2× SSC at the temperature of the hybridization. We listed the used oligonucleotide probes in Table rr in S2 File. Results and discussion Small RNA sequencing and data pre-processing In this study, we used high-throughput sequencing method to analyze the tissue and development-specific expression of sRNAs in the fruit development of sweet pepper (Capsicum annuum cv. FeheÂroÈzoÈn). We collected samples at 28 and 40 DAA (Fig 1A). At these stages, the phase of rapid cell division is completed, and the development is already in the cell expansion phase. Since cell expansion accounts for the largest increase in fruit volume, the regulation of these developmental stages are extremely important from an economic point of view. To obtain tissue-specific sRNA expression profiles, we dissected the pepper fruit. Seed, placenta and flesh tissues were collected separately (Fig 1B and Figure A in S1 File) and processed for sRNA sequencing in two independent biological repeats. A total of more than 135 million reads were generated, in the range of 4.4±17.5 million reads per individual samples (Table 1). After removing adaptor sequences and filtering out low-quality tags, a total of 131.19 million clean reads were obtained, which were 16±28-nt in length. Next, we removed low abundance (< 2) reads, and those that aligned to known noncoding RNAs (rRNA, snRNA, snoRNA, tRNA), which resulted in 74.8 million filtered reads. More than 85% of these filtered reads matched the genome (version 1.55) of pepper cultivar CM334 [ 27 ] perfectly. The most critical filtering step was the removal of sequences with low abundance. The small numbers of sequences with an invalid length, low-quality tags, or low complexity indicate that the quality of our libraries was appropriate. Size distribution of the small RNA sequences The size distribution of all the filtered sequences also confirms the high-quality of sRNA sequencing data since more than 94% of the filtered sequences fall into the 21±24-nt size range, while the 16±20-nt and the 25±28-nt range represent only 4.8% and 1%, respectively (Fig 2). The 46.73% of the total sRNAs were 24-nt-long while 24.05% were 21-nt-long. The 22 and 23-nt-long sRNAs were the less abundant classes, 14%, and 10%, respectively (Fig 2). The analysis of the individual samples revealed that the 24-nt-long sRNA species was the dominant size class in the majority of the samples [ 13, 43 ]. However, we observed an increased ratio of 21-ntlong sRNAs in the flesh at 28 DAA and even more prominently at 40 DAA (Fig 2) indicating that the proportion of different size classes of sRNAs can vary significantly during the pepper fruit development. To analyze this phenomenon, we calculated the 21/24-nt ratio of the different samples and observed differences between various tissue types and developmental stages. The 21/24 ratio of the placenta samples increased only moderately, while in developing seeds it decreased from 55% (seed 28 DAA) to about 25% (seed 40 DAA). In the flesh, however, it increased from 75% (flesh 28 DAA) to about 200% (flesh 40 DAA). We found the ratio of the 21-nt-long sRNAs increasing during the expansion of the fruit except for the seed. The dynamic changes in the relative accumulation of various sRNA species with the advance of time indicate that sRNAs may play a pivotal role during fruit expansion. Moreover, the high percentage of the 21-nt-long sRNA species in the flesh suggests that miRNAs, phasiRNAs, or 5 / 20 Fig 1. Developmental stages of the sweet pepper cv. FeheÂroÈzoÈn investigated in this study. (A) The two selected fruit developmental stages at 28 and 40 days after anthesis (DAA). Grid = 1 cm. (B) Vertical and horizontal sections of expanding fruit showing the composing tissues as indicated. The seeds, the placenta, and the flesh were isolated and used for RNA extractions separately. other 21-nt-long sRNA molecules play a crucial role in the development of the economically most important tissue of the pepper fruit. 6 / 20 Filter 1: removal of reads with low abundance (1±2). Filter 2: removal of tRNA, rRNA, snRNA and snoRNA and low complexity sequences. Aligned: reads that matched the pepper genome. The numbers denote millions of reads. Fig 2. Size distribution of small RNAs in the sequencing libraries. The diagram shows the small RNA size composition of each library. Biological replicates are displayed next to each other. Numbers at the top indicate the length of sRNAs. PLOS ONE | https://doi.org/10.1371/journal.pone.0200207 7 / 20 Conserved miRNAs in the developing pepper fruit Two approaches have been applied to identify known miRNAs in the pepper fruit-specific sRNA libraries. MirProf has been used to profile known miRNAs by allowing two mismatches. MirProf uses Patman to align reads against miRBase [ 44 ]. MiRBase currently does not contain any pepper miRNAs, so we aligned our reads to the recently published conserved pepper miRNAs (can-mir) using Patman with the same parameters [ 25, 27 ]. We discarded sRNA reads, which were not present at least in two libraries (213 different sequences remained). However, this dataset contains some 18-nt-long sequences, which we did not consider for further analysis. Our final set of conserved miRNAs contained 193 non-redundant sequences, from 38 conserved miRNA families (Table A in S2 File). The majority of these miRNAs are in the mirBase database with two mismatches: 143 in Solanum tuberosum, 132 in Solanum lycopersicum, and 110 in Nicotiana tabacum. A recent comprehensive survey of plant miRNAs from 31 vascular plants species specified 82 known miRNA families that can be sorted into eight groups depending on their distribution across lineages and species [ 45 ]. The first group contains miRNA families that are widely present and highly expressed in all terrestrial species. In our samples, all of these miRNA families expressed at high level according to the RPM data (miR156: 249 RPM, miR166: 466210 RPM, miR167: 95 RPM, miR168: 4444 RPM), except for the miR172 family, which expressed at a lower level (17 RPM). The group 2 miRNA families are represented in all taxonomic lineages but can be absent or present with low abundance in some species (miR158, miR159, miR160, miR162, miR164, miR169, miR171, miR319, miR390, miR393, miR394, miR396, miR397, miR529, miR535 and miR4414). Most of the group 2 miRNA families were also present in our data sets except for the miR158, miR529, and miR535 families. Seven families (miR159, miR162, miR169, miR171, miR319, miR393, miR396) were abundant (more than 75 RPM) at least in one library. Other six families (miR160, miR164, miR390, miR394, miR397, miR4414) were moderately abundant (around 5±50 RPM). MiR158, miR529, and miR535 were absent from all the investigated pepper libraries. In a recent study investigating Nicotiana benthamiana sRNAs, miR529 and miR535 were absent from sRNA libraries [ 46 ]. One of the other six groups, group 8, was enriched in Solanaceae species, like Petunia, Nicotiana tabacum, Solanum lycopersicum, Solanum tuberosum and Capsicum annuum. From the ten miRNA families present in this group (miR1919, miR4376, miR5300, miR5301, miR6022, miR6023, miR6024, miR6025, miR6026, miR6027, miR6149), we identified 5 (miR4376, miR5300, miR5301, miR6024, miR6027) in our libraries. We checked the missing miRNA families in the libraries analyzed by Chavez and his colleagues, and some of them expressed weakly in their pepper samples, too. One example is mir1919, which is missing from all of their pepper libraries. Known pepper-specific miRNAs in the developing pepper fruit Recent studies described new pepper-specific miRNAs [ 25, 27, 28 ]. Based on published pepper miRNA sequences, we prepared a merged list of known pepper-specific miRNA sequences. The latest paper on pepper fruit miRNAs [ 28 ] described four sequences (Novel_miR265, Novel_miR218, Novel_miR23, Novel_miR279) from our set of new miRNAs, so we transferred these miRNAs to the known pepper-specific group, with the IDs given in that paper. The other published two sets of new miRNA sequences have only two common elements: can-miR-n010/can-miRC73p.1 and can-miR-n003a/can-miRC13-5p. We created a non-redundant list of 68 sequences and we found 50 of these sequences (including the two common ones) in our libraries (Table B in S2 File). We found more than 83% of the miRNAs described by Hwang et al. and 33% of the miRNAs described by Qin et al. [ 25, 27 ]. We did not detect some of the miRNAs, such as can-mirn05 and can-mir-n027, possibly because they are expressed in the root and the stem but not in 8 / 20 the fruit. However, we detected some known miRNAs which appeared to be flower- or leaf-specific in the studies mentioned above (i.e., can-mir-n019a-d and can-mir-n030). New miRNAs in developing pepper fruit To reveal still unknown miRNAs, we constructed a novel miRNA identification pipeline (Figure B in S1 File). 74.82 million filtered reads were used to identify new miRNAs by paralleled utilization of two computational plant miRNA prediction tools, miRDeep-P, and miRCat. Both of these tools map the reads to the genome and extract reference sequences by extending mapped reads in the flanking regions based on their ability to fold predicted miRNA precursor RNA stem-loop structures. These programs use the Vienna RNA package [ 47 ] to determine the secondary structure of the predicted precursors and apply a scoring system to rank the obtained potential miRNAs precursor structures. Although these programs use similar computational approaches, they can generate significantly different results. To enhance the reliability of our miRNA predictions, we used both programs to identify the potentially new pepper miRNAs and accepted only the overlapping set of different predictions. The next filtering steps were the elimination of candidate reads which did not reach the abundance level of 10 RPM and appeared at least in two sRNA libraries. Utilization of this method allowed us to describe 73 different, highly confident novel pepper miRNA sequences (Table C in S2 File). In all cases, we detected the miRNA sequence in at least one of our libraries; also their predicted precursors have a stable secondary structure (Table D in S2 File). We were able to cluster some of the novel miRNAs into six miRNA families, but the majority of them are singletons (Figure C in S1 File). We also checked the similarity between our newly described and the previously described miRNAs and found that five (can-mir-f11, can-mir-f14, can-mirf27, can-mir-f46, can-mir-f57) of our novel miRNAs were similar to previously described miRNAs [ 28 ]. These family members differ mainly in their sequence composition indicating that they probably derive from different precursors. We show the structures of some selected novel miRNA precursors in Fig 3. Validation of the miRNA expression profiles We validated the results by sRNA Northern blot analysis. We used DNA, or LNA oligonucleotide probes to detect the expression patterns of selected conserved miRNAs (miR171, miR172, miR396, miR159, and miR167) (Fig 4, green bars). The result of the Northern analysis was in a good agreement with the sequencing data, although in some cases we observed discrepancies. This could be attributed to the fact that the hybridization-based technology measures the abundance of heterogeneous miRNA population depending on the parameters of the hybridization procedure. Next, we analyzed the expression patterns of selected new pepper-specific miRNAs (Fig 4, red bars). Similarly to the expression of the conserved miRNAs, the experimental analyses revealed highly tissue-specific expression of the novel miRNAs. Altogether, these data indicate that the gained sRNA sequencing data are suitable for differential expression analysis of miRNAs during the developmental process. Differential miRNA expression in the pepper fruit To gain insights into the putative roles of the identified miRNAs, we conducted a differential expression analysis between the tissues and developmental stages using DESeq2. First, we determined the differentially expressed miRNAs between the tissues at 28 DAA (Fig 5). We detected 20, potentially influential conserved miRNAs from nine families, four new, and one known pepper-specific differentially expressed miRNAs in the placenta relative to the flesh and the seed (Tables E-F, N-O, and W-X in S2 File). The comparison of the placenta and the 9 / 20 Fig 3. Novel miRNA representatives. (A) Read alignment diagram for the selected miRNA precursors. Mature strand is highlighted with green and star strand with red. (B) Predicted secondary structure of pre-miRNAs (Mfold). Mature miRNA and miRNA sequences are highlighted. seed samples revealed less conservative (13), but more known pepper-specific (9) and new (13) differentially expressed miRNAs. We detected an increased number of changes between the flesh and the seed samples, where 44 miRNAs displayed significantly altered expression (Tables G, P, and Y in S2 File). Altogether, the majority of the differentially expressed miRNAs were 10 / 20 Fig 4. Small RNA Northern blot analyses of selected miRNAs. Expression of the conservative (green) and the newly predicted (red) miRNAs. The total RNA samples from the seed, the placenta and the flesh at 28 and 40 days after anthesis (DAA) were run on agarose gels, transferred to nylon membranes, and hybridized with probes detecting miRNAs as indicated. We show the NGS read counts (read/million) in the upper panels, the results of small RNA northern blots in the middle panels, and the rRNA levels as loading controls in the lower panels. We used the same membrane for the hybridization of miR171, miR396, and miR159, can-mir-f13 specific probes. conservative, and they show both up and downregulation in all comparisons. At this stage, differential expression of known and new pepper-specific miRNAs mostly related to their altered expression levels in the seed since fewer changes were detected between the placenta and flesh samples. At 40 DAA, we observed extensive changes in expression profiles of all miRNA classes between every tissue (Fig 5, Figure E in S1 File, Tables H-J, Q-S, and Z-bb in S2 File). The majority of miRNAs with altered accumulation expressed differentially in more than one comparisons. We found miRNAs in all classes, which showed significantly different levels of expression in all 11 / 20 Fig 5. Heat map showing the expression levels of the differentially expressed miRNAs. The differential expression analysis between the seed (S), the placenta (P), and the flesh (F) samples at (A) 28 DAA and (B) 40 DAA was conducted using DESeq2. The results were filtered to keep only those miRNAs whose mean normalized expression level was at least 10, and the expression change was at least six-fold. We colored the conservative (black), known (blue) and novel (red) miRNAs. We found 21 differentially expressed miRNAs between different stages (28, 40) in the flesh (C), 28 miRNAs in the seed (D), and none in the placenta. Arrows represent up- (") and down-regulation (#) compared to the two other tissues (for example S" means, the expression of the marked miRNA are upregulated in seed compared to both, placenta and flesh. We also listed those changes, which occurred only between two tissues. For example, we marked a miRNA up-regulation between the placenta and the seed like this: P!S". tested tissues (can-mir172a-5p, can-mir172f, can-miR319a, can-mir319bc, can-miR394b, canmiR-n035, can-mir-f04, and can-mir-f25). The differentially expressed conservative miRNAs at this stage (40 different mature sequences) are mainly downregulated in comparison of placenta 12 / 20 and flesh samples while in the other comparisons up-, and downregulation of miRNAs are mostly balanced. The known pepper specific miRNAs, similarly to the 28 DAA stage, are mostly upregulated specifically in the seed. The majority of the changes in the expression level of new miRNAs are also related to seed, but we identified expression differences between flesh and placenta samples as well. We also compared the miRNA expression profiles between 28 and 40 DAA in the different tissues (Fig 5C and 5D, and Figure F in S1 File). We observed significant differences in the case of the seed and the flesh samples (Tables L-M, U-V, and dd-ee in S2 File). The majority of the differentially expressed sequences were conservative miRNAs, and only limited changes were found in the expression of pepper-specific known and novel miRNAs. The placenta, however, showed practically no changes in miRNA expression profiles (Tables K, T, and cc in S2 File). The fact that most of these changes affected mostly the conservative miRNAs suggests evolutionary conservation of the miRNA-mediated regulation of the fruit tissue development during the expansion phase. Altogether, these data suggest that the investigated developmental stages are associated with specialized miRNA expression patterns and that the miRNA-mediated regulation is prevalent during the expansion phase. Expression of AGO1 protein in the pepper fruit tissues We found that the pepper fruit expansion is associated with profound temporal changes in the miRNA expression profile which can indicate the central role of miRNA-mediated regulation in this developmental process. AGO1 is the most important AGO protein in the miRNA pathway which is responsible for the cleavage or translational inhibition of target mRNAs determined by the loaded miRNA. AGO1 homeostasis is controlled by the action of miR168 on AGO1 mRNA [ 6 ]. Next, we investigated the accumulation of AGO1 protein, as the key executor component of RNAi, postulating that the accumulation level of AGO1 reflects the activity RNAi in the particular tissue types. First, we analyzed the accumulation of the AGO1 key regulator miR168 by sRNA Northern blot. We found, in agreement with the NGS data, that miR168 expresses abundantly in general with the highest level in the flesh samples at both time points (Fig 6). Using protein extracts representing the same samples we analyzed the expression of AGO1 with western blot. At 28 DAA no significant AGO1 accumulation was detected in the seed and the placenta samples at the detection level of the used method while in the flesh sample AGO1 accumulated at very high level. At 40 DAA we experienced drastic and moderate AGO1 induction in the seed and the placenta samples and importantly a continuously high AGO1 expression in the flesh. In accordance with the differential expression results, these data indicate that extensive miRNA expression changes are associated with the expansion of the pericarp. The importance of miRNAs in the expansion of the flesh is further supported by the observation that the proportions of the 21-nt-long sRNAs are highest in the flesh samples, especially at 40 DAA, where the 21/24-ratio becomes higher than one (Fig 2). In the placenta, the AGO1 level remains relatively low suggesting that temporal changes in the miRNA expressions are not required for the proper formation of this tissue. Indeed, the differential expression of the pattern of the placenta samples at the two investigated time points showed no changes (Figures D-F in S1 File). In the case of the seed, a very low AGO1 expression is detected at the early time point, however, in the later time point, the amount of the AGO1 protein drastically increased probably due to the increasing size of the embryo (Figure A in S1 File). We also found that the levels of miR168 in the samples do not correlate with the amount of AGO1 protein suggesting a potentially unusual regulation. This phenomenon requires further investigation. 13 / 20 Fig 6. Accumulation of AGO1 protein and miR168 in pepper fruit tissues. RNA and protein samples of seed, placenta, and flesh, representing developmental stages at 28 and 40 days after anthesis (DAA), were used for small RNA northern blot analyses and Western blotting. Upper panel is the schematic representation of miR168 NGS read counts (read/million). Middle panels show the results miR168 northern blots. EtBr stained rRNAs serve as loading control. Bottom panels show Western blot analyses of total protein extracts for AGO1 accumulation. Equal loading was verified by Ponceau staining of the membrane after western blotting. Expression of 24-nt-long hetsiRNA clusters The most abundant class of endogenous siRNAs are 24-nt-long and are predominantly produced from transposons and repeated sequence elements, especially heterochromatic regions [ 48 ]. Genome associated siRNAs are engaged in various nuclear processes such as developmental gene and transposon regulation, heterochromatin formation, and genome stability [ 49 ]. It has also been shown that 24-nt-long siRNAs play a crucial role in transgenerational epigenetic inheritance [ 50 ]. To analyze the expression pattern of the 24-nt siRNA-producing clusters during pepper fruit expansion, we investigated our datasets using ShortStack [ 41 ]. We identified 1,012,092 sRNA-producing clusters altogether. Most of the de novo identified siRNA clusters associated with the dominant sequence length of 24-nt (850,801; 84.1%). Of these clusters, 426,309 (42.1%) appeared to be associated with sequence repeats or mobile genetic elements. Others are localized to the promoters or gene bodies of protein-coding genes, probably affecting the methylation state of the DNA and consequently the transcriptional activity. We performed a differential expression analysis on the non-miRNA producing sRNA clusters using DESeq2 in the same way as in the case of miRNAs. Pairwise expression changes were graphically represented on MA-plots to visualize the distribution of expression values and differences (Figure G in S1 File). The associated information provided by ShortStack and DESeq2 for the significantly changed clusters are listed in Tables ff-nn in S2 File. We can observe a 14 / 20 similar trend in changes between tissues and time points as in the case of miRNAs. Most importantly, there are only a small number of significantly changed clusters in the placenta samples between the two time points, a moderate number of changing clusters between the placenta and the flesh samples, whereas the greatest changes can be observed between the seed and other tissues. Also, it is noteworthy that the expression levels of most of the sRNA clusters are higher in the seed samples at 40 DAA than at 28 DAA. This can be explained by the growth and differentiation of the embryonic tissue in the seed (Figure A in S1 File), suggesting that the 24-nt-long hetsiRNA-producing clusters are more active in this tissue to prevent genome instability during embryonic development. The detailed analyses of the top 50 most abundant 24-nt hetsiRNA-producing clusters revealed tissue-specific and temporally highly regulated expression patterns (Fig 7A, Table oo in S2 File). Some siRNA-producing clusters expressed at a very high level in the seed at 28 DAA, drastically reduced at 40 DAA, and almost completely silent in other tissues. To verify this remarkable tissue- and developmental stage-specific expression of the named clusters we investigated the expression of two representative clusters by sRNA Northern blot analyses (Fig 7C). Cluster_297979 is pepper-specific, covering approximately a 5 Kb-long intergenic region and expresses at a very high level. The locus is stranded but cannot be folded into a simple hairpin structure. A dominant 24-nt-long siRNA species is produced from the negative strand. Cluster_461037, approximately 1 Kb long, exhibits an opposite pattern of expression compared to Cluster_297979, showing enhanced expression at 40 DAA in the seed (Fig 7A). This locus was described earlier as a miRNA gene (can-miR-n008) [ 25 ], but in our dataset, it looks more like an inverted repeat-associated siRNA locus. The dominant 24-ntlong sequence contains one mismatch compared to the can-miR-n008 mature sequence. Furthermore, no star sequence was found in any of the libraries. Under our experimental conditions, an abundant 21-nt-long siRNA species is also produced. In line with the NGS data, we Fig 7. Expression pattern and sequence length composition of 24-nt-long hetsiRNA-producing clusters. Heat map showing (A) the log2 transformed expression levels (Read Per Million) and (B) sequence length composition of the 50 most abundant clusters, and (C) Northern blot analyses of the expression of Cluster 297979 and Cluster 461037. 15 / 20 detected the enhanced expression of both abundant sequence species in the seed at 28 DAA which declined at 40 DAA using DNA probe specific to the dominant siRNA species (Fig 7C). At this sensitivity level, in agreement with the NGS data, both probes detected no signals in the flesh and the placenta samples. Altogether, these data reveal a dynamic spatiotemporal regulation of the expression of hetsiRNAs in the developing pepper fruit. These observations are in accordance with previous Arabidopsis results demonstrating the dynamic expression of Pol IV-dependent siRNAs during silique development [ 51 ]. It has been shown that the Pol IV-dependent siRNA accumulation is initiated in the maternal gametophyte and continues in the seed, predominantly in the developing endosperm. Our observation that siRNA clusters are activated in a developmental phase-dependent and tissue-specific manner suggests that hetsiRNAs are potent regulatory components of the pepper fruit development. Expression of phased siRNA and secondary siRNA-producing clusters Phased secondary siRNAs are produced from Pol II-dependent transcripts that are cleaved by a predominantly 22-nt-long phase-initiating miRNA [ 52 ]. One exception is TAS3, which is cleaved by a 21-nt-long miR390 by a different mechanism [ 53 ]. This cleavage then initiates the synthesis of a complementary RNA by an RNA-dependent RNA polymerase (RDR6 in Arabidopsis), resulting in a double-stranded intermediate that is diced by the DCL4 endonuclease at every 21-nt starting from the miRNA cleavage site. ShortStack identified 93, and 948 potential phased clusters of 21-nt and 24 nt, respectively, the 50 most abundant are shown in Figures H and I in S1 File (the associated data are in Tables pp and qq in S2 File). The majority of the 21-nt PHAS loci overlap with genes coding for leucine-rich repeat-containing receptors, which a play a role in pathogen resistance and fruit ripening as well [ 54 ]. They have already been described as phasiRNA-producing genes in other species. The phasing is usually initialized by miRNAs belonging to the miR482/ miR2118 family or by other phasiRNAs. Surprisingly, we could not identify TAS homolog genes by BLAST alignment of Arabidopsis or other known plant TAS sequences to the Zunla RNA sequences, possibly because only the functional tasiRNA sequences are conserved, not the full noncoding transcripts. We searched for short matching patterns to the known tasiRNA sequences [ 55 ] using Patman [ 36 ]. This way we managed to identify three possible TAS3 homologue genes (ARF2-like, ARF3, ARF4), two TAS4 homologous genes (Potassium transporter 5 and Aquaporin PIP2-1-like) with microhomology only to the conserved tasiRNA sequences, and two TAS1-like genes (coding for AT4G04775-like uncharacterized protein and ATP-dependent DNA helicase Q-like 5) that can only be found in apple and peach, according to the tasiRNA database [ 55 ]. We did not find these potential TAS loci to produce phased siRNAs, probably because their known phase initializing miRNAs are hardly expressed in all of our libraries. Beside the above mentioned large protein family, there are single phasiRNA-producing genes, which are more specific to certain plants. It is worth mentioning that a DICER-LIKE PROTEIN 2-like gene is also a phased siRNA producing locus under our experimental conditions, which was also found to be a PHAS locus in Glycine max and Medicago truncatula [ 8 ] [ 56 ]. Moreover, this locus was differentially expressed between 28 and 40 DAA, in all tissues. DCL2 is a member of the sRNA processing machinery that plays a role in endogenous and virus-derived siRNA production. Some 24-nt PHAS loci have been described earlier [ 57, 58 ] and appear to be developmentally regulated. It is noteworthy that the most abundant 24-nt PHAS loci we identified are also developmental stage-specific, they are much more abundant in the seed sample at 28 DAA (Figure I in S1 File). According to the previous findings, these loci are associated with transposons, and their function is, similarly to the PIWI-associated siRNA in animals, to maintain genomic stability in reproductive tissues. 16 / 20 Conclusion In conclusion, our work provides comprehensive tissue-specific miRNA and siRNA expression landscape linked to expansion of pepper fruit. The separated analyses of fruit tissue components allowed us to dissect and list the sRNA components dominantly expressing in different tissue types in the expanding pepper fruit at two time points. Our analyses revealed several abundantly expressing tissue- and pepper-specific small regulatory RNA species implicating their regulatory role. Altogether, the presented data extend the potential perspectives of sRNA-mediated regulatory processes associated with the complex program of fruit development and provide fundamental data sets for further experiments. Supporting information S1 File. Supplementary Figures. Figures A-I. (PDF) (XLSX) S2 File. Supplementary Tables. Tables A-Z and Tables aa-rr. Acknowledgments helpful discussions. Author Contributions Conceptualization: ZoltaÂn Havelda. Data curation: PeÂter Gyula, Tibor Nagy. We thank ErzseÂbet PoldaÂn for her valuable assistance in the laboratory and EÂva VaÂrallyay for Formal analysis: DeÂnes Taller, PeÂter Gyula, Tibor Nagy, Endre Barta. Funding acquisition: ZoltaÂn Havelda. Investigation: PeÂter Gyula. Methodology: DeÂnes Taller, Jeannette BaÂlint, Tibor Nagy, Ivett Baksa, GyoÈrgy Szittya. Project administration: ZoltaÂn Havelda. Supervision: GyoÈrgy Szittya, JaÂnos Taller, ZoltaÂn Havelda. Validation: Jeannette BaÂlint. Visualization: DeÂnes Taller, PeÂter Gyula. Writing ± original draft: DeÂnes Taller, Jeannette BaÂlint, PeÂter Gyula, ZoltaÂn Havelda. Writing ± review & editing: Endre Barta, JaÂnos Taller. 17 / 20 18 / 20 2008; 18(10):1602±9. Epub 2008/07/26. https://doi.org/10.1101/gr.080127.108 PMID: 18653800; PubMed Central PMCID: PMC2556272. 19 / 20 1. Martinez de Alba AE , Elvira-Matelot E , Vaucheret H . Gene silencing in plants: a diversity of pathways . Biochimica et biophysica acta . 2013 ; 1829 (12): 1300 ± 8 . Epub 2013/11/05. https://doi.org/10.1016/j. bbagrm. 2013 . 10 .005 PMID: 24185199 . 2. Rogers K , Chen X . Biogenesis, turnover, and mode of action of plant microRNAs. The Plant cell . 2013 ; 25 ( 7 ): 2383 ± 99 . Epub 2013/07/25. https://doi.org/10.1105/tpc.113.113159 PMID: 23881412; PubMed Central PMCID : PMC3753372 . 3. Borges F , Martienssen RA . The expanding world of small RNAs in plants . Nature reviews Molecular cell biology . 2015 ; 16 ( 12 ): 727 ± 41 . Epub 2015/11/05. https://doi.org/10.1038/nrm4085 PMID: 26530390; PubMed Central PMCID : PMC4948178 . 4. Bologna NG , Voinnet O. The diversity, biogenesis, and activities of endogenous silencing small RNAs in Arabidopsis. Annual review of plant biology . 2014 ; 65 : 473 ± 503 . Epub 2014/03/04. https://doi.org/10. 1146/annurev-arplant- 050213 -035728 PMID: 24579988 . 5. Zhang H , Xia R , Meyers BC , Walbot V . Evolution, functions, and mysteries of plant ARGONAUTE proteins . Current opinion in plant biology . 2015 ; 27 : 84 ± 90 . Epub 2015/07/21. https://doi.org/10.1016/j.pbi. 2015 . 06 .011 PMID: 26190741 . 6. Vaucheret H. AGO1 homeostasis involves differential production of 21-nt and 22-nt miR168 species by MIR168a and MIR168b . PloS one . 2009 ; 4 ( 7 ): e6442 . Epub 2009 /08/04. https://doi.org/10.1371/journal. pone.0006442 PMID: 19649244; PubMed Central PMCID : PMC2714465 . 7. Chen HM , Chen LT , Patel K , Li YH , Baulcombe DC , Wu SH . 22 - Nucleotide RNAs trigger secondary siRNA biogenesis in plants . Proceedings of the National Academy of Sciences of the United States of America . 2010 ; 107 ( 34 ): 15269 ± 74 . Epub 2010/07/21. https://doi.org/10.1073/pnas.1001738107 PMID: 20643946; PubMed Central PMCID : PMC2930544 . 8. Fei Q , Xia R , Meyers BC . Phased, secondary, small interfering RNAs in posttranscriptional regulatory networks . The Plant cell . 2013 ; 25 ( 7 ): 2400 ± 15 . Epub 2013/07/25. https://doi.org/10.1105/tpc.113. 114652 PMID: 23881411; PubMed Central PMCID : PMC3753373 . 9. Dotto MC , Petsch KA , Aukerman MJ , Beatty M , Hammell M , Timmermans MC . Genome-wide analysis of leafbladeless1-regulated and phased small RNAs underscores the importance of the TAS3 ta-siRNA pathway to maize development . PLoS genetics . 2014 ; 10 ( 12 ): e1004826 . Epub 2014 /12/17. https://doi. org/10.1371/journal.pgen.1004826 PMID: 25503246; PubMed Central PMCID : PMC4263373 . 10. Kumar R . Role of microRNAs in biotic and abiotic stress responses in crop plants . Applied biochemistry and biotechnology . 2014 ; 174 ( 1 ): 93 ± 115 . Epub 2014/05/30. https://doi.org/10.1007/s12010-014-0914- 2 PMID: 24869742 . 11. Zhang B. MicroRNA: a new target for improving plant tolerance to abiotic stress . Journal of experimental botany . 2015 ; 66 ( 7 ): 1749 ± 61 . Epub 2015/02/24. https://doi.org/10.1093/jxb/erv013 PMID: 25697792; PubMed Central PMCID : PMC4669559 . 12. Li C , Zhang B . MicroRNAs in Control of Plant Development . Journal of cellular physiology . 2016 ; 231 ( 2 ): 303 ± 13 . Epub 2015/08/08. https://doi.org/10.1002/jcp.25125 PMID: 26248304 . 13. Axtell MJ . Classification and comparison of small RNAs from plants . Annual review of plant biology . 2013 ; 64 : 137 ± 59 . Epub 2013/01/22. https://doi.org/10.1146/annurev-arplant- 050312 -120043 PMID: 23330790 . 14. Qin Z , Li C , Mao L , Wu L . Novel insights from non-conserved microRNAs in plants . Frontiers in plant science . 2014 ; 5 : 586 . Epub 2014/11/13. https://doi.org/10.3389/fpls. 2014 .00586 PMID: 25389431; PubMed Central PMCID : PMC4211545 . 15. Kamthan A , Chaudhuri A , Kamthan M , Datta A . Small RNAs in plants: recent development and application for crop improvement . Frontiers in plant science . 2015 ; 6 : 208 . Epub 2015/04/18. https://doi.org/10. 3389/fpls. 2015 .00208 PMID: 25883599; PubMed Central PMCID : PMC4382981 . 16. Xu Q , Liu Y , Zhu A , Wu X , Ye J , Yu K , et al. Discovery and comparative profiling of microRNAs in a sweet orange red-flesh mutant and its wild type . BMC genomics . 2010 ; 11 : 246 . Epub 2010/04/20. https://doi.org/10.1186/ 1471 -2164-11-246 PMID: 20398412; PubMed Central PMCID : PMC2864249 . 17. Li H , Mao W , Liu W , Dai H , Liu Y , Ma Y , et al. Deep sequencing discovery of novel and conserved microRNAs in wild type and a white-flesh mutant strawberry . Planta . 2013 ; 238 ( 4 ): 695 ± 713 . Epub 2013/06/ 29. https://doi.org/10.1007/s00425-013 -1917-x PMID : 23807373 . 18. Ge A , Shangguan L , Zhang X , Dong Q , Han J , Liu H , et al. Deep sequencing discovery of novel and conserved microRNAs in strawberry (Fragariaxananassa) . Physiologia plantarum . 2013 ; 148 ( 3 ): 387 ± 96 . Epub 2012/10/16. https://doi.org/10.1111/j.1399- 3054 . 2012 . 01713 . x PMID : 23061771 . 19. Yanik H , Turktas M , Dundar E , Hernandez P , Dorado G , Unver T. Genome-wide identification of alternate bearing-associated microRNAs (miRNAs) in olive (Olea europaea L.) . BMC plant biology . 2013 ; 13 : 10 . Epub 2013/01/17. https://doi.org/10.1186/ 1471 -2229-13-10 PMID: 23320600; PubMed Central PMCID : PMC3564680 . 20. Mohorianu I , Schwach F , Jing R , Lopez-Gomollon S , Moxon S , Szittya G , et al. Profiling of short RNAs during fleshy fruit development reveals stage-specific sRNAome expression patterns . The Plant journal: for cell and molecular biology . 2011 ; 67 ( 2 ): 232 ± 46 . Epub 2011/03/30. https://doi.org/10.1111/j. 1365 - 313X . 2011 . 04586 . x PMID : 21443685 . 21. Moxon S , Jing R , Szittya G , Schwach F , Rusholme Pilcher RL , Moulton V , et al. Deep sequencing of tomato short RNAs identifies microRNAs targeting genes involved in fruit ripening . Genome research. 22. Zuo J , Zhu B , Fu D , Zhu Y , Ma Y , Chi L , et al. Sculpting the maturation, softening and ethylene pathway: the influences of microRNAs on tomato fruits . BMC genomics . 2012 ; 13 : 7 . Epub 2012/01/11. https://doi. org/10.1186/ 1471 -2164-13-7 PMID: 22230737; PubMed Central PMCID : PMC3266637 . 23. Karlova R , van Haarst JC , Maliepaard C , van de Geest H , Bovy AG , Lammers M , et al. Identification of microRNA targets in tomato fruit development using high-throughput sequencing and degradome analysis . Journal of experimental botany . 2013 ; 64 ( 7 ): 1863 ± 78 . Epub 2013/03/15. https://doi.org/10.1093/ jxb/ert049 PMID: 23487304; PubMed Central PMCID : PMC3638818 . 24. Gao C , Ju Z , Cao D , Zhai B , Qin G , Zhu H , et al. MicroRNA profiling analysis throughout tomato fruit development and ripening reveals potential regulatory role of RIN on microRNAs accumulation . Plant biotechnology journal . 2015 ; 13 ( 3 ): 370 ± 82 . Epub 2014/12/18. https://doi.org/10.1111/pbi.12297 PMID: 25516062 . 25. Hwang DG , Park JH , Lim JY , Kim D , Choi Y , Kim S , et al. The hot pepper (Capsicum annuum) microRNA transcriptome reveals novel and conserved targets: a foundation for understanding MicroRNA functional roles in hot pepper . PloS one . 2013 ; 8 ( 5 ): e64238 . Epub 2013 /06/06. https://doi.org/10.1371/ journal.pone.0064238 PMID: 23737975; PubMed Central PMCID : PMC3667847 . 26. Kim S , Park M , Yeom SI , Kim YM , Lee JM , Lee HA , et al. Genome sequence of the hot pepper provides insights into the evolution of pungency in Capsicum species . Nature genetics . 2014 ; 46 ( 3 ): 270 ± 8 . Epub 2014/01/21. https://doi.org/10.1038/ng.2877 PMID: 24441736 . 27. Qin C , Yu C , Shen Y , Fang X , Chen L , Min J , et al. Whole-genome sequencing of cultivated and wild peppers provides insights into Capsicum domestication and specialization . Proceedings of the National Academy of Sciences of the United States of America . 2014 ; 111 ( 14 ): 5135 ± 40 . Epub 2014/03/05. https://doi.org/10.1073/pnas.1400975111 PMID: 24591624; PubMed Central PMCID : PMC3986200 . 28. Liu Z , Zhang Y , Ou L , Kang L , Liu Y , Lv J , et al. Identification and characterization of novel microRNAs for fruit development and quality in hot pepper (Capsicum annuum L.) . Gene . 2017 ; 608 : 66 ± 72 . Epub 2017/01/26. https://doi.org/10.1016/j.gene. 2017 . 01 .020 PMID: 28122266 . 29. Nag A , Jack T. Sculpting the flower; the role of microRNAs in flower development . Current topics in developmental biology . 2010 ; 91 : 349 ± 78 . Epub 2010/08/14. https://doi.org/10.1016/S0070- 2153 ( 10 ) 91012 - 0 PMID: 20705188 . 30. Garcia D. A miRacle in plant development: role of microRNAs in cell differentiation and patterning . Seminars in cell & developmental biology . 2008 ; 19 ( 6 ): 586 ± 95 . Epub 2008/08/19. https://doi.org/10.1016/j. semcdb. 2008 . 07 .013 PMID: 18708151 . 31. Valoczi A , Varallyay E , Kauppinen S , Burgyan J , Havelda Z. Spatio-temporal accumulation of microRNAs is highly coordinated in developing plant tissues . The Plant journal: for cell and molecular biology . 2006 ; 47 ( 1 ): 140 ± 51 . Epub 2006/07/11. https://doi.org/10.1111/j. 1365 - 313X . 2006 . 02766 . x PMID : 16824182 . 32. Martin M. Cutadapt removes adapter sequences from high-throughput sequencing reads . EMBnetjournal . 2011 ; 17 :pp. 10 ± 2 . https://doi.org/10.14806/ej.17.1. 200 33. Andrews S. FastQC: a quality control tool for high throughput sequence data 2010 [updated 2010] . 34. Stocks MB , Moxon S , Mapleson D , Woolfenden HC , Mohorianu I , Folkes L , et al. The UEA sRNA workbench: a suite of tools for analysing and visualizing next generation sequencing microRNA and small RNA datasets . Bioinformatics . 2012 ; 28 ( 15 ): 2059 ± 61 . Epub 2012/05/26. https://doi.org/10.1093/ bioinformatics/bts311 PMID: 22628521; PubMed Central PMCID : PMC3400958 . 35. Yang X , Li L . miRDeep-P: a computational tool for analyzing the microRNA transcriptome in plants . Bioinformatics . 2011 ; 27 ( 18 ): 2614 ± 5 . Epub 2011/07/22. https://doi.org/10.1093/bioinformatics/btr430 PMID: 21775303 . 36. Prufer K , Stenzel U , Dannemann M , Green RE , Lachmann M , Kelso J. PatMaN : rapid alignment of short sequences to large databases . Bioinformatics . 2008 ; 24 ( 13 ): 1530 ± 1 . Epub 2008/05/10. https:// doi.org/10.1093/bioinformatics/btn223 PMID: 18467344; PubMed Central PMCID : PMC2718670 . 37. Love MI , Huber W , Anders S. Moderated estimation of fold change and dispersion for RNA-seq data with DESeq2 . Genome Biol . 2014 ; 15 ( 12 ): 550 . Epub 2014/12/18. https://doi.org/10.1186/s13059-014- 0550-8 PMID: 25516281; PubMed Central PMCID : PMC4302049 . 38. Huang Y , Niu B , Gao Y , Fu L , Li W. CD-HIT Suite: a web server for clustering and comparing biological sequences . Bioinformatics . 2010 ; 26 ( 5 ): 680 ± 2 . Epub 2010/01/08. https://doi.org/10.1093/ bioinformatics/btq003 PMID: 20053844; PubMed Central PMCID : PMC2828112 . 39. Stothard P. The sequence manipulation suite: JavaScript programs for analyzing and formatting protein and DNA sequences . BioTechniques . 2000 ; 28 ( 6 ): 1102 , 4 . Epub 2000/06/27. PMID: 10868275 . 40. Crooks GE , Hon G , Chandonia JM , Brenner SE . WebLogo: a sequence logo generator . Genome research . 2004 ; 14 ( 6 ): 1188 ± 90 . Epub 2004/06/03. https://doi.org/10.1101/gr.849004 PMID: 15173120; PubMed Central PMCID : PMC419797 . 41. Axtell MJ . ShortStack: comprehensive annotation and quantification of small RNA genes . Rna . 2013 ; 19 ( 6 ): 740 ± 51 . Epub 2013/04/24. https://doi.org/10.1261/rna.035279.112 PMID: 23610128; PubMed Central PMCID : PMC3683909 . 42. Varallyay E , Burgyan J , Havelda Z. MicroRNA detection by northern blotting using locked nucleic acid probes . Nature protocols . 2008 ; 3 ( 2 ): 190 ± 6 . Epub 2008/02/16. https://doi.org/10.1038/nprot. 2007 .528 PMID: 18274520 . 43. Parent JS , Martinez de Alba AE , Vaucheret H. The origin and effect of small RNA signaling in plants. Frontiers in plant science . 2012 ; 3 : 179 . Epub 2012/08/22. https://doi.org/10.3389/fpls. 2012 .00179 PMID: 22908024; PubMed Central PMCID : PMC3414853 . 44. Kozomara A , Griffiths-Jones S. miRBase: annotating high confidence microRNAs using deep sequencing data . Nucleic acids research . 2014 ; 42 (Database issue): D68±73. Epub 2013 /11/28. https://doi.org/ 10.1093/nar/gkt1181 PMID: 24275495; PubMed Central PMCID : PMC3965103 . 45. Chavez Montes RA , de Fatima Rosas-Cardenas F , De Paoli E , Accerbi M , Rymarquis LA , Mahalingam G , et al. Sample sequencing of vascular plants demonstrates widespread conservation and divergence of microRNAs . Nature communications . 2014 ; 5 : 3722 . Epub 2014/04/25. https://doi.org/10.1038/ ncomms4722 PMID: 24759728 . 46. Baksa I , Nagy T , Barta E , Havelda Z , Varallyay E , Silhavy D , et al. Identification of Nicotiana benthamiana microRNAs and their targets using high throughput sequencing and degradome analysis . BMC genomics . 2015 ; 16 : 1025 . Epub 2015/12/03. https://doi.org/10.1186/s12864-015-2209-6 PMID: 26626050; PubMed Central PMCID : PMC4667520 . 47. Lorenz R , Bernhart SH , Siederdissen CHZ , Tafer H , Flamm C , Stadler PF , et al. ViennaRNA Package 2.0. Algorithm Mol Biol . 2011 ; 6 . https://doi.org/Artn 2610.1186/ 1748 -7188-6-26 PMID: ISI: 000302425300001 . 48. Kasschau KD , Fahlgren N , Chapman EJ , Sullivan CM , Cumbie JS , Givan SA , et al. Genome-wide profiling and analysis of Arabidopsis siRNAs . PLoS biology . 2007 ; 5 ( 3 ): e57 . Epub 2007 /02/15. https://doi. org/10.1371/journal.pbio.0050057 PMID: 17298187; PubMed Central PMCID : PMC1820830 . 49. Castel SE , Martienssen RA . RNA interference in the nucleus: roles for small RNAs in transcription, epigenetics and beyond . Nature reviews Genetics . 2013 ; 14 ( 2 ): 100 ± 12 . Epub 2013/01/19. https://doi.org/ 10.1038/nrg3355 PMID: 23329111; PubMed Central PMCID : PMC4205957 . 50. Bond DM , Baulcombe DC . Small RNAs and heritable epigenetic variation in plants . Trends in cell biology . 2014 ; 24 ( 2 ): 100 ± 7 . Epub 2013/09/10. https://doi.org/10.1016/j.tcb. 2013 . 08 .001 PMID: 24012194 . 51. Mosher RA , Melnyk CW , Kelly KA , Dunn RM , Studholme DJ , Baulcombe DC . Uniparental expression of PolIV-dependent siRNAs in developing endosperm of Arabidopsis . Nature . 2009 ; 460 ( 7252 ): 283 ± 6 . Epub 2009/06/06. https://doi.org/10.1038/nature08084 PMID: 19494814 . 52. Cuperus JT , Carbonell A , Fahlgren N , Garcia-Ruiz H , Burke RT , Takeda A , et al. Unique functionality of 22-nt miRNAs in triggering RDR6-dependent siRNA biogenesis from target transcripts in Arabidopsis . Nature structural & molecular biology . 2010 ; 17 ( 8 ): 997 ± 1003 . Epub 2010/06/22. https://doi.org/10. 1038/nsmb.1866 PMID: 20562854; PubMed Central PMCID : PMC2916640 . 53. Montgomery TA , Howell MD , Cuperus JT , Li D , Hansen JE , Alexander AL , et al. Specificity of ARGONAUTE7-miR390 interaction and dual functionality in TAS3 trans-acting siRNA formation . Cell . 2008 ; 133 ( 1 ): 128 ± 41 . Epub 2008/03/18. https://doi.org/10.1016/j.cell. 2008 . 02 .033 PMID: 18342362 . 54. Li B , Yan J , Jia W. FERONIA/FER-like receptor kinases integrate and modulate multiple signaling pathways in fruit development and ripening . Plant signaling & behavior . 2017 ; 12 ( 12 ): e1366397 . Epub 2017 / 12/08. https://doi.org/10.1080/15592324. 2017 .1366397 PMID: 29215944 . 55. Zhang C , Li G , Zhu S , Zhang S , Fang J. tasiRNAdb: a database of ta-siRNA regulatory pathways . Bioinformatics . 2014 ; 30 ( 7 ): 1045 ± 6 . Epub 2013/12/29. https://doi.org/10.1093/bioinformatics/btt746 PMID: 24371150 . 56. Fei Q , Li P , Teng C , Meyers BC . Secondary siRNAs from Medicago NB-LRRs modulated via miRNAtarget interactions and their abundances . The Plant journal: for cell and molecular biology . 2015 ; 83 ( 3 ): 451 ± 65 . Epub 2015/06/05. https://doi.org/10.1111/tpj.12900 PMID: 26042408 . 57. Song X , Li P , Zhai J , Zhou M , Ma L , Liu B , et al. Roles of DCL4 and DCL3b in rice phased small RNA biogenesis . The Plant journal: for cell and molecular biology . 2012 ; 69 ( 3 ): 462 ± 74 . Epub 2011/10/07. https://doi.org/10.1111/j. 1365 - 313X . 2011 . 04805 . x PMID : 21973320 . 58. Zhai J , Zhang H , Arikit S , Huang K , Nan GL , Walbot V , et al. Spatiotemporally dynamic, cell-type-dependent premeiotic and meiotic phasiRNAs in maize anthers . Proceedings of the National Academy of Sciences of the United States of America . 2015 ; 112 ( 10 ): 3146 ± 51 . Epub 2015/02/26. https://doi.org/10. 1073/pnas.1418918112 PMID: 25713378; PubMed Central PMCID : PMC4364226 .

This is a preview of a remote PDF: http://journals.plos.org/plosone/article/file?id=10.1371/journal.pone.0200207&type=printable

Dénes Taller, Jeannette Bálint, Péter Gyula, Tibor Nagy, Endre Barta, Ivett Baksa, György Szittya, János Taller, Zoltán Havelda. Expansion of Capsicum annum fruit is linked to dynamic tissue-specific differential expression of miRNA and siRNA profiles, PLOS ONE, 2018, DOI: 10.1371/journal.pone.0200207