Skip to main content
Download PDF

Response of microRNAs to cold treatment in the young spikes of common wheat

BMC Genomics volume 18, Article number: 212 (2017) Cite this article

Abstract

Background

MicroRNAs (miRNAs) are a class of small non-coding RNAs that play important roles in biotic and abiotic stresses by regulating their target genes. For common wheat, spring frost damage frequently occurs, especially when low temperature coincides with plants at early floral organ differentiation, which may result in significant yield loss. Up to date, the role of miRNAs in wheat response to frost stress is not well understood.

Results

We report here the sequencing of small RNA transcriptomes from the young spikes that were treated with cold stress and the comparative analysis with those of the control. A total of 192 conserved miRNAs from 105 families and nine novel miRNAs were identified. Among them, 34 conserved and five novel miRNAs were differentially expressed between the cold-stressed samples and the controls. The expression patterns of 18 miRNAs were further validated by quantitative real time polymerase chain reaction (qRT-PCR). Moreover, nearly half of the miRNAs were cross inducible by biotic and abiotic stresses when compared with previously published work. Target genes were predicted and validated by degradome sequencing. Gene Ontology (GO) enrichment analysis showed that the target genes of differentially expressed miRNAs were enriched for response to the stimulus, regulation of transcription, and ion transport functions. Since many targets of differentially expressed miRNAs were transcription factors that are associated with floral development such as ARF, SPB (Squamosa Promoter Binding like protein), MADS-box (MCM1, AG, DEFA and SRF), MYB, SPX (SYG1, Pho81 and XPR1), TCP (TEOSINTE BRANCHED, Cycloidea and PCF), and PPR (PentatricoPeptide Repeat) genes, cold-altered miRNA expression may cause abnormal reproductive organ development.

Conclusion

Analysis of small RNA transcriptomes and their target genes provide new insight into miRNA regulation in developing wheat inflorescences under cold stress. MiRNAs provide another layer of gene regulation in cold stress response that can be genetically manipulated to reduce yield loss in wheat.

Background

Low temperature is one of the most important environmental stimuli that affect plant growth and development. Common wheat, as a classical vernalization-required plant, is tolerant to lower temperatures at the vegetative development stage. However, it will suffer frost injury when low temperatures coincide with the reproductive stage of wheat. In spring, the low temperature occurring at the stage when pistil and stamen differentiate into anthers will cause anther sterility and floret abortion, leading to significant yield loss in wheat [1]. In the vegetative developmental stage, vernalization genes VRN1, VRN2, and VRN3 regulate the transition between the vegetative and reproductive stages and interconnect with the cold acclimation regulation locus FR2 including the CBF and COR genes. These genes have been characterized to respond to freezing tolerance [2, 3]. For the spring freeze, stress-related, photosynthesis, and plastid-associated genes are demonstrated as responsive genes by transcriptomic and proteomic analyses [4]. For example, a set of spring freeze-responsive genes were identified by using the Affymetrix GeneChip in barley [5], while using a similar approach a variety of genes are reported to be involved in wheat spring freeze stimulus including WCOR413, LEA, glycine-rich RNA-binding protein, ferritin, aquaporin 2, and a pathogen-induced protein, the ice recrystallization protein, cold-related proteins, CBF transcription factors, calcium-dependent protein kinases, Na+/H+ antiporters, aquaporins, and many metabolic enzymes [6]. Recently, RNA-seq and digital gene expression analysis also identified a series of protein-coding genes from cold-treated young spikes in wheat [7]. In general, frost resistances at the vegetative stage and the reproductive stage seem to confer overlapping regulatory networks such as the one mediated by the CBF gene. Despite this, gene networks for frost resistance at the reproductive stage appear to be more complicated and involve more physiological processes.

MicroRNAs (miRNAs) are classical small non-coding RNAs that guide post-transcriptional gene regulation. MiRNAs have been shown to play key roles in various biological processes, including development, hormone responses and stress adaptation [8,9,10,11]. Recently, many studies have demonstrated the role of plant miRNAs in cold stress response. In Arabidopsis, 16 miRNAs, including miR156, miR159, miR164, miR165, miR168, miR169, miR172, miR319, miR389, miR393, miR396, miR397, miR398, miR400, miR402, and miR408, were identified by RNA gel blot analysis [12], microarray analysis [13], and a computation-based approach to be related to cold response [14]. In Populus, 19 cold stress-responsive miRNAs were identified by miRNA microarray [15], among which miR156, miR164, miR168, miR169, miR393, and miR396 were overlapped with those in Arabidopsis. In addition, 25 species-specific miRNAs were identified as cold regulated in Brachypodium by high-throughput sequencing [16]. These works indicate conserved mechanisms for cold-responsive miRNAs as well as the diverged species-specific regulation. In common wheat (Triticum aestivum), many conserved, Triticum-specific, and wheat-specific miRNAs have been identified, including the wheat A genome progenitor, the wheat D genome progenitor, and the AABBDD hexaploid wheat [17,18,19,20,21,22,23,24,25,26]. MiRNAs have been identified to play important roles in various stresses in wheat [9], including abiotic stresses such as salt [27, 28], drought [28, 29], dehydration [30], phosphorus [31], heat [32] and biotic stresses such as the fungal infection of powdery mildew [32] and Puccinia striiformis f. sp. Tritici [33]. So far, six miRNAs, miR167c, tae-miR167d, tae-miR172a, tae-miR393, tae-miR396a, and tae-miR444c.1 were considered to be cold responsive in a thermosensitive genic male sterile (TGMS) lines of wheat [34]. In addition, seven miRNAs, miR159, miR164, miR169, miR319, miR398, miR1029, and miR1126 were also identified to be cold-responsive in the seedling of wheat [9].

Despite these efforts, our knowledge about the roles of miRNAs in wheat response to cold stress is still limited. Here, we sequenced small RNAs of cold treated young spikes at the pistil and stamen differentiation stage and controls. A diverse set of wheat small RNAs were identified, some of which were differentially expressed between cold-treated and control samples. MiRNA targets were then validated by degradome sequencing. By comparing with miRNAs reported previously under other stress and developmental conditions, we provided a comprehensive picture of miRNA functions in cold treated wheat inflorescences that may assist cold tolerance improvement in wheat.

Methods

Plant materials

The seeds of common wheat cultivar JM22 were germinated and grown in a 4 °C chamber for vernalization. After 40 days, the seedlings were transplanted to an artificial climate chamber with a relative humidity of 75% and 22 °C/16 °C day/night temperatures with a light intensity of 3000 lx. Approximately 20 days after transplanting, most of the young spikes progressed into the pistil and stamen differentiation stage, and then the seedlings were divided into two groups. One group was moved to a chamber and grew under cold treatment at 0 °C for 48 h, while the other was used as the control sample and remained in the 22 °C/16 °C growth chamber.

RNA isolation and small RNA and degradome library construction

Young spikes were pooled for the control samples or cold treatment samples after cold treatment for 48 h at 0 °C with two replications, respectively, and then were stored in liquid nitrogen for RNA extraction. Total RNA was isolated using the TRIzol reagent according to the manufacturer’s instructions (Invitrogen, Carlsbad, CA, USA). RNA quality was tested using a 2100 Bioanalyzer RNA Nanochip (Agilent, Santa Clara, CA, USA). The RNA concentration was quantified using a NanoDrop ND-1000 Spectrophotometer (Nano-Drop, Wilmington, DE, USA). Two degradome libraries for the control and cold treatment samples were constructed following the method by German et al [35].

Bioinformatics analysis of small RNA sequencing data

Small RNA and degradome reads were generated from Illumina HiSeqTM analysis. Using the Fastx-toolkit pipeline (http://hannonlab.cshl.edu/fastx_toolkit/), the raw sequencing data were pre-processed to remove low-quality reads, reads smaller than 18 nucleotides (nt), trim adaptor sequences, and contamination formed by adaptor-adaptor ligation. Regarding the small RNA sequencing data, clean reads ranging from 18 to 30 nt were aligned to the Chinese Spring genome from the IWGSC (version 2) website and the transcriptome dataset for young spikes of the common wheat downloaded from the GEO dataset as SRX375489 by the bowtie program [7]. The reads were categorized as the exon, intron, and repeat-associated small RNAs with a perfect match to the reference genome. Perfectly matching sequences were identified for further analysis. The reads were annotated as rRNA, tRNA, small nuclear RNA (snRNA), and small nucleolar RNA (snoRNA) by alignments against the Rfam 12.0 database (http://rfam.xfam.org/). The reads were identified as known miRNAs by the BLAST search against the miRBase 21 database with at most two mismatches (http://www.mirbase.org/), and the wheat microRNAs from PNRD (http://structuralbiology.cau.edu.cn/PNRD), PmiRExAt (http://pmirexat.nabi.res.in/), and IWGSC miRNA dataset [36] with the perfect match. After the removal of the known small RNAs and repeat-associated small RNAs, we used the miRdeep-P program with default parameters to predict the novel miRNAs using a plant-specific scoring system and filtering criteria based on a probabilistic model of miRNA biogenesis [37].

Bioinformatics analysis of degradome sequencing data

For the degradome sequencing data, the clean reads were generated using the same method for the small RNA sequencing analysis. Next, the clean reads were used to identify the cleavage sites of the miRNA-target interaction by the program CleaveLand [38]. The transcripts included the CDS, UTR of Chinese Spring (version 2.2) downloaded from the IWGSC website, and the assembled transcriptome from the Zhang et al downloaded from the GEO dataset as SRX375489 [7]. The cleavage sites were classified as the 0, 1, 2, 3, and 4 categories according to the abundance of reads in the slice sites along the whole transcripts with a P-value less than 0.05.

Quantification of miRNA expression level by qRT-PCR

qRT-PCR was performed to determine the validity of miRNA by deep sequencing for expression profile analysis. Total RNA was extracted as described above. According to previous reports [20, 39], the miRNA abundance was detected. Briefly, we first polyadenylated the total RNA (3 μg) including miRNAs, and then used the Mir-XTM miRNA First-Strand Synthesis kit (Clontech, Inc., Terra Bella, USA) to reverse transcribe poly (T) adapters into cDNA. Using the internal reference gene UBQ, we normalized the cDNA products [20]. These products were used as templates for qRT-PCR. qRT-PCR was performed on a LightCycler® 480 (Roche, Switzerland) using LightCycler 480 SYBR Green I (Roche, Switzerland). Along the entire miRNA sequence and adapter sequence provided by miRNA First-Strand synthesis kit (Clontech, Inc., Terra Bella, USA), the miRNA-specific forward primer for each miRNA and the universal reverse primer were designed. The primers were shown in the Additional file 1: Table S1. qPCR was performed according to the following protocol: 95 °C for 2 min, 40 cycles of 95 °C for 5 s and 60 °C for 20 s, followed by a thermal denaturing step to generate the melt curves. All of the reactions were performed in triplicate, including non-template controls. Statistical analysis was performed using the 2-ΔΔct method [40].

Quantification of miRNA targets by qRT-PCR

Total RNA was extracted from samples according to the manufacturer’s instructions of the RNAiso plus kit (TaKaRa, Dalian, China). They were reverse transcribed with poly (T) adapters into cDNA using a PrimeScriptTM RT reagent Kit with gDNA Eraser kit (TaKaRa, Dalian, China). The cDNA products of each tissue were normalized using Actin (AB181991) as the internal reference gene [7]. These products were used as templates for qRT-PCR. qRT-PCR was performed on a LightCycler® 480 (Roche, Switzerland) using LightCycler 480 SYBR Green I (Roche, Switzerland). The qRT-PCR programmes were performed according to the following protocol: 95 °C for 2 min, 40 cycles of 95 °C for 5 s and 56 °C for 20 s, according to a thermal denaturing step to generate the melt curves. All reactions were run in triplicates, including the non-template controls. Statistical analysis was performed using the 2-ΔΔct method [40]. The primers were shown in the Additional file 1: Table S1.

Results

Sample collection and small RNA sequencing

To investigate the response of miRNAs to cold stress, four libraries of small RNAs were constructed from young spikes at the pistil and stamen differentiation stage with and without cold stress, respectively (Fig. 1), and with two biological replications. After the removal of contaminating reads, 20–22 million clean data representing 6.4–8.0 million distinct sequences were obtained (Table 1). The high correlation between the two replicates, with R2 = 0.96 for the control samples and R2 = 0.98 for the cold-stress samples, indicated the high quality of the small RNA libraries (Fig. 2a & b). Next, we mapped these sequences to the genome of Chinese Spring. Approximately 81–82% of the small RNAs could be perfectly mapped to the genome, and there was no bias distribution among the three sets of subgenomes A, B, and D (Table 1 & Additional file 1: Table S2). Most of the small RNAs were located in the intergenic and repeat element regions, while only 3 to 6% small RNAs were mapped to exon and intron regions. Intriguingly, in the exon regions, the small RNAs located in the subgenome B were significantly lower than those in the other two subgenomes, A and D (Additional file 1: Table S2 & Additional file 2: Figure S1), suggesting possible subgenome dominance of small RNA regulation functions. The lengths of the small RNAs ranged from 14 to 43 nucleotides (nts). Most of the small RNAs ranged from 20 to 24 nt and were used for subsequent analysis. For distinct small RNAs, the 24-nt category was the most abundant, representing 71–73% among unique reads and 45–53% for redundant ones, followed by 21-nt small RNAs (Fig. 2c & d). As shown in other plants, 24-nt siRNAs may play a role in maintaining genome integrity and stabilization by heterochromatin formation [41].

Fig. 1
figure 1

Young spikes at the pistil and stamen differentiation stage

Full size image
Table 1 Small RNA mapping information using the Chinese Spring genome as a referencea
Full size table
Fig. 2
figure 2

a-b Correlation of the small RNA libraries of the control samples between the two replicates (a) and cold-stressed samples between the two replicates (b) in young spikes at the differentiation stage of pistil and stamen. c-d Length distribution of distinct (c) and abundant (d) small RNAs mapped on the Chinese Spring genome. Red, coral red, cyan and blue colours represent the control samples with two replicates and cold-stressed libraries with two replicates in young spikes at the differentiation stage of pistil and stamen

Full size image

Identification of miRNAs in wheat young spikes

The functions of many miRNAs are conserved in plants, and these miRNAs play important roles in development and stress responses [8]. MiRNAs were generated across the evolution of plants, the miRNAs with important functions were retained as the plant-conserved miRNAs, and some miRNAs as monocot-specific or eudicot-specific miRNAs came into being after the divergence of monocot and eudicot, and some miRNAs as species-specific miRNAs were generated after the species divergence [42]. To identify the plant-conserved miRNAs in wheat, we aligned the cleaned small RNAs to known plant miRNAs downloaded from the miRBase (release 21.0), PNRD, PmiRExAt, and IWGSC miRNA database using the BLASTN program. In total, 192 known miRNAs from 105 families were obtained. Among them, 25 conserved miRNA families covered almost all of the plant-conserved miRNA families. Twenty-four grass-specific miRNA families were identified, including those in rice, maize, sorghum, and Brachypodium. Recently, using the next-generation sequencing, many miRNAs were identified in the bread wheat, their progenitor T. Urartu and A. tauschii, emmer wheat T. turgidum, and barley [43]. These miRNAs were defined as Triticum-specific miRNAs. Here, 13 Triticum-specific miRNAs, as well as 43 wheat-specific miRNAs, were identified (Additional file 1: Table S3). After the removal of these known miRNAs and Rfam-associated small RNAs, we aligned the remaining small RNAs to the genome for novel miRNA discovery. To date, many wheat-specific miRNAs have been deposited in the miRBase across different development stages and different stress responses. However, based on the special feature of the canonical stem-loop regions of longer RNA precursors, we identified additional nine novel miRNAs in the four small RNA samples that met the previous reported criteria by the miRDeep-P program [37]. Seven of them presented the miRNA* sequences indicating bona fide miRNAs (Additional file 1: Table S4). Five of them were 24 nt in length, and four had 5’ terminal adenine nucleotides. The adenine preferences in the 5’ terminus were processed by the AGO4 protein. Thus, these newly discovered novel miRNAs should increase the number of known miRNA families.

The abundance of miRNAs in wheat young spikes

To predict the potential roles of novel miRNAs, we investigated their expression profiles. We found that 14 miRNAs from nine families had reads more than ten thousands in at least one sample. Most highly abundant miRNAs were conserved miRNAs such as miR168, followed by miR166, miR167, and miR172, which also highly expressed in other plants [34, 44]. Intriguingly, grass-specific miR5062, Triticum-specific miRNAs including tae-miR9672 and miR9863, as well as wheat-specific tae-miR9663 and miR5168, also showed quite high expression levels (Additional file 1: Table S3). There were 24 miRNAs from 16 families with more than one thousand reads in at least one sample. They were four plant conserved miRNAs, four grass-specific miRNAs, eight Triticum-specific miRNAs and four wheat-specific miRNAs. The high expression levels of these miRNAs may suggest their important biology functions. By contrast, some miRNA families, including plant-conserved miRNAs, were expressed at a lower abundance (Additional file 1: Table S3). For the novel miRNAs, miR1002 and miR1003 were intermediately expressed with reads of more than one thousand, while the other miRNAs were lowly expressed with only about one hundred reads (Additional file 1: Table S4). Thus, miRNAs were transcribed at various levels in young spikes during the differentiation of the pistil and stamen in wheat.

Differential expression of miRNAs under cold treatment

As have been shown in previous studies, miRNAs are often associated with cold response in many plants such as Arabidopsis, Populus, Brachypodium, and Medicago [45]. Here, we performed differential expression analysis of miRNAs in young spikes with and without cold treatment. After normalizing reads of each miRNA as ‘transcripts per million’ (TPM) for each sample, we characterized 39 differentially expressed miRNAs (FDR < 0.001) from 28 miRNA families using the DEGseq program (Fig. 3a) [46]. Twelve of these differentially expressed miRNA families were conserved in plants, two were grass-specific miRNAs, six Triticum-specific, three wheat-specific miRNAs, and five were newly identified novel miRNAs. Among them, 23 were down-regulated, and 16 were up-regulated with differential expression fold changes ranging from 1.2 to 30 times (Fig. 3a).

Fig. 3
figure 3

Differential expression profiles of miRNAs. a Fold change of differentially expressed miRNAs between the control and cold-stressed libraries in young spikes at the pistil and stamen differentiation stage. Red and blue colours represent up-regulated and down-regulated expression. b Correlation of the differentially expressed miRNAs validated by qRT-PCR compared to the Illumina sequencing. c Venn diagram of the microRNAs response to cold stress among wheat, Brachypodium, Medicago, Populus, and Arabidopsis

Full size image

To confirm the sequencing data, we selected 18 differentially expressed miRNAs for qRT-PCR detection. The expression means of cold responsive miRNAs as detected as by qRT-PCR was very similar to the fold change of the sequencing reads with a high correlation coefficients (R2 = 0.832; P-value = 5.344e-08) (Fig. 3b). For example, miR156, miR169, miR171a, miR444c, miR827, and miR5048a were identified by sequencing as down-regulated. Such patterns were confirmed by qRT-PCR. Similarly, miR159a, miR160b, miR165a, miR171c, miR319a, miR2118a, miR1130c, and the novel miRNAs miR10008, miR10004, and miR10009 were also validated as being up-regulated by qRT-PCR, consistent with the results of Illumina sequencing (Additional file 1: Table S5). These data demonstrated the effectiveness and accuracy of high-throughput sequencing.

Among the 23 differentially expressed known miRNA families, nine were overlapped with cold-stress miRNAs from other species, such as Arabidopsis [12,13,14], Medicago [45], Populus [15], and Brachypodium [16] (Fig. 3c). MiR169, miR172, miR156, miR319, miR159, and miR396 showed the cold-stress response in at least three species, and miR160, miR165, miR167, and miR171 overlapped with Populus, Arabidopsis, and Medicago (Additional file 1: Table S6). However, some miRNAs were found to be species-specific (Additional file 1: Table S6), suggesting divergent regulatory pathways in wheat.

MiRNA target prediction

To understand miRNA functions, we performed target prediction of miRNAs by the TargetFinder program [47] and targets with scores less than four were further analysed. We found that targets of many known miRNAs were mostly associated with development and stress (Additional file 1: Table S7). For example, miR156 targets an SPB gene while miR172 targets an AP2-like ethylene-responsive transcription factor (AP2), both of which mediate the transition of the vegetative stage to reproductive stage [48,49,50]. Other cold responsive miRNAs, MiR319 and miR444, targeted TCP and MADS-box genes were involved flower development [51, 52]. Additionally, miR390, miR5071, miR2118, miR9863 and miR7757 were predicted to target the Leucine-rich Repeat Receptor-like protein kinase family (LRR) that are involved in disease resistance [53, 54], as well as miR1120, miR1127, miR1130, miR1137, miR1439, miR5049, miR5062, and miR9673 that regulated the WD domain gene for flower development and the immune system [55,56,57] (Additional file 1: Table S7).

For novel miRNAs, a total of eight miRNA-target interaction were predicted. As shown in Additional file 1: Table S8, the predicted targets include gene encoding zipper proteins, binding proteins, protein kinases, enzymes and transporters that are potentially involved in multiple cellular processes. Some target genes have been shown to be also involved in stress responses. For instances, the homeobox-leucine zipper protein family regulated by tae-miR1002 and tae-miR1007 has been shown to be involved in the water stress and abscisic acid treatment in Arabidopsis [58], while the CASP-like protein, which is targeted by tae-miR1005, mediates the environmental stress response in plants except for the pathogen-induced membrane protein gene in bacterial disease resistance and oomycete disease susceptibility in pepper [59, 60]. In addition, the LRR, which participates in plant defence, regulated of the developmental process and the sensing or transduction of hormone signals in plants [61]. These genes were predicted to be regulated by tae-miR1004 and tae-miR1005. Also, the sugar transporter gene that was targeted by the tae-miR1007 plays a central role in pathogen resistance in rice and cassava [62].

To confirm the coordinated expressions of target genes with those of cold responsive miRNAs, 12 predicted target genes were selected to perform qRT-PCR assays. The results showed that six target genes displayed coordinated expression patterns with corresponding miRNAs while other six showed expression patterns opposite to the expected (Additional file 1: Table S5), which could be caused by more complicated regulation of the target genes such as by other transcription factors.

MiRNA target validation by degradome sequencing

To validate that miRNAs found there can indeed cleave predicted targets, we performed degradome sequencing using the same wheat inflorescence samples. After removing contaminating reads, 18,427,417 and 20,253,810 clean reads were obtained in the two degradome libraries, respectively, one control and one cold treated. 99.3% and 99.4% of them were 20 and 21 nts in length, respectively (Additional file 2: Figure S2a). A total of 23,253,400 reads (60.12%) were overlapped between the two libraries. Among them, 15,879,947 (86.17%) and 17,063,879 (84.25%) reads were aligned to the CDS sequences comprising the Chinese Spring and transcriptome datasets, respectively [7]. Using the CleaveLand program, 109 known miRNAs from 62 families and 112 from 65 families were confirmed to cleave the targets with a P-value less than 0.05 using the degradome sequences in the control and cold-stressed libraries, respectively (Additional file 1: Table S9). One hundred and twenty two miRNAs from 70 miRNA families were verified to cleave the targets using the degradome sequences in at least one library. Totally, 99 miRNAs from 57 families were demonstrated to interact with targets in the two libraries. In total, 457 and 484 transcript cleaved sites were confirmed in degradome libraries, and more than half of them were either defined as the ‘0 Category’ with the maximum tags or only one peak at the cleavage site on the transcript (Additional file 2: Figure S2b). Examples for the ‘0 Categories’ are shown in Fig. 4a-d. In addition, 312 miRNA-transcript pairs were shared between the two libraries, and 145 and 172 were characterized specifically in the control and cold stress libraries, respectively. Among transcripts detected by degradome sequencing, 344 transcripts were predicted as targets with a score less than four using the TargetFinder program. For the novel miRNAs, all nine miRNAs were identified to regulate the 32 transcripts in the two libraries. There were six miRNAs with 20 interaction transcripts overlapped between the two libraries, in which three miRNAs with three interaction transcripts specifically were identified in the control libraries and four with nine targeting transcripts were specifically identified in the cold stress libraries (Additional file 1: Table S10). Moreover, as shown in Table 2, a significant number of transcripts detected by degradome sequencing have been reported previously such as miR156-SBP, miR172-AP2, miR160-ARF (Auxin Response Factor), miR169-NFYA (Nuclear transcription Factor Y subunit A), miR319-TCP, and MiR9863-NBS, indicating conserved miRNA-target interaction in plants.

Fig. 4
figure 4

Target plots (T-plots) of miR156 (a), miR159 (b), miR169 (c), and miR5028 (d) characterized by degradome sequencing. The abundance of the signature tags was plotted along the indicated transcript. The red dots indicate the predicted cleavage sites on the x axis; the black line indicates the signatures produced by miRNA-directed cleavage

Full size image
Table 2 MiRNA targets validated by degradome sequencinga
Full size table

Gene enrichment (GO) analysis for target genes of differentially expressed miRNAs

To better understand the function of these differentially expressed miRNAs, enrichment analysis for 685 target genes of significantly differentially expressed miRNAs was analysed using the BINGO program [63]. The GO term annotation of wheat genes was performed according to the Chinese Spring genome release available from the IWGSC website. The results showed that 30, 19 and 5 GO terms were significantly enriched (FDR < 0.05) in the biological process, molecular function, and cellular component categories, respectively. For the biological process, most of the enriched GO terms were classified as the biological regulation such as regulation of transcription, RNA metabolic process, and macromolecule biosynthesis process, suggesting that changes in gene transcription regulation occur in response to the cold stimulus in wheat. For the response stimulus pathway, GO terms for the response to the endogenous stimulus, organic substance, and hormone stimulus were enriched. The analysis also showed that plants may have induced the salvage pathway in response to cold stress via the enriched GO terms such nucleotide salvage, purine nucleotide salvage and IMP salvage (Fig. 5a & Additional file 1: Table S11). Although most of differentially expressed miRNAs are involved in stress responses, some of them may also affect ion transport, such as the presence of the enriched GO term sulphate transport. Regarding the molecular function category, the transcript regulation activity pathway was enriched, such as the DNA binding, ADP binding, and transcript factor activity GO terms. The GO terms related to the sulphate transmembrane transport activity were both enriched in molecular function and biological process (Fig. 5b & Additional file 1: Table S11). For the cellular component category, the membrane-bounded organelle, intracellular membrane-bounded organelle, and nucleus GO terms were enriched, suggesting that the membrane component was affected under cold stress, leading to the change in membrane permeability that may influence ion transport. Interestingly, most response genes were found to be located in the nucleus (Fig. 5c & Additional file 1: Table S11). Overall, our work showed that differentially expressed miRNAs may respond to cold stress by regulating their target genes.

Fig. 5
figure 5

GO enrichment analysis for the targets of significantly differentially expressed miRNAs. Cytoscape views of enriched GO terms were shown in biological process (a), molecular function (b), and cellular component (c) for the targets of differentially expressed miRNAs by cold stress. The coloured nodes and white nodes represent the GO terms enriched and non-enriched in the target genes, respectively

Full size image

Discussion

MiRNAs demonstrate important regulatory functions in adaptive responses to biotic and abiotic stresses in plants [9, 11, 12, 43]. Budak and colleagues had summarized the stress-associated miRNAs and their related regulatory mechanisms in cereals and Triticum [9,10,11, 43]. MiRNAs as the post-transcript regulators, had been detected as multi-responsive to multi-environmental conditions individually and/or together with their various miRNA partners [9, 11]. Some of miRNA behaved the regulation roles with the great potential for the stress responses such as drought stress in the wild emmer wheat [64]. In Triticum, several studies have shown miRNAs responding to abiotic stresses such as salt, phosphorous, drought, and heat stress, and biotic stresses such as powdery mildew and stripe rust infection [27, 28, 31, 32, 43, 64, 65]. Here, we compared cold responding miRNAs with those responsive to various other stresses and developmental conditions and present an interactive analysis of miRNAs in multiple regulatory pathways, which resulted from the high-through sequencing of small RNAs (Fig. 6).

Fig. 6
figure 6

The interactive pathways of stress-responsive miRNAs and their targets. The cold-stressed microRNAs include information summarized from Xin et al (chartreuse oval) and this study (forest green oval). The forest green link indicated the unique regulation between the miRNAs and targets in response to cold stress only in our study, and the dark blue link indicated the target genes involved in the flower development pathway

Full size image

MiRNAs are extensively involved in cold responses in plants

The sensitive response of miRNA expression levels during cold treatment indicates that miRNAs are important regulatory nodes for plant cold response. Although different species confer different set of miRNAs responding to cold, there a set of core miRNAs that are shared by most species. In wheat, 39 miRNAs from 28 families were differentially expressed under cold treatment. Among them, 12 (31%) were conserved in plants, i.e., conserved in Arabidopsis, populus, Brachypodium, and Medicago [41]. The remaining appears to be temperate grass-specific. It is not clear if plants growing in cold weather develop more complicated cold response mechanisms. Judging from the functions of their targets, we found that genes regulated by cold are often involved in other stress response such as heat, drought, and salt stress. Many conserved miRNAs are also closely involved in plant development. Despite this, some miRNAs seem to be functioning mainly in stress responses. For instance, MiR394 is a cold response miRNA in wheat as shown in this study. MiR394 targets the LEAF CURLING RESPONSIVENESS (LCR) gene which is involved in salt and drought stress responses, and Fe-responsive in Arabidopsis [9, 66]. Similar, miR394 was found to respond to drought stress in wheat as well [29], while grass-specific miRNA miR5169 targets YELLOW STRIPE LIKE 6 (YSL6) which together with YSL4 controls iron release from chloroplasts [67]. The differential expression of miR5169 under cold stress may affect the iron transporter as suggested by Gene Ontology analysis. In addition, miR169 was identified to respond to cold stress by the regulation of NFYA in Arabidopsis [12, 13], Populus [15], and Brachypodium [16]. MiR169 as a big miRNA family in plants were detected with various expression in different tissues [9, 42], which also response to drought, salinity, heat and powdery mildew resistance [11, 43]. However, only two miRNA members miR169g and miR169n(o) with the drought responsive element and ABA-responsive element were showed to the response to the drought and salinity stresses through ABA-dependent and ABA-independent stress-inducible gene expression pathways [11]. Therefore, the post-transcriptional regulatory mode of miRNAs to their targets probably provides more efficient regulation on multiple genes that function in different biological processes.

Co-induction of miRNAs by cold and other stresses

Although cold response miRNAs seem to be prevalent in plants, they are also responsive to other stress conditions. Up to date, most cold-responsive miRNAs also respond to additional abiotic and biotic stresses [9]. In wheat, for example, miR160 and miR167 are involved in the response to multiple biotic and abiotic stresses [27, 28, 31, 32, 34, 65]. Both miR160 and miR167 targets auxin-response factors (ARFs), but they displayed opposite expression patterns. MiR160 was up-regulated and miR167 down-regulated in this study, which is similar to their expression profiles in Medicago in response to cold and freezing stresses [45]. In addition, they were both heat-responsive miRNAs [9, 11]. In wheat and Aegilops tauschii, miR167 exhibited the reverse responses to drought stress. MiR167 is up regulated after drought stress in Aegilops tauschii, while down regulated in wheat [68]. For the nutrient stresses, miR160 had a role in the N homeostasis [9], and miR167 was observed to responsive to Pi deprivation in wheat [9]. In the biotic stress, miR167 were also involved in the defence responses to the infection of Puccinia striiformis f. sp. Tritici in wheat [43], while miR160 were responsive to the Magnaporthe oryzae infection in rice [11]. Thus, miRNA may regulate different auxin regulation factors involved in different stress response pathways. In wheat, miR159 is inducible by abiotic stresses, including heat, cold in this study, salt, drought, dehydration and nutrient stress such as N deficiency, Cd Stress, and phosphorus stresses, as well as biotic stress such as powdery mildew infection and Puccinia striiformis stress [9,31,, 11, 27, 28, 30–32, 34, 43, 65, 69], whereas it showed transient or mild regulation under cold stress in Arabidopsis [13]. Differential expression of the MYB-miR159 module was detected as Aluminum toxicity-responsive in soybean and Fe-responsive in Arabidopsis [9]. Since over-expression of miR159 leads to male sterility, cold stress may be one important factor to inhibit anther development via miRNAs [70], which could be as promising candidates to be used in development of higher yielding in wheat [9, 11].

The correlation of cold responsive miRNAs with biotic response was exemplified by miR9672 which targets BTB/POB, a disease-resistance domain protein that is involved in the immune system responses, including those to biotic stress such as powdery mildew infection [32, 55]. Similarly, miR9863 targets disease resistance with the classical domain NB-ARC (Nucleotide-Binding, Apaf-1, R-proteins, and CED-4) that triggers isolate-specific immune responses against the powdery mildew fungus in Triticum [54]. The cross response of miRNA to biotic and abiotic stress was also observed. MiR827, for example, responds to powdery mildew infection, heat, dehydration and cold stress via targeting the SPX family [30, 32]. The SPX family is associated with sensitivity to cold stress in rice seedlings [71] and semi-male sterility, resulting in the reduction of grain yield in rice [72]. These versatile regulatory capabilities demonstrate miRNAs as plural potent regulator for plant stress responses.

Development related miRNAs are sensitive to cold

Among cold responsive miRNAs, miR156 and miR172 are most important developing related miRNAs. In Arabidopsis, the interaction of miR156 and miR172 controls developmental timing by the regulation of SPB and AP2 in Arabidopsis [48]. Furthermore, miR156 acts in concert to secure male fertility by targeting the SPB transcript factors in Arabidopsis [73]. Thus, overexpressed miRNA156 will increase the biomass in cereals such as the energy crop switchgrass. In wheat, miR156 responds to heat, drought, dehydration, and Cd stresses, and fungus stresses in wheat [28, 30, 32, 64], while Tae-miR172 expresses in the tapetum and microsporocytes at the anther development stage [34]. The response of miR156 and miR172 to cold, heat, and drought are conserved in Arabidopsis, rice, and Brachypodium [13, 74]. MiR171 responds to powdery mildew infection [32], dehydration [30], drought [29, 64], and salinity stresses [9] in wheat, which is also shown responding to cold stress in this study. It targeted the GRAS (GAI, RGA, and SCR) family. In GRAS, the member of GS6, plays an important role in the control of grain development in rice [75]. MiR319, which responds to both cold stress and drought stress in wheat [29], was also characterized in Arabidopsis and rice for cold stress [12,13,14, 76]. It targets members of the TCP gene family, which is specifically expressed in the flower organs [77], which will contribute to the higher yield in crops under stress conditions. MiR444 regulates the MADs-box genes, especially in monocot species, and has been demonstrated to regulate flower development [78]. In wheat, it was observed to respond to powdery mildew infection [32], dehydration [30], drought [29] and cold stress in our study. MADS28 was demonstrated to specifically regulate the floral organ number, filament length and pollen release in soybean [79]. Other less conserved miRNAs that respond to cold treatment are also involved in developmental processes. For instance, miR2118 triggers the generation of phased siRNAs, which are involved in panicle development [80, 81], while miR5062 and miR9674 targets the PPR protein gene, which regulates inflorescence branch development in rice [82] and can restore fertility in the cytoplasmic male-sterile line in Brassica napus [83]. MiR5048 targets the U-box. The U-Box/ARM was discovered to regulate flowering time in Arabidopsis [84]. Therefore, there are a number of targets of cold responsive miRNAs that also work in plant development, especially flower development, male sterility, and floret and spike development. Such an observation is interesting and may be the result of the tissue and treatment we sued in this study because the real cold tolerance genes such as CBF, FR2, COR, and LEA seem not to be the target of any miRNAs. In other words, the dramatic change in miRNA expression levels in wheat young inflorescences will cause sensitive reaction of their targets which will regulate proper flower organ development and hence the correct timing in the wild. Such mechanisms may be explored to manipulate genetically to enhance cold tolerance in wheat contributing to yield increase.

Conclusions

This work provides the first small RNA expression profiles of young spikes at the differentiation stage of the pistil and stamen in common wheat. Our small RNA sequencing data and degradome dataset revealed regulatory roles of miRNAs during wheat cold response. The cross-response of miRNAs to multi-biotic and abiotic stresses indicated that wheat has evolved sophisticated miRNA-mediated pathways to cope with ever changing environments. Further study of these mechanisms should help our understanding of plant cold response and further secure grain yield in wheat.

Abbreviations

AP2:

AP2-like ethylene-responsive transcription factor

ARF:

Auxin response factor

GO:

Gene ontology

GRAS:

GAINSENSITIVE (GAI), REPRESSOR OF gal-3 (RGA), and SCARECROW (SCR)

LRR:

Leucine-rich repeat receptor-like protein kinase family

MADS:

MCM1, AG, DEFA and SRF

MiRNA:

MicroRNAs

NB-ARC:

Nucleotide-binding, Apaf-1, R-proteins, and CED-4

NFYA:

Nuclear transcription factor Y subunit A

PPR:

PentatricoPeptide repeat

qRT-PCR:

Quantitative real time polymerase chain reaction

SPB:

Squamosa promoter binding like protein

SPX:

SYG1, Pho81 and XPR1

TCP:

TEOSINTE BRANCHED, Cycloidea and PCF

TGMS:

Thermosensitive genic male sterile

TPM:

Transcripts per million

YSL6:

YELLOW STRIPE LIKE 6

References

  1. Paulsen GM, Heyne EG. Grain production of winter wheat after spring freeze injury. Agron J. 1983;75(4):705–7.

    Article  Google Scholar 

  2. Zhu J, Pearce S, Burke A, See DR, Skinner DZ, Dubcovsky J, Garland-Campbell K. Copy number and haplotype variation at the VRN-A1 and central FR-A2 loci are associated with frost tolerance in hexaploid wheat. Theor Appl Genet. 2014;127(5):1183–97.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  3. Knox AK, Dhillon T, Cheng H, Tondelli A, Pecchioni N, Stockinger EJ. CBF gene copy number variation at Frost Resistance-2 is associated with levels of freezing tolerance in temperate-climate cereals. Theor Appl Genet. 2010;121(1):21–35.

    Article  PubMed  Google Scholar 

  4. Herman EM, Rotter K, Premakumar R, Elwinger G, Bae H, Ehler-King L, Chen S, Livingston 3rd DP. Additional freeze hardiness in wheat acquired by exposure to -3 degreesC is associated with extensive physiological, morphological, and molecular changes. J Exp Bot. 2006;57(14):3601–18.

    Article  CAS  PubMed  Google Scholar 

  5. Koo BC, Bushman BS, Mott IW. Transcripts associated with non-acclimated freezing response in two barley cultivars. Plant Genome. 2008;1(1):21–32.

    Article  CAS  Google Scholar 

  6. Kang G, Li G, Yang W, Han Q, Ma H, Wang Y, Ren J, Zhu Y, Guo T. Transcriptional profile of the spring freeze response in the leaves of bread wheat (Triticum aestivum L.). Acta Physiol Plant. 2013;35(2):575–87.

    Article  CAS  Google Scholar 

  7. Zhang S, Song G, Gao J, Li Y, Guo D, Fan Q, Sui X, Chu X, Huang C, Liu J, et al. Transcriptome characterization and differential expression analysis of cold-responsive genes in young spikes of common wheat. J Biotechnol. 2014;189:48–57.

    Article  CAS  PubMed  Google Scholar 

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

    Article  CAS  PubMed  Google Scholar 

  9. Alptekin B, Langridge P, Budak H. Abiotic stress miRNomes in the Triticeae. Funct Integr Genomics. 2016.

  10. Budak H, Kantar M. Harnessing NGS and big data optimally: comparison of miRNA prediction from assembled versus non-assembled sequencing data—the case of the grass Aegilops tauschii complex genome. OMICS. 2015;19(7):407–15.

    Article  CAS  PubMed  Google Scholar 

  11. Budak H, Kantar M, Bulut R, Akpinar BA. Stress responsive miRNAs and isomiRs in cereals. Plant Sci. 2015;235:1–13.

    Article  CAS  PubMed  Google Scholar 

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

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

  14. Zhou X, Wang G, Sutoh K, Zhu JK, Zhang W. Identification of cold-inducible microRNAs in plants by transcriptome analysis. Biochim Biophys Acta. 2008;1779(11):780–8.

    Article  CAS  PubMed  Google Scholar 

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

    Article  CAS  PubMed  Google Scholar 

  16. Zhang J, Xu Y, Huan 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 

  17. Yao Y, Guo G, Ni Z, Sunkar R, Du J, Zhu JK, Sun Q. Cloning and characterization of microRNAs from wheat (Triticum aestivum L.). Genome Biol. 2007;8(6):R96.

    Article  PubMed  PubMed Central  Google Scholar 

  18. Wei B, Cai T, Zhang R, Li A, Huo N, Li S, Gu YQ, Vogel J, Jia J, Qi Y, et al. Novel microRNAs uncovered by deep sequencing of small RNA transcriptomes in bread wheat (Triticum aestivum L.) and Brachypodium distachyon (L.) Beauv. Funct Integr Genomics. 2009;9(4):499–511.

    Article  CAS  PubMed  Google Scholar 

  19. Sun F, Guo G, Du J, Guo W, Peng H, Ni Z, Sun Q, Yao Y. Whole-genome discovery of miRNAs and their targets in wheat (Triticum aestivum L.). BMC Plant Biol. 2014;14:142.

    Article  PubMed  PubMed Central  Google Scholar 

  20. Han R, Jian C, Lv J, Yan Y, Chi Q, Li Z, Wang Q, Zhang J, Liu X, Zhao H. Identification and characterization of microRNAs in the flag leaf and developing seed of wheat (Triticum aestivum L.). BMC Genomics. 2014;15:289.

    Article  PubMed  PubMed Central  Google Scholar 

  21. Alptekin B, Budak H. Wheat miRNA ancestors: evident by transcriptome analysis of A, B, and D genome donors. Funct Integr Genomics. 2016.

  22. Kantar M, Akpinar BA, Valarik M, Lucas SJ, Dolezel J, Hernandez P, Budak H. Subgenomic analysis of microRNAs in polyploid wheat. Funct Integr Genomics. 2012;12(3):465–79.

    Article  CAS  PubMed  Google Scholar 

  23. Lucas SJ, Budak H. Sorting the wheat from the chaff: identifying miRNAs in genomic survey sequences of Triticum aestivum chromosome 1AL. PLoS One. 2012;7(7):e40859.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  24. International Wheat Genome Sequencing Consortium. A chromosome-based draft sequence of the hexaploid bread wheat (Triticum aestivum) genome. Science. 2014;345(6194):1251788.

    Article  Google Scholar 

  25. Ling HQ, Zhao S, Liu D, Wang J, Sun H, Zhang C, Fan H, Li D, Dong L, Tao Y, et al. Draft genome of the wheat A-genome progenitor Triticum urartu. Nature. 2013;496(7443):87–90.

    Article  CAS  PubMed  Google Scholar 

  26. Jia J, Zhao S, Kong X, Li Y, Zhao G, He W, Appels R, Pfeifer M, Tao Y, Zhang X, et al. Aegilops tauschii draft genome sequence reveals a gene repertoire for wheat adaptation. Nature. 2013;496(7443):91–5.

    Article  CAS  PubMed  Google Scholar 

  27. Lu W, Li J, Liu F, Gu J, Guo C, Xu L, Zhang H, Xiao K. Expression pattern of wheat miRNAs under salinity stress and prediction of salt-inducible miRNAs targets. Front Agric China. 2011;5(4):433–22.

    Article  Google Scholar 

  28. Pandey R, Joshi G, Bhardwaj AR, Agarwal M, Katiyar-Agarwal S. A comprehensive genome-wide study on tissue-specific and abiotic stress-specific miRNAs in Triticum aestivum. PLoS One. 2014;9(4):e95800.

    Article  PubMed  PubMed Central  Google Scholar 

  29. Akpinar BA, Kantar M, Budak H. Root precursors of microRNAs in wild emmer and modern wheats show major differences in response to drought stress. Funct Integr Genomics. 2015;15(5):587–98.

    Article  CAS  PubMed  Google Scholar 

  30. Ma X, Xin Z, Wang Z, Yang Q, Guo S, Guo X, Cao L, Lin T. Identification and comparative analysis of differentially expressed miRNAs in leaves of two wheat (Triticum aestivum L.) genotypes during dehydration stress. BMC Plant Biol. 2015;15:21.

    Article  PubMed  PubMed Central  Google Scholar 

  31. Zhao X, Liu X, Guo C, Gu J, Xiao K. Identification and characterization of microRNAs from wheat (Triticum aestivum L.) under phosphorus deprivation. J Plant Biochem Biotechnol. 2013;22(1):113–23.

    Article  CAS  Google Scholar 

  32. Xin M, Wang Y, Yao Y, Xie C, Peng H, Ni Z, Sun Q. Diverse set of microRNAs are responsive to powdery mildew infection and heat stress in wheat (Triticum aestivum L.). BMC Plant Biol. 2010;10:123.

    Article  PubMed  PubMed Central  Google Scholar 

  33. Feng H, Wang T, Feng C, Zhang Q, Zhang X, Huang L, Wang X, Kang Z. Identification of microRNAs and their corresponding targets involved in the susceptibility interaction of wheat response to Puccinia striiformis f. sp. tritici. Physiol Plant. 2016;157(1):95–107.

    Article  CAS  PubMed  Google Scholar 

  34. Tang Z, Zhang L, Xu C, Yuan S, Zhang F, Zheng Y, Zhao C. Uncovering small RNA-mediated responses to cold stress in a wheat thermosensitive genic male-sterile line by deep sequencing. Plant Physiol. 2012;159(2):721–38.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  35. German MA, Luo S, Schroth G, Meyers BC, Green PJ. Construction of Parallel Analysis of RNA Ends (PARE) libraries for the study of cleaved miRNA targets and the RNA degradome. Nat Protoc. 2009;4(3):356–62.

    Article  CAS  PubMed  Google Scholar 

  36. Consortium IWGS. A chromosome-based draft sequence of the hexaploid bread wheat (Triticum aestivum) genome. Science. 2014;345(6194):1251788.

    Article  Google Scholar 

  37. Yang X, Li L. miRDeep-P: a computational tool for analyzing the microRNA transcriptome in plants. Bioinformatics. 2011;27(18):2614–5.

    CAS  PubMed  Google Scholar 

  38. Addo-Quaye C, Miller W, Axtell MJ. CleaveLand: a pipeline for using degradome data to find cleaved small RNA targets. Bioinformatics. 2009;25(1):130–1.

    Article  CAS  PubMed  Google Scholar 

  39. Shi R, Chiang VL. Facile means for quantifying microRNA expression by real-time PCR. Biotechniques. 2005;39(4):519–25.

    Article  CAS  PubMed  Google Scholar 

  40. 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(4):402–8.

    Article  CAS  PubMed  Google Scholar 

  41. Li A, Liu D, Wu J, Zhao X, Hao M, Geng S, Yan J, Jiang X, Zhang L, Wu J, et al. mRNA and small RNA transcriptomes reveal insights into dynamic homoeolog regulation of allopolyploid heterosis in nascent hexaploid wheat. Plant Cell. 2014;26(5):1878–900.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  42. Budak H, Akpinar BA. Plant miRNAs: biogenesis, organization and origins. Funct Integr Genomics. 2015;15(5):523–31.

    Article  CAS  PubMed  Google Scholar 

  43. Budak H, Khan Z, Kantar M. History and current status of wheat miRNAs using next-generation sequencing and their roles in development and stress. Brief Funct Genomics. 2015;14(3):189–98.

    Article  PubMed  Google Scholar 

  44. Yang X, Wang L, Yuan D, Lindsey K, Zhang X. Small RNA and degradome sequencing reveal complex miRNA regulation during cotton somatic embryogenesis. J Exp Bot. 2013;64(6):1521–36.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  45. Shu Y, Liu Y, Li W, Song L, Zhang J, Guo CP. Genome-wide investigation of MicroRNAs and their targets in response to freezing stress in Medicago sativa L, based on high-throughput sequencing. G3 (Bethesda). 2016;6(3):755–63.

    Article  Google Scholar 

  46. Wang L, Feng Z, Wang X, Wang X, Zhang X. DEGseq: an R package for identifying differentially expressed genes from RNA-seq data. Bioinformatics. 2010;26(1):136–8.

    Article  PubMed  Google Scholar 

  47. Bo X, Wang S. TargetFinder: a software for antisense oligonucleotide target site selection based on MAST and secondary structures of target mRNA. Bioinformatics. 2005;21(8):1401–2.

    Article  CAS  PubMed  Google Scholar 

  48. Wu G, Park MY, Conway SR, Wang JW, Weigel D, Poethig RS. The sequential action of miR156 and miR172 regulates developmental timing in Arabidopsis. Cell. 2009;138(4):750–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  49. Zhu QH, Helliwell CA. Regulation of flowering time and floral patterning by miR172. J Exp Bot. 2011;62(2):487–95.

    Article  CAS  PubMed  Google Scholar 

  50. Wang JW, Czech B, Weigel D. miR156-regulated SPL transcription factors define an endogenous flowering pathway in Arabidopsis thaliana. Cell. 2009;138(4):738–49.

    Article  CAS  PubMed  Google Scholar 

  51. Schommer C, Debernardi JM, Bresso EG, Rodriguez RE, Palatnik JF. Repression of cell proliferation by miR319-regulated TCP4. Mol Plant. 2014;7(10):1533–44.

    Article  CAS  PubMed  Google Scholar 

  52. Nag A, King S, Jack T. miR319a targeting of TCP4 is critical for petal growth and development in Arabidopsis. Proc Natl Acad Sci U S A. 2009;106(52):22534–9.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  53. Zhai J, Jeong DH, De Paoli E, Park S, Rosen BD, Li Y, Gonzalez AJ, Yan Z, Kitto SL, Grusak MA, et al. MicroRNAs as master regulators of the plant NB-LRR defense gene family via the production of phased, trans-acting siRNAs. Genes Dev. 2011;25(23):2540–53.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  54. Liu J, Cheng X, Liu D, Xu W, Wise R, Shen QH. The miR9863 family regulates distinct Mla alleles in barley to attenuate NLR receptor-triggered disease resistance and cell-death signaling. PLoS Genet. 2014;10(12):e1004755.

    Article  PubMed  PubMed Central  Google Scholar 

  55. Islas-Flores T, Rahman A, Ullah H, Villanueva MA. The receptor for activated C kinase in plant signaling: tale of a promiscuous little molecule. Front Plant Sci. 2015;6:1090.

    Article  PubMed  PubMed Central  Google Scholar 

  56. Miller JC, Chezem WR, Clay NK. Ternary WD40 repeat-containing protein complexes: evolution, composition and roles in plant immunity. Front Plant Sci. 2015;6:1108.

    PubMed  Google Scholar 

  57. Kumar S, Jordan MC, Datla R, Cloutier S. The LuWD40-1 gene encoding WD repeat protein regulates growth and pollen viability in flax (Linum Usitatissimum L.). PLoS One. 2013;8(7):e69124.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  58. Lee YH, Chun JY. A new homeodomain-leucine zipper gene from Arabidopsis thaliana induced by water stress and abscisic acid treatment. Plant Mol Biol. 1998;37(2):377–84.

    Article  CAS  PubMed  Google Scholar 

  59. Hong JK, Hwang BK. The promoter of the pepper pathogen-induced membrane protein gene CaPIMP1 mediates environmental stress responses in plants. Planta. 2009;229(2):249–59.

    Article  CAS  PubMed  Google Scholar 

  60. Hong JK, Choi DS, Kim SH, Yi SY, Kim YJ, Hwang BK. Distinct roles of the pepper pathogen-induced membrane protein gene CaPIMP1 in bacterial disease resistance and oomycete disease susceptibility. Planta. 2008;228(3):485–97.

    Article  CAS  PubMed  Google Scholar 

  61. van der Knaap E, Song WY, Ruan DL, Sauter M, Ronald PC, Kende H. Expression of a gibberellin-induced leucine-rich repeat receptor-like protein kinase in deepwater rice and its interaction with kinase-associated protein phosphatase. Plant Physiol. 1999;120(2):559–70.

    Article  PubMed  PubMed Central  Google Scholar 

  62. Eom JS, Chen LQ, Sosso D, Julius BT, Lin IW, Qu XQ, Braun DM, Frommer WB. SWEETs, transporters for intracellular and intercellular sugar translocation. Curr Opin Plant Biol. 2015;25:53–62.

    Article  CAS  PubMed  Google Scholar 

  63. Maere S, Heymans K, Kuiper M. BiNGO: a Cytoscape plugin to assess overrepresentation of gene ontology categories in biological networks. Bioinformatics. 2005;21(16):3448–9.

    Article  CAS  PubMed  Google Scholar 

  64. Kantar M, Lucas SJ, Budak H. miRNA expression patterns of Triticum dicoccoides in response to shock drought stress. Planta. 2011;233(3):471–84.

    Article  CAS  PubMed  Google Scholar 

  65. Kumar RR, Pathak H, Sharma SK, Kala YK, Nirjal MK, Singh GP, Goswami S, Rai RD. Novel and conserved heat-responsive microRNAs in wheat (Triticum aestivum L.). Funct Integr Genomics. 2015;15(3):323–48.

    Article  CAS  PubMed  Google Scholar 

  66. Song JB, Gao S, Sun D, Li H, Shu XX, Yang ZM. miR394 and LCR are involved in Arabidopsis salt and drought stress responses in an abscisic acid-dependent manner. BMC Plant Biol. 2013;13:210.

    Article  PubMed  PubMed Central  Google Scholar 

  67. Divol F, Couch D, Conejero G, Roschzttardtz H, Mari S, Curie C. The Arabidopsis YELLOW STRIPE LIKE4 and 6 transporters control iron release from the chloroplast. Plant Cell. 2013;25(3):1040–55.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  68. Akpinar BA, Budak H. Dissecting miRNAs in wheat D genome progenitor, Aegilops tauschii. Front Plant Sci. 2016;7:606.

    Article  PubMed  PubMed Central  Google Scholar 

  69. Zhu JK. Salt and drought stress signal transduction in plants. Annu Rev Plant Biol. 2002;53:247–73.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  70. 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(3):705–21.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  71. Wang C, Wei Q, Zhang K, Wang L, Liu F, Zhao L, Tan Y, Di C, Yan H, Yu J, et al. Down-regulation of OsSPX1 causes high sensitivity to cold and oxidative stresses in rice seedlings. PLoS One. 2013;8(12):e81849.

    Article  PubMed  PubMed Central  Google Scholar 

  72. Zhang K, Song Q, Wei Q, Wang C, Zhang L, Xu W, Su Z. Down-regulation of OsSPX1 caused semi-male sterility, resulting in reduction of grain yield in rice. Plant Biotechnol J. 2016.

  73. Xing S, Salinas M, Hohmann 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(12):3935–50.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  74. Zhou L, Liu Y, Liu Z, Kong D, Duan M, Luo L. 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 

  75. Sun L, Li X, Fu Y, Zhu Z, Tan L, Liu F, Sun X, Sun C. GS6, a member of the GRAS gene family, negatively regulates grain size in rice. J Integr Plant Biol. 2013;55(10):938–49.

    CAS  PubMed  Google Scholar 

  76. Yang C, Li D, Mao D, Liu X, Ji C, Li X, Zhao X, Cheng Z, Chen C, Zhu L. Overexpression of microRNA319 impacts leaf morphogenesis and leads to enhanced cold tolerance in rice (Oryza sativa L.). Plant Cell Environ. 2013;36(12):2207–18.

    Article  CAS  PubMed  Google Scholar 

  77. De Paolo S, Gaudio L, Aceto S. Analysis of the TCP genes expressed in the inflorescence of the orchid Orchis italica. Sci Rep. 2015;5:16265.

    Article  PubMed  PubMed Central  Google Scholar 

  78. Lu C, Jeong DH, Kulkarni K, Pillay M, Nobuta K, German R, Thatcher SR, Maher C, Zhang L, Ware D, et al. Genome-wide analysis for discovery of rice microRNAs reveals natural antisense microRNAs (nat-miRNAs). Proc Natl Acad Sci U S A. 2008;105(12):4951–6.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  79. Huang F, Xu G, Chi Y, Liu H, Xue Q, Zhao T, Gai J, Yu D. A soybean MADS-box protein modulates floral organ numbers, petal identity and sterility. BMC Plant Biol. 2014;14:89.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  80. Ta KN, Sabot F, Adam H, Vigouroux Y, De Mita S, Ghesquiere A, Do NV, Gantet P, Jouannic S. miR2118-triggered phased siRNAs are differentially expressed during the panicle development of wild and domesticated African rice species. Rice (N Y). 2016;9(1):10.

    CAS  PubMed Central  Google Scholar 

  81. Song X, Li P, Zhai J, Zhou M, Ma L, Liu B, Jeong DH, Nakano M, Cao S, Liu C, et al. Roles of DCL4 and DCL3b in rice phased small RNA biogenesis. Plant J. 2012;69(3):462–74.

    Article  CAS  PubMed  Google Scholar 

  82. Jiang G, Xiang Y, Zhao J, Yin D, Zhao X, Zhu L, Zhai W. Regulation of inflorescence branch development in rice through a novel pathway involving the pentatricopeptide repeat protein sped1-D. Genetics. 2014;197(4):1395–407.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

  83. Koizuka N, Imai R, Fujimoto H, Hayakawa T, Kimura Y, Kohno-Murase J, Sakai T, Kawasaki S, Imamura J. Genetic characterization of a pentatricopeptide repeat protein gene, orf687, that restores fertility in the cytoplasmic male-sterile Kosena radish. Plant J. 2003;34(4):407–15.

    Article  CAS  PubMed  Google Scholar 

  84. Li W, Ahn IP, Ning Y, Park CH, Zeng L, Whitehill JG, Lu H, Zhao Q, Ding B, Xie Q, et al. The U-Box/ARM E3 ligase PUB13 regulates cell death, defense, and flowering time in Arabidopsis. Plant Physiol. 2012;159(1):239–50.

    Article  CAS  PubMed  PubMed Central  Google Scholar 

Download references

Acknowledgements

We are thankful for the helpful modifications and comments on the manuscript from Professor Aili Li, Institute of Crop Science, CAAS.

Funding

This work was supported by the Shandong Science and Technology Development Plan (2014GSF121001), the National Natural Science Foundation of China (31501312, 31401378 and 31601301), the Natural Science Foundation of Shandong Province (ZR2014CM006), the Agricultural Science and Technology Innovation Project of Shandong Academy of Agricultural Sciences (CXGC2016C09), the Youth Talent Program of Shandong Academy of Agricultural Sciences, and the Youth Foundation of Shandong Academy of Agricultural Sciences (2014QNZ02).

Availability of data and materials

The four small RNA datasets and two degradome datasets have been deposited in the NCBI SRA database with the accessions SRP076764 and SRP076763, respectively. The publicly available genome sequences of the Chinese Spring and plant miRNA sequences supporting the results of this article are available at the IWGSC website (http://www.wheatgenome.org), miRBase website (http://www.mirbase.org), PNRD (http://structuralbiology.cau.edu.cn/PNRD), PmiRExat (http://pmirexat.nabi.res.in) and miRNAs identified by Mayers et al [24]. The datasets supporting the conclusions of this article are included within the article and its additional files.

Authors’ contributions

GYL and GQS conceived and supervised this study. RZZ and GQS analysed the data. GQS, YLL, SJZ, JG, XDH, MLC, JW and WL performed the experiments. RZZ, GQS, and GYL interpreted the results and wrote the manuscript. All of the authors read and approved this manuscript.

Competing interests

The authors declare that they have no competing interests.

Consent for publication

Not applicable.

Ethics approval and consent to participate

Not applicable.

Author information

Author notes

    Authors and Affiliations

    1. Crop Research Institute, Shandong Academy of Agricultural Sciences, Jinan, 250100, China

      Guoqi Song, Rongzhi Zhang, Shujuan Zhang, Yulian Li, Jie Gao, Xiaodong Han, Mingli Chen, Jiao Wang, Wei Li & Genying Li

    2. Key Laboratory of Wheat Biology & Genetic Improvement on North Yellow & Huai River Valley, Ministry of Agriculture, Jinan, 250100, China

      Guoqi Song, Rongzhi Zhang, Shujuan Zhang, Yulian Li, Jie Gao, Xiaodong Han, Mingli Chen, Jiao Wang, Wei Li & Genying Li

    3. National Engineering Laboratory for Wheat & Maize, Jinan, 250100, China

      Guoqi Song, Rongzhi Zhang, Shujuan Zhang, Yulian Li, Jie Gao, Xiaodong Han, Mingli Chen, Jiao Wang, Wei Li & Genying Li

    Authors
    1. Guoqi Song
      View author publications

      You can also search for this author in PubMed Google Scholar

    2. Rongzhi Zhang
      View author publications

      You can also search for this author in PubMed Google Scholar

    3. Shujuan Zhang
      View author publications

      You can also search for this author in PubMed Google Scholar

    4. Yulian Li
      View author publications

      You can also search for this author in PubMed Google Scholar

    5. Jie Gao
      View author publications

      You can also search for this author in PubMed Google Scholar

    6. Xiaodong Han
      View author publications

      You can also search for this author in PubMed Google Scholar

    7. Mingli Chen
      View author publications

      You can also search for this author in PubMed Google Scholar

    8. Jiao Wang
      View author publications

      You can also search for this author in PubMed Google Scholar

    9. Wei Li
      View author publications

      You can also search for this author in PubMed Google Scholar

    10. Genying Li
      View author publications

      You can also search for this author in PubMed Google Scholar

    Corresponding author

    Correspondence to Genying Li.

    Additional files

    Additional file 1

    Table S1. Nucleotide sequences of primers for qRT-PCR of miRNAs and their targets. Table S2. Small RNA statistics and genome mapping information referred to the Chinese Spring genome. Table S3. Expression level of miRNAs in the control and cold-treated libraries. Table S4. Novel miRNA information, including mature, star, and precursor sequences, referred to the Chinese Spring genome. Table S5. Fold change of differentially expressed miRNAs and their targets by qRT-PCR. Table S6. Overlapped miRNAs in response to cold stress among wheat, Brachypodium, Medicago, Populus, and Arabidopsis. Table S7. Predicted targets of conserved miRNAs by the TargetFinder program. Table S8. Predicted targets of novel miRNAs by the TargetFinder program. Table S9. Targets of known miRNAs validated by degradome sequencing in the control and cold stress libraries. Table S10. Targets of novel miRNAs validated by degradome sequencing in the control and cold stress libraries. Table S11. Gene Ontology enrichment, including biological process, molecular function, and cellular components for the target genes of differentially expressed miRNAs after cold stress, referred as Chinese Spring. (ZIP 643 kb)

    Additional file 2

    Figure S1. Percentage of small RNAs mapped to the exon regions in control and cold stress libraries. White, grey and black rectangles represent the subgenome A, subgenome B, and subgenome D, respectively. Figure S2. Length distribution of tags (a) and distribution for the number of transcripts in the categories of cleavage (b) for the degradome sequencing between the control and cold stress libraries. (ZIP 10 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

    Song, G., Zhang, R., Zhang, S. et al. Response of microRNAs to cold treatment in the young spikes of common wheat. BMC Genomics 18, 212 (2017). https://doi.org/10.1186/s12864-017-3556-2

    Download citation

    • Received:

    • Accepted:

    • Published:

    • DOI: https://doi.org/10.1186/s12864-017-3556-2

    Keywords

    BMC Genomics

    ISSN: 1471-2164

    Contact us