António M. Baptista

Constant-pH molecular dynamics using stochastic titration

A new method is proposed for performing constant-pH molecular dynamics (MD) simulations, that is, MD simulations where pH is one of the external thermodynamic parameters, like the temperature or the pressure. The protonation state of each titrable site in the solute is allowed to change during a molecular mechanics (MM) MD simulation, the new states being obtained from a combination of continuum electrostatics (CE) calculations and Monte Carlo (MC) simulation of protonation equilibrium. The coupling between the MM/MD and CE/MC algorithms is done in a way that ensures a proper Markov chain, sampling from the intended semigrand canonical distribution. This stochastic titration method is applied to succinic acid, aimed at illustrating the method and examining the choice of its adjustable parameters. The complete titration of succinic acid, using constant-pH MD simulations at different pH values, gives a clear picture of the coupling between the trans/gauche isomerization and the protonation process, making i...

Progress in the prediction of pKa values in proteins

The pKa‐cooperative aims to provide a forum for experimental and theoretical researchers interested in protein pKa values and protein electrostatics in general. The first round of the pKa‐cooperative, which challenged computational labs to carry out blind predictions against pKas experimentally determined in the laboratory of Bertrand Garcia‐Moreno, was completed and results discussed at the Telluride meeting (July 6–10, 2009). This article serves as an introduction to the reports submitted by the blind prediction participants that will be published in a special issue of PROTEINS: Structure, Function and Bioinformatics. Here, we briefly outline existing approaches for pKa calculations, emphasizing methods that were used by the participants in calculating the blind pKa values in the first round of the cooperative. We then point out some of the difficulties encountered by the participating groups in making their blind predictions, and finally try to provide some insights for future developments aimed at improving the accuracy of pKa calculations. Proteins 2011;

Simulation of protein conformational freedom as a function of pH: constant-pH molecular dynamics using implicit titration

Solution pH is a determinant parameter on protein function and stability, and its inclusion in molecular dynamics simulations is attractive for studies at the molecular level. Current molecular dynamics simulations can consider pH only in a very limited way, through a somewhat arbitrary choice of a set of fixed charges on the titrable sites. Conversely, continuum electrostatic methods that explicitly treat pH effects assume a single protein conformation whose choice is not clearly defined. In this paper we describe a general method that combines both titration and conformational freedom. The method is based on a potential of mean force for implicit titration and combines both usual molecular dynamics and pH‐dependent calculations based on continuum methods. A simple implementation of the method, using a mean field approximation, is presented and applied to the bovine pancreatic trypsin inhibitor. We believe that this constant‐pH molecular dynamics method, by correctly sampling both charges and conformation, can become a valuable help in the understanding of the dependence of protein function and stability on pH.

Simulation of Electron-Proton Coupling with a Monte Carlo Method: Application to Cytochrome c3 Using Continuum Electrostatics

A new method is presented for simulating the simultaneous binding equilibrium of electrons and protons on protein molecules, which makes it possible to study the full equilibrium thermodynamics of redox and protonation processes, including electron-proton coupling. The simulations using this method reflect directly the pH and electrostatic potential of the environment, thus providing a much closer and realistic connection with experimental parameters than do usual methods. By ignoring the full binding equilibrium, calculations usually overlook the twofold effect that binding fluctuations have on the behavior of redox proteins: first, they affect the energy of the system by creating partially occupied sites; second, they affect its entropy by introducing an additional empty/occupied site disorder (here named occupational entropy). The proposed method is applied to cytochrome c3 of Desulfovibrio vulgaris Hildenborough to study its redox properties and electron-proton coupling (redox-Bohr effect), using a continuum electrostatic method based on the linear Poisson-Boltzmann equation. Unlike previous studies using other methods, the full reduction order of the four hemes at physiological pH is successfully predicted. The sites more strongly involved in the redox-Bohr effect are identified by analysis of their titration curves/surfaces and the shifts of their midpoint redox potentials and pKa values. Site-site couplings are analyzed using statistical correlations, a method much more realistic than the usual analysis based on direct interactions. The site found to be more strongly involved in the redox-Bohr effect is propionate D of heme I, in agreement with previous studies; other likely candidates are His67, the N-terminus, and propionate D of heme IV. Even though the present study is limited to equilibrium conditions, the possible role of binding fluctuations in the concerted transfer of protons and electrons under nonequilibrium conditions is also discussed. The occupational entropy contributions to midpoint redox potentials and pKa values are computed and shown to be significant.

Improved modeling of side-chains in proteins with rotamer-based methods: A flexible rotamer model

Side‐chain modeling has a widespread application in many current methods for protein tertiary structure determination, prediction, and design. Of the existing side‐chain modeling methods, rotamer‐based methods are the fastest and most efficient. Classically, a rotamer is conceived as a single, rigid conformation of an amino acid side‐chain. Here, we present a flexible rotamer model in which a rotamer is a continuous ensemble of conformations that cluster around the classic rigid rotamer. We have developed a thermodynamically based method for calculating effective energies for the flexible rotamer. These energies have a one‐to‐one correspondence with the potential energies of the rigid rotamer. Therefore, the flexible rotamer model is completely general and may be used with any rotamer‐based method in substitution of the rigid rotamer model. We have compared the performance of the flexible and rigid rotamer models with one side‐chain modeling method in particular (the self‐consistent mean field theory method) on a set of 20 high quality crystallographic protein structures. For the flexible rotamer model, we obtained average predictions of 85.8% for χ 1, 76.5% for χ 1+2 and 1.34 Å for root‐mean‐square deviation (RMSD); the corresponding values for core residues were 93.0%, 87.7% and 0.70 Å, respectively. These values represent improvements of 7.3% for χ 1, 8.1% for χ 1+2 and 0.23 Å for RMSD over the predictions obtained with the rigid rotamer model under otherwise identical conditions; the corresponding improvements for core residues were 6.9%, 10.5% and 0.43 Å, respectively. We found that the predictions obtained with the flexible rotamer model were also significantly better than those obtained for the same set of proteins with another state‐of‐the‐art side‐chain placement method in the literature, especially for core residues. The flexible rotamer model represents a considerable improvement over the classic rigid rotamer model. It can, therefore, be used with considerable advantage in all rotamer‐based methods commonly applied to protein tertiary structure determination, prediction, and design and also in predictions of free energies in mutational studies. Proteins 1999;37:530–543. ©1999 Wiley‐Liss, Inc.

Constant-pH molecular dynamics simulations reveal a β-rich form of the human prion protein.

The misfolding of the prion protein (PrP) into a pathogenic β-rich form (PrP(Sc)) has been suggested to occur in the endocytic pathway, triggered by low pH. In this work we performed several constant-pH molecular dynamics simulations of human PrP 90-231 in the pH range 2-7, totaling more than 2 μs. We observed a strong conformational pH dependence where on average the helix content decreased and the β content increased toward acidic pH. Unlike some proposed models, the flexible N-terminus region did not gain stable structure at low pH. Rather, the main structural changes occurred on the helix-rich C-terminus core, as proposed in other models, namely, in the regions around 135-155 and 185-200. The protonation of His187 is found to be associated with a loss of interaction between two PrP subdomains, potentially playing a major role in the misfolding process. In one of the simulations at pH 2, a stable β-rich structure was formed that may be an intermediate of PrP(Sc) formation, indicating that misfolding may precede dimerization.

Acidic range titration of HEWL using a constant‐pH molecular dynamics method

In this work, we present the first application to a protein of the stochastic constant‐pH molecular dynamics (MD) method with the inclusion of proton tautomerism. The acidic titration of HEWL was performed under different conditions. Both generalized reaction field (GRF) and particle mesh Ewald (PME) methods were used in the treatment of the long range electrostatics and, even though the PME simulations revealed to be more stable, the better results were obtained using GRF (pKa RMSD of 0.82 for GRF and 1.13 for PME). The results using PME at different dielectric constants (2, 4, and 8) also revealed that there was no significant improvement in pKas prediction upon increasing the dielectric constant. The secondary structure analysis of HEWL revealed a remarkably stable protein in the acidic pH range. The β‐sheet strands (unlike the α‐helices) seem to be destabilized upon pH decrease, suggesting that the β‐domain is less stable than the α‐domain. The four principal α‐helices were also ordered according to their stability in the acidic pH range and the results (4 < 1 < 2 ≈ 3) were consistent with the ones obtained in thermal denaturation studies. Proteins 2008.

Lipases and esterases: a review of their sequences, structure and evolution

This chapter aims to provide a brief review on the enzyme family of lipases and esterases. The sequences, 3D structures and pH dependent electrostatic signatures are presented and analyzed. Since the family comprises more than 100 sequences, we have tried to focus on the most interesting features from our perspective, which translates into finding similarities and differences between members of this family, in particular in and around the active sites, and to identify residues that are partially or totally conserved. Such residues we believe are either important for maintaining the structural scaf-fold of the protein or to maintain activity or specificity. The structure function relationship for these proteins is therefore of central interest. Can we uniquely identify a protein from this large family of sequences--and if so, what is the identifier? The protein family displays some highly complex features: many of the proteins are interfacially activated, i.e. they need to be in physical contact with the aggregated substrate. Access to the active site is blocked with either a loop fragment or an alpha-helical fragment in the absence of interfacial contact. Although the number of known, relevant protein 3D structures is growing steadily, we are nevertheless faced with a virtual explosion in the number of known or deduced amino acid sequences. It is therefore unrealistic to expect that all protein sequences within the foreseeable future will have their 3D structure determined by X-ray diffractional analysis or through other methods. When feasible the gene and/or the amino acid sequences will be analyzed from an evolutionary perspective. As the 3D folds are often remarkably similar, both among the triglyceride lipases as well as among the esterases, the functional diversities (e.g. specificity) must originate in differences in surface residue utilization, in particular of charged residues. The pH variations in the isopotential surfaces of some of the most interesting lipases are presented and a qualitative interpretation proposed. Finally we illustrate that NMR has potential for becoming an important tool in the study of lipases, esterases and their kinetics.

Conformational changes induced by ATP‐hydrolysis in an ABC transporter: A molecular dynamics study of the Sav1866 exporter

ATP‐Binding Cassette (ABC) transporters are ubiquitous membrane proteins that use energy from ATP binding or/and hydrolysis to actively transport allocrites across membranes. In this study, we identify ATP‐hydrolysis induced conformational changes in a complete ABC exporter (Sav1866) from Staphylococcus aureaus, using molecular dynamics (MD) simulations. By performing MD simulations on the ATP and ADP+IP bound states, we identify the conformational consequences of hydrolysis, showing that the major rearrangements are not restricted to the NBDs, but extend to the transmembrane domains (TMDs) external regions. For the first time, to our knowledge, we see, within the context of a complete transporter, NBD dimer opening in the ADP+IP state in contrast with all ATP‐bound states. This opening results from the dissociation of the ABC signature motif from the nucleotide. In addition, in both states, we observe the opening of a gate entrance in the intracellular loop region leading to the exposure of the TMDs internal cavity to the cytoplasm. To see if this opening was large enough to allow allocrite transport, the adiabatic energy profile for doxorubicin passage was determined. For both states, this profile, although an approximation, is overall downhill from the cytoplasmatic to the extracellular side, and the local energy barriers along the TMDs are relatively small, evidencing the exporter nature of Sav1866. The major difference between states is an energy barrier located in the cytoplasmic gate region, which becomes reduced upon hydrolysis, suggesting that allocrite passage is facilitated, and evidencing a possible molecular mechanism for the active transport in these proteins. Proteins 2011;

Proton pathways in a [NiFe]‐hydrogenase: A theoretical study

We present here a theoretical study to investigate possible proton pathways in the [NiFe]‐hydrogenase from Desulfovibrio gigas. The approach used in this study consists of a combination of Poisson–Boltzmann and Monte Carlo simulations together with a distance‐based network analysis to find possible groups involved in the proton transfer. Results obtained at different pH values show a reasonable number of proton active residues distributed by the protein interior and surface, with a concentration around the metal centres. The electrostatic interactions in this protein are strong, as shown by the unusual shape of the titration curves of several sites. Some residue pairs show strongly correlated protonations, indicating the sharing and probably exchange of a proton between them. The conjugation of the PB and MC simulations with the distance‐based analysis allows a detailed characterization of the possible proton pathways. We discuss previous suggestions and propose a new complete pathway for the proton transfer between the active site and the surface. This pathway is mainly composed of histidines and glutamic acid residues. Proteins 2008.