Insights and Challenges in Correcting Force Field Based Solvation Free Energies Using a Neural Network Potential

Abstract

We present a comprehensive study investigating the potential gain in accuracy for calculating absolute solvation free energies (ASFE) using a neural network potential to describe the intramolecular energy of the solute. We calculated the ASFE for most compounds from the FreeSolv database using the Open Force Field (OpenFF) and compared them to earlier results obtained with the CHARMM General Force Field (CGenFF). By applying a nonequilibrium (NEQ) switching approach between the molecular mechanics (MM) description (either OpenFF or CGenFF) and the neural net potential (NNP)/MM level of theory (using ANI-2x as the NNP potential), we attempted to improve the accuracy of the calculated ASFEs. The predictive performance of the results did not change when this approach was applied to all 589 small molecules in the FreeSolv database that ANI-2x can describe. When selecting a subset of 156 molecules, focusing on compounds where the force fields performed poorly, we saw a slight improvement in the root-mean-square error (RMSE) and mean absolute error (MAE). The majority of our calculations utilized unidirectional NEQ protocols based on Jarzynski's equation. Additionally, we conducted bidirectional NEQ switching for a subset of 156 solutes. Notably, only a small fraction (10 out of 156) exhibited statistically significant discrepancies between unidirectional and bidirectional NEQ switching free energy estimates.

Figure 1

Free energy calculations between high and low levels of theory can be used to correct alchemical free energy calculations performed at a low level of theory. Left: The thermodynamic cycle used to compute an ASFE ΔG_low^solv at the MM level of theory using annihilation of the solute’s nonbonded interactions. Right: The indirect free energy cycle to correct ΔG_low^solv (in this work using either the CGenFF or the OpenFF force field) is from the MM to the NNP/MM level of theory.

Figure 2

Absolute solvation free energy calculations using the OpenFF force field. ΔG_calc are the calculated values using protocol EXS, while experimental values (ΔG_exp) are taken from the literature. The dark and light gray areas depict the ±1 and ±2 kcal/mol confidence interval.

Figure 3

The KDE of δΔG for the NNP corrected free energy estimate for the 589 compounds of the FreeSolv database. Overlay of the KDE for δ ΔG (eq 5) of the MM (blue line, ΔG_OpenFF) and the NNP/MM (green line, ΔG_OpenFF/ANI2x) free energy estimates. Additionally, the plot shows dots indicating the individual deviations from the experimental values for the MM ASFE results (blue dots) and the NNP/MM corrected ASFE results using the unidirectional correction according to protocol EXS (green dots). The gray line displays the ideal behavior, δΔG = 0 with the two dotted gray lines indicating a deviation of ±1 kcal/mol.

Figure 4

The KDE of δ ΔG (eq 5) for the 156 compound subset. Top: KDE of the errors obtained with OpenFF 2.0 (blue line) and for the NNP corrected results (green line). Bottom: KDE of the errors for CGenFF (blue) and NNP corrected results (green). Results for the individual molecules are shown as dots: blue for the MM, and green for the NNP corrected values.

Figure 5

Panel A: Box plot of the unidirectional NNP correction (ΔG_MM→NNP/MM^corr) as a function of the number of rotatable bonds for the 156 compound subset using protocol UVIE (an analogous plot for EXS is shown in Figure S2 in the Supporting Information). Molecules that are outliers in terms of ΔG_MM→NNP/MM^corr are indicated as diamonds and labeled by numbers, starting with 1 for the compound having the highest correction value. Compounds that are also outliers when using protocol EXS are colored in purple. The bars at the top of the plot indicate how many molecules have this number of rotatable bonds. Panel B: Molecular structures of the outliers 1–12. Molecules in the purple box are also outliers when using protocol EXS. Panel C: Density of the indicated dihedral angle of compounds 4 (left side) and 5 (right side), respectively, in the gas phase and in aqueous solution when using CGenFF (blue) and ANI-2x (green).