Monte carlo simulations of proteins at constant pH with generalized born solvent, flexible sidechains, and an effective dielectric boundary

Abstract

Titratable residues determine the acid/base behavior of proteins, strongly influencing their function; in addition, proton binding is a valuable reporter on electrostatic interactions. We describe a method for pK(a) calculations, using constant-pH Monte Carlo (MC) simulations to explore the space of sidechain conformations and protonation states, with an efficient and accurate generalized Born model (GB) for the solvent effects. To overcome the many-body dependency of the GB model, we use a "Native Environment" approximation, whose accuracy is shown to be good. It allows the precalculation and storage of interactions between all sidechain pairs, a strategy borrowed from computational protein design, which makes the MC simulations themselves very fast. The method is tested for 12 proteins and 167 titratable sidechains. It gives an rms error of 1.1 pH units, similar to the trivial "Null" model. The only adjustable parameter is the protein dielectric constant. The best accuracy is achieved for values between 4 and 8, a range that is physically plausible for a protein interior. For sidechains with large pKa shifts, ≥2, the rms error is 1.6, compared to 2.5 with the Null model and 1.5 with the empirical PROPKA method.

Continuum Electrostatics Approaches to Calculating pKas and Ems in Proteins

“…The high dimensionality of this sum makes it impractical to evaluate for most protein systems. Instead of direct evaluation, ρ β is often calculated through limited conformational sampling; e.g., via Monte Carlo (MC) simulations (Song, Mao & Gunner, 2009; Polydorides & Simonson, 2013). Conformational degrees of freedom can range from sampling side chain rotameric states to relatively inexpensive optimization of steric clashes (Song, 2011) and simple enumeration of different tautomeric forms for the hydrogen position on protonated side chains.…”

“…Monte Carlo methods can be used to incorporate side chain conformer sampling on a rigid protein backbone (Rabenstein, Ullmann & Knapp, 1998; Song et al, 2009; Polydorides & Simonson, 2013). Such sampling attempts to explicitly evaluate the ensemble average described above and thus incorporates a significant amount of side chain response to charge state changes.…”

“…This sets a challenging metric for evaluating the performance of more sophisticated titration prediction methods. For example, the RMSD using modern Monte Carlo methods with continuum electrostatics force field, an ε solute between 4 and 8 and addition of conformational sampling or modification of the dielectric boundary and distribution can be between 0.9 and 1.1 (Song et al, 2009; Polydorides & Simonson, 2013; Wang et al, 2015). However, informatics-based methods such as PROPKA3 can do much better while sacrifice the underlying physical interactions for knowledge-based potentials (Olsson, 2011).…”

Comparing pairwise‐additive and many‐body generalized Born models for acid/base calculations and protein design

“…Acid/base calculations are related to CPD, since proton binding/unbinding to a side chain can be treated formally as a kind of mutation. Indeed, p K a calculations have been used recently to test CPD model and software packages …”

Model for screened, charge-regulated electrostatics of an eye lens protein: Bovine gammaB-crystallin

“…In so doing it is important to recognize that other factors can contribute that we do not incorporate, including changes in conformation that are important in allosteric effects and in calculations that use more microscopic representations of dielectric properties [31,37–42], hydration and the hydrophobic effect [43], hydrogen-bonding [44,45], static dipole potentials [45–47], and ion binding [48], each of which can also be expected to produce changes in local charge patterns. We note that γ B-crystallin is believed to have a fairly robust internal structure; for example, circular dichroism measurements [49] showed no significant spectroscopic changes between −20 °C and 60 °C, though this does not rule out the possible role of conformational flexibility in affecting the present model.…”