Donald Bashford

New protein parameters are reported for the all-atom empirical energy function in the CHARMM program. The parameter evaluation was based on a self-consistent approach designed to achieve a balance between the internal (bonding) and interaction (nonbonding) terms of the force field and among the solvent-solvent, solvent-solute, and solute-solute interactions. Optimization of the internal parameters used experimental gas-phase geometries, vibrational spectra, and torsional energy surfaces supplemented with ab initio results. The peptide backbone bonding parameters were optimized with respect to data for N-methylacetamide and the alanine dipeptide. The interaction parameters, particularly the atomic charges, were determined by fitting ab initio interaction energies and geometries of complexes between water and model compounds that represented the backbone and the various side chains. In addition, dipole moments, experimental heats and free energies of vaporization, solvation and sublimation, molecular volumes, and crystal pressures and structures were used in the optimization. The resulting protein parameters were tested by applying them to noncyclic tripeptide crystals, cyclic peptide crystals, and the proteins crambin, bovine pancreatic trypsin inhibitor, and carbonmonoxy myoglobin in vacuo and in crystals. A detailed analysis of the relationship between the alanine dipeptide potential energy surface and calculated protein φ, χ angles was made and used in optimizing the peptide group torsional parameters. The results demonstrate that use of ab initio structural and energetic data by themselves are not sufficient to obtain an adequate backbone representation for peptides and proteins in solution and in crystals. Extensive comparisons between molecular dynamics simulations and experimental data for polypeptides and proteins were performed for both structural and dynamic properties. Energy minimization and dynamics simulations for crystals demonstrate that the latter are needed to obtain meaningful comparisons with experimental crystal structures. The presented parameters, in combination with the previously published CHARMM all-atom parameters for nucleic acids and lipids, provide a consistent set for condensed-phase simulations of a wide variety of molecules of biological interest.

Implicit solvation models provide, for many applications, a reasonably accurate and computationally effective way to describe the electrostatics of aqueous solvation. Here, a popular analytical Generalized Born (GB) solvation model is modified to improve its accuracy in calculating the solvent polarization part of free energy changes in large-scale conformational transitions, such as protein folding. In contrast to an earlier GB model (implemented in the AMBER-6 program), the improved version does not overstabilize the native structures relative to the finite-difference Poisson-Boltzmann continuum treatment. In addition to improving the energy balance between folded and unfolded conformers, the algorithm (available in the AMBER-7 and NAB molecular modeling packages) is shown to perform well in more than 50 ns of native-state molecular dynamics (MD) simulations of thioredoxin, protein-A, and ubiquitin, as well as in a simulation of Barnase/Barstar complex formation. For thioredoxin, various combinations of input parameters have been explored, such as the underlying gas-phase force fields and the atomic radii. The best performance is achieved with a previously proposed modification to the torsional potential in the Amber ff99 force field, which yields stable native trajectories for all of the tested proteins, with backbone root-mean-square deviations from the native structures being approximately 1.5 A after 6 ns of simulation time. The structure of Barnase/Barstar complex is regenerated, starting from an unbound state, to within 1.9 A relative to the crystal structure of the complex.

It would often be useful in computer simulations to use a simple description of solvation effects, instead of explicitly representing the individual solvent molecules. Continuum dielectric models often work well in describing the thermodynamic aspects of aqueous solvation, and approximations to such models that avoid the need to solve the Poisson equation are attractive because of their computational efficiency. Here we give an overview of one such approximation, the generalized Born model, which is simple and fast enough to be used for molecular dynamics simulations of proteins and nucleic acids. We discuss its strengths and weaknesses, both for its fidelity to the underlying continuum model and for its ability to replace explicit consideration of solvent molecules in macromolecular simulations. We focus particularly on versions of the generalized Born model that have a pair-wise analytical form, and therefore fit most naturally into conventional molecular mechanics calculations.

A macroscopic electrostatic model is used to calculate the pKa values of the titratable groups in lysozyme. The model makes use of detailed structural information and treats solvation self-energies and interactions arising from permanent partial charges and titratable charges. Both the tetragonal and triclinic crystal structures are analyzed. Half of the experimentally observed pKa shifts (11 out of 21) are well reproduced by calculations for both structures; this includes the unusually high pKa of Glu 35 in the active site. For more than half the pKa's (13 out of 21), there is a large difference (1-3.3 pK units) between the results from the two structures. Many of these correspond to the titrating groups for which the calculations are in error. Since for an ionic strength of 0.1 M the Debye screening between titratable groups leads to a very high effective dielectric constant (the average value for all pairs of titrating groups is approximately 900), near-neighbor interactions dominate the pKa perturbations. Thus, the pKa values are very sensitive to the details of the local protein conformation, and it is likely that side-chain mobility has an important role in determining the observed pKa shifts.

The analytic generalized Born approximation is an efficient electrostatic model that describes molecules in solution. Here it is modified to permit a more accurate description of large macromolecules, while its established performance on small compounds is nearly unaffected. The modified model is also adapted to describe molecules with an interior dielectric constant not equal to unity. The model is tested by computations of pK shifts for a number of titratable residues in lysozyme, myoglobin, and bacteriorhodopsin. In general, except for some deeply buried residues of bacteriorhodopsin, the results show reasonable agreement with both experimental data and calculations based on numerical solution of the Poisson−Boltzmann equation. A very close agreement between the two models is obtained in an application to the prediction of the pK shifts associated with conformational change. The calculations based on this version of the generalized Born approximation are much faster than finite difference solutions of the Poisson−Boltzmann equation, which makes the present method useful for a variety of other applications where computational time is a critical factor. The model may also be integrated into molecular dynamics programs to replace explicit solvent simulations which are particularly time-consuming for large molecules.