A three-dimensional protein-ligand activity prediction method based on attention mechanism

Abstract

This invention belongs to the field of artificial intelligence-assisted drug design technology, and specifically discloses a three-dimensional protein-ligand activity prediction method based on attention mechanism. This invention constructs a molecular three-dimensional information extraction network for molecular characterization, and extracts the three-dimensional information of protein ligand binding pose by combining spherical information transmission method with physically meaningful features; it then constructs a three-dimensional graph neural network model based on attention mechanism to predict protein-ligand interactions, including conformation prediction and binding affinity prediction, and finally outputs the prediction results.

Description

Technical Field

[0001] This invention belongs to the field of artificial intelligence-assisted drug design technology, and specifically relates to a three-dimensional protein-ligand activity prediction method based on attention mechanism. Background Technology

[0002] With the development of deep learning technology, researchers have applied convolutional neural networks to predict protein-ligand affinity in the drug discovery process. One approach to make convolutional neural networks suitable for this application is to regularize the molecular structure. For example, the Pafnucy model trims the complex structure to fit the geometric center of the ligand molecule. Cube size, in each The resolution is used to mesh the molecule, and the chemical information of each atom within it is described by a 19-dimensional vector, resulting in a 4D tensor representation of the protein ligand molecule. OnionNet computes this by arranging the ligand atom around a radius of... The protein structure is divided into an onion-like hierarchical structure, and the occurrence frequency of specific atomic pairs in each layer is calculated. These methods, in order to adapt to the data regularity requirements of convolutional neural networks, all perform molecular structure segmentation, which to some extent results in the loss of some important information.

[0003] Another approach is to treat molecules as graph data structures. AGL-Score proposes a multi-scale weighted label algebraic subgraph method to represent molecules. After using the subgraph to calculate the eigenvalues ​​and eigenvectors of the Laplacian matrix and adjacency matrix to describe the molecule, these are input into a gradient boosting tree model to train the final model. C-SGEN is a molecular property prediction method based on spatial graph embedding, including the prediction of protein-ligand binding affinity. This method first represents the molecule as a graph embedding representation, uses an adjacency matrix to represent the inter-atomic connections, and then inputs the adjacency matrix and graph embedding representation into a convolutional neural network for one-dimensional convolution.

[0004] In recent years, graph convolutional neural networks have been used to predict protein-ligand binding affinity. PotentialNet is also a graph neural network method proposed for predicting molecular properties. When defining the edges in the graph, it considers the Euclidean distance between atoms in addition to chemical bonds. The network structure uses gated recurrent units to update atomic features. GraphDelta uses the BPS symmetric function to establish a single-atom environment vector to represent the molecule. It directly concatenates the vectors of the atom, its neighboring atoms, and the bonds between them, and inputs them into a fully connected layer. The output vectors of the neighboring atoms are summed to obtain molecular information, which is then directly concatenated with the vector of the central atom and input into a second fully connected layer for information updates. The information updates of all atoms of the ligand are combined, and after fuzzing, they are input into a third fully connected layer to obtain the final result.

[0005] An IGN scoring function model based on graph representation learning discloses intramolecular graph convolution and intermolecular graph convolution modules to sequentially learn intramolecular and intermolecular interactions in protein-ligand complexes. The learned intermolecular interactions are then applied to downstream task predictions, including protein-ligand binding affinity prediction, large-scale virtual screening, and small molecule binding conformation prediction. An invention of a protein-ligand affinity evaluation method based on graph attention networks addresses the issue of small dataset size by simulating a semi-flexible docking process. A data augmentation method based on molecular flexibility is designed to expand the data scale, improve model quality, and ensure the scientific validity and rationality of the method. A molecular feature extraction method based on a graph attention mechanism is designed to extract effective molecular features, improving the accuracy and performance of the scoring function.

[0006] Traditional methods and traditional machine learning methods still have some shortcomings in the prediction of receptor and ligand activity. (1) The huge computational cost limits their application in high-throughput screening. Classic molecular docking methods, including DOCK, AutoDock, GOLD, etc., have been used to predict binding affinity. Although these methods speed up the computation through principled parameter fitting, their computational accuracy is not satisfactory. Among them, the machine learning-based method ΔvinaRF 20 It excels in predictive accuracy. However, traditional machine learning methods rely on manual feature engineering, with models heavily dependent on feature construction and selection. Currently, it's difficult to find a feature engineering method universally applicable to the entire small molecule chemical space. This limitation can be overcome in deep learning models. Therefore, deep learning methods are increasingly favored by researchers.

[0007] (2) In most studies on molecular representation using scoring functions, numerous molecular representation learning methods have been proposed. Most of these methods use SMILES strings as input. While these models are powerful, the fact that SMILES represents a one-dimensional linear representation of molecular structure makes it difficult for neural networks to learn the original structural information of molecules from SMILES. Therefore, methods based on handcrafted descriptors, one-dimensional protein sequences, and / or two-dimensional graphical representations limit their ability to learn generalized molecular interactions in three-dimensional space.

[0008] (3) In the protein-ligand interaction prediction model, deep learning methods provide a new approach for this invention through end-to-end learning of complex neural networks. The emerging GNN has been introduced into the structure-based protein-ligand interaction prediction workflow. Graph Neural Networks (GNNs) apply deep learning to graph structures, and the Graph Convolutional Network (GCN) within it can perform convolution operations on the graph. However, GCNs have some drawbacks: they rely on the Laplacian matrix and cannot be directly used for directed graphs; model training depends on the entire graph structure and cannot be used for dynamic graphs; they cannot assign different weights to neighboring nodes during convolution; training is a fully incremental update, making it difficult to scale to large-scale networks, and convergence is slow. Summary of the Invention

[0009] To address the shortcomings of existing technologies, this invention provides a three-dimensional protein-ligand activity prediction method based on an attention mechanism. It uses a three-dimensional model for molecular characterization and combines spherical information transmission with physically meaningful features derived from the Schrödinger equation as a feature extraction network to effectively extract the three-dimensional information of protein-ligand binding pose. This invention constructs an attention-based prediction model to predict protein-ligand interactions. The attention mechanism can assign different weights to different nodes and relies on paired adjacent nodes during training, rather than specific network structures. A distance-aware attention algorithm is designed to distinguish various types of intermolecular interactions. Furthermore, a variant of a graph neural network suitable for learning protein-ligand interactions is introduced, allowing the model to focus on intermolecular interactions.

[0010] To solve the above-mentioned technical problems, the technical solution adopted by the present invention is as follows:

[0011] A three-dimensional protein-ligand activity prediction method based on an attention mechanism includes the following steps:

[0012] Step 1: Construct a dataset by taking the three-dimensional structural information of proteins and ligands as input to predict binding affinity, build a virtual screening model, and simulate a virtual screening scenario.

[0013] Step 2: Construct a molecular 3D information extraction network for molecular characterization and extraction of 3D information of protein-ligand complexes. The 3D information extraction network includes an input module, a message passing module, multiple interaction modules, and an output module. The input module is used to construct initial messages for the edges k of the graphical structure. The message passing module uses a spherical message passing method and multiple physical representations to uniquely represent the graphical structure with 3D information and merges the relative position information of atoms into a spherical coordinate system. The interaction module updates the initial messages by merging all physical representations. The output module aggregates all incoming messages to update the features of the nodes.

[0014] Step 3: Construct a three-dimensional graph neural network model based on an attention mechanism to predict protein-ligand interactions, including conformation prediction and binding affinity prediction; wherein the three-dimensional graph neural network model based on the attention mechanism adopts a graph attention model, including a graph attention layer, a connection layer, a graph pooling layer, and a fully connected layer, and the following method is used for prediction:

[0015] Step 3.1: Construct two adjacency matrices and embed the extracted 3D information of the protein-ligand complex into the two adjacency matrices to learn the influence of protein-ligand interactions on the node features of each atom.

[0016] Step 3.2: Design a distance-aware attention algorithm to distinguish various types of intermolecular interactions; use a graph attention layer to learn protein-ligand interactions, so that the model focuses on intermolecular interactions, and finally obtain the node features updated by the graph attention layer gating enhancement.

[0017] Step 3.3: After combining the features obtained in Step 3.2 to obtain aggregate features, input them into the graph pooling layer to learn the pairwise non-bonded interactions of atoms to generate intermolecular interactions of the complex.

[0018] Step 3.4: Use fully connected layers enhanced with batch normalization layers as the task network for the final protein-ligand interaction prediction to generate prediction results.

[0019] Furthermore, the dataset constructed in step 1 includes four datasets: PDBbind, DUD-E, DEKOIS2.0, and LIT-PCBA. The PDBbind dataset is used to combine affinity prediction, DUD-E and DEKOIS2.0 are multi-objective datasets used to build target-free virtual screening models, and the LIT-PCBA dataset is used to simulate real virtual screening scenarios.

[0020] Furthermore, the dataset processing method is as follows:

[0021] Step 1.1: For each sample in the PDBbind dataset, perform docking calculations to generate the binding pose of the protein-ligand complex. If the root mean square error (RMSD) of the generated pose compared to the experimentally validated binding structure is less than 100%, then the binding pose is considered valid. Then mark it as a positive sample if the RMSD is greater than 1. Then it is marked as a negative sample, with RMSD between 2 and Samples between them were omitted;

[0022] Step 1.2: Divide the PDBbind dataset into training and testing sets based on proteins to ensure that the training and testing sets do not share the same proteins; in addition, remove PDBbind samples containing proteins from the DUD-E dataset from the training and testing sets.

[0023] Furthermore, the specific methods for molecular characterization in step 2 are as follows:

[0024] Step 2.1: Use the relative position information of atoms as input to the three-dimensional information extraction network to identify the structure of three-dimensional molecules in spherical coordinates;

[0025] For any point in a spherical coordinate system, its position is determined by a set of triples. Specify, where d, θ and Let represent radial distance, polar angle, and azimuth angle, respectively; when modeling a 3D molecular diagram in a spherical coordinate system, any atom i can be the origin of the local spherical coordinate system, where d, θ, and These are respectively bond length, bond angle, and torsion angle; therefore, the relative positions of each adjacent atom i can be represented by the corresponding tuples. designated;

[0026] Step 2.2: Define a three-dimensional molecular diagram G, represented as a 4-tuple G = (u, V, E, P). It is the global eigenvector of the molecular graph G, where d u It is the dimension of the global features; V = {v i} i＝1:n These are atomic features, each of which Let be the eigenvector of atom i, and n be the total number of atoms. E = {(e_i)} k ,r k ,s k )} k＝1:m It is a set of edges, where each r is the edge eigenvector of edge k. k It is the receiver atom, s k Let $k$ be the sending atom of edge $k$, and $m$ be the total number of edges. $P = {r$ k} k＝1:n It is a set of three-dimensional Cartesian coordinates, which contains the three-dimensional spatial information of each atom;

[0027] Step 2.3: Perform spherical message passing in spherical coordinates and update the message;

[0028] make N represents the set of edges pointing to atom i. i This represents the index of the input node of atom i; the output after the message passing process includes the updated global feature vector. Updated atomic features V′={v′ i} i＝1:n and the updated edge set E′={(e′ k ,r k ,s k )} k＝1:m , where v′ i Let e′ represent the updated atomic eigenvector. k Let represent the updated edge feature vector. Spherical message passing is defined in spherical coordinates as shown in equation (1):

[0029]

[0030] Where φ e ,φ v and φ u These are three information update functions for edges, atoms, and the entire molecular graph, respectively. The receiver atom r of edge k k eigenvectors, The sender atom s of edge k k eigenvectors, The sender atom s that points to edge k k edge set, r h This indicates the accepting atom of edge h. The sender atom s that points to edge k k The index of the input node, E′ i This represents the updated set of edges pointing to atom i; the aggregation function ρ e→v and ρ v→u Aggregating information across different geometric types, particularly in spherical messaging, the three-dimensional information in P is transformed and merged to update each message; therefore, spherical messaging uses another position aggregation function ρ. p→e Update the message;

[0031] Step 2.4: Select basis functions to convert triples. It is then converted into a physics-based representation for molecular characterization.

[0032] Furthermore, in step 2.4, spherical Bessel functions are chosen as the basis functions for encoding distance, and spherical harmonic functions are chosen as the basis functions for encoding angle; the final molecular characterization is expressed as follows: One representation is shown in formula (2):

[0033]

[0034] Where j l (·) is a spherical Bessel function of order l. It is a spherical harmonic function of order m and l, where c represents the cutoff radius and β lnIt is the nth root of the l-th order Bessel function; where l∈[0,…,L] H -1], m∈[-l, ...,l] and n∈[1, ...,N B ], L H and N B represents the highest order of the spherical harmonic function and the spherical Bessel function, respectively.

[0035] Furthermore, in step 3.1, the adjacency matrix is ​​represented as A, which is an N matrix ×N matrix A matrix, where if the i-th and j-th atomic nodes are connected, then A ij >0, otherwise A ij =0, N matrix This represents the number of atomic nodes in the graph; the two adjacency matrices A are constructed. 1 and A 2 A 1 The adjacency matrix A represents pure covalent interactions. 2 The adjacency matrix A represents the interaction between covalent and non-covalent molecules. 1 and A 2 The construction is shown in formulas (3) and (4):

[0036]

[0037]

[0038] Where d ij It is the distance between the i-th and j-th atoms, and μ and σ are learnable parameters, as shown in formula (4). This reflects that intermolecular bonds are weaker than covalent bonds, and their strength decreases as the bond distance increases.

[0039] Furthermore, step 3.2 designs a distance-aware attention algorithm, where the input to the graph attention layer is the adjacency matrix A and node features. Where x∈R F N atom F represents the number of nodes, i.e., the number of atoms, and F is the dimensionality of the nodes; the graph attention layer generates a new set of node features. Where x∈R F ;

[0040] The specific method is as follows:

[0041] Step 3.2.1: Calculate the attention coefficient of each node;

[0042] First, in order to obtain sufficient expressive power, each node feature By a learnable weight matrix W∈R F×FThe transformation is represented as Then, the attention coefficient is obtained, as shown in formula (5):

[0043] e ij =x i ′T Ex′ j +x ′T Ex′ i (5)

[0044] Where E∈R F×F It is also a learnable matrix, with attention coefficient e ij This indicates the importance of the feature of the j-th node to the feature of the i-th node; by analyzing x... i ′T Ex′ j and The summation is used to force e ij =e ji To reflect the structure of the graph, the attention coefficient e ij Only calculate j∈N i , where N i A is a neighbor node of the i-th node. ij >0;

[0045] Then, in order to manipulate the scale of the attention coefficients of each node, the attention coefficients are normalized among their neighbors, and A is... ij Multiply by the normalized attention coefficient, normalized attention coefficient a ij Given by formula (6):

[0046]

[0047] Step 3.2.2: Update node features;

[0048] After obtaining the normalized attention coefficient a ij Then, each node feature is updated to a linear combination of the node features of its neighboring nodes and the normalized attention coefficient, as shown in Equation (7).

[0049]

[0050] Step 3.2.3: Introduce the gating enhancement algorithm;

[0051] A gated augmentation algorithm is introduced to directly pass the feature information of the previous node to the next layer, and the output of the graph attention layer is then used. Implemented as x i 'and x i The linear combination of "" is shown in formulas (8) and (9).

[0052] x i out=z i x′ i +(1-z i )x″ i (8)

[0053] z i =σ(U(x) i ′||x″ i (9)

[0054] The output of the graph attention layer that combines all nodes Form a new node feature set Where U∈R 2F×1 Let z be a learnable vector, b be a learnable scalar value, σ represent the sigmoid activation function, (·||·) is the concatenation of two vectors, and z is a learnable vector. i This indicates how much information about the features of the input node will be directly passed to the next layer.

[0055] Furthermore, in step 3.2, for the two adjacency matrices A 1 and A 2 The graph attention layer generates new node features x. out1 and x out2 x out1 =GAT(x in A 1 ), x out2 =GAT(x in A 2 ), where GAT represents a gated augmented graph attention layer, used in calculating x. out1 and x out2 At the same time, they share a gate-enhanced graph attention layer; the output node features x out It is through x out2 Subtract x from the middle out1 The result is shown in formula (10).

[0056] x out =x out2 -x out1 (10)

[0057] After updating the node feature vectors through multiple graph attention layers, the feature vectors of the ligand atoms are summarized into a vector representing the protein-ligand complex graph, as shown in Equation (11):

[0058]

[0059] Finally, the multilayer perceptron was applied to x graphThe activity of protein-ligand complexes or binding postures is classified, and the Leaky ReLU activation function is used between layers, with the sigmoid function used after the last layer.

[0060] Compared with the prior art, the advantages of this invention are:

[0061] 1. Three-dimensional structural information is crucial for drug discovery, but a principled message passing framework is lacking to address it. This invention uses spherical message passing as a unified and efficient solution, enabling near-complete molecular representations without increasing computational budget. Using meaningful three-dimensional information can significantly improve the performance of prediction tasks.

[0062] 2. This invention proposes a three-dimensional graph neural network model based on an attention mechanism to predict protein-ligand interactions, enabling the model to distinguish the contribution of each interaction to binding affinity.

[0063] The model employs an attention mechanism, which allows different weights to be assigned to different nodes. During training, it relies on pairs of adjacent nodes rather than on the specific network structure. The model with the added attention mechanism can distinguish key active site residues and interactions through attention coefficients. Using active site information in the form of active site maps and fingerprints helps to effectively learn the interaction between target proteins and drug molecules.

[0064] Protein-ligand interaction prediction is performed using graphical features obtained by subtracting each feature of the target protein and a given ligand from the graphical features of their complexes. These strategies allow the model to learn key factors for accurately predicting protein-ligand interactions by focusing on intermolecular interactions rather than simply memorizing certain patterns of ligand molecules. Furthermore, it outperforms previous deep learning models in virtual screening, pose prediction, and docking.

[0065] Furthermore, the model of this invention can reproduce the natural population distribution of active and inactive molecules. Attached Figure Description

[0066] To more clearly illustrate the technical solutions of the embodiments of the present invention, the drawings used in the description of the embodiments will be briefly introduced below. Obviously, the drawings described below are only some embodiments of the present invention. For those skilled in the art, other drawings can be obtained based on these drawings without creative effort.

[0067] Figure 1 This is a flowchart of the method of the present invention;

[0068] Figure 2 This is a diagram of the three-dimensional graph neural network model architecture based on the attention mechanism of the present invention;

[0069] Figure 3 This is a schematic diagram of the molecular three-dimensional information extraction model of the present invention. Detailed Implementation

[0070] The present invention will be further described below with reference to the accompanying drawings and specific embodiments.

[0071] Combination Figures 1-3 The present invention provides a three-dimensional protein-ligand activity prediction method based on an attention mechanism, comprising the following steps:

[0072] Step 1: Construct a dataset. Use all .mol2 / .sdf files containing the three-dimensional structure information of proteins and ligands as the basic input format. Use RDKit to convert the basic input format into the .pkl format that can be recognized by the deep learning model. This data will be used to predict binding affinity, build virtual screening models, and simulate virtual screening scenarios.

[0073] This embodiment constructs four datasets: PDBBind V2020 (2020 version), DUD-E, DEKOIS2.0, and LIT-PCBA. For affinity prediction, due to the high data requirements of deep learning models, both the general set from PDBBind V2020 and the core set from PDBBind V2016 are used. Recent research indicates that various deep learning models with good predictive ability in affinity prediction perform poorly in large-scale structure-based virtual screening. Therefore, this embodiment tests the model's performance on structure-based virtual screening using the large-scale DUD-E, DEKOIS2.0, and LIT-PCBA datasets. More specifically, DUD-E and DEKOIS2.0 are multi-objective datasets used to construct target-free virtual screening models. LIT-PCBA is a more challenging dataset with a Decoys / Actives ratio exceeding 1000, thus it was used to simulate realistic virtual screening scenarios. LIT-PCBA is a target-specific dataset, and this embodiment uses it to construct a target-specific virtual screening model.

[0074] As a preferred implementation method, the data processing method in step 1 is as follows:

[0075] Step 1.1: For each sample in the PDBbind dataset, perform docking calculations to generate possible binding poses of the protein-ligand complex. If the root mean square error (RMSD) of a generated pose compared to the experimentally validated binding structure is less than 100%, then the possible binding pose is considered a valid binding pose. Then mark it as a positive sample if the RMSD is greater than 1. Then it is marked as a negative sample, with RMSD between 2 and Samples between [the specified samples] were omitted.

[0076] Step 1.2: Divide the PDBbind dataset into training and testing sets based on proteins to ensure that the training and testing sets do not share the same proteins; in addition, remove PDBbind samples containing proteins from the DUD-E dataset from the training and testing sets.

[0077] Step 2: Construct a molecular three-dimensional information extraction network, perform molecular characterization, and extract the three-dimensional information of the protein-ligand complex.

[0078] Combination Figure 3 As shown in the diagram, LB2 represents a linear block with two linear layers, σ(LB) represents a linear layer followed by an activation function, || represents concatenation, and ⊙ represents element-wise multiplication. The 3D information extraction network includes an input module, a message passing module, multiple interaction modules, and an output module. The input module constructs an initial message for the edge k of the graph structure; the message passing module uses a spherical message passing method and multiple physical representations to uniquely represent the graph structure with 3D information and merges the relative position information of atoms into a spherical coordinate system; the interaction module updates the initial message by merging all physical representations; and the output module aggregates all incoming messages to update the node features.

[0079] As a preferred embodiment, the specific method for molecular characterization in step 2 is as follows:

[0080] Step 2.1: Use the relative position information of atoms as input to the three-dimensional information extraction method to identify the structure of three-dimensional molecules in the spherical coordinate system.

[0081] The relative three-dimensional information can be distance or angle, which has inherent invariance to the translation and rotation of the input molecule. For any point in the spherical coordinate system, its position is determined by a set of triplets. Specify, where d, θ and Let represent radial distance, polar angle, and azimuth angle, respectively; when modeling a 3D molecular diagram in a spherical coordinate system, any atom i can be the origin of the local spherical coordinate system, where d, θ, and These are respectively bond length, bond angle, and torsion angle; therefore, the relative positions of each adjacent atom i can be represented by the corresponding tuples. specified.

[0082] Step 2.2: Define a three-dimensional molecular diagram G, represented as a 4-tuple G = (u, V, E, P). It is the global eigenvector of the molecular graph G, where d u It is the dimension of the global features; V = {v i} i＝1:n These are atomic features, each of which Let be the eigenvector of atom i, and n be the total number of atoms. E = {(e_i)} k ,r k ,s k )} k＝1:m It is a set of edges, where each r is the edge eigenvector of edge k. k It is the receiver atom, s k Let $k$ be the sending atom of edge $k$, and $m$ be the total number of edges. $P = {r$ k} k＝1:n It is a set of three-dimensional Cartesian coordinates, which contains the three-dimensional spatial information of each atom.

[0083] Step 2.3: Perform spherical message passing in spherical coordinates and update the message;

[0084] make N represents the set of edges pointing to atom i. i This represents the index of the input node of atom i; the output after the message passing process includes the updated global feature vector. Updated atomic features V′={v′ i} i＝1:n and the updated edge set E′={(e′ k ,r k ,s k )} k＝1:m , where v i e′ represents the updated atomic eigenvector, e′ k Let represent the updated edge feature vector. Spherical message passing is defined in spherical coordinates as shown in equation (1):

[0085]

[0086] Where φ e ,φ v and φ u These are three information update functions for edges, atoms, and the entire molecular graph, respectively. The receiver atom r of edge k k eigenvectors, The sender atom s of edge k k eigenvectors, The sender atom s that points to edge k k edge set, rh This represents the receiving atom of edge h. The sender atom s that points to edge k k The index of the incoming node, h is r k s k , The union of sets, E′ i This represents the updated set of edges pointing to atom i. The aggregation function ρ e→v and ρ v→u Aggregating information across different geometric types, particularly in spherical messaging, the three-dimensional information in P is transformed and merged to update each message; therefore, spherical messaging uses another position aggregation function ρ. p→e Update the message;

[0087] Step 2.4: Select basis functions to convert triples. It is then converted into a physics-based representation for molecular characterization.

[0088] Molecular characterization design needs to follow physical laws. An important aspect is the selection of appropriate basis functions. In this embodiment, spherical Bessel functions are chosen as the basis functions for encoding distance, and spherical harmonic functions are chosen as the basis functions for encoding angle. The final molecular characterization is expressed as follows: One representation is shown in formula (2):

[0089]

[0090] Where j l (·) is a spherical Bessel function of order l. It is a spherical harmonic function of order m and l, where c represents the cutoff radius and β ln It is the nth root of the l-th order Bessel function; where l∈[0,…,L] H -1], m∈[-l, ...,l] and n∈[1, ...,N B ], L H and N B and represent the highest order of the spherical harmonic function and the spherical Bessel function, respectively; these are hyperparameters in the experimental environment. Furthermore, this embodiment can also be derived from... Two simplified representations, Ψ(d) and Ψ(d,θ), are derived.

[0091] Step 3: Construct a three-dimensional graph neural network model based on the attention mechanism to predict protein-ligand interactions.

[0092] The attention-based 3D graph neural network model adopts the GAT graph attention model, which includes a graph attention layer, a connection layer, a graph pooling layer, and a fully connected layer.

[0093] A graph can be defined by (V, E, A), where V is a set of nodes, E is a set of edges, and A is an adjacency matrix. In a property graph, the attributes of each node are typically represented by a vector. The adjacency matrix A is an N... matrix ×N matrix A matrix, where if the i-th and j-th atomic nodes are connected, then A ij >0, otherwise A ij =0, N matrix This represents the number of atomic nodes in the graph. Step 3 mainly includes three stages: (i) updating node features, (ii) aggregating node features and processing graph features, and (iii) predicting the graph label. In the first stage, the node features representing the attributes of the i-th atomic node are updated in several message passes between neighboring nodes. This stage aims to obtain a high-level representation of the node features. Then, the updated node features are aggregated to generate graph features. Here, the result of the aggregation must be invariant to changes in the node order.

[0094] As a preferred implementation, the following method is used for prediction:

[0095] Step 3.1: Construct two adjacency matrices and embed the three-dimensional structural information of the protein-ligand complex into the two adjacency matrices to learn the influence of protein-ligand interactions on the node features of each atom.

[0096] The two adjacency matrices A are constructed 1 and A 2 A 1 The adjacency matrix A represents pure covalent interactions. 2 The adjacency matrix A represents the interaction between covalent and non-covalent molecules. 1 and A 2 The construction is shown in formulas (3) and (4):

[0097]

[0098]

[0099] Where d ij It is the distance between the i-th and j-th atoms, and μ and σ are learnable parameters, as shown in formula (4). This reflects that intermolecular bonds are weaker than covalent bonds, and their strength decreases as the bond distance increases.

[0100] Step 3.2: Design a distance-aware attention algorithm to distinguish various types of intermolecular interactions; use a graph attention layer to learn protein-ligand interactions, so that the model focuses on intermolecular interactions, and finally obtain the node features updated by attention layer gating enhancement.

[0101] The input to the graph attention layer is the adjacency matrix A and the node features. Where x∈R F N atom F represents the number of nodes, i.e., the number of atoms, and F is the dimensionality of the nodes; the graph attention layer generates a new set of node features. Where x∈R F .

[0102] The specific method is as follows:

[0103] Step 3.2.1: Calculate the attention coefficient of each node;

[0104] First, in order to obtain sufficient expressive power, each node feature By a learnable weight matrix W∈R F×F The transformation is represented as Then, the attention coefficient is obtained, as shown in formula (5):

[0105] e ij =x i ′T Ex′ j +x ′T Ex′ i (5)

[0106] Where E∈R F×F It is also a learnable matrix, with attention coefficient e ij This indicates the importance of the feature of the j-th node to the feature of the i-th node; by analyzing x... i ′T Ex′ j and The summation is used to force e ij =e ji To reflect the structure of the graph, the attention coefficient e ij Only calculate j∈N i , where N i A is a neighbor node of the i-th node. ij >0, because the adjacency matrix reflects both connectivity and normalized distance.

[0107] Then, in order to manipulate the scale of the attention coefficients of each node, the attention coefficients are normalized among their neighbors. Furthermore, this invention will... ij Multiplying by a normalized attention coefficient reflects the likelihood that nodes with shorter Euclidean distances are more important than others. This can be considered a form of inductive bias. Therefore, the normalized attention coefficient a ij Given by formula (6):

[0108]

[0109] Step 3.2.2: Update node features;

[0110] After obtaining the normalized attention coefficient a ij Then, each node feature is updated to a linear combination of the node features of its neighboring nodes and the normalized attention coefficient, as shown in Equation (7).

[0111]

[0112] Step 3.2.3: Introduce the gating enhancement algorithm;

[0113] A gated augmentation algorithm can significantly improve model performance. This invention introduces a gated augmentation algorithm mechanism that directly passes the feature information of the previous node to the next layer, and passes the output of the graph attention layer... Implemented as x i 'and x i The linear combination of "" is shown in formulas (8) and (9).

[0114] x i out =z i x i ′+(1-z i )x i " (8)

[0115] z i =σ(U(x) i ′||x i ")+b) (9)

[0116] The output of the graph attention layer that combines all nodes Form a new node feature set Where U∈R 2F×1 Let z be a learnable vector, b be a learnable scalar value, σ represent the sigmoid activation function, (·||·) is the concatenation of two vectors, and z is a learnable vector. i This indicates how much information about the features of the input node will be directly passed to the next layer.

[0117] This invention, through the aforementioned method, enables the model to focus on intermolecular interactions rather than merely memorizing certain patterns of ligand molecules, thereby learning key factors for accurately predicting protein-ligand interactions. The input to the neural network of this invention is the node feature x. in Adjacency matrix A 1 and A 2 After the above steps, the graph attention layer generates new node features x. out1 and x out2 x out1 =GAT(x in A1 ), x out2 =GAT(x in A 2 ), where GAT represents the gated augmented graph attention layer. It's important to note that when calculating x... out1 and x out2 At the same time, they share a gate-enhanced graph attention layer.

[0118] Output node features x out It is through x out2 Subtract x from the middle out1 The result is shown in formula (10).

[0119] x out =x out2 -x out1 (10)

[0120] Then, by subtracting the features of the two nodes, the model of this invention learns the difference between the structure under the binding posture and the structure after separation. After updating the node feature vectors through multiple graph attention layers, the feature vectors of the ligand atoms are summarized into a vector representing the protein-ligand complex graph, as shown in Equation (11):

[0121]

[0122] Finally, the multilayer perceptron was applied to x graph The activity of protein-ligand complexes or binding postures is classified, and the Leaky ReLU activation function is used between layers, with the sigmoid function used after the last layer.

[0123] Step 3.3: After combining the features obtained in Step 3.2 to obtain aggregate features, input them into the graph pooling layer to learn the pairwise non-bonded interactions of atoms to generate intermolecular interactions of the complex.

[0124] In the classical molecular force field, the mapping from pairwise atomic interactions to total nonbonded interactions is a sum, which is implemented by a neural network, where the specific functional form is automatically learned by the model, as shown in Equation (12).

[0125]

[0126] Among them, graph pooling layers are only applied to B. new This can be viewed as the intermolecular interactions of paired atom nonbonding interactions to generate complexes; B ij new Represents graph G pl The edge; is a random parameter; tanh is a nonlinear activation; E plIs it Figure G? pl The edge set in the MLP; D″ is the output dimension of the MLP; MAXPOOL is the maximum element-wise value of all edge features; G * This is the output of the graph pooling layer.

[0127] The output G of the graph pooling layer * Weighted and pooled and max pooling This approach combines information from multiple views, where max pooling highlights edge features with the highest information intensity, and weighted pooling utilizes the total information intensity of the edge features.

[0128] Step 3.4: Use fully connected layers enhanced with batch normalization layers as the task network for the final protein-ligand interaction prediction to generate prediction results.

[0129] This invention can perform various prediction tasks, such as virtual screening and attitude prediction. The loss calculation during prediction is as follows:

[0130] (1) The training objective is to minimize the loss function. Given N... sample Given a dataset of samples, the predicted value and the true value of a given sample are respectively... and y i The loss function (MSE loss) for the affinity prediction task (regression task) is defined as follows, as shown in Equation (13):

[0131]

[0132] (2) The Focal loss for structure-based virtual screening and pose prediction tasks (classification tasks) is defined as follows, as shown in Equation (14):

[0133]

[0134] Focal loss is an improved function based on cross-entropy to handle class imbalance. It reduces the loss contribution of easily classified samples, forcing the model to focus on hard-to-classify samples during training. Here, α is a weighting factor balancing the importance of the active and the decoy, set as (number of decoys) / (number of all complexes). γ is a pooling parameter that adjusts the reduction rate of easily classified samples, set to 2.0 in the experiments.

[0135] The training, performance verification, and evaluation of the model of this invention are described below.

[0136] 1. Prepare the dataset.

[0137] (1) For the affinity prediction task (regression task), the entire core set of PDBbind V2016 was used as the external test set, 20% of the general set was randomly selected as the validation set, and the remainder was used as the training set. To further test the robustness of the model, another dataset (PDBBind V2013 core set) was also used as the external test set. Since some complexes in the PDBbind V2013 core set were also included in the general set of PDBbind V2020, these complexes were directly removed from the original PDBbind V2013 core set, resulting in this external test set containing only 95 complexes.

[0138] (2) For structure-based virtual screening tasks (classification tasks), two models (target-free and target-specific virtual screening models) are trained. In training the target-free virtual screening model, the entire DUD-E dataset is randomly divided into training and validation sets using stratified sampling at a 4:1 ratio. Given the potential hidden bias in the DUD-E dataset, this invention does not use any test complexes from the DUD-E dataset, as the statistics may show varying degrees of over-optimism. Instead, the entire DEKOIS2.0 dataset, with different sources and targets, serves as an independent test set. In training the target-specific virtual screening model, the LIT-PCBA dataset is extracted from three representative targets with a Decoys / Actives ratio greater than 1000 (FEN1: structure-specific nuclease, KAT2A: acetyltransferase, and PKM2: M2-type pyruvate kinase) to better simulate real virtual screening scenarios. For each target, the corresponding test set is directly taken from the original dataset, and the remaining datasets are stratified at a 4:1 ratio to divide them into corresponding training and validation sets.

[0139] (3) For the conformation prediction task, the combined conformations of the training, validation, and test sets of the affinity prediction task were used as the training, validation, and test conformations. To further verify the generalization ability of the model, the model trained on the DUD-E dataset was also tested on three challenging objectives in the LIT-PCBA dataset.

[0140] The primary evaluation metrics are the root mean square error (RMSE) and the Pearson correlation coefficient (Pearson R) to assess the quality of the regression model (in conjunction with the affinity prediction task). For more diverse evaluation methods, the coefficient of determination (R²) is also provided. 2 The calculation methods for these indicators are as follows:

[0141]

[0142]

[0143]

[0144]

[0145]

[0146]

[0147]

[0148] Among them, y i ′ is the predicted binding affinity, y i It is a practical combination of affinity, y and These refer to the average of the predicted value and the experimental value, respectively.

[0149] For structure-based virtual screening tasks, five widely used metrics are reported, including AUC_ROC, AUC_PRC, LogAUC, BEDROC, and EFs at different thresholds (0.1, 0.5, 1, and 5%). For docking conformation prediction tasks, three main metrics are considered: AUC_ROC, AUC_PRC, and Top-1 success rate.

[0150] 2. Using the prepared dataset, the model of this invention is trained, validated, and evaluated through pre-planned implementation, training, and evaluation methods.

[0151] (1) Binding affinity prediction evaluation using the PDBBind dataset is a robust benchmark for evaluating the performance of protein-ligand interaction prediction models. The model of this invention was trained using 8298 complexes from PDBBind V2020 and validated using the remaining 2067 complexes.

[0152] (2) This invention uses two independent external test sets (the core sets of PDBBind V2016 and V2013) to check the generalization ability of the model.

[0153] (3) To evaluate the model in large-scale virtual screening, this invention systematically evaluates the model’s virtual screening capability using the DUD-E, DEKOIS2.0 and LIT-PCBA datasets.

[0154] In the first evaluation phase, a target-free model was trained and validated using the DUD-E dataset, and then tested using the DEKOIS2.0 dataset. Unlike many previous studies, this invention did not use the DUD-E dataset as an external test set because of its potential hidden biases.

[0155] In the second evaluation phase, this invention selected three more challenging and representative targets from the LIT-PCBA dataset, with a decoy-to-active ratio exceeding 1000, to better simulate realistic virtual screening scenarios. The models trained in this phase were target-specific and externally tested using the corresponding test set provided in the original paper.

[0156] In the third evaluation phase, the model trained on the DUD-E dataset was also tested on three challenging targets in the LIT-PCBA dataset to further validate the model's generalization ability.

[0157] The application of this invention is illustrated below by way of example:

[0158] Binding affinity prediction using the PDBBind dataset is a perfect benchmark for evaluating the performance of protein-ligand interaction prediction models. This invention trained the model using 8298 complexes from PDBBind V2020 and validated it using the remaining 2067 complexes. The model's generalization ability was examined using two independent external test sets (core sets from PDBBind V2016 and V2013). The statistics for the model, determined by the minimum validation error criterion, are shown in Table 1.

[0159] Table 1. Different datasets are provided to validate the model's minimum loss performance metric when incorporating affinity prediction.

[0160]

[0161] It can be observed that the model achieved good performance metrics, R0 p The mean square error (RMSE) is >0.750, and the root mean square error (RMSE) on the validation set and the two external test sets is 1.294, 1.220 and 1.372, respectively.

[0162] To further explore the predictive capabilities of the models, this invention analyzes several representative machine learning-based protein-ligand interaction prediction models and lists the corresponding statistics for the PDBBind V2016 and 2013 core sets in Table 2. As shown in Table 2, most reported models are based on CNN algorithms, and some traditional machine learning methods such as Random Forest (RF) and Gradient Boosting Tree (GBT) are also highly favored, but GNN algorithms have received less attention. For the PDBBind V2016 core set, the model's performance in Table 2 is Rp = 0.837 and RMSE = 1.220.

[0163] Table 2. Performance comparison of PDBbind V2016 and V2013 core sets in predicting protein-ligand interactions.

[0164]

[0165] In summary, when constructing the prediction model, this invention (1) directly incorporates the three-dimensional structural information of the protein-ligand binding pose without using heuristic chemical rules to handle non-covalent interactions; it achieves an approximately complete representation of the molecule without increasing the computational budget, combining features based on spherical message passing and physically meaningful features derived from the Schrödinger equation as a unified and effective scheme for extracting the three-dimensional structural information of the protein-ligand binding pose. (2) A distance-aware graph attention mechanism is designed to enable the model to distinguish the contribution of each interaction to the binding affinity. (3) This invention utilizes graphical features obtained by subtracting each feature of the target protein and the given ligand from the graphical features of their complex for protein-ligand interaction prediction. These strategies allow the model to learn key factors for accurately predicting protein-ligand interactions by focusing on intermolecular interactions rather than simply memorizing certain patterns of ligand molecules. (4) This invention improves the performance of the model by employing a gated skip connection mechanism.

[0166] Of course, the above description is not intended to limit the present invention, and the present invention is not limited to the examples given above. Any changes, modifications, additions or substitutions made by those skilled in the art within the scope of the present invention should be protected by the present invention.

Claims

1. A method for predicting three-dimensional protein-ligand activity based on an attention mechanism, characterized by, The method comprises the following steps: Step 1, constructing a data set, taking the three-dimensional structure information of proteins and ligands as input to predict binding affinity, construct a virtual screening model, and simulate a virtual screening scenario; Step 2, constructing a three-dimensional information extraction network for molecular characterization to extract three-dimensional information of protein-ligand complexes; the three-dimensional information extraction network comprises an input module, a message passing module, a plurality of interaction modules, and an output module, the input module is used to construct an initial message for the edge k of the graph structure; The message passing module uses a spherical message passing method and a plurality of physical representations to uniquely represent the molecular graph structure from three-dimensional information, and combines the relative position information of atoms into a spherical coordinate system; the interaction module updates the initial message by combining all physical representations; The output module aggregates all incoming messages to update the features of the nodes; Step 3, constructing a three-dimensional graph neural network model based on an attention mechanism to predict protein-ligand interactions, including conformation prediction and binding affinity prediction; The three-dimensional graph neural network model based on the attention mechanism comprises a graph attention layer, a connection layer, a graph pooling layer, and a fully connected layer, and the following methods are used during prediction: Step 3.1, constructing two adjacency matrices, embedding the extracted three-dimensional information of protein-ligand complexes into the two adjacency matrices, and learning the influence of protein-ligand interactions on the node features of each atom; Step 3.2, designing a distance-aware attention algorithm to distinguish various types of intermolecular interactions; using the graph attention layer to learn protein-ligand interactions, so that the model focuses on intermolecular interactions, and finally obtains node features updated by the gate of the graph attention layer; Step 3.3, after combining the features obtained in step 3.2 to obtain aggregated features, inputting the pair-wise atomic non-bonding interactions learned by the graph pooling layer to generate intermolecular interactions of the complex; Step 3.4, using a batch normalization layer enhanced fully connected layer as a task network for final protein-ligand interaction prediction to generate a prediction result.

2. The method of claim 1, wherein the method is a method of predicting three- dimensional protein-ligand activity based on an attention mechanism. The data set constructed in step 1 includes four data sets: PDBbind, DUD-E, DEKOIS2.0 and LIT-PCBA, the PDBbind data set is used for binding affinity prediction, DUD-E and DEKOIS2.0 are multi-target data sets, which are used to construct a target-free virtual screening model, and the LIT-PCBA data set is used to simulate a real virtual screening scenario.

3. The method of claim 2, wherein the method is based on an attention mechanism. The processing method of the data set is as follows: Step 1.1: For each sample in the PDBbind dataset, perform docking calculations to generate the binding pose of the protein-ligand complex. If the root mean square error (RMSD) of the generated pose compared to the experimentally validated binding structure is less than 100%, then the binding pose is considered valid. Then mark it as a positive sample if the RMSD is greater than 1. Then it is marked as a negative sample, with RMSD between 2 and Samples between them were omitted; Step 1.2, according to the protein, the PDBbind data set is divided into a training set and a test set, so that the training set and the test set do not share the same protein; in addition, the PDBbind samples containing the proteins in the DUD-E data set are deleted from the training set and the test set.

4. The method of claim 1, wherein the method is based on an attention mechanism. The specific method of molecular characterization in step 2 is as follows: Step 2.1, taking the relative position information of atoms as the input of the three-dimensional information extraction network, identifying the structure of the three-dimensional molecule in the spherical coordinate system; For any point in the spherical coordinate system, its position is specified by the triplet where d, θ, and respectively represent the radial distance, the polar angle, and the azimuthal angle; when modeling a three-dimensional molecular graph in the spherical coordinate system, any atom i can be the origin of a local spherical coordinate system, in which case d, θ, and become the bond length, the bond angle, and the torsion angle, respectively; Thus, the relative position of each neighboring atom of atom i can be represented by the respective tuple designated; Step 2.

2. Define a three-dimensional molecular graph G, denoted as a 4-tuple G = (u, V, E, P), is the global feature vector of the molecular graph G, where d u is the dimension of the global feature; V = {v i} i＝1:n is the atomic feature, where each is the feature vector of atom i, and n is the total number of atoms; E = {(e k , r k , s k} k＝1:m is the set of edges, where each is the edge feature vector of edge k, r k is the receiver atom, s k is the sender atom of edge k, and m is the total number of edges; P = {r k} k＝1:n is a set of three-dimensional Cartesian coordinates containing three-dimensional spatial information for each atom; Step 2.3, performing spherical message passing in the spherical coordinate system to update the message; Let denote the edge set pointing to atom i, N i denote the index of incoming node of atom i; the output after the message passing procedure includes the updated global feature vector updated atom feature V' = {v' i} i＝1:n and updated edge set E' = {(e' k ,r k ,s k )} k＝1:m where v' i denotes the updated atom feature vector, e' k denotes the updated edge feature vector, the spherical message passing is defined in the spherical coordinate system as shown in equation (1): where φ e ,φ v and φ u are three information update functions on edges, atoms and the whole molecule graph respectively, denotes the feature vector of the receiver atom r k of edge k, denotes the feature vector of the sender atom s k of edge k, denotes the set of edges pointing to the sender atom s k of edge k, r h denotes the receiver atom of edge h, denotes the index of the incoming node pointing to the sender atom s k of edge k, E i ' denotes the updated set of edges pointing to atom i; the aggregation functions p e→v and p v→u aggregate information between different geometric types, in spherical message passing, the three-dimensional information in P is transformed and merged to update each message; spherical message passing uses another position aggregation function p p→e for message update; Step 2.

4. Select basis functions, convert triplets to physics-based representations, perform molecular characterization.

5. The method of claim 4, wherein the method is based on an attention mechanism. In step 2.4, the spherical Bessel function is chosen as the basis function for the encoding distance, and the spherical harmonic function is chosen as the basis function for the encoding angle; the final molecular representation is expressed as One representation is expressed as shown in equation (2): where is a spherical Bessel function of order is a spherical harmonic of order c denotes the cut-off radius, is the nth root of the Bessel function of order and n e [1,..., N B ], L H and N B denote the highest order of the spherical harmonic and the spherical Bessel function, respectively.

6. The attention-mechanism-based three-dimensional protein-ligand activity prediction method according to claim 1, wherein, In step 3.1, the adjacency matrix, denoted as A, is an N matrix x N matrix matrix, where A ij > 0 if the ithand jthatomic node is connected, and A ij = 0 otherwise, N matrix denotes the number of atomic nodes in the graph; two adjacency matrices A 1 and A 2 are constructed, A 1 represents the adjacency matrix of pure covalent interactions, and A 2 represents the adjacency matrix of covalent interactions and non-covalent intermolecular interactions, the construction of A 1 and A 2 is shown in equations (3), (4): where d ij is the distance between the ith and jth atoms, μ and σ are learnable parameters, and the formula in equation (4) reflects that intermolecular bonds are weaker than covalent bonds, and their strength decreases as the bond distance increases.

7. The attention-mechanism-based three-dimensional protein-ligand activity prediction method according to claim 6, wherein, In step 3.2, a distance-aware attention algorithm is designed, wherein the input of the graph attention layer is the adjacency matrix A and the node feature set where x e R F , N atom is the number of nodes, i.e. the number of atoms, and F is the dimension of the node features; the graph attention layer generates a new set of node features where x e R F ; The specific method is as follows: Step 3.2.1, calculate the attention coefficient of each node; First, in order to obtain sufficient expressiveness, each node feature is transformed by a learnable weight matrix W ∈ R F×F is transformed by a learnable weight matrix W ∈ R Then, the attention coefficient is obtained as shown in equation (5): where E ∈ R F×F is also a learnable matrix, attention coefficient e ij represents the importance of the jth node feature to the ith node feature; by summing and , e ij is forced to be ji ; in order to reflect the structure of the graph, attention coefficient e ij is only calculated for j ∩ N i , where N i is the neighbor node of the ith node, A ij > 0; Then, to manipulate the scale of the attention coefficients of each node, the attention coefficients are normalized among neighbors, and A ij is multiplied by the normalized attention coefficients a ij is given by equation (6): Step 3.2.2, update the node feature; After obtaining the normalized attention coefficients a ij Each node feature is then updated as a linear combination of the node features of the neighboring nodes and the normalized attention coefficients, as shown in equation (7) Step 3.2.3, introduce the gating enhancement algorithm; A gating enhancement algorithm is introduced to directly pass the feature information of the previous node to the next layer, and the output of the graph attention layer is implemented as x i and x i The linear combination of and x is shown in formulas (8) and (9). x i out = z i x' i + (1 - z i )x" i (8) z i =σ(U(x′ i ||x″ i )+b) (9) combining the outputs of the graph attention layers for all nodes composing a new set of node features where U ∈ R 2F×1 is a learnable vector, b is a learnable scalar value, σ denotes a sigmoid activation function, (· || ·) is concatenation of two vectors, z i denotes how much information about the input node features will be directly passed to the next layer.

8. The method of claim 7, wherein the method is based on an attention mechanism. In step 3.2, for two adjacency matrices A 1 and A 2 , the graph attention layer generates new node features x out1 and x in , x 1 = GAT(x out2 , A in ), x 2 = GAT(x out1 , A out2 ), where GAT represents a gated graph attention layer, and one gated graph attention layer is shared when computing x out and x out2 ; Output node feature x out is obtained by subtracting x out2 from x out1 , as shown in equation (10), x out = x out2 - x out1 (10) After updating the node feature vector through multiple graph attention layers, the feature vector of the ligand atom is summarized into a vector representing the protein-ligand complex graph, as shown in equation (11): Finally, a multilayer perceptron is applied to the x graph The classification whether a protein-ligand complex or binding pose is active uses Leaky ReLU activation functions between the layers and a sigmoid function after the last layer.

Images