Advertisement
Browse Subject Areas
?

Click through the PLOS taxonomy to find articles in your field.

For more information about PLOS Subject Areas, click here.

  • Loading metrics

Open Access

Peer-reviewed

Research Article

Identification of microRNAs and their targets in inflorescences of an Ogura-type cytoplasmic male-sterile line and its maintainer fertile line of turnip (Brassica rapa ssp. rapifera) via high-throughput sequencing and degradome analysis

Abstract

Cytoplasmic male sterility (CMS) is a widely used trait in angiosperms caused by perturbations in nucleus-mitochondrion interactions that suppress the production of functional pollen. MicroRNAs (miRNAs) are small non-coding RNAs that act as regulatory molecules of transcriptional or post-transcriptional gene silencing in plants. The discovery of miRNAs and their possible implications in CMS induction provides clues for the intricacies and complexity of this phenomenon. Previously, we characterized an Ogura-CMS line of turnip (Brassica rapa ssp. rapifera) that displays distinct impaired anther development with defective microspore production and premature tapetum degeneration. In the present study, high-throughput sequencing was employed for a genome-wide investigation of miRNAs. Six small RNA libraries of inflorescences collected from the Ogura-CMS line and its maintainer fertile (MF) line of turnip were constructed. A total of 120 pre-miRNAs corresponding to 89 mature miRNAs were identified, including 87 conversed miRNAs and 33 novel miRNAs. Among these miRNAs, the expression of 10 differentially expressed mature miRNAs originating from 12 pre-miRNAs was shown to have changed by more than two-fold between inflorescences of the Ogura-CMS line and inflorescences of the MF line, including 8 down- and 2 up-regulated miRNAs. The expression profiles of the differentially expressed miRNAs were confirmed by stem-loop quantitative real-time PCR. In addition, to identify the targets of the identified miRNAs, a degradome analysis was performed. A total of 22 targets of 25 miRNAs and 17 targets of 28 miRNAs were identified as being involved in the reproductive development for Ogura-CMS and MF lines of turnip, respectively. Negative correlations of expression patterns between partial miRNAs and their targets were detected. Some of these identified targets, such as squamosa promoter-binding-like transcription factor family proteins, auxin response factors and pentatricopeptide repeat-containing proteins, were previously reported to be involved in reproductive development in plants. Taken together, our results can help improve the understanding of miRNA-mediated regulatory pathways that might be involved in CMS occurrence in turnip.

Citation: Lin S, Su S, Jin L, Peng R, Sun D, Ji H, et al. (2020) Identification of microRNAs and their targets in inflorescences of an Ogura-type cytoplasmic male-sterile line and its maintainer fertile line of turnip (Brassica rapa ssp. rapifera) via high-throughput sequencing and degradome analysis. PLoS ONE 15(7): e0236829. https://doi.org/10.1371/journal.pone.0236829

Editor: Yun Zheng, Kunming University of Science and Technology, CHINA

Received: April 14, 2020; Accepted: July 14, 2020; Published: July 30, 2020

Copyright: © 2020 Lin et al. This is an open access article distributed under the terms of the Creative Commons Attribution License, which permits unrestricted use, distribution, and reproduction in any medium, provided the original author and source are credited.

Data Availability: All of the small RNA sequencing data and degradome sequencing data were submitted to the Sequence Read Archive of the NCBI under accession number PRJNA552762 (DOI: http://www.ncbi.nlm.nih.gov/sra/PRJNA552762).

Funding: This work was jointly supported by funds from the Special Project on Science and Technology Innovation of Seed and Seedling of Wenzhou[Z20160008] (SL), and the National Natural Science Foundation of China [31972418] (SL).

Competing interests: The authors have declared that no competing interests exist.

Abbreviations: A, adenine; AP2, APETALA 2; ARF, auxin response factor; CMS, cytoplasmic male sterility; COG, Clusters of Orthologous Group; CUC1, CUP-SHAPED COTYLEDON 1; eggNOG, orthologous groups of genes; FC, fold change; FDR, false discovery rate; FPKMs, fragments per kilobase of transcript per million mapped reads; GMUCT, genome-wide mapping of uncapped and cleaved transcripts; GRF, growth-regulating factor; KEGG, Kyoto Encyclopedia of Genes and Genomes; MF, maintainer fertile; MFEI, minimum folding energy indices; MFEs, minimum free energies; miRNAs, microRNAs; NCBI, National Center for Biotechnology Information; ncRNA, non-coding RNA; Nr, NCBI non-redundant protein; nt, nucleotide; Pfam, Homologous protein family; PPR-containing protein, pentatricopeptide repeat-containing protein; qRT-PCR, quantitative real-time PCR; rRNA, ribosomal RNA; SBP, squamosal-promoter binding protein; SD, standard deviation; snoRNA, small nucleolar RNA; snRNA, small nuclear RNA; TPM, transcripts per million; tRNA, transfer RNA; U, uracil

Introduction

Cytoplasmic male sterility (CMS) is a maternally inherited trait that is common across angiosperms and is caused by the interaction of a nuclear fertility restorer gene and a mitochondrial CMS gene, resulting in the inability to produce functional pollen [1]. CMS halts pollen development at almost all developmental stages (before, during and after meiosis of pollen mother cells) in higher plants [2]. To date, several CMS systems have been characterized, such as Ogura, Polima, Kos and NWB systems, among which there are many intra- and inter-specific variations within each cytoplasmic genotype [3–8]. Regardless of the type, the mechanism of CMS is believed to arise as a consequence of nucleus-cytoplasm incompatibility [1,2].

In most previous studies on CMS, the focus has been mainly on CMS-inducing mitochondrial genes and fertility restorer genes [9–15]. Mitochondrial dysfunction is considered to be attributed to unusual open reading frames (orfs) in the mitochondrial genome [10]. Many CMS-associated orfs have been identified and characterized in various plant species. For instance, orf138 is associated with Ogura-CMS in Brassiceae, the co-transcription of orf355 and orf77 is responsible for maize CMS-S, and orf224 co-transcribes together with atp6 and down-regulates pollen generation in Polima-CMS of Brassica napus [16]. The effects of mitochondrial genes governing CMS can be suppressed by nuclear restorer-of-fertility (Rf) genes, which mostly encode pentatricopeptide repeat (PPR) proteins [12]. To date, numerous and variable Rfs or Rf candidates have been identified [17]. It is well established that the CMS phenotype is due to incompatibility resulting from the combination of mitochondrial dysfunction and a lack of Rf genes [1,2]. However, one fact that cannot be ignored is that nuclear genes affect the expression and function of cytoplasmic genes, as evidenced by the presence of nuclear fertility restorers, and mutant mitochondrial gene expression also contributes to retrograde regulation that fine-tunes nuclear gene expression, thus affecting pollen development [2,17–20]. Mitochondrial retrograde regulation especially manifests with variation in the stages at which pollen/anther abortion occur and the inconsistency in microsporogenesis and tapetum development processes with specific mitochondrial genes in different nuclear backgrounds [21–24]. All these findings have directed attention in recent years to studies of nucleus-cytoplasm interactions at the whole-genome level.

MicroRNAs (miRNAs) are a distinct class of single-stranded, short (~20–24 nucleotides in length), endogenously expressed non-coding RNAs (ncRNAs) that are widespread throughout the plant kingdom [25]. These small molecules are processed by Dicer-like proteins from longer RNA precursors that can form stem-loop regions [26]. Mature miRNAs are incorporated into the RNA-induced silencing complex to target their respective complementary or nearly complementary mRNAs and then act as inhibitory signals that direct mRNA cleavage or trigger translational repression [27]. MiRNAs emerge as small regulatory molecules of vital plant developmental processes, from vegetative growth to reproduction and stress responses [28–34].

Advancements in miRNA arrays and sequencing technology have led to the identification of an increasing number of miRNAs and precursor miRNAs in pollen [35–40]. In Arabidopsis, unexpectedly diverse miRNA populations belonging to 33 different miRNA families were detected in mature pollen grains, most of which displayed an enriched expression pattern in pollen [37]. There is an increasing amount of evidence confirming that some miRNAs, such as miR159 and miR167, are master modulators of plant male sterility [41–44]. High-throughput sequencing and degradome analysis have highlighted the differential expression of miRNAs and their respective targets between CMS and fertile lines in many species [45–51]. In B. juncea, 47 miRNAs were differentially expressed between CMS and fertile lines, 101 in soybean, 42 in cybrid pummelo (Citrus grandis), and 87 in Chinese cabbage [33,45,48]. Some of these differentially expressed miRNAs were predicted to target transcription factor family proteins or functional proteins with potential roles in male gametophyte development. For instance, two novel miRNAs (novel-miR-448 and novel-miR-335) highly expressed in CMS buds of Chinese cabbage were confirmed to significantly suppress the expression of sucrose transporter SUC1 and H+-ATPase 6, which perform essential roles in pollen development [48]. The discovery of miRNAs and their roles in the regulation of gene expression during pollen development shed light on the possible connection between miRNA action and CMS phenotypes. However, the regulatory network of CMS occurrence, especially the understanding of the involvement of miRNAs in CMS, is still limited.

Previously, an Ogura-CMS line ‘BY10-2A’ and its maintainer fertile (MF) line ‘BY10-2B’ of turnip (B. rapa ssp. rapifera) were characterized [52]. Mutation of mitochondrial orf138 retro-regulates the expression of nuclear genes, and interactions between them are responsible for male sterility in Ogura-CMS turnip. The first sign of disintegration shown by the anthers of Ogura-CMS line is that tapetum swells at the center of the locule during the transition from the microspore mother cells to tetrads, leading to failure of microspore development and thus complete male sterility. Using RNA sequencing analysis and bioinformatics, a large number of differentially expressed genes have been identified, which make good candidates for CMS-related genes. In this study, we further investigated the expression of miRNAs and their targets in inflorescences of the Ogura-CMS line and its MF line of turnip by high-throughput sequencing and degradome analysis. The current impact of our study is that the biogenesis of miRNAs could be regulated during retrograde signaling involved in CMS, revealing potential roles of miRNAs and their targets in regulating anther development in turnip.

Materials and methods

Plant materials, sample collection and total RNA isolation

An Ogura-CMS line of turnip was developed previously through consecutive back-crossing and inter-specific hybridization with B. rapa ssp. chinensis constituting the Ogura-CMS cytoplasm donor and the fertile turnip serving as the recurrent parent. The stable Ogura-CMS turnip line ‘BY10-2A’ and its corresponding MF line ‘BY10-2B’ were planted at the experimental farm of Wenzhou Vocational College of Science and Technology, Wenzhou, Zhejiang, China. All floral buds of an inflorescence from the Ogura-CMS and MF lines of turnip were collected after flowering. In each case, samples were harvested from ten individual plants and pooled, with transcriptome profiles representing ‘f’ differences. Three biological replicates were performed. The samples were then immediately frozen in liquid nitrogen and stored at -70°C until RNA isolation. Total RNA was extracted using Trizol reagent (Invitrogen, CA, USA) according to the manufacturer’s protocol. A NanoPhotometer spectrophotometer (IMPLEN, CA, USA), a Qubit RNA Assay Kit in conjunction with a Qubit2.0 Fluorometer (Life Technologies, CA, USA), and an Agilent Bioanalyzer 2100 system (Agilent Technologies, CA, USA) were used to determine the RNA purity, RNA concentration, and RNA integrity, respectively, to ensure equivalent amounts of RNA samples were used for sequencing.

Small RNA library preparation and sequencing

A total amount of 1.5 μg of RNA per sample was used as input material for small RNA library preparations using a NEBNext Ultra RNA Library Prep Kit for Illumina (NEB, CA, USA) according to the manufacturer’s recommendations. Briefly, small RNA was first ligated to 3’ and 5’ RNA adapters. Second, the adapter-ligated RNAs were transcribed into cDNA, after which PCR was performed. Afterward, 18–30 nucleotide (nt) fragments were selected and screened via PAGE. Last, the PCR products were purified and the library was prepared. The library preparations were then used for cluster generation on a cBot Cluster Generation System of Illumina and subsequently sequenced on an Illumina HiSeq 2500 platform by Biomarker Biotechnology Corporation (Beijing, China).

Identification of known and novel miRNAs

Clean data were obtained by removing reads containing adapters, reads containing poly-N sequences, reads of low-quality, and reads smaller than 18 nt or longer than 30 nt from the raw data. All downstream analyses were performed based on high-quality clean data as assessed by the Q20 value, Q30 value, GC-content and sequence duplication level. Ribosomal RNA (rRNA), transfer RNA (tRNA), small nuclear RNA (snRNA), small nucleolar RNA (snoRNA), other ncRNAs and repeats were then filtered and removed from the clean reads using Bowtie software [53] in conjunction with the Silva database, GtRNAdb, the Rfam database and the Repbase database. The remaining unannotated reads were then matched with assembled mRNA sequences uploaded to the Sequence Read Archive of the National Center for Biotechnology Information (NCBI) (https://www.ncbi.nlm.nih.gov/sra/PRJNA505114; accession number PRJNA505114). The matched sequences were used to detect miRNAs predicted by comparisons with previously known pre-miRNA sequences of all plants in the miRBase database (version 21.0) by miRDeep2 software [54,55]. Complete alignment of the sequences and hits with zero mismatches were considered as candidate conserved miRNAs, while the others were reserved as candidate novel miRNAs. Randfold software (version 2.1.7) and the RNAfold web server (http://rna.tbi.univie.ac.at/cgi-bin/RNAWebSuite/RNAfold.cgi) were used for secondary structure predictions of putative pre-miRNAs [54,56].

Differential expression analysis of miRNAs

To estimate miRNA expression levels, the small RNAs were mapped back to the precursor sequences, and the read count for each miRNA was obtained from the mapping results. The frequency of miRNAs in the six libraries was normalized to the expression of transcripts per million (TPM; normalized expression = Read count × 1,000,000/Mapped reads) as proposed previously [57]. Differentially expressed miRNAs in the inflorescences of the Ogura-CMS line and its MF line were screened via the DESeq (2010) R Package. MiRNAs with a Benjamini-Hochberg false discovery rate (FDR) ≤ 0.05 and a |log2 fold change (FC)| ≥ 1 were considered as differentially expressed.

Degradome sequencing and analysis

Six degradome libraries were constructed from the inflorescences of the Ogura-CMS line and its MF line based on the genome-wide mapping of uncapped and cleaved transcripts (GMUCT) version 2.0 method as described previously [58]. In brief, poly-A-enriched mRNA was ligated with a 5’-RNA adapter, and the 5’-ligated RNA was then purified and separated from the unligated adapter by performing a second poly-A selection. Reverse transcription was performed using a primer that was a random hexamer fused to the 3’-adapter, allowing for the adapter to be added, followed by PCR amplification. The amplified products were subsequently gel purified and subjected to deep sequencing on an Illumina HiSeq 2500 platform by Biomarker Biotechnology Corporation (Beijing, China). Clean data were obtained by trimming 5’ and 3’ adapters and filtering low-quality reads. The clean tags were mapped to assembled mRNA sequences (https://www.ncbi.nlm.nih.gov/sra/PRJNA505114; accession number PRJNA505114) and then annotated with rRNA, tRNA, snRNA, snoRNA, and other ncRNA from the Rfam database to obtain the unannotated tags that were used to predict subsequent degradation sites. The statistical data of the cleaved sites and target plots (t-plots) were analyzed by CleaveLand 3.0 (http://axtell-lab-psu.weebly.com/cleaveland.html) [59]. The transcript abundance was plotted for each transcript, and the heights of degradome peaks at each sliced-target transcript position were grouped into five categories according to the relative abundance of tags at the target sites [45,50]. MiRNA target sequences were searched against the NCBI non-redundant protein (Nr) database, the Swiss-Prot database and the evolutionary genealogy of genes: Non-supervised Orthologous Groups (eggNOG) database, and aligned to sequences within the Clusters of Orthologous Group (COG) database, Kyoto Encyclopedia of Genes and Genomes (KEGG) database, and Homologous protein family (Pfam) database to predict and classify their functions.

Stem-loop and regular quantitative real-time PCR (qRT-PCR) of miRNAs and their targets

Fragments per kilobase of transcript per million mapped reads (FPKMs) expression data of miRNA targets were generated from our previous study using RNA sequencing [52]. The expression levels of miRNAs and their targets were validated by stem-loop qRT-PCR [60] and qRT-PCR, respectively. First, miRNAs were extracted from the residual plant samples for small RNA library preparations using a PureLink miRNA Isolation Kit (Thermo Fisher Scientific). For each miRNA, 1 μg of the RNA sample was reverse-transcribed with SuperScript III Reverse Transcriptase (Thermo Fisher Scientific) using a specific stem-loop reverse transcription primer (S1 Table). The reactions were then incubated for 5 min at 25°C, followed by 50°C for 15 min, and a final incubation at 85°C for 5 min. The cDNA template for the miRNA targets was reverse-transcribed using SuperScript III First-Strand Synthesis SuperMix (Invitrogen, CA, USA). qRT-PCR was subsequently performed using PowerUp SYBR Green Master Mix (Applied Biosystems) on a Bio-Rad CFX384 instrument (Bio-Rad, http://www.bio-rad.com). 5S rRNA was used as an internal standard for miRNA analysis, and ACT7 was used as an endogenous control for miRNA target analysis [61]. The melting curve was determined for the amplification specificity, and the amplified products were confirmed by sequencing. Three biological replicates and three independent technical replicates were assessed. The 2-ΔΔCt method was applied to calculate the relative expression levels of the miRNAs and their targets [62]. All the primers used are listed in S1 and S2 Tables.

Results

Analysis of small RNA library data sets from inflorescences of Ogura-CMS line and its MF line of turnip

To identify miRNAs related to Ogura-CMS in turnip, six independent small RNA libraries from inflorescences collected from the Ogura-CMS line and its MF line were constructed. On average, 30,898,875 and 26,650,640 raw reads from the Ogura-CMS line and its MF line, respectively, were generated by high-throughput sequencing (S3 Table). All of the raw reads of the six libraries were uploaded to the Sequence Read Archive of the NCBI under accession number PRJNA552762 (http://www.ncbi.nlm.nih.gov/sra/PRJNA552762). After removing reads containing adapters, reads containing poly-N sequences, reads of low quality, averages of 27,637,247 and 24,111,786 clean reads from the Ogura-CMS line and its MF line, respectively, were obtained, with lengths ranging from 18 nt to 30 nt, accounting for approximately 90% of the total reads (S3 Table). After filtering the rRNA, tRNA, snRNA, snoRNA, other ncRNA and repeats from the clean reads, totals of 21,524,682 and 18,774,516 unannotated reads on average, respectively, were obtained (S3 Table). Among these reads, averages of 4,471,741 and 3,673,198 small RNA sequences, respectively, were mappable, accounting for approximately 15% of the clean reads from the two lines (S3 Table). In the length distribution analysis, the majority of the small RNAs in the six libraries were 21–24 nt, with the 21 and 24 nt lengths being dominant (Fig 1), which is consistent with results of other Brassica species, such as Chinese cabbage (B. rapa ssp. pekinensis) and B. juncea [45,46,48]. The highest sequence redundancy was observed in the 21 nt long fraction of Ogura-CMS line libraries and the 24 nt long fraction of MF line libraries. In addition, the results showed that the length distribution of small RNAs was similar between the Ogura-CMS line and its MF line (Fig 1).

thumbnail
Download:
Fig 1. Length distribution of small RNAs of inflorescences from the Ogura-CMS line and MF line of turnip.

https://doi.org/10.1371/journal.pone.0236829.g001

Identification of known miRNAs in inflorescences of turnip

Because the whole genome of turnip is not publicly available, the mRNA transcriptome of inflorescences of turnip uploaded to the Sequence Read Archive of the NCBI (https://www.ncbi.nlm.nih.gov/sra/PRJNA505114; accession number PRJNA505114) was used as a reference sequence. To identify miRNAs involved in reproductive development in turnip, all mappable small RNA sequences aligned to the unigenes for precursor identification were compared with previously known plant miRNAs in the miRBase database (miRBase 21.0, http://www.mirbase.org/). In the present study, 120 pre-miRNAs corresponding to 89 mature 5’- or 3’-miRNAs were detected (S4 and S5 Tables).

A total of 87 pre-miRNAs corresponding to 59 mature miRNAs that have the same sequences as known pre-miRNAs in miRBase with no mismatches allowed were identified; their read numbers among the six small RNA libraries are listed in S4 Table. The length of most known miRNAs was 21 nt, followed by 24 nt (S1A Fig). Here, there were 59 pre-miRNAs corresponding to 34 mature miRNAs originating from the Arabidopsis genome, and six mature miRNAs originating from seven pre-miRNAs came from the Brassica genome. Another 19 mature miRNAs derived from 21 pre-miRNAs identified were mapped to the genomes of other plant species. Among these 87 conserved miRNAs, ath-miR159a, ath-miR319a and ath-miR165a-3p showed very high expression levels in both the Ogura-CMS line and its MF line (S4 Table).

Among these known miRNAs, 38 miRNAs originating from 59 pre-miRNAs belonging to 30 families were detected (S6 Table). For the MIR159 family, two family members belonging to ath-miR159a were identified using deep sequencing. They shared the same mature sequences but were derived from different precursors, i.e., they originated from different loci in the turnip genome. This was also the case for several other miRNAs; a single mature miRNA sequence may originate from two to five different precursors, such as bra-miR5712, ath-miR171b-5p, ath-miR165a-3p and ath-miR394a (S4 Table). This may be attributed to the idea that some highly similar miRNA genes might be produced by a replication event from one origin sequence to another, resulting in more copies of the miRNA group. For other miRNAs, one miRNA member was generated from a single precursor (S4 Table).

Identification of novel miRNAs in inflorescences of turnip

To identify novel miRNAs in turnip, miRDeep2 software [54,55] and RNAfold were used to explore the secondary structures, Dicer cleavage sites and minimum free energies (MFEs) of the unannotated small RNA sequences that could be mapped to the unigenes assembled in the mRNA transcriptome sequences of inflorescences of turnip [52]. The small RNAs that mapped exactly to the assembled sequences but not to the previously known plant miRNAs in the miRBase database were classified as candidate novel miRNAs. As a result, a total of 30 new miRNAs derived from 33 different precursors belonging to 29 families (Tables 1 and 2) were identified as being specific to turnip in this study. Among these miRNAs, 15 new miRNA members belonging to 13 known miRNA families were discovered in six independent small RNA libraries (Table 1). In addition, a total of 18 candidate novel miRNAs belonging to 16 novel families (Table 2) were also identified. These 33 miRNAs have not been previously reported as bra-miRNAs in miRBase. The lengths of most miRNAs were 21 nt and 24 nt (S1B Fig), which was in agreement with common characteristics of plant miRNAs. Precursors of these novel miRNAs were identified via miRDeep2 software and varied from 79 to 250 nt in length, with MFE values ranging from -13.00 to -111.60 kcal/mol, and the minimum folding energy indices (MFEIs) ranging from 0.1 to 4.9 (Tables 1 and 2; S5 Table). Most of the new mature miRNA sequences presented a uracil (U) (42.4%) or an adenine (A) (39.4%) as the first nucleotide. Like the known miRNAs, some of the novel miRNAs, such as two new miRNA members (N-bra-miR4415a-3p and N-bra-miR4415b-3p) that belonged to the MIR4415 family, shared the same mature sequence, but their precursors originated from different loci of the turnip genome. These types of miRNAs were referred to as sub-members. However, unlike the known miRNAs, the expression levels of most novel miRNAs were very low in the inflorescences of turnip (S5 Table). The secondary hairpin structures of the representative miRNA precursors to N-bra-miR8007-3p and N-bra-miRn14-3p were selected as an example that is shown in Fig 2, and others are listed in S2 Fig.

thumbnail
Download:
Fig 2. Predicted secondary structure of two novel miRNAs from inflorescences of Ogura-CMS line and MF line.

(A) N-bra-miR8007-3p. (B) N-bra-miRn14-3p. The red shaded areas indicate mature miRNAs; the blue shaded areas indicate star miRNAs.

https://doi.org/10.1371/journal.pone.0236829.g002

thumbnail
Download:
Table 1. Identification of new miRNA members of known miRNA families in Ogura-CMS and MF inflorescences.

https://doi.org/10.1371/journal.pone.0236829.t001

thumbnail
Download:
Table 2. Novel miRNAs of new miRNA families identified in inflorescences of Ogura-CMS line and MF line.

https://doi.org/10.1371/journal.pone.0236829.t002

Differential expression profiling of miRNAs in inflorescences of the Ogura-CMS line and its MF line of turnip

To ensure the reliability of the data, the count of miRNAs in the six libraries was normalized to the expression of TPM [57] and then log10 transformed. The correlations of the normalized data between three biological replicates for each of the two samples were then determined. As expected, the sequencing data of miRNAs in the three biological replicates showed good repeatability (all correlation coefficients ≥ 0.942) (S3 Fig), demonstrating that the results were reliable. The miRNAs were then subjected to differential expression analysis between the inflorescences of the Ogura-CMS line and the inflorescences of the MF line of turnip. The expression abundance of all known and novel miRNAs was compared based on filter parameters of a relative fold change ≥ 2 and a FDR ≤ 0.05 after being subjected to TPM normalization. As a result, most of the miRNAs were equally expressed; a total of 10 mature miRNAs originating from 12 pre-miRNAs showed differential expression between the inflorescences of the Ogura-CMS line and the inflorescences of the MF line according to high-throughput sequencing (Fig 3 and Table 3). Among these differentially expressed mature miRNAs, 8 miRNAs belonging to 8 miRNA families were significantly down-regulated in the Ogura-CMS line compared with the MF line, whereas 2 miRNAs belonging to the MIR319 family were up-regulated, according to the normalized sequence reads (Table 3). Notably, among the 8 down-regulated miRNAs, N-bra-miR1222-3p and N-bra-miRn9-5p were expressed specifically in the MF line.

thumbnail
Download:
Fig 3. Cluster heat map of differentially expressed miRNAs in Ogura-CMS and MF inflorescences of turnip.

Differentially expressed miRNAs were filtered according to a relative fold change ≥ 2 and a Benjamini-Hochberg false discovery rate (FDR) ≤ 0.05 after being normalized to the expression of transcripts per million (TPM). The green color indicates a down-regulated pattern, and the red color indicates an up-regulated pattern.

https://doi.org/10.1371/journal.pone.0236829.g003

thumbnail
Download:
Table 3. List of differential expressed miRNAs between Ogura-CMS and MF inflorescences of turnip.

https://doi.org/10.1371/journal.pone.0236829.t003

Stem-loop qRT-PCR was conducted to validate the expression profiles of differentially expressed mature miRNAs. The results for most of the miRNAs were in agreement with those of the sequencing data (Fig 4; S4 and S5 Tables). However, the expression levels for ath-miR319a were inconsistent between the stem-loop qRT-PCR and small RNA sequencing technology. The expression of ath-miR319a was found to not significant differ between the inflorescences of the Ogura-CMS line and the inflorescences of the MF line according to stem-loop qRT-PCR, but was upregulated in the Ogura-CMS line with a |log2 FC| = 1.285 based on high-throughput sequencing analysis. This inconsistent trend possibly occurred because of the difference in sensitivity between the stem-loop qRT-PCR and high-throughput sequencing technology. Additionally, the expression patterns for each biological replicate showed high reproducibility, indicating the high reliability of the results of our deep sequencing.

thumbnail
Download:
Fig 4. Relative expression detection of miRNAs in Ogura-CMS inflorescences and MF inflorescences using stem-loop qRT-PCR.

The columns indicate the relative expression levels of miRNAs. 5S rRNA was used as an endogenous control. The results were generated from three biological replicates and three independent technical replicates, and the error bars indicate the mean ± SD (standard deviation). The lines show the normalized transcripts per million (TPM) expression data of miRNAs generated by high-throughput sequencing.

https://doi.org/10.1371/journal.pone.0236829.g004

Identification of miRNA target genes in inflorescences of the Ogura-CMS line and its MF line of turnip using degradome analysis

A degradome sequencing analysis was performed to validate the miRNA targets as miRNAs function by regulating their target genes and especially by degrading their target mRNAs in plants [48]. Six independent degradome libraries derived from inflorescences of Ogura-CMS line and its MF line were constructed. After removing the adapter sequences and low-quality reads from the raw reads, averages of 12,056,606 and 11,839,122 clean reads from the Ogura-CMS line and its MF line, respectively, were obtained (Table 4). These clean reads included 4,444,400 and 4,257,755 unique reads, respectively, of which 2,600,784 (58.4%) and 2,288,874 (53.6%), respectively, were perfectly matched to the unigenes assembled in the mRNA transcriptome sequencing of the inflorescences of turnip. After BLAST searches of the Rfam database and filtering out a small proportion of the hits that were annotated as ncRNAs, such as rRNA, tRNA, snoRNA, and snRNA, the remaining reads were further analyzed to identify miRNA targets.

thumbnail
Download:
Table 4. Analysis of degradome sequences from inflorescences of Ogura-CMS line and MF line in turnip.

https://doi.org/10.1371/journal.pone.0236829.t004

A total of 22 targets for 25 miRNAs (23 known miRNAs and 2 novel miRNAs) and 17 targets for 28 miRNAs (26 known miRNAs and 2 novel miRNAs) were identified as being involved in the reproductive development of the Ogura-CMS and MF lines of turnip, respectively (Table 5). Based on the distribution of cleavage sites near raw sequence tags, the cleavage sites of identified miRNA targets could be divided into two categories: category ‘0’ and category ‘1’ (Table 5), among which category ‘0’ is described as > 1 raw read at the position with abundance at the position equal to the maximum on the transcript and only one maximum on the transcript, and category ‘1’ is defined as > 1 raw read at the position with abundance at the position equal to the maximum on the transcript and more than one maximum position on the transcript [45,50]. Six representative miRNAs and their corresponding targets were selected, and t-plots were constructed, in which the red line indicates the cleavage site of each transcript (Fig 5). Among these targets, most were identified for known miRNAs, and only two were for two novel miRNAs (N-bra-miR4415a-3p and N-bra-miR4415b-3p) that accumulated to a relatively high expression level (Table 5 and S5 Table). However, the target genes of most identified miRNAs could not be detected in the present degradome analysis. Unfortunately, for the vast majority of the novel miRNAs, their targets could not be identified.

thumbnail
Download:
Fig 5. T-plots of representative identified miRNA targets in inflorescences of Ogura-CMS line and MF line.

The X axis indicates the site position of miRNA target sequences, and the Y axis indicates the abundance of raw tags. The red peak indicates the predicted cleavage site and p ≤ 0.05.

https://doi.org/10.1371/journal.pone.0236829.g005

thumbnail
Download:
Table 5. Target genes identified by degradome sequencing in inflorescences of Ogura-CMS line and MF line.

https://doi.org/10.1371/journal.pone.0236829.t005

According to the annotation analysis of miRNA targets, the miRNAs targeted different genes involved in a wide variety of biological functions, including several transcription factors that play important roles in gene regulation, such as homeobox-leucine zipper proteins, squamosal-promoter binding proteins (SBPs), auxin response factors (ARFs), APETALA 2 (AP2) and growth-regulating factors (GRFs) (Table 5). Some of these transcription factor targets, such as AP2, TOE3, CUP-SHAPED COTYLEDON 1 (CUC1), may be related to CMS or anther development. Some non-transcription factor targets were further identified as miRNA targets in our study, such as PPR-containing proteins and disease resistance proteins. Additionally, some target genes have unknown functions (Table 5).

Expression of miRNA targets in inflorescences of the Ogura-CMS line and its MF line of turnip

To test the causal relationship between the level of miRNAs and their target expression, four selective miRNAs and their corresponding cleavage targets identified in the degradome sequencing analysis were satisfactorily verified by stem-loop qRT-PCR and regular qRT-PCR, respectively. As expected, the expression levels of the targets of ath-miR156a-5p (TRINITY_DN18236_c0_g1_2308), two sub-members of ath-miR157a-5p (TRINITY_DN23691_c0_g1_12492 and TRINITY_DN40654_c0_g1_45037) and ath-miR5654-5p (TRINITY_DN24251_c3_g2_15972) were higher in the CMS line than in the MF line, which contrasted with their corresponding miRNA data (Fig 6). ath-miR156a-5p targets SBPs (TRINITY_DN23417_c1_g4, TRINITY_DN23417_c1_g13 and TRINITY_DN23946_c2_g2), two sub-members of ath-miR157a-5p also target SBPs (TRINITY_DN23417_c1_g4 and TRINITY_DN23417_c1_g13), and ath-miR5654-5p targets TRINITY_DN25443_c0_g1 (PPR-containing protein) (Table 5). Taken together, the results showed that there are negative correlations between the expression patterns of miRNAs and their targets.

thumbnail
Download:
Fig 6. Expression of miRNAs and their corresponding targets in inflorescences of Ogura-CMS line and MF line.

The normalized transcripts per million (TPM) expression data of miRNAs were generated from small RNA sequencing. The fragments per kilobase of transcript per million mapped reads (FPKMs) expression data of miRNA targets were generated from RNA-sequencing [52]. 5S rRNA was used as an internal standard for miRNA analysis, and ACT7 was used as an endogenous control for miRNA target analysis. The results were generated from three biological replicates and three independent technical replicates, and the error bars indicate the mean ± SD (standard deviation).

https://doi.org/10.1371/journal.pone.0236829.g006

Discussion

The predominant role of miRNAs in plant biological processes and responses to environmental stresses has directed attention in recent years to plant miRNA research [31,32]. With the development of in silico tools and innovative techniques, such as small RNA deep sequencing coupled with degradome analysis and miRNA arrays, a substantial number of miRNAs have been identified in Arabidopsis, rice and many other plant species [63,64]. However, the focus of the research community is still concentrated on model and key plant species [34]. Although the vast majority of miRNA families are conserved in angiosperms, miRNA members and their expression levels vary by plant species [27]. MiRNA identification should be expanded to non-model and lesser-explored plant species. Turnip is one of the most important local root vegetables in China and has been consumed for thousands of years [65]. In the present study, using high-throughput sequencing, a total of 120 pre-miRNAs corresponding to 89 mature miRNAs were identified in inflorescences of the Ogura-CMS line and its MF line of turnip (S4 and S5 Tables). Among these identified miRNAs, 33 miRNAs were novel (Tables 1 and 2). Based on the miRNA family analysis, 38 known miRNAs originating from 59 pre-miRNAs were subjected to 30 families (S6 Table). Some of the identified miRNAs are highly conserved among diverse angiosperms, such as miR156, miR160, miR164, miR172 and miR396, while some are species specific, such as miR5718, which is specific to B. campestris [27,66].

The participation of miRNAs in pollen development is becoming increasingly evident [66]. In the case of CMS, the regulatory network responsible for pollen generation refers to the rearrangements of mitochondrial CMS-inducing genes and their interactions with nuclear genes [10]. Although the exact nature of miRNAs with respect to CMS in plants remains elusive, discovery of miRNAs with possible implications for CMS induction is intriguing in this respect. The recent explosion of high-throughput sequencing has enabled genome-wide discovery of these regulatory molecules and their targets, implying their potential involvement in CMS [67]. Using this method combined with bioinformatic analysis, researchers have successively identified differentially expressed miRNAs and their targets between CMS lines and fertile lines of rice, B. juncea, cybrid pummelo (Citrus grandis) and many other species [45–51]. In the present study, 10 differentially expressed mature miRNAs originating from 12 pre-miRNAs in turnip were identified with more than two-fold changes in expression between inflorescences of the Ogura-CMS line and inflorescences of the MF line (Fig 3 and Table 3). Although direct evidence substantiating miRNAs as causative agents for CMS phenotypes remains unavailable, the differences in miRNA biogenesis in the Ogura-CMS line in our study may partly reflect mitochondrial retrograde regulation of the nuclear transcriptome.

Using high-throughput sequencing coupled with degradome analysis, we identified a total of 22 targets for 25 miRNAs and 17 targets for 28 miRNAs involved in reproductive development for the Ogura-CMS line and its MF line of turnip, respectively (Table 5). Plant miRNAs have been implicated in the negative regulation of the stability of their target genes by RNA cleavage and translational repression, of which directing target cleavage is the best-known and predominant mode of miRNA action [27]. The most common targets for CMS-related miRNAs are genes encoding transcription factors. MiRNA-mediated regulation of transcription factors affects various downstream related genes, suggesting a core role of the involvement of miRNAs in gene regulatory networks [68]. In the present study, several transcription factor targets were identified, such as homeobox-leucine zipper proteins, SBPs, ARFs, AP2s, GRFs and MYB domain proteins (Table 5). Some of these transcription factor targets resemble those of previous studies with respect to their defined or inferred roles in regulating reproductive development processes. For instance, miR156 and miR157 are well known for controlling anther development by targeting SBP-like genes [69,70]. miR396 targets GRF genes, including GRF9 [68]. In addition, miR172-mediated regulation of AP2 and AP2-like genes such as TOE3 is essential for floral organ identity in Arabidopsis [71–73]. Perturbing miR164-directed repression of CUC1 disrupts floral development [74]. Although conclusive evidence ascertaining the linkage between ARF10 and reproductive development remains obscure, overexpression of two other ARF genes (ARF16 and ARF17) free from miR160 regulation was reported to result in aberrant flowers with floral organ defects and reduced fertility [75–77].

Apart from transcription factors, some functional protein targets were also identified for the miRNAs in our study, indicating a role for miRNAs in regulating protein stability. PPR-containing proteins constitute an important class of fertility restorer proteins that are required for generating a functional male gametophyte in CMS plants [14,19]. PPR-containing proteins have been previously confirmed as targets of miR158, miR161, miR475 and miR476 [69,78,79]. Here, miR5654 and miR1507 were also predicted to target PPR-containing proteins (Table 5). In addition, other protein targets were also identified, such as disease resistance proteins and phosphoribosyltransferase (Table 5). Together, all of these transcription factor and functional protein target prediction analysis results indicate that some of the miRNAs are involved in anther development and male sterility, leading to a possible connection between miRNA action and Ogura-CMS phenotype.

In the present study, the discovery of miRNAs and their targets helps to improve the understanding of the contribution of miRNAs in CMS and the influence of mitochondrial CMS-inducing genes on global transcriptome changes, although, the effects of miRNAs on CMS occurrence and the exact roles of miRNA-mediated regulation of transcription factors or other downstream related genes in CMS induction is still demanding complementing genetic technologies to substantiate the findings of our present study.

Conclusion

In this study, high-throughput sequencing and degradome analysis were employed to enable genome-wide discovery of conserved and novel miRNAs and their targets possibly involved in anther development in turnip. Comprehensive analysis of miRNAs and their targets in inflorescences of the Ogura-CMS line and its MF line of turnip suggested that mutation of mitochondrial orf138 leads to the fine-tuned expression of miRNAs, which may further participate in the regulatory network of CMS occurrence by regulating their target genes.

Supporting information

S1 Fig. Length distribution of known and novel miRNAs in inflorescences from Ogura-CMS line and MF line.

(A) Known miRNAs. (B) Novel miRNAs.

https://doi.org/10.1371/journal.pone.0236829.s001

(TIF)

S2 Fig. Predicted secondary structures of novel miRNAs from inflorescences of Ogura-CMS line and MF line.

The red shaded areas indicate mature miRNAs; the blue shaded areas indicate star miRNAs.

https://doi.org/10.1371/journal.pone.0236829.s002

(TIF)

S3 Fig. Correlation coefficiencies between three biological replicates for the sequencing data of miRNAs.

https://doi.org/10.1371/journal.pone.0236829.s003

(TIF)

S1 Table. Specific primers used for stem-loop quantitative real-time PCR (qRT-PCR) validation of miRNAs in turnip.

https://doi.org/10.1371/journal.pone.0236829.s004

(XLSX)

S2 Table. Primers designed for quantitative real-time PCR (qRT-PCR) validation of miRNA targets in turnip.

https://doi.org/10.1371/journal.pone.0236829.s005

(XLSX)

S3 Table. Analysis of small RNA sequences from inflorescences of Ogura-CMS line and MF line of turnip.

https://doi.org/10.1371/journal.pone.0236829.s006

(XLSX)

S4 Table. Known miRNAs identified in inflorescences of Ogura-CMS line and MF line of turnip.

https://doi.org/10.1371/journal.pone.0236829.s007

(XLSX)

S5 Table. Novel miRNAs identified in inflorescences of Ogura-CMS line and MF line of turnip.

https://doi.org/10.1371/journal.pone.0236829.s008

(XLSX)

S6 Table. Identification of known miRNA members of known miRNA families in Ogura-CMS and MF inflorescences.

https://doi.org/10.1371/journal.pone.0236829.s009

(XLSX)

Acknowledgments

We thank Dr. Heng Dong (Hangzhou Normal University, China) for the critical comments and editing on this manuscript.

References

  1. 1. Yamagishi H, Bhat SR. Cytoplasmic male sterility in Brassicaceae crops. Breed Sci. 2014;64: 38–47. pmid:24987289
  2. 2. Hanson MR, Bentolila S. Interactions of mitochondrial and nuclear genes that affect male gametophyte development. Plant Cell. 2004;16: S154–S169. pmid:15131248
  3. 3. Ogura H. Studies on the new male-sterility in Japanese radish, with special reference to the utilization of this sterility. Mem Fac Agric Kagoshima Univ. 1968;6: 39–78.
  4. 4. Fu T, Yang G, Yang X. Studies on “three line” Polima cytoplasmic male sterility developed in Brassica napus L. Plant Breed. 1990;104: 115–120.
  5. 5. Iwabuchi M, Koizuka N, Fujimoto H, Sakai T, Imamura J. Identification and expression of the kosena radish (Raphanus sativus cv. Kosena) homologue of the ogura radish CMS-associated gene, orf138. Plant Mol Biol. 1999;39: 183–188. pmid:10080721
  6. 6. Yamagishi H, Terachi T. Intra- and inter-specific variations in the mitochondrial gene orf138 of Ogura-type male-sterile cytoplasm from Raphanus sativus and Raphanus raphanistrum. Theor Appl Genet. 2001;103: 725–732.
  7. 7. Nahm SH, Lee HJ, Lee SW, Joo GY, Harn CH, Yang SG, et al. Development of a molecular marker specific to a novel CMS line in radish (Raphanus sativus L.). Theor Appl Genet. 2005;111: 1191–1200. pmid:16142466
  8. 8. Lee YP, Park S, Lim C, Kim H, Lim H, Ahn Y, et al. Discovery of a novel cytoplasmic male-sterility and its restorer lines in radish (Raphanus sativus L.). Theor Appl Genet. 2008;117: 905–913. pmid:18597066
  9. 9. Dahan J, Mireau H. The Rf and Rf-like PPR in higher plants, a fast-evolving subclass of PPR genes. RNA Biol. 2013;10: 1469–1476. pmid:23872480
  10. 10. Horn R, Gupta KJ, Colombo N. Mitochondrion role in molecular basis of cytoplasmic male sterility. Mitochondrion. 2014;19: 198–205. pmid:24732436
  11. 11. Hu J, Huang W, Huang Q, Qin X, Yu C, Wang L, et al. Mitochondria and cytoplasmic male sterility in plants. Mitochondrion. 2014;19: 282–288. pmid:24566371
  12. 12. Gaborieau L, Brown GG, Mireau H. The propensity of pentatricopeptide repeat genes to evolve into restorers of cytoplasmic male sterility. Front Plant Sci. 2016;7: 1816. pmid:27999582
  13. 13. Chen Z, Zhao N, Li S, Grover CE, Nie H, Wendel JF, et al. Plant mitochondrial genome evolution and cytoplasmic male sterility. Crit Rev Plant Sci. 2017;36: 55–69.
  14. 14. Liu Z, Dong F, Wang X, Wang T, Su R, Hong D, et al. A pentatricopeptide repeat protein restores nap cytoplasmic male sterility in Brassica napus. J Exp Bot. 2017;68: 4115–4123. pmid:28922764
  15. 15. Rizzolatti C, Bury P, Tatara E, Pin PA, Rodde N, Bergès H, et al. Map-based cloning of the fertility restoration locus Rfm1 in cultivated barley (Hordeum vulgare). Euphytica. 2017;213: 276.
  16. 16. Singh S, Dey SS, Bhatia R, Kumar R, Behera TK. Current understanding of male sterility systems in vegetable Brassicas and their exploitation in hybrid breeding. Plant Reprod. 2019;32: 231–256. pmid:31053901
  17. 17. Kubo T, Arakawa T, Honma Y, Kitazaki K. What does the molecular genetics of different types of restorer-of-fertility genes imply? Plants. 2020;9: 361.
  18. 18. Lee B, Lee H, Xiong L, Zhu JK. A mitochondrial complex I defect impairs cold-regulated nuclear gene expression. Plant Cell. 2002;14: 1235–1251. pmid:12084824
  19. 19. Schmitz-Linneweber C, Small I. Pentatricopeptide repeat proteins: a socket set for organelle gene expression. Trends Plant Sci. 2008;13: 663–670. pmid:19004664
  20. 20. Dong X, Kim WK, Lim YP, Kim YK, Hur Y. Ogura-CMS in Chinese cabbage (Brassica rapa ssp. pekinensis) causes delayed expression of many nuclear genes. Plant Sci. 2013;199: 7–17. pmid:23265314
  21. 21. Laser KD, Lersten NR. Anatomy and cytology of microsporogenesis in cytoplasmic male sterile angiosperms. Bot Rev. 1972;38: 425–454.
  22. 22. Gottschalk W, Kaul MLH. The genetic control of microsporogenesis in higher plants. Nucleus. 1974;17: 133–166.
  23. 23. Colhoun CW, Steer MW. Microsporogenesis and the mechanism of cytoplasmic male sterility in maize. Ann Bot. 1981;48: 417–424.
  24. 24. Chen L, Liu L, Jin P, Gong Y, Sun X, Ma E. Cytological and molecular identification of cytoplasm in two male sterile lines in radish. Mol Plant Breed. 2009;7: 757–762.
  25. 25. Axtell MJ. Classification and comparison of small RNAs from plants. Annu Rev Plant Biol. 2013;64: 137–159. pmid:23330790
  26. 26. Bartel DP. MicroRNAs: genomics, biogenesis, mechanism, and function. Cell. 2004;116: 281–297. pmid:14744438
  27. 27. Jones-Rhoades MW, Bartel DP, Bartel B. MicroRNAs and their regulatory roles in plants. Annu Rev Plant Biol. 2006;57: 19–53. pmid:16669754
  28. 28. Teotia S, Tang G. To bloom or not to bloom: role of microRNAs in plant flowering. Mol Plant. 2015;8: 359–377. pmid:25737467
  29. 29. Li C, Zhang B. MicroRNAs in control of plant development. J Cell Physiol. 2016;231: 303–313. pmid:26248304
  30. 30. D’Ario M, Griffiths-Jones S, Kim M. Small RNAs: big impact on plant development. Trends Plant Sci. 2017;22: 1056–1068. pmid:29032035
  31. 31. Singh A, Gautam V, Singh S, Sarkar Das S, Verma S, Mishra V, et al. Plant small RNAs: advancement in the understanding of biogenesis and role in plant development. Planta. 2018;248: 545–558. pmid:29968061
  32. 32. Xu J, Hou QM, Khare T, Verma SK, Kumar V. Exploring miRNAs for developing climate-resilient crops: A perspective review. Sci Total Environ. 2019;653: 91–104. pmid:30408672
  33. 33. Ding Y, Ding L, Xia Y, Wang F, Zhu C. Emerging roles of microRNAs in plant heavy metal tolerance and homeostasis. J Agric Food Chem. 2020;68: 1958–1965. pmid:32003983
  34. 34. Zhou X, Khare T, Kumar V. Recent trends and advances in identification and functional characterization of plant miRNAs. Acta Physiol Plant. 2020;42: 25.
  35. 35. Chambers C, Shuai B. Profiling microRNA expression in Arabidopsis pollen using microRNA array and real-time PCR. BMC Plant Biol. 2009;9: 87. pmid:19591667
  36. 36. Grant-Downton R, Hafidh S, Twell D, Dickinson HG. Small RNA pathways are present and functional in the angiosperm male gametophyte. Mol Plant. 2009;2: 500–512. pmid:19825633
  37. 37. Grant-Downton R, Le Trionnaire G, Schmid R, Rodriguez-Enriquez J, Hafidh S, Mehdi S, et al. MicroRNA and tasiRNA diversity in mature pollen of Arabidopsis thaliana. BMC Genomics. 2009;10: 643. pmid:20042113
  38. 38. Wei LQ, Yan LF, Wang T. Deep sequencing on genome-wide scale reveals the unique composition and expression patterns of microRNAs in developing pollen of Oryza sativa. Genome Biol. 2011;12: R53. pmid:21679406
  39. 39. Peng H, Chun J, Ai T, Tong Y, Zhang R, Zhao M, et al. MicroRNA profiles and their control of male gametophyte development in rice. Plant Mol Biol. 2012;80: 85–102. pmid:22403030
  40. 40. He H, Yang T, Wu W, Zheng B. Small RNAs in pollen. Sci China Life Sci. 2015;58: 246–252. pmid:25634522
  41. 41. Achard P, Herr A, Baulcombe DC, Harberd NP. Modulation of floral development by a gibberellin-regulated microRNA. Development. 2004;131: 3357–3365. pmid:15226253
  42. 42. Millar AA, Gubler F. The Arabidopsis GAMYB-like genes, MYB33 and MYB65, are microRNA-regulated genes that redundantly facilitate anther development. Plant Cell. 2005;17: 705–721. pmid:15722475
  43. 43. Ru P, Xu L, Ma H, Huang H. Plant fertility defects induced by the enhanced expression of microRNA167. Cell Res. 2006;16: 457–465. pmid:16699541
  44. 44. Wu MF, Tian Q, Reed JW. Arabidopsis microRNA167 controls patterns of ARF6 and ARF8 expression, and regulates both female and male reproduction. Development. 2006;133: 4211–4218. pmid:17021043
  45. 45. Yang J, Liu X, Xu B, Zhao N, Yang X, Zhang M. Identification of miRNAs and their targets using high-throughput sequencing and degradome analysis in cytoplasmic male-sterile and its maintainer fertile lines of Brassica juncea. BMC Genomics. 2013;14: 9. pmid:23324572
  46. 46. Jiang J, Lv M, Liang Y, Ma Z, Cao J. Identification of novel and conserved miRNAs involved in pollen development in Brassica campestris ssp. chinensis by high-throughput sequencing and degradome analysis. BMC Genomics. 2014;15: 146. pmid:24559317
  47. 47. Omidvar V, Mohorianu I, Dalmay T, Fellner M. Identification of miRNAs with potential roles in regulation of anther development and male-sterility in 7B-1 male-sterile tomato mutant. BMC Genomics. 2015;16: 878. pmid:26511108
  48. 48. Wei X, Zhang X, Yao Q, Yuan Y, Li X, Wei F, et al. The miRNAs and their regulatory networks responsible for pollen abortion in Ogura-CMS Chinese cabbage revealed by high-throughput sequencing of miRNAs, degradomes, and transcriptomes. Front Plant Sci. 2015;6: 894. pmid:26557132
  49. 49. Yan J, Zhang H, Zheng Y, Ding Y. Comparative expression profiling of miRNAs between the cytoplasmic male sterile line MeixiangA and its maintainer line MeixiangB during rice anther development. Planta. 2015;241: 109–123. pmid:25228384
  50. 50. Ding X, Li J, Zhang H, He T, Han S, Li Y, et al. Identification of miRNAs and their targets by high-throughput sequencing and degradome analysis in cytoplasmic male-sterile line NJCMS1A and its maintainer NJCMS1B of soybean. BMC Genomics. 2016;17: 24. pmid:26729289
  51. 51. Fang YN, Zheng BB, Wang L, Yang W, Wu XM, Xu Q, et al. High-throughput sequencing and degradome analysis reveal altered expression of miRNAs and their targets in a male-sterile cybrid pummelo (Citrus grandis). BMC Genomics. 2016;17: 591. pmid:27506907
  52. 52. Lin S, Miao Y, Su S, Xu J, Jin L, Sun D, et al. Comprehensive analysis of Ogura cytoplasmic male sterility-related genes in turnip (Brassica rapa ssp. rapifera) using RNA sequencing analysis and bioinformatics. PLoS ONE. 2019;14: e0218029. pmid:31199816
  53. 53. Langmead B, Trapnell C, Pop M, Salzberg SL. Ultrafast and memory-efficient alignment of short DNA sequences to the human genome. Genome Biol. 2009;10: R25. pmid:19261174
  54. 54. Friedländer MR, Mackowiak SD, Li N, Chen W, Rajewsky N. miRDeep2 accurately identifies known and hundreds of novel microRNA genes in seven animal clades. Nucleic Acids Res. 2012;40: 37–52. pmid:21911355
  55. 55. Zhang Z, Jiang L, Wang J, Gu P, Chen M. MTide: an integrated tool for the identification of miRNA-target interaction in plants. Bioinformatics. 2015;31: 290–291. pmid:25256573
  56. 56. Bonnet E, Wuyts J, Rouzé P, Van de Peer Y. Evidence that microRNA precursors, unlike other non-coding RNAs, have lower folding free energies than random sequences. Bioinformatics. 2004;20: 2911–2917. pmid:15217813
  57. 57. Kasschau KD, Fahlgren N, Chapman EJ, Sullivan CM, Cumbie JS, Givan SA, et al. Genome-wide profiling and analysis of Arabidopsis siRNAs. PLoS Biol. 2007;5: e57. pmid:17298187
  58. 58. Willmann MR, Berkowitz ND, Gregory BD. Improved genome-wide mapping of uncapped and cleaved transcripts in eukaryotes—GMUCT 2.0. Methods. 2014;67: 64–73. pmid:23867340
  59. 59. Addo-Quaye C, Miller W, Axtell MJ. CleaveLand: a pipeline for using degradome data to find cleaved small RNA targets. Bioinformatics. 2009;25: 130–131. pmid:19017659
  60. 60. Varkonyi-Gasic E, Hellens RP. Quantitative stem-loop RT-PCR for detection of microRNAs. In: Kodama H, Komamine A, editors. RNAi and Plant Gene Function Analysis. Totowa, NJ: Humana Press; 2011. pp. 145–157.
    • 61. Wang C, Cui HM, Huang TH, Liu TK, Hou XL, Li Y. Identification and validation of reference genes for RT-qPCR analysis in non-heading Chinese cabbage flowers. Front Plant Sci. 2016;7: 811. pmid:27375663
    • 62. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2(-Delta Delta C(T)) method. Methods. 2001;25: 402–408. pmid:11846609
    • 63. Banerjee A, Roychoudhury A, Krishnamoorthi S. Emerging techniques to decipher microRNAs (miRNAs) and their regulatory role in conferring abiotic stress tolerance of plants. Plant Biotechnol Rep. 2016;10: 185–205.
    • 64. Slatko BE, Gardner AF, Ausubel FM. Overview of next‐generation sequencing technologies. Curr Protoc Mol Biol. 2018;122: e59. pmid:29851291
    • 65. Zhang H, Schonhof I, Krumbein A, Gutezeit B, Li L, Stützel H, et al. Water supply and growing season influence glucosinolate concentration and composition in turnip root (Brassica rapa ssp. rapifera L.). J Plant Nutr Soil Sci. 2010;171: 255–265.
    • 66. Mishra A, Bohra A. Non-coding RNAs and plant male sterility: current knowledge and future prospects. Plant Cell Rep. 2018;37: 177–191. pmid:29332167
    • 67. Štorchová H. The role of non-coding RNAs in cytoplasmic male sterility in flowering plants. Int J Mol Sci. 2017;18: 2429.
    • 68. Jones-Rhoades MW, Bartel DP. Computational identification of plant microRNAs and their targets, including a stress-induced miRNA. Mol Cell. 2004;14: 787–799. pmid:15200956
    • 69. Rhoades MW, Reinhart BJ, Lim LP, Burge CB, Bartel B, Bartel DP. Prediction of plant microRNA targets. Cell. 2002;110: 513–520. pmid:12202040
    • 70. Xing S, Salinas M, Höhmann S, Berndtgen R, Huijser P. miR156-targeted and nontargeted SBP-box transcription factors act in concert to secure male fertility in Arabidopsis. Plant Cell. 2010;22: 3935–3950. pmid:21177480
    • 71. Aukerman MJ, Sakai H. Regulation of flowering time and floral organ identity by a microRNA and its APETALA2-like target genes. Plant Cell. 2003;15: 2730–2741. pmid:14555699
    • 72. Chen X. A microRNA as a translational repressor of APETALA2 in Arabidopsis flower development. Science. 2004;303: 2022–2025. pmid:12893888
    • 73. Schwab R, Palatnik JF, Riester M, Schommer C, Schmid M, Weigel D. Specific effects of microRNAs on the plant transcriptome. Dev Cell. 2005;8: 517–527. pmid:15809034
    • 74. Mallory AC, Dugas DV, Bartel DP, Bartel B. MicroRNA regulation of NAC-domain targets is required for proper formation and separation of adjacent embryonic, vegetative, and floral organs. Curr Biol. 2004;14: 1035–1046. pmid:15202996
    • 75. Kasschau KD, Xie Z, Allen E, Llave C, Chapman EJ, Krizan KA, et al. P1/HC-Pro, a viral suppressor of RNA silencing, interferes with Arabidopsis development and miRNA function. Dev Cell. 2003;4: 205–217. pmid:12586064
    • 76. Mallory AC, Bartel DP, Bartel B. MicroRNA-directed regulation of Arabidopsis AUXIN RESPONSE FACTOR17 is essential for proper development and modulates expression of early auxin response genes. Plant Cell. 2005;17: 1360–1375. pmid:15829600
    • 77. Wang JW, Wang LJ, Mao YB, Cai WJ, Xue HW, Chen XY. Control of root cap formation by microRNA-targeted auxin response factors in Arabidopsis. Plant Cell. 2005;17: 2204–2216. pmid:16006581
    • 78. Lu S, Sun YH, Shi R, Clark C, Li L, Chiang VL. Novel and mechanical stress–responsive microRNAs in Populus trichocarpa that are absent from Arabidopsis. Plant Cell. 2005;17: 2186–2203. pmid:15994906
    • 79. Ma Z, Jiang J, Hu Z, Lyu T, Yang Y, Jiang J, et al. Over-expression of miR158 causes pollen abortion in Brassica campestris ssp. chinensis. Plant Mol Biol. 2017;93: 313–326. pmid:27909970