Molecular dynamics simulations of RNA: an in silico single molecule approach

Abstract

RNA molecules are now known to be involved in the processing of genetic information at all levels, taking on a wide variety of central roles in the cell. Understanding how RNA molecules carry out their biological functions will require an understanding of structure and dynamics at the atomistic level, which can be significantly improved by combining computational simulation with experiment. This review provides a critical survey of the state of molecular dynamics (MD) simulations of RNA, including a discussion of important current limitations of the technique and examples of its successful application. Several types of simulations are discussed in detail, including those of structured RNA molecules and their interactions with the surrounding solvent and ions, catalytic RNAs, and RNA-small molecule and RNA-protein complexes. Increased cooperation between theorists and experimentalists will allow expanded judicious use of MD simulations to complement conceptually related single molecule experiments. Such cooperation will open the door to a fundamental understanding of the structure-function relationships in diverse and complex RNA molecules. .

FIGURE 1

Limited accuracy of a nonpolarizable force field for cation–solute interactions. The dependence of the interaction energy between O6(G) and K⁺ in a G-DNA like geometry. Black, reference QM data; blue and red, standard CHARMM (radius 1.76375 A and well depth 0.087 kcal/mol⁻¹) and AMBER (radius 2.6580 A and well depth 0.000328 kcal/mol⁻¹) parameters, respectively. The cation in this particular interaction appears too large while the binding energy is clearly underestimated. Note that the discrepancy originates from the lack of polarization in the force field and cannot be fully overcome within the framework of pair additive (nonpolarizable) potentials. For example, reduction of the cation radius could worsen the calculated ion solvation energy.

FIGURE 2

Stability of a 500 ps simulation of a fully neutralized and solvated tRNA molecule. (A) Yeast tRNA^Asp secondary structure. (B) Stereo view of snapshots of the tRNA backbone, taken at 20 ps intervals. (C) Time course of RMS deviations from the starting structure. Reproduced with permission from Auffinger, P.; Louise-May, S.; Westhof, E. Biophys J 1999, 76, 50–64, © Biophysical Society.

FIGURE 3

Structural and dynamic signatures of 5S ribosomal (r)RNA loop E. (A) Base pairing pattern (using classification from Ref. 119) of (left) the bacterial Loop E X-ray structure and (right) the spinach chloroplast Loop E architecture predicted from isostericity rules and MD simulations. The mutations from one sequence to the other are highlighted by red boxes. (B) MD simulations of Helix IV Loop E 5S rRNA in complex with ribosomal protein L25 reveal unique hydration sites with single water molecules bound throughout 25 ns trajectories. (C) The deep major groove of Loop E forms a unique pocket of negative electrostatic potential heavily occupied by monovalent cations (green balls). In summary, MD characterizes Loop E as a rigid RNA segment with unique shape and prominent ion-binding properties.

FIGURE 4

V-shaped Kink-turns (K-turns) are among the most recurrent RNA motifs. MD simulations predict (top) that K-turns are uniquely flexible elbow-like RNA building blocks, where subtle local conformational changes in the kink area propagate as large scale motions towards the attached helical stems. The local dynamics associated with K-turn flexibility are due to dynamical insertion of long-residency waters between the C and A nucleotides of the A-minor type I interaction between the two helical arms (middle). Large-scale elbow-like dynamics are observed in a simulation of the Helix 42–44 RNA portion of the large ribosomal subunit upon dynamical water insertion into the A-minor interaction in the universally conserved K-turn 42 of Helix 42 (bottom). The MD results appear to agree with cryo-EM data that show large-scale dynamics of this rRNA segment (the RNA part of the L7/L12 stalk) during the elongation cycle. Reproduced with permission from Razga, F.; Zacharias, M.; Reblova, K.; Koca, J.; Sponer, J. Structure 2006, 14, 825–835, © Cell Press.

FIGURE 5

The central pocket of the HIV-1 DIS RNA kissing complex is permanently occupied by 2–3 delocalized monovalent cations (green), smoothly exchanging with bulk solvent on a timescale of a few nanoseconds.

FIGURE 6

An example trajectory showing diverse conformational sampling in the folding of a small RNA hairpin. Na⁺ ions are shown in green, blue and red arrows indicate native and non-native base pairing. The black arrow in the last frame shows a hydrated ion bound to a site also seen in simulations of the native structure. Reproduced with permission from Sorin, E. J.; Rhee, Y. M.; Pande, V. S. Biophys J 2005, 88, 2516–2524, © Biophysical Society.

FIGURE 7

MD simulations of pre- and post-cleavage HDV ribozyme constructs showing hydrogen bonding patterns consistent with a general base mechanism. (A) Sequence and secondary structure of the simulated genomic HDV ribozyme with structural elements color-coded. The product form lacks U-1. The open arrow indicates the cleavage site. (B) Overlap of the crystal structures of the precursor (color-coded as in A) and product (silver) with key nucleotides indicated. Bottom panels: Overlay of representative averaged structures of the catalytic pocket from precursor simulations (color-coded as in A; cyan spheres, Na⁺ ions; broken lines, hydrogen bonds and inner-sphere ion contacts) with the precursor crystal structure (gray; yellow sphere, crystallographically resolved Mg ion). (C) Simulation with an unprotonated C75, where the U-1(O2′)–C75(N3) hydrogen bond necessary for general base activity forms. (D) Simulation with a protonated C75, where the C75H⁺(N3)–G1(O5′) hydrogen bond necessary for general acid activity does not form. Reproduced with permission from Krasovska, M. V.; Sefcikova, J.; Spackova, N.; Sponer, J.; Walter, N. G. J Mol Biol 2005, 351, 731–748, © Academic Press.

FIGURE 8

Possible catalytic role of intracavity water in the hairpin ribozyme. (a) The initial active site geometry is that of the crystal structure (gray), but over the course of the MD simulations, the 2′OH shifts from bulk solvent to the site of chemistry and forms a hydrogen bond with the asterisked water molecule (color). (b) The electrostatic potential map shows a minimum of −51 kT/e near the asterisked water molecule. Reproduced with permission from Rhodes, M.M., Reblova, K., Sponer, J., Walter, N.G. Proc Natl Acad Sci USA 103, 2006, 13381–13385, © National Academy of Sciences.

FIGURE 9

(A) Experimental and calculated hydration sites in the binding of paromomycin to the rRNA A-site. Red sites are found in both the crystal structure and the simulation, while gray sites are found only in the simulation. Sites are numbered according to the peak height of the calculated density, with 1 as the highest peak. (B) Stereo view of the experimental hydration sites (red) and the closest corresponding sites from the simulation (gray). Reproduced with permission from Vaiana, A. C.; Westhof, E.; Auffinger, P. Biochimie 2006, 88, 1061–1073, © Editions Scientifiques Elsevier.

FIGURE 10

Structures from MD simulations investigating the role of conserved aromatic amino acids in U1A–RNA binding. These simulations provide a model to explain the stabilization of a phenylalanine to alanine mutant protein by a modified base (A-4CPh). This study predicts an extended conformation for the modified base when bound to the WT protein, and a folded conformation when bound to the mutant protein. (A) WT protein and WT RNA show stacking. (B) WT protein with a folded modified base interrupts stacking. (C) WT protein with an extended modified base preserves stacking. (D) Mutant protein with a folded modified RNA preserves stacking. (E) Mutant protein with an extended modified base interrupts stacking. Reproduced with permission from Zhao, Y.; Kormos, B. L.; Beveridge, D. L.; Baranger, A. M. Biopolymers 2006, 81, 256–269, © Wiley Interscience.

FIGURE 11

Four stages of accommodation during tRNA movement from the A/T to the A/A site in a targeted all-atom MD simulation. (a) tRNA interaction regions colored by accommodation stage. (b) Time evolution of parameters describing tRNA deformation. (c) Snapshots from each stage of accommodation. Reproduced with permission from Sanbonmatsu, K. Y.; Joseph, S.; Tung, C. S. Proc Natl Acad Sci USA 2005, 102, 15854–15859, © National Academy of Sciences.