Explainable chemical artificial intelligence from accurate machine learning of real-space chemical descriptors

Abstract

Machine-learned computational chemistry has led to a paradoxical situation in which molecular properties can be accurately predicted, but they are difficult to interpret. Explainable AI (XAI) tools can be used to analyze complex models, but they are highly dependent on the AI technique and the origin of the reference data. Alternatively, interpretable real-space tools can be employed directly, but they are often expensive to compute. To address this dilemma between explainability and accuracy, we developed SchNet4AIM, a SchNet-based architecture capable of dealing with local one-body (atomic) and two-body (interatomic) descriptors. The performance of SchNet4AIM is tested by predicting a wide collection of real-space quantities ranging from atomic charges and delocalization indices to pairwise interaction energies. The accuracy and speed of SchNet4AIM breaks the bottleneck that has prevented the use of real-space chemical descriptors in complex systems. We show that the group delocalization indices, arising from our physically rigorous atomistic predictions, provide reliable indicators of supramolecular binding events, thus contributing to the development of Explainable Chemical Artificial Intelligence (XCAI) models.

Fig. 1. Architecture of SchNet4AIM (modified SchNetPack toolbox).

Schematics of the SchNetPack toolbox targeted for the prediction of one-body (1P) and two-body (2P) terms, showing the main representation (left) and prediction (right) blocks. The atomic environment (AE) tensor, created by the representation block, describes an atom or atomic pair in a molecule. This is parsed to the prediction block which outputs a molecular property (P) decomposed into a collection of locally defined terms (K). The labels Z and R denote atomic numbers and positions of an M atom system, while r_ij signifies the distance between atoms i and j. AW, EAW and EPAW refer to AIMwise, ElementalAIMwise, and ElementalPairAIMwise output modules, respectively, used to construct universal and particle-specific models.

Fig. 2. Dispersion plots for T^O and EinterO−H.

Dispersion plots for the kinetic energy of the O atoms, T^O, and O–H interaction energies, ${\mathrm{E}}_{\mathrm{inter}}^{\mathrm{O}-\mathrm{H}}$, of the water cluster database. The training and testing data points are shown in blue and red, respectively. All values are reported in atomic units (a.u.), whereas the testing error metrics are given in kcal mol⁻¹. The mean absolute error (MAE) and root mean squared error (RMSE) were employed as measures of accuracy. The labels U and PS are used to refer to universal (AIMwise) and particle-specific (ElementalAIMwise or ElementalPairAIMwise) models, respectively. Source data are provided as a Source Data file.

Fig. 3. 1,3-dipolar cycloaddition.

Schematic representation of the addition of acetonitrile oxide to ethylene to yield a 5-membered heterocycle, showing the transition state (TS) structure involved in the transformation. The numbers and labels of the main atoms involved in the reaction are shown. The remaining, and with the exception of the C atom in the CH₃ moiety of acetonitrile oxide, are H atoms.

Fig. 4. SchNet4AIM predictions in extrapolation regimes: electronics of a 1,3-dipolar cycloaddition.

Evolution of the atomic charges of the main atoms involved in the reaction (a), the total molecular charge of each of the reactants (b), and the electron delocalization between the terminals atoms (1,3 and 4,7) of both species (c). Atomic or molecular charges are denoted as Q, whereas the term δ is used to refer to the delocalized electron counts. Solid and dashed lines are used to represent the predicted (SchNet4AIM) and observed data, respectively. The SchNet4AIM models are tested in never-seen (extrapolation) regions of the chemical space. All values are reported in electrons (e⁻) relative to the progress of the reaction coordinate (χ). The atomic labels correspond to the numbering shown in Fig. 3. d The transition state of the reaction, showing the emergent C–C and C–O bonds that will be formed, is also shown. Source data are provided as a Source Data file.

Fig. 5. 13P-CO₂ complexation and guest release.

a 13P-CO₂ system at 500 (i), 950 (ii), 1318 (iii), and 1600 (iv) fs throughout the 300  K simulation, corresponding to some of the local maxima in the electron delocalization (δ) between CO₂ and NH₂, δ(CO₂,NH₂). b Correlation maps between the electronic and geometrical descriptors (i.e., the distance between the centers of mass, ∣R∣) for the CO₂-Ph⁴ (j), CO₂-${\mathrm{NH}}_2^1$ (jj), and ${\mathrm{NH}}_2^3$-${\mathrm{NH}}_2^4$ (jjj) contacts throughout the 300 K simulation. Blue dots indicate those binding events exclusively predicted by SchNet4AIM, whereas magenta ones show contacts simultaneously (±10 fs) estimated by the geometrical and electronic metrics. Each tick corresponds to 1000 fs. c 13P-CO₂ system at 0, 100, 150, and 250 fs throughout the 900 K simulation, along with the SchNet4AIM computed electron delocalization between neighboring NH₂ moieties, δ(NH₂–NH₂). For the electronic metrics, reported in electrons (e⁻), the raw and bin-averaged data are shown. Time is reported in femtoseconds (fs). Source data are provided as a Source Data file.

Fig. 6. Chemical interactions between CO₂ and 13P after the cage cleavage.

Evolution of the SchNet4AIM predicted electron delocalization (δ) between CO₂ and the phenyl (Ph), oxo (O), and amino (NH₂) moieties of 13P throughout the 13P-CO₂ 900 K simulation and after the cage rupture event. Time is reported in femtoseconds (fs). The labels R¹, R², R³, and R⁴ denote the position of the NH₂, O, and Ph moieties in 13P using the same numbering and color code as that shown in Fig. 5A. Source data are provided as a Source Data file.