Calculation of absolute protein-ligand binding free energy from computer simulations

Abstract

A general methodology for calculating the equilibrium binding constant of a flexible ligand to a protein receptor is formulated on the basis of potentials of mean force. The overall process is decomposed into several stages that can be computed separately: the free ligand in the bulk is first restrained into the conformation it adopts in the bound state, position, and orientation by applying biasing potentials, then it is translated into the binding site, where it is released completely. The conformational restraining potential is based on the root-mean-square deviation of the peptide coordinates relative to its average conformation in the bound complex. Free energy contributions from each stage are calculated by means of free energy perturbation potential of mean force techniques by using appropriate order parameters. The present approach avoids the need to decouple the ligand from its surrounding (bulk solvent and receptor protein) as is traditionally performed in the double-decoupling scheme. It is believed that the present formulation will be particularly useful when the solvation free energy of the ligand is very large. As an application, the equilibrium binding constant of the phosphotyrosine peptide pYEEI to the Src homology 2 domain of human Lck has been calculated. The results are in good agreement with experimental values.

Fig. 1.

Crystallographic structure of the p56^lck SH2 domain in complex with a pYEEI peptide (21).

Fig. 2.

Schematic representation of the local reference frame used to define the position and orientation of the ligand relative to the receptor protein and construct the restraining potentials. The spherical coordinate system establishing the position of the ligand relative to the protein are the P3–L1 distance r₁, the P2–P1–L1 angle θ₁, and the P3–P2–P1–L1 dihedral angle φ₁. The Euler angles needed to define the orientation of the ligand relative to the protein are the P3–L1–L2 angle Θ₁, the P2–P1–L1–L2 dihedral angle Φ₁, and the P1–L1–L2–L3 dihedral angle ψ₁.

Fig. 3.

Calculated PMF w(ξ) of the ligand conformational degrees of freedom as function of the rmsd ξ in the bound (A) and bulk (B) states.

Fig. 4.

Calculated PMF W(r₁) as a function of the radial distance r₁ between the ligand center of mass relative to the protein.