Phosphorylation is one of the most important post-translational modifications in living cells and the phosphate groups are present in many of its biomolecules. Despite the fact that the phosphate group is ubiquitous, their role and the environment to which is is exposed can change dramatically. Therefore, in computational methods, the phosphate parameters are probably the among the least transferable, e.g. a phosphate group in phospholipids is significantly different from a phosphoryl tyrosine, in terms of electrostatics. In current forcefields, phosphate groups have been parameterized in the context of lipids and in connection with nucleotides, while parameters for phosphorylated amino acids are still scarce. Most attempts made at these parameters, were based on RESP calculations on the ESP generated from QM calculations [1]. The phosphate group is pH-active, with the second a pKa value around 6, which is problematic to be studied using conventional MD simulations. With constant-pH MD (CpHMD) simulations [2], we can capture the coupling between protonation and conformation for the phosphate group in their different environments.
We have devised several strategies to obtain the GROMOS 54A7 charge parameters for phosphoryl tyrosine (pTyr), serine (pSet), and threonine (pThr), namely: RESP calculations performed on the ESP generated either by Hartree-Fock or DFT calculations; directly adapting the charge set from Wojciechowski et al. [1]; and by manually curating the best charge set obtained taking in consideration the experimental data available. With all charge sets it was possible to run CpHMD simulations on simple pentapeptides to calibrate the pKa value against experimental data. We tested the parameters with simple systems like the Gly-Gly-X-Ala tetrapeptides and with more complex systems, like the phosphorylated dodecapeptide, free or complexed with a phospholipase [3]. We observed that for the simple systems, all charge sets, after calibration, are able to predict their unshifted pKa values. However, all recipe-based charge parameterization processes, failed (sometimes even qualitatively) in predicting the pKa shift for the complex, which could only be corrected by manually curating the charges.
[1] Wojciechowski, M., Grycuk, T., Antosiewicz, J.M., Lesyng, B., Biophys. J., 2003, 84, 750-6.
[2] Machuqueiro, M., Baptista, A. M., J. Phys Chem. B, 2006, 110, 2927-2933.
[3] Singer, A.U., Forman-Kay, J.D., Prot. Sci., 1997, 6, 1910-9.
