Mazhar Iqbal and
Syed Tarique Moin*
Third World Center for Science and Technology, H. E. J. Research Institute of Chemistry, International Center for Chemical and Biological Sciences, University of Karachi, Karachi-75270, Pakistan. E-mail: tarique.syed@iccs.edu; Fax: +92-21-348-19018; Tel: +92-21-99261774
First published on 13th June 2022
Riboswitch is a regulatory segment of messenger RNA (mRNA), which by binding to various cellular metabolites regulates the activity of mRNA via modulating transcription, translation, alternative splicing, and stability of the mRNA. yybP–ykoY riboswitch of Lactococcus lactis, which is present upstream of the yoaB gene, functions as a Mn2+-specific genetic ON-switch, and modulates expression of proteins which are significant for Mn2+ homeostasis. The P1.1 switch helix of the aptamer domain of the riboswitch contains an intrinsic transcription terminator structure, which gets stabilized with Mn2+ binding and causes disruption of terminator structure and allows the continuation of transcription. The current research work involved the evaluation of structural and dynamical properties of the yybP-ykoY riboswitch of L. lactis in its Mn2+-free, Mn2+-bound (wild-type), and Mn2+-bound mutant (A41U) states by applying molecular dynamics simulations. Based on the simulations, the effects of Mn2+ absence and A41U mutation were evaluated on the structure and dynamics of the riboswitches followed by the computation of the free energy of metal binding in the wild-type and the mutant riboswitches. The simulation results provided insights into the properties of the riboswitch with the focus on the dynamics of the P1.1 switch helix, and the manganese binding site designated as MB site, as well as the relative stability of the wild-type and the mutant riboswitches, which helped to understand the structural and dynamical role of the metal ion involved in the function of Mn2+-sensing riboswitch.
yybP-ykoY riboswitch that is mostly present in plants and human bacterial pathogens for instance Lactococcus lactis, Bacillus subtilis, Escherichia coli, Xanthomonas oryzae, etc., has high sensitivity for manganese (Mn2+) ion. It senses charge, geometry and Lewis acid hardness of Mn2+ by making inner coordination sphere contacts directly with five phosphoryl oxygen atoms and one nitrogen atom prototyped as OP1G9, OP1U44, OP2C45, OP2U39, OP2C40, and N7A41 in L. lactis.10 It also regulates distinct genes and is involved in Mn2+ homeostasis.10 There are two metal binding sites designated as MA and MB in the riboswitch of L. lactis. Magnesium (Mg2+) ion, owing to its more abundance than Mn2+ at physiological conditions, and having high similarity in ionic radius, charge and octahedral coordination scheme with Mn2+, occupies the MA site, whereas the MB site has high selectivity for Mn2+ due to its coordination with nitrogen atom (N7) of invariable adenosine (A41). The hairpin like structure of yybP–ykoY riboswitch contains four helices which are designated as P1, P2, P3 and P4. Two super helices formed by the co-axial stacking of the four helices (P1–P4) are docked at the conserved L1 and L3 loops, and this docked interface is stabilized by two highly coordinated metal ions (Mg2+ and Mn2+) illustrated in Fig. 1. In vitro as well as in vivo studies revealed that yybP–ykoY riboswitch at upstream of yoaB gene functions as Mn2+-specific genetic ON-switch, and modulates expression of proteins significant for Mn2+ homeostasis.11 The P1.1 helix of the aptamer domain contains intrinsic transcription terminator. In the absence of Mn2+ ion, RNA polymerase constructs predominantly terminated transcript, however with Mn2+ binding, the aptamer domain gets stabilized and causes disruption of transcription-terminator structure and allows continuation of transcription.11,12
One experimental study was reported on the crystallization of L. lactis yybP–ykoY riboswitch in its ligand-bound form, and its comparison with that of the E. coli riboswitch in Mn2+-free form.13 Unlike Mn2+-bound crystal structure of yybP–ykoY riboswitch of L. lactis, the crystal structure of the Mn2+-free riboswitch has not been reported yet. Furthermore, to the best of our knowledge, no computational study has also been carried out which could describe complete dynamics of Mn2+-free riboswitch of L. lactis, thus leaving behind a question on how yybP–ykoY riboswitch of L. lactis changes its structural and dynamical properties in the absence of Mn2+ from the MB site. Although, molecular dynamics simulations were carried out for yybP–ykoY riboswitches of L. lactis and X. oryzae in which divalent metal ions present in the metal-binding sites (MA and MB) were replaced with monovalent K+ ions and their effect on the structural flexibility was determined,14 however, a detail study has not been reported yet, comparing the dynamics of Mn2+-bound yybP–ykoY riboswitch regarded as the wild-type, which contains Mg2+ and Mn2+ at the MA and the MB site, respectively, A41U Mn2+-bound – a mutant in which an invariable Mn2+-sensing adenosine is replaced with uracil, and Mn2+-free with Mg2+ at MA site only of yybP–ykoY riboswitches of L. lactis. The aim of this study was to gain deep understanding of the structural and dynamical properties of the Mn2+-free, Mn2+-bound and A41U Mn2+-bound yybP-ykoY riboswitches of L. lactis at the atomistic details.
In one experimental study, single molecule kinetic analysis of RNA transient structure (SiM-KARTS) assay showed that the anti-terminating P1.1 switch helix in yybP–ykoY riboswitch gets stabilized as Mn2+ binds to MB site via signal transduction through metal ion sensing core, thus restraining the formation of transcription-terminator hairpin. This stabilizing effect is lost upon mutation of the invariable Mn2+-sensing adenosine (A41) with uridine.14 Through molecular dynamics simulations, we elaborated the metal binding and the mutation effect in L. lactis yybP–ykoY riboswitch as well as found how the stability of P1.1 switch helix is affected in the absence of Mn2+. To achieve the objective, we ran 1 μs all-atom MD simulations of the Mn2+-free, Mn2+-bound and A41U Mn2+-bound riboswitches and came up with the detailed description of the structural and dynamical properties with the focus on the dynamic flexibility of the P1.1 switch helix.
Initial configuration of the simulation systems was obtained by using xLEaP module of AMBER20 in which ff99bsc0χOL3 AMBER force field16–18 was employed for modeling the ribonucleotides, and non-bonded model approach was adopted to treat the metal ions in which the metal center was described as simple Lennard-Jones sphere with an integer charge considering only electrostatic and van der Waals (vdW) terms which were used to describe the interaction between the metal ions and its surroundings. Here, the metal–ligand interactions were described through electrostatic and van der Waals parameters, with σ and ε values of 1.40 Å and 0.0168 kcal mol−1, and 1.36 Å and 0.0102 kcal mol−1, for Mn2+ and Mg2+, respectively.19 The implementation of the force field parameters to these three states of the riboswitch was followed by protonation, neutralization with Na+ as counter ions and solvation with 12 Å buffer of TIP3P water model20 in a truncated octahedral box with periodic boundary conditions.
The relative stability of the wild-type (Mn2+-bound) and the mutant (A41U Mn2+-bound) riboswitches were estimated in terms of relative binding free energy (ΔΔGbind) that was calculated as the free energy (ΔGAB) difference of the transformation of wild-type form into mutant (A41U) form with Mn2+ ion at MB site in process 1 and without Mn2+ ion at MB site in process 2 employing thermodynamic integration (TI) approach (Scheme 1). TI based free energy estimation method computes the free energy difference between the two states A and B, in this case wild-type and mutant riboswitches, by coupling them via a parameter λ that serves as an additional, non-spatial coordinate, which allows the free energy difference between the states to be computed as:
![]() | (1) |
![]() | ||
Scheme 1 : Thermodynamic cycle employed for the computation of relative binding free energy difference between wild type and mutant yybP–ykoY riboswitches. |
The Hamiltonian, H that indicates the extent of change came about between the wild-type (state A) and mutant (state B), was calculated as a function of coupling parameter λ. TI simulations were performed for each process employing one-step protocol that is based on disappearing one unique residue and appearing other unique residue simultaneously. Atoms of both unique residues were in softcore region and were individually specified during calculation, and both van der Waals (vdW) and charge interactions of disappearing or appearing unique atoms with neighboring atoms were illustrated by softcore potentials. We used dual topology approach in which both the wild-type and the mutant molecules were considered, simultaneously and atoms of common residues in both states were having same coordinates. To remove redundant bonding term calculations, both the wild-type and the mutant molecules were merged in such a way that the new topology contained only the common atoms and the softcore atoms which were atoms of both unique residues. TI simulations were run for each λ value and gradient of Hamiltonian H with respect to the coupling parameter, λ, was averaged over numerous equilibrated configurations generated during simulations that was then numerically integrated to get free energy difference (ΔΔGbind) of the two states. In the current study, the free energy simulation protocol employed 6 λ points with an array of 0.2 from 0.0 to 1.0, and 1 ns TI simulation was performed for each λ point at the isobaric-isothermal ensemble to obtain ΔΔGbind. The simulation systems were minimized at each λ point under periodic boundary conditions, employing steepest-descent algorithm prior to TI simulation. Afterward, the systems were heated up to ∼300 K over 20 ps under NVT conditions using weak-coupling algorithm with constraints applied to bonds involving hydrogen atoms using SHAKE algorithm, and positional restraint with force constant of 5 kcal mol−1 Å−2 was also applied to all atoms except water molecules and hydrogen atoms. Finally, each simulation system was simulated in the isobaric–isothermal (NPT) ensemble at each λ point for 1 ns with time-step of 2 fs employing CPU version of pmemd module implemented in AMBER 20;26 the stretching of bonds involving hydrogen atoms was constrained with SHAKE algorithm, temperature at 300 K and pressure at 1 atm were regulated by Langevin thermostat with collision frequency of 2 ps−1 and Berendsen barostat with pressure relaxation time of 2 ps, respectively. Softcore potential was employed, and α and β parameters of the softcore potential were set to 0.5 and 12 Å2, respectively. All essential ΔG values were obtained by numerically integrating ∂H/∂λ values collected during the free energy simulations.
The simulation trajectories were primarily analyzed in terms of root mean square deviation (RMSD) of the backbone atoms (P, O5′, C5′, C4′, C3′ & O3′) and heavy atoms belonging to nucleotides of the MB site employing crystallographic structure as a reference via CppTraj,27 an analysis tool implemented in AMBERTOOLS 20. Moreover, radius of gyration (Rg), root mean square fluctuation (RMSF), and dynamic cross correlation matrix (DCCM) were also evaluated. A correlation plot between Rg and RMSD was obtained with the help of pytraj, a python version of CppTraj. Essential dynamics were obtained via principal component analysis (PCA) using the gromacs analysis tools.28–32 The trajectories were visualized employing VMD program and the images were captured with the help of CHIMERA and PYMOL,33 along with VMD.34
Fig. 2c illustrated the RMSD profiles of the MB site nucleotides in Mn2+-free, wild-type and mutant riboswitches. In case of Mn2+-free riboswitch, the RMSD profile exhibited more fluctuation compared to the wild-type case, even though the averaged RMSD of the Mn2+-free riboswitch (0.59 ± 0.09 Å) was close to the wild-type case (0.56 ± 0.06 Å) but lower than that of the case of mutant riboswitch (0.78 ± 0.08 Å). The variation in RMSDs pointed toward a contrast dynamical behaviour of the two ions, however, both were reportedly having very similar structural properties.35 Moreover, the variation in the backbone RMSD was lower in both cases (Mn2+-free and wild-type riboswitches) compared to that of the mutant case thus suggesting that the overall dynamics of the MB site in Mn2+-free riboswitch was lower than that of the mutant riboswitch, however it was conversely observed for the MB site RMSD as the mutant riboswitch case demonstrated high fluctuations than that of the wild-type thus pointing toward the lower stability of the metal-binding site due to mutation. The contrasting behaviour was attributed to the mutation in the riboswitch that increased dynamic flexibility of the metal-binding site but caused the overall structure of the riboswitch less flexible. The RMSD data thus reflected that the dynamics of yybP–ykoY riboswitch was altered with the mutation and absence of manganese ion, which therefore led to further evaluate the effect of mutation and Mn2+ absence.
Time evolution radius of gyration (Rg) of the backbone atoms of the three riboswitches was computed to obtain information of the overall dynamics in the absence and presence of the metal ion, and due to the mutation shown in Fig. 3a. Rg gives estimation of a system's compactness as the higher Rg value indicates lower compactness.36 The average Rg value for the Mn2+-free riboswitch (24.82 ± 0.35 Å) was very close to that of the wild-type (24.64 ± 0.35 Å) but higher than the mutant riboswitch (23.94 ± 0.27 Å), which indicated that the volume of the Mn2+-free and the wild-type riboswitches was indifferent demonstrating that Mn2+ absence has little effect on the overall structure of the riboswitch, however, their volume was slightly higher compared to the mutant riboswitch. Fig. 3b displayed the Rg vs. RMSD correlation plots; in case of Mn2+-free riboswitch, distribution of the sampled conformers was located at RMSD from 4.05 to 8.05 Å corresponding to Rg value between 24.0 and 25.7 Å, whereas for the wild-type case, the corresponding distribution showed similar RMSD from 4.01 to 8.29 Å and Rg from 25.69 to 23.8 Å. In the case of the mutant riboswitch, RMSD appeared between 2.99 and 5.60 Å which corresponded to Rg value from 23.29 to 24.63 Å. The correlation plots exhibited that the dynamics of the Mn2+-free and wild-type riboswitches yielded to be similar, whereas the mutant case demonstrated distinct dynamics than the other two riboswitches.37
Effect of Mn2+ absence and A41U mutation on the dynamics of the riboswitch, in particular the P1.1 switch helix, was analyzed by root mean square fluctuation (RMSF), which facilitated to understand the flexibility of different nucleotide regions during simulations.38 Fig. 4a showed the RMSF profiles of the Mn2+-free, the wild-type and the mutant riboswitch; in the case of Mn2+-free riboswitch, residual fluctuations were higher than that of the wild-type riboswitch except for U66, G67, U68 and U86 nucleotides, and the MB site nucleotides (G9, U39, C40, A41, U44 and C45) exhibited high fluctuations depicted in Fig. 4d–f. The higher flexibility of different regions in Mn2+-free riboswitch compared to that in the wild-type was attributed to the absence of the manganese ion which caused the riboswitch more flexible. A similar fluctuation trend was observed for the P1.1 switch helix comprising of A5, G6, G7, G8, C96, C97, U98, U99, U100 and C101 nucleotides in the Mn2+-free riboswitch demonstrating high flexibility as deduced from the high intensity fluctuations in the RMSF profiles shown in Fig. 4b and c. This was interpreted as manganese and magnesium ions reportedly have distinct dynamical properties therefore, the binding of Mg2+ ion to the Mn2+-free riboswitch was assumed altered thereby destabilizing the MB site and the P1.1 switch helix. The stabilization of these regions can only be achieved only if the suitable metal ion like Mn2+ essentially binds to the metal-binding site as in the wild-type riboswitch, since according to an experimental study, the binding of Mn2+ ion disrupts transcription terminator hairpin structure and thus stabilizes anti-terminating P1.1 switch helix leading to transcription continuance.11,20 The RMSF profile of the mutant riboswitch showed mixed residual fluctuations in comparison to the wild-type case as shown in Fig. 4a, and the MB site residues in particular G9 and U41 nucleotides showed higher fluctuations, and U39 and C40 nucleotides showed lower fluctuations, whereas U44 and C45 nucleotides possessed similar fluctuations as observed in case of wild-type riboswitch (Fig. 4d–f). Nucleotides of the P1.1 switch helix, in the mutant riboswitch showed distinct fluctuation pattern compared to the wild-type, since A3 and A4 nucleotides produced low fluctuations compared to G6, G7, U100, C101 and C102 nucleotides having high fluctuations, whereas A5, C96 C97 U98 and U99 nucleotides possessed similar fluctuations like those of the wild-type riboswitch (Fig. 4b and c). This implied that the mutation had mixed effect on the overall structure flexibility of the riboswitch, including the MB site and the P1.1 helix.
It was therefore observed that the mutation perturbed the MB site despite an overall structural stabilization occurred in the mutant riboswitch. To evaluate the contrast in the structure and dynamics of the wild-type riboswitch in comparison to other two riboswitches, comparative dynamic cross-correlation matrix (DCCM) analysis was carried out for Mn2+-free, wild-type and mutant riboswitch to figure out the correlation or anti-correlation between nucleotides of the MB site and other regions, in particular the P1.1 switch helix (Fig. 5). DCCM plots depicted correlated movements between residues belonging to different regions of the riboswitches which were observed throughout 1 μs trajectories of each riboswitch. In DCCM plot, correlation range was from +1 to −1 representing with variation in the color intensity; red-colored contours correspond to positively correlated movements, white-colored regions represent zero correlation or uncorrelated movements and blue-colored contours signify negatively or anti-correlated movements. The characteristic features in the DCCM plots were depicted as boxed regions labelled as 1 to 16 showing effects of Mn2+ absence and A41U mutation on the structure and dynamics of the riboswitch. In the case of Mn2+-free riboswitch, region 8, and 3 and 5 corresponded to positively correlated and zero correlated motions, respectively (Fig. 5a), whereas the presence of Mn2+ ion in the wild-type riboswitch caused the riboswitch to attain mostly positively correlated motions (Fig. 5b). The A41U mutation was demonstrated to affect the atomic motions causing the mutant riboswitch to adopt distinct motions for instance regions 12, 13, 15 and 16 corresponded to zero correlation (Fig. 5c). Moreover, a contrast was observed between Mn2+-free and wild-type cases, as regions 1 and 7 in the Mn2+-free case depicted positively correlated motions compared to the corresponding region with zero correlated motions. Similarily, a constrast was observed for the case of mutant riboswitch when comparing with the wild-type case; regions 7, and 2, 9, 14 and 11 in the mutant case corresponded to positively and negatively correlated motions, respectively, whereas Mn2+-free case had mixed correlations consisting of both positive and negative correlations for the atomic motions corresponding to respective regions 1 and 7, and 2, 4, 6, 9 and 10. Nevertheless, the difference in the DCCM plots was too small to observe the obvious effect of Mn2+-absence and A41U mutation on the correlation dynamics among different regions, particularly, between MB site nucleotides (G9, U44, A41, C45, U39 and C40) and nucleotides of P1.1 switch helix (A3, A4, A5, G6, G7, C96, C97, U98, U99, U100, C101 and C102) of the riboswitches, which provoked to perform essential dynamics (ED) or principal component analysis (PCA) in order to evaluate distinct dynamics of the riboswitches.
PCA reduces the dimensionality of huge data set and helps to identify significant motions which are incumbent for the biological activity of a molecule.39 As the essential motions are extracted from sampled conformations, the method of applying PCA to a simulated trajectory is called essential dynamics (ED).40–42 With the removal of roto-translation motions of all atoms as of the nucleic acid, conformations are collected throughout the MD trajectory. Consequently, covariance matrix is set up that upon diagonalization of system's atomic coordinates targets the significant motions of conformation which are accountable for the mass variance in atomic positions as eigenvalues. Eigenvectors (essential modes) are determined by ranking eigenvalues of principal components (PCs) represented by scree plot shown in Fig. 6. Afterward, essential conformational space is determined by projecting the original configuration of the system on the PCs.39 Essential subspace comprises a few degree of freedom, in which positional fluctuations occur an-harmonically, therefore, essential degree of freedom points out the motions relevant to the function of biomolecules.42 Here, eigenvalues from covariance matrix were constructed from all-atom coordinates of the Mn2+-free, the wild-type and the mutant riboswitches which are shown in descending order. Eigenvalues of the Mn2+-free case showed steeper decrease than those of the wild-type and the mutant cases, and the wild-type case also showed steeper decrease than the mutant case. Trace of the covariance matrix of the Mn2+-free case happened to be 332.71 nm2 whereas for the wild-type and mutant cases, the values yielded to be 203.66 and 205.39 nm2, respectively. Since the high values of trace of the covariance matrix generally pointed toward more flexible system,39 therefore, high value of trace of covariance in case of Mn2+-free riboswitch demonstrated high dynamic flexibility due to Mn2+ absence in comparison to other two riboswitches, and the mutant case with slightly high trace of covariance matrix was assessed to have slightly more flexible structure compared to the wild-type riboswitch. Moreover, large eigenvalues for the case of Mn2+-free riboswitch suggested to be less compact structure compared to the wild-type and the mutant riboswitches, thus giving a more clear description of the volume of Mn2+-free and wild-type riboswitches as described earlier by Rg calculations (see Fig. 3a). The first two eigenvectors were selected for the evaluation of essential dynamics of the riboswitches since they corresponded to the highest eigenvalues demonstrated by visible kink in the scree plot (Fig. 6).
![]() | ||
Fig. 6 Screen plot depicting first 20 eigenvalues of covariance matrix obtained during principal component analysis of Mn2+-free (black), wild-type (red) & mutant (blue) yybP–ykoY riboswitches. |
Conformational sampling of Mn2+-free, wild-type and mutant riboswitches in the essential subspace was depicted in Fig. 7, illustrating the conformations along the two eigenvectors (EV1 and EV2) projected by all-atoms. It was observed that Mn2+-free riboswitch occupied larger and distinct conformational subspace than wild-type riboswitch, however, atomic positions in the Mn2+-free case was also overlapping with that of the wild-type case showing some differences in the sampling of Mn2+-free and wild-type riboswitches, thus pointing toward the change in the conformational space of the riboswitch due to Mn2+ ion (Fig. 7a). On the other hand, the mutant riboswitch occupied smaller conformational subspace than the wild-type riboswitch and atomic positions nearly overlapped with that of the wild-type riboswitch demonstrating that the mutation stabilized the overall structure of riboswitch and also decreased conformational subspace occupied by the mutant riboswitch (Fig. 7b). To further ascertain the mutation effect on the riboswitch dynamics in particular the MB site and the P1.1 switch helix, atomic fluctuations of nucleotides, and dynamic correlation between the MB site and the P1.1 switch helix, accountable for the essential subspace dynamics were analyzed via RMSF calculations along EV1 during PCA, since it reduced the functionally irrelevant motions and unmasked functionally significant motions (Fig. 8). The PCA-based RMSF helped to gain further insight into the essential dynamics of the riboswitch that was not achieved earlier by simple RMSF calculations and DCCM analysis (see Fig. 4 and 5).
Fig. 8a showed that the mutation decreased the flexibility of the overall structure of the riboswitch, as low fluctuations were observed in the RMSF profile of the mutant riboswitch thus revealing that the mutation stabilized the overall structure which also correlated the backbone RMSD of the mutant riboswitch showing enhanced structural stability (see Fig. 2a). PCA based RMSF analysis also demonstrated low fluctuations of the MB site nucleotides in the mutant riboswitch such as G9, U39, U44 and C45 except C40 and U41 depicted in Fig. 8d–f. The RMSF profile also revealed that the mutation affected the overall dynamics of the riboswitch but not the flexibility of the MB site nucleotides rather than make the metal site more rigid. On the other hand, all nucleotides (A4, A5, C96, C97, U98, U99, U100, C101 and C102) except A3 in the P1.1 switch helix of the mutant riboswitch showed relatively high fluctuations depicted in Fig. 8b and c, demonstrating that the mutation increased the helix flexibility. This finding complemented an experimental study reporting Mn2+ binding at the MB site that stabilized P1.1 switch helix in yybP-ykoY riboswitch, and this stabilizing effect was lost upon Mn2+ sensing adenosine mutation.14
The first two principal components, PC1 and PC2 of each Mn2+-free, wild-type and mutant riboswitch were employed to obtain information of the global minimum-energy (stable) conformations responsible for the essential dynamics, which were depicted in the form of Gibbs free energy landscapes (FEL) shown in Fig. 9. The lowest energy states or stable conformations represented by blue colored maps whereas red colored maps corresponded to high-energy conformations. Each riboswitch produced different FEL pattern; the Mn2+-free case showed two free energy basins (valleys) with pointed and flat ends, which were less broad than the single energy basin observed in the wild-type case (Fig. 9a and b). The minima with flat ends corresponded to clusters of stable conformations, the minima with pointed end represented a stable conformation, and energy barriers were large between basins unlike between the minima in basins. Since narrow energy valleys indicated small number of conformations, the Mn2+-free case yielded less number of stable conformations in comparison to the wild-type riboswitch. FEL of the mutant riboswitch depicted broader minima with little energy barrier than that of the wild-type riboswitch; minima with flat ends indicated clustering of the stable or low energy conformations, and a canonical end suggested the presence of a stable conformation. In the case of mutant riboswitch, large number of stable conformations were produced than the wild-type riboswitch as reflected by FEL pattern (Fig. 9c). The stable conformations indicated by the minima were responsible for essential dynamics in the riboswitch, and as the number of stable conformations decreases flexibility of the system increases and vice versa. Consequently, FEL pattern served to figure out the effect of the mutation and Mn2+ absence on the sampled essential subspace.
The coordination sphere of the metal binding MB site was also analyzed by evaluating dynamics of its geometric configuration in terms of bond lengths and bond angles. Fig. 10 illustrated time evolution coordination distances between the metal ion and the coordinating atoms of the nucleotides in the Mn2+-free, the wild-type and the mutant riboswitches. Based on the distance plots, the averaged coordination distances were computed for the metal-nucleotide atomic pairs in each riboswitch and listed in Table 1. In case of Mn2+-free riboswitch, distances between Mg2+ and the MB site residues were averaged to be less than 2 Å throughout the simulation time except Mg–N7A41 distance (2.31 ± 0.11 Å) that was same as the Mn–N7A41 distance (2.30 ± 0.09 Å) in the wild-type riboswitch, thus demonstrating the similar structure behavior of these two metal ions. Nevertheless, based on the distance plots, the average coordination distances were computed between the metal ions and nucleotides of the MB site in each riboswitch. In case of Mn2+-free riboswitch, distances between Mg2+ ion and the MB site residues were averaged to be less than 2 Å throughout the simulation time except the Mg–N7A41 distance (2.31 ± 0.11 Å) that was same as the Mn–N7A41 distance (2.30 ± 0.09 Å) in the wild-type riboswitch. However, the amplitude of these two distance plots had a characteristic constrast as the distance in Mn2+-free riboswitch was slightly larger than the wild-type riboswitch (Fig. 10a), reflecting less stability of the MB site in the absence of essential Mn2+ ion. Likewise, a relatively large amplitude was observed for Mg–OP1G9 distance (1.96 ± 0.05 Å) compared to Mn–OP1G9 distance (2.02 ± 0.05 Å) indicating less stability of the metal coordination site even though the distance in the Mn2+-free riboswitch was low than the wild-type case (Fig. 10b). Other metal-nucleotide distances in Mn2+-free riboswitch were also smaller than those of the wild-type riboswitch but the Mg–nucleotide distances had similar amplitudes (Fig. 10c), suggesting its rigidity than those in the wild-type riboswitch. Based on the distance plots, it was evaluated as the Mg2+ ion unlike Mn2+ ion showed constant variations while making interactions with the MB site nucleotides, whereas Mn2+ being an essential constituent of the riboswitch as of the wild-type, forms dynamically stable metal coordination site which could be attributed to the function of the Mn2+-sensing riboswitches. The rigidity of the metal-nucleotide distances in Mn2+-free riboswitch also corroborated the RMSF data, demonstrating structural flexibility of the riboswitch (see Fig. 4d–f), however, irrespective of very similar RMSD values for the MB site nucleotides in Mn2+-free and wild-type riboswitches. The average distances between Mn2+ and MB site residues in the wild-type riboswitch were very close as that of the mutant riboswitch as listed in Table 1, except Mn–N7A41 distance (2.30 ± 0.09 Å) in the wild-type that was larger than Mn–O4U41 distance (2.12 ± 0.09 Å) in the mutant form shown in Fig. 10d, in which U41 nucleotide was not having N7 atom. This indicated that the mutation led to shorter and rigid contacts between Mn2+ and MB site residues which also correlated the PCA results signifying that the flexibility of MB site decreases with mutation (see Fig. 8d–f). Fig. 10d depicted a fragile coordination distance of Mn2+ with O4 atom of uracil in the mutant riboswitch that was formed at the beginning of the simulation time as deduced from the distance plot for the first 5 ns, which upon equilibration at ∼1 ns got decreased thus causing the coordination shell to be stabilized, and this finding was found to be correlating with an experimental study11 (Fig. 10f). According to the experimental study, the A41U mutation led Mn2+ coordination scheme to adopt certain changes in such a way that the coordination scheme of MA site nucleotides and G10, U39, C40, and C45 of MB site remained unaffected, whereas the backbone phosphate atoms of U44 was shifted away and replaced by a water molecule as resolved by 2.2 Å electron density.11 Conversely, the simulation study showed that the contact between OP1 atom of U44 nucleotide and Mn2+ remained stable and was not substituted by water as demonstrated by the distance plots illustrated in Fig. 10e and f.
Bond distance | MD simulation (Å) | X-ray (Å)11 |
---|---|---|
Mn2+-free riboswitch | ||
Mg–OP1G9 | 1.96 ± 0.05 Å | — |
Mg–OP1U44 | 1.95 ± 0.05 Å | — |
Mg–N7A41 | 2.31 ± 0.11 Å | — |
Mg–OP2C45 | 1.97 ± 0.05 Å | — |
Mg–OP2U39 | 1.95 ± 0.05 Å | — |
Mg–OP2C40 | 1.95 ± 0.05 Å | — |
![]() |
||
Wild-type riboswitch | ||
Mn–OP1G9 | 2.02 ± 0.05 Å | 2.21 |
Mn–OP1U44 | 2.02 ± 0.05 Å | 2.23 |
Mn–N7A41 | 2.30 ± 0.09 Å | 2.26 |
Mn–OP2C45 | 2.02 ± 0.04 Å | 2.27 |
Mn–OP2U39 | 2.02 ± 0.05 Å | 2.30 |
Mn–OP2C40 | 2.01 ± 0.05 Å | 2.33 |
![]() |
||
Mutant riboswitch | ||
Mn–OP1G9 | 2.02 ± 0.05 Å | — |
Mn–OP1U44 | 2.02 ± 0.05 Å | — |
Mn–O4U41 | 2.12 ± 0.09 Å | — |
Mn–OP2C45 | 2.02 ± 0.05 Å | — |
Mn–OP2U39 | 2.01 ± 0.05 Å | — |
Mn–OP2C40 | 2.02 ± 0.05 Å | — |
After having analyzed the effects of Mn2+ absence and the mutation in terms of structure and dynamics of yybP–ykoY riboswitch, and in particular, the flexibility (stability) of the P1.1 switch helix and the MB site, the A41U mutation effect was quantified as the binding free energy of the mutant riboswitch. The free energy was estimated with respect to the Mn2+ ion bound to the wild-type and the mutant riboswitches through relative binding free energy difference (ΔΔGbind) via free energy simulations based on TI approach,43 in order to evaluate the relative stability of the two Mn2+-bound riboswitches. Fig. 11 illustrated the plot of relative free energy for each interval of λ that is between adjacent Hamiltonian for both riboswitches, where process 1 illustrated binding free energies in the presence of Mn2+, computed as a function of coupling parameter, λ, whereas binding energy calculation (without Mn2+) was carried out in process 2 to obtain useful information employing non-physical thermodynamic cycle (Scheme 1) for the relative binding free energy difference based on the following relation,
ΔΔGbind = ΔGBBind − ΔGABind = ΔG![]() ![]() | (2) |
Table 2 listed the free energies computed during process 1 (ΔGB1) and process 2 (ΔG
B2), and the relative binding free energy difference (ΔΔGbind) revealed that the mutant riboswitch was more stable than the wild-type riboswitch with relative binding energy difference of 23.75 kcal mol−1. The increase in the binding free energy for the mutant riboswitch was attributed to the more rigid contact between O4U41 and Mn2+ compared to the Mn–N7A41 coordination bond in the wild-type riboswitch (Table 1 and Fig. 10). Based on the free energy data, it was obvious that the A41U mutation stabilizes A41U yybP–ykoY riboswitch of L. lactis.
Process | Free energy (kcal mol−1) |
---|---|
ΔGA→B1 (process 1) | 45.09 ± 0.09 |
ΔGA→B2 (process 1) | 21.34 ± 0.04 |
ΔΔGbind | 23.75 ± 0.05 |
This journal is © The Royal Society of Chemistry 2022 |