Standard binding free energies from computer simulations: What is the best strategy?

Abstract

Accurate prediction of standard binding free energies describing protein:ligand association remains a daunting computational endeavor. This challenge is rooted to a large extent in the considerable changes in conformational, translational and rotational entropies underlying the binding process that atomistic simulations cannot easily capture. In spite of significant methodological advances, reflected in a continuously improving agreement with experiment, a characterization of alternate strategies aimed at measuring binding affinities, notably their respective advantages and drawbacks, is somewhat lacking. Here, two distinct avenues to determine the standard binding free energy are compared in the case of a short, proline-rich peptide associating to the Src homology domain 3 of tyrosine kinase Abl. These avenues - one relying upon alchemical transformations and the other on potentials of mean force (PMFs) - invoke a series of geometrical restraints acting on collective variables designed to alleviate sampling limitations inherent to classical molecular dynamics simulations. The experimental binding free energy of ΔG_bind = -7.99 kcal/mol is well reproduced by the two strategies developed herein, with ΔG_bind = -7.7 for the alchemical route and ΔG_bind = -7.8 kcal/mol for the alternate PMF-based route. In detailing the underpinnings of these numerical strategies devised for the accurate determination of standard binding free energies, many practical elements of the proposed rigorous, conceptual framework are clarified, thereby paving way to tackle virtually any recognition and association phenomenon.

Figure 1

Definition of a frame of reference utilized to characterize the binding of peptide APSYSPPPPP (licorice-like representation colored by atom type) to the SH3 domain of the Abl kinase (secondary-structure representation). For both the protein and the ligand, three groups of atoms are arbitrarily selected, forming two triplets {P₁, P₂, P₃} and {L₁, L₂, L}, respectively. The position of APSYSPPPPP with respect to the SH3 domain is defined by means of the set of spherical coordinates {r; θ ϕ}, where r is the P1-L1 distance, θ is the L1-P1-P2 angle, and ϕ is the P1-P1-P2-P3 angle. The relative orientation of the substrate is expressed using the three Euler angles {Θ (P1-L1-L2), Φ (P1-L1-L2-L3), Ψ(P2-P1-L1-L2)}.

Figure 2

Potentials of mean force for (A) RMSD restraints and (B) extraction of the ligand from the binding site. (A) PMF for RMSD of p41 in the bound state (solid line) and in the unbound state (dashed line). (B) Separation PMF, W(r), calculated from REMD-US (solid line) and ABF (dashed line).