A Combination of Docking, QM/MM Methods, and MD Simulation for Binding Affinity Estimation of Metalloprotein Ligands

Abstract

To alleviate the problems in the receptor-based design of metalloprotein ligands due to inadequacies in the force-field description of coordination bonds, a four-tier approach was devised. Representative ligand-metalloprotein interaction energies are obtained by subsequent application of (1) docking with metal-binding-guided selection of modes; (2) optimization of the ligand-metalloprotein complex geometry by combined quantum mechanics and molecular mechanics (QM/MM) methods; (3) conformational sampling of the complex with constrained metal bonds by force-field based molecular dynamics (MD); and (4) a single point QM/MM energy calculation for the time-averaged structures. The QM/MM interaction energies are, in a linear combination with the desolvation-characterizing changes in the solvent-accessible surface areas, correlated with experimental data. The approach was applied to structural correlation of published binding free energies of a diverse set of 28 hydroxamate inhibitors to zinc-dependent matrix metalloproteinase 9 (MMP-9). Inclusion of step 3 and step 4 significantly improved both correlation and prediction. The two descriptors explained 90% of variance in inhibition constants of all 28 inhibitors, ranging from 0.08 to 349 nM, with the average unassigned error of 0.318 log units. The structural and energetic information obtained from the timeaveraged MD simulation results helped understand the differences in binding modes of related compounds.

Improved estimation of ligand–macromolecule binding affinities by linear response approach using a combination of multi-mode MD simulation and QM/MM methods

“…Conformational sampling using 200-ps MD simulations did not improve the statistics significantly (r 2 = 0.589). We showed previously [18,32] that these correlations can be amended by a careful selection of the used intervals of MD simulations.…”

“…The QM/MM optimization of the best FlexX poses showed that the hydroxamate groups of all 28 derivatives created, along with the histidine triad of MMP-9, rather similar trigonal bipyramidal coordination spheres around the catalytic zinc [18] (Figure 1). Glu402 formed the hydrogen bonds with the nonionized hydroxyl group of the hydroxamates.…”

“…(1). We have recently shown that this modification of the LR approach with electrostatic and van der Waals energies improves both descriptive and predictive abilities of the correlation [18].…”

“…The computation details were described previously [18,22]. Briefly, the protocol consisted of: (1) FlexX docking of the ligands to the processed structure of MMP-9 taken from the Protein Data Bank [30] (file 1GKC [31]), with the pose selected primarily on the basis of acceptable distances between the catalytic zinc and the hydroxamate atoms; (2) QM/MM minimization of the selected poses, with the QM region including the ligand, tha catalytic zinc, Glu402 and the His triad as shown in Figure 1; (3) (Table S1).…”

“…The LR decomposition of the free energy of binding, ΔG b , reads: (1) where K is the association constant, R is the universal gas constant, T is temperature, E denotes the energy (van der Waals energy and electrostatic energy, scaled by separate coefficients α and β [8,17] or overall QM/MM energy scaled by α [18]), PSA and NPSA represent polar and nonpolar solvent-accessible surface areas [17], and γ, δ, and κ are adjustable coefficients. The angle brackets denote the ensemble averages and Δ indicates the difference between the ensemble averages in the bound and free ligand states.…”