A fast and accurate method for predicting pKa of residues in proteins

Accurate Computation of Aqueous pKas of Biologically Relevant Organic Acids: Overcoming the Challenges Posed by Multiple Conformers, Tautomeric Equilibria, and Disparate Functional Groups with the Fully Black-Box pK-Yay Method

The use of static electronic structure calculations to compute solution-phase pKas offers a great advantage in that a macroscopic bulk property could be computed via microscopic computations involving very few molecules. There are various sources of errors in the quantum chemical calculations though. Overcoming these errors to accurately compute pKas of a plethora of acids is an active area of research in physical chemistry pursued by both computational as well as experimental chemists. We recently developed the pK-Yay method in our attempt to accurately compute aqueous pKas of strong and weak acids. The method is fully black-box, computationally inexpensive, and is very easy for even a nonexpert to use. However, the method was thus far tested on very few molecules (only 16 in all). Herein, in order to assess the future applicability of pK-Yay, we study the effect of multiple conformers, the presence of tautomers under equilibrium, and the impact of a wide variety of functional groups (derivatives of acetic acid with substituents at various positions, dicarboxylic acids, aromatic carboxylic acids, amines and amides, phenols and thiols, and fluorine bearing organic acids). Starting with more than 1000 conformers and tautomers, this study establishes that overall errors of ∼ 1.0 pKa units are routinely obtained for a majority of the molecules. Larger errors are noted in cases where multiple charges, intramolecular hydrogen bonding, and several ionizable functional groups are simultaneously present. An important conclusion to emerge from this work is that, the computed pKas are insensitive (difference <0.5) to whether we consider multiple conformers/tautomers or only choose the most stable conformer/tautomer. Further, pK-Yay captures the stereoelectronic effects arising due to differing axial vs equatorial pattern, and is useful to predict the dominant acid-base equilibrium in a system featuring several equilibria. Overall, pK-Yay may be employed in several chemical applications featuring organic molecules and biomonomers.

3D Interaction Homology: Computational Titration of Aspartic Acid, Glutamic Acid and Histidine Can Create pH-Tunable Hydropathic Environment Maps

Aspartic acid, glutamic acid and histidine are ionizable residues occupying various protein environments and perform many different functions in structures. Their roles are tied to their acid/base equilibria, solvent exposure, and backbone conformations. We propose that the number of unique environments for ASP, GLU and HIS is quite limited. We generated maps of these residue's environments using a hydropathic scoring function to record the type and magnitude of interactions for each residue in a 2703-protein structural dataset. These maps are backbone-dependent and suggest the existence of new structural motifs for each residue type. Additionally, we developed an algorithm for tuning these maps to any pH, a potentially useful element for protein design and structure building. Here, we elucidate the complex interplay between secondary structure, relative solvent accessibility, and residue ionization states: the degree of protonation for ionizable residues increases with solvent accessibility, which in turn is notably dependent on backbone structure.

Quantum mechanics/molecular mechanics restrained electrostatic potential fitting

We present a quantum mechanics/molecular mechanics (QM/MM) method to evaluate the partial charges of amino acid residues for use in MM potentials based on their protein environment. For each residue of interest, the nearby residues are included in the QM system while the rest of the protein is treated at the MM level of theory. After a short structural optimization, the partial charges of the central residue are fit to the electrostatic potential using the restrained electrostatic potential (RESP) method. The resulting charges and electrostatic potential account for the individual environment of the residue, although they lack the transferable nature of library partial charges. To evaluate the quality of the QM/MM RESP charges, thermodynamic integration is used to measure the pKa shift of the aspartic acid residues in three different proteins, turkey egg lysozyme, beta-cryptogein, and Thioredoxin. Compared to the AMBER ff99SB library values, the QM/MM RESP charges show better agreement between the calculated and experimental pK(a) values for almost all of the residues considered.

The multiple roles of histidine in protein interactions

Background

Among the 20 natural amino acids histidine is the most active and versatile member that plays the multiple roles in protein interactions, often the key residue in enzyme catalytic reactions. A theoretical and comprehensive study on the structural features and interaction properties of histidine is certainly helpful.

Results

Four interaction types of histidine are quantitatively calculated, including: (1) Cation-π interactions, in which the histidine acts as the aromatic π-motif in neutral form (His), or plays the cation role in protonated form (His⁺); (2) π-π stacking interactions between histidine and other aromatic amino acids; (3) Hydrogen-π interactions between histidine and other aromatic amino acids; (4) Coordinate interactions between histidine and metallic cations. The energies of π-π stacking interactions and hydrogen-π interactions are calculated using CCSD/6-31+G(d,p). The energies of cation-π interactions and coordinate interactions are calculated using B3LYP/6-31+G(d,p) method and adjusted by empirical method for dispersion energy.

Conclusions

The coordinate interactions between histidine and metallic cations are the strongest one acting in broad range, followed by the cation-π, hydrogen-π, and π-π stacking interactions. When the histidine is in neutral form, the cation-π interactions are attractive; when it is protonated (His⁺), the interactions turn to repulsive. The two protonation forms (and pK_a values) of histidine are reversibly switched by the attractive and repulsive cation-π interactions. In proteins the π-π stacking interaction between neutral histidine and aromatic amino acids (Phe, Tyr, Trp) are in the range from -3.0 to -4.0 kcal/mol, significantly larger than the van der Waals energies.

Rapid calculation of protein pKa values using Rosetta

We developed a Rosetta-based Monte Carlo method to calculate the pK(a) values of protein residues that commonly exhibit variable protonation states (Asp, Glu, Lys, His, and Tyr). We tested the technique by calculating pK(a) values for 264 residues from 34 proteins. The standard Rosetta score function, which is independent of any environmental conditions, failed to capture pK(a) shifts. After incorporating a Coulomb electrostatic potential and optimizing the solvation reference energies for pK(a) calculations, we employed a method that allowed side-chain flexibility and achieved a root mean-square deviation (RMSD) of 0.83 from experimental values (0.68 after discounting 11 predictions with an error over 2 pH units). Additional degrees of side-chain conformational freedom for the proximal residues facilitated the capture of charge-charge interactions in a few cases, resulting in an overall RMSD of 0.85 pH units. The addition of backbone flexibility increased the overall RMSD to 0.93 pH units but improved relative pK(a) predictions for proximal catalytic residues. The method also captures large pK(a) shifts of lysine and some glutamate point mutations in staphylococcal nuclease. Thus, a simple and fast method based on the Rosetta score function and limited conformational sampling produces pK(a) values that will be useful when rapid estimation is essential, such as in docking, design, and folding.

Bluues: a program for the analysis of the electrostatic properties of proteins based on generalized Born radii

BACKGROUND

The Poisson-Boltzmann (PB) equation and its linear approximation have been widely used to describe biomolecular electrostatics. Generalized Born (GB) models offer a convenient computational approximation for the more fundamental approach based on the Poisson-Boltzmann equation, and allows estimation of pairwise contributions to electrostatic effects in the molecular context.

RESULTS

We have implemented in a single program most common analyses of the electrostatic properties of proteins. The program first computes generalized Born radii, via a surface integral and then it uses generalized Born radii (using a finite radius test particle) to perform electrostaic analyses. In particular the ouput of the program entails, depending on user's requirement: 1) the generalized Born radius of each atom; 2) the electrostatic solvation free energy; 3) the electrostatic forces on each atom (currently in a developmental stage); 4) the pH-dependent properties (total charge and pH-dependent free energy of folding in the pH range -2 to 18; 5) the pKa of all ionizable groups; 6) the electrostatic potential at the surface of the molecule; 7) the electrostatic potential in a volume surrounding the molecule;

CONCLUSIONS

Although at the expense of limited flexibility the program provides most common analyses with requirement of a single input file in PQR format. The results obtained are comparable to those obtained using state-of-the-art Poisson-Boltzmann solvers. A Linux executable with example input and output files is provided as supplementary material.