Skip to main content

LncRNA TCONS_00021861 is functionally associated with drought tolerance in rice (Oryza sativa L.) via competing endogenous RNA regulation



Water deficit is an abiotic stress that retards plant growth and destabilizes crop production. Long non coding RNAs (lncRNAs) are a class of non-coding endogenous RNAs that participate in diverse cellular processes and stress responses in plants. lncRNAs could function as competing endogenous RNAs (ceRNA) and represent a novel layer of gene regulation. However, the regulatory mechanism of lncRNAs as ceRNA in drought stress response is yet unclear.


In this study, we performed transcriptome-wide identification of drought-responsive lncRNAs in rice. Thereafter, we constructed a lncRNA-mediated ceRNA network by analyzing competing relationships between mRNAs and lncRNAs based on ceRNA hypothesis. A drought responsive ceRNA network with 40 lncRNAs, 23 miRNAs and 103 mRNAs was obtained. Network analysis revealed TCONS_00021861/miR528-3p/YUCCA7 regulatory axis as a hub involved in drought response. The miRNA-target expression and interaction were validated by RT-qPCR and RLM-5’RACE. TCONS_00021861 showed significant positive correlation (r = 0.7102) with YUCCA7 and negative correlation with miR528-3p (= -0.7483). Overexpression of TCONS_00021861 attenuated the repression of miR528-3p on YUCCA7, leading to increased IAA (Indole-3-acetic acid) content and auxin overproduction phenotypes.


TCONS_00021861 could regulate YUCCA7 by sponging miR528-3p, which in turn activates IAA biosynthetic pathway and confer resistance to drought stress. Our findings provide a new perspective of the regulatory roles of lncRNAs as ceRNAs in drought resistance of rice.

Peer Review reports


Rice (Oryza sativa L.) is a major staple crop worldwide and is highly susceptible to water deficit. Drought stress is a limiting factor in rice production and has caused severe losses on the rice yields over the past decades. Plants can adapt to stresses by modulating multiple signaling pathways, among which noncoding RNAs play essential roles [1].

Long non-coding RNAs (lncRNAs) are a diverse class of ncRNAs larger than 200 bp in length without apparent coding potential. LncRNAs play widespread roles in virtually every biological process in plants [2]. RNA sequencing technologies have provided powerful tools for transcriptomic research and shed light into the lncRNAs of both model and crop plants, e.g. Arabidopsis [3, 4], wheat [5], maize [6, 7], rice [8, 9] and Medicago truncatula [10]. In particular, many stress-responsive lncRNAs were revealed in rice. lncRNAs were found to mediate anti-pathogen immunity [1, 7] as well as the response to abiotic stresses such as drought [11], salinity [12], heavy metal [13, 14] and nutrient deficiency [15].

Competing endogenous RNA (ceRNA) is a widespread type of gene regulation for lncRNA. LncRNAs acting as potential ceRNAs harbor microRNA response elements (MREs), thereby they compete with mRNAs for shared miRNA and thus regulate miRNA-mediated gene silencing. The first ceRNA was reported in Arabidopsis and the study identified an endogenous lncRNA, IPS1 (Induced by Phosphate Starvation 1). IPS1 acts as target mimic of miR399 and inhibits miR399-induced cleavage of PHO2, which modulates phosphate accumulation [16].

Many lncRNA acting as endogenous miRNA target mimics (eTMs) have been identified in plants. lnc_253 and lnc_973 may function as eTMs for miR156e under salt stress in cotton [17]. In tomato, lncRNA23468 acts as a ceRNA to modulate NBS-LRR genes by decoying miR482b under Phytophthora infestans [18]. Feng et al. profiled lncRNAs in Cd-stressed rapeseed and identified 67 lncRNAs as competing target mimics for 36 Cd responsive miRNAs [19]. In pigeonpea, lncRNA_1231 during flowering sequesters miR156b through the target mimicry manner leading to increased expression of flower-specific SPL-12 transcription factor and regulates flower development [20]. lncRNA39026 acts as eTMs of miR168a and may decoy miRNAs to upreglate SlAGO1, and enhances tolerance to Phytophthora infestans [21]. In Cassava, linRNA159 and linRNA340 were predicted to be target mimics of miR164 and miR169 under cold stress [22]. In winter wheat, lncRNAs were found to competitively bind with miR398 and regulate CSD1 expression, which improves cold resistance [5]. However, only a little proportion of these lncRNAs have been experimentally and functionally characterized.

In the current study, we systematically analyzed the transcriptional landscape of drought-responsive lncRNAs in rice to identify differentially expressed lncRNAs, miRNAs and mRNAs. A deregulated lncRNA-associated ceRNA network was constructed according to the ceRNA hypothesis. A ceRNA regulatory triplet were identified, validated and functionally characterized. Finally, a hypothetic model based on lncRNA-miRNA-mRNA-IAA signaling was established to reveal the role of lncRNAs in drought response. Our findings may provide insights into the mechanism of lncRNAs in drought response regulation.


Transcriptome sequencing of rice in response to drought stress

High-throughput paired-end sequencing (Illumina) was performed in control and PEG treated rice seedlings. A total of 578.64 million raw reads were obtained from 6 samples. After adapter trimming and filtration, 564.11 million clean reads were obtained and 93.79 % were mapped to the rice reference sequence. A total of 301,926 transcripts were reconstructed with cufflink 2.0. We obtained 542 high-confidence putative lncRNAs, including 48 intronic lncRNAs, 331 antisense lncRNAs and 163 intergenic lncRNAs. Most of lncRNAs lie between FPKM value 1–10. Box plots (Fig. 1a and c) illustrate the overall distribution of the expression level of lncRNA, mRNA and miRNA respectively. The lncRNA density distribution was given in Fig. 1d. As shown in Fig. 1e, most lncRNAs(93.81 %)were 201–1900 bp in length and 93.47 % of lncRNAs had no more than 3 exons (Fig. 1f). The short transcript length and low exon number were two common features of lncRNAs.

Fig. 1

Characteristics of lncRNAs expression. a: Expression levels of lncRNAs. b: Expression levels of mRNAs. c: Expression levels of miRNAs. d: Expression density of lncRNAs. e: Lengths distributions of lncRNAs. f: Exon numbers distribution of lncRNAs

Expression Profiles of lncRNAs, miRNAs and mRNAs

Deferential expression was analyzed with the cutoff log2FC > 1 and adjusted P-value < 0.05. DE lncRNAs, miRNAs and mRNAs were displayed in volcano plot (Fig. 2a) and heatmap (Fig. 2b). A total of 98 DE-lncRNAs (64 upregulated and 34 downregulated), 57 DE-miRNAs (34 upregulated and 23 downregulated) and 7020 DE-mRNAs (3222 upregulated and 3798 downregulated) were identified in the PEG treated group as against the control, respectively.

Fig. 2

Analysis of differentially expressed of lncRNAs, miRNAs and mRNAs. a: volcano plot of the differential lncRNAs, miRNAs and mRNAs. Red and green dots indicate the upregulation and downregulation, respectively. Grey dots represent non-differential expression. b: Heatmap for the differential lncRNAs, miRNAs and mRNAs. The red indicates upregulation, whereas blue represents downregulation

ceRNA Network was constructed in rice leaves

To evaluate the lncRNA-associated ceRNA interaction landscape involved in drought stressed rice, we constructed the lncRNA-mRNA-miRNA network according to the “ceRNA hypothesis” as follows:

First, miRNA-target interaction was predicted via nucleotide sequence and coexpression values, respectively. In total, we obtained 4739 miRNA-mRNA and 381 miRNA-lncRNA interactions. The lncRNA‐mRNA interaction was identified based on the Pearson correlation coefficient of the expression value.

Next, ceRNAs were filtered using ceRNA score, which was defined as the number of MREs between shared mRNAs and lncRNAs versus the total number of lncRNA MREs. For a given lncRNA-miRNA-mRNA triplet, both mRNAs and lncRNAs are targeted and negatively co-expressed with a common miRNA. As a result, 2727 lncRNA-miRNA-mRNA triplets were identified, consisting of 72 lncRNAs, 44 miRNAs and 257 mRNAs.

Finally, to make the networks more concise and robust, we further filtered the 103 key DE mRNAs in drought response through literature mining and functional annotation. Based on the selected mRNAs implicated in drought response, we constructed a network of 180 lncRNA-miRNA-mRNA triplets, including 40 lncRNAs, 23 miRNAs and 103 mRNAs (Fig. 3).

Fig. 3

lncRNA-based ceRNA network. Red, blue and green nodes represent lncRNAs, miRNAs and mRNAs, respectively. The edges represented the competing interactions among them

Hub nodes with high connectivity are important to a network. In our analysis, the hub was defined as the node with highest degree. miR528-3p with the highest degree of 15 was found to be the core component of the network. Among the miR528-3p involved mRNA-lncRNA coexpression, TCONS_00021861 relative to YUCCA7 (LOC4331709, Indole-3-pyruvate monooxygenase) had most significant positive correlation (r = 0.89, p = 0.016). Therefore, the TCONS_00021861/miR528-3p/YUCCA7 triplet was chosen for further functional investigation.

RT-qPCR confirmed the differential expression of drought-responsive lncRNAs

RT-qPCR was performed to validate TCONS_00021861, miR528-3p and YUCCA7 expression levels. TCONS_00021861 and YUCCA7 were downregulated, whereas miR528-3p was upregulated in the samples obtained from PEG‐treated rice compared with the control samples (Fig. 4a). The results were consistent with lncRNA sequencing. We then assessed the correlations among TCONS_00021861, miR528-3p and YUCCA7. As illustrated in Fig. 4b, the significant positive correlation between TCONS_00021861 and YUCCA7, and their negative correlation to miR528-3p expression levels support the ceRNA hypothesis.

Fig. 4

Verification of competitive interaction between TCONS_00021861/miR528-3p/YUCCA7. a: The expression levels of TCONS_00021861, YUCCA7 and miR528-3p were examined by RT-qPCR in drought stressed vs. control. b: Expression correlation between TCONS_00021861/miR528-3p/YUCCA7 in stressed leaf tissues. c: Potential miR528-3p cleavage sites on TCONS_00021861 and YUCCA7 identified by RLM-5′ RACE. d: YUCCA7 expression in response to overexpression of TCONS_00021861 and miR528-3p. e: IAA content of mesophyll cells determined by HPLC. f: Pearson correlation coefficients between TCONS_00021861/miR528-3p/YUCCA7. EV: empty vector; LNC: TCONS_00021861 overexpression; MIR: miR528-3p overexpression; LNC + MIR: co-overexpression of TCONS_00021861 and miR528-3p; C: the control without any treatment; D: treatment by 20 % PEG 6000; *P < 0.05, **P < 0.01, ***P < 0.001

RLM-5’ RACE confirmed cleavage of TCONS_00021861 and YUCCA7 by miR528-3p

Interaction between miR528-3p and its targets, TCONS_00021861 and YUCCA7, were experimentally confirmed by RLM-5’ RACE. After sequencing the RLM-5’ RACE PCR products, it’s found that both TCONS_00021861 and YUCCA7 had specific cleavage sites in the complementary miR528-3p sequence (Fig. ​4c).

TCONS_00021861 specifically sequesters miR528-3p to upregulate YUCCA7

To test whether TCONS_00021861 could regulate YUCCA7 by decoying miR528-3p, we overexpressed TCONS_00021861/miR528-3p and both to obtain three independent transgenic lines. qPCR was performed to evaluate the transfection effect. As shown in Fig. 4d, overexpressed TCONS_00021861 induced the YUCCA7 levels. In contrast, the overexpression of miR528-3p reduced YUCCA7 level, and the reduction could be rescued by overexpressing TCONS_00021861 simultaneously. In summary, the data demonstrated that TCONS_00021861 is a positive regulator of YUCCA7 by sequestering miR528-3p.

Interaction between TCONS_00021861 and miR528-3p regulated IAA pathway

Since YUCCA7 is known to regulate IAA synthesis, we inferred that TCONS_00021861 could participate in IAA biosynthesis pathway via ceRNA. To test this hypothesis, transgenic lines overexpressing TCONS_00021861/miR528-3p and both were used for the IAA content analysis. IAA content was detected by HPLC. As shown in Fig. 4e, IAA level significantly increased in transgenic plants overexpressing TCONS_00021861 whereas decreased obviously in plants overexpressing miR528-3p. However, no obvious IAA change was observed in the transgenic rice overexpressing both TCONS_00021861 and miR528-3p.

Correlation analysis between TCONS_00021861/miR528-3p/YUCCA7 found that IAA content was positively correlated with TCONS_00021861 (r = 0.83), YUCCA7 (r = 0.93) and negatively correlated with miR528-3p (r = -0.70) (Fig. 4f). Therefore, it’s assumed that TCONS_00021861 may function under drought treatment to upregulate the YUCCA7 concentration, which further increased the IAA level.

Overexpression of TCONS_00021861 alleviated drought-induced growth reduction and ROS accumulation

Plant growth indices were measured to investigate the impact of TCONS_00021861 and miR528-3p on the resistance of rice to drought stress. Transgenic lines overexpressing TCONS_00021861/miR528-3p and both were used for the growth observation. As shown in Fig. 5a-e, plant growth was affected by drought stress in terms of leaf and root length. Compared with WT, the overexpression of TCONS_00021861 increased root dry weight, leaf dry weight, root length and leaf length by 10.37 %, 4.38 %, 12.19 and 9.88 % respectively. Indices in the miR528-3p overexpression groups invariably decreased (from 7.85 to 18.25 %). However, upon co-overexpression of TCONS_00021861 and miR528-3p, growth indices were rescued (from 4.05 to 18.76 %).

Fig. 5

Effects of overexpressed TCONS_00021861 and miR528-3p on growth indices and ROS accumulation of drought-stressed rice seedlings. a: Seedlings were photographed and measured 1 week after PEG simulated drought stress. b: root dry weight. c: leaf dry weight. d: root length. e: leaf length. f: H2O2 content. g: O2 content. All the data represent the mean ± standard deviation (SD) of 20 rice seedlings

Two main ROS species, H2O2 and O2, were quantitatively determined. As shown in Fig. 5f, before PEG treatment, H2O2 and O2 accumulations were relatively low and no obvious difference was observed between wildtype and transgenic lines. However, drought stress significantly increased the accumulation of H2O2, especially in miR528-3p overexpression lines. Co-overexpression of TCONS_00021861 and miR528-3p had a smaller increase range in H2O2. In TCONS_00021861 overexpression line, no obvious increase in H2O2 accumulation was detected. Similar result was observed in O2 accumulation (Fig. 5g). The results showed that TCONS_00021861 overexpression alleviates the accumulation of ROS species.

Interestingly, enhanced IAA content in TCONS_00021861 overexpression plants is consistent with the enhanced growth index and reduced ROS content (Fig. 4e), indicating that TCONS_00021861 could relieve drought-induced growth retardation and ROS accumulation via IAA signaling.

Ultrastructure by Transmission Electron Microscopy (TEM)

The ultrastructure of mesophyll cells was observed using TEM. As indicated in Fig. 6a, in control group, the morphology of mesophyll cells appeared to be normal. The chloroplasts were elliptical. The basal thylakoid sheets were neatly stacked and connected with the interstitial thylakoid, and the envelope was intact. The shape of mitochondria was regular with clear double membranes and tubular spines. The nucleus was elliptical with clear and complete nucleolus.

Fig. 6

Ultrastructure of rice mesophyll cells. a: control group. b: drought stressed group. c: TCONS_00021861 overexpression group. d: miR528-3p overexpression group; e: TCONS_00021861 and miR528-3p co-overexpression group. Red arrows indicate the damaging effect of drought stress. Ch: chloroplast, CM: cell membrane, CW: cell wall, G: grana, M: mitochondrion, N: nucleus, O: osmiophilic granules

In drought stressed group however, many typical signs of damage were observed (Fig. 6b) especially in chloroplasts, mitochondria, including the disrupted cellular and chloroplast membranes, mitochondrial swelling, globular chloroplasts, disorder of thylakoid arrangement, and a reduced thickness of granal stacking. The nucleus is relatively stable with slight deformation. Part of the nucleus melts and becomes an empty nucleus.

Compared with drought stressed group, the cellular organelles were less affected in the TCONS_00021861 overexpression group (Fig. 6c). The shape of chloroplast appeared slightly changed and the granal stacking remained intact, suggesting alleviated cell damage after TCONS_00021861 overexpression. In the miR528-3p overexpression group (Fig. 6d), the cellular and chloroplast membranes were distorted. The grana lamellae became disintegrated and osmiophilic granules increased, which were the signs of damage to chloroplast ultrastructure. When TCONS_00021861 and miR528-3p were coexpressed, as shown in Fig. 6e, the damage effect was somewhat similar to that of the control with drought stress, including mitochondrial swelling and disrupted chloroplast membranes. Compared with miR528-3p overexpression lines however, the number of osmiophilic granules reduced and thylakoid arrangement was relatively normal, indicating that the degree of damage was relatively minor.


Dwindling water resource is a main challenge for the quality and quantity of rice production. lncRNAs are essential regulatory components in response to abiotic stresses. Previous studies have well-characterized lncRNAs in a set of plant species by high-throughput sequencing [3,4,5,6,7,8, 10].

Sequestration of miRNAs via target mimicry represents an important post-transcriptional regulation function of lncRNAs [22]. Several lncRNA acting as eTMs have been reported in plants under various abiotic stress [5, 17,18,19,20,21,22,23,24]. Nevertheless, most of them were merely identified by expression correlation and bioinformatic prediction. The detailed regulatory mechanisms of lncRNAs remain largely unknown.

This study provided a systematic perspective on the regulatory function of lncRNAs in the drought tolerance of Oryza Sativa. A lncRNA-miRNA-mRNA interaction network was comprehensively integrated via high throughput screening and bioinformatics analysis.

Using RNA sequencing technology, we identified 542 unique lncRNAs based on strict criteria, which constitute a credible group of rice lncRNAs. Among them, we identified 98 lncRNAs that were differentially expressed under drought stress.

To investigate the possible functions of the lncRNAs as ceRNAs, we constructed the potential regulatory network of lncRNA-miRNA-mRNA considering the similarity of the miRNA regulatory patterns and the similarity in expression. We next filtered the key DE mRNAs based on the literature and functional annotation. As a result, a drought-responsive ceRNA network with 40 lncRNAs, 23 miRNAs and 103 mRNAs was built. TCONS_00021861/miR528-3p/YUCCA7 triplet is located at the hub of the ceRNA network, with the highest degree and the most significant positive correlation. We then explored the comprehensive functional landscape of the TCONS_00021861/miR528-3p/YUCCA7 triplet.

There have been several reports about miR528 and members of the YUCCA (YUC) family in drought resistance. miR528 was found to be deregulated and participate in the regulation of defense response to drought stress in Sugarcane [25] and Brachy podium [26]. The YUCCA gene family, which encodes a flavin monooxygenase-like enzyme, mediates auxin overproduction that endows plants with drought resistance. For example, activation of YUCCA7 enhances drought tolerance in Arabidopsis [27]. Transgenic poplar plants [28] and potato [29] expressing YUCCA6 exhibits high-auxin induced phenotypes and increased tolerance to water deficit.

In line with previous findings, this study found that lncRNA TCONS_00021861 and YUCCA7 were co-expressed and dramatically down-regulated, whereas miR528-3p expression was upregulated in the PEG-treated rice compared with the controlled group. Their differential expression and coregulation were verified experimentally by RT-qPCR (Fig. 4a). Further, TCONS_00021861 showed strong positive correlation (r = 0.7102) with YUCCA7 and strong negative correlation with miR528-3p (r = − 0.7483) (Fig. 4b). RLM-5’ RACE results indicated that both TCONS_00021861 and YUCCA7 had adequate sequence complementarity with miR528-3p (Fig. 4c), and might compete binding with miR-528‐3p via shared MREs.

Credibility of the TCONS_00021861/miR528-3p/YUCCA7 regulation was further characterized in transgenic plants. TCONS_00021861 is a positive regulator of YUCCA7 in vivo because transcript level of YUCCA7 increased in the transgenic plant overexpressing TCONS_00021861 relative to wild type. In contrast, the overexpression of miR528-3p reduced YUCCA7 level and this reduction was rescued by simultaneous overexpression of TCONS_00021861 (Fig. 4d). These data indicate that TCONS_00021861 could regulate the expression of YUCCA7 by sponging miR528-3p.

IAA is an important plant auxin that is synthesized by YUCCA family members. IAA plays major roles in many plant growth and developmental processes, which can be synthesized via Trp-dependent or Trp-independent routes [30]. Genomic and functional studies in Arabidopsis have shown that de novo IAA synthesis was dominated by Trp-dependent pathway. Mounting evidence shows the involvement of YUCCA7 in Trp-dependent IAA synthesis. YUCCA7 catalyzes the conversion of TAM to NHT, a rate-limiting step in auxin synthesis in multiple plants. Functional assay revealed that TCONS_00021861 could induce characteristic auxin overproduction phenotypes by increasing the expression level of IAA, which could be reversed by miR528-3p overexpression (Fig. 4e). These observations signify the probable role for TCONS_00021861 in auxin biosynthesis. Consistent with this notion is the alleviated mesophyll cell damage observed by TEM in the transgenic plants overexpressing the TCONS_00021861 (Fig. 6).

A further correlation analysis was conducted between the IAA content, TCONS_00021861, miR528-3p and YUCCA7 expression level to confirm their regulatory relationship. As illustrated in Fig. 4f, there was a significant positive Pearson correlation between IAA amounts and TCONS_00021861 as well as YUCCA7 expression level, whereas miR528-3p was negatively correlated, suggesting that TCONS_00021861 mainly regulates IAA biosynthesis via decoying miR528-3p and subsequently increasing the level of YUCCA7.

Based on the results above, we proposed a model for the involvement of TCONS_00021861/miR528-3p/YUCCA7 in the regulation of Trp-dependent pathway for IAA synthesis (Fig. 7). In this model, drought induces the overexpression of TCONS_00021861, which acts as a specific decoy for miR528-3p and upregulates YUCCA7. YUCCA7 is involved in Trp-dependent IAA biosynthesis as a rate-limiting step converting TAM to NHT. Upregulation of YUCCA7 thus rendered a greater accumulation of IAA in rice under water deficient conditions, promoted the auxin overproduction phenotypes and ultimately endowed rice with drought resistance.

Fig. 7

Putative model for the function of TCONS_00021861 in IAA biosynthesis under drought stress


In summary, this study constructed a ceRNA regulatory network for lncRNA-miRNA-mRNA interactions and found TCONS_00021861 as positive regulator of YUCCA7 gene through interaction with miR528-3p. Functional elucidation of TCONS_00021861 highlights its important roles in the regulation IAA biosynthesis in response to water deficiency, which provides a foundation for understanding of the mechanism of lncRNAs in drought response.


Rice culture and drought stress

The rice cultivar (Oryza sativa L. subsp. Japonica cv. Suxiang 3) used for this study were purchased from Zhongnongke company, Jiangsu, China. The seeds were grown in the greenhouse under 16/8 h photoperiod, 26 °C /18°C (day/night), 55 % RH, and light intensity of 250 µmol m− 2 s− 1. Six-leaf stage seedlings were exposed to 20 % PEG 6000 for 48 h to simulate drought stress of -0.5 MPa osmotic potential. The control group was cultured under well irrigation.

Library construction and sequencing

Total RNA was extracted from rice leaves using TRIzol (Invitrogen) following the manufacturer’s protocols. For each sample, 1.5 µg RNA was used for subsequent analysis. Ribo-Zero™ Magnetic Kit was used to remove rRNAs. The strand-specific sequencing library was generated for 150 bp paired-end run on Illumina HiSeq 4000 platform, at a depth of 120 million reads per library. The stranded RNA-seq data have been submitted to the NCBI SRA with accession number SRP316304.

Sequencing data analysis

After removing adapter contaminations or poly-N and low-quality sequences (Q20 < 20) from the raw reads, clean reads were mapped to the rice reference genome with Tophat2 [31] and assembled using Cufflinks (v2.1.1) [32] to construct transcriptomes. CPAT [33], CNCI [34] and CPC [35] and were used to filter transcripts with coding potential. Coding transcripts with CPAT=“yes”, CNCI scores > 0 and CNCI scores > 0 were removed.

Expression levels of lncRNA and mRNA were determined with RPKM, and miRNA expression level was calculated with transcript per million (TPM). The differential expresion between groups was determined by DEGseq (adjusted P-value < 0.05 and log2FC > 1).

Construction of the ceRNA network

Nucleotide sequence-based miRNA-target interaction was predicted by miRanda tools with default parameters. Coexpression-based miRNA‐target interaction was determined by calculating Pearson correlation coefficient (PCC) for miRNA-target pairs. Negatively correlated pairs PCC<-0.7 and p < 0.01were selected for further analysis. Overlapping results of miRanda and coexpression thus represent a candidate miRNA‐target pair.

Based on the ceRNA principle, we calculate the ceRNA score of an lncRNA-mRNA pair targeted by common miRNAs according to the following formula:

$$\text{c}\text{e}\text{R}\text{N}\text{A} \text{s}\text{c}\text{o}\text{r}\text{e}=\frac{\text{M}\text{R}\text{E} \text{f}\text{o}\text{r} \text{s}\text{h}\text{a}\text{r}\text{e}\text{d} \text{m}\text{R}\text{N}\text{A}\_\text{l}\text{n}\text{c}\text{R}\text{N}\text{A}}{\text{M}\text{R}\text{E} \text{f}\text{o}\text{r} \text{l}\text{n}\text{c}\text{R}\text{N}\text{A}}$$

The lncRNA-mRNAs correlation was analyzed using the PCC from matched lncRNA and mRNA expression profiles. The lncRNA-mRNA pairs with PCC > 0.7 and p < 0.01 were considered to have potential relevance. The shared pairs from ceRNA scoring and coexpression were considered as the true ceRNAs, and were then mapped using Cytoscape software (V.3.5.1) [36] to visualize the lncRNA-miRNA‐mRNA network.


miRNA, lncRNA and mRNA expression levels were analyzed by qRT-PCR. The actin gene was used as an internal control and the 2 − ΔΔCt method was used to quantify the expression level. All PCR reactions were performed by three technical replicates and three biological replicates.

RNA Ligase-Mediated 5’ RACE(RLM-5’ RACE)

RLM-5’ RACE was performed using 5’-Full RACE Kit (Takara, Japan). Briefly, 5 µg total RNA was ligated to RNA adapter followed by reverse transcription using the Oligo (dT) primer and superscript reverse transcriptase. PCR was performed with 3’ gene-specific primers and 5’ adaptor primers. The PCR product was cloned into pEASY-T1 vector and plasmid DNA was sequenced.

Vector construction and plant transformation

PCR product of miR528-3p was cloned into the KpnI and PstI sites of plant expression vector p1301-35 S-NOS. The full-length sequences of lncRNA were amplified and inserted into the pEASY-Blunt vector. After Sanger sequencing, the correct fragments were digested with XhoI and KpnI and cloned into the PMV2 vector. These plasmids were transformed into rice plants by the floral dip method [37]. The whole pMV2 plasmid was transfected as a negative control. T0 transgenic plants were grown to maturity. T1 seeds were harvested and germinated on half MS medium containing 25 mg/L hygromycin. T1 progenies were analyzed by PCR for the segregation patterns using TCONS_00021861 and miR528-3p specific primers. T2 seeds from independent T1 lines showing 3:1 Mendelian segregation ratios were harvested separately and were germinated on hygromycin-containing rooting medium. T2 seedlings with absence segregation patterns were considered to be homozygous lines and were then screened and confirmed by RT-PCR analysis. Among the positive T2 transformants screened, the homozygous lines with highest expression level of the target transgene (named as LNC5-1-4, MIR9-1-2, LNC + MIR4-1-1) were used for subsequent functional analyses.

Measurement of growth indices, ROS, IAA and TEM

Growth was evaluated in terms of elongation and quantity of leaf and root 48 h after stress. The H2O2 and O2 content were determined following previous research methods [38]. For IAA content quantification, 5 g leaves was homogenized and extracted with methanol and ethyl acetate. Extract of leaves was then subjected to HPLC on Agilent Zorbax Eclipse XDB-C18 column (250 × 4.6 mm, 5 μm). TEM images of mesophyll cells were obtained (FEI Tecnai Spirit) according to Ju et al. [39].

Statistical analysis

All data were analyzed using SPSS v22.0. Student’s t-test was performed to compare significant differences between the two groups of data. Data were presented as SD calculated from three replicates.

Availability of data and materials

The RNA-seq data have been submitted to the NCBI SRA database with accession number SRP316304 ( Other datasets supporting the conclusions of this article are included within the article and its additional files.



Competing endogenous RNA


endogenous miRNA target mimics


Indole-3-acetic acid


Long non coding RNA


microRNA response element


RNA ligase-mediated 5′ rapid amplification of cDNA ends


  1. 1.

    Yu Y, Zhou YF, Feng YZ, He H, Lian JP, Yang YW, et al. Transcriptional landscape of pathogen-responsive lncRNAs in rice unveils the role of ALEX1 in jasmonate pathway and disease resistance. Plant Biotechnol J. 2020;18(3):679–90.

    CAS  Article  Google Scholar 

  2. 2.

    Wang Y, Luo X, Sun F, Hu J, Zha X, Su W, et al. Overexpressing lncRNA LAIR increases grain yield and regulates neighbouring gene cluster expression in rice. Nat Commun. 2018;9(1):3516.

    Article  Google Scholar 

  3. 3.

    Lu Z, Xia X, Jiang B, Ma K, Zhu L, Wang L, et al. Identification and characterization of novel lncRNAs in Arabidopsis thaliana. Biochem Biophys Res Commun. 2017;488(2):348–54.

    CAS  Article  Google Scholar 

  4. 4.

    Severing E, Faino L, Jamge S, Busscher M, Kuijer-Zhang Y, Bellinazzo F, et al. Arabidopsis thaliana ambient temperature responsive lncRNAs. BMC Plant Biol. 2018;18(1):145.

    Article  Google Scholar 

  5. 5.

    Lu Q, Guo F, Xu Q, Cang J. LncRNA improves cold resistance of winter wheat by interacting with miR398. Funct Plant Biol. 2020;47(6):544–57.

    CAS  Article  Google Scholar 

  6. 6.

    Li L, Eichten SR, Shimizu R, Petsch K, Yeh CT, Wu W, et al. Genome-wide discovery and characterization of maize long non-coding RNAs. Genome Biol. 2014;15(2):R40.

    Article  Google Scholar 

  7. 7.

    Zhu M, Zhang M, Xing L, Li W, Jiang H, Wang L, et al. Transcriptomic Analysis of Long Non-Coding RNAs and Coding Genes Uncovers a Complex Regulatory Network That Is Involved in Maize Seed Development. Genes (Basel). 2017;8(10):274.

  8. 8.

    Huang X, Zhang H, Wang Q, Guo R, Wei L, Song H, et al. Genome-wide identification and characterization of long non-coding RNAs involved in flag leaf senescence of rice. Plant Mol Biol. 2021;105(6):655–84.

    CAS  Article  Google Scholar 

  9. 9.

    Qi W, Chen H, Yang Z, Hu B, Luo X, Ai B, et al. Systematic characterization of long non-coding rnas and their responses to drought stress in Dongxiang wild rice. Rice Science. 2020;27(1):21–31.

    Article  Google Scholar 

  10. 10.

    Pecrix Y, Staton SE, Sallet E, Lelandais-Briere C, Moreau S, Carrere S, et al. Whole-genome landscape of Medicago truncatula symbiotic genes. Nat Plants. 2018;4(12):1017–25.

    CAS  Article  Google Scholar 

  11. 11.

    Li P, Yang H, Wang L, Liu H, Huo H, Zhang C, et al. Physiological and transcriptome analyses reveal short-term responses and formation of memory under drought stress in rice. Front Genet. 2019;10:55.

    CAS  Article  Google Scholar 

  12. 12.

    Singh U, Khemka N, Rajkumar MS, Garg R, Jain M. PLncPRO for prediction of long non-coding RNAs (lncRNAs) in plants and its application for discovery of abiotic stress-responsive lncRNAs in rice and chickpea. Nucleic Acids Res. 2017;45(22):e183.

    CAS  Article  Google Scholar 

  13. 13.

    Chen L, Shi S, Jiang N, Khanzada H, Wassan GM, Zhu C, et al. Genome-wide analysis of long non-coding RNAs affecting roots development at an early stage in the rice response to cadmium stress. BMC Genomics. 2018;19(1):460.

    Article  Google Scholar 

  14. 14.

    Tang Z, Xu M, Ito H, Cai J, Ma X, Qin J, et al. Deciphering the non-coding RNA-level response to arsenic stress in rice (Oryza sativa). Plant Signal Behav. 2019;14(9):1629268.

    Article  Google Scholar 

  15. 15.

    Shin SY, Jeong JS, Lim JY, Kim T, Park JH, Kim JK, et al. Transcriptomic analyses of rice (Oryza sativa) genes and non-coding RNAs under nitrogen starvation using multiple omics technologies. BMC Genomics. 2018;19(1):532.

    Article  Google Scholar 

  16. 16.

    Franco-Zorrilla JM, Valli A, Todesco M, Mateos I, Puga MI, Rubio-Somoza I, et al. Target mimicry provides a new mechanism for regulation of microRNA activity. Nat Genet. 2007;39(8):1033–7.

    CAS  Article  Google Scholar 

  17. 17.

    Deng F, Zhang X, Wang W, Yuan R, Shen F. Identification of Gossypium hirsutum long non-coding RNAs (lncRNAs) under salt stress. BMC Plant Biol. 2018;18(1):23.

    Article  Google Scholar 

  18. 18.

    Jiang N, Cui J, Shi Y, Yang G, Zhou X, Hou X, et al. Tomato lncRNA23468 functions as a competing endogenous RNA to modulate NBS-LRR genes by decoying miR482b in the tomato-Phytophthora infestans interaction. Hortic Res. 2019;6:28.

    Article  Google Scholar 

  19. 19.

    Feng SJ, Zhang XD, Liu XS, Tan SK, Chu SS, Meng JG, et al. Characterization of long non-coding RNAs involved in cadmium toxic response in Brassica napus. RSC Adv. 2016;6(85):82157–73.

  20. 20.

    Antara D, Swati S, Kuldeep K, Tribhuvan KU, Singh NK, Kishor G. Non-coding RNAs having strong positive interaction with mRNAs reveal their regulatory nature during flowering in a wild relative of pigeonpea (Cajanus scarabaeoides). Mol Biol Rep. 2020;47(5):3305–17.

    Article  Google Scholar 

  21. 21.

    Cui J, Jiang N, Hou X, Wu S, Zhang Q, Meng J, et al. Genome-Wide Identification of lncRNAs and Analysis of ceRNA Networks During Tomato Resistance to Phytophthora infestans. Phytopathology. 2020;110(2):456–64.

    CAS  Article  Google Scholar 

  22. 22.

    Li S, Yu X, Lei N, Cheng Z, Zhao P, He Y, et al. Genome-wide identification and functional prediction of cold and/or drought-responsive lncRNAs in cassava. Sci Rep. 2017;7:45981.

    CAS  Article  Google Scholar 

  23. 23.

    Zhou D, Du Q, Chen J, Wang Q, Zhang D. Identification and allelic dissection uncover roles of lncRNAs in secondary growth of Populus tomentosa. DNA Res. 2017;24(5):473–86.

    CAS  Article  Google Scholar 

  24. 24.

    Wang M, Wu HJ, Fang J, Chu CC, Wang XJ. A long noncoding RNA involved in rice reproductive development by negatively regulating osa-miR160. Sci Bull. 2017;62:470–5.

  25. 25.

    Ferreira TH, Gentile A, Vilela RD, Costa GG, Dias LI, Endres L, et al. microRNAs associated with drought response in the bioenergy crop sugarcane (Saccharum spp.). PLoS One. 2012;7(10):e46703.

    CAS  Article  Google Scholar 

  26. 26.

    Budak H, Akpinar A. Dehydration stress-responsive miRNA in Brachypodium distachyon: evident by genome-wide screening of microRNAs expression. OMICS. 2011;15(11):791–9.

    CAS  Article  Google Scholar 

  27. 27.

    Lee M, Jung JH, Han DY, Seo PJ, Park WJ, Park CM. Activation of a flavin monooxygenase gene YUCCA7 enhances drought resistance in Arabidopsis. Planta. 2012;235(5):923–38.

    CAS  Article  Google Scholar 

  28. 28.

    Ke Q, Wang Z, Ji CY, Jeong JC, Lee HS, Li H, et al. Transgenic poplar expressing Arabidopsis YUCCA6 exhibits auxin-overproduction phenotypes and increased tolerance to abiotic stress. Plant Physiol Biochem. 2015;94:19–27.

    CAS  Article  Google Scholar 

  29. 29.

    Kim JI, Baek D, Park HC, Chun HJ, Oh DH, Lee MK, et al. Overexpression of Arabidopsis YUCCA6 in potato results in high-auxin developmental phenotypes and enhanced resistance to water deficit. Mol Plant. 2013;6(2):337–49.

    CAS  Article  Google Scholar 

  30. 30.

    Abu-Zaitoon YM, Bennett K, Normanly J, Nonhebel HM. A large increase in IAA during development of rice grains correlates with the expression of tryptophan aminotransferase OsTAR1 and a grain-specific YUCCA. Physiol Plant. 2012;146(4):487–99.

    CAS  Article  Google Scholar 

  31. 31.

    Kim D, Pertea G, Trapnell C, Pimentel H, Kelley R, Salzberg SL. TopHat2: accurate alignment of transcriptomes in the presence of insertions, deletions and gene fusions. Genome Biol. 2013;14(4):R36.

    Article  Google Scholar 

  32. 32.

    Trapnell C, Williams BA, Pertea G, Mortazavi A, Kwan G, van Baren MJ, et al. Transcript assembly and quantification by RNA-Seq reveals unannotated transcripts and isoform switching during cell differentiation. Nat Biotechnol. 2010;28(5):511–5.

    CAS  Article  Google Scholar 

  33. 33.

    Wang L, Park HJ, Dasari S, Wang S, Kocher JP, Li W. CPAT: Coding-Potential Assessment Tool using an alignment-free logistic regression model. Nucleic Acids Res. 2013;41(6):e74.

    CAS  Article  Google Scholar 

  34. 34.

    Sun L, Luo H, Bu D, Zhao G, Yu K, Zhang C, et al. Utilizing sequence intrinsic composition to classify protein-coding and long non-coding transcripts. Nucleic Acids Res. 2013;41(17):e166.

    CAS  Article  Google Scholar 

  35. 35.

    Kong L, Yong Z, Ye ZQ, Liu XQ, Ge GJNAR. CPC: assess the protein-coding potential of transcripts using sequence features and support vector machine. Nucleic Acids Res. 2007;35(Web Server issue):W345-9.

  36. 36.

    Shannon P, Markiel A, Ozier O, Baliga NS, Wang JT, Ramage D, et al. Cytoscape: a software environment for integrated models of biomolecular interaction networks. Genome Res. 2003;13(11):2498–504.

    CAS  Article  Google Scholar 

  37. 37.

    Clough SJ, Bent AF. Floral dip: a simplified method for Agrobacterium-mediated transformation of Arabidopsis thaliana. Plant J. 1998;16(6):735–43.

    CAS  Article  Google Scholar 

  38. 38.

    Kong F, Deng Y, Zhou B, Wang G, Wang Y, Meng Q. A chloroplast-targeted DnaJ protein contributes to maintenance of photosystem II under chilling stress. J Exp Bot. 2014;65(1):143–58.

    CAS  Article  Google Scholar 

  39. 39.

    Ju S, Wang L, Zhang C, Yin T, Shao S. Alleviatory effects of silicon on the foliar micromorphology and anatomy of rice (Oryza sativa L.) seedlings under simulated acid rain. PLoS One. 2017;12(10):e0187021.

    Article  Google Scholar 

Download references


Not applicable.


This work was supported by the National Natural Science Foundation of China (Grant 31770903) which has supported RNA-Seq analysis, Technology R&D Program of Suzhou (Grant SNG201905) which has supported RLM-5’ RACE analysis and Qinglan Project of Jiangsu Province, which has supported TEM analysis.

Author information




JC designed the experiments and wrote the manuscript. YZ performed the experiments. XQ edited the manuscript. All authors read and approved its content.

Corresponding author

Correspondence to Jiajia Chen.

Ethics declarations

Ethics approval and consent to participate

The rice cultivar used for this study were purchased from Zhongnongke company, Jiangsu, China. The experimental research on plants performed in this study complies with institutional, national and international guidelines.

Consent for publication

Not applicable.

Competing interests

The authors declare that they have no competing interests.

Additional information

Publisher’s Note

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

Supplementary Information

Additional file 1.

Statistics of raw, clean and mapped reads of all 6 samples.

Additional file 2.

Full list of lncRNA-miRNA-mRNA interactions in drought responsive ceRNA network.

Additional file 3.

Oligonucleotide primers used in RT-qPCR.

Additional file 4.

The specific primer sequences for RLM-5’ RACE of target genes.

Additional file 5.

Genetic analysis of the T1 progeny of transgenic rice lines.

Rights and permissions

Open Access This article is licensed under a Creative Commons Attribution 4.0 International License, which permits use, sharing, adaptation, distribution and reproduction in any medium or format, as long as you give appropriate credit to the original author(s) and the source, provide a link to the Creative Commons licence, and indicate if changes were made. The images or other third party material in this article are included in the article's Creative Commons licence, unless indicated otherwise in a credit line to the material. If material is not included in the article's Creative Commons licence and your intended use is not permitted by statutory regulation or exceeds the permitted use, you will need to obtain permission directly from the copyright holder. To view a copy of this licence, visit The Creative Commons Public Domain Dedication waiver ( applies to the data made available in this article, unless otherwise stated in a credit line to the data.

Reprints and Permissions

About this article

Verify currency and authenticity via CrossMark

Cite this article

Chen, J., Zhong, Y. & Qi, X. LncRNA TCONS_00021861 is functionally associated with drought tolerance in rice (Oryza sativa L.) via competing endogenous RNA regulation. BMC Plant Biol 21, 410 (2021).

Download citation


  • Drought stress
  • lncRNA
  • Network
  • ceRNA
  • IAA