Ray Luo

The Amber biomolecular simulation programs.

We describe the development, current features, and some directions for future development of the Amber package of computer programs. This package evolved from a program that was constructed in the late 1970s to do Assisted Model Building with Energy Refinement, and now contains a group of programs embodying a number of powerful tools of modern computational chemistry, focused on molecular dynamics and free energy calculations of proteins, nucleic acids, and carbohydrates.

A Point-Charge Force Field for Molecular Mechanics Simulations of Proteins Based on Condensed-Phase Quantum Mechanical Calculations

Molecular mechanics models have been applied extensively to study the dynamics of proteins and nucleic acids. Here we report the development of a third‐generation point‐charge all‐atom force field for proteins. Following the earlier approach of Cornell et al., the charge set was obtained by fitting to the electrostatic potentials of dipeptides calculated using B3LYP/cc‐pVTZ//HF/6‐31G** quantum mechanical methods. The main‐chain torsion parameters were obtained by fitting to the energy profiles of Ace‐Ala‐Nme and Ace‐Gly‐Nme di‐peptides calculated using MP2/cc‐pVTZ//HF/6‐31G** quantum mechanical methods. All other parameters were taken from the existing AMBER data base. The major departure from previous force fields is that all quantum mechanical calculations were done in the condensed phase with continuum solvent models and an effective dielectric constant of ε = 4. We anticipate that this force field parameter set will address certain critical short comings of previous force fields in condensed‐phase simulations of proteins. Initial tests on peptides demonstrated a high‐degree of similarity between the calculated and the statistically measured Ramanchandran maps for both Ace‐Gly‐Nme and Ace‐Ala‐Nme di‐peptides. Some highlights of our results include (1) well‐preserved balance between the extended and helical region distributions, and (2) favorable type‐II poly‐proline helical region in agreement with recent experiments. Backward compatibility between the new and Cornell et al. charge sets, as judged by overall agreement between dipole moments, allows a smooth transition to the new force field in the area of ligand‐binding calculations. Test simulations on a large set of proteins are also discussed.

Accelerated Poisson–Boltzmann calculations for static and dynamic systems

We report here an efficient implementation of the finite difference Poisson–Boltzmann solvent model based on the Modified Incomplete Cholsky Conjugate Gradient algorithm, which gives rather impressive performance for both static and dynamic systems. This is achieved by implementing the algorithm with Eisenstats two optimizations, utilizing the electrostatic update in simulations, and applying prudent approximations, including: relaxing the convergence criterion, not updating Poisson–Boltzmann‐related forces every step, and using electrostatic focusing. It is also possible to markedly accelerate the supporting routines that are used to set up the calculations and to obtain energies and forces. The resulting finite difference Poisson–Boltzmann method delivers efficiency comparable to the distance‐dependent dielectric model for a system tested, HIV Protease, making it a strong candidate for solution‐phase molecular dynamics simulations. Further, the finite difference method includes all intrasolute electrostatic interactions, whereas the distance dependent dielectric calculations use a 15‐Å cutoff. The speed of our numerical finite difference method is comparable to that of the pair‐wise Generalized Born approximation to the Poisson–Boltzmann method.

A Poisson–Boltzmann dynamics method with nonperiodic boundary condition

We have developed a well-behaved and efficient finite difference Poisson–Boltzmann dynamics method with a nonperiodic boundary condition. This is made possible, in part, by a rather fine grid spacing used for the finite difference treatment of the reaction field interaction. The stability is also made possible by a new dielectric model that is smooth both over time and over space, an important issue in the application of implicit solvents. In addition, the electrostatic focusing technique facilitates the use of an accurate yet efficient nonperiodic boundary condition: boundary grid potentials computed by the sum of potentials from individual grid charges. Finally, the particle–particle particle–mesh technique is adopted in the computation of the Coulombic interaction to balance accuracy and efficiency in simulations of large biomolecules. Preliminary testing shows that the nonperiodic Poisson–Boltzmann dynamics method is numerically stable in trajectories at least 4 ns long. The new model is also fairly e...

Comparison of generalized born and poisson models: Energetics and dynamics of HIV protease

This study characterizes the accuracy of energies and forces computed with a generalized Born (GB) model and the distance‐dependent dielectric (DDD) model with respect to detailed finite solutions of the Poisson equation (FDPE). Tests are done for a small molecule in solution and for HIV‐1 protease with inhibitor, KNI‐272. GB agrees well with FDPE for the small molecule, but less well for the protein system. The correlation between GB and FDPE energies is poorest in calculations of changes upon binding. Also, forces computed with the GB model are less accurate than energies. The DDD model is far less accurate than GB. Nanosecond stochastic dynamics simulations of HIV‐1 protease with an empty active site are used to examine the consequence of the models for the conformational preferences of the active site. Interestingly, the active site flaps remain near their starting conformations in the FDPE and GB simulations but collapse into the active site in the DDD simulation.

Virtual screening using molecular simulations

Effective virtual screening relies on our ability to make accurate prediction of protein‐ligand binding, which remains a great challenge. In this work, utilizing the molecular‐mechanics Poisson‐Boltzmann (or Generalized Born) surface area approach, we have evaluated the binding affinity of a set of 156 ligands to seven families of proteins, trypsin β, thrombin α, cyclin‐dependent kinase (CDK), cAMP‐dependent kinase (PKA), urokinase‐type plasminogen activator, β‐glucosidase A, and coagulation factor Xa. The effect of protein dielectric constant in the implicit‐solvent model on the binding free energy calculation is shown to be important. The statistical correlations between the binding energy calculated from the implicit‐solvent approach and experimental free energy are in the range of 0.56–0.79 across all the families. This performance is better than that of typical docking programs especially given that the latter is directly trained using known binding data whereas the molecular mechanics is based on general physical parameters. Estimation of entropic contribution remains the barrier to accurate free energy calculation. We show that the traditional rigid rotor harmonic oscillator approximation is unable to improve the binding free energy prediction. Inclusion of conformational restriction seems to be promising but requires further investigation. On the other hand, our preliminary study suggests that implicit‐solvent based alchemical perturbation, which offers explicit sampling of configuration entropy, can be a viable approach to significantly improve the prediction of binding free energy. Overall, the molecular mechanics approach has the potential for medium to high‐throughput computational drug discovery. Proteins 2011;

Development of polarizable models for molecular mechanical calculations I: parameterization of atomic polarizability.

In this work, four types of polarizable models have been developed for calculating interactions between atomic charges and induced point dipoles. These include the Applequist, Thole linear, Thole exponential model, and the Thole Tinker-like. The polarizability models have been optimized to reproduce the experimental static molecular polarizabilities obtained from the molecular refraction measurements on a set of 420 molecules reported by Bosque and Sales. We grouped the models into five sets depending on the interaction types, that is, whether the interactions of two atoms that form the bond, bond angle, and dihedral angle are turned off or scaled down. When 1-2 (bonded) and 1-3 (separated by two bonds) interactions are turned off, 1-4 (separated by three bonds) interactions are scaled down, or both, all models including the Applequist model achieved similar performance: the average percentage error (APE) ranges from 1.15 to 1.23%, and the average unsigned error (AUE) ranges from 0.143 to 0.158 Å(3). When the short-range 1-2, 1-3, and full 1-4 terms are taken into account (set D models), the APE ranges from 1.30 to 1.58% for the three Thole models, whereas the Applequist model (DA) has a significantly larger APE (3.82%). The AUE ranges from 0.166 to 0.196 Å(3) for the three Thole models, compared with 0.446 Å(3) for the Applequist model. Further assessment using the 70-molecule van Duijnen and Swart data set clearly showed that the developed models are both accurate and highly transferable and are in fact have smaller errors than the models developed using this particular data set (set E models). The fact that A, B, and C model sets are notably more accurate than both D and E model sets strongly suggests that the inclusion of 1-2 and 1-3 interactions reduces the transferability and accuracy.

Development of Polarizable Models for Molecular Mechanical Calculations II: Induced Dipole Models Significantly Improve Accuracy of Intermolecular Interaction Energies

In the companion paper, we presented a set of induced dipole interaction models using four types of screening functions, which include the Applequist (no screening), the Thole linear, the Thole exponential model, and the Thole Tinker-like (another form of exponential screening function) functions. In this work, we evaluate the performance of polarizability models using a large set of amino acid analog pairs in conformations that are frequently observed in protein structures as a benchmark. For each amino acid pair, we calculated quantum mechanical interaction energies at the MP2/aug-cc-pVTZ//MP2/6-311++G(d,p) level with the basis set superposition error (BSSE) correction and compared them with molecular mechanics results. Encouragingly, all polarizable models significantly outperform the additive F94 and F03 models (mimicking AMBER ff94/ff99 and ff03 force fields, respectively) in reproducing the BSSE-corrected quantum mechanical interaction energies. In particular, the root-mean-square errors (RMSEs) for three Thole models in Set A (where the 1-2 and 1-3 interactions are turned off and all 1-4 interactions are included) are 1.456, 1.417, and 1.406 kcal/mol for model AL (Thole Linear), model AE (Thole exponential), and model AT (Thole Tinker-like), respectively. In contrast, the RMSEs are 3.729 and 3.433 kcal/mol for F94 and F03 models, respectively. A similar trend was observed for the average unsigned errors (AUEs), which are 1.057, 1.025, 1.011, 2.219, and 2.070 kcal/mol for AL, AE, AT, F94/ff99, and F03, respectively. Analyses based on the trend line slopes indicate that the two fixed charge models substantially underestimate the relative strengths of noncharge-charge interactions by 24 (F03) and 35% (F94), respectively, whereas the four polarizable models overestimate the relative strengths by 5 (AT), 3 (AL, AE), and 13% (AA), respectively. Agreement was further improved by adjusting the van der Waals parameters. Judging from the notably improved accuracy in comparison with the fixed charge models, the polarizable models are expected to form the foundation for the development of high quality polarizable force fields for protein and nucleic acid simulations.

The Physical Basis of Nucleic Acid Base Stacking in Water

It has been argued that the stacking of adenyl groups in water must be driven primarily by electrostatic interactions, based upon NMR data showing stacking for two adenyl groups joined by a 3-atom linker but not for two naphthyl groups joined by the same linker. In contrast, theoretical work has suggested that adenine stacking is driven primarily by nonelectrostatic forces, and that electrostatic interactions actually produce a net repulsion between adenines stacking in water. The present study provides evidence that the experimental data for the 3-atom-linked bis-adenyl and bis-naphthyl compounds are consistent with the theory indicating that nonelectrostatic interactions drive adenine stacking. First, a theoretical conformational analysis is found to reproduce the observed ranking of the stacking tendencies of the compounds studied experimentally. A geometric analysis identifies two possible reasons, other than stronger electrostatic interactions, why the 3-atom-linked bis-adenyl compounds should stack more than the bis-naphthyl compounds. First, stacked naphthyl groups tend to lie further apart than stacked adenyl groups, based upon both quantum calculations and crystal structures. This may prevent the bis-naphthyl compound from stacking as extensively as the bis-adenyl compound. Second, geometric analysis shows that more stacked conformations are sterically accessible to the bis-adenyl compound than to the bis-naphthyl compound because the linker is attached to the sides of the adenyl groups, but to the ends of the naphthyl groups. Finally, ab initio quantum mechanics calculations and energy decompositions for relevant conformations of adenine and naphthalene dimers support the view that stacking in these compounds is driven primarily by nonelectrostatic interactions. The present analysis illustrates the importance of considering all aspects of a molecular system when interpreting experimental data, and the value of computer models as an adjunct to chemical intuition.

Assessment of Linear Finite-Difference Poisson-Boltzmann Solvers

CPU time and memory usage are two vital issues that any numerical solvers for the Poisson–Boltzmann equation have to face in biomolecular applications. In this study, we systematically analyzed the CPU time and memory usage of five commonly used finite‐difference solvers with a large and diversified set of biomolecular structures. Our comparative analysis shows that modified incomplete Cholesky conjugate gradient and geometric multigrid are the most efficient in the diversified test set. For the two efficient solvers, our test shows that their CPU times increase approximately linearly with the numbers of grids. Their CPU times also increase almost linearly with the negative logarithm of the convergence criterion at very similar rate. Our comparison further shows that geometric multigrid performs better in the large set of tested biomolecules. However, modified incomplete Cholesky conjugate gradient is superior to geometric multigrid in molecular dynamics simulations of tested molecules. We also investigated other significant components in numerical solutions of the Poisson–Boltzmann equation. It turns out that the time‐limiting step is the free boundary condition setup for the linear systems for the selected proteins if the electrostatic focusing is not used. Thus, development of future numerical solvers for the Poisson–Boltzmann equation should balance all aspects of the numerical procedures in realistic biomolecular applications.