Source-linked AI summary
Non-covalent interactions across organic and biological subsets of chemical space: Physics-based potentials parametrized from machine learning
Tristan Bereau, Robert A. DiStasio, Alexandre Tkatchenko, O. Anatole von Lilienfeld
TL;DR
Classical intermolecular potentials are difficult to parametrize consistently for new molecules, limiting their use across chemical space. IPML combines physics-based interaction models with machine-learned local atomic properties and achieves accurate results across diverse molecular dimers while extending to denser systems.
Problem
Classical potentials are limited to narrow molecular sets because parametrization is tedious, non-systematic, time consuming, and difficult to automate.
Method
IPML uses machine learning to predict atom-in-molecule properties that parameterize physics-based electrostatics, polarization, repulsion, and dispersion models.
Results
0.4 kcal/mol MAE is obtained for 528 angularly displaced molecular dimers, while the benzene crystal cohesive energy is −14.3 kcal/mol at equilibrium density.
Takeaways & Limitations
IPML provides transferable intermolecular potentials for small molecules using local properties predicted by machine learning and global parameters optimized across chemical space.
Abstract
from arXiv · showhide
Classical intermolecular potentials typically require an extensive parametrization procedure for any new compound considered. To do away with prior parametrization, we propose a combination of physics-based potentials with machine learning (ML), coined IPML, which is transferable across small neutral organic and biologically-relevant molecules. ML models provide on-the-fly predictions for environment-dependent local atomic properties: electrostatic multipole coefficients (significant error reduction compared to previously reported), the population and decay rate of valence atomic densities, and polarizabilities across conformations and chemical compositions of H, C, N, and O atoms. These parameters enable accurate calculations of intermolecular contributions---electrostatics, charge penetration, repulsion, induction/polarization, and many-body dispersion. Unlike other potentials, this model is transferable in its ability to handle new molecules and conformations without explicit prior parametrization: All local atomic properties are predicted from ML, leaving only eight global parameters---optimized once and for all across compounds. We validate IPML on various gas-phase dimers at and away from equilibrium separation, where we obtain mean absolute errors between 0.4 and 0.7 kcal/mol for several chemically and conformationally diverse datasets representative of non-covalent interactions in biologically-relevant molecules. We further focus on hydrogen-bonded complexes---essential but challenging due to their directional nature---where datasets of DNA base pairs and amino acids yield an extremely encouraging 1.4 kcal/mol error. Finally, and as a first look, we consider IPML in denser systems: water clusters, supramolecular host-guest complexes, and the benzene crystal.
I. INTRODUCTION
Classical potentials encode interaction physics but remain difficult to parametrize systematically for new molecules. IPML combines physics-based functional forms with ML predictions of local atomic properties to improve transferability while reducing molecule-specific reference calculations.
- Motivation: Classical potentials are limited to narrow molecular and material sets because their parametrization is tedious and non-systematic.Reproducing conformational and thermodynamic properties consistently across molecules remains challenging, time consuming, and difficult to automate.
- Motivation: Machine-learning potentials can reproduce reference energies accurately but typically interpolate only within the interactions represented in their training samples.A model trained on water clusters can describe liquid-water properties accurately while remaining specific to water interactions.
- IPML strategy: IPML combines physics-based models with ML-predicted parameters to preserve physical interaction forms and reduce reference calculations for new molecules.The approach is designed to leverage known symmetries and functional forms while alleviating molecule-specific optimization.
- IPML strategy: The model predicts distributed multipoles, atomic polarizabilities, and valence-density populations and decay rates for electrostatics, polarization, repulsion, and dispersion.These local atom-in-molecule properties replace much of the explicit parametrization normally required for each compound.
- IPML strategy: The electrostatic model improves multipole learning by using MBIS partitioning, neighbor-defined local axes, and the aSLATM representation.The representation encodes chemical elements, dispersion-scaled pair distances, and three-body configurations while remaining atom-index invariant.
2. Atomic-density overlap
The atomic-density-overlap model represents short-range interactions through environment-dependent valence-density populations and decay rates. ML predicts these quantities from molecular environments, reducing the number of fitted short-range parameters.
- Atomic-density overlap: Exchange-repulsion and other short-range interactions are modeled from overlapping electron densities.The approach uses Slater-type valence atomic densities to describe the overlap contribution.
- Atomic-density overlap: The decay rate for an atom pair is defined as σij = √σiσj from the individual atomic-density decay rates.Earlier models obtained decay rates from reference DFT calculations and fitted atom-type-dependent prefactors to short-range energies.
- Atomic-density overlap: Including valence-density populations reduces unknown dimer prefactors to one value for repulsion and short-range polarization, with no free penetration parameter.The populations are volume integrals of the valence atomic densities.
- ML parametrization: IPML trains ML models to predict atomic-density population N and decay rate σ from atom-in-molecule environments.Reference N and σ values were computed for 1,102 molecules, providing 16,945 atom-in-molecule properties.
- ML parametrization: Hirshfeld ratios describe environment-induced changes in atomic volumes and are used to estimate atomic polarizabilities from free-atom polarizabilities.Kernel-ridge regression predicts ratios for H, C, O, and N atoms from molecular geometry.
B. Intermolecular interactions from physics-based models
The physics-based interaction model combines ML-derived atomic properties with multipole electrostatics and short-range corrections. Its electrostatic treatment uses distributed multipole coefficients and addresses charge penetration when molecular electron densities overlap.
- Interaction-energy components: The interaction-energy model uses ML-derived local properties to construct the different intermolecular energy contributions.The section introduces how those properties enter the interaction-energy terms.
- Electrostatics: Distributed multipole electrostatics expands each atom’s electrostatic potential into multipole coefficients and contracts them through an interaction matrix.The multipole vector contains charges, dipoles, and higher moments, while the matrix contains derivatives of 1/r_ij.
- Electrostatics: The multipole coefficients used for electrostatics are supplied by the ML model developed and improved for this work.This connects the electrostatic energy expression to the learned atom-in-molecule properties.
- Charge penetration: At short intermolecular distances, wavefunction overlap invalidates the no-overlap assumption of the multipole expansion and produces charge-penetration effects.Penetration is linked to charge-density overlap and modeled by separating point charges into core and damped valence contributions.
- Charge penetration: The penetration treatment contains core-core, core–smeared-density damping, and smeared-density overlap terms, but does not correct penetration from higher multipoles.The higher-multipole limitation makes the separation between core and smeared contributions conceptually unclear.
3. Repulsion
The repulsive energy is modeled from overlaps of valence atomic densities, with element-dependent prefactors and multiplicative mixing rules.
- Valence atomic-density overlap provides the basis for parametrizing the repulsive energy.
- Element-dependent prefactors enter the repulsion model and have units of (energy)1/2.
- A multiplicative mixing rule combines atomic contributions to produce U rep.
5. Many-body dispersion
Many-body dispersion is treated with a quantum-harmonic-oscillator formulation, within an IPML potential combining five intermolecular contributions and only eight global parameters.
- 5. Many-body dispersion: Many-body dispersion uses an efficient random-phase-approximation formulation cast as a system of quantum harmonic oscillators.
- 5. Many-body dispersion: The intermolecular IPML model combines electrostatics, charge penetration, repulsion, induction/polarization, and many-body dispersion.
- 5. Many-body dispersion: Eight global parameters are optimized simultaneously across different compounds to assess transferability.
- 5. Many-body dispersion: The work provides a Python implementation using kernel ridge regression and makes the property-training datasets available.
- 5. Many-body dispersion: The potential terms are parametrized against reference energies on part of S22x5 and validated on additional intermolecular datasets.
C. Training of Hirshfeld ratios
The IPML workflow predicts local atomic properties for intermolecular potentials and evaluates these predictions through polarizabilities and dimer energies. Hirshfeld-ratio predictions correlate strongly with reference values, while energy accuracy varies with separation and depends on multipole quality.
- Training of Hirshfeld ratios: R2 = 99.5% and MAE = 0.006 for predicted Hirshfeld ratios on a separate test set.The model was trained on 12,300 atoms from 1,000 small organic molecules and tested on 17,100 atoms.
- Training of Hirshfeld ratios: The ML polarizability predictions show excellent agreement with experiment, with a MARE of 8.6%, virtually identical to Tkatchenko–Scheffler calculations after SCS.Fractional anisotropy tends to be underestimated, although overall agreement remains reasonable.
- Intermolecular energies: The model combines ML-predicted local properties with a small set of optimized global parameters fitted against representative dimer and host–guest datasets.Model 1 uses S22x5 at 0.9x and 1.0x distances, whereas Model 2 adds S12L host–guest complexes at 1.0x.
- Intermolecular energies: 0.7 kcal/mol is the overall MAE across S22x5 distance factors, decreasing from 1.0 kcal/mol at 0.9x to 0.2 kcal/mol at 2.0x.The distance-dependent trend indicates robust asymptotic behavior, while strongly hydrogen-bonded complexes produce outliers beyond ±1 kcal/mol.
- Intermolecular energies: Accurate intermolecular energies depend critically on multipole predictions because the few global parameters leave little room for error compensation.A poorer multipole model produced artifacts in hydrogen-cyanide partial charges and artificially strong hydrogen polarization.
- Intermolecular energies: The fitted polarization contribution is small, suggesting that short-range polarization energy may be absorbed into repulsion terms.The authors note that related AIM- and physics-based force fields use the same overlap model for repulsion and short-range polarization.
IV. PERFORMANCE OF THE IPML MODEL
IPML generalizes beyond its training dimers to angularly displaced, non-equilibrium molecular complexes. It achieves low error on the larger S66a8 benchmark and is competitive with a related overlap-based model, although comparisons use different datasets and error measures.
- S66a8 performance: 0.4 kcal/mol MAE was obtained across the 528 angularly displaced S66a8 dimers at the CCSD(T)/CBS reference level.The benchmark contains 66 dimers with eight geometries each and is described as a larger representative set than S22.
- Comparison with MEDFF: MEDFF reports an RMSE of 0.36 kcal/mol for dispersion-dominated S66 complexes at equilibrium, but the comparison is limited by different datasets and error measurements.The authors note that hydrogen-bonded complexes are typically more challenging and suggest their model compares favorably under that qualification.
B. Amino-acid side chains (SSI dataset)
IPML reproduces neutral amino-acid and biomolecular intermolecular energies across increasingly complex systems, with strong dimer performance but substantial errors for larger clusters and host–guest complexes. The benzene-crystal results improve when host–guest complexes inform global-parameter optimization, although systematic condensed-phase extension remains limited.
- Amino-acid side chains: 2,216 neutral amino-acid side-chain dimers show excellent agreement with CCSD(T)/CBS reference energies.Charged and sulfur-containing dimers were excluded; the dataset includes only neutral HCON compounds.
- Hydrogen-bonded complexes: 1.4 kcal/mol MAE is obtained for 127 neutral DNA-base and amino-acid dimers dominated by strong hydrogen-bonded complexes.The result is described as encouraging despite the directional difficulty of hydrogen-bond modeling.
- Water clusters: 8.1 kcal/mol MAE arises for water clusters containing 2–10 molecules as errors compound and progressively overstabilize larger clusters.The model still correlates highly with CCSD(T)/CBS reference energies.
- Supramolecular complexes: 9.7 kcal/mol MAE is found for host–guest complexes, despite high correlation with diffusion Monte Carlo reference energies.Including larger complexes in the global-parameter fit substantially improves Model 2, but one multiply hydrogen-bonded outlier remains overstabilized by 8 kcal/mol.
- Benzene crystal: −14.3 kcal/mol is Model 2’s benzene-crystal cohesive energy at equilibrium density, only 2 kcal/mol from experiment.Model 1 gives −17.2 kcal/mol versus the experimental −12.2 kcal/mol; Model 2 better reproduces the crystal despite understabilizing the two dimer configurations.
V. CONCLUSIONS AND FUTURE OUTLOOK
The conclusions present IPML as a transferable intermolecular-potential framework that predicts local atomic properties with ML while retaining physics-based interaction models. The outlook identifies improved multipoles, richer interactions, and computational-cost management as necessary for broader force-field and condensed-phase use.
- Conclusions: IPML predicts distributed multipoles, Hirshfeld ratios, valence-density parameters, and related local properties to model electrostatics, polarization, repulsion, and many-body dispersion.Its global parameters are optimized across H, C, N, and O chemical space rather than reoptimized for each new compound.
- Conclusions: The framework combines physics-based interaction forms with ML-predicted parameters instead of learning intermolecular energies purely from data.This preserves physical constraints and reduces the need for reference calculations when new molecules are encountered.
- Future outlook: Analytical derivatives could extend IPML toward a conformationally dependent force field, but balancing added accuracy against computational overhead remains unresolved.Derivatives are straightforward or already available for several terms, while force-field optimization itself is left open.
- Future outlook: 1–100 s evaluation times for 10–100-atom systems reflect the cost of explicit polarization, many-body dispersion, and especially multipole prediction.Approximately 90% of the reported time is spent predicting multipoles with the aSLATM representation.
- Future outlook: More accurate multipoles and advanced interactions such as anisotropic or many-body repulsion are identified as routes toward more transferable models.The conclusions emphasize that multipole accuracy is critical for reliable molecular and condensed-phase energies.
Appendix A: Many-body dispersion
The many-body dispersion implementation represents atoms as coupled quantum harmonic oscillators whose interactions depend on ML-derived polarizabilities and characteristic frequencies. Range separation and chemistry-independent parameters control the dipole coupling and resulting MBD energy.
- MBD construction: Atomic polarizabilities and characteristic frequencies provide the inputs for a coupled quantum harmonic-oscillator model of many-body dispersion.The model contains N atoms and uses frequency-dependent atomic polarizabilities to obtain dispersion coefficients.
- Dipole coupling: The dipole interaction tensor T_pq is built from derivatives of a modified Coulomb potential and is range-separated by a Fermi-function scaling.The range separation uses β and van der Waals radii scaled by a chemistry-independent fitting parameter.
- MBD energy: The eigenvalues of the coupled interaction matrix provide the many-body dispersion energy.The implementation follows the stated MBD formulation after constructing the range-separated coupling matrix.
- Parameters: β, γ, and d are chemistry-independent parameters in the MBD methodology.These parameters control the shared dispersion model rather than being fitted separately for each chemical system.
Appendix B: Covariant kernels
The appendix develops covariant kernels for predicting vector and second-rank tensor quantities while preserving rotational behavior. The construction represents atomic environments with Gaussian functions and analytically integrates the kernel over all 3D rotations.
- The covariant kernel is introduced for predicting dipoles from samples related by rotations.
- Atomic environments are encoded using atom-centered Gaussian functions.
- The kernel construction analytically integrates over all 3D rotations to enforce rotational covariance.
- The construction is extended to quadrupole moments by adapting the procedure to second-rank tensors.
- The kernel expression uses rotation matrices to align environment vectors with the z axis, and the outer product combines tensor components.