Learning protein-ligand binding affinity with atomic environment vectors

Scoring functions for the prediction of protein-ligand binding affinity have seen renewed interest in recent years when novel machine learning and deep learning methods started to consistently outperform classical scoring functions. Here we explore the use of atomic environment vectors (AEVs) and feed-forward neural networks, the building blocks of several neural network potentials, for the prediction of protein-ligand binding affinity. The AEV-based scoring function, which we term AEScore, is shown to perform as well or better than other state-of-the-art scoring functions on binding affinity prediction, with an RMSE of 1.22 pK units and a Pearson’s correlation coefficient of 0.83 for the CASF-2016 benchmark. However, AEScore does not perform as well in docking and virtual screening tasks, for which it has not been explicitly trained. Therefore, we show that the model can be combined with the classical scoring function AutoDock Vina in the context of \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta$$\end{document}Δ-learning, where corrections to the AutoDock Vina scoring function are learned instead of the protein-ligand binding affinity itself. Combined with AutoDock Vina, \documentclass[12pt]{minimal} \usepackage{amsmath} \usepackage{wasysym} \usepackage{amsfonts} \usepackage{amssymb} \usepackage{amsbsy} \usepackage{mathrsfs} \usepackage{upgreek} \setlength{\oddsidemargin}{-69pt} \begin{document}$$\Delta$$\end{document}Δ-AEScore has an RMSE of 1.32 pK units and a Pearson’s correlation coefficient of 0.80 on the CASF-2016 benchmark, while retaining the docking and screening power of the underlying classical scoring function. Supplementary Information The online version contains supplementary material available at 10.1186/s13321-021-00536-w.


Example of Atomic Environment Vector Computation Water
As a simple example of how atom-centred symmetry functions (ACFSs) are used to construct atomic environment vectors (AEVs), let us consider a simple water molecule with the fictitious coordinates of Table 1  If this is the only system we want to describe, we have only two elements (N e = 2) and we need to compute three AEVs, one for each atom.
Radial symmetry functions are parametrised by η R and R s ; for simplicity we only consider here η R = 1 and R s = {0, 1}. Angular symmetry functions are parametrised by η A , R s , θ s and ζ; for simplicity we only consider here η A = 1, R s = 0, θ s = 0 and ζ = 1.
The atomic environment vector for the oxygen atom (for N e = 2) has the following form: the cutoff distance R c (which we consider here large enough to include all atoms of the system). The atomic environment vector for the oxygen atom therefore reduces to Explicitly, the non-zero ACSFs composing the AEV for the oxygen atom are If we consider R c to be large enough so that f c (R) ≈ 1 and we use the geometry defined in Table 1, we can perform an explicit calculation for the particular configuration considered here (where R 3,1 = R 3,2 = 1 and θ 3,1,2 = π/2): The AEV for oxygen encodes its atomic environment and it is, by construction, rotationally and translationally invariant.
The same procedure can be repeated for every atom of the system, so that all atoms are described by their own AEV, so that we can describe the whole system with a matrix of AEVs of dimension N atoms × N AEVs (where N AEVs depends on the number of elements n e as well as the number of different values for the parameters R s , η R , . . . ).

Ammonia
Let us consider another simple example: ammonia. Again, we only have two elements (N e = 2) and we need to compute four AEVs, one for each atom. If we consider the same parametrisation for radial and angular symmetry functions described above, we have the following atomic environment vector for nitrogen: since there are no other nitrogen atoms in the system, the AEV for the only nitrogen atom reduces to Explicitly, denoting d N H the nitrogen-hydrogen distance, we have Using d N H = 1 and θ N ;HH = 109.5 we have the following atomic environment vector for nitrogen G N = 3e −1 , 3, 0, 0; 3 (1 + cos(109.5)) e −1 , 0, 0 . Discrepancies with the analytical calculation above come from the fact that the radial part is multiplied by a factor of 1/4 in the TorchANI implementation (see TorchANI code) and that the angle between two vectors is computed as

Gradients
The model, consisting of the AEVComputer and a collection of NNs, can be described as a function of the atomic coordinates f ( R), returning the binding affinity a: The loss function (MSE loss) is defined as the mean square difference between predicted and experimental affinities. For a single prediction a and the corresponding experimental affinity A, the loss function is: Given that a = f ( R), the gradient of the loss function with respect to the atomic coordinates can be computed using the chain rule The negative gradient of the loss function indicates the directions where the atoms can be moved to minimise the loss function, i.e. bring the predicted binding affinity a closer to the expected binding affinity A.
In some cases, such as classification of actives and decoys, it is clear what the desired outcome of the model is (i.e. an active molecule) and therefore the loss can be computed with respect to the desired output. In other cases, such as the prediction of the binding affinity, the desired outcome A (experimental binding affinity) is usually not known, and therefore it is not possible to compute the loss and its gradient. However, it remains possible to compute the gradient of the output with respect to the atomic coordinates.
The gradient of the predicted binding affinity a can be computed by differentiating f 8 with respect to the atomic coordinates R: This gradient can be computed with back-propagation and only requires a forward pass within the network. The positive gradient of the input indicates the directions where the atoms can be moved in order to increase (maximise) the predicted binding affinity.