Predicting 3D structure and stability of RNA pseudoknots in monovalent and divalent ion solutions

Abstract

RNA pseudoknots are a kind of minimal RNA tertiary structural motifs, and their three-dimensional (3D) structures and stability play essential roles in a variety of biological functions. Therefore, to predict 3D structures and stability of RNA pseudoknots is essential for understanding their functions. In the work, we employed our previously developed coarse-grained model with implicit salt to make extensive predictions and comprehensive analyses on the 3D structures and stability for RNA pseudoknots in monovalent/divalent ion solutions. The comparisons with available experimental data show that our model can successfully predict the 3D structures of RNA pseudoknots from their sequences, and can also make reliable predictions for the stability of RNA pseudoknots with different lengths and sequences over a wide range of monovalent/divalent ion concentrations. Furthermore, we made comprehensive analyses on the unfolding pathway for various RNA pseudoknots in ion solutions. Our analyses for extensive pseudokonts and the wide range of monovalent/divalent ion concentrations verify that the unfolding pathway of RNA pseudoknots is mainly dependent on the relative stability of unfolded intermediate states, and show that the unfolding pathway of RNA pseudoknots can be significantly modulated by their sequences and solution ion conditions.

About the Authors:

Ya-Zhou Shi

Contributed equally to this work with: Ya-Zhou Shi, Lei Jin

Roles Conceptualization, Data curation, Formal analysis, Funding acquisition, Investigation, Methodology, Validation, Writing – original draft, Writing – review & editing

Affiliations Research Center of Nonlinear Science, School of Mathematics and Computer Science, Wuhan Textile University, Wuhan, China, Department of Physics and Key Laboratory of Artificial Micro- and Nano-structures of Ministry of Education, School of Physics and Technology, Wuhan University, Wuhan, China

Lei Jin

Contributed equally to this work with: Ya-Zhou Shi, Lei Jin

Roles Conceptualization, Methodology, Writing – review & editing

Affiliation: Department of Physics and Key Laboratory of Artificial Micro- and Nano-structures of Ministry of Education, School of Physics and Technology, Wuhan University, Wuhan, China

Chen-Jie Feng

Roles Data curation, Investigation, Methodology

Affiliation: Department of Physics and Key Laboratory of Artificial Micro- and Nano-structures of Ministry of Education, School of Physics and Technology, Wuhan University, Wuhan, China

Ya-Lan Tan

Roles Data curation, Formal analysis, Validation

Affiliation: Department of Physics and Key Laboratory of Artificial Micro- and Nano-structures of Ministry of Education, School of Physics and Technology, Wuhan University, Wuhan, China

Zhi-Jie Tan

Roles Conceptualization, Data curation, Formal analysis, Funding acquisition, Methodology, Project administration, Resources, Supervision, Validation, Writing – review & editing

* E-mail: [email protected]

Affiliation: Department of Physics and Key Laboratory of Artificial Micro- and Nano-structures of Ministry of Education, School of Physics and Technology, Wuhan University, Wuhan, China

ORCID logo http://orcid.org/0000-0002-7367-0182

Introduction

RNAs can fold into complex three-dimensional (3D) structures to carry out their various biological functions [1]. An RNA pseudoknot represents a very common structure motif, which is not only one of the fundamental structure elements in various classes of RNAs such as human telomerase RNA, self-splicing introns of ribozyme and S-adenosylmethionine-responsive riboswitches, but also involved in many biological functions, including regulation and catalysis [2,3]. For instance, an RNA pseudoknot can be present within the coding regions of an mRNA, where it stimulates programmed -1 ribosomal frameshifting to control the relative expression levels of proteins [2–4]. Generally, an RNA pseudoknot is formed when a sequence of nucleotides within a single-stranded loop region forms base pairs with a complementary sequence outside that loop [2,3,5,6]. Many experiments have shown that this special 3D topology is key to realize the various functions of RNA pseudoknots [2–4]. In addition, the stability of RNA pseudoknots can also play important roles in modulating their biological functions, and structure changes of RNA pseudoknots could cause diseases such as dyskeratosis [3,7,8]. Thus, to determine 3D structures and quantify stability of RNA pseudoknots is essential to unveil the mechanisms of their functions and to further aid the related drug design [5,9].

There have been several successful experimental methods to obtain 3D structures of RNAs, such as X-ray crystallography, nuclear magnetic resonance spectroscopy, and newly developed cryo-electron microscopy [9–12]. However, it is still very time-consuming and expensive to derive high-resolution 3D structures of RNAs and the RNA structures deposited in Protein Data Bank (PDB) are still limited [9,12]. To complement experimental measurements, some computational models have been developed to predict 3D structures for RNAs [13–22]. The knowledge-based models [23–34] such as MC-Fold/MC-Sym pipeline [24], FARNA [25], 3dRNA [29,35,36], RNAComposer [30] and pk3D [31] are rather successful and efficient in constructing 3D structures for RNA pseudoknots through fragments assembly based on limited experimental structures/fragments or reliable secondary structures, while it is still a problem to exactly predict secondary structures of RNA pseudoknots [11,20]. Furthermore, most of the above methods cannot give reliable predictions for the thermodynamic properties of RNA pseudoknots from their sequences [9–11].

Simultaneously, some coarse-grained (CG) models have been developed to predict the thermodynamic stability of RNAs including pseudoknots [37–46]. The Vfold model enables predictions for the structure, stability, and the free energy landscape for RNA pseudoknots from sequences through enumerating loop conformations on a diamond lattice [37,38]. The model is applicable to secondary structure folding while the 3D structures need to be built through fragment assembly based on secondary structures [47]. Several other CG models such as the iFoldRNA [39], the HiRE-RNA [40] and the oxRNA [42] have been used to predict 3D structure and stability for a few RNA pseudoknots, but the parameters of these models may need further validation for quantifying RNA thermodynamics to accord with experiments. In addition, due to the polyanionic nature of RNAs, metal ions (e.g., Na+ and Mg2+) in solutions can play an essential role in RNA folding [48–53], and Mg2+ can play a more special role in stabilizing the compact folded structures of RNA pseudoknots [54–57]. However, the above structure prediction models seldom consider the conditions departing from the high salt (e.g., 1M NaCl). Although all-atomic molecular dynamics simulations can be used to probe ion-RNA interactions, it is still difficult to simulate RNA structure folding at present due to the huge computation cost [56,57]. In simplified CG models, the effect of ions (especially Mg2+) is seldom properly involved due to the interplay between ion binding and structure deformation as well as the particularly efficient role of Mg2+ beyond mean-field description [51–53]. Recently, a Gö-like CG model has been introduced to reproduce the folding thermodynamics of several RNA pseudoknots in the presence of monovalent ions [46,58,59], and another structure-based model can well capture the ion atmosphere around RNAs with an explicit treatment of divalent ions [60]. However, the two structure-based models could not be used to predict 3D structures for RNA pseudoknots solely from the sequences [11,20,46,60]. Therefore, it still remains an important problem to predict 3D structures and thermodynamic stability for RNA pseudoknots especially in monovalent/divalent ion solutions only from the sequences.

In this study, we focused on predicting 3D structures and stability for extensive RNA pseudoknots in monovalent and divalent ion solutions from their sequences through our previously developed three-bead CG model [61,62]. In the following, we first revisited the key features of our CG model such as the CG representation and the implicit-solvent/salt force field for RNAs. We then employed the model to predict 3D structures for various RNA pseudoknots from their respective sequences. Afterward, we made the prediction for the stability of typical pseudoknots with different lengths and sequences over a wide range of monovalent/divalent ion concentrations. Finally, we made the comprehensive analyses on the unfolding pathway for various RNA pseudoknots in ion solutions and examined the effect of monovalent/divalent ions on the unfolding pathway of RNA pseudoknots. Throughout the article, we have made the comparisons between the predictions and the extensive experimental data as well as the comparisons with the existing models.

Materials and methods

Coarse-grained structure representation for RNAs

In our model, an RNA is represented as a chain of nucleotides, where each nucleotide is reduced to three beads retaining the key structure features of an RNA chain [46,47,61,62]. As shown in Fig 1A, the backbone phosphate bead (P) and sugar bead (C) coincide with the phosphate and C4’ atoms of a nucleotide, and the base beads (N) are placed at the base atoms linked to the sugar, that is N1 atom for pyrimidine or the N9 atom for purine [61,62]. The P, C and N beads are treated as spheres with van der Waals radii of 1.9Å, 1.7Å and 2.2 Å, respectively, and each P bead has a charge of–e on its center [61,63].

[Figure omitted. See PDF.]

Fig 1. The representation of all-atom and CG model for an RNA pseudoknot, and 3D structure prediction for a paradigm pseudoknot in the present model.

(a) The 3D structure of MMTV pseudoknot (PDB code: 1rnk) in all-atomistic (left) and our CG representation (right) as well as the schematic representation for one nucleotide in the present CG model (middle). (b) The secondary structure for MMTV pseudoknot consisting of two stems and three loops. The corresponding secondary structure elements in (a) and (b) are in same colors: Stem 1 (red), Stem 2 (green), Loop 1 (blue), Loop 2 (cyan), Loop 3 (magenta) and the 3' dangling nucleotide U (orange). (c, d) The time-evolution of the energy (top panel in (c)), the number of base pairs (middle panel in (c)), the RMSDs between predicted structures and the native structure in PDB (bottom panel in (c)), and the typical 3D structures (d) during the Monte Carlo simulated annealing simulation of the Aquifex aeolicus tmRNA pseudoknot PK1 (PDB code: 2g1w). The inset in bottom panel in (c) shows the zoomed portion of the figure in the interval of [1.6✕107, 2✕107]. The RMSDs are calculated over C beads from the corresponding C4’ atoms in native structure, and the structures in (d) are shown with the PyMol (http://www.pymol.org).

https://doi.org/10.1371/journal.pcbi.1006222.g001

Coarse-grained force field and simulation procedure

In the CG model, the effective potential energy of an RNA conformation is given by [61,62](1)where bond length energy Ub, bond angle energy Ua and dihedral energy Ud account for chain connectivity and angular rotation for an RNA chain, and Uexc represents for excluded volume interactions between two CG beads. Ubp and Ubs are the base-pairing and base-stacking interactions, and Ucs is the coaxial stacking interaction between two neighbor stems. The last term Uel corresponds to electrostatic interactions between phosphate groups, which are ignored by most of the existing predictive models for RNA 3D structures [11,20].

The detailed description of the potentials in Eq 1 and the determination of the potential parameters have been described in S1 Text and also in Refs. [61,62]. Briefly, two sets of parameters of the bonded potentials (Ub, Ua and Ud), Paranonhelical used in RNA folding process and Parahelical used only in structure refinement for helical stems, are derived respectively from single strands/loops and stems in the PDB [12,61,64]. The sequence-dependent strength of base-staking energy is derived from the combination of the experimental thermodynamic parameters [65–67]. In most occurring pseudoknots with interhelix loop length ≤ 1nt, two helical stems can be often coaxially stacked to form a quasi-continuous double helix (Fig 1), and the strength of Ucs depends on sequences of two interfaced base pairs [65]. The coaxial stacking could stimulate high levels of -1 frameshifting [3,4], and consequently, could be important for stabilizing functional structures of RNA pseudoknots. The electrostatic interaction Uel is taken into account through the combination of the Debye-Hückel approximation and the concept of counterion condensation (CC) [68]. Notably, based on the tightly bound ion (TBI) model [51,52,69], the competition between monovalent and divalent ions was also taken into account in Uel to enable the CG model to simulate RNA pseudoknot folding in mixed monovalent/divalent ion solutions [61,62]. Although the present model has been described by us in Refs. 61 and 62, the model is still not employed for 3D structure predictions of extensive RNA pseudoknots and it has never been used to predict the stability of RNA pseudoknot in ion solutions, especially in the presence of divalent ions [61,62]. Here, the model will be tested by extensive RNA pseudoknots on 3D structure prediction, and be further used to predict thermodynamic stability and the unfolding pathway for various RNA pseudoknots over the wide range of monovalent/divalent ion conditions.

Based on the CG force field, the Monte Carlo (MC) simulations with simulated annealing algorithm are used to predict 3D structures of RNA pseudoknot [43,61,62], where an initial simulation is started at a high temperature and a given solution condition from a totally random chain configuration generated from an RNA sequence. The system is then gradually cooled in steps, and the ion condition is fixed during the cooling process. At each temperature, RNA conformational changes are accomplished via the pivot moves which have been demonstrated to be rather efficient in sampling conformations of polymers [63], and the changes are accepted or rejected according to the standard Metropolis algorithm [43,61]. The final structures obtained at the lowest target temperature (e.g., room/body temperature) are the folded conformations of the RNA predicted by the CG model. Notably, the recorded trajectories at different temperatures during the cooling process allow us to analyze the stability of the RNAs [61,62].

Results

In this section, first, we employed the present CG model to predict the 3D structures for extensive RNA pseudoknots. Afterwards, the CG model was used to predict the stability of various RNA pseudoknots and the effects of monovalent and divalent ions. Finally, we made the comprehensive analyses on the unfolding pathway of RNA pseudoknots and the ion effect. Our predictions were compared with the available experimental data and existing models.

Predicting 3D structures of RNA pseudoknots

Folding process and structure refinement.

For each RNA pseudoknot, a random chain is generated from its sequence only based on the potentials of Ub and Uexc in Eq 1. Afterwards, for the random configuration, the MC simulation with simulated annealing algorithm is performed from high temperature to the target temperature (e.g., 298 K) with the use of Paranonhelical parameters. As an example, Fig 1C shows the folding process of a small RNA pseudoknot (22nt; PDB code: 2g1w; sequence: 5’-GGGGUGGCUCCCCUAACAGCCG-3’) in the present model. As temperature is gradually decreased from 130°C to 25°C, the energy of the RNA chain reduces with the formation of base pairs (Fig 1C), and the initial random chain folds into its native-like pseudoknotted structures; see Fig 1D. Following that, another MC simulation (e.g., final 5✕106 steps in Fig 1C) is performed at target temperature based on the final structure predicted by the preceding annealing process, and the two sets of bonded potential parameters Paranonhelical and Parahelical are employed respectively for the single strands/loops and base-pairing regions to better capture the geometry of helical part [61]. As a result, an ensemble of refined 3D structures are obtained, and can be evaluated by their root-mean-square deviation (RMSD) values calculated over C beads from the corresponding C4’ atoms in the native structure in PDB [70]; see Fig 1C. The mean RMSD (the averaged value over the structure ensemble in the refinement process) and minimum RMSD (corresponding to the structure closest to the native one in the refinement process) are used to evaluate the reliability of our predictions on 3D structures. As shown in the inset of the bottom panel of Fig 1C, the mean and minimum RMSDs of the paradigm RNA pseudoknot (PDB code: 2g1w) between predicted structures and its native structure are 4.8 Å and 3.3 Å, respectively, and the corresponding predicted 3D structures as well as the native one are also shown in Fig 2A.

[Figure omitted. See PDF.]

Fig 2. Comparisons of RMSDs between the present model and other models.

(a) The predicted 3D structures (ball-stick) with the mean RMSDs (top) and the minimum RMSDs (bottom) for four sample RNA pseudoknots (PDB codes: 2g1w, 2tpk, 1e95, and 2lc8) from their native structures (cartoon). The mean (minimum) RMSDs for three pseudoknots are 4.8 Å (3.3 Å), 4.4 Å (2.8 Å), 5.4 Å (3.5 Å) and 7.4 Å (5.1 Å), respectively, and the 3D structures are shown with the PyMol (http://www.pymol.org). (b) The predictions for the 3D structures of 17 RNA pseudoknots from the present model, from the MC-Fold/MC-Sym pipeline and from the RNAComposer. The RMSDs of predicted structures for 17 RNA pseudoknots are calculated over C beads from the corresponding C4’ atoms in native structures.

https://doi.org/10.1371/journal.pcbi.1006222.g002

Structure prediction and comparisons with previous models.

To examine the ability of the model on predicting 3D structures of RNA pseudoknots, 17 common pseudoknots (≤ 56nt) which have been determined by experiments as individual molecules were used in our 3D structure prediction. The detailed descriptions of these pseudoknots are listed in Table B in S1 Text. Due to lack of the ion conditions for the experimental structures determined by X-ray crystallography, here we only predicted the 3D structures for all RNA pseudoknots at 1M [Na+].

Fig 2A shows the predicted 3D structures (ball-stick) with the mean and minimum RMSDs and the experimental structures (cartoon) for four typical RNA pseudoknots with different lengths and sequences. As shown in Fig 2, the present model can effectively capture the 3D shapes of RNA pseudoknots, in which Loop 1 crosses the deep major groove of the lower Stem 2, while Loop 2 generally crosses the minor groove side of Stem 1 [2,3]. The mean and minimum RMSDs for the 17 pseudoknots are shown in Fig 2B, and for most of pseudoknots, the mean and minimum RMSDs are less than 6 Å and 4 Å, respectively, which suggest that the model can make reliable predictions for 3D structures of RNA pseudoknots.

Furthermore, we also made comparisons with the MC-Fold/MC-Sym pipeline [24] and RNAComposer [30], which are well-established web services (http://www.major.iric.ca/MC-Fold/ and http://rnacomposer.ibch.poznan.pl/, respectively) for predicting tertiary structures of RNAs including pseudoknots with high accuracy [10,13]. For the 17 pseudoknots used here (Table B in S1 Text), we first employed the MC-Fold (option: consider H-type pseudoknots, return the best 1000 structures, and explore the best 50% sub-optimal structures) to predict their secondary structures in MC-Sym format, the best ones of which are further submitted (or edited and then submitted) to MC-Sym (option: model_limit = 9999 and time_limit = none) for tertiary structure predictions. The RMSDs of best structures (top 1) predicted by the MC-Fold/MC-Sym pipeline are calculated over C4’ atoms from the corresponding atoms in the experimental structures in PDB; see Fig 2B. It should be noted that although we chose the best option of MC-Fold/MC-Sym, the pipeline still fails to predict the 3D structures for two RNA pseudoknots (PDB codes: 2m8k and 2lc8), even though the experimental secondary structures are taken as input. For the 15 pseudoknots except for 2m8k and 2lc8, the overall mean RMSD from the present model is ~5.4 Å, which is smaller than that of 6.3 Å from the MC-Fold/MC-Sym pipeline, suggesting that the present model gives slightly better predictions for tested sequences. Notably, for 2m8k of 48nt and 2lc8 of 56nt, the present model also gives good predictions with mean RMSDs of 8.9 Å and 7.4Å and minimum RMSDs of 6.0 Å and 5.1 Å, respectively. Second, we further employed the RNAComposer (in interactive mode) to predict tertiary structures for the 17 pseudoknots by entering their sequences and experimental secondary structures, and the RMSD between the predicted and experimental structures is also calculated over C4’ atoms. As shown in Fig 2B, the average prediction accuracy (overall mean RMSD = ~5.6 Å and overall minimum RMSD = ~3.9 Å) of the present model for the 17 pseudoknots is slightly better than that of the RNAComposer (mean RMSD = 7.7 Å). Therefore, the comparisons with the MC-Fold/MC-Sym pipeline and RNAComposer show that the present model can make reliable predictions for 3D structures of RNA pseudoknots from sequences.

To clarify the contributions of the coaxial stacking potential, we further made the additional predictions on 3D structures for 17 RNA pseudoknots using the present model without involving the coaxial stacking potential. As shown in Table C in S1 Text, for the RNA pseudoknots except for the five ones (PDB codes: 2a43, 1l2x, 1yg4, 437d, 2ap5), the present model with coaxial stacking can make the better predictions with lower RMSDs compared to those without involving the coaxial stacking potential, which suggests that the inclusion of the coaxial stacking in the model can generally improve the 3D structure prediction for RNA pseudoknots. However, for the RNA pseudoknots such as 2a43 and 1l2x from plant luteovirus which do not have coaxial stacking interactions [4,71–74], the involvement of the coaxial stacking potential gives slightly worse predictions than those without the coaxial stacking. This suggests that the coaxial stacking potential may need to be further developed in more details for various RNA pseudoknots.

Predicting stability of RNA pseudoknots in monovalent/divalent ion solutions

Beyond 3D structure prediction, the present model was also employed to predict the stability of RNA pseudoknots in monovalent and divalent salt solutions.

Predicting RNA pseudoknot stability.

First, we used the MMTV frameshifting pseudoknot [75] as an example to show how to examine RNA pseudoknot stability with the use of our CG model. MMTV pseudoknot is an H-type RNA pseudoknot containing an unpaired adenosine at the junction of the two helical stems; see Fig 1B. The stability of MMTV pseudoknot directly affects the efficiency of frameshifting activity [75]. Beyond the 3D structure predictions for MMTV pseudoknot (PDB code: 1rnk) (Fig 2), we further employed the present model to examine the stability of MMTV pseudoknot. Fig 3A shows that the number of formed base pairs changes at different temperatures, and there are mainly three states over temperature: the fully folded pseudoknot state (F) at low temperatures (e.g., <~40°C); the unfolded coil state (U) at high temperatures (e.g., >~100°C); and partially unfolded intermediate hairpin states (I) at medium temperature which is coincident with experiments [75–77]. Then, the fractions of the three states at each temperature can be calculated (Fig 3B), and the fractions of folded and unfolded states (fF(T) and fU(T)) can be fitted to a two-state model through the following equations [61,65]:(2)(3)Here, Tm1 and Tm2 are two melting temperatures of the corresponding transitions (F→I and I→U), respectively. dT1 and dT2 are corresponding adjustable parameters. As shown in Fig 3, for MMTV pseudoknot at 50mM [K+], the predicted Tm1 and Tm2 are 49.0°C and 81.9°C, which agree well with the experimental data (50.2°C and 83.0°C) [75].

[Figure omitted. See PDF.]

Fig 3. The stability prediction for a sample RNA pseudoknot in the present model.

(a) The time-evolution of the number of base pairs for MMTV pseudoknot (shown in Fig 1A) at different temperatures (100°C, 80°C, 60°C, 40°C from top to bottom, respectively) in 50mM KCl solution. (b) The fractions of folded state (F, black), unfolded state (U, blue), and intermediate state (I, red) as a function of temperature for MMTV pseudoknot at 50mM [K+]. The dotted lines are fitted to the predicted data through Eqs. 2 and 3. (c) The fraction of denatured base pairs f as a function of temperature for MMTV pseudoknot at 50mM [K+]; symbols: from the present model; line: from Eq 4. Ball-stick: the typical 3D structures predicted at different temperatures shown with the PyMol (http://www.pymol.org).

https://doi.org/10.1371/journal.pcbi.1006222.g003

Furthermore, based on the fF(T) and fU(T), the fraction of number of denatured base pairs f can be calculated through the following equation [76](4)where fI is the fraction of number of denatured base pairs when the fraction of intermediate state is maximum; see Fig 3C. As shown in Figs 4A and 5A, df/dT (the first derivative of f calculated by Eq 4 with respect to temperature) profile for MMTV pseudoknot is in good agreement with previous differential scanning calorimetry profile [75]. This suggests that the present model can give reliable predictions for stability of the RNA pseudoknot.

[Figure omitted. See PDF.]

Fig 4. The comparisons between predictions (lines) and experiments (symbols) for six RNA pseudoknots with various sequences at high salt concentrations.

(a) MMTV pseudoknot at 1000mM [K+] [75]; (b) T2 pseudoknot at 1000mM [K+] [77]; (c) PEMV-1 pseudoknot at 500mM [K+] [73]; (d) BWYV pseudoknot at 500mM [K+] [54,72]; (e) PLRV pseudoknot at 500mM [K+] [73]; and (f) ScYLV pseudoknot at 500mM [K+] [74]. Lines: df/dT, the first derivative of f with the temperature. Symbols: the heat capacity Cp or dA/dT, the first derivative of absorbance with respect to temperature.

https://doi.org/10.1371/journal.pcbi.1006222.g004

[Figure omitted. See PDF.]

Fig 5. The comparisons between predictions (lines) and experiments (symbols) for MMTV and T2 pseudoknots in monovalent ion solutions.

(a) MMTV pseudoknot at 1000mM [K+] (green square) and 50mM [K+] (red triangle), respectively. Lines: df/dT, the first derivative of f with respect to temperature. Symbols: the heat capacity Cp. (c) T2 pseudoknot at 1000mM [K+] (green square) and 100mM [K+] (red triangle), respectively. Symbols: dA/dT, the first derivative of absorbance with respect to temperature. (b, d) The melting temperatures Tm1 and Tm2 of two transitions (F→I and I→U) as functions of [K+] for MMTV (b) and T2 (d) pseudoknots. Symbols: experimental Tm1 (blue square) and Tm2 (red triangle) [75,77]. Lines: corresponding predictions.

https://doi.org/10.1371/journal.pcbi.1006222.g005

In addition to MMTV pseudoknot, five other pseudoknots with different sequences from gene 32 mRNA [77] and plant luteoviruses [71–74] are also examined by the present model, and the sequences and predicted secondary structures for the pseudoknots can be found in S1 Fig. As shown in Table 1, our predictions on melting temperatures (Tm1 and Tm2) for the six RNA pseudoknots at high salt concentrations (1000mM [K+] for MMTV and T2, and 500mM [K+] for the other four pseudoknots) agree well with the experimental data with the mean deviations of ~2.5°C for Tm1 and ~1.8°C for Tm2. In addition, Fig 4 shows the comparisons between the calculated and the experimental thermal unfolding curves [71–77] for the six RNA pseudoknots at high salt, indicating that the predicted melting curves are also in good accordance with the experiments over different sequences. These suggest that the present model with stacking interactions parametrized using a set of sequence-specific thermodynamics parameters can quantitatively predict the stability of RNA pseudoknots with various sequences. Nevertheless, the predicted Tm1’s for BWYV and PEMV-1 are slightly lower than the corresponding experimental data; see Table 1. One possible reason is that both pseudoknots contain a collection of loop-stem tertiary structural interactions, which could make contributions to the stability of the pseudoknotted structures [4,72–74].

[Figure omitted. See PDF.]

Table 1. The melting temperatures (Tm1 and Tm2)a of six RNA pseudoknots at high salt concentrations.

https://doi.org/10.1371/journal.pcbi.1006222.t001

To further evaluate the present model, we made the comparisons with the available predictions from the existing models such as the model of Denesyuk and Thirumalai [46,58,59] and the Vfold model with monovalent salt-corrected thermodynamic parameters and the fitting parameters of loop-stem tertiary contacts [37,38,78], for the stability of MMTV, BWYV and PEMV-1 pseudoknots in monovalent ion solutions. As shown in S2 Fig, for BWYV and PEMV-1 pseudoknots at 500mM [K+], the mean deviation between Tm’s (Tm1+Tm2) from the Vfold model and the experiments is ~8.4°C [78], which is slightly higher than that from the present model (~6.4°C). For MMTV pseudoknot at 1000mM [K+], the deviation between Tm’s (Tm1+Tm2) from the model of Denesyuk and Thirumalai and the experiment is ~8.2°C [46], a higher deviation than that from the present model (~3.1°C). Such deviation from the model of Denesyuk and Thirumalai and the present model becomes rather small for MMTV at 50mM [K+]. Thus, the present model can be reliable in predicting thermal stability for RNA pseudoknots in monovalent ion solutions, and it is noted that the present model can also provide 3D structures for RNA pseudoknots at different temperatures from the sequences; see Fig 3C and S3 Fig. Recently, the HiRE-RNA model and the oxRNA model have been proposed and both can predict the presence of two peaks in the melting curves of several RNA pseudoknots, however, there is still lack of the quantitative comparisons between the two models and the experiments for extensive RNA pseudoknots [40,42].

Furthermore, we made the additional calculations for the stability of the six RNA pseudoknots (Table 1 and Table D in S1 Text) using the present model without involving the coaxial stacking potential. As shown in Table 1, S3 Fig and Table D in S1 Text, Tm1’s from the present model with coaxial stacking are generally higher than those without coaxial stacking and appear closer to the experimental values [72–77], suggesting that the involvement of coaxial stacking enhances the stability of RNA pseudoknotted structures and improves our predictions on RNA pseudoknot stability. It is also noted that, in some pseudoknots such as BWYV pseudoknot, no coaxial stacking is found between two stems, while triple/quadruple base interactions are observed at the junctions in their experimental structures [71–74]. It is interesting that the present model without considering triple/quadruple interactions still make good predictions on the stability for this kind of pseudoknots. This should be attributed to that the involvement of coaxial stacking in the present model partially compensates the lack of triple/quadruple base interactions in the model.

Monovalent ion effect on RNA pseudoknot stability.

Due to the high density of negative charges on its backbone, RNA pseudoknot stability is sensitive to ionic condition of solution [49,51–54,77,79], while the effect of salts is generally ignored in the existing RNA structure prediction models [11,20]. For MMTV pseudoknot, we further performed the simulations at different temperatures over a broad range of [K+], and calculated the melting temperatures based on the data from the simulations. As shown in Fig 5A, the unfolding of MMTV pseudoknot at different [K+]’s has the similar two transitions in accordance with the available experiments [75], and the predicted melting temperatures (Tm1 and Tm2) agree well with the experimental data with the maximum deviation of ~1.8°C over different K+ concentrations; see Fig 5B.

In addition, we examined the stability of the bacteriophage T2 gene 32 mRNA pseudoknot in solutions of different [K+]’s. The sequence and predicted secondary structure of T2 pseudoknot are shown in S1 Fig, and the predicted and experimental 3D structures (PDB code: 2tpk) are shown in Fig 2. As shown in Fig 5C, in contrast to MMTV pseudoknot, there is one visible peak in experimental unfolding curves for T2 pseudoknot at 100mM [K+] and 1000mM [K+] due to the smaller difference between two transition temperatures Tm1 and Tm2 for T2 pseudoknot [75,77], i.e, at 1M [K+], Δ = |Tm1-Tm2| ≈ 12°C for T2 pseudoknot and ~24°C for MMTV pseudoknot. In addition, the predicted Tm1 and Tm2 for T2 pseudoknot are also in good accordance with the experimental data with the mean deviation of ~1.1°C over the wide range of [K+] [77]; see Fig 5D.

As shown in Fig 5B and 5D, the increase of [K+] enhances the stability of RNA pseudoknots as well as the stability of intermediate hairpins, and Tm1 of the transition from folded state to intermediate state is more sensitive to [K+] than Tm2 of the transition from intermediate state to fully unfolded chain. The phenomenon is interesting and is reasonable. Generally, the formation of a pseudoknot structure involves higher charge buildup than the formation of intermediate hairpin states, and hence, the stability of a pseudoknot is more dependent on ions than that of a hairpin [51–54].

Divalent ion effect on RNA pseudoknot stability.

Previous studies have shown that divalent ions such as Mg2+ are especially effective in stabilizing RNA tertiary structure [48,52–54,79]. As described in Section of Material and methods, the effect of divalent ions has been implicitly accounted for in the present model by combining the CC theory [68] and the results from the TBI model [51,69]. Here, beyond previous computational models for RNA pseudoknot stability in monovalent ion solutions, the present model is examined by predicting the stability for two typical RNA pseudoknots (MMTV and T2) in mixed monovalent/divalent ion solutions. As shown in Fig 6, the thermal unfolding curves and Tm’s predicted by the present model for MMTV and T2 pseudoknots are in accordance with the experiments [75,77] over a wide range of [Mg2+] with fixed [K+]’s (50mM for MMTV pseudoknot and 100mM for T2 pseudoknot). For MMTV pseudoknot, the mean deviations of Tm1 and Tm2 between predictions and experiments are ~3.9°C and ~2.3°C, respectively, and for T2 pseudoknot, the corresponding mean deviations are ~2.2°C for Tm1 and ~1.5°C for Tm2. This suggests that the present model can nearly make quantitative predictions for the stability of RNA pseudoknots in mixed ion solutions from their sequences, even though the ion effect is involved implicitly in the present model.

[Figure omitted. See PDF.]

Fig 6. The comparisons between predictions (lines) and experiments (symbols) for MMTV and T2 pseudoknots in divalent ion solutions.

(a) MMTV pseudoknot at 50mM [K+] with 0.5mM [Mg2+] and 3.0mM [Mg2+], respectively. Symbols: the heat capacity Cp for MMTV pseudoknot at 50mM [K+] with 0.5mM [Mg2+]. (c) T2 pseudoknot at 100mM [K+] with 0.1mM [Mg2+] and 1.0mM [Mg2+], respectively. Symbols: dA/dT, the first derivative of absorbance with respect to temperature. (b, d) The melting temperatures Tm1 and Tm2 of two transitions (F→I and I→U) as functions of [Mg2+] for MMTV pseudoknot in the presence of 50mM [K+] (b) and T2 pseudoknot in the presence of 100mM [K+] (d). Symbols: experimental Tm1 (blue square) and Tm2 (red triangle) [75,77]. Lines: corresponding predictions.

https://doi.org/10.1371/journal.pcbi.1006222.g006

Fig 6 also shows that the competition between K+ and Mg2+ on RNA pseudoknot stability is captured by the present model: (i) when [Mg2+] is very low, K+ ions dominate the stability of pseudoknots and the values of Tm including Tm1 and Tm2 are close to those of pure K+ solutions; (ii) the increase of [Mg2+] significantly enhances the stability of RNA pseudoknots against monovalent K+ and such effect would become saturated at very high [Mg2+] due to the strong electrostatic neutralization; see Fig 6B and 6D. This is attributed to the anticooperative binding of K+ and Mg2+ and more efficient role of Mg2+ binding [49–52]. Nonetheless, the predicted melting temperatures, e.g., Tm1 of the first transition, are slightly lower than the corresponding experimental values at high [Mg2+]; see Fig 6B and 6D. This may be attributed to the fact that the implicit Mg2+ treatment in the present model might slightly underestimate the role of Mg2+ in stabilizing compact pseudoknot structure [53].

Thermally unfolding pathway of RNA pseudoknots in ion solutions

Since intermediate states of RNAs can be important to their biological functions [5,76,80–85], unfolding pathway of RNAs including some pseudoknots has been studied through theoretical modeling and experiments [75–77,81–88]. To examine the unfolding pathway of RNA pseudoknots, we made comprehensive analyses for six RNA pseudoknots; see Fig 7 and S4 Fig. Based on the simulations for each pseudoknot at a given solution condition, beyond the fractions of states F and U, the fractions of different intermediate hairpin states (named as S1 and S2 for intermediate states reserving one of Stem 1 and Stem 2, respectively) at different temperatures can also be calculated; see Figs 7 and 8 and S4 and S6 Figs. Furthermore, we employed the model to predict the unfolding pathway for various RNA pseudoknots in monovalent/divalent ion solutions and examined the effect of monovalent/divalent ions on the unfolding pathway of RNA pseudoknots, which was seldom covered in previous studies since the effect of divalent ions is generally difficult to be involved.

[Figure omitted. See PDF.]

Fig 7. The predicted thermal unfolding of MMTV, T2 and T2 variant pseudoknots at 1M [K+].

(a-c) The fractions of F, S1, S2 and U states as functions of temperature in thermal unfolding of MMTV (a), T2 (b), and T2 variant (c) pseudoknots at 1M [K+]. F stands for fully folded RNA; S1, hairpin intermediate with Stem 1; S2, hairpin intermediate with Stem 2; U, fully unfolded RNA. (d) Schematic diagram shows the dominating structural transitions of MMTV pseudoknot along the unfolding pathway inferred from the fraction of states shown in panel (a). (e) Schematic diagram shows the structural transitions of T2 pseudoknot and the variant of T2 pseudoknot along the unfolding pathway inferred from their respective fraction of states shown in panels (b) and (c).

https://doi.org/10.1371/journal.pcbi.1006222.g007

Unfolding pathway of RNA pseudoknot varies with its sequence.

As shown in Fig 7A for MMTV pseudoknot at 1M [K+], at a low temperature (e.g., <~40°C), the RNA is completely in folded pseudoknot state with two structural motifs of Stem 1 and Stem 2. As temperature is increased (e.g., 40°C-80°C), the fraction of F state decreases gradually, and simultaneously, the fraction of S1 state increases gradually, which indicates that the Stem 2 in the pseudoknot melts first with increasing temperature. When temperature is increased to higher level (e.g., >~80℃), the fraction of completely unfolded conformations increases accompanied with the decrease of the fraction of S1 state. Although the S2 state can also be found in this process, the fraction of S2 state relative to all intermediate states is relatively small (~18%). Thus, with the increase of temperature, the dominating unfolding pathway of MMTV pseudoknot is F→S1→U overwhelming the pathway of F→S2→U; see Fig 7. But for the four (PEMV-1, BWYV, PLRV and ScYLV) pseudoknots from plant luteoviruses, the unfolding processes almost undergo the only pathway F→S1→U through the intermediate state S1 and the other pathway of F→S2→U appears negligibly, as shown in S4 Fig. Interestingly, for T2 pseudoknot at 1M [K+], there are two comparable unfolding pathways: F→S1→U and F→S2→U, where the population of S1 state is only slightly higher than that of S2 state; see Fig 7B. Our calculations are in accordance with the experiments [75–77,80–83] as well as the recent theoretical studies [59,84–88]. The above analyses show that RNA pseudoknots of different sequences can have very different unfolding pathways.

What dominates the unfolding pathway of RNA pseudoknots?

It is interesting that different RNA pseudoknots can have apparently different unfolding pathways. Thirumalai and his coworkers have proposed that the unfolding pathway of several pseudoknots is largely determined by the stabilities of constituent secondary structures [59,84,86]. To examine that and to further understand the dominant factor on the unfolding pathway for extensive RNA pseudoknots, we calculated the free energies of the two intermediate states (hairpins) for each of these pseudoknots at different temperatures and 1M [Na+] using the Mfold algorithm [89]; see S5 Fig. For the pseudoknots except for T2, the two intermediate hairpin states are sufficiently different in stability, and the unstable hairpin (S2 state) containing Stem 2 with apparently higher free energy would melt first when temperature is increased. For example, the relative free energy ΔΔG = ΔGS1 − ΔGS2 between the intermediate states S1 and S2 of MMTV pseudoknot is ~-0.6kcal/mol at T~85°C around which fractions of S1 and S2 have maximum values, while for the other four (PEMV-1, BWYV, PLRV and ScYLV) pseudoknots, ΔΔG is ~-4.0kcal/mol, ~-3.6kcal/mol, ~-2.7kcal/mol and ~-3.2kcal/mol in the range of 40–90°C, respectively. That is why the intermediate S2 state occurs with a lower and visible population in unfolding process of MMTV pseudoknot and it almost never appears in unfolding processes of other four (PEMV-1, BWYV, PLRV and ScYLV) pseudoknots. Unlike the above described pseudoknots, the two intermediate hairpins of T2 pseudoknot are very similar in stability with the relative free energy ΔΔG of ~0kcal/mol (in the range of 70–80°C), which would lead to two parallel unfolding pathways with similar probabilities through the intermediate states of S1 and S2, respectively. To further examine the role of stability of intermediate state on unfolding pathway, we also made the analysis for a variant of T2 pseudoknot, in which the A3-U16 base pair in Stem 1 is substituted by a more stable G-C base pair; see Fig 7C and S5 Fig. Due to the significant enhancements of stability for S1 state containing Stem 1 with ΔΔG~-3.4kcal/mol (in the range of 70–80°C), the intermediate state of S2 nearly disappears, as shown in Fig 7C. As a result, the unfolding pathway of the variant of T2 pseudoknot is F→S1→U, which is distinctively different from that of T2 pseudoknot.

Therefore, the above comprehensive analyses for six wild-type pseudoknots and a variant of T2 pseudoknot show that the unfolding pathway of an RNA pseudoknot is mainly dependent on the stability of intermediate states, which is in consistent with the conclusion from Thirumalai and his coworkers [59,84,86]. Furthermore, despite the same stability of the two intermediate hairpins in T2 pseudoknot, the experiments show that Stem 2 can be slightly easier to melt in unfolding process [77,82,83], which indicates that other factors may affect the unfolding of pseudoknots. For example, Wang, Zhang and their coworkers have recently proposed that the contribution of the noncanonical interactions between helices and loops may make contribution to the unfolding of pseudoknots [85,88], while such interaction is not involved in the present model.

Unfolding pathway of RNA pseudoknots can be modulated by ions.

Since ions can significantly affect the stability of intermediate states [52,53,69], we further examined the effect of monovalent/divalent ions on the unfolding pathway of MMTV and T2 pseudoknots, beyond the previous analyses on unfolding pathway of RNA pseudoknots. As shown in Fig 8 and S6 Fig, as K+ concentration decreases, the structure melting transitions can get easier at low temperatures due to the enhancement of electrostatic repulsion in pseudoknots [52,53], and correspondingly, the unfolding pathway of the two pseudoknots changes due to the change of ion condition. For example, for MMTV pseudoknot, S2 state nearly never appears at low salt (e.g., 10mM K+), while gradually appears as the [K+] increase, i.e., the fraction of S2 state increases from ~0% at 10mM [K+] to ~18.0% at 1M [K+]; see Fig 8. Such prediction is in good accordance with the very recent experimental results (Roca, Hori, Velmurugu, Narayanan, Narayanan, Thirumalai, and Ansari, arXiv: 1710.0695). This is because the increase of [K+] stabilizes S2 state more pronouncedly than S1 state since S2 state has a large hairpin than S1 state (excluding dangling tails), i.e., the relative free energy ΔΔG calculated from the salt extension from the TBI model is ~-2.5kcal/mol at 10mM [K+] (at T~45°C) and becomes ~-0.6kcal/mol at 1M [K+] (at T~85°C) [55,69]. For T2 pseudoknot, as K+ concentration is decreased from 1M, the probabilities of the two parallel pathways change visibly, i.e., the fraction of S2 decreases from ~44.0% at 1M [K+] to ~18.5% at 10mM [K+]. This is also attributed to the larger hairpin of S2 state and the corresponding stronger ion effect in structure stabilization, i.e., the relative free energy ΔΔG is ~-0.1kcal/mol at 1M [K+] (at T~75°C) and becomes ~-1.6kcal/mol at 10mM [K+] (at T~40°C) [55,69]. Similarly, as shown in S6 Fig, Mg2+ can also significantly affect the unfolding pathway of MMTV and T2 pseudoknots although they are in the buffers containing 50mM and 100mM K+, respectively.

[Figure omitted. See PDF.]

Fig 8.

The fractions of F, S1, S2 and U states as functions of temperature in thermal unfolding of MMTV (a-c) and T2 (d-f) pseudoknots at different [K+]’s. F stands for fully folded RNA; S1, hairpin intermediate with Stem1; S2, hairpin intermediate with Stem2; U, fully unfolded RNA.

https://doi.org/10.1371/journal.pcbi.1006222.g008

The above analyses for MMTV and T2 pseudoknots indicate that the change of ion conditions can apparently modulate the unfolding pathway of RNA pseudoknots through changing the relative stability between the two unfolding intermediate states at different ion conditions.

Discussion

It is important to predict 3D structures and stability of RNA pseudoknots in monovalent/divalent ion solutions from their sequences. In this work, we employed our previously developed model to address this problem. Beyond mainly focusing on reproducing structures, as many previous structure prediction models have done, the present model enables us to predict and analyze 3D structure stability for RNA pseudoknots in different monovalent/divalent ion solutions. The following are the major conclusions:

1. The present model predicts the native-like 3D structures for RNA pseudoknots with an overall mean RMSD of 5.6 Å and an overall minimum RMSD of 3.9 Å from experimental structures, and the overall prediction accuracy of our model is slightly higher than previous models.

2. The present model successfully predicts the stability of RNA pseudoknots with different lengths and sequences over a wide range of monovalent/divalent ion concentrations, and the predicted melting temperatures for the two unfolding transitions are in good accordance with extensive experiment data.

3. Our comprehensive analyses show that the unfolding pathway of RNA pseudoknots is mainly determined by the stabilities of intermediate states which can be significantly modulated by the sequences and solution ion conditions.

Despite the extensive agreements between our predictions and experiments, the present model has several limitations that should be overcome in future model development. First, the present model does not treat possible noncanonical interactions such as base triple interactions between loops and stems, self-stacking in loop nucleotides and special hydrogen bonds involving phosphates and sugars, which could be important for some more complex pseudoknotted structures [7,17,38]. Beyond the common H-type pseudoknots (≤ 56nt) used in this work, larger RNAs with complex structures should be incorporated in to further improve the present model [19,90–94]. Second, the effect of monovalent/divalent salts is implicitly accounted for in the present model by the combination of CC theory and the TBI model. Such implicit-salt treatment may be responsible for the underestimation on the stability of RNA pseudoknots at high [Mg2+]. Mg2+ can play an efficient and special role in stabilizing compact RNA structures [51–54,79], and further development may need to involve Mg2+ explicitly in our model. Third, in this work, we mainly focused on the 3D structures and thermodynamic stability of RNA pseudoknots, and did not involve the stability under mechanical force. Mechanical forces can be not only considered as a useful probe for RNA stability, but also important for the functions of some RNA pseudoknots [3,81,86,95–97]. For example, the frameshifting efficiency may be affected by the magnitude of unfolding force for RNA pseudoknots [3,81,96]. Fortunately, the present model can be extended to study the mechanical stability of RNA pseudoknots by including external force in the energy functions of the model [67,86,95]. Finally, the 3D structure predicted by the present model is at the CG level, and it is still necessary to develop the model to reconstruct all-atomistic structures based on the CG structures for further practical applications. Nevertheless, the present model could be a reliable predictive model for predicting 3D structures and stability of RNA pseudoknots in ion solutions from their sequences and the analyses can be helpful to understand the physical mechanism for the unfolding pathway of RNA structures.

Supporting information

S1 Text. The force field of the present model, the RNA pseudoknots for 3D structure prediction used in this work, and the 3D structures and stability of RNA pseudoknots predicted by the present model with/without the coaxial stacking potential.

https://doi.org/10.1371/journal.pcbi.1006222.s001

(PDF)

S1 Fig. The predicted secondary structure of the six pseudoknots for stability prediction used in this work.

(a) MMTV pseudoknot; (b) T2 pseudoknot; (c) PEMV-1 pseudoknot; (d) BWYV pseudoknot; (e) PLRV pseudoknot; (f) ScYLV pseudoknot.

https://doi.org/10.1371/journal.pcbi.1006222.s002

(TIF)

S2 Fig. The comparisons between predictions from the present model (solid lines) and other models (dotted lines) for several pseudoknots.

(a,b) BWYV (a) and PEMV-1 (b) pseudoknots at 500mM [K+], respectively. Solid lines: df/dT, the first derivative of f with respect to temperature from the present model. Dotted lines: the heat capacity Cp from Ref. 78. Symbols: the heat capacity Cp from experiments [72,73]. (c,d) MMTV pseudoknot at 1000mM [K+] (c) and 50mM [K+] (d), respectively. Solid lines: df/dT, the first derivative of f with respect to temperature from the present model. Dotted lines: the heat capacity Cp from Ref. 46. Symbols: the heat capacity Cp from experiments [75].

https://doi.org/10.1371/journal.pcbi.1006222.s003

(TIF)

S3 Fig. The comparisons between stability of two typical RNA pseudoknots predicted by the present model with and without the coaxial stacking potential.

(a) MMTV pseudoknot at 1000mM [K+]; (b) BWYV pseudoknot at 500mM [K+]. Lines: df/dT, the first derivative of f with the temperature; red: predictions from the model with the coaxial stacking potential; green: predictions from the model without the coaxial stacking potential. Cartoon: the predicted 3D structures of the two pseudoknots at different temperatures. Red/black arrow: predictions from the model with coaxial stacking potential; green arrow: predictions from the model without the coaxial stacking potential.

https://doi.org/10.1371/journal.pcbi.1006222.s004

(TIF)

S4 Fig. The fractions of F, S1, S2 and U states as functions of temperature in thermal unfolding of pseudoknots at 500mM [Na+].

(a) PEMV-1, (b) BWYV, (c) PLRV, and (d) ScYLV pseudoknots. F stands for fully folded RNA; S1, hairpin intermediate with Stem1; S2, hairpin intermediate with Stem2; U, fully unfolded RNA.

https://doi.org/10.1371/journal.pcbi.1006222.s005

(TIF)

S5 Fig. The folding free energies of the constructs associated with the stems of six wild-type RNA pseudoknots and the variant of T2 pseudoknots at 1M [Na+] as a function of temperature.

(a) MMTV pseudoknot; (b) T2 and T2 variant pseudoknots; (c) PEMV-1 pseudoknot; (d) BWYV pseudoknot; (e) PLRV pseudoknot; and (f) ScYLV pseudoknot. Here, the free energies are computed using Mfold (http://unafold.rna.albany.edu/) [89].

https://doi.org/10.1371/journal.pcbi.1006222.s006

(TIF)

S6 Fig. The fractions of F, S1, S2 and U states as functions of temperature in thermal unfolding of MMTV and T2 pseudoknots in divalent ion solutions.

(a-c) MMTV pseudoknot at 50mM [K+] and different [Mg2+]’s: (a) 0.1mM [Mg2+], (b) 1mM [Mg2+], and (c) 10mM [Mg2+]. (d-f) T2 pseudoknot at 100mM [K+] and different [Mg2+]’s: (d) 0.1mM [Mg2+], (e) 1mM [Mg2+], and (f) 10mM [Mg2+].

https://doi.org/10.1371/journal.pcbi.1006222.s007

(TIF)

Acknowledgments

Parts of the numerical calculation in this work are performed on the super computing system in the Super Computing Center of Wuhan University.

Citation: Shi Y-Z, Jin L, Feng C-J, Tan Y-L, Tan Z-J (2018) Predicting 3D structure and stability of RNA pseudoknots in monovalent and divalent ion solutions. PLoS Comput Biol 14(6): e1006222. https://doi.org/10.1371/journal.pcbi.1006222

References

1. Atkins JF, Gesteland RF, Cech TR. RNA worlds: From life’s origins to diversity in gene regulation. Cold Spring Harbor Laboratory Press Cold Spring Harbor, NY; 2011.

2. Staple DW, Butcher SE. Pseudoknots: RNA structures with diverse functions. PLoS Biol. 2005; 3: 956–959.

3. Giedroc DP, Cornish PV. Frameshifting RNA pseudoknots: Structure and mechanism. Virus Res. 2009; 139: 193–208. pmid:18621088

4. Su L, Chen L, Egli M, Berger JM. Rich A. Minor groove RNA triplex in the crystal structure of a ribosomal frameshifting viral pseudoknot. Nat. Struct. Biol. 1999; 6: 285–292. pmid:10074948

5. Tinoco I, Bustamante C. How RNA folds. J. Mol. Biol. 1999; 293: 271–281. pmid:10550208

6. Liu B, Mathews DH, Turner DH. RNA pseudoknots: folding and finding. F1000 Biol. Rep. 2010; 2: 8. pmid:20495679

7. Yingling YG, Shapiro BA. The impact of dyskeratosis congenita mutations on the structure and dynamics of the human telomerase RNA pseudoknot domain. J. Biomol. Struct. Dyn. 2007; 24: 303–319. pmid:17206847

8. Gong S, Wang Y, Zhang W. Kinetic regulation mechanism of pbuE riboswitch. J. Chem. Phys. 2015; 142: 015103. pmid:25573585

9. Schlick T, Pyle AM. Opportunities and challenges in RNA structural modeling and design. Biophys. J. 2017; 113: 225–234. pmid:28162235

10. Hajdin CE, Ding F, Dokholyan NV, Weeks KM. On the significance of an RNA tertiary structure prediction. RNA 2010; 16: 1340–1349. pmid:20498460

11. Pyle AM, Schlick T. Challenges in RNA structural modeling and design. J. Mol. Biol. 2016; 428:733–735. pmid:26876599

12. Berman HM, Westbrook J, Feng Z, Gilliland G, Bhat TN, Weissig H, Shindyalov IN, Bourne PE. The Protein Data Bank. Nucleic Acids Res. 2000; 28: 235–242. pmid:10592235

13. Miao Z, Adamiak RW, Antczak M, Batey RT, Becka AJ, Biesiada M, et al. RNA-Puzzles Round III: 3D RNA structure prediction of five riboswitches and one ribozyme. RNA 2017; 23: 655–672. pmid:28138060

14. Miao Z, Westhof E. RNA structure: advances and assessment of 3D structure prediction. Annu. Rev. Biophys. 2017; 46: 483–503. pmid:28375730

15. Seetin MJ, Mathews DH. Automated RNA tertiary structure prediction from secondary structure and low-resolution restraints. J. Comput. Chem. 2011; 32: 2232–2244. pmid:21509787

16. Shapiro BA, Yingling YG, Kasprzak W, Bindewald E. Bridging the gap in RNA structure prediction. Curr. Opin. Struct. Biol. 2007; 17: 157–165. pmid:17383172

17. Sim AY, Minary P, Levitt M. Modeling nucleic acids. Curr. Opin. Struct. Biol. 2012; 22: 273–278. pmid:22538125

18. Dufour D, Marti-Renom MA. Software for predicting the 3D structure of RNA molecules. WIREs Comput. Mol. Sci. 2014;

19. Somarowthu S. Progress and current challenges in modeling large RNAs. J. Mol. Biol. 2016; 428: 736–747. pmid:26585404

20. Shi YZ, Wu YY, Wang FH, Tan ZJ. RNA structure prediction: progress and perspective. Chin. Phys. B 2014; 23: 078701.

21. Cragnolini T, Derreumaux P, Pasquali S. Ab initio RNA folding. J. Phys. Condens. Matt. 2015; 27: 233102.

22. Bida JP, Maher LJ. Improved prediction of RNA tertiary structure with insights into native state dynamics. RNA 2012; 18: 385–393. pmid:22279150

23. Jonikas MA, Radmer RJ, Laederach A, Das R, Pearlman S, Herschlag D, Altman RB. Coarse-grained modeling of large RNA molecules with knowledge-based potentials and structural filters. RNA 2009; 15:189–199. pmid:19144906

24. Parisien M, Major F. The MC-Fold and MC-Sym pipeline infers RNA structure from sequence data. Nature 2008; 452: 51–55. pmid:18322526

25. Das R, Baker D. Automated de novo prediction of native-like RNA tertiary structures. Proc. Natl. Acad. Sci. USA 2007; 104: 14664–14669. pmid:17726102

26. Kim N, Liang C, Elmetwaly S, Jung S, Curuksu J, Schlick T. Graph-based sampling for approximating global helical topologies of RNA. Proc. Nalt. Acad. Sci. USA 2014; 111: 4079–4084.

27. Li J, Zhang J, Wang J, Li W, Wang W. Structure prediction of RNA loops with a probabilistic approach. PLoS Comput. Biol. 2016; 12: e1005032. pmid:27494763

28. Boniecki M, Lach G, Dawson WK, Tomala K, Lukasz P, Soltysinski T, Rother KM, Bujnicki JM. SimRNA: a coarse-grained method for RNA folding simulations and 3D structure prediction. Nucleic Acids Res. 2016; 44: e63. pmid:26687716

29. Zhao Y, Huang Y, Gong Z, Wang Y, Man J, Xiao Y. Automated and fast building of three-dimensional RNA structures. Sci. Rep. 2012; 2: 734. pmid:23071898

30. Popenda M, Szachniuk M, Antczak M, Purzycka KJ, Lukasiak P, Bartol N, Blazewicz J, Adamiak RW. Automated 3D structure composition for large RNAs. Nucleic Acids Res. 2012; 40: e112. pmid:22539264

31. Zhang J, Dundas J, Lin M, Chen M, Wang W, Liang J. Prediction of geometrically feasible three-dimensional structures of pseudoknotted RNA through free energy estimation. RNA 2009; 15: 2248–2263. pmid:19864433

32. Yao J, Reinharz V, Major F, Waldispuhl J. RNA-MoIP: prediction of RNA secondary structure and local 3D motifs from sequence data. Nucleic Acids Res. 2017; 45: W440–W444. pmid:28525607

33. Bernauer J, Huang X, Sim AYL, Levitt M. Fully differentiable coarse-grained and all-atom knowledge-based potentials for RNA structure evaluation. RNA 2011; 17: 1066–1075. pmid:21521828

34. Capriotti E, Norambuena T, Marti-Renom MA, Melo F. All-atom knowledge-based potential for RNA structure prediction and assessment. Bioinfomatics 2011; 27: 1086–1093.

35. Wang J, Mao K, Zhao Y, Zeng C, Xiang J, Xiao Y. Optimization of RNA 3D structure prediction using evolutionary restraints of nucleotide-nucleotide interactions from direct coupling analysis. Nucleic Acids Res. 2017; 45: 6299–6309. pmid:28482022

36. Wang J, Zhao Y, Wang J, Xiao Y. Computational study of stability of an H-H-type pseudoknot motif. Phys. Rev. E 2015; 92: 062705.

37. Cao S, Chen SJ. Predicting RNA pseudoknot folding thermodynamics. Nucleic Acids Res. 2006; 34: 2634–3652. pmid:16709732

38. Cao S, Chen SJ. Predicting structures and stabilities for H-type pseudoknots with interhelix loops. RNA 2009; 15: 696–706. pmid:19237463

39. Ding F, Sharma S, Chalasani P, Demidov VV, Broude NE, Dokholyan NV. Ab initio RNA folding by discrete molecular dynamics: from structure prediction to folding mechanisms. RNA 2008; 14: 1164–1173. pmid:18456842

40. Cragnolini T, Laurin Y, Derreumaux P, Pasquali S. Coarse-grained HiRE-RNA model for ab initio RNA folding beyond simple molecules, including noncanonical and multiple base pairings. J. Chem. Theory and Comput. 2015; 11: 3510–3522.

41. Uusitalo JJ, Ingolfsson HI, Marrink SJ, Faustino I. Martini coarse-grained force field: extension to RNA. Biophys. J. 2017; 113: 246–256. pmid:28633759

42. Sulc P, Romano F, Ouldridge TE, Doye JPK, Louis AA. A nucleotide-level coarse-grained model of RNA. J. Chem. Phys. 2011; 140: 235102.

43. Bell DR, Cheng SY, Salazar H, Ren P. Capturing RNA folding free energy with coarse-grained molecular dynamics simulations. Sci. Rep. 2017; 7: 45812. pmid:28393861

44. Xia Z, Bell DR, Shi Y, Ren, P. RNA 3D structure prediction by using a coarse-grained model and experimental data. J. Phys. Chem. B 2013; 117: 3135–3144. pmid:23438338

45. Zhou HX. Theoretical frameworks for multiscale modeling and simulation. Curr. Opin. Struct. Biol. 2014; 25: 67–76. pmid:24492203

46. Denesyuk NA, Thirumalai D. Coarse-grained model for predicting RNA folding thermodynamics. J. Phys. Chem. B 2013; 117: 4901–4911. pmid:23527587

47. Cao S, Chen SJ. Physics-based de novo prediction of RNA 3D structures. J. Phys. Chem. B 2011; 115: 4216–4226. pmid:21413701

48. Woodson SA. Metal ions and RNA folding: a highly charged topic with a dynamic future. Curr. Opin. Struct. Biol. 2005; 9: 104–109.

49. Pabit SA, Sutton JL, Chen H, Pollack L. Role of ion valence in the submillisecond collapse and folding of a small RNA domain. Biochemistry 2013; 52: 1539–1546. pmid:23398396

50. Sun LZ, Zhang D, Chen SJ. Theory and modeling of RNA structure and interactions with metal ions and small molecules. Annu. Rev. Biophys. 2017; 46: 227–246. pmid:28301768

51. Tan ZJ, Chen SJ. Predicting ion binding properties for RNA tertiary structures. Biophys. J. 2010; 99: 1565–1576. pmid:20816069

52. Tan ZJ, Chen SJ. Salt contribution to RNA tertiary structure folding stability. Biophys. J. 2011; 101: 176–187. pmid:21723828

53. Lipfert J, Doniach S, Das R, Herschlag D. Understanding nucleic acid-ion interactions. Annu. Rev. Biochem. 2014; 83: 19.1–19.29.

54. Soto AM, Misra V, Draper DE. Tertiary structure of an RNA pseudoknot is stabilized by “diffuse” Mg2+ ions. Biochemistry 2007; 46: 2973–2983. pmid:17315982

55. Tan ZJ, Chen SJ. RNA helix stability in mixed Na+/Mg2+ solution. Biophys. J. 2007; 92: 3615–3632. pmid:17325014

56. Wu YY, Zhang ZL, Zhang JS, Zhu XL, Tan ZJ. Multivalent ion-mediated nucleic acid helix-helix interactions: RNA versus DNA. Nucleic Acids Res. 2015; 43: 6156–6165. pmid:26019178

57. Zhang ZL, Wu YY, Xi K, Sang JP, Tan ZJ. Divalent ion-mediated DNA-DNA interactions: a comparative study of triplex and duplex. Biophys. J. 2017; 113: 517–528. pmid:28793207

58. Denesyuk NA, Thirumalai D. Crowding promotes the switch from hairpin to pseudoknot conformation in human telomerase RNA. J. Am. Chem. Soc. 2011; 133: 11858–11861. pmid:21736319

59. Cho SS, Pincus DL, Thirumalai D. Assembly mechanisms of RNA pseudoknots are determined by the stabilities of constituent secondary structures. Proc. Nat. Acad. Sci. USA 2009; 106: 17349–17354. pmid:19805055

60. Hayes RL, Noel JK, Mandic A, Whitford PC, Sanbonmatsu KY, Mohanty U, Onuchic JN. Generalized Manning condensation model captures the RNA ion atmosphere. Phys. Rev. Lett. 2015; 114: 258105. pmid:26197147

61. Shi YZ, Wang FH, Wu YY, Tan ZJ. A coarse-grained model with implicit salt for RNAs: predicting 3D structure, stability and salt effect. J. Chem. Phys. 2014; 141: 105102. pmid:25217954

62. Shi YZ, Jin L, Wang FH, Zhu XL, Tan ZJ. Predicting 3D structure, flexibility, and stability of RNA hairpins in monovalent and divalent ion solutions. Biophys. J. 2015; 109: 2654–2665. pmid:26682822

63. Wang FH, Wu YY, Tan ZJ. Salt contribution to the flexibility of single-stranded nucleic acid of finite length. Biopolymers 2013; 99: 370–381. pmid:23529689

64. Zhang X, Bao L, Wu YY, Zhu XL, Tan ZJ. Radial distribution function of semiflexible oligomers with stretching flexibility. J. Chem. Phys. 2017; 147: 054901. pmid:28789545

65. Xia T, SantaLucia J, Burkand ME, Kierzek R, Schroeder SJ, Jiao X, Cox C, Turner DH. Thermodynamic parameters for an expanded nearest-neighbor model for formation of RNA duplexes with Watson-Crick base pairs. Biochemistry 1998; 37: 14719–14735. pmid:9778347

66. Wang Y, Gong S, Wang Z, Zhang W. The thermodynamics and kinetics of a nucleotide base pair. J. Chem. Phys. 2016; 144: 115101. pmid:27004898

67. Zhang Y, Zhou H, Ouyang Z. Stretching single-stranded DNA: Interplay of electrostatic, base-pairing, and base-pair stacking interactions. Biophys. J. 2001; 81: 1133–1143. pmid:11463654

68. Manning GS. The molecular theory of polyelectrolyte solutions with applications to the electrostatic properties of polynucleotides. Q. Rev. Biophys. 1978; 11: 179–246. pmid:353876

69. Tan ZJ, Chen SJ. Salt dependence of nucleic acid hairpin stability. Biophys. J. 2008; 95: 738–752. pmid:18424500

70. Parisien M, Cruz JA, Westhof E, Major F. New metrics for comparing and assessing discrepancies between RNA 3D structures and models. RNA 2009; 15: 1875–1885. pmid:19710185

71. Kim YG, Su L, Maas S, O’Neill A, Rich A. Specific mutations in a viral RNA pseudoknot drastically change ribosomal frameshifting efficiency. Proc. Nalt. Acad. Sci. USA 1999; 96: 14234–14239.

72. Nixon PL, Giedroc DP. Energetics of a strongly pH dependent RNA tertiary structure in a frameshifting pseudoknot. J. Mol. Biol. 2000; 296: 659–671. pmid:10669615

73. Nixon PL, Cornish PV, Suram SV, Giedroc DP. Thermodynamic analysis of conserved loop-stem interactions in P1-P2 frameshifting RNA pseudoknots from plant luteoviridae. Biochemistry 2002; 41: 10665–10674. pmid:12186552

74. Cornish PV, Henning M, Giedroc DP. A loop 2 cytidine-stem 1 minor groove interaction as a positive determinant for pseudoknot-stimulated -1 ribosomal frameshifting. Proc. Natl. Acad. Sci. USA 2005; 102: 12694–12699. pmid:16123125

75. Theimer CA, Giedroc DP. Contribution of the intercalated adenosine at the helical junction to the stability of the gag-pro frameshifting pseudoknot from mouse mammary tumor virus. RNA 2000; 6: 409–421. pmid:10744025

76. Narayanan R, Velmurugu Y, Kuznetsov SV, Ansari A. Fast folding of RNA pseudoknots initiated by laser temperature-jump. J. Am. Chem. Soc. 2011; 133: 18767–18774. pmid:21958201

77. Nixon PJ, Giedroc DP. Equilibrium unfolding (folding) pathway of a model H-type pseudoknotted RNA: the role of magnesium ions in stability. Biochemistry 1998; 37: 16116–16129. pmid:9819204

78. Cao S, Giedroc DP, Chen SJ. Predicting loop-helix tertiary structural contacts in RNA pseudoknots. RNA 2010; 16: 538–552. pmid:20100813

79. Leipply D, Draper DE. Effects of Mg2+ on the free energy landscape for folding a purine riboswitch RNA. Biochemistry 2011; 50: 2790–2799. pmid:21361309

80. Chadalavada DM, Senchak SE, Bevilacqua PC. The folding pathway of the genomic hepatitis delta virus ribozyme is dominated by slow folding of the pseudoknots. J. Mol. Biol. 2002; 4: 559–575.

81. White KH, Orzechowski M, Fourmy D, Visscher K. Mechanical unfolding of the beet western yellow virus -1 frameshift signal. J. Am. Chen. Soc. 2011; 133: 9775–9782.

82. Zhang X, Zhang D, Zhao C, Tian K, Shi R, Du X, Burcke AJ, Wang J, Chen SJ, Gu LQ. Nanopore electric snapshots of an RNA tertiary folding pathway. Nat. Commun. 2017; 8: 1458. pmid:29133841

83. Zhang X, Xu X, Yang Z, Burcke AJ, Gates KS, Chen SJ, Gu LQ. Mimicking ribosomal unfolding of RNA pseudoknot in a protein channel. J. Am. Chem. Soc. 2015; 137: 15742–15752. pmid:26595106

84. Biyun S, Cho SS, Thirumalai D. Folding of human telomerase RNA pseudoknot using ion-jump and temperature-quench simulations. J. Am. Chem. Soc. 2011; 133: 20634–20643. pmid:22082261

85. Bian Y, Zhang J, Wang J, Wang J, Wang W. Free energy landscape and multiple folding pathways of an H-Type RNA pseudoknot. PLoS One 2015; 10: e0129089. pmid:26030098

86. Hori N, Denesyuk NA, Thirumalai D. Salt effects on the thermodynamics of a frameshifting RNA pseudoknot under tension. J. Mol. Biol. 2016; 428: 2847–2859. pmid:27315694

87. Zhang L, Bao P, Leibowitz MJ, Zhang Y. Slow formation of a pseudoknot structure is rate limiting in the productive co-transcriptional folding of the self-splicing Candida intron. RNA 2009; 15: 1986–1992. pmid:19710184

88. Zhang Y, Zhang J, Wang W. Atomistic analysis of pseudoknotted RNA unfolding. J. Am. Chen. Soc. 2011; 133: 6882–6885.

89. Zuker M. Mfold web server for nucleic acid folding and hybridization prediction. Nucleic Acids Res. 2003; 31: 3406–3415. pmid:12824337

90. Bao L, Zhang X, Shi YZ, Wu YY, Tan ZJ. Understanding the relative flexibility of RNA and DNA duplexes: stretching and twist-stretch coupling. Biophys. J. 2017; 112: 1094–1104. pmid:28355538

91. Bailor MH, Mustoe AM, Brooks CL, Al-Hashimi HM. Topological constraints: using RNA secondary structure to model 3D conformation, folding pathways, and dynamic adaptation. Curr. Opin. Struct. Biol. 2011; 21: 296–305. pmid:21497083

92. Williams B, Zhao B, Tandon A, Ding F, Weeks KM, Zhang Q, Dokholyan NV. Structure modeling of RNA using sparse NMR constraints. Nucleic Acids Res. 2017; 45: 12638–12647. pmid:29165648

93. Magnus M, Matelska D, Lach G, Chojnowski G, Boniecki M, Purta E, et al. Computational modeling of RNA 3D structures with the aid of experimental restraints. RNA Biol. 2014; 11: 522–536. pmid:24785264

94. Sim AYL, Levitt M. Clustering to identify RNA conformations constrained by secondary structure. Proc. Nalt. Acad. Sci. USA 2011; 108: 3590–3595.

95. Hyeon C, Thirumalai D. Force-unfolding and force-quench refolding of RNA hairpins. Biophys. J. 2006; 90: 3410–3427. pmid:16473903

96. Chen G, Chang KY, Chou MY, Bustamante C, Tinoco I. Triplex structures in an RNA pseudoknot enhance mechanical stability and increase efficiency of –1 ribosomal frameshifting. Proc. Nalt. Acad. Sci. USA 2009; 106: 12706–12711.

97. Green L, Kim CH, Bustamante C, Tinoco I. Characterization of the mechanical unfolding of RNA pseudoknots. J. Mol. Biol. 2008; 375: 511–528. pmid:18021801