Binding free energy decomposition and multiple unbinding paths of buried ligands in a PreQ1 riboswitch

Abstract

Riboswitches are naturally occurring RNA elements that control bacterial gene expression by binding to specific small molecules. They serve as important models for RNA-small molecule recognition and have also become a novel class of targets for developing antibiotics. Here, we carried out conventional and enhanced-sampling molecular dynamics (MD) simulations, totaling 153.5 μs, to characterize the determinants of binding free energies and unbinding paths for the cognate and synthetic ligands of a PreQ1 riboswitch. Binding free energy analysis showed that two triplets of nucleotides, U6-C15-A29 and G5-G11-C16, contribute the most to the binding of the cognate ligands, by hydrogen bonding and by base stacking, respectively. Mg2+ ions are essential in stabilizing the binding pocket. For the synthetic ligands, the hydrogen-bonding contributions of the U6-C15-A29 triplet are significantly compromised, and the bound state resembles the apo state in several respects, including the disengagement of the C15-A14-A13 and A32-G33 base stacks. The bulkier synthetic ligands lead to significantly loosening of the binding pocket, including extrusion of the C15 nucleobase and a widening of the C15-C30 groove. Enhanced-sampling simulations further revealed that the cognate and synthetic ligands unbind in almost opposite directions. Our work offers new insight for designing riboswitch ligands.

Fig 1. Structures of the PreQ₁-bound aptamer domain from the Tte PreQ₁ riboswitch and of the cognate and synthetic ligands.

(A) Left: sequence and secondary structure of the aptamer; middle: three-dimensional structure of the PreQ₁-aptamer complex, with the aptamer shown in both cartoon and surface representations. Nucleotides in the sequence and in the structure are color-matched. Top right: zoomed version showing PreQ₁ in an oblique view to highlight the base stacking with G11 above and with G5 and C16 below. Bottom right: top view highlighting the in-plane hydrogen bonding with U6, C15, and A29. This structure was prepared using coordinates from the PDB entry 6E1W [23], with missing nucleotides copied from PDB entry 3Q50 [24]. (B) Chemical structures of PreQ₁ and PreQ₀. (C) Chemical structures of three synthetic ligands, L₁ to L₃.

Fig 2. Interactions of cognate and synthetic ligands with the PreQ₁ aptamer in cMD simulations with Mg²⁺.

(A) Contributions of individual nucleotides to the binding free energies. A dashed horizontal line is drawn at -1.5 kcal/mol, which separates the pocket-lining nucleotides from the rest of the sequence. Inset: a table listing the correlation coefficients between the individual contributions of any two complexes. (B) Ligand-nucleotide base stacking statistics. Top: illustration of when a nucleobase is (“in”) or is not (“out”) in a stacking position with the ligand. A rectangle is drawn around the ligand rings atoms, with a minimum of 0.5 Å separation (shown by red lines). A cytosine is in an “in” position, with vertical distance drawn as a solid line, as the projection of its center is inside the rectangle; a guanine is in an “out” position, with vertical distance drawn as a dashed line, as the projection of its center is outside the rectangle. Bottom: in-fractions of three nucleobases and their average vertical distances from the ligand rings. (C)-(H) In-plane hydrogen bonds between ligands and nucleobases, shown as dashed lines, in representative structures from cMD simulations.

Fig 3. Distributions and effects of Mg²⁺.

(A) Density contours of Mg²⁺ ions in the Q₁-bound complex, shown as wireframe. Three Mn²⁺ ions in PDB entry 6VUI are shown as ochre spheres; the corresponding Mg²⁺ sites are labeled as M1, M2, and M3. Phosphate groups in G4 and A13 are shown in stick representation, to highlight the bridging role of M2. (B) Two views showing the coordination of the Mg²⁺ ion at the M2 site. Top: coordination of Mg²⁺ by six water molecules and the latter’s hydrogen bonding with the G3 and G4 nucleobases and A14 phosphate. Bottom: two distances, d_4-13 (between G4 O6 and A13 OP1) and d_14-16 (between A14 P and C16 N4), introduced to characterize the effects of the Mg²⁺ ion. (C) Probability densities of d_4-13 and d_14-16, in cMD simulations with and without Mg²⁺ ions. (D) Effect of Mg²⁺ ions on the separation of the L2 loop from the S1 helix in the Q₁-bound form. Two representative structures are superimposed, with the aptamer in the presence of Mg²⁺ shown in the same multi-color scheme as in Fig 1A and the aptamer in the absence of Mg²⁺ shown in a uniform cyan color. In the bottom view, the C15 nucleotides in the two structures are shown in a stick representation.

Fig 4. Contrast between Q₁- and L₁-bound complexes.

(A) Conformations of the L2 loop (U12-A13-A14-C15). Left: conformations adopted in the apo (yellow) and Q₁-bound (orange) forms, from PDB entries 6VUH and 3Q50, respectively. Right: conformations in representative structures of the L₁-bound (yellow) and Q₁-bound (orange) forms from cMD simulations. (B) Comparison between representative Mg²⁺ sites from cMD simulations and Mn²⁺ positions from crystal structures (PDB entries 6VUI and 6VUH). Left: similarity between cMD Mg²⁺ (green) and crystal Mn²⁺ (ochre) positions in the Q₁-bound form. Right: similarity between cMD Mg²⁺ (green) positions in the L₁-bound from and crystal Mn²⁺ (ochre) positions in the apo form. The L₁-bound structure has a larger separation between the ligand and the 3’ end. (C) The orthogonal orientations of the C15-A14-A13 nucleobases in the Q₁- and L₁-bound forms in cMD simulations, leading to one stack or two separate stacks, respectively, with the A32-G33 nucleobases. (D) Density maps in the space of two collective variables: the distance from the center of the binding pocket to the center of the A32 and G33 nucleobases, and the average angle (Θ) between the A13, A14, and C15 nucleobases and the A32 and G33 nucleobases.

Fig 5. Loosening in the back of the binding pocket when bound with synthetic ligands, in cMD simulations with Mg²⁺.

(A)-(C) Density contours of ligands, shown as wireframe in reference to a representative Q₁-bound structure, with the six pocket-lining nucleotides displayed in stick representation. Panel (A) shows contours in red for Q₁ and in green for Q₀; panel (B) shows contours in red for L₁ and in green for L₂; and panel (C) shows contours in red for L_3m and in green for L_3M. The L_3m and L_3M contours were calculated only from snapshots where hydrogen bonding with A29 or U6 was present. (D)-(F) Representative conformations of the Q₁-, L₁-, and L_3m-bound forms, respectively. Ligands are shown in both stick representation and as dot surface. The C1’ atoms of C15 and C30 are connected to define the distance d_15-30 and to illustrate the back door. (G) The probability densities of d_15-30 in the apo form and the five liganded forms. (H) Two views into Q₁ in the binding pocket. The aptamer is shown as gray surface while the ligand is shown as green spheres. Left: front view showing Q₁ exposure; right: back view showing buried Q₁. (I) Corresponding presentation for L₁, except that this ligand is exposed both on the front and on the back.

Fig 6. Unbinding and rebinding pathways of ligands.

(A) The trajectories of ligand centers shown as dots colored according to the MD simulation time. The aptamer and bound ligands are shown in cartoon and stick representations, respectively. Left: Q₁; right: L₃. A plane in purple bisects the binding pocket into the front half and the back half. Two nucleotides defining the front door in the Q₁-bound complex are labeled in red; two nucleotides defining the back door in the L₃-bound complex are labeled in blue. (B) The front and back unbinding paths. (C) The center-to-center distance r between the ligand and the binding pocket and the z coordinate of the ligand center along the normal of the pocket-bisecting plane. The left-most panel illustrates r and z, and presents interpretations of arrow directions and colors that appear in the right three panels, which show time traces of r and z in three metadynamics simulations. Red dashed and green dotted horizontal lines are drawn at r = 1 and 9 Å, respectively, to indicate the times of entrance to and exit from the binding pocket.