MD simulation studies to investigate iso-energetic conformational behaviour of modified nucleosides m2G and m2 2G present in tRNA

Modified nucleic acid bases are most commonly found in tRNA. These may contain modifications from simple methylation to addition of bulky groups. Methylation of the four canonical nucleotide bases at a wide variety of positions is particularly prominent among the known modification. Methylation of N2 group of guanine is a relatively common modification in tRNA and rRNA. N2-methylguanosine (m2G) is the second most often encountered nucleoside in E. coli tRNAs. N2, N2- dimethylguanosine (m2 2G) is found in the majority of eukaryotic tRNAs and involved in forming base pair interactions with adjacent bases. Hence, in order to understand the structural significance of these methylated nucleic acid bases we have carried out molecular dynamics simulation to see the salvation effect. The results obtained shows iso-energetic conformational behaviors for m2G and m2 2G. The simulation trajectory of m2G shows regular periodical fluctuations suggesting that m2G is equally stable as either s-cis or s-trans rotamers. The two rotamers of m2G may interact canonically or non-canonically with opposite base as s-trans m2G26:C/A/U44 and s-cis m2G26:A/U44. The free rotations around the C-N bond could be the possible reason for these iso-energetic conformations. Dimethylation of G has almost no influence on base pairing with either A or U. Thus, these results reveal that modified nucleosides m2G and m2 2G may play an important role to prevent tRNA from adopting the unusual mitochondrial like conformation.


Introduction
RNA molecules undergo extensive post-transcriptional modifications that are important for their biological activities. Posttranscriptional modifications have been known as a natural mechanism to provide structural stability across the wide range of temperature in archaea as well as bacteria [1]. Transfer RNAs have the largest number and the greatest diversity of modifications: base or ribose methylation, base isomerization, base reduction, base thiolation and more complex hypermodifications [2,3]. An important characteristic of tRNA is the presence of high content of modified nucleosides of which methylation represents the principle posttranscriptional modification during its maturation. In the maturation process of tRNA, transfer of methyl group occurs at polynucleotide level through an S-adenosyl-L-methionine donor, resulting in modification of heterocyclic base, the ribose moiety, or both [4]. The family of structurally related nucleosides m 2 G, m 2 2G, m 2 Gm and m 2 2Gm, are from known archaeal tRNA sequences. These modified nucleosides are conserved at only two locations, position 10 first base in the proximal position of the dihydrouridine (DHU) arm and at position 26, junction between the D-stem and the anticodon stem, where they play crucial roles in the control and stabilization of the tertiary L fold structure of the tRNA [5,6]. The m 2 G and m 2 2G modifications in tRNA are found not only at position 26 but also at positions 6, 7, 9, 10, 18 and 27 in various organisms [7].
Experimentally, it has been found that level of certain modified nucleosides in archaeal thermophiles play major stabilizing role beyond the effects of magnesium ion binding and G-C content of tRNA [8]. Earlier study involving three-dimensional models of yeast tRNA Phe derived from X-ray crystallographic data implies that m 2 2G26 functions as a molecular hinge.
This hinge adjusts the angular position of the D-stem and the anticodon stem during protein synthesis, thus maintaining a certain rigidity/flexibility in this part of tRNA [9]. Nuclear magnetic resonance studies on the resonance of the methyl proton in yeast tRNA Phe also provide evidence to support the notation that m 2 2G26 has a significant role in regulating the stacking and conformational dynamics of this region of tRNA molecule [10]. The yeast tRNA (m 2 2G26) methyltransferase is dependent on the D-stem sequence and size of variable loop for the synthesis of N 2 -N 2 dimethyl guanosine at 26 th position [11]. Mutations were introduced in both the D-stem and the variable loop of tRNA Asp to obtain dimethylation of the normally unmodified G26 by the yeast N 2 , N 2 -dimethyl G26methyltransferase [12]. The presence of m 2 2G26 in cytosolic tRNA may avert the molecule from adopting an unusual mitochondrial tRNA pattern folding and instead, allow it to fold into the canonical cloverleaf model. Through screening of the tRNA sequence and gene database it was revealed that some cytosolic tRNAs have the potential to fold into alternate structures. It was further noted that when a tRNA had the potential for this alternate folding, m 2 2G was found at position 10 and 26 presumably to block the formation of this nonstandard folding pattern [13]. The methylated guanosine from 26 position of tRNA may have role in regulating the stacking interactions and the conformational dynamics [14].
N 2 -methylguanosine is found in both helical and loop regions of RNA secondary structure [15,16] and it can exist in either s-cis or strans rotamers [18]. Incorporation of m 2 G was found to be isoenergetic with G in the duplex context as well as in GNRA (N = any CSBJ Abstract: Modified nucleic acid bases are most commonly found in tRNA. These may contain modifications from simple methylation to addition of bulky groups. Methylation of the four canonical nucleotide bases at a wide variety of positions is particularly prominent among the known modification. Methylation of N2 group of guanine is a relatively common modification in tRNA and rRNA. N 2 -methylguanosine (m 2 G) is the second most often encountered nucleoside in E. coli tRNAs. N 2 , N 2dimethylguanosine (m 2 2G) is found in the majority of eukaryotic tRNAs and involved in forming base pair interactions with adjacent bases. Hence, in order to understand the structural significance of these methylated nucleic acid bases we have carried out molecular dynamics simulation to see the salvation effect. The results obtained shows iso-energetic conformational behaviors for m 2 G and m 2 2G. The simulation trajectory of m 2 G shows regular periodical fluctuations suggesting that m 2 G is equally stable as either s-cis or s-trans rotamers. The two rotamers of m 2 G may interact canonically or non-canonically with opposite base as s-trans m 2 G26:C/A/U44 and s-cis m 2 G26:A/U44. The free rotations around the C-N bond could be the possible reason for these isoenergetic conformations. Dimethylation of G has almost no influence on base pairing with either A or U. Thus, these results reveal that modified nucleosides m 2 G and m 2 2G may play an important role to prevent tRNA from adopting the unusual mitochondrial like conformation.
nucleotide and R = purines) tetra loops [17]. The two rotamers of m 2 G, s-cis and s-trans have been found equally stable in RNA duplex [17] and in tRNA [18]. Free rotations around C-N bond aids m 2 G to get energetically two stable conformations and to form base pair interactions such as s-cis m 2 G:A/U or s-trans m 2 G:C/A/U [18].
This phenomenon has not been studied in detail at atomic level. Hence, present attempt has been made to investigate the dynamic behavior of iso-energetic conformations of m 2 G/m 2 2G using MD simulation technique. The results clearly show that the iso-energetic nature of s-trans and s-cis conformations could be because of free rotations around the C-N bond. Before this various computational techniques have been used to understand the conformational behavior and dynamics of many complex modified nucleosides [18][19][20][21][22][23]. However, all-atom molecular dynamics simulations of entire solvated ribosome, mRNA and tRNA complex have been studied to find out motion of tRNA from the A/T state into the A site [24].

Computational details
Molecular dynamics (MD) simulations were performed using Amber 10 simulation suite on (HP ProLiant-ML150G6) server in order to highlight the influence of explicit solvation on the conformation of modified nucleosides N 2 -methylguanosine (m 2 G) and N 2 -N 2 dimethyl guanosine (m 2 2G). PCILO predicted preferred conformations of m 2 G and m 2 2G [18] were used as starting geometries for MD simulation studies. Antechamber software was used to calculate charges. Each nucleoside was solvated by 729 SPC/E water molecules filling a 34.20 × 31.05 × 31.85 Å 3 rectilinear box with 1.0 water density [25]. Simulations were performed under periodic boundary conditions by employing the Particle Mesh Ewald [26] method to calculate long-range electrostatic interactions. MD trajectories were propagated at 2.0 fs time step using the shake algorithm [27] to all hydrogen atoms with non-bonded cutoff of 10 Å. The non-bonded pair list was updated by every 10 steps. The trajectories were calculated by maintaining constant temperature (300 K) and constant pressure (1atm) at 2 fs time step according to Berendsen coupling algorithm [28].
An equilibration protocol similar to the earlier molecular dynamics simulation study of nucleic acids was applied [29,30]. The equilibration protocol consisted of 5000 steps of steepest descent minimization followed by 5 ps of MD at 300 K applied to relaxation of initial strain present between water molecules and N 2 -methyl derivatives of guanosine. In the next step N 2 -methyl derivatives were fixed while water molecules were allowed to relax at 100 K (1 ps), 200 K (1 ps), and finally at 300 K for 198 ps, thus equilibration protocol was completed at 200 ps.
Equilibrated system was further subjected to 5000 steps of steepest descent minimization to remove bad contacts between water molecules and nucleic acid bases. In further steps of MD simulation, no positional constrains were applied to the system and the temperature was progressively increased to 300 K in steps of 50 K with 1 ps at each step. Finally system was subjected to production MD of 20 ns at 300 K temperature and constant pressure (1atm) with fully solvated and neutralized system. PTRAJ module of Amber Tool 10 was used for analysis of average structures [31].

Results and Discussion
The preferred conformation of N 2 -methylguanosine [18] ( Fig.  2A) has been used as a starting geometry for 20 ns molecular dynamics simulation study. In order to confirm the iso-energetic conformational behavior of m 2 G we have analyzed four different average structures particularly at 0 to 1 ns (Fig. 2B), 3 to 4.5 ns ( The average structure obtained at 0 to 1 ns (Fig. 2B) maintains the initial geometry ( Fig. 2A) [18] by preserving s-trans or "proximal" conformation for N 2 -methyl substituent of guanosine (26 th ), which is stabilized by hydrogen bonding interaction between N(3)…HC(10) (Fig. 2B  . A large deviation around the torsion angle β is due to rotations around C-N bond. Next average structure taken at 3 to 4.5 ns (Fig. 2C) is also stabilized by hydrogen bonding between N(3)…HC(10) along with this, interaction between N(3)…HC2′ ( Fig. 2C and Table 1) provides an additional structural stability to this average structure (Fig. 2C), as observed in our earlier conformational study of m 2 G [18]. The torsion angle β shows large deviation (120°) whereas α maintains starting value as compared with initial structure of m 2 G [18].
Second snapshot structure ( Fig. 2E and Table 1)  The geometrical parameters for torsion angles and hydrogen bonding interactions analyzed from average and snapshot structures are given in table 1. The snapshot structure for m 2 G (Fig. 3A) taken at 2 ns prefers s-cis conformation due to change in α torsion angle which deviates from 180° to 29°, while other torsion angles β and χ diverges to great extent from initial structure and are in close agreement with crystal conformer 1OB5.pdb [34]. This structure is stabilized by the hydrogen bond between N(3)…HC2′ (Table 1) as found in earlier conformational study of m 2 G [18]. Average structure (Fig. 3B) chosen for the period 5 to 11 ns when α torsion angle flipped by 179° as compared to preferred structure of m 2 G ( Fig. 2A).
The obtained average structure maintains distal conformation for methyl substituent of guanosine and gets stabilized by N(3)…HC1′ interaction which was not observed in s-trans conformer of m 2 G. Last average structure (Fig. 3C) was taken within the range of 13 to 19 ns, showing s-cis conformation of m 2 G. The s-cis conformation is obtained due to change in α torsion angle from 180° to 355°. This average structure (Fig. 3C) shows deviations for torsion angle β by 116° and χ by 73°. Obtained average structure (Fig. 3C) shows similar values for torsion angle α and χ as compared with crystal structure 1OB5.pdb [34]. A large deviation around the torsion angle β is due to fluctuations from s-trans to s-cis conformation by rotating C-N bond of methyl group.
The average structure obtained from 2 to 3 ns (Fig. 5B) prefers distal conformation for m 2 2G and prevent Watson-Crick base pairing with C, instead it would prefer non-canonical Watson-Crick interactions to pair with A/U44. Compared with crystal conformer (1EHZ.pdb) [33], average structure retains quite similar torsion angle values for α, β and χ. The average structure selected at last 1000 ps (19 to 20 ns) does not show much difference as compared with earlier average structure. The only difference between these two average structures is variation around α torsion angle which is positioned to 99° (Fig. 5C) from its preferred value (Fig. 5A). Due to this small change in conformational property of m 2 2G, dimethylation of guanosine has almost no influence on pairing with either A or U, because the N 2 position of guanosine has no impact on these base pairing interactions. The m 2 2G-A pair can be formed with little hindrance, because, even though the methyl groups of m 2 2G26 are in the plane of the base, as in case of yeast tRNA Phe , they take part in a propeller-type arrangement with the base. The m 2 2G-U pair would not be affected by any conformational arrangement of the methyl groups similarly as discussed in [18].  The α torsion angle fluctuates periodically in between ±180° or ± 60° (Fig. 6A and Table 1) over the 20 ns molecular dynamics simulation period suggesting free rotation around C-N bond.

Conclusion
The regular periodical fluctuations around the bond C(2)-N(2) of m 2 G was observed throughout the 20 ns molecular dynamics simulation, which confirms the existence of iso-energetic s-cis or strans rotamers of m 2 G. These iso-energetic rotamers interconvert  easily during the simulation period. These results are in favor with preferred and alternative conformations of m 2 G obtained by our earlier conformational energy calculations [18] as well as crystal structure (1EHZ.pdb [33] and 1OB5.pdb [34]). The periodical fluctuations of s-trans to s-cis and vice versa could be possible due to free rotations around the C-N bond of methyl group. According to tRNA sequence analysis [18] and this MD simulation results we would like to say that m 2 G26 can form three different canonical as well as non-canonical Watson-Crick base pairing interactions with other bases. Such base pairing may be summarized as i) an usual Watson-Crick base pairing of m 2 G26-C44 where the methyl substituent must be in s-trans orientation, ii) non Watson-Crick base pairing between m 2 G26-A where the methyl substituent is likely to be s-cis orientation, and iii) non Watson-Crick m 2 G26-U base pairing where the methyl group can adopt one of them, i.e. s-cis or s-trans conformation. These results reveal that m 2 G is equally stable as either the s-cis or s-trans rotamers and the rotational preference of methyl group may be specific to the sequence context reliant upon which face of the base contributes in hydrogen bonding. Thus, MD simulation results confirm that the N 2 -methyl group of m 2 G26 may prefer energetically two stable rotamers, i.e., s-trans m 2 G26:C/A/U44 and s-cis m 2 G26:A/U44 as found in earlier results [18].
Similarly, the presence of two methyl groups unlike in case of single methyl in m 2 G virtually eliminates the possibility of pairing with C and, indeed, m 2 2G26 pairs exclusively with A or U at position 44 and is flanked by C27:G43 on one side and the m 2 G10-C25-G45 triple on the other [35]. Hence, these results suggest that the modified nucleosides m 2 G26 and m 2 2G26 play an important role in tRNA folding and may prevent tRNA from adopting the unusual mitochondrial like conformation.