Thermodynamic linkage between the binding of protons and inhibitors to HIV-1 protease

The aspartyl dyad of free HIV-1 protease has apparent pK(a)s of approximately 3 and approximately 6, but recent NMR studies indicate that the aspartyl dyad is fixed in the doubly protonated form over a wide pH range when cyclic urea inhibitors are bound, and in the monoprotonated form when the inhibitor KNI-272 is bound. We present computations and measurements related to these changes in protonation and to the thermodynamic linkage between protonation and inhibition. The Poisson-Boltzmann model of electrostatics is used to compute the apparent pK(a)s of the aspartyl dyad in the free enzyme and in complexes with four different inhibitors. The calculations are done with two parameter sets. One assigns epsilon = 4 to the solute interior and uses a detailed model of ionization; the other uses epsilon = 20 for the solute interior and a simplified representation of ionization. For the free enzyme, both parameter sets agree well with previously measured apparent pK(a)s of approximately 3 and approximately 6. However, the calculations with an internal dielectric constant of 4 reproduce the large pKa shifts upon binding of inhibitors, but the calculations with an internal dielectric constant of 20 do not. This observation has implications for the accurate calculation of pK(a)s in complex protein environments. Because binding of a cyclic urea inhibitor shifts the pK(a)s of the aspartyl dyad, changing the pH is expected to change its apparent binding affinity. However, we find experimentally that the affinity is independent of pH from 5.5 to 7.0. Possible explanations for this discrepancy are discussed.     

Impact of ligand and protein desolvation on ligand binding to the S1 pocket of thrombin

In the present study, we investigate the impact of a tightly bound water molecule on ligand binding in the S1 pocket of thrombin. The S1 pocket contains a deeply buried deprotonated aspartate residue (Asp189) that is, due to its charged state, well hydrated in the uncomplexed state. We systematically studied the importance of this water molecule by evaluating a series of ligands that contains pyridine-type P1 side chains that could potentially alter the binding properties of this water molecule. All of the pyridine derivatives retain the original hydration state albeit sometimes with a slight perturbance. In order to prevent a direct H-bond formation with Asp189, and to create a permanent positive charge on the P1 side chain that is positioned adjacent to the Asp189 carboxylate anion, we methylated the pyridine nitrogen. This methylation resulted in displacement of water but was accompanied by a loss in binding affinity. Quantum chemical calculations of the ligand solvation free energy showed that the positively charged methylpyridinium derivatives suffer a large penalty of desolvation upon binding. Consequently, they have a substantially less favorable enthalpy of binding. In addition to the ligand desolvation penalty, the hydration shell around Asp189 has to be overcome, which is achieved in nearly all pyridinium derivatives. Only for the ortho derivative is a partial population of a water next to Asp189 found. Possibly, the gain of electrostatic interactions between the charged P1 side chain and Asp189 helps to compensate for the desolvation penalty. In all uncharged pyridine derivatives, the solvation shell remains next to Asp189, partly mediating interactions between ligand and protein. In the case of the para-pyridine derivative, a strongly disordered cluster of water sites is observed between ligand and Asp189.     

Population density analysis for determining the protonation state of the catalytic dyad in BACE1-tertiary carbinamine-based inhibitor complex

BACE1 is an aspartyl protease with a very relevant role in medicinal chemistry related to Alzheimer Disease since it has demonstrated to be a promising therapeutic target for inhibition and possible control for the progress of the peptide accumulation characteristic of this pathology. The enzymatic activity of this protein is given by the aspartic dyad, Asp93 and Asp289, which can adopt several protonation states depending on the chemical nature of its inhibitors, this is, monoprotonated, diprotonated and di-deprotonated states. In the present study, the analysis of the population density, for a series of protein-inhibitor molecular dynamics simulations, was carried out to identify the most feasible protonation state adopted by the catalytic dyad in the presence of tertiary carbinamine (TC) transition state analog inhibitors. The results revealed that the monoprotonated Asp289i state, in which the Asp93 and Asp289 residue side chains are deprotonated and protonated on the inner oxygen, respectively, is the most preferred in the presence of TC family inhibitors. This result was obtained after evaluating, for all 9 possible protonation state configurations, the individual and combined population densities of a set of parameters sensitive to protonation state of the Aspartic dyad, using an X-ray experimental BACE1/TC crystallographic structure as reference. This case study demonstrates again the usefulness of the concept of population density as a quantitative tool to establish the most stable system settings, among all possible, by measuring the level of occurrence of simultaneous events obtained from a sampling over time. These results will help to clear the phenomena related to the TCs inhibitory pathway, as well as assist in the design of better TC inhibitors against Alzheimer’s protease.     

pH‐dependent conformational dynamics of beta‐secretase 1: A molecular dynamics study

Beta‐secretase 1 (BACE‐1) is an aspartyl protease implicated in the overproduction of β‐amyloid fibrils responsible for Alzheimer disease. The process of β‐amyloid genesis is known to be pH dependent, with an activity peak between solution pH of 3.5 and 5.5. We have studied the pH‐dependent dynamics of BACE‐1 to better understand the pH dependent mechanism. We have implemented support for graphics processor unit (GPU) accelerated constant pH molecular dynamics within the AMBER molecular dynamics software package and employed this to determine the relative population of different aspartyl dyad protonation states in the pH range of greatest β‐amyloid production, followed by conventional molecular dynamics to explore the differences among the various aspartyl dyad protonation states. We observed a difference in dynamics between double‐protonated, mono‐protonated, and double‐deprotonated states over the known pH range of higher activity. These differences include Tyr 71‐aspartyl dyad proximity and active water lifetime. This work indicates that Tyr 71 stabilizes catalytic water in the aspartyl dyad active site, enabling BACE‐1 activity.     

Thermodynamic dissection of the binding energetics of KNI-272, a potent HIV-1 protease inhibitor

KNI-272 is a powerful HIV-1 protease inhibitor with a reported inhibition constant in the picomolar range. In this paper, a complete experimental dissection of the thermodynamic forces that define the binding affinity of this inhibitor to the wild-type and drug-resistant mutant V82F/184V is presented. Unlike other protease inhibitors, KNI-272 binds to the protease with a favorable binding enthalpy. The origin of the favorable binding enthalpy has been traced to the coupling of the binding reaction to the burial of six water molecules. These bound water molecules, previously identified by NMR studies, optimize the atomic packing at the inhibitor/protein interface enhancing van der Waals and other favorable interactions. These interactions offset the unfavorable enthalpy usually associated with the binding of hydrophobic molecules. The association constant to the drug resistant mutant is 100-500 times weaker. The decrease in binding affinity corresponds to an increase in the Gibbs energy of binding of 3-3.5 kcal/mol, which originates from less favorable enthalpy (1.7 kcal/mol more positive) and entropy changes. Calorimetric binding experiments performed as a function of pH and utilizing buffers with different ionization enthalpies have permitted the dissection of proton linkage effects. According to these experiments, the binding of the inhibitor is linked to the protonation/deprotonation of two groups. In the uncomplexed form these groups have pKs of 6.0 and 4.8, and become 6.6 and 2.9 in the complex. These groups have been identified as one of the aspartates in the catalytic aspartyl dyad in the protease and the isoquinoline nitrogen in the inhibitor molecule. The binding affinity is maximal between pH 5 and pH 6. At those pH values the affinity is close to 6 x 10(10) M(-1) (Kd = 16 pM). Global analysis of the data yield a buffer- and pH-independent binding enthalpy of -6.3 kcal/mol. Under conditions in which the exchange of protons is zero, the Gibbs energy of binding is -14.7 kcal/mol from which a binding entropy of 28 cal/K mol is obtained. Thus, the binding of KNI-272 is both enthalpically and entropically favorable. The structure-based thermodynamic analysis indicates that the allophenylnorstatine nucleus of KNI-272 provides an important scaffold for the design of inhibitors that are less susceptible to resistant mutations.     

DAPI binding to the DNA minor groove: a continuum solvent analysis

A continuum solvent model based on the generalized Born (GB) or finite-difference Poisson-Boltzmann (FDPB) approaches has been employed to compare the binding of 4'-6-diamidine-2-phenyl indole (DAPI) to the minor groove of various DNA sequences. Qualitative agreement between the results of GB and FDPB approaches as well as between calculated and experimentally observed trends regarding the sequence specificity of DAPI binding to B-DNA was obtained. Calculated binding energies were decomposed into various contributions to solvation and DNA-ligand interaction. DNA conformational adaptation was found to make a favorable contribution to the calculated total interaction energy but did not change the DAPI binding affinity ranking of different DNA sequences. The calculations indicate that closed complex formation is mainly driven by nonpolar contributions and was found to be disfavored electrostatically due to a desolvation penalty that outbalances the attractive Coulomb interaction. The calculated penalty was larger for DAPI binding to GC-rich sequences compared with AT-rich target sequences and generally larger for the FDPB vs the GB continuum model. A radial interaction profile for DAPI at different distances from the DNA minor groove revealed an electrostatic energy minimum a few Angstroms farther away from the closed binding geometry. The calculated electrostatic interaction up to this distance is attractive and it may stabilize a nonspecific binding arrangement.     

Electrostatic interactions in the GCN4 leucine zipper: substantial contributions arise from intramolecular interactions enhanced on binding.

The GCN4 leucine zipper is a peptide homodimer that has been the subject of a number of experimental and theoretical investigations into the determinants of affinity and specificity. Here, we utilize this model system to investigate electrostatic effects in protein binding using continuum calculations. A particularly novel feature of the computations made here is that they provide an interaction-by-interaction breakdown of the electrostatic contributions to the free energy of docking that includes changes in the interaction of each functional group with solvent and changes in interactions between all pairs of functional groups on binding. The results show that (1) electrostatic effects disfavor binding by roughly 15 kcal/mol due to desolvation effects that are incompletely compensated in the bound state, (2) while no groups strongly stabilize binding, the groups that are most destabilizing are charged and polar side chains at the interface that have been implicated in determining binding specificity, and (3) attractive intramolecular interactions (e.g., backbone hydrogen bonds) that are enhanced on binding due to reduced solvent screening in the bound state contribute significantly to affinity and are likely to be a general effect in other complexes. A comparison is made between the results obtained in an electrostatic analysis carried out calculationally and simulated results corresponding to idealized data from a scanning mutagenesis experiment. It is shown that scanning experiments provide incomplete information on interactions and, if overinterpreted, tend to overestimate the energetic effect of individual side chains that make attractive interactions. Finally, a comparison is made between the results available from a continuum electrostatic model and from a simpler surface-area dependent solvation model. In this case, although the simpler model neglects certain interactions, on average it performs rather well.     

Conformational flexibility of the catalytic Asp dyad in HIV-1 protease: An ab initio study on the free enzyme

The enzyme protease from the human immunodeficiency virus type 1 (HIV-1 PR) is one of the main targets for therapeutic intervention in AIDS. Computer modeling is useful for probing the binding of novel ligands, yet empirical force field-based methods have encountered problems in adequately describing interactions of the catalytic aspartyl pair. In this work we use ab initio dynamic methods to study the molecular interactions and the conformational flexibility of the Asp dyad in the free enzyme. Calculations are performed on model complexes that include, besides the Asp dyad, the conserved Thr26 and Gly27 residues and water molecules present in the active site channel. Our calculations provide proton location and binding mode of the active-site water molecule, which turn out to be different from those of the eukariotic isoenzyme. Furthermore, the calculations reproduce well the structural features of the aspartyl dyad in the protein. Finally, they allow the identification of both dipole/charge interactions and a low-barrier hydrogen bond as important stabilizing factors for the peculiar conformation of the active site. These findings are consistent with site-directed mutagenesis experiments on the 27, 27; positions (Bagossi et al., Protein Eng 1996;9:997-1003). The electric field of the protein frame (included in some of the calculations) does not affect significantly the chemical bonding at the cleavage site. Proteins 2000;39:26-36.     

Structural determinants of trypsin affinity and specificity for cationic inhibitors

The binding free energies of four inhibitors to bovine beta-trypsin are calculated. The inhibitors use either ornithine, lysine, or arginine to bind to the S1 specificity site. The electrostatic contribution to binding free energy is calculated by solving the finite difference Poisson-Boltzmann equation, the contribution of nonpolar interactions is calculated using a free energy-surface area relationship and the loss of conformational entropy is estimated both for trypsin and ligand side chains. Binding free energy values are of a reasonable magnitude and the relative affinity of the four inhibitors for trypsin is correctly predicted. Electrostatic interactions are found to oppose binding in all cases. However, in the case of ornithine- and lysine-based inhibitors, the salt bridge formed between their charged group and the partially buried carboxylate of Asp189 is found to stabilize the complex. Our analysis reveals how the molecular architecture of the trypsin binding site results in highly specific recognition of substrates and inhibitors. Specifically, partially burying Asp189 in the inhibitor-free enzyme decreases the penalty for desolvation of this group upon complexation. Water molecules trapped in the binding interface further stabilize the buried ion pair, resulting in a favorable electrostatic contribution of the ion pair formed with ornithine and lysine side chains. Moreover, all side chains that form the trypsin specificity site are partially buried, and hence, relatively immobile in the inhibitor-free state, thus reducing the entropic cost of complexation. The implications of the results for the general problem of recognition and binding are considered. A novel finding in this regard is that like charged molecules can have electrostatic contributions to binding that are more favorable than oppositely charged molecules due to enhanced interactions with the solvent in the highly charged complex that is formed.     

Experimental and 'in silico' analysis of the effect of pH on HIV-1 protease inhibitor affinity: implications for the charge state of the protein ionogenic groups

The pH dependence of the HIV-1 protease inhibitor affinity was studied by determining the interaction kinetics of a series of inhibitors at three pH values by surface plasmon resonance (SPR) biosensor analysis. The results were rationalized by molecular mechanics based protocols that have as a starting point the structures of the HIV-1 protease inhibitor complexes differing in the protonation states as predicted by our calculations. The SPR experiments indicate a variety of binding affinity pH dependencies which are rather well reproduced by our simulations. Moreover, our calculations are able to pinpoint the possible changes in the charged state of the protein binding site and of the inhibitor that underlie the observed effects of the pH on binding affinity. The combination of SPR and molecular mechanics calculations has afforded novel insights into the pH dependence of inhibitor interactions with their target. This work raises the possibility of designing inhibitors with different pH binding affinity profiles to the ones described here.     

Determinants of ligand binding to cAMP-dependent protein kinase

Protein kinases are essential for the regulation of cellular growth and metabolism. Since their dysfunction leads to debilitating diseases, they represent key targets for pharmaceutical research. The rational design of kinase inhibitors requires an understanding of the determinants of ligand binding to these proteins. In the present study, a theoretical model based on continuum electrostatics and a surface-area-dependent nonpolar term is used to calculate binding affinities of balanol derivatives, H-series inhibitors, and ATP analogues toward the catalytic subunit of cAMP-dependent protein kinase (cAPK or protein kinase A). The calculations reproduce most of the experimental trends and provide insight into the driving forces responsible for binding. Nonpolar interactions are found to govern protein-ligand affinity. Hydrogen bonds represent a negligible contribution, because hydrogen bond formation in the complex requires the desolvation of the interacting partners. However, the binding affinity is decreased if hydrogen-bonding groups of the ligand remain unsatisfied in the complex. The disposition of hydrogen-bonding groups in the ligand is therefore crucial for binding specificity. These observations should be valuable guides in the design of potent and specific kinase inhibitors.     

Optimization of binding electrostatics: charge complementarity in the barnase-barstar protein complex

Theoretical and experimental studies have shown that the large desolvation penalty required for polar and charged groups frequently precludes their involvement in electrostatic interactions that contribute strongly to net stability in the folding or binding of proteins in aqueous solution near room temperature. We have previously developed a theoretical framework for computing optimized electrostatic interactions and illustrated use of the algorithm with simplified geometries. Given a receptor and model assumptions, the method computes the ligand-charge distribution that provides the most favorable balance of desolvation and interaction effects on binding. In this paper the method has been extended to treat complexes using actual molecular shapes. The barnase-barstar protein complex was investigated with barnase treated as a target receptor. The atomic point charges of barstar were varied to optimize the electrostatic binding free energy. Barnase and natural barstar form a tight complex (K(d) approximately 10(-14) M) with many charged and polar groups near the interface that make this a particularly relevant system for investigating the role of electrostatic effects on binding. The results show that sets of barstar charges (resulting from optimization with different constraints) can be found that give rise to relatively large predicted improvements in electrostatic binding free energy. Principles for enhancing the effect of electrostatic interactions in molecular binding in aqueous environments are discussed in light of the optima. Our findings suggest that, in general, the enhancements in electrostatic binding free energy resulting from modification of polar and charged groups can be substantial. Moreover, a recently proposed definition of electrostatic complementarity is shown to be a useful tool for examining binding interfaces. Finally, calculational results suggest that wild-type barstar is closer to being affinity optimized than is barnase for their mutual binding, consistent with the known roles of these proteins.     

Determination of the protonation state of the Asp dyad: conventional molecular dynamics versus thermodynamic integration

The protonation state of the Asp dyad is important as it can reveal enzymatic mechanisms, and the information this provides can be used in the development of drugs for proteins such as memapsin 2 (BACE-1), HIV-1 protease, and rennin. Conventional molecular dynamics (MD) simulations have been successfully used to determine the preferred protonation state of the Asp dyad. In the present work, we demonstrate that the results obtained from conventional MD simulations can be greatly influenced by the particular force field applied or the values used for control parameters. In principle, free-energy changes between possible protonation states can be used to determine the protonation state. We show that protonation state prediction by the thermodynamic integration (TI) method is insensitive to force field version or to the cutoff for calculating nonbonded interactions (a control parameter). In the present study, the protonation state of the Asp dyad predicted by TI calculations was the same regardless of the force field and cutoff value applied. Contrary to the intuition that conventional MD is more efficient, our results clearly show that the TI method is actually more efficient and more reliable for determining the protonation state of the Asp dyad.     

Binding free energy calculations between bovine β‐lactoglobulin and four fatty acids using the MMGBSA method

The bovine dairy protein β‐lactoglobulin (βlg) is a promiscuous protein that has the ability to bind several hydrophobic ligands. In this study, based on known experimental data, the dynamic interaction mechanism between bovine βlg and four fatty acids was investigated by a protocol combining molecular dynamics (MD) simulations and molecular mechanics generalized Born surface area (MMGBSA) binding free energy calculations. Energetic analyses revealed binding free energy trends that corroborated known experimental findings; larger ligand size corresponded to greater binding affinity. Finally, binding free energy decomposition provided detailed information about the key residues stabilizing the complex. © 2014 Wiley Periodicals, Inc. Biopolymers 101: 1010–1018, 2014.     

Tricyclic ureas: A new class of HIV-1 protease inhibitors

A new class of tricyclic ureas containing a conformationally constrained proline was designed with the aid of molecular modeling. Efficient stereoselective intermolecular pinacol coupling represented the highlight of the synthesis. These rigid cyclic ureas are active towards HIV-1 protease, with 9 being the most potent compound (Ki = 9 nM) despite interacting with only three side chain binding pockets of HIV protease.     

Insights into the functional role of protonation states in the HIV-1 protease-BEA369 complex: molecular dynamics simulations and free energy calculations

The molecular mechanics Poisson-Boltzmann surface area (MM-PBSA) method combined with molecular dynamics (MD) simulations were used to investigate the functional role of protonation in human immunodeficiency virus type 1 (HIV-1) protease complexed with the inhibitor BEA369. Our results demonstrate that protonation of two aspartic acids (Asp25/Asp25′) has a strong influence on the dynamics behavior of the complex, the binding free energy of BEA369, and inhibitor–residue interactions. Relative binding free energies calculated using the MM-PBSA method show that protonation of Asp25 results in the strongest binding of BEA369 to HIV-1 protease. Inhibitor–residue interactions computed by the theory of free energy decomposition also indicate that protonation of Asp25 has the most favorable effect on binding of BEA369. In addition, hydrogen-bond analysis based on the trajectories of the MD simulations shows that protonation of Asp25 strongly influences the water-mediated link of a conserved water molecule, Wat301. We expect that the results of this study will contribute significantly to binding calculations for BEA369, and to the design of high affinity inhibitors.     

The activation of the nicotinic acetylcholine receptor by the transmitter

Experimental evidence has been published from isolated guinea pig muscle in vitro, and from direct ligand binding to receptors from T. californica, indicating that two agonist ions react with the nicotinic receptor by exchanging for one magnesium ion. It is the basis of the ion exchange receptor pair model, in which two acetylcholine ions exchange for one magnesium ion in contact with and between a pair of negatively charged receptor groups about 4 A apart. In the resting state the electrostatic attraction between the negatively charged receptor groups and the Mg2+ ion exerts a binding force. This binding force is opposed by the quantum mechanical repulsions of the electron clouds of the charged groups and ions in contact, together with the mutual repulsion of the pair of receptor oxyanions. When the Mg2+ ion is replaced by two acetylcholine ions the quaternary heads of the latter are positioned so that they form two mutually repelling ACh+ receptor group dipoles. As the Mg2+ ion leaves, its rehydration energy contributes to the sum of the electron cloud repulsions and the ACh+ receptor group dipole repulsions, causing the receptor groups to be forced apart activating the receptor macromolecule. The subsequent decrease in ACh+ concentration results in the reestablishment of the resting state. The coulombic electrostatic energy, the Born repulsion energy, the London attraction energy and the oxyanion ACh+ dipole repulsion energies have been calculated and shown to be consistent with the model. The displacement of the Mg2+ by two ACh+ ions makes several hundred kcals of energy available for receptor group separation and receptor activation.     

A hypothesis concerning the structure of cAMP-and cGMP-dependent protein kinases.

cAMP-and cGMP-dependent protein kinases have been purified. Each enzyme demonstrates high specificity and affinity for the cyclic nucleotide with binding of two moles of nucleotide per holoenzyme and each enzyme is an ATP: phosphotransferase. The holoenzymes have similar molecular weights and demonstrate similar molecular asymmetry. A structural model relating the two enzymes is proposed. cGMP-dependent protein kinase is proposed to be a dimer composed of two identical protomers in isologous association with the chains arranged in anti-parallel fashion. cAMP-dependent protein kinase is proposed to have a similar structure with a dyad axis of symmetry but with a discontinuity in each chain. These structures account for the differing mechanisms of cyclic nucleotide activation of the two enzymes.     

Protein molecular dynamics with the generalized Born/ACE solvent model

Implicit solvent models are increasingly important for the study of proteins in aqueous solution. Here, the generalized Born (GB) solvent polarization model as implemented in the analytical ACE potential [Schaefer and Karplus (1996) J Phys Chem 100:1578] is used to perform molecular dynamics simulations of two small, homologous proteins: the immunoglobulin-binding domain of streptococcal protein G and the Ras binding domain of Raf. Several model parameterizations are compared through more than 60 ns of simulation. Results are compared with two simpler solvent models-an accessible surface area model and a distant-dependent dielectric model, with finite-difference Poisson calculations, with existing explicit solvent simulations, and with experimental data. The simpler models yield stable but distorted structures. The best GB/ACE implementation uses a set of atomic Voronoi volumes reported recently, obtained by averaging over a large database of crystallographic protein structures. A 20% reduction is applied to the volumes, compensating in an average sense for an excessive de-screening of individual charges inherent in the ACE self-energy and for an undersolvation of dipolar groups inherent in the GB screening function. This GB/ACE parameterization yields stable trajectories on the 0.5-1-ns time scale that deviate moderately (approximately 1.5-2.5 A) from the X-ray structure, reproduce approximately the surface distribution of charged, polar, and hydrophobic groups, and reproduce accurately backbone flexibility as measured by amide NMR-order parameters. Over longer time scales (1.5-3 ns), some of the protein G runs escape from the native energy basin and deviate strongly (3 A) from the native structure. The conformations sampled during the transition out of the native energy basin are overstabilized by the GB/ACE solvation model, as compared with a numerical treatment of the full dielectric continuum model.     

Conformational changes in HIV-1 proteinase: effect of protonation of the active center on conformation of HIV-1 proteinase in water

At weak acidic pH, where HIV-1 proteinase is most stable and active, its catalytic Asp 25/25' dyad shares one proton. At a physiological pH the dyad is deprotonated, however, 2 ns molecular dynamics simulations of the HIV-1 protease with monoprotonated and deprotonated Asp25/25' dyad is performed, in order to investigate the influence of Asp25/25' protonation state on the proteinase dynamics. For net charge neutralization the 4 Cl- ions were included. In case of deprotonated active site the significant tertiary structure deviation of HIV-1 PR structure from crystal structure is observed, while in the monoprotonated one the tertiary structure fluctuates near starting structure. Possible mechanism of the influence of the Asp25/25' protonation state on proteinase dynamics is discussed.