Computing Binding Enthalpy Values

(Note: These tutorials are meant to provide illustrative examples of how to use the AMBER software suite to carry out simulations that can be run on a simple workstation in a reasonable period of time. They do not necessarily provide the optimal choice of parameters or methods for the particular application area.)
Copyright Andrew T. Fenley & Michael K. Gilson 2014

Section 5: Potential Pitfalls in Computing Binding Enthalpies

While the direct method of computing binding enthalpies is simple in concept and execution, there are many subtle details that can lead to erroneous results. Here we will touch upon these pitfalls and refer the reader to the Supporting Information of Fenley et al. for additional information.

Identical Simulation Parameters

One, somewhat obvious detail is making sure all of the simulation parameters are identical for the production simulations of each system. In particular, one must use the same ensemble settings, nonbonded cutoff (cut) value, and PME settings. If the parameters are not identical, then differences in the ensemble average of the absolute potential energies calculated throughout the simulations will include nonphysical contributions directly related to the mismatch in parameters.

Proper Balance of Atoms

Care must be taken in setting up the systems such that the total number of solute, explicit solvent, and ionic species are equal between the bound (complex and solvent) and unbound (free host and free guest) states. Furthermore, the presence of charged solutes creates the additional constraint that the number of explicit solvent molecules used to solvate the charged solutes is equal across the simulations. For example, let's assume the charge of the host is neutral but the charge of the guest is +1, then the number of solvent molecules in the free guest system should match the number of solvent molecules in the complex system. Otherwise, if the imbalance of solvent molecules is large enough, the binding enthalpy estimate will contain a potentially non-negligible heat of dilution contribution.

NPT versus NVT

The use of the NVT ensemble for the production simulations has a few notable benefits over the use of the NTP ensemble. The primary benefit is computational throughput, particularly when compared to the Berendsen barostat as implemented in AMBER. Also, when recomputing potential energies from a trajectory file, an MPI version of Sander (sander.MPI) can be used if the box volume does not change, i.e. NVT. Otherwise, the single thread implementation of Sander must be used. Note, this is only a concern when one is interested in computing molecular contributions, e.g. solvent-solvent, to the total binding enthalpy. Finally, the long range, vdW correction (vdwmeth=1) is a constant value for every frame if the box volume remains fixed, and thus, easy to ascertain its contribution to the binding enthalpy.

The major drawback to the use of the NVT ensemble is the requirement of a suitable amount of equilibration in the NPT ensemble such that an accurate estimate of the mean box volume is obtainable, as discussed in Section 2 of this tutorial. A failure to do this step properly can result in the binding enthalpy estimate to be off by hundreds of kcal/mol! Unfortunately, the suitable amount of simulation time during the NPT equilibration stage is difficult to determine a priori. Also, only use volume measurements of the simulation box after the density of the system has stabilized.

A new Monte Carlo barostat was introduced in AMBER 14 with a throughput similar to NVT. However, at this time, we have not sufficiently tested this barostat with computing binding enthalpies and thus cannot endorse its use for such computations at this time.

System Size

Finally, we caution regarding the use of the direct approach for calculating binding enthalpies of large systems, e.g. protein-ligand. Converging the total potential energy of these systems to the same precision as in this tutorial could take a significant amount of simulation time; perhaps tens or hundreds of microseconds or more! That said, we think that current computer hardware is capable of reaching that level of simulation so long as the researcher(s) is(are) willing to dedicate a significant amount of compute resources and time to the problem.

Click here to go back to the Introduction

Click here to go back to Section 4