Protein pKa's from Adaptive Landscape Flattening ... - ACS Publications

Subscriber access provided by Kaohsiung Medical University

Biomolecular Systems

Protein pKa's from adaptive landscape flattening instead of constant-pH simulations Francesco Villa, and Thomas Simonson J. Chem. Theory Comput., Just Accepted Manuscript • DOI: 10.1021/acs.jctc.8b00970 • Publication Date (Web): 15 Nov 2018 Downloaded from http://pubs.acs.org on November 18, 2018

Just Accepted “Just Accepted” manuscripts have been peer-reviewed and accepted for publication. They are posted online prior to technical editing, formatting for publication and author proofing. The American Chemical Society provides “Just Accepted” as a service to the research community to expedite the dissemination of scientific material as soon as possible after acceptance. “Just Accepted” manuscripts appear in full in PDF format accompanied by an HTML abstract. “Just Accepted” manuscripts have been fully peer reviewed, but should not be considered the official version of record. They are citable by the Digital Object Identifier (DOI®). “Just Accepted” is an optional service offered to authors. Therefore, the “Just Accepted” Web site may not include all articles that will be published in the journal. After a manuscript is technically edited and formatted, it will be removed from the “Just Accepted” Web site and published as an ASAP article. Note that technical editing may introduce minor changes to the manuscript text and/or graphics which could affect content, and all legal disclaimers and ethical guidelines that apply to the journal pertain. ACS cannot be held responsible for errors or consequences arising from the use of information contained in these “Just Accepted” manuscripts.

is published by the American Chemical Society. 1155 Sixteenth Street N.W., Washington, DC 20036 Published by American Chemical Society. Copyright © American Chemical Society. However, no copyright claim is made to original U.S. Government works, or works produced by employees of any Commonwealth realm Crown government in the course of their duties.

Page 1 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

Journal of Chemical Theory and Computation

Protein pKa’s from adaptive landscape flattening instead of constant-pH simulations Francesco Villa and Thomas Simonson∗ Laboratoire de Biochimie (CNRS UMR7654), Ecole Polytechnique, Palaiseau, France E-mail: [email protected]

Abstract Protein acid/base constants, or pKa ’s are often computed from Monte Carlo or molecular dynamics simulations at a series of constant pH values. Instead, we propose to adaptively flatten the free energy landscape in the space of protonation states. The flattening is achieved by a Wang-Landau Monte Carlo, where a bias potential is constructed adaptively during an initial phase, such that all protonation states achieve comparable probabilities. Biased ensembles of states are then reweighted by subtracting out the bias and adding a pH-dependent free energy term. Titration curves constructed for three test proteins agreed, within the small numerical uncertainty, with those obtained earlier from the constant-pH approach.

Running title: Protein pKa ’s with adaptive landscape flattening Key words: Monte Carlo simulation, Proteus software, Generalized Born solvent, free energy simulation

1 ACS Paragon Plus Environment

Journal of Chemical Theory and Computation 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60

1

Introduction

Acid/base reactions play an important role in protein chemistry, influencing or participating in protein folding and stability, enzyme catalysis, and ligand or protein binding. 1,2 Conversely, acid/base constants, or pKa ’s, can report on the structural, dynamic, electrostatic, and dielectric properties of proteins. 3,4 Computer simulations are an increasingly powerful tool to obtain pKa ’s, complementary to experiments. 5 Importantly, simulations can provide mechanistic understanding, especially when they employ atomistic Monte Carlo (MC) or molecular dynamics (MD) simulations. 6,7 Thus, simulations can reveal not only the equilibrium constants but which structural rearrangements occur when the protonation state of a particular protein side chain changes. One popular method performs MD or MC under constant-pH conditions. 6,8–11 This is fairly straightforward with MC. With MD, technical difficulties have had to be overcome, arising from the discrete nature of protonation states. All methods then proceed by scanning the solution pH from a low to a high value (or vice versa), in a series of simulations that mimic a titration experiment. Here, we propose a different approach, based on an adaptive flattening of the energy surface as a function of protonation states. The flattening is achieved by constructing an appropriate bias potential during an initial, adaptive simulation phase. 12,13 The bias adaptation rule is analogous to one used in the well-tempered metadynamics method. 13–15 In a second, production phase, thanks to the flattened landscape, all protonation states are (approximately) equally populated. The equilibrium population of any protonation state at any pH can then be recovered a posteriori by substracting out the bias and adding back a pH term. The method is implemented using MC and an implicit solvent model, within the Proteus software, 16 but in principle it is more general. Indeed, landscape flattening was used recently to compute


Page 2 of 28

Page 3 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


ligand binding in the context of the lambda dynamics MD approach. 17,18 Our implementation is an instance of Wang-Landau MC, 12 whereas an MD implementation would be a form of metadynamics. 13 There is also a similarity to multicanonical simulations, which can be thought of as exploring a range of temperatures in a single simulation; 19,20 here, we effectively explore a range of pH values at once. Our approach is of interest for two reasons. First, in some cases, it is more efficient than the usual series of constant-pH simulations (see examples below). Second, it is conceptually interesting and could lead to more complex protocols where not only the acid/base landscape is flattened, but other landscape coordinates as well, related to ligand binding or side chain mutations, for example. Since the method is implemented in the Proteus protein design software, which allows mutations and ligand binding, such applications are in principle already feasible. Indeed, a related method was used recently to design PDZ-binding peptide sequences. 21,22 To test and illustrate the method, we considered three small proteins, BPTI, OMTKY3, and barnase. All three were treated previously with a constant-pH titration method, using the same physical model. 11 This includes a molecular mechanics force field, Amber ff99SB for the protein, a Generalized Born (GB) implicit solvent model, 23,24 the use of an average, effective, protein–solvent dielectric boundary, a fixed protein backbone, and protein side chains that explore a discrete library of conformers or rotamers. Thus, while backbone motions were described implicitly (through the protein dielectric constant), side chain rearrangements due to protonation events were modeled explicitly. These ingredients are fairly typical of constant-pH MC methods. 6 The goal here was not to test the physical model, which has already been benchmarked for pKa ’s. 11 The goal was to illustrate the adaptive method, test its efficiency, and show that it gives the same results as constant-pH MC, within the simulation uncertainty. The test proteins included 12, 13, 3 ACS Paragon Plus Environment


Page 4 of 28

and 14 titratable groups, respectively. The new, adaptive calculations used exactly the same physical model as the earlier, constant-pH calculations, including the same GBSA energy, the same fixed backbone conformations, and the same set of side chain rotamers. Using a single MC simulation of 200 million steps, we collected populations that could be post-processed to produce titration curves spanning a pH range from 0 to 14. The number of MC stes (and CPU time) was about 1/3 the amount used earlier for the constant-pH method. Both the pKa ’s and the full titration curves were the same as those obtained earlier with the earlier method.

2 2.1

Theoretical methods Acid/base Monte Carlo simulation

We recall briefly the standard framework for MC with acid/base reactions. 6,25 We describe proton binding using a classical, molecular mechanics framework. 26 For a deprotonation reaction, say AH → A− + H+ , the equilibrium constant is Ka = [A− ][H + ]/[AH] and the standard reaction free energy is

∆G0 = −kT ln Ka = ckT pKa

(1)

where c = ln 10 and pKa = -log10 Ka . At non-standard concentrations, the reaction free energy is

∆G = ∆G0 + kT ln

[A− ] [A− ][H + ] = ckT (pKa − pH) + kT ln [AH] [AH]


(2)

Page 5 of 28 1 2 3 4 5 6 7 8 9 10 11 12 13 14 15 16 17 18 19 20 21 22 23 24 25 26 27 28 29 30 31 32 33 34 35 36 37 38 39 40 41 42 43 44 45 46 47 48 49 50 51 52 53 54 55 56 57 58 59 60


We use a procedure that compares the deprotonation of an amino acid side chain in the protein environment to that of an analog, or model compound in solution. For each titratable amino acid side chain R, its model compound X(R) is composed of the side chain R with its backbone capped between two methyl groups (N-acetyl-R-N-methylamide). The experimental pKa of X(R), denoted pKmod , is known. a We consider an event where a proton is removed from a particular protein side chain R, with the protein conformation kept fixed. The specific conformation is designated C, which indicates both a choice of side chain conformers and a certain number and location of protons bound to the titratable side chains. Upon deprotonation, the implicit solvent adjusts to the new state, while the proton is released to a proton bath at a given pH. The corresponding free energy change ∆G(C) can be considered a potential of mean force or PMF, since the solute conformation is fixed while the solvent degrees of freedom have been integrated out (through the implicit solvent model). It can be written

∆G(C) = ∆W (C) − µ◦ − ckT pH

(3)

∆W (C) represents the contribution of the protein to the PMF. In practice, ∆W (C) includes a molecular mechanics energy and a contribution from the GB solvent. The deprotonation event can be seen to occur at equal concentrations of the two protein forms, so there is no term involving

[A− ] . [AH]

The standard chemical potential of the proton,

µ◦ , appears separately (it is not included in W ). Its value can be obtained by writing the protonation free energy of the model compound X(R) in the same form:

∆Gmod = ∆W mod − µ◦ − ckT pH


(4)


Page 6 of 28

Applying Eq. (1) to X, we have

∆G(C) = ∆W (C) − ∆W mod − ckT pKamod − pH

(5)

Notice that that X denotes the protonated form of the model compound. The MC simulation should populate the protonation states according to a Boltzmann distribution governed by G(C). This is accomplished by applying the standard, Metropolis-Hastings MC algorithm, 27,28 where deprotonation/protonation moves are associated with the energy change ±∆G(C). Additional MC moves are used to sample side chain conformers (see Numerical methods, below). For a move that deprotonates a side chain of type X, the appropriate value ∆W mod was tabulated earlier. 29 It includes a molecular mechanics and a GB energy contribution, averaged over the side chain rotamers available to X (in each of its protonation states). It also depends on the experimental pKamod value of X.

2.2

Adaptive landscape flattening in protonation space

The protein includes p titratable side chains that can change their protonation states, such as Asp or His. The protonation state of the titratable side chain i at time t is denoted si (t). We make use of a bias potential E b that varies over time and depends on the protonation states of individual positions i and of pairs of positions i, j:

E b (s1 (t), s2 (t), ..., sp (t); t) =

X

Eib (si (t); t) +

i

X

Eijb (si (t), sj (t); t)

(6)

i

Protein pKa's from Adaptive Landscape Flattening ... - ACS Publications

Recommend Documents