Re: [AMBER] Reference energies in CpHMD (constant pH MD) for HIP

From: Brian Radak <>
Date: Thu, 28 Feb 2019 12:14:41 -0500

Just my 2 cents:

I don't have much firsthand experience with the specific constant-pH MD
code in AMBER, but this is something I also had to grapple with while
implementing a similar approach in NAMD. Our paper, as well as another from
the GROMACS community, deals with the nuances in excruciating detail. The
even older work from Shen and Brooks probably also addresses this, but I
haven't looked at it recently.

Jason (I believe) is speaking to the issue of microscopic versus
macroscopic values. In general, the reference free energies are microscopic
(they describe a single, specific transition) while the reference pKa
values are macroscopic (they describe groups of states that are usually
experimentally indistinguishible). For identical sites (like carboxylates)
the "extra entropy" is exactly the log ratio of identical species on either
side of the transition (1 deprotonated state versus 4 protonated states, so
the shift is roughly log(4)). Any deviations from the sites not being
equivalent are then completely predicted by the force field. For the
carboxylate syn/anti ratio this is probably a good approximation, since you
are essentially asking that the PMF of proton rotation be "correct". For
histidine, this is probably a horrible approximation, because the force
field probably doesn't accurately account for the different bond energies
of HID vs HIE, which is why the separate values of 6.5 and 7.1 are used.

The titration curves of HIS are also themselves tricky. As was pointed out,
if you just track the HID population, the high pH population will not
saturate at 1.0 bc some of those deprotonated states will be HIE instead.
The ratio of the two saturated populations should be P(HIE)/P(HID) =
10**(pKa(HID) - pKa(HIE)) ~= 0.3. You can also work out the exact
saturation values (and the midpoint pH values), but the math is a bit


On Thu, Feb 28, 2019 at 1:08 AM Jason Swails <> wrote:

> On Thu, Feb 21, 2019 at 9:11 PM Cruzeiro,Vinicius Wilian D <
>> wrote:
> > The reference energies can also be obtained by performing TI
> calculations.
> >
> A word of warning here: TI calculations will always underestimate the
> reference energy of protonated species for residues with multiple
> protomers. For instance, carboxylate residues have, broadly speaking, 4
> different protonation states possible -- they can be in the syn- or anti-
> position of either carboxylate oxygen. These are two sets of degenerate
> states, with the anti- positions being higher energy than the syn-
> positions (but the two anti- and two syn- positions are degenerate with
> each other). So Amber defines 5 different states for carboxylate residues
> (deprotonated, syn-O1, anti-O1, syn-O2, anti-O2).
> By contrast, residues like Tyrosine that have only a single protonable site
> (or histidine if you treat each site independently) should have TI energies
> that basically exactly match the reference energies (before correcting for
> the pKa offset).
> As a result, you need to increase the energy gap between the deprotonated
> and protonated states in order to compensate for the fact that MC attempts
> on protonated states are 4x more likely to happen than MC attempts on the
> deprotonated state.
> Mongan gives a good description of this in his original 2004 paper
> introducing the method.
> Good luck,
> Jason
> --
> Jason M. Swails
> _______________________________________________
> AMBER mailing list
AMBER mailing list
Received on Thu Feb 28 2019 - 09:30:02 PST
Custom Search