- Research article
- Open Access
Genome-wide identification and functional analysis of lincRNAs acting as miRNA targets or decoys in maize
BMC Genomicsvolume 16, Article number: 793 (2015)
Long intergenic noncoding RNAs (lincRNAs) are endogenous non-coding RNAs (ncRNAs) that are transcribed from ‘intergenic’ regions of the genome and may play critical roles in regulating gene expression through multiple RNA-mediated mechanisms. MicroRNAs (miRNAs) are single-stranded small ncRNAs of approximately 21–24 nucleotide (nt) that are involved in transcriptional and post-transcriptional gene regulation. While miRNAs functioning as mRNA repressors have been studied in detail, the influence of miRNAs on lincRNAs has seldom been investigated in plants.
LincRNAs as miRNA targets or decoys were predicted via GSTAr.pl script with a set of rules, and lincRNAs as miRNA targets were validated by degradome data. Conservation analysis of lincRNAs as miRNA targets or decoys were conducted using BLASTN and MAFFT. The function of lincRNAs as miRNA targets were predicted via a lincRNA-mRNA co-expression network, and the function of lincRNAs as miRNA decoys were predicted according to the competing endogenous RNA (ceRNA) hypothesis.
In this work, we developed a computational method and systematically predicted 466 lincRNAs as 165 miRNA targets and 86 lincRNAs as 58 miRNA decoys in maize (Zea mays L.). Furthermore, 34 lincRNAs predicted as 33 miRNA targets were validated based on degradome data. We found that lincRNAs acting as miRNA targets or decoys are a common phenomenon, which indicates that the regulated networks of miRNAs also involve lincRNAs. To elucidate the function of lincRNAs, we reconstructed a miRNA-regulated network involving 78 miRNAs, 117 lincRNAs and 8834 mRNAs. Based on the lincRNA-mRNA co-expression network and the competing endogenous RNA hypothesis, we predicted that 34 lincRNAs that function as miRNA targets and 86 lincRNAs that function as miRNA decoys participate in cellular and metabolic processes, and play role in catalytic activity and molecular binding functions.
This work provides a comprehensive view of miRNA-regulated networks and indicates that lincRNAs can participate in a layer of regulatory interactions as miRNA targets or decoys in plants, which will enable in-depth functional analysis of lincRNAs.
Long noncoding RNAs (lncRNAs) are generally long transcripts of more than 200 nucleotide (nt) that lack a coding sequence (CDS) or open reading frame (ORF) [1, 2]. Despite exhibiting lower expression levels compared with mRNAs, lncRNAs can regulate gene expression at the transcriptional and post-transcriptional levels [3–7]. As one of the largest classes of lncRNAs, long intergenic noncoding RNAs (lincRNAs) are endogenous lncRNAs that are transcribed from ‘intergenic’ regions of the genome. They play critical roles in regulating multiple important biological processes in humans and other animals, including cell cycle regulation, immune surveillance and embryonic stem cell differentiation [8–14], while they primarily participate in the environmental stimulus response, vernalization and nodulation in plants, including Arabidopsis thaliana , Triticum aestivum , Cucumis sativus , Setaria italic , Populus trichocarpa  and Zea mays [20–22]. However, compared with animal lincRNAs, the functions of plant lincRNAs and their regulatory roles remain largely undiscovered.
Unlike lincRNAs, plant microRNAs (miRNAs) are approximately 21–24 nucleotide (nt) single-stranded, small non-coding RNAs that typically form near-perfect duplexes with their targets and mediate cleavage or translation repression at the post-transcriptional level [23, 24]. They play vital roles in regulating a broad range of biological metabolic processes, including roles in plant development, flowering time, leaf morphogenesis, hormone signaling and responses to environmental stresses, such as phosphate or/and sulfate stress [25–30]. miRNAs usually regulate the expression of their mRNA targets through cleavage in plants [31, 32]. However, recent studies suggest that miRNAs function in a more sophisticated way than was initially assumed. In addition to protein-coding RNAs acting as miRNA targets, lincRNAs can also be directly targeted by miRNAs for cleavage [19, 33–35].
More interestingly, lincRNAs can also serve as miRNA decoys, miRNA sponges, target mimicry, or target mimics to interfere with the miRNA-mediated regulation of their mRNA targets. Similar to the sequence-dependent interactions of miRNAs with their mRNA targets, miRNA decoys also rely on the sequence-dependent interaction of miRNAs with lincRNAs, except for the bulges in the middle of miRNA-lincRNA duplexes. If lincRNAs acting as miRNA decoys and mRNAs acting as miRNA targets can be bound by the same miRNAs, then lincRNAs could function as competing endogenous RNAs (ceRNAs); they could directly interact with the specific miRNA and sequester it in a type of target mimicry to protect target mRNAs from repression, which is known as the “ceRNA hypothesis” [36, 37]. In animals, the long noncoding RNA linc-MD1 can act as a miR133 and miR135 sponge and up-regulate muscle-specific expression of the respective miR133 and miR135 targets MAML1 and MEF2C . In plants, the classical example of an miRNA decoy is IPS1, which is a long non-coding RNA that contains an ath-miR399 decoy site and can serve as an miRNA decoy to inactivate ath-miR399 and up-regulate the expression of the ath-miR399 primary target PHO2 . In rice, it has been reported that two lincRNAs that act as decoys of miR160 and miR164 can regulate floral and/or seed development . Recently, 25 miRNA decoys from Arabidopsis and 94 miRNA decoys from rice were identified; overexpressing the decoys of miR160 and miRNA166 can alter plant development, indicating that ncRNAs, short ORF encoding genes and intergenic sequences acting as miRNA decoys are functional in plants .
Maize (Zea mays L.) is one of the most important crops worldwide. It serves as a food source for people around the world and as a model organism in genetics research . With the release of the maize genome, increasing amounts of transcriptome data; degradome data; and specific data on miRNAs, lincRNAs and mRNAs have been accumulated. It is now possible for us to investigate the function of lincRNAs as miRNA targets or decoys in maize. Here, lincRNAs acting as miRNA targets were initially identified based on degradome data, and lincRNAs that may act as miRNA decoys were subsequently predicted. To explore the function of lincRNAs acting as miRNA targets or decoys, a genome-scale network among miRNAs, lincRNAs acting as miRNA targets, lincRNAs acting as miRNA decoys, and mRNAs was first constructed. Then, the functions of lincRNAs acting as miRNA targets were predicted and annotated via a co-expression network between lincRNAs and mRNAs, and the functions of lincRNAs acting as miRNA decoys were predicted and annotated according to the ceRNA hypothesis. Our research demonstrates that lincRNAs can act as miRNA targets or decoys to mediate the regulation of gene expression, and the annotation of lincRNA functions will facilitate the validation of the lincRNA functions in the future.
LincRNA and cDNA data
Primary data on lincRNAs were first integrated from three published studies on maize, consisting of 1704 high-confidence lncRNAs, 439 lincRNAs, and 664 putative maize lncRNAs [20–22]. Then, lncRNAs that were not located in intergenic regions and lincRNAs that were small RNA precursors were filtered out, and a total of 1831 lincRNAs were obtained and used in further analyses (Additional file 1). To distinguish the lincRNAs from these three data sources, the first authors’ names were added to the IDs of the lincRNAs. Maize cDNA data were downloaded from MaizeGDB ftp://ftp.ensemblgenomes.org/pub/plants/release-22/fasta/zea_mays/.)
Data on mature miRNAs were downloaded from miRBase (version 21: June 2014, http://www.mirbase.org/) [43, 44], and 321 maize miRNAs were extracted. A total of 203 unique miRNAs were obtained after merging these sequences with different miRNA IDs.
The degradome data from maize were downloaded from NCBI’s Gene Expression Omnibus (GEO) with the accession numbers of SRX222260, SRX222262, SRX222264, SRX222266 (http://www.ncbi.nlm.nih.gov/sra/?term=SRP018376) [45–47]. The raw reads from the above data were first processed using the FASTAX-Toolkit to trim adapter sequences with many “N” and ignored reads that were less than 18 nt. Then, the redundant reads were merged and 3268059, 4106567, 2682186 and 2965163 unique reads were obtained from the optional nitrate root tip, low nitrate root tip, low nitrate leaf and optional nitrate leaf, respectively (Additional file 2).
Prediction of miRNA targets
The miRNA targets of lincRNAs or cDNAs were predicted using GSTAr.pl script, and the minimum free energy (MFE) of miRNA-lincRNA or miRNA-cDNA duplexes was calculated with the RNAhybrid program [48–50]. Then, a modified version of the CleaveLand4 program was used to identify the potential cleavage sites of miRNAs in the corresponding targets based on degradome data http://sites.psu.edu/axtell/software/cleaveland4/) . To obtain high-quality lincRNAs acting as miRNA targets and to distinguish those lincRNAs acting as miRNA decoys, the following rules were used: at most, one mismatch or indel was allowed between the 9th and 12th positions of the 5′ end of miRNA sequences, the total number of bulges or mismatches in the other regions was not allowed to exceed 4 nt, and no continuous mismatches were allowed [41, 51]. In addition, target plots indicating the abundance of each distinct read for the lincRNAs acting as miRNA targets were generated.
Prediction of miRNA decoys
LincRNAs potentially acting as miRNA decoys were predicted based on Wu’s methods with a slight modification [41, 52]. Generally, the following set of rules was used: (1) the number of mismatches or indels should be larger than 1 and less than 6 between the 9th and 12th positions of the 5′ end of the miRNA sequences; (2) perfect nucleotide pairing was required between the 2nd and 8th positions of the 5′ end of miRNA sequences; and (3) the number of mismatches and indels should be no more than 4 in other regions. These rules were implemented using in-house Perl scripts.
Conservation analysis of lincRNAs acting as miRNA targets or decoys
To investigate the conservation of lincRNAs acting as miRNA targets or decoys, five genomes of other monocotyledons (monocots) (Sorghum bicolor, Setaria italica, Panicum virgatum, Oryza sativa and Brachypodium distachyon) were downloaded from Phytozome (v9.1) (http://www.phytozome.net/) , and the lincRNA regions that paired with miRNA targets or decoys were searched against the 5 monocot genomes using BLASTN with a cutoff threshold of an E-value less than 1e-1 . Then, the significantly matched regions plus their flanking regions (100 bp in total) were obtained . Finally, multiple sequence alignment was conducted with MAFFT v6.864b, using parameter settings of maxiterate 1000 and localpair . If the identities between the conserved sites were greater than 80%, then the conserved sites were highlighted.
Construction of miRNA-lincRNA-mRNA networks
To infer the function of lincRNAs, networks were constructed based on the complementary pairs between miRNAs and lincRNAs and between miRNAs and mRNAs. The nodes in the networks consisted of miRNAs, lincRNAs acting as miRNA targets, lincRNAs acting as miRNA decoys, mRNAs acting as miRNA targets, and mRNAs acting as miRNA decoys. The miRNA-lincRNA-mRNA networks were visualized with Cytoscape 3.1.1 .
Functional prediction of lincRNAs acting as miRNA targets based on the lincRNA-mRNA co-expression networks
Fifty-four datasets, including 30 RNA-seq experiments performed in 13 different tissues (leaf, immature ear, immature tassel, seed, endosperm, embryo, embryo sac, anther, ovule, pollen, silk, root and shoot apical tissues), were applied to construct a co-expression network between lincRNAs acting as miRNA targets and mRNA genes [58–63]. The construction method was similar to that of Liao  and Hao . In general, the pipeline for constructing the co-expression network was as follows: (1) genes, including mRNAs and lincRNAs, whose variances ranked in the top 75 % of the expression profiles were retained; (2) the p-value of Pearson’s correlation coefficient (Pcc) was calculated for each pair of genes using Fisher’s asymptotic test in the WGCNA library of R , and these values were adjusted using the Bonferroni correction method; and (3) co-expression relationships showing adjusted p-values of less than 0.05 and ranking in the top 5 % and bottom 5 % of Pcc were selected for further analysis. The Bonferroni multiples test was executed using the multtest package from R. The co-expression networks were also visualized using Cytoscape .
Based on the co-expression network between lincRNAs acting as miRNA targets and mRNAs, we used the AgriGO toolkit and input the list of mRNA genes to predict the function of these lincRNAs .
Functional prediction of lincRNAs acting as miRNA decoys based on miRNA-lincRNA-mRNA networks
Based on the ceRNA hypothesis and gene ontology (GO) analysis, the function of lincRNAs acting as miRNA decoys can be speculated based on the miRNA-lincRNA-mRNA networks. AgriGO, an integrated web-based GO analysis toolkit, was employed for the functional annotation and enrichment analysis . The IDs of all of the listed mRNAs connected with lincRNAs acting as miRNA decoys were submitted for GO analysis, and the overrepresented GO terms in the “biological process”, “cellular component” and “molecular function” categories were obtained using Fisher’s exact test and the Bonferroni multiples test (P-value < 0.05).
Identification of lincRNAs as putative miRNA targets
Previous research has suggested that miRNAs play roles in regulating the expression of mRNAs, but the comprehensive patterns of miRNA regulation of lincRNAs remain unknown. To systematically investigate the miRNA-mediated regulatory mechanism of lincRNAs, a method for predicting miRNA targets among lincRNAs was applied (see Materials and Methods). The results revealed 789 miRNA-lincRNA interactions (Additional file 3). In total, 466 lincRNA targets were predicted for 165 miRNAs in Zea mays.
To eliminate potential false-positive lincRNAs predicted as miRNA targets, we applied degradome reads to validate miRNA targets using a modified version of the CleaveLand pipeline . The results showed that 42 miRNA-lincRNA duplexes were supported by the degradome reads, which were formed by 33 miRNAs and 34 lincRNAs (Table 1, Additional file 4). When the degradome reads were mapped on each lincRNA, the abundance of the degradome reads at each position of the lincRNAs and the cleaved positions in each lincRNA could be obtained. For example, the abundance of degradome reads and cleavage sites in the lincRNA Boerner_Z27kG1_17085, which can act as a target of zma-miR166h-5p, is shown in Fig. 1.
Conservation of lincRNAs as miRNA targets between six monocotyledons
To investigate the conservation of lincRNAs as miRNA targets, the lincRNA regions that paired with miRNAs were searched against the 5 genomes of monocots, and the significant matched regions plus their flanking regions were obtained. Conservation analysis was performed, and 12 of 33 miRNAs were found to show conserved target regions in lincRNAs among maize and three to five other species. For example, the sequence logo and multiple sequence alignment of zma-miR166n-5p targets in lincRNAs provide a precise description of the conservation of these target regions (Fig. 2). However, the lincRNA regions outside of the predicted miRNA binding sites were not conserved, except for the lincRNAs targeted by zma-miR160b/g-3p and zma-miR169l-3p, which were conserved among 4 and 5 species, respectively (Additional file 5). In summary, lincRNAs acting as miRNA targets are a common phenomenon among monocots.
Identification of lincRNAs acting as miRNA decoys
Previous studies have shown that the duplexes formed by miRNAs and miRNA decoys usually contain bulges or mismatches in the middle of the miRNA binding sites, which is thought to block the interaction between miRNAs and their specific mRNA targets [41, 52]. GSTAr.pl can efficiently identify sites with large bulges in the alignments between miRNAs and lincRNAs. Therefore, we used a computational pipeline to identify lincRNAs acting as miRNA decoys in maize. In total, we found that 86 lincRNAs that may act as miRNA decoys could be bound by 58 miRNAs and formed 104 miRNA-lincRNA duplexes (Table 2, Additional file 6).
Conservation of lincRNAs as miRNA decoys between six monocotyledons
Similar to the analysis of the conservation of lincRNAs as miRNA targets, a conservation analysis of lincRNAs acting as miRNA decoys was also performed between the lincRNA regions that paired with miRNAs. Altogether, 10 of 58 miRNAs showed conserved decoy regions in lincRNAs among four to six species. For example, the sequence logo and multiple sequence alignment of zma-miR171f-5p decoys provide a precise description of the conservation of decoy regions (Fig. 3). Except for the zma-miR159e-3p and zma-miR482-3p decoys, other lincRNAs as miRNA decoy sites were conserved, but all of the surrounding regions were non-conserved (Additional file 7).
LincRNAs may participate in miRNA-lincRNA-mRNA networks
Previous research has demonstrated that engineered miRNA decoys can affect the regulation of miRNAs in plants [39, 52, 67]. To investigate the function of lincRNAs acting as miRNA targets or miRNA decoys, comprehensive genome-wide networks mediated by miRNAs were constructed. The networks were composed of 9402 nodes and 10,529 edges, and the nodes included 78 miRNAs, 117 lincRNAs (lincRNAs acting as miRNA targets, lincRNAs acting as miRNA decoys) and 8834 mRNAs (mRNAs acting as miRNA targets, mRNAs acting as miRNA decoys) (Fig. 4, Additional file 8). There were 42 interactions between miRNAs and lincRNAs acting as miRNA targets, which included 33 miRNAs and 34 lincRNAs, and 104 interactions between miRNAs and lincRNAs acting as miRNA decoys, which included 58 miRNAs and 86 lincRNAs. Moreover, 3714 mRNAs as 78 miRNA targets and 5490 mRNAs as 78 miRNA decoys are also shown. Interestingly, we found that the majority of nodes participated in other miRNA-regulated networks, but only three miRNAs including zma-miR529-5p, zma-miR399g-5p and zma-miR393c-5p:zma-miR393a-5p, formed separate sub-networks.
To further investigate the patterns of the miRNA-lincRNA-mRNA networks, we compared the number of four types of RNAs, including lincRNAs acting as miRNA targets, lincRNAs acting as miRNA decoys, mRNAs acting as miRNA targets and mRNAs acting as miRNA decoys, and found that the numbers of the four types were unevenly distributed for each miRNA. Additionally, the number of miRNA decoys was often greater than that of miRNA targets in most sub-networks, and only a small number of sub-networks had more miRNA targets than decoys (Fig. 5).
Furthermore, we found that miRNAs could bind to one or more lincRNAs (Fig. 6, Additional file 9). For example, Boerner_Z27kG1_07658 and Boerner_Z27kG1_17312 acted as decoys of zma-miR166n-5p, and Boerner_Z27kG1_01291 acted as a target of this miRNA. We also found that some lincRNAs could be bound by miRNAs from the same or different miRNA families. For example, Li_TCONS_00096947 and Li_TCONS_00064018 could be bound by zma-miR169n-3p and zma-miR169q-3p, and Boerner_Z27kG1_01046 could be bound by zma-miR408b-3p:zma-miR408a, zma-miR528a-3p:zma-miR528b-3p and zma-miR164b-3p (Fig. 6). Amazingly, the same lincRNA could be used as both a miRNA target and decoy using different binding sites in the lincRNAs. For example, Boerner_Z27kG1_08283 could be a target of zma-miR160b-3p:zma-miR160g-3p and zma-miR482-3p, and it could act as a decoy for zma-miR399e-5p (Fig. 6).
Functional prediction of lincRNAs acting as miRNA targets based on the lincRNA-mRNA co-expression network
To speculate on the functions of the 34 validated lincRNAs acting as miRNA targets, a co-expression network between lincRNAs and mRNAs was first constructed and then visualized (see materials and methods). The lincRNA-mRNA co-expression network was composed of 32 lincRNA nodes, 9043 mRNA nodes and 17968 edges (Fig. 7, Additional file 10), and we were able to infer that 32 lincRNAs could be co-expressed with 9043 mRNAs. In the network, we could see that one or more mRNAs were centered around lincRNAs and were connected to lincRNAs based on the Pearson correlation coefficient. Therefore, we could infer the function of each lincRNA based on the function of connected mRNAs. Through GO enrichment and functional analysis of the mRNAs that were co-expressed with lincRNAs, we found that lincRNAs mainly participate in cellular, metabolic and other biological processes, such as regulation of biological processes, metabolic processes, cellular processes, as well as in the response to stress (Fig. 8, Additional files 11 and 12). These lincRNAs were also highly enriched in cellular component terms including thylakoid and photosynthetic membrane (Additional files 11 and 13). Moreover, we found that the GO terms “hydrolase activity, acting on acid anhydrides”, “tetrapyrrole binding”, “iron ion binding” and “heme binding” were enriched in the “molecular function” category (Additional files 11 and 14).
Functional prediction of lincRNAs acting as miRNA decoys based on miRNA-lincRNA-mRNA networks
Based on the ceRNA hypothesis, which suggests that when lincRNAs acting as miRNA decoys and mRNAs are targeted by the same miRNAs, the function of the lincRNAs acting as miRNA decoys can be inferred from the mRNAs, we speculated on the function of 86 lincRNAs acting as 58 miRNA decoys. After using the AgriGO toolkit to perform GO analysis of mRNAs that could be targeted by the same miRNAs acting on lincRNAs , we found that lincRNAs acting as miRNA decoys were involved in multiple biological processes, participated in the formation of many cellular components, and influenced the activities of molecular functions (Fig. 9, Additional file 15). They were mainly involved in cellular and metabolic processes, and the molecular functions of lincRNAs acting as miRNA decoys were focused on catalytic activity and binding functions (Fig. 9).
To obtain a global function of lincRNAs acting as miRNA decoys, we performed enrichment analysis again, and found that these 86 lincRNAs may participate in diverse biological processes, such as cellular component organization; cellular component biogenesis; cellular processes; and metabolic processes, including macromolecular complex subunit organization, nucleosome assembly, DNA packing, superoxide metabolic processes, ribosome biogenesis, oxidation reduction, biosynthetic processes and translation (Additional file 16). They could also be involved in the formation of cells, macromolecular complexes, cell projections, cytoplasm, microtubules and protein-DNA complexes (Additional file 17). Moreover, these lincRNAs might modulate the effects of multiple molecular functions, including binding, structural molecular activity, transporter activity, catalytic activity, and electron carrier activity, and they may exhibit translation elongation factor activity, unfolded protein binding activity, monooxygenase activity, ammonia-lyase activity and GTPase activity (Additional file 18).
lincRNAs can be direct miRNA targets in maize
With their importance in regulating gene expression, lincRNAs have garnered significant attention in the life science field. Although increasing lincRNAs have been predicted and identified in plants [15–22], the relationship between miRNAs and lincRNAs have seldom been investigated by comparing the mRNAs as miRNA targets [19, 33, 34]. Recently, 51 lincRNAs were identified as putative targets of 30 miRNAs in Populus trichocarpa , but the evidence of lincRNAs acting as miRNA targets in plants are still lacking.
In plants, degradome sequencing is a new technology to identify and validate targets of miRNAs [68–72], and it has been used to directly validate miRNA targets in plants. Currently, only mRNAs as miRNA targets, but not lincRNAs as miRNA targets, have been validated by degradome data. Thus, using degradome data, we validated 34 lincRNAs as 33 miRNA targets, which indicates that, similar to mRNAs acting as miRNA targets, lincRNAs can also directly act as miRNA targets.
LincRNAs can also be miRNA decoys in maize
Functional target mimicries (miRNA decoys) were first studied in Arabidopsis ; consequently, computational methods have been used to identify miRNA decoys, but most of the identified miRNA decoys were protein-coding genes [52, 73, 74]. Only few studies were performed in ncRNAs as putative miRNA decoys [19, 41, 55], and no lincRNAs as miRNA decoys had previously been investigated in maize.
In our study, we found that a portion of lincRNAs could not be directly cleaved by the miRNA-associated silencing complex due to the existence of mismatches or large bulges at the 9th to 12th nucleotide positions of the miRNA-lincRNA pairing site. Using bioinformatics, we identified 86 lincRNAs acting as 58 miRNA decoys in maize and found that the miRNA decoy sites were conserved; however, most of the flanking regions of the miRNA decoy sites were not conserved. Our results indicate that lincRNAs acting as miRNA decoys widely exist in plants, which supports previously published data that lincRNAs as miRNA decoys could also be regulators of miRNA [19, 41].
The potential function of lincRNAs as miRNA targets or decoys
To investigate the function of maize lincRNAs acting as miRNA targets or decoys, two methods were used in this study: a co-expression network and the ceRNA hypothesis. The co-expression network, which is commonly used to predict gene function [64, 75, 76], was used to predict the function of lincRNAs as miRNA targets. By using the co-expression network, we predicted the function of 32 maize lincRNAs, and these lincRNAs were enriched in signaling processes, the regulation of biological processes, multicellular organismal processes, metabolic processes and immune system processes. Interestingly, these lincRNAs were enriched in multiple molecular functions, mainly in the catalytic activity and binding categories. Furthermore, when comparing with drought response lincRNAs previously reported, we found that three lincRNAs as miRNA targets in stress category were differentially expressed between the control and drought-stressed leaves (Additional file 19), which indicated that lincRNAs as miRNA targets may be involved in drought-stress .
The ceRNA hypothesis implies a network relationship between miRNAs, lincRNAs as miRNA decoys, and mRNA as miRNA targets; in these networks, lincRNAs could act as miRNA decoys, sequestering miRNAs and thereby favoring the expression of repressed mRNA targets [36, 77], and such networks can be used to predict the function of lincRNAs as miRNA decoys. Here, the functions of 86 lincRNAs acting as 56 miRNA decoys were predicted, and it was found that they can inhibit miRNA functions in a spatial- or temporal-specific manner, thus contributing to the regulation of transcript complexity in maize. Furthermore, when comparing the lincRNAs as miRNA decoys in the stress category using the previously reported drought response lincRNAs, 7 lincRNAs as miRNA decoys had been investigated previously and were differentially expressed between the control and drought-stressed leaves (Additional file 19), which indicated that lincRNAs associated with drought stress could potentially regulate miRNAs through lincRNAs as miRNA decoys.
Of the 1831 identified lincRNAs in maize, the number of lincRNAs that had the inferred function (34 lincRNAs as miRNA target, 86 lincRNAs as miRNA decoys) was still limited, which is consistent with the diverse mechanism of action of lincRNAs [15, 22]. We think that the lincRNAs as miRNA targets or miRNA decoys are just one type of lincRNAs, and we hope to investigate the function of other types of lincRNAs by using other methods, such as lincRNA-protein interaction prediction. In summary, our study lays a solid foundation for elucidating the regulatory mechanisms of miRNAs in maize and provides a source for exploring the function of lincRNAs in the future.
This study employed a computational pipeline for the systematic analysis of putative miRNA-lincRNA duplexes to better understand the role of lincRNAs. We found that 42 miRNA-lincRNA duplexes remained after filtering based on degradome evidence, and they were composed of 33 miRNAs and 34 lincRNAs that may be directly cleaved by miRNAs. Furthermore, 32 of the 34 lincRNAs could be co-expressed with mRNAs, and 86 lincRNAs were predicted as miRNA decoys that may competitively bind to miRNAs. According to the obtained co-expression networks and the ceRNA hypothesis, we effectively predicted the function of lincRNAs as miRNA targets or decoys. Future experimental studies are required to elucidate the mechanisms of miRNA-lincRNA duplexes and to reveal the functions of these lincRNAs in plants.
Availability of supporting data
The datasets supporting the results of this article are included within the article and its additional files.
Long intergenic noncoding RNAs
Long non-coding RNAs
Open reading frame
Competing endogenous RNAs
Gene expression omnibus
Minimum free energy
Ponjavic J, Ponting CP, Lunter G. Functionality or transcriptional noise? Evidence for selection within long noncoding RNAs. Genome Res. 2007;17(5):556–65.
Ponting CP, Oliver PL, Reik W. Evolution and functions of long noncoding RNAs. Cell. 2009;136(4):629–41.
Mercer TR, Dinger ME, Mattick JS. Long non-coding RNAs: insights into functions. Nat Rev Genet. 2009;10(3):155–9.
Ng SY, Johnson R, Stanton LW. Human long non-coding RNAs promote pluripotency and neuronal differentiation by association with chromatin modifiers and transcription factors. EMBO J. 2012;31(3):522–33.
Hawkins PG, Morris KV. Transcriptional regulation of Oct4 by a long non-coding RNA antisense to Oct4-pseudogene 5. Transcription. 2010;1(3):165–75.
Zhu QH, Wang MB. Molecular functions of long non-coding RNAs in plants. Genes (Basel). 2012;3(1):176–90.
Rinn JL, Chang HY. Genome regulation by long noncoding RNAs. Annu Rev Biochem. 2012;81:145–66.
Guttman M, Donaghey J, Carey BW, Garber M, Grenier JK, Munson G, et al. lincRNAs act in the circuitry controlling pluripotency and differentiation. Nature. 2011;477(7364):295–300.
Guttman M, Rinn JL. Modular regulatory principles of large non-coding RNAs. Nature. 2012;482(7385):339–46.
Khalil AM, Guttman M, Huarte M, Garber M, Raj A, Rivea Morales D, et al. Many human large intergenic noncoding RNAs associate with chromatin-modifying complexes and affect gene expression. Proc Natl Acad Sci U S A. 2009;106(28):11667–72.
Moran VA, Perera RJ, Khalil AM. Emerging functional and mechanistic paradigms of mammalian long non-coding RNAs. Nucleic Acids Res. 2012;40(14):6391–400.
Cabili MN, Trapnell C, Goff L, Koziol M, Tazon-Vega B, Regev A, et al. Integrative annotation of human large intergenic noncoding RNAs reveals global properties and specific subclasses. Genes Dev. 2011;25(18):1915–27.
Guttman M, Amit I, Garber M, French C, Lin MF, Feldser D, et al. Chromatin signature reveals over a thousand highly conserved large non-coding RNAs in mammals. Nature. 2009;458(7235):223–7.
Li T, Wang S, Wu R, Zhou X, Zhu D, Zhang Y. Identification of long non-protein coding RNAs in chicken skeletal muscle using next generation sequencing. Genomics. 2012;99(5):292–8.
Liu J, Jung C, Xu J, Wang H, Deng S, Bernad L, et al. Genome-wide analysis uncovers regulation of long intergenic noncoding RNAs in Arabidopsis. Plant Cell. 2012;24(11):4333–45.
Xin M, Wang Y, Yao Y, Song N, Hu Z, Qin D, et al. Identification and characterization of wheat long non-protein coding RNAs responsive to powdery mildew infection and heat stress by using microarray analysis and SBS sequencing. BMC Plant Biol. 2011;11:61.
Hao Z, Fan C, Cheng T, Su Y, Wei Q, Li G. Genome-wide identification, characterization and evolutionary analysis of long intergenic noncoding RNAs in cucumber. PLoS One. 2015;10(3):e0121800.
Qi X, Xie S, Liu Y, Yi F, Yu J. Genome-wide annotation of genes and noncoding RNAs of foxtail millet in response to simulated drought stress by deep sequencing. Plant Mol Biol. 2013;83(4–5):459–73.
Shuai P, Liang D, Tang S, Zhang Z, Ye CY, Su Y, et al. Genome-wide identification and functional prediction of novel and drought-responsive lincRNAs in Populus trichocarpa. J Exp Bot. 2014;65(17):4975–83.
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.
Boerner S, McGinnis KM. Computational identification and functional predictions of long noncoding RNA in Zea mays. PLoS One. 2012;7(8):e43047.
Zhang W, Han Z, Guo Q, Liu Y, Zheng Y, Wu F, et al. Identification of maize long non-coding RNAs responsive to drought stress. PLoS One. 2014;9(6):e98958.
Bartel DP. MicroRNAs: target recognition and regulatory functions. Cell. 2009;136(2):215–33.
Seitz H. Redefining microRNA targets. Curr Biol. 2009;19(10):870–3.
Ding Y, Tao Y, Zhu C. Emerging roles of microRNAs in the mediation of drought stress response in plants. J Exp Bot. 2013;64(11):3077–86.
Chiou TJ. The role of microRNAs in sensing nutrient stress. Plant Cell Environ. 2007;30(3):323–32.
Sunkar R, Chinnusamy V, Zhu J, Zhu JK. Small RNAs as big players in plant abiotic stress responses and nutrient deprivation. Trends Plant Sci. 2007;12(7):301–9.
Aukerman MJ, Sakai H. Regulation of flowering time and floral organ identity by a MicroRNA and its APETALA2-like target genes. Plant Cell. 2003;15(11):2730–41.
Palatnik JF, Allen E, Wu X, Schommer C, Schwab R, Carrington JC, et al. Control of leaf morphogenesis by microRNAs. Nature. 2003;425(6955):257–63.
Fujii H, Chiou TJ, Lin SI, Aung K, Zhu JK. A miRNA involved in phosphate-starvation response in Arabidopsis. Curr Biol. 2005;15(22):2038–43.
Bartel DP. MicroRNAs: genomics, biogenesis, mechanism, and function. Cell. 2004;116(2):281–97.
Rhoades MW, Reinhart BJ, Lim LP, Burge CB, Bartel B, Bartel DP. Prediction of plant microRNA targets. Cell. 2002;110(4):513–20.
Jalali S, Bhartiya D, Lalwani MK, Sivasubbu S, Scaria V. Systematic transcriptome wide analysis of lncRNA-miRNA interactions. PLoS One. 2013;8(2):e53823.
Juan L, Wang G, Radovich M, Schneider BP, Clare SE, Wang Y, et al. Potential roles of microRNAs in regulating long intergenic noncoding RNAs. BMC Med Genomics. 2013;6 Suppl 1:S7.
Liang H, Zhang J, Zen K, Zhang CY, Chen X. Nuclear microRNAs and their unconventional role in regulating non-coding RNAs. Protein Cell. 2013;4(5):325–30.
Salmena L, Poliseno L, Tay Y, Kats L, Pandolfi PP. A ceRNA hypothesis: the Rosetta Stone of a hidden RNA language? Cell. 2011;146(3):353–8.
Rubio-Somoza I, Weigel D, Franco-Zorilla JM, Garcia JA, Paz-Ares J. ceRNAs: miRNA target mimic mimics. Cell. 2011;147(7):1431–2.
Cesana M, Cacchiarelli D, Legnini I, Santini T, Sthandier O, Chinappi M, et al. A long noncoding RNA controls muscle differentiation by functioning as a competing endogenous RNA. Cell. 2011;147(2):358–69.
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.
Zhang YC, Liao JY, Li ZY, Yu Y, Zhang JP, Li QF, et al. Genome-wide screening and functional analysis identify a large number of long noncoding RNAs involved in the sexual reproduction of rice. Genome Biol. 2014;15(12):512.
Wu HJ, Wang ZM, Wang M, Wang XJ. Widespread long noncoding RNAs as endogenous target mimics for microRNAs in plants. Plant Physiol. 2013;161(4):1875–84.
Bennetzen JL, Hake S. Handbook of maize. New York: Springer;2009.
Griffiths-Jones S, Saini HK, van Dongen S, Enright AJ. miRBase: tools for microRNA genomics. Nucleic Acids Res. 2008;36(Database issue):D154–158.
Kozomara A, Griffiths-Jones S. miRBase: integrating microRNA annotation and deep-sequencing data. Nucleic Acids Res. 2011;39(Database issue):D152–157.
Zhao Y, Xu Z, Mo Q, Zou C, Li W, Xu Y, et al. Combined small RNA and degradome sequencing reveals novel miRNAs and their targets in response to low nitrate availability in maize. Ann Bot. 2013;112(3):633–42.
Edgar R, Domrachev M, Lash AE. Gene expression omnibus: NCBI gene expression and hybridization array data repository. Nucleic Acids Res. 2002;30(1):207–10.
Barrett T, Troup DB, Wilhite SE, Ledoux P, Evangelista C, Kim IF, et al. NCBI GEO: archive for functional genomics data sets--10 years on. Nucleic Acids Res. 2011;39(Database issue):D1005–1010.
Kruger J, Rehmsmeier M. RNAhybrid: microRNA target prediction easy, fast and flexible. Nucleic Acids Res. 2006;34(Web Server issue):W451–454.
Rehmsmeier M, Steffen P, Hochsmann M, Giegerich R. Fast and effective prediction of microRNA/target duplexes. RNA. 2004;10(10):1507–17.
Tafer H, Hofacker IL. RNAplex: a fast tool for RNA-RNA interaction search. Bioinformatics. 2008;24(22):2657–63.
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.
Ivashuta S, Banks IR, Wiggins BE, Zhang Y, Ziegler TE, Roberts JK, et al. Regulation of gene expression in plants through miRNA inactivation. PLoS One. 2011;6(6):e21330.
Goodstein DM, Shu S, Howson R, Neupane R, Hayes RD, Fazo J, et al. Phytozome: a comparative platform for green plant genomics. Nucleic Acids Res. 2012;40(Database issue):D1178–1186.
Mount DW. Using the Basic Local Alignment Search Tool (BLAST). CSH Protoc. 2007;2007:pdb.top17.
Ye CY, Xu H, Shen E, Liu Y, Wang Y, Shen Y, et al. Genome-wide identification of non-coding RNAs interacted with microRNAs in soybean. Front Plant Sci. 2014;5:743.
Katoh K, Misawa K, Kuma K, Miyata T. MAFFT: a novel method for rapid multiple sequence alignment based on fast Fourier transform. Nucleic Acids Res. 2002;30(14):3059–66.
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.
Li P, Ponnala L, Gandotra N, Wang L, Si Y, Tausta SL, et al. The developmental dynamics of the maize leaf transcriptome. Nat Genet. 2010;42(12):1060–7.
Davidson RM, Hansey CN, Gowda M, Childs KL, Lin H, Vaillancourt B, et al. Utility of RNA sequencing for analysis of maize reproductive transcriptomes. Plant Genome. 2011;4(3):191–203.
Chang YM, Liu WY, Shih AC, Shen MN, Lu CH, Lu MY, et al. Characterizing regulatory and functional differentiation between maize mesophyll and bundle sheath cells by transcriptomic analysis. Plant Physiol. 2012;160(1):165–77.
Bolduc N, Yilmaz A, Mejia-Guerra MK, Morohashi K, O’Connor D, Grotewold E, et al. Unraveling the KNOTTED1 regulatory network in maize meristems. Genes Dev. 2012;26(15):1685–90.
Paschold A, Jia Y, Marcon C, Lund S, Larson NB, Yeh CT, et al. Complementation contributes to transcriptome complexity in maize (Zea mays L.) hybrids relative to their inbred parents. Genome Res. 2012;22(12):2445–54.
Waters AJ, Makarevitch I, Eichten SR, Swanson-Wagner RA, Yeh CT, Xu W, et al. Parent-of-origin effects on gene expression and DNA methylation in the maize endosperm. Plant Cell. 2011;23(12):4221–33.
Liao Q, Liu C, Yuan X, Kang S, Miao R, Xiao H, et al. Large-scale prediction of long non-coding RNA functions in a coding-non-coding gene co-expression network. Nucleic Acids Res. 2011;39(9):3864–78.
Langfelder P, Horvath S. WGCNA: an R package for weighted correlation network analysis. BMC Bioinformatics. 2008;9:559.
Du Z, Zhou X, Ling Y, Zhang Z, Su Z. agriGO: a GO analysis toolkit for the agricultural community. Nucleic Acids Res. 2010;38(Web Server issue):W64–70.
Todesco M, Rubio-Somoza I, Paz-Ares J, Weigel D. A collection of target mimics for comprehensive analysis of microRNA function in Arabidopsis thaliana. PLoS Genet. 2010;6(7):e1001031.
Addo-Quaye C, Eshoo TW, Bartel DP, Axtell MJ. Endogenous siRNA and miRNA targets identified by sequencing of the Arabidopsis degradome. Curr Biol. 2008;18(10):758–62.
Guo W, Zhang Y, Wang Q, Zhan Y, Zhu G, Yu Q, et al. High-throughput sequencing and degradome analysis reveal neutral evolution of Cercis gigantea microRNAs and their targets. Planta. 2015; doi:10.1007/s00425-015-2389-y.
Li YF, Zheng Y, Addo-Quaye C, Zhang L, Saini A, Jagadeeswaran G, et al. Transcriptome-wide identification of microRNA targets in rice. Plant J. 2010;62(5):742–59.
Xu T, Wang Y, Liu X, Lv S, Feng C, Qi M, et al. Small RNA and degradome sequencing reveals microRNAs and their targets involved in tomato pedicel abscission. Planta. 2015;242(4):963–84.
Liu H, Qin C, Chen Z, Zuo T, Yang X, Zhou H, et al. Identification of miRNAs and their target genes in developing maize ears by combined small RNA and degradome sequencing. BMC Genomics. 2014;15:25.
Meng Y, Shao C, Wang H, Jin Y. Target mimics: an embedded layer of microRNA-involved gene regulatory networks in plants. BMC Genomics. 2012;13:197.
Banks IR, Zhang Y, Wiggins BE, Heck GR, Ivashuta S. RNA decoys: an emerging component of plant regulatory networks? Plant Signal Behav. 2012;7(9):1188–93.
Liao Q, Shen J, Liu J, Sun X, Zhao G, Chang Y, et al. Genome-wide identification and functional annotation of Plasmodium falciparum long noncoding RNAs from RNA-seq data. Parasitol Res. 2014;113(4):1269–81.
Zhao Y, Luo H, Chen X, Xiao Y, Chen R. Computational methods to predict long noncoding RNA functions based on co-expression network. Methods Mol Biol. 2014;1182:209–18.
Kartha RV, Subramanian S. Competing endogenous RNAs (ceRNAs): new entrants to the intricacies of gene regulation. Front Genet. 2014;5:8.
This work was supported by grants from the National Natural Science Foundation of China (Grant No.31070256, No.31370329), the Program for New Century Excellent Talents in University [NCET-12-0896], the Co-Innovation Center for Qinba Region’s Sustainable Development (CIC-QBRSD), the Natural Science Basic Research Plan of Shaanxi Province, China (Program No.2014JM3074), and the Fundamental Research Funds for the Central Universities (No. GK201403004). The funding agencies had no role in the study, its design, the data collection and analysis, the decision to publish, or the preparation of the manuscript.
The authors declare that they have no competing interests.
CYF and GLL designed the pipeline and drafted the manuscript. CYF, ZQH, and JHY performed the data analysis. All of the authors critically revised and provided final approval of this manuscript.
LincRNA information derived from three articles. (XLS 20 kb)
Summary of degradome reads in Zea mays . (XLS 20 kb)
LincRNAs predicted as putative miRNA targets. (XLS 344 kb)
The number of degradome reads that are perfectly matched with lincRNAs identified as miRNA targets. (XLS 2627 kb)
The sequence logos of the 12 conserved lincRNAs as miRNA targets. (ZIP 3605 kb)
LincRNAs predicted as miRNA decoys. (XLS 71 kb)
The sequence logos of the 10 conserved lincRNA as miRNA decoys. (ZIP 1503 kb)
The nodes and edges information of miRNA-lincRNA-mRNA networks. (XLS 755 kb)
The nodes and edges information of miRNA-lincRNA networks. (XLS 31 kb)
The nodes and edges information of lincRNA-mRNA co-expression networks. (XLS 1090 kb)
The potential function of lincRNAs as miRNA targets. (XLS 5756 kb)
GO enrichment analysis of maize lincRNAs functioning as miRNA targets associated with “biological processes”. (PDF 153 kb)
GO enrichment analysis of maize lincRNAs functioning as miRNA targets associated with “cellular components”. (PDF 50 kb)
GO enrichment analysis of maize lincRNAs functioning as miRNA targets associated with “molecular functions”. (PDF 65 kb)
The potential function of lincRNAs as miRNA decoys. (XLS 1841 kb)
GO enrichment analysis of maize lincRNAs functioning as miRNA decoys associated with “biological processes”. (PDF 1350 kb)
GO enrichment analysis of maize lincRNAs functioning as miRNA decoys associated with “cellular components”. (PDF 773 kb)
GO enrichment analysis of maize lincRNAs functioning as miRNA decoys associated with “molecular functions”. (PDF 1498 kb)
The expression of drought-stressed lincRNAs as miRNA targets or decoys. (XLS 22 kb)