Molecular dynamics simulation study on Bacillus subtilis EngA: the presence of Mg2+ at the active-sites promotes the functionally important conformation

Abstract EngA, a GTPase contains two GTP binding domains [GD1, GD2], and the C-terminal KH domain shown to be involved in the later stages of ribosome maturation. Association of EngA to the ribosomal subunit in the intermediate stage of maturation is essential for complete ribosome maturation. However, this association was shown to be dependent on the nucleotide bound combinations. This nucleotide dependent association tendency is attributed to the conformational changes that occur among different nucleotide bound combinations. Therefore, to explore the conformational changes, all-atom molecular dynamics simulations for Bacillus subtilis EngA in different nucleotide bound combinations along with the presence or absence of Mg2+ in the active-sites were carried out. The presence of Mg2+ along with the bound nucleotide at the GD2 active-site dictates the GD2-Sw-II mobility, but the GD1-Sw-II mobility has not shown any nucleotide or Mg2+ dependent movement. However, the GD1-Sw-II secondary conformations are shown to be influenced by the GD2 nucleotide bound state. This allosteric connection between the GD2 active-site and the GD1-Sw-II is also observed through the dynamic network analysis. Further, the exploration of the GD1-KH interface interactions exhibited a more attractive tendency when GD1 is bound to GTP-Mg2+. In addition, the presence of Mg2+ stabilizes active-site water and also increases the distances between the α- and γ- phosphates of the bound GTP. Curiously, three water molecules in the GD1 active-site and only one water molecule in the GD2 active-site are stabilized. This indicates that the probability of GTP hydrolysis is more in GD1 compared to GD2. Communicated by Ramaswamy H. Sarma


Introduction
GTPases (Guanosine tri-phosphatases) are also called GTP binding proteins or G-proteins which bind and hydrolyze GTP into GDP and inorganic c-phosphate.GTPases constitute an important class of regulatory proteins that regulate many biological processes such as signal transduction, cell multiplication, protein synthesis, ribosome biogenesis, etc., These regulations happen through conformational changes of GTPases between GTP-bound (On/active) and GDP-bound (Off/inactive) states (Bourne et al., 1991;Britton, 2009;Leipe et al., 2002;Vetter & Wittinghofer, 2001).
GTPases bind GTP within the GTP binding domain (Gdomain) through the conserved G-motifs (G1 [GxxxxGK(S/T)], G2 [T], G3 [DxxG], G4 [(N/T)(K/Q)x(D/E)] and a weakly conserved G5 [SAK] motif) (Bourne et al., 1990;Wittinghofer & Vetter, 2011).The G1 motif facilitates interaction with the phosphate moieties and hence the loop of the G1 motif is referred to as P-loop (phosphate-binding loop).The G2 and G3 motifs tend to interact with the c-phosphate and the regions around the G2 and G3 motifs are termed switch-I (Sw-I) and switch-II (Sw-II) respectively.Since these regions tend to undergo significant conformational changes and switch the neighboring molecules to perform their functions.
Therefore, GTPases are also called molecular switches.The extent of conformational changes varies from one GTPase to another and this may also occur in the regions which are other than loops of G2 and G3 (Upendra & Krishnaveni, 2022;Vetter & Wittinghofer, 2001).The G4 and G5 motifs interact and stabilize the guanine base of the bound nucleotide.
GTPases involved in ribosome biogenesis are referred to as RA-GTPases (ribosome assembly GTPases) (Britton, 2009).Ribosome is a complex macromolecular machine that synthesizes protein, it consists of two subunits (30S and 50S for 70S prokaryotic ribosome; and 40S and 60S for 80S eukaryotic ribosome) which are made up of rRNAs (ribosomal RNAs) and rproteins (ribosomal proteins).Assembling these constituents into a functional ribosome is referred to as ribosome biogenesis/ribosome maturation/ribosome assembly.The process of ribosome maturation requires a larger number of enzymatic factors.GTPases from TRAFAC (translation factor) class are a few among them and they are (i) RsgA, Era, and YqeH -which are involved in the maturation of the 30S ribosomal subunit; (ii) RbgA, ObgE, YsxC, and YphCwhich are involved in the maturation of the 50S ribosomal subunit; (iii) EngA, YihA, and Obg -which are involved in the maturation of both the 50S and 30S ribosomal subunits (Britton, 2009;Goto et al., 2013;Karbstein, 2007).
EngA (Essential Neisserial GTPase A) is a unique GTPase in which it contains two contiguous G-domains [GD1, GD2] at the N-terminal.This is followed by the C-terminal KH (K-homology) domain (Figure 1A) (Foucher et al., 2012;Hwang & Inouye, 2001;Muench et al., 2006;Robinson et al., 2002).EngA is shown to be involved in the later stages of ribosome maturation, depletion or mutation of EngA results in the incomplete formation of ribosome.Therefore, EngA has to associate with the ribosome in the intermediate stage of the maturation process (ribosomal intermediate subunits) for complete maturation (Bharat et al., 2006;Hwang & Inouye, 2006;Tomar et al., 2009).EngA has been shown to associate with both the ribosomal subunits (30S and 50S).For some species, this association of EngA to a particular subunit depends on the combination of the bound nucleotides.In [GTP, GTP] bound combination (both G-domains bound to GTP), EngA favors associating with the 50S subunit or its intermediate.In [GDP, GTP] bound combination, EngA is shown to associate with the 30S subunit (Lamb et al., 2007;Majumdar et al., 2017;Schaefer et al., 2006;Tomar et al., 2009).
This nucleotide dependent association with a particular subunit is attributed to the conformational changes of EngA in different nucleotide bound states.This attribution is further supported by the comparison of the crystal structures of Thermotoga maritima EngA (TmDer) in [GTP, GDP] bound state and Bacillus subtilis EngA (YphC) in [GDP, GDP] bound state, in which a significant conformational change is observed with a larger GD1 domain motion of 39 Å between the centroids of the TmDer-GD1 and the YphC-GD1 (Figure 1B) (Muench et al., 2006).Further, the comparison of G-domains of these YphC and TmDer shows a significant Sw-II helical shift for TmDer-GD1 compared to others (Figure 1C).Based on this, Muench et al. hypothesized that the exchange of nucleotide in the Gdomain setup the local fluctuations which eventually triggers to undergo larger domain motion (Muench et al., 2006).
In contrast, the comparison of the later determined crystal structures of YphC bound to different nucleotides shows neither the larger domain motion (Figure 1B) nor the local Sw-II helical shift (upon the comparison of mere G-domains) (Figure 1C).However, the SAXS profiles of YphC in various nucleotide bound combinations exhibit different patterns, suggesting the attainment of various conformations (da Silveira Tom� e et al., 2018;Foucher et al., 2012;Zhang et al., 2014).In addition, the comparison of the other structures of EngA from different species bound to various nucleotides (Escherichia coli EngA[GMPPNP, GMPPNP], Coxiella burnetii EngA[GDP, GDP], and Neisseria gonorrhoeae EngA[GDP, GDP]) exhibit totally different conformation from YphC and TmDer, in which GD1 domains are moved to various different positions, whereas GD2 and KH domains are superposed over each other (Figure 1B).
These variations in the conformational arrangement of EngA from species to species in various nucleotide bound combinations motivated us to explore conformational variations of EngA in different nucleotide bound combinations through molecular dynamics (MD) simulations.Our previous MD simulation study on TmDer in different nucleotide bound combinations show no larger GD1 domain motion but exhibited that the GD2 nucleotide bound state influences the GD1-KH interface interactions (Upendra & Krishnaveni, 2020).In continuation with this, for the present study, we have carried out 1000 ns all-atom MD simulations for YphC in different nucleotide bound combinations along with the presence or absence of Mg 2þ at the active-sites such as [GDP, GDP], [GDP, GTP], [GDP, GTP-Mg 2þ ], [GTP, GDP], [GTP-Mg 2þ , GDP], [GTP, GTP], and [GTP-Mg 2þ , GTP-Mg 2þ ].Basic conformational analyses through root-mean-square deviation (RMSD), radius of gyration (Rg), root-mean-square fluctuation (RMSF), and distances calculations exhibited the active-site dependent functionally important motions.The secondary structural analyses show an allosteric connection between the Sw-II regions of GD1 and GD2, which is further confirmed through the network analyses.Also, visualization of the trajectory revealed that the presence of Mg 2þ stabilizes active-site water molecules, the same is shown in distance calculations.

System preparation and MD simulations
The structure of YphC (PDB: 2HJG) was retrieved from PDB (Muench et al., 2006).In this structure of YphC (PDB:2HJG) a few amino acids from different regions and a few atoms from different residues are missing.The details of the missing residues and atoms are given in Tables S1 and S2.These missing amino acid residues and atoms were modeled using the loop modeling method of the Modeller-9.20 through the UCSF Chimera interface ( � Sali & Blundell, 1993;Webb & Sali, 2014;Yang et al., 2012).
Here the structure of YphC (PDB: 2HJG) is used as a template structure and residues are fixed based on the geometrical restraints.The complexes of YphC-guanine nucleotide combinations were prepared using Chimera by matching the existing nucleotide in the modeled structure with the required guanine nucleotide (Pettersen et al., 2004).The validity of the modeled complexes was assured based on the hydrogen bond interactions of G-motifs with the different parts of bound nucleotide as per the GTPase conventional interaction with the nucleotide (which is described in the introduction).All modeled complexes were simulated using the GROMACS-18.6 package (Abraham et al., 2015;Berendsen et al., 1995;Hess et al., 2008;Lindahl et al., 2001;Pronk et al., 2013;Van Der Spoel et al., 2005) for 1000 ns using all-atom MD simulations.Initially, the topologies for protein were generated using the AMBER99SB-ILDN force field (Lindorff-Larsen et al., 2010) and the topologies for guanine nucleotides were produced using the ANTECHAMBER module of the AMBERTOOLS20 package (Salomon-Ferrer et al., 2013).YphC-nucleotide complexes were placed in a triclinic box with a minimum distance of 1.2 nm from the box edges.The periodic boundary conditions (PBC) were used to avoid edge effects.Subsequently, the solvent was filled using the TIP3P water model (Jorgensen et al., 1983).To make the system to be neutral, some of the solvent molecules were replaced by Na þ and Cl -.Further, energy minimization was performed for 50000 steps using the steepest descent algorithm with an energy tolerance of 1000 kJ mol À 1 nm À 1 .In order to attain the temperature of 300 K and pressure of 1 bar, NVT and NPT equilibration simulations were performed for 2 ns using a V-rescale thermostat (Bussi et al., 2007) and Parrinello-Rahman barostat (Parrinello & Rahman, 1981).Finally, 1000 ns production simulations were conducted using an NPT ensemble.The cut-off of 1.2 nm was used for short-range interactions.The particle-mesh Ewald (PME) scheme was used for handling long-range interactions (Darden et al., 1993;Essmann et al., 1995).The LINCS algorithm was used to constrain the bonds (Hess et al., 1997).

Analyses
Preliminary conformational analyses such as root-meansquare deviation (RMSD), radius of gyration (Rg), root-meansquare fluctuation (RMSF) and distance calculations were carried out using built-in GROMACS tools, visual molecular dynamics (VMD) and PLUMED-2.5.2 (Bonomi et al., 2009).Secondary structure analyses were performed using the DSSP (Dictionary of protein secondary structure) (Kabsch & Sander, 1983) module integrated within GROMACS.Solvent accessible surface area (SASA) was calculated with a spherical probe radius of 0.14 nm using the SASA module of the GROMACS (Eisenhaber et al., 1995).Visualizations and image rendering were made with VMD (Humphrey et al., 1996) and Chimera package.Gnuplot was used for drawing graphs.

Extraction of non-bonded interaction energy
The non-bonded interaction energies within the cut-off of 1.2 nm were extracted between the required regions of the protein using the following procedure.First, the index groups for the required region have to be appended in the index file, subsequently, those appended index groups were used as energy index groups in the molecular dynamics parameters file (mdp file).Finally, the simulations have to be re-run with the new mdp file along with the original trajectory file.This produces the energy file which contains the information of the non-bonded interaction energy between the required regions.This can be extracted from the GROMACS energy module.

Dynamical network construction and analyses
Dynamical protein contact networks (PCNs) have been constructed for the simulated trajectories using the network view plugin of the VMD (Eargle & Luthey-Schulten, 2012;Sethi et al., 2009).Here each amino acid is treated as the node with Ca-atom as the representative atom and edges are drawn between the nodes when any of the two heavy atoms (non-hydrogen) from two different residues are within the 4.5 Å for 75% of the total simulation time.We have also assigned two nodes for each bound nucleotide, one for atoms of the guanine base and the other for the atoms of the phosphate moieties and ribose sugar.Covalently bonded neighbors were not considered for edges.Edges were weighted based on the correlations between the residues, where W ij and C ij are edge weight and correlation value between the residues i and j: The correlation calculation were performed with the Carma software (Glykos, 2006).These weight factors are used to calculate the path length between the distant nodes.Path length ðD ij Þ between the nodes i and j is the sum of the edge weights between the consecutive nodes ðk, lÞ along the path ði, jÞ is given by, D ij ¼ P k, l W kl : The shortest path or optimal path between the functionally important distant residues was determined based on the Floyd-Warshall algorithm (Floyd, 1962).The densely interconnected nodes of the network are identified and grouped into communities based on the Girvan-Newman algorithm (Girvan & Newman, 2002).

Results and discussion
In order to explore the conformational variations of YphC, all-atom MD simulations were performed upon different nucleotide bound combinations along with the presence or absence of Mg 2þ such as [GDP, GDP], [GDP, GTP], [GDP, GTP-Mg 2þ ], [GTP, GDP], [GTP-Mg 2þ , GDP], [GTP, GTP], and [GTP-Mg 2þ , GTP-Mg 2þ ].The choice of Mg 2þ to position in the active-site for the present study (instead of other ions like K þ ) is based on the survey of EngA structures available in the RCSB-PDB (research collaboratory for structural bioinformatics -protein data bank), where we found a structure of a G-domain of EngA (PDB: 5X4B) (this data is yet to be published) in which Mg 2þ is present in the active-site and three K þ ions are at the peripheral region.The details of the simulated systems are given in Table 1.

Structural stability and flexibility
Overall conformational variations were measured through RMSD and Rg using Ca atoms of YphC (Figure S1) over the simulation period for all nucleotide bound combinations.All simulated systems exhibit stable conformations with small variations.To explore the residue-wise conformational variations among different nucleotide bound combinations, RMSF using Ca atoms of YphC were obtained (Figure 2).RMSF plots exhibit significant fluctuations at all loop regions of YphC, the extent of fluctuations depends on the length and position of the loop.The switch regions -GD1-Sw-I (residues: 30-40), GD1-Sw-II (residues: 60-82), GD2-Sw-I (residues: 202-212), GD2-Sw-II (residues: 232-256), and the GD1 interface (GD1-IF) region (residues: 125-144) exhibit larger mobility as well as variation in mobilities among different nucleotide bound combinations.This prompted us to explore the dynamics of these highly mobile regions.The Sw-I regions of both GD1 and GD2 are not located near the bound nucleotide but are in the peripheral regions, where it can move freely.Visualization of trajectory (data not shown) exhibited random fluctuations, this suggests that the Sw-I regions of both GD1 and GD2 may stabilize through interaction with the ribosomal constituents (rRNA or rprotein) or with the other ribosome maturation factors (RbgA, ObgE, etc.,).

GD2-Sw-II mobility is influenced by the presence of Mg 21 at the active-site
The Sw-II regions of both GD1 and GD2 are located nearer to their bound nucleotide and exhibit larger fluctuations.In addition, variations of Sw-II fluctuations among different nucleotide bound combinations are relatively larger for the GD2 domain and small for GD1.This prompted us to explore the movements of Sw-II regions of both GD1 and GD2 with respect to their bound nucleotide during the simulation period.Therefore, the distances between the terminal atoms of side-chains of the charged residues (Asp, Glu, Lys, and Arg) of the Sw-II loop (Figure 3) with the terminal phosphates of the bound nucleotide (PG/PB) were obtained (Figures S2 and S3).Here charged residues are considered because they are more sensitive/influential to the variation in the charge state of the active-site environment as compared to other neutral or polar amino acid residues.
The charged residues of GD1-Sw-II (Asp62, Asp65, and Glu66) have been shown to stay away from the bound nucleotide and their movements are not much influenced by the bound nucleotide as well as the Mg 2þ (Figure S2).
However, the movement of charged residues of GD2-Sw-II (Arg234, Lys235, Lys236, Lys238, Glu241, and Glu244) have been shown to be influenced by the presence or absence of Mg 2þ along with their bound nucleotide.Here, one or the other or a few residues approach the bound nucleotide (Figure S3).In the [GDP, GDP] bound system (depicted with brown color), the positively charged residue Lys235 approach the negatively charged phosphate moieties of the bound nucleotide and stay near to it over the entire simulation period (Figure S3B, brown).In the [GDP, GTP] bound system (blue), the residues Lys236 and Lys238 stay near the bound nucleotide around 0-200 ns and 200-400 ns respectively (Figure S3C and D, blue).In the [GDP, GTP-Mg 2þ ] bound system (red), due to the presence of positively charged Mg 2þ , negatively charged Glu244 tend to approach the active-site for a significant course of the total simulation time (Figure S3F, red).In the [GTP, GDP] bound system (green), the residues Lys236 and Lys238 approach the bound nucleotide around 0-400 ns and 600-900 ns respectively (Figure S3C and D, green).In the [GTP-Mg 2þ , GDP] bound system (orange), the residue Lys235 stays near to the bound nucleotide for almost the total simulation period (Figure S3B, orange).In [GTP, GTP] bound system (cyan), the residues Lys235 and Lys238, stay near to the bound nucleotide for almost the entire simulation period (Figure S3B and D).In the [GTP-Mg 2þ , GTP-Mg 2þ ] bound system (magenta), no charged residues of the Sw-II loop have approached the bound nucleotide but distances seem to be stable (Figure S3, magenta).These distance measurements (Figures S2 and S3) provide the reasons for observed fluctuations in the RMSF plot (Figure 2).Larger fluctuations are due to the larger jumps (fluctuations) in the distances between the charged residues of Sw-II and the phosphate moieties of the bound nucleotide.

The GD2 nucleotide bound state influences the GD1-Sw-II secondary conformations
Exploration of Sw-II movements with respect to the bound nucleotide through distances shows well equilibration (stable distances) for some nucleotide bound combinations and random shifts (shoots up/down in distances at various time points) for some other nucleotide bound systems (Figures S2  and S3).This led us to obtain the secondary structural conformations of the Sw-II regions of both GD1 and GD2 over the simulation period (Figure 4).The secondary structure analyses of GD1-Sw-II exhibited the formation of loosely coiled p-helix (5-helix) for systems having GTP at GD2, in particular, the occupation of the 5helix form is longer when both G-domains possess GTP (Figure 4B, C, F, & G), whereas, in other systems, regular a-helix form is occupied almost the entire simulation period (Figure 4A, D, & E).The constituent residues (residues: 75-80) of the formed p-helix are located at the edge of the GD1-Sw-II, which is at the interface of the linker, this is at the interface of GD2-Sw-II (Figure 1A).This may indicate that the GD2 nucleotide bound state influences the GD1-Sw-II.This allosteric connection between the GD2 nucleotide bound state and the GD1-Sw-II was also observed in our previous simulation work on TmDer (Upendra & Krishnaveni, 2020).In that, the GD2 nucleotide bound state influences the GD1-KH interface interactions, and the GD1 interface is nothing but the Sw-II helical region of GD1.
The secondary structure analyses of GD2-Sw-II have exhibited a pattern of tightly packed 3 10 -helix (residues: 240-243) and p-helix (residues: 245-250) for the [GTP-Mg 2þ , GTP-Mg 2þ ] bound system (Figure 4N).The region where p-helix formed is at the interface to the linker, and the linker is embedded between the Sw-II regions of both GD1 and GD2 (Figure 1A).This suggests that the Sw-II of both domains tends to approach when both G-domains have GTP-Mg 2þ .The region where 3 10 -helix formed is in the peripheral region, indicating the stability of that region.This pattern of 3 10 -helix is rarely formed in the other nucleotide bound combinations.But, In the [GDP, GTP-Mg 2þ ] bound system, a-helix (at residues: 242-246) is formed for a significant period of the total simulation time after 500 ns (Figure 4J).This stabilizes the GD2-Sw-II peripheral region (residues: 240-246).Thus, when GD2 possesses GTP-Mg 2þ , the GD2-Sw-II peripheral region stabilizes and probably has more chances to interact with rRNA or rprotein.

Communication between the GD2 active-site and the GD1-Sw-II loop is stronger for systems that possess GTP at both G-domains
Secondary structure analyses of Sw-II regions of both GD1 and GD2 have indicated some sort of allosteric connection between the GD2 active-site and the GD1-Sw-II.This prompted us to construct the dynamic PCNs and obtained their associated network metrics (Figures S4, S5, 5, and S6) (Tables 2 and 3) in order to extract the connection between the GD2 active-site and the GD1-Sw-II.
The communication pathways between the source the Ser189 of the GD2 active-site and the targets the Glu82 and Gly64 of the GD1-Sw-II were extracted by means of the shortest and their associated suboptimal paths within the length offset of 5 in different nucleotide bound combinations (Figures 5 and S6) (Tables 2 and 3).The Ser189 is chosen as a source point, since it is from the P-loop of the GD2 and interacts with the gamma-phosphate.Thus, the Ser189 immediately experiences the environmental changes in the active-site and communicates the effect to the functionally important distant residue through interconnected residues.The residues Glu82 and Gly64 are chosen as target points since Glu82 lies at the edge of the GD1-Sw-II and the Gly64 is in the loop of the GD1-Sw-II and shows higher mobility in RMSF (Figure 2).The [GTP-Mg 2þ , GTP-Mg 2þ ] bound system shows the smaller shortest path length compared to other nucleotide bound systems that possess Mg 2þ (Table 2).This indicates that communication between the chosen source (Ser189) and the target (Glu82) is stronger for the respective system as anticipated from the secondary structure analyses (Figure 4).The stronger communication strength of the [GTP-Mg 2þ , GTP-Mg 2þ ] bound system can also be evidenced by the larger number of sub-optimal paths and these sub-optimal paths are almost degenerate.This degeneracy of sub-optimal paths can be observed in Figure 5 with the larger size of the communicating edges.However, the [GTP, GTP] bound system exhibits a smaller shortest path length compared to the [GTP-Mg 2þ , GTP-Mg 2þ ] bound system but the number of sub-optimal paths is smaller.Thus, the communication between the Ser189 and the target Glu82 is stronger and highly correlated for the [GTP-Mg 2þ , GTP-Mg 2þ ] bound system.Similarly, stronger communication is observed between the source (Ser189) and target (Gly64) for the [GTP-Mg 2þ , GTP-Mg 2þ ] bound system (Table 3 and Figure S6) compared to other nucleotide bound systems.
The enhancement of the communication between the GD2 active-site and Gly64 of the GD1-Sw-II loop suggests that the region around Gly64 has more chances to interact with the ribosome as it is at the periphery.
Further, the extraction of the total non-bonded interaction energies (combination of Coulomb and Lennard-Jones interaction energies) between the Sw-II regions of the GD1 and Table 2. Communication between the Ser189 of the GD2 P-loop and the Glu82 of the GD1-SwII in the Dynamic PCN for different nucleotide bound systems.

System
Residues involved in the shortest Path Length of the shortest Path No. of suboptimal paths within the length offset of 5 The notations 2GDP438, 2GTP438, and 2GDP439 are used for the bound guanine nucleotides (GDP/GTP) of the GD2.Since we also have assigned nodes for each bound nucleotide along with the amino acid residues in the constructed network, as a result, these nucleotides appeared in some shortest paths as one of the communicating nodes.Here the number '2' prefix to 'GTP/GDP' represent the second G-domain.Whereas the numbers 438 and 439 suffix to 'GTP/ GDP' represents the assigned residue numbers, as YphC contains 436 amino acid residues, therefore the extra added ions or nucleotides will be assigned to new residue numbers.GD2 exhibited a relatively more attractive tendency for the [GTP-Mg 2þ , GTP-Mg 2þ ] bound system compared to others (Figure S7) as anticipated from the secondary structure analyses (Figure 4).
Thus, in YphC, the GD1-KH interface interaction tends to be more for systems that possess GTP-Mg 2þ at GD1.In contrast, our previous MD simulation study on TmDer (Upendra & Krishnaveni, 2020) shows that the GD1-KH interface interaction energies are larger when GD2 possesses GTP.That could be due to different interface regions in TmDer owing to the different conformational arrangement compared to YphC.Further, the MD simulations of TmDer were not performed in the presence of Mg 2þ at the active-sites.This suggests that further simulation studies for TmDer in the presence of Mg 2þ at the active-site need to be performed in order to have a clear picture of GD1-KH interface interactions.
Further, the solvent accessible surface area (SASA) for the entire region of YphC was obtained to check the effect of this transition at the GD1-KH interface (Figure S10).The [GTP-Mg 2þ , GTP-Mg 2þ ] bound system has slightly more SASA compared to other nucleotide bound systems, this probably enhances the interaction with the ribosome.

The presence of Mg 21 stabilizes the activesite water
The presence or absence of Mg 2þ at the active-site influences the Sw-II fluctuations as well as GD1-KH interface interactions (Figures 3 and 7).Further, to explore the role of activesite Mg 2þ on the GTPase activity, the stability of the solvated water molecules at the active-site (active-site water molecule/s) was monitored through the distances between the active-site water molecules and the terminal phosphate (PG/ PB) of the bound nucleotide over the simulation period (Figures 8 and 9).
In the presence of Mg 2þ , the active-site water molecule/s stay in the active-site over the complete simulation period through the coordination of Mg 2þ with the oxygen atom of water molecules (Figures 8,9,and S11).The presence of Mg 2þ also increases the distance between the a-and cphosphates of the bound GTP, this probably favors the GTP hydrolysis process (Figures 10 and 11).In contrast, in the absence of Mg 2þ , the active-site water molecules are displaced immediately from the active-site and mixed with the bulk solvent, this probably reduces the GTPase activity.
Curiously, three water molecules reside at the GD1 activesite (Figure 8E and G) and only one water molecule resides  at the GD2 active-site (Figure 9C and G) over the complete simulation period.This suggests that the probability of GTP hydrolysis is more in GD1 than in GD2.Interestingly, this is in accordance with the experimentally reported data (Foucher et al., 2012), in which they have shown that the GTPase activity of GD1 is 150 fold higher than the GTPase activity of GD2.This difference in the stabilization of three and one water molecules for the respective GD1 and GD2 may be due to the minute sequential and structural difference.However, the exact role of these water molecules in the GTP hydrolysis mechanism has to be studied in future studies using QM/MM simulations.
The water molecules stay at the active-site due to the intrinsic nature of Mg 2þ to coordinate with six oxygen atoms.However, only three oxygen atoms from water molecules in GD1 active-site and only one oxygen atom from a water molecule in GD2 active-site are coordinated with the Mg 2þ .The remaining coordinated atoms are given in Figure S12.For GD1-Mg 2þ , the coordinated atoms are oxygen atoms from phosphate moieties of bound GTP (O2B, O2G, O3G) along with the oxygen atoms from three water molecules (Figure S12 B and C) .For GD2-Mg 2þ , the coordinated atoms are oxygen from Ser189(OG), carbonyl oxygen from Met233(O), and oxygen atoms from phosphate moieties of bound GTP (O2B, O2G, O3G) along with the oxygen atom of a water molecule (Figure S12 A and D).

Conclusion
1000 ns MD simulations were carried out for YphC in different nucleotide bound combinations.The stability of the simulations was monitored through RMSD and Rg of Ca atoms of the entire protein.Exploration of the residue-wise mobility through RMSF exhibited larger fluctuation at the switch regions as well as the GD1-IF region.This prompted us to explore the dynamics of these highly mobile regions.Visualization of the Sw-I movements of both the G-domains exhibited random fluctuations and independent of the bound nucleotide owing to its position at the peripheral region.This Sw-I dynamics probably stabilize through the interaction with the ribosome constituents or ribosome maturation factors.Exploration of the dynamics of Sw-II regions of both the G-domains through the distances between the charged residues of the Sw-II loops with the terminal phosphates (PB/PG) of the bound nucleotide exhibited that the GD2-Sw-II dynamics are influenced by the presence or absence of Mg 2þ along with the bound nucleotide.This tendency is not observed in the GD1-Sw-II, however, the GD1-Sw-II secondary conformations were shown to be influenced by the GD2 nucleotide bound state.This allosteric connection between the GD2 active-site and the GD1-Sw-II was confirmed through the dynamical PCN.
Exploration of the dynamics of the GD1-IF through RMSD exhibited a transition tendency for the [GTP-Mg 2þ , GTP-Mg 2þ ] bound system, this transition is directed towards the KH domain as observed from the total non-bonded interaction energy.
Further, the role of Mg 2þ on the GTPase activity were investigated.The presence of Mg 2þ at the active-sites stabilizes the active-site water molecules as well as increases the distances between the a-and c-phosphates of the bound GTP, this probably favors the GTPase activity.Surprisingly, it is found that three water molecules in the GD1 active-site and only one water molecule in the GD2 active-site are stabilized, this suggests that the probability of GTP hydrolysis is more in GD1 compared to GD2.Interestingly, this is also in accordance with experimentally reported data.However, the exact role of these water molecules in the GTP hydrolysis mechanism has to be investigated in future studies through QM/MM simulations.

Figure 2 .
Figure 2. Root-mean-square fluctuation (RMSF) of Ca atoms of YphC with respect to the average structure for different nucleotide bound systems.

Figure 3 .
Figure 3. Depiction of the positions of the charged amino acid residues of the Sw-II loops of the (A) GD1 domain and (B) GD2 domain.

Figure 5 .
Figure5.The shortest paths (red) and their associated suboptimal paths (blue) between the Ser189 of the GD2 active-site and the Glu82 of the GD1-Sw-II in different nucleotide bound combinations.The size of the edge is proportional to the number of paths passing through that edge.

Figure 6 .
Figure 6.Root-mean-square deviation (RMSD) of Ca atoms of GD1 interface during the simulation period with respect to the starting structure for the different nucleotide bound combinations.

Figure 7 .
Figure 7.The total non-bonded interaction energy between the GD1 and KH interface regions for different nucleotide bound combinations during the simulation period.

Figure 8 .
Figure 8. Distances between the terminal phosphate and the active-site water molecules of the GD1 domain, emphasizing the stability of the active-site water molecules in the presence of Mg 2þ in the active-site.

Figure 9 .
Figure 9. Distances between the terminal phosphate and the active-site water molecules of the GD2 domain, emphasizing the stability of the active-site water molecules in the presence of Mg 2þ in the active-site.

Figure 10 .
Figure 10.(A and B) Distances between the PA and PG of the bound GTP of the GD1 domain for different nucleotide bound systems and (C) their corresponding histogram plots.(D) Structures of GTP extracted from the simulation frames of different nucleotide bound systems, depicting more PA-PG distance for the GD1 domains having Mg 2þ .

Figure 11 .
Figure 11.(A and B) Distances between the PA and PG of the bound GTP of the GD2 domain for different nucleotide bound systems and (C) their corresponding histogram plots.(D) Structures of GTP extracted from the simulation frames of different nucleotide bound systems, depicting more PA-PG distance for the GD2 domains having Mg 2þ .

Table 1 .
Details of the simulated systems -YphC bound to different nucleotides.

Table 3 .
Communication between the Ser189 of the GD2 P-loop and the GLy64 of the GD1-SwII in the Dynamic PCN for different nucleotide bound systems.