Improved model of hydrated calcium ion for molecular dynamics simulations using classical biomolecular force fields

Abstract

Calcium ions (Ca2+) play key roles in various fundamental biological processes such as cell signaling and brain function. Molecular dynamics (MD) simulations have been used to study such interactions, however, the accuracy of the Ca2+ models provided by the standard MD force fields has not been rigorously tested. Here, we assess the performance of the Ca2+ models from the most popular classical force fields AMBER and CHARMM by computing the osmotic pressure of model compounds and the free energy of DNA–DNA interactions. In the simulations performed using the two standard models, Ca2+ ions are seen to form artificial clusters with chloride, acetate, and phosphate species; the osmotic pressure of CaAc2 and CaCl2 solutions is a small fraction of the experimental values for both force fields. Using the standard parameterization of Ca2+ ions in the simulations of Ca2+-mediated DNA–DNA interactions leads to qualitatively wrong outcomes: both AMBER and CHARMM simulations suggest strong inter-DNA attraction whereas, in experiment, DNA molecules repel one another. The artificial attraction of Ca2+ to DNA phosphate is strong enough to affect the direction of the electric field-driven translocation of DNA through a solid-state nanopore. To address these shortcomings of the standard Ca2+ model, we introduce a custom model of a hydrated Ca2+ ion and show that using our model brings the results of the above MD simulations in quantitative agreement with experiment. Our improved model of Ca2+ can be readily applied to MD simulations of various biomolecular systems, including nucleic acids, proteins and lipid bilayer membranes.

Molecular electrometer and binding of cations to phospholipid bilayers

“…(We note that while this is the most likely scenario for the models that overestimated changes in both order parameters, for CaCl 2 it is possible also that the headgroup response is oversensitive to bound cations, see ESI † .) In CHARMM36 with the heptahydrated Ca 2+ by Yoo et al 76 , ∆S β CH was overestimated but ∆S α CH underestimated (Fig. 2), in line with the ∆S β CH /∆S α CH ratio in CHARMM36 being larger than in experiments (Fig.…”

“…As a function of CaCl 2 concentration, all models but one (CHARMM36 with the recent heptahydrated Ca 2+ by Yoo et al 76 ) overestimated the order parameter decrease (Fig. 2), which according to the molecular electrometer indicates too strong Ca 2+ binding.…”

“…4). Further, compensation of missing electronic polarizability by scaling the ion charge 112,114 76 significantly reduced Ca 2+ binding in CHARMM36 (Fig. 6), however, the model must be further analysed to fully interpret the results.…”

“…5). It should also be noted that the low binding affinity of Na + in CHARMM36 is due to the additional repulsion (NBFIX 68 ) added between the sodium ions and lipid oxygens (ESI † ), and that in the Ca 2+ model by Yoo et al 76 the calcium is forced to be solvated solely by water. Altogether, our results indicate that probably both, lipid and ion force field parameters, need improvement to correctly predict the cation binding affinity, and the associated structural changes.…”

Development of Nonbonded Models for Metal Cations Using the Electronic Continuum Correction

“…These data are not in agreement with the experimental one of 8 provided by Marcus. However, the exact number of water molecules coordinating Ca 2+ is not known in the literature, [53] and theoretical simulations suggest a CN number in the range of 6-8. [54,55] Changes in the CN are detected also for Sr 2+ , for which an average CN of 7.3 has been calculated: although not in agreement with the CN of 8 suggested by Marcus, this value is within the range of 6-8 defined by several computational results provided by the literature.…”