A General Picture of Cucurbit[8]uril Host-Guest Binding: Recalibrating Bonded Interactions

Abstract

Atomic-level understanding of the dynamic feature of host-guest interactions remains a central challenge in supramolecular chemistry. The remarkable guest binding behavior of the Cucurbiturils family of supramolecular containers makes them promising drug carriers. Among Cucurbit[n]urils, Cucurbit[8]uril (CB8) has an intermediate portal size and cavity volume. It can exploit almost all host-guest recognition motifs formed by this host family. In our previous work, an extensive computational investigation of the binding of seven commonly abused and structurally diverse drugs to the CB8 host was performed, and a general dynamic binding picture of CB8-guest interactions was obtained. Further, two widely used fixed-charge models for drug-like molecules were investigated and compared in great detail, aiming at providing guidelines in choosing an appropriate charge scheme in host-guest modelling. Iterative refitting of atomic charges leads to improved binding thermodynamics and the best root-mean-squared deviation from the experimental reference is 2.6 kcal/mol. In this work, we focus on a thorough evaluation of the remaining parts of classical force fields, i.e., the bonded interactions. The widely used general Amber force fields are assessed and refitted with generalized force-matching to improve the intra-molecular conformational preference, and thus the description of inter-molecular host-guest interactions. The interaction pattern and binding thermodynamics show a significant dependence on the modelling parameters. The refitted system-specific parameter set improves the consistency of the modelling results and the experimental reference significantly. Finally, combining the previous charge-scheme comparison and the current force-field refitting, we provide general guidelines for the theoretical modelling of host-guest binding.

Keywords: Cucurbit[8]uril; abused drugs; binding mode; force fields; host–guest interaction.

Figure 1

(a) The 3D structure of the CB8 host and the 2D chemical structures of its 7 guests. For the guest G5 Ketamine, due to the small difference between the experimental and ChemAxon predicted pKa and the pH condition that the binding affinity is measured, both the protonated and deprotonated forms of the guest G5 are considered in our modelling. (b) An illustration of the 3D spherical coordinates CV used to bias the simulation.

Figure 1

(a) The 3D structure of the CB8 host and the 2D chemical structures of its 7 guests. For the guest G5 Ketamine, due to the small difference between the experimental and ChemAxon predicted pKa and the pH condition that the binding affinity is measured, both the protonated and deprotonated forms of the guest G5 are considered in our modelling. (b) An illustration of the 3D spherical coordinates CV used to bias the simulation.

Figure 2

(a) The correlations between the MM (i.e., GAFF2 and FM-PM6) and PM6-D3H4X energetics calculated from 25 ns trajectories generated at 300 K in vacuo for the host, (b) those between energetics produced by MM force fields and a different SQM Hamiltonian of RM1 instead of the fitting target PM6-D3H4X. (c,d) The time series of force errors (i.e., Frobenius norm of $\Delta F$) produced by GAFF2 and the refitted FM-PM6 force field calculated with configurations generated at 300 K for the host CB8. Red dots for the force errors larger than 50 kcal/(mol·Å), green for force errors larger than 30 kcal/(mol·Å), blue for errors larger than 10 kcal/(mol·Å), and white for the other small-error points. The first 96 atoms are non-hydrogen heavy atoms, the force errors of which are larger than the other hydrogen atoms. The errors of (e) the energetics and (f) atomic forces produced by GAFF2 and the refitted FM-PM6 force field calculated with configurations generated at 300 K for all host and guest molecules under investigation.

Figure 2

(a) The correlations between the MM (i.e., GAFF2 and FM-PM6) and PM6-D3H4X energetics calculated from 25 ns trajectories generated at 300 K in vacuo for the host, (b) those between energetics produced by MM force fields and a different SQM Hamiltonian of RM1 instead of the fitting target PM6-D3H4X. (c,d) The time series of force errors (i.e., Frobenius norm of $\Delta F$) produced by GAFF2 and the refitted FM-PM6 force field calculated with configurations generated at 300 K for the host CB8. Red dots for the force errors larger than 50 kcal/(mol·Å), green for force errors larger than 30 kcal/(mol·Å), blue for errors larger than 10 kcal/(mol·Å), and white for the other small-error points. The first 96 atoms are non-hydrogen heavy atoms, the force errors of which are larger than the other hydrogen atoms. The errors of (e) the energetics and (f) atomic forces produced by GAFF2 and the refitted FM-PM6 force field calculated with configurations generated at 300 K for all host and guest molecules under investigation.

Figure 3

(a) The correlations between the MM (i.e., GAFF2 and FM-BLYP) and BLYP-D4 energetics calculated from 10 ns trajectories generated at 300 K in vacuo for the host, (b) those between the energetics produced by MM and a different reference level of PM6-D3H4X (i.e., the previous target in FM-PM6 refitting) instead of the target level BLYP-D4. Note that for this reference level, the FM-BLYP parameter set produces larger errors compared with the previous FM-PM6, which is fitted directly to the data at this reference level. (c,d) The time series of errors of atomic forces produced by GAFF2 and the refitted FM-BLYP force field calculated with configurations generated at 300 K for the host CB8. Red dots for the force errors larger than 50 kcal/(mol·Å), green for force errors larger than 30 kcal/(mol·Å), blue for errors larger than 10 kcal/(mol·Å), and white for the other small-error points. The first 96 atoms are non-hydrogen heavy atoms, the force errors of which are larger than the other hydrogen atoms. The errors of (e) energetics and (f) atomic forces produced by GAFF2 and the refitted FM-BLYP force field calculated with configurations generated at 300 K for all host and guest molecules under investigation.

Figure 3

(a) The correlations between the MM (i.e., GAFF2 and FM-BLYP) and BLYP-D4 energetics calculated from 10 ns trajectories generated at 300 K in vacuo for the host, (b) those between the energetics produced by MM and a different reference level of PM6-D3H4X (i.e., the previous target in FM-PM6 refitting) instead of the target level BLYP-D4. Note that for this reference level, the FM-BLYP parameter set produces larger errors compared with the previous FM-PM6, which is fitted directly to the data at this reference level. (c,d) The time series of errors of atomic forces produced by GAFF2 and the refitted FM-BLYP force field calculated with configurations generated at 300 K for the host CB8. Red dots for the force errors larger than 50 kcal/(mol·Å), green for force errors larger than 30 kcal/(mol·Å), blue for errors larger than 10 kcal/(mol·Å), and white for the other small-error points. The first 96 atoms are non-hydrogen heavy atoms, the force errors of which are larger than the other hydrogen atoms. The errors of (e) energetics and (f) atomic forces produced by GAFF2 and the refitted FM-BLYP force field calculated with configurations generated at 300 K for all host and guest molecules under investigation.

Figure 4

(a) The dihedral term that shows the most significant difference between GAFF derivatives describing the conformational preference (stiffness) of the host. The superposition of the host configurations during 20 ns unbiased simulations in explicit solvent with bonded parameters from (b) GAFF, (c) GAFF2, (d) the FM-PM6 parameters initiated from GAFF, (e) the FM-PM6 parameters initiated from GAFF2, and (f) FM-BLYP initiated from GAFF2. (g) The time series of the radius of gyration of the CB8 ring under the GAFF2 and the refitted parameter sets.

Figure 4

(a) The dihedral term that shows the most significant difference between GAFF derivatives describing the conformational preference (stiffness) of the host. The superposition of the host configurations during 20 ns unbiased simulations in explicit solvent with bonded parameters from (b) GAFF, (c) GAFF2, (d) the FM-PM6 parameters initiated from GAFF, (e) the FM-PM6 parameters initiated from GAFF2, and (f) FM-BLYP initiated from GAFF2. (g) The time series of the radius of gyration of the CB8 ring under the GAFF2 and the refitted parameter sets.

Figure 5

(a) The number of contacts between all atoms of the host CB8 and the guest deprotonated G5 and (b) its decomposition by each atom of the host during 1000 ns enhanced sampling simulations with the newly refitted FM-PM6 parameters. The y-axis represents the serial number of host atoms. The first 96 atoms of the host are heavy atoms, and the others are hydrogen atoms. Red dots denote contacts larger than 10, green dots represent contact numbers between 5 and 10, blue ones are those larger than 1, and the other are represented by white dots. The cyan oval provides an example of the typical center-binding pose explored during enhanced sampling simulations, while in the magenta one there is a side-binding mode.

Figure 6

2D $C-C_{Ph}$ free energy surfaces in kcal/mol obtained under the GAFF2 force field.

Figure 6

2D $C-C_{Ph}$ free energy surfaces in kcal/mol obtained under the GAFF2 force field.

Figure 7

2D $C-C_{Ph}$ free energy surfaces in kcal/mol obtained under the newly refitted force field FM-PM6.

Figure 7

2D $C-C_{Ph}$ free energy surfaces in kcal/mol obtained under the newly refitted force field FM-PM6.

Figure 8

2D $C-C_{Ph}$ free energy surfaces in kcal/mol obtained under the newly refitted force field FM-BLYP.

Figure 8

2D $C-C_{Ph}$ free energy surfaces in kcal/mol obtained under the newly refitted force field FM-BLYP.

Figure 9

Correlation between the binding affinities obtained from our computational modeling and experimental reference for CB8–guest systems. The exact values of the binding affinities are presented in Table 1.