Identification and characterization of CircRNAs involved in the regulation of wheat root length

Background Recent studies indicate that circular RNAs (circRNAs) may play important roles in the regulation of plant growth and development. Plant roots are the main organs of nutrient and water uptake. However, whether circRNAs involved in the regulation of plant root growth remains to be elucidated. Methods LH9, XN979 and YN29 are three Chinese wheat varieties with contrasting root lengths. Here, the root circRNA expression profiles of LH9, XN979 and YN29 were examined by using high-throughput sequencing technology. Results Thirty-three and twenty-two differentially expressed circRNAs (DECs) were identified in the YN29-LH9 comparison and YN29-XN979 comparison, respectively. Among them, ten DECs coexisted in both comparisons. As the roots of both LH9 and XN979 were significantly larger and deeper than YN29, the ten DECs coexisting in the two comparisons were highly likely to be involved in the regulation of wheat root length. Moreover, three of the ten DECs have potential miRNA binding sites. Real-time PCR analysis showed that the expression levels of the potential binding miRNAs exhibited significant differences between the long root plants and the short root plants. Conclusions The expression levels of some circRNAs exhibited significant differences in wheat varieties with contrasting root phenotypes. Ten DECs involved in the regulation of wheat root length were successfully identified in which three of them have potential miRNAs binding sites. The expression levels of putative circRNA-binding miRNAs were correlated with their corresponding circRNAs. Our results provide new clues for studying the potential roles of circRNAs in the regulation of wheat root length. Electronic supplementary material The online version of this article (10.1186/s40659-019-0228-5) contains supplementary material, which is available to authorized users.


Background
Plant roots are the main organs of nutrient and water uptake from soil. Genetic improvement of root traits is of vital importance for improving crop nutrient and water use efficiencies [1][2][3][4][5]. Understanding the molecular mechanisms controlling plant root growth and development is beneficial to molecular breeding aimed at improving root traits and resource utilization efficiency of crops. Previous studies have shown that plant root growth and development are controlled by phytohormones and sustained by the root apical meristem (RAM) [6][7][8][9][10][11][12]. Besides, redox regulation and the balance of ROS also have important roles in maintaining RAM activity [13][14][15]. Therefore, the growth and development of plant roots are regulated by complex hormone signals and other related pathways.
Wheat (Triticum aestivum L.) is one of the most important food crops in the world. In wheat, circR-NAs have been shown to be involved in response to drought stress. Wang et al. isolated 88 circRNAs and found that 62 circRNAs were differentially expressed under drought stress conditions compared with control [38]. In a recent study, Ren et al. identified six circRNAs involved in the common response to nitrogen deficiency stress and 23 circRNAs involved in the regulation of low nitrogen-promoted root growth in wheat [41]. These studies show that circRNAs may play important roles in responding to abiotic stresses. However, it is unclear whether circRNAs are involved in the regulation of wheat root growth under normal growth conditions. To explore this question, the root circRNA expression profiles of three wheat varieties (including two long root varieties and one short root variety) were obtained using high-throughput sequencing technology. Differentially expressed circRNAs (DECs) were identified and further validated using real-time PCR technology. Target miRNAs of the DECs were predicted and the expression levels of these miRNAs in the roots of LH9, XN979 and YN29 were also examined. This is the first report on the identification of differentially expressed circRNAs in wheat varieties with contrasting root length. The results are helpful for further investigate of the potential roles of circRNAs in regulating wheat root length.

Plant materials
XN979, LH9 and YN29 are three Chinese wheat varieties. The roots of LH9 and XN979 are significantly longer than those of YN29. Here, XN979, LH9 and YN29 were selected as materials to identify circRNAs involved in the regulation of wheat root length.

Plant growth conditions and evaluation of root phenotype
Seed sterilization, germination and the growth conditions of wheat plants were conducted as previously described by Ren et al. [42]. Plants were randomly placed and grown in a greenhouse with six replications each. The maximum root length (MRL) and total root length (TRL) of XN979, LH9 and YN29 were analyzed by using WinRHIZO software (Regent Instruments, Canada) after 15 days of transferring. The developmental stages of wheat plants were Zadoks growth scale 13 [43]. The roots of XN979, LH9 and YN29 were fast-frozen using liquid nitrogen for RNA extraction.

Libraries construction and sequencing
Root total RNA was extracted using Trizol reagent following the manufacturer's procedure. The concentration and purity of total RNA were measured by NanoDrop ND-1000 spectrophotometer (NanoDrop Technologies, Wilmington, DE, USA). Approximately 10 µg total RNA was used to deplete the ribosomal RNA (rRNA) following the instructions of the Epicentre Ribo-Zero Gold Kit (Illumina, San Diego, USA). The remaining RNAs were used as templates for the construction of cDNA libraries following the instructions of the RNA-Seq sample preparation kit (Illumina, San Diego, USA) [41]. Three biological replicates were analyzed. The samples were named as XN979-1, -2, -3; LH9-1, -2, -3; YN29-1, -2, -3, respectively. Paired-end sequencing was performed on an Illumina Hiseq 2500 platform (Hangzhou Shangyi biotechnology company, Hangzhou, China).

Identification of circRNAs
The obtained clean reads were tempted to align with wheat reference genome (Triticum aestivum TGACv1.0) using the bowtie2 (bowtie2-2.2.2) alignment method. The reads of the linear RNA can be mapped appropriately to wheat reference genome, while the reads at the loop forming junctions of circRNAs cannot be directly aligned to the wheat reference genome [41]. Then the find_circ software was employed to detect head-to-tail splicing (back-spliced) of the remaining unmapped RNAseq reads (default setting). The detected back-spliced reads were further filtered to predict circRNAs following the recommended setting rules (GU/AG appears on both sides of the splice site, clear breakpoint can be detected; ≤ 2 mismatches; The length of the circRNA junctions ≤ 100 kb) [41,44].

Differential expression analysis of circRNAs
Differential expression analysis of circRNAs between different varieties was performed using the DEseq R package [41]. Only the circRNAs with p values ≤ 0.05 and |log2 (foldchange)| ≥ 1 were defined as DECs [45,46]. psRNATarget software was used to predict circRNA-miRNA interactions of the DECs [47].

Real-time PCR analysis
To confirm and quantify the predicted circRNAs, divergent primers of five randomly selected DECs were designed based on the flanking sequences of the headto-tail splicing site of each circRNA (Additional file 1: Table S1). The cDNA samples were used as templates and mixed with primers and SYBR Green PCR Real Master Mix (Tiangen, China) for real-time PCR analyses. In the quantitative analysis of selected miRNAs, miRcute Plus miRNA First-Strand cDNA Synthesis Kit (Tiangen, China) was used for cDNA synthesis and miRcute Plus miRNA qPCR Detection Kit (Tiangen) was used for quantitative analysis according to the kit instructions. The reverse primers were provided in the miRcute Plus miRNA qPCR Detection Kit and the forward ones were designed according to the instructions of the kit (Additional file 1: Table S1). Real-time PCRs were performed on a Thermal Cycler CFX96 Real-Time System (Bio-Rad, USA). The program settings are as follows: 95 °C for 5 min, then 40 cycles of 95 °C for 15 s, 60 °C for 15 s, 72 °C for 15 s. TaActin was used as an internal reference gene to normalize the expression levels of the investigated DECs and miRNAs. SPSS 21.0 software was used to analyze the statistical significance of the data.

Root phenotypes of different wheat varieties
The root phenotypes of LH9, XN979 and YN29 are shown in Fig. 1. There are significant differences between the long root varieties (LH9, XN979) and the short root variety (YN29). The MRL of LH9 and XN979 were 13.7 cm and 9.2 cm longer than that of YN29 (Fig. 1a), and the TRL of LH9 and XN979 were 302.4 cm and 211.9 cm longer than that of YN29 respectively (Fig. 1b). Overall, the roots of LH9 and XN979 were significantly larger and deeper than those of YN29 (Fig. 1c).

CircRNAs identification
To determine circRNAs involved in the regulation of root length, nine cDNA libraries from wheat roots of LH9, XN979 and YN29 were constructed and sequenced. Over 100 million raw reads were generated in each library (among 100.9 to 104.7 million raw reads in each library). About 76% to 85% reads were successfully mapped to wheat genome (TGAC v1.0). Find_circ software was employed to detect head-to-tail splicing (back-spliced) of the remaining 15-24% unmapped RNA-seq reads. The number of circRNAs identified in each sample ranged from 285 to 478. Among them, more than 70% of cir-cRNAs are exonic circRNAs (Fig. 2a). The proportion of intronic circRNAs is the lowest, no more than 6.5% in each sample. The proportion of intergenic circRNAs is between the two, accounting for 18.7-24.9% (Fig. 2a).
Since the roots of both LH9 and XN979 are significantly larger and deeper than those of YN29 (Fig. 1c), the ten DECs coexisting in the two comparisons are highly likely to be involved in the regulation of wheat root length.

Real-time PCR verification of DECs
Five randomly selected DECs were selected for data verification of RNA-seq by using real-time PCR technology. Results showed that the expression levels of the five DECs detected using real-time PCR technology matched well with the results of RNA-seq (Fig. 3), indicating that the results of RNA-seq are reliable.

The regulation of miRNAs levels by circRNAs
It has been reported that circRNAs can act as ceRNAs (competitive endogenous RNA) of miRNAs, inhibiting their binding to mRNA molecules, thereby regulating gene expression [21,44]. To detect whether the ten DECs coexisted in both comparisons can perform the function, psRNATarget software was employed to predict potential miRNA binding sites of the circRNAs. Results showed that three out of the ten DECs were predicted to have two to nine corresponding miRNAs binding sites (Table 1). CircRNA23 have putative miRNA binding sites for tae-miR1134 and tae-miR5085. CircRNA3645 have putative miRNA binding sites for tae-miR1133, tae-miR1137a and tae-miR9773. According to the psRNATarget-based prediction, circRNA2473 can bind to nine different miRNAs. To investigate whether these circRNAs regulate the expression of the corresponding miRNAs, eight potential circRNA-binding miRNAs were selected to check their expression levels in the long root plants (LH9 and XN979) and the short root plants (YN29). Results showed that six of them exhibited significant differences between the long root plants (LH9 and XN979) and the short root plants (YN29) (Fig. 4).

Discussion
Deeper rooting is beneficial for plant production and survival under water and nitrogen deficient conditions [48][49][50]. Therefore, genetic improvement of root traits is essential for improving crop nutrient and water use efficiencies. Recent studies have shown that circRNAs are widespread in plants and play important roles in regulating plant growth and development, and stress tolerance [23][24][25][26][27][28][29][30][31][32][33][34][35][36][37][38][39][40][41]. However, whether circRNAs, as a new type of non-coding RNA molecules, participate in the regulation of root growth is not clear. Here, an integrated comparative root transcriptomic study of long root plants (LH9 and XN979) and short root plants (YN29) was conducted to explore this question. The number of circRNAs identified in each sample ranged from 285 to 478. More than 70% circRNAs belong to exonic circRNAs category. The intronic circRNAs only account for 3.9-6.5% of the total number of circRNAs in different samples (Fig. 2a). The proportion of intronic circRNAs is similar to the results of previous reports in Arabidopsis (3.8% intronic circR-NAs), tomato (3.6%) and wheat (less than 6.5%) [16,19,41]. In the YN29-LH9 comparison and YN29-XN979 comparison, 33 and 22 DECs were identified, respectively (Fig. 2b). In total, ten circRNAs existed in both the YN29-LH9 comparison and the YN29-XN979 comparison. As has been mentioned above, both the roots of LH9 and XN979 were significantly larger and deeper than those of YN29 (Fig. 1c), therefore, the ten coexisted DECs in both comparisons should play more critical roles in the regulation of root length. Previous studies in animal and human have confirmed that circRNAs can act as miRNA sponges to capture tae-miR1133 tae-miR9773 tae-miR164 tae-miR1125 tae-miR1134 tae-miR1118 tae-miR1137a tae-miR5085 Relative expression LH9 XN979 YN29 Fig. 4 Relative expression analysis of potential binding miRNAs of the DECs in the roots of LH9, XN979 and YN29. Significant difference at p < 0.05 is indicated by different letters above the columns miRNAs from their target genes via ceRNA (competing endogenous RNAs) networks [51,52]. Similar reports have recently been made in plants [26,34,38,41]. To reveal whether the ten DECs can target miRNAs and participate in transcriptional regulation of genes, psR-NATarget software was used to identify potential miRNA binding sites. As expected, three of the ten circRNAs had two to nine putative miRNA-binding sites (Table 1), which was consist with another study in wheat [38]. Moreover, the results of real-time PCR analysis showed that the expression levels of six out of the eight investigated potential binding miRNAs exhibited significant differences between the long root plants (LH9 and XN979) and the short root plants (YN29) (Fig. 4). Therefore, the expression levels of these miRNAs were correlated with their corresponding circRNAs. These results indicate that the three circRNAs may bind to the target miRNAs and modulate the transcript levels of these miRNAs. Interestingly, it has been reported that some putative target miRNAs are involved in the regulation of plant growth and development [53][54][55][56][57]. For example, miR164 has been reported controlling root development in Arabidopsis, maize and potato [53][54][55]. In Arabidopsis, one of miR164 target genes is NAC1, which transduces auxin signals for lateral root emergence [58,59]. In mir164 mutants, plants express less miR164 and more NAC1 mRNA, and produced more lateral roots. Moreover, the mutant phenotypes can be complemented by expression miR164a and miR164b genomic sequences [53]. Besides, the target genes of the potential circRNA binding miRNAs such as miR1122, miR1125, miR1134 and miR1133 are also involved in plant growth, development, metabolism and stress response [56,57]. Therefore, these circRNAs may regulate wheat root length by modulating target miR-NAs levels. We also noticed that different circRNAs may have the same miRNA binding site, and one circRNA may have several corresponding miRNA binding sites (Table 1). Since the mechanism of how circRNAs regulate their target miRNA are not clear in plants so far and these miRNA binding sites are predicted by psRNATarget software, whether circRNA can bind to each of them still requires further experimental validation. Furthermore, how do the circRNAs that do not have miRNAs binding sites work? What are the mechanisms by which different circRNAs coordinately regulate the length of wheat roots? There is still a lot of work to be done to clarify the regulatory roles of these circRNAs.

Conclusions
Our study revealed that the expression levels of some cir-cRNAs in roots exhibited significant differences between the long root plants (LH9 and XN979) and the short root plants (YN29). Ten DECs involved in the regulation of wheat root length were successfully identified, three of which have potential miRNAs binding sites. The expression levels of putative circRNA-binding miRNAs were correlated with their corresponding circRNAs. These results provide new clues for investigating the functions of circRNAs in the regulation of wheat root length.

Additional files
Additional file 1: Table S1. Primers used for real-time PCR.