Skip to main content
  • Research article
  • Open access
  • Published:

Identification of cold stress responsive microRNAs in two winter turnip rape (Brassica rapa L.) by high throughput sequencing

Abstract

Background

Low temperature is a major abiotic stress affecting the production of rapeseed in China by impeding plant growth and development. A comprehensive knowledge of small-RNA expression pattern in Brassica rapa under cold stress could improve our knowledge of microRNA-mediated stress responses.

Results

A total of 353 cold-responsive miRNAs, 84 putative novel and 269 conserved miRNAs, were identified from the leaves and roots of two winter turnip rape varieties ‘Longyou 7’ (cold-tolerant) and ‘Tianyou 4’ (cold-sensitive), which were stressed under − 4 °C for 8 h. Eight conserved (miR166h-3p-1, miR398b-3p, miR398b-3p-1, miR408d, miR156a-5p, miR396h, miR845a-1, miR166u) and two novel miRNAs (Bra-novel-miR3153-5p and Bra-novel-miR3172-5p) were differentially expressed in leaves of ‘Longyou 7’ under cold stress. Bra-novel-miR3936-5p was up-regulated in roots of ‘Longyou 7’ under cold stress. Four and five conserved miRNAs were differentially expressed in leaves and roots of ‘Tianyou 4’ after cold stress. Besides, we found two conserved miRNAs (miR319e and miR166m-2) were down-regulated in non-stressed roots of ‘Longyou 7’ compared with ‘Tianyou 4’. After cold stress, we found two and eight miRNAs were differentially expressed in leaves and roots of ‘Longyou 7’ compared with ‘Tianyou 4’. The differentially expressed miRNAs between two cultivars under cold stress include novel miRNAs and the members of the miR166 and miR319 families. A total of 211 target genes for 15 known miRNAs and two novel miRNAs were predicted by bioinformatic analysis, mainly involved in metabolic processes and stress responses. Five differentially expressed miRNAs and predicted target genes were confirmed by quantitative reverse transcription PCR, and the expressional changes of target genes were negatively correlated to differentially expressed miRNAs. Our data indicated that some candidate miRNAs (e.g., miR166e, miR319, and Bra-novel-miR3936-5p) may play important roles in plant response to cold stress.

Conclusions

Our work indicates that miRNA and putative target genes mediated metabolic processes and stress responses are significant to cold tolerance in B. rapa.

Background

Cold tolerance is a survival strategy of some overwintering plants in frigid winter. Cold-tolerant plants respond and adapt to cold stress through a number of biochemical changes, such as accumulation of metabolites [1, 2], synthesis of unsaturated fatty acids [3], production of antifreeze proteins [4], and increase of heat shock proteins [5]. In response to cold stress, plants have established stress-responsive molecular regulation mechanisms at the transcriptional, translational, and post-translational levels [6]. Plant cells perceive cold stress signals by activating specific transcription factors (TFs) and a series of complex signal transductions. Subsequently, the TFs combine with cis-acting elements to promote and regulate the transcription and expression of target genes, thereby improving the cold-tolerance capability of plants [7, 8].

MicroRNAs have been recently reported as new regulator in plant adaptation to environmental stresses [9, 10], which are small (20–24 nt), endogenous, and noncoding RNAs originate from 70 to 80 nt pre-miRNAs with stem-loop structures. Pre-miRNAs are sequentially processed by Dicer-like proteins in cell nucleus [11] and then transferred to cytoplasm before assembled to RNA-induced silencing complexes (RISCs), which combine with Argonaute proteins. RISCs negatively regulate target gene expression by undergoing base pairing with the 3′ UTR and coding regions of mRNA [12, 13], inducing the degradation of target mRNA or inhibition of translation [14,15,16]. Therefore, miRNAs play a vital role in regulating gene expression [17, 18]. Numerous studies concentrated on their role in plant response to various stresses have been reported [7, 8]. To date, stress-responsive miRNAs in a diverse number of plants under adverse conditions were reported [19, 20]. The first plant defense miRNA discovered was miR393 in Arabidopsis [21, 22], which could increase the amount of cold-induced proteins by depressing its target gene (ubiquitin E3 ligase) [23, 24]. Furthermore, many miRNAs related to cold tolerance, such as miR156, miR166, miR172,miR319, miR396, and miR397, are identified in plants with capability of regulating cold stress responses [24,25,26,27]. Presently, miRNA related to cold stress is in need to be identified in Brassica, which is of great economic value in edible vegetable and oil production.

Winter turnip rape (B. rapa L.) is a valuable oil crop for its high-quality edible oil in northwestern China, which is also capable of conserving soil and water in winter and spring [28]. However, extremely low temperature negatively affects the growth and development of winter turnip rape, which consequently fails to overwinter and propagate in northwestern China. The aboveground parts of winter turnip rape are withered during winter, whereas the roots are capable of overwintering, indicating the root tolerance to cold is important for its survival in winter [29]. Unlike Arabidopsis and other plants, cold-tolerance in winter turnip rape is rarely studied. ‘Longyou 7’ is the first cultivated cold-tolerant winter turnip rape variety that can survive at extremely low temperatures (− 32 °C, overwinter survival rate is more than 90%). And ‘Tianyou 4’ is a cold-sensitive winter turnip rape variety. Recently, our group has investigated the physiological, biochemical, and molecular mechanisms of cold-tolerance in winter turnip rape [28,29,30,31]. Whereas the miRNAs and target genes with putative functions in cold tolerance of B. rapa has not been reported yet. In the present study, we constructed and sequenced small RNA (sRNA) libraries from leaves and roots of two winter turnip rape varieties (‘Longyou 7’ and ‘Tianyou 4’) treated with cold stress (− 4 °C), to acknowledge the expressional differences of miRNAs related to cold stress in B. rapa. The miRNAs and target genes identified in this study could facilitate the understanding of molecular mechanisms of cold stress signaling in B. rapa.

Methods

Plant materials, cold stress treatment and physiological analysis

Winter turnip rape ‘Longyou 7’ (cold-tolerant cultivar) and ‘Tianyou 4’ (cold-sensitive cultivar) were provided by the College of Agriculture, Gansu Agricultural University (Lanzhou, China). Plants were grown in plastic pots until six true leaves in a ventilated greenhouse (20 °C, 16 h light/8 h dark cycle), and then moved into growth chamber (Safu, Ningbo, China) for treatment. Cold treatment was started by gradually decreasing the temperature from 20 °C to − 4 °C. The plants were grown at 20 °C for 48 h and then induced at 10 °C and 4 °C for 48 h, respectively. Finally, they were stressed at − 4 °C for 8 h. The leaves and roots of plants grown at 20 °C (as control, CK) and − 4 °C (as treatments, TR) were collected and frozen in liquid nitrogen for RNA extraction and physiological analysis. Generally, eight samples were collected and named 4LCK (leaf of ‘Tianyou 4’ at 20 °C), 4LTR (leaf of ‘Tianyou 4’ treated at − 4 °C), 4RCK (root of ‘Tianyou 4’ at 20 °C), 4RTR (root of ‘Tianyou 4’ treated at − 4 °C), 7LCK (leaf of ‘Longyou 7’ at 20 °C), 7LTR (leaf of ‘Longyou 7’ treated at − 4 °C), 7RCK (root of ‘Longyou 7’ at 20 °C), and 7RTR (root of ‘Longyou 7’ treated at − 4 °C), respectively. Each sample were polled from three plants and two biological replicates were included. Malondialdehyde (MDA) content and peroxidase (POD) activities was analyzed according to Campos et al. [32] and Quiroga et al. [33].

RNA isolation, library construction and sequencing

Total RNA was extracted using TRNzol Universal Reagent (Tiangen, China) following the manufacturer’s instructions. The quality and concentration of RNA were tested with Agilent 2100 bio-analyzer. The sRNA libraries were constructed as previously described [34, 35]. In brief, the 18–30 nt sRNAs of each sample were gel purified from the total RNA. Two adaptors (5′ and 3′ adaptors) were ligated to the purified sRNAs successively and then reversely transcribed and amplified for 15 cycles by PCR. Then sRNA libraries (including two biological replicates) of the eight samples were constructed. The quality and production of the libraries were further tested by Agilent 2100 Bioanalyzer and ABI StepOnePlus Real-Time PCR system, and then sequenced on an Illumina Genome Analyzer by the Beijing Genomics Institute in Shenzhen (http://www.genomics.cn/en/index).

Analysis of small RNA sequencing data

Raw data from high-throughput sequencing were filtered as previously described to obtain clean reads [35, 36]. The process included the removal of adapters and low-quality sequences containing poly A and less than 18 nucleotides. Next, the clean sequences were mapped onto the B. rapa genome and Rfam12.2 (http://rfam.xfam.org/) to remove noncoding RNAs, such as rRNAs, tRNAs, snRNAs, and snoRNA. The obtained sequences after screening were searched in the miRBase 21.0 database (http://www.mirbase.org) to identify known mature miRNAs. The sRNAs matched to known plant miRNAs with fewer than three mismatches, and could form typical stem-loop secondary structures, were identified as conserved miRNAs [37]. Sequences that were unmapped to any miRNA in the miRBase and had less than two mismatches were taken as novel miRNAs by miRDeep2 (https://www.mdc-berlin.de/8551903/en).

Differential expression analysis of miRNAs

miRNA expression level was calculated by Transcripts Per Kilobase Million (TPM), which can eliminate sequencing discrepancy [38]. Differential expression of miRNA between two samples was analyzed by DESeq2 [39, 40]. miRNAs were considered up-regulated when the fold change was greater than 2 or down-regulated when the Padj value (adjusted p-value) was less than 0.01. TPM and fold change were calculated as previously described [41].

Target gene predication and enrichment analysis

Target gene prediction of miRNAs was analyzed by psRobot (http://omicslab.genetics.ac.cn/psRobot) and TargetFinder (https://github.com/carringtonlab/TargetFinder) according to the related references [42, 43]. The targets of known and novel miRNAs with different expression between two samples were subjected to Gene Ontology (GO) annotations using the WEGO package [44].

Expression analysis of miRNAs and predicted target genes using qRT-PCR

Stem-loop quantitative reverse transcription PCR (qRT-PCR) was performed to analyze miRNA expression. According to the mature miRNA sequence, the specific stem-loop reverse transcription (RT) primers and forward primers of seven selected miRNAs were designed following the method described by Chen et al. [45]. The universal reverse primer was used for all qRT-PCR reactions and U6 snRNA was used as an endogenous reference gene. For qRT-PCR validation of target gene, first-strand cDNA was generated using a SuperScript®III RT Reagent Kit (Invitrogen, China, Beijing). The Actin of B. rapa was used as an endogenous reference gene for normalization. The primer sequences for validation of miRNA and target gene expression are listed in Additional file 1: Table S1 and Additional file 2: Table S2, respectively. qRT-PCR analysis was carried out on 7900HT Fast Real-Time PCR System (Applied Biosystems, USA) using SYBR qPCR Mix (Invitrogen). The PCR conditions were 95 °C for 2 min, followed by 40 cycles of 94 °C for 20 s, 60 °C for 20 s, and 72 °C for 30 s. Triplicate technical replicates were run for qPCR, and amplification products were further confirmed by analyzing the melting curve and gel electrophoresis. The relative expression of miRNAs was calculated using the 2−ΔΔCt method [46].

Results

Physiological responses to cold stress

The leaf and root tissues of ‘Longyou 7’ and ‘Tianyou 4’ were used to study cold responses. The degree of cold tolerance was estimated by measuring MDA content and POD activity, which were significantly increased in the leaves and roots of the two varieties under cold stress (Fig. 1). We found the MDA content and POD activity with no significant difference between two cultivars grown at 20 °C. However, the MDA content of 7LTR was higher than that of 4LTR. No significant difference of MDA content was found between the cold stressed roots. POD activity in the leaves and roots of ‘Longyou 7’ were higher than ‘Tianyou 4’ under cold stress.

Fig. 1
figure 1

The responses of MDA content and POD activity to cold stress

Overview of sRNA sequencing in winter turnip rapes

To identify the miRNAs responding to cold stress, 16 sRNA libraries from eight samples (two biological replicates for each sample) were constructed. These libraries yielded a total of 228,117,863 raw reads by high-throughput sequencing (approximately 14.3 million raw reads per library). After removing the low-quality reads, sequences with adapter and polyA, sequences less than 18 nucleotides and more than 30, we obtained 214,293,155 clean reads and an average of 13,393,322 reads (93.97%) per library (Table 1). Of these reads, 181,942,519 clean reads were mapped to the database of miRBase [47] and Rfam (Rfam 12.2). The mapped sRNAs were annotated and classified into structural RNA (rRNA, tRNA, snRNA, snoRNA), known miRNAs (miRNAs in miRBase 21.0), precursor RNA, intergenic RNA, and unknown sRNAs (Table 2). Notably, approximately 2% of the reads were mapped to mature miRNAs in each sample (Table 2). The proportion of known miRNAs increased in 7RTR compared with 7RCK, implying that a number of miRNAs were induced in the roots of ‘Longyou 7’ under cold stress. The distribution of various size of sRNAs showed that 24 nt sRNAs are the most abundant, accounting for 23.9% (4LTR) to 40.6% (7RTR), following by 23 nt sRNAs (Fig. 2). The distribution is consistent with previous reports in various plants [48, 49]. To verify the results, we analyzed the Pearson correlation coefficients of the miRNA expression between the biological duplicates (Additional file 3: Figure S1), and the Pearson r value of miRNA expression were about 0.8.

Table 1 The statistics of sRNA-seq reads and mapped reads
Table 2 Distribution of small RNAs among different categories in winter turnip rape
Fig. 2
figure 2

Sequence length distribution of sRNAs in winter turnip rape

Identification of conserved and novel miRNAs

Conserved miRNAs were identified by aligning against the miRBase 21.0 database. We obtained 269 conserved miRNAs belonging to 90 miRNA families from 16 libraries (Additional file 4: Table S3). And 96 members of them belonged to the known Bra-miRNAs, and 173 were similar to known miRNAs from other plant species. The result demonstrated that some conserved miRNAs in evolution also mediated the regulatory pathways of winter turnip rape. The expression profiles of the conserved miRNAs showed great variation in the amount of reads, which varied from 10,000 to less than 10 reads (Additional file 4: Table S3). Six miRNAs (miR158-3p, miR159a-1, miR162-3p, miR162a-3p, miR166, miR166a) presented high expression levels in all the samples. miR159a-1 was identified with the highest expression level, covering about 10,000 reads in 4LCK and 7LCK. In contrast, about 70 miRNAs were identified with less than 10 reads in all libraries. Besides, miR156 was identified as the largest miRNA family including 22 members, following by miR172 and miR319. Moreover, 222 and 226 known miRNAs were common between control and treatment from leaves of ‘Tianyou 4’ and ‘Longyou 7’, respectively (Fig. 3a). And 231 and 229 known miRNAs were common between control and treatment from roots of ‘Tianyou 4’ and ‘Longyou 7’, respectively (Fig. 3b), implying that known miRNAs are stable in different cultivars even under cold stress. The lengths of most of these miRNAs was 21 nt long, similar to those observed in previous reports [50, 51]. Characteristic information of all the conserved miRNA sequences was summarized in Additional file 4: Table S3.

Fig. 3
figure 3

Analysis of identified conserve microRNAs in leaves and roots from two winter turnip rape varieties

For the identification of novel miRNA candidates, the unannotated sRNAs matched to the B. rapa genome sequence were analyzed by miRDeep2. All the novel miRNA sequences were named in the form of Bra-novel-miR-number. Generally, we identified 3538 novel miRNAs (Additional file 5: Table S4), and 84 of them have star sequences (Additional file 6: Table S5). The lengths of novel miRNAs with star sequences varied from 20 nt to 24 nt, with the majority of 21 nt in length. The mean length of these miRNA precursors was 148 nt, ranging from 65 to 373 nt (Additional file 6: Table S5). Most of the novel miRNAs were identified with relatively low read counts compared with conserved miRNAs. Only six novel miRNAs were identified with more than 100 reads, and most of the remaining novel miRNAs had less than 10 reads (Additional file 6: Table S5). To confirm whether these sRNAs were true miRNAs of winter turnip rape, we analyzed their hairpin structures with Mfold or RNAfold [52]. The minimum free energies of these pre-miRNA hairpin structures varied from − 170.6 kcal·mol− 1 to − 19.1 kcal·mol− 1 with an average of − 58.9 kcal·mol− 1 similar to those of other plant miRNA precursors [48]. The secondary structures of the three most abundant novel miRNAs (Bra-novel-miR0167-3p, Bra-novel-miR0837-3p, and Bra-novel-miR3610-5p) were predicted, indicating they can form typical stem-loop hairpins, and the folding free energy was − 52.99, − 54.80, and − 32.40 kcal/mol, respectively (Fig. 4).

Fig. 4
figure 4

Secondary structures of differentially expressed novel miRNA precursors

Differentially expressed miRNAs under cold stress

The heat map presenting the expressional pattern of conserved and novel miRNAs in stressed and non-stressed leaves and roots of two B. rapa cultivars were shown in Additional file 7: Figure S2 and Additional file 8: Figure S3. Eleven differentially expressed miRNAs with more than two fold change were identified in ‘Longyou 7’ after cold stress, and nine differentially expressed miRNAs were identified in ‘Tianyou 4’ under cold stress (Additional file 9: Table S6). In detail, seven miRNAs (five known and two novel miRNAs) were up-regulated and three conserved miRNAs were down-regulated in leaves of ‘Longyou 7’ under cold stress. And Bra-novel-miR3936-5p was up-regulated in roots of ‘Longyou 7’ under cold stress. miR166e-3p was up regulated in leaves of ‘Tianyou 4’ under cold stress, while miR166a, miR167h and miR396-3p-3 were down-regulated in leaves of ‘Tianyou 4’ after cold stress. Three miRNAs (miR166e-3p, miR408-5p-2 and miR319–2) and two miRNAs (miR396a-3p-3 and miR166u) were up-regulated and down-regulated in stressed roots of ‘Tianyou 4’ compared with non-stressed roots of ‘Tianyou 4’, respectively (Fig. 5). We found two differential miRNAs (up-regulated miR166e-3p and down-regulated miR396a-3p-3) both in the leaves and roots of ‘Tianyou 4’ under cold stress.

Fig. 5
figure 5

Differential expression of miRNAs under cold stress

Comparative analysis of differentially expressed miRNAs between two varieties

The differentially expressed miRNAs between two varieties were compared. No significant difference was observed between non-stressed leaves of ‘Tianyou 4’ and ‘Longyou 7’, whereas miR166e-3p was significantly down-regulated and miR166h-3p was up-regulated in stressed leaves of ‘Longyou 7’compared with stressed leaves of ‘Tianyou 4′. We also identified two down-regulated miRNAs (miR319e and miR166m-2) between non-stressed roots of ‘Tianyou 4′ and ‘Longyou 7′, and eight differentially expressed miRNAs (Bra-novel-miR0040-3p, miR166e-3p, miR319e-1, miR319–2, and miR319a were down-regulated; Bra-novel-miR3153-5p, Bra-novel-m0894-3p and Bra-novel-m3936-5p were up-regulated) in stressed roots of ‘Longyou 7′ compared with stressed roots of ‘Tianyou 4′. We detected one co-expressed differential miRNA (miR166e-3p), which was down-regulated in both cold stressed leaves and roots of ‘Longyou 7′ compared with stressed organs of ‘Tianyou 4′ (Additional file 10: Table S7). By further analyzing these results, we found that all the differentially expressed known miRNAs of the two varieties belonged to miR166 and miR319 families. The identified members of miR319 were down-regulated in cold stressed leaf and root of ‘Longyou 7′ compared with ‘Tianyou 4′. However, members of the miR166 family showed differences in respect to their responses to cold stress. For example, miR166e-3p was down-regulated in the leaves and roots of ‘Longyou 7′ under cold stress, but miR166h-3p was up-regulated in the leaves. This result indicate that different expression pattern of miR166 members in different library might be related to their roles in the regulatory networks of cold stress.

Target prediction of differentially expressed miRNAs and gene enrichment analysis

To acknowledge the potential functions of differentially expressed miRNAs, we used psRobot and TargetFinder to predict the target genes [53, 54]. All the target genes were subjected to GO enrichment analysis. A total of 211 genes for fifteen known and two novel miRNAs were predicted as potential miRNA targets (Additional file 11: Table S8). GO annotation of these target genes showed that some stress-related GO terms were enriched, including response to stimulus and signaling. The most enriched GO terms were metabolic process, cellular process, cell, cell part, and binding (Fig. 6). The target genes of the differentially expressed miRNAs (Bra-novel-miR3936-5p) in 7RCK and 7RTR were involved in the regulation of metabolic processes and binding (Fig. 6). The diverse functions of the predicted target genes suggested that the differentially expressed miRNAs may play important roles during plant response to cold stress.

Fig. 6
figure 6

Gene Ontology (GO) categories and sub-categories of differentially expressed miRNA targets between the following samples: 7LCK versus7LTR, 4LCK versus 4LTR (A); 4RCK versus 4RTR, 7RCK versus 7RTR (B)

Validation of expression patterns of miRNAs and predicted target genes by qRT-PCR

For the verification of the sequencing results, five miRNAs (miR166e-3p, miR396a-3p-3, miR845a-1, miR319e-1, and Bra-novel-miR3936-5p) and the corresponding target genes (BraA06001949, BraA09000083, BraA09002428, BraA01004006, and BraA09006072) were selected for qRT-PCR analysis. The results showed that the qRT-PCR data of the differentially expressed miRNAs were in accordance with the sequencing results, and the expressional changes of target genes were negatively correlated to differentially expressed miRNAs (Fig. 7).

Fig. 7
figure 7

qRT-PCR analysis of differentially expressed miRNAs and the target genes

Discussion

Cold is a common stress facing by many plant species grown in northwest regions of China, which can induce cellular dehydration and membrane injury. B. rapa is a vegetative plant broadly grown in many regions, the selection and breeding of cold tolerance cultivars is of great importance to the vegetative supports in northwest regions. Plants enduring cold stress could adapt themselves with defense mechanisms, including physiological responses, adaptation in metabolic and signal transduction pathways, as well as gene expressional changes [55]. MDA is a product of lipid peroxidation, and the degree of membrane damage can be reflected by MDA level in the stressed tissue [56]. POD can maintain the balance of formation and elimination of ROS by scavenging excess ROS, which helps plants to endure stress conditions [29, 57]. Therefore, MDA and POD are important physiological indexes of stress responses. Under cold stress, the MDA content in leaves and POD activity in roots of ‘Longyou 7’ were higher than ‘Tianyou 4’, indicating cold stress led to more serious membrane damage in leaves of ‘Longyou 7’, but cold tolerance of ‘Longyou 7’ roots was stronger than ‘Tianyou 4’ under cold stress. These physiological changes under cold stress agreed to the previous reports [28, 58].

miRNAs as systemic regulators are involved in plant response to various abiotic stresses. Cold stress severely limits the geographic distribution of plant species by suppressing growth and development. Plants could adapt to cold stress by regulating gene expression. Moreover, 24 nt siRNAs can induce transcriptional gene silencing by DNA and histone modifications [59,60,61,62]. The relationship between miRNAs and cold response was reported in many plants [24, 25, 63, 64]. In winter turnip rape, the genome-wide expression profile of miRNAs in response to cold stress remains unreported. In the present study, millions of short reads were obtained from the sRNA libraries of two winter turnip rape cultivars under 20 °C and − 4 °C treatments. Similar to other plants, winter turnip rapes was enriched with 24 nt sRNAs [65, 66].The identified 269 conserved miRNAs (Additional file 4: Table S3) and 84 novel miRNAs (Additional file 6: Table S5) showed great variations in the amount of reads (from 1 to 94,237). The expression levels of novel miRNAs was relatively lower than that of conserved miRNAs, similar to the previous reports [67, 68]. The novel miRNAs could enrich the B. rapa miRNA database, which may also provide a foundation for cold-tolerant regulation of miRNA.

High-throughput sequencing platforms are widely used to identify miRNAs at a genome-wide scale. Our study was the first to identify differences between the miRNAs of the same variety under control and cold-stressed treatments. Cold stress triggered 10 differentially expressed miRNAs in the leaves of ‘Longyou 7’, but only one miRNA was up-regulated in the roots of ‘Longyou 7’. As previously reported, root tolerance to cold is very important for survival of ‘Longyou 7’ in winter [28]. Thus, we suspect the differentially expressed miRNA induced in roots of ‘Longyou 7’ might be valuable to its better survival rates compared with ‘Tianyou 4’. The aboveground part of ‘Longyou 7’ becomes more sensitive to cold stress, exhibiting poor growth and senescence [28]. Therefore, the number of differentially expressed miRNAs in the leaves of ‘Longyou 7’ was higher, which might be related to the solutions of ‘Longyou 7’ to alleviate cold induced injuries such as cytoplasmic damages caused by ice [55]. According to previous studies, minimal changes in miRNA expression can cause significant alterations in mRNA levels and related biological processes [69]. Thus, Bra-novel-miR3936-5p expressed in the roots may be a pivotal regulatory factor for the strong tolerance of ‘Longyou 7’ to cold stress. Furthermore, cold stress triggers four differentially expressed miRNAs in the leaves and five differentially expressed miRNAs in the roots of ‘Tianyou 4’. Of these miRNAs, two (miR166e-3p and miR396a-3p-3) were common in the leaves and roots of ‘Tianyou 4’, indicating that they were the most basic regulatory factors of cold stress in ‘Tianyou 4’. In summary, the differential expression of miRNAs implied the different response to cold stress in the two varieties.

In both varieties, no significant change was observed in the miRNAs expression of leaves, except for two miRNAs down-regulated in their roots under 20 °C. However, the number of differential miRNAs was increased under cold stress. Two differentially expressed miRNAs were found in the leaves and eight in the roots. These results indicated that the different expression levels of miRNAs under cold stress account for the different cold tolerance between two varieties. miRNAs with different expression pattern may provide a basis for the identification of molecular mechanism of cold tolerance in root of winter turnip rape. As reported, miRNAs regulate gene expression at the post-transcriptional level [70, 71]. The different expression levels of miRNAs between two varieties could attribute to explain the mechanisms of different cold tolerance. The root is the only overwintering organ of winter turnip rape, and its tolerance to cold stress is critical to their survival during winter. Under cold stress, eight miRNAs was differentially expressed between roots of ‘Longyou 7’ and roots of ‘Tianyou 4’, which might be related to the different root tolerance between two cultivars. We detected the expression of miR166e-3p was lower in the leaves and roots of ‘Longyou 7’ than in ‘Tianyou 4’, but its family member miR166h-3p was highly expressed in the leaves of ‘Longyou 7’ than in ‘Tianyou 4’, indicating that miR166 plays important role in the regulation of plant response to cold.

miRNAs participate in many biological processes through regulating the expression of target genes. Thus, identify the putative function of the target genes is important to acknowledge the miRNA-mediated mechanisms. We predicted 211 target genes with diverse functions for 15 known miRNAs and two novel miRNAs using bioinformatics approach (Additional file 11: Table S8). The majority of the predicted target genes encoded transcription factors, such as ATHB, MYB, GAMYB, TCP, bHLH transcription factor, Argonaute, and SQUAMOSA promoter binding proteins. Furthermore, other target genes, encoding glucosidase, glycylpeptide transferase, laccase, acyl-lipid desaturase, and galactinol synthase, are involved in various physiological and metabolic processes. Many other biological processes, including development, signal transduction, metabolism, and responses to stress are also regulated by miRNA.

In our study, the members of miR156, miR166 miR396 and miR845 family showed different expression pattern under cold stress (Additional file 9: Table S6). Squamosa-promoter binding protein-like (SPL) is the target gene of miR156, which has been reported with function in regulating leaf development [72, 73]. We found miRNA156a-5p was up-regulated in leaves of ‘Longyou 7’ under cold stress, whereas no difference of miR156a expression was observed in roots (Additional file 9: Table S6). The same results has been reported in barley under dehydration stress [74]. miR156 was up-regulated in Arabidopsis seedlings under salt stress [24]. The expressional changes of miR156a in leaves of ‘Longyou 7’ indicate it might be responsive to cold stress. Freeze injury can lead to plant dehydration and damage intracellular membrane system [75, 76]. Thus, we speculate that up-regulation of miRNA156a-5p might be resulted from the leaf dehydration and intracellular membrane damages in ‘Longyou 7’ under cold stress. And the aboveground of ‘Longyou 7’ could be influenced by the expressional changes of miR156 and the target genes, whereas there is no effect on roots.

miR845 has been reported with functions in regulating abiotic stress [77], but the mechanism of its role remains unclear. In our study, miR845a-1 targets a gene encoding probable protein S-acyltransferase 14 (PAT14) that involved in leaf senescence, and was over-expressed in senescent leaves of Arabidopsis [78]. Leaf senescence is a critical process in plant development and stress responses [79]. miR845a-1 was down-regulated in leaves of ‘Longyou 7’ under cold stress, and the target gene PAT14 was up-regulated (Fig. 7), which might be responsible for the earlier leaf senescence in ‘Longyou 7’ under cold stress.

Transcription factor TCP4-like (TCP4) is the main target gene of miR319. A recent research found miR319-mediated TCP4 controlled secondary cell wall formation by activating the expression of a downstream VND7 transcriptional network in A. thaliana. And the overexpression of TCP4 could significantly increase cell wall thickness, promote the formation of xylem vessel and deposition of lignin and cellulose [80]. Increase of cell wall thickness and lignification could increase cold tolerance in plants [81]. Under cold stress, we found three miR319 members (miR319e-1, miR319a and miR319–2) was down-regulated in roots of ‘Longyou 7’ compared with ‘Tianyou 4’, and expression of TCP4 was increased in roots of ‘Longyou 7’ (Fig. 7). These changes might lead to thicker cell wall in roots of ‘Longyou 7’ than that of ‘Tianyou 4’, which could be helpful to the cold tolerance of ‘Longyou 7’.

Homeobox-leucine zipper protein ATHB15, encoded by the target gene of miR166e-3p, is important to the synthesis of plant cell wall and cellulose [82]. In our study, the down-regulation of miR166e-3p and up-regulation of ATHB-15 in roots of ‘Longyou 7’ compared with ‘Tianyou 4’ under cold stress was identified (Fig. 7), which might improve cold tolerance in ‘Longyou7’ by increasing cell wall thickness and cellulose accumulation.

Most galactinol synthase (GolS) genes were reported to be induced by abiotic stresses, which can increase amounts of galactinol and raffinose to improve stress tolerance of plants [83,84,85,86]. In this study, we found Bra-novel-miR3936-5p targets GolS4 (BraA09006072), which could encode galactinol synthase and take part in the biosynthesis of raffinose oligosaccharides (ROFs) (Additional file 11: Table S8). Interestingly, Bra-novel-miR3936-5p was up-regulated in roots of ‘Longyou 7’ under cold stress compared with ‘Tianyou 4’, which led to the down-regulation of Gols4 (BraA09006072). In Populus trichocarpa, PtrGolS4 was down-regulated under cold stress, whereas the expression of some GolS homologous was induced, indicating PtrGolS4 with no contribution to the accumulation of galactinol and raffinose under cold stress [87]. Feng et al. (2017) found that VvGolS4 could participate in H2S-mediated pathways that was responsive to temperature changes [88]. The expressional changes of Bra-novel-miR3936-5p and GolS4 in ‘Longyou 7’ might be related to cold stress responses involving some metabolic pathways. Bra-novel-miR3936-5p was confirmed as a cold stress-induced up-regulated miRNA in roots of ‘Longyou 7’, this is valuable for the dissection of cold tolerance in ‘Longyou 7’ compared with ‘Tianyou 4’.

Overall, we identified and compared cold-responsive miRNAs in the leaves and roots of two winter turnip rape varieties using high-throughput sequencing, which were selected for stem-loop qRT-PCR validation. The data could lay a base for the identification of miRNAs related to cold tolerance. Prediction and functional annotation of target genes were carried out to known the putative function of miRNAs in cold stress responses. However, whether these differentially expressed miRNAs are true regulators in cold tolerance of ‘Longyou 7’, especially in the roots, need more experimental validations. For instance, genetic transformation of precursors of these candidate miRNAs into ‘Tianyou 4’ to validate its molecular function, accompany cold-tolerant phenotypes analysis in a segregating population with different miRNA expression levels, would help to confirm the molecular function of these miRNAs in regulating cold stress responses.

Conclusions

In the present study, small RNA libraries of cold stressed and non-stressed leaves and roots of two winter turnip rape varieties were constructed separately, including two replicates for each sample. As a result, a set of conserved and novel miRNAs were identified by high throughput sequencing. The differentially expressed miRNAs were screened to find stress responsive miRNAs by comparing different organs and cultivars, as well as different treatments (cold stressed and non-stressed). The results indicated many candidate miRNAs and target genes with potential functions in regulating cold tolerance in ‘Longyou 7’, which are valuable for dissecting post-transcriptional and transcriptional regulations of cold responses, and B. rapa improvement in cold regions in the world.

Abbreviations

ATHB15:

Homeobox-leucine zipper protein 15

CK:

Control

GO:

Gene ontology

GolS:

Galactinol synthase

MDA:

Malondialdehyde

miRNA:

microRNA

nt:

Nucleotide

PAT14:

Probable protein S-acyltransferase14

PCR:

Polymerase chain reaction

POD:

Peroxidase

qRT-PCR:

Quantitative real time polymerase chain reaction

RISCs:

RNA-induced silencing complexes

ROS:

Reactive oxygen species

rRNA:

Ribosomal RNA

snoRNA:

Small nucleolar RNA

snRNA:

Small nuclear RNA

sRNA:

Small RNA

TCP4:

Transcription factor TCP4-like

TF:

Transcription factor

TPM:

Transcripts per kilobase million

TR:

Treatment

tRNA:

Transfer RNA

UTR:

Untranslated region

WEGO:

Web gene ontology annotation plot

References

  1. Arbaoui M, Balko C, Link W. Study of faba bean (Viciafaba L.) winter-hardiness and development of screening methods. Field Crops Res. 2008;106(1):60–7.

    Article  Google Scholar 

  2. Kerepesi I, Bányai-Stefanovits E, Galiba G. Cold acclimation and abscisic acid induced alterations in carbohydrate content in call of wheat genotypes differing in frost tolerance. J Plant Physiol. 2004;161(1):131–3.

    Article  CAS  PubMed  Google Scholar 

  3. Frentzen M, Nishida I, Murata N. Properties of the plastidial acyl-(acyl-carrier-protein): glycerol-3-phosphate acyltransferase from the chilling sensitive plant squash (Cucurbita moschata). Plant Cell Physiol. 1987;28(7):1195–201.

    CAS  Google Scholar 

  4. Yeh S, Moffatt BA, Griffith M, Xiong F, Daniel SC, Yang SB, Sarhan F, Danyluk J, Xue YQ, Hew CL, et al. Chitinase genes responsive to cold encode antifreeze proteins in winter cereals. Plant Physiol. 2000;124(3):1251–63.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  5. Medina J, Catalá R, Salinas J. Developmental and stress regulation of RCI2A and RCI2B, two cold-inducible genes of Arabidopsis encoding highly conserved hydrophobic proteins. Plant Physiol. 2001;125(4):1655–66.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  6. Krutzfeldt J, Stoffel M. MicroRNAs: a new class of regulatory gene affecting metabolism. Cell Metab. 2006;4(1):9–12.

    Article  CAS  PubMed  Google Scholar 

  7. Bartel B, Bartel DP. MicroRNAs: at the root of plant development? Plant Physiol. 2003;132(2):709–17.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  8. 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(5):1360–75.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  9. Biggar KK, Storey KB. The emerging roles of microRNAs in the molecular responses of metabolic rate depression. J Mol Cell Biol. 2010;3(3):167–75.

    Article  PubMed  Google Scholar 

  10. Biggar KK, Storey KB. Insight into post-transcriptional gene regulation: stress-responsive microRNAs and their role in the environmental stress survival of tolerant animals. J Exp Biol. 2015;218(Pt9):1281–9.

    Article  PubMed  Google Scholar 

  11. Bartel DP. MicroRNAs: genomics, biogenesis, mechanism, and function. Cell. 2004;116(2):281–97.

    Article  CAS  PubMed  Google Scholar 

  12. Cho Y, Lim JH, Jeong SC, Song MK, Ryu JC. Hexanal-induced changes in miRNA-mRNA interactions in A549 human alveolar epithelial cells. Toxicol Environ Heal Sci. 2015;7(2):143–59.

    Article  Google Scholar 

  13. Son GW, Yang H, Park HR, Lee SE, Jin YH, Park CS, Park YS. Analysis of miRNA expression profiling in melatonin-exposured endothelial cells. Mol Cell Toxicol. 2016;12(1):73–81.

    Article  CAS  Google Scholar 

  14. Bartel DP. MicroRNAs: target recognition and regulatory functions. Cell. 2009;136(2):215–33.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  15. Voinnet O. Origin, biogenesis, and activity of plant microRNAs. Cell. 2009;136(4):669–87.

    Article  CAS  PubMed  Google Scholar 

  16. Siomi H, Siomi MC. Posttranscriptional regulation of microRNA biogenesis in animals. Mol Cell. 2010;38(3):323–32.

    Article  CAS  PubMed  Google Scholar 

  17. Morgan CP, Bale TL. Sex differences in microRNA regulation of gene expression: no smoke, just miRs. Biol Sex Differ. 2012;3(1):22.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  18. Riviere G, Wu GC, Fellous A, Goux D, Sourdaine P, Favrel P. DNA methylation is crucial for the early development in the Oyster C. gigas. Mar Biotechnol. 2013;15(6):739–53.

    Article  CAS  PubMed  Google Scholar 

  19. Ding YF, Zhu C. The role of microRNAs in copper and cadmium homeostasis. Biochem Bioph Res Co. 2009;386(1):6–10.

    Article  CAS  Google Scholar 

  20. Wei LY, Zhang DF, Xiang F, Zhang ZX. Differentially expressed miRNAs potentially involved in the regulation of defense mechanism to drought stress in maize seedlings. Int J Plant Sci. 2009;170(8):979–89.

    Article  CAS  Google Scholar 

  21. Navarro L, Dunoyer P, Jay F, Arnold B, Dharmasiri N, Estelle M, Voinnet O, Jones JD. A plant miRNA contributes to antibacterial resistance by repressing auxin signaling. Science. 2006;312:436–9.

    Article  CAS  PubMed  Google Scholar 

  22. Sunkar R, Zhu JK. Novel and stress-regulated microRNAs and other small RNAs from Arabidopsis. Plant Cell. 2004;16(8):2001–19.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  23. Samuel MA, Mudgil Y, Salt JN, Delmas F, Ramachandran S, Chilelli A, Goring DR. Interactions between the S-domain receptor kinases and at PUB-ARM E3 ubiquitin ligases suggest a conserved signaling pathway in Arabidopsis. Plant Physiol. 2008;147(4):2084–95.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. Liu HH, Tian X, Li YJ, Wu CA, Zheng CC. Microarray-based analysis of stress regulated microRNAs in Arabidopsis thaliana. RNA. 2008;14(5):836–43.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  25. Zhang JY, Xu YY, Hua Q, Chong K. Deep sequencing of Brachypodium small RNAs at the global genome level identifies microRNAs involved in cold stress response. BMC Genomics. 2009;10:449.

    Article  PubMed  PubMed Central  Google Scholar 

  26. Lv DK, Bai X, Li Y, Ding XD, Ge Y, Cai H, Ji W, Wu N, Zhu YM. Profiling of cold-stress-responsive miRNAs in rice by microarrays. Gene 2010; 459(1–2):39–47.

    Article  CAS  PubMed  Google Scholar 

  27. Huang R, Cang J, Yu J, Lu BW, Liu LJ, Wang JF, Guo RM, Xu C. Solexa sequencing and bioinformatics analysis of small RNA in winter wheat. Chin Bull Bot. 2014;49(1):8–18.

    Article  CAS  Google Scholar 

  28. Sun WC, Wu JY, Fang Y, Liu Q, Yang RY, Ma WG, Wang XF, Ye J, Zeng J, Zhang YH, et al. Growth and development characteristics of winter rapeseed northern-extended from the cold and arid regions in China. Acta Agron Sin. 2010;36(12):2124–34.

    Article  Google Scholar 

  29. Liu ZG, Sun WC, Yang NN, Wang Y, He L, Zhao CX, Shi PF, Yang G, Li XC, Wu JY, et al. Morphology and physiological characteristics of cultivars with different levels of cold-resistance in winter rapeseed (Brassica campestris L.) during cold acclimation. Sci Agric Sin. 2013;46(22):4679–87.

    CAS  Google Scholar 

  30. Zeng XC, Liu ZG, Shi PH, Xu YZ, Sun J, Fang Y, Yang G, Wu JY, Kong DJ, Sun WC. Cloning and expression analysis of copper and zinc superoxide dismutase (Cu/Zn-SOD) gene from Brassica campestris L. Acta Agron Sin. 2014;40(4):636–43.

    Article  CAS  Google Scholar 

  31. Zeng XC, Sun WC, Sun J, Xu YZ, Fang Y, Shi PH, Yang G, Kong DJ, Wu JY, Liu ZG. Cloning and expression analysis of fe superoxide dismutase (Fe-SOD) gene from winter turnip rape (Brassica rapa L.). Sci Agric Sin. 2013;46(21):4603–11.

    CAS  Google Scholar 

  32. Campos PS, Quartin V, Ramalho JC, Nunes MA. Electrolyte leakage and lipid degradation account for cold sensitivity in leaves of Coffea sp. plants. Plant Physiol. 2003;160(3):283–92.

    Article  CAS  Google Scholar 

  33. Quiroga M, Guerrero C, Botella MA, Barceló A, Amaya I, Medina MI, Alonso FJ, Forchetti SM, Tigier H, Valpuesta V. A tomato peroxidase involved in the synthesis of lignin and suberin. Plant Physiol. 2000;122(4):1119–27.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  34. Chen HY, Chen XL, Chen D, Li JF, Zhang Y, Wang A. A comparison of the low temperature transcriptomes of two tomato genotypes that differ in freezing tolerance: Solanum lycopersicum and Solanum habrochaites. BMC Plant Biol. 2015;15:132–47.

    Article  PubMed  PubMed Central  Google Scholar 

  35. Li H, Hu B, Wang W, Zhang ZH, Liang Y, Gao XK, Li P, Liu YQ, Zhang LH, Chu CC. Identification of microRNAs in rice root in response to nitrate and ammonium. J Genet Genomics. 2016;43(11):651–61.

    Article  PubMed  Google Scholar 

  36. Yi SK, Gao ZX, Zhao HH, Zeng C, Luo W, Chen BX, Wang WM. Identification and characterization of microRNAs involved in growth of blunt snout bream (Megalobrama amblycephala) by Solexa sequencing. BMC Genomics. 2013;5(14):754–65.

    Article  Google Scholar 

  37. Kozomara A, Griffiths-Jones S. miRBase: integrating microRNA annotation and deep-sequencing data. Nucleic Acids Res. 2011;39(suppl_1):152–7.

    Article  Google Scholar 

  38. t Hoen PA, Ariyurek Y, Thygesen HH, Vreugdenhil E, Vossen RH, de Menezes RX, Boer JM, van Ommen GJ, den Dunnen JT. Deep sequencing-based expression analysis shows major advances in robustness, resolution and inter-lab portability over five microarray platforms. Nucleic Acids Res. 2008;36(21):e141.

    Article  Google Scholar 

  39. 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.

    Article  PubMed  PubMed Central  Google Scholar 

  40. Wang LK, Feng ZX, Wang X, Wang XW, Zhang XG. DEGseq: an R package for identifying differentially expressed genes from RNA-seq data. Bioinformatics. 2010;26(1):136–8.

    Article  PubMed  Google Scholar 

  41. Yu HL, Cong L, Zhu ZX, Wang CY, Zou JQ, Tao CG, Shi ZS, Lu XC. Identification of differentially expressed microRNA in the stems and leaves during sugar accumulation in sweet sorghum. Gene. 2015;571(2):221–30.

    Article  CAS  PubMed  Google Scholar 

  42. Allen E, Xie ZX, Gustafson AM, Carrington JC. MicroRNA-directed phasing during trans-acting siRNA biogenesis in plants. Cell. 2005;121(2):207–21.

    Article  CAS  PubMed  Google Scholar 

  43. Schwab R, Palatnik JF, Riester M, Schommer C, Schmid M, Weigel D. Specific effects of microRNAs on the plant transcriptome. Cell. 2005;8(4):517–27.

    CAS  Google Scholar 

  44. Ye J, Fang L, Zheng H, Zhang Y, Chen J, Zhang Z, Wang J, Li S, Li R, Bolund L, Wang J, et al. WEGO: a web tool for plotting GO annotations. Nucleic Acids Res. 2006;34:W293–7.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  45. Chen C, Ridzon DA, Broomer AJ, Zhou Z, Lee DH, Nguyen JT, Barbisin M, Xu NL, Mahuvakar VR, Andersen MR, et al. Real-time quantification of microRNAs by stem-loop RT-PCR. Nucleic Acids Res. 2005;33(20):179–87.

    Article  Google Scholar 

  46. Livak KJ, Schmittgen TD. Analysis of relative gene expression data using real-time quantitative PCR and the 2 −ΔΔC T method. Methods. 2001;25(4):402–8.

    CAS  PubMed  Google Scholar 

  47. Kozomara A, Griffiths-Jones S. MiRBase: annotating high confidence microRNAs using deep sequencing data. Nucleic Acids Res. 2014;42(D1):D68–73.

    Article  CAS  PubMed  Google Scholar 

  48. Wei XC, Xu W, Yuan YX, Yao QJ, Zhao YY, Wang ZY, Jiang WS, Zhang XW. Genome-wide investigation of microRNAs and their targets in Brassica rapa ssp. pekinensis root with Plasmodiophora brassicae infection. Horticultural Plant J. 2016;2(4):209–16.

    Article  Google Scholar 

  49. Liu ZB, Zhang YP, Ou LP, Kang LY, Liu YH, Lv JH, Wei G, Yang BZ, Yang S, Chen WC, et al. Identification and characterization of novel microRNAs for fruit development and quality in hot pepper (Capsicum annuum L.). Gene. 2017;608:66–72.

    Article  CAS  PubMed  Google Scholar 

  50. Ye XL, Song TF, Liu C, Feng H, Liu ZY. Identification of fruit related microRNAs in cucumber (Cucumis sativus L.) using high-throughput sequencing technology. Hereditas. 2014;151(6):220–8.

    Article  PubMed  Google Scholar 

  51. Lelandais-Briere C, Naya L, Sallet E, Calenge F, Frugier F, Hartmann C, Gouzy J, Crespia M. Genome-wide medicago truncatula small RNA analysis revealed novel microRNAs and isoforms differentially regulated in roots and nodules. Plant Cell. 2009;21(9):2780–96.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  52. Zuker M. Mfold web server for nucleic acid folding and hybridization prediction. Nucleic Acids Res. 2003;31(13):3406–15.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Wu HJ, Ma YK, Chen T, Wang M, Wang XJ. PsRobot: a web-based plant small RNA meta-analysis toolbox. Nucleic Acids Res. 2012;40(W1):W22–8.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Fahlgren N, Carrington JC. MiRNA target prediction in plants. Methods Mol Biol. 2010;592:51–7.

    Article  CAS  PubMed  Google Scholar 

  55. Knight MR, Knight H. Low-temperature perception leading to gene expression and cold tolerance in higher plants. New Phytol. 2012;195:737–51.

    Article  CAS  PubMed  Google Scholar 

  56. Upadhyaya H, Khan MP. Hydrogen peroxide induces oxidative stress in detached leaves of Oryza sativa L. Gen Appl Plant Physiol. 2007;33:83–95.

    CAS  Google Scholar 

  57. Zhang MC, Duan LS, Tian XL, He ZP, Li JM, Wang BM, Li ZH. Uniconazole-induced tolerance of soybean to water deficit stress in relation to changes in photosynthesis, hormones and antioxidant system. Plant Physiol. 2007;164:709–17.

    Article  CAS  Google Scholar 

  58. Sun WC, Wu JY, Zeng J, Zhu HX, Liu YL, Zhang YH. Primary evaluation of cold tolerance among eight winter Brassica rapa. J Hunan Agric Univ. 2007;33:151–5.

    Google Scholar 

  59. Lu S, Sun Y, Chiang VL. Stress-responsive microRNAs in Populus. Plant J. 2008;55(1):131–51.

    Article  CAS  PubMed  Google Scholar 

  60. Melnikova NV, Dmitriev AA, Belenikin MS, Speranskaya AS, Krinitsina AA, Rachinskaia OA, Lakunina VA, Krasnov GS, Snezhkina AV, Sadritdinova AF, et al. Excess fertilizer responsive miRNAs revealed in Linum usitatissimum L. Biochimie. 2015;109:36–41.

    Article  CAS  PubMed  Google Scholar 

  61. Hamilton A, Voinnet O, Chappell L, Baulcombe D. Two classes of short interfering RNA in RNA silencing. EMBO J. 2002;21:4671–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  62. Zilberman D, Cao XF, Jacobsen SE. ARGONAUTE4 control of locus-specific siRNA accumulation and DNA and histone methylation. Science. 2003;299(5607):716–9.

    Article  CAS  PubMed  Google Scholar 

  63. Leviatan N, Alkan N, Leshkowitz D, Fluhr R. Genome-wide survey of cold stress regulated alternative splicing in Arabidopsis thaliana with tiling microarray. PLoS One. 2013;8(6):e66511.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  64. Sunkar R, Li YF, Jagadeeswaran G. Functions of microRNAs in plant stress responses. Trends Plant Sci. 2012;17(4):196–203.

    Article  CAS  PubMed  Google Scholar 

  65. Fahlgren N, Howell MD, Kasschau KD, Chapman EJ, Sullivan CM, Cumbie JS, Givan SA, Law TF, Grant SR, Dangl JL, et al. High-throughput sequencing of Arabidopsis microRNAs: evidence for frequent birth and death of MIRNA genes. PLoS One. 2007;2(2):e219.

    Article  PubMed  PubMed Central  Google Scholar 

  66. Martinez G, Forment J, Llave C, Pallás V, Gómez G. High-throughput sequencing, characterization and detection of new and conserved cucumber miRNAs. PLoS One. 2011;6(5):e19523.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  67. Mao WH, Li ZY, Xj X, Li YD, Yu JQ. A combined approach of high-throughput sequencing and degradome analysis reveals tissue specific expression of microRNAs and their targets in cucumber. PLoS One. 2012;7(3):e33040.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  68. Chi XY, Yang QL, Chen XP, Wang JY, Pan LJ, Chen MN, Yang Z, He Y, Liang XQ, Yu SL. Identification and characterization of microRNAs from peanut (Arachis hypogaea L.) by high-throughput sequencing. PLoS One. 2011;6(11):e 27530.

    Article  CAS  Google Scholar 

  69. Farh KK-H, Grimson A, Jan C, Lewis BP, Johnston WK, Lim LP, Burge CB, Bartel DP. The widespread impact of mammalian microRNAs on mRNA repression and evolution. Science. 2005;310(5755):1817–21.

    Article  CAS  PubMed  Google Scholar 

  70. Cho Y, Lim JH, Jeong SC, Song MK, Ryu JC. Hexanal-induced changes in miRNA-mRNA interactions in A549 human alveolar epithelial cells. Toxicol Environ Heal Sci. 2015;7:143–59.

    Article  Google Scholar 

  71. Son GW, Yang H, Park HR, Lee SE, Jin YH, Park CS, Park YS. Analysis of miRNA expression profiling in melatonin-exposured endothelial cells. Mol Cell Toxicol. 2016;12:73–81.

    Article  CAS  Google Scholar 

  72. Wu G, Poethig R. Temporal regulation of shoot development in Arabidopsis thaliana by miR156 and its target SPL3. Development. 2006;133:3539–47.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  73. Bai SL, Saito T, Ito A, Tuan PA, Xu Y, Teng YW, Moriguchi T. Small RNA and PARE sequencing in flower bud reveal the involvement of sRNAs in endodormancy release of Japanese pear( Pyrus pyrifolia Kosui). BMC Genomics. 2016;17:230.

    Article  PubMed  PubMed Central  Google Scholar 

  74. Kantar M, Unver T, Budak H. Regulation of barley miRNAs upon dehydration stress correlated with target gene expression. Funct Integr Genomics. 2010;10:493–507.

    Article  CAS  PubMed  Google Scholar 

  75. Steponkus PL. Role of plasma membrane in freezing injury and cold acclimation. Annu Rev Plant Physiol. 1984;35:543–84.

    Article  CAS  Google Scholar 

  76. Steponkus PL, Uemura M, Webb MS. Membrane destabilization during freezing-induced dehydration. Curr Topics Plant Physiol. 1993;10:37–47.

    CAS  Google Scholar 

  77. Zhou LG, Liu YH, Liu ZC, Kong D, Duan M, Luo LJ. Genome-wide identification and analysis of drought-responsive microRNAs in Oryza sativa. J Exp Bot. 2010;61(15):4157–68.

    Article  CAS  PubMed  Google Scholar 

  78. Lai JB, Boya Y, Cao ZD, Chen YM, Wu Q, Huang JY, Yang CW. Two homologous protein S-acyltransferases, PAT13 and PAT14, cooperatively regulate leaf senescence in Arabidopsis. J Exp Bot. 2015;66(20):6345–53.

    Article  CAS  PubMed  Google Scholar 

  79. Lim PO, Kim HJ, Nam HG. Leaf senescence. Annu Rev Plant Biol. 2007;58:115–36.

    Article  CAS  PubMed  Google Scholar 

  80. Sun XD, Wang CD, Xiang N, Li X, Yang SH, Du JC, Yang YP, Yang YQ. Activation of secondary cell wall biosynthesis by miR319-targeted TCP4 transcription factor. Plant Biotech J. 2017;15:1284–94.

    Article  CAS  Google Scholar 

  81. Rajashekar CB, Lafta A. Cell-wall changes and cell tension in response to cold acclimation and exogenous abscisic acid in leaves and cell cultures. Plant Physiol. 1996;111:605–12.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  82. Kim J, Jung JH, Reyes JL, Kim YS, Kim SY, Kim JA CKS, Lee M, Lee Y, Kim VN, et al. MicroRNA directed cleavage of ATHB15 mRNA regulates vascular development in Arabidopsis inflorescence stems. Plant J. 2005;42(1):84–94.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  83. Taji T, Ohsumi C, Iuchi S, Seki M, Kasuga M, Kobayashi M, Yamaguchi-Shinozaki K, Shinozaki K. Important roles of drought- and cold-inducible genes for galactinol synthase in stress tolerance in Arabidopsis thaliana. Plant J. 2002;29:417–26.

    Article  CAS  PubMed  Google Scholar 

  84. Takahashi R, Joshee N, Kitagawa Y. Induction of chilling resistance by water stress, and cDNA sequence analysis and expression of water stress regulated genes in rice. Plant Mol Biol. 1994;26:339–52.

    Article  CAS  PubMed  Google Scholar 

  85. Sprenger N, Keller F. Allocation of raffinose family oligosaccharides to transport and storage pools in Ajuga reptans: the roles of two distinct galactinol synthases. Plant J. 2000;21:249–58.

    Article  CAS  PubMed  Google Scholar 

  86. Cunningham SM, Nadeau P, Castonguay Y, Laberge S, Volenec JJ. Raffinose and stachyose accumulation, galactinol synthase expression, and winter injury of contrasting alfalfa germplasms. Crop Sci. 2003;43:562–70.

    Article  CAS  Google Scholar 

  87. Zhou J, Yang Y, Yu J, Wang LK, Yu X, Ohtani M, Kusano M, Saito K, Demura T, Zhuge Q. Responses of Populus trichocarpa galactinol synthase genes to abiotic stresses. J Plant Res. 2014;127:347–58.

    Article  CAS  PubMed  Google Scholar 

  88. Feng RJ, Hou LX, Gou Y, Ma Q, Liu X. Gene cloning and expression analysis of galactinol synthase gene VvGolS2-4 in grape. Acta Agriculturae boreali-sinica. 2017;32(2):87–95.

    Google Scholar 

Download references

Acknowledgements

We sincerely thank the College of Agronomy, Gansu Agricultural University for kindly providing seeds of two turnip rapa.

Funding

This work was supported by the National Natural Science Foundation of China (Grant No. 31560397, 31660401, 31460356), National Key Basic Research Program of China (2015CB150201), the Natural Science Foundation of Gansu Province (1506RJZG051, 145RJZG050), and the Priority Academic Program Development of Jiangsu Higher Education Institutions.

Availability of data materials

All the data pertaining to the present study has been included in the tables and figures of the manuscript, and the raw data of sRNA-seq has been submitted to NCBI (Accession: PRJNA412963). The authors are pleased to share all the data and plant materials upon reasonable request.

Author information

Authors and Affiliations

Authors

Contributions

YW and WS conceived and designed the study. XZ, YX, FZ and LM participated in the experiments. JJ and DW analyzed the data. All the authors drafted and approved the final manuscript.

Corresponding authors

Correspondence to Youping Wang or Wancang Sun.

Ethics declarations

Ethics approval and consent to participate

Not applicable.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Publisher’s Note

Springer Nature remains neutral with regard to jurisdictional claims in published maps and institutional affiliations.

Additional files

Additional file 1:

Table S1. Primers for qRT-PCR analysis of differentially expressed miRNAs (DOCX 14 kb)

Additional file 2:

Table S2. Primers for qRT-PCR analysis of target genes of differentially expressed miRNAs (DOC 32 kb)

Additional file 3:

Figure S1. The correlation analysis of miRNA expression between two biological duplicates (JPEG 3923 kb)

Additional file 4:

Table S3. The conserved miRNAs identified in winter turnip rape (XLS 76 kb)

Additional file 5:

Table S4. The novel miRNAs identified in winter turnip rape (XLS 1404 kb)

Additional file 6:

Table S5. The novel miRNAs identified in winter turnip rape with star sequences (XLS 58 kb)

Additional file 7:

Figure S2. Heat map of conserved miRNAs identified in winter turnip rape (JPEG 1334 kb)

Additional file 8:

Figure S3. Heat map of novel miRNAs identified in winter turnip rape (JPEG 526 kb)

Additional file 9:

Table S6. Differentially expressed miRNAs between cold-stressed and non-stressed samples of two winter turnip rape varieties (DOC 53 kb)

Additional file 10:

Table S7. Differentially expressed miRNAs between two varieties (DOC 41 kb)

Additional file 11:

Table S8. Targets prediction of all the sequenced miRNAs (XLS 57 kb)

Rights and permissions

Open Access This article is distributed under the terms of the Creative Commons Attribution 4.0 International License (http://creativecommons.org/licenses/by/4.0/), which permits unrestricted use, distribution, and reproduction in any medium, provided you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons license, and indicate if changes were made. The Creative Commons Public Domain Dedication waiver (http://creativecommons.org/publicdomain/zero/1.0/) applies to the data made available in this article, unless otherwise stated.

Reprints and permissions

About this article

Check for updates. Verify currency and authenticity via CrossMark

Cite this article

Zeng, X., Xu, Y., Jiang, J. et al. Identification of cold stress responsive microRNAs in two winter turnip rape (Brassica rapa L.) by high throughput sequencing. BMC Plant Biol 18, 52 (2018). https://doi.org/10.1186/s12870-018-1242-4

Download citation

  • Received:

  • Accepted:

  • Published:

  • DOI: https://doi.org/10.1186/s12870-018-1242-4

Keywords