Subject: Summary: How to calculate pKa values in large proteins?...

From: Vladislav Vassiliev <[email-masked]>
Date: Mon, 12 May 1997 17:59:46 -0400 (EDT)

Dear CCL members,

Recently, I sent out the following post: 

+++++++++++ BEGIN ORIGINAL MESSAGE ++++++++++++++

Let's imagine I have coordinates of a real protein (for example, from
PDB) and now I would like to know the pKa values of all the residues
(Asp's, Glu's, Lys's, Arg's, His's) at a given pH. What kind of programs
could predict these pKa values?

I would also be interested in any references concerning the approaches of
this kind of prediction.

++++++++++++ END OF ORIGINAL MESSAGE ++++++++++++++++

Here is a summary of the responses.  Many thanks to all who responded. 

++++++++++++++++++++++++++++++++++++++++++++++++++++++
From: Paul Beroza <[email-masked]>

Several research groups have investigated this problem, 
commonly referred to as the "multiple site titration problem".  

The most widespread approach is to assume that the pKa of a titrating 
residue is perturbed from its solution value by the electrostatic 
environment it experiences in the protein.  The electrostatic energies 
involved can be calculated from continuum models using finite difference,
and other, numerical methods.

I've appended a list of references dealing with the subject
(liberally sprinkled with references to my own work, of course).
There are many more.

I think the Bashford and Karplus paper is a good place to start,
though more recent work uses more sophisticated models.

Alternatively, Arieh Warshel has a more microscopic PDLD model 
(protein dipoles langevin dipoles) which has also been applied to pKa
calculation in proteins.

Several available computer programs can do these sorts of calculations:

        MEAD    (from Don Bashford ftp.scripps.edu:pub/bashford)
        pep     (my finite difference code ftp.scripps.edu:pub/beroza)
        DELPHI  (from Barry Honig's group)
        UHBD    (from Andy McCammon's group)
        polaris (from Arieh Warshel's group)

---------------------------------------------------------------

@ARTICLE{Russel,
        AUTHOR="S. T. Russel and A. Warshel",
        JOURNAL="J. Mol. Biol.",
        TITLE="The Energetics of Ionized Groups in Bovine
                Pancreatic Trypsin Inhibitor",
        VOLUME="185",
        PAGES="389-404",
        YEAR="1985"
        }

@ARTICLE{Bashford90,
        AUTHOR="D. Bashford and M. Karplus",
        JOURNAL="Biochemistry",
        TITLE=" {${\rm pK_a}$'s}
                of Ionizable Groups in Proteins:
                Atomic Detail from a Continuum Model",
        VOLUME="29",
        PAGES="10219-10225",
        YEAR="1990"
        }

AUTHOR="P. Beroza and D. R. Fredkin and M. Y. Okamura and G.  Feher",
        JOURNAL="Proc. Natl. Aca. Sci.",
        TITLE="Protonation of Interacting Residues in a
                Protein by a {M}onte {C}arlo
                Method: {A}pplication to Lysozyme and the
                Photosynthetic Reaction Center of
                {{\it Rhodobacter sphaeroides}}",
        VOLUME="88",
        PAGES="5804-5808",
        YEAR="1991"
        }

@ARTICLE{Oberoi93,
        AUTHOR="H. Oberoi and N. M. Allewell",
        JOURNAL="Biophys. J.",
        TITLE="Multigrid Solution of the Nonlinear Poisson-Boltzmann
                Equation and Calculation of Titration Curves",
        VOLUME="65",
        PAGES="48-55",
        YEAR="1993"
        }

@ARTICLE{Yang93,
        AUTHOR="A.-S. Yang and M. R. Gunner and R. Sampogna
                and K. Sharp and B. Honig",
        JOURNAL="Proteins",
        TITLE="On the Calculation of
                {${\rm pK_a}$'s} in Proteins",
        VOLUME="15",
        PAGES="252-265",
        YEAR="1993"
        }

@ARTICLE{Antosiewicz94,
        AUTHOR="J. Antosiewicz and J.A. McCammon and M.K. Gilson",
        TITLE="Prediction of pH-dependent Properties of Proteins",
        JOURNAL="J. Mol. Biol.",
        VOLUME="238",
        PAGES="415-436",
        YEAR="1994"
        }

@ARTICLE{You95,
        AUTHOR="T.J. You and D. Bashford",
        TITLE="Conformation and hydrogen ion titration of proteins:
                 a continuum electrostatic model with conformational
                flexibility",
        JOURNAL="Biophys. J.",
        VOLUME="69",
        PAGES="1721-1733",
        YEAR="1995"
        }

@ARTICLE{Beroza95,
        AUTHOR="P. Beroza and D. R. Fredkin and M. Y. Okamura and 
                G. Feher",
        JOURNAL="Biophys. J.",
        TITLE="Electrostatic Calculations of Amino Acid Titration and
        Electron Transfer, {${\rm Q_A^- Q_B \rightarrow  Q_A Q_B^-}$},
        in the Reaction Center",
        VOLUME="68",
        PAGES="2233-2250",
        YEAR="1995"
        }

@ARTICLE{Beroza96,
        AUTHOR= "P. Beroza and D. R. Fredkin",
        TITLE="Calculation of Amino Acid {${\rm pK_a}$s} in a Protein
        from a Continuum Electrostatic Model: Method and Sensitivity
        Analysis",
        JOURNAL="J. Comp. Chem.",
        VOLUME="17",
        PAGES="1229-1244",
        YEAR="1996"
        }

@ARTICLE{Beroza96,
        AUTHOR= "P. Beroza and D. A. Case",
        TITLE="Including Side Chain Flexibility in
        Continuum Electrostatic Calculations of Protein
        Titration",
        JOURNAL="J. Phys. Chem.",
        VOLUME="100",
        PAGES="20156-20163",
        YEAR="1996"
        }

+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

From: "Andrey Bliznyuk" <[email-masked]>

Two recent references are:

1. J. Comp. Chem. 1996, v.17, N14, p.1633-1644
2. J. Phys. Chem. 1995, v.99, pp.7180-7187

see also a refs inside.

Unfortunately, there is no reliable method to do the predictions. 
Even the best methods [ref.1] may routinely give an error of about 1 
pK unit. The reason, I believe, is that current MM forcefield do not 
have a good charge model, i.e. charges are located at the atomic 
centers and are not polarizable. While it seems to be not very 
limiting while doing interaction energy calculations, the pKa 
calculations fails. I tried several FF (Amber, OPLS, Charm, Gromos, 
CFF91) and could not reproduce the pKa shift going from CH3COOH to 
C2H5COOH. Hence, I would not trust the enzyme pKa calculations.

+++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

From: Stefan Grzybek <[email-masked]>

There are different approaches available for this kind of calculation. The
main differences are how the electrostatic potentials and interaction energies
are calculated and how the interacting charges are treated to calculate
ionization states.

To the first point: It is possible to calculate the electrostatic
interactions simply using Coulomb's law with an appropriately choosen
dielectric constant or to use solutions of the Poisson-Boltzmann equation
to calculate interaction energies. The most programs I know of, use
solutions to the Poisson-Boltzmann equation for calculating interaction
energies.

The second point: First is possible to use the Tanford-Roxby approximation
to solve the problem of interactiong ionization states; however this breaks
down if there are strongly interacting residues present. Better is to
consider the real statistical mechanical average over all possible ionization
states. This however has the drawback that calculating the ionization states
of a small protein with only about 25 to 30 ionizable residues already becomes
problematic. Therefore people have invented the 'reduced site' approximation,
monte carlo methods or mixed tanford-roxby/statistical mechanics methods.

Recent developments in this field focus on the incorporation of protein dynamics
or sidechain flexibility into the treatment of pK-calculations.

Some programs I know of:
MacroDox by Scott H. Northrup (uses the Tanford-Roxby approximation)
  http://pirn.chem.tntech.edu/macrodox.html

MEAD by Donald Bashford ([email-masked])
  ftp://ftp.scripps.edu/pub/electrostatics/

There has also been developed a pk-calculation program suite in the
lab of Barry Honig by Andreas Windemuth. Mail to [email-masked]
or [email-masked] for further information.

Since there is a large amount of literature in this field, I don't append a
list here. If you are interested, please contact me directly.

++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++++

From: Mohamad Aman Jairajpuri <[email-masked]>

I got your forwarded letter enquiring about the calculation of pKa of 
ionizable group in proteins through one of the profs in my department.
I have been involoved in such type of pKa calculations for the last six 
months using finite difference poission boltzmann algorithim developed by 
Dr James Warwicker. These are a complex set of programms in FORTRAN and 
requires VMX/VAX environment for the operation. I have used this programm 
extensively for the calculation of pKa of arginines in several proteins. 
The programm is not commercialized and hence is in an extremely crude 
form and very very difficult to impliment. Fortunately I was able to use 
it after much effort and it finally paved the way for my thesis, (where I 
was able to match experimentally determined pKas with that of the 
theoritical) which I would be writing soon. This programm can be down 
loaded from FTP site (FTP calvados.ukc.ac.uk) If you need any help I would 
be willing to help you.

You can also obtain a programm developed by Dr Gilson and Dr McCammon for 
the calculation of the pKa in proteins. This is a commercial programm 
which is routinely released to academic institution on very nominal costs 
($400 approx). You can contact Dr J. Andrew McCammon at Department of 
Chemistry and Biochemistry, university of California at San Diego, La 
Jolla, California 92093-0365.