PresRAT: a server for identification of bacterial small-RNA sequences and their targets with probable binding region

ABSTRACT

Bacterial small-RNA (sRNA) sequences are functional RNAs, which play an important role in regulating the expression of a diverse class of genes. It is thus critical to identify such sRNA sequences and their probable mRNA targets. Here, we discuss new procedures to identify and characterize sRNA and their targets via the introduction of an integrated online platform ‘PresRAT’. PresRAT uses the primary and secondary structural attributes of sRNA sequences to predict sRNA from a given sequence or bacterial genome. PresRAT also finds probable target mRNAs of sRNA sequences from a given bacterial chromosome and further concentrates on the identification of the probable sRNA-mRNA binding regions. Using PresRAT, we have identified a total of 66,209 potential sRNA sequences from 292 bacterial genomes and 2247 potential targets from 13 bacterial genomes. We have also implemented a protocol to build and refine 3D models of sRNA and sRNA-mRNA duplex regions and generated 3D models of 50 known sRNAs and 81 sRNA-mRNA duplexes using this platform. Along with the server part, PresRAT also contains a database section, which enlists the predicted sRNA sequences, sRNA targets, and their corresponding 3D models with structural dynamics information.

KEYWORDS:

• sRNA identification
• sRNA target prediction
• RNA structure
• RNA-protein interaction
• RNA-RNA interaction

Introduction

Small-RNA (sRNA) based regulation of mRNA sequences play a central role in various bacterial gene expressions in response to environmental changes [1,2]. These sRNA sequences are extensively studied in context with different bacterial regulatory processes such as quorum sensing [3], the formation of biofilms [4], expression of diverse virulence factors [5], iron uptake [6], etc. In the prokaryotic system, these functional RNA regulators exist as small transcripts about 30–600 nucleotides long and modulate the translation and stability of mRNAs by binding with them [7]. Binding regions of these sRNA sequences with their cognate mRNAs are not fixed. Some sRNA sequences primarily bind at the ribosome binding site (RBS) of their target mRNA sequences [1] while other sRNA sequences modulate translation of target genes by attaching at distant locations [1,7]. The antisense mechanisms exhibited by the sRNA sequences to modulate their targets are the key to understand the function of these RNA regulators [1,2,7]. The significance of these sRNA sequences and their targeted genes in bacteria has driven the research community to identify and more of these sRNA sequences and their targets concerning various biological functions. Different experimental methods such as RNA-seq, tiling arrays, sRNA pulse expression, immunoprecipitation are available to detect sRNA sequences and their targets. Still, they are tedious, expensive, and the outputs of these experiments also need exhaustive filtering [8–10]. Hence, computational approaches can augment the experimental methods [11] in the detection of sRNA sequences and their targets.

Comparative genomics are a standard computational practice to find sRNA sequences within a bacterial genome [12]. QRNA [13] and Intergenic Sequence Inspector [14] programs implemented a comparative genomics approach to find sRNAs. Programs like RNAz [15] and sRNApredict [16] use the conserved RNA secondary structure information to calculate the thermodynamic energy for sRNA prediction. GLASSgo [17] combined an iterative BLAST strategy with pairwise identity filtering and structure-based clustering to find sRNA homologs from scratch. Some of the programs also utilize the promoter and terminator sequence signal annotations for sRNA prediction. sRNAscanner [18] applies a position weighted matrix of sRNA specific transcriptional signals, and sRNAfinder [19] utilizes a generalized probabilistic model that integrates primary sequence data, transcript expression data, microarray experiments, and conserved RNA structure information to predict sRNA genes. Diverse sets of methods also utilize different machine learning approaches [20,21], antisense RNA transcript expression profiling [22], and motif-based discovery [23] algorithms to predict sRNA sequences.

Identification of sRNA target genes is also important and represents a challenging task. Experimental procedures, including genetic screens [24], knockouts studies [25], microarray analysis [26], and other systematic proteomic investigations [27], are implemented to identify sRNA target genes. However, these methods are expensive, time-consuming, and extremely tedious. Therefore, the development of sensitive and accurate computational approaches would aid and complement the experimental analysis. The first step in sRNA target prediction is to search for a complementary region between the sRNA and its target mRNA sequence. Comparing to microRNA and its target mRNA interaction, which requires a generally conserved 6–8 nucleotide long seed region for initiation of binding, the sRNA-mRNA binding region is relatively more extended and more evolutionary flexible [28]. Non-Watson-Crick base-pairing plays a key role in sRNA-mRNA target interaction and programs such as TargetRNA [26], GUUGle [29] utilizes both Watson-Crick and non-Watson-Crick base pairing models to find sRNA targets and their interacting regions. RNAcofold [30] computes the hybridization energy and base-pairing arrangement between a pair of RNA molecules, while RNAplex [31] calculates the same hybridization energy using modified thermodynamic energy models. RNAup [32] utilizes the thermodynamics of RNA-RNA interactions and the hybridization energy, while IntaRNA [33] incorporates binding site accessibility information for searching the best interaction with minimum extended hybridization energy. SPOT [34] incorporated existing computational tools to search for sRNA binding sites and utilizes experimental data for filtering.

Despite the availability of various sRNA and its target prediction algorithms, the sensitivity remains quite low (20–49%) [11]. Hence, it is fair to assume that many more sRNA sequences and their target genes remain undiscovered, particularly in less studied prokaryotic organisms. Here, in this article, we introduce an integrated platform with an alternative approach to identify bacterial small-RNA and their target genes from a given nucleotide sequence and a whole bacterial chromosome. This online platform, PresRAT (Prediction of sRNA and their Targets) integrates two separate but related procedures of sRNA and their target identification into a simple and effective single dais.

Apart from identifying the sRNA sequences and their targets, a key feature of this platform is to provide sRNA and sRNA-mRNA duplex structures generated via three dimensional (3D) modelling and further refined through molecular dynamics simulation. Predicting the accurate 3D structure of RNA from the primary sequence is still a fundamental and challenging problem in computational biology [35]. There are few notable programs such as iFoldRNA [36], FARNA/FARFAR [37,38], MC-Sym [39], MANIP [40], 3DRNA [41], RNA2D3D [42], NAST [43], and ASSEMBLE [44] try to tackle the existing bottlenecks of RNA structure prediction and improve the algorithms for RNA 3D structure generation. As per the latest RCSB protein databank [45] statistics, there are 1528 experimentally solved 3D structures of RNA representing only ~0.01% of the whole deposited structures. However, there are no experimentally solved 3D structures available for sRNA only apart from a few conjugated with chaperone proteins [46–50]. PresRAT uses the RNA2D3D [42] program to generate an initial RNA model followed by extensive refinement of the structure using GROMACS v4.6.3 molecular dynamics simulation package [51].

The validation and benchmarking results of the PresRAT for sRNA, sRNA-target, and sRNA-mRNA binding regions suggest improvements compared to those achieved by other available programs. In this current study, we have studied 292 bacterial genomes to identify a total of 21,286 and 44,923 sRNA sequences from genic and non-genic regions of the respective genomes. We also have identified 2447 potential targets for 50 unique sRNA sequences from 13 bacterial genomes. Further, we have shown the structural details of 50 sRNA and 81 sRNA-mRNA duplex 3D models through molecular dynamics simulation studies. All these procedures and information are automated and integrated into the PresRAT platform freely available at http://www.hpppi.iicb.res.in/presrat.

Results

PresRAT: the webserver interfaces and available programs

PresRAT webserver and the database sections are built using CGI-Perl scripts. Integrated CGI-Perl scripts allow the users to utilize this webserver efficiently. The webserver is freely available at http://www.hpppi.iicb.res.in/presrat. The PresRAT webserver section is divided into four parts to accomplish different tasks. The first two parts are called Pred^sRNA and Pred^GsRNA. Pred^sRNA interface can be used to identify sRNA sequences one at a time while the Pred^GsRNA interface directs the user to identify potential sRNA sequences from different bacterial sequences. The next two parts of the PresRAT webserver compromises the sRNA target identification protocols. The Pred^Tar Interface helps to identify interacting nucleotides between a pair of sRNA and mRNA sequence, and Pred^WG performs the function of sRNA target identification from a given bacterial genome. The following subsections briefly describe the methodology and the results of each of these parts of the PresRAT server.

Pred^sRNA: prediction of sRNA from its sequence and secondary structural attributes

Pred^sRNA uses sequence and secondary structural information of existing sRNA and non-sRNA sequences to calculate a combined score to predict novel sRNA sequences. In sequence analysis, a directional (5ʹ->3ʹ) dinucleotide Sequence_score is first calculated for the input nucleotide sequence using Log Odds (LOD) ratio matrices (Please refer the detailed sRNA prediction methodology in Supplementary File 1 and Figure S1). We also calculate a U-richness score that reflects the presence of uracil (U) nucleotides at the 3ʹ end of the sRNA sequence. PresRAT effectively uses the secondary structure information by first calculating a set of local minima conformations of the provided RNA sequence, followed by estimating the two different scores, namely ABE_score and ALE_score. ABE_score and ALE_score indicate the Average Base Energy and Average Loop Energy of local minima secondary structure conformations. The final combined score i.e. sRNA_score, is calculated as follows:

(1) $sRN{A}_{score}=(Sequenc{e}_{score}+Uric{h}_{score}-AB{E}_{score}-AL{E}_{score})$(1)

The combined sRNA_score (details of the scoring schemes can be found in Supplementary File 1) is used as the final scoring scheme for filtering the true positive sRNA instances from the true negative ones. Equally weighted bacterial small-RNA (sRNA) feature set of Sequence_score, Urich_score, ABE_score, and ALE_score described in the Supplementary File 1 offers a simple but powerful discriminating scoring scheme to filter out sRNA sequences from non-sRNA sequences. All the individual feature sets can differentiate an sRNA from a non-sRNA sequence, as demonstrated in Figure S2.

For benchmarking and validation of sRNA prediction by PresRAT, a set of non-redundant 1174 sRNA (Additional file1), 5870 genic, and 5870 non-genic non-sRNA sequences were randomly collected from 54 different bacterial species (Additional file 2). Receiver operating characteristic (ROC) statistics were calculated by comparing the sRNA to non-sRNA sRNA_score in a 1:5 ratio. The ROC statistics generated using the sRNA_score distribution demonstrate sensitivity values of 27.7%, 49.6%, 65.0% and 77.5% at 1.0% 5.0%, 10.0%, and 20.0% error rates, against non-sRNA sequences collected from genic regions whereas sensitivities of 60.3%, 64.6%, 70.1%, and 78.7% were obtained at 1.0% 5.0%, 10.0%, and 20.0% error rates against non-sRNA sequences collected from non-genic parts of the bacterial chromosomes (Fig. 1A and Additional file 3), respectively. A comparative study against RNAz [15] program with Pred^sRNA also showed better performance in predicting sRNA sequences as compared to genic non-sRNA sequences (Fig. 1B-C). However, it is also important to mention that RNAz requires evolutionary information to predict functional RNA sequences, and thus it is restricted only to homologous RNA sequences in bacterial chromosomes. Therefore, the final comparison was performed on a limited set of commonly predicted sRNA (334 sequences), genic (1670), and non-genic (1670) non-sRNA sequences by the two programs.

Figure 1. Receiver operating characteristic curve for sRNA prediction. Panel A shows the ROC curve for sRNA prediction performance by the Pred^sRNA program in genic and non-genic category. A total of 1174 sRNA, 5870 genic, and 5870 non-genic non-sRNA sequences were considered to measure the performance. Panel B and C shows the performance of sRNA prediction by Pred^sRNA and RNAz programs in genic and non-genic categories, respectively. A total of common 334 sRNA, 1670 genic, and 1670 non-genic non-sRNA sequences were considered to measure the performance

Pred^sRNA page accepts query nucleotide sequence and based on its sRNA_score and predicts whether the given sequence is likely to be sRNA or not above a user-defined confidence level. Pred^sRNA result page provides a minimum free energy secondary structure of the probable sRNA sequence. Pred^sRNA also provides a useful option of building the three-dimensional (3D) structural model of the predicted sRNA sequence via the Pred^sRNA3D module, which generates 3D sRNA models using RNA2D3D [42] and RNAfold [525354] programs, respectively. Following the generation of a 3D sRNA structure, a molecular dynamics (MD) simulation is performed on the RNA structure using the GROMACS [51] MD simulation package (details of the structural analysis protocol can be found at Supplementary File 1). The initial and final structures are provided from the MD simulation trajectory file along with other necessary information.

Pred^GsRNA: prediction of sRNA sequences within the bacterial chromosome

Pred^GsRNA aims to predict potential sRNA sequences within a bacterial genome. Pred^GsRNA applies a heuristic approach to find potential sRNA sequences from the bacterial genome. Briefly, Pred^GsRNA first extracts the 300 nucleotides starting 5ʹ and ending 3ʹ regions of both the genic and non-genic sequences (Figure S3 in Supplementary File 1), followed by a filtering step to exclude sequences containing an exclusive set of 8mer nucleotide sequences (Figure S4 in Supplementary File 1). The retained sequences are then subjected to sRNA_scorecalculation. The result page of Pred^GsRNA contains the sRNA sequences predicted from the genic sequences and non-genic sequences of the bacterial genome, respectively. The corresponding result page of this program consists of the top 100 predicted sRNA sequences along with their start and end loci, sRNA_score, and p-values. Similarly, a pictorial representation of the predicted sRNA sequences on a genome browser view is also provided for better visualization (please refer to the Supplementary File 1 for details of this approach).

Pred^GsRNA is used to predict new probable sRNA sequences from the 292 bacterial genomes (Additional file 4) to identify a total of 21,286 and 44,923 hypothetical sRNA sequences from genic and non-genic parts of the bacterial genomes (Additional file 5). sRNAscanner [18], a separate program for bacterial chromosomal sRNA detection, is used for comparative analysis. sRNAscanner uses a transcriptional signal based computational method for sRNA detection and detected only 114 unique sRNA sequences. In comparison, our protocol identified a total of 219 sRNA sequences out of 1174 known sRNA sequences (Additional file 6 and Figure S5). Hence, Pred^GsRNA demonstrates ~2 fold higher sensitivity over existing program in identification of known sRNA sequences from the whole bacterial genome. However, Pred^GsRNA is a greedy and heuristic approach to search sRNA sequences from bacterial genomes. The method is restricted to search the sRNA sequences within the 300bps of genic and non-genic feature. The users are encouraged to use other methods like sRNAscanner to complement sRNA sequence searching from a bacterial genome.

Pred^TAR: prediction of sRNA-mRNA target binding regions

Pred^TAR aims to predict potential binding regions between bacterial sRNA and their target mRNA binding regions. In short, Pred^TAR combines two different approaches to find the sRNA-mRNA target binding region. The first approach (Component I) uses the base un-pairing probability of both the sRNA and its target mRNA sequence (Fig. 2A). The second approach (Component II) (Fig. 2B) utilizes the local minima energy structures of sRNA and the minimum free energy structure of the target mRNA to find the potential binding region (Please refer to Supplementary File 1 and Figure S6 for details). The webserver output page provides the alignment of the given sRNA and mRNA sequences along with scores and duplex energy values.

Figure 2. Prediction of sRNA targets. Panel A demonstrates the sRNA-mRNA target binding region identification protocol (Component I) of Pred^TAR. The main aim of this protocol is to identify the un-pairing probability of each base from both sRNA and mRNA sequences and utilize this probability information to identify the binding region. Panel B describes the Component II protocol used by PresRAT to identify the sRNA-mRNA target binding region. In this protocol, the energetic value of both sRNA and mRNA secondary structures are taken into account to identify binding regions

Experimentally known binding regions in 88 out of 91 sRNA and their respective target mRNA sequences (Additional file 7) were used to validate the Component I and II algorithms. The selection of these 88 binding regions was made based on the frequency distribution of the sRNA binding site on mRNA sequences relative to its translation initiation site (Fig. 3A). The frequency distribution of the 91 sRNA binding sites on mRNA sequences reveals that in 97% of the cases, the sRNA binding site is located within the −250 to +100 window region of a target mRNA (‘0’ being the translation start site). The Component I and Component II could correctly identify the binding region for 49% and 67% of the true sRNA-mRNA pair within the window mentioned above. The Venn diagram in Fig. 3B shows the overlap of correctly predicted sRNA-mRNA binding regions by the two algorithms. A 74% (65 out of 88) sensitivity was achieved when a region with the lowest free energy of sRNA-mRNA binding by either of the algorithms was considered, while 82% (72 out of 88) sensitivity (Fig. 3C) was achieved when both the common and uniquely predicted binding regions were considered. A comparative analysis (Fig. 3C) with other published methods (RNAplex [31], RNAup [32], and IntaRNA [33], respectively) shows a much better performance of Pred^TAR in predicting sRNA-mRNA targets (a complete list of predicted sRNA-mRNA binding regions is provided in Additional file 8).

Figure 3. sRNA target finding performance. Panel A shows the frequency distribution of sRNA binding sites in the 91 mRNA sequences. In 98% of the cases, it is found that the sRNA binding site is present within the window of −250 to +100 nucleotides (0 being the translation initiation site) of the target gene. Panel B shows the overlap of correctly predicted target by the PresRAT Component I and II protocols, respectively. Panel C shows the PresRAT sensitivity in predicting the correct sRNA-mRNA target binding regions in 88 samples compared to other programs. Panel D shows the number of true sRNA target gene identification through whole-genome search by PresRAT compared to other existing programs

Pred^TAR also provides a useful option of building the three-dimensional structural (3D) model of the sRNA-mRNA bound duplex region via the Pred^dup3D module, which generates the duplex model using RNA2D3D [42] program for the secondary structure information followed by molecular dynamics (MD) simulation. Details of the structural analysis can be found in Supplementary File 1. The initial and final structures are provided from the MD simulation trajectory file along with other necessary information.

Pred^WG: prediction of sRNA target genes and binding regions within a bacterial genome

Pred^WG aims to predict potential bacterial sRNA targets within a bacterial genome. The details about this approach are documented in Supplementary File 1. The result page of Pred^WG contains the list of the predicted target genes along with individual alignment and the corresponding duplex energy value of the binding region. Similarly, a pictorial representation of the predicted sRNA sequences on a genome browser view is also provided for better visualization.

A performance evaluation using 91 known sRNA-mRNA target pairs, Pred^WG achieved 86% sensitivity in predicting the correct target genes (within the top 100 predictions/genome) within bacterial genomes. The other competing programs, like IntaRNA [33], RNAPredator [55], TargetRNA2 [56], and RNAplex [31], achieved sensitivity values of 38%, 23%, 22%, and 14%, respectively (Fig. 3D) (complete list of targets found by the various programs are provided in Additional file 9). For a fair comparison, a similar search window (−250 to +100 of mRNA translation start site) length was maintained for all the programs. Within the top 100 predictions for individual sRNA sequences, Pred^WG predicted 2447 novel mRNA targets with a similar binding energy distribution profile (p values < 0.05) (Figure S7) as that of the known 91 sRNA-mRNA interactions. The complete list of the identified target mRNAs for 50 sRNA sequences along with their corresponding functions, is provided in the Additional file 10. We also performed a Gene Ontology (GO) [57] based molecular functional network of E. coli known sRNA sequences (18 sequences) with their known (27 unique targets) and predicted mRNA target sequences (192) (Figure S8A). The molecular functional network suggests the involvement of functions such as binding, catalysis, and transporter activities. Additionally, a smaller subset of the predicted target genes also offers unique molecular functions like electron transport, guanyl-nucleotide exchange, and other regulator activities (Figure S8B).

PresRAT: the database interface

PresRAT database section enlists the predicted sRNA sequences, sRNA targets, and the corresponding 3D models of selected sRNAs and sRNA-mRNA duplex structures. Predicted sRNA sequences from 292 bacterial genomes are accessible to users from a drop-down selection menu provided on the database page. Similarly, the list of predicted target genes for a given sRNA within a bacterial genome is readily accessible in a similar manner. The database also harbours mRNA targets of 50 known sRNA sequences from 13 bacterial genomes along with their details.

Exhaustive structural analysis was performed to predict the 3D structures of the known sRNA sequences and their corresponding binding regions with the target genes. Three-dimensional modelling of sRNA sequences and the target binding regions using their current secondary structure information provides the opportunity to gain structural insight, which in turn may become quite valuable in RNA based therapeutics development. There are 50 sRNA models in the PresRAT database. All the 50 sRNA and 81 unique sRNA-mRNA bound duplex models are explicitly refined using the GROMACS package through one-nanosecond of molecular simulation runs. The final models with other structural details are available to users in the database section. Further methodical details and validation results can be found in Supplementary File 1.

Discussion

Here, we introduce a web-based platform PresRAT to identify and analyse bacterial small-RNAs using both individual and genomic sequence level searches. The platform integrates a combination of scoring schemes developed based on primary and secondary structural attributes of small-RNAs. Using PresRAT, we have scanned a total of 292 bacterial genomes to identify 21,286 and 44,923 potential sRNA sequences. One of the key features of this platform is the generation of sRNA three-dimensional models and their dynamic analysis. For in-depth structural analysis, the availability of energy minimized and refined three-dimensional models of sRNAs are essential. PresRAT provides an option for automated model building and molecular dynamics simulation of sRNA in the presence of water. Using this automated procedure, we have successfully built 50 different small-RNA structures and analysed some of their basic properties calculated from their dynamics data (Figure S9 in Supplementary File 1). Comparison of our server derived model of DsrA and the NMR structure of the SL1 of the DsrA (pdb code: 5WQ1) [46] shows high similarities (RMSD: 2.75 Å) between the structures, indicating reliability of our 3D model (Figure S10 in Supplementary File 1). However, the NMR structure covers only a short fraction of the DsrA sequence.

A set of 81 experimentally known sRNA-mRNA duplex binding regions was modelled using the Pred^dup3D module and the resultant models were further stabilized via MD simulations (Figure S11 in Supplementary File 1). We also investigated the influence of Hfq on four different sRNA sequences using this platform (please refer to the Supplementary File 1 and Figure S12 for more details). In corroboration with previous experimental results, our molecular dynamics simulations demonstrated a major influence of Hfq binding in the regulation of the flexibility and conformations of the bound and unbound sRNA sequences. PresRAT enables a user to search for probable target mRNAs of sRNA sequences from a given bacterial chromosome and further concentrate on the identification of the probable mRNA target binding regions. For identification of target binding regions, PresRAT uses base un-pairing probability information of sRNA-mRNA sequences and local minima of sRNA secondary structures. Both the approaches combined with other additional steps provide a simple but effective way of identification of target mRNA sequences from the whole bacterial chromosomes. In this way, we have scanned 13 bacterial genomes and identified 2447 novel mRNA target sequences for existing 50 sRNA sequences. The identified target mRNAs for a given sRNA are provided in the database section of PresRAT along with their interacting nucleotide region and their related functional information. Existing model generation and molecular dynamics simulation procedure are extended to build the energy minimized and refined three-dimensional models of the sRNA-mRNA duplex regions, along with their structural characterization. Our study revealed that in general, these duplex regions are thermodynamically unfavourable, which indicates the involvement of other critical factors for the generation of stable sRNA-mRNA complexes.

Author contributions

A.C. collected and organized the data, developed the program, webserver and drafted the manuscript. K.K collected data, worked on the web server and drafted the manuscript. S. C. drafted the manuscript and coordinated the project. All authors read and approved the final manuscript.

Acknowledgments

The Council of Scientific and Industrial Research (CSIR), Government of India supported the study through funding grant BSC0121. K.K. acknowledges Department of Biotechnology (DBT), Goverment of India for fellowship.

Disclosure statement

The author(s) declare no competing financial interests.

Supplementary material

Supplemental data for this article can be accessed here.

Additional information

Funding

This work was supported by the Department of Biotechnology, Ministry of Science and Technology [GAP357]; Department of Science and Technology, Ministry of Science and Technology [HRR/2016/000093].